Asked by Marlon Saveri Silva
on 31 Mar 2018

Hello,

I need to solve the following differential equation

odes=diff(T) == A*T + B;

[VF,Subs] = odeToVectorField(odes);

Sys = matlabFunction(VF,'Vars',{'t','Y'});

[T,Y] = ode45(Sys, [0 5], T0);

Therefore, I need to declare T(t). It could be done as below

syms T1(t) T2(t) T3(t) T4(t) T5(t)

However, the the length of T vector is not necessarily 5, but an user defined number. So, I tried the following:

syms t; T(t)=sym('T',[Num,1]);

But T was not considered as a t function and. When "odes" differentiates it, dT/dt is considered zero and I can't solve the problem.

So, how to create time dependent symbolic variables in an array with Num lines?

Answer by Marlon Saveri Silva
on 31 Mar 2018

I got a solution from another topic, but I'm not sure if it's the best, since some warnings appears on screen:

T=sym([]);

for i=1:N

syms(sprintf('T%d(t)', i)) %declare each element in the array as a single symbolic function

TT = symfun(sym(sprintf('T%d(t)', i)), t); %declare each element to a symbolic "handle"

T = [T;TT]; %paste the symbolic "handle" into an array

end

Q
on 25 Sep 2018

Walter Roberson
on 4 Oct 2018

T=cell(N,1);

for I=1:N

syms(sprintf('T%d(t)', i)) %declare each element in the array as a single symbolic function

TT = symfun(sym(sprintf('T%d(t)', i)), t); %declare each element to a symbolic "handle"

T{I} = TT; %paste the symbolic "handle" into an array

end

Walter Roberson
on 4 Oct 2018

For R2018a and later,

T = cell(N,1);

for I=1:N

TT = symfun(str2sym(sprintf('T%d(t)', i)), t); %declare each element to a symbolic "handle"

T{I} = TT; %paste the symbolic "handle" into an array

end

Sign in to comment.

Answer by Sven
on 4 Oct 2018

Edited by Sven
on 4 Oct 2018

Hi Marlon,

if you substitute sym with eval, defining T, then there is no warning anymore. Also, @Q please see below where I used a cell array.

for i=1:5

syms(sprintf('T%d(t)', i))

TT = symfun(eval(sprintf('T%d(t)', i)), t); % use eval instead of sym

T{i} = TT; % use cell array

end

Sign in to comment.

Answer by ahmad mahmoodi
on 13 Oct 2018

Edited by Walter Roberson
on 13 Oct 2018

Hello all,

I used the same method. However I can not multiply my cell array of coefficients (AB) into the cell array of symbolic functions (T, Tm). In the following is the abstract of my code:

T=cell(4*n,1);

Tm=cell(6*n,1);

AB=cell(4*n,6*n);

GF=cell(4*n,1);

for i=1:n

AB, Tm an GF (using your comments) and the relation between them are defined here

end

odes=diff(T)==AB*Tm+GF

or

(odes=diff(T)==AB.*Tm+GF;)

errors:

Undefined operator '*' for input arguments of type 'cell'.

(Undefined operator '.*' for input arguments of type 'cell'.)

Also I tried to use matrix instead as follow:

T=sym(zeros(4*n,1));

Tm=sym(zeros(6*n,1));

for i=1:n

AB, Tm an GF using your comments as matrices and the relation between them are defined here

end

odes=diff(T)==AB*Tm+GF;

[VF,Subs] = odeToVectorField(odes);

Sys = matlabFunction(VF,'Vars',{'t','Y'});

Y0(1:4*n,1)=14.7;

[T,Y] = ode45(Sys, [0 5], Y0);

and here is the error:

Error using mupadengine/feval (line 163)

Cannot convert the initial value problem to an equivalent dynamical system. Either the differential equations cannot be solved for the highest derivatives or

inappropriate initial conditions were specified.

Error in odeToVectorField>mupadOdeToVectorField (line 177)

T = feval(symengine,'symobj::odeToVectorField',sys,x,stringInput);

Error in odeToVectorField (line 126)

sol = mupadOdeToVectorField(varargin);

Error in Fb1 (line 154)

[VF,Subs] = odeToVectorField(odes);

Thanks in advance for your reply!

Walter Roberson
on 13 Oct 2018

It seems unlikely to me that Mathworks will ever define a meaning for the * operator between a cell array and anything else.

