Archive5
- To: mathgroup at yoda.ncsa.uiuc.edu
- Subject: Archive5
- From: uunet!lakisis.umd.edu!mbp
- Date: Mon, 22 Jan 90 12:45:10 EST
Below is Archive 5. There are a few short programs which people
submitted with their article, and I didn't bother to extract them
so that they can be separately included in the archives, since they
are not packages.
--Mark Phillips
[Thanks to Mark for his efforts. --smc]
---cut here------cut here------cut here------cut here------cut here---
Mathematica Notes from sci.math.symbolic, sci.math, and comp.sys.next
Archive 5: 20 November 1989 --- 22 January 1990
Editor: Mark Phillips
mbp at lakisis.umd.edu
Subjects:
Mathematica Syntax
Postscript output from mathematica
Mathematica question
Re: Mathematica question
Re: Mathematica question
Mathematica: Are They Crazy Or Am I?
Re: Mathematica: Are They Crazy Or Am I?
Re: Mathematica: Are They Crazy Or Am I
Re: Mathematica: Are They Crazy Or Am I?
Mathematica Wanted / MAC-SE/020/882
Re: Mathematica Wanted / MAC-SE/020/882
parser for Mathematica
3x+1 package for Mathematica
"Do YOU still do these integrals by hand?"
Re: non-commutative algebra in Mathematica?
Questions about Mathematica Special Input Forms
Re: non-commutative algebra in Mathematica?
Re: non-commutative algebra in Mathematica?
Simplification in Mathematica
Making a postscript file from a composite plot in Mathematica
factorizing polynomials
Fast RLE routines
Re: A request
Re: non-commutative algebra in Mathematica?
Re: A request
Re: A request
RunEncode competition
(*======================================================================*)
From: fateman at renoir.Berkeley.EDU (Richard Fateman)
Subject: Mathematica Syntax
Summary: How is Mathematica syntax defined?
Date: 27 Dec 89 17:49:52 GMT
Newsgroups: sci.math.symbolic
I find Mathematica grammar to be a hodge-podge, and I wonder
if anyone else has found this to be bothersome. Because neither
the manual (and Roman Maeder's additional book on programming
in Mathematica) have a formal definition of the language
(e.g. in Backus-Naur Form), it is difficult to predict the
acceptability of certain (admittedly marginal) forms.
It is sometimes difficult to predict the structure of those
forms accepted by the system, as well.
My question really is,
Does anyone (presumably at WRI) actually have a complete formal
grammar definition for Mathematica? (and can it be made public?)
In case you think you think the answer is already in the manual:
The list of "precedences" in the manual are inadequate because
each operator must have left- and right- binding powers, and
the notion that there are only "expressions" and "symbols"
in the language is unsupported by the behavior of the system.
Furthermore, there are just plain errors and omissions in the table.
Reserving for a future time comments on the unpredictability
of the semantics, let me give just a few examples of the grammar..
Input Output in FullForm -- comments
1 . 2 Dot[1,2]
1. 2 2.0
1 .2 0.2
1.2.3 --Error-- not Dot[1.2,3]
a.!b.c Dot[a,Not[Dot[b,c]] -- not? Dot[a,Not[b],c]
a**!b**c --Error-- not NonCommutativeMultiply[a,Not[NCM[b,c]]]
a>b< c Inequality[a,Greater,b,Less,c] ... but
a==b>c Equal[a,Greater[b,c]] not Inequality[a,Equal,b,Greater,c]
a!*!b Times[Factorial[a],Not[b]] ... but
a! !b Times [Factorial[Factorial[a]],b]
a!!b Times [Factorial2[a],b] -- there are lots of oddities with "!"
(a,b,c) --almost impossible to print out-- it is Segment[a,b,c]
16^^1a 26 ...but
16^^1 a a
16^^1x --Error-- not 1*x or x
f[!b] F[Not[b]] .. but
!b -- attempts to execute b as shell command
There are many other oddities, some may be merely a lack of conventionality
for someone familiar with Pascal or C, and some are fairly esoteric...
So, to repeat the question: Does anyone have a BNF for Mathematica?
(probably also required: a definition of the lexical structure...)
If so, is it available for public consumption?
By the way, I found some of the problems in the language in the process of
writing a Mathematic scanner/parser, in Common Lisp.
Though I am still tinkering with the program, my hope is to be able to
parse Mathematica programs and libraries independently of Mathematica.
(and to allow others to do so as well)
Richard J. Fateman
fateman at berkeley.edu
(*----------------------------------------------------------------------*)
From: jeremy at math.lsa.umich.edu (Jeremy Teitelbaum)
Subject: Postscript output from mathematica
Date: 29 Dec 89 20:06:50 GMT
Newsgroups: sci.math.symbolic
I would like to include postscript output from Mathematica into
TeX files, using the \special command in TeX. Has anyone had
any luck with this?
--
jeremy at math.lsa.umich.edu
Jeremy Teitelbaum
Math Dept.
U. of Michigan
Ann Arbor, MI 48109
(*----------------------------------------------------------------------*)
From: jesperse at amelia.nas.nasa.gov (Dennis C. Jespersen)
Subject: Mathematica question
Summary: NIntegrate of function defined with Which fails?
Date: 30 Dec 89 00:02:27 GMT
Newsgroups: sci.math.symbolic
This is Mathematica NeXT Release 1.0, Kernel version 1.2, Front end
version 1.2, on a NeXT.
With psi[x_] := Which[0<x && x<=1, 1, True, 0]
I find Plot[psi[x], {x,-1,2}] works fine, but
NIntegrate[psi[x], {x,-1,2}] returns 0;
NIntegrate[psi[x], {x,0,1}] also returns 0.
With PSI[x_] := Which[x<0,0,x>1,0,True,1] both Plot and NIntegrate
are OK.
With phooey[x_] := Which[x<0,0,x<1,x^2,x<2,-2*x^2+6*x-3,x<3,(3-x)^2,x>3,0],
Plot[phooey[x[], {x,0,3}] is OK but
NIntegrate[phooey[x], {x,0,3}] gives
Integrand Null is not numerical at {2.94043}.
So how should one go about defining a function in piecewise form so that
NIntegrate can be successful?
Dennis Jespersen
MS 202A-1
NASA Ames Research Center
Moffett Field, CA 94035
(415)694-6742
jesperse at wk18.nas.nasa.gov
jesperse at amelia.nas.nasa.gov
(*----------------------------------------------------------------------*)
From: moler at simplicity.Stanford.EDU (Cleve Moler)
Subject: Re: Mathematica question
Date: 30 Dec 89 05:43:30 GMT
Newsgroups: sci.math.symbolic
Dennis Jespersen asks about Mathematica:
> So how should one go about defining a function in piecewise form so that
> NIntegrate can be successful?
One of his examples is:
> phooey[x_] := Which[x<0,0,x<1,x^2,x<2,-2*x^2+6*x-3,x<3,(3-x)^2,x>3,0],
The numerical integrators in MATLAB, "quad" and "quad8", can integrate
this function reliably, but their behavior is interesting. Consider
integration over the interval [-s, 3+s] for some s >= 0. This
contains the support interval [0, 3] and so the result should be the
same for all s . MATLAB does get the right answer (2.0) for any s,
but the amount of work it does can vary considerably. The two
integrators use adaptive methods based on Simpson's rule and the
8-panel Newton rule. They begin by dividing the interval
into subintervals. If the subdivision points happen to fall at
the break points of Jespersen's piecewise quadratic example,
the rules get the exact answer immediately and terminate quickly.
This happens, for example, with s = 0.5. But for "most" values
of s, the break points fall within the subintervals, the
discontinuities in the second derivative are "felt" and many
levels of recursion and subdivision are required to get a
specified accuracy.
By the way, I expressed the function in a manner that allows
it to be evaluated at many points in a vector x :
function y = phi(x)
y = (x.^2) .* (0<x) .* (x<=1) + ...
(-2*x.^2+6*x-3) .* (1<x) .* (x<=2) + ...
(3-x).^2 .* (2<x) .* (x<=3);
The .* and .^ are MATLAB's notation for element-by-element operations
on vectors. Although this is still a little cryptic, I find it
easier to read than the Which notation.
-- Cleve Moler
The MathWorks
moler at mathworks.com
or, na.moler at na-net.stanford.edu
(*----------------------------------------------------------------------*)
From: rivin at Gang-of-Four.Stanford.EDU (Igor Rivin)
Subject: Re: Mathematica question
Date: 30 Dec 89 19:09:58 GMT
Newsgroups: sci.math.symbolic
In article <4365 at amelia.nas.nasa.gov> jesperse at amelia.nas.nasa.gov (Dennis C. Jespersen) writes:
>This is Mathematica NeXT Release 1.0, Kernel version 1.2, Front end
>version 1.2, on a NeXT.
>
>With psi[x_] := Which[0<x && x<=1, 1, True, 0]
>
>I find Plot[psi[x], {x,-1,2}] works fine, but
>
>NIntegrate[psi[x], {x,-1,2}] returns 0;
>
>NIntegrate[psi[x], {x,0,1}] also returns 0.
>
>With PSI[x_] := Which[x<0,0,x>1,0,True,1] both Plot and NIntegrate
>are OK.
>
>With phooey[x_] := Which[x<0,0,x<1,x^2,x<2,-2*x^2+6*x-3,x<3,(3-x)^2,x>3,0],
>Plot[phooey[x[], {x,0,3}] is OK but
>NIntegrate[phooey[x], {x,0,3}] gives
> Integrand Null is not numerical at {2.94043}.
>
>So how should one go about defining a function in piecewise form so that
>NIntegrate can be successful?
>
>Dennis Jespersen
>MS 202A-1
>NASA Ames Research Center
>Moffett Field, CA 94035
>(415)694-6742
>jesperse at wk18.nas.nasa.gov
>jesperse at amelia.nas.nasa.gov
Here is a reply courtesy of Jerry Keiper (the designer of NIntegrate):
The reason
phooey[x_] := Which[x<0,0,x<1,x^2,x<2,-2*x^2+6*x-3,x<3,(3-x)^2,x>3,0]
cannot be integrated with NIntegrate[] is that Which[] always evaluates. The
way to fix this behavior is to allow phooey[] to evaluate only when x_ is
a number. e.g.
In[1]:= phooey[x_?NumberQ] := Which[ x<0, 0,
x<1, x^2,
x<2, -2*x^2+6*x-3,
x<3, (3-x)^2,
x>3, 0]
Since it is clear where the singularities occur, we can tell NIntegrate[]
to work between them rather than trying to resolve them:
In[2]:= NIntegrate[phooey[x], {x,-5,0,1,2,3,7}]
Out[2]= 2.
In[3]:= %-2
-16
Out[3]= -8.88178 10
--
Igor Rivin Wolfram Research, Inc.
rivin at Gang-of-Four.Stanford.EDU or
rivin at wri.com
(*----------------------------------------------------------------------*)
From: iphwk at TERRA.OSCS.MONTANA.EDU (Bill Kinnersley)
Subject: Mathematica: Are They Crazy Or Am I?
Date: 1 Jan 90 13:56:44 GMT
Newsgroups: sci.math.symbolic
What am I doing wrong here?
---------------
Mathematica 1.2 (August 30, 1989) [With pre-loaded data]
by S. Wolfram, D. Grayson, R. Maeder, H. Cejtin,
S. Omohundro, D. Ballman and J. Keiper
with I. Rivin and D. Withoff
Copyright 1988,1989 Wolfram Research Inc.
In[1]:= m*z^2/(p^2+z^2)^(3/2) + z/(p^2+z^2)
2
m z z
Out[1]= ------------ + -------
2 2 3/2 2 2
(p + z ) p + z
In[2]:= Together[%]
2 2 3
p z + m z + z
Out[2]= ----------------
2 2 3/2
(p + z )
--
--Bill Kinnersley
Physics Department Montana State University Bozeman, MT 59717
INTERNET: iphwk at terra.oscs.montana.edu BITNET: IPHWK at MTSUNIX1
226 Transfer complete.
(*----------------------------------------------------------------------*)
From: rivin at Gang-of-Four.Stanford.EDU (Igor Rivin)
Subject: Re: Mathematica: Are They Crazy Or Am I?
Date: 2 Jan 90 20:35:22 GMT
Newsgroups: sci.math.symbolic
In article <9001011356.AA13770 at terra.oscs.montana.edu> iphwk at TERRA.OSCS.MONTANA.EDU (Bill Kinnersley) writes:
>What am I doing wrong here?
>
>---------------
>
>Mathematica 1.2 (August 30, 1989) [With pre-loaded data]
>by S. Wolfram, D. Grayson, R. Maeder, H. Cejtin,
> S. Omohundro, D. Ballman and J. Keiper
>with I. Rivin and D. Withoff
>Copyright 1988,1989 Wolfram Research Inc.
>
>In[1]:= m*z^2/(p^2+z^2)^(3/2) + z/(p^2+z^2)
>
> 2
> m z z
>Out[1]= ------------ + -------
> 2 2 3/2 2 2
> (p + z ) p + z
>
>In[2]:= Together[%]
>
> 2 2 3
> p z + m z + z
>Out[2]= ----------------
> 2 2 3/2
> (p + z )
>
>
>--
>--Bill Kinnersley
> Physics Department Montana State University Bozeman, MT 59717
> INTERNET: iphwk at terra.oscs.montana.edu BITNET: IPHWK at MTSUNIX1
>226 Transfer complete.
You are not doing anything wrong. This is a known bug in Mathematica
1.2. Together command.
There was a maintenance upgrade (which you, apparently, have not
received) which fixed this problem around two months ago. If you send
me (or support at wri.com) the relevant information (your machine type,
registration #, etc) we will send you the upgrade. This also applies to
any other Mathematica user reading this who can duplicate the behavior
exhibited in the included message.
--
Igor Rivin Wolfram Research, Inc.
rivin at Gang-of-Four.Stanford.EDU or
rivin at wri.com
(*----------------------------------------------------------------------*)
From: chrstnsn at ux1.cso.uiuc.edu
Subject: Re: Mathematica: Are They Crazy Or Am I
Date: 3 Jan 90 01:40:52 GMT
Newsgroups: sci.math.symbolic
Bill:
I believe that you have run into a bug in an early version of
1.2 that was known and is now fixed:
Mathematica (sun4) 1.2 (November 6, 1989) [With pre-loaded data]
by S. Wolfram, D. Grayson, R. Maeder, H. Cejtin,
S. Omohundro, D. Ballman and J. Keiper
with I. Rivin and D. Withoff
Copyright 1988,1989 Wolfram Research Inc.
-- SunView graphics initialized --
In[1]:= m*z^2/(p^2+z^2)^(3/2) + z/(p^2+z^2)
2
m z z
Out[1]= ------------ + -------
2 2 3/2 2 2
(p + z ) p + z
In[2]:= Together[%]
2 2 2
m z + z Sqrt[p + z ]
Out[2]= ----------------------
2 2 3/2
(p + z )
Steve Christensen
University of Illinois
steve at ncsa.uiuc.edu
(*----------------------------------------------------------------------*)
From: flynn at pixel.cps.msu.edu (Patrick J. Flynn)
Subject: Re: Mathematica: Are They Crazy Or Am I?
Date: 3 Jan 90 12:32:54 GMT
Newsgroups: sci.math.symbolic
In article <1990Jan2.203522.15261 at Neon.Stanford.EDU> rivin at Gang-of-Four.Stanford.EDU (Igor Rivin) writes:
+In article <9001011356.AA13770 at terra.oscs.montana.edu> iphwk at TERRA.OSCS.MONTANA.EDU (Bill Kinnersley) writes:
+>What am I doing wrong here?
+>
+> [description of bug in 1.2 version of Together[] omitted]
+
+You are not doing anything wrong. This is a known bug in Mathematica
+1.2. Together command.
+There was a maintenance upgrade (which you, apparently, have not
+received) which fixed this problem around two months ago. If you send
+me (or support at wri.com) the relevant information (your machine type,
+registration #, etc) we will send you the upgrade. This also applies to
+any other Mathematica user reading this who can duplicate the behavior
+exhibited in the included message.
Does this upgrade policy apply to the 1.2 Mathematica shipped with the
NeXT 1.0 system software? It does exhibit the bug. The Release notes
for Mathematica on the NeXT state that registering your copy of
Mathematica "... does not entitle you to free technical support,
updates, or other services from Wolfram Research", and suggests that we
contact NeXT about these issues.
So-
How is this upgrade handled? And by whom (NeXT or Wolfram)?
--
Patrick Flynn, CS, Mich. State U., flynn at cps.msu.edu
(*----------------------------------------------------------------------*)
From: bannon at andromeda.rutgers.edu.rutgers.edu (Ron Bannon)
Subject: Mathematica Wanted / MAC-SE/020/882
Date: 4 Jan 90 17:39:23 GMT
Newsgroups: sci.math.symbolic
MATHEMATICA
^^^^^^^^^^^
WANTED
^^^^^^
I'm thinking about purchasing Mathematica for my MAC but the price is a
little steep. Is it worth the $700, or is there any way to get it cheaper?
I plan to use Mathematica on my MAC-SE which has a 68020/68882 Radius
accelerator and 4 MB of RAM. The 4 MB is the maximum configuration and I'm a
little concerned that it may not be enough especially under MultiFinder.
Suggestions will be greatly appreciated.
Thanks
Ron Bannon
bannon at andromeda.rutgers.edu
bannon at math.rutgers.edu
Ron Bannon
bannon at andromeda.rutgers.edu
::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
Is George Bush a traitor? Read "October Surprise" by Honegger. Send for details.
::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
(*----------------------------------------------------------------------*)
From: fateman at renoir.Berkeley.EDU (Richard Fateman)
Subject: parser for Mathematica
Date: 6 Jan 90 00:22:26 GMT
Newsgroups: sci.math.symbolic
OK, its ready for a few brave souls. If you want to try
a Mathematica (tm) parser written in Common Lisp, I can make
it available to those who can ftp from some computer hereabouts.
I'm not going to send it by tape anywhere.
Respond to me (not to the whole net), and I'll tell you how to
get it. There are 2 parts: a LaTeX article describing it,
and a Lisp text file, about 30,000 characters. It would have been
smaller, but I decided to expand out most of the macros so
I could fiddle with borderline cases.
It parses about 50 lines/CPU second on an IBM RT.
Note that this is NOT a interpreter or a pattern matcher
for the internal form that is produced.
Example of input/output:
f[x_]:=If[x==0,1,x f[x-1]]
(SetDelayed (f (Pattern x (Blank)))
(If (Equal x (Integer 0)) (Integer 1)
(Times x (f (Plus x (Times (Integer -1) (Integer 1)))))))
Here's another Mathematica puzzle. What do you get if you type in
f(0,x)
and why is it always equal to f(x,0)? If you have a Mathematica,
try FullForm[Hold[f(x,0)]] for a hint.
Richard Fateman
fateman at renoir.berkeley.edu
(*----------------------------------------------------------------------*)
From: ilan at Gang-of-Four.Stanford.EDU (Ilan Vardi)
Subject: 3x+1 package for Mathematica
Date: 6 Jan 90 07:47:26 GMT
Newsgroups: sci.math.symbolic
(*
This Package is a new version of CollatzProblem.m
(the "3x+1 problem"). I have used a slightly different
definition of Collatz[n] than in the original
CollatzProblem.m but it follows the discovery by R. Terras
and C. Everett that you should really look at (3x + 1)/2
(see J. Lagarias, American Math. Monthly, January 1985).
I have also allowed negative values of n. You should be able
to check the conjecture for any value of n, n <= 2^(10^5).
After that I started to run into memory problems. Also, the programs
runs at least 4 times faster on the development version than in
version 1.2. This is probably due to division of large numbers
(1000 digits) by small numbers (3 digits) being speeded up.
You have to create a table so it may take about 5 minutes
to read in the file.
-Ilan Vardi, January 5, 1990. *)
Collatz[1] := {1}
Collatz[-1] := {-1}; Collatz[-5] := {-5}; Collatz[-17] := {-17};
Collatz[n_Integer] := Prepend[Collatz[(3 n +1) / 2], n] /; OddQ[n]
Collatz[n_Integer] := Prepend[Collatz[n / 2], n] /; EvenQ[n]
(* The following is a faster version of Collatz[n] *)
CollatzFaster[n_Integer] := Collatz[n] /; n < CollatzSmallRange
CollatzFaster[n_Integer] :=
Block[{cr},
cr = Flatten[Map[NestList[If[EvenQ[#], #/2, (3# + 1)/2]&,
#, CollatzIterateConstant - 1]&,
CollatzIterate[n]]];
Join[cr, Rest[Collatz[Last[cr]]]]]
(* CollatzFaster uses this function which lists the k'th Collatz
iterates, where k = CollatzIterateConstant (default = 10). *)
ColllatzIterate[n_Integer] := {} /; n < CollatzSmallRange
CollatzIterate[n_Integer] :=
Prepend[CollatzIterate[{n, 1} .
CollatzCoefficientListConstant [[1 +
Mod[n, CollatzListLength]]]], n]
(* The following is a more efficient method of computing how
many steps it takes to reach 1 (total stopping time) *)
CollatzTotalStoppingTime[n_Integer] := Length[Collatz[n]] - 1/;
n < CollatzSmallRange
CollatzTotalStoppingTime[n_Integer] :=
CollatzIterateConstant + CollatzTotalStoppingTime[{n, 1} .
CollatzCoefficientListConstant [[1 + Mod[n, CollatzListLength]]]]
CollatzTotalStoppingTime[n_Integer] :=
CollatzLargeIterateConstant + CollatzTotalStoppingTime[
NestList[Block[
{a = CollatzCoefficientListConstant[[1 +
Mod[#[[2]], CollatzListLength] ]]},
{SemiProduct[#[[1]], a],
Mod[{#[[2]], 1} . a, CollatzLargeListLength]}] &,
{{1, 0}, Mod[n, CollatzLargeListLength]},
CollatzLargeFactor] [[-1, 1]] . {n, 1}] /;
n > CollatzLargeRange
SemiProduct[{a_, b_}, {c_, d_}] := {a c, b c + d}
(* *******************************************************************
The above programs require a table that can be precomputed using
the following program.
**********************************************************************)
CollatzCoefficientList[k_Integer] :=
RotateRight[
Map[Function[x,
{3^Apply[Plus, x] / 2^Length[x],
Reverse[x] .
Accumulate[#1 #2 &, 1/2 Prepend[3^Reverse[Rest[x]],1]]}],
Mod[Map[NestList[If[EvenQ[#], #/2, (3 # + 1)/2]&, #, k-1] &,
Range[1, 2^k]], 2]]]
(***********************************************************************
Constants: Default for iteration speedup is 10. This seems optimal
so far.
************************************************************************)
CollatzIterateConstant = 10
CollatzListLength = 2^CollatzIterateConstant
CollatzCoefficientListConstant =
CollatzCoefficientList[CollatzIterateConstant]
CollatzSmallRange = 2^(CollatzIterateConstant + 1) + 2
CollatzLargeFactor = 100
CollatzLargeIterateConstant = CollatzIterateConstant CollatzLargeFactor
CollatzLargeListLength = 2^(CollatzIterateConstant CollatzLargeFactor)
CollatzLargeRange = Max[2^(CollatzIterateConstant CollatzLargeFactor +2),
2^1000]
Off[General::recursion]
$RecursionLimit = Infinity;
(* Running Times of above on Dec 3100
lasso [11] collatz> math13
Mathematica (DECstation 3100) 1.3.beta (January 3, 1990) [With pre-loaded data]
by S. Wolfram, D. Grayson, R. Maeder, H. Cejtin,
S. Omohundro, D. Ballman, and J. Keiper
with I. Rivin, D. Withoff, B. K. Smith, and H.-W. Tam
Copyright 1988,1989,1990 Wolfram Research, Inc.
-- Terminal graphics initialized --
In[1]:= Timing[<<CollatzProblemFaster.m]
Out[1]= {55.5667 Second, Null}
In[2]:= CollatzTotalStoppingTime[2^50 -1]//Timing
Out[2]= {0.383333 Second, 383}
In[3]:= CollatzTotalStoppingTime[2^500 - 1]//Timing
Out[3]= {1.86667 Second, 4331}
In[4]:= CollatzTotalStoppingTime[2^500 + 1]//Timing
Out[4]= {0.983333 Second, 2204}
In[5]:= CollatzTotalStoppingTime[2^5000 - 1]//Timing
Out[5]= {44.8667 Second, 43247}
In[6]:= CollatzTotalStoppingTime[2^50000 - 1]//Timing
Out[6]= {622.767 Second, 428838}
In[7]:= CollatzTotalStoppingTime[2^(10^5) - 1]//Timing
Out[7]= {1563.97 Second, 863323}
*)
(*----------------------------------------------------------------------*)
From: t19 at nikhefh.nikhef.nl (Geert J v Oldenborgh)
Subject: "Do YOU still do these integrals by hand?"
Date: 8 Jan 90 17:56:47 GMT
Newsgroups: sci.math.symbolic
"Do YOU still do these intergals by hand?"
it said in large letters on the Maple and Mathematica adds in our
(physics) lab. I decided not to be caught old-fashioned, got myself
an account on a machine with maple and mathematica and tried
a few integrals which I needed and which take 20 min to 1 hr by hand.
(plus some time to check numerically..)
k = sqrt(x^2-a^2);
int( log((p0*x - p*k)/(p0*x + p*k))/x,x=a..b,continuous);
(I really only need lim(a->0))
Result: Maple thinks a long time, and finally gives back the integral.
Mathematica the same. (The answer, in the limit a->0, is a set of 5
simple dilogs plus a divergent logarithm).
OK, try again:
int(int(1/(pDp*x^2 + pDq*x*y + qDq*y^2 + pDs*x + qDs*y + sDs - I*eps),
y=0..x),x=0..1);
in which the limit eps 0-> 0, eps > 0 is understood. Nothing. OK, let's
put eps=0 and not care about the imaginary part (which is simple).
Maple: nothing. Mathematica hangs up the machine. (The answer is a sum
of six fairly simple dilogs). These are the kind of integrals I do
daily.
"Do YOU still do these integrals by hand?" YES!
Geert Jan van Oldenborgh
working at, but not representing, NIKHEF-H, Amsterdam
(*----------------------------------------------------------------------*)
From: baez at x.ucr.edu (john baez)
Subject: Re: non-commutative algebra in Mathematica?
Date: 9 Jan 90 00:17:41 GMT
Newsgroups: sci.math.symbolic
In article <PALMIERI.89Dec14155157 at m2-225-6.athena.mit.edu> palmieri at m2-225-6.athena.mit.edu (John H Palmieri) writes:
>Does anyone out there know a good way to do non-commutative algebraic
>calculations with, say, Mathematica (on a Sun, in case that matters)?
>(I am particularly interested in the Steenrod algebra, by the way.)
>Thanks.
It's too bad I didn't run into you while I was at MIT over Xmas
(helping Segal with his book on quantum field theory). I've
used Mathematica to do some calculations in the Weyl algebra
over a 2-dimensional symplectic space; the main trick is
to use the ** (NonCommutativeMultiply) operator, and to
put in relations in such a way that expressions will be
reduced in the desired direction (taking care to avoid loops).
I never got that far in algebraic topology but if the Steenrod
algebra is associative and easily thought of in terms of
generators and relations this should work. Here's all that
was needed in my case:
Unprotect[NonCommutativeMultiply]
q1 ** p1 -> p1 ** q1 + I
q2 ** p2 -> p2 ** q2 + I
q1 ** p2 -> p2 ** q1
q2 ** p1 -> p1 ** q2
p2 ** p1 -> p1 ** p2
q2 ** q1 -> q1 ** q2
a_ ** (b_ + c_) -> (a ** b) + (a ** c)
(a_ + b_) ** c_ -> (a ** c) + (b ** c)
a_ ** b_ :> b a /; NumberQ[b]
a_ ** b_ :> a b /; NumberQ[a]
a_ ** (b_ c_) :> b (a ** c) /; NumberQ[b]
(a_ b_) ** c_ :> a (b ** c) /; NumberQ[a]
a_ ** (b_ + c_) -> (a ** b) + (a ** c)
(a_ + b_) ** c_ -> (a ** c) + (b ** c)
I could do calculations that were far too grungy to (want to) do by
hand, but ran into long run times before I could do all I
wanted to, probably because the above rules, while conceptually
simple, ate up time unnecessarily. I was too lazy to figure
out how to optimize things; I suspect calculations would go
much faster in a well-written C program, but the point is that
like many mathematicians I find such programming a tiresome
distraction when trying to prove theorems.
(*----------------------------------------------------------------------*)
From: stef at zweig.sun.com (Stephane Payrard)
Subject: Questions about Mathematica Special Input Forms
Date: 9 Jan 90 03:21:59 GMT
Newsgroups: sci.math.symbolic
Is there a function which returns the precedence of a
given operator (in its input form)?
Is there a way for a user to extend the set of Input Forms?
Disclaimers: I am a newcomer to this group, to Mathematica; and English
is not ny mother tongue.
--
--
Stephane Payrard -- stef at sun.com -- (415) 336 3726
Sun Microsystems -- 2550 Garcia Avenue -- M/S 16-40 -- Mountain View CA 94043
(*----------------------------------------------------------------------*)
From: t19 at nikhefh.nikhef.nl (Geert J v Oldenborgh)
Newsgroups: sci.math.symbolic
Subject: Re: non-commutative algebra in Mathematica?
Date: 9 Jan 90 15:31:19 GMT
In article <PALMIERI.89Dec14155157 at m2-225-6.athena.mit.edu> palmieri at m2-225-6.athena.mit.edu (John H Palmieri) writes:
>>Does anyone out there know a good way to do non-commutative algebraic
>>...
>I've used Mathematica to do some calculations in the Weyl algebra
>over a 2-dimensional symplectic space; ...
You might want to try FORM, which is very good at these kind of things.
Noncommuting objects and the linearity properties are built in and you can
sweep out commutators with a loop like
repeat;
id, p*q = q*p + whatever;
endrepeat;
The program is free and the executables for a dozen computers (including
sun4) and docs can be obtained by anonymous ftp from nikhefh.nikhef.nl.
It is also quite a bit faster than anything else I have used.
Geert Jan van Oldenborgh
NIKHEF-H, Amsterdam.
(*----------------------------------------------------------------------*)
From: t68 at nikhefh.nikhef.nl (Jos Vermaseren)
Newsgroups: sci.math.symbolic
Subject: Re: non-commutative algebra in Mathematica?
Date: 10 Jan 90 10:37:45 GMT
The program FORM can do most of these noncummutative things as
multiplication of functions is default noncummutative.
(You can define commuting functions though). So your objects
happen to be functions without arguments and everything drops
in place. You can make procedures that reduce your expressions
when it suits you, or put some smart reductions at the proper
locations (general reductions have a tendency to blow you out
of the computer).
You can pick up FORM (and a .dvi file of its manual) by means
of anonymous ftp from nikhefh.nikhef.nl (it works already for
many places and the rest should follow rather soon). There are
the following executables:
Apollo DN3000, Apollo DN10000, SUN4/Sparc, Gould 9080, Gould NP1,
DECstation 3100(ultrix), VMS5, ultrix (on regular VAX), Atari-ST,
IBM-PC (had to be pruned a little), MAC, Alliant.
There will be more.
Read the instructions before you print the manual.
Jos Vermaseren
(*----------------------------------------------------------------------*)
From: desouza at math.berkeley.edu (Paulo de Souza)
Newsgroups: sci.math.symbolic
Subject: Simplification in Mathematica
Date: 11 Jan 90 08:55:54 GMT
A question on simplification:
How can one make Mathematica perform simplifications of the type t^(2r) = -1,
on a series of polynomials computed, where the value of "r" changes during the
execution ?
What I have done is to write the rule in a package, but everytime the value "r"
changes I have to reread it, in order for it to accept the new value in the
simplifications. I am doing this inside a package where the value of "r" is an
input of one of the functions.
Just for completeness here is the essential part of the code I am using:
BeginPackage["KirMelInvariant`"]
.................................
Unprotect[Power]
Power/: t^(p_) := t^(Mod[p,4r]) /; p >= 4r
Power/: t^(p_) := t^(-Mod[-p,4r]) /; p <= -4r
Power/: t^(p_) := -t^(Mod[p,2r]) /; p >= 2r
Power/: t^(p_) := -t^(-Mod[-p,2r]) /; p <= -2r
Protect[Power]
....................................
Begin["`Private`"]
KM[p_,q_,r_] :=
Block[ { a, n, k, BigSum = 0 },
a = ContFract[p, q];
n = Length[a];
k = Table[1, {n}];
LoopFun[1,a,n,r,t];
BigSum = Expand[ BigSum ];
BigSum = Factor[ BigSum ];
BigSum
]
End[]
EndPackage[]
Paulo Ney desouza at math.berkeley.edu
(*----------------------------------------------------------------------*)
From: yip at cs.yale.edu (Ken Yip)
Subject: Making a postscript file from a composite plot in Mathematica
Date: 18 Jan 90 07:47:38 GMT
Newsgroups: sci.math.symbolic
On the sun4 version of Mathematica, it is easy to get a postscript
file from a plot by DISPLAY and psfix. But what about a plot
that consists of several plots superposed together?
One reason I want to do this is to annotate a plot with several text strings
(not just labels on axes). Maybe there are better ways to do this?
Any help??
(*----------------------------------------------------------------------*)
From: arf at maths.warwick.ac.uk (Anthony Iano-Fletcher)
Subject: factorizing polynomials
Date: 17 Jan 90 12:06:05 GMT
Newsgroups: sci.math.symbolic
Can anyone suggest good algorthms for
(i) factorising polynomials (with integer or
real coeficients)
(ii) finding the hcf of 2 such polynomials.
At the moment I am using Euclids algorthm for part (ii)
but this is a little slow.
Algorthms for (i) must exists since packages like Maple,
Mathematica etc can do it.
Also, does anyone know if Maple and Mathematica use efficient
algorthms for handling multiple precision integers.
Thanks in advance,
Anthony Iano-Fletcher.
--
UUCP: ...!ukc!warwick!arf Internet: arf at maths.warwick.ac.uk
JANET: arf at uk.ac.warwick.maths PHONE: +44 203 523523 x2677
(*----------------------------------------------------------------------*)
From: perry at phoenix.Princeton.EDU (Kevin R. Perry)
Subject: Fast RLE routines
Summary: find the fastest mathematica RLE encoder!
Date: 19 Jan 90 16:40:09 GMT
Newsgroups: sci.math.symbolic
At the recent MATHEMATICA conference in beautiful Redwood City, CA,
there was a programming competition; the challenge (for those of
you who weren't there) was to construct "the best" Run-Length encoding
and decoding routines in mathematica. For example, the RunEncode
function, when given an argument of { 1, 1, 2, 2, 2, 3, 3, 4 }
should return { {1,2}, {2,3}, {3,2}, {4,1} }. The competition
announcement stated "Entries will be judged on the basis of elegance,
simplicity, and efficiency." The winning entry certainly deserves
to be noted for it's elegant simplicity; as best I can remember, it was:
RunEncode[{rest___, w:(x_)..}] :=
Append[ RunEncode[{rest}], {x, Length[{w}]} ]
RunEncode[{}] := {}
But this leaves something to be desired in efficiency, which
is really the most important practical consideration to many of us.
I have also re-constructed a few other competition entries, and decided
to run some timing tests of my own. I found the results interesting.
I have duplicated these tests on each of SGI 4D/20, NeXT, and MacII,
and I used the built-in Timing[] function to determine speeds, on a
variety of datasets from 1 to 10,000 datapoints in length.
It turns out that the above function is the SLOWEST one I have yet
tried! So much for the "efficiency" part of the contest! On
*extremely* easy problems it does a little bit better than the most
inane procedural definition I could think up, but it does horrendously
worse on any reasonably realistic problem, presumably because of the
awesome number of recursive calls it ends up doing. The speediest
version I have been able to create so far is (oh happy surprise :-) the
function I submitted as an entry to the contest. [And I didn't even
get an honorable mention! :-( ]
Indeed, though mine is not particularly *simple*, I think
there is a certain practical elegance in it, and it IS extremely
fast - checking it on data of increasing length and complexity,
it quickly attains order-of-magnitude improvements over any
other entry I've been able to re-create from my conference notes.
There were some contenders that I haven't managed to test, and it's
quite possible one of them is the fastest. I'd like to find out what the
fastest possible routine is, so if you have any hints about methods of
doing it faster, or if anybody out there had a good competition entry,
send me a note so I can try it. Or compare it yourself first...
For the record, here's my "Most Efficient" run-length encoder:
RunEncode[list_] :=
( Transpose[ { list[[#]], # - Drop[ Prepend[#, 0], -1] } ])& [
Flatten[ Append[ Position[
((Apply[Unequal,#])& /@ Partition[list,2,1]), True
], Length[list] ] ] ]
RunEncode[{}] := {}
This guy processes the list Flatten[Table[Table[k, {100}], {k, 100}]]
in about 9 seconds on a 4D/20, compared to 1000 seconds or more for
most methods. Which says a lot about the care one must take in
designing efficient routines in MATHEMATICA! A highly functional
approach would appear to often work best, rather than procedural or
even rule-based or syntactic approaches. This seems to mean applying a
small number of functions directly to an entire list of data, often through
judicious use of Map or Apply, rather than using formulas that iterate
or recurse at a higher level.
"Long Live MATHEMATICA!"
Kevin Perry
CIT-ICGL
Princeton University
<perry at princeton.edu>
(*----------------------------------------------------------------------*)
From: wilkins at jarthur.Claremont.EDU (Mark Wilkins)
Subject: Re: A request
Date: 20 Jan 90 08:24:30 GMT
Newsgroups: sci.math.symbolic
In article <34892 at mips.mips.COM> mark at mips.COM (Mark G. Johnson) writes:
>I would dearly like to know the first 75 or so nonzero
>digits of e to the -2500 power.
1.83567266916215689307725305187314756746974669626209549936356569265433541
* 10^(-1086)
I got this using Mathematica 1.2 on a Macintosh II, using the command
In[1]:=
N[E^(-2500),75]
-- Mark Wilkins
wilkins at jarthur.claremont.edu
(*----------------------------------------------------------------------*)
From: seiford at ecs.umass.edu
Subject: Re: non-commutative algebra in Mathematica?
Date: 18 Jan 90 13:01:25 GMT
Newsgroups: sci.math.symbolic
In article <716 at nikhefh.nikhef.nl>, t19 at nikhefh.nikhef.nl (Geert J v Oldenborgh) writes:
> In article <PALMIERI.89Dec14155157 at m2-225-6.athena.mit.edu> palmieri at m2-225-6.athena.mit.edu (John H Palmieri) writes:
>>>Does anyone out there know a good way to do non-commutative algebraic
>>>...
>
>>I've used Mathematica to do some calculations in the Weyl algebra
>>over a 2-dimensional symplectic space; ...
>
> You might want to try FORM, which is very good at these kind of things.
> Noncommuting objects and the linearity properties are built in and you can
You should contact Eran Yehudai at Stanford (eran at slacvm.bitnet) who has
implemented several non-commutative algebras in Mathematica and has developed
some general tools.
larry
(*----------------------------------------------------------------------*)
From: stan at valley.UUCP (Stanley L. Kameny)
Subject: Re: A request
Date: 20 Jan 90 08:16:14 GMT
Newsgroups: sci.math.symbolic
In article <3859 at jarthur.Claremont.EDU> wilkins at jarthur.Claremont.EDU
(Mark Wilkins) writes:
> In article <34892 at mips.mips.COM> mark at mips.COM (Mark G. Johnson) writes:
>
> >I would dearly like to know the first 75 or so nonzero
> >digits of e to the -2500 power.
>
> 1.83567266916215689307725305187314756746974669626209549936356569265433541
> * 10^(-1086)
>
> I got this using Mathematica 1.2 on a Macintosh II, using the command
>
> In[1]:=
> N[E^(-2500),75]
>
> -- Mark Wilkins
> wilkins at jarthur.claremont.edu
Using REDUCE on a PC AT, and inputting precision 75; e^(-2500);
results in
1.835 67266 91621 56893 07725 30518 73147 56746 97466 96262 09549 936
35 65692 65433 54141 4 E - 1086
which is clearly 75 places. (I wonder where Mathematica lost the last 3
places?)
Stan Kameny
stan_kameny at rand.org
(*----------------------------------------------------------------------*)
From: wilkins at jarthur.Claremont.EDU (Mark Wilkins)
Subject: Re: A request
Date: 21 Jan 90 00:31:49 GMT
Newsgroups: sci.math.symbolic
In article <101.UUL1.2#239 at valley.UUCP> stan at valley.UUCP (Stanley L. Kameny) writes:
>
>which is clearly 75 places. (I wonder where Mathematica lost the last 3
>places?)
User error in cut-and-paste on my Mac. Oh, well. Those three digits may
cause World War III or something but I hope not...
-- Mark Wilkins
wilkins at jarthur.claremont.edu
(*----------------------------------------------------------------------*)
From: ilan at Gang-of-Four.Stanford.EDU (Ilan Vardi)
Subject: RunEncode competition
Date: 22 Jan 90 13:04:57 GMT
Newsgroups: sci.math.symbolic
I believe that I have discovered the simplest and most efficient
RunEncode program:
Attributes[times] = {Flat}
times[{a_, m_}, {a_, n_}] := times[{a, m + n}]
RunEncodeIlan[x_] := List @@ times @@ ({#, 1}& /@ x)
The idea of the program is that, algebraically, RunEncode is a
model of a free semigroup (it seemes that that the best approach is
to convert the problem into an algebraic question). For example,
if *** is a semigroup product then when you are not given any
relations (i.e., it's a free semigroup)
a *** b *** b *** a *** a => a *** b^2 *** a^2.
So all you need to do is to program the associative law, but this is
done for you by setting Attributes => {Flat}. I believe that most of
the work done by the other programs is essentially spent reprogramming
the associative law. Note that the "times" on the right hand side of
the second line is need to get the right answer for input like {a, a},
i.e., of length > 1 but with only one element.
You can use similar ideas to write weird programs that almost work:
Unprotect[Times]; Attributes[Times] := {Flat, Listable}
RunEncodeStrange[x_] := {1, -1} List @@ # & /@ List @@ Times @@ (1/x)
Here is a comparison of running times on a Sun 3 of my program with
the winning program and an unsubtle "C type" program
RunEncodeWin[{rest___, w:(x_)..}] :=
Append[RunEncodeWin[{rest}], {x, Length[{w}]}]
RunEncodeWin[{}] := {}
RunEncodeC[list_List] :=
Block[{t = {}, i = j = 1, k = Length[list]},
While[i <= k,
While[j <= k && list[[i]] == list[[j]], j++];
AppendTo[t, {list[[i]], j - i}]; i += j - i];
t]
timing[x_] := Block[{t},
t = {Timing[RunEncodeWin[x]],
Timing[RunEncodeIlan[x]],
Timing[RunEncodeC[x]]};
{First /@ t, SameQ @@ Last /@ t}]
Mathematica (sun3.68881) 1.2 (August 3, 1989) [With pre-loaded data]
by S. Wolfram, D. Grayson, R. Maeder, H. Cejtin,
S. Omohundro, D. Ballman and J. Keiper
with I. Rivin and D. Withoff
Copyright 1988,1989 Wolfram Research Inc.
Temporary license: expiration date is 1 Jun 1990.
-- Terminal graphics initialized --
In[1]:= <<RunEncode.m
In[2]:= x = {a,a,b,b,a,a,b,c,c,a,a}
Out[2]= {a, a, b, b, a, a, b, c, c, a, a}
In[3]:= RunEncodeWin[x]
Out[3]= {{a, 2}, {b, 2}, {a, 2}, {b, 1}, {c, 2}, {a, 2}}
In[4]:= RunEncodeIlan[x]
Out[4]= {{a, 2}, {b, 2}, {a, 2}, {b, 1}, {c, 2}, {a, 2}}
In[5]:= RunEncodeStrange[x]
Out[5]= {{a, 2}, {b, 2}, {a, 2}, {b, 1}, {c, 2}, {a, 2}}
In[6]:= x = Digits[3^20, 2]
Out[6]= {1, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 1,
> 1, 1, 0, 0, 1, 0, 0, 0, 1}
In[7]:= timing[x]
Out[7]= {{0.283333 Second, 0.133333 Second, 0.25 Second}, True}
In[8]:= timing[Digits[3^100, 2]]
Out[8]= {{10.9 Second, 2.13333 Second, 2.4 Second}, True}
In[9]:= timing[Flatten[{#,#,#}& /@ Range[100]]]
Out[9]= {{37.0667 Second, 6.85 Second, 4.65 Second}, True}
In[10]:= timing[Range[100]]
Out[10]= {{6.5 Second, 0.116667 Second, 2.2 Second}, True}
(*----------------------------------------------------------------------*)