1 Introduction
Assignment 3
A Small Numerical Library
Prof. Darrell Long CSE 13SFall 2019
Due: October 20th at 11:59 pm
Just look at the graceful way that he lectures, one hand waves while the other conjectures.
As we know, computers are simple machines that carry out a sequence of very simple steps, albeit very quickly. Unless you have a specialpurpose processor, a computer can only compute addition, subtrac tion, multiplication, and division. If you think about it, you will see that the functions that might interest you when dealing with real or complex numbers can be built up from those four operations. We use them for nearly every program that we write, so we ought to understand how they are created.
If you recall from your calculus class, with some conditions a function fx can be represented by its Taylor series expansion near some point fa:
xakk
fxfa
k! f
a.
k1
If you have forgotten or never taken calculus, do not despair. Go to a laboratory section for review: the concepts required here are just derivatives.
Since we cannot compute an infinite series, we must be content to calculate a finite number of terms. In general, the more terms that we compute, the more accurate our approximation. For example, if we expand to 10 terms we get:
fxfaxafa 1xa2fa 1f3axa3 26
1 f4axa41 f5axa51 f6axa6
24
2019 Darrell LongRegents of the University of California 1
120
f8 axa8
362880 Oxa .
Taylor series, named after Brook Taylor, requires that we pick a point a where we will center the
f7 axa75040
720
f9 axa9 10
40320
approximation. In the case a0, then it is called a Maclaurin series. Often we choose 0, but the closer
to the value of x the better we will approximate the function. For example, lets consider ex centered around 0:
x x2 x3 x4 x5 x6 x7 x8 x9 10
e 1x 26 24120720504040320362880Ox .
This is one of the simplest series when centered at 0, since e01. Consider the general case: exeaeaxa1eaxa21eaxa3 1 eaxa4 1 eaxa5 1 eaxa6
2 6 24 120 720 eaxa7 eaxa8 eaxa9 eaxa10 11
5040403203628803628800 Oxa .
Since d exex the exponential function does not drop out as it does for a0, leaving us with our
dx
original problem. If we knew ea for ax then we could use a small number of terms. However, we do
1 2 3 4
not know it so we must use a0.
WhatistheOxa11term? Thatistheerrortermthatisontheorderofthevalueinparentheses.
This is different from the bigO that we will discuss with regard to algorithm analysis. 2 Your Task
Programming is one of the most difficult branches of applied mathematics; the poorer mathematicians had better remain pure mathematicians.
Edsger Dijkstra
For this assignment, you will be creating a small numerical library. Our goal is for you to have some idea of what must be done to implement functions that you use all of the time.
You will be writing and implementing sin, cos, tan, and exp using the Taylor series approximations for exp and Pade approximants for sin, cos, and tan. You will then run these functions and compare them to the standard library math.h implementations and output the results into a table similar to what is shown in Figures 1 and 2.
x Sin Library Difference
6.2832 0.00000000 0.00000000 0.0000000000 6.0868 0.19509032 0.19509032 0.0000000000
Figure 1: Example of program output for sine.
From left to right, the columns represent the input number, your programs cosine value from the input number, the actual math librarys value from the input number and lastly, the difference between your value and the librarys value.
You will test sine and cosine from 2 to 2 with steps of 16. Tangent will be tested from 20.001 to 20.001 with steps of 16. For the exponential function, ex will be from 0 to 10 with steps of 0.1.
2019 Darrell LongRegents of the University of California 2
1 2 3 4 5 6 7 8 9
10 11 12
Figure 2: Example of program output for cosine.
Each implementation will be a separate function. You must label the functions Exp, Sin, Cos, and Tan. Since the math library uses exp, sin, cos, and tan, you will not be able to use the same names. To use math.h in your program you must be sure to link the math library at compile time using the lm flag.
The following is an example function that implements Newtons method of computing square roots thatdoesntconflictwithsqrtfoundinmath.h.NotethatthefunctionisnamedSqrtdouble x.
define EPSILON 0.00001 defineABSx x0?x:x
double Sqrtdouble x
double y1.0;
double old0.0;
while ABSyoldEPSILON
oldy;
y0.5yxy;
1 x Cos
2
3 6.2832 1.00000000 4 6.0868 0.98078528
Library
1.00000000 0.98078528
Difference
0.0000000000 0.0000000000
return y;
2.1 Sine and Cosine
Computing x using Newtons method.
The domain of sin and cos is 2, 2, and so centering them around 0 makes sense. Since the domain is restricted, you can reduce any parameter to 2, 2 making your approximation better. The Taylor series for sinx centered about 0 is:
x3 x5 x7 x9 x11 x13 x15 16 sinxx612050403628803991680062270208001307674368000O x
and the corresponding series for cosx is:
x2 x4 x6 x8 x10 x12 x14 16
cosx1 22472040320362880047900160087178291200 Ox .
We can use what is called a Pade Approximant. Its beyond the scope of this course to go into com puting them, but essentially it is the ratio of two polynomials that conform to certain properties. It is often easier to compute and more accurate than a truncated series. The Pade approximant for a 15 term series for sinx centered around 0 is:
sinx 72x 62077121×67479608290×4242757835320×21727021696400 . 1768969×8366075360×643268148000×43245696213760×2124345562140800
2019 Darrell LongRegents of the University of California
3
It is a lot easier to square a number than to raise it to a power, so we can simplify it by putting the formula into Horner normal form, by factoring out x as much as possible:
sinxx x2 5385317968804469552712×2 17478564143040 x2124345562140800 . 1768969×2366075360 x243268148000 x23245696213760 x2124345562140800
If you want to be clever you can compute x2 once, store it in a variable, using that result again and save a few instructions. Does this matter? A good compiler will recognize the common subexpression and do it for you behind the scenes, but numerical code tends to be heavily used so every little bit helps.
Consider the corresponding approximant for cosx centered around 0 written in Horner normal form. Why Horner normal form? It has the fewest multiplications.
cosxx2 24694099×25773801680 x23440482214405595410620800 x211797350364800 . 209794×237502280 x24124237040 x2303264561600 x211797350364800
Restricting the domain is important. As we move away from the center of the series, our approxima tion gets worse. Consider Figure 3, and you can see that when we get much beyond 2,2 the graphs diverge not just a little, but wildly.
sin x versus Pade approximant 4
2
105 510
2
4
Figure 3: Comparing sinx with an order 10 Pade approximant.
2.2 Tangent
You will recall that tanxsinxcosx but that it is undefined when cosx0, that is when x is a multiple of 2. You could just do the division:
1 double tandouble xreturn sinxcosx;
2019 Darrell LongRegents of the University of California 4
While doing the division is correct in the real numbers, remember that we are working with approxi mations and so you will be more accurate with a formula for directly computing your approximation for tanx. A Pade approximant for tanx is:
tanxx x2 x2 2574055×2 2837835 x291891800654729075 . x2 x2 x21485 x231531518918900 x2310134825654729075
2.3 ex
The other important function is ex called the exponential function. While you might hope for a Pade approximant, you are bound for disappointment. The domain of ex is unrestricted so we will not be able to use any fixed formula. As we get further away from where we center our series, the more terms that will need to get good results. Fortunately, the series for ex converges very quickly.
Lets consider how well this works. In Figure 4, we will use our expansion to 10 terms and plot for e0,…,e10. After about 7, the approximation starts to diverge significantly. What this tells us is that 10 terms is not enough.
ex
8000
6000
4000
2000
Taylor Approximation
2 4 6 8 10
x
Figure 4: Comparing ex with its Taylor approximation centered at zero.
If we are naive about computing the terms of the series we can quickly get into troublethe values
of k! get large very quickly. We can do better if we observe that: xn x xn1
n!nn1! .
At first, that looks like a recursive definition and in fact, you could write it that way, but it would be wasteful. As we progress through the computation, assume that we know the previous result. We then just have to compute the next term and multiply it by the previous term. At each step we just need to
2019 Darrell LongRegents of the University of California 5
25
20
15
10
5
include stdio.h
2 4 6 8 10
Figure 5: Comparing xnn! for x2,3,4,5.
compute xn, starting with n0!1 and multiply it by the previous value and add it into the total. It turns into a simple for or while loop.
We can use anepsilon to halt the computation since xnn! for a sufficiently large n. Consider Figure 5: briefly, xn dominates but is quickly overwhelmed by n! and so the ratio rapidly approaches zero.
2.4 Commandline Arguments
Your program will use commandline arguments to select which one of your implemented functions to run. You will use getopt to parse commandline arguments, which is included under getopt.h. In most C programs, you will see two parameters in the main function: int argc and char argv. The parameter argc refers to the argument counter, or how many arguments were supplied on the commandline. Arguments are delimited by white space, which includes spaces and tabs. The parameter argv refers to the argument values, and is interpreted as an array of strings, where argv0 corresponds to the first argument on the commandline: the executable itself.
Try this code, and make sure that you understand it:
1
2
3
4
5
6
7 return 0; 8
Commandline options must be defined in order for getopt to parse them. These options are de fined in a string, where each character in the string corresponds to an option character that can be speci
2019 Darrell LongRegents of the University of California 6
int mainint argc, char argv
for int i0; iargc; i1
printfargvdsn, i, argvi;
1 2 3 4 5 6 7 8 9
10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38
fied on the on the commandline. Upon running the executable, getopt scans through the command line arguments, checking for option characters.
As an example, the following program supports two commandline options, p and i. Note that the option character i in the defined option string OPTIONS has a colon following it. The colon signifies that, when the i option is enabled on the commandline using i, getopt is looking for an argument to be supplied following it. An error is thrown by getopt if an argument for a flag requiring one is not supplied.
include getopt.h include stdio.h include stdbool.h
define OPTIONS pi:
int mainint argc, char argvint c0;
bool printfalse;
char infileNULL;
while cgetoptargc, argv, OPTIONS ! 1switch c
case p:Print option.
printtrue;
break;
case i:Input file option.
infileoptarg;A pointer to the next element in argv.
break;
if argc1
putsError: no arguments supplied!; return 1;
if !infile
putsError: input file required!; return 1;
if printputsinfile;
return 0;
2019 Darrell LongRegents of the University of California 7
What this program does is check for the p option to print and the i option to specify some input file name. The getopt defined variable optarg points at the following argvelement, which in this case should be the input file name. The condition argc1 checks if no other commandline arguments were supplied other than the executable itself and errors out if so. The program also errors out in the event that an input file name was not supplied. The supplied input file name is only printed if the print option was specified on the commandline.
Example usage of program to specify input file and print its name assume the executable name is filename:
1 .filename p i input file name
Your program must support the following commandline options:
1. s:runssintests 2. c:runscostests 3. t:runstantests 4. e:runsexptests 5. a:runsalltests
These options are mutually exclusive; only one may be chosen at a time. Is a bool the best choice, or would an enum be better?
3 Deliverables
You will need to turn in:
Thinking doesnt guarantee that we wont make mistakes. But not thinking guarantees that we will.
Leslie Lamport
1. math.c: This is your main program that determines whether exp, sin, cos, or tan will be run as well as contain your implementations. getopt must be used to handle commandline argu ments dictating which tests to run. The getopt options you must support are s to run sin tests, c to run cos tests, t to run tan tests, e to run exp tests, and a to run all tests. Your output for each function must match the form shown in Figures 1 and 2. Your compiled program must be called math.
2. Makefile: This is a file that will allow the grader to type make to compile your program. Typing make must build your program and running .math alone as well as with option flags must run your program.
CFLAGSWall Wextra Werror Wpedanticmustbeincluded.CCclangmustbespecified.
make cleanmustremoveallfilesthatarecompilergenerated.
2019 Darrell LongRegents of the University of California 8
make infermustruninferonyourprogram.Complaintsgeneratedbyinfershouldeither be fixed or explained in your README.
makeshouldbuildyourprogram,asshouldmake all.
3. README.md:Thismustbeinmarkdown.ThismustdescribehowtouseyourprogramandMakefile.
This is also where you put any explanations for complaints generated by infer.
4. DESIGN.pdf:ThismustbeaPDF.Thedesigndocumentshoulddescribeyourdesignforyourpro gram with enough detail that a sufficiently knowledgeable programmer would be able to replicate your implementation. This does not mean copying your entire program in verbatim. You should instead describe how your program works with supporting pseudocode. For this assignment pay extra attention to how you describe your representation of a Taylor series approximation in C.
5. WRITEUP.pdf: This must be a PDF. Your WRITEUP should be a discussion of the results for your experiments. This means analyzing the differences in the output of your implementations versus those in the math.h library. Include possible reasons for the differences between your implemen tation and the math.h version. Graphs can be especially useful in showing the differences and backing up your arguments.
4 Submission
We in science are spoiled by the success of mathematics. Mathematics is the study of problems so simple that they have good solutions.
Whitfield Diffie
To submit your assignment, refer back to assignment0 for the steps on how to submit your assign ment through git. Remember: add, commit, and push!
Your assignment is turned in only after you have pushed. If you forget to push, you have not turned in your assignment and you will get a zero. I forgot to push is not a valid excuse. It is highly recommended to commit and push your changes often.
2019 Darrell LongRegents of the University of California 9
Reviews
There are no reviews yet.