Riassunto esame Elementi di calcolo numerico, Prof. Antonini Mario, libro consigliato Numerical methods for engineers, Chapra, S., Canale, R

Introduction

Numerical analysis: find a numerical solution to a mathematical problem. This process divides into 5 steps:

• Real problem
• Model
• Mathematical model
• Solve the mathematical model
• Interpret the solution →

Sometimes the exact solution might not be available we need an approximation by using numerical method Other times we could have a system of equations or equation of a high grade, which both require high computation to solve.

In the formulation of a real problem into a numerical method we introduce an error, which can have different sources/natures:

• Modelling errors
• Physical measurement errors
• Machine representation and arithmetic errors (computers are finite machines and can’t represent infinite numbers)
• Mathematical approximation errors

The point is to guarantee that the error is under control and the different errors are of the −16∼ 10 same order (→ this applies for example in the case in which an error is and the −4∼ 10 other is. The effort required for the first one is wasted because the second one is huge compared to it).

Notation

Absolute error on a scalar: = | − | this does not take into account the order of magnitude of the value under examination −3 −3 | | → | 1. 000 − 1. 001 | = 0. 001 2 · 10 − 1 · 10 = 0. 001 (compare and in the first case the error is negligible, while in the second one it is not negligible).

ε True relative error on a scalar: | −| ε = | | this, compared to the absolute error, takes into account the order of magnitude of the value under examination.

Mathematical approximation errors: → → In most cases the true value is not known the error is not computable we approximate the error: | | ε = | | If we apply a numerical method, it will produce better and better approximations: | − | ε = | | and the iterative numerical method will stop once the approximation goes under a certain threshold: ε < ℎℎ

Computer number representation

The fundamental unit is a word which consists of bits. How many digits are available defines → a representation system only a limited range of integers can be represented. E.g. in a 16 bit system the range of the integers that can we store is -32768 to 32767.

In general, in an n-bit system the range is −1 −1 [ ] − 2, 2 − 1.

Integers are stored by specifying the word size (how many bits are available). Then the first bit indicates the sign:

• 0 = positive
• 1 = negative

Then the remaining bits are used to store the absolute value of the number. E.g.: 7 5 3 2 0 ( 10101101 ) = 2 + 2 + 2 + 2 + 2 = 128 + 32 + 8 + 4 + 1 = (173) 2 10.

Numerical analysis

Find a numerical solution to a mathematical problem.

Floating point representation

