NEW COMPUTATIONAL TECHNOLOGIES for RELIABLE COMPUTER SIMULATIONS - PowerPoint PPT Presentation

✬ ✩ NEW COMPUTATIONAL TECHNOLOGIES for RELIABLE COMPUTER SIMULATIONS Sergey Korotov ⋆, † ⋆ Institute of Mathematics Helsinki University of Technology, Finland † Academy of Finland ✫ ✪ 1

✬ ✩ Main Problem in Mathematical Modelling and Numerical Analysis: RELIABLE VERIFICATION OF ACCURACY OF APPROXIMATE SOLUTIONS OBTAINED IN COMPUTER SIMULATIONS • What are the sources of various errors affecting the reliability of numerical solutions then ? ✫ ✪ 2

✬ ✩ U Physical process ⇓ − → Modelling error ε 1 ⇓ Differential model Au = f u ⇓ − → Approximation error ε 2 ⇓ u h Discrete model A h u h = f h ⇓ − → Computational error ε 3 ⇓ ✫ ✪ u ǫ A ǫ h u ǫ h = f ǫ Numerical solution h h 3

✬ ✩ Two Principal Relations • Computations on the base of a reliable (certified) model. Here ε 1 is assumed to be small and computed u ǫ h gives a desired information on U � U − u ǫ h � ≤ ε 1 + ε 2 + ε 3 (1) • Verification of a mathematical model. Here physical data U and u ǫ h are compared to judge on the quality of mathematical model � ε 1 � ≤ � U − u ǫ h � + ε 2 + ε 3 (2) ✫ ✪ 4

✬ ✩ Thus, two major problems of mathematical modelling: • reliable computer simulation • verification of mathematical models by comparing physical and mathematical experiments require efficient methods able to provide COMPUTABLE AND REALISTIC estimates of ε 2 + ε 3 ✫ ✪ 5

✬ ✩ No Error Control - No Reliability 2 1 1.5 0.9 1 0.8 0.5 0.7 0 0.6 −0.5 0.5 1 0.4 0.3 0.5 1 0.2 0.8 0.6 0.1 0.4 0.2 0 0 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 exact solution mesh with 741 nodes error = 0.09212 1 1 0.9 0.9 0.8 0.8 0.7 0.7 0.6 0.6 0.5 0.5 0.4 0.4 0.3 0.3 0.2 0.2 0.1 0.1 0 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 mesh with 785 nodes mesh with 855 nodes error = 0.11562 error = 0.44662 ✫ ✪ 6

✬ ✩ NEW TECHNOLOGIES OF ERROR CONTROL Elliptic Type BVPs Model elliptic BVP: Find a function u such that − div( A ∇ u ) = f in Ω & u = 0 on ∂ Ω (3) Ω is bounded domain in R d with Lipschitz boundary ∂ Ω, f ∈ L 2 (Ω), matrix of coefficients A is symmetric, with entries a ij ∈ L ∞ (Ω), i, j = 1 , . . . , d, and is such that c 2 | ξ | 2 ≥ A ( x ) ξ · ξ ≥ c 1 | ξ | 2 ∀ ξ ∈ R d ∀ x ∈ Ω (4) ✫ ✪ 7

✬ ✩ Weak formulation: Find u ∈ H 1 0 (Ω) such that � � ∀ w ∈ H 1 A ∇ u · ∇ w dx = fw dx 0 (Ω) (5) Ω Ω u be any function from H 1 • Let ¯ 0 (Ω) (e.g., computed by some numerical method) considered as approximation of u • We are interested in diverse and reliable control of deviation (or error) e := u − ¯ u during computational process, depending on the final goal(s) of calculations • So far two possibilities for the error control are well realized and used in mathematical and engineering communities ✫ ✪ 8