The point of using cell arrays was that the original poster was trying to define a vector of symbolic functions. MATLAB does not have any problem with vectors of symbolic expressions, but when you have vectors of symbolic functions, then the meaning of

TheVector(Value)

becomes ambiguous. Does it mean to index the vector at the given value, pulling out one of the functions? Or does it mean to execute all of the functions in turn, passing in the Value to each one? MATLAB decided that it should mean executing all of them in turn. Because of that, you cannot index into a vector of symbolic functions. A vector of symbolic functions becomes the same as if you had created a single symbolic function whose body was the vector.

When you are working with differential equations, chances are that what you need is a vector of symbolic expressions, rather than a vector of functions.

... Basically, you should probably be avoiding creating functions dynamically.

Cannot convert the initial value problem to an equivalent dynamical system. Either the differential equations cannot be solved for the highest derivatives or

inappropriate initial conditions were specified.

We will probably need your actual code (and the inputs) to analyze that.

ahmad mahmoodi
on 13 Oct 2018

Hello Walter. Thank you for your quick response.

My problem consists of layers of differential equations. For the layers between each layer has 4 differential equations and 4 variables. The number of equations and the coefficients in the first and the last layer differ from the ones between (3 equations and 4 variable for the first layer and 4 equations 3 variables for the last layer). So totally I have n equations and n variables. Also two of the variables from each that layer are used in the next layer. So I have a for loop that I need to use if loop inside it. All of my variables are time dependant.

I solved the problem with just two layers by entering the equations manually using ode45 and the results match. But I need my code to calculate when I have multiple layers . Here is the code:

% parameters

c_p_grout=1500;

rho_grout=1460;

lambda_grout=1.5;

c_p_fluid=4183;

rho_fluid=997.6;

lambda_soil=2.51225;

lambda_fluid=0.5913;

m_dot=0.46;

L=100;

dl=50;

n=L/dl+1;

d_b=0.2;

d_a=0.04;

d_i=d_a-2*0.0037;

d_s1=0.22;

d_z1=sqrt((d_b^2+d_s1^2)/2);

R_conv=0.00670;

R_cond_1=0.08568;

x=0.70275;

R_a=0.23701;

R_b=0.073;

R_g=4*R_b-R_conv-R_cond_1;

R_gg=1/((1/(R_a-R_conv-R_cond_1-x*R_g))-(1/(1-x)/R_g))*4;

R_fg=R_conv+R_cond_1+x*R_g;

R_gb=(1-x)*R_g;

R_gg1=R_gg;

R_gg2=R_gg1;

R_fgred=R_fg/2;

R_gbred_wo_soil=R_gb/2;

R_ggred_wo_soil=1/((2/R_gg1)+(2/R_gg2));

R_bz1=log(d_z1/d_b)/(2*pi*lambda_soil);

R_gbred=R_gbred_wo_soil+2*R_bz1;

R_ggred=1/((1/R_ggred_wo_soil)+(1/2/R_gbred_wo_soil)-(1/2/R_gbred));

C_fluid=2*rho_fluid*c_p_fluid*pi/4*d_i^2*dl/2; %C_fluid at layer 1

C_g=rho_grout*(pi/4)*((d_b^2/2)-(2*d_a^2))*c_p_grout*dl/2; % In EES 2*d_a^2 instead of d_a^2

T_b=15;

Tinlet=90;

A=(-dl/2)*((1/C_g*R_ggred)+(1/C_g*R_gbred)+(1/C_g*R_fgred));

B=(dl/2)*(1/C_g*R_fgred);

C=(dl/2)*(1/C_g*R_ggred);

D=(-dl/2)*((1/C_fluid*R_fgred))-(m_dot*c_p_fluid/C_fluid);

E=(m_dot*c_p_fluid/C_fluid);

F=(dl/2)*(1/C_fluid*R_fgred);

G=(dl/2)*(T_b/C_g*R_gbred);

A1=2*A;

B1=2*B;

C1=2*C;

F1=2*F;

Gi=2*G;

D1=(-dl)*((1/C_fluid*R_fgred))-(m_dot*c_p_fluid/C_fluid);

T=sym(zeros(4*n,1));

Tm=sym(zeros(6*n,1));

for i=1:n

syms(sprintf('T%d(t)', 4*i-3))

TT=symfun(sym(sprintf('T%d(t)', 4*i-3)), t);

T(4*i-3,:)=TT;

syms(sprintf('T%d(t)', 4*i-2))

TT1=symfun(sym(sprintf('T%d(t)', 4*i-2)), t);

