oscar-system / Singular.jl

Julia package for the Singular library
Other
31 stars 33 forks source link

Making `Singular.LibTropical.tropicalVariety` usable #807

Open YueRen opened 4 months ago

YueRen commented 4 months ago

I am trying to make Singular.LibTropical.tropicalVariety usable for small examples. Can anybody point me in the right direction where to start?

In Singular, the function takes as input an ideal and an optional number and returns a fan:

                     SINGULAR                                 /  Development
 A Computer Algebra System for Polynomial Computations       /   version 4.4.0
                                                           0<
 by: W. Decker, G.-M. Greuel, G. Pfister, H. Schoenemann     \   Apr 2024
FB Mathematik der Universitaet, D-67653 Kaiserslautern        \
> LIB "tropical.lib";
[...]
> ring r = 0,(a,b,c,d,e,f,g,h,i,j),dp;
> ideal I = bf-ah-ce, bg-ai-de, cg-aj-df, ci-bj-dh, fi-ej-gh;
> tropicalVariety(I); # function in tropical.lib
_application PolyhedralFan
_version 2.2
_type PolyhedralFan

AMBIENT_DIM
10

DIM
7

LINEALITY_DIM
5

RAYS
-3 1 1 1 1 1 1 -1 -1 -1 # 0
-1 -1 1 1 -1 1 1 1 1 -3 # 1
-1 1 -1 1 1 -1 1 1 -3 1 # 2
-1 1 1 -1 1 1 -1 -3 1 1 # 3
1 -3 1 1 1 -1 -1 1 1 -1 # 4
1 -1 -1 1 1 1 -3 -1 1 1 # 5
1 -1 1 -1 1 -3 1 1 -1 1 # 6
1 1 -3 1 -1 1 -1 1 -1 1 # 7
1 1 -1 -1 -3 1 1 1 1 -1 # 8
1 1 1 -3 -1 -1 1 -1 1 1 # 9

[...]

MAXIMAL_CONES
{0 1}   # Dimension 7
{0 2}
{0 3}
{1 4}
{1 8}
{2 6}
{3 5}
{2 7}
{3 9}
{4 5}
{4 6}
{5 7}
{6 9}
{7 8}
{8 9}

Two things need to be done:

  1. make Singular.LibTropical.tropicalVariety callable from julia
  2. make OSCAR understand the return value

Currently I am still stuck at 1:

julia> R,(a,b,c,d,e,f,g,h,i,j) = QQ["a","b","c","d","e","f","g","h","i","j"];

julia> I = ideal([b*f-a*h-c*e, b*g-a*i-d*e, c*g-a*j-d*f, c*i-b*j-d*h, f*i-e*j-g*h]);

julia> Singular.LibTropical.tropicalVariety(I)
ERROR: unrecognized argument Ideal (-a*h + b*f - c*e, -a*i + b*g - d*e, -a*j + c*g - d*f, -b*j + c*i - d*h, -e*j + f*i - g*h)
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] prepare_argument(x::MPolyIdeal{QQMPolyRingElem})
   @ Singular ~/.julia/packages/Singular/C3EMh/src/caller.jl:358
 [3] (::Singular.var"#226#227")(i::MPolyIdeal{QQMPolyRingElem})
   @ Singular ./none:0
 [4] iterate
   @ ./generator.jl:47 [inlined]
 [5] collect(itr::Base.Generator{Tuple{MPolyIdeal{QQMPolyRingElem}}, Singular.var"#226#227"})
   @ Base ./array.jl:834
 [6] low_level_caller(lib::String, name::String, args::MPolyIdeal{QQMPolyRingElem})
   @ Singular ~/.julia/packages/Singular/C3EMh/src/caller.jl:389
 [7] tropicalVariety(args::MPolyIdeal{QQMPolyRingElem})
   @ Singular.LibTropical ~/.julia/packages/Singular/C3EMh/src/Meta.jl:44
 [8] top-level scope
   @ REPL[6]:1

Any help would be much appreciated!

YueRen commented 4 months ago

Also, I'm happy to drop the optional number. Just being able to call it with an ideal would already be great!

thofma commented 4 months ago

From slack:

Singular.jl interpreter functionality only understandards Singular.jl types

See also See also https://discourse.julialang.org/t/singular-jl-polynomial-ring-over-a-polynomial-ring/110839/2

YueRen commented 4 months ago

Thanks, that seems to have worked, I know get an expected error in the return:

julia> R,(a,b,c,d,e,f,g,h,i,j) = Singular.polynomial_ring(QQ,["a","b","c","d","e","f","g","h","i","j"]);

