Also at ]Oxford Suzhou Centre for Advanced Research

Also at ]Oxford Suzhou Centre for Advanced Research

# Tracking the vortex motion by using Brownian fluid particles

Zhongmin Qian<EMAIL_ADDRESS>Mathematical Institute, University of
Oxford, OX2 6GG, England  [  Youchun Qiu  Institut de Mathématiques de
Toulouse, UMR 5219, Université de Toulouse, CNRS, UPS, F-31062, Toulouse Cedex
9, France.  Yihuang Zhang  Mathematical Institute, University of Oxford, OX2
6GG, England  [

###### Abstract

In this paper we propose a simple yet powerful vortex method to numerically
approximate the dynamics of an incompressible flow. The idea is to sample the
distribution of the initial vortices of the fluid flow in question then follow
vortex dynamics along Taylor’s Brownian fluid particles. The weak convergences
of this approximation scheme are obtained for 2D and 3D cases, though only for
small time in 3D case. Based on our method, the simulation results are quite
attracting.

††preprint: AIP/123-QED

##  I Introduction

Numerical methods have become important components in the study of fluid
dynamics, in particular for modeling turbulence flows. With the advance of
computational power DNS (Direct Numerical Simulation), LES (Large Eddy
Simulation) and other new technologies have been developed in recent years for
solving the Navier-Stokes equations numerically. Among them, various vortex
methods which are based on the vorticity transport equation have become an
attractive approach for simulating fluid flows in particular turbulent flows.
One may find a comprehensive account in the monographs Cottet and
KoumoutsakosCottet and Koumoutsakos (2000), Majda and BertozziMajda and
Bertozzi (2002), SaffmanSaffman (1995), Ting and KnioTing and Knio (2007), and
the recent review Mimeau and MortazaviMimeau and Mortazavi (2021).

The idea of vortex methods was introduced in ChorinChorin (1973). Several
vortex approximation procedures and convergence results have been established
for 2D flows, see for example Anderson et al.Anderson and Greengard (1985);
Beale and Majda (1982); Hald and del Prete (1978); Hald (1979, 1987) and the
literature therein. For 3D vortex methods of inviscid fluid flows, the
convergence with Lagrangian stretching was proved in Beale and MajdaBeale and
Majda (1982). A different approach which updates the vorticity through the
velocity field was proposed and the corresponding convergence result was shown
in BealeBeale (1986). For viscous fluid flows, the vortex dynamics is replaced
by a random dynamical system in which a Brownian motion term is added to the
equation of motion of the fluid particles, and the fluid particles become
Brownian fluid particles. Itô’s stochastic differential equations (SDEs) take
place of ordinary differential equations (ODEs). A remarkable convergence
result for 2D random vortex method has been established in Long Long (1988).
For 3D viscous fluid flows however, to the best of authors’ knowledge, there
is no satisfactory solution yet so far.

In the present paper, we propose a simple method of tracking the vortex
dynamics of an incompressible fluid which gives surprisingly satisfactory
simulations for both inviscid and viscous fluid flows. Our method is to
introduce the vorticity evolution directly in terms of Brownian fluid
particles specified by a set of stochastic differential equations. We do not
evolve the Taylor diffusion by mollifying the Biot-Savart kernel as in the
traditional vortex methods, but instead we sample the distribution of the
initial vortices and develop the initial distribution according to the SDEs
determined by the vorticity equation. Let us describe this approach in more
detail and at the same time establish the notations we will use throughout the
paper.

Let $u=(u^{1},u^{2},u^{3})$ denote the velocity of an incompressible fluid
flow moving in a range without boundary constraint. Hence the velocity
$u(x,t)$ satisfies the equations of motion, the Navier-Stokes equations

$\frac{\partial}{\partial t}u^{i}+u^{j}\frac{\partial}{\partial
x^{j}}u^{i}=\nu\Delta u^{i}-\frac{\partial}{\partial
x^{i}}p,\quad\frac{\partial}{\partial x^{k}}u^{k}=0$ (1)

in $\mathbb{R}^{3}$, where $i=1,2,3$, $\nu\geq 0$ is the kinetic viscosity and
$p(x,t)$ is the pressure which is uniquely determined by $u(x,t)$ up to a
constant at every $t$. Einstein’s convention that the term with a pair of
repeated indices are summed over from 1 to 3 has been applied. The case where
$\nu=0$ corresponds to inviscid fluid flows, and the Navier-Stokes equations
are reduced to the Euler equations.

The vorticity of $u$ is denoted by
$\omega=(\omega^{1},\omega^{2},\omega^{3})$, although, which will be
specified, we will depart from this convention in order to present some
results in a general setting. By definition
$\omega^{i}=\varepsilon^{ijk}\frac{\partial}{\partial x^{j}}u^{k}$ for
$i=1,2,3$. The vorticity equations, which are the equations of vortex motion,
play a dominated role in vortex methods, and are obtained by differentiating
the Navier-Stokes equations

$\frac{\partial}{\partial t}\omega^{i}+u^{j}\frac{\partial}{\partial
x^{j}}\omega^{i}=\nu\Delta\omega^{i}+\omega^{j}\frac{\partial}{\partial
x^{j}}u^{i}$ (2)

for $i=1,2,3$. The non-linear term $u\cdot\nabla\omega$ on the left-hand side
is the convection of $\omega$ which appears for both 2D and 3D flows. The
second non-linear term $\omega\cdot\nabla u$, representing the stretching of
the vorticity, appears only in 3D flows. This fact makes substantial
difference between 2D flows and 3D flows. The study of some turbulence
problems (see for example SaffmanSaffman (1995), Pullin and SaffmanPullin and
Saffman (1998)) can be formulated in terms of Cauchy’s initial value problem
to the vorticity equations (2) together with the equation that
$\omega=\nabla\wedge u$, subject to the initial vorticity
$\omega_{0}=\omega(\cdot,0)$. Vortex methods aim to provide numerical schemes
to the initial value problem.

The velocity $u(x,t)$, under the assumption that both $u$ and $\omega$ decay
sufficiently fast at the infinity, may be recovered from reading the vorticity
via the Biot-Savart law

$u^{i}(x,t)=\int_{\mathbb{R}^{3}}\varepsilon^{ijk}G^{j}(x-y)\omega^{k}(y,t)dy$
(3)

where $G=(G^{1},G^{2},G^{3})$ and

$G(z)=-\frac{1}{4\pi}\frac{z}{|z|^{3}}$

is the Biot-Savart singular integral kernel.

Now we are in a position to describe our simple random vortex dynamics. Two
(random) vector fields $V(x,t)$ and $W(x,t)$ will be defined below which do
not necessarily satisfy the Navier-Stokes equations nor the relation that
$W=\nabla\wedge V$. Our goal is in fact to construct approximation solutions
to the vorticity equation (2) directly, in the spirit which is quite like
Feynman’s functional integration for Schrödinger’s equations. Hence $V$ will
be the approximate velocity of $u$, and $W$ the approximation of $\omega$.

Following the general ideas in the vortex methods, we propose the following
dynamics scheme for the vortex motion of an incompressible fluid flow. At the
initial time $t=0$, we sample a (finite) collection of locations
$x_{n}\in\mathbb{R}^{3}$, where $n$ runs through a finite index set, at which
the major vortex motion may be demonstrated. For example $x_{n}$ can be the
centers of vortex rings. Suppose the initial vorticity of the fluid flow is
given by

$W^{i}(x,0)=\sum_{n}A_{n}^{i}(0)\varphi(x-x_{n})$ (4)

where $A_{n}(0)$, vectors with components $A_{n}^{i}(0)$ (for $i=1,2,3$) are
the initial vortices, and $\varphi$ is a wavelet type function which should be
close to the Dirac delta (at $0$) function. Hence at the initial stage, the
vortices of the fluid flow are distributed among $x_{n}$ so that

$W(x,0)\simeq\sum_{n}A_{n}(0)\delta(x-x_{n}),$

see for example Cottet and KoumoutsakosCottet and Koumoutsakos (2000).
Therefore (4) provides us with the distribution of the initial vortices in the
fluid flow. The key idea is that both the loci $x_{n}$ and the initial
vortices $A_{n}(0)$ are sampled so that they are representative for the vortex
motion at the initial stage. The initial vorticity $W(x,0)$ will be
transported along the Brownian fluid particles. More precisely the vorticity
$A_{n}(0)$ at time $t>0$ is transported to a new location $X_{n}(t)$ with a
new vorticity $A_{n}(t)=(A_{n}^{1}(t),A_{n}^{2}(t),A_{n}^{3}(t))$, so that the
distribution of the vortices at time $t>0$ is given by

$W^{i}(x,t)=\sum_{n}A_{n}^{i}(t)\varphi(x-X_{n}(t))$ (5)

for $i=1,2,3$. The velocity $V(x,t)$ at time $t>0$ is defined in terms of the
Biot-Savart law (3):

$V^{i}(x,t)=\int_{\mathbb{R}^{3}}\varepsilon^{ijk}G^{j}(x-y)W^{k}(y,t)dy,$ (6)

so that the relation between $W$ and $V$ is maintained at least partly. It
remains to determine the dynamics of $A_{n}$ and $X_{n}$. By initiative, the
velocity of $X_{n}=(X_{n}^{1},X_{n}^{2},X_{n}^{3})$ should be the velocity of
the fluid flow, so that

$dX_{n}^{i}(t)=V^{i}(X_{n}(t),t)dt+\sqrt{2\nu}dB^{i}(t),\quad X_{n}(0)=x_{n}$
(7)

for $i=1,2,3$ and $B=(B^{1},B^{2},B^{3})$ is a standard 3D Brownian motion on
a probability space $(\Omega,\cal{F},\mathbb{P})$. That is, $X_{n}$ are
Taylor’s diffusions initiated from $x_{n}$, and therefore $X_{n}$ are the
Brownian fluid particles started at $x_{n}$. The dynamics of the vortices
$A_{n}=(A_{n}^{1},A_{n}^{2},A_{n}^{3})$ are determined by the following
ordinary differential equations

$dA_{n}^{j}(t)=A_{n}^{k}(t)\frac{\partial V^{j}}{\partial
x^{k}}(X_{n}(t),t)dt,\quad A_{n}(0)=A_{n}(0)$ (8)

where $i=1,2,3$, which are responsible for the vortex stretching of the fluid
flow. For 2D fluid flows, since there is no vorticity stretching, so $A_{n}$
stay as constant vectors along the fluid particle trajectories, see the
section about 2D flows below.

In the next section, we will translate the previous system (7), (8), (5) and
(6) together into a closed system of stochastic differential equations which
thus define our random vortex method.

The main contribution of the present paper is to show that the previous simple
vortex dynamics gives rise to good approximations to the motions of vortex
dynamics. We demonstrate this by showing several theoretical results about the
approximation solutions to the vorticity equations constructed in terms of
$X_{n}$ and $A_{n}$, and also by simulations based on this simple vortex
method.

The paper is organized as the following. In the next Section 2, we describe
the system of stochastic differential equations which implement a simple
vortex method. This system of SDEs allows us to employ the Monte-Carlo
simulation to the study of incompressible fluid flows, which will be
demonstrated in Section 7. In Section 3, we derive the approximation vorticity
equation, which takes a form of a simple stochastic partial differential
equation. In Section 4 and 5 we discuss the weak convergence results with
respect to certain Sobolev norms, which justify our simple vortex method. In
Section 6, we discuss the 2D incompressible flows, and not surprisingly we
show that weak convergence holds for all time.

##  II Description of the random vortex method

In this section we propose, following the approach outlined in the
Introduction, a random vortex dynamics system in a slightly general setting,
otherwise we will maintain the notations established in the previous section.

The main change will be made for the relation between the vorticity field $W$
and the velocity $V$ which may be not given by the vector identity
$W=\nabla\wedge V$, instead it will be given in terms of a vector integral
kernel $K=(K^{1},K^{2},K^{3})$, where $K^{i}$ are locally integrable functions
on $\mathbb{R}^{3}$. The reason to work with a slightly general kernel $K$
than the Biot-Savart kernel $G$ is the following. For solving the initial
value problem to the vorticity equations (2) where $\omega=\nabla\wedge u$ and
$\textrm{div}u=0$, which are equivalent to the system (2) and (3),
mathematical difficulty arises since the Biot-Savart kernel $G$ is singular,
so it is natural to replace $G$ by its smooth approximations. The simple way
to create an approximation is to mollify the Biot-Savart kernel $G$. That is,
choosing a smooth function $\psi\geq 0$ with a compact support in
$(-\frac{1}{2},\frac{1}{2})^{3}$ and $\int_{\mathbb{R}^{3}}\psi(x)dx=1$. For
each $\delta>0$, set $\psi_{\delta}(x)=\delta^{-3}\psi(\delta^{-1}x)$ and
$G_{\delta}(x)=G\star\psi_{\delta}(x)$, where $\star$ denotes the convolution,
i.e.

$G_{\delta}(x)=\int_{\mathbb{R}^{3}}G(y)\psi_{\delta}(x-y)dy.$ (9)

Then $G_{\delta}$ is smooth and

$\left\|D^{k}G_{\delta}\right\|\leq C_{k}\frac{1}{\delta^{2+k}}$ (10)

for some constant $C_{k}>0$, for $k=1,2,\cdots$, depending only on the
regularization $\psi$. Moreover $G_{\delta}\rightarrow G$ as $\delta\downarrow
0$ in distribution sense. Therefore the initial value problem to the following
system

$\frac{\partial}{\partial t}\omega^{i}+u^{j}\frac{\partial}{\partial
x^{j}}\omega^{i}=\nu\Delta\omega^{i}+\omega^{j}\frac{\partial}{\partial
x^{j}}u^{i}$ (11)

with the initial data $\omega(x,0)=\omega_{0}(x)$, and

$u^{i}(x,t)=\int_{\mathbb{R}^{3}}\varepsilon^{ijk}G_{\delta}^{j}(x-y)\omega^{k}(y,t)dy$
(12)

for $\delta>0$ gives rise to approximation solutions, denoted by
$\omega^{\delta}$ and $u^{\delta}$. In fact as long as the initial data
$\omega_{0}$ is smooth with a compact support, then
$\omega^{\delta}\rightarrow\omega$ as $\delta\downarrow 0$ with respect to
some Sobolev norm at least for small time. Therefore, from the computational
view-point, we only need to develop numerical schemes for the approximation
solutions $\omega^{\delta}$ and $u^{\delta}$, thus it is important to work
with singular kernels $K$ such as $G_{\delta}$. While we would like to point
out that the procedure for going to the approximation equations (11) and (12)
seems unnecessary for implementing the random vortex method below, rather than
for the technical reason that a priori estimates are not available for 3D
Navier-Stokes equations, see Lemma 1 below.

Let us now define the random vortex system. In defining the distribution of
the initial vortices (4) the wavelet type function $\varphi$ is assumed to be
smooth with a compact support about the original $0$. Let us assume that the
support of $\varphi\geq 0$ lies inside the ball at $0$ with radius
$r_{\varphi}>0$ which will be chosen to be small, and the total mass
$\int_{\mathbb{R}^{3}}\varphi(x)dx=1$. These are the structure data for our
random vortex scheme, which are fixed if we are not care about the convergence
issue.

Let $W(x,0)$ and $W(x,t)$ (for $t>0$) be given by (4) and (5) respectively,
where $x_{n}$ and $A_{n}(0)$ are sampled initially, so they are the fixed data
too. We then define the vector field $V(x,t)$ by the convolution (together
with the wedge product)

$V^{i}(x,t)=\int_{\mathbb{R}^{3}}\varepsilon^{ijk}K^{j}(x-y)W^{k}(y,t)dy$ (13)

which coincides with (3) if $K$ is the Biot-Savart kernel $G$. The dynamics of
$(X_{n},A_{n})$ are still defined by (7) and (8). Together with (5) we deduce
that

$V^{i}(x,t)=-\sum_{n}\varepsilon^{ikj}A_{n}^{k}(t)K_{\varphi}^{j}(x-X_{n}(t))$
(14)

for $i=1,2,3$, $t\geq 0$ and $x\in\mathbb{R}^{3}$, where

$K_{\varphi}^{j}(x)=\int_{\mathbb{R}^{3}}K^{j}(x-y)\varphi(y)dy$

(where $j=1,2,3$) turns out to be the mollification of $K^{j}$ by $\varphi$.
This is a nice feature in this scheme. By utilizing equation (14) and
substituting it into (7) the dynamic system for the Brownian fluid particles
$X_{n}$ can be reformulated as the following SDEs

$\displaystyle
dX_{n}^{i}(t)=-\sum_{m}\varepsilon^{ikj}A_{m}^{k}(t)K_{\varphi}^{j}(X_{n}(t)-X_{m}(t))dt+\sqrt{2\nu}dB^{i}(t),$
$\displaystyle X_{n}(0)=x_{n},$ (15)

and similarly the dynamics system (8) for $A_{n}$ may be written as the
following ODEs

$\displaystyle
dA_{n}^{i}(t)=A_{n}^{l}(t)\sum_{m}\varepsilon^{ijk}\frac{\partial
K_{\varphi}^{j}}{\partial x^{l}}(X_{n}(t)-X_{m}(t))A_{m}^{k}(t)dt,$
$\displaystyle A_{n}(0)=A_{n}(0),$ (16)

where $i=1,2,3$, $m$ and $n$ run through the finite range of the initial
locations $x_{n}$.

The system of SDEs (15, 16) is closed and depends only on the structure data
$K$ and $\varphi$. We observe that $K_{\varphi}$ is smooth if $K$ is locally
integrable under our assumption on $\varphi$.

###### Theorem 1

Given a finite collection of loci $x_{n}\in\mathbb{R}^{3}$ and a family of
initial vortices $A_{n}(0)$, there is a unique maximal strong solution
$(X_{n},A_{n})$ to SDEs (15, 16) up to the explosion time $\tau>0$. The vector
fields defined by

$W^{i}(x,t)=\sum_{n}A_{n}^{i}(t)\varphi(x-X_{n}(t))$ (17)

and

$V^{i}(x,t)=-\sum_{n}\varepsilon^{ikj}A_{n}^{k}(t)K_{\varphi}^{j}(x-X_{n}(t))$
(18)

are smooth in $x$ for $t<\tau$, where $i=1,2,3$. Moreover the explosion time
$\tau\geq\tau_{A(0),\varphi}$, where

$\tau_{A(0),\varphi}=\frac{1}{\left\|A(0)\right\|\left\|DK_{\varphi}\right\|_{\infty}},$
(19)

$\left\|A\right\|=\sum_{n}|A_{n}|$ and $\left\|f\right\|_{\infty}$ denotes the
$L^{\infty}$-norm of a function on $\mathbb{R}^{3}$. We note that
$\tau_{A(0),\varphi}>0$ is deterministic.

Proof. The proof follows from the standard result in Itô’s theory of
stochastic differential equations, see for example Theorem 2.3 on page 173 in
Ikeda and WatanabeIkeda and Watanabe (1989). The coefficients in defining SDEs
(15) and (16) are in general not globally Lipschitz, but nevertheless locally
Lipschitz continuous. Therefore the explosion time $\tau$ may be finite, and
the maximal strong solution $(X_{n},A_{n})$ is unique for a given 3D Brownian
motion $B=(B^{1},B^{2},B^{3})$ on a probability space
$(\varOmega,\mathcal{F},\mathbb{P})$. To derive the estimate (19), we utilize
a truncation technique. For each $L>0$, let $\theta_{L}(s)$ be the cutting-off
function which equals $s$ if $|s|\leq L$ and equals $L$ if $|s|>L$. The unique
strong solution $(X^{L},A^{L})$ to the truncated SDEs

$dX_{n}^{i}(t)=-\sum_{m}\varepsilon^{ikj}\theta_{L}(A_{m}^{k}(t))K_{\varphi}^{j}(X_{n}(t)-X_{m}(t))dt+\sqrt{2\nu}dB^{i}(t)$

and

$dA_{n}^{i}(t)=\theta_{L}(A_{n}^{l}(t))\sum_{m}\varepsilon^{ijk}\frac{\partial
K_{\varphi}^{j}}{\partial x^{l}}(X_{n}(t)-X_{m}(t))\theta_{L}(A_{m}^{k}(t))dt$

subject to the same initial data exist for all $t$ (Theorem 2.4, page 177 in
Ikeda and WatanabeIkeda and Watanabe (1989)). It follows from the second
equation that

$\left\|A_{n}^{L}(t)\right\|\leq\left\|A_{n}^{L}(0)\right\|+\int_{0}^{t}\left\|A_{n}^{L}(s)\right\|\sum_{m}\left\|DK_{\varphi}\right\|_{\infty}\left\|A_{m}^{L}(s)\right\|ds,$
(20)

by adding the inequalities where $n$ runs through the index set, to obtain
that

$\left\|A^{L}(t)\right\|\leq\left\|A(0)\right\|+\left\|DK_{\varphi}\right\|_{\infty}\int_{0}^{t}\left\|A^{L}(s)\right\|^{2}ds$
(21)

for $t>0$. By Gronwall’s inequality

$\left\|A^{L}(t)\right\|\leq\frac{\left\|A(0)\right\|}{1-\left\|A(0)\right\|\left\|DK_{\varphi}\right\|_{\infty}t}$
(22)

for all $t<\tau_{A(0),\varphi}$. The estimate (19) now follows immediately by
sending $L\uparrow\infty$.

###### Remark 1

From (22), it is clear that

$\lVert A(t)\rVert\leq 2\lVert A(0)\rVert$ (23)

for $t<\frac{1}{2}\tau_{A(0),\varphi}$.

By Gronwall’s inequality, it follows from (20) that

$\left\|A_{n}(t)\right\|\leq\left\|A_{n}(0)\right\|e^{\left\|DK_{\varphi}\right\|_{\infty}\int_{0}^{t}\left\|A(s)\right\|ds}$
(24)

for $t<\tau_{A(0),\varphi}$, and therefore, by (23)

$\left\|A_{n}(t)\right\|\leq\left\|A_{n}(0)\right\|e^{2\left\|DK_{\varphi}\right\|_{\infty}\left\|A(0)\right\|t}$
(25)

for $t<\frac{1}{2}\tau_{A(0),\varphi}$.

In the next section, we demonstrate that $(W,V)$ defined in Theorem 1 is an
approximation solution to the vorticity equation (2) in certain sense.

##  III Approximating the vorticity equation

In this section, we assume that the structure data $K,\varphi$, $x_{n}$ and
$A_{n}(0)$ are given as in the previous section. $(X_{n},A_{n})$ is the unique
maximal solution pair to the SDEs (15) and (16), and $V(x,t)$ and $W(x,t)$ are
defined in Theorem 1.

###### Theorem 2

For each $x\in\mathbb{R}^{3}$, $W(x,t)$ and $V(x,t)$ are continuous semi-
martingales for $t\in[0,\tau_{A(0),\varphi})$, and

$\displaystyle dW^{i}$ $\displaystyle=\left(W^{j}\frac{\partial
V^{i}}{\partial x^{j}}-V^{j}\frac{\partial W^{i}}{\partial x^{j}}+\nu\Delta
W^{i}\right)dt$ $\displaystyle-\sqrt{2\nu}\frac{\partial W^{i}}{\partial
x^{j}}dB^{j}(t)+\left(F^{i}+G^{i}\right)dt$ (26)

for $t\in[0,\tau_{A(0),\varphi}]$, where the error terms

$F^{i}(x,t)=\sum_{n}A_{n}^{i}(t)\left(V^{j}(x,t)-V^{j}(X_{n}(t),t)\right)\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))$ (27)

