这是用户在 2024-3-28 11:00 为 https://app.immersivetranslate.com/pdf-pro/ed8ea26c-faf2-48de-bf67-d00d403798c4 保存的双语快照页面,由 沉浸式翻译 提供双语支持。了解如何保存?
2024_03_28_7eed7b69f809da89bcfeg

Implicit-explicit Runge-Kutta methods for time-dependent partial differential equations
时变偏微分方程的隐含-显式 Runge-Kutta 方法

Uri M. Ascher , Steven J. Ruuth , Raymond J. Spiteri a Institute of Applied Mathematics and Department of Computer Science, University of British Columbia, Vancouver, BC,
a 不列颠哥伦比亚大学应用数学研究所和计算机科学系,不列颠哥伦比亚省温哥华市、
Canada V6T  加拿大 V6T Department of Mathematics, University of California, Los Angeles, CA, USA
美国加州大学洛杉矶分校数学系
Institute of Applied Mathematics and Department of Mathematics, University of British Columbia, Vancouver, BC,
不列颠哥伦比亚大学应用数学研究所和数学系,不列颠哥伦比亚省温哥华市、
Canada V6T 1Z2 加拿大 V6T 1Z2

Abstract 摘要

Implicit-explicit (IMEX) linear multistep time-discretization schemes for partial differential equations have proved useful in many applications. However, they tend to have undesirable time-step restrictions when applied to convection-diffusion problems, unless diffusion strongly dominates and an appropriate BDF-based scheme is selected (Ascher et al., 1995).
在许多应用中,偏微分方程的隐式-显式(IMEX)线性多步时间离散化方案已被证明非常有用。然而,当它们应用于对流扩散问题时,往往会有不理想的时间步长限制,除非扩散占主导地位并选择了适当的基于 BDF 的方案(Ascher 等人,1995 年)。

In this paper, we develop Runge-Kutta-based IMEX schemes that have better stability regions than the best known IMEX multistep schemes over a wide parameter range. (c) 1997 Elsevier Science B.V.
在本文中,我们开发了基于 Runge-Kutta 的 IMEX 方案,与已知的最佳 IMEX 多步方案相比,该方案在较宽的参数范围内具有更好的稳定性区域。(c) 1997 Elsevier Science B.V. 版权所有。

1. Introduction 1.导言

