3 visualizzazioni (ultimi 30 giorni)
Ananthakrishnan il 6 Mag 2024
Spostato: Torsten il 6 Mag 2024
Tried to solve a system of PDEs
Then Weq is depended on Td by equation :
where ϕ(i)=94.92211.*q(i)^4 - 76.5155.*q(i)^3 + 19.08577.*q(i)^2 + 0.03858.*q(1)+ 0.14278 ,ϕ(phi , relative humidity varying by input of q)
My code : its running but answers not coming as expected , probably cause information flow is not right. please show help me with it
clc
clear
close all
a = 0.05;%m
b = 0.025;%m
Lp = 2*b + 2*(sqrt(b^2 + (a*pi)^2))*((3 + (2*b/a*pi)^2)/(4 + (2*b/a*pi)^2));
A = 2*a*b;
Pa = 101325; % Pa
Ta = 30; % Celsius
D = (2.256/Pa)*((Ta+273.13)/256)^1.81;
V = 2;
rhoa = 1.15; % kg/m^3
Sha = 2.45;
Ky = rhoa*Sha*D*Lp/4*A;
Nua = 2.45;
ka = 0.0321; % W/(m K)
h = Nua*ka*Lp/4*A;
fd = 0.914; % kg/m
fm = 0.101; % kg/m
cd = 0.9; % J/kg K
cw = 4.19; % J/kg K
cm = 0.9; % J/kg K
ca = 1.005; % J/kg K
cv = 2.028; % J/kg K
w10=0.0133; % intital values
Ta10=303;
Td10=303;
RH10=0.5 ;
q10=0.01;
%increment for lenght
Nz=100;
z=linspace(0,1,Nz);
dz=z(2)-z(1);
%time
t=0:3:299;%secs
%inital conditions all zeros
ICi=zeros(1,Nz);
IC=[ICi,ICi,ICi,ICi].';
%4 variables os IC must be 4 c1-c6
%odesolver
%re calculation
%define value
q=y(:,1:Nz);
w=y(:,Nz+1:2*Nz);
Td=y(:,2*Nz+1:3*Nz);
Ta=y(:,3*Nz+1:4*Nz);
%bounday cond
q(:,1)=q10; %kg/kgda
w(:,1)=w10; %celcius
Td(:,1)=Td10; %celcius
Ta(:,1)=Ta10;
q(:,end)=w(:,end-1);
w(:,end)=w(:,end-1);
Ta(:,end)=Ta(:,end-1);
Td(:,end)=Td(:,end-1);
%plot
plot(Td,Ta)
%PDE-ODe
dydt=zeros(length(y),1);
dqdt=zeros(Nz,1);
dwdt=zeros(Nz,1);
dTddt=zeros(Nz,1);
%define value
q=y(1:Nz);
w=y(Nz+1:2*Nz);
Td=y(2*Nz+1:3*Nz);
Ta=y(3*Nz+1:4*Nz);
%bounday condition
q(1)=q10;%kg/kg
w(1)=w10;%30 celcius
Td(1)=Td10;%30celcius
Ta(1)=Ta10;
q(end)=q(end-1);
w(end)=0.09;
Ta(end)=297;
Td(end)=Td(end-1);
%interior
for i=2:Nz-1
pws(i)=exp(23.196-(3816.44/(Td(i)-46.13)));
RH(i)=94.92211.*q(i)^4 - 76.5155.*q(i)^3 + 19.08577.*q(i)^2 + 0.03858.*q(1)+ 0.14278;
weq(i)=0.62188*RH(i)./((Pa./pws(i))-RH(i));
dqdt(i)=-2*Ky*Pa/fd*(weq(i)-w(i)); % one of the sysetm equations but its an ODE
dqdz(i)=(q(i+1)-q(i-1))./2./dz;
dwdz(i)=(w(i+1)-w(i-1))./2./dz;
dTddz(i)=(Td(i+1)-Td(i-1))./2./dz;
dwdt(i)=-V.*dwdz(i)-(fd/2*A*rhoa).*dqdt(i);
end
end
##### 0 CommentiMostra -2 commenti meno recentiNascondi -2 commenti meno recenti

Accedi per commentare.

### Risposte (1)

Torsten il 6 Mag 2024
Spostato: Torsten il 6 Mag 2024
For Ta and w, you can only prescribe a value at z = 0, not at z = 1.
Further use upwind differencing for the first-order derivatives, not central differencing for Ta and w. Central differencing will introduce oscillations for the solutions.
The equations for Td and q are ordinary differential equations for which you are not allowed to prescribe any boundary values.
##### 0 CommentiMostra -2 commenti meno recentiNascondi -2 commenti meno recenti

Accedi per commentare.

### Categorie

Scopri di più su Structural Mechanics in Help Center e File Exchange

R2024a

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by