and

$G^{i}(x,t)=\sum_{n}A_{n}^{j}(t)\left(\frac{\partial V^{i}}{\partial
x^{j}}(X_{n}(t),t)-\frac{\partial V^{i}}{\partial
x^{j}}(x,t)\right)\varphi(x-X_{n}(t)).$ (28)

Moreover, according to (18) and (17), $V=K\star W$ where the convolution is
made with the wedge product of two vectors.

Proof. Since $A_{n}(x,t)$ are continuous semi-martingales with finite
variations (in $t$), by applying Itô’s formula to $W^{i}$ define in (17), we
have

$\displaystyle dW^{i}(x,t)$ $\displaystyle=$
$\displaystyle\sum_{n}\varphi(x-X_{n}(t))dA_{n}^{i}(t)-\sum_{n}A_{n}^{i}(t)\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))dX_{n}^{j}(t)$ $\displaystyle+$
$\displaystyle\nu\sum_{n}A_{n}^{i}(t)\Delta\varphi(x-X_{n}(t))dt.$

Next we substitute $dX_{n}$ and $dA_{n}$ by using the SDEs (15) and (16) to
obtain that

$\displaystyle dW^{i}(x,t)$ $\displaystyle=$
$\displaystyle\sum_{n}\sum_{m}\varepsilon^{iqk}A_{n}^{l}(t)A_{m}^{k}(t)\varphi(x-X_{n}(t))\frac{\partial
K_{\varphi}^{q}}{\partial x^{l}}(X_{n}(t)-X_{m}(t))dt$ $\displaystyle+$
$\displaystyle\sum_{n}\sum_{m}\varepsilon^{jkq}A_{m}^{k}(t)A_{n}^{i}(t)K_{\varphi}^{q}(X_{n}(t)-X_{m}(t))\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))dt$ $\displaystyle+$
$\displaystyle\nu\sum_{n}A_{n}^{i}(t)\Delta\varphi(x-X_{n}(t))dt$
$\displaystyle-$
$\displaystyle\sqrt{2\nu}\sum_{n}A_{n}^{i}(t)\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))dB^{j}(t).$ (29)

On the other hand, according to the construction (17), (18) we have

