How to define this variable on MATLAB

8 views (last 30 days)
Yakoob on 19 Dec 2022
Commented: Paul on 12 Jan 2023
I have this positive equilibrium (x,l,y,v,M) where
x = (q/b)*(eta+(eps*M)/(1+M))*(1+mu+((rho*M)/(1+M)))
l = (q/(alp*sig))*(zeta*M-omega0)
y = (1/sig)*(zeta*M-omega0)
v = (b/sig)*(1/(eta+(eps*M)/(1+M)))*(zeta*M-omega0)
M is a positive root of the following equation
Does anyone have an idea of how to enter these informaion on MATLAB, espically the last three lines?

Answers (2)

Suvansh Arora
Suvansh Arora on 22 Dec 2022
In order to solve this equation with 5 unknows, follow the MATLAB answers article mentioned below:
Steven Lord
Steven Lord on 22 Dec 2022
If you want to assign a sym object into an array and operate symbolically, preallocate the array as a sym array rather than as a double. If you then wanted to pass that into a function that only accepts double arrays you may need to convert at that point (which may involve using subs to substitute values for the symbolic variables.)
A = zeros(3, 4, 'sym')
A = 
syms x
A(2, 2) = x^2 % this works
A = 
Compare with:
B = zeros(3, 4)
B = 3×4
0 0 0 0 0 0 0 0 0 0 0 0
B(2, 2) = x^2 % this throws an error
Unable to perform assignment because value of type 'sym' is not convertible to 'double'.

Caused by:
Error using symengine
Unable to convert expression containing symbolic variables into double array. Apply 'subs' function first to substitute values for variables.

Sign in to comment.

Paul on 22 Dec 2022
We can use the Symbolic Math Toolbox as one option
Set up the equations
syms M positive
x = (q/b)*(eta+(eps*M)/(1+M))*(1+mu+((rho*M)/(1+M)));
l = (q/(alp*sig))*(zeta*M-omega0);
y = (1/sig)*(zeta*M-omega0);
v = (b/sig)*(1/(eta+(eps*M)/(1+M)))*(zeta*M-omega0);
eqn = r0*(1+k1*(mu+(rho*M)/(1+M))*(1-k2*(mu+(rho*M)/(1+M))))*(eta+(eps*M)/(1+M))*(K-alp*sig*q*(eta+(eps*M)/(1+M))*(1+mu+(rho*M)/(1+M))-b*(q+alp)*(zeta*M-omega0))-b^2*alp*K*(zeta*M-omega0)==0
eqn = 
eqn = simplify(lhs(eqn)) == 0
eqn = 
Try solve
Msol1 = solve(eqn,M)
Msol1 = 
Use vpa to get the root of the polynomial
ans = 
Verify the solution
ans = 
Alternatively, we can extract the numerator of the LHS and use roots, there is only ony positive solution:
Msol2 = roots(sym2poly(numden(lhs(eqn))))
Msol2 = 6×1
2.1266 -1.2184 -1.0161 -0.9851 -0.7685 -0.2448
If the equation to solve was more complicated then just a sixth order polynmial, we could use vpasolve
Msol3 = vpasolve(eqn,M,[0 inf])
Msol3 = 
Or if desired to just get a numerical solutiion, use fsolve. The function to solve could have been defined directly without going through the symbolic stuff, but since we already have the symbolic stuff, just use it
fun = matlabFunction(lhs(eqn));
Msol4 = fsolve(fun,1)
Equation solved. fsolve completed because the vector of function values is near zero as measured by the value of the function tolerance, and the problem appears regular as measured by the gradient.
Msol4 = 2.1285
ans = -1.2185e-05
Paul on 12 Jan 2023
Unless solve can find an exact answer, any other solution, whether in double precision or variable precision arithmetic, is going to be an approximation. It's just a matter of how accurate that approximation is and what you're willing to accept. I guess I'm saying the same thing as Steven.

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!

Translated by