Structure preserving numerical schemes for complex dissipative/conservative nonlinear systems Jie Shen Purdue University Collaborators: Jie Xu, Jiang Yang, Qing Cheng ICERM Workshop on Numerical Methods and New Perspectives for Extended Liquid Crystalline Systems, Dec 9-13, 2019 Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Outline Motivation The scalar auxiliary variable (SAV) approach Application to a Q-tensor model SAV approach for Navier-Stokes equations and more general dissipative/conservative system The SAV approach with global constraints Concluding remarks Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Motivation Many physical problems can be modeled by PDEs that take the form of gradient flows or Hamiltonian systems. Examples include Allen-Cahn, Cahn-Hilliard, phase-field models, liquid crystals, super conductivity, image processing, ...; nonlinear Schr¨ odinger equation, Sine-Gordon equation, ... Gradient flows or Hamiltonian systems are dynamically driven by a free energy or Hamiltonian E ( φ ), and take the form: ∂φ ∂ t = −G δ E ( φ ) , δφ where G is a positive operator (gradient flows) or a skew-symmetric operator (Hamiltonian systems), and satisfy a dissipative or conservative energy law: dt E ( φ ) = − ( G δ E ( φ ) d , δ E ( φ ) ) . δφ δφ Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Some examples 2 |∇ φ | 2 and G = I ; 1 � heat equation: E ( φ ) = Ω 2 |∇ φ | 2 + 4 ǫ 2 ( φ 2 − 1) 2 ) and G = I ; Ω ( 1 1 � Allen-Cahn: E ( φ ) = 2 |∇ φ | 2 + 1 4 ǫ ( φ 2 − 1) 2 ) and Ω ( ǫ � Cahn-Hilliard: E ( φ ) = G = − ∆; 4 φ 4 + α 2 φ 2 − |∇ φ | 2 + 1 Ω ( 1 � 2 | ∆ φ | 2 ) Phase-field crystal: E ( φ ) = and G = − ∆; L 1 minimization: E ( φ ) = � Ω |∇ φ | and G = I ; Nonlinear Schr¨ odinger equation: 2 |∇ φ | 2 + 1 Ω ( 1 � 2 F ( | φ | 2 )) and G = i ; E ( φ ) = 2 | ∂ x φ | 2 + φ 3 ), G = ∂ x . Ω ( 1 � KDV equation: E ( φ ) = Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Existing approaches for time discretization Full implicit schemes: Can be unconditionally energy stable, but have to solve nonlinear equations at each step and may need ∆ t sufficiently small to have a unique solution. stabilized linearly implicit schemes and ETD schemes: With some conditions, first-order schemes can be unconditionally energy stable; efficient but difficult to get to second- or higher-order. Convex splitting schemes: First-order schemes always unconditionally energy stable; difficult to get to second- or higher-order; still need to solve (much easier) nonlinear equations at each step. (invariant) energy quadratization (IEQ/EQ) schemes: Linear, second-order, unconditionally energy stable; work with a large class of gradient flows; need to solve coupled systems with variable coefficients. Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
The scalar auxiliary variable (SAV) approach Ω [ 1 � Let E ( φ ) = 2 φ L φ + F ( φ )] dx and assume that � E 1 ( φ ) := Ω F ( φ ) dx is bounded from below, i.e., E 1 ( φ ) > − C 0 for some C 0 > 0. We introduce one scalar auxiliary variable (SAV): � 1 � 2 φ L φ dx + r 2 ( t ) − C 0 . r ( t ) = E 1 ( φ ) + C 0 so that E ( φ ) = Ω Then, the original system ∂φ ∂ t = −G δ E ( φ ) can be recast as: δφ ∂φ ∂ t = −G µ r ( t ) F ′ ( φ ) µ = L φ + � E 1 ( φ ) + C 0 1 � r t = F ′ ( φ ) φ t dx . � 2 E 1 ( φ ) + C 0 Ω Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
The SAV scheme based on Crank-Nicolson Second-order CN-AB scheme: φ n +1 − φ n = −G µ n +1 / 2 , ∆ t r n +1 / 2 µ n +1 / 2 = L φ n +1 / 2 + F ′ (˜ φ n +1 / 2 ) , � E 1 [˜ φ n +1 / 2 ] + C 0 r n +1 − r n φ n +1 − φ n F ′ (˜ φ n +1 / 2 ) � = dx , ∆ t � ∆ t E 1 [˜ Ω φ n +1 / 2 ] + C 0 2 φ n +1 / 2 ) := 3 g ( φ n ) − g ( φ n − 1 ) where ψ n +1 / 2 := ψ n +1 + ψ n and g (˜ . 2 2 Taking the inner products of the above with µ n +1 / 2 , − φ n +1 − φ n and ∆ t 2 r n +1 / 2 , respectively, one derives the following energy law: 1 2( φ n +1 , L φ n +1 )+ | r n +1 | 2 − 1 2( φ n , L φ n ) −| r n | 2 = − ( G µ n +1 / 2 , µ n +1 / 2 ) , i.e., it maintains the dissipation rate for dissipative systems and is energy conserving for conservative systems. Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Difference with a semi-implicit scheme Eliminating r n +1 from the above scheme, we find φ n +1 − φ n = −G µ n +1 / 2 , ∆ t r n µ n +1 / 2 = L φ n +1 / 2 + F ′ (˜ φ n +1 / 2 ) � E 1 [˜ φ n +1 / 2 ] + C 0 F ′ (˜ φ n +1 / 2 ) � φ n +1 / 2 )( φ n +1 − φ n ) dx F ′ (˜ + . 4( E 1 [˜ φ n +1 / 2 ] + C 0 ) Ω Without the last term, the scheme is a semi-implicit scheme. The last term serves as an extra stabilizing term! Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Efficient implementation We can write the schemes as a matrix system φ n +1 c 1 I G 0 = ¯ µ n +1 b n , −L c 2 I ∗ r n +1 ∗ 0 c 3 So we can solve r n +1 with a block Gaussian elimination, which requires solving a system with constant coefficients of the form � c 1 I � � φ � G = ¯ b . −L µ c 2 I With r n +1 known, we can obtain ( φ n +1 , µ n +1 ) by solving one more equation in the above form. So the cost is essentially twice the cost of a semi-implicit scheme, but it enjoys many additional benefits. Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Convergence and error analysis The SAV schemes are semi-implicit schemes. Previous stability and error analysis on semi-implicit schemes usually assume the uniform Lipschitz condition on the second derivative of the nonlinear free energy, which is not satisfied by even the double-well potential. Thanks to the unconditional energy stability of the SAV schemes, we can derive H 2 bounds for the numerical solution under mild conditions on the free energy. The H 2 bounds on the numerical solution will enable us to establish the convergence, and with additional smoothness assumption, the error estimates. Error estimates for the SAV approach have been established for (i) semi-discrete case (S. & J. Xu, SINUM ’18); (ii) fully discrete with finite-elements (H. Chen & S.); (iii) fully discrete with finite-differences (X. Li., & H. Rui & S., Math Comp). Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
SAV approach with stabilization If the nonlinear term is too ”strong”, the SAV approach may require restrictive time steps for accuracy. However, this situation can be easily improved with a stabilization. Given ǫ ≪ 1. Consider the free energy � 1 2 |∇ φ | 2 + 1 E ( φ ) = ǫ 2 F ( φ ) . 1 � Then, the SAV approach with E 1 ( φ ) = ǫ 2 F ( φ ) will require small time steps to get accurate results. Choose S > 0, and split the free energy as follows: � (1 2 |∇ φ | 2 + S ǫ 2 φ 2 ) + 1 ǫ 2 ( F ( φ ) − S φ 2 ) . E ( φ ) = 1 � ǫ 2 ( F ( φ ) − S φ 2 ) > − C 0 , so SAV can be We still have applied with this splitting, and leads to much improved results. Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
small time step large time step 1.5 1.5 SAVnoSTA SAVnoSTA SAVwithSTA SAVwithSTA 1 1 0.5 0.5 0 0 -0.5 -0.5 -1 -1 -1.5 -1.5 0 1 2 3 4 5 6 0 1 2 3 4 5 6 Figure: (Effect of stabilization) The solution at T = 0 . 1. Left: ∆ t = 10 − 4 ; Right: ∆ t = 4 × 10 − 3 . The red dashed lines represent solutions with stabilization, while the black solid lines represent solutions without stabilization. Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Adaptive time stepping Thanks to its unconditionally energy stability, one can (and should) couple the scheme with an adaptive time stepping strategy. Note that | r n +1 − 1 | serves as a simple and effective criterion. Figure: Numerical comparisons among small time steps, adaptive time steps, and large time steps Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
SAV approach for the Q-tensor model Consider the Landau-de Gennes free energy written as E [ Q ] = E b + E e with � � [ a 2tr Q 2 − b 3tr Q 3 + c 4(tr Q 2 ) 2 ] d x , E b ( Q ) = f b ( Q ) d x = Ω Ω 3 3 � [ L 1 2 |∇ Q | 2 + L 2 ∂ i Q ik ∂ j Q jk + L 3 � � E e ( Q ) = ∂ i Q jk ∂ j Q ik ] d x , 2 2 Ω k =1 k =1 where L 1 , L 1 + L 2 + L 3 > 0 so that E b ≥ 0, and c > 0 such that E b ( Q ) + C 0 > 0 for some C 0 > 0. Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
The SAV formulation Let L Q = δ E e δ Q , and introduce a SAV � r ( t ) = E b + C 0 . We can rewrite the system as: ∂ Q ∂ t = − µ, r ( t ) δ E b µ = L Q + √ E b + C 0 δ Q ; d r 1 ( δ E b δ Q , ∂ Q d t = 2 √ E b + C 0 ∂ t ) . Jie Shen Structure preserving numerical schemes for complex dissipative/conservative
Recommend
More recommend
