bug?: difficulty to extract imaginary part

34 views
Skip to first unread message

Ralf Hemmecke

unread,
May 14, 2024, 9:04:18 AM5/14/24
to fricas-devel
Admittedly, it might be difficult to extract the imaginary part of a
radical expression. But that seems to look like a bug.

%%% (412) -> xx :=
(sqrt((((-12669586846893008563685878644359325*sqrt(3)+15974693204716576905254784724655502)*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)^2+(-676159865925870354666914169665824080551250000*sqrt(3)+8217402890479177097807248979603419220411058375726750000)*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)-45931279578807095346392460907109478385796132783696494926118495077505668392625000000000)*sqrt(((-5401066809391479730461765*sqrt(3)+6810039372933885991985941)*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)^2+3503093079743488011787977136025306152751500000*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)+19580584011451989397556080279474514673471759017621573422247980153750000000000)/138)+(75458077201551937891754164660043236993041255703201346421110021520000000*sqrt(3)-95142773619263899472829762733153800309338692262462996235137171088000000)*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)^2-48941566061703353092027102734653684371144117066589419560010587508867726730881752000000000000*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)+547119031158247887135235714889320242561133043433631241769109111538886308888911873326867936830640527579160000000000000000000)/32996199040131120862458002308760594)-7644000*sqrt(((-5401066809391479730461765*sqrt(3)+6810039372933885991985941)*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)^2+3503093079743488011787977136025306152751500000*nthRoot(137198597903998437385921091448494025437519087343750000000000000*sqrt(3)+172989501261663064201623395407144143954141328843750000000000000,
3)+19580584011451989397556080279474514673471759017621573422247980153750000000000)/138)+91052936304600216411782349862348428000000000)/133412830008771839271564000000;

Type:
AlgebraicNumber
%%% (413) -> ee := xx::Expression(INT);

Type:
Expression(Integer)
%%% (414) -> imag(ee) $ TrigonometricManipulations(ZZ, EX(ZZ))

(414) 0
Type:
Expression(Integer)
%%% (415) -> imag ee

(415) 0
Type:
Expression(Integer)
%%% (416) -> ee::Complex(Float)

(416) - 69137.1165280576_4221 + 123509.3125610854_8141 %i
Type:
Complex(Float)

In fact, the value xx is one of radicalRoots(pp) where

pp :=
x^4-2729960418308000*x^3-395258439243352250000*x^2-55499520947716391500000000*x-345363656226658026765625000000

Interestingly, when I put xx into Mathematica, I get a much nicer
expressions.

In[15]:= p1 = Root[pp, 1] // ToRadicals

Out[15]= 250 (2729960418308 + 1930373524352 Sqrt[2] -
23569 Sqrt[2 (13416226688183641 + 9486704869150589 Sqrt[2])])

In[25]:= p3 = Root[pp, 3] // ToRadicals

Out[25]= 250 (2729960418308 - 1930373524352 Sqrt[2] -
23569 I Sqrt[2 (-13416226688183641 + 9486704869150589 Sqrt[2])])

Can I somehow "convince" FriCAS to return similarly "simple" radical
expresssions?

Thank you
Ralf

Ralf Hemmecke

unread,
May 14, 2024, 9:42:21 AM5/14/24
to fricas...@googlegroups.com
Trying

pp :=
x^4-2729960418308000*x^3-395258439243352250000*x^2-55499520947716391500000000*x-345363656226658026765625000000
xxx := radicalRoots(pp,'x)$RadicalSolvePackage(INT)
x4 := xxx.4
)set mess bot on

in a fresh session gives

%%% (14) -> imag(x4)

Function Selection for imag
Arguments: EXPR(INT)

[1] signature: EXPR(INT) -> EXPR(INT)
implemented: slot (Expression (Integer))(Expression (Integer))
from TRIGMNIP(INT,EXPR(INT))
[2] signature: EXPR(COMPLEX(INT)) -> EXPR(INT)
implemented: slot (Expression (Integer))(Expression (Complex
(Integer))) from CTRIGMNP(INT,EXPR(COMPLEX(INT)))


(14) 0
Type:
Expression(Integer)

Interestingly, in the session where I originaly started with, I now get
... (see below). Isn't that weird?

Unfortunately, I cannot reproduce what actually caused that different
behaviour for the same thing.

Oh... the modifying behaviour comes from

%%% (20) -> digits 1000

After that, I also get a non-zero expression for

imag(x4)

Very weird.

Ralf

%%% (476) -> imag(x4)

Function Selection for imag
Arguments: EXPR(INT)

[1] signature: EXPR(INT) -> EXPR(INT)
implemented: slot (Expression (Integer))(Expression (Integer))
from TRIGMNIP(INT,EXPR(INT))
[2] signature: EXPR(COMPLEX(INT)) -> EXPR(INT)
implemented: slot (Expression (Integer))(Expression (Complex
(Integer))) from CTRIGMNP(INT,EXPR(COMPLEX(INT)))


(476)
ROOT
9
*
ROOT

