Error solving nonlinear system of equations

[iils-jdinkelacker]

Solver fails to solve the following nonlinear system of equations (4th order terms)

th = 2400
tk = 320
Solve({th^4 - tw^4 -1.77E-6*qw == 0,
		tw^2 - t1^2 +2454*(tw - t1) - 0.2319 *qw ==0,
		t1^4 - t2^4 +4.601E6*(t1^2 - t2^2) + 1.129E10*(t1 - t2) - 1.6229E8 * qw == 0,
		t2^2 - t3^2 + 2454 *(t2 -t3) -0.2319 *qw ==0,
		tw -t3 -0.0549*qs == 0,
		t3^4 - tb^4 -4.141E5*(qs + qw) == 0,
		tb - tk -8.164E-8*(qs +qw) == 0}, {tw,t1,t2,t3,tb,qw,qs})

Error result:

    Error: Rationalize: No converging fraction found for value 0.6574477170497257.Solve: The system cannot be solved with the methods available to Solve.

This system of equations might to be solvable with the current alvailable algorithms in Symja. Trust-region algorithms with start values (fsolve) in Matlab and Pyhton SciPy get good approximate results. I could not find any such functions in Symja, but they might be a good addition to the functionality.

[axkr]

Not a solution only a comment for Symja syntax usage.
The expression contains float numbers like "1.129E10" which should be rewritten in Symja syntax as 1.129*^10.

Solve({th==2400,tk==320,th^4 - tw^4 -1.77*^-6*qw == 0,
        tw^2 - t1^2 +2454*(tw - t1) - 0.2319 *qw ==0,
        t1^4 - t2^4 +4.601*^6*(t1^2 - t2^2) + 1.129*^10*(t1 - t2) - 1.6229*^8 * qw == 0,
        t2^2 - t3^2 + 2454 *(t2 -t3) -0.2319 *qw ==0,
        tw -t3 -0.0549*qs == 0,
        t3^4 - tb^4 -4.141*^5*(qs + qw) == 0,
        tb - tk -8.164*^-8*(qs +qw) == 0}, {tw,t1,t2,t3,tb,qw,qs,th,tk})

Internally I get division by 0 problems for this expression.

[iils-hwellmann]

Trust-region algorithms with start values (fsolve) in Matlab and Pyhton SciPy get good approximate results. I could not find any such functions in Symja, but they might be a good addition to the functionality.

@iils-jdinkelacker can you list the guess-values that lead to these good results?

[iils-jdinkelacker]

Thanks for the syntax tip. Now I also get a lot of errors liks "Power: Indeterminate expression (-1.7710^-6qw+th^4-tw^4)/(-1.13589*10^-24) encountered.", which I assume are the division by 0 problems.

The start values for scipy fsolve are (th = 2400, tw = 2000, t1 = 1400, t2 = 1200, t3 = 1000, tb = 800, tk = 320, qw = 400000, qs = 30000) and get results (2400.000000000453, 2399.9993912855166, 2385.212934513328, 708.4971515532554, 680.6456294208236, 320.04024139792955, 319.99999999963455, 461594.8472077818, 31317.919159643687)

[axkr]

The start values for scipy fsolve are (th = 2400, tw = 2000, t1 = 1400, t2 = 1200, t3 = 1000, tb = 800, tk = 320, qw = 400000, qs = 30000) and get results (2400.000000000453, 2399.9993912855166, 2385.212934513328, 708.4971515532554, 680.6456294208236, 320.04024139792955, 319.99999999963455, 461594.8472077818, 31317.919159643687)

If I'm doing a "cross-check" with your scipy results by inserting the variable values into the expressions:

{th^4 - tw^4 -1.77*^-6*qw,
 tw^2 - t1^2 +2454*(tw - t1) - 0.2319 *qw,
 t1^4 - t2^4 +4.601*^6*(t1^2 - t2^2) + 1.129*^10*(t1 - t2) - 1.6229*^8 * qw,
 t2^2 - t3^2 + 2454 *(t2 -t3) -0.2319 *qw,
 tw -t3 -0.0549*qs,
 t3^4 - tb^4 -4.141*^5*(qs + qw),
 tb - tk -8.164*^-8*(qs +qw)} /. {th->2400.000000000453,tw->2399.9993912855166,t1->2385.212934513328,t2->708.4971515532554,t3->680.6456294208236,tb->320.04024139792955,tk->319.99999999963455,qw->461594.8472077818,qs->31317.919159643687}

I get:

{33659487.499802814,-1.544947365003024,22487008.039317191,-6.4688450895313243,0.0000000002545837,20731337.636887173,0.00000000004876339}

