Проблема с использованием ode45 для решения дифференциального уравнения (ошибка с вектором начальных условий?) - PullRequest
Новая
       11

Проблема с использованием ode45 для решения дифференциального уравнения (ошибка с вектором начальных условий?)

0 голосов
/ 07 января 2019

Здравствуйте, я пытаюсь настроить ode45 для решения этого дифференциального уравнения

То, чего я пытаюсь достичь, это дать векторную переменную столбца скорости ветра и интенсивности света наружной температуры. Я хочу, чтобы ode45 дал мне внутреннюю температуру для соответствующего значения скорости ветра наружной температуры и т. Д., Т.е. я хочу олово (1 ) из ode45, который будет использовать первое значение времени из моего tspan при первой внешней температуре и т. д.

То, как я кодировал его, вводя точку останова внутри функции, я заметил, что с первого выполнения моего кода вектор внутренней температуры уже имеет размер столбца 297x1 такой же длины, что и мой tspan., проходящий после точки останова результатов dT. в ошибке, что вектор начальных условий должен иметь то же количество элементов, что и возвращаемый вектор. Вопрос в том, есть ли способ достичь желаемого результата, как я уже говорил выше для ode45, использовать одно значение из моих векторов для получения одного результата? Извините за мой плохой английский, это не мой родной язык, и спасибо заранее.

error and sizes of my vectors

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;
...