17298950126166306420162339540714414395414132884375_
0000000000000
*
+-+
\|3
+

4115957937119953121577632743454820763125572620312500_
00000000000
/
+-+
27 \|3
,
3
^
2
+
-
33537077489620857808043115000000
*
ROOT

1729895012616630642016233954071441439541413288437_
50000000000000
*
+-+
\|3
+

411595793711995312157763274345482076312557262031250_
000000000000
/
+-+
27 \|3
,
3
+
- 298305252478017928101288965487125000000000
*
ROOT
9
*
ROOT

172989501261663064201623395407144143954141328_
843750000000000000
*
+-+
\|3
+

41159579371199531215776327434548207631255726203_
1250000000000000
/
+-+
27 \|3
,
3
^
2
+
16768538744810428904021557500000
*
ROOT

17298950126166306420162339540714414395414132884_
3750000000000000
*
+-+
\|3
+

4115957937119953121577632743454820763125572620312_
50000000000000
/
+-+
27 \|3
,
3
+
- 298305252478017928101288965487125000000000
/
9
*
ROOT

1729895012616630642016233954071441439541413288437_
50000000000000
*
+-+
\|3
+

411595793711995312157763274345482076312557262031250_
000000000000
/
+-+
27 \|3
,
3
+
45777447049433703438392223807463148736000000000
*
ROOT

17298950126166306420162339540714414395414132884375000000_
0000000
*
+-+
\|3
+

4115957937119953121577632743454820763125572620312500000000_
00000
/
+-+
27 \|3
,
3
/
ROOT

1729895012616630642016233954071441439541413288437500000000_
00000
*
+-+
\|3
+

411595793711995312157763274345482076312557262031250000000000000
/
+-+
27 \|3
,
3
*
ROOT
9
*
ROOT

17298950126166306420162339540714414395414132884_
3750000000000000
*
+-+
\|3
+

4115957937119953121577632743454820763125572620312_
50000000000000
/
+-+
27 \|3
,
3
^
2
+
16768538744810428904021557500000
*
ROOT

1729895012616630642016233954071441439541413288437_
50000000000000
*
+-+
\|3
+

411595793711995312157763274345482076312557262031250_
000000000000
/
+-+
27 \|3
,
3
+
- 298305252478017928101288965487125000000000
/
9
*
ROOT

172989501261663064201623395407144143954141328843750_
000000000000
*
+-+
\|3
+

41159579371199531215776327434548207631255726203125000_
0000000000
/
+-+
27 \|3
,
3
/
6
Type:
Expression(Integer)

Ralf Hemmecke

unread,
May 14, 2024, 10:40:33 AM5/14/24
to fricas...@googlegroups.com
Some more experiments with the imaginary part...

Can someone reproduce/explain the strange behaviour of the attached script?

Please modify the path to efstruc.spad accordingly.

Ralf
foo.input

Qian Yun

unread,
May 14, 2024, 7:22:05 PM5/14/24
to fricas...@googlegroups.com
One possible explanation:

)compile a spad file clears internal Kernel cache.

So the following computation can have messed up Kernel order.

So after compiling a spad file, do not use kernels created before that.

- Qian

Ralf Hemmecke

unread,
May 15, 2024, 1:19:11 AM5/15/24
to fricas...@googlegroups.com
On 5/15/24 01:21, Qian Yun wrote:
> One possible explanation:
>
> )compile a spad file clears internal Kernel cache.
>
> So the following computation can have messed up Kernel order.
>
> So after compiling a spad file, do not use kernels created before that.


OK, that explains some of the oddity, but not why "digits 20" makes it
to fail and "digits 21" sometimes helps to make it work and sometimes not.

Ralf

Waldek Hebisch

unread,
May 15, 2024, 5:42:49 AM5/15/24
to fricas...@googlegroups.com
TrigonometricManipulations (which implements imag) is using
sign$ElementaryFunctionSign to determine if a root is complex.
'sign' may return "failed" and 'imag' treats this the same as
nonnegative, so assumes that root is real. 'sign' for numbers
tries interval approximation. If precision it too low, then
this fails and 'sign' returns "failed"...

Mathematically determining sign of numbers (which is needed to choose
form of roots) requires numerical computation potentially of much
higher precision than input data. So there will be a cutoff or
we would run out of memory (or time) in some cases. So one question
is how hard FriCAS should try, currently FriCAS just uses current
numeric preciion and do not try to escalate. Another question
is what to do in case of failure. Honestly FriCAS should return
"failed" or signal error in such case. 'sign' returns "failed"
which IMO is OK. Probably 'imag' should signal error. I can
guess that possibly users call 'imag' and related functions
without real need, so lying (as done now) could be a pragmatic
choice. Or simply this is piece of research code that never
really get finished.

--
Waldek Hebisch

Ralf Hemmecke

unread,
May 15, 2024, 5:58:22 AM5/15/24
to fricas...@googlegroups.com
Thank you, Waldek,