julia> I = Singular.Ideal(R,[b*f-a*h-c*e, b*g-a*i-d*e, c*g-a*j-d*f, c*i-b*j-d*h, f*i-e*j-g*h]);

julia> Singular.LibTropical.tropicalVariety(I)
ERROR: unrecognized object with Singular type number 545
Note that if Singular.jl cannot interpret the type, it is doubtful that a interpreter procedure returning such a type can be useful to Singular.jl.
Note also that Singular.jl does not support the attributes used by the interpreter.
Finally, Singular.lookup_library_symbol can be used to fetch the current value of global variables stored in the interpreter.
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] (::Singular.var"#216#217"{Int64})()
   @ Singular ~/.julia/packages/Singular/C3EMh/src/caller.jl:142
 [3] get
   @ ./dict.jl:547 [inlined]
 [4] convert_normal_value(valueptr::Ptr{Nothing}, typ::Int64, R::Singular.PolyRing{Singular.n_FieldElem{…}})
   @ Singular ~/.julia/packages/Singular/C3EMh/src/caller.jl:141
 [5] convert_return(value::Vector{Any}, R::Singular.PolyRing{Singular.n_FieldElem{Singular.FieldElemWrapper{…}}})
   @ Singular ~/.julia/packages/Singular/C3EMh/src/caller.jl:229
 [6] low_level_caller(lib::String, name::String, args::Singular.sideal{Singular.spoly{Singular.n_FieldElem{…}}})
   @ Singular ~/.julia/packages/Singular/C3EMh/src/caller.jl:403
 [7] tropicalVariety(args::Singular.sideal{Singular.spoly{Singular.n_FieldElem{Singular.FieldElemWrapper{…}}}})
   @ Singular.LibTropical ~/.julia/packages/Singular/C3EMh/src/Meta.jl:44
 [8] top-level scope
   @ REPL[13]:1
Some type information was truncated. Use `show(err)` to see complete types.

Is there a way to teach Singular.jl how to deal with this "unrecognized object with Singular type number 545"? I guess cannot just make it into a PolyhedralComplex, since the latter is part of OSCAR and not Singular. But I could ask Singular to return:

That way, I can wrap Singular.LibTropical.tropicalVariety in OSCAR and have the wrapper return a proper polyhedral complex.

thofma commented 4 months ago

Is there a way to teach Singular.jl how to deal with this "unrecognized object with Singular type number 545"?

I don't want to stop you from any work, but this might be highly non-trivial.

I guess cannot just make it into a PolyhedralComplex, since the latter is part of OSCAR and not Singular. But I could ask Singular to return:

  • a matrix with the ray generators
  • a matrix with the lineality space generators
  • a matrix with the maximal cone incidences

That way, I can wrap Singular.LibTropical.tropicalVariety in OSCAR and have the wrapper return a proper polyhedral complex.

Yes, I think you can execute arbitrary things with https://oscar-system.github.io/Singular.jl/dev/caller/#Global-Interpreter-Variables and set the global variables to something that you could then digest on the julia side.

YueRen commented 4 months ago

I just realised that Singular has no function that returns the incidence matrix of a fan, meaning which maximal cones contain which rays (this is essential for reconstructing the fan in OSCAR).

So the easiest option would probably be parsing the string of a Singular fan, which does contain said information:

julia> Singular.call_interpreter("ring r=0,(a,b,c,d,e,f,g,h,i,j),dp; ideal I = b*f-a*h-c*e, b*g-a*i-d*e, c*g-a*j-d*f, c*i-b*j-d*h, f*i-e*j-g*h; fan TropI = tropicalVariety(I);  string TropIString = string(TropI);");

