Rating: **4.1**/5.0 (7 Votes)

Category: Coursework

As far as I know there are no functions of the type you are searching for in the standard library. Here is one implementation of the:

For a fixed function f(x) to be integrated between fixed limits a and b. one can double the number of intervals in the extended trapezoidal rule without losing the benefit of previous work. The coarsest implementation of the trapezoidal rule is to average the function at its endpoints a and b. The first stage of refinement is to add to this average the value of the function at the halfway point. The second stage of refinement is to add the values at the 1/4 and 3/4 points.

A number of elementary quadrature algorithms involve adding successive stages of refinement. It is convenient to encapsulate this feature in a Quadrature structure:

Then the Trapzd structure is derived from this as follows:

The Trapzd structure could be used in several ways. The simplest and crudest is to integrate a function by the extended trapezoidal rule where you know in advance the number of steps you want. If you want 2^M + 1. you can accomplish this by the fragment:

with the answer returned as val. Here Ftor is a functor containing the function to be integrated.

Much better, of course, is to refine the trapezoidal rule until some specified degree of accuracy has been achieved. A function for this is:

Suppose $x_0,\ldots,x_n$ gives the partition of interval $[a,b]$ into $n$ evenly spaced subintervals, i.e. suppose $x_i = a+i\Delta x$ for $i=0,\ldots,n$ with $\Delta x = \frac

The code below computes the trapezoidal approximations to $\int_0^3 x^2\,dx$ first using the partition with 20 evenly spaced subintervals, then using 100 subintervals.

>> a=0; b=3; n=20; >> h=(b-a)/n; >> x=a:h:b; >> y=x.^2; >> h*trapz(y) ans = 9.0113 >> n=100; >> h=(b-a)/n; >> x=a:h:b; >> y=x.^2; >> h*trapz(y) ans = 9.0004

Not surprisingly using more subintervals yields a result closer to the exact value 9.