that explains the behaviour much better. I suspected that there must be
some numerics involved, but when I had put debugging output near places
where I saw Float-s, that never appeared, so I got a bit confused and
didn't want to waste more time on this.

Potentially, I would be in favour of letting imag or real fail if the
sign cannot be determined correctly. In that case the error message
should contain a hint that it might be helpful to increase the Float digits.

On the other hand, it looks like FriCAS is anyway spitting out radical
expressions that are more complicated than needed. I tried rsimp on them
but that failed, unfortunately. Seemingly Mathematica is trying harder
in ToRadicals.

https://reference.wolfram.com/language/ref/ToRadicals.html

Side remark. I am actually not so happy with the name rsimp or rootSimp.
Cannot FriCAS stick to the convention that we have full words and no
abbreviations in function names?

Ralf

Qian Yun

unread,
May 15, 2024, 6:47:52 AM5/15/24
to fricas...@googlegroups.com
You said before AlgebraicNumber should represent any branch
of the root, but here we are talking about the sign of radical.

Does that mean we are taking a specific branch interpretation here?

- Qian

Waldek Hebisch

unread,
May 15, 2024, 7:12:58 AM5/15/24
to fricas...@googlegroups.com
On Wed, May 15, 2024 at 06:47:45PM +0800, Qian Yun wrote:
> You said before AlgebraicNumber should represent any branch
> of the root, but here we are talking about the sign of radical.
>
> Does that mean we are taking a specific branch interpretation here?

Yes. Othewise there is problem with meaning of 'imag'. Also
note that this is in extra package, not in core domains.

--
Waldek Hebisch

Qian Yun

unread,
May 15, 2024, 7:28:54 AM5/15/24
to fricas...@googlegroups.com
OK, so for radicals, there should be a much more effective
method to determine its sign. Seems like current method is
mostly focused on polynomials and pretty limited for radicals.

- Qian

Waldek Hebisch

unread,
May 15, 2024, 7:34:45 AM5/15/24
to fricas...@googlegroups.com
On Wed, May 15, 2024 at 07:28:50PM +0800, Qian Yun wrote:
> OK, so for radicals, there should be a much more effective
> method to determine its sign. Seems like current method is
> mostly focused on polynomials and pretty limited for radicals.

We are looking at sign of argument to a radical. This argument
can be a general expression, so we 'sign' needs to handle
general expressions. Possibly we could add special case in
'sign' to handle expressions in radicals, I am not sure how
much this would help.

--
Waldek Hebisch

Waldek Hebisch

unread,
May 16, 2024, 9:12:24 AM5/16/24
to fricas...@googlegroups.com
On Tue, May 14, 2024 at 03:04:15PM +0200, Ralf Hemmecke wrote:
>
> In fact, the value xx is one of radicalRoots(pp) where
>
> pp := x^4-2729960418308000*x^3-395258439243352250000*x^2-55499520947716391500000000*x-345363656226658026765625000000
>
> Interestingly, when I put xx into Mathematica, I get a much nicer
> expressions.
>
> In[15]:= p1 = Root[pp, 1] // ToRadicals
>
> Out[15]= 250 (2729960418308 + 1930373524352 Sqrt[2] -
> 23569 Sqrt[2 (13416226688183641 + 9486704869150589 Sqrt[2])])
>
> In[25]:= p3 = Root[pp, 3] // ToRadicals
>
> Out[25]= 250 (2729960418308 - 1930373524352 Sqrt[2] -
> 23569 I Sqrt[2 (-13416226688183641 + 9486704869150589 Sqrt[2])])
>
> Can I somehow "convince" FriCAS to return similarly "simple" radical
> expresssions?

With the attached patch I get:

(2) -> radical_solve(univariate(pp))

(2)
[
+--------------------------------------------------------------------+
| +-+
\|658730414260118770166403625000 \|2 + 931585485794307216540975125000
+
+-+
482593381088000 \|2 + 682490104577000
,

+--------------------------------------------------------------------+
| +-+
- \|658730414260118770166403625000 \|2 + 931585485794307216540975125000
+
+-+
482593381088000 \|2 + 682490104577000
,

+---+
\|- 1
*
+--------------------------------------------------------------------+
| +-+
\|658730414260118770166403625000 \|2 - 931585485794307216540975125000
+
+-+
- 482593381088000 \|2 + 682490104577000
,

-
+---+
\|- 1
*
ROOT
+-+
658730414260118770166403625000 \|2
+
- 931585485794307216540975125000
+
+-+
- 482593381088000 \|2 + 682490104577000
]
Type: Union(List(Expression(Integer)),...)


The result looks more complicated, but this is because FriCAS does
not pull squares outside square roots.

I would be easy to hook this into 'radicalSolve' so that such result
is obtainde by default.

--

Waldek Hebisch
ROOTUT.diff

Ralf Hemmecke

unread,
May 16, 2024, 12:44:15 PM5/16/24
to fricas...@googlegroups.com
Waldek,

you are amazing!!! A big thanks for the radicalSolve patch.

Ralf
Reply all
Reply to author
Forward
0 new messages