# Hello, I have problem with odeToVectorField in ODE second order system of equations.

13 vues (au cours des 30 derniers jours)
Vahid le 18 Août 2023
% clear all;
g= 9.81;
m= 0.468;
Ix= 4.856e-3;
Iy= 4.856e-3;
Iz= 8.801e-3;
Ax= 0.25;
Ay= 0.25;
Az= 0.25;
l= 0.225;
k= 2.98e-6;
b= 1.14e-7;
Im= 3.357e-5;
f= 0; tu0= 0; tu1=0; tu2= 0;
Tuo= [tu0 tu1 tu2];
syms x0(t) x1(t) x2(t) t0(t) t1(t) t2(t) t
X= [x0(t) x1(t) x2(t)];
T= [t0(t) t1(t) t2(t)];
dX= diff(X,t);
dT= diff(T,t);
% d2X= diff(dX,t);
% d2T= diff(dT,t);
s1= sin(t0) ;s2=sin(t1) ;s3=sin(t2);
c1= cos(t0) ;c2=cos(t1) ;c3=cos(t2);
t0d= dT(1); t1d=dT(2) ; t2d=dT(3) ;
% C Matrix
c11= 0;
c12= (Iy-Iz)*(t1d*c1*s1+ t2d*s1*s1*c2)+(Iz-Iy)*(t2d*c1*c1*c2)- Ix*t2d*c2;
c13= (Iz-Iy)*t2d*c1*s1*c2*c2;
c21= (Iz-Iy)*(t1d*c1*s1+t2d*s1*c2)+(Iy-Iz)*t2d*c1*c1*c2+Ix*t2d*c2;
c22= (Iz-Iy)*t0d*c1*s1 ;
c23= -Ix*t2d*s2*c2+Iy*t2d*s1*s1*s2*c2+Iz*t2d*c1*c1*s2*c2;
c31= (Iy-Iz)*t2d*c2*c2*s1*c1-Ix*t1d*c2;
c32= (Iz-Iy)*(t1d*c1*s1*s2+t0d*s1*s1*c2)+(Iy-Iz)*t0d*c1*c1*c2+...
Ix*t2d*s2*c2-Iy*t2d*s1*s1*s2*c2-Iz*t2d*c1*c1*s2*c2;
c33= (Iy-Iz)*t0d*c1*s1*c2*c2-Iy*t1d*s1*s1*c2*s2-Iz*t1d*c1*c1*c2*s2+Ix*t1d*c2*s2;
C= [c11 c12 c13; c21 c22 c23; c31 c32 c33;];
% J matrix -invers
j11= Ix;
j12= 0;
j13= -Ix*s2;
j21= 0 ;
j22= Iy*c1*c1+Iz*s1*s1;
j23= (Iy-Iz)*c1*s1*c2;
j31= -Ix*s2;
j32= (Iy-Iz)*c1*s1*c2;
j33= Ix*s2*s2+Iy*s1*s1*c2*c2+Iz*c1*c1*c2*c2;
J = ([j11 j12 j13; j21 j22 j23; j31 j32 j33]);
DOM00= [c3*s2*c1+s3*s1, s3*s2*c1-c3*s1, c2*c1];
Ar= [Ax 0 0; 0 Ay 0; 0 0 Az];
% Equation of Motion
DOM01= J\(Tuo'- C*dT');
DOM02= [0 0 -g]' + (f/m)*DOM00'-(1.0/m)*Ar*dX';
Eq1 = diff(t0,2)== DOM01(1);
Eq2 = diff(t1,2)== DOM01(2);
Eq3 = diff(t2,2)== DOM01(3);
Eq4 = diff(x0,2)== DOM02(1);
Eq5 = diff(x1,2)== DOM02(2);
Eq6 = diff(x2,2)== DOM02(3);
vars = [t0(t); t1(t); t2(t); x0(t); x1(t); x2(t)];
[V,S] = odeToVectorField([Eq1 Eq2 Eq3 Eq4 Eq5 Eq6]);
Unable to 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.

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

Error in odeToVectorField (line 119)
M = matlabFunction(V,'vars', {'t'});
##### 3 commentairesAfficher 1 commentaire plus ancienMasquer 1 commentaire plus ancien
Vahid le 18 Août 2023
correct, it is not linear ...
Roshana Sandadini le 7 Nov 2023
Can I see the original Problem?

Connectez-vous pour commenter.

