Error using mupadmex Error in MuPAD command: symbolic:s​ym:isAlway​s:LiteralC​ompare|0 < root(z^20 - (555455*z^19)/28224 + (21279533*z^18)/112896 - (524070713*z^17)/451584 + (331

2 次查看(过去 30 天)
syms r;
eq = (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^4 * sqrt(((1/2) + 1)^2 - r^2)^4 * sqrt(((1/2) - 1)^2 - r^2)^4 ...
- 2 * (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^2 * (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^2 ...
* sqrt(((1/2) - 1)^2 - r^2)^2 * sqrt(((1/2) + 1)^2 - r^2)^4 * r^4 ...
- 2 * (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^2 * (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^2 ...
* sqrt(((1/2) - 1)^2 - r^2)^4 * sqrt(((1/2) + 1)^2 - r^2)^2 * r^4 ...
+ (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^4 * sqrt(((1/2) + 1)^2 - r^2)^4 * r^8 ...
- 2 * (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^2 * (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^2 ...
* sqrt(((1/2) + 1)^2 - r^2)^2 * sqrt(((1/2) - 1)^2 - r^2)^2 * r^8 ...
+ (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^4 * sqrt(((1/2) - 1)^2 - r^2)^4 * r^8;
sol = solve(eq, r);
positive_solutions = sol(sol > 0);
Error using symengine
Unable to prove '0 < root(z^20 - (555455*z^19)/28224 + (21279533*z^18)/112896 - (524070713*z^17)/451584 + (331982393245*z^16)/65028096 - (1096679671049*z^15)/65028096 + (22231132110163*z^14)/520224768 - (86269923640667*z^13)/1040449536 + (659908638483491*z^12)/5549064192 - (614010908943583*z^11)/5549064192 + (1542961696333139*z^10)/66588770304 + (3684172462923767*z^9)/29595009024 - (86167575482678911*z^8)/355140108288 + (125425715847634333*z^7)/532710162432 - (428976223694222489*z^6)/4261681299456 - (115352974498514309*z^5)/2130840649728 + (2032084483563487*z^4)/16911433728 - (1527197305978447*z^3)/16911433728 + (778786651869991*z^2)/21743271936 - (134768480539*z)/18874368 + 2325457729/4194304, z, 1)^(1/2)' literally. Use 'isAlways' to test the statement mathematically.

Error in indexing (line 920)
X = find(mupadmex('symobj::logical',A.s,9)) - 1;

Error in indexing (line 968)
R_tilde = builtin('subsref',L_tilde,Idx);
disp(['q: ', num2str(length(positive_solutions))]);

采纳的回答

Dyuman Joshi
Dyuman Joshi 2024-1-2
编辑:Dyuman Joshi 2024-1-2
solve() is unable to find an explicit solution to the given equation, which is expected given the degree of the polynomial - on simplifying the equation, you can observe that the highest power is 40.
There is no known method to find explicit solution for a polynomial of that degree.
You can use vpasolve -
Also, you should update the condition to check for non-complex roots as well.
syms r;
eq = (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^4 * sqrt(((1/2) + 1)^2 - r^2)^4 * sqrt(((1/2) - 1)^2 - r^2)^4 ...
- 2 * (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^2 * (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^2 ...
* sqrt(((1/2) - 1)^2 - r^2)^2 * sqrt(((1/2) + 1)^2 - r^2)^4 * r^4 ...
- 2 * (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^2 * (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^2 ...
* sqrt(((1/2) - 1)^2 - r^2)^4 * sqrt(((1/2) + 1)^2 - r^2)^2 * r^4 ...
+ (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^4 * sqrt(((1/2) + 1)^2 - r^2)^4 * r^8 ...
- 2 * (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^2 * (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^2 ...
* sqrt(((1/2) + 1)^2 - r^2)^2 * sqrt(((1/2) - 1)^2 - r^2)^2 * r^8 ...
+ (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^4 * sqrt(((1/2) - 1)^2 - r^2)^4 * r^8;
disp(simplify(eq))
sol = vpasolve(eq, r)
sol = 
%Update the condition
positive_solutions = sol(real(sol) > 0 & imag(sol) == 0 );
disp(['q: ', num2str(length(positive_solutions))]);
q: 2

更多回答(1 个)

Walter Roberson
Walter Roberson 2024-1-2
Your equality is equivalent to a polynomial of degree 20.
MATLAB is not able to find a closed form solution for the roots (which is not surprising -- theory says that degree 4 is the largest degree that can routinely be solved algebraically.) So MATLAB returns a form that "stands in" for the polynomial roots.
You then ask to reduce the solutions to the positive numbers. That requires being able to prove that particular solutions are positive... which is difficult when you cannot express the solutions explicitly.
  2 个评论
Walter Roberson
Walter Roberson 2024-1-2
syms r;
eq = (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^4 * sqrt(((1/2) + 1)^2 - r^2)^4 * sqrt(((1/2) - 1)^2 - r^2)^4 ...
- 2 * (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^2 * (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^2 ...
* sqrt(((1/2) - 1)^2 - r^2)^2 * sqrt(((1/2) + 1)^2 - r^2)^4 * r^4 ...
- 2 * (3040*r^6 - 12432*r^4 + 18544*r^2 - 18592)^2 * (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^2 ...
* sqrt(((1/2) - 1)^2 - r^2)^4 * sqrt(((1/2) + 1)^2 - r^2)^2 * r^4 ...
+ (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^4 * sqrt(((1/2) + 1)^2 - r^2)^4 * r^8 ...
- 2 * (768*r^8 - 768*r^6 + 288*r^4 - 304*r^2 + 69)^2 * (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^2 ...
* sqrt(((1/2) + 1)^2 - r^2)^2 * sqrt(((1/2) - 1)^2 - r^2)^2 * r^8 ...
+ (768*r^8 - 6912*r^6 + 23328*r^4 - 35248*r^2 + 33081)^4 * sqrt(((1/2) - 1)^2 - r^2)^4 * r^8;
sol = solve(eq, r);
dsol = double(sol);
positive_solutions = sol(imag(dsol)==0 & real(dsol) > 0)
positive_solutions = 
disp(['q: ', num2str(length(positive_solutions))]);
q: 2
disp(double(positive_solutions))
0.5000 0.5000

请先登录,再进行评论。

Community Treasure Hunt

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

Start Hunting!

Translated by