T(4*i-2,:)=TT1;

syms(sprintf('T%d(t)', 4*i-1))

TT2=symfun(sym(sprintf('T%d(t)', 4*i-1)), t);

T(4*i-1,:)=TT2;

syms(sprintf('T%d(t)', 4*i))

TT3=symfun(sym(sprintf('T%d(t)', 4*i)), t);

T(4*i,:)=TT3;

Tm(6*i-5,:)=TT;

Tm(6*i-4,:)=TT1;

% Tm(6*i-3,:) inside if loop

Tm(6*i-2,:)=TT2;

Tm(6*i-1,:)=TT3;

% Tm(6*i,:) inside if loop

if i==1

%Tg1(i)=A*Tg1(i)+B*T1(i)+C*Tg2(i)+G;

%Tg2(i)=A*Tg2(i)+B*T2(i)+C*Tg1(i)+G;

%dT2(i)=D*T2(i)+E*T2(i+1)+F*Tg2(i);

T(4*i-1,:)=Tinlet;

Tm(6*i-3,:)=0;

syms(sprintf('T%d(t)', 4*i+4))

TT6=symfun(sym(sprintf('T%d(t)', 4*i+4)), t);

Tm(6*i,:)=TT6;

AB1=[A C 0 B 0 0; C A 0 0 B 0 ; 0 0 0 0 0 0; 0 F 0 0 D E];

G1=[G; G; 0; 0];

AB(4*(i-1)+1:4*(i-1)+4,6*(i-1)+1:6*(i-1)+6)=AB1;

GF(4*(i-1)+1:4*(i-1)+4,1)=G1;

elseif i==n

% Tg1(i)=A*Tg1(i)+B*T1(i)+C*Tg2(i)+G;

% Tg2(i)=A*Tg2(i)+B*T2(i)+C*Tg1(i)+G;

% T1(i)=D*T1(i)+E2*T1(i-1)+F*Tg1(i);

Tm(6*i,:)=0;

syms(sprintf('T%d(t)', 4*i-5))

TT7=symfun(sym(sprintf('T%d(t)', 4*i-5)), t);

Tm(6*i-3,:)=TT7;

Tm(6*i-1,:)=Tm(6*i-2,:);

T(4*i-1,:)=T(4*i,:);

AB3=[A C 0 B 0 0; C A 0 0 B 0 ; F 0 E D 0 0; 0 0 0 0 0 0];

G3=[G; G; 0; 0];

AB(4*(i-1)+1:4*(i-1)+4,6*(i-1)+1:6*(i-1)+6)=AB3;

GF(4*(i-1)+1:4*(i-1)+4,1)=G3;

else

% Tg1(i)=A*Tg1(i)+B*T1(i)+C*Tg2(i)+G;

% Tg2(i)=A*Tg2(i)+B*T2(i)+C*Tg1(i)+G;

% T1(i)=D*T1(i)+E2*T1(i-1)+F*Tg1(i);

% T2(i)=D*T2(i)+E2*T2(i+1)+F*Tg2(i)

syms(sprintf('T%d(t)', 4*i-5))

TT8=symfun(sym(sprintf('T%d(t)', 4*i-5)), t);

Tm(6*i-3,:)=TT8;

syms(sprintf('T%d(t)', 4*i+4))

TT9=symfun(sym(sprintf('T%d(t)', 4*i+4)), t);

Tm(6*i,:)=TT9;

AB2=[A1 C1 0 B1 0 0; C1 A1 0 0 B1 0 ; F1 0 E D1 0 0; 0 F1 0 0 D E];

G2=[Gi; Gi; 0; 0];

AB(4*(i-1)+1:4*(i-1)+4,6*(i-1)+1:6*(i-1)+6)=AB2;

GF(4*(i-1)+1:4*(i-1)+4,1)=G2;

end

end

odes=diff(T)==AB*Tm+GF;

[VF,Subs] = odeToVectorField(odes);

Sys = matlabFunction(VF,'Vars',{'t','Y'});

Y0(1:4*n,1)=14.7;

[T,Y] = ode45(Sys, [0 5], Y0);

I also tried to solve using

function dT = SRBZ(t, T)

the previous loops and coeficients here...

end

T0(1:4*n,1)=14.7;

[t,T] = ode45(@SRBZ,[0 2],T0);

But it gives error again and can not solve

Sign in to comment.

Opportunities for recent engineering grads.

Apply Today
## 0 Comments

Sign in to comment.