When a time-dependent partial differential equation (PDE) involves terms of different types, it is a natural idea to employ different discretizations for them. Implicit-explicit (IMEX) time-discretization schemes are an example of such a strategy. Linear multistep IMEX schemes have been used by many researchers, especially in conjunction with spectral methods . Some schemes of this type were proposed and analyzed as far back as the late 1970s [5,15]. Instances of these methods have been successfully applied to the incompressible Navier-Stokes equations [9] and in environmental modeling studies [16]. A systematic, comparative study for PDEs of convection-diffusion type was carried out in [2], and a corresponding study for reaction-diffusion problems arising in morphology is reported in [11].
当一个与时间相关的偏微分方程(PDE)涉及不同类型的项时,对它们采用不同的离散方法是一个自然的想法。隐式-显式(IMEX)时间离散化方案就是这种策略的一个例子。线性多步 IMEX 方案已被许多研究人员采用,特别是与频谱方法结合使用 。早在 20 世纪 70 年代末,就有人提出并分析了这类方案[5,15]。这些方法已成功应用于不可压缩 Navier-Stokes 方程 [9] 和环境建模研究 [16]。文献[2]对对流-扩散类型的 PDE 进行了系统的比较研究,文献[11]对形态学中出现的反应-扩散问题进行了相应的研究。
In this work, we consider problems of convection-diffusion type (or hyperbolic-parabolic equations), e.g.,
在这项工作中,我们考虑的是对流-扩散型问题(或双曲-抛物方程),例如:
Discretization of the spatial derivatives (e.g., by finite-difference, finite-element, finite-volume, or spectral methods) yields a very large ODE system in time,
空间导数的离散化(如采用有限差分法、有限元法、有限体积法或频谱法)会在时间上产生一个非常大的 ODE 系统、
where corresponds to the convection (hyperbolic) term, , and corresponds to the diffusion (parabolic) term, . An IMEX scheme consists of applying an implicit discretization for and an explicit one for . This is natural because the system (1.2) with is generally stiff and linear, whereas with the system is not too stiff and it is often nonlinear. Moreover, if one insists on using an implicit discretization for the latter case, the construction of iterative solvers is made challenging by the properties of the matrix to be inverted.
其中 对应对流(双曲)项 对应扩散(抛物线)项 。IMEX 方案包括对 采用隐式离散化,对 采用显式离散化。这是很自然的,因为带有 的系统 (1.2) 通常是刚性和线性的,而带有 的系统则不太刚性,通常是非线性的。此外,如果坚持在后一种情况下使用隐式离散化,那么迭代求解器的构造就会因需要反演的矩阵的特性而面临挑战。
In [2] we have analyzed and experimented with linear multistep schemes for (1.2). For diffusiondominated problems, schemes utilizing Backward Differentiation Formulae (BDF) for have optimal damping properties. The corresponding IMEX scheme is a semi-explicit BDF (SBDF) [5,9,15]. However, when diffusion is not dominant, the choice of method is less obvious. Ideally, one would like a dissipative scheme for the hyperbolic term, so that the resulting IMEX scheme would have good stability and smoothing properties, independent of . Yet, it is well known that this necessitates looking for schemes of order (or more accurately, backward steps) at least 3 (cf. [2]). The main difficulty plaguing multistep IMEX schemes is the possibility that relatively small time steps become necessary due to stability restrictions (e.g., [14]). Such restrictions become worse when working with higher-order multistep schemes, as demonstrated in [2]. On the other hand, the stability region of explicit Runge-Kutta schemes actually increases slightly when three- and four-stage schemes are contemplated. Thus, we are led to investigate the development and performance of IMEX Runge-Kutta schemes.
在 [2] 中,我们分析并试验了 (1.2) 的线性多步方案。对于以扩散为主的问题,利用 的后向微分公式(BDF)的方案具有最佳阻尼特性。相应的 IMEX 方案是半显式 BDF(SBDF)[5,9,15]。然而,当扩散不占主导地位时,方法的选择就不那么明显了。理想情况下,我们希望双曲项采用耗散方案,这样得到的 IMEX 方案将具有良好的稳定性和平滑特性,与 无关。然而,众所周知,这就需要寻找阶数(或更准确地说,后退步数)至少为 3 的方案(参见 [2])。困扰多步 IMEX 方案的主要困难是,由于稳定性的限制,可能需要相对较小的时间步长(例如 [14])。正如文献[2]所示,在使用高阶多步方案时,这种限制会变得更加严重。另一方面,当考虑三阶和四阶方案时,显式 Runge-Kutta 方案的稳定区域实际上会略有增加。因此,我们开始研究 IMEX Runge-Kutta 方案的发展和性能。
When the system (1.2) arises from a PDE in more than one spatial variable, the simplicity and efficiency of solving the algebraic equations corresponding to the implicit part of the discretization at each time step is of paramount importance. In [2], we show that, under appropriate circumstances and for certain IMEX schemes, not much more than one multigrid W-cycle per time step is needed. Here, we aim at such a performance per Runge-Kutta stage. Hence, we concentrate on diagonally-implicit Runge-Kutta (DIRK) schemes for , paying particular attention to their attenuation properties. We construct a number of IMEX schemes of this sort, and investigate their properties. Since the effective order of DIRK schemes drops to 1 for very stiff ODEs [8], we do not expect the Runge-Kutta schemes to be competitive with SBDF in the limit where is heavily dominant and very stiff. However, our Runge-Kutta schemes are shown to have excellent properties for other parameter ranges.
当系统 (1.2) 由一个以上空间变量的 PDE 引起时,求解每个时间步离散化的隐含部分对应的代数方程的简便性和效率至关重要。在 [2] 中,我们证明了在适当的情况下,对于某些 IMEX 方案,每个时间步所需的多网格 W 循环不超过一个。在这里,我们的目标是实现每个 Runge-Kutta 阶段的性能。因此,我们专注于 的对角隐式 Runge-Kutta (DIRK) 方案,尤其关注其衰减特性。我们构建了许多此类 IMEX 方案,并研究了它们的特性。由于对于非常僵硬的 ODEs,DIRK 方案的有效阶数降到了 1 [8],因此我们预计在 占主导地位且非常僵硬的极限情况下,Runge-Kutta 方案与 SBDF 相比并不具有竞争力。不过,我们的 Runge-Kutta 方案在其他参数范围内具有出色的特性。
It is well known that, despite its recent rise in popularity, the use of Runge-Kutta schemes for the integration of time-dependent PDEs is not without drawbacks. One problem is a loss of accuracy order when time-dependent inflow boundary conditions are prescribed [4]. A remedy is proposed in [1]. Note also that, in the (usual) case where the accuracy of the solution at the end of the step is higher than at internal stages, the accuracy to which these internal stages must be calculated (say by an iterative method, as in Section 5) cannot be lowered.
众所周知,尽管 Runge-Kutta 方案近来大受欢迎,但使用 Runge-Kutta 方案对依赖时间的 PDEs 进行积分并非没有缺点。其中一个问题是,当规定了随时间变化的流入边界条件时,精确度会下降 [4]。[1]中提出了一种补救方法。还需注意的是,在(通常)步末解的精度高于内部阶段解的情况下,这些内部阶段的计算精度(如第 5 节中的迭代法)不能降低。
The system (1.2) can be cast as a partitioned system,
系统 (1.2) 可以看作一个分区系统、
with and . The IMEX Runge-Kutta schemes devised below are then partitioned Runge-Kutta schemes. Their orders are discussed in [7, Section II.15]. They can also be seen as a particular class of splitting methods. We note that for (1.2), and need not be known individually (indeed we may not have initial conditions for them). All we need to have is their smooth, local existence.
。下面设计的 IMEX Runge-Kutta 方案就是分区 Runge-Kutta 方案。它们的阶次在 [7, 第 II.15 节] 中讨论。它们也可以看作是一类特殊的分割方法。我们注意到,对于 (1.2), 不需要单独知道(事实上,我们可能没有它们的初始条件)。我们只需要知道它们平稳的局部存在。
In Section 2, we develop a number of IMEX Runge-Kutta schemes of up to four stages, which are up to third-order accurate. The schemes divide naturally into two classes-those whose last internal stage is identified with the solution at the next time instance (i.e., at the end of the time step), and those where an additional quadrature is used at the end of the step. The first class is particularly good for highly-stiff problems, but it is the second class which seems to yield some of the more promising variants for a wide range of parameters, particularly the scheme identified as ( 3 internal stages for the implicit formula, 4 stages for the explicit and a combined accuracy of order 3 ).
在第 2 节中,我们开发了多达四级的 IMEX Runge-Kutta 方案,其精度可达三阶。这些方案自然分为两类--一类是最后一个内部级数与下一个时间实例(即时间步结束时)的解相一致,另一类是在步结束时使用额外的正交。第一类方案尤其适用于高刚度问题,但第二类方案似乎能产生一些更有前途的变量,适用于各种参数,尤其是确定为 的方案(隐式为 3 个内部阶段,显式为 4 个阶段,综合精度为 3 阶)。
In Section 3, we investigate the performance of the various schemes for a simple test equation, which arises from a von Neumann analysis of an advection-diffusion equation. Stability regions are plotted in Figs. 1 and 2. Results of test runs for a linear advection-diffusion equation are reported in Section 4.1 and the Burgers equation is experimented with in Section 4.2. In Section 5, we investigate the performance of our various schemes in conjunction with employing a multigrid method for resolving the implicitness at each stage.
在第 3 节中,我们研究了各种方案在一个简单测试方程中的性能,该方程来自对一个平流扩散方程的 von Neumann 分析。稳定区域如图 1 和图 2 所示。第 4.1 节报告了线性平流-扩散方程的测试运行结果,第 4.2 节对布尔格斯方程进行了试验。在第 5 节中,我们将结合使用多网格方法解决每个阶段的隐含性问题,研究各种方案的性能。

