Numerical solution of complexSolve (decimal digits displayed)

[Grégory Vanuxem]

Hello,

Any idea why I obtain a huge numerical precision using complexSolve here, at least, displayed?

(1) -> complexSolve(x^3+3*x^2+x+17=0, 0.0000000001)

   (1)
   [x = - 3.8744001380_296332337,
    x = 0.4372000690_1481618814_9406395091_6108308765_2840423528 - 2.0485683086_313182244 %i,

       x
     =
         0.4372000690_1481618814_9406395091_6108308765_2840423528
       +
         2.0485683086_3131822439_9819882819_4290399551_391601562 %i
     ]
                                                              Type: List(Equation(Polynomial(Complex(Float))))
(2) -> precision()$Float

   (2)  68
                                                                                         Type: PositiveInteger

(3) -> digits()$Float

   (3)  20

Regards,

Greg

[Qian Yun]

It is happening in numsolve.spad, function "findGenZeros"
that precision is increased by 100.

I think in the attached patch, by reordering the transform
to be after resetting precision, this problem is solved.

- Qian

> --
> You received this message because you are subscribed to the Google
> Groups "FriCAS - computer algebra system" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to fricas-devel...@googlegroups.com <mailto:fricas-
> devel+un...@googlegroups.com>.
> To view this discussion visit https://groups.google.com/d/msgid/fricas-
> devel/CAHnU2dZfeWi0kGYXSMQWzi1m-
> FmT_giYQUqU9qt2rpUMZ%3D9CJg%40mail.gmail.com <https://groups.google.com/
> d/msgid/fricas-devel/CAHnU2dZfeWi0kGYXSMQWzi1m-
> FmT_giYQUqU9qt2rpUMZ%3D9CJg%40mail.gmail.com?
> utm_medium=email&utm_source=footer>.

[Waldek Hebisch]

By default FriCAS displays all computed digits. If that bothers
you, you can do:

(1) -> outputGeneral(20)
Type: Void
(2) -> complexSolve(x^3+3*x^2+x+17=0, 0.0000000001)

(2)
[x = - 3.8744001380_296332337,
x = 0.4372000690_1481618815 - 2.0485683086_313182244 %i,
x = 0.4372000690_1481618815 + 2.0485683086_313182244 %i]
Type: List(Equation(Polynomial(Complex(Float))))

--
Waldek Hebisch

[Waldek Hebisch]

On Wed, Jun 03, 2026 at 05:34:55PM +0800, Qian Yun wrote:
> It is happening in numsolve.spad, function "findGenZeros"
> that precision is increased by 100.
>
> I think in the attached patch, by reordering the transform
> to be after resetting precision, this problem is solved.

I am not sure if there is a problem. Namely, to compute result
at requested accuracy we need to increase intermediate precision.
Since final rounding may add to total error computing "exact
20 digits" (or similar) is rather complex process. So, instead,
normal policy in FriCAS is to keep results at whatever precision
they were computed. For the stated problem, that is display we
have controls like 'outputGeneral' which allow setting number
of displayed digits without affecting future use in computations.

Patch like below are tricky, because in 99.999999% of cases they
work fine. But there could be some pesky combination of input
data were added rounding goes outside requested bounds. To
prove that this does not happen requires effort. And if you
aim at "exact digits" then usually this requires much more
complex approach, basically estimating error and using error
estimate to decide if rounding could change displayed digits.
If not, then one rounds, otherwise repeats computation at
higher precision.

It is not clear to me if it is worth going in such direcition.
Current approach allows getting needed accuracy with reasonable
effort and output controls are enough for nicer printouts.

> To unsubscribe from this group and stop receiving emails from it, send an email to fricas-devel...@googlegroups.com.
> To view this discussion visit https://groups.google.com/d/msgid/fricas-devel/887ce321-5acf-4b62-bc36-6610034a0c21%40gmail.com.

> diff --git a/src/algebra/numsolve.spad b/src/algebra/numsolve.spad
> index 2400b653..5db6c0f5 100644
> --- a/src/algebra/numsolve.spad
> +++ b/src/algebra/numsolve.spad
> @@ -246,9 +246,9 @@ InnerNumericFloatSolvePackage(K, F, Par) : Cat == Cap where
> nfeps := (1/etol)*nfeps
> lz := innerSolve1(f, neps)
> ok := true
> - sol : L L F := []
> + sol : L L CI := []
> for z in lz while ok repeat
> - sol1 : L F := [CI_to_F(F_to_CI1(z, feps))]
> + sol1 : L CI := [F_to_CI1(z, feps)]
> for pol in rlp for xvar in rest(rlvar) repeat
> pp := ieval(pol, xvar, zvar, z, nfeps)
> pp case "failed" =>
> @@ -258,11 +258,12 @@ InnerNumericFloatSolvePackage(K, F, Par) : Cat == Cap where
> width(real(ppi)) > feps or width(imag(ppi)) > feps =>
> ok := false
> break
> - sol1 := cons(CI_to_F(ppi), sol1)
> + sol1 := cons(ppi, sol1)
> sol := cons(sol1, sol)
> ok =>
> bits(obits)
> - return reverse(sol)
> + res := [[CI_to_F(x) for x in l] for l in sol]
> + return reverse!(res)
> etol := etol^2
> ebits := 2*ebits
>


--
Waldek Hebisch

[Grégory Vanuxem]

