`FindRoot [ ]`::jsing: Encountered a singular Jacobian at a point...WHY
Fixing the typo 1/\[Rho] t1
--> 1/\[Rho]t1
, I can reproduce the error. I think it is just that the Jacobian is too poorly conditioned at the initial points:
LinearSolve[
jacob /. parameters /. Rule @@@ initialVal]["ConditionNumber"]
LinearSolve::luc: Result for LinearSolve of badly conditioned matrix {{0,-1.37143,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},<<48>>} may contain significant numerical errors.
(* 4.75451*10^16 <-- condition number > 10.^MachinePrecision *)
Use a working precision of roughly around $\log_{10}(\hbox{condition number})$ plus the PrecisionGoal
:
sol = FindRoot[SetPrecision[system /. parameters, 24],
SetPrecision[givenPoint, 24], PrecisionGoal -> 8,
WorkingPrecision -> 24] // N
(*
{EG -> 116.757, EXR -> 1.57378, FSAV -> 10.052, IADJ -> 1.51676,
PA1 -> 0.796624, PA2 -> 1.09585, PD1 -> 0.505239, PD2 -> 1.09585,
PE1 -> 1.57378, PM2 -> 2.15832, PQ1 -> 0.52553, PQ2 -> 1.19291,
PVA1 -> 0.424617, PVA2 -> 0.722743, PX1 -> 0.796624, PX2 -> 1.09585,
QA1 -> 330.749, QA2 -> 337.713, QD1 -> 361.553, QD2 -> 337.713,
QE1 -> 51.3489, QF11 -> 104.764, QF12 -> 72.2356, QF21 -> 70.705,
QF22 -> 137.295, QH11 -> 61.6308, QH12 -> 97.7302, QH21 -> 149.331,
QH22 -> 80.8366, QINT11 -> 99.5802, QINT12 -> 47.1428,
QINT21 -> 59.2738, QINT22 -> 84.8569, QINV1 -> 42.4692,
QINV2 -> 128.924, QM2 -> 132.401, QQ1 -> 361.553, QQ2 -> 570.222,
QX1 -> 330.749, QX2 -> 337.713, WALRAS -> -5.72072*10^-14,
WF1 -> 0.936218, YF11 -> 94.7432, YF12 -> 125., YF21 -> 81.7784,
YF22 -> 82.9999, YG -> 159.551, YH1 -> 307.694, YH2 -> 194.959}
*)
The solution is accurate at machine precision:
system /. parameters /. sol
(*
{0., 0., 0., 0., 0., 0., 0., 0., 0., 2.22045*10^-16, -5.68434*10^-14, 5.68434*10^-14,
1.11022*10^-16, 0., 0., -2.22045*10^-16, 0., 0., 0.,
-1.42109*10^-14, 0., 0., 1.13687*10^-13, 1.11022*10^-16, 0., -5.68434*10^-14,
2.77556*10^-17, 0., -1.42109*10^-14, 0., -1.42109*10^-14, 0., 0.,
1.42109*10^-14, 1.42109*10^-14, 0., 0., 0., 0., 0., 0., 0., 0., -2.84217*10^-14, 0.,
5.68434*10^-14, 0., 0., 0.}
*)
Numerical analysis remark: Why is the error message "singular" and not "badly conditioned"?
By Gastinel's Theorem (which may be found here), a matrix differs from a singular one by a relative error of $1/K$, where $K$ is the condition number (with the size of the error and condition number measured with respect to the same operator norm). When that relative error is less than the error to be expected from ordinary rounding, then the matrix might have been obtained by rounding a singular matrix to the working precision. I don't know what tolerance FindRoot
uses to determine whether a matrix is effectively singular, but given that in this case $1/K$ is less than $MachineEpsilon
, which is the smallest expected nonzero relative rounding error in a single entry, it is not surprising that the Jacobian is considered singular. Artificially increasing the precision of the data in parameters
and initialVal
does not really change this fact. But since FindRoot
is an iterative method, a bad first step, whether or not there was one, could be corrected in subsequent iterations. Further, FindRoot
uses a "damped" method that shortens extremely large steps.
When you fix the typo you can use the new in 12.0 affine covariant Newton solver:
sol = FindRoot[system /. parameters, givenPoint,
Method -> {"AffineCovariantNewton"}];
MinMax[system /. parameters /. sol]
{-6.821210263296962`*^-13, 2.50716114535976`*^-13}
More documentation is here.
sol
{EG -> 116.757, EXR -> 1.57378, FSAV -> 10.052, IADJ -> 1.51676,
PA1 -> 0.796624, PA2 -> 1.09585, PD1 -> 0.505239, PD2 -> 1.09585,
PE1 -> 1.57378, PM2 -> 2.15832, PQ1 -> 0.52553, PQ2 -> 1.19291,
PVA1 -> 0.424617, PVA2 -> 0.722743, PX1 -> 0.796624, PX2 -> 1.09585,
QA1 -> 330.749, QA2 -> 337.713, QD1 -> 361.553, QD2 -> 337.713,
QE1 -> 51.3489, QF11 -> 104.764, QF12 -> 72.2356, QF21 -> 70.705,
QF22 -> 137.295, QH11 -> 61.6308, QH12 -> 97.7302, QH21 -> 149.331,
QH22 -> 80.8366, QINT11 -> 99.5802, QINT12 -> 47.1428,
QINT21 -> 59.2738, QINT22 -> 84.8569, QINV1 -> 42.4692,
QINV2 -> 128.924, QM2 -> 132.401, QQ1 -> 361.553, QQ2 -> 570.222,
QX1 -> 330.749, QX2 -> 337.713, WALRAS -> -4.07116*10^-14,
WF1 -> 0.936218, YF11 -> 94.7432, YF12 -> 125., YF21 -> 81.7784,
YF22 -> 82.9999, YG -> 159.551, YH1 -> 307.694, YH2 -> 194.959}