Transcript Numerical Methods for Generalized Zakharov System

MA5233: Computational Mathematics
Weizhu Bao
Department of Mathematics
& Center for Computational Science and Engineering
National University of Singapore
Email: [email protected]
URL: http://www.math.nus.edu.sg/~bao
Computational Science
Third paradigm for
– Discovery in Science
– Solving scientific &
engineering problems
Interdisciplinary
–
–
–
–
Problem-driven
Mathematical models
Numerical methods
Algorithmic aspects—
computer science
– Programming
– Software
– Applications, ……
Dynamics of soliton in quantum physics
Wave interaction in plasma physics
Wave interaction in particle physics
Vortex-pair dynamics in superfluidity
Vortex-dipole dynamics in superfluidity
Vortex lattice dynamics in superfluidity
Vortex lattice dynamics in BEC
Computational Science
Computational Mathematics – Scientific computing/numerical
analysis
Computational Physics
Computational Chemistry
Computational Biology
Computational Fluid Dynamics
Computational Enginnering
Computational Materials Sciences
……...
Steps for solving a practical problems
Physical or engineering problems
Mathematical model – physical laws
Analytical methods – existence, regularity, solution, …
Numerical methods – discretization
Programming -- coding
Results -- computing
Compare with experimental results
Computational Mathematics
Numerical analysis/Scientific computing
A branch of mathematics interested in constructive
methods
Obtain numerically the solution of mathematical
problems
Theory or foundation of computational science
– Develop new numerical methods
– Computational error analysis:
• Stability
• Convergence
• Convergence rate or order of accuracy,….
History
Numerical analysis can be traced back to a
symposium with the title ``Problems for the
Numerical Analysis of the Future, UCLA, July 2931, 1948.
Volume 15 in Applied Mathematics Series,
National Bureau of Standards
Boom of research and applications: Fluid flow,
weather prediction, semiconductor, physics, ……
Milestone Algorithms
1901: Runge-Kutta methods for ODEs
1903: Cholesky decomposition
1926: Aitken acceleration process
( Sn1  Sn )2
Tn  Sn 
Sn2  2Sn1  Sn
1946: Monte Carlo method
1947: The simplex algorithm
1955: Romberg method
1956: The finite element method
if
Sn1  S
lim
 a 1
n S  S
n
Milestone algorithms
1957: The Fortran optimizing compiler
1959: QR algorithm
1960: Multigrid method
1965: Fast Fourier transform (FFT)
1969: Fast matrix manipulations
1976: High Performance computing & packages:
LAPACK, LINPACK – Matlab
1982: Wavelets
1982: Fast Multipole method
Top 10 Algorithms
1946: Monte Carlo method
1947: Simplex method for linear programming
1950: Krylov subspace iterative methods
1951: Decompositional approach for matrix computation
1957: Fortran optimizing compiler
1959-61: QR algorithms
1962: Quicksort
1965: Fast Fourier Transform (FFT)
1977: Integer relation detection algorithm
1982: Fast multipole algorithm
http://amath.colorado.edu/resources/archive/topten.pdf
Contents
Basic numerical methods
– Round-off error
– Function approximation and interpolation
– Numerical integration and differentiation
Numerical linear algebra
– Linear system solvers
– Eigenvalue probems
Numerical ODE
Nonlinear equations solvers & optimization
Contents
Numerical PDE
–
–
–
–
Finite difference method (FDM)
Finite element method (FEM)
Finite volume method (FVM)
Spectral method
Problem driven research:
– Computational Fluid dynamics (CFD)
– Computational physics
– Computational biology, ……
How to do it well
Three key factors
– Master all kinds of different numerical methods
– Know and aware the progress in the applied science
– Know and aware the progress in PDE or ODE
Ability for a graduate student
–
–
–
–
Solve problem correctly
Write your results neatly
Speak your results well and clear – presentation
Find good problems to solve
Numerical error
Example 1:
1 2  3
no error!!
Example 2:
1    1  3.14159 4.14159
round  off error!!!
Example 3:
cos(0.1)  1  0.12 / 2  0.995 cos(x)  1  x2 / 2 Truncationerror!!
Example 4:
cos(1 / 3)  1  0.3332 / 2  1  0.110/ 2  0.945
 cos(x )  1  x 2 / 2 T runcation round- off error!!
Numerical error
Truncation error or error of the method
cos(x)  P( x)  1  x2 / 2 Truncationerror|R(x)| O(|x|4 )
Round-off error: due to finite digits of numbers in computer
  3.14159  R    3.14159  0.0000026 ....
Numerical errors for practical problems
– Truncation error
– Round-off error
– Model error & observation error & empirical error etc.
Absolute error
Absolute error:
x  x *  e*  x *  x
  3.14  e*    3.14  0.00159....
1 / 3  0.33  e*  1 / 3  0.33  0.003333....
Absolute error bound (not unique!!):
| e* || x *  x |  *  x *  *  x  x *  *  x  x *  *
measure length with ruler with minimummm   *  0.5mm
  3.14   *  0.0016or 0.002
Relative error
An example: x  10  1   x*  1 y  1000 5   *y  5
which approximation is better?
Relative error:
e * x * x
e 

, in practical
x
x
*
r
Relative error bound:
*
 
*
r
1
 10%
| x* | 10
*
5
y  1000 5   r* 

 0.5%
| y* | 1000
x  10  1   r* 

e * x * x
e 

