Complete the following problems with a single script and custom functions as appropriate. Important lines of code should include descriptive comments.All graphs/figures/sketches should include a grid, labels, legend (if necessary), and the appropriate fontsize/linewidth/markersize.
AERO 300 Laboratory 6 Numerical Differentiation and Newton-Cotes Integration (Quadrature) The pre-lab assignment (Section 4) is due at the beginning of this lab (can be hand-written or typed). The lab (Section 5) will be due before your next lab section (.zip file submitted to PolyLearn). 1 Objectives This lab covers concepts and methods used for approximating the derivative of a function using finite difference equations. In addition, this lab covers numerical techniques for integrating a function with Newton-Cotes equations. In both cases, this lab assumes the closed form of the function is unknown; the function is only known at discrete points. In this lab, you will: • Develop code to approximate the first and second derivative of data representing the displacement of a falling object. • Use the Rectangle, Trapezoid and Simpson’s rule methods of quadrature to integrate data from an accelerometer to calculate speed and displacement. • In both cases, the data sets include measurement noise. 2 Introduction Quadrature is the computation of a univariate definite integral. It can refer to either numerical or analytic techniques; one must gather from context which is meant. The term refers to the geometric origin of integration in determining the area of a plane figure by approximating it with squares. There are several reasons for carrying out numerical integration. First, the integrand f(x) may be known only at certain points, such as obtained by sampling. Some embedded systems and other computer applications may need numerical integration for this reason. Even if a formula for the integrand may be known, it may be difficult or impossible to find an antiderivative which is an elementary function. An example of such an integrand is ?(?) = ?−? 2 , the antiderivative of which cannot be written in elementary form. For some functions, it may be possible to find an antiderivative symbolically, but it may be easier to compute a numerical approximation than to compute the antiderivative. That may be the case if the antiderivative is given as an infinite series or product, or if its evaluation requires a special function which is not available. Numerical Differentiation is the “inverse” process of numerical integration. All the reasons to take a numerical derivative carry over from those for numerical integration. In addition, many equations which govern the motion of fluids and structure are most easily solved numerically using finite difference equations in multiple dimensions (space and time). We will return to this use of finite difference equations later in the course. In the meantime, we will use finite difference schemes to perform straight ahead estimate of the derivatives of data. 2.1 Numerical Differentiation As discussed in class, numerical differentiation is accomplished by using the Taylor series expansion to develop a finite difference equation. With respect to errors in the approximation of the derivative, different finite difference formulas were also developed. Table 1 below summarizes those formulas. Two-Point Forward-Difference ?′(?) ≈ ?(? + ℎ) − ?(?) ℎ Three-Point Centered Difference ?′(?) ≈ ?(? + ℎ) − ?(? − ℎ) 2ℎ Three-Point Centered-Difference (2nd derivative) ? ′′(?) ≈ ?(? − ℎ) − 2?(?) + ?(? + ℎ) ℎ2 Oftentimes you are given data and you want to take the derivative of the data. In this case, you may know what the function should look like, but the data contains noise. In this lab, you will investigate how the presence of noise impact the numerical differentiation of data. 2.2 Numerical Integration – Quadrature As with numerical differentiation, in some instances, you have measurement data and you want to calculate the integral of the data. Again, you may know what the function you are integrating looks like, but the data contains noise. In addition to calculating the total integral using composite Newton-Cotes methods, we can compute the integral at each time step, store that information so we can then plot the integral of a function versus time. For example, using an accelerometer, you can measure the acceleration of a falling object every 0.1 seconds. If you integrate the measurement, you would then have the speed of the object. If you integrate again, you would then have the displacement. We call this “propagation” of the object. In this case, we know from physics what the speed and displacement of the object should be as a function of time. We can use Newton-Cotes quadrature methods to propagate the speed and displacement by applying any of the Newton-Cotes formulas to integrate the speed or displacement over one time step and add the calculated integral to the previous value. For example, if the speed and displacement of an object are 0 ? ? and 3 ? at ? = 0, and we measure the acceleration as ?(0) = 9.85 ? ?2 and ?(0.1) = 9.75 ? ?2 , we can then estimate the speed and displacement using the trapezoid rule as: ?(0.1) = ?(0) + ℎ 2 (?(0) + ?(0.1)) = 0 + 0.1 2 (9.85 + 9.75) = 0.98 ? ? ?(0.1) = ?(0) + ℎ 2 (?(0) + ?(0.1)) = 3 + 0.1 2 (0 + 0.98) = 0.98 ? ? = 3.049 ? In general, we could continue this process for each time step we measure the acceleration and use the following propagation equations: ??+1 = ?? + ℎ 2 (?? + ??+1) ??+1 = ?? + ℎ 2 (?? + ??+1) In these equations, the ? represents the time step, and ℎ is the time between samples. As described below, we can use the same process using any of the quadrature methods. Rectangle Rule In integral calculus, the rectangle method uses an approximation to a definite integral, made by finding the area of a series of rectangles. Either the left or right corners, or top middle of the boxes lie on the graph of a function, with the bases run along the x-axis. Each of these schemes is shown below in figure 1. Figure 1. Left, midpoint, and right rectangular Integration schemes To propagate a measurement using the left-rectangular rule, we can use the following formula: ? ?+1 = ? ? + ℎ?(? ? ) (1) where we have taken ∫ ?(?)?? ??+1 ?? ≈ ℎ?(??) Trapezoid Rule The trapezoid rule uses trapezoids, as opposed to rectangles to approximate an integral. An advantage of the trapezoid rule is that the sign of the error of the approximation is easily known. An integral approximated with this rule on a concave-up function will be an overestimate because the trapezoids include all the area under the curve and extend over it. Using this method on a concave-down function yields an underestimate because area is unaccounted for under the curve, but none is counted above. If the interval of the integral being approximated includes an inflection point, then the error is harder to identify. Figure 2 – Trapezoid rule scheme To propagate a measurement using the trapezoid rule, we can use the following formula: ??+1 = ?? + ℎ 2 (?(??) + ?(??+1)) Simpson’s Rule Simpson's rule is another method for approximating definite integrals. It works in a similar way to the trapezoidal rule except that the integrand is approximated to be a quadratic rather than a straight line within each subinterval. To propagate a measurement using Simpson’s rule, we can use the following formula: ??+1 = ??−1 + ℎ 3 (?(??+1) + 4?(??) + ?(??−1)) 2.3 Accuracy Each method has an accuracy associated with it. The rectangle rule, for example, is accurate for a zeroth order polynomial. Using the approximation on any higher order polynomial results in an error; if the rectangle rule was used to approximate a first order polynomial the error would be on the order of h. In other words, the rectangle rule is a first order accurate method. A similar argument can be made for the trapezoid and Simpson's rules. 3 Help Potentially useful terms for this lab: • for • linspace • ones • plot • grid • hold • legend • title • xlabel • ylabel 4 Pre-lab Assignment NOTE: This must be done before lab. If it is not completed, you will not be allowed into lab. Figure 5: Simpson's approximation Andy G Cristales Use the concept of left-rectangular propagation from equation (1) to show that the total integral from ? to ? is given as: ∫ ?(?)?? ??=? ?1=? ≈ ℎ?(?) + ℎ ∑ ?(??) ?−1 ?=2 5 Lab Assignment Complete the following problems with a single script and custom functions as appropriate. Important lines of code should include descriptive comments. You are encouraged to work in groups, however the code you submit must be your own. All graphs/figures/sketches should include a grid, labels, legend (if necessary), and the appropriate fontsize/linewidth/markersize. 1. Numerical Differentiation a. Load the data file dispData.mat. The data comes from measuring the displacement of a falling object versus time. b. Use the two-point forward-difference and the three-point centered-different schemes to estimate the derivative of the given data. c. Use the three-point centered-difference scheme to determine the 2nd derivative of the data. d. Assuming the initial displacement and speed of the object are zero meters and zero meters per second, plot the following (plot the displacements, speeds and accelerations on separate subplots): i. The expected displacement, speed and acceleration versus time. ii. The measured displacement and estimated speed and acceleration versus time. e. What is the effect of taking a numerical derivative of noisy data? Can you think of a better method for estimating the acceleration of gravity from this data? 2. Numerical Integration – Quadrature a. Load the data file accelData.mat. The data comes from measure tin the acceleration of a falling object versus time. b. Use the rectangular, trapezoid and Simpson’s rules to integrate the data twice, hence propagating the speed and displacement of the object. c. Assuming the initial displacement and speed of the object are zero meters and zero meters per second, plot the following (plot the