### Réponse acceptée

Torsten le 18 Août 2023
Modifié(e) : Torsten le 18 Août 2023
dydt = fun(0.1,ones(12,1)).'
dydt = 1×12
1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 3.7411 -0.6185 3.2759 -0.5342 -0.5342 -10.3442
function dydt = fun(t,y)
g= 9.81;
m= 0.468;
Ix= 4.856e-3;
Iy= 4.856e-3;
Iz= 8.801e-3;
Ax= 0.25;
Ay= 0.25;
Az= 0.25;
l= 0.225;
k= 2.98e-6;
b= 1.14e-7;
Im= 3.357e-5;
f= 0; tu0= 0; tu1=0; tu2= 0;
Tuo= [tu0 tu1 tu2];
t0 = y(1);
t1 = y(2);
t2 = y(3);
x0 = y(4);
x1 = y(5);
x2 = y(6);
t0d = y(7);
t1d = y(8);
t2d = y(9);
x0d = y(10);
x1d = y(11);
x2d = y(12);
X = [x0 x1 x2];
T = [t0 t1 t2];
dX = [x0d x1d x2d];
dT = [t0d t1d t2d];
s1= sin(t0) ;s2=sin(t1) ;s3=sin(t2);
c1= cos(t0) ;c2=cos(t1) ;c3=cos(t2);
% C Matrix
c11= 0;
c12= (Iy-Iz)*(t1d*c1*s1+ t2d*s1*s1*c2)+(Iz-Iy)*(t2d*c1*c1*c2)- Ix*t2d*c2;
c13= (Iz-Iy)*t2d*c1*s1*c2*c2;
c21= (Iz-Iy)*(t1d*c1*s1+t2d*s1*c2)+(Iy-Iz)*t2d*c1*c1*c2+Ix*t2d*c2;
c22= (Iz-Iy)*t0d*c1*s1 ;
c23= -Ix*t2d*s2*c2+Iy*t2d*s1*s1*s2*c2+Iz*t2d*c1*c1*s2*c2;
c31= (Iy-Iz)*t2d*c2*c2*s1*c1-Ix*t1d*c2;
c32= (Iz-Iy)*(t1d*c1*s1*s2+t0d*s1*s1*c2)+(Iy-Iz)*t0d*c1*c1*c2+...
Ix*t2d*s2*c2-Iy*t2d*s1*s1*s2*c2-Iz*t2d*c1*c1*s2*c2;
c33= (Iy-Iz)*t0d*c1*s1*c2*c2-Iy*t1d*s1*s1*c2*s2-Iz*t1d*c1*c1*c2*s2+Ix*t1d*c2*s2;
C= [c11 c12 c13; c21 c22 c23; c31 c32 c33;];
% J matrix -invers
j11= Ix;
j12= 0;
j13= -Ix*s2;
j21= 0 ;
j22= Iy*c1*c1+Iz*s1*s1;
j23= (Iy-Iz)*c1*s1*c2;
j31= -Ix*s2;
j32= (Iy-Iz)*c1*s1*c2;
j33= Ix*s2*s2+Iy*s1*s1*c2*c2+Iz*c1*c1*c2*c2;
J = ([j11 j12 j13; j21 j22 j23; j31 j32 j33]);
DOM00= [c3*s2*c1+s3*s1, s3*s2*c1-c3*s1, c2*c1];
Ar= [Ax 0 0; 0 Ay 0; 0 0 Az];
% Equation of Motion
DOM01= J\(Tuo.'- C*dT.');
DOM02= [0 0 -g].' + (f/m)*DOM00.'-(1.0/m)*Ar*dX.';
dydt = [y(7:12);DOM01;DOM02];
end
##### 0 commentairesAfficher -2 commentaires plus anciensMasquer -2 commentaires plus anciens

Connectez-vous pour commenter.

### Plus de réponses (2)

ProblemSolver le 18 Août 2023
The error you are getting indicates that the system of equations defined by DOM01 is non-linear. The odeToVectorField function in MATLAB requires the system of equations to be linear for it to work properly. But in your case, DOM01 contains trigonometric functions like sin() and cos() which makes it a nonlinear system.
Some options to resolve this:
1. Try linearizing the equations around an operating point using small angle approximations. This would eliminate the trig terms and result in a linear system that odeToVectorField can handle.
2. Use an alternate ODE solver in MATLAB that supports nonlinear systems, like ode45:
[T,X] = ode45(@ODEfunc,tspan,X0);
function dXdt = ODEfunc(t,X)
% Implement RHS of equations
dXdt = ...
end
3. Convert the equations to state-space form with nonlinear dynamics:
xdot = f(x,u)
y = g(x,u)
And use Simulink to simulate the state-space model instead of solving analytically with odeToVectorField.
##### 0 commentairesAfficher -2 commentaires plus anciensMasquer -2 commentaires plus anciens

