symengine
Symbolic solution of linear equation system returns wrong results
Related to #442. Solving a linear system Ax=b symbolically via A.LUsolve(b) returns a wrong result. Both when using .subs or Lambdify. The result contains values of up to e+100. When solving the linear system numerically I get the correct result, with no value exceeding e-02.
What causes this difference? In this comment (https://github.com/sympy/sympy/issues/25208#issuecomment-1581122169) one contributor of sympy also mentions the accuracy difference. Where does this come from?
Small update: By converting every number in A and b from python's float to sympy's (and later symengine's) Float and specifying a precision via the dps argument, I got quite some improvement in the result. However, the result is still quite off and not yet usable. Even setting dps=1000 does not help.
With symbolic values, you can't do pivoting correctly. LU decomposition using Gaussian elimination can be numerically unstable when small values are chosen as pivots and roundoff errors become larger.
