Slove function return empty solutions
古いコメントを表示
Hello, I'm trying to solve the attached syntax, but the aolve function return empty solutions. Please help.
syms V_1 V_2 x_1 x_2 r
pi1 = (V_1) * (x_1^r/(x_1^r+x_2^r)) - x_1
pi2 = (V_2) * (x_2^r/(x_1^r+x_2^r)) - x_2
dpi1dx = diff(pi1, x_1)
dpi2dx = diff(pi2, x_2)
s = solve(dpi1dx==0, dpi2dx==0, x_1, x_2)
回答 (2 件)
Walter Roberson
2023 年 3 月 16 日
1 投票
Use dsolve for differential equations
20 件のコメント
Walter Roberson
2023 年 3 月 16 日
Note that if the derivative of a function is identical to 0 then the function must be constant with respect to that variable.
If you are working with initial conditions you would test things such as subs(dpi1dx, x1, 0)==0 (but remember to include the defining equations as well.) See the dsolve documentation
Roy
2023 年 3 月 17 日
Walter Roberson
2023 年 3 月 17 日
Sorry, I did misunderstand what you are trying to do.
solve() is for trying to find indefinitely-precise closed form solutions. The derivative of pi1 involves x_1^(2*r+1) . Even if we assume that r is a positive integer, there is no general solution for polynomials with degree greater than 4.
Therefore there are only a small number of values of r that have closed form solutions.
Roy
2023 年 3 月 18 日
Walter Roberson
2023 年 3 月 18 日
If you show your paper calculations for r = 3 and for r = 3/2 then perhaps someone would be able to come up with something. Use V_1 = 1 to make the calculation easier.
syms V_1 V_2 x_1 x_2 r
pi1 = (V_1) * (x_1^r/(x_1^r+x_2^r)) - x_1
pi2 = (V_2) * (x_2^r/(x_1^r+x_2^r)) - x_2
dpi1dx = diff(pi1, x_1)
dpi2dx = diff(pi2, x_2)
eqn = subs([dpi1dx, dpi2dx], V_2, V_1)
simplify(subs(eqn, [x_1, x_2], [r*V_1/2, r*V_1/2]), 'steps', 20)
These are not 0, so x_1 == x_2 == r*V/2 is not a root of the derivatives, and therefore is not a critical point.
I made some further tests. You can solve dpi1dx for x_2 or you can solve dpi2dx for x_1 but you cannot get anywhere on the next steps, which involve exp(i*pi*ANGLE) times something. You can rewrite to get (sin(ANGLE) + 1i*cos(ANGLE)) times something, but doing that does not help.
@Roy, really? Is r*V/4 really the solution?
syms V_1 V_2 x_1 x_2 r V
pi1 = (V_1) * (x_1^r/(x_1^r+x_2^r)) - x_1;
pi2 = (V_2) * (x_2^r/(x_1^r+x_2^r)) - x_2;
subs(pi1,[x_1,x_2],[r*V/4,r*V/4])
subs(pi2,[x_1,x_2],[r*V/4,r*V/4])
If it was, I would have thought the result, from substituting your claimed "solution" directly back into those equations, it would yield 0. I might be getting old though.
Even if you now claim that, oh, you made a mistake, and you have r==2 exactly, it still fails, unless also V_1==V_2. But then why do you have two different variables?
Roy
2023 年 3 月 19 日
Under the simplification of x_1 and x_2 being equal:
syms x_1 x_2
syms r V positive
pi1 = (V) * (x_1^r/(x_1^r+x_2^r)) - x_1
pi2 = (V) * (x_2^r/(x_1^r+x_2^r)) - x_2
dpi1dx = diff(pi1, x_1)
dpi2dx = diff(pi2, x_2)
eqns = subs([dpi1dx==0, dpi2dx==0], [x_2], [x_1])
seqns = simplify(eqns)
simplify(solve(seqns(2), x_1), 'steps', 20)
Roy
2023 年 3 月 19 日
Roy
2023 年 3 月 21 日
If you set x_2 to be a constant multiple, c, of x_1, then
syms x_1 x_2
syms c r V_1 V_2 positive
pi1 = (V_1) * (x_1^r/(x_1^r+x_2^r)) - x_1
pi2 = (V_2) * (x_2^r/(x_1^r+x_2^r)) - x_2
dpi1dx = diff(pi1, x_1)
dpi2dx = diff(pi2, x_2)
eqns = subs([dpi1dx==0, dpi2dx==0], [x_2], [c*x_1])
seqns = simplify(eqns)
sol = solve(seqns(2), x_1, 'returnconditions', true)
simplify(sol.conditions)
simplify(sol.x_1, 'steps', 20)
seqns2 = simplify(subs(seqns(1), x_1, sol.x_1), 'steps', 20)
which is to say that if x_1 and x_2 have a particular ratio, then V_1 and V_2 must have the same ratio.
Or you could read this as saying that if you know the ratio of V_1 and V_2 then x_1 and x_2 must have the same ratio, and x_1 is as given in ans above.
Notice that with V_2 == V_1 * c then V_2 / c is V_1 so ans could be rewritten as V_1 * c^r * r / (c^r+1)^2
Walter Roberson
2023 年 3 月 21 日
There is no mathematics that is able to find a closed form solution for that. It requires finding the set of Z such that
V_1*r*Z^(r - 1)*x_2^r - Z^(2*r) - 2*Z^r*x_2^r - x_2^(2*r)
is 0. But as I discussed earlier, it has been proven by Able and Rosser that there does not exist a general "algebraic" solution for the case where the power (2*r) > 4 -- so r = 1 and r = 2 are the positive integer cases that can be solved .
There are also solutions for r = 1/5, r = 1/3, r = 1/2, r = 3/5, and possibly some others.
Roy
2023 年 3 月 21 日
Walter Roberson
2023 年 3 月 21 日
The problem is not solveable for most r .
For example for r = 3/2 then the solutions are
RootOf(4*Z^3*x_2^(3/2) + 2*Z^6 - 3*Z*x_2^(3/2)*V_1 + 2*x_2^3,Z)^2
which is the set of Z such that the expression 4*etc becomes 0. But notice the Z^6 part -- so you would need the closed-form solution for a degree 6 polynomial, and such solutions only exist if the expression can be factored into polynomials of degree 4 or lower.
If r = N/4 for odd integer N, then you need to solve something of degree either 2*N+4 (for small N) or degree 2*N (starting at N = 5). r = 1/5 and r = 3/5 are tractable (but long!!), the other N/5 are not tractable.
Roy
2023 年 3 月 21 日
Roy
2023 年 3 月 21 日
0 投票
3 件のコメント
Are you sure the two expressions below turn out to be 0 ?
syms V_1 V_2 x_1 x_2 r
pi1 = (V_1) * (x_1^r/(x_1^r+x_2^r)) - x_1
pi2 = (V_2) * (x_2^r/(x_1^r+x_2^r)) - x_2
dpi1dx = diff(pi1, x_1)
dpi2dx = diff(pi2, x_2)
simplify(subs(dpi1dx,[x_1 x_2],[V_2*(r*(V_2/V_1)^(r-1))/(1+(V_2/V_1)^r)^2,V_1*(r*(V_1/V_2)^(r-1))/(1+(V_1/V_2)^r)^2]))
simplify(subs(dpi2dx,[x_1 x_2],[V_2*(r*(V_2/V_1)^(r-1))/(1+(V_2/V_1)^r)^2,V_1*(r*(V_1/V_2)^(r-1))/(1+(V_1/V_2)^r)^2]))
If you add the assumption of positive then they do resolve to 0
syms V_1 V_2 x_1 x_2 r positive
pi1 = (V_1) * (x_1^r/(x_1^r+x_2^r)) - x_1
pi2 = (V_2) * (x_2^r/(x_1^r+x_2^r)) - x_2
dpi1dx = diff(pi1, x_1)
dpi2dx = diff(pi2, x_2)
simplify(subs(dpi1dx,[x_1 x_2],[V_2*(r*(V_2/V_1)^(r-1))/(1+(V_2/V_1)^r)^2,V_1*(r*(V_1/V_2)^(r-1))/(1+(V_1/V_2)^r)^2]))
simplify(subs(dpi2dx,[x_1 x_2],[V_2*(r*(V_2/V_1)^(r-1))/(1+(V_2/V_1)^r)^2,V_1*(r*(V_1/V_2)^(r-1))/(1+(V_1/V_2)^r)^2]))
カテゴリ
ヘルプ センター および File Exchange で Linear Algebra についてさらに検索
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!










