代写Math 104A、Matlab/python程序

日期：2022-11-19 04:03

Math 104A Computation Homework # 2

1 Main problem

We have studied some quadrature rules and their theoretical rates of convergence (or error terms). In this

homework, you compute and report actual rates of convergence of some quadrature rules when they are

implemented for selected integrands. Also, you interpret the result in your own language.

(Quadrature rules) Your program computes the following three quadrature rules for integrals on [?1, 1.1]

with respect to constant weight w(x) ≡ 1.

1. composite trapezoidal rules (“Trapezoidal” below) (textbook p. 481)

2. composite Simpson’s rule (“Simpson” below) (textbook p. 484)

3. composite Gaussian quadrature with three nodes (“Gauss3” below), i.e., n = 2 so we have Gaussian

nodes x0, x1, and x2 for each subinterval (textbook p. 496 provides a version on a single interval

[?1, 1]. You will need to do some work to obtain a composite version from it. In particular, you

will need to apply change of interval.)

(Integrands) Implement these quadrature rules for the following three integrands: consine, a C1– func-

tion, and (a scaled) normal distribution density.

f (x) = cos(pix), g(x) =x2 (x ≥ 0)?x2 (x < 0) , h(x) = exp(?x22 ). (1)

(Experiment) Compute the errors caused by each quadrature rule for each integrand (i.e., the abso-

lute value of the difference between approximation and the true value) implemented with a number of

different resolutions.

– Here, the resolution means the number of subintervals involved in composite rules. (2i for i =

0, 1, 2, · · · , 9 is a good option.) For example, if you decide to use 10 different resolutions, you will

have total 3 × 3 × 10 errors from each combination of a method, an integrand, and a resolution.

– To measure the error, you will need true values of the integrals. Think about how you can obtain

those. If you have no idea, see the Section 2 for tips.

(Report) Report a table that contains the rate of convergence of each quadrature rule for each integrand

as the resolution gets one step higher.

– For example, if you decide to use 10 different resolutions, you will have 9 convergence rate for

each combination of a quadrature rule and an integrand because you have 9 steps of refinement. If

you are not sure what is going on, read ‘how to compute the rate of convergence’ below first, and

possibly all other instructions, then come back here.

1

– How to compute the rate of convergence? We will discuss this in class in a bit more detail. But

the summary is the same. Suppose you have implemented a numerical method using two different

resolutions N1 and N2, where N2 is the higher resolution than N1, say N2 = 2N1. If you have

obtained errors e1 and e2 using the resolutions N1 and N2 respectively, the rate of convergence

based on these two results is given by

r =

log(e1/e2)

log(N2/N1)

.

See the class notes for details.

(Summary)

data

– Target integral information (interval [?1, 1.1]), Three integrands listed above

parameters to choose

– Resolutions (i.e., an array of the number of subintervals for composite rules; 2i for i = 0, 1, 2, · · · , 9

recommended)

output

– MethodIntegrand.csv: A csv file (comma separated values) that contains the rates of conver-

gence so that one can easily compares them across integrands for each method used. Ideally, we

want something like this:

#Subintervals Trapezoidal Simpson Gauss3

consine C1 function normal dist consine C1 function normal dist consine C1 function normal dist

However, csv files work the best when the data is complete and in the simplest rectangular form

like the following. Thus, feel free to submit in this form:

Method Trapezoidal Trapezoidal Trapezoidal Simpson Simpson Simpson Gauss3 Gauss3 Gauss3

Integrand consine C1 function normal dist consine C1 function normal dist consine C1 function normal dist

– IntegrandMethod.csv: A csv file (comma separated values) that contains the rates of conver-

gence so that one can easily compares them across methods for each integrand.

For the same reason as above about csv files, feel free to report csv-friendly table:

2

#Subintervals consine C1 function normal dist

Trapezoidal Simpson Gauss3 Trapezoidal Simpson Gauss3 Trapezoidal Simpson Gauss3

1

2

4

