26 views (last 30 days)
Show older comments
Álvaro about 8 hours ago
-
-
Link
Direct link to this question
https://au.mathworks.com/matlabcentral/answers/2122626-i-have-a-nonlinear-equation-with-a-symbolic-variable-and-cant-solve-it
Edited: John D'Errico about 4 hours ago
Accepted Answer: John D'Errico
Open in MATLAB Online
Hi, I have to plot a phase diagram of an implicit equation, but in there lies a function p that depends nonlinearly on one of the variables (alpha), so I am trying to solve the nonlinear equation in terms of my variable alpha. However, I only get a solution in terms of a z2 variable. Can anybody help? Here is the code
n = 2.1;
alpha0 = 0;
syms p alpha
eqn = p^(n+1) - alpha0*p^n + p -(alpha0 + alpha) == 0;
sol = solve(eqn,p);
Warning: Solutions are parameterized by the symbols: z2. To include parameters and conditions in the solution, specify the 'ReturnConditions' value as 'true'.
Warning: Solutions are only valid under certain conditions. To include parameters and conditions in the solution, specify the 'ReturnConditions' value as 'true'.
For n = 2 works just fine, but not for n = 2.1, which prompts my problem.
0 Comments Show -2 older commentsHide -2 older comments
Show -2 older commentsHide -2 older comments
Sign in to comment.
Sign in to answer this question.
Accepted Answer
John D'Errico about 6 hours ago
⋮
Edited: John D'Errico about 4 hours ago
Open in MATLAB Online
Look carefully at your problem.
n = 2.1;
alpha0 = 0;
syms p alpha
eqn = p^(n+1) - alpha0*p^n + p -(alpha0 + alpha) == 0
eqn=
So you have sort of, almost kind of, a cubic polynomial, but not really. Raising numbers to non-integrer powers causes problems. There is NO analytical solution. Effectively, we can think of this as a very high order polynomial problem, of degree 31 in fact.) And polynomials higher than order 4 have no solution in algebraic form. This has been known for well over 100 years.
Is there something you could do? If you have the value of alpha, then you could simple use fzero. For example, you might do this:
Psol = @(alpha) fzero(@(p) p^(n+1) - alpha0*p^n + p -(alpha0 + alpha),1)
Psol = function_handle with value:
@(alpha)fzero(@(p)p^(n+1)-alpha0*p^n+p-(alpha0+alpha),1)
Psol(3)
ans = 1.2072
So plug in any value of alpha, and fzero will (hopefully) return a solution. You could also use vpasolve, which will survive as long as alpha is given in numerical form. Or fsolve, etc. Another way of looking at this is in the form of a plot. (If we recognize that whenever p is less than zero, all sorts of hell will break lose in terms of complex variables, then we can see that p should be positive. So write the problem by isolating alpha from the rest.
alpha = p + p^3.1
Now we can plot alpha, as a function of p.
fplot(@(p) p + p.^(3.1),[0,3])
xlabel p
ylabel alpha
grid on
Given any value of alpha on the vertical axis, you wish to solve for p. Again, there is NO algebraic solution to this. Yes, if you wanted, you could approximate the inverse relationship. What you want to see is that same relation, but with the axes flipped. And you could easily enough do that. For example...
pint = linspace(0,10);
alphaint = pint + pint.^3.1;
plot(alphaint,pint,'-')
xlabel alpha
ylabel p
grid on
You could now easily enough fit those points with a spline. For example...
pspl = spline(alphaint,pint);
fnval(pspl,5)
ans = 1.4982
So when alpha is 5, p would now directly be predicted as 1.4982. Of course, this is just an interpolation. And again, there is no algebraic function you can write down. Sorry.
3 Comments Show 1 older commentHide 1 older comment
Show 1 older commentHide 1 older comment
Álvaro about 4 hours ago
Direct link to this comment
https://au.mathworks.com/matlabcentral/answers/2122626-i-have-a-nonlinear-equation-with-a-symbolic-variable-and-cant-solve-it#comment_3171911
Open in MATLAB Online
The span of alpha for the figure I want to plot is [1 10^4], so there exists the posibility of running the equation for all possible values of alpha, let's say for x = linspace(1,10^4,10^5). However, when running this in a "for" loop, it almost instantaneously goes crazy as the initial condition is probably off (I think, though it says it is due to encountering complex solutions). The code is like this
n = 2.1;
alpha0 = 0;
x = linspace(1,10^4,10^5);
A = zeros(10^5,1);
Psol = @(alpha) fzero(@(p) p^(n+1) - alpha0*p^n + p -(alpha0 + alpha),1);
for i = 1:length(x)
A(i,1) = Psol(i);
end
Incorrect number or types of inputs or outputs for function solve.
Error in solution>@(alpha)solve(@(p)p^(n+1)-alpha0*p^n+p-(alpha0+alpha),1) (line 5)
Psol = @(alpha) solve(@(p) p^(n+1) - alpha0*p^n + p -(alpha0 + alpha),1);
Álvaro about 5 hours ago
Direct link to this comment
https://au.mathworks.com/matlabcentral/answers/2122626-i-have-a-nonlinear-equation-with-a-symbolic-variable-and-cant-solve-it#comment_3171926
Open in MATLAB Online
I have also tried with the vpasolve, which has finally worked with the following code
n = 2.1;
alpha0 = 0;
x = linspace(1,10^4,10^5);
A = zeros(10^5,1);
syms p
for i = 1:length(x)
eq = p^(n+1) - alpha0*p^n + p -(alpha0 + i) == 0;
A(i,1) = vpasolve(eq,p);
end
My computer has been working for 30 minutes but I finally have an array for all the values of p from alpha = 1 to alpha = 10^4. I do not know whether the values are correct or how to implement this result onto a fimplicit equation, but I will manage. Thanks a lot!
John D'Errico about 2 hours ago
Direct link to this comment
https://au.mathworks.com/matlabcentral/answers/2122626-i-have-a-nonlinear-equation-with-a-symbolic-variable-and-cant-solve-it#comment_3171946
Edited: John D'Errico about 2 hours ago
Open in MATLAB Online
Honestly, you are going the wrong way. First, how far out does p go, if alpha is 1e4?
n = 2.1;
alpha0 = 0;
x = linspace(1,10^4,10^5);
A = zeros(10^5,1);
Note that I used fzero slightly different here, giving it a bracket. As long as the upper end is high enough, it will have no problem.
Psol = @(alpha) fzero(@(p) p^(n+1) - alpha0*p^n + p -(alpha0 + alpha),[0,100]);
Psol(1e4)
ans = 19.5007
Ok, so we want to go as high as 19.5007, or so.
pint = linspace(0,20,10000); % a nice fine spacing
alphaint = pint.^(n+1) - alpha0*pint.^n + pint -alpha0;
Now we have a set of values that relate alpha and p.
plot(alphaint,pint,'.')
xlabel alpha
ylabel p
grid on
And the nice thing is, this was virtually instantaneous. Now you can use any interpolant to recover p, as a function of alpha.
Perhaps even better yet, see this relation is fairly well-behaved if we plot it in loglog form
loglog(alphaint,pint,'-')
xlabel alpha
ylabel p
grid on
And that suggests you have a simple scheme that will work very nicely. Interpolate log(p), as a function of log(alpha). Even a linear interpolant will sufficent there. Then exponentiate the result. something like this:
ppred = @(alpha) exp(interp1(log(alphaint(2:end)),log(pint(2:end)),log(alpha)));
ppred(3000)
ans = 13.2140
Note that I dropped the first point in those vectors, since it will be zero, and log(0) is a problem.
Sign in to comment.
More Answers (0)
Sign in to answer this question.
See Also
Tags
- nonlinear
- nonlinear equation
- symbolic variable
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
An Error Occurred
Unable to complete the action because of changes made to the page. Reload the page to see its updated state.
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
You can also select a web site from the following list
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- Deutsch
- English
- Français
- United Kingdom(English)
Asia Pacific
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)
Contact your local office