$\displaystyle V^{j}\frac{\partial W^{i}}{\partial x^{j}}(x,t)$
$\displaystyle=$
$\displaystyle-\sum_{m}\sum_{n}\varepsilon^{jkq}A_{m}^{k}(t)A_{n}^{i}(t)K_{\varphi}^{q}(x-X_{m}(t))\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))$ $\displaystyle=$
$\displaystyle-\sum_{m}\sum_{n}\varepsilon^{jkq}A_{m}^{k}(t)A_{n}^{i}(t)K_{\varphi}^{q}(X_{n}(t)-X_{m}(t))\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))$ $\displaystyle+$
$\displaystyle\sum_{n}A_{n}^{i}(t)\left(V^{j}(x,t)-V^{j}(X_{n}(t),t)\right)\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t)),$ $\displaystyle W^{j}\frac{\partial V^{i}}{\partial
x^{j}}(x,t)$ $\displaystyle=$
$\displaystyle-\sum_{n}\sum_{m}\varepsilon^{ikq}A_{m}^{j}(t)A_{n}^{k}(t)\varphi(x-X_{m}(t))\frac{\partial
K_{\varphi}^{q}}{\partial x^{j}}(x-X_{n}(t))$ $\displaystyle=$
$\displaystyle\sum_{n}\sum_{m}\varepsilon^{iqk}A_{n}^{l}(t)A_{m}^{k}(t)\varphi(x-X_{n}(t))\frac{\partial
K_{\varphi}^{q}}{\partial x^{l}}(x-X_{m}(t))$ $\displaystyle=$
$\displaystyle\sum_{n}\sum_{m}\varepsilon^{iqk}A_{n}^{l}(t)A_{m}^{k}(t)\varphi(x-X_{n}(t))\frac{\partial
K_{\varphi}^{q}}{\partial x^{l}}(X_{n}(t)-X_{m}(t))$ $\displaystyle-$
$\displaystyle\sum_{n}A_{n}^{l}(t)\left(\frac{\partial V^{i}}{\partial
x^{l}}(X_{n}(t),t)-\frac{\partial V^{i}}{\partial
x^{l}}(x,t)\right)\varphi(x-X_{n}(t))$

and

$\Delta W^{i}(x,t)=\sum_{n}A_{n}^{i}(t)\Delta\varphi(x-X_{n}(t)).$

Substituting these equations into (29), we obtain (26).

We next show that $W$ and $V$ are approximation solutions to the vorticity
equation (2) with the initial vorticity $W(x,0)$.

###### Lemma 1

Suppose the support of $\varphi$ is contained in the ball centered at $0$ with
radius $r_{\varphi}>0$.

1) It holds that

$|F(x,t)|\leq
2\left\|A(t)\right\|^{2}\left\|DK_{\varphi}\right\|_{\infty}\left\|D\varphi\right\|_{\infty}r_{\varphi}$

and

$|G(x,t)|\leq
2\left\|A(t)\right\|^{2}\left\|D^{2}K_{\varphi}\right\|_{\infty}\left\|\varphi\right\|_{\infty}r_{\varphi}$

for $t<\tau_{A(0),\varphi}$ and $x\in\mathbb{R}^{3}$.

2) There is a positive constant $C$ depending only on $\left\|A(0)\right\|$,
$\left\|DK_{\varphi}\right\|_{\infty}$,
$\left\|D^{2}K_{\varphi}\right\|_{\infty}$, $\left\|D\varphi\right\|_{\infty}$
and $\left\|\varphi\right\|_{\infty}$ such that

$|F(x,t)|\leq Cr_{\varphi}\quad\textrm{ and }|G(x,t)|\leq Cr_{\varphi}$

for all $t\in\left[0,\frac{1}{2}\tau_{A(0),\varphi}\right]$ and
$x\in\mathbb{R}^{3}$.

Proof. By definition (18)

$V^{i}(x,t)=-\sum_{n}\varepsilon^{ikj}A_{n}^{k}(t)K_{\varphi}^{j}(x-X_{n}(t)),$
(30)

we have

$|V(x,t)-V(y,t)|\leq 2\left\|A(t)\right\|\left\|DK_{\varphi}\right\||x-y|$
(31)

and

$|DV(x,t)-DV(y,t)|\leq
2\left\|A(t)\right\|\left\|D^{2}K_{\varphi}\right\||x-y|.$ (32)

Therefore 1) follows immediately and 2) follows from 1) and (23).

###### Corollary 1

Suppose the support of the wavelet function $\varphi$ is contained in the ball
centered at $0$ with radius $r_{\varphi}>0$. Let

$\displaystyle U^{i}(x,t)$ $\displaystyle=$ $\displaystyle
W^{i}(x,t)-W^{i}(x,0)$ $\displaystyle-$
$\displaystyle\int_{0}^{t}\left(W^{j}\frac{\partial V^{i}}{\partial
x^{j}}V^{j}\frac{\partial W^{i}}{\partial x^{j}}+\nu\Delta
W^{i}\right)(x,s)ds$ (33)

for $i=1,2,3$ and $t<\tau_{A(0),\varphi}$ and $x\in\mathbb{R}^{3}$. Then

$\left|\mathbb{E}\left[U^{i}(x,t)\right]\right|\leq
C_{1}\tau_{A(0),\varphi}r_{\varphi}$ (34)

for all $t\in\left[0,\frac{1}{2}\tau_{A(0),\varphi}\right]$.

Proof. By (26) we have

$U^{i}(x,t)=-\sqrt{2\nu}\int_{0}^{t}\frac{\partial W^{i}}{\partial
x^{j}}dB^{j}(s)+\int_{0}^{t}\left(F^{i}+G^{i}\right)ds.$

Since under the assumptions on $\varphi$, $\int_{0}^{t}\frac{\partial
W^{i}}{\partial x^{j}}dB^{j}(s)$ are martingales for
$t\in\left[0,\frac{1}{2}\tau_{A(0),\varphi}\right]$ so that

$\mathbb{E}\left[U^{i}(x,t)\right]=\mathbb{E}\int_{0}^{t}\left(F^{i}+G^{i}\right)ds$

for $i=1,2,3$ and (34) follows from Lemma 1, 2).

###### Remark 2

If $K=G$ is the Biot-Savart kernel, then $\textrm{div}V=0$ according to the
formula (18) and

$V^{i}(x,t)=\int_{\mathbb{R}^{3}}\varepsilon^{ijk}G^{j}(x-y)W^{k}(y,t)dy$ (35)

which implies that $\Delta V=-\nabla\wedge W$. Therefore $W=\nabla\wedge
V+\nabla f$ for some scalar function $f$, and $V$ and $W$ satisfy the
vorticity equation approximately in the sense stated in Corollary 1. The
equation that $W=\nabla\wedge V$ may fail in general due to the fact that
$W(x,0)$ may be not divergence-free. The reason why the approximation $V$ and
$W$ still do the job (see the simulations below) nicely is that in practice
$\varphi$ is close to an indicator function of a small ball, and therefore
$W(x,t)$ is nearly divergence-free, and the relation $W\simeq\nabla\wedge V$
may be restored approximately.

##  IV Weak convergence

A simple procedure of sampling the initial vortices, for the sake of
theoretical study, may be described as the following. Suppose $\omega_{0}$ is
the initial vorticity, which is smooth with a compact support $S$, to the
initial value problem to the vorticity equations (2). Let $\psi$ be described
at the beginning of the previous section which defines $G_{\delta}$ for every
$\delta>0$.

###### Lemma 2

Let $p\geq 1$, and the error terms in the random vortex scheme, $F(x,t)$ and
$G(x,t)$ be defined by (27) and (28).

1) The following two estimates hold:

$\left\|F(\cdot,t)\right\|_{W^{-1,p}}\leq
6\left\|A(t)\right\|^{2}\left(\left\|DK_{\varphi}\right\|\left\|\varphi\right\|_{L^{p}}r_{\varphi}+\left\|\varphi\nabla\wedge
K_{\varphi}\right\|_{L^{p}}\right)$ (36)

and

$\left\|G(x,t)\right\|_{L^{p}}\leq
2\left\|A(t)\right\|^{2}\left\|D^{2}K_{\varphi}\right\|_{\infty}\left\|\varphi\right\|_{L^{p}}r_{\varphi}$
(37)

for all $t\in[0,\frac{1}{2}\tau_{A(0),\varphi}]$, where
$\left\|\cdot\right\|_{W^{s,p}}$ denote the Sobolev norms, see Adams (1975)

2) If the singular kernel $K=G_{\delta}$, where $\delta>0$, given in (9), then
there are universal constant $C_{1},C_{2}$ such that

$\left\|F(\cdot,t)\right\|_{W^{-1,p}}\leq
C_{2}\frac{r_{\varphi}}{\delta^{3}}e^{C_{1}\frac{2}{\delta^{3}}\left\|A(0)\right\|}\left\|A(0)\right\|^{2}\left\|\varphi\right\|_{L^{p}}$

and

$\left\|G(x,t)\right\|_{L^{p}}\leq
C_{2}\frac{r_{\varphi}}{\delta^{4}}e^{C_{1}\frac{2}{\delta^{3}}\left\|A(0)\right\|}\left\|A(0)\right\|^{2}\left\|\varphi\right\|_{L^{p}}$

for $t\in\left[0,\frac{1}{2}\tau_{A(0),\varphi}\wedge 1\right]$.

Proof. It is clear that from definition

$\left\|D^{k}K_{\varphi}\right\|_{\infty}\leq\left\|D^{k}K\right\|_{\infty}.$

In particular, if $K=G_{\delta}$, then

$\left\|D^{k}K_{\varphi}\right\|_{\infty}\leq C_{k}\frac{1}{\delta^{2+k}}$

for some universal constants $C_{1},C_{2},\ldots$. Consider the linear
functional $f_{i}(h)=\int_{\mathbb{R}^{3}}F^{i}(x,t)h(x)dx$ where $h$ defined
on $\mathbb{R}^{3}$is smooth with a compact support, i.e. a test function.
Since

$F^{i}(x,t)=\sum_{n}A_{n}^{i}(t)\left(V^{j}(x,t)-V^{j}(X_{n}(t),t)\right)\frac{\partial\varphi}{\partial
x^{j}}(x-X_{n}(t))$ (38)

so that

$\displaystyle f_{i}(h)$ $\displaystyle=$
$\displaystyle-\sum_{n}A_{n}^{i}(t)\int_{\mathbb{R}^{3}}\left(V^{j}(x,t)-V^{j}(X_{n}(t),t)\right)\varphi(x-X_{n}(t))\frac{\partial
h(x)}{\partial x^{j}}dx$ $\displaystyle-$
$\displaystyle\int_{\mathbb{R}^{3}}\textrm{div}V(x,t)h(x)\varphi(x-X_{n}(t))\sum_{n}A_{n}^{i}(t)dx.$

By the definition for $V$ (see (18)) one has the elementary estimate (see the
proof of Lemma 1)

$|V^{j}(x,t)-V^{j}(y,t)|\leq
2\left\|A(t)\right\|\left\|DK_{\varphi}\right\||x-y|.$

Since

$\displaystyle\textrm{div}V(x,t)$
$\displaystyle=-\sum_{n}\varepsilon^{ikj}A_{n}^{k}(t)\frac{\partial}{\partial
x^{i}}K_{\varphi}^{j}(x-X_{n}(t))$
$\displaystyle=\sum_{n}A_{n}^{k}(t)\nabla\wedge K_{\varphi}(x-X_{n}(t))$

so that

$|\textrm{div}V(x,t)|\leq\left\|A(t)\right\||\nabla\wedge
K_{\varphi}(x-X_{n}(t))|.$

By using these estimates and the assumption that the support of $\varphi$ lies
in the ball centered at $0$ with radius $r_{\varphi}$, we therefore deduce
that

$\displaystyle\left|f_{i}(h)\right|$ $\displaystyle\leq
2\left\|A(t)\right\|^{2}\left\|DK_{\varphi}\right\|r_{\varphi}\int_{\mathbb{R}^{3}}\varphi(x-X_{n}(t))|Dh(x)|dx$
$\displaystyle+\left\|A(t)\right\|^{2}\int_{\mathbb{R}^{3}}|\nabla\wedge
K_{\varphi}(x-X_{n}(t))||\varphi(x-X_{n}(t))|h(x)|dx$ $\displaystyle\leq
2\left\|A(t)\right\|^{2}\left\|DK_{\varphi}\right\|\left\|\varphi\right\|_{L^{p}}r_{\varphi}\left\|Dh\right\|_{L^{q}}$
$\displaystyle+\left\|A(t)\right\|^{2}\left\|\varphi\nabla\wedge
K_{\varphi}\right\|_{L^{p}}\left\|h\right\|_{L^{q}}$

where $p\geq 1$ and $\frac{1}{p}+\frac{1}{q}=1$, for every test function $h$
and $i=1,2,3$. Therefore

$\left\|F(\cdot,t)\right\|_{W^{-1,p}}\leq
6\left\|A(t)\right\|^{2}\left(\left\|DK_{\varphi}\right\|\left\|\varphi\right\|_{L^{p}}r_{\varphi}+\left\|\varphi\nabla\wedge
K_{\varphi}\right\|_{L^{p}}\right).$

