Substituing symbolic variable with an element of a matrix

Asked by VincentLec on 27 Jun 2016

Latest activity Edited by VincentLec on 28 Jun 2016

1 view (last 30 days)

So, I have this ODE system (named Y):

Y=
2.09e-4*thetaxpp*cos(4.0*sin(3.0*t)) - 9.16e-4*thetaypp*sin(4.0*sin(3.0*t)) - 0.00125*thetayp*cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125*thetaxp*cos(3.0*t)*sin(4.0*sin(3.0*t))
2.09e-4*thetaypp*cos(4.0*sin(3.0*t)) + 9.16e-4*thetaxpp*sin(4.0*sin(3.0*t)) + 0.00125*thetaxp*cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125*thetayp*cos(3.0*t)*sin(4.0*sin(3.0*t))
                  3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))*thetaxp^2 + 3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))*thetayp^2 + 2.09e-4*thetazpp - 1.11*sin(3.0*t)

Now, I want to use the ode45i function on it, which means that I want to do:

f2=@(t,y,yp)(Y,[0 10],y0,yp0)

To do so, I need to perform these substitutions (substitute a symbolic variable with an element of a matrix) before implementing Y in the function ode15i:

thetaxpp=yp(1)
thetaypp=yp(2)
thetazpp=yp(3)
thetaxp=y(1)
thetayp=y(2)
thetazp=y(3)

(my calculations are mostly inspired by this)

Now, I've tried to do:

f2=@(t,y,yp)(subs(Y,[thetaxpp, thetaypp, thetazpp, thetaxp, thetayp, thetazp],[yp(1), yp(2) , yp(3), y(1), y(2), y(3)]))

but I get the message

Error using odearguments (line 101)
Inputs must be floats, namely single or double.

I also tried isolating the equations of Y using fprintf (to transform the equations into text before doing the substitution and to transform them again into equations after) like this:

formatSpec = '\n %4.2f \n %4.3f \n %4.4f\n';
fprintf(formatSpec,Y(1),Y(2),Y(3))

but I got the message:

Error using fprintf
Function is not defined for 'sym' inputs.

Do you guys have any idea what I could do to solve my problem?

Thanks a bunch

0 Comments

1 Answer

Link

Answer by Star Strider on 27 Jun 2016

You don’t need to perform the substitutions! Let MATLAB do it! Use the matlabFunction function, and create an anonymous function to use with your ODE.

See if this works for you:

syms thetaxpp thetaypp thetaxp thetayp thetazpp thetazp t
Y = [2.09e-4*thetaxpp*cos(4.0*sin(3.0*t)) - 9.16e-4*thetaypp*sin(4.0*sin(3.0*t)) - 0.00125*thetayp*cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125*thetaxp*cos(3.0*t)*sin(4.0*sin(3.0*t))
2.09e-4*thetaypp*cos(4.0*sin(3.0*t)) + 9.16e-4*thetaxpp*sin(4.0*sin(3.0*t)) + 0.00125*thetaxp*cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125*thetayp*cos(3.0*t)*sin(4.0*sin(3.0*t))
                  3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))*thetaxp^2 + 3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))*thetayp^2 + 2.09e-4*thetazpp - 1.11*sin(3.0*t)];
Yfcn = matlabFunction(Y);

You will have to add ‘y’ as the second argument in the argument list of ‘Yfcn’ to create:

Yfcn = @(t,y,thetaxp,thetayp,thetaxpp,thetaypp,thetazpp)[thetaxpp.*cos(sin(t.*3.0).*4.0).*2.09e-4-thetaypp.*sin(sin(t.*3.0).*4.0).*9.16e-4-thetayp.*cos(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2)-thetaxp.*sin(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2);thetaypp.*cos(sin(t.*3.0).*4.0).*2.09e-4+thetaxpp.*sin(sin(t.*3.0).*4.0).*9.16e-4+thetaxp.*cos(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2)-thetayp.*sin(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2);thetazpp.*2.09e-4-sin(t.*3.0).*(1.11e2./1.0e2)+thetaxp.^2.*cos(sin(t.*3.0).*4.0).*sin(sin(t.*3.0).*4.0).*3.54e-4+thetayp.^2.*cos(sin(t.*3.0).*4.0).*sin(sin(t.*3.0).*4.0).*3.54e-4];
tspan = [ ... ];
initcond = [0; 0; 0];
[t,y] = ode15i(@(t,y) Yfcn(t,y,thetaxp,thetayp,thetaxpp,thetaypp,thetazpp), tspan, initcond);