Integrand consine consine consine C1 function C1 function C1 function normal dist normal dist normal dist

Method Trapezoidal Simpson Gauss3 Trapezoidal Simpson Gauss3 Trapezoidal Simpson Gauss3

(Submission): Please, submit your result as indicated below. Part of the credit is assigned to submitting

things as instructed.

– Source code: (.m, .py, or .ipynb; in particular, NO .pdf)

– two output files: MethodIntegrand.csv, IntegrandMethod.csv

– Interpretation: interpretation.txt (A text file that describes how you interpret the result)

– (Optional): ExtraExperiment.txt – Interpretation on an additional experiment: Try the same

experiment with for integrals from -1 to 1. (You must change the true integrals accordingly.) What

do you observe? What is your interpretation?

2 Tips for this homework

2.1 Obtain true integrals

Although there should be many different ways, using Wolfram Alpha is one good way. Mathematica is

case-sensitive. So be careful of capitalization. Also, use brackets and curly brackets correctly.

Type in the Wolfram Alpha: N[Integrate[Exp[-x*x/2], {x, -1, 1.1}], 20]

This tells Wolfram Alpha that find

∫ 1.1

1 exp(?x2/2)dx up to 20 decimal digits

20 decimal digits will do its job as the true value.

Do similar things for other integrands.

2.2 Tools for reporting

As I showed in class, it is very good to know data processing tools when you report the results. However,

it takes time to learn those tools. Thus, consider learning and using them if you can budget some time. If

you bring this up in office hours, I am happy to discuss it.

(Python) Use the pandas package (along with numpy). Among many great tools that belong to it,

the most relevant things include:

multiindex, pivot, stack, unstack, sort_values, to_csv, to_excel, etc.

3

(Matlab) Use table object. Note that the way Matlab approaches tables are quite different from

how pandas treats data. I believe they are developed based on different philosophy. The most

relevant tools include:

splitvars, mergevars, stack, unstack, writetable, etc.

3 Tips for programming

3.1 Possible time savers

Here I share some issues/tricks you may want to know. Of course, your code may not experience any of

these. But I am trying to save your time in case you decide to use the same ideas/tools and/or have the

same issues.

Common

You can use for loop to conduct many different experiments in a single line. To do that, you can

collect functions in a list (Python) or a cell array (Matlab). (See below for a bit more details for

each language.) For example, if you use Python (Matlab), you can create a list (cell array) named

fn of functions and a list (cell array) mthd of quadrature rules. Then, you can use for loop to

conduct all experiments: “for i from 0 to 2 (from 1 to 3 for Matlab), and for j from 0 to 2 (from 1

to 3 for Matlab), compute the numerical integration of fn[j] using mthd[i].

Python

Depending on the version of Python, 5/3 may “equal” 1 instead of 1.66666 · · · . You can avoid this

by typing “5./3”.

Depending on how you construct the resolutions, they may be considered float64 to Python

rather than int (integer), which may cause error when Python carries out composite quadrature

rules. Consider data type casting by NN = NN.astype(int), where NN is the name of the numpy

array containing resolutions.

(List of functions) To collect functions, you simply gather them by creating a Python list of the

function names. For example, you may want something like the following:

fn = [f, g, h]

mthd = [Trapez, Simpson, Gauss3].

If needed, Google “list of functions.”

Matlab

(cell array) You will want to use cell array for this homework. If you don’t know yet, study bare

minimum of cell arrays. (Again, Google is a good friend). a quick piece of basic is that you have

to use curly brackets to access its entry: e.g., fn{1}. if you use fn(1), Matlab will throw an error.

(cell array of functions) To collect functions, you have to collect so-called function handles of

the Matlab functions. Anonymous functions’ names themselves are function handles (things you

define using the syntax f = @(x) 2*x; ). But the functions you define using function output

= f(x) .... end, you have to refer to it by @f. And you have to collect them in a cell array.

4

For example, if you have defined the integrands using anonymous functions and quadrature rules

using regular functions, you may need something like the following: