# Non-linear Algebraic 36 equations unsloved

48 visualizzazioni (ultimi 30 giorni)
MINATI PATRA il 11 Mag 2024
Commentato: MINATI PATRA il 17 Mag 2024 alle 1:00
Hi forum
I have 36 Non-linear Algebraic equations with 36 unknowns to find out.
'solve' comand failed.
'fsolve' comand requires initial choices which I dont have.
##### 1 CommentoMostra -1 commenti meno recentiNascondi -1 commenti meno recenti
Torsten il 11 Mag 2024
Modificato: Torsten il 11 Mag 2024
An idea ? Without the equations ?
Usually, solving 36 equations without no idea for a good initial guess is hopeless. But why not trying with x0 = ones(36,1) and "fsolve" ?

Accedi per commentare.

### Risposte (3)

John D'Errico il 11 Mag 2024
Actually, a starting value set like ones(36,1) is often a dangerous thing to do. Making all of the starting values the same can sometimes put the solver in a state where it cannot proceed, due to symmetries in the problem. A better choice of simple random numbers is a far better idea, perhaps rand(36,1), or randn(36,1).
Is this likely to be sufficient? Possibly, but that is impossible to know wothout even a clue at your actual equations. It MAY yield a solution, but you need to recognize that without a decent start, a solver can often get lost, and diverge to some bad place.
Worse, even if you use a random start point, this may be insufficient if some of the parameters have very different scalings than the others.
Yes, it is absolutely no surprise an analytical solver like solve failed. In fact, that is far more the expected result than to find a solution from solve on such a problem. (Though I once recall a nonlinear problem with 15 unknowns I posed, where solve found the correct solution, but that itself was a great surprise to me.) As such, you almost certainly will need to use a numerical solver like fsolve. The problem then lies in how to choose a starting set. And if you understand nothing about the parameters in the model to be able to provide an intelligent solution, well, expect random crapola much of the time.
Think of such a problem as if you were to put a blind person down on the face of the earth, and ask them to find the point of lowest altitude, given nothing more than a cane, and an altimeter. (Yes, an altimeter that works underwater, so they may also want scuba gear.) If you start them out in some completely random spot, what are the odds they will find the actual minimum point? Far more likely they may get stuck in a locally sub-optimal point, where there is no direction they can move that will allow them to go downhill anymore.
And I am sorry to say, there is no magical way to automatically choose an intelligent set of starting values. What can you do? You may be left with two alternatives.
1. Use a multi-start method. That has you start a solver like fsolve off with a randomly generated point. Save the result you get. Then repeat. Some of the time, you will get crap, non-convergence, whatever. But save the set of solutions where fsolve thinks it did converge, and choose from among that set of "results", taking the one that makes you happy.
2. Use a solver like GA, or one of the many other tools from the particle swarm family, or a simulated annealing tool, etc. These are tools that are not so heavily tied down to a start point. I'd start with GA, a genetic algorithm, as my personal suggestion.
##### 1 CommentoMostra -1 commenti meno recentiNascondi -1 commenti meno recenti
MINATI PATRA il 16 Mag 2024 alle 9:43
@ Dear John D'Errico
have a try

Accedi per commentare.

