Solving Ill-Posed Cauchy Problems using Rational Krylov Methods - PowerPoint PPT Presentation

Solving Ill-Posed Cauchy Problems using Rational Krylov Methods Lars Eldén and Valeria Simoncini GAMM 2008 Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 1 / 60

Outline Cauchy Problem 1 Regularization Error Estimate Rational Krylov method 2 Krylov/Regularization Numerical Examples 3 Example 1 Animation Example 2 Conclusions 4 Future work Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 2 / 60

Motivating example: Ilmenite iron melting furnace Thermocouple Level K to level D Electrode thermocouples with level K highest Under Between Center electrodes electrodes The furnace material properties are temperature dependent. Problem: find the inner shape of the furnace. Nonlinear, and (rather) complex geometry PhD thesis: I M Skaar, Monitoring the Lining of a Melting Furnace, NTNU, Trondheim, 2001 Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 3 / 60

Inverse Heat Conduction Problem Steady state heat conduction problem: The upper boundary is unavailable for measurements Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 4 / 60

Ill-Posed Cauchy Problem Ω : connected in R 2 with smooth boundary ∂ Ω L : linear, self-adjoint, positive definite elliptic in Ω . u zz − Lu = 0 , ( x , y ) ∈ Ω , z ∈ [ 0 , z 1 ] , u ( x , y , z ) = 0 , ( x , y ) ∈ ∂ Ω , z ∈ [ 0 , z 1 ] , u ( x , y , 0 ) = g ( x , y ) , ( x , y ) ∈ Ω , u z ( x , y , 0 ) = 0 , ( x , y ) ∈ Ω . Sought: f ( x , y ) = u ( x , y , z 1 ) , ( x , y ) ∈ Ω . Formal solution √ u ( x , y , z ) = cosh ( z L ) g BUT: L is a pos. def. unbounded operator ⇒ ILL-POSED! Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 5 / 60

Standard Iterative Procedure Guess f ( 1 ) for k = 1 , 2 , . . . until convergence 1 Solve 1 u zz − Lu = 0 , ( x , y ) ∈ Ω , z ∈ [ 0 , z 1 ] , u ( x , y , z ) = 0 , ( x , y ) ∈ ∂ Ω , z ∈ [ 0 , z 1 ] , u ( x , y , z 1 ) = f ( k ) , ( x , y ) ∈ Ω , u z ( x , y , 0 ) = 0 , ( x , y ) ∈ Ω . giving u ( k ) Evaluate � g ( · , · ) − u ( k ) ( · , · , 0 ) � and adjust f ( k ) − → f ( k + 1 ) 2 In every iteration: Solve a 3D well-posed problem Often slow convergence. Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 6 / 60

Other possible methods? Tikhonov regularization? Impossible, because we do not know the integral operator for equations with variable coefficients and/or complicated geometry. Replace unbounded L by a bounded approximation? Possible in connection with finite difference approximation, but more difficult with finite elements. BUT: Krylov method! Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 7 / 60

Regularization Formal solution √ u ( x , y , z ) = cosh ( z L ) g BUT: L is a pos. def. unbounded operator ⇒ ILL-POSED! High frequency perturbations in g are blown up Regularization Replace unbounded operator by a bounded one! Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 8 / 60