2. IMEX Runge-Kutta schemes
第二届 IMEX Runge-Kutta 计划

We now develop some IMEX Runge-Kutta schemes. For , we consider an implicit -stage DIRK scheme [8] with coefficients , in the usual Butcher notation. Let . For , we consider an -stage explicit scheme with the abscissae and coefficients . To cast the DIRK scheme as an -stage scheme as well, we can pad the -stage scheme with zeroes, obtaining the tableau
现在我们开发一些 IMEX Runge-Kutta 方案。对于 ,我们考虑系数为 的隐式 - 阶段 DIRK 方案 [8],采用通常的 Butcher 符号。设 。对于 ,我们考虑一个 级的显式方案,其主轴为 ,系数为 。为了把 "迪尔克 "方案也转换成一个 级方案,我们可以在 级方案中加入零,得到表格
Referring to the coefficients of this padded tableau as , we see that . This property simplifies significantly the form of the order conditions on the method's coefficients due to coupling between the individual schemes.
把这个填充表的系数称为 ,我们可以看到 。由于各个方案之间的耦合,这一特性大大简化了方法系数的阶次条件形式。
Set
One step from to of the IMEX scheme is given as follows:
IMEX 方案从 的一个步骤如下:
For do: 对于 这样做:
  • Solve for : 求出
where 其中
  • Evaluate 评估
Finally, evaluate 最后,评估
We consider two special cases, which lead to two sub-families of methods: those which satisfy (2.2) and those which satisfy (2.3)-(2.4):
我们考虑了两种特殊情况,由此产生了两个方法子系列:满足 (2.2) 的方法和满足 (2.3)-(2.4) 的方法:
(1) In the case that (and in particular ), we have in place of (2.1e)
(1) 在 (尤其是 )的情况下,我们可以用(2.1e)来代替
(2) In the case that need not be evaluated. Furthermore, if
(2) 在 无需求值的情况下。此外,如果
(implying that the implicit scheme is stiffly accurate and ), then substituting (2.1e) into the expression for , we see that
(意味着隐式方案是刚性精确的,且 ),然后将 (2.1e) 代入 的表达式,我们可以看到
This is useful for very stiff problems. We note also that the explicit scheme can now be generated from a general explicit -stage scheme based on the abscissae ,
这对于非常僵硬的问题非常有用。我们还注意到,显式方案现在可以从基于缺省值 的一般显式 级方案生成、
In (2.1b), there are nonlinear systems to solve involving alone, and the operators to be inverted are
在 (2.1b) 中,有 个仅涉及 的非线性系统需要求解,需要反演的算子是
Assuming that the Jacobian is evaluated only once (at ), this operator is the same for each if is independent of , i.e., for SDIRK schemes. However, if an iterative method like 1-cycle multigrid is used for the implicit scheme, then perhaps it is less critical for the overall efficiency to insist on SDIRK, and a more general DIRK scheme can be contemplated. See [8, Section IV.6], for stifflyaccurate SDIRK schemes. In our search for particularly accurate and stable schemes, we have allowed the diagonal elements of to differ from one another. Despite this, the schemes we recommend are all SDIRK, since there were sufficient excess degrees of freedom in the solution process to allow for this computationally-efficient choice.
假设 Jacobian 只求值一次(在 处),那么如果 无关,即对于 SDIRK 方案,该算子对每个 都是相同的。然而,如果隐式方案使用了像 1 周期多网格这样的迭代法,那么坚持使用 SDIRK 也许对整体效率并不那么重要,可以考虑使用更通用的 DIRK 方案。关于刚性精确的 SDIRK 方案,请参见 [8,第 IV.6 节]。在寻找特别精确和稳定的方案时,我们允许 的对角元素彼此不同。尽管如此,我们推荐的方案都是 SDIRK,因为在求解过程中有足够多的自由度,允许这种计算效率高的选择。
Example 1. Consider the case of (1.2) where
例 1.考虑 (1.2) 的情况,其中
Thus, we apply the implicit scheme to advance and the explicit one to advance . Observe that if is influenced by a control function, for instance, then at each stage , this influence is propagated to from . This does not happen when using a "normal" explicit scheme (cf. [13]). Moreover, if is independent of then
因此,我们采用隐式方案推进 ,采用显式方案推进 。请注意,如果 受到控制函数的影响,那么在每个阶段 ,这种影响会从 传播到 。在使用 "正常 "显式方案时,这种情况不会发生(参见 [13])。此外,如果 无关,那么
and the scheme becomes explicit!
,该方案就变得清晰明了!
In particular, for Hamiltonian systems we can set . If is independent of , then we obtain an explicit scheme which still retains some advantages of the implicit one.
特别是,对于哈密顿系统,我们可以设置 。 如果 无关,那么我们得到的显式方案仍然保留了隐式方案的一些优点。
We now proceed to construct some instances of our IMEX RK family of schemes. We will use the triplet to identify a scheme (2.1), where is the number of stages of the implicit scheme, is the number of explicit scheme stages (so for case (1), i.e., (2.2), and for case (2), i.e., (2.3)-(2.4) above), and is the combined order of the scheme.
现在,我们开始构建 IMEX RK 系列方案的一些实例。我们将使用三元组 来标识方案 (2.1),其中 是隐式方案的级数, 是显式方案的级数(因此 表示情况 (1),即 (2.2), 表示情况 (2),即上文 (2.3)-(2.4)), 是方案的组合阶数。

2.1. Forward-backward Euler
2.1.前后欧拉

The pair of backward and forward Euler schemes
一对后向和前向欧拉方案
can be padded to read
可以填充为
This yields the linear one-step IMEX (cf. [2]):
这就产生了线性一步 IMEX(参见 [2]):
Note that this scheme satisfies (2.3).
请注意,该方案符合 (2.3)。