x*
x*
*
r
*
| x* |
Absolute error bounds for basic operations
x  x * with  x*,
Suppose
Error bounds
y  y * with  *y
x  y  x *  y*   x* y   x*   *y
*
x y  x * y*   xy
| x* |  *y  | y* |  x*
y / x  y * / x*  
*
y/x

| x* |  *y  | y* |  x*
| x* |2
f ( x )  f ( x*)   *f ( x ) | f ' ( x*) |  x*
*



| a | 1
a x  a x*   ax* | a |  x*   x*
  x | a | 1
x/a  x*/a 
*
x/a
  x* | a | 1
  / | a |  *
  x | a | 1
*
x
Significant digits
  3.1415926....
An example
  3.14  3 significant digits
  3.14159 6 significant digits
Definition: n significant digits
x  0.a1 a2 an   10m with a1  0
x  x*  0.a1 a2 an  10m such that |x  x* | 0.5  10n 1m
Method:
– Write in the standard form
– Count the number of digits after decimal
Error spreading: An example
1
 x n1

0  I n  e  x e dx  e e x  n  x n 1e x dx  1  n I n 1 , n  1,2,3
0
0
0


1
I0  e
1
1
1
n x
1
x
1
e
dx

1

e

0
I0  1  e1, I n  1  n I n1, n  1,2,3,
Algorithm 1:
– Use 4 significant digits for practical computation
~
~
~
I0  0.6321, I n  1  n I n1, n  1,2,3,
– Results n
0
1
2
3
~
In
n
~
In
4
0.6321 0.3679 0.2642 0.2074 0.1704
5
6
7
8
9
0.1480 0.1120 0.2160  0.7280 7.552
0*  0.5 104 ,  n*  n  n*1  n! 0*, 8*  8! 0*  2!!!
Error spreading: An example
Algorithm 2
– Result
n
Iˆn
n
Iˆn
e 1
1
 e 1  x 9 dx  I 9  e 1  x 9e x dx  e 1  x 9e dx 
10
10
0
0
0
1
1
1
1


1
1
e
1
  0.0684, Iˆn 1  (1  Iˆn ), n  9,8,0
Iˆ9   
2  10 10 
n
0
1
2
3
4
0.6321 0.3679 0.2642 0.2073 0.1708
5
6
7
8
9
0.1455 0.1268 0.1121 0.1035 0.0684
– Truncation error analysis

*
n1
  / n,    / 9!  !!!!
*
n
*
0
*
9
*
9
Convergence and its rate
Numerical integration

I   f ( x ) dx with f ( x )  x[2  sin(k x )] for someintegerk
0
Exact solution
sin( k )  cos(  k )
I  

2
k
k
2
Numerical methods
h
 0
N
,
xi  i h,
xi 1 / 2  (i  1 / 2)h, i  0,1,2, N
Composite midpoint rule
I  I Mh  h f ( x1/ 2 )  f ( x11/ 2 )  f ( x21/ 2 )   f ( xN 1/ 2 )
Composite Simpson’s rule
2
1
2
1
2
1
I  I Sh  h  f ( x0 )  f ( x1/2 )  f ( x1 )  f ( x11/2 )  f ( x2 )  f ( x21/2 ) 
3
3
3
3
3
6

2
1

f ( x N 1/2 )  f ( x N ) 
3
6

Composite trapezoidal rule
1
1

I  ITh  h  f ( x0 )  f ( x1 )  f ( x2 )    f ( xN 1 )  f ( xN )
2
2

Error estimate
| I  I Mh | O(h2 ), | I  ITh | O(h2 ),
| I  I Sh | O(h4 )
Numerical results
Numerical errors
Observations
Before h0
– Truncation error is too large !!
After h1
– Round-off error is dominated!!
Between h0 and h1
– Clear order of accuracy is observed for the method
We can observe clear convergence rate for
proper region of the mesh size!!!
Numerical Differentiation
Numerical differentiation
~
~
f ( x  h)  f ( x  h) f ( x  h)  f ( x  h)
f ' ( x) 

2h
2h
~
~
f ( x  h )  f ( x  h )  1 , f ( x  h )  f ( x  h )   2
The total error
h2
h3
f ( x  h)  f ( x)  h f ' ( x) 
f ' ' ( x) 
f ' ' ' ( x)  
2
6
~
~
f ( x  h)  f ( x  h)
f ( x  h )  1  f ( x  h )   2
| e* |
 f ' ( x) 
 f ' ( x)
2h
2h


 2  1
2h

h
h2
 2  1 h 2

f ' ' ' ( x)   

f ' ' ' ( x)  
6
2h
6
 A h 2 : E ( h )
Numerical Differentiation
Numerical Differentiation
Total error depends
O( h 2 )
– Truncation error:
16


10
in Matlab
– Round-off error:
– Minimizer of E(h):


E ' (h )  
 2 A h  0  h*  (
h2
E ' ' (h*)  2 A  0
2A
)1 / 3
– Double precision:
  10 , A  O(1)  h*  10
16
5
– Clear region to observe truncation error:
5
h*  10  h  h0  O(1)
How to determine order of accuracy
Numerical approximation or method
b
I   f ( x )dx  I h
| I  I h | O(h p )
a
How to determine p and C??
– By plot log E(h) vs log h
E (h) :| I  I h | C h p
How to determine order of accuracy
– By quotation
| I  Ih |
Ch p
p


2
| I  I h / 2 | C (h / 2) p