您好,我正在尝试设置一个 ode45 来求解这个微分方程 ρaCaVdTin/dT= Qshort - Qconv,cond - Qinfilt 其中 Qshort Qconv 和 Qinflit 都是描述温室中热量损失或增益的方程。
我想要实现的是通过给出外部温度风速和光强度的列向量变量,我希望 ode45 为我提供外部温度风速等对应值的内部温度,即我希望 Tin(1) 来来自 ode45,它将使用我的 tspan 中的第一个时间值、第一个外部温度等。
我通过在函数中输入断点对其进行编码的方式我注意到,从我的代码第一次执行开始,内部温度向量已经是一个 297x1 列向量,长度与我的 tspan 相同。越过 dT 的断点会导致错误初始条件向量必须与返回的向量具有相同数量的元素。问题是,有没有一种方法可以实现我的预期结果,正如我上面所说的 ode45 使用我的向量的一个值来获得一个结果?抱歉我的英语不好,这不是我的第一语言,提前感谢。
L=4;
W=3.7;
H=3;
volup=1*(1*4)/2; %ogkos trigwnikhs orofhs
voldown=W*L*2; %ogkos katw or8ogwniou parellhlepipedou
V=voldown+volup; %Sunolikos ogkos
taxythtaaera = importdata('taxythta aera.txt');
exwterikh=importdata('metablhtes.txt');
ilios=importdata('hlios.txt');
S=L*W;
Vw=taxythtaaera*0.27777777777778;
Tac=15; %8ermokrasia se kelsioy
Ta=exwterikh+273.15; %se kelvin kai 15 se kelsiou kserw gw
Ta(265)=[];
Ta(143)=[];
Ta(53)=[];
Pa=1.137; %air density inside
Ca=1005; %specific heat of air
sta=Pa*Ca*V;
ac=0.1; %cover absorbtivity
tc=0.85; %cover transmittance
I=ilios; %kanonika pinakas olikh iliakh aktinobolia w/m^2
Sc=H*W;
Tsky= 0.0552.*(Ta.^1.5);
Tskyc=Tsky-300;
R=0.75; %air changes per hour APO FUSIKO AERISMO
Qshort=ac*tc*S*I;
stat=Qshort/sta;
stat2=S/sta;
stat3=5.2*(R/((Sc*L))^1/2);
stat4=R/V;
ho=2.8+1.2.*Vw;
ho(132)=[];
ho(53)=[];
ho(274)=[];
stat5=((ho*S)/sta)*(1-tc);
Kc=0.028; %cover K w*m^-1*k^-1
Lc=0.03; %cover thickness
%dT=Qshort/pa*ca*V-U*(S/Pa*Ca*V)(T-Tout)-(Pa*Ca*R/Pa*Ca*V)*(T-Ti)/3600-(ho*S*(1-tc)/(Pa*Ca*V)(T-Tsky)
t0=0;
To=286.15;
tend=35640;
n=297; %giati allazei toso to apotelesma analoga me ta deigmata? sta 300 fainetai arketa sta8ero
tspan=linspace(t0,tend,n);
[tsol,Tsol]=ode45(@(t,T) funcODE(t,T,ho,Ta,I),tspan,To);
plot(tsol,Tsol-273.15)
x=Tsol-273.15;
这是我的主要程序
这是我的功能
function dT= funcODE(t,T,ho,Ta,I)
Tout=Ta;
ho=ho;
I=I;
L=4;
W=3.7;
H=3;
volup=1*(1*4)/2; %ogkos trigwnikhs orofhs
voldown=W*L*2; %ogkos katw or8ogwniou parellhlepipedou
V=voldown+volup; %Sunolikos ogkos
S=L*W;
i=1;
Pa=1.137; %air density inside
Ca=1005; %specific heat of air
sta=Pa*Ca*V;
ac=0.1; %cover absorbtivity
tc=0.85; %cover transmittance
Sc=H*W;
Tsky= 0.0552.*(Tout.^1.5);
R=0.75; %air changes per hour APO FUSIKO AERISMO
Qshort=ac*tc*S*I;
stat=Qshort/sta;
stat2=S/sta;
stat3=5.2*(R/((Sc*L)^1/2));
stat4=R/V;
stat5=((ho*S)/sta)*(1-tc);
Kc=0.028; %cover K w*m^-1*k^-1
Lc=0.03; %cover thickness
U=(1.0714+(1/ho)+(1/(1.52*(T-Tout).^1/3+stat3))).^-1;
dT=stat -(((R/(V*3600))*(T-Tout(i)))+(((ho*S*(1-tc))/sta)*(T-Tsky(i)))+(S*(T-Tout(i))*U(i)));
i=i+1;