julia> Singular.lookup_library_symbol("Top", "TropIString")
"_application PolyhedralFan\n_version 2.2\n_type PolyhedralFan\n\nAMBIENT_DIM\n10\n\nDIM\n7\n\nLINEALITY_DIM\n5\n\nRAYS\n-3 1 1 1 1 1 1 -1 -1 -1\t# 0\n-1 -1 1 1 -1 1 1 1 1 -3\t# 1\n-1 1 -1 1 1 -1 1 1 -3 1\t# 2\n-1 1 1 -1 1 1 -1 -3 1 1\t# 3\n1 -3 1 1 1 -1 -1 1 1 -1\t# 4\n1 -1 -1 1 1 1 -3 -1 1 1\t# 5\n1 -1 1 -1 1 -3 1 1 -1 1\t# 6\n1 1 -3 1 -1 1 -1 1 -1 1\t# 7\n1 1 -1 -1 -3 1 1 1 1 -1\t# 8\n1 1 1 -3 -1 -1 1 -1 1 1\t# 9\n\nN_RAYS\n10\n\nLINEAL" ⋯ 229 bytes ⋯ " 0 0\t# 2\n-1 0 0 1 1 0 0 0 -1 0\t# 3\n-1 -1 1 1 1 0 0 0 0 -1\t# 4\n\nF_VECTOR\n1 10 15\n\nSIMPLICIAL\n1\n\nPURE\n1\n\nCONES\n{}\t# Dimension 5\n{0}\t# Dimension 6\n{1}\n{2}\n{3}\n{4}\n{5}\n{6}\n{7}\n{8}\n{9}\n{0 1}\t# Dimension 7\n{0 2}\n{0 3}\n{1 4}\n{1 8}\n{2 6}\n{3 5}\n{2 7}\n{3 9}\n{4 5}\n{4 6}\n{5 7}\n{6 9}\n{7 8}\n{8 9}\n\nMAXIMAL_CONES\n{0 1}\t# Dimension 7\n{0 2}\n{0 3}\n{1 4}\n{1 8}\n{2 6}\n{3 5}\n{2 7}\n{3 9}\n{4 5}\n{4 6}\n{5 7}\n{6 9}\n{7 8}\n{8 9}\n"
thofma commented 4 months ago

This sounds horrible

YueRen commented 4 months ago

It's a good band-aid. The only other two alternatives are patching gfanlib or reimplementing the algorithm in Oscar.

fingolfin commented 4 months ago

Why can't you write a little Singular glue code which returns the desired data and then invoke that from OSCAR?

YueRen commented 4 months ago

Because I would have to change gfanlib itself, Singular's interface to gfanlib, and make OSCAR use the latest master version of Singular (the latter of which I don't know how). And this issue is time-critical.

fingolfin commented 4 months ago

Why would you have to change gfanlib? You surely can write a Singular function calling gfanlib and converting the output into a format that Singular.jl can read. Then just load that code and call that function in Oscar.

YueRen commented 4 months ago

Because the class gfan::ZFan (which is what Singular uses for fans) lacks a function for returning its maximal cone incidences. Though I guess you are right that is not really necessary, I can just access the members directly.

I will still have to make some changes to Singular's interface to gfanlib to make that functionality available in the interpreter, and make OSCAR depend on the latest master version of Singular.

fingolfin commented 3 weeks ago

I just rediscovered this issue.

@hannes14 there doesn't seem to be a FAN_CMD or similar. Is there any chance we could determine during run- or compiletime that "545" is the Singular "type" id for fans; and then use that to provide a proper interface for these fan objects?

Looking at the Singular source code, there seems to be a C global fanID which stores this value.

And diving a bit deeper, this seems to be the C++ code which computes a "tropical variety" just given an ideal (with multiple variants of this code in Singular/dyn_modules/gfanlib/tropicalVariety.cc for different kinds of input)

      poly g = I->m[0];
...
        try
        {
          tropicalStrategy currentStrategy(I,currRing);
          std::set<gfan::ZCone> maxCones = tropicalVariety(g,currRing,&currentStrategy);
          res->rtyp = fanID;
          res->data = (char*) toZFan(maxCones,currentStrategy.getExpectedAmbientDimension());
          return FALSE;
        }
        catch (const std::exception& ex)
        {
          Werror("ERROR: %s",ex.what());
          return TRUE;
        }

where

static gfan::ZFan* toZFan(std::set<gfan::ZCone> maxCones, int d)
{
  gfan::ZFan* zf = new gfan::ZFan(d);
  for (std::set<gfan::ZCone>::iterator sigma = maxCones.begin(); sigma!=maxCones.end(); sigma++)
    zf->insert(*sigma);
  return zf;
}

So now I wonder if it would make even more sense to directly expose the gfan::ZFan type via libsingular_julia (resp. gfan::ZCone )

YueRen commented 3 weeks ago

So now I wonder if it would make even more sense to directly expose the gfan::ZFan type via libsingular_julia (resp. gfan::ZCone)

That largely depends on how difficult it is to implement it. I personally don't have any long-term plans on using Singular's tropicalVariety in OSCAR, I'm already working towards a better implementation in OSCAR with @VictoriaSchleis and @APMS04.

fingolfin commented 3 weeks ago

Ah OK! If you plan to rewrite all that that code in the coming months anyway, then I don't think it makes sense to worry about a better way to call it.