✬ ✩ Control of error in global norm: We are usually interested in estimation of gradient of deviation in energy norm 1 / 2   � | | |∇ ( u − ¯ u ) | | | := A ∇ ( u − ¯ u ) · ∇ ( u − ¯ u ) dx (6)   Ω • Such estimation gives a general presentation on the quality of ¯ u • Such type of control is in focus since [Babuˇ ska, Rheinboldt, 1978] Control of local errors: Another trend in a posteriori error estimation is based on concept of local error control in addition to classical control in energy norm • This approach is motivated by practical needs: analysts are often interested not only in the value of the overall error, but also in controlling errors over certain parts of solution domain, or relative ✫ ✪ to some interesting characteristics (“quantities of interest”) 9

✬ ✩ • Common way for performing such type control is to introduce linear functional ℓ associated with subdomain of interest (and/or with a “quantity of interest”) and to obtain a computable estimate for ℓ ( u − ¯ u ) • For example, one can be interested in estimation of � ℓ ( e ) = ℓ ( u − ¯ u ) = ϕ ( u − ¯ u ) dx (7) Ω where ϕ ∈ L 2 (Ω) and supp ϕ = ω ⊆ Ω, which provides us with info on error e locally, in subdomain ω • Estimates for ℓ ( u − ¯ u ) can also be used for estimation of unknown “quantity of interest” ℓ ( u ) since ℓ ( u ) = ℓ (¯ u ) + ℓ ( u − ¯ u ) • Certain “quantities of interest” (e.g., in linear elasticity) are often more interesting to practitioners than solution itself ✫ ✪ 10