Be certain to define:

thetaxp,thetayp,thetaxpp,thetaypp,thetazpp

in your workspace first.

NOTE — This is UNTESTED CODE but it should work if I understand your Question correctly.

1 Comment

VincentLec on 28 Jun 2016

Converting the symbolic function into a Matlab function is actually a really good idea (not sure why I haven't thought about it). However, when I tried using your code (and adding the initial conditions of yp considering that when using ode15i, you need two vectors of initial conditions), I got the following message:

Error using @(t,y)Yfcn(t,y,thetaxp,thetayp,thetaxpp,thetaypp,thetazpp)
Too many input arguments.
Error in odearguments (line 87)
f0 = feval(ode,t0,y0,args{:});   % ODE15I sets args{1} to yp0.
Error in ode15i (line 116)
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, ~, odeFcn, ...
Error in BoucleForGrosseMatrice (line 237)
[t,y] = ode15i(@(t,y) Yfcn(t,y,thetaxp,thetayp,thetaxpp,thetaypp,thetazpp), tspan, initcond,initcond2);

Which then lead me to adding yp as a variable, such as

Yfcn = @(t,y,yp,thetaxp,thetayp,thetaxpp,thetaypp,thetazpp)[thetaxpp.*cos(sin(t.*3.0).*4.0).*2.09e-4-thetaypp.*sin(sin(t.*3.0).*4.0).*9.16e-4-thetayp.*cos(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2)-thetaxp.*sin(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2);thetaypp.*cos(sin(t.*3.0).*4.0).*2.09e-4+thetaxpp.*sin(sin(t.*3.0).*4.0).*9.16e-4+thetaxp.*cos(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2)-thetayp.*sin(sin(t.*3.0).*4.0).*cos(t.*3.0).*(1.0./8.0e2);thetazpp.*2.09e-4-sin(t.*3.0).*(1.11e2./1.0e2)+thetaxp.^2.*cos(sin(t.*3.0).*4.0).*sin(sin(t.*3.0).*4.0).*3.54e-4+thetayp.^2.*cos(sin(t.*3.0).*4.0).*sin(sin(t.*3.0).*4.0).*3.54e-4];

and

[t,y] = ode15i(@(t,y,yp) Yfcn(t,y,yp,thetaxp,thetayp,thetaxpp,thetaypp,thetazpp), tspan, initcond,initcond2);

However, this gave me the following message

 Error using odearguments (line 101)
Inputs must be floats, namely single or double.

Also, to clarify, what I wish to do is to transform this symbolic matrix

 Y = [2.09e-4*thetaxpp*cos(4.0*sin(3.0*t)) - 9.16e-4*thetaypp*sin(4.0*sin(3.0*t)) - 0.00125*thetayp*cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125*thetaxp*cos(3.0*t)*sin(4.0*sin(3.0*t))
2.09e-4*thetaypp*cos(4.0*sin(3.0*t)) + 9.16e-4*thetaxpp*sin(4.0*sin(3.0*t)) + 0.00125*thetaxp*cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125*thetayp*cos(3.0*t)*sin(4.0*sin(3.0*t))
                  3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))*thetaxp^2 + 3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))*thetayp^2 + 2.09e-4*thetazpp - 1.11*sin(3.0*t)];

Into this Matlab function

 Y = [2.09e-4* yp(1) *cos(4.0*sin(3.0*t)) - 9.16e-4* yp(2) *sin(4.0*sin(3.0*t)) - 0.00125* y(2) *cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125* y(1) *cos(3.0*t)*sin(4.0*sin(3.0*t))
2.09e-4* yp(2) *cos(4.0*sin(3.0*t)) + 9.16e-4* yp(1) *sin(4.0*sin(3.0*t)) + 0.00125* y(1) *cos(3.0*t)*cos(4.0*sin(3.0*t)) - 0.00125* y(2) *cos(3.0*t)*sin(4.0*sin(3.0*t))
                  3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))* y(1)^2  + 3.54e-4*cos(4.0*sin(3.0*t))*sin(4.0*sin(3.0*t))* y(2)^2  + 2.09e-4* yp(3) - 1.11*sin(3.0*t)];

so it can be used with the function ode15i.

You are now following this question

VincentLec (view profile)