diff --git a/M2/Macaulay2/d/ffi.d b/M2/Macaulay2/d/ffi.d index aec240824df..65be1e85765 100644 --- a/M2/Macaulay2/d/ffi.d +++ b/M2/Macaulay2/d/ffi.d @@ -207,13 +207,18 @@ setupconst("ffiVoidType", toExpr(Ccode(voidPointer, "&ffi_type_void"))); -- integer types -- ------------------- +-- returns pointer to ffi_type object for given integer type +-- input: (n:ZZ, signed:Boolean) +-- n = number of bits (or 0 for mpz_t) ffiIntegerType(e:Expr):Expr := ( when e is a:Sequence do ( if length(a) == 2 then ( when a.0 is n:ZZcell do ( - when a.1 + if isZero(n.v) + then toExpr(Ccode(voidPointer, "&ffi_type_pointer")) + else when a.1 is signed:Boolean do ( bits := toInt(n); if signed.v then ( @@ -243,23 +248,24 @@ ffiIntegerType(e:Expr):Expr := ( else WrongNumArgs(2)); setupfun("ffiIntegerType", ffiIntegerType); +-- returns pointer to integer object with given value +-- inputs (x:ZZ, y:ZZ, signed:Boolean) +-- x = value +-- y = number of bits (or 0 for mpz_t) ffiIntegerAddress(e:Expr):Expr := ( when e - is x:ZZcell do ( - y := newZZmutable(); - set(y, x.v); - ptr := getMem(pointerSize); - Ccode(void, "*(void **)", ptr, " = ", y); - Ccode(void, "GC_REGISTER_FINALIZER(", ptr, - ", (GC_finalization_proc)mpz_clear, ", y, ", 0, 0)"); - toExpr(ptr)) is a:Sequence do ( if length(a) == 3 then ( when a.0 is x:ZZcell do ( when a.1 is y:ZZcell do ( - when a.2 + if isZero(y.v) then ( + z := copy(x.v); + ptr := getMem(pointerSize); + Ccode(void, "*(mpz_srcptr *)", ptr, " = ", z); + toExpr(ptr)) + else when a.2 is signed:Boolean do ( bits := toInt(y); ptr := getMemAtomic(bits / 8); @@ -302,16 +308,21 @@ ffiIntegerAddress(e:Expr):Expr := ( else WrongNumArgs(3)); setupfun("ffiIntegerAddress", ffiIntegerAddress); +-- returns value (as ZZ) of integer object given its address +-- inputs: (x:Pointer, y:ZZ, signed:ZZ) +-- x = pointer to integer object +-- y = number of bits (or 0 for mpz_t) ffiIntegerValue(e:Expr):Expr := ( when e - is x:pointerCell do toExpr(moveToZZ(Ccode(ZZmutable, "*(mpz_ptr*)", x.v))) is a:Sequence do ( if length(a) == 3 then ( when a.0 is x:pointerCell do ( when a.1 is y:ZZcell do ( - when a.2 + if isZero(y.v) + then toExpr(moveToZZ(Ccode(ZZmutable, "*(mpz_ptr*)", x.v))) + else when a.2 is signed:Boolean do ( bits := toInt(y); if signed.v then ( @@ -347,16 +358,23 @@ setupfun("ffiIntegerValue", ffiIntegerValue); -- real types -- ---------------- +-- returns pointer to ffi_type object for given real number type +-- input: n:ZZ (0 = mpfr_t, 32 = float, 64 = double) ffiRealType(e:Expr):Expr := ( when e is n:ZZcell do ( bits := toInt(n); - if bits == 32 then toExpr(Ccode(voidPointer, "&ffi_type_float")) + if bits == 0 then toExpr(Ccode(voidPointer, "&ffi_type_pointer")) + else if bits == 32 then toExpr(Ccode(voidPointer, "&ffi_type_float")) else if bits == 64 then toExpr(Ccode(voidPointer, "&ffi_type_double")) else buildErrorPacket("expected 32 or 64 bits")) else WrongArgZZ()); setupfun("ffiRealType", ffiRealType); +-- returns pointer to real number object with given value +-- inputs: (x:RR, y:ZZ) +-- x = value +-- y = type (0 = mpfr_t, 32 = float, 64 = double) ffiRealAddress(e:Expr):Expr := ( when e is a:Sequence do ( @@ -366,23 +384,33 @@ ffiRealAddress(e:Expr):Expr := ( when a.1 is y:ZZcell do ( bits := toInt(y); - ptr := getMemAtomic(bits / 8); - if bits == 32 - then ( - z := toFloat(x); - Ccode(void, "*(float *)", ptr, " = ", z)) - else if bits == 64 - then ( - z := toDouble(x); - Ccode(void, "*(double *)", ptr, " = ", z)) - else return buildErrorPacket("expected 32 or 64 bits"); - toExpr(ptr)) + if bits == 0 then ( + z := copy(x.v); + ptr := getMem(pointerSize); + Ccode(void, "*(mpfr_srcptr *)", ptr, " = ", z); + toExpr(ptr)) + else if bits == 32 || bits == 64 then ( + ptr := getMemAtomic(bits / 8); + if bits == 32 + then ( + z := toFloat(x); + Ccode(void, "*(float *)", ptr, " = ", z)) + else if bits == 64 + then ( + z := toDouble(x); + Ccode(void, "*(double *)", ptr, " = ", z)); + toExpr(ptr)) + else return buildErrorPacket("expected 0, 32, or 64")) else WrongArgZZ(2)) else WrongArgRR(1)) else WrongNumArgs(2)) else WrongNumArgs(2)); setupfun("ffiRealAddress", ffiRealAddress); +-- returns value (as RR) of real number object given its address +-- inputs: (x:Pointer, y:ZZ) +-- x = pointer to real number object +-- y = type (0 = mpfr_t, 32 = float, 64 = double) ffiRealValue(e:Expr):Expr := ( when e is a:Sequence do ( @@ -392,11 +420,13 @@ ffiRealValue(e:Expr):Expr := ( when a.1 is y:ZZcell do ( bits := toInt(y); - if bits == 32 + if bits == 0 + then toExpr(moveToRR(Ccode(RRmutable, "*(mpfr_ptr*)", x.v))) + else if bits == 32 then toExpr(Ccode(float, "*(float *)", x.v)) else if bits == 64 then toExpr(Ccode(double, "*(double *)", x.v)) - else buildErrorPacket("expected 32 or 64 bits")) + else buildErrorPacket("expected 0, 32, or 64")) else WrongArgZZ(2)) else WrongArgPointer(1)) else WrongNumArgs(2)) diff --git a/M2/Macaulay2/d/gmp.d b/M2/Macaulay2/d/gmp.d index c6513041d8f..d84c77fd29e 100644 --- a/M2/Macaulay2/d/gmp.d +++ b/M2/Macaulay2/d/gmp.d @@ -203,7 +203,7 @@ clear(x:RRimutable) ::= Ccode( void, "mpfi_clear(", x, ")" ); clear(z:CCmutable):void := ( clear(z.re); clear(z.im); ); -export moveToZZ(z:ZZmutable):ZZ := ( +export copy(z:ZZ):ZZ := ( y := GCmalloc(ZZmutable); Ccode(void, " int s = z->_mp_size, ss = s>=0 ? s : -s; @@ -212,13 +212,14 @@ export moveToZZ(z:ZZmutable):ZZ := ( ",y,"->_mp_alloc = ss, ",y,"->_mp_size = s, ",y,"->_mp_d = p; "); Ccode(ZZ,y)); +export moveToZZ(z:ZZmutable):ZZ := copy(Ccode(ZZ, z)); export moveToZZandclear(z:ZZmutable):ZZ := ( w := moveToZZ(z); clear(z); w); -export moveToRR(z:RRmutable):RR := ( +export copy(z:RR):RR := ( y := GCmalloc(RRmutable); Ccode(void, " int limb_size = (",z,"->_mpfr_prec - 1) / GMP_NUMB_BITS + 1; @@ -231,6 +232,7 @@ export moveToRR(z:RRmutable):RR := ( "); Ccode(RR,y) ); +export moveToRR(z:RRmutable):RR := copy(Ccode(RR, z)); export moveToRRi(z:RRimutable):RRi := ( y := GCmalloc(RRimutable); diff --git a/M2/Macaulay2/packages/ForeignFunctions.m2 b/M2/Macaulay2/packages/ForeignFunctions.m2 index 2a22f759f65..2ef4e745433 100644 --- a/M2/Macaulay2/packages/ForeignFunctions.m2 +++ b/M2/Macaulay2/packages/ForeignFunctions.m2 @@ -1,7 +1,7 @@ newPackage("ForeignFunctions", Headline => "foreign function interface", Version => "0.3", - Date => "January 28, 2023", + Date => "February 8, 2024", Authors => {{ Name => "Doug Torrance", Email => "dtorrance@piedmont.edu", @@ -15,9 +15,9 @@ newPackage("ForeignFunctions", -* -0.3 (2024-01-28, M2 1.23) +0.3 (2024-02-08, M2 1.23) * add subscripted assignment for various pointer types -* add support for GMP integers +* add support for GMP integers and MPFR reals * add support for describe, expression, toExternalString, and toString * use null coalescing operator @@ -77,6 +77,7 @@ export { "mpzT", "float", "double", + "mpfrT", "voidstar", "charstar", "voidstarstar", @@ -258,12 +259,7 @@ if version#"pointer size" == 4 then ( ) else ( long = int64; ulong = uint64) - -mpzT = new ForeignIntegerType; -mpzT.Name = "mpz_t"; -mpzT.Address = ffiPointerType -new mpzT from ZZ := (T, n) -> new T from {Address => ffiIntegerAddress n} -value mpzT := ffiIntegerValue @@ address +mpzT = foreignIntegerType("mpz_t", 0, true) ForeignIntegerType Number := ForeignIntegerType Constant := (T, x) -> new T from truncate x @@ -289,12 +285,13 @@ foreignRealType(String, ZZ) := (name, bits) -> ( float = foreignRealType("float", 32) double = foreignRealType("double", 64) +mpfrT = foreignRealType("mpfr_t", 0) ForeignRealType Number := ForeignRealType Constant := (T, x) -> new T from realPart numeric x ForeignRealType RRi := (T, x) -> T toRR x -isAtomic ForeignRealType := T -> true +isAtomic ForeignRealType := T -> if T === mpfrT then false else true -------------------------- -- foreign pointer type -- @@ -957,6 +954,29 @@ doc /// double /// +doc /// + Key + mpfrT + (NewFromMethod, mpfrT, RR) + (value, mpfrT) + Headline + MPFR multiple-precision floating-point type + Description + Text + Macaulay2's native @TO RR@ real number type wraps around @TT + "mpfr_t"@ from @HREF{"https://mpfr.org/", "MPFR"}@. This type + (which is an instance of @TO ForeignRealType@) allows for + conversion between Macaulay2 reals and MPFR reals without loss + of precision. + Example + mpfrT numeric(100, pi) + value oo + mpfrAdd = foreignFunction("mpfr_add", void, {mpfrT, mpfrT, mpfrT, int}) + x = mpfrT 0p100 + mpfrAdd(x, numeric(100, pi), exp 1p100, 0) + x +/// + doc /// Key (symbol SPACE, ForeignRealType, Number) @@ -1962,6 +1982,7 @@ assert Equation(value mpzT 10^100, 10^100) -- real types assert Equation(value float 3.14159, 3.14159p24) assert Equation(value double 3.14159, 3.14159p53) +assert Equation(value mpfrT 3.14159p100, 3.14159p100) -- pointer types ptr = address int 3