Alex Sha il 11 Mag 2024
@MINATI PATRA please giving out the code you used, or the detail description of your 36 equations
##### 1 CommentoMostra -1 commenti meno recentiNascondi -1 commenti meno recenti
MINATI PATRA il 15 Mag 2024 alle 12:45
%% Sorry for the late reply, Actually my Laptop was out of order, Please have a try with this code
tic
syms x N1 N2 N3 Re Peh Pem
N1 = 0.1;N2 = 0.1; N3 = 0.1; Re = 0.1; Peh = 0.1; Pem = 0.1;
A = sym('a', [1 10]);B = sym('b', [1 10]);C = sym('c', [1 8]);D = sym('d', [1 8]); j = 1:10; k = 1:8;
f = sum(A.*x.^(j-1));g = sum(B.*x.^(j-1)); t = sum(C.*x.^(k-1));p = sum(D.*x.^(k-1));
Df = diff(f); D2f = diff(Df); D3f = diff(D2f); D4f = diff(D3f); D1g = diff(g); D2g = diff(D1g); D1p = diff(p); D2p = diff(D1p); D1t = diff(t); D2t = diff(D1t);
fP = subs(f,x,1); fN = subs(f,x,-1); fdP = subs(Df,x,1); fdN = subs(Df,x,-1);gP = subs(g,x,1); gN = subs(g,x,-1);tP = subs(t,x,1); tN = subs(t,x,-1);pP = subs(p,x,1); pN = subs(p,x,-1);
F = (1+N1)*D4f - N1*D2g - Re*(f*D3f - Df*D2f); G = N2*D2g - N1*(D2f - 2*g) - N3*Re*(f*D2g - Df*g);T = D2t + Peh*(Df*t - f*D1t); P = D2p + Pem*(Df*p - f*D1p);
DF = diff(F); D2F = diff(DF); D3F = diff(D2F); D4F = diff(D3F); DG = diff(G); D2G = diff(DG);D3G = diff(D2G); DP = diff(P); D2P = diff(DP); DT = diff(T); D2T = diff(DT);
FP = subs(F,x,1); FN = subs(F,x,-1); FdP = subs(DF,x,1); FdN = subs(DF,x,-1); Fd2P = subs(D2F,x,1); Fd2N = subs(D2F,x,-1);
GP = subs(G,x,1); GN = subs(G,x,-1); GdP = subs(DG,x,1); GdN = subs(DG,x,-1); Gd2P = subs(D2G,x,1); Gd2N = subs(D2G,x,-1);Gd3P = subs(D3G,x,1); Gd3N = subs(D3G,x,-1);
TP = subs(T,x,1); TN = subs(T,x,-1); TdP = subs(DT,x,1); TdN = subs(DT,x,-1); Td2P = subs(D2T,x,1); Td2N = subs(D2T,x,-1);
PP = subs(P,x,1); PN = subs(P,x,-1); PdP = subs(DP,x,1); PdN = subs(DP,x,-1); Pd2P = subs(D2P,x,1); Pd2N = subs(D2P,x,-1);
e1 = fN == 0; e2 = fP == 0; e3 = fdN == 0; e4 = fdP == -1; e5 = gN == 0; e6 = gP == 1; e7 = pN == 1; e8 = pP == 0; e9 = tN == 1; e10 = tP == 0;
e11 = FN == 0; e12 = FP == 0; e13 = FdN == 0; e14 = FdP == 0;e15 = Fd2N == 0; e16 = Fd2P == 0;
e17 = GN == 0; e18 = GP == 0; e19 = GdN == 0; e20 = GdP == 0;e21 = Gd2N == 0; e22 = Gd2P == 0;e23 = Gd3N == 0; e24 = Gd3P == 0;
e25 = TN == 0; e26 = TP == 0; e27 = TdN == 0; e28 = TdP == 0; e29 = Td2N == 0; e30 = Td2P == 0;
e31 = PN == 0; e32 = PP == 0; e33 = PdN == 0; e34 = PdP == 0;e35 = Pd2N == 0; e36 = Pd2P == 0;
EQN = [e1 e2 e3 e4 e5 e6 e7 e8 e9 e10 e11 e12 e13 e14 e15 e16 e17 e18 e19 e20 e21 e22 e23 e24 e25 e26 e27 e28 e29 e30 e31 e32 e33 e34 e35 e36];
VARS = [A B C D];
S = solve(EQN);
toc
S.a1

Accedi per commentare.

