prime? of UP

1 view
Skip to first unread message

Serge D. Mechveliani

unread,
Feb 28, 2012, 12:24:27 PM2/28/12
to fricas...@googlegroups.com
For f := (x^2 - 1) :: UP(x,INT)

factor f
yields
(12) (x - 1)(x + 1)
Type: Factored(UnivariatePolynomial(x,Integer)),

and prime?(f) $UP(x,INT)

reports in the interpreter that no operation named prime? fits here.
And in Spad, the compiler reports
...
****** level 3 ******
$x:= ((elt (UnivariatePolynomial x (Integer)) prime?) f)
$m:= $EmptyMode
$f:=
((((|f| #) (|fA| #) (|fInput| #) (|expInput| #) ...)
((|sumInts| #) (* #) (+ #) (- #) ...)))

>> Apparent user error:
not known that (UnivariatePolynomial x (Integer)) has
(has (Integer) (PolynomialFactorizationExplicit))

It also does not know `irreducible?', `reducible?'.
As it can factor in UP(x,INT) as in an UniqueFactorizationDomain,
why cannot it find there value of prime?(f) ?
As I recall, prime? is of UniqueFactorizationDomain.

Regards,

------
Sergei
mec...@botik.ru


Ralf Hemmecke

unread,
Feb 28, 2012, 4:10:24 PM2/28/12
to fricas...@googlegroups.com
On 02/28/2012 06:24 PM, Serge D. Mechveliani wrote:
> For f := (x^2 - 1) :: UP(x,INT)
>
> factor f
> yields
> (12) (x - 1)(x + 1)
> Type: Factored(UnivariatePolynomial(x,Integer)),
>
> and prime?(f) $UP(x,INT)

That's funny. HyperDoc (and the code) says that factor only exists if
the coefficient domain has PolynomialFactorizationExplicit. One can
easily check that *no* domain in FriCAS exports this category
unconditionally. Strange that FriCAS finds a factorization routine at all.

That the compiler cannot factor is then no surprise, but how to find out
what the interpreter actually did?

Aha!

(1) -> )set message selection on
(1) -> factor f

Function Selection for factor
Arguments: UP(x,INT)

[1] signature: UP(x,INT) -> FR(UP(x,INT))
implemented: slot (Factored (UnivariatePolynomial x
(Integer)))(UnivariatePolynomial x (Integer)) from
MULTFACT(SAOS,NNI,INT,UP(x,INT))
[2] signature: SUP(UP(x,INT)) -> FR(SUP(UP(x,INT)))
implemented: slot (Factored (SparseUnivariatePolynomial
(UnivariatePolynomial x (Integer))))(SparseUnivariatePolynomial
(UnivariatePolynomial x (Integer))) from MULTFACT(SAOS,NNI,INT,UP(x,INT))


(1) (x - 1)(x + 1)
Type:
Factored(UnivariatePolynomial(x,Integer))

Looking at MULTFACT reveals that this only applies to very special
situation.

https://github.com/hemmecke/fricas-svn/blob/master/src/algebra/multfact.spad.pamphlet#L465
(copy attached)

It's quite interesting how much special code is hidden in the
interpreter so that a function of this package can be called. :-(
Lot's of work ahead.

Ralf

==============================================================
MultivariateFactorize(OV,E,R,P) : C == T
where
R : Join(EuclideanDomain, CharacteristicZero)
-- with factor on R[x]
OV : OrderedSet
E : OrderedAbelianMonoidSup
P : PolynomialCategory(R,E,OV)
Z ==> Integer
MParFact ==> Record(irr:P,pow:Z)
USP ==> SparseUnivariatePolynomial P
SUParFact ==> Record(irr:USP,pow:Z)
SUPFinalFact ==> Record(contp:R,factors:List SUParFact)
MFinalFact ==> Record(contp:R,factors:List MParFact)

-- contp = content,
-- factors = List of irreducible factors with exponent
L ==> List

C == with
factor : P -> Factored P
++ factor(p) factors the multivariate polynomial p over its
coefficient
++ domain
factor : USP -> Factored USP
++ factor(p) factors the multivariate polynomial p over its
coefficient
++ domain where p is represented as a univariate polynomial with
++ multivariate coefficients
T == add
factor(p:P) : Factored P ==
R is Fraction Integer =>
factor(p)$MRationalFactorize(E,OV,Integer,P)
R is Fraction Complex Integer =>
factor(p)$MRationalFactorize(E,OV,Complex Integer,P)
R is Fraction Polynomial Integer and OV has convert: % -> Symbol =>
factor(p)$MPolyCatRationalFunctionFactorizer(E,OV,Integer,P)
factor(p,factor$GenUFactorize(R))$InnerMultFact(OV,E,R,P)

factor(up:USP) : Factored USP ==
factor(up,factor$GenUFactorize(R))$InnerMultFact(OV,E,R,P)

Serge D. Mechveliani

unread,
Feb 29, 2012, 5:13:47 AM2/29/12
to fricas...@googlegroups.com
On Tue, Feb 28, 2012 at 10:10:24PM +0100, Ralf Hemmecke wrote:
> On 02/28/2012 06:24 PM, Serge D. Mechveliani wrote:
>> For f := (x^2 - 1) :: UP(x,INT)
>>
>> factor f
>> yields
>> (12) (x - 1)(x + 1)
>> Type: Factored(UnivariatePolynomial(x,Integer)),
>>
>> and prime?(f) $UP(x,INT)
>
> That's funny. HyperDoc (and the code) says that factor only exists if
> the coefficient domain has PolynomialFactorizationExplicit. One can
> easily check that *no* domain in FriCAS exports this category
> unconditionally. Strange that FriCAS finds a factorization routine at
> all.

This is somow against classical algorithmic algebra.
1) Classical algebra says that
Z[x] is an unique factorization domain with an algorithmic factorization.
2) Axiom announces that it has polynomial factorization, over Integer too.
3) UniqueFactorizationDomain has `factor' and `prime?' among its operations.

Hence it is natural for Axiom to provide `factor' and `prime?' for Z[x]
for all its representations given in Axiom for UP, SUP, POLY, etc.
Because as it is done for POLY(..INT..), it is easy to spread them also
to UP(x,INT) too.
Probably, it does this, but it is difficult to find.


> Looking at MULTFACT reveals that this only applies to very special
> situation.

> [..]

> ==============================================================
> MultivariateFactorize(OV,E,R,P) : C == T
> where
> R : Join(EuclideanDomain, CharacteristicZero)
> -- with factor on R[x]
> OV : OrderedSet
> E : OrderedAbelianMonoidSup
> P : PolynomialCategory(R,E,OV)
>

> [..]


Thank you.
This looks like some hint.
But to tell it simply: how to factor in UP(x,INT) ?

There are many people who know that a polynomial f with integer
coefficients can be factored and do not know the meaning of the scripts
like


> R : Join(EuclideanDomain, CharacteristicZero)
> -- with factor on R[x]
> OV : OrderedSet
> E : OrderedAbelianMonoidSup

They just have f := x^3 -1 : UP(x,INT),
and nothing else.
Now: do they need to call in Spad
factor(f) $ MultivariateFactorize(OV,E,R,P)

and to cleverly decide what to put for OV, E and P ?
Also may be, to consider convertion from INT to Fraction INT.
This invesigation will, probably, take a week.

I suspect that there exists a simple way for this in Spad.
I am looking into
"Browse -- factor -- Operations -- Origins"

(in the same time as e-mailing this question)
and trying to choose ...
May be, GenUFactorize(INT)
fits?
I try in Spad:
f() : Factored SUP ==
X := monomial(1 ::INT, 1 ::NNI) $UPol
y := X pretend SUP
f := (4*y^2 - 1) :: SUP
factor(f) $GenUFactorize(INT)
It works!
Is this a reasonable approach?

And what about `prime?' ? We need to search a package ...

Generally, it is very difficult to find the needed package and also
often difficult to put a correct parameter values.

May be, there is possible import AllFactorization
?
(in order to apply `factor' without specification,
and let the compiler give a warning:
such and such variants for `factor' fit this particular piece of code).

And there is something strange in this total necessity for the user to
search cleverly among tons of packages.
The compiler knows itself that SUP INT has UniqueFactorizationDomain,
and that `factor' and `prime?' are operations from
UniqueFactorizationDomain. The user program needs not to mention this
knowledge. Hence, for
factor(f :: SUP INT) :: Factored SUP INT

it must know itself what `factor' signature variant fits this piece of
code. Often such is unique. For example, the compiler knows that `factor'
of @Boolean does not fit the above expression. In this case, the thing
is compiled.
If more than one implementation signatures fit, the compiler could give
a warning: choose among implementations in these packages ...

In Haskell + glasgow-ext (DoCon), it is much simpler.
In most cases the choice is unique. In other cases (of several class
instances that overlap at this type of f) the compiler needs to choose
an instance of an operation (say, `factor') among a finite set of intances
for this particular type of the argument, say, f :: UPol R.
The compiler chooses the instance the most specific by the type parameter
substitution. If such does not exist, reports the error.
For example,
If Integer has Foo, then the type UPol Integer
is obtained by a substitution from Foo rR => UPol rR.
Hence the instance of `factor' is taken from the implementatiion given
for UPol Integer, it is more specific than the genetal one of
Foo rR => UPol rR, and it overrides. Such is an agreement for the
language (this can be generalized).
And it occurs that the user needs only to specify the type of f, and
the compiler finds all the rest
-- with respect to the program _modules_ imported by the compiled module
(modules are like packages in Spad, only without parameters).

The instances for the same type may be given in many modules.
There is possible nested import. And DoCon allows

import DExport -- import all open items DoCon --
-- all items exported from all DoCon modules.

Axiom is larger, so this may be several large blocks, say,
import AllFactorization AllPlottingGraphis

Thanks,
------
Sergei
mec...@botik.ru


Waldek Hebisch

unread,
Feb 29, 2012, 1:04:48 PM2/29/12
to fricas...@googlegroups.com
Serge D. Mechveliani wrote:
>
> On Tue, Feb 28, 2012 at 10:10:24PM +0100, Ralf Hemmecke wrote:
> > On 02/28/2012 06:24 PM, Serge D. Mechveliani wrote:
> >> For f := (x^2 - 1) :: UP(x,INT)
> >>
> >> factor f
> >> yields
> >> (12) (x - 1)(x + 1)
> >> Type: Factored(UnivariatePolynomial(x,Integer)),
> >>
> >> and prime?(f) $UP(x,INT)
> >
> > That's funny. HyperDoc (and the code) says that factor only exists if
> > the coefficient domain has PolynomialFactorizationExplicit. One can
> > easily check that *no* domain in FriCAS exports this category
> > unconditionally. Strange that FriCAS finds a factorization routine at
> > all.
>
> This is somow against classical algorithmic algebra.
> 1) Classical algebra says that
> Z[x] is an unique factorization domain with an algorithmic factorization.
> 2) Axiom announces that it has polynomial factorization, over Integer too.
> 3) UniqueFactorizationDomain has `factor' and `prime?' among its operations.

Actually situation is more tricky than you think: it is not enough
to know that Z[x] is an unique factorization domain, but you
also need an algorithm to factorize in Z[x]. Typicaly Axiom/Fricas
builds algorithm in Z[x] using algorithm in Z as building block
(where Z now stands for arbitrary Ring). However, IIRC there is
a theorem that essentially says that it is impossible to build
factorization algorithm in R[x] just from arithmetic operations
and factorization algorithm in R. This is quite different then
say gcd, where fraction free version of Euclidean algorithm
uses only arithmetic in R and together with gcd in R is enough
to get gcd in R[x]. So you need more than unique factorization
domain to have factorization. Theoretically Axiom developers
solved the problem -- they invented stronger property namely
PolynomialFactorizationExplicit which correctly propagetes trough
towers of extensions. And they wrote code to handle it.
So, why we do not use this code? The problem is performance,
in cases handled by current code the general code is a lot
slower.

So, currently we have code which handles several special
cases. Mathematically those cases cover the same area
as "general" code. The difference is that currently
factor function has to solve problem similar to what you
(Sergei) is doing in Haskell interface: analyze type and
dispatch to correct factorization routine. This dispatch
only works for some types and some which could work are
not handled.


> Hence it is natural for Axiom to provide `factor' and `prime?' for Z[x]
> for all its representations given in Axiom for UP, SUP, POLY, etc.
> Because as it is done for POLY(..INT..), it is easy to spread them also
> to UP(x,INT) too.
> Probably, it does this, but it is difficult to find.

Orignal Axiom developers wanted to treat 'factor' and `prime?' via
PolynomialFactorizationExplicit. But this is not practical
for performance reasons. For specifice domains, like Z[x] we
could solve the problem in ad hoc manner. But I am not sure
if possible gain is worth the effort: interpreter finds 'factor'
when needed and in Spad one freqently uses more general domains
and need to call factor from appropriate package (possibly
converting polynomial to supported type).

>
> > Looking at MULTFACT reveals that this only applies to very special
> > situation.
> > [..]
> > ==============================================================
> > MultivariateFactorize(OV,E,R,P) : C == T
> > where
> > R : Join(EuclideanDomain, CharacteristicZero)
> > -- with factor on R[x]
> > OV : OrderedSet
> > E : OrderedAbelianMonoidSup
> > P : PolynomialCategory(R,E,OV)
> >
> > [..]
>
>
> Thank you.
> This looks like some hint.
> But to tell it simply: how to factor in UP(x,INT) ?

<snip>


> I suspect that there exists a simple way for this in Spad.
> I am looking into
> "Browse -- factor -- Operations -- Origins"
>
> (in the same time as e-mailing this question)
> and trying to choose ...
> May be, GenUFactorize(INT)
> fits?
> I try in Spad:
> f() : Factored SUP ==
> X := monomial(1 ::INT, 1 ::NNI) $UPol
> y := X pretend SUP
> f := (4*y^2 - 1) :: SUP
> factor(f) $GenUFactorize(INT)
> It works!
> Is this a reasonable approach?

Yes. For factorization over integers the actual factor routine
is in 'UnivatiateFactorize'. GenUFactorize analyses its argument
and dispatches to UnivatiateFactorize or another specific
package.

> And what about `prime?' ? We need to search a package ...

Sorry, unimplemented. But implementation is trival: first
use 'factor' and then check the result for a non-trivial
factors.

> Generally, it is very difficult to find the needed package and also
> often difficult to put a correct parameter values.
>
> May be, there is possible import AllFactorization
> ?
> (in order to apply `factor' without specification,
> and let the compiler give a warning:
> such and such variants for `factor' fit this particular piece of code).
>
> And there is something strange in this total necessity for the user to
> search cleverly among tons of packages.
> The compiler knows itself that SUP INT has UniqueFactorizationDomain,
> and that `factor' and `prime?' are operations from
> UniqueFactorizationDomain.

No. Currently SUP INT does not have UniqueFactorizationDomain.
And we can not make it into UniqueFactorizationDomain without
adding factor to it, because then interpreter would ignore
factor from factorization packages and try to use factor from
UniqueFactorizationDomain (and fail calling non-existent function).
And I wrote, for Integer it would be small amount of work to
add a special case. But then we probably would like to handle
Fraction(Integer) and then tons of other cases. So simple
special casing does not scale.

I know that current situation with 'factor' is not nice and
I am working on solution. But I believe that we need general
solution (like PolynomialFactorizationExplicit, but without
performance problems) not patching to handle few special
cases. General solution takes time to develop, and as long
as it is not finished it does not work at all.


--
Waldek Hebisch
heb...@math.uni.wroc.pl

Serge D. Mechveliani

unread,
Mar 1, 2012, 1:51:15 PM3/1/12
to fricas...@googlegroups.com
On Wed, Feb 29, 2012 at 07:04:48PM +0100, Waldek Hebisch wrote:
> Serge D. Mechveliani wrote:
> >
> > On Tue, Feb 28, 2012 at 10:10:24PM +0100, Ralf Hemmecke wrote:
> > > On 02/28/2012 06:24 PM, Serge D. Mechveliani wrote:
> > >> For f := (x^2 - 1) :: UP(x,INT)
> > >>
> > >> factor f
> > >> yields
> > >> (12) (x - 1)(x + 1)
> > >> Type: Factored(UnivariatePolynomial(x,Integer)),
> > >>
> > >> and prime?(f) $UP(x,INT)
> > >
> > > That's funny. HyperDoc (and the code) says that factor only exists if
> > > the coefficient domain has PolynomialFactorizationExplicit. One can
> > > easily check that *no* domain in FriCAS exports this category
> > > unconditionally. Strange that FriCAS finds a factorization routine at
> > > all.
> >
> > This is somow against classical algorithmic algebra.
> > 1) Classical algebra says that
> > Z[x] is an unique factorization domain with an algorithmic factorization.
> > 2) Axiom announces that it has polynomial factorization, over Integer too.
> > 3) UniqueFactorizationDomain has `factor' and `prime?' among its operations.
>
> Actually situation is more tricky than you think: it is not enough
> to know that Z[x] is an unique factorization domain, but you
> also need an algorithm to factorize in Z[x].