✬ ✩ Error Control in Energy Norm u ∈ H 1 Upper estimate: Since u − ¯ 0 (Ω) we have � � | 2 = | | |∇ ( u − ¯ u ) | | f ( u − ¯ u ) dx − A ∇ ¯ u · ∇ ( u − ¯ u ) dx = Ω Ω � � � y ∗ ·∇ ( u − ¯ u − y ∗ ) ·∇ ( u − ¯ = f ( u − ¯ u ) dx − ( A ∇ ¯ u ) dx − u ) dx = Ω Ω Ω � � u − A − 1 y ∗ ) · ∇ ( u − ¯ ( f + div y ∗ )( u − ¯ = u ) dx − A ( ∇ ¯ u ) dx ≤ Ω Ω u − A − 1 y ∗ | ≤ � f + div y ∗ � � u − ¯ u � + | | |∇ ¯ | | | | |∇ ( u − ¯ u ) | | | where y ∗ ∈ H (div; Ω) and � · � is a standard norm in L 2 (Ω) • In above Green formula and CBS inequality are used ✫ ✪ 11

✬ ✩ From Friedrichs inequality ∀ w ∈ H 1 � w � ≤ c Ω �∇ w � 0 (Ω) (8) and conditions on A we observe that u ) � 2 ≤ c 2 u � 2 ≤ c 2 Ω | 2 � u − ¯ Ω �∇ ( u − ¯ | | |∇ ( u − ¯ u ) | | (9) c 1 Thus, we immediately get | ≤ c Ω u − A − 1 y ∗ | � f + div y ∗ � + | | | |∇ ( u − ¯ u ) | | | |∇ ¯ | | (10) √ c 1 After squaring and using Young inequality, we get � c 2 � 1 + 1 | 2 ≤ � f + div y ∗ � 2 + (1 + β ) | Ω u − A − 1 y ∗ | | 2 | | |∇ ( u − ¯ u ) | | | |∇ ¯ | β c 1 (11) valid for any y ∗ ∈ H (div; Ω) and any β > 0 ✫ ✪ 12

✬ ✩ • Upper bound (11) was first obtained in [Repin, 1997] in a general form using quite complicated tools of duality theory • More general problems with mixed Dirichlet/Neumann/Robin boundary conditions are treated in a just presented (and more simple) way in [Korotov, 2005], where (11) follows as a particular case Lower estimate: It is easy to show that | 2 = 2( J ( v ) − J ( u )) ∀ v ∈ H 1 | | |∇ ( u − v ) | | 0 (Ω) (12) where J ( w ) = 1 � � A ∇ w · ∇ w dx − fw dx (energy functional) 2 Ω Ω ∀ w ∈ H 1 • Function u minimizes J = ⇒ J ( u ) ≤ J ( w ) 0 (Ω) = ⇒ | 2 ≥ 2( J (¯ | | |∇ e | | u ) − J ( w )) (13) where w is any function from H 1 0 (Ω) ✫ ✪ 13

✬ ✩ Practical Realization • Assume that computations are performed on T h 1 , T h 2 , T h 3 , . . . , where h = h 1 > h 2 > h 3 > . . . , and, thus, we always have several successive approximations u h 1 , u h 2 , u h 3 , . . . On computation of constant c Ω : • Upper bound contains the only unknown constant c Ω 1 • c Ω = √ λ Ω , where λ Ω is smallest eigenvalue of Laplacian for Ω • Other estimation techniques involve many unknown constants: • which are very hard to estimate • which have to be always recomputed if computational mesh has changed ✫ • But c Ω remains the same under any changes of meshes ✪ 14

✬ ✩ On minimization of upper bound: now u ≡ u h ¯ � c 2 � 1 + 1 Ω � f +div y ∗ � 2 +(1 + β ) | |∇ u h − A − 1 y ∗ | | 2 M ⊕ ( u h , β, y ∗ ) := | | β c 1 • Coarse upper bounds can be computed fast using, e.g., y ∗ = G µ ( ∇ u µ ) ∈ H (div , Ω), µ = h 1 , h 2 , h 3 , . . . , and G µ is some gradient averaging operator • Sharp estimates require a real minimization of M ⊕ ( u h , β, y ∗ ) with respect to “free variables” y ∗ and β On computation of lower bound: We usually try to have J ( u h 1 ) > J ( u h 2 ) > J ( u h 3 ) > . . . , which suggests lower bounds | 2 = | | 2 ≥ 2( J ( u h ) − J ( u µ )) > 0 | | |∇ e | | | |∇ ( u − u h ) | | (14) ✫ where µ = h 2 , h 3 , . . . ✪ 15

✬ ✩ Test No. 1 1 5 4.5 4 2 1.5 3.5 1 0 3 0.5 0 2.5 −1 −1 −0.5 2 −0.5 1.7514 0 0 1.5 0.5 0.5 −1 1 1 1 −1 0 1 0 1 2 3 4 • A ≡ I, f ≡ 10 , T h with 92 nodes • Several successive meshes are used in computations of bounds • Upper bound decreases from 4.6426 to 2.4443 • Lower bound grows from 1.2191 to 1.7469 | 2 ≈ 1 . 7514 ✫ • | | |∇ e | | ✪ 16

✬ ✩ Error Control Via Linear Functionals • Second technology is developed for controlling local errors Consider e = u − ¯ u measured in terms of linear functional ℓ � Example: ℓ ( e ) = ϕ ( u − ¯ u ) dx (15) Ω where ϕ ∈ L 2 (Ω) and supp ϕ = ω ⊆ Ω Adjoint problem: Find v ∈ H 1 0 (Ω) such that � ∀ w ∈ H 1 A ∇ v · ∇ w dx = ℓ ( w ) 0 (Ω) (16) Ω • Adjoint problem is uniquely solvable ✫ ✪ 17

✬ ✩ • Usually we only have some approximation ¯ v (e.g., computed by FEM on adjoint mesh T τ ) • It can be shown that ℓ ( u − ¯ u ) = E 0 (¯ u, ¯ v ) + E 1 ( e, e ϕ ) (17) where � � E 0 (¯ u, ¯ v ) = f ¯ v dx − A ∇ ¯ v · ∇ ¯ (18) u dx Ω Ω and � E 1 ( e, e ϕ ) = A ∇ e · ∇ e ϕ dx (19) Ω • Deviation e ϕ := v − ¯ v is computational error in adjoint problem ✫ ✪ 18