Unformatted text preview:

Computer Algebra Systems: Numerics“Symbolic Computation” includes numeric as a subset“Symbolic Computation” vs...What is the added value for Symb.+Num?Numerics tend to be misunderstoodSquare root of 9 by Newton IterationIt looks like it was getting worse, and then got betterMathematica has gotten more elaborateOther possibilities: IEEE binary FP stdStart with a numeric programming systemExplicitly add numeric libraries to CASRewrite all the code in lispNon-functional vs functional: the Fortran versionThe functional versionHow functional losesRepairing the functional versionRepairing using “registers”So the Problems can be fixed at some inconvenience.Even if CAS has bignums, link to outside..Why might GMP be faster?The size of k is criticalWhat about MPFUN, ARPREC ?Any other clever ideas?Why restrict outside libraries to floats?Richard Fateman CS 282 Lecture 17 1Computer Algebra Systems: NumericsLecture 17Richard Fateman CS 282 Lecture 17 2“Symbolic Computation” includes numeric as a subset•Why do CAS not entirely replace numeric programming environments?Richard Fateman CS 282 Lecture 17 3“Symbolic Computation” vs...•Purely Numeric Systems prosper . Why?–loss in efficiency is not tolerated–unless something sophisticated is going on, the symbolic system adds more complexity than necessary. (learning curve)–CAS systems are “extra cost”•Other reasons.–People are successful in the first approach they learned. They don’t change.–How else to explain FortranRichard Fateman CS 282 Lecture 17 4What is the added value for Symb.+Num?•SENAC-like systems (Computer Algebra, front end help systems)•Code-generation systems (GENTRAN)•integrated visualization, interaction, plotting•exact integer and rational arithmetic•extra precision (seamlessly)•interval arithmetic–explicit endpoints (range in Maple, Interval in MMa)–implicit intervals (significance arithmetic)Richard Fateman CS 282 Lecture 17 5Numerics tend to be misunderstood•Insufficient explanation about what is going on•Peculiar user expectations. Is 3.000 more accurate than 3.0? Is it more precise? •Why is sum(0.001,i=1,1000) only 0.99994?•Mathematica default makes simple convergent processes diverge.Richard Fateman CS 282 Lecture 17 6Square root of 9 by Newton Iteration•s[x_]:= x- (x2-9)/(2*x);•Nest[s,2,5]  (11641532182693481445313/3880510727564493815104... differs from 3 by•1/3880510727564493815104•Nest[s,2.0,5]-3 = “0.0” ;start interation at 2.•Nest[s,2.000000000000000,5]-3 = 0.x10-18•Nest[s,2.000000000000000,50]-3 = 0.x10-5•r=Nest[s,2.000000000000000,70]-3 = 0. Nest[s,2.00000000000000000000000000,88]-2 =0. //umh, you mean the iteration also converges to 2??Richard Fateman CS 282 Lecture 17 7It looks like it was getting worse, and then got better•InputForm[r] is 0``-0.4771•furthermore, r+1 prints as 0.Richard Fateman CS 282 Lecture 17 8Mathematica has gotten more elaborate•AccuracyGoal•WorkingPrecision•SetPrecision–beyond simple characterization•Claims (v 3) to run all routines to enough accuracy to provide (conservatively) as many digits correct as requested. [subsequently retracted?]•Decisions (e.g. sin(tan(x))< tan(sin(x)) for x near zero) can be tricky. Taylor series of difference starts as –x7/30+ ... )Richard Fateman CS 282 Lecture 17 9Other possibilities: IEEE binary FP std•Start with standard (IEEE float) and extend toward symbolic. IEEE 754, 854 (any radix). •Problematical: there are symbols like +/- infinity, not-a-number, signed 0, in IEEE, which take on some of the properties of symbols. What to do? In particular….•Is NaN a way to represent a symbol z? (a symbol is a number that is not a number?)•Rounding modes (etc) in software are time consuming when implemented poorly.Richard Fateman CS 282 Lecture 17 10Start with a numeric programming system•Matlab: add a “Maple Toolbox”. Allow symbols or expressions as strings in a matrix. •Limited integration of facilities.•Excel: add functionalities (again, using strings) as patches to a spreadsheet program.Richard Fateman CS 282 Lecture 17 11Explicitly add numeric libraries to CAS•Treat (say) numeric matrices as a special case: transfer to ordinary double-precision floats to do numerics.•Put all the work into good interfaces so that the CAS can guide the computation.•From lisp systems, “foreign function” calls?Richard Fateman CS 282 Lecture 17 12Rewrite all the code in lisp•How hard would it be to compile C or Fortran into Lisp, and then compile it from Lisp into binary code?•A program: f2CL exists. Major efforts to pound on it have improved it (credits: Kevin Broughan, Raymond Toy, me..)•How does this compare to FF?Richard Fateman CS 282 Lecture 17 13Non-functional vs functional: the Fortran version •x = x+1 in Fortran–load value of x from location L into a register Ra–add 1 into Ra [ignore overflow?]–store Ra into location L–Three assembler instructions. No memory.Richard Fateman CS 282 Lecture 17 14The functional version•(setf x (+ x 1)) in Lisp [or other functional style languages] –Load pointer to value of x from location L into register Rx–Load value of x into register Ra–Add 1 into Ra –Check for overflow: jump to bignumber routine–Check for a HEAP location for the answer: L2•If no space available, do garbage collection–Store L2 in heap and store (pointer to L2) in LRichard Fateman CS 282 Lecture 17 15How functional loses•a loop like this:do 100 times: x x+1can use up 100 cells of memory (heap)Richard Fateman CS 282 Lecture 17 16Repairing the functional version•(dsetv x (+ x 1)) in Lisp [or other functional style languages] // macro defining destructive version, generates (e.g. in GMP)(gmpz_add_ui (inside x) (inside x) 1)………………………… TARGET addend addendRichard Fateman CS 282 Lecture 17 17Repairing using “registers”How to generate temporary spaces/ registers/ at compile time?(let ((temp1 #.(runtime-allocated-temp)) (temp2 #.(runtime-allocated-temp))….) ….) (…hairy arithmetic needing temporaries temp1, temp2, ..)If compiled nicely, “temp2” might even be allocated on a stack, and the loop might use 1 (or zero) cells of memory...Richard Fateman CS 282 Lecture 17 18So the Problems can be fixed at some inconvenience.Superfast GCVery clever compiler (stack allocate vars etc.)Special encoding for likely inner-loop


View Full Document
Download Lecture Notes
Our administrator received your request to download this document. We will send you the file to your email shortly.
Loading Unlocking...
Login

Join to view Lecture Notes and access 3M+ class-specific study document.

or
We will never post anything without your permission.
Don't have an account?
Sign Up

Join to view Lecture Notes 2 2 and access 3M+ class-specific study document.

or

By creating an account you agree to our Privacy Policy and Terms Of Use

Already a member?