How to remove "root(f(z), z, n)" from solve solution?

4 次查看(过去 30 天)
Wondering if it is possible to remove the expression "root(f(z), z, n)" from the solution presented by the solve function. In this case, my code is listed below and the outputs of "steadystate = solve(eq5,eq6,eq7,eq8, L,S,B,A)" include the expression "root(f(z), z, n)." I understand it to be the case that this stands in for the nth root of the polynomial f(z), but I'm wondering if it is possible to force a solution from MatLab which does not feature this expression.
clear
syms f B r L k c g t S u v h m w J a n p V D A mu q m1 m2 m3
assume(f,'positive')
assume(r,'positive')
assume(k,'positive')
assume(c,'positive')
assume(g,'positive')
assume(t,'positive')
assume(u,'positive')
assume(v,'positive')
assume(h,'positive')
assume(m,'positive')
assume(w,'positive')
assume(J,'positive')
assume(a,'positive')
assume(n,'positive')
assume(p,'positive')
assume(V,'positive')
assume(D,'positive')
assume(mu,'positive')
assume(q,'positive')
assume(A,'positive')
assume(m1,'positive')
assume(m2,'positive')
assume(m3,'positive')
eq1 = (f*B + r*L - k*L - c*L - m1*A*L)
eq2 = (g*B + t*S - u*S - v*S - m2*A*S)
eq3 = (h*S + m*L + w*B - a*B - n*B - p*B - m3*A*B)
eq4 = (q*(A+S+B)-mu*A)
eq5 = (eq1 == 0)
eq6 = (eq2 == 0)
eq7 = (eq3 == 0)
eq8 = (eq4 == 0)
steadystate = solve(eq5,eq6,eq7,eq8, L,S,B,A)
  2 个评论
Torsten
Torsten 2022-4-6
Maybe
steadystate = solve(eq5,eq6,eq7,eq8, L,S,B,A);
steadystate = vpa(steadystate)
Tristen Jackson
Tristen Jackson 2022-4-6
Unfortunately, using vpa(steadystate) is not effective, nor is applying this function to any of the components of 'steadystate' individually.

请先登录,再进行评论。

采纳的回答

Walter Roberson
Walter Roberson 2022-4-6
syms f B r L k c g t S u v h m w J a n p V D A mu q m1 m2 m3
assume(f,'positive')
assume(r,'positive')
assume(k,'positive')
assume(c,'positive')
assume(g,'positive')
assume(t,'positive')
assume(u,'positive')
assume(v,'positive')
assume(h,'positive')
assume(m,'positive')
assume(w,'positive')
assume(J,'positive')
assume(a,'positive')
assume(n,'positive')
assume(p,'positive')
assume(V,'positive')
assume(D,'positive')
assume(mu,'positive')
assume(q,'positive')
assume(A,'positive')
assume(m1,'positive')
assume(m2,'positive')
assume(m3,'positive')
eq1 = (f*B + r*L - k*L - c*L - m1*A*L)
eq1 = 
eq2 = (g*B + t*S - u*S - v*S - m2*A*S)
eq2 = 
eq3 = (h*S + m*L + w*B - a*B - n*B - p*B - m3*A*B)
eq3 = 
eq4 = (q*(A+S+B)-mu*A)
eq4 = 
eq5 = (eq1 == 0)
eq5 = 
eq6 = (eq2 == 0)
eq6 = 
eq7 = (eq3 == 0)
eq7 = 
eq8 = (eq4 == 0)
eq8 = 
steadystate = solve(eq5,eq6,eq7,eq8, [L,S,B,A], 'maxdegree', 3)
Warning: Solutions are only valid under certain conditions. To include parameters and conditions in the solution, specify the 'ReturnConditions' value as 'true'.
steadystate = struct with fields:
L: [3×1 sym] S: [3×1 sym] B: [3×1 sym] A: [3×1 sym]
steadystate.L
ans = 
  3 个评论
Torsten
Torsten 2022-4-6
Is there any way to "reconstruct" how MATLAB came to the final expressions for the variables ?
Walter Roberson
Walter Roberson 2022-4-6
syms A0 A1 A2 A3 x
sols = solve(A3*x^3 + A2*x^2 + A1*x^1 + A0*x^0, x, 'MaxDegree', 3);
string(sols(1))
ans = "(((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3) - (A1/(3*A3) - A2^2/(9*A3^2))/(((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3) - A2/(3*A3)"
string(sols(2))
ans = "(A1/(3*A3) - A2^2/(9*A3^2))/(2*(((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3)) - A2/(3*A3) - (3^(1/2)*((A1/(3*A3) - A2^2/(9*A3^2))/(((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3) + (((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3))*1i)/2 - (((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3)/2"
string(sols(3))
ans = "(3^(1/2)*((A1/(3*A3) - A2^2/(9*A3^2))/(((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3) + (((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3))*1i)/2 - A2/(3*A3) + (A1/(3*A3) - A2^2/(9*A3^2))/(2*(((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3)) - (((A0/(2*A3) + A2^3/(27*A3^3) - (A1*A2)/(6*A3^2))^2 + (A1/(3*A3) - A2^2/(9*A3^2))^3)^(1/2) - A0/(2*A3) - A2^3/(27*A3^3) + (A1*A2)/(6*A3^2))^(1/3)/2"
This shows the pattern; it is the classic solution for roots of a cubic. Everything beyond that is a matter of MATLAB having used coeffs() internally to figure out what A0, A1, A2, A3 match to in the original equation, then substituting.
Well, except that if you are working with polynomials with non-trivial coefficients, it can often be significantly more efficient to ask for the solution to the model polynomial and then subs() in coefficients.

请先登录,再进行评论。

更多回答(1 个)

Torsten
Torsten 2022-4-6
编辑:Torsten 2022-4-6
syms f B r L k c g t S u v h m w J a n p V D A mu q m1 m2 m3
eq1 = (f*B + r*L - k*L - c*L - m1*A*L)
eq2 = (g*B + t*S - u*S - v*S - m2*A*S)
eq3 = (h*S + m*L + w*B - a*B - n*B - p*B - m3*A*B)
eq4 = (q*(A+S+B)-mu*A)
eq5 = (eq1 == 0)
eq6 = (eq2 == 0)
eq7 = (eq3 == 0)
eq8 = (eq4 == 0)
vL = solve(eq5,L)
vS = solve(eq6,S)
vB = subs(eq8,S,vS)
vB = solve(vB,B)
vLL = subs(vL,B,vB)
vSS = subs(vS,B,vB)
vA = subs(eq3,[S,L,B],[vSS,vLL,vB])
simplify(vA)
If you run this code, you will see that you end up with product of poynomials in A of degree 1, 3 and 1.
So the maximum degree is 3 - you can solve analytically for A.
The expressions for B, S and L now follow easily by substituting the expressions for A in vLL, vSS and vB.

类别

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

标签

产品


版本

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by