10.2 Interface boundary conditions

Interface boundary conditions are needed when there are multiple grids in the computational domain, as discussed in Section 11 below. This applies to complex geometries, when multiple patches are used for computational efficiency, to adapt the domain decomposition to the wave zone, mesh refinement, or a combination of these.

A simple approach for exchanging the information between two or more grids is a combination of interpolation and extrapolation of the whole state vector at the intersection of the grids. This is the method of choice in mesh refinement, for example, and works well in practice. In the case of curvilinear grids and very–high-order FD or spectral methods, though, it is in general not only difficult to prove numerical stability, but even to find a scheme that exhibits stability from a practical point of view.

10.2.1 Penalty conditions

The penalty method discussed above for outer boundary conditions can also be used for multi-domain interface ones, including those present in complex geometries [115, 235, 118, 311, 312, 140]. It is simple to implement, robust, leads to stability for a very large class of problems and preserves the accuracy of arbitrary high-order FD and spectral methods.

Finite differences.
As an example, consider again an advection equation but now on the whole real line

ut = λux, − ∞ < x < ∞, t ≥ 0, u (0,x) = f(x), − ∞ < x < ∞,
where we assume that the initial data f is C ∞ smooth and has compact support. The domain of the real line is chosen just for simplicity, to focus on the interface procedure at x = 0. In the realistic case of compact domains, outer boundaries are also present and these can be treated by any of the methods discussed in Section 10.1.

Consider two grids, left and right, covering the intervals (− ∞, 0], and [0,+ ∞ ), respectively:

xl = iΔxl, i = 0,− 1,− 2,..., ir r x i = iΔx , i = 0, 1,2,...,
where the gridspacings need not agree, Δxl ⁄= Δxr, and the difference operators Dl and Dr, which are not necessarily equal to each other either, and do not even need to be of the same order of accuracy, satisfy SBP with respect to scalar products given by the weights l r σ ,σ on their individual grids:
∑0 +∑∞ ⟨vl,vl⟩ = Δxl σl vlvl, ⟨vr,vr⟩ = Δxr σr vrvr. (10.37 ) Σl ij ij Σr ij i j i,j= −∞ i,j=0

For the time, being assume that both SBP scalar products are diagonal, though. Then, the SAT semi-discrete approximation to the problem is

l d-vli = λDlvli +-δi,0S--(vr0 − vl0), i = 0,− 1,− 2,..., (10.38 ) dt Δxl σl00 d r r r δi,0Sr l r --vi = λD vi + ---r-r-(v0 − v0), i = 0,1,2,.... (10.39 ) dt Δx σ00
Notice that in this approach the numerical solution at any fixed resolution is bi-valued at the interface boundary x = 0, and in the same spirit as the penalty approach for outer boundary conditions; continuity of the fields at the interface is not enforced strongly but weakly.

Defining the energy

Ed := ⟨vl,vl⟩Σ + ⟨vr,vr⟩Σr , (10.40 ) l
and using the approximation (10.38View Equation, 10.39View Equation) and the SBP property of the difference operators, its time derivative is
d l l 2 r r 2 l r l r dtEd = (λ − 2S )(u0) + (− λ − 2S )(u 0) + 2(S + S )u 0u 0. (10.41 )
Then an estimate follows if two conditions are satisfied. One of them is λ + Sr − Sl = 0. The other one imposes an additional constraint on the values of S l and S r:

The following remarks are in order:

Spectral methods.
The standard procedure for interface spectral methods is to penalize each incoming characteristic variable, exactly as in the outer boundary condition case. Namely, as in Eq. (10.33View Equation) with lower bounds for the penalty strengths given by Eqs. (10.34View Equation) and (10.36View Equation) for Legendre and Chebyshev polynomials, respectively. We know from the FD analysis above, though, that in general this does not imply an energy estimate, and in order to achieve one, outgoing modes also need to be penalized, with strengths that are coupled to the penalty for incoming modes. However, the procedure of penalizing just incoming modes at interfaces appears to works well in practice, so we analyze this in some detail.

Figure 6View Image shows the spectrum of the Chebyshev penalty method as described, for an advection equation,

ut = ux, − 1 ≤ x ≤ 1, t ≥ 0, (10.47 )
where there is an interface boundary at x = 0. In more detail: the figure shows the maximum real component in the spectrum for N = 20 collocation points as a function of the penalty strength S. There are no real positive values, and this remains true for different values of N, in agreement with the fact that penalizing just the incoming mode works well in practice. The figure also supports the possibility that S ≳ 0.3 might actually be enough for stability, as in the outer boundary case.

The figure also shows the spectral radius as a function of the penalty strength. Beyond S = 1 it grows very quickly, and even though as mentioned the timestep is usually determined by keeping the time integration error below that one due to spatial discretization, that might not be the case if the spectral radius is too large. Thus, it is probably good to keep S ≲ 1.

View Image

Figure 6: Maximum real component (left) and spectral radii (right) versus penalty strength S, for the Chebyshev penalty method for the advection equation with two domains (N = 20 collocation points).

  Go to previous page Go up Go to next page