The estimate for the error term $G^{i}$ is trivial. In fact by (28) and (32)
to obtain

$\displaystyle\left\|G^{i}(x,t)\right\|_{L^{p}}$ $\displaystyle\leq$
$\displaystyle\left\|\sum_{n}A_{n}^{j}(t)\left(\frac{\partial V^{i}}{\partial
x^{j}}(X_{n}(t),t)-\frac{\partial V^{i}}{\partial
x^{j}}(x,t)\right)\varphi(x-X_{n}(t))\right\|_{L^{p}}$ $\displaystyle\leq$
$\displaystyle
2\left\|A(t)\right\|^{2}\left\|D^{2}K_{\varphi}\right\|_{\infty}\left\|\varphi\right\|_{L^{p}}r_{\varphi}$

In particular, if $K=G_{\delta}$, then $\nabla\wedge K_{\varphi}=0$, so that

$\displaystyle\left\|F(\cdot,t)\right\|_{W^{-1,p}}$ $\displaystyle\leq
6\left\|A(t)\right\|^{2}\left\|DK_{\varphi}\right\|\left\|\varphi\right\|_{L^{p}}r_{\varphi}$
$\displaystyle\leq
6C_{1}\frac{r_{\varphi}}{\delta^{3}}\left\|A(t)\right\|^{2}\left\|\varphi\right\|_{L^{p}}$

and

$\left\|G(x,t)\right\|_{L^{p}}\leq
6C_{2}\frac{r_{\varphi}}{\delta^{4}}\left\|A(t)\right\|^{2}\left\|\varphi\right\|_{L^{p}}$

The other conclusions follow from the following estimate: if $K=G_{\delta}$,
then $\left\|DK_{\varphi}\right\|\leq C_{1}\frac{1}{\delta^{3}}$, thus by
estimate (25) we obtain that

$\displaystyle\left\|A(t)\right\|$
$\displaystyle\leq\left\|A(0)\right\|^{2}e^{2\left\|DK_{\varphi}\right\|_{\infty}\left\|A(0)\right\|t}$
$\displaystyle\leq\left\|A(0)\right\|^{2}e^{2C_{1}\frac{1}{\delta^{3}}\left\|A(0)\right\|t}$
$\displaystyle\leq\left\|A(0)\right\|^{2}e^{C_{1}\frac{1}{\delta^{3}}\left\|A(0)\right\|}$

for all $t\in\left[0,\frac{1}{2}\tau_{A(0),\varphi}\wedge 1\right]$. We are
now in a position to state a weak convergence theorem. Choose $\delta>0$ so
the solutions $\omega^{\delta}$ and $u^{\delta}$ solving (11) and (12) are
good approximations to the vorticity equations (2). Let $h>0$. Let $x_{n}$
where $n=(n^{1},n^{2},n^{3})\in\mathbb{Z}^{3}$ the center of the lattice
(open) box $B_{n}$ with size $h$ whose lower-left corner is $hn$. Then
$\omega_{0}$ may be approximated by $\sum_{n}\omega_{0}(x_{n})1_{B_{n}}$ in
$L^{p}$ space for $p\geq 1$. Let
$W^{h}(x,0)=\sum_{n}A_{n,h}(0)\phi_{h}(x-x_{n})$ with
$A_{n,h}(0)=\omega_{0}(x_{n})h^{3}$, where $n$ runs over all $n$ such that
$B_{n}\cap S\neq\textrm{\~{A}}$, and
$\phi_{h}(x)=\frac{8}{h^{3}}\phi\left(\frac{2}{h}x\right)$, where $\phi\geq 0$
is chosen so that $\phi$ is an approximation of $1_{B}$ in some Sobolev space,
$\phi$ is smooth with a compact support in $B$, where $B=(-1/2,1/2)^{3}$. Then

$\left\|A_{,h}(0)\right\|=\sum_{n}|\omega_{0}(x_{n})|h^{3}\leq
h|S|+\left\|\omega_{0}\right\|_{L^{1}}\leq|S|+\left\|\omega_{0}\right\|_{L^{1}}$

which is independent of $h\in(0,1)$ and the lattice size. With this
$W^{h}(x,0)$ as the initial sampling distribution, according to Corollary 1,
$V^{\delta,h},W^{\delta,h}$ defined by (18) and (17) with $K=G_{\delta}$ and
$W(x,0)=W^{h}(x,0)$, are approximations in mean of the vorticity equations (2)
with initial vorticity $\omega_{0}$. Notice that

$\displaystyle K_{\phi_{h}}^{j}(x)$
$\displaystyle=\frac{8}{h^{3}}\int_{\mathbb{R}^{3}}K^{j}(y)\phi_{h}\left(\frac{2}{h}(x-y)\right)dy$
$\displaystyle=\frac{8}{h^{3}}\int_{\mathbb{R}^{3}}K^{j}(x-y)\phi_{h}\left(\frac{2}{h}y\right)dy$
$\displaystyle=\int_{|y|_{\infty}<\frac{1}{2}}K^{j}(x-\frac{h}{2}y)\phi\left(y\right)dy$

which yields that

$\left\|D^{k}K_{\phi_{h}}\right\|_{\infty}\leq\left\|D^{k}K\right\|_{\infty}.$

Suppose $K=G_{\delta}$, then

$\left\|D^{k}K_{\phi_{h}}\right\|_{\infty}\leq C_{k}\frac{1}{\delta^{2+k}}.$

###### Theorem 3

Let $\delta>0$ and $h>0$. Let $\psi$, which is used to define $K=G_{\delta}$
for every $\delta>0$, and $\phi$, which is used to define the initial
$W^{h}(x,0)$ for every $h>0$, be two non-negative smooth functions with
supports lying inside the box $(-\frac{1}{2},\frac{1}{2})^{3}$ with
$\int_{\mathbb{R}^{3}}\psi(x)dx=\int_{\mathbb{R}^{3}}\phi(x)dx=1$ as above.
Let

$\displaystyle U^{\delta,h}(x,t)$ $\displaystyle=$ $\displaystyle
W^{\delta,h}(x,t)-W^{h}(x,0)$ $\displaystyle-$
$\displaystyle\int_{0}^{t}\left(W^{\delta,h}\cdot\nabla
V^{\delta,h}-V^{\delta,h}\cdot\nabla W^{\delta,h}+\nu\Delta
W^{\delta,h}\right)ds.$ (39)

Then there is $T$ depending only on $|S|+\left\|\omega_{0}\right\|_{L^{1}}$
such that

$\left\|\mathbb{E}\left[U^{\delta,h}(\cdot,t)\right]\right\|_{W^{-1,1}}\rightarrow
0$

for all $t\in[0,T]$, as $h\downarrow 0$ for any fixed $\delta>0$.

Proof. This follows from the above estimates and Lemma 2 item 2). While from
the proof we can see that the error term tends to zero as $\delta\downarrow 0$
and $h\downarrow 0$ such that $h\ll\delta^{4}e^{-C/\delta^{3}}$ (with some
positive constant $C$). This means that we need to choose the lattice size $h$
much smaller than the regularization $\delta$ in order to ensure the
convergence result. We should point out that this estimate is very crude due
to the lack of a priori estimates for solutions to the Navier-Stokes
equations. The simulations in section 7 demonstrate that the convergence rate
of our vortex method is much fast even for the case where $\delta=0$, but its
proof is beyond the reach of the current mathematical analysis for these
partial differential equations.

##  V Modified random vortex dynamics

In the previous section, we have shown our random vortex system (15, 16)
converges to the vorticity equations in mean, while if the viscosity $\nu>0$
is not small, then the random perturbation term, i.e. the noise part involving
Brownian motion, appearing in (26) may be not small although its mean always
stays zero. To make the noise as small as possible, we may split each particle
into $N$ copies of the same particle and apply independent copies of Brownian
motion for these particles. Hence we propose the following SDEs

$\displaystyle dX_{n,a}^{i}(t)$
$\displaystyle=-\sum_{m,b}\varepsilon^{ikj}A_{m,a}^{k}(t)K_{\varphi}^{j}(X_{n,a}(t)-X_{m,b}(t))dt$
$\displaystyle+\sqrt{2\nu}dB_{n,a}^{i}(t),$ $\displaystyle X_{n,a}(0)$
$\displaystyle=x_{n},$ (40)

and

$\displaystyle
dA_{n,a}^{i}(t)=A_{n,a}^{l}(t)\sum_{m,b}\varepsilon^{ijk}\frac{\partial
K_{\varphi}^{j}}{\partial x^{l}}(X_{n,a}(t)-X_{m,b}(t))A_{m,b}^{k}(t)dt,$
$\displaystyle A_{n,a}(0)=\frac{1}{N}A_{n}(0),$ (41)

where $a,b$ runs from $1$ up to $N$, where $N$ is a fixed natural number, and
$B_{n,a}$ are independent copies of 3D Brownian motion on some probability
space $(\varOmega,\mathcal{F},\mathbb{P})$. Thus we have split the
distribution of the initial vortices as

$\displaystyle W(x,0)$ $\displaystyle=\sum_{n}A_{n}(0)\varphi(x-x_{n})$
$\displaystyle=\sum_{n,a}A_{n,a}(0)\varphi(x-x_{n})$

and we have used independent Brownian motions for different locations to
reduce the eventual noise in the approximation vorticity equation (26).
$W(x,t)$ and $V(x,t)$ are still defined in terms of (18) and (17), where the
only modifications one has to make are the sums over $n,a$ rather than $n$.
The approximation voticity equations (26) has the same form but with different
noise term:

$dW^{i}=\left(W^{j}\frac{\partial V^{i}}{\partial x^{j}}-V^{j}\frac{\partial
W^{i}}{\partial x^{j}}+\nu\Delta
W^{i}\right)dt+H^{i}+\left(F^{i}+G^{i}\right)dt$

where $i=1,2,3$ and

$H(x,t)=-\sqrt{2\nu}\sum_{n,a}A_{n,a}(t)\cdot\nabla\varphi(x-X_{n,a}(t))dB_{n,a}(t)$

which can be verified by using Itô’s formula too. Hence

$\displaystyle\mathbb{E}\left|\int_{0}^{t}dH(x,s)\right|^{2}$ $\displaystyle=$
$\displaystyle
2\nu\frac{1}{N}\sum_{n}\int_{0}^{t}|A_{n}(s)|^{2}|\nabla\varphi(x-X_{n,a}(s))|^{2}ds$
$\displaystyle\leq$
$\displaystyle\frac{4\nu}{N}\left\|\nabla\varphi\right\|_{\infty}^{2}\left\|A(0)\right\|\max_{n}|A_{n}(0)|e^{2\left\|DK_{\varphi}\right\|_{\infty}\left\|A(0)\right\|t}$

for $t\leq\frac{1}{2}\tau_{A(0),\varphi}$, which goes to zero as
$N\rightarrow\infty$. The error term estimates for $F$ and $G$ remain the same
which are independent of $N$. Therefore

$\left\|\mathbb{E}\left[|U^{\delta,h}(\cdot,t)|^{2}\right]\right\|_{W^{-1,1}}\rightarrow
0$

for all $t\in[0,T]$, as $h\downarrow 0$ for any fixed $\delta>0$.

##  VI 2D case

We retain the same notations as used in the previous sections, but we work
with the two dimensional vorticity equation instead. 2D case is special as we
have indicated – there is no non-linear stretching term in the vorticity
equation. The vorticity $\omega(x,t)$ of an incompressible fluid flow with
velocity $u=(u^{1},u^{2})$ can be identified with the scalar function
$\frac{\partial}{\partial x^{1}}u^{2}-\frac{\partial}{\partial x^{2}}u^{1}$,
and the vorticity equation

$\frac{\partial}{\partial t}\omega+u^{j}\frac{\partial}{\partial
x^{j}}\omega=\nu\Delta\omega$ (42)

appears as a “linear” parabolic equation. Hence the initial vortices are
transported along the motion of Brownian fluid particles. Since
$\textrm{div}u=0$, so that

$\Delta u^{1}=-\frac{\partial}{\partial x^{2}}\omega,\quad\Delta
u^{2}=\frac{\partial}{\partial x^{1}}\omega$ (43)

and

$u(x,t)=\int_{\mathbb{R}^{2}}G(x-y)\omega(y,t)dy$

where in 2D case, the Biot-Savart kernel
$G(x)=\frac{1}{2\pi}(-\frac{x_{2}}{|x|^{2}},\frac{x_{1}}{|x|^{2}})$ which is
singular at the original $0$. We therefore need to treat it as for the 3D case
to introduce $G_{\delta}(x)=G\ast\psi_{\delta}(x)$, where
$\psi_{\delta}(x)=\delta^{-2}\psi(\delta^{-1}x)$, and $\psi(x)$ is a smooth
function with total integral $1$ and a compact support lying inside
$(-1/2,1/2)^{2}$. Then we have $\lVert DG_{\delta}(x)\rVert\leq
C\frac{1}{\delta^{2}}$, where $C$ depends only on the regularization $\psi$.
We then implement the vortex method to the approximation system (42) together
with the integral relation