For example, you write CommutativeRing R
and presume that it is algorithmic, and do not recall of
"also need an algorithm" for arithmetics.

Also Browse -- UniqueFactorizationDomain -- Categories shows:

"A constructive unique factorization domain. ...
can constructively factor members into ...
"
-- which means exactly AlgorithmicalUniqueFactorizationDomain.

This doc (and my DoCon too) removes the prefix Algorithmical (or Constructive)
in this case because
1) in practice, this is a rare case to have R : UniqueFactorizationDomain
which is not constructive,
2) it is better to solve dynamically, by extracting an attribute from a
domain, whether it has `Constructive'.

In the library, this question can be solved simply:
hasFactorizationAlgorithm R = true
means that `factor' is really a function representing a true factorization
in R, and in this case R has ConstructiveUniqueFactorizationDomain
(in this particular library).
hasFactorizationAlgorithm R just extracts the attribute value from the
domain description.
This is not sufficient to apply factor f, one needs to check
hasFactorizationAlgorithm R before this.
`false' means that `factor' may not present a correct factorization.
This all must be explained in doc,
also in the 'false' case, it may be set factor f = error "...",
for more safety.

For example, suppose Axiom has a correct factorOverInteger for
R = UP(x,INT). So, UP(x,INT) : UniqueFactorizationDomain,
and the HasFactorizationAlgorithm attribute is put true.
And in the implementation in this domain R it is put
factor(f : UP(x,INT)) == factorOverInteger f,

