Time marcher
Integrating-factor technique
I would like to integrate the following equations in time:
\[\der{\wav{\ux}{\ix \iy}}{t}
+
\frac{1}{Re} \left( \mkx^2 + \mky^2 \right) \wav{\ux}{\ix \iy}
=
\wav{f_x}{\ix \iy},\]
\[\der{\wav{\uy}{\ix \iy}}{t}
+
\frac{1}{Re} \left( \mkx^2 + \mky^2 \right) \wav{\uy}{\ix \iy}
=
\wav{f_y}{\ix \iy},\]
\[\der{\wav{T}{\ix \iy}}{t}
+
\frac{1}{Re Sc} \left( \mkx^2 + \mky^2 \right) \wav{T}{\ix \iy}
=
\wav{g}{\ix \iy},\]
where \(f_x\), \(f_y\), and \(g\) are the other terms such as the non-linear terms.
In general, these equations have the following form:
\[\der{p}{t}
+
C p
=
q,\]
where \(C\) is a constant.
Here I introduce an integrating factor
\[\exp{
\left( C t \right)
}\]
and multiply it to the equation:
\[\exp{
\left( C t \right)
}
\left(
\der{p}{t}
+
C p
\right)
=
\exp{
\left( C t \right)
}
q.\]
The left-hand-side terms are re-arranged as
\[\der{}{t}
\left\{
p
\exp{
\left( C t \right)
}
\right\},\]
and thus I obtain the resulting differential equation:
\[\der{P}{t}
=
Q,\]
where
\[P
\equiv
p
\exp{
\left( C t \right)
}\]
and
\[Q
\equiv
q
\exp{
\left( C t \right)
}.\]
Up to this point, no approximations have been employed, ensuring that all relations are analytically accurate.
From this point onward, I will introduce an approximation for the right-hand-side integral in order to handle it numerically.
Specifically, I will utilise the classical fourth-order Runge-Kutta scheme in conjunction with the integrating-factor technique.
General explicit Runge-Kutta schemes
Here I aim to numerically solve
\[\der{P}{t}
=
Q,\]
where
\[p
\equiv
P
\exp{
\left( C t \right)
}\]
and
\[Q
\equiv
q
\exp{
\left( C t \right)
},\]
where \(\exp{\left( C t \right)}\) is the integrating factor.
A general explicit Runge-Kutta scheme reads
\[\begin{split}& P^1 = P^n + Q^0 a_{10} \Delta t \\
& P^2 = P^n + Q^0 a_{20} \Delta t + Q^1 a_{21} \Delta t \\
& P^3 = P^n + Q^0 a_{30} \Delta t + Q^1 a_{31} \Delta t + Q^2 a_{32} \Delta t \\
& \vdots\end{split}\]
The \(k\)-th row is written as
\[P^k = P^n + \sum_{l = 0}^{k - 1} Q^l a_{kl} \Delta t.\]
Note that
\[t^0 = t^n,
P^0 = P^n.\]
Integrating factor
For notational simplicity, I define
\[E \left( x \right)
\equiv
\exp{
\left(
C x
\right)
}.\]
Assigning \(P\) and \(Q\) to the above relation yields
\[p^k
E \left(
t^n
+
c_k \Delta t
\right)
=
p^n
E \left(
t^n
\right)
+
\sum_{l = 0}^{k - 1}
q^l
E \left(
t^n
+
c_l \Delta t
\right)
a_{kl} \Delta t.\]
Dividing the equation by \(E \left( t^n \right)\) leads to
\[p^k
E \left(
c_k \Delta t
\right)
=
p^n
+
\sum_{l = 0}^{k - 1}
q^l
E \left(
c_l \Delta t
\right)
a_{kl} \Delta t.\]
Fourth-order Runge-Kutta scheme
The Butcher tableau of the classical fourth-order scheme
\[\begin{split}\begin{array}{c|cccc}
c_0 & & & & \\
c_1 & a_{10} & & & \\
c_2 & & a_{21} & & \\
c_3 & & & a_{32} & \\
\hline
c_4 & a_{40} & a_{41} & a_{42} & a_{43} \\
\end{array}\end{split}\]
is
\[\begin{split}\begin{array}{c|cccc}
0 & & & & \\
1/2 & 1/2 & & & \\
1/2 & & 1/2 & & \\
1 & & & 1 & \\
\hline
1 & 1/6 & 1/3 & 1/3 & 1/6 \\
\end{array}\end{split}\]
The whole process to update a field from
\(p^n\) to \(p^{n+1}\) is as follows:
\[p^0 = p^n.\]
\[\newcommand{\va}{1}
\newcommand{\vb}{0}
p^{\va}
E \left(
c_{\va} \Delta t
\right)
=
p^n
+
q^{\vb}
E \left(
c_{\vb} \Delta t
\right)
a_{\va \vb} \Delta t.\]
\[\newcommand{\va}{2}
\newcommand{\vb}{1}
p^{\va}
E \left(
c_{\va} \Delta t
\right)
=
p^n
+
q^{\vb}
E \left(
c_{\vb} \Delta t
\right)
a_{\va \vb} \Delta t.\]
\[\newcommand{\va}{3}
\newcommand{\vb}{2}
p^{\va}
E \left(
c_{\va} \Delta t
\right)
=
p^n
+
q^{\vb}
E \left(
c_{\vb} \Delta t
\right)
a_{\va \vb} \Delta t.\]
\[\begin{split}p^{4}
E \left(
c_{4} \Delta t
\right)
& =
p^n \\
& +
q^{0}
E \left(
c_{0} \Delta t
\right)
a_{40} \Delta t \\
& +
q^{1}
E \left(
c_{1} \Delta t
\right)
a_{41} \Delta t \\
& +
q^{2}
E \left(
c_{2} \Delta t
\right)
a_{42} \Delta t \\
& +
q^{3}
E \left(
c_{3} \Delta t
\right)
a_{43} \Delta t.\end{split}\]
\[p^{n+1}
=
p^4.\]
Here, \(p^k\) is stored to a buffer which is allocated to store the intermediate field.
Since \(a_{kl}\) is all zero except \(k - l = 1\), only one additional buffer is needed for this purpose.
On the other hand, the derivatives \(q^k\) are all to be stored since the classical scheme is not a low-storage scheme.