$u(x,t)=\int_{\mathbb{R}^{2}}G_{\delta}(x-y)\omega(y,t)dy$ (44)

for each fixed $\delta>0$. Since $\omega(x,t)$ is scalar, so that there is a
simplification when we sample the distribution of initial vortices. As for the
3D case, $(x_{n})$ is a finite collection of points in $\mathbb{R}^{2}$ and
the initial vortices $A_{n}(0)$ are scalars. The dynamic equation for
$A_{n}(t)$ is no longer needed as $A_{n}(t)=A_{n}(0)=A_{n}$ are independent of
$t$, and the dynamic equation for $X_{n}(t)$ is reduced to the following SDE

$\displaystyle
dX_{n}^{i}(t)=\sum_{m}A_{m}K_{\varphi}^{i}(X_{n}(t)-X_{m}(t)))dt+\sqrt{2\nu}dB^{i}(t),$
$\displaystyle X_{n}(0)=x_{n},$ (45)

whose coefficients are globally Lipschitz, so that it has a unique strong
solution $X_{n}(t)$ defined for all $t$. Let

$W(x,t)=\sum_{n}A_{n}\varphi(x-X_{n}(t))$ (46)

and

$V^{i}(x,t)=\sum_{n}A_{n}K_{\varphi}^{i}(x-X_{n}(t))$ (47)

By using Itô’s formula one may deduce that

$\displaystyle dW(x,t)$ $\displaystyle=$ $\displaystyle-V(x,t)\cdot\nabla
W(x,t)dt+\nu\Delta W(x,t)dt$ $\displaystyle-$
$\displaystyle\sum_{n}A_{n}\nabla\varphi(x-X_{n}(t))\cdot\sqrt{2\nu}dB(t)$
$\displaystyle-$
$\displaystyle\sum_{n}A_{n}\nabla\varphi(x-X_{n}(t))\cdot(V(X_{n}(t),t)-V(x,t))dt.$

The non-martingale error term

$F(x,t)=\sum_{n}A_{n}\nabla\varphi(x-X_{n}(t))\cdot(V(X_{n}(t),t)-V(x,t)).$

Following the same procedure as in the proof of Lemma 2, we may deduce the
following global estimates for 2D case.

###### Lemma 3

1) For any $p\geq 1$

$\lVert F\rVert_{W^{-1,p}}\leq\lVert A\rVert^{2}\bigg{(}\lVert
DK_{\varphi}\rVert\lVert\varphi\rVert_{L^{p}}r_{\varphi}+\lVert\varphi\text{\text{div}}K_{\varphi}\rVert_{L^{p}}\bigg{)}.$

2) If the singular kernel $K=G_{\delta}$, where $\delta>0$, then

$\lVert F(\cdot,t)\rVert_{W^{-1,p}}\leq C_{5}\frac{1}{\delta^{2}}\lVert
A\rVert^{2}\lVert\varphi\rVert_{L^{p}}r_{\varphi}$

for all $t>0$ where $\lVert A\rVert=\sum_{n}|A_{n}|$.

Now we state the weak convergence theorem. Let $h>0$, and let $x_{n}$ where
$n=(n_{1},n_{2})\in\mathbb{Z}^{2}$ the center of the lattice box $B_{n}$ with
size $h$ whose lower left corner is $hn$. Then the initial vorticity
$\omega_{0}$ may be approximated by $\sum_{n}\omega_{0}(x_{n})1_{B_{n}}$ for
example in an $L^{p}$-space. Let
$W^{h}(x,0)=\sum_{n}A_{n,h}\varphi_{h}(x-x_{n})$ with
$A_{n,h}=\omega_{0}(x_{n})h^{2}$, where $n$ runs over all $n$ such that
$B_{n}\cap S\neq\tilde{A}$, and
$\varphi_{h}(x)=\frac{1}{h^{2}}\phi(\frac{x}{h})$, where $\phi$ is an
approximation of $1_{B}$ where $B$ is the unit square $(-1/2,1/2)^{2}$ in
$L^{p}$ space for $p>1$. Then

$\lVert A_{,h}\rVert=\sum_{n}|\omega_{0}(x_{n})|h^{2}\leq
h|S|+\lVert\omega_{0}\rVert_{L^{1}}\leq|S|+\lVert\omega_{0}\rVert_{L^{1}},$

which is independent of $h\in(0,1)$ and the lattice size. With this
$W^{h}(x,0)$ as the initial sampling distribution, let
$V^{\delta,h},W^{\delta,h}$ be defined as in equations (47) and (46) with
$K=G_{\delta}$. Therefore, by Lemma 3, item 2), we have the following theorem:

###### Theorem 4

Let

$\displaystyle U^{\delta,h}(x,t)$ $\displaystyle=W^{\delta,h}(x,t)-W^{h}(x,0)$
$\displaystyle-\int_{0}^{t}\left(-V^{\delta,h}\cdot\nabla
W^{\delta,h}+\nu\Delta W^{\delta,h}\right)ds.$ (48)

Then

$\left\|\mathbb{E}\left[U^{\delta,h}(\cdot,t)\right]\right\|_{W^{-1,1}}\rightarrow
0$

for all $t>0$ as long as $h\ll\delta^{2}$, and for any $T>0$, this convergence
is uniform for all $(x,t)\in\mathbb{R}^{2}\times(0,T]$.

##  VII Simulation results

First we show some figures of the velocity field by our approximation scheme
for 3D flows, which illustrates the velocity field projected at the plane
$z=0$. The flow starts with two points at $(-1,0,0)$ and $(1,1,0)$ with
vorticity $(0,0,10)$ and $(10,10,10)$ respectively. i.e
$x_{1}=(-1,0,0),x_{2}=(1,1,0)$, $,A_{1}(0)=(0,0,10),A_{2}(0)=(10,10,10)$.
Throughout the simulation in 3D flows, We choose the Biot-Savart kernel $G$,
and we let $\varphi$ be $c(|x|^{2}-1)^{2}$, where $c$ is a constant such that
the total integral of $\varphi=1$ (We find a polynomial helps improve the
computational speed of numerical integration in the packages we used in
program). We let $h=0.5$. Figure. 1 shows the case when $\nu=0$.

(a) $t=0$

(b) $t=1$

(c) $t=3$

(d) $t=5$

Figure 1: An inviscid Flow

(a)

(b)

(c)

(d)

Figure 2: Split each particle into $5$ parts

(a)

(b)

(c)

(d)

Figure 3: Split each particle into $50$ parts

When $\nu\neq 0$, we are going to split the particles as in the previous
section, and we will see that more splitting helps minimise the stochastic
error. We keep the choice of the kernel, $\varphi$, $h$, the initial condition
as in the inviscid flow. The only difference here is that we consider the
viscosity $\nu=1$ in this case, which means we have to solve a system of SDEs
with the appropriate noise. We are looking at the projection of velocity field
when $t=0.5$. (Figure. 2 and Figure. 3). When we split each particle into $50$
parts, the result from different randomisation is more consistent.

For 2D flows simulation, we use the 2D Biot-Savart kernel $H$, and
$\varphi=c(|x|^{2}-1)^{2}$, $h=0.5$. we start the flow with randomized
$(x_{n},A_{n})$ $n=1,2,\dots,10$. $x_{n}$ are sampled in $[-1,1]\times[-1,1]$
by uniform distribution independently, and $A_{n}$ are sampled in $[-10,-10]$
by uniform distribution. The particles inside the flow are denoted by
different colors, and then we look at their locations at different time.
(Figure. 4)

(a) $t=0$

(b) $t=0.3$

(c) $t=0.6$

(d) $t=1$

Figure 4: The dynamics of an inviscid flow

Figure. 5 is the vorticity picture for a non-inviscid flow with viscosity
$\nu=1$. We use the same kernel, $\varphi$, $h$ as in the previous simulation.
We let $(x_{1},x_{2},x_{3},x_{4})$ be sampled in $\mathbb{R}^{2}$ by $N(0,I)$,
with vorticity $(10,10,-10,-10)$ respectively. We split each particle into
$50$ parts to minimise the stochastic error. The vorticity value at each point
is represented by its color, and the relations between the color and vorticity
are shown in the bar on the right of each figure. The blue color represents
positive vorticity, and the red color represents negative vorticity.

Figure 5: Vorticity of a non-inviscid flow

## Data Availability Statement

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

###### Acknowledgements.

This publication is based on work partially supported by the EPSRC Centre for
Doctoral Training in Mathematics of Random Systems: Analysis, Modelling and
Simulation (EP/S023925/1)

## References

  * Cottet and Koumoutsakos (2000) G. H. Cottet and P. D. Koumoutsakos,  _Vortex methods: theory and practice_ (Cambridge University Press, 2000).
  * Majda and Bertozzi (2002) A. J. Majda and A. L. Bertozzi,  _Vorticity and incompressible flow_ (Cambridge University Press, 2002).
  * Saffman (1995) P. G. Saffman,  _Vortex dynamics_ (Cambridge University Press, 1995).
  * Ting and Knio (2007) R. Ting, L.; Klein and O. M. Knio,  _Vortex dominated flows: analysis and computation for multiple scale phenomena_ (Springer-verlag, 2007).
  * Mimeau and Mortazavi (2021) C. Mimeau and I. Mortazavi, “A review of vortex methods and their applications: from creation to recent advances,” Fluids 6, 68 (2021).
  * Chorin (1973) A. J. Chorin, “Numerical study of slightly viscous flow,” J. Fluid Mech 57, 785–796 (1973).
  * Anderson and Greengard (1985) C. Anderson and C. Greengard, “On vortex methods,” SIAM J. Numer. Anal. 22, 413–440 (1985).
  * Beale and Majda (1982) J. T. Beale and A. Majda, “Vortex methods. i. convergence in three dimensions,” Mathematics of Computation 39, 1–27 (1982).
  * Hald and del Prete (1978) O. Hald and V. M. del Prete, “Convergence of vortex methods for euler’s equations,” Mathematics of Computation 32, 791–809 (1978).
  * Hald (1979) O. Hald, “Convergence of vortex methods for euler’s equations, ii,” SIAM Journal on Numerical Analysis 16, 726–755 (1979).
  * Hald (1987) O. Hald, “Convergence of vortex methods for euler’s equations, iii,” SIAM Journal on Numerical Analysis 24, 538–582 (1987).
  * Beale (1986) J. T. Beale, “A convergent 3-d vortex method with grid-free stretching,” Mathematics of Computation 46, 401–424 (1986).
  * Long (1988) D. G. Long, “Convergence of the random vortex method in two dimensions.” Journal of the American Mathematical Society 1, 779–804 (1988).
  * Pullin and Saffman (1998) D. I. Pullin and P. G. Saffman, “Vortex dynamics in turbulence,” Annu. Rev. Fluid Mech. 30, 31–51 (1998).
  * Ikeda and Watanabe (1989) N. Ikeda and S. Watanabe,  _Stochastic differential equations and diffusion processes_ , 2nd ed. (North-Holland Publishing Company, 1989).
  * Adams (1975) R. A. Adams,  _Sobolev spaces_ (Academic Press, 1975).


so that $\tilde{D}_{i,0}:=D_{i,0}+\frac{1}{p}f\in\hat{\mathbf{N}}_{p}$. Also
recall that $D_{i,0}=-\sum_{j=1}^{p-1}D_{i,j}$ and
$\sum_{j=0}^{p-1}\tilde{D}_{i,j}=f$. Extend the $C_{p}$-action to
$\mathbf{N}_{p}\oplus\tfrac{1}{p}\mathbf{Z}$ by letting it act trivially on
$\tfrac{1}{p}\mathbf{Z}$. Since $\varrho$ is an element of
$p\mathbf{L}_{p}^{\vee}$, and $\sigma-1$ takes $\mathbf{L}_{p}^{\vee}$ (and
hence $\mathbf{N}_{p}$) to $\mathbf{L}_{p}$, it follows that $i(\sigma-1)$
takes its values in $\mathbf{Z}$. In other words, $\hat{\mathbf{N}}_{p}$ is
$C_{p}$-invariant.

Suppose we are given primitive embedding
$\hat{\mathbf{N}}_{p}\hookrightarrow\mathbf{\Lambda}$. If $M$ is endowed with
a complex structure $J$ for which the image of this embedding is the Picard
lattice (i.e., the orthogonal complement of $H^{2,0}(M,J)$), then the theory
of K3 surfaces ensures that after possibly multiplying this embedding with
$-1$ (in order to make it compatible with the spinor orientation), $f$ is the
class of an elliptic fibration. This elliptic fibration has for each
$i=1,\dots,\nu$ a singular fiber of type $I_{p}$ whose irreducible components
have as classes the images of $\\{\tilde{D}_{i,j}\\}_{j=0}^{p-1}$, and $\nu$
other singular fibers of type $I_{1}$. The Torelli theorem ensures that
$C_{p}$ uniquely lifts to the automorphism group of $(M,J)$ as in Figure (1)
so that $\mathbf{L}_{C_{p}}$ is then just the image of $\mathbf{L}_{p}$.