and this is the fastest kwnown special method to factor in UP(x,INT).
Then, one implements the method to factor UP over any
D : SuchAndSuchCleverRing,
and INT : SuchAndSuchCleverRing.
In the implementation for UP(x,D) it must be set

factor f == if D has IntegerLike then factorOverInteger f
else genericFactorUP(D, f)

General algorithms overlap with special ones, and a more special overrides.
And this must be expressed in the Spad code in the library.
And I expect that in Axiom it is arranged a similar approach, may be, by
some different means.

Has Axiom a category like AlgorithmicalUniqueFactorizationDomain ?
Is it PolynomialFactorizationExplicit ?

> Typicaly Axiom/Fricas builds algorithm in Z[x] using algorithm in Z as
> building block
> (where Z now stands for arbitrary Ring). However, IIRC there is
> a theorem that essentially says that it is impossible to build
> factorization algorithm in R[x] just from arithmetic operations
> and factorization algorithm in R.

I am aware of this.
And this should not reject the definition of UniqueFactorizationDomain given
in Browse (see above) -- for the reasons I write above.

> Theoretically Axiom developers solved the problem -- they invented
> stronger property namely PolynomialFactorizationExplicit which correctly
> propagetes trough towers of extensions. And they wrote code to handle it.

Lenstra, Grigoriev, Gianni, Davenport, Trinks, and others (and some of them
even wrote a code for Axiom), developed the methods and algorithms
(not necesssarily programs) which, in particular, describe a wide class of
commutative rings R for which it is built a correct factorization algorithm
in R[x1..xn]. I know of this for a long time.
I expect that the Axiom library somehow follows the main stream of this whole
theory.

