5. Convergence Analysis

In this chapter, we will prove convergence of the vorticity redistribution method for the Stokes or heat equation (3.21). The redistribution fractions $f_{ij}^n$ are assumed to satisfy the redistribution equations (4.8) and following, to satisfy the positivity constraint (4.13), and to be restricted to vortices within a mutual distance (4.2), with $R>1$. Since we are enforcing consistency in the $L_2$ norm, using the Fourier transform, while we have stability in the $l_1$ norm, and the redistribution fractions are only partly determined, the conventional convergence arguments need some modifications.

We will show convergence in the $L_2$ norm by showing convergence of the Fourier transform of the numerical solution,

$$\widehat\omega^n = \widehat\phi(k \delta) \sum_i \Gamma_i^n e^{-{\rm i}\vec k\vec x_i}$$ (5.1)

to the Fourier transform of the exact solution,

$$\widehat\omega(t) = \widehat\omega_0 e^{- k^2 \nu t}$$ (5.2)

in the $L_2$ norm. Here $\omega_0$ is the given initial vorticity and we do not explicitly show the dependence on $\vec k$.

The total error consists of the error induced by discretizing the initial data and the error induced by the redistribution method itself:

$$\Vert\widehat\omega_0 e^{-k^2\nu t} -\widehat\omega^n\Vert \... ...Vert\widehat\omega^0 e^{-k^2\nu t} - \widehat\omega^n\Vert \ .$$ (5.3)

The first error due to discretizing the initial data can be important if the initial data have only limited smoothness or if a low-order smoothing function is used. It depends on how the initial discretization is performed. Typically the initial vortices are given a uniform spacing $h = O(h_v)$ and the initial vortex strength is taken as $\Gamma_i^0=h^2 \omega_0({\vec x}_i)$. Since the initial vorticity field is evaluated only at the vortices, some information is lost; aliasing makes $\omega_0$ indistinguishable from the Fourier interpolant $\omega_h$ through the vorticity values. The total error due to discretization of the initial data can be written:

$$\Vert(\widehat\omega_0 - \widehat\omega^0)e^{-k^2\nu t}\Vert... ...ert(\widehat\omega_h - \widehat\omega^0)e^{-k^2\nu t}\Vert \ .$$ (5.4)

The magnitude of the first of these two errors depends on the number of square integrable derivatives of the initial vorticity. It may be shown that if $\sigma$ derivatives are square integrable, this error is of order $h^\sigma$ ([219] pp. 198-206). In two dimensions $\sigma$ has to be greater than one, but fractional values are allowed.

The second error is due to the vortex core. Assuming $\widehat\phi_\delta$ to be bounded, for nonzero times the order of this error is simply the order of accuracy of the vortex core. Thus, if the core is accurate $O(h_v^M)$, the overall accuracy of the computation is not affected by the core.

It follows that for sufficiently accurate smoothing function and smooth initial data, the only important error will be that due to the redistribution process. To estimate this error, we first define the local error in the Fourier transform at time-level $n$ to be the difference between the redistribution solution and the exactly diffused solution from the previous time-step:

$$\widehat \epsilon^n \equiv \widehat\omega^{n+1} - \widehat\omega^n e^{-k^2 \nu \Delta t}\ .$$ (5.5)

By repeated application of this definition, the error in the Fourier transform due to redistribution can be bounded by

$$\vert\widehat\omega^n - \widehat\omega^0 e^{-k^2\nu t}\vert ... ...ac{\widehat\epsilon}{k^2 h_v^2} \left(1 + k^2 h_v^2\right) \ ,$$ (5.6)

where $h_v=\sqrt{\nu\Delta t}$ and $\widehat\epsilon=\max_{i=0}^{n-1}\{\vert\widehat\epsilon^i\vert\}$.

To estimate $\widehat\epsilon$, recall from chapter 4 that the redistribution equations (4.8) and following ensure vanishing of the first few powers of $\Delta t$ in the error in (5.5). The Taylor series remainder theorem can be used to express the remaining difference $\widehat\epsilon^n$. That expression is shown in Appendix A; it can be bounded as

$$\widehat\epsilon \le \vert\widehat\phi(k \delta)\vert \max_n(\Gamma^n F^n) R^{M+2} (k h_v)^{M+2} (1 + k^2 h_v^2)\ ,$$ (5.7)

$$\qquad F^n = \max_i \sum_j \vert f^n_{ij}\vert\ , \qquad \Gamma^n = \sum_i \vert\Gamma_i^n\vert\ .$$ (5.8)

For the assumed positivity of the redistribution fractions (4.13), $F^n=1$ and $\Gamma^n$ cannot increase. Thus the total error in the Fourier transform is bounded by

$$\vert\widehat\omega^n - \widehat\omega^0 e^{-k^2\nu t}\vert ... ...t\phi(k \delta)\vert 4 \Gamma^0 R^{M+2} \min\{(k h_v)^M,1\}\ ,$$ (5.9)

where the second bound comes from the bound $\vert\phi_\delta\vert\Gamma^n$ to (5.1).

In this work, we will assume that $\Gamma^0$, the absolute circulation of the discretized initial data, is finite. Note that this is a restriction on the $l_1$ norm of the initial discrete vortex strengths, rather than on the $L_2$ norm of the initial vorticity distribution. However, the Cauchy inequality applied to

$$\sum_i h \vert \omega_i\vert (1 + x_i^2 + y_i^2)^{(1+\alpha)/2} \cdot h (1 + x_i^2 + y_i^2)^{-(1+\alpha)/2}\ ,$$ (5.10)

with $\alpha$ an arbitrary positive constant, readily shows that $\Gamma$ can be bounded in terms of the $L_2$ norms of the initial vorticity and aliasing error provided that the initial vorticity is restricted to a finite region or at least decays sufficiently rapidly at large distances. For example, it would suffice that $\omega_0 = O(x^{-2-\alpha})$ for $x \to \infty$ for some $\alpha>0$.

The final $L_2$ error in the vorticity is found from square integration of (5.9) over all wavenumbers. Thus the error due to redistribution is found to be:

$$\Vert\widehat\omega^n - \widehat\omega^0 e^{-k^2\nu t}\Vert ... ...vert\widehat\phi(k)\vert^2 k^{2M+1} {\rm d} k\right)^{1/2} \ .$$ (5.11)

To minimize this error, a relatively large core size is desirable. If we take the core size proportional to some small power $\alpha$ of $h_v$, the error will be $O(h_v^{M-\alpha(M+1)})$. Since we can take $\alpha$ as any positive number, we can obtain any order of accuracy arbitrarily close to $O(h_v^M)$. Note however that for a core with a finite order of accuracy, the first error in (5.3), due to discretizing the initial data, limits the maximum size of $\delta$. We will discuss using a smoothing function to evaluate the vorticity further in the next chapter.

This completes the discussion of convergence for the Stokes equations. It is interesting to note that the true stability conditions are that $F^n$ and $\Gamma^n$ are bounded. Next we need to address how to find the redistribution fractions $f^n_{ij}$ in an actual application of the scheme. We will address this in the next chapter on the numerical implementation of the convection and diffusion steps described in section 3.2.