Connectez-vous pour commenter.

Paul le 18 Août 2023
Modifié(e) : Paul le 19 Août 2023
The problem is confusion between sym expressions and symfun objects.
g= 9.81;
m= 0.468;
Ix= 4.856e-3;
Iy= 4.856e-3;
Iz= 8.801e-3;
Ax= 0.25;
Ay= 0.25;
Az= 0.25;
l= 0.225;
k= 2.98e-6;
b= 1.14e-7;
Im= 3.357e-5;
f= 0; tu0= 0; tu1=0; tu2= 0;
Tuo= [tu0 tu1 tu2];
syms x0(t) x1(t) x2(t) t0(t) t1(t) t2(t) t
X= [x0(t) x1(t) x2(t)];
T= [t0(t) t1(t) t2(t)];
dX= diff(X,t);
dT= diff(T,t);
% d2X= diff(dX,t);
% d2T= diff(dT,t);
s1= sin(t0) ;s2=sin(t1) ;s3=sin(t2);
c1= cos(t0) ;c2=cos(t1) ;c3=cos(t2);
t0d= dT(1); t1d=dT(2) ; t2d=dT(3) ;
% C Matrix
c11= 0;
c12= (Iy-Iz)*(t1d*c1*s1+ t2d*s1*s1*c2)+(Iz-Iy)*(t2d*c1*c1*c2)- Ix*t2d*c2;
c13= (Iz-Iy)*t2d*c1*s1*c2*c2;
c21= (Iz-Iy)*(t1d*c1*s1+t2d*s1*c2)+(Iy-Iz)*t2d*c1*c1*c2+Ix*t2d*c2;
c22= (Iz-Iy)*t0d*c1*s1 ;
c23= -Ix*t2d*s2*c2+Iy*t2d*s1*s1*s2*c2+Iz*t2d*c1*c1*s2*c2;
c31= (Iy-Iz)*t2d*c2*c2*s1*c1-Ix*t1d*c2;
c32= (Iz-Iy)*(t1d*c1*s1*s2+t0d*s1*s1*c2)+(Iy-Iz)*t0d*c1*c1*c2+...
Ix*t2d*s2*c2-Iy*t2d*s1*s1*s2*c2-Iz*t2d*c1*c1*s2*c2;
c33= (Iy-Iz)*t0d*c1*s1*c2*c2-Iy*t1d*s1*s1*c2*s2-Iz*t1d*c1*c1*c2*s2+Ix*t1d*c2*s2;
C= [c11 c12 c13; c21 c22 c23; c31 c32 c33;];
% J matrix -invers
j11= Ix;
j12= 0;
j13= -Ix*s2;
j21= 0 ;
j22= Iy*c1*c1+Iz*s1*s1;
j23= (Iy-Iz)*c1*s1*c2;
j31= -Ix*s2;
j32= (Iy-Iz)*c1*s1*c2;
j33= Ix*s2*s2+Iy*s1*s1*c2*c2+Iz*c1*c1*c2*c2;
J = ([j11 j12 j13; j21 j22 j23; j31 j32 j33]);
DOM00= [c3*s2*c1+s3*s1, s3*s2*c1-c3*s1, c2*c1];
Ar= [Ax 0 0; 0 Ay 0; 0 0 Az];
% Equation of Motion
Using ctranspose here will result in a lot of unnecessary conj() calls later on. Use transpose.
%DOM01= J\(Tuo'- C*dT');
%DOM02= [0 0 -g]' + (f/m)*DOM00'-(1.0/m)*Ar*dX';
DOM01= J\(Tuo.'- C*dT.');
DOM02= [0 0 -g].' + (f/m)*DOM00.'-(1.0/m)*Ar*dX.';
At this point, DOM01 is scalar symfun that takes in a scalar value and returns a 3 x 1 sym expression
whos DOM01
Name Size Bytes Class Attributes DOM01 1x1 8 symfun
temp = DOM01(0);
whos temp
Name Size Bytes Class Attributes temp 3x1 8 sym
So the following assignment is evaluating DOM01 at t = 1 and returning the 3 x 1 result and setting that equal to the 3x1 expansion of the lhs. The rhs is too painful to look at, but if you do you'll see various terms being evaluated at t = 1.
Eq1 = diff(t0,2)== DOM01(1);
lhs(Eq1)
ans(t) =
I suspect that you really want to pick off the elements of DOM01, like so
Eq1 = diff(t0,2) == [1 0 0]*DOM01(t)
Eq1(t) =
You could do that for all of the equations. Or just do it all at once:
%{
Eq1 = diff(t0,2)== DOM01(1);
Eq2 = diff(t1,2)== DOM01(2);
Eq3 = diff(t2,2)== DOM01(3);
Eq4 = diff(x0,2)== DOM02(1);
Eq5 = diff(x1,2)== DOM02(2);
Eq6 = diff(x2,2)== DOM02(3);
vars = [t0(t); t1(t); t2(t); x0(t); x1(t); x2(t)];
[V,S] = odeToVectorField([Eq1 Eq2 Eq3 Eq4 Eq5 Eq6]);
%}
Eq = [diff(t0,2);diff(t1,2);diff(t2,2);diff(x0,2);diff(x1,2);diff(x2,2)] == [DOM01(t);DOM02(t)];
[V,S] = odeToVectorField(Eq); % no errors
odefun = matlabFunction(V,'vars',{'t','Y'})
odefun = function_handle with value:
@(t,Y)[Y(2);((cos(Y(3)).*sin(Y(3)).*Y(2).*Y(4).*-2.546391437701197e+31-cos(Y(1)).*cos(Y(3)).^2.*sin(Y(1)).*Y(6).^2.*5.680808882942518e+31+cos(Y(1)).*cos(Y(3)).^4.*sin(Y(1)).*Y(6).^2.*1.029588117355377e+32-cos(Y(1)).*sin(Y(1)).*sin(Y(3)).^2.*Y(6).^2.*3.134417445241321e+31+cos(Y(1)).*sin(Y(1)).*sin(Y(3)).^4.*Y(6).^2.*3.134417445241321e+31+cos(Y(1)).*cos(Y(3)).^2.*Y(4).*Y(6).*5.680808882942518e+31-cos(Y(1)).*cos(Y(3)).^4.*Y(4).*Y(6).*4.615072290611251e+31+cos(Y(3)).*sin(Y(3)).^3.*Y(2).*Y(4).*7.161463728312448e+31+cos(Y(3)).^3.*sin(Y(3)).*Y(2).*Y(4).*7.161463728312448e+31+cos(Y(1)).*sin(Y(3)).^2.*Y(4).*Y(6).*3.134417445241321e+31+cos(Y(1)).*sin(Y(3)).^3.*Y(4).*Y(6).*2.546391437701197e+31+cos(Y(3)).*sin(Y(1)).*sin(Y(3)).*Y(2).*Y(6).*2.546391437701197e+31+cos(Y(1)).*cos(Y(3)).^2.*sin(Y(3)).*Y(4).*Y(6).*4.615072290611251e+31-cos(Y(3)).*sin(Y(1)).*sin(Y(3)).^3.*Y(2).*Y(6).*7.161463728312448e+31-cos(Y(3)).^3.*sin(Y(1)).*sin(Y(3)).*Y(2).*Y(6).*7.161463728312448e+31+cos(Y(1)).*cos(Y(3)).^2.*sin(Y(1)).*sin(Y(3)).^2.*Y(6).^2.*1.343029861879509e+32-cos(Y(1)).*cos(Y(3)).^2.*sin(Y(3)).^2.*Y(4).*Y(6).*6.683753143521304e+31).*(-1.760312695965974e-32))./(cos(Y(3)).^2.*sin(Y(3)).^2.*2.0+cos(Y(3)).^4+sin(Y(3)).^4);Y(4);(cos(Y(1)).*cos(Y(3)).*sin(Y(3)).^5.*Y(2).^2.*4.615072290611251e+31+cos(Y(1)).*cos(Y(3)).^5.*sin(Y(3)).*Y(2).^2.*4.615072290611251e+31+cos(Y(3)).^2.*sin(Y(1)).*Y(2).*Y(4).*3.134417445241321e+31+cos(Y(3)).^4.*sin(Y(1)).*Y(2).*Y(4).*2.546391437701197e+31+sin(Y(1)).*sin(Y(3)).^2.*Y(2).*Y(4).*5.680808882942518e+31-sin(Y(1)).*sin(Y(3)).^4.*Y(2).*Y(4).*4.615072290611251e+31+cos(Y(1)).*cos(Y(3)).^3.*sin(Y(3)).^3.*Y(2).^2.*9.230144581222501e+31-cos(Y(1)).^3.*cos(Y(3)).*sin(Y(3)).^5.*Y(6).^2.*4.615072290611251e+31-cos(Y(1)).^3.*cos(Y(3)).^5.*sin(Y(3)).*Y(6).^2.*4.615072290611251e+31+cos(Y(1)).^2.*cos(Y(3)).^4.*Y(2).*Y(6).*5.680808882942518e+31-cos(Y(1)).^2.*cos(Y(3)).^6.*Y(2).*Y(6).*4.615072290611251e+31-cos(Y(3)).^2.*sin(Y(1)).^2.*Y(2).*Y(6).*3.134417445241321e+31+cos(Y(1)).^2.*sin(Y(3)).^4.*Y(2).*Y(6).*5.680808882942518e+31+cos(Y(3)).^4.*sin(Y(1)).^2.*Y(2).*Y(6).*8.81522632818384e+31+cos(Y(1)).^2.*sin(Y(3)).^6.*Y(2).*Y(6).*4.615072290611251e+31-sin(Y(1)).^2.*sin(Y(3)).^2.*Y(2).*Y(6).*5.680808882942518e+31+sin(Y(1)).^2.*sin(Y(3)).^4.*Y(2).*Y(6).*1.597669005649629e+32-cos(Y(1)).^3.*cos(Y(3)).^3.*sin(Y(3)).^3.*Y(6).^2.*9.230144581222501e+31-cos(Y(3)).^2.*sin(Y(1)).*sin(Y(3)).^2.*Y(2).*Y(4).*2.068680852910054e+31+cos(Y(1)).^2.*cos(Y(3)).^2.*sin(Y(3)).^2.*Y(2).*Y(6).*1.136161776588504e+32+cos(Y(1)).^2.*cos(Y(3)).^2.*sin(Y(3)).^4.*Y(2).*Y(6).*4.615072290611251e+31-cos(Y(1)).^2.*cos(Y(3)).^4.*sin(Y(3)).^2.*Y(2).*Y(6).*4.615072290611251e+31+cos(Y(3)).^2.*sin(Y(1)).^2.*sin(Y(3)).^2.*Y(2).*Y(6).*2.479191638468013e+32+cos(Y(1)).*cos(Y(3)).*sin(Y(1)).^2.*sin(Y(3)).*Y(6).^2.*2.546391437701197e+31-cos(Y(1)).*cos(Y(3)).*sin(Y(1)).^2.*sin(Y(3)).^3.*Y(6).^2.*7.161463728312448e+31-cos(Y(1)).*cos(Y(3)).^3.*sin(Y(1)).^2.*sin(Y(3)).*Y(6).^2.*7.161463728312448e+31-cos(Y(1)).*cos(Y(3)).*sin(Y(1)).*sin(Y(3)).*Y(4).*Y(6).*2.546391437701197e+31-cos(Y(1)).*cos(Y(3)).*sin(Y(1)).*sin(Y(3)).^2.*Y(4).*Y(6).*2.068680852910054e+31+cos(Y(1)).*cos(Y(3)).*sin(Y(1)).*sin(Y(3)).^3.*Y(4).*Y(6).*9.230144581222501e+31+cos(Y(1)).*cos(Y(3)).^3.*sin(Y(1)).*sin(Y(3)).*Y(4).*Y(6).*7.161463728312448e+31)./(cos(Y(1)).*(cos(Y(3)).^2.*sin(Y(3)).^2.*2.0+cos(Y(3)).^4+sin(Y(3)).^4).*5.680808882942518e+31);Y(6);(cos(Y(3)).^2.*Y(2).*Y(4).*3.134417445241321e+31+cos(Y(3)).^4.*Y(2).*Y(4).*2.546391437701197e+31+sin(Y(3)).^2.*Y(2).*Y(4).*5.680808882942518e+31-sin(Y(3)).^4.*Y(2).*Y(4).*4.615072290611251e+31-cos(Y(3)).^2.*sin(Y(1)).*Y(2).*Y(6).*3.134417445241321e+31+cos(Y(3)).^4.*sin(Y(1)).*Y(2).*Y(6).*8.81522632818384e+31-sin(Y(1)).*sin(Y(3)).^2.*Y(2).*Y(6).*5.680808882942518e+31+sin(Y(1)).*sin(Y(3)).^4.*Y(2).*Y(6).*1.597669005649629e+32-cos(Y(3)).^2.*sin(Y(3)).^2.*Y(2).*Y(4).*2.068680852910054e+31+cos(Y(3)).^2.*sin(Y(1)).*sin(Y(3)).^2.*Y(2).*Y(6).*2.479191638468013e+32+cos(Y(1)).*cos(Y(3)).*sin(Y(1)).*sin(Y(3)).*Y(6).^2.*2.546391437701197e+31-cos(Y(1)).*cos(Y(3)).*sin(Y(3)).*Y(4).*Y(6).*2.546391437701197e+31-cos(Y(1)).*cos(Y(3)).*sin(Y(1)).*sin(Y(3)).^3.*Y(6).^2.*7.161463728312448e+31-cos(Y(1)).*cos(Y(3)).^3.*sin(Y(1)).*sin(Y(3)).*Y(6).^2.*7.161463728312448e+31-cos(Y(1)).*cos(Y(3)).*sin(Y(3)).^2.*Y(4).*Y(6).*2.068680852910054e+31+cos(Y(1)).*cos(Y(3)).*sin(Y(3)).^3.*Y(4).*Y(6).*9.230144581222501e+31+cos(Y(1)).*cos(Y(3)).^3.*sin(Y(3)).*Y(4).*Y(6).*7.161463728312448e+31)./(cos(Y(1)).*(cos(Y(3)).^2.*sin(Y(3)).^2.*2.0+cos(Y(3)).^4+sin(Y(3)).^4).*5.680808882942518e+31);Y(8);Y(8).*(-1.25e+2./2.34e+2);Y(10);Y(10).*(-1.25e+2./2.34e+2);Y(12);Y(12).*(-1.25e+2./2.34e+2)-9.81e+2./1.0e+2]
Also, you may want consider defining all those constants as sym constants up front, like
g = sym(9.81);
In general, it's best (IMO) to stay in sym world as long as possible to allow the engine to find cancellations and simplifications, should there be any.
##### 6 commentairesAfficher 4 commentaires plus anciensMasquer 4 commentaires plus anciens
Paul le 19 Août 2023
To your point, is it possible that using symbolic constants might not always work. That is if we do this
syms a b
temp = a*b;
c = temp/a
c =
b
the 'a's will always cancel.
If a and b are assigned symbolic constants
a = sym(1.23456);
b = sym(5.43210);
temp = a*b
temp =
c = temp/a
c =
then in this case c is returned with the value of b. Will this result always obtain because the symbolic constants are always rational numbers and so cancelling the common factors in computations is straightforward? Or can there be a case, say where c is result of a sequence of expressions involving who knows what kind of operations, where something gets lost in the shuffle, so to speak, and the result for c would be not quite what would be expected?
Walter Roberson le 19 Août 2023

sym() of a floating point number can return multiples of π and square roots when the default conversion is used. The default conversion also approximates... so for example a floating point number sufficiently close to 1/3 will be converted to 1/3, not just the single double precision number that most closely represents 1/3. π ± about 999 eps is converted to π

But once you have converted to rational or sqrt or pi, the symbolic engine sticks with the representation and calculates with it until you have a double() or vpa() or vpaintegral(), or you have a relational operation in a context that forces a decision to be made. "if A<B" can be difficult to answer when A and B are symbolic formula. The engine cannot rely on being able to come up with mathematical proofs so sometimes it has to resort to numeric evaluation. But it holds off until then.

For example sum(sin(sym(linspace(0,pi))) will give a mostly unevaluated sum of hundreds of sin terms: it does not approximate the sin of rationals.

Now if you construct a symbolic floating-point number then the engine will create approximations as it goes. for example sym('1.23') or str2sym('1.23') or sym() of a floating point number with particular options to the call. Once you are working with symbolic floating-point including if you vpa() then you can run into the same kind of issues as you can with double precision

Connectez-vous pour commenter.

### Catégories

En savoir plus sur Numbers and Precision dans Help Center et File Exchange

### Community Treasure Hunt

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

Start Hunting!

Translated by