> So, why we do not use this code? The problem is performance,
> in cases handled by current code the general code is a lot slower.

??
The method of resovling an overlap of a special algorithm and a general one
may be, for example, as I describe above.
Axiom may apply a different method of resolving, but I believe that the
Spad language is expressive enough to have a good resolution method for this.
Why cannot Axiom arrange this clearly?


> > Hence it is natural for Axiom to provide `factor' and `prime?' for Z[x]
> > for all its representations given in Axiom for UP, SUP, POLY, etc.
> > Because as it is done for POLY(..INT..), it is easy to spread them also
> > to UP(x,INT) too.
> > Probably, it does this, but it is difficult to find.
>
> Orignal Axiom developers wanted to treat 'factor' and `prime?' via
> PolynomialFactorizationExplicit. But this is not practical
> for performance reasons.

> [..]

> For specifice domains, like Z[x] we
> could solve the problem in ad hoc manner. But I am not sure
> if possible gain is worth the effort: interpreter finds 'factor'
> when needed and in Spad one freqently uses more general domains
> and need to call factor from appropriate package (possibly
> converting polynomial to supported type).

I do not understand why Axiom cannot efficiently resolve this overlap.
Probably, it can, but I cannot find the decription of this dispatch method.
Currently, I even do not need to observe it. It is only desirable to have a
clear rule for an user of how to choose the most efficient `factor'
implementation.

Example: for determinant of a SquareMatrix over rR, DoCon has

(1) det for CommutativeRing rR => SquareMatrix rR,

(2) det_euc for Euclidean rR => SquareMatrix rR,

(3) det_upol_finField ::
Field k => [ResidueE (UPol k)] -> [[UPol k]] -> UPol k
-- ext mM
-- A degree cost method to compute det M over k[x] for a finite
-- field k.

On average, (2) is faster than (1), (3) is faster than (2).
DoCon provides three different functions, for some reasons.
If it introduced class WithDet, may be, it could have a unique operation
`det', whose instances are efficienly resolved for these three functions
(I am not sure, needs recollection, trying).

But the Spad language is more expressive, and it is likely to have an easier
resolution.

Regards,

------
Sergei
mec...@botik.ru

Reply all
Reply to author
Forward
0 new messages