Torsten il 16 Mag 2024 alle 10:18
Spostato: Torsten il 16 Mag 2024 alle 12:36
Seems your system has several solutions.
syms x N1 N2 N3 Re Peh Pem
N1 = 0.1;N2 = 0.1; N3 = 0.1; Re = 0.1; Peh = 0.1; Pem = 0.1;
A = sym('a', [1 10]);B = sym('b', [1 10]);C = sym('c', [1 8]);D = sym('d', [1 8]); j = 1:10; k = 1:8;
f = sum(A.*x.^(j-1));g = sum(B.*x.^(j-1)); t = sum(C.*x.^(k-1));p = sum(D.*x.^(k-1));
Df = diff(f); D2f = diff(Df); D3f = diff(D2f); D4f = diff(D3f); D1g = diff(g); D2g = diff(D1g); D1p = diff(p); D2p = diff(D1p); D1t = diff(t); D2t = diff(D1t);
fP = subs(f,x,1); fN = subs(f,x,-1); fdP = subs(Df,x,1); fdN = subs(Df,x,-1);gP = subs(g,x,1); gN = subs(g,x,-1);tP = subs(t,x,1); tN = subs(t,x,-1);pP = subs(p,x,1); pN = subs(p,x,-1);
F = (1+N1)*D4f - N1*D2g - Re*(f*D3f - Df*D2f); G = N2*D2g - N1*(D2f - 2*g) - N3*Re*(f*D2g - Df*g);T = D2t + Peh*(Df*t - f*D1t); P = D2p + Pem*(Df*p - f*D1p);
DF = diff(F); D2F = diff(DF); D3F = diff(D2F); D4F = diff(D3F); DG = diff(G); D2G = diff(DG);D3G = diff(D2G); DP = diff(P); D2P = diff(DP); DT = diff(T); D2T = diff(DT);
FP = subs(F,x,1); FN = subs(F,x,-1); FdP = subs(DF,x,1); FdN = subs(DF,x,-1); Fd2P = subs(D2F,x,1); Fd2N = subs(D2F,x,-1);
GP = subs(G,x,1); GN = subs(G,x,-1); GdP = subs(DG,x,1); GdN = subs(DG,x,-1); Gd2P = subs(D2G,x,1); Gd2N = subs(D2G,x,-1);Gd3P = subs(D3G,x,1); Gd3N = subs(D3G,x,-1);
TP = subs(T,x,1); TN = subs(T,x,-1); TdP = subs(DT,x,1); TdN = subs(DT,x,-1); Td2P = subs(D2T,x,1); Td2N = subs(D2T,x,-1);
PP = subs(P,x,1); PN = subs(P,x,-1); PdP = subs(DP,x,1); PdN = subs(DP,x,-1); Pd2P = subs(D2P,x,1); Pd2N = subs(D2P,x,-1);
e1 = fN == 0; e2 = fP == 0; e3 = fdN == 0; e4 = fdP == -1; e5 = gN == 0; e6 = gP == 1; e7 = pN == 1; e8 = pP == 0; e9 = tN == 1; e10 = tP == 0;
e11 = FN == 0; e12 = FP == 0; e13 = FdN == 0; e14 = FdP == 0;e15 = Fd2N == 0; e16 = Fd2P == 0;
e17 = GN == 0; e18 = GP == 0; e19 = GdN == 0; e20 = GdP == 0;e21 = Gd2N == 0; e22 = Gd2P == 0;e23 = Gd3N == 0; e24 = Gd3P == 0;
e25 = TN == 0; e26 = TP == 0; e27 = TdN == 0; e28 = TdP == 0; e29 = Td2N == 0; e30 = Td2P == 0;
e31 = PN == 0; e32 = PP == 0; e33 = PdN == 0; e34 = PdP == 0;e35 = Pd2N == 0; e36 = Pd2P == 0;
EQNS = [e1 e2 e3 e4 e5 e6 e7 e8 e9 e10 e11 e12 e13 e14 e15 e16 e17 e18 e19 e20 e21 e22 e23 e24 e25 e26 e27 e28 e29 e30 e31 e32 e33 e34 e35 e36];
EQNS = lhs(EQNS)-rhs(EQNS);
VARS = [A B C D];
EQNS = matlabFunction(EQNS,"Vars",{VARS});
VARS0 = rand(1,36);
VARS = fsolve(EQNS,VARS0)
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.
VARS = 1x36
3.6725 0.0906 -7.4061 0.1528 4.0630 -0.3481 -0.3478 0.1256 0.0184 -0.0208 -520.1053 79.3286 666.5327 -148.6066 -177.2206 97.7998 35.9712 -31.6665 -4.6780 3.6447 0.6472 -0.5924 -0.1550 0.1378 -0.0042 -0.0492 0.0120 0.0037 0.6472 -0.5924
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
A = VARS(1:10)
A = 1x10
3.6725 0.0906 -7.4061 0.1528 4.0630 -0.3481 -0.3478 0.1256 0.0184 -0.0208
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
B = VARS(11:20)
B = 1x10
-520.1053 79.3286 666.5327 -148.6066 -177.2206 97.7998 35.9712 -31.6665 -4.6780 3.6447
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
C = VARS(21:28)
C = 1x8
0.6472 -0.5924 -0.1550 0.1378 -0.0042 -0.0492 0.0120 0.0037
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
D = VARS(29:36)
D = 1x8
0.6472 -0.5924 -0.1550 0.1378 -0.0042 -0.0492 0.0120 0.0037
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
norm(EQNS(VARS))
ans = 4.1370e-11
##### 3 CommentiMostra 1 commento meno recenteNascondi 1 commento meno recente
Torsten il 16 Mag 2024 alle 17:13
Note that each run with initial conditions changed seems to give a different solution. I don't know if this is to be expected for your underlying problem.
MINATI PATRA il 17 Mag 2024 alle 1:00
Ok, I will take care of that. Thanks for your concern.

Accedi per commentare.

### Categorie

Scopri di più su Programming in Help Center e File Exchange

### Community Treasure Hunt

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

Start Hunting!

Translated by