If x stores a *non* -evenly spaced partitition of $[a,b]$, and y holds the corresponding values of $f$, then trapz(x,y) computes the trapezodial approximation to $<\displaystyle \int_a^b f(x)\,dx>$ for the given partition. (This syntax can be used for evenly spaced partitions as well, but it's slightly less efficient than our first approach for that common special case.)

>> x=[0 .05 .1 .3 .32 .4 .55 .7 1 1.3 1.4 1.45 1.6 1.8 1.9 2 2.2 2.3 2.6 2.8 2.91 3]; >> y=x.^2; >> trapz(x,y) ans = 9.0217

It's not difficult to write your own implementation of the trapezoidal method in MATLAB. Feel free to do so as an exercise.

Although (as of this writing) MATLAB does not offer a builtin implementation of Simpson's rule in its most basic form, you can download one from the MATLAB Central file exchange. The code in simp.m is a barebones implementation of Simpson's rule.

function area = simp( fcn, a, b, n ) % Uses Simpson's rule to approximate the integral of the function % with handle fcn from a to b using n *pairs* of evenly spaced % subintervals. % Author: R Smyth % Version 1.0 8/22/2013 h = (b-a)/(2*n); % Subinterval spacing x = a:h:b; % Partition points y = fcn(x); % Function values at partition points i = 0:2*n; coeffs = 3+(-1).^(i+1); % Kludge to generate pattern of alternating 2s and 4s coeffs(1)=1; coeffs(end)=1; area = h/3 * sum(coeffs.*y); % Simpson's rule formula

octave:1> simp( @(x)x.^2, 0, 3, 1 ) ans = 9 octave:2> format long octave:3> simp( @(x)1./x, 1, exp(1), 3 ) ans = 1.00018687674140 octave:4> simp( @(x)1./x, 1, exp(1), 22 ) ans = 1.00000007582682

Notice that Simpson's rule produced the exact value for $\int_0^3 x^2\,dx$ using just one pair of subintervals. Is this a surprise? Shouldn't be. Simpson's rule uses quadratic approximations and thus produces exact values for quadratic integrands.

MATLAB has several builtin methods for numerically approximating integrals including the integral method (introduced earlier ) which uses global adaptive quadrature and the quadgk method (also introduced earlier ) which uses adaptive Gauss-Kronrad quadrature. Follow the links to the MATLAB docs for more information about these methods.

Mathematica's numerical integration routine NIntegrate[] (introduced earlier ) can be made to use trapezoidal approximations by setting the Method option to "TrapezoidalRule".

If you wish to control the partition you can implement the trapezoidal method as below.

Ordinarily you would use the numerical_integral() command (introduced earlier ) to compute numerical approximations to definite integrals in Sage. If you wish to implement the trapezoidal method there are of course many ways to do so, one of which is illustrated below.

To explore complex systems, physicists, engineers, financiers and mathematicians require computational methods since mathematical models are only rarely solvable algebraically. Numerical methods, based upon sound computational mathematics, are the basic algorithms underpinning computer predictions in modern systems science. Such methods include techniques for simple optimisation, interpolation from the known to the unknown, linear algebra underlying systems of equations, ordinary differential equations to simulate systems, and stochastic simulation under random influences. Topics covered are: the mathematical and computational foundations of the numerical approximation and solution of scientific problems; simple optimisation; vectorisation; clustering; polynomial and spline interpolation; pattern recognition; integration and differentiation; solution of large scale systems of linear and nonlinear equations; modelling and solution with sparse equations; explicit schemes to solve ordinary differential equations; random numbers; stochastic system simulation.

School of Mathematical Sciences

North Terrace Campus

Up to 3.5 hours per week

Available for Study Abroad and Exchange

MATHS 1012 and (COMP SCI 1012, 1101, MECH ENG 1100, 1102, 1103, 1104, 1105, C&ENVENG 1012)

MATHS 2102 or MATHS 2201

To explore complex systems, physicists, engineers, financiers and mathematicians require computational methods since mathematical models are only rarely solvable algebraically. Numerical methods, based upon sound computational mathematics, are the basic algorithms underpinning computer predictions in modern systems science. Such methods include techniques for simple optimisation, interpolation from the known to the unknown, linear algebra underlying systems of equations, ordinary differential equations to simulate systems, and stochastic simulation under random influences. Topics covered are: the mathematical and computational foundations of the numerical approximation and solution of scientific problems; simple optimisation; vectorisation; clustering; polynomial and spline interpolation; pattern recognition; integration and differentiation; solution of large scale systems of linear and nonlinear equations; modelling and solution with sparse equations; explicit schemes to solve ordinary differential equations; random numbers; stochastic system simulation.

The full timetable of all activities for this course can be accessed from Course Planner .

Demonstrate understanding of common numerical methods and how they are used to obtain approximate solutions to otherwise intractable mathematical problems.

Apply numerical methods to obtain approximate solutions to mathematical problems.

Derive numerical methods for various mathematical operations and tasks, such as interpolation, differentiation, integration, the solution of linear and nonlinear equations, and the solution of differential equations.

Analyse and evaluate the accuracy of common numerical methods.

E. Kreyszig, Advanced engineering mathematics, 9th edition, Wiley, 2006.

A. Greenbaum & T. P. Chartier, Numerical methods, Princeton University Press, 2012.

W. Cheney & D. Kincaid, Numerical mathematics and computing, Thomson, 2004.

D. P. O'Leary, Scientific computing with case studies, SIAM, 2008.

D. M. Etter, Engineering problem solving with Matlab, Prentice-Hall, 1993.

W. H. Press et al, Numerical recipes in [C, Fortran. ], Cambridge University Press, c1996-1999.

Lecture recordings and screencasts, MapleTA exercises, partial lecture notes, assignments, tutorial exercises, and course announcements will be posted on MyUni.

This course uses a variety of methods for delivery of the course material.

Some lecture material is delivered using online screencasts together with interactive Maple TA exercises and quizzes. Other lecture material is delivered in traditional face-to-face lecture format.

Tutorials are held fortnightly. In these classes, you will complete short quizzes and work on tutorial problems that aim to enhance your understanding of the lecture material and ability to solve theoretical problems. You are encouraged to attempt the problems before the tutorial and to complete all the remaining problems afterwards.

Practicals are held fortnightly, alternating with tutorials. In these classes, you will use Matlab to implement numerical algorithms developed in lectures. Practical work must be submitted to show that you have completed the session.

Assignments are set fortnightly. In the assignments, you are usually asked to write a Matlab program to solve a mathematical problem and present your results in a written report. Questions about theoretical aspects of the problem may also be asked.

The information below is provided as a guide to assist students in engaging appropriately with the course requirements.

Tutorial quizzes and MapleTA exercises will be set throughout the course. They are of equal weight.

You will need to submit both electronic and hardcopy components for each assignment. The electronic component must be submitted according to the assignment instructions. It will be marked electronically and the result added to your hardcopy mark. The hardcopy component must be submitted to the designated hand-in boxes within the School of Mathematical Sciences with a signed cover sheet attached.

Late assignments will not be accepted. Students may be excused from an assignment for medical or compassionate reasons. Documentation is required and the lecturer must be notified as soon as possible.

Assignments will have a two week turn-around time for feedback to students.

Grades for your performance in this course will be awarded in accordance with the following scheme:

M10 (Coursework Mark Scheme)

Further details of the grades/results can be obtained from Examinations .

Grade Descriptors are available which provide a general guide to the standard of work that is expected at each grade level. More information at Assessment for Coursework Programs .

Final results for this course will be made available through Access Adelaide .

The University places a high priority on approaches to learning and teaching that enhance the student experience. Feedback is sought from students in a variety of ways including on-going engagement with staff, the use of online discussion boards and the use of Student Experience of Learning and Teaching (SELT) surveys as well as GOS surveys and Program reviews.

SELTs are an important source of information to inform individual teaching practice, decisions about teaching duties, and course and program curriculum design. They enable the University to assess how effectively its learning environments and teaching practices facilitate student engagement and learning outcomes. Under the current SELT Policy (http://www.adelaide.edu.au/policies/101/) course SELTs are mandated and must be conducted at the conclusion of each term/semester/trimester for every course offering. Feedback on issues raised through course SELT surveys is made available to enrolled students through various resources (e.g. MyUni). In addition aggregated course SELT data is available.

This section contains links to relevant assessment-related policies and guidelines. all university policies can be obtained from: http://www.adelaide.edu.au/policies/

Students are reminded that in order to maintain the academic integrity of all programs and courses, the university has a zero-tolerance approach to students offering money or significant value goods or services to any staff member who is involved in their teaching or assessment. Students offering lecturers or tutors or professional staff anything more than a small token of appreciation is totally unacceptable, in any circumstances. Staff members are obliged to report all such incidents to their supervisor/manager, who will refer them for action under the university's student’s disciplinary procedures.

The University of Adelaide is committed to regular reviews of the courses and programs it offers to students. The University of Adelaide therefore reserves the right to discontinue or vary programs and courses without notice. Please read the important information contained in the disclaimer .

· Midpoint Rule (Rectangular Rule)

· Euler-Maclaurin Formula

· Richardson Extrapolation

· Trapezoidal Rule

· Simpson�s Rule

There are a number of numerical methods that can be used to approximate an integral. Some of them are already experienced in the course notes: midpoint rule, Euler-Maclaurin, trapezoidal and Simpson rule. In this section we consider again in a more broader way the trapezoidal rule, the Richardson extrapolation, Simpson�s rules and many other examples not covered in the course. Numerical integration methods are also called quadrature formulas. The term quadrature means the process of finding a square with the same area as the area enclosed by an arbitrary closed curve.

An interactive Matlab code for the *Midpoint rule* is given in this application. The algorithm is applied again on the *f306* example function. The student has the chance to experience the *Midpoint* rule for different functions (by changing the function equation in the f306.m file), different number of strips and different abscissas extremes. The Matlab code of this algorithm is given below:

*%function I=midpnt(a,N,b)**% Use Midpoint rule to integrate f(x) between a and b.**% a: lower bound, b: upper bound, N: number of strips*

*while 1**fprintf( ' \n Give the number of strips')**fprintf( ' \n N = '); N = input('');**if N > 2; break; end**fprintf( ' \n N should be greater than 2! ')**fprintf( ' \n Repeat your input for N! ')**end*

*while 1**fprintf( ' \n The absciss extremes a and b (a < b)')**fprintf( ' \n a = '); a = input('');**fprintf( ' b = '); b = input('');**if a < b; break; end**fprintf( ' \n a should be less than b! ')**fprintf( ' \n Repeat your input for a and b! ')**end*

*Below, the algorithm was run with the following input values: N = 1000 strips, a = -2 and b = 6.7. The calculated value of the integral approximation is I = 19.13.*

*Give the number of strips**N = 1000*

*The absciss extremes a and b (a < b)**a = -2**b = 6.7*

The ploted curve of the function looked like the one in the figure below:

With a *midpnt* function like the one given in the course notes (midpnt.m ) the Richardson extrapolation can be executed like in midpnt_richard.m and midpnt_richard.m (see class notes).

Next the *Trapezoidal* rule is approached.

A graphical approach of the trapezoidal integration rule is given in this example. The integral value is calculated at the same time with hatching the area which gives the integral value. The code is given below:

% function I = trapez_g(f_name, a, b, n)

% function trapez_g: same as trpez_n except this

% will plot the function.

*while 1**fprintf( ' \n Give the number of strips')**fprintf( ' \n N = '); n = input('');**if n > 2; break; end**fprintf( ' \n N should be greater than 2! ')**fprintf( ' \n Repeat your input for N! ')**end*

*while 1**fprintf( ' \n The absciss extremes a and b (a < b)')**fprintf( ' \n a = '); a = input('');**fprintf( ' b = '); b = input('');**if a < b; break; end**fprintf( ' \n a should be less than b! ')**fprintf( ' \n Repeat your input for a and b! ')**end*

For the f(x) = x-cos(x) given in f306.m the following results and graphics are generated:

Give the number of strips

N = 1000

The absciss extremes a and b (a < b)

a = -6

b = 9.8

The initial function proposed by Nakamura trapez_n.m has been made interactive. The function which implements the trapezoidal rule without any graphical feature is given in trapez_g.m while a vector variables approach is written in trapez_v.m

Rctrapez_gen.m Recursive trapezoidal integration of order *n* extended from the function given in Rctrapez.m proposed in "NUMERICAL METHODS Using MATLAB" by John H. Mathews and Kurtis D. Fink. The initial recursive function is given below:

*function T=rctrap(f,a,b,n)**% recursive trapezoidal integration**%Input - f is the integrand input as a string 'f'**% - a and b are upper and lower limits of integration**% - n is the number of times for recursion**%Output - T is the recursive trapezoidal rule list*

*% NUMERICAL METHODS: MATLAB Programs**%(c) 1999 by John H. Mathews and Kurtis D. Fink**%To accompany the textbook:**%NUMERICAL METHODS Using MATLAB,**%by John H. Mathews and Kurtis D. Fink**%ISBN 0-13-270042-5, (c) 1999**%PRENTICE HALL, INC.**%Upper Saddle River, NJ 07458*

Next this function is invoked in the Rctrapez_gen.m interactive shell and gives the following output for N = 10 recursions, a = -2 and b = 6.7 input values.

Give the number of strips

N = 10

The absciss extremes a and b (a < b)

a = -2

b = 6.7

Columns 1 through 7

18.2777 22.4181 19.6950 19.2630 19.1634 19.1390 19.1329

Columns 8 through 11

19.1314 19.1310 19.1309 19.1309

Next the *Simpson* rule is experienced for students.

The interactive shell given in mat_simp1.m invokes the function written in simp1.m. This function executed *Simpson* rule using a vectorwise approach. Again the algorithm is run for the test function given in f(x) = x-cos(x) given in f306.m. Below are given the code of this implementation and then the input variables and the output results. The *time* value gives the time spent on runing this algorithm, while *flops* variable gives the total number of floating point operations.

*while 1**fprintf( ' \n Give the number of strips')**fprintf( ' \n N = '); N = input('');**if N > 2; break; end**fprintf( ' \n N should be greater than 2! ')**fprintf( ' \n Repeat your input for N! ')**end*

*while 1**fprintf( ' \n The absciss extremes a and b (a < b)')**fprintf( ' \n a = '); a = input('');**fprintf( ' b = '); b = input('');**if a < b; break; end**fprintf( ' \n a should be less than b! ')**fprintf( ' \n Repeat your input for a and b! ')**end*

And the simp1.m function:

*function q = simp1(func,a,b,m)**% Implements Simpson's rule using vectors.**% Example call: q = simp1(func,a,b,m)**% Integrates user defined function func from a to b, using m divisions*

And the inputs/outputs of this algorithm:

Give the number of strips

N = 100

The absciss extremes a and b (a < b)

a = -2

b = 6.7

==============RESULTS================ n integral value

2 23.798293193

4 18.787258916

8 19.119000418

16 19.130191252

32 19.130812410

64 19.130850154

time= 0.02 secs flops= 870

Then the interactive shell given in mat_simp2.m invokes the function written in simp2.m. This function executed *Simpson* rule using a for loop approach. The algorithm is run for the test function given in f(x) = x-cos(x) given in f306.m. Below are given the code of this implementation and then the input variables and the output results. The *time* value gives the time spent on runing this algorithm, while *flops* variable gives the total number of floating point operations.

*while 1**fprintf( ' \n Give the number of strips')**fprintf( ' \n N = '); N = input('');**if N > 2; break; end**fprintf( ' \n N should be greater than 2! ')**fprintf( ' \n Repeat your input for N! ')**end*

*while 1**fprintf( ' \n The absciss extremes a and b (a < b)')**fprintf( ' \n a = '); a = input('');**fprintf( ' b = '); b = input('');**if a < b; break; end**fprintf( ' \n a should be less than b! ')**fprintf( ' \n Repeat your input for a and b! ')**end*

And the simp2.m function:

*function q = simp2(func,a,b,m)**% Implements Simpson's rule using for loop.**% Example call: q = simp2(func,a,b,m)**% Integrates user defined function**% func from a to b, using m divisions*

And the inputs/outputs of this algorithm:

Give the number of strips

N = 100

The absciss extremes a and b (a < b)

a = -2

b = 6.7

==============RESULTS================ n integral value

2 23.798293193

4 18.787258916

8 19.119000418

16 19.130191252

32 19.130812410

64 19.130850154

time= 0.04 secs flops= 897

A *Simpson* rule approach for two variable functions is given in mat_simp2var.m. which uses the function written in simp2v.m.

Another *Simpson* rulefunction proposed in "NUMERICAL METHODS Using MATLAB" by John H. Mathews and Kurtis D. Fink is given in simprl.m

Another *Extended Simpson* rulefunction is given in simp_n.m and a vector approach of the same algorithm is given in simp_v.m

Other integration method examples of numerical integration algorithms that are at the moment beyond course's scope are:

· Gaussian integration (Gauss quadrature of order n) fgauss_integr.m. gauss_q_integr.m. gaussint_mat4.m

· Gaussian integration for a two variable function (Gauss quadrature of order n) gauss2v.m

· Newton closed integration formula of order n newt_integr.m

· Algorithms which use already implemented Matlab functions based on numerical methods:

- Matlab *gaherm* function that implements Gauss-Hermite integration method: mat_gaherm_integr.m

- Matlab *galag* function that implements Gauss-Laguerre integration method: mat_galag_integr.m

- Matlab *quad* function that numerically evaluates the integral (low order method): mat_quadrature.m

The student is free to experience more of these methods. A reminder that the algorithms given as functions can be run in the command line mode in Matlab's command window. For example for the fgauss_integr.m example the algorithm can be run for function given in f306.m with the following command line:

» fgauss_integr('f306',-5, 2, 16)

and as a result the following output is given in the same command window:

All the above algorithms give the same result for this case example.

Valentin Muresan, Dublin City University, muresanv@eeng.dcu.ie