Cut Off High Frequencies Eigenvalues of L : λ 2 j , j = 1 , 2 , . . . and λ j → + ∞ as j → ∞ General approach: Compute the k eigenvalues of smallest modulus: LX k = X k D k , where X k holds orthonormal eigenvectors Approximate by projection √ √ � L ) X k X ⊤ D k ) X ⊤ cosh ( z L ) g ≈ cosh ( z k g = X k cosh ( z k g Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 9 / 60

Eigenvalues of L L is large and sparse (of the order 10 4 − 10 5 , say) Compute the smallest eigenvalues ⇓ Operate with L − 1 ⇓ Solve many standard 2D elliptic problems − Lw = v Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 10 / 60

Error Estimate L 2 (Ω) setting, u is an “exact solution” v is an approximate solution with perturbed data g m Theorem Assume that � u ( · , · , 1 ) � ≤ M and that the data perturbation satisfies � g − g m � ≤ ǫ . Then if v is computed by projection using the eigenvalues satisfying λ j ≤ λ c , where λ c = ( 1 / z 1 ) log ( M /ǫ ) then � u ( · , · , z ) − v ( · , · , z ) � ≤ 3 ǫ 1 − z / z 1 M z / z 1 , 0 ≤ z ≤ z 1 . Optimal error bound Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 11 / 60

Eigenvalue method Is it necessary to compute the eigenvalues and eigenvectors accurately? Do we need all the information that we get in the eigenvalues? Can we take advantage of the fact that we want to compute an approximation of √ cosh ( z L ) g for this particular vector? Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 12 / 60

Eigenvalue method Is it necessary to compute the eigenvalues and eigenvectors accurately? NO! Do we need all the information that we get in the eigenvalues? NO! Can we take advantage of the fact that we want to compute an approximation of √ cosh ( z L ) g for this particular vector? YES! Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 13 / 60

Lanczos tridiagonalization Choose q 1 and iterate L − 1 q k = q k − 1 β k − 1 + q k α k + q k + 1 β k , k = 1 , 2 , . . . , with α k = q ⊤ k L − 1 q k and β k = q ⊤ k + 1 L − 1 q k ; One matrix-vector multiply L − 1 q k per step One standard 2D elliptic solve (black box) per step − Lw = q k Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 14 / 60

Lanczos properties Initial convergence influenced by starting vector v 1 . Choose q 1 = 1 /β g m , β = � g m � Faster for largest eigenvalues of L − 1 ⇒ Fast convergence for some of the smallest eigenvalues of L Optional: To get faster convergence for eigenvalues in [ 0 , λ c ] operate with ( L − τ I ) − 1 , where τ = λ c / 2 Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 15 / 60

Lanczos reduction L − 1 Q k = Q k T k + β k + 1 q k + 1 e ⊤ k + 1 ≈ Q k T k Approximation √ L ) g ≈ Q k cosh ( zT − 1 / 2 ) Q ⊤ cosh ( z k g k Problem: We cannot prevent T k from approximating large eigenvalues! Regularize T k : cut off large eigenvalues 1 Solution: 1 Krylov + regularization: O’Leary & Simmons (1981), Björck, Grimme & Van Dooren (1994) Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 16 / 60

Projected and Truncated Approximation Let (( θ ( k ) ) 2 , y ( k ) ) , j = 1 , . . . , k j j be the eigenpairs of T − 1 k Define F ( z , λ ) = cosh ( z λ 1 / 2 ) and S k = T − 1 k Truncated approximation: Q k F ( z , S c k ) Q ⊤ v k ( z ) = k g m y ( k ) cosh ( z θ ( k ) )( y ( k ) � ) ⊤ e 1 � g m � . := Q k j j j θ ( k ) ≤ λ c j Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 17 / 60

Error Estimate Recall F ( z , λ ) = cosh ( z λ 1 / 2 ) and S k = T − 1 k Theorem Let u be the “exact solution” and v k ( z ) = Q k F ( z , S c k ) Q ⊤ k g m Under the same hypotheses as earlier, 3 ǫ 1 − z / z 1 M z / z 1 � u ( z ) − v k ( z ) � ≤ 2 � [ F ( z , L c ) − Q k F ( z , S c k ) Q ⊤ + k ] g � . Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 18 / 60

Krylov/Regularization, First version Starting vector q 1 = 1 /β g m for k = 2 , 3 , . . . until “stable” [ Q k , T k ] = krylovstep ( L − 1 , Q k − 1 , T k − 1 ) Compute v k ( z ) = Q k F ( z , S c k ) Q ⊤ k g m end Check residual � Kv k − g m � < ǫ Kv k is the solution of the 3D problem with u = v k at the upper boundary and u z = 0 at the lower. Expensive! Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 19 / 60

Residual: � Kv k − g m � < ǫ Solve 3D problem (denote solution u k ) u zz − Lu = 0 , ( x , y ) ∈ Ω , z ∈ [ 0 , z 1 ] , u ( x , y , z ) = 0 , ( x , y ) ∈ ∂ Ω , z ∈ [ 0 , z 1 ] , u ( x , y , 1 ) = v k ( x , y ) , ( x , y ) ∈ Ω , u z ( x , y , 0 ) = 0 , ( x , y ) ∈ Ω . Well-posed but expensive! Kv k = u k ( x , y , 0 ) We only want to compute this when we are sure that � Kv k − g m � < ǫ Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 20 / 60

Krylov/Regularization: “stable” Starting vector q 1 = 1 /β g m for k = 1 , 2 , . . . until “stable” [ Q k , T k ] = krylovstep ( L − 1 , Q k − 1 , T k − 1 ) Compute v k ( z ) = Q k F ( z , S c k ) Q ⊤ k g m end Check residual � Kv k − g m � < ǫ How can we quantify “stable”? Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 21 / 60

Error Estimate for Krylov Procedure Recall F ( z , λ ) = cosh ( z λ 1 / 2 ) and S k = T − 1 k 3 ǫ 1 − z / z 1 M z / z 1 � u ( z ) − v k ( z ) � ≤ 2 � [ F ( z , L c ) − Q k F ( z , S c k ) Q ⊤ + k ] g � . Second term: Krylov approximation error for low frequency operator applied to the exact right hand side! Heuristic When the approximate solution does not change much between successive Krylov steps, then this error is small. Lars Eldén and Valeria Simoncini () Rational Krylov Method GAMM 2008 22 / 60