Howard Elman University of Maryland
• Vicki Howles • John Shadid • David Kay • David Silvester • Daniel Loghin • Ray Tuminaro • Alison Ramage • Andy Wathen 2
Incompressible Navier-Stokes Equations α − ν ∇ + ⋅ + = 2 ( grad ) grad u u u u p f t − = div 0 u α =0 ! steady state problem α =1 ! evolutionary problem Discretization and linearization Matrix equation x=b T u f F B = − 0 p B C Goal: Robust general solution algorithms Easy to implement Derived from subsidiary building blocks Adaptible to a variety of scenarios (steady / evolutionary / Stokes / Boussinesq ) 3
Outline 1. General approach preconditioning for saddle point problems 2. Relation to traditional approaches projection methods SIMPLE 3. Details for Navier - Stokes equations 4. Analytic / experimental results 5. Potential for more general problems 4
Preliminary: Steady Stokes Equations − ∇ + = 2 grad Rusten & Winther, 1992 u p f − = div 0 Silvester & Wathen , 1993 u Algebraic equation: A = discrete vector T u f A B = , 0 Laplacian − B C p Symmetric indefinite ! MINRES algorithm is applicable 0 Preconditioning operator: A 0 Q S 0 T Generalized eigenvalue problem: A B u A u = λ − 0 p Q p B C S Au+ B T p = λ Au , λ ≠ 1 ) u= 1/( λ - 1) A -1 B T p Case C=0: Bu = λ Q S p BA -1 B T p= λ ( λ - 1) Q S p
BA -1 B T p= µ Q S p, µ = λ ( λ - 1) Verf ü rth , 1984: For Q S = pressure mass matrix, µ 2 [a s ,b s ] independent of discretization parameter h µ λ =1 § (1+4 µ ) 1/2 , Under mapping λ ½ [ ] [ ] 0 d -a -b c / 2 k − 1 ( ) /( ) bc ad ≤ Convergence bound for MINRES : || || 2 || || r r 0 k + 1 ( ) /( ) bc ad − 1 0 Computational requirements, for times a vector A 0 Q S Poisson solve: can be approximated, e.g. with multigrid Mass matrix solve: cheap 6
Generalize to Navier-Stokes Equations Linearization via Picard iteration (slightly inaccurate notation) : − ν ∇ + ⋅ + = α + + + + ( 1 ) 2 ( 1 ) ( ) ( 1 ) ( 1 ) ( ) m m ( m grad ) m grad m m u u u u p f ∆ t − = + ( 1 ) div m 0 u T u f F B = Discretization 0 0 p B 0 Analogue of Stokes strategy: preconditioner F 0 Q S BF -1 B T p= µ Q S p , µ = λ ( λ - 1) Same analysis ! ) seek approximation Q S to Schur complement N.B. Same question arises for other strategies for linearization 7
Suppose Q S ¼ BF -1 B T so that eigenvalues of BF -1 B T p = µ Q S p are tightly clustered. Under mapping µ λ =1 § (1+4 µ ) 1/2 , eigenvalues λ are clustered in two regions, one on each side of imaginary axis Can improve this: Observation: symmetry is important for Stokes solver MINRES: optimal with fixed cost per step ) need block diagonal preconditioner For Navier - Stokes: don’t have symmetry need Krylov subspace method for nonsymmetric matrices (e.g. GMRES) 8
Alternative: block triangular preconditioner T F B − 0 Q S ! Generalized eigenvalue problem T T u F B u F B = µ − 0 0 B p Q p S Fu+ B T p = µ Fu , µ ≠ 1 ) u= -F -1 B T p Bu = - µ Q S p BF -1 B T p= µ Q S p µ ½ {1} [ Theorem (Fischer, Ramage, Silvester, Wathen ): For preconditioned GMRES iteration, let p 0 be arbitrary and u 0 = F -1 (f -B T u 0 ) ( ) r 0 =(0,w 0 )) (u k ,p k ) T be generated with block triangular preconditioner , (u k ,p k ) D be generated with block diagonal preconditioner . Then (u 2k ,p 2k ) D = (u 2k+1 ,p 2k+1 ) D = ( u k ,p k ) T . 9
Computational requirements, to implement block triangular preconditioner : − 1 T w F B v = Compute − 0 r q Q S Solve Q s r = - q, then solve Fw = v -B T r The only difference from block diagonal solve: matrix - vector product B T r (negligible) For second step: convection - diffusion solve: can be approximated, e.g. with multigrid For first step: something new needed: Q S 10
One more interpretation: 0 T T F B I F B = − − + − − 1 1 0 ( T ) B C BF I BF B C − 1 T T 0 I F B F B = ) − − + − − 1 1 0 ( T ) BF I B C BF B C − 1 ¼ T Q B F − 0 Q S Shows what is needed for stabilized discretization : Q S ¼ BF -1 B T + C 11
Relation to Projection Methods “Classical” O( ∆ t ) projection method ( Chorin 1967, Temam 1969): − (*) ( ) m u u Step 1: − ν ∇ + ⋅ = 2 * ( ) ( ) ( m grad ) m u u u f ( ) ∆ t − ν ∇ = − ⋅ + 1 1 2 * ( ) ( ) ( ) ( m grad ) m m I u f u u u ∆ ∆ In matrix form: ( t ) t +ν = − + 1 1 * ( ) ( ) m m M A u f Nu Mu ∆ ∆ t t ∇ + ( 1 ) m 1 1 * I u Step 2: u = ∆ ∆ t t − ∇ + ( 1 ) m 0 0 p + ( 1 ) In matrix form: T m 1 * Mu 1 M B u = ∆ ∆ t t + ( 1 ) m 0 0 B p Performed via pressure - Poisson solve 12
Substitute u * from Step 1 into Step 2 : ( )( ) − − + + + ν + ν 1 ( 1 ) m 1 1 1 u f Nu ( ) Mu ( ) T m m 1 M A M A M B = ∆ ∆ ∆ ∆ t t t t + ( 1 ) m 0 p 0 B ( ) + ν − ( ) 1 0 1 M A 1 T I M B ∆ ∆ t t Perot, 1993 − − 1 1 T B B M B 0 I ∆ t Contrast: update derived purely from linearization & discretization : − + M ν + + ( 1 ) T m 1 1 A B u f Nu ( ) Mu ( ) m m = ∆ ∆ t t + ( 1 ) 0 m 0 B p ( ) + ν − ( ) + ν 1 0 1 M A 1 T I M A B ∆ ∆ t t − − + ν 1 1 T B B M A B 0 I ( )( ) ∆ t − − + ν 1 = − ∆ ν = ∆ − 1 1 1 Error : T T ( ) T ( ) B M A M B t M AB O t ∆ ∆ t t 13
For higher order accuracy in time and related approaches: Dukowicz & Dvinsky 1992 Perot 1993 Quarteroni, Saleri & Veneziani 2000 Henriksen & Holman 2002 14
Relation to SIMPLE Patankar & Spaulding, 1972 − 0 1 T T F B F I F B = − − 1 T 0 0 B B BF B I ˆ ¼ − 0 Q 1 T I F B F − ˆ − 1 T B B F B 0 I Q F : approximate convection -diffusion solve ^ F: diagonal part of F Many variants 15
Perspective of New Approach • Take on Schur complement directly • Separate time discretization from algebraic algorithm • Enable flexible treatment of time discretization, linearization Allow choice of linearization Allow large time steps / CFL numbers if circumstances warrant 16
Kay, Loghin , Approximation for the Schur Complement (I) Wathen 2002 Suppose the gradient and convection - diffusion operators approximately commute (w=u (m) ): ∇ − ν ∇ + ⋅ ∇ ≈ − ν ∇ + ⋅ ∇ ∇ 2 2 ( ) ( ) w w p u Require pressure convection - diffusion operator − − − − = 1 1 1 1 T T Discrete analogue: M B M F M F M B u p p u u ⇒ − = − − 1 1 1 T T BF B BM B F M u p p A p ≡ − 1 Q A F M In practice: don’t have equality, instead S p p p − 1 Requirements: Poisson solve A − 1 p for Q − 1 S M Mass matrix solve p − 1 + Convection -diffusion solve F 17
Recommend
More recommend