Floating point representation: ± · where: - it’s called mantissa or significant. It may be positive or negative …) - it’s the base (2, 8, 10, 16, - it’s the integer exponent of the base used.

NB: this representation is normalized: there’s only one digit at the left of the dot: −3 −6 0. 00345 · 10 = 3. 45 · 10.

The floating-point representation in computers is usually used with binary base: → - each digit (bit) is either 0 or 1 since the representation is normalized the first bit in is always 1: ± (1 + ) · 2 → - is the fractional part of the mantissa the first bit of is always 0 −5 (1101. 1) · 10.

E.g.: normalize the binary number: 2 −5 −2 (1101. 1) · 2 = 1. 1011 · 22 NB: we only have to store and.

Most computers adopt IEEE double-precision, eight bytes (64 bits) are used to represent floating-point numbers:

• 1 bit for the sign
• 52 bits for the fractional part of the mantissa
• 11 bits for the exponent as an integer (ranging from -1023 to 1024)

Example with normalized floating-point system

Example: = 10 assume you have a normalized floating-point system with base using 5 digits (1 for the sign, 2 for the exponent and 2 for the mantissa). For simplicity, assume that one of the exponent digits is used for its sign. A general representation of the number following normalization would be 0 0, · 101 1 2.

• What’s the largest possible positive number that can be represented?

99. 9 · 10 = 0 = 0 with both and that means both and are positive 1 0 1 0.

• What’s the smallest possible positive number that can be represented?

−9 1. 0 · 10 = 0 = 1 with and which means that is positive and is negative 1 0 1 0.

In order to obtain the respective negatives you just change in the first case to obtain the 19 − 9. 9 · 10 smallest possible negative number ( ) and in the second case to obtain the 1 −9 − 1. 0 · 10 largest possible negative number ( ).

∈ () For any number we denote its representation in the floating-point system. Assume that digits are available for the fractional part of the mantissa. There are 2 options to store the number:

• Chopping: → () = ±, ... · we simply stop after we are out of 1 2 3 +1 digits available
• Rounding: < () = ±, ... · - if: 2+2 1 2 3 +1 ( ) ≥ () = ±, ... + 1 · - if: 2+2 1 2 3 +1

The upper bound of the relative error is called machine precision and it can be shown that it is: 1 − |−()| ≤ : = ε for chopping || 1 − |−()| 1 ≤ : = ε for rounding || 2 and this means that rounding has slightly higher machine precision. Note that it is an upper bound for the machine error, but still, rounding is overall better. E.g. keep 4 digits: 1. 5636 → 1. 563 chopping: 1. 5636 → 1. 564 rounding:

Arithmetic operations

Arithmetic operations:

• Sum of 2 real numbers:
• The numbers are expressed so that they have the same exponent
• The mantissas are added together
• The number is normalized
• The float of the result is computed
• Multiplication and division of 2 real numbers:
• Do the product or division of the mantissas of the two numbers
• Sum or subtract the exponents
• The number is normalized
• The float of the result is computed

Main sources of errors in arithmetic operations:

• Subtracting cancellation: subtracting 2 nearly equal numbers can rise errors
• Large computations: an individual round-off error could be small, but the cumulative effect over the course of many operations can be not negligible
• Loss of significant digits: adding a very large number to a very small one can lead into a loss of precision

E.g.: 8 1 = 1. 53 · 10 = 2. 01 · 10 Sum and 8 8 8 ( ) ( ) + = 1. 53 · 10 + 0. 000000201 · 10 = 1. 530000201 · 10.

Conditioning and stability

A problem is said to be well-conditioned if to small variations on the data correspond small variations on the results. On the contrary, if the small variations on the data lead to big variations on the results, the problem is ill-conditioned.

E.g.: 4 *( − 1 ) = 0 ⇒ = 14 −8 * −2 ( − 1 ) = 10 ⇒ = 1 + 10 0.

In the second one a perturbation on the right hand side is introduced (instead of we have −8 10 1). This generates a perturbation also on the result, which is not anymore but −2 1 + 10. There are 6 orders of magnitude in difference between the perturbation introduced and the perturbation on the result. This means that the problem is ill-conditioned.

An algorithm is said to be stable if it amplifies a little the roundoffs errors introduced during the computation.

E.g.: consider the problem of calculating the series { } = φ, = 0, 1, 2,... of integer powers of the golden mean 5 −1 φ = 2.

This can be done in 2 ways:

• Set = 1 and calculate 0 = φ · = 1, 2,... for −1
• Set = 1 = φ and compute 0 1 = − = 2, 3,... for −2 −1

Both methods are mathematically equivalent, but the second one is said to be unstable, since it gives entirely wrong results for large values of. This source of error in the second → method is introduced from the loss of significant digits adding/subtracting large and small numbers together.

Computational complexity

The speed of a computer is usually measured in FLOPS (FLoating-point OPerations per Second).

The computational complexity is the rate at which the total number of FLOPS required to () solve a problem of size grows as grows. Denoting the number of operations required to solve a problem, we have that:

• The more efficient the algorithm, the slower () grows as grows
• The more inefficient the algorithm, the faster () grows as grows

Remember that we are always interested in the asymptotic bounds of the behavior of → ∞ as.

We have that () ∼ (()) ↔ | () | ≤ · (), for some positive () (()) which means that if is proportionate to it’s said to be superiorly limited by () some scalar multiplied by the function.

Linear systems and factorizations

Linear systems and factorizations → Matrices are the fundamental element in the computations intro matrices.

Matrix operations:

The sum of two matrices and is a matrix defined as × × × = + × × × where every element of is = + = 1,..., = 1,..., with and.

The product of two matrices and is a matrix defined as × × × = × × × × where every element of is = ∑ = 1,..., = 1,..., with and =1.

The transpose of a matrix, called, is obtained by turning every row into a column × and vice versa.

The determinant of a square matrix is recursively defined as × (1+) ( )() = ∑ ( − 1 ) · · > 2 if 1=1 () = − = 2 if 11 22 12 21 where is obtained by removing the -th row and the -th column from the matrix ×.

() ≠ 0 If a square matrix has it’s said to be non-singular. In this case, the inverse × −1 of is the matrix defined as −1 −1 = =.

Matrices allow a compact notation for representing systems of linear equations: + +... + = 11 1 12 2 1 1 + +... + = 21 1 22 2 2 2 … + +... + = 1 1 2 2 as = where:

• Is the matrix containing the coefficients ×
• Is the vector of unknowns ×1
• Is the vector of solutions ×1

→ = > If the system has more equations (rows) than unknowns (columns) the system is said to be overdetermined, while if the system has less equations that unknowns → < is said to be underdetermined.

In principle, if is non-singular, we could obtain the solution in the following way: −1 −1 −1 = → = → = −1 but this is NOT an efficient way to solve a system of equations since computing usually is very expensive and unstable.

Solving linear systems

Cramer’s Rule: each unknown variable in a system of linear algebraic equations may be expressed as a fraction of 2 determinants:

• () The denominator is the determinant of the matrix:
• The numerator is the determinant of the matrix where the column of the coefficients of the considered unknown variable is substituted with the known values: ( ℎ 1 ) = () 1

Note that Cramer’s Rule is inefficient as the number of equations grows.

What do we do then? Strategy: elimination of the unknowns

• Multiply the equations by constants
• Combine two equations so that one of the unknowns will be eliminated
• Solve the resulting single equation for the remaining unknown
• Substitute this value into either of the original equations to compute the other variable

E.g.: Apply the eliminate the unknown methods to the following equations: + = 11 1 12 2 1 + = 21 1 22 2 2.

Step 1: multiply eq.1 by and eq.2 by: 21 11 + = 21 11 1 21 12 2 21 1 + = 11 21 1 11 22 2 11 2.

Step 2: subtract eq.3 from eq.4 eliminating: 1 − = − 11 22 2 21 12 2 11 2 21 1.

Step 3: solve the eq. for: 2 − 11 2 21 1 = − 2 11 22 21 12.

Step 4: substitute the value into the original eq. to get the value of the other unknown: 1 − 22 1 12 2 = − 1 11 22 21 12.

This of course gets really hard with many equations and unknowns.

General rule: Naive Gaussian Elimination

The elimination of the unknowns technique can be extended to large sets of equations by → developing a systematic scheme or algorithm to eliminate unknowns and back-substitute Gaussian Elimination. We need some modifications to obtain a reliable algorithm, in particular we have to avoid divisions by zero. The “naive” version does not issue this problem. Nonetheless the approach is designed to solve a general set of equations: + +... + = 11 1 12 2 1 1 + +... + = 21 1 22 2 2 2 ... + +... + = 1 1 2 2.

Step 1: Forward elimination of the unknown: the goal is to get an upper triangular matrix. We need to eliminate the corresponding unknowns accordingly: + +... + = 11 1 12 2 1 1 (1) (1) (1) +... + = 22 2 2 2 ... (−1) (−1) = which gives as a result: (−1) = (−1) − 1.

Then we substitute into the equation to get, then substitute into the −1 −1 − 2 equation, and so on, proceeding backwards.

Substep 1: → Eliminate the first unknown from the second through the -th equations 1 21 21 + +... + = - 21 1 12 2 1 111 11 + +... + = 21 1 22 2 2 2 ---------------------------------------------- ( ) ( ) 21 21 21 − +... + − = − 22 12 2 2 1 211 11 11 where we can distinct: ( ) (1) 21 − = - 22 12 2211 ... - ( ) (1) 21 − = - 2 1 211 (1) 21 − = - 2 211.

We go on from the second to the third, from the third to the fourth, and we reach: + + +... + = 11 1 12 2 13 3 1 1 (1) (1) (1) (1) + +... + = 22 2 23 3 2 2 (1) (1) (1) (1) +