Error using mupadengin​e/feval_in​ternal

I'm trying to run this simple code in Matlab but I'm getting this error:

" Error using mupadengine/feval_internal

System contains a nonlinear equation in 'diff(y(x), x, x)'. The system must be quasi-linear: highest

derivatives must enter the differential equations linearly.

Error in odeToVectorField>mupadOdeToVectorField (line 171)

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

Error in odeToVectorField (line 119)

sol = mupadOdeToVectorField(varargin);

Error in Untitled (line 11)

[VF,Subs] = odeToVectorField(ode);"

Could anyone please help me with this problem? I'm new to Matlab.

My code is :

syms y(x) x Y 
r=0.05;
m=0.01;
p=1;
s=1;
a=0.25;
b=0.25;
Dy = diff(y);
D2y = diff(y,2);
ode= y-1/x*1/(r-((p*Dy*x*(s^2))/y)+((D2y*x*(s^2))/(2*y))+((Dy*(m-(s^2)))/y))^(1/(1-a-b));
[VF,Subs] = odeToVectorField(ode);

Error using mupadengine/feval_internal
System contains a nonlinear equation in 'diff(y(x), x, x)'. The system must be quasi-linear: highest derivatives must enter the differential equations linearly.

Error in odeToVectorField>mupadOdeToVectorField (line 171)
T = feval_internal(symengine,'symobj::odeToVectorField',sys,x,stringInput);

Error in odeToVectorField (line 119)
sol = mupadOdeToVectorField(varargin);

odefcn = matlabFunction(VF, 'Vars',{x,Y});

tspan = [20 80];                                     
ic = [1 0];                                      
[x,y] = ode45(odefcn, tspan, ic);

figure
plot(x, y)
grid
legend(string(Subs), 'loc','best')

Answer 1

[Walter Roberson]

syms y(x) x Y 
r=0.05;
m=0.01;
p=1;
s=1;
a=0.25;
b=0.25;
Dy = diff(y);
D2y = diff(y,2);
ode= y-1/x*1/(r-((p*Dy*x*(s^2))/y)+((D2y*x*(s^2))/(2*y))+((Dy*(m-(s^2)))/y))^(1/(1-a-b))

ode(x) = 

Notice that the denominator has the second derivative of y(x), and notice that the entire term is squared -- so the square of the second derivative of y(x) is used in the ODE.

The error is saying that whatever the highest degree of derivative is used, that derivative must be used linearly -- it can be multiplied by constants but it cannot be raised to any power (other than 0, which would remove the term, or 1, which would leave the term unchanged.)

[Fatemeh]

Thank you for your help. Is it still possible to solve this equation using another way?

[Walter Roberson]

At the moment, I do not know of any way to automatically process it into a function handle.

[Walter Roberson]

syms y(x) x Y 
r=0.05;
m=0.01;
p=1;
s=1;
a=0.25;
b=0.25;
Dy = diff(y);
D2y = diff(y,2);
ode= y-1/x*1/(r-((p*Dy*x*(s^2))/y)+((D2y*x*(s^2))/(2*y))+((Dy*(m-(s^2)))/y))^(1/(1-a-b))

ode(x) = 

syms dy d2y
temp = subs(ode, {D2y, Dy}, {d2y, dy})

temp(x) = 

rewritten = subs(isolate(temp==0, d2y), {d2y, dy}, {D2y, Dy})

rewritten = 

[eqs,vars] = reduceDifferentialOrder(rewritten,y(x))

eqs = 

vars = 

[M,F] = massMatrixForm(eqs,vars)

M = 

F = 

f = M\F

f = 

odefun = odeFunction(f,vars)

odefun = function_handle with value:
    @(x,in2)[in2(2,:);(in2(2,:).*9.9e+1-in2(1,:).*5.0+in2(2,:).*x.*1.0e+2+in2(1,:).*sqrt(1.0./(x.*in2(1,:))).*1.0e+2)./(x.*5.0e+1)]

tspan = [20 80];                                     
ic = [1 0];

[x, y] = ode15s(odefun, tspan, ic);
plot(x, y)
grid
legend(string(f), 'location','best')