2.2. Forward-backward Euler
2.2.前后欧拉

The above form of Euler's method does not satisfy . The variant that does is
上述形式的欧拉方法不满足 。满足要求的变体是
i.e., , where . This scheme involves an additional evaluation of per step, compared to the previous one.
,其中 。与前一种方案相比,这种方案每一步需要额外评估

2.3. Implicit-explicit midpoint
2.3.隐式-显式中点

The Euler pairs are first-order accurate and the scheme has well-known disadvantages when (cf. Section 3 and [2]). The following pair of implicit-explicit schemes:
欧拉对是一阶精确的,而 方案在 时有众所周知的缺点(参见第 3 节和 [2])。下面是一对隐式-显式方案:
correspond to applying explicit midpoint for and implicit midpoint for . It is second-order accurate because the two schemes from which it is composed are each second-order accurate and [7]. This scheme performs usually comparably to the unreasonably popular CNAB considered in [2], and has better symmetry properties.
它是二阶精确的,因为由它组成的两个方案分别是二阶精确的和 [7]。这种方案的性能通常可与 [2] 中考虑的不合理流行的 CNAB 方案相媲美,并且具有更好的对称性。
Example 2. Applying the scheme to the separable Hamiltonian system
例 2.将 方案应用于可分离的哈密顿系统
and defining 并界定
we obtain 我们得到
We note that this scheme is explicit, and it can also be shown to be symplectic. It is identical to the leap-frog/Verlet scheme [12]
我们注意到这个方案是显式的,也可以证明它是交映的。它与跃迁/韦莱方案[12]相同

2.4. A third-order combination
2.4.三阶组合

The two-stage, third-order DIRK scheme with the best damping properties turns out to be the DIRK scheme [7, p. 207]
具有最佳阻尼特性的两级三阶 DIRK 方案是 DIRK 方案 [7, 第 207 页]。
with . For the test equation , if we let , then as . This is a marked improvement over the midpoint scheme which has no attenuation at the stiffness limit .
。对于测试方程 ,如果我们让 ,则为 。这比中点方案 有了明显的改进,因为中点方案在刚度极限 没有衰减。
The corresponding third-order explicit Runge-Kutta scheme (ERK) is
相应的三阶显式 Runge-Kutta 方案(ERK)为
The resulting IMEX combination is third-order accurate, has some dissipation near the imaginary axis (like all three-stage third-order ERK methods) and has some attenuation of the stability function at .
由此产生的 IMEX 组合具有三阶精度,在虚轴附近有一些耗散(与所有三级三阶 ERK 方法一样),在 处的稳定函数有一些衰减。

2.5. L-stable, two-stage, second-order DIRK
2.5.L 稳定、两阶段、二阶 DIRK

One may be concerned that the attenuation that the above scheme has may not be sufficient for some problems. A two-stage, second-order DIRK scheme which is stiffly accurate is [8, p. 106]
有人可能会担心,上述方案的衰减可能不足以解决某些问题。一个两级二阶 DIRK 方案是刚性精确的 [8,第 106 页]。
where . The corresponding three-stage second-order ERK is
其中 。相应的三级二阶 ERK 为
Varying to get a dissipative stability region, we match the terms of the exponential in the stability function up to third order. This yields the stability region of a three-stage, third-order explicit RK scheme, with . The resulting IMEX combination is second-order accurate.
改变 以获得耗散稳定区域,我们匹配稳定函数中指数项,直到三阶。这样就得到了一个三阶、三阶显式 RK 方案的稳定区域, 。由此得到的 IMEX 组合具有二阶精度。

2.6. L-stable, two-stage, second-order DIRK
2.6.L 稳定、两阶段、二阶 DIRK

We use the same form (with the same but with not yet specified) as in the previous IMEX
我们使用与之前的 IMEX 相同的形式( 相同,但 尚未指定),即
scheme, except for requiring that (2.3) hold instead of (2.2), i.e., . This gives a second-order scheme
除了要求 (2.3) 成立而不是 (2.2),即 。这就给出了一个二阶方案
with .  .

2.7. L-stable, three-stage, third-order DIRK
2.7.L 稳定、三级、三阶 DIRK

We would now like to exploit the larger dissipativity region, which four-stage, fourth-order ERK schemes have near the imaginary axis. A three-stage, third-order DIRK scheme which is stiffly accurate is , p. 106]
现在,我们想利用四级四阶 ERK 方案在虚轴附近较大的耗散区域。刚性精确的三阶 DIRK 方案是 , 第 106 页]。
0 0
0
1
,
where is the middle root of ,
其中 的中根、
Numerically, this evaluates to
从数值上看,其值为
0.4358665215 0.4358665215 0 0
0.7179332608 0.2820667392 0.4358665215 0
1 1.208496649 -0.644363171 0.4358665215
1.208496649 -0.644363171 0.4358665215
.
The corresponding four-stage, third-order ERK, constructed such that it has the same (desirable) stability region as all four-stage, fourth-order ERK schemes, is
相应的四级三阶 ERK 与所有四级四阶 ERK 方案具有相同的(理想的)稳定区域,即
0 0 0 0 0
0 0 0
0 0
1 0
0
.
The third-order conditions indicate that this is a two-parameter family of schemes. Taking our degrees of freedom to be , the remaining expressions are
三阶条件表明,这是一个双参数方案族。假设我们的自由度为 ,其余表达式为
Selecting particular values for and so as to match the terms of the exponential in the stability function up to fourth order, we obtain the scheme
选择 的特定值,以匹配稳定函数中指数项的四阶以下值,我们得到了以下方案
0 0 0 0 0
0.4358665215 0.4358665215 0 0 0
0.7179332608 0.3212788860 0.3966543747 0 0
1 -0.105858296 0.5529291479 0.5529291479 0
0 1.208496649 -0.644363171 0.4358665215
The resulting IMEX combination is third-order accurate.
由此产生的 IMEX 组合具有三阶精度。