######  Remark 4.15.

If the singularities of the quotient $X/G$ are resolved in a minimal manner
(as in § 3.5) we obtain a K3 surface $X^{\prime}$ that is also elliptically
fibered, but one for which the $G$-orbit space of each $I_{p}$-fiber resp.
$I_{1}$-fibers of $X$ yields an $I_{1}$-fiber resp. resolves to an
$I_{p}$-fiber of $X^{\prime}$. So $X^{\prime}$ is of the same type as $X$. In
each smooth fiber $F$ the group $G$ is contained in its translation group
$T(F)$ (an elliptic curve for which $F$ is a principal homogeneous space); it
is in fact a subgroup of order $p$ of the $p$-torsion group $T(F)[p]$. The
$G$-orbit space $F^{\prime}$ of $F$ has translation group $T(F)/G$ and the
image of $T(F)[p]$ in $T(F)/G$ is a group that is via the Weil pairing (or
equivalently, the symplectic form on $H_{1}(F)$) naturally identified with
$G^{\prime}:=\operatorname{Hom}(G,\mathbb{C}^{\times})$. One can show that
$G^{\prime}$ in fact acts on $X^{\prime}$ and that repeating this construction
with $(X^{\prime},G^{\prime})$ returns a copy of $(X,G)$.

Elliptic fibrations with $C_{p}$-action permuting the fibers. The advantage of
the above model is that it helps us to understand the geometry of the
$C_{p}$-action. But a K3 surface $X$ with an elliptic fibration and an
automorphism of order $p$ is not necessarily of the above type, as we will see
below this can happen when $p=2,3,5$ (the algebraic counterpart of this
observation is somewhat implicit in §5 of Garbagnati-Sarti [GS]), but not for
$p=7$. We here use that for such $p$ there exists an embedding
$i:\mathbf{L}_{p}\hookrightarrow\mathbf{\Lambda}$ such that
$i(\mathbf{L}_{p})^{\perp}$ contains a copy of $\mathbf{U}$ (as is shown in
[GS]). Let $f$ be one of the isotropic generators of this copy so that
$i(\mathbf{L}_{p})\oplus\mathbb{Z}f$ is a primitive sublattice of
$\mathbf{\Lambda}$. The Torelli theorem guarantees the existence of a complex
structure $J$ on $M$ such that $X:=(M,J)$ has
$i(\mathbf{L}_{p})\oplus\mathbb{Z}f$ as its Picard lattice. This K3 surface
comes with a faithful action of $C_{p}$ and $\pm f$ is the class of an
elliptic fibration $X\to B\cong\mathbb{P}^{1}$ that will be left invariant
under the $C_{p}$-action. We want to understand how $C_{p}$ acts on this
fibration.

Since $X$ has no $(-2)$-curves, the Kodaira classification of singular fibers
implies that these must be all irreducible: either nodal or cuspidal. So if a
fiber $X_{t}$ is singular and $C_{p}$-invariant, then the number of fixed
points of $C_{p}$ it contains is equal to its Euler characteristic: one when
$X_{t}$ is nodal (the fixed point is its singular point) and two when $X_{t}$
is cuspidal (the singular point and a smooth point). This implies that $C_{p}$
cannot preserve each fiber, for then the number of fixed points would at least
be equal to the Euler characteristic of $X$, which is $24$ and this is $>\nu$.
So $C_{p}$ acts nontrivially on the base $B\cong\mathbb{P}^{1}$ of the pencil.
In terms of an affine coordinate $t$ on $B$ this action will be as
multiplication with a primitive $p$th root of unity, so that the fibers
$X_{0}$ (over $t=0$) and $X_{\infty}$ (over $t=\infty$), contain the fixed
points. Each of these comes with an automorphism of order $p$ and the number
of the fixed points must add up to $\nu$. If $X_{0}$ (or $X_{\infty}$) is
smooth, then for $p=2$ it can have four fixed points (minus the identity with
respect to a fixed point) and for $p=3$ it can have three (but only when
$X_{0}$ has $j$-invariant zero). Singular fibers different from $0,\infty$ lie
over free $C_{p}$-orbits in $\mathbb{P}^{1}\smallsetminus\\{0,\infty\\}$. This
then leads to the following possibilities only.

For $(p,\nu)=(2,8)$, $X_{0}$ and $X_{\infty}$ must be _smooth_ and each must
contain four fixed points, the involution being the standard one described
above. So there will be in general $12$ free $C_{2}$-orbits of singular
fibers. So we will have in general 24 nodal curves. The 24 points in
$\mathbb{P}^{1}\smallsetminus\\{0,\infty\\}$ will consist of $12$ antipodal
pairs.

For $(p,\nu)=(3,6)$, $X_{0}$ and $X_{\infty}$ must be smooth and each contains
three fixed points (so they have $j$-invariant zero) and there will be in
general eight $C_{3}$-orbits of singular fibers.

For $(p,\nu)=(5,4)$, $X_{0}$ and $X_{\infty}$ must each be cuspidal and
contain two fixed points and there will be in general four $C_{5}$-orbits of
singular fibers.

This type of reasoning also shows that the case $(p,\nu)=(7,3)$ cannot occur.
For then the three fixed points are contained in the union of the invariant
fibers $X_{0}$ and $X_{\infty}$, so that one them (say $X_{0}$) must be
cuspidal (the other is then nodal). At the cusp, we then choose complex-
analytic coordinates $(z_{1},z_{2})$ so that $t=z_{1}^{3}+z_{2}^{2}$ and the
$C_{p}$-action is linearized. It is then clear that $C_{p}$ must act on the
tangent space with characters $(\chi^{2},\chi^{3})$ with $\chi$ primitive.
This however contradicts the fact that it acts through
$\textrm{SL}_{2}(\mathbb{C})$. The impossibility of this case answers a
question raised by Garbagnati-Sarti at the end of their paper [GS].

The lattice of an elliptic fibration with a section. We only discuss the more
interesting case, when $C_{p}$ preserves the fibers. Elliptic fibrations
endowed with a section are obtained by enlarging $\hat{\mathbf{N}}_{p}$ to a
lattice $\mathbf{K}_{p}=\hat{\mathbf{N}}_{p}+\mathbb{Z}s$ whose quadratic form
is determined by stipulating that $s$ is perpendicular to $i(\mathbf{N}_{p})$,
$s\cdot f=1$ and $s\cdot s=-2$. It is then clear that the orthogonal
complement of $i(\mathbf{N}_{p})$ in $\mathbf{K}_{p}$ is spanned by $f$ and
$s$. The basis $\\{(e:=s+f,f)\\}$ identifies this orthogonal complement with
$\mathbf{U}$, so that

$\mathbf{K}_{p}=i(\mathbf{N}_{p})\oplus\mathbf{U}\cong\mathbf{N}_{p}\oplus\mathbf{U}.$

Let us also determine the orthogonal complement of $\mathbf{L}_{p}$ in
$\mathbf{K}_{p}$. This clearly contains $f$. Since for $x\in\mathbf{N}_{p}$,

$0=s\cdot i(x)=s\cdot(x+\tfrac{1}{p}(\varrho\cdot x)f)=s\cdot
x+\tfrac{1}{p}\varrho\cdot x=(s+\tfrac{1}{p}\varrho)\cdot x,$

we see that $ps+\varrho\perp\mathbf{L}_{p}$. This also shows that
$s\cdot\varrho=-\tfrac{1}{p}\varrho\cdot\varrho=2(p-1)$ and so

$(ps+\varrho)\cdot(ps+\varrho)=-2p^{2}+4p(p-1)-2(p-1)p=-2p.$

Since $s\cdot f=1$, it follows that $e^{\prime}:=ps+\varrho+f$ is isotropic
and that $\mathbf{L}_{p}^{\perp}=\mathbb{Z}e^{\prime}+\mathbb{Z}f$ is a copy
of $\mathbf{U}(p)$. It also follows that the $C_{p}$-action on
$\mathbf{L}_{p}$ extends to $\mathbf{K}_{p}$ by letting it fix $e^{\prime}$
and $f$, for then $(\sigma-1)(s)=(\sigma-1)(-\varrho/p)$ and we know that
$\sigma-1$ takes $\varrho/p\in\mathbf{L}_{p}^{\vee}$ to $\mathbf{L}_{p}$.

The algebro-geometric interpretation is that for a primitive embedding of
$\mathbf{K}_{p}\hookrightarrow\mathbf{\Lambda}$ that is compatible with the
spinor orientation there exist complex structures $J$ on $M$ for which the
image of $\mathbf{K}_{p}$ is the Picard lattice, $f$ the class of an elliptic
fibration, and $s$ the class of a section. The group $C_{p}$ then lifts
uniquely to $(M,J)$ so that every $C_{p}$-translate of this section is also
one and the image of $\mathbf{L}_{p}$ is $\mathbf{L}_{C_{p}}$.

## References

  * [AB] M. Atiyah and R. Bott: A Lefschetz fixed point formula for elliptic complexes: II, Annals of Math. (2), 88 (1968), 451–491.
  * [BK] D. Baraglia, H. Konno: A note on the Nielsen realization problem for K3 surfaces, https://arxiv.org/pdf/1908.03970.pdf
  * [BPHV] W. Barth, K. Hulek, C. Peters, A. Van de Ven: Compact complex surfaces (2nd ed.) Ergebnisse der Mathematik und ihrer Grenzgebiete (3), 4 Springer-Verlag, Berlin (2004).
  * [Bo] C. Borcea: Diffeomorphisms of a K3 surface, Math. Ann. 275 (1986), 1–4.
  * [Bry] J. Bryan: Seiberg-Witten theory and $\mathbb{Z}/2p$-actions on spin $4$-manifolds, Math. Res. Lett. 5 (1998), 165–183.
  * [Car] R.W. Carter: Conjugacy classes in the Weyl group, Compositio Math. 25 (1972), 1–59.
  * [GKT] J. Giansiracusa, A. Kupers, B. Tshishiku: Characteristic classes of bundles of K3 manifolds and the Nielsen realization problem, Tunesian J. of Math. 3 (2021), 75–92.
  * [CK] W. Chen, S. Kwasik: Symplectic symmetries of 4-manifolds, Topology 46 (2007), 103–128.
  * [D-or] S.K. Donaldson: The orientation of Yang-Mills moduli spaces and 4-manifold topology, J. Differential Geom. 26 (1987), 397–428.
  * [Ed] A. Edmonds: Aspects of group actions on four-manifolds, Topology and its Applications 31 (1989), 109–124.
  * [Fr] M.H. Freedman: The topology of four-dimensional manifolds, J. Differential Geometry 17 (1982), 357–453.
  * [FM] R. Friedman, J.W. Morgan: Smooth four-manifolds and complex surfaces, Ergebnisse der Mathematik und ihrer Grenzgebiete (3), 27, Springer-Verlag, Berlin (1994).
  * [GS] A. Garbagnati, A. Sarti: Symplectic automorphisms of prime order on K3 surfaces, J. of Algebra 318 (2007), 323-350.
  * [vGS] B. van Geemen, A. Sarti: Nikulin involutions on K3 surfaces, Math. Z. 255 (2007), 731–753.
  * [Ha] K. Hashimoto, Finite symplectic actions on the K3 lattice, Nagoya Math. J. 206 (2012), 99–153.
  * [HZ] F. Hirzebruch, D. Zagier: The Atiyah-Singer Theorem and Elementary Number Theory, Math. Lect. Series 3, Publish or Perish, Inc. (1974).
  * [HKLR] Hitchin, N.; Karlhede, A.; Lindström, U.; Roček, M.: Hyper-Kähler metrics and supersymmetry, Comm. Math. Phys. 108 (1987), 535–589.
  * [I] M. Itoh: Conformal geometry of Ricci flat 4-manifolds, Kodai Math. J. 17 (1994), 179–200.
  * [Ke] S. Kerckhoff: The Nielsen realization problem, Annals of Math. (2) 117 (1983), 235–265.
  * [KM] P.B. Kronheimer, T.S. Mrowka: Gauge theory for embedded surfaces. II, Topology 34 (1995), 37–97.
  * [Li] T-J. Li: Smoothly Embedded Spheres in Symplectic 4-Manifolds, Proc. Amer. Math. Soc. 127 (1999), 609–613.
  * [Lo] E. Looijenga: Teichmüller spaces and Torelli theorems for hyperkähler manifolds, Math. Zeit. 298 (2021), 261–279.
  * [Ma1] T. Matumoto: On diffeomorphisms of a K3 surface, in Algebraic and topological theories (Kinosaki, 1984), 616–621, Kinokuniya, Tokyo (1986).
  * [Ma2] T. Matumoto, Homologically trivial smooth involutions on K3 surfaces, in Aspects of low-dimensional manifolds, Adv. Stud. Pure Math., 20, Kinokuniya, 365–376 (1992).
  * [MS] J.W. Milnor, J.D. Stasheff: Characteristic classes, Annals of Mathematics Studies 76, Princeton University Press, Princeton, N. J (1974).
  * [Mor] D.R. Morrison, On K3 surfaces with large Picard number, Invent. Math. 75 (1984), 105–121.
  * [Mu] S. Mukai: Finite groups of automorphisms of K3 surfaces and the Mathieu group, Invent. Math. 94 (1988),183–221.
  * [Nik1] V.V. Nikulin: Finite groups of automorphisms of Kählerian K3 surfaces, Trudy Moskov. Mat. Obshch. 38 (1979), 75–137, translated as: _Finite automorphism groups of Kähler K3 surfaces,_ Trans. Moscow Math. Soc. 38 (1980) 71–135.
  * [Nik2] V.V. Nikulin: Integer symmetric bilinear forms and some of their geometric applications, Izv. Akad. Nauk SSSR Ser. Mat. 43 (1979), 111–177.
  * [Ru] D. Ruberman: Involutions on spin $4$-manifolds, Proc. Amer. Math. Soc. 123 (1995), 593–596.
  * [Se] P. Seidel: Lectures on four-dimensional Dehn twists, In: Symplectic 4-manifolds and algebraic surfaces, 231–267, Lecture Notes in Math. 1938, Springer, Berlin (2008).
  * [Serre] J.-P. Serre: Cours d’arithmétique, Collection SUP: Le Mathématicien 2, Presses Universitaires de France, Paris (1970).
  * [Siu] Y.-T. Siu: Every K3 surface is Kähler, Invent. Math. 73 (1983), 139–150.
  * [Yau] S.-T. Yau: On the Ricci curvature of a compact Kähler manifold and the complex Monge-Ampère equation. I, Comm. Pure Appl. Math. 31 (1978), 339–411.

