How to declare time dependent symbolic array?

39 次查看(过去 30 天)
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?
This question complements another.

回答(4 个)

Marlon Saveri Silva
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
  4 个评论
Walter Roberson
Walter Roberson 2018-10-4
编辑:Walter Roberson 2018-10-4
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
Walter Roberson 2018-10-4
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

请先登录,再进行评论。


Sven
Sven 2018-10-4
编辑:Sven 2018-10-4
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

ahmad mahmoodi
ahmad mahmoodi 2018-10-13
编辑:Walter Roberson 2018-10-13
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!
  2 个评论
Walter Roberson
Walter Roberson 2018-10-13
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
ahmad mahmoodi 2018-10-13
编辑:Walter Roberson 2018-10-14
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

请先登录,再进行评论。


Andrés Arciniegas
Andrés Arciniegas 2019-6-16
Greetings. This simple loop generates a vector of time dependent variables. Works smoothly and might be useful to someone :).
Addressing the problem above, the script generates the vector: .
% Script to create time dependent variables
Tt =[] % Array
syms(var_ind)
for i=1:n
var_str = strcat('T',num2str(i),'(t)'); % Builds the symbolic variable string
syms(var_str);
Tt = [Tt str2sym(var_str)]; % store the variable in the array
end

类别

Help CenterFile Exchange 中查找有关 Assumptions 的更多信息

产品

Community Treasure Hunt

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

Start Hunting!

Translated by