2.8. A four-stage, third-order combination
2.8.四级三阶组合

No pair consisting of a three-stage, L-stable DIRK and a four-stage ERK satisfying (2.3) with a combined third-order accuracy could be found. In order to obtain a third-order formula satisfying (2.3), we need to go to a four-stage, third-order, L-stable DIRK, coupled with a four-stage ERK in the manner indicated in case (2) of Section 2. After having satisfied the order conditions and imposed L-stability, we are left with a few families of schemes.
我们找不到一个由三阶 L 稳定 DIRK 和四阶 ERK 组成的组合,能满足(2.3)的三阶精度。为了得到满足 (2.3) 的三阶公式,我们需要按照第 2 节例 (2) 所示的方式,将四阶、三阶、L 稳定 DIRK 与四阶 ERK 结合起来。在满足了阶次条件并施加了 L 稳定性之后,我们只剩下几组方案了。
Experimenting with the test equation of Section 3, we have selected a scheme which has rational coefficients which are not too large, and a diagonal entry which is close to that of the four-stage, fourth-order, L-stable DIRK. In particular, we have chosen , accompanied by the choices , and . The resulting scheme is
通过对第 3 节测试方程的实验,我们选择了一个有理系数不太大的方案,其对角线条目与四阶 L 稳定 DIRK 的对角线条目相近。特别是,我们选择了 ,并同时选择了 。由此得出的方案是
{
}
0 0 0 0 0 0 0 0 0
0 0 0 0
2 0 0
2
2
2
11
1 0
0 0 0 U
3 3 1 1 0 0
1 0
0

3. Test equation for convection-diffusion
3.对流-扩散测试方程

As is shown in [2], using a centered discretization scheme for the spatial derivatives of advectiondiffusion PDEs, a von Neumann analysis yields that a simple yet appropriate test equation is the scalar ODE (1.2) with
如文献[2]所示,通过对平流扩散 PDE 的空间导数采用居中离散化方案,冯-诺依曼分析得出一个简单而又合适的测试方程是标量 ODE (1.2),其计算公式为
where are real constants, usually , and
其中 是实常数,通常为 ,而
For a given time step , we can define
对于给定的时间步长 ,我们可以定义
and write an IMEX step (2.1)-(2.2) as
并将 IMEX 步骤 (2.1)-(2.2) 写为
where 其中
For schemes where , i.e., where (2.3) holds, an expression like (3.1) need not be used because .
对于 ,即 (2.3) 成立的方案,无需使用类似 (3.1) 的表达式,因为 .
Let us consider the two Euler schemes first. For the usual pair (1,1,1), Eq. (3.2) implies
让我们先考虑两个欧拉方案。对于通常的一对 (1,1,1),公式 (3.2) 意味着
Fig. 1. Time step stability restrictions for various values of and , for IMEX schemes satisfying (2.3). Note that all these schemes perform well in the high-attenuation limit.
图 1.满足 (2.3) 的 IMEX 方案在不同 值下的时间步长稳定性限制。请注意,所有这些方案在高衰减极限时都表现良好。
On the imaginary axis, , and the scheme is unconditionally unstable for any (cf. [2]). As and with bounded, however, the scheme is unconditionally stable and , i.e., attenuation compatible with that of the exact solution is realized. For the other pair, , we get from (3.1)-(3.2)
在虚轴上, ,对于任何 ,方案都是无条件不稳定的(参见 [2])。然而,当 有界时,方案是无条件稳定的,并且 ,即实现了与精确解兼容的衰减。对于另一对,即 ,我们可以从 (3.1)-(3.2) 得到
Setting initially, we obtain that, along the imaginary axis,
设为初始值,我们可以得到沿虚轴方向的结果、
so stability is achieved provided that
因此,只要
This corresponds to the CFL condition. On the other hand, letting , we also obtain the restriction in the high-attenuation limit. Thus, we see that the variant has stability characteristics which are preferable over the usual backward-forward Euler scheme in the hyperbolic limit and worse in the parabolic limit.
这相当于 CFL 条件。另一方面,让 ,我们还可以得到高衰减极限下的限制条件 。由此可见,变式 在双曲极限时的稳定性优于通常的后向欧拉方案,而在抛物线极限时则较差。
Now, we can ask, given a fixed ratio , what is the maximum for which the scheme is stable, i.e., ? This gives the stability restriction on the step size for a given problem. In Fig. 1, we plot the resulting curves for the IMEX scheme presented earlier of the second type, i.e., those satisfying (2.3). Corresponding curves for the schemes satisfying (2.2) are plotted in Fig. 2.
现在,我们可以问,在给定比率 的情况下,方案稳定的最大值 是多少,即 ?这就给出了给定问题的步长稳定性限制。在图 1 中,我们绘制了前面介绍的第二类 IMEX 方案(即满足 (2.3) 的方案)的结果曲线。满足(2.2)的方案的相应曲线见图 2。
Common to all the schemes in Fig. 2 is that, as , there is a time step restriction involving , typically . In contrast, no such restriction occurs for schemes satisfying (2.3), and hence (2.4), in Fig. 1. The behavior of these latter schemes is qualitatively similar in the high-attenuation limit to that of the SBDF schemes [2]. On the other hand, the scheme in Fig. 2 suggests a potential advantage, in terms of allowable time steps, over a large range of parameter values.
图 2 中所有方案的共同点是,当 时,有一个涉及 的时间步长限制,通常是 。 相比之下,图 1 中满足 (2.3) 的方案没有这种限制,因此也没有 (2.4)。在高衰减极限下,后一种方案的行为与 SBDF 方案[2]在性质上相似。另一方面,图 2 中的 方案表明,在允许的时间步长方面,它在很大的参数值范围内具有潜在优势。
Fig. 2. Time step stability restrictions for various values of and , for IMEX schemes satisfying (2.2). Note that some of these schemes, especially , perform well when is not too large in magnitude, but they all have time step restrictions in the high-attenuation limit.
图 2.满足 (2.2) 的 IMEX 方案在不同 值下的时间步长稳定性限制。请注意,其中一些方案,尤其是 ,在 的量级不太大时表现良好,但它们在高衰减极限时都有时间步长限制。