Is the result really sufficient for your purposes?

[iils-jdinkelacker]

No, you are right. It seems, that this solution is pretty far of. The full output of fsolve gives the same picture for the errors [ 6.36292317e+10, 1.26735730e+04, -2.89330158e+11, 1.12445091e+03, -1.42911858e-07, -1.04461625e+10, -2.97211720e-08, -2.73585556e-07, 2.03255468e-07]
It still gave the message that it converged and no warning. With very different start values fsolve still converges towards this result. It didn't raise suspicions, because the results for Q are consistent with data from a document, which this is system of equations is derived from.

So there seems to be something wrong with the set up of this system of equation.

[iils-jdinkelacker]

I had one more look at the equation system:

th = 2400
tk = 320
Solve({th^4 - tw^4 - 1.77*^-6 * qw == 0,
		tw^2 - t1^2 + 2454 * (tw - t1) - 0.2319 * qw == 0,
		t1^4 - t2^4 + 4.601*^6 * (t1^2 - t2^2) + 1.129*^10 * (t1 - t2) - 1.6229*^8 * qw == 0,
		t2^2 - t3^2 + 2454 *(t2 -t3) -0.2319 * qw == 0,
		tw - t3 - 0.0549 * qs == 0,
		t3^4 - tb^4 - 4.141*^5 * (qs + qw) == 0,
		 tb - tk - 8.164*^-6 * (qs + qw) == 0}, {tw,t1,t2,t3,tb,qw,qs})

(the last equation the wrong factor 8.164*^-8 --> 8.164*^-6)

Yields the following Error:
Power: Indeterminate expression (-1.77*10^-6*qw+th^4-tw^4)/(-1.13589*10^-24) encountered.FindRoot: maximal count (100) exceeded.Power: Indeterminate expression (-8.164*e-6.0*qs-63492.88*t2-25.87322*t2^2+63492.88*t3+25.87322*t3^2+tb-tk)/(-3.01204*10^-16) encountered.Power: Indeterminate expression (-9.08426*10^-6*qs^4-0.0187304*t2-7.6326*10^-6*t2^2+0.0187304*t3-0.000661877*qs^3*t3+7.6326*10^-6*t3^2-0.0180841*qs^2*t3^2-0.2196*qs*t3^3-t3^4+th^4)/(-3.51765*10^-26) encountered.Power: Indeterminate expression (-9.08426*10^-6*qs^4-0.0187304*t2-7.6326*10^-6*t2^2+0.0187304*t3-0.000661877*qs^3*t3+7.6326*10^-6*t3^2-0.0180841*qs^2*t3^2-0.2196*qs*t3^3-t3^4)/(-3.51765*10^-26) encountered.InverseFunction: Inverse functions are being used. Values may be lost for multivalued inverses.Solve: The system cannot be solved with the methods available to Solve.Set: Symbol Times is Protected.$RecursionLimit: Recursion depth of 256 exceeded during evaluation of th=2400*tk=320*Solve({-1.77*10^-6*qw+th^4-tw^4==0,-0.2319*qw+2454*(-t1+tw)+-t1**2+tw**2==0,-1.6229*10^8*qw+1.129*<<SHORT>>t3^2==0,-0.0549*qs-t3+tw==0,-414100.0*(qs+qw)+t3**4-tb**4==0,-8.164*e-6*(qs+qw)+tb-tk==0},{tw,t1,t2,t3,tb,qw,qs}).

with scipy.optimize least_squares, I get good results with initial values for th, tw, t1, t2, t3, tb, tk, qw, qs initial_guess = [2400.0, 2000.0, 1400.0, 1200.0, 1000.0, 800.0, 320.0, 400000.0, 30000.0]

result:

th = 2400.0
tw = 2399.9999999999854
t1 = 2385.214748463352
t2 = 708.8789196244769
t3 = 681.0339568421458
tb = 324.02372187906377
tk = 320.0
qw = 461550.7157556547
qs = 31310.85688812094

with small errors:

eq1: -0.012257266887508722
eq2: 1.3533281162381172e-09
eq3: -0.015625
eq4: 5.820766091346741e-11
eq5: 0.0
eq6: 0.0
eq7: -1.3322676295501878e-14
eq8: 0.0
eq9: 0.0

Maybe least squares method with initial values would be a good addition to the functionality.

[axkr]

Maybe least squares method with initial values would be a good addition to the functionality.

or FindFit could be extended for such cases:

• https://github.com/axkr/symja_android_library/blob/master/symja_android_library/doc/functions/FindFit.md