Le mer. 3 juin 2026 à 12:21, Waldek Hebisch <de...@fricas.org> a écrit :

On Wed, Jun 03, 2026 at 08:11:01AM +0200, Grégory Vanuxem wrote:
> Hello,
>
> Any idea why I obtain a huge numerical precision using complexSolve here,
> at least, displayed?
>
> (1) -> complexSolve(x^3+3*x^2+x+17=0, 0.0000000001)
>
>    (1)
>    [x = - 3.8744001380_296332337,
>     x = 0.4372000690_1481618814_9406395091_6108308765_2840423528 -
> 2.0485683086_313182244 %i,
>
>        x
>      =
>          0.4372000690_1481618814_9406395091_6108308765_2840423528
>        +
>          2.0485683086_3131822439_9819882819_4290399551_391601562 %i
>      ]
>                                                               Type:
> List(Equation(Polynomial(Complex(Float))))
> (2) -> precision()$Float
>
>    (2)  68
>
>              Type: PositiveInteger
>
> (3) -> digits()$Float
>
>    (3)  20
>

By default FriCAS displays all computed digits.  If that bothers
you, you can do:

 

Ah ok, I was surprised. I didn't know or do not remember. I better understand this function now. I remember when you modified the output* function for Float and that was not very clear to me this function. In fact, even the documentation tells almost nothing about this feature ("outputGeneral(n) sets the output mode to general notation with n significant digits displayed.").

Good,

Greg

PS: BTW sorry for the bad pull request this afternoon, I was experimenting with Julia related category/domain/package GitHub generated documentation pages (a personal Python script to update github.io).

 

(1) -> outputGeneral(20)
                                                                   Type: Void
(2) -> complexSolve(x^3+3*x^2+x+17=0, 0.0000000001)

   (2)
   [x = - 3.8744001380_296332337,
    x = 0.4372000690_1481618815 - 2.0485683086_313182244 %i,
    x = 0.4372000690_1481618815 + 2.0485683086_313182244 %i]
                             Type: List(Equation(Polynomial(Complex(Float))))

--
                              Waldek Hebisch

--
You received this message because you are subscribed to the Google Groups "FriCAS - computer algebra system" group.
To unsubscribe from this group and stop receiving emails from it, send an email to fricas-devel...@googlegroups.com.

To view this discussion visit https://groups.google.com/d/msgid/fricas-devel/aiAANrVv-Xu71PPN%40fricas.org.

[Grégory Vanuxem]

Le mer. 3 juin 2026 à 11:35, Qian Yun <oldk...@gmail.com> a écrit :

It is happening in numsolve.spad, function "findGenZeros"
that precision is increased by 100.

Wow

 

[Qian Yun]

On 6/3/26 6:38 PM, Waldek Hebisch wrote:
> On Wed, Jun 03, 2026 at 05:34:55PM +0800, Qian Yun wrote:
>> It is happening in numsolve.spad, function "findGenZeros"
>> that precision is increased by 100.
>>
>> I think in the attached patch, by reordering the transform
>> to be after resetting precision, this problem is solved.
>
> I am not sure if there is a problem. Namely, to compute result
> at requested accuracy we need to increase intermediate precision.
> Since final rounding may add to total error computing "exact
> 20 digits" (or similar) is rather complex process. So, instead,
> normal policy in FriCAS is to keep results at whatever precision
> they were computed.

But the increased precision is temporary, and not all bits are
accurate. Also "complexSolve" comes with a parameter "eps",
which dictates the accuracy requirement from user.

- Qian

[Waldek Hebisch]

Consider first some basics. Using n-bit integer mantissa we
have at most 1/2 of error due to rounding. Mantissa itself
is at least 2^(n - 1), which leads to 2^(-n) of releative error.
This error is due to representation, so in a sense unavidable.
Our basic arithmetic operations provide results with such an
accuracy.

Now, in general I would prefer to deliver result with relative
error at most 2^(-n). This is not always reasonable, but I would
like to have such accuracy when it is relatively easy to attain.
Since such an error can be produced by rounding alone, we either
would need to do some heroic work to deliver very accurate
result which even after rounding has required accuracy or
we do not round and return more bits. The second is much
easier.

To put it a bit differently, unless we internaly compute with
much higher accuracy some of printed digits will be inaccurate.
I find it more reasonable to print more digits and have
guaranteed error bound than to print less digits and less
accurate results. However, if user prefers the second, then
'outputGeneral' still allows printing with prescribed number
of digits. Changing computational routines to always round
does not give choice: once accuracy is lost it can not be
easily restored.

Maybe a little theoretical example. Consider say 10 bit
accuracy and exact reasult which in binary is 1000000000.1001.
"Correct" rounding would give 1000000001. Computation may give
something like 1000000000.01110. After rounding it becomes
1000000000. So the relative error is bigger than desired
bound. OTOH if we keep few extra bits the result is within
bound.

--
Waldek Hebisch

[Qian Yun]

My point is, for this specific "complexSolve" function,
say user requests 10 digits accuracy, FLOAT domain has
20 digits accuracy by default, the result of "complexSolve"
is accurate up to 15 digits (a random guess), but it is
a 50 digits float.

If digits()$Float is x, then user requests accuracy greater
than x, we should signal an error.
If the eps accuracy is smaller than x, then we can return
results with x digits perfectly fine.

- Qian