4. Finite-difference approximations in 1D
4.一维有限差分近似值

4.1. A linear problem
4.1.线性问题

Consider the one-dimensional variable-coefficient problem
考虑一维变系数问题
subject to periodic boundary conditions on the interval and initial condition
取决于区间 上的周期性边界条件和初始条件
We use centered, second-order differences for the spatial derivatives. Note that the solution is smooth for all .
我们使用居中的二阶差分求空间导数。请注意,对于所有 ,解都是平滑的。

4.1.1. Large-viscosity examples
4.1.1.大粘度示例

To examine the behavior of IMEX RK schemes for large viscosities, the model problem (4.1) was approximated using the time step . For the spatial discretization, we chose a grid spacing of or . Utilizing several schemes, computations to time were performed for viscosities, , in the range . For , these values correspond to mesh Reynolds numbers, , in the range . For smaller is correspondingly smaller.
为了检验 IMEX RK 方案在大粘度情况下的表现,我们使用 时间步长 对模型问题(4.1)进行了近似。对于空间离散化,我们选择了 的网格间距。利用几种方案,对时间 的粘度 进行了计算,粘度范围为 。对于 ,这些值对应的网格雷诺数 范围内。对于更小的 则相应更小。
The relative errors measured in the maximum norm are plotted against in Fig. 3. In Fig. 3(a), the errors for the three- and four-stage schemes are omitted since they are always within a factor of 2 of the two-stage, second-order scheme . As expected from the theory of the previous section, the
图 3 中是以最大常模测量的相对误差与 的对比图。在图 3(a) 中,省略了三阶和四阶方案的误差,因为它们与两阶二阶方案 的误差始终在 2 倍以内。正如前一节的理论所预期的那样,三阶和四阶方案的误差都在 2 倍以内。
(a)
(c)
(b)
(d)
Fig. 3. Large-viscosity behavior for .
图 3. 时的大粘度行为。
three-stage DIRK and the schemes satisfying property (2.3) allow for the largest stable time steps when is large.
较大时,三阶段 DIRK 和满足特性(2.3)的方案可以实现最大的稳定时间步长。
When we decrease , keeping fixed, the methods and become unstable over the entire depicted range of . Among the stable methods, the higher order of and becomes more apparent as the spatial error ceases to dominate the temporal error. For , only the methods satisfying (2.3) remain stable. No effective error reduction is detected for these values of and .
当我们减小 时,在保持 不变的情况下,方法 的整个描述范围内变得不稳定。在稳定的方法中,随着空间误差不再主导时间误差, 的高阶方法变得更加明显。对于 ,只有满足 (2.3) 的方法保持稳定。对于 这两个值,没有检测到有效的误差减小。
Fig. 4. Large-viscosity behavior for .
图 4. 时的大粘度行为。
Doubling the time step while keeping also causes the two-stage methods and to become unstable for almost the entire -interval depicted. Nonetheless, the three- and four-stage methods remain stable over the entire interval and the scheme is stable over most of the interval. Even upon tripling the time step (see Fig. 4), the three- and four-stage schemes are appropriate for approximating strongly-damped flows. Indeed, a comparison with the results of [2] indicates that the three-stage method is more accurate and stable in this example than the popular CNAB scheme with one-third the step size, even for this relatively coarse step size .
将时间步长加倍,同时保持 也会导致两阶段方法 在几乎整个 -区间内变得不稳定。不过,三阶段和四阶段方法在整个区间内保持稳定,方案 在大部分区间内保持稳定。即使将时间步长增加两倍(见图 4),三阶段和四阶段方案也适合用于近似强阻尼流。事实上,与文献[2]的结果比较表明,在本例中,三阶段方法 比流行的 CNAB 方案(步长为其三分之一)更精确、更稳定,即使在步长相对较粗的情况下 也是如此。

4.1.2. Small-viscosity examples
4.1.2.小粘度示例

To examine the behavior of IMEX RK methods for large mesh Reynolds numbers, example (4.1) was approximated using discretization step sizes in space and in time. Using several IMEX schemes, computations to time are performed for viscosities, , in the range . These values correspond to mesh Reynolds numbers, , in the range .
为了检验 IMEX RK 方法在大网格雷诺数下的表现,使用空间离散步长 和时间离散步长 对示例 (4.1) 进行了近似计算。使用几种 IMEX 方案,对 范围内的粘度 进行了时间 的计算。这些值对应的网格雷诺数 范围内。
The relative errors in maximum norm for these computations are plotted against in Fig. 5. Here, we find that the errors for the third-order schemes and the two-stage, second-order DIRK nearly coincide since they arise predominantly from the spatial discretization. It is especially noteworthy that the errors for these four schemes are always smaller than those arising for the two-step IMEX schemes considered in [2] even though twice the time step is applied.
图 5 是这些计算的最大常模相对误差与 的对比图。在这里,我们发现三阶方案和两阶二阶 DIRK 的误差几乎是一致的,因为它们主要产生于空间离散化。特别值得注意的是,这四种方案的误差始终小于 [2] 中考虑的两步 IMEX 方案的误差,即使采用了两倍的时间步长。

4.2. Burgers equation 4.2.伯格斯方程

Qualitatively similar results to those obtained above for the linear problem (4.1) are obtained for the Burgers equation
布尔格斯方程的计算结果与上述线性问题 (4.1) 的计算结果基本相似
Fig. 5. Small-viscosity example for .
图 5. 的小粘度示例。
with which we have experimented, subject to homogeneous Dirichlet boundary conditions and the initial conditions
我们已对其进行了试验,受均相 Dirichlet 边界条件和初始条件的限制
For small , a boundary layer develops near . This necessitates the use of upwinding. Note that upwinding can be viewed as an addition of a discrete diffusion term (where is the spatial step-size) to a centered-difference approximation of .
较小的情况下,边界层会在 附近形成。这就需要使用上卷法。请注意,上卷可以看作是在 的中心差分近似值上增加了一个离散的 扩散项(其中 是空间步长)。