Dept. of Mathematics, University of Chicago
E-mail<EMAIL_ADDRESS>

Yau Mathematical Sciences Center, Tsinghua U. (Beijing), U. of Chicago, U.
Utrecht
E-mail<EMAIL_ADDRESS>


An analogue of the ElGamal scheme based on the Markovski algorithm

Nadeghda Malyutina, Victor Shcherbacov

###### Abstract

We give an analogue of the ElGamal encryption system based on the Markovski
algorithm [4; 5].

AMS Classification 20N05, 05B15, 94A60

Keywords: quasigroup, ElGamal’s scheme, Markovski algorithm, encryption,
decryption

In memory of academician M. M. Choban

##  1 Introduction

Usually the classical Taher ElGamal encryption system is formulated in the
language of number theory using multiplication modulo a prime [1].

ElGamal’s scheme is a public key cryptosystem based on the difficulty of
computing discrete logarithms in a finite field. The cryptosystem includes an
encryption algorithm and a digital signature algorithm. ElGamal Scheme
Underlies US Former Electronic Digital Signature Standards (DSA) and Russia
(GOST R 34.10-94).

The scheme was proposed by Taher ElGamal in 1985. ElGamal developed one of the
variants of the Diffie-Hellman algorithm. He improved the Diffie-Hellman
system and obtained two algorithms that were used for encryption and for
authentication. Unlike RSA, the ElGamal algorithm was not patented and,
therefore, became a cheaper alternative, since it did not require payment of
license fees [2].

The sender of messages and their recipient can be individuals, organizations,
or technical systems. These may be subscribers of a network, users of a
computer system, or abstract “parties” involved in information interaction.
But more often participants are identified with people and replaced with the
formal designations A and B by Alice and Bob. It is assumed that messages are
transmitted through the so-called “open” communication channel, available for
listening to some other persons.

In cryptography, it is usually assumed that a person sending messages or
receiving them has some opponent E and this opponent can intercept messages
transmitted over an open channel. The enemy is considered as a certain person
named Eve, who has at her disposal powerful computing equipment and owns
cryptanalysis methods. Naturally, Alice and Bob want their messages to be
incomprehensible to Eve, and use special ciphers for this.

Before sending a message over an open communication channel from A to B, A
encrypts the message, and B, having received the encrypted message, decrypts
it, restoring the original text. The important thing is that Alice and Bob can
agree on the cipher they use not on an open channel, but on a special ”closed”
channel, inaccessible for listening to the enemy. It should be borne in mind
that usually the organization of such a closed channel and the transmission of
messages through it is too expensive compared to an open channel or a closed
channel cannot be used at any time. Each attempt to break the cipher is called
an attack on the cipher. In cryptography, it is generally accepted that the
adversary can know the encryption algorithm used, the nature of the
transmitted messages and the intercepted ciphertext, but does not know the
secret key.

Developers of modern cryptosystems strive to make attacks on known and
selected text invulnerable. Significant successes have been achieved along
this path.

##  2 ElGamal’s scheme

Suppose there are subscribers $A,B,C,\dots$ who want to transmit encrypted
messages to each other without having any secure communication channels. We
will consider the code proposed by ElGamal, which solves this problem, using,
in contrast to the Shamir code, only one message forwarding. In fact, the
Diffie -Hellman scheme is used here to form a common secret key for two
subscribers transmitting a message to each other, and then the message is
encrypted by multiplying it by this key. For each subsequent message, the
secret key is recalculated. A large prime number is selected $p$ and number
$g$, such that different degrees of $g$ are different numbers modulo p. The
numbers p and g are transmitted to subscribers in the clear.

Then each subscriber of the group selects his secret number $c_{i}$,
$1<c_{i}<p-1$, and calculates the corresponding open number $d_{i}$,

$d_{i}=g^{c_{i}}\pmod{p}$ (1)

Table 1. (User keys in the ElGamal system)

$\begin{array}[]{|l|l|l|}\hline\cr$\text{S}ubscriber $&$\text{P}rivate \, key
$&$\text{P}ublic \> key $\\\ \hline\cr A&c_{A}&d_{A}\\\ \hline\cr
B&c_{B}&d_{B}\\\ \hline\cr C&c_{C}&d_{C}\\\ \hline\cr\end{array}$

We show now how A sends message m to the subscriber B. We will assume, that
the message is presented as a number $m<p$.

Step 1. A forms a random number $k,1\geq k\geq p-2$, calculates numbers:

$r=g^{k}\mod p$ (2)

$e=m\cdot{d_{B}}^{k}\mod p$ (3)

and passes a couple of numbers $(r,e)$ to the subscriber $B$.

Step 2 B, getting $(r,e)$, calculates

$m^{\prime}=e\cdot r^{p-1-c_{B}}\mod p$ (4)

Statement 1 (properties of the ElGamal cipher).

(1) The subscriber B received a message, i.e. $m^{\prime}=m$;

(2) the adversary, knowing $p,g,d_{B},r$ and $e$, cannot calculate $m$.

######  Example 2.1.

Consider the transmission of message m = 15 from A to B.

We take p = 23, g = 5. Let subscriber B choose for himself a secret number
$c_{B}=13$ and calculate (1): $d_{B}=5^{13}\pmod{23}=21$.

Subscriber A randomly selects the number $k$, for example, $k=7$, and
calculates from (2), (3): $r=5^{7}\pmod{23}=17$, $e=15\cdot
21^{7}\pmod{23}=15\cdot 10\pmod{23}=12.$

Now $A$ sends to $B$ an encrypted message in the form of a pair of numbers
$(17,12)$ and $B$ calculates: $m^{\prime}=12\cdot
17^{\,23-1-13}\pmod{23}=12\cdot 17^{9}\pmod{23}=12\cdot 7\pmod{23}=15.$ So $B$
was able to decrypt the transmitted message.

By a similar scheme, all subscribers in the network can send messages.
Moreover, any subscriber who knows the public key of subscriber B can send him
messages encrypted using the public key $d_{B}$. But only subscriber B, and no
one else, can decrypt these messages using the secret key $c_{B}$ known only
to him.

The Shamir cipher completely solves the problem of exchanging messages that
are closed for reading, in the case when subscribers can use only open
communication lines.

However, this message is sent three times from one subscriber to another,
which is a drawback. The ElGamal cipher allows you to solve the same problem
in one data transfer, but the amount of transmitted ciphertext is twice the
size of the message.

It is easy to see that this system can also be formulated in terms of a
residue ring modulo $p$ or, equivalently, using the language of the Galois
field $GF(p)$.

In addition, we can use the concept of the action of a group of automorphisms
of a cyclic group $(Z_{p},+)$ on this group. Let $(Z_{p},+)$ be a cyclic group
of residues of large simple order with respect to addition of residues and
element $a$ be the generator of the group $(Z_{(p-1)},\cdot)\cong
Aut(Z_{p},+)(gcd(a,p-1)=1)$.

Alice’s keys are the following: Public key $p$, $a$ and $a^{m}$,
$m\in\mathbb{N}$. Private key $m$.

Encryption. To send a message $b\in(Z_{(p-1)},\cdot)$, Bob is calculating
$a^{r}$ and $a^{mr}$ for random $r\in\mathbb{N}$ (sometimes the number $r$ is
called an ephemeral key [3]).

Ciphertext: $(a^{r};a^{mr}\cdot b)$.

Decryption. Alice knows $m$, so if she gets the ciphertext $(a^{r};a^{mr}\cdot
b)$, she will calculate $a^{mr}$ from $a^{r}$ and then $a^{(-mr)}$ and then
from $a^{mr}\cdot b$ calculate b.

######  Example 2.2.

Example Alice picks $p=107,a=2,m=67,$ and calculates $a^{m}=2^{67}\equiv
94\pmod{107}$.

Her public key $(p,a^{m})=(107,94)$, and her private key is $m=67$. Bob wants
to send a message “B ”to Alice. He selects a random integer $r=45$ and
encrypts $B=66$ like $(a^{m})^{r}\cdot B$.

Bob gets: $(2^{45},94^{45}\cdot 66)\equiv(28,5\cdot
66)\equiv(28,9)\pmod{107}$.

He sends an encrypted message $(28,9)$ to Alice. Alice receives this message
and using her private key $m=67$ she decrypts as follows: $28^{(-67)}\cdot
9=28^{(106-67)}\cdot 9\equiv 28^{39}\cdot 9\equiv 43\cdot 9\equiv
66\pmod{107}.$

The complexity of this system is based on the complexity of the discrete
logarithm problem. ElGamal’s encryption system is not secure according to the
selected attack ciphertext [3]. ElGamal cryptosystems are usually used in a
hybrid cryptosystem, i.e. the message itself is encrypted using a symmetric
cryptosystem and ElGamal also uses a symmetric cryptosystem to encrypt the
key.

##  3 An analogue of the ElGamal scheme based on the Markovski algorithm

We give an analogue of the ElGamal encryption system based on the Markovski
algorithm [4; 5].

Let $(Q,f)$ be a binary quasigroup and $T=(\alpha,\beta,\gamma)$ its isotopy.

Alice’s keys are as follows: The public key is $(Q,f)$, $T$,
$T^{(m,n,k)}=(\alpha^{m},\beta^{n},\gamma^{k})$ , $m,n,k\in\mathbb{N}$, and
the Markovski algorithm.

Private key $m,n,k$.

Encryption. To send a message $b\in(Q,f)$, Bob calculated $T^{(r,s,t)}$,
$T^{(mr,ns,kt)}$ for random $r,s,t\in\mathbb{N}$ and $(T^{(mr,ns,kt)}(Q,f))$.

The ciphertext is $(T^{(r,s,t)},T^{(mr,ns,kt)}(Q,f)b).$

To obtain $(T^{(mr,ns,kt)}(Q,f)b)$, Bob uses the Markovski algorithm which is
known to Alice.

Decryption Alice knows $m,n,k,$ so if she gets the ciphertext
$(T^{(r,s,t)},(T^{(mr,ns,kt)}(Q,f))b)$, she will calculate
$(T^{(mr,ns,kt)}(Q,f))^{(-1)}$ using $T^{(r,s,t)}$ and finally she will
calculate $b$.

######  Example 3.1.

Let $(Q,f)$ be a binary quasigroup defined by the following Cayley table:

Table 2.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr f&0&1&2&3&4&5&6\\\ \hline\cr
0&5&2&6&4&0&3&1\\\ \hline\cr 1&1&6&5&3&4&2&0\\\ \hline\cr 2&0&5&4&6&3&1&2\\\
\hline\cr 3&4&1&3&0&2&6&5\\\ \hline\cr 4&2&4&0&1&6&5&3\\\ \hline\cr
5&6&3&1&2&5&0&4\\\ \hline\cr 6&3&0&2&5&1&4&6\\\ \hline\cr\end{array}$

and $T=(\alpha,\beta,\gamma)$ its isotopy, where: $\alpha=(234)(0516)$
corresponds to a permutation of rows of a quasigroup table $Q$;
$\beta=(0321)(56)$ corresponds to a permutation of the columns of a quasigroup
table $Q$ obtained after application $\alpha$; $\gamma=(1236054)$ substitution
applied to the table obtained after application $\beta$.

And for $\gamma$ we have the inverse $\gamma^{(-1)}$ the following kind:
$\gamma^{(-1)}=(1450632)$.

Cayley tables of these permutations are of the form:

Table 3.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\alpha&0&1&2&3&4&5&6\\\ \hline\cr
0&6&3&1&2&5&0&4\\\ 1&3&0&2&5&1&4&6\\\ 2&4&1&3&0&2&6&5\\\ 3&2&4&0&1&6&5&3\\\
4&0&5&4&6&3&1&2\\\ 5&1&6&5&3&4&2&0\\\ 6&5&2&6&4&0&3&1\\\ \hline\cr\end{array}$

Table 4.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\beta&0&1&2&3&4&5&6\\\ \hline\cr
0&2&6&3&1&5&4&0\\\ 1&5&3&0&2&1&6&4\\\ 2&0&4&1&3&2&5&6\\\ 3&1&2&4&0&6&3&5\\\
4&6&0&5&4&3&2&1\\\ 5&3&1&6&5&4&0&2\\\ 6&4&5&2&6&0&1&3\\\ \hline\cr\end{array}$

Table 5.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\gamma^{-1}&0&1&2&3&4&5&6\\\
\hline\cr 0&1&3&2&4&0&5&6\\\ 1&0&2&6&1&4&3&5\\\ 2&6&5&4&2&1&0&3\\\
3&4&1&5&6&3&2&0\\\ 4&3&6&0&5&2&1&4\\\ 5&2&4&3&0&5&6&1\\\ 6&5&0&1&3&6&4&2\\\
\hline\cr\end{array}$

Then Alice’s keys are as follows: The private key: $m=3,n=6,k=5$. The public
key is $(Q,f),T,T^{(3,6,5)}=(\alpha^{3},\beta^{6},\gamma^{5})$ and the
Markovski algorithm, where:
$\alpha^{3}=(0615);\beta^{6}=(02)(13);\gamma^{5}=(0315624),\gamma^{-5}=(0426513)$.

As a result, we get the following Cayley tables:

Table 6.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\alpha^{3}&0&1&2&3&4&5&6\\\
\hline\cr 0&3&0&2&5&1&4&6\\\ 1&6&3&1&2&5&0&4\\\ 2&0&5&4&6&3&1&2\\\
3&4&1&3&0&2&6&5\\\ 4&2&4&0&1&6&5&3\\\ 5&5&2&6&4&0&3&1\\\ 6&1&6&5&3&4&2&0\\\
\hline\cr\end{array}$

Table 7.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\beta^{\,6}&0&1&2&3&4&5&6\\\
\hline\cr 0&2&5&3&0&1&4&6\\\ 1&1&2&6&3&5&0&4\\\ 2&4&6&0&5&3&1&2\\\
3&3&0&4&1&2&6&5\\\ 4&0&1&2&4&6&5&3\\\ 5&6&4&5&2&0&3&1\\\ 6&5&3&1&6&4&2&0\\\
\hline\cr\end{array}$

Table 8.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr(\gamma^{5})^{(-1)}&0&1&2&3&4&5&6\\\
\hline\cr 0&6&1&0&4&3&2&5\\\ 1&3&6&5&0&1&4&2\\\ 2&2&5&4&1&0&3&6\\\
3&0&4&2&3&6&5&1\\\ 4&4&3&6&2&5&1&0\\\ 5&5&2&1&6&4&0&3\\\ 6&1&0&3&5&2&6&4\\\
\hline\cr\end{array}$

Encryption. To send a message $b=630512403$, Bob computes from the known
$T=(\alpha,\beta,\gamma)$: $\alpha=(234)(0516)$; $\beta=(0321)(56)$;
$\gamma=(1236054)$, calculates isotopy $T^{(r,s,t)}$ for random numbers
$r=5,s=3,t=6$, i.e. $T^{(5,3,6)}$:

$\alpha^{5}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 5&6&4&2&3&1&0\\\ \end{pmatrix}$

$\beta^{3}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 1&2&3&0&4&6&5\\\ \end{pmatrix}$

$\gamma^{6}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 6&4&1&2&5&0&3\\\ \end{pmatrix}$

In our example $T^{(5,3,6)}$ we get: $\alpha^{5}=(0516)(243)$;
$\beta^{3}=(0123)(56)$; $\gamma^{6}=(0632145)$.

Then he calculates $T^{(mr,ns,kt)}$ using the public key:

$T^{(m,n,k)}=(\alpha^{m},\beta^{n},\gamma^{k})=(\alpha^{*},\beta^{*},\gamma^{*})$:

$\alpha^{*}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 6&5&2&3&4&0&1\\\ \end{pmatrix}$

$\beta^{*}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 2&3&0&1&4&5&6\\\ \end{pmatrix}$
$\gamma^{*}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 3&5&4&1&0&6&2\\\ \end{pmatrix}$

Then he raises these permutations, respectively, to the power $r=5,s=3,t=6$
and gets:

$\alpha^{*5}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 6&5&2&3&4&0&1\\\ \end{pmatrix}$
$\beta^{*3}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 2&3&0&1&4&5&6\\\ \end{pmatrix}$

$\gamma^{*6}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 4&3&6&0&2&1&5\\\ \end{pmatrix}$

$\alpha^{5m}=\alpha^{5}=(0615)$; $\beta^{3n}=\beta^{3}=(02)(13)$;
$\gamma^{6k}=(0426513)$, $(\gamma^{6k})^{(-1)}=(0315624)$.

As a result of the application of the new isotopy $T^{(5m,3n,6k)}$ to the
quasigroup $(Q,f)$ we obtain:

Table 9.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\alpha^{5m}&0&1&2&3&4&5&6\\\
\hline\cr 0&3&0&2&5&1&4&6\\\ 1&6&3&1&2&5&0&4\\\ 2&0&5&4&6&3&1&2\\\
3&4&1&3&0&2&6&5\\\ 4&2&4&0&1&6&5&3\\\ 5&5&2&6&4&0&3&1\\\ 6&1&6&5&3&4&2&0\\\
\hline\cr\end{array}$

Table 10.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\beta^{3n}&0&1&2&3&4&5&6\\\
\hline\cr 0&2&5&3&0&1&4&6\\\ 1&1&2&6&3&5&0&4\\\ 2&4&6&0&5&3&1&2\\\
3&3&0&4&1&2&6&5\\\ 4&0&1&2&4&6&5&3\\\ 5&6&4&5&2&0&3&1\\\ 6&5&3&1&6&4&2&0\\\
\hline\cr\end{array}$

Table 11.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr(\gamma^{6k})^{(-1)}&0&1&2&3&4&5&6\\\
\hline\cr 0&4&6&1&3&5&0&2\\\ 1&5&4&2&1&6&3&0\\\ 2&0&2&3&6&1&5&4\\\
3&1&3&0&5&4&2&6\\\ 4&3&5&4&0&2&6&1\\\ 5&2&0&6&4&3&1&5\\\ 6&6&1&5&2&0&4&3\\\
\hline\cr\end{array}$

To obtain $(T^{(mr,ns,kt)}(Q,f)b$, Bob uses the Markovski algorithm known to
Alice, with the known leader value $l=3$, then the ciphertext for
$b=6\,3\,0\,5\,1\,2\,4\,0\,3$ will look like: $v_{1}=3\cdot 6=6,v_{2}=6\cdot
3=2,v_{3}=2\cdot 0=0,v_{4}=0\cdot 5=0,v_{5}=0\cdot 1=6,v_{6}=6\cdot
2=5,v_{7}=5\cdot 4=3,v_{8}=3\cdot 0=1,v_{9}=1\cdot 3=1,b^{\prime}=620065311$.

Decryption. Alice knows $m=3,n=6,k=5$, so if she gets an isotopy $T^{(r,s,t)}$
and ciphertext $(T^{(mr,ns,kt)}(Q,f))b)=620065311$, she will calculate the
isotopy first $T^{(mr,ns,kt)}$ using $T^{(r,s,t)}=T^{(**,**,**)}$ :

$\alpha^{\ast\ast}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 5&6&4&2&3&1&0\\\
\end{pmatrix}$

$\beta^{\ast\ast}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 1&2&3&0&4&6&5\\\
\end{pmatrix}$

$\gamma^{\ast\ast}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 6&4&1&2&5&0&3\\\
\end{pmatrix}$

She calculates $T^{(mr,ns,kt)}$:

$\alpha^{**3}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 6&5&2&3&4&0&1\\\ \end{pmatrix}$

$\beta^{**6}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 2&3&0&1&4&5&6\\\ \end{pmatrix}$

$\gamma^{**5}=\begin{pmatrix}0&1&2&3&4&5&6\\\ 4&3&6&0&2&1&5\\\ \end{pmatrix}$

As a result, she receives the same table as Bob received in the encryption
process. For Table $(\gamma^{6k})^{-1}$ Alice builds a parastrophe $(23)$ used
in the Markovski algorithm for decryption:

Table 12.

$\begin{array}[]{|l|l|l|l|l|l|l|l|}\hline\cr\backslash&0&1&2&3&4&5&6\\\
\hline\cr 0&5&2&6&3&0&4&1\\\ 1&6&3&2&5&1&0&4\\\ 2&0&4&1&2&6&5&3\\\
3&2&0&5&1&4&3&6\\\ 4&3&6&4&0&2&1&5\\\ 5&1&5&0&4&3&6&2\\\ 6&4&1&3&6&5&2&0\\\
\hline\cr\end{array}$

and finally, using this table, she calculates b: $u_{1}=3\backslash
6=6,u_{2}=6\backslash 2=3,u_{3}=2\backslash 0=0,u_{4}=0\backslash
0=5,u_{5}=0\backslash 6=1,u_{6}=6\backslash 5=2,u_{7}=5\backslash
3=4,u_{8}=3\backslash 1=0,u_{9}=1\backslash 1=3.$

Therefore, $b=630512403$.

In this algorithm, isostrophy [6] can also be used instead of isotopy, the
modified algorithm instead of the Markovski algorithm and n-ary ($n>2$)
quasigroups [7; 8] instead of binary quasigroups.

A generalization of the Diffie-Hellman scheme of the open key distribution is
given in [9].

The generalization is based on the concepts of the left and right powers of
the elements of some non-associative groupoids.

For medial quasigroups, this approach is implemented in [10]. The protocol of
the elaboration of a common secret key based on Moufang loops is given in
[10].

This protocol is a generalization of the results from [11]. Generalizations of
the ElGamal scheme based on Moufang loops are given in [10].

In [12], the discrete logarithmic problem with Moufang loops is reduced to the
same problem over finite simple fields. Another generalization of the ElGamal
scheme based on quasi-automorphisms of quasigroups is presented in [10].

##  4 Conclusion

Today, different points of view on the same mathematical idea lead to
different generalizations. We considered in our work an analogue of the
ElGamal encryption system based on the Markovski algorithm. This algorithm is
under improvement and its other modifications are planned.

## References

  * [1] T. ElGamal. A public key cryptosystem and a signature scheme based on discrete logarithms. - IEEE Transactions on Information Theory, 31(4): p.469-472, 1985.
  * [2] B.Ya. Ryabko, A.N. Fionov. Cryptographic methods of information protection: a training manual. - M. Hotline-Telecom, 2005, p. 12-34. ISBN 5-89176-233-1.
  * [3] Wikipedia. ElGamal encryption, 2014. http:// en.wikipedia.org/wiki/ ElGamal-encryption.
  * [4] N.A. Moldovyan, A.V. Shcherbacov, and V.A. Shcherbacov. On some applications of quasigroups in cryptology. In Workshop on Foundations of Informatics, August 24-29, 2015, Chisinau, Proceedings, pages 331-341.
  * [5] V.A. Shcherbacov. On generalisation of Markovski cryptoalgorithm. In Workshop on General Algebra, February 26-March 1, 2015, Technische Universitat at Dresden, Technical Report, Technische Universitat at Dresden, Dresden, 36-37, 2015.
  * [6] V.A. Shcherbacov. On the structure of left and right F-, SM- and E-quasigroups. J. Gen. Lie Theory Appl., 3(3): p. 197-259, 2009.
  * [7] V.D. Belousov. n-Ary Quasigroups. Stiintsa, Kishinev, 1971 (in Russian).
  * [8] V.A. Shcherbacov. Quasigroups in cryptology. Comput. Sci. J. Moldova, 17(2): p. 193-228, 2009.
  * [9] S. Yu. Katyshev, V.T. Markov, and A.A. Nechaev. Utilization of nonassociative groupoids for the realization of an open key-distribution procedure. Diskret. Mat., 26: p. 45- 64, 2014 (in Russian).
  * [10] A.V. Gribov. Algebraic Non-Associative Structures and Its Applications in Cryptology. PhD thesis, Moscov State University, 2015 (in Russian).
  * [11] E. Stickel. A new method for exchanging secret keys. In Proceedings of the Third International Conference on Information Technology and Applications, volume 2, pages 426-430, 2005.
  * [12] G. Maze. Algebraic Methods For Constructing One-Way Trapdoor Functions. PhD thesis, University of Notre Dame, 2003.

Nadeghda Malyutina1

1 Ph.D. Student/Moldova State University

E–mail<EMAIL_ADDRESS>
Victor Shcherbacov2

2 Researcher/Vladimir Andrunachievici Institute of Mathematics and Computer
Science

E–mail<EMAIL_ADDRESS>