# /* This file is part of msolve. # * # * msolve is free software: you can redistribute it and/or modify # * it under the terms of the GNU General Public License as published by # * the Free Software Foundation, either version 2 of the License, or # * (at your option) any later version. # * # * msolve is distributed in the hope that it will be useful, # * but WITHOUT ANY WARRANTY; without even the implied warranty of # * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # * GNU General Public License for more details. # * # * You should have received a copy of the GNU General Public License # * along with msolve. If not, see # * # * Authors: # * Christian Eder # * Jorge Garcia Fontan # * Huu Phuoc Le # * Mohab Safey El Din */ #WARNING: this file uses a sed command which does not work for mac users #Mac users should use the msolve-to-maple-file-interface-macos.mpl file FormatOutputMSolve:=proc(ll, _Z) local dim, deg, degquot, params, nvars, elim, den, cfs, i, varstr, linearform: dim:=ll[1]: if dim > 0 then return -1, [], []; fi; nvars:=ll[2]: degquot:=ll[3]: varstr:=ll[4]: linearform:=ll[5]: deg:=ll[6][1]: if nops(ll[6]) > 0 then elim:=PolynomialTools[FromCoefficientList](ll[6][2],_Z); else #computation failed return -2, [], []; fi: if nops(ll[5]) > 0 then den:=PolynomialTools[FromCoefficientList](ll[7][2],_Z); fi: params:=[]: cfs:=[1]: if degquot > 0 then for i from 1 to nvars-1 do params:=[op(params), PolynomialTools[FromCoefficientList](ll[8][i][1][2],_Z)]: cfs:=[op(cfs), ll[8][i][2]]: od: fi: return varstr, linearform, elim, den, [diff(elim, _Z), op(params)], cfs; end: GetRootsFromMSolve:=proc(l) local _Z, e, d, p, c, v, lf, sols, vals, i, sols_and_vals, realroots, OldDigits, mysols; local exist_gcd_p_e_nontrivial, length_p, gcd_p_e; v, lf, e, d, p, c := FormatOutputMSolve(l, _Z); if degree(e) = 0 then #positive dimension (-1) or computation failed (-2) return e, []; else return 0, [v,lf,e,d,p,c]; fi: end proc: ToMSolve:=proc(F, char, vars, fname) local i, fd, F2, str; fd:=fopen(fname, WRITE): for i from 1 to nops(vars)-1 do fprintf(fd, "%a, ", vars[i]): end; fprintf(fd, "%a ", vars[nops(vars)]): fprintf(fd,"\n"); fprintf(fd,"%d\n", char); if char = 0 then F2:=map(f->sort(expand(numer(f)), order=tdeg(op(vars))), F): F2:=remove(member, F2, [0]): for i from 1 to nops(F2)-1 do fprintf(fd, "%a,\n", F2[i]): od: fprintf(fd, "%a\n", F2[nops(F2)]): else F2:=map(f->sort(expand(numer(f)), order=tdeg(op(vars))) mod char, F): F2:=remove(member, F2, [0]): for i from 1 to nops(F2)-1 do fprintf(fd, "%a,\n", F2[i]): od: fprintf(fd, "%a\n", F2[nops(F2)]): fi: fclose(fd): str := cat("sed -i -e ':a' -e 'N' -e '$!ba' -e 's/\\\\\\n//g' ", fname): # str := cat("sed -i '' -e ':a' -e 'N' -e '$!ba' -e 's/\\\\\\n//g' ", fname): system(str): end proc: GetOptions:=proc(opts) local str, msolve_path, fname1, fname2, file_dir, verb, param, nthreads, output, gb; str:=subs(opts,"verb"); if type(str, integer) then verb:=str; else verb:=0: end if; str:=subs(opts,"leadmons"); if type(str, integer) and str > 0 then gb:=str; else gb:=2: end if; str:=subs(opts,"elim"); if type(str, integer) and str > 0 then elim:=str; else elim:=0: end if; str:=subs(opts,"output"); if type(str, integer) then output:=str; else output:=0: end if; str:=subs(opts,"nthreads"); if type(str, integer) then nthreads:=str; else nthreads:=1: end if; if member("mspath", map(lhs, opts)) then str:=subs(opts, "mspath"); if type(str, string) then msolve_path:=str; else printf("Error in format options"); end if; else msolve_path := "../msolve"; end if; if member("file_dir", map(lhs, opts)) then str:=subs(opts, "file_dir"); if type(str, string) then file_dir:=str; else printf("Error in format options"); end if; else file_dir := "/tmp/"; end if; if member("file_in", map(lhs, opts)) then str:=subs(opts, "file_in"); if type(str, string) then fname1:=cat(file_dir,str); else printf("Error in format options"); end if; else fname1 := cat(file_dir, RandomTools[Generate](string(8,alpha)), ".ms");; end if; if member("file_out", map(lhs, opts)) then str:=subs(opts, "file_out"); if type(str, string) then fname2:=cat(file_dir,str); else printf("Error in format options"); end if; else fname2 := cat(file_dir, RandomTools[Generate](string(8,alpha)), ".ms");; end if; param:=0; msolve_path, fname1, fname2, file_dir, verb, param, nthreads, output, gb, elim; end proc; #Input. #F: list of polynomials with rational coefficients #fc: field characteristic #vars: list of variables #Optional argument: # {"mspath"=, # "verb"=, # "file_dir"=, # "file_in"=, # "file_out"=, # "leadmons"= <1 to get leading monomials only> # "elim"=} #Output. #[] -> an error occured during the computation #else returns a Groebner basis of the ideal generated by F for the # grevlex ordering over monomials involving variables in vars # with vars[1] > ... > vars[n] # if "leadmons"=1 is part of the third optional argument, then only # the leading monomials are returned MSolveGroebner:=proc(F, fc, vars, opts:={}) local results, dim, fname1, fname2, verb, param, msolve_path, file_dir, field_char, lsols, nl, i, gb, output, nthreads, str, elim; if type(F, list(polynom(rational))) = false then printf("First argument is not a list of polynomials with rational coefficients\n"); end if; if type(fc, integer)=true then if fc < 0 then error "Field characteristic cannot be negative"; end if; if fc > 0 and isprime(fc)=false then error "Field characteristic should be a prime number"; end if; if fc > 2^31 then error "Field characteristic is too large to be supported"; end if; if fc > 0 then field_char := fc; end if; else printf("Second argument should be a prime integer < 2^31\n"); end if; if not(indets(F) subset indets(vars)) then printf("Given variables do not match the variables in the input polynomials\n"); end if; msolve_path:="../msolve"; file_dir:="/tmp/"; fname1:=cat(file_dir, RandomTools[Generate](string(8,alpha)), ".ms"); fname2:=cat(file_dir, RandomTools[Generate](string(8,alpha)), ".ms"); while evalb(fname1=fname2) do fname2:=cat(RandomTools[Generate](string(8,alpha)), ".ms"); od: verb:=0; param:=0: nthreads:=1; output:=0; gb:=2; if nops(opts) > 0 then msolve_path, fname1, fname2, file_dir, verb, param, nthreads, output, gb, elim := GetOptions(opts); fi; ToMSolve(F, field_char, vars, fname1); str := cat(msolve_path, " -v ", verb, " -g ", gb, " -e ", elim, " -P ", param, " -t ", nthreads, " -f ", fname1, " -o ", fname2): try system(str): read(fname2): catch: printf("There has been an issue in msolve computation\n"); return []; end: results:=%: system(cat("rm ", fname2)); system(cat("rm ", fname1)); return results; end proc: #Input. #F: list of polynomials with rational coefficients #vars: list of variables #Optional argument: # {"mspath"=, # "verb"=, # "file_dir"=, # "file_in"=, # "file_out"=, # "output"= <1 to get mid points of isolating boxes>} #Output. #[] -> an error occured during the computation #[1] -> input system has infinitely many complex solutions #[-1, []] -> input polynomial system has no real solution #[0, [sols]] -> input polynomial system has finitely many complex solutions # sols is the list of real solutions given with isolating boxes with the following format # [vars[1] = [a1, b1], ..., vars[n] = [an,bn]] #if "output"=1 is part of the third (optional) argument output format is # [ vars[1] = (a1+b1)/2, ..., vars[n] = (an+bn)/2 ] MSolveRealRoots:=proc(F, vars, opts:={}) local results, dim, fname1, fname2, verb, param, msolve_path, file_dir, lsols, nl, i, j, gb, output, nthreads, str, sols, prec; if type(F, list(polynom(rational))) = false then printf("First argument is not a list of polynomials with rational coefficients\n"); end if; if not(indets(F) subset indets(vars)) then printf("Given variables do not match the variables in the input polynomials\n"); end if; msolve_path:="../msolve"; file_dir:="/tmp/"; fname1:=cat(file_dir, RandomTools[Generate](string(8,alpha)), ".ms"); fname2:=cat(file_dir, RandomTools[Generate](string(8,alpha)), ".ms"); while evalb(fname1=fname2) do fname2:=cat(RandomTools[Generate](string(8,alpha)), ".ms"); od: verb:=0; param:=0: nthreads:=1; output:=0; gb:=0; if nops(opts) > 0 then msolve_path, fname1, fname2, file_dir, verb, param, nthreads, output, gb := GetOptions(opts); fi; ToMSolve(F, 0, vars, fname1); if Digits <= 10 then prec:=128: else prec:= iquo(Digits/10)*128: fi: str := cat(msolve_path, " -v ", verb, " -P ", param, " -t ", nthreads, " -f ", fname1, " -o ", fname2): gb:=0; #Needed to avoid the user stops GB comp once a prime computation is done param:=0; try system(str): read(fname2): catch: printf("There has been an issue in msolve computation\n"); return []; end: results:=%: system(cat("rm ", fname2)); system(cat("rm ", fname1)); if nops(results) = 0 then printf("There has been an issue in msolve computation\n"); return []; end if; dim := results[1]; if dim = -1 then if verb >= 1 then printf("Algebraic set is empty\n"); end if; return [-1, []]; end if; if dim > 0 then if verb >= 1 then printf("Algebraic set has positive dimension\n"); end if; return [1]; end if; if dim = 0 then lsols := results[2]; nl := lsols[1]: sols:=[]: for i from 1 to nl do sols:=[op(sols), op(map(_p->[seq(vars[j] = _p[j], j=1..nops(vars))], lsols[i+1]))]; od: if output=1 then sols := map(_p->map(_c->lhs(_c)=(rhs(_c)[1]+rhs(_c)[2])/2, _p), sols); end if; return [0, sols]; end if; return results; end proc: # #Example: Katsura-6 # #Input data # F:=[x1+2*x2+2*x3+2*x4+2*x5+2*x6-1, # x1^2+2*x2^2+2*x3^2+2*x4^2+2*x5^2+2*x6^2-x1, # 2*x1*x2+2*x2*x3+2*x3*x4+2*x4*x5+2*x5*x6-x2, # x2^2+2*x1*x3+2*x2*x4+2*x3*x5+2*x4*x6-x3, # 2*x2*x3+2*x1*x4+2*x2*x5+2*x3*x6-x4, # x3^2+2*x2*x4+2*x1*x5+2*x2*x6-x5]; # vars:=[x1,x2,x3,x4,x5,x6]; # #Usage # #with rational parametrization # param, sols:=MSolveRealRoots(F,vars,"../binary/msolve","/tmp/in.ms","/tmp/out.ms",1); # #or with solutions only # sols:=MSolveRealRoots(F,vars);