5. Time-dependent multigrid in 2D
5.二维随时间变化的多网格

We consider the two-dimensional convection-diffusion problem,
我们考虑二维对流扩散问题、
where . We carry out our computations on the square and consider periodic boundary conditions and the initial conditions
其中 。我们在正方形 上进行计算,并考虑周期性边界条件和初始条件
For the spatial discretization, we use standard second-order centered differences. Time stepping is carried out using a variety of IMEX RK schemes (the convection term, , is handled explicitly and the diffusion term, , is handled implicitly). This treatment yields a positive-definite, symmetric, sparse, linear system to be solved at each stage. Such systems are solved efficiently using a multigrid algorithm, the components of which are outlined in [2].
在空间离散化方面,我们使用标准的二阶中心差分。时间步进采用多种 IMEX RK 方案(对流项 是显式处理,扩散项 是隐式处理)。这种处理方法产生了一个正有限、对称、稀疏的线性系统,需要在每个阶段求解。这种系统可以使用多网格算法高效求解,其组成部分在 [2] 中作了概述。
The model problem (5.1) was approximated using a mesh size and a residual tolerance of TOL at each stage. The time step was selected to be to achieve stable
对模型问题(5.1)进行了近似处理,网格大小为 ,每个阶段的残差公差为 TOL 。时间步长选择为 ,以实现稳定的
Table 1 表 1
Viscosity,
Scheme 0.01 0.02 0.03 0.04 0.05
Forward-backward Euler 1.94 1.90 1.71 1.29 1.19
Forward-backward Euler 1.97 1.94 1.48 1.19 1.19
Implicit-explicit midpoint 2.23 2.87 3.13 3.29 3.68
Third-order combination 4.32 4.42 3.68 3.71 3.81
Griepentrog's scheme 4.58 5.23 5.39 5.61 5.90
L-stable, 2-stage, 2nd-order DIRK 3.81 3.35 2.74 2.58 2.61
L-stable, 2-stage, 2nd-order DIRK 3.90 3.35 2.81 2.58 2.58
L-stable, 3-stage, DIRK 5.87 5.35 3.81 3.61 3.58
L-stable, 4-stage, 3rd-order scheme 7.84 7.35 5.16 4.81 4.77
results. For several IMEX RK schemes, the average number of fine-grid iterations at each time step was computed. The results for -cycles are provided in Table 1 . In addition to the schemes developed in Section 2, we include results for Griepentrog's scheme [6,7], which is not of the form (2.1).
结果。对于几种 IMEX RK 方案,计算了每个时间步的细网格平均迭代次数。表 1 列出了 循环的结果。除了第 2 节中的方案外,我们还包括 Griepentrog 方案[6,7]的结果,该方案的形式不是 (2.1)。
From the table, we see that the strongly-damping L-stable schemes require the fewest fine-grid iterations per stage for large-viscosity problems. Weakly-damping schemes require far more effort to solve the implicit equations accurately, because lingering high-frequency modes necessitate more work on the finest grid. This is especially evident for the implicit-explicit midpoint scheme since more than three iterations were required per stage when . Griepentrog's scheme also requires extra fine-grid iterations, presumably to damp out the high-frequency components which arise during the initial stage.
从表中可以看出,对于大粘度问题,强阻尼 L 稳定方案每阶段所需的细网格迭代次数最少。弱阻尼方案需要更多的努力才能准确求解隐式方程,因为残留的高频模式需要在最细网格上做更多的工作。这一点在隐式-显式中点方案 中尤为明显,因为当 时,每个阶段需要进行三次以上的迭代。Griepentrog 的方案也需要额外的细网格迭代,可能是为了抑制初始阶段出现的高频成分。

References 参考资料

