Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Brent Baccala
committed
Apr 27, 2018
1 parent
8dd9f87
commit 285c14a
Showing
6 changed files
with
789 additions
and
0 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,153 @@ | ||
Comparison of open source computer algebra systems | ||
|
||
Axiom | ||
doesn't use readline | ||
builtin hyperdoc is way behind HTML | ||
couldn't figure how to parse: | ||
integrand := x*((x^2*exp(2*x^2)-log(x+1)^2)^2) / ((x+1)*(log(x+1)^2 - x^2*exp(2*x^2))^2) | ||
|
||
|
||
Maxima | ||
uses LISP-style function definition syntax | ||
|
||
MAXIMA function tex_displ(MLABEL %O14 EXPR) | ||
- fetch tex-environment - either "" (strings) OR tex-environment property on car OR *tex-environment-default* | ||
- print tex1(LABEL %O14) (minus 14 chars because of trailing \mathbb{false}) | ||
- then print env around tex1(EXPR) | ||
|
||
MAXIMA function print(EXPR) calls displa for output | ||
|
||
|
||
Sage | ||
struggling to do polynomial long division | ||
calls Maxima or Singular to do this | ||
implements algebraic closure of Q (QQbar) | ||
|
||
+ fix latex printing on power series | ||
- subs should raise NotImplementedError if element not callable | ||
- subs should work recursively on complicated rings | ||
- option to specify how QQbar elements print | ||
- can't type convert between two slightly different NumberRings: | ||
Number Field in a with defining polynomial a^7 - 18 | ||
Number Field in a with defining polynomial y^7 - 18 | ||
- add map_coefficients method to Laurent series (Taylor, too?) | ||
|
||
|
||
Macaulay 2 | ||
doesn't implemente fraction fields over complicated (non-Z, non-Q) constant fields | ||
|
||
Singular | ||
can't compute Riemann-Roch spaces in characteristic p | ||
can compute Hamber-Noether expansions (can this be converted to a Puiseux expansion?) | ||
|
||
CoCoA | ||
C++ library | ||
|
||
|
||
Macaulay 2 - Riemann-Roch space calculator (Noether-Brill algorithm ??) | ||
|
||
completeLinearSystemOnNodalPlaneCurve=method() | ||
completeLinearSystemOnNodalPlaneCurve(Ideal,List):=(J,D)->( | ||
singJ:=saturate(ideal jacobian J+J); | ||
-- adjoint ideal | ||
H:=ideal (mingens ideal(gens intersect(singJ,D_0)%J))_(0,0); | ||
-- a curve passing through singJ and D_0 | ||
E0:=((J+H):D_0):(singJ^2); -- residual divisor | ||
if not(degree J *degree H - degree D_0 -2*degree singJ==degree E0) | ||
then error"residual divisor of has wrong degree"; | ||
L1:=mingens ideal (gens truncate(degree H, intersect(E0,D_1,singJ)))%J; | ||
h0D:=(tally degrees source L1)_{degree H}; -- h^0 O(D) | ||
L:=L1_{0..h0D-1}; -- matrix of homogeneous forms, L/H =L(D) subset K(C) | ||
(L,(gens H)_(0,0))) | ||
|
||
http://www2.macaulay2.com/Macaulay2/doc/Macaulay2-1.11/share/doc/Macaulay2/Macaulay2Doc/html/___Tutorial_co_sp__Divisors.html | ||
|
||
|
||
|
||
|
||
|
||
a=QQ['a'].0 | ||
aRing = NumberField(a^2 + 1, 'a') | ||
|
||
b=QQ['b'].0 | ||
bRing = NumberField(b^2 + 1, 'a') | ||
|
||
aRing is bRing | ||
|
||
aa=aRing.0 | ||
bb=bRing.0 | ||
|
||
bRing(aa) | ||
|
||
def ringmap(x): | ||
# First map x from the polynomial ring down to its coefficent ring (a NumberField) | ||
x = x.parent()(x) | ||
# Then convert it to a polynomial and map to into numfield | ||
return x.polynomial()(numfield.gens()[0]) | ||
|
||
|
||
HOW TO CREATE A MAP BETWEEN TWO POLYNOMIAL RINGS INDUCED BY A MAP ON THEIR COEFFICIENTS | ||
|
||
R.<a> = QQ[] | ||
aRing = NumberField(a^2 + 1, 'a') | ||
R.<b> = QQ[] | ||
bRing = NumberField(b^2 + 1, 'a') | ||
|
||
aRing2 = aRing['x'] | ||
bRing2 = bRing['x'] | ||
|
||
h = aRing.hom(bRing.gens()) | ||
|
||
h2 = aRing2.hom(h, bRing2) | ||
|
||
DOESN'T WORK FOR LAURENT SERIES | ||
|
||
|
||
|
||
BUG REPORT TO SHOW HOW SINGULAR RETURNS RESULTS IN WRONG FIELD | ||
|
||
R.<b> = QQ[] | ||
S = NumberField(b^7 - 18, 'a') | ||
|
||
R.<x,y> = S[] | ||
f = x + y | ||
|
||
singular.lib('hnoether.lib') | ||
hne = f._singular_().hnexpansion().sage() | ||
|
||
hne[0][0].base_ring() is f.parent() | ||
|
||
|
||
|
||
SAGE - SINGULAR IDEALS | ||
|
||
rings/polynomial/multi_polynomial_ideal.py | ||
lines 558-592 | ||
create Singular ideal | ||
lift generators into parent ring before str-ifying them for Singular | ||
|
||
interfaces/singular.py | ||
lines 1957-1959 | ||
converts Singular ideal to sage ideal | ||
|
||
1958 R = R or self.sage_global_ring() | ||
1959 return R.ideal([p.sage_poly(R) for p in self]) | ||
|
||
if R is a quotient ring, | ||
p.sage_poly(R) can't be lifted into the base ring | ||
|
||
sage_poly: | ||
|
||
1779 p = MPolynomial_polydict(R,PolyDict(sage_repr,force_int_exponents=False,force_etuples=False)) | ||
1780 if isinstance(R, MPolynomialRing_polydict): | ||
1781 return p | ||
1782 else: | ||
1783 return QuotientRingElement(R,p,reduce=False) | ||
|
||
what we get back is a QuotientRingElement containing an MPolynomial_polydict | ||
|
||
self.sage_global_ring() will return a quotient ring different from the one we started with | ||
its generators will be named 'x', 'y', 'z', instead of 'xbar', 'ybar', 'zbar' | ||
|
||
subs doesn't work because _call_ isn't defined on QuotientRingElement's | ||
(subs in RingElement returns self if _call_ isn't defined) |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,34 @@ | ||
#!/bin/sh | ||
# | ||
# gnuplot script to plot roots of an equation | ||
|
||
equation=$1 | ||
if [ "x$equation" = "x" ]; then equation="z**6-1"; fi | ||
|
||
/usr/bin/gnuplot -geometry 500x500 -persist <<EOF | ||
#cat > cyclo.out <<EOF | ||
unset ztics | ||
set ticslevel 0 | ||
set grid | ||
set zeroaxis lt 3 | ||
min(a,b) = a>b ? b : a | ||
set isosamples 500,500 | ||
set view 0,0 | ||
set key off | ||
set size 1.4,1.4 | ||
set origin -0.2,-0.2 | ||
strength = 1.0 | ||
poly(z) = (strength * ($equation)) | ||
splot [-2:2] [-2:2] [5:100] min(abs(1/poly(x+{0,1}*y)),10) | ||
bind Up "strength = strength * 2; replot" | ||
bind Down "strength = strength / 2; replot" | ||
EOF |
Binary file not shown.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,157 @@ | ||
A bunch of responses: | ||
|
||
1. Thanks for providing a brief way of producing this error message. | ||
2. If you really want an answer, you could try | ||
|
||
integrate_use_rootof:true; | ||
|
||
integrate(1/(x^(5/2)+3*x^(1/3)),x); | ||
|
||
though I'm not sure the answer returned is correct or | ||
even sensible. | ||
|
||
It returns | ||
6*lsum(log(x^(1/6)-%r1)/(13*%r1^9),%r1,rootsof(x^(13/6)+3)) | ||
|
||
Mathematica produces an answer, but one which is quite large. | ||
Wolfram Alpha produces an answer that is different but is | ||
also large, and runs out of time, it seems, in analyzing the | ||
expression. Incidentally, it apparently allows the Maxima syntax :) | ||
|
||
I'm not sure the answer is correct or even sensible. The | ||
Mathematica system 10 could not simplify the derivative back | ||
to the input, at least not before I lost patience. | ||
|
||
|
||
|
||
3. My understanding of algebraic:true is that it uses | ||
only roots of integers (absent in your expression) | ||
%i (also absent), and tellrat() info (also absent). | ||
So setting algebraic to true is maybe not the right thing to do. | ||
|
||
Introducing both sqrt(x) and x^(1/3) by a single algebraic | ||
extension of degree 6 is possible via tellrat. I'm not | ||
sure this gets you to where you want. | ||
|
||
Or even if you really want this integral or are just poking at Maxima. | ||
|
||
4. If you really want to try to patch this bug in the handling of | ||
algebraic:true, | ||
+ integrate + ratsimp, not by fixing the bug, but by catching it in a | ||
different | ||
place, you can look at code in src/rat3a.lisp where you can see that the | ||
polynomial division in the function pquotient, at some point says, in | ||
effect, | ||
HEY! you want me to divide exactly -- that's what pquotient does -- but | ||
the remainder is not zero. There's no bug in pquotient, but there is a bug | ||
in someone who called it on bad argumemtns. (This is probably because | ||
the GCD routine in use fails when there are algebraically dependent | ||
indeterminates, but this handling of algebraics is touchy. maybe because | ||
of comment 1. (thanks) someone can find this??? Anyway, this | ||
error-catch of the throw in pquotient is what you would have to change, so | ||
that | ||
it doesn't just vomit all the way to the top, but to some other err-catch. | ||
Maybe you can figure out enough from the comments and the code to | ||
do this. The main body of code (but not all the comments) are almost | ||
50 years old, circa 1966 I think. | ||
|
||
5. I am not sure what you mean by "throw to the outside" and how this makes | ||
life easier. Perhaps you are resisting the "easy" way -- which is to use | ||
the Maxima read-meval-print loop to call your functionality -- and | ||
calling | ||
Maxima via some other language at the end of a pipe? | ||
(Why easy? load into Maxima some lisp function from a file foo, by | ||
load("c:/..../foo.lisp"); | ||
|
||
and in the file foo.lisp | ||
|
||
(defun $sherfgen(x) (call_whatever_you_want_next x)) | ||
|
||
;; call_whatever .... can link to fortran, python, jscript, java, most any | ||
windows dll, typical unix .o files, libraries etc | ||
|
||
And the load() command can be put into a maxima initialization file so | ||
it is loaded up without user attention. | ||
|
||
|
||
Anyway. I hope this info is useful, and you take the opportunity to learn | ||
more lisp and | ||
Maxima stuff. | ||
|
||
RJF | ||
|
||
|
||
On Mon, Sep 7, 2015 at 12:54 AM, David Scherfgen <d.scherfgen at googlemail.com | ||
> wrote: | ||
|
||
> Hello, | ||
> | ||
> the following integral results in an error "quotient by zero": | ||
> | ||
> integrate(ev(ratsimp(1/(x^(5/2)+3*x^(1/3))),algebraic),x); | ||
> | ||
> Unfortunately, I can't catch the error using errcatch. It goes right | ||
> through it, i.e. the whole evaluation is aborted, and I don't get any | ||
> result (usually, errcatch would give me an empty list if it caught an | ||
> error). | ||
> | ||
> I've had a similar problem with questions being asked by Maxima - there is | ||
> no official way to intercept them. This is really bad if you run Maxima | ||
> "automatedly" - without anybody being there to answer questions. | ||
> | ||
> But luckily, there is Robert Dodier's noninteractive package that | ||
> re-defines the Maxima functions that are responsible for asking those | ||
> questions. The re-defined functions don't ask the questions, but throw | ||
> them, so that it can be caught from the outside. | ||
> | ||
> I tried to do something similar with the function that's responsible for | ||
> generating the error, the function merror in merror.lisp. | ||
> | ||
> This is my attempt at it: | ||
> | ||
> ; If this is enabled, errors will be intercepted and thrown. | ||
> (defmvar $no_errors t) | ||
> | ||
> ; A wrapper around MERROR in merror.lisp. | ||
> (defvar *real-merror* (symbol-function 'merror)) | ||
> (defun merror (sstring &rest l) | ||
> (if $no_errors | ||
> (meval `(($throw) '(($merror) ,sstring ,@ l))) | ||
> (apply *real-merror* sstring l))) | ||
> | ||
> In a simple test case, it seems to work: | ||
> | ||
> (%i1) catch(1/0); | ||
> (%o1) merror("expt: undefined: 0 to a negative exponent.") | ||
> | ||
> It throws the error instead of aborting the evaluation, just as I want. | ||
> | ||
> But in my original problem, it doesn't work, but crashes Maxima instead: | ||
> | ||
> (%i1) catch(integrate(ev(ratsimp(1/(x^(5/2)+3*x^(1/3))),algebraic),x)); | ||
> Maxima encountered a Lisp error: | ||
> Condition in MEVAL [or a callee]: INTERNAL-SIMPLE-ERROR: Bind stack | ||
> overflow. | ||
> Automatically continuing. | ||
> To enable the Lisp debugger set *debugger-hook* to nil. | ||
> | ||
> Unfortunately, my LISP knowledge is basically nonexistent. | ||
> | ||
> Does anyone of you know why this doesn't work? Why can't errcatch catch | ||
> the error in the first place? Is there something wrong with my merror | ||
> wrapper? | ||
> | ||
> Thanks for any help! | ||
> | ||
> Best regards, | ||
> David | ||
> | ||
> | ||
> ------------------------------------------------------------------------------ | ||
> | ||
> _______________________________________________ | ||
> Maxima-discuss mailing list | ||
> Maxima-discuss at lists.sourceforge.net | ||
> https://lists.sourceforge.net/lists/listinfo/maxima-discuss | ||
> | ||
> |
Oops, something went wrong.