The purpose of this paper is to investigate the second-order partitioned time stepping method for a coupled system of heat equations with linear coupling condition.Our motivation is to consider the numerical simulations for the models of atmosphereocean interactions. Many numerical methods were developed for such problems, for example, operator-splitting and Lagrange multiplier domain decomposition methods were presented by Bresch and Koko[1]for two coupled Navier-Stokes fluids; Burman and Hansbo[2]gave an interior penalty stabilized method for an elliptic interface problem by treating the interface data as a Lagrange multiplier. However, solving the monolithic and coupled problem via global discretizations may preclude the usage of highly optimized black box subdomain solvers and limit the computational efficiency.Alternatively, the partitioned time stepping method provides a convenient decoupling strategy, the basic idea is based on the implicit-explicit (IMEX) approach, in which the action across the interface is lagged. It means that the subdomain solvers can be solved individually as black boxes. Connors and his co-workers developed the partitioned time stepping method for the atmosphere-ocean coupling[3-5]. Besides, these approaches have also been applied in decoupling the Stokes-Darcy model[6-13].
Figure 1 Two subdomains coupled by an interface I
For i=1,2, we introduce two Sobolev spaces
and the corresponding product space X = X1×X2and L2(Ω) = L2(Ω1)×L2(Ω2).Besides, let (·,·)Ωidenote the standard L2inner product on Ωi. For u,v ∈X with u=[u1,u2]Tand v=[v1,v2]T, ui,vi∈Xi, we define the L2and H1inner products in X as follows
and the induced L2and H1norms are ‖u‖=(u,u)and ‖u‖X=(u,u)
A natural subdomain variational formulation for (1)-(4), obtained by the general variation process, is to find (for i=1,2, i ̸=j) ui:[0,T]→Xisatisfying
For u ∈X,we define the operators A,B :X →X′via the Riesz representation theorem as follows
where [·] denotes the jump of the indicated quantity across the interface I. Thus, the coupled or monolithic variational formulation for (1)-(4) is obtained by summing (5)over i,j =1,2 and i ̸=j and is to find u:[0,T]→X satisfying
where f = [f1,f2]T. From[3], we know that the monolithic problem (8) has a global energy that is exactly conserved.
Let Tibe a triangulation of Ωiand Th= T1∪T2. We denote Xi,h⊂Xias the conforming finite element spaces with i = 1,2, and define Xh= X1,h×X2,h. The discrete operators Ah,Bh: Xh→X′h= Xhare defined analogously by restricting (6)and (7) to Xh. With these notations, the coupled finite element method for (8) can be written as: find u ∈Xhsatisfying
for any v ∈Xhwith the initial condition u(x,0)=u0.
In this section, we propose two partitioned time stepping methods for (1)-(4). In both schemes, the coupling terms on the interface conditions are treated explicitly so that only two decoupled diffusion equations are solved at each time step. Therefore,subproblems can be implemented in parallel and the legacy code for each one can be utilized. Here, we denote the time step size by △t.
The first scheme, we discretize in time via a second-order BDF, whereas the interface term is treated via a second-order explicit Gear’s extrapolation formula. The BDF2 scheme states as below.
For the second scheme,we combine the second-order implicit Adams-Moulton treatment of symmetric terms and the second-order explicit Adams-Bashforth treatment of the interface term to propose the following second-order scheme.
To prove the unconditional stabilities of two second-order schemes proposed in section 3, we give some basic facts and notation first. The G-matrix associated with the classical second-order BDF is given by
for any w ∈X2, define G-norm by |w|2G= 〈w,Gw〉. It is easy to verify that, for any vi∈X, i=0,1,2, we have
where w0= [v0,v1]Tand w1= [v1,v2]T. This G-norm is an equivalent norm on(L2(Ω))2in the sense that there exist Cl,Cu>0 such that
Besides, we also recall the following three basic inequalities:
Theorem 1(Unconditional stability of BDF2) Let T >0 be any fixed time,then Algorithm 1 is unconditionally by stable on (0,T].
Proof For Step I in Algorithm 1, we set v=u1in (10), it gives that
From Young’s and trace inequalities, we have
For Step II in Algorithm 1, by setting v=un+1in (11), we have
From (14), we have
where wn=[un+1,un]Tand δun+1=un+1−2un+un−1. Note that
Thus, by combining with (17), the unconditional stability of BDF2 is proved.
Next, to analyze the stability of AMB2 scheme, we introduce the following parameters
Substituting (30)-(32) into (29) yields
Define the energy
Then, by adding
to both sides, we have
In this section,we study the convergence results of both BDF2 and AMB2 schemes.We assume that the mesh is regular and the parameter h denotes the grid size. We use continuous piecewise polynomial of degree l for both finite element spaces X1,hand X2,h.
Definition 1 For any u ∈X, define a projection Phu ∈Xhsatisfying
It is easy to verify that if u ∈(Hl+1(Ω1))d×(Hl+1(Ω2))d,we have the following property
To analyze the error estimate, we define the error at t=tnas
Theorem 3(Convergence of BDF2) Assume that the exact solution of the couping problem(1)-(4)is sufficient regular in the sense of u ∈H3(0,T;H1)∩H2(0,T;Hl+1),and the time-step restriction
holds. Then, the solution of the BDF2 scheme satisfies the following error estimate
Proof By subtracting (11) from (9) at time tn, we derive the following error equation
From the definition of projection (38), (43) can be rewritten as
By setting vh=θn+1in (44), we have
Denote ϖn=[θn+1,θn]T, we discard the positive term Bh(θn+1,θn+1), it gives that
For the term Bh(δθn+1,θn+1),by using Cauchy-Schwarz inequality and trace inequality,we have
The terms on the RHS side of (47) can be bounded by using Young’s inequalities as
The desired error estimate follows from (58) and the interpolation error (39).
Theorem 4(Convergence of AMB2) Assume that the solution of the coupling problem (1)-(4) is sufficient regular in the sense of u ∈H3(0,T;H1)∩H1(0,T;Hl+1).Then the solution of AMB2 scheme satisfies the following error estimate
Proof By subtracting (13) from (9) at time tn+12, we derive the following error equation
It can be rewritten as
where we use the definition of projection
By setting vh=θn+1in (61), we derive
From Cauchy-Schwarz inequality, we have
and
For the interface term, there exists a constant C1, which is the same as that in (31)such that
By combining these inequalities with (63), we obtain
and discard the second positive term on the LHS of (66), we have
For the terms on the RHS side of (68), we have
The same as (57), we have
From Taylor’s theorem with the integral form of the remainder, we have
Similarly, we have
By combining (69)-(73) with (68) and discarding the positive terms on the LHS, we have
By recursion, we have
The desired error estimate follows from (75) and the interpolation error (39).
In this section, we carry out the numerical experiments for BDF2 and AMB2 schemes. We focus on the convergent rates of both schemes. Assume that Ω1=[0,1]×[0,1] and Ω2=[0,1]×[−1,0], the interface I is the portion of the x−axis from 0 to 1. Then ˆn1= [0,−1]Tand ˆn2= [0,1]T. The forcing term f is chosen to ensure that the exact solutions are as follows[3]
u1(t,x,y)=ax(1 −x)(1 −y)e−t, u2(t,x,y)=ax(1 −x)(c1+c2y+c3y2)e−t,
with
Computational results comparing the performance of two schemes are listed for two test problems:
Test problem 1: a=ν1=ν2=κ=1;
Test problem 2: a=4, ν1=5, ν2=10, κ=1/4.
For test problem 1, by setting △t=h with h=1/16, 1/32, 1/64 successively, we present the errors and convergent orders in Table 1 for both BDF2 and AMB2 with P1 finite element (here and later, we fix α = 0.8 for AMB2). The results illustrate the second-order in time accuracy for ‖un−unh‖. Besides, we notice that BDF2 has a significantly smaller error than AMB2. In Table 2, we set △t2= h3with h =1/8, 1/16, 1/32 and P2 finite element is chosen, the results illustrate the second-order in time accuracy and three-order in space accuracy for ‖un−unh‖. In this case, we can also find that BDF2 has a little better accuracy than AMB2. These results verify our theoretical results given in Theorem 3 and Theorem 4.
In the same way, in Table 3 and Table 4, we implement test problem 2 for both P1 and P2 finite element spaces,respectively. The expected convergence rates are obtained for BDF2 and AMB2 schemes.
Table 1 L2−error for BDF2 and AMB2 with P1, △t=h
Table 2 L2−errors for BDF2 and AMB2 with P2, △t2 =h3
Table 3 L2−errors for BDF2 and AMB2 with P1, △t=h
Table 4 L2−errors for BDF2 and AMB2 with P2, △t2 =h3
We proposed and investigated two second-order partitioned time stepping methods for a parabolic two domain problem. We have shown that our schemes are unconditionally stable and optimally convergent. The second-order partitioned methods for the fully nonlinear fluid-fluid problem is a subject of our future research.