[1] S. Abarbanel, D. Gottlieb and M. Carpenter, On the removal of boundary errors caused by Runge-Kutta integration of nonlinear partial differential equations, SIAM J. Sci. Comput. 17 (1996) 777-782.
[1] S. Abarbanel, D. Gottlieb and M. Carpenter, On the removal of boundary errors caused by Runge-Kutta integration of nonlinear partial differential equations, SIAM J. Sci. Comput.17 (1996) 777-782.
[2] U. Ascher, S. Ruuth and B. Wetton, Implicit-explicit methods for time-dependent PDE's, SIAM J. Numer. Anal. 32 (1995) 797-823.
[2] U. Ascher, S. Ruuth and B. Wetton, Implicit-explicit methods for time-dependent PDE's, SIAM J. Numer.Anal.32 (1995) 797-823.
[3] C. Canuto, M.Y. Hussaini, A. Quarteroni and T.A. Zang, Spectral Methods in Fluid Dynamics (Springer, 1987).
[3] C. Canuto、M.Y. Hussaini、A. Quarteroni 和 T.A. Zang,《流体力学光谱方法》(施普林格出版社,1987 年)。
[4] M. Carpenter, D. Gottlieb, S. Abarbanel and W. Don, The theoretical accuracy of Runge-Kutta time discretizations for the initial-boundary value problem: A study of the boundary error, SIAM J. Sci. Comput. 16 (1995) 1241-1252.
[4] M. Carpenter, D. Gottlieb, S. Abarbanel and W. Don, The theoretical accuracy of Runge-Kutta time discretizations for the initial-boundary value problem: A study of the boundary error, SIAM J. Sci.16 (1995) 1241-1252.
[5] M. Crouziex, Une méthode multipas implicite-explicite pour l'approximation des équations d'évolution paraboliques, Numer. Math. 35 (1980) 257-276.
[5] M. Crouziex, Une méthode multipas implicite-explicite pour l'approximation des équations d'évolution paraboliques, Numer.35 (1980) 257-276.
[6] E. Griepentrog, Gemischte Runge-Kutta-verfahren für steife systeme, in: Seminarbereicht Sekt. Math. 11 (Humboldt-Universität, Berlin, 1978) 19-29.
[6] E. Griepentrog, Mixed Runge-Kutta methods for stiff systems, in: Seminarbereicht Sekt.11 (Humboldt University, Berlin, 1978) 19-29.
[7] E. Hairer, S.P. Nørsett and G. Wanner, Solving Ordinary Differential Equations I: Nonstiff Problems (Springer, New York, 1993).
[7] E. Hairer, S.P. Nørsett and G. Wanner, Solving Ordinary Differential Equations I. Nonstiff Problems (Springer, New York, 1993):Nonstiff Problems (Springer, New York, 1993).
[8] E. Hairer and G. Wanner, Solving Ordinary Differential Equations II: Stiff and Differential-Algebraic Problems (Springer, New York, 1991).
[9] G.E. Karniadakis, M. Israeli and S.A. Orszag, High-order splitting methods for the incompressible NavierStokes equations, J. Comput. Phys. 97 (1991) 414-443.
[9] G.E. Karniadakis, M. Israeli and S.A. Orszag, High-order splitting methods for the incompressible NavierStokes equations, J. Comput.97 (1991) 414-443.
[10] J. Kim and P. Moin, Application of a fractional-step method to incompressible Navier-Stokes equations, J. Comput. Phys. 59 (1985) 308-323.
[10] J. Kim and P. Moin, Application of a fractional-step method to incompressible Navier-Stokes equations, J. Comput.59 (1985) 308-323.
[11] S. Ruuth, Implicit-explicit methods for reaction-diffusion problems in pattern formation, J. Math. Biol. 34 (2) (1995) 148-176.
[11] S. Ruuth, Implicit-explicit methods for reaction-diffusion problems in pattern formation, J. Math. Biol.34 (2) (1995) 148-176.
[12] R.D. Skeel, G. Zhang and T. Schlick, A family of symplectic integrators: stability, accuracy, and molecular dynamics applications, SIAM J. Sci. Comput. 18 (1) (1997) 203-222.
[12] R.D. Skeel, G. Zhang and T. Schlick, A family of symplectic integrators: stability, accuracy, and molecular dynamics applications, SIAM J. Sci. Comput.18 (1) (1997) 203-222.
[13] R. Spiteri, U. Ascher and D. Pai, Numerical solution of differential systems with algebraic inequalities arising in robot programming, in: Proceedings of the IEEE Conference on Robotics and Automation (1995).
[13] R. Spiteri, U. Ascher and D. Pai, Numerical solution of differential systems with algebraic inequalities arising in robot programming, in:电气和电子工程师学会机器人与自动化会议论文集(1995 年)。
[14] S. Turek, A comparative study of time stepping techniques for the incompressible Navier-Stokes equations: from fully implicit non-linear schemes to semi-explicit projection methods, Internat. J. Numer. Methods Fluids 22 (1996) 987-1011.
[14] S. Turek, A comparative study of time stepping techniques for the incompressible Navier-Stokes equations: from fully implicit non-linear schemes to semi-explicit projection methods, Internat.Numer.Methods Fluids 22 (1996) 987-1011.
[15] J.M. Varah, Stability restrictions on second order, three level finite difference schemes for parabolic equations, SIAM J. Numer. Anal. 17 (2) (1980) 300-309.
[15] J.M. Varah, Stability restrictions on second order, three level finite difference schemes for parabolic equations, SIAM J. Numer. Anal.Anal.17 (2) (1980) 300-309.
[16] J.G. Verwer, J.G. Blom and W. Hundsdorfer, An implicit-explicit approach for atmospheric transportchemistry problems, Appl. Numer. Math. 20 (1996) 191-209.
[16] J.G. Verwer, J.G. Blom and W. Hundsdorfer, An implicit-explicit approach for atmospheric transportchemistry problems, Appl. Numer.Math.20 (1996) 191-209.

  1. Corresponding author. E-mail: ascher@cs.ubc.ca.
    通讯作者。电子邮件:ascher@cs.ubc.ca.
    The work of this author was partially supported under NSERC Canada Grant OGP0004306.
    作者的部分研究工作得到了加拿大国家科学研究中心(NSERC Canada)OGP0004306 号基金的资助。
    The work of this author was partially supported by an NSERC Postdoctoral Scholarship and NSF DMS94-04942. E-mail: ruuth@math.ucla.edu.
    本文作者的研究工作得到了美国国家科学研究中心博士后奖学金和国家自然科学基金 DMS94-04942 的部分资助。电子邮件:ruuth@math.ucla.edu.
    The work of this author was partially supported under NSERC Canada Grant OGP0004306. E-mail: spiteri@math.ubc.ca. 0168-9274/97/$17.00 (C) 1997 Elsevier Science B.V. All rights reserved.
    作者的部分研究工作得到了加拿大国家科学研究中心(NSERC Canada)OGP0004306 号基金的资助。电子邮件:spiteri@math.ubc.ca.0168-9274/97/$17.00 (C) 1997 Elsevier Science B.V. 版权所有。保留所有权利。
  2. W. Hunsdorfer (Workshop on Innovative Time Integrators, CWI, Amsterdam, 1996) has indicated that a more general test equation where and are allowed to be complex is useful, but we stay with the simpler test equations and lend further support to the analysis with more general numerical experiments in the following sections.
    W. Hunsdorfer(创新时间积分器研讨会,CWI,阿姆斯特丹,1996 年)指出,允许 复杂的更一般的测试方程是有用的,但我们仍然使用较简单的测试方程,并在下面的章节中通过更一般的数值实验为分析提供进一步的支持。
  3. By large viscosity, we mean that the mesh Reynolds number is small ; see, e.g., [2]). An estimate for the mesh Reynolds number is given by where represents viscosity, represents characteristic speed, and represents grid spacing.
    粘度大,意味着网格雷诺数小 ;见,例如,[2])。网格雷诺数的估计值为 ,其中 代表粘度, 代表特征速度, 代表网格间距。