# MARKOV JUMP PROCESSES APPROXIMATING A NON-SYMMETRIC GENERALIZED DIFFUSION:
NUMERICS EXPLAINED TO PROBABILISTS††thanks: Supported by grant 0037014 of the
Ministry of Science, Higher Education and Sports, Croatia.

Nedžad Limić  Dept. of Mathematics, University of Zagreb, Bijenička 30, 10002
Zagreb, Croatia, e–mail<EMAIL_ADDRESS>
###### Abstract

Consider a non-symmetric generalized diffusion $X(\cdot)$ in ${{\hbox{\msbm
R}}}^{d}$ determined by the differential operator
$A(\mbox{\boldmath{$x$}})=-\sum_{ij}\partial_{i}a_{ij}(\mbox{\boldmath{$x$}})\partial_{j}+\sum_{i}b_{i}(\mbox{\boldmath{$x$}})\partial_{i}$.
In this paper the diffusion process is approximated by Markov jump processes
$X_{n}(\cdot)$, in homogeneous and isotropic grids $G_{n}\subset{{\hbox{\msbm
R}}}^{d}$, which converge in distribution to the diffusion $X(\cdot)$. The
generators of $X_{n}(\cdot)$ are constructed explicitly. Due to the
homogeneity and isotropy of grids, the proposed method for $d\geq 3$ can be
applied to processes for which the diffusion tensor
$\\{a_{ij}(\mbox{\boldmath{$x$}})\\}_{11}^{dd}$ fulfills an additional
condition. The proposed construction offers a simple method for simulation of
sample paths of non-symmetric generalized diffusion. Simulations are carried
out in terms of jump processes $X_{n}(\cdot)$. For $d=2$ the construction can
be easily implemented into a computer code.

AMS subject classification: 60H35, 60J22, 60J27, 60J60, 65C30

Key words: Symmetric diffusion, Approximation of diffusion, Simulation of
diffusion, Divergence form operators

##  1 INTRODUCTION

A symmetric tensor valued function $\mbox{\boldmath{$x$}}\mapsto
a(\mbox{\boldmath{$x$}})=\\{a_{ij}(\mbox{\boldmath{$x$}})\\}_{11}^{dd}$ which
is measurable, bounded and strictly positive definite on ${{\hbox{\msbm
R}}}^{d}$ defines a second order differential operator in divergence form on
${{\hbox{\msbm R}}}^{d}$,
$A_{0}(\mbox{\boldmath{$x$}})=-\sum_{ij}\partial_{i}a_{ij}(\mbox{\boldmath{$x$}})\partial_{j}$.
Each $A_{0}(\mbox{\boldmath{$x$}})$ determines a symmetric diffusion
$X(\cdot)$ in ${{\hbox{\msbm R}}}^{d}$. In [SZ] the process $X(\cdot)$ is
approximated by a sequence of Markov jump processes $X_{n}(\cdot)$ (MJP) in
grids, i.e. processes in continuous time and discrete state spaces. Coordinate
systems in in [SZ] depend locally on the structure of the tensor valued
function $a$. The object of the present analysis is a non-symmetric diffusion
determined by
$A(\mbox{\boldmath{$x$}})=A_{0}(\mbox{\boldmath{$x$}})+B(\mbox{\boldmath{$x$}})$,
where
$B(\mbox{\boldmath{$x$}})=\sum_{i}b_{i}(\mbox{\boldmath{$x$}})\partial_{i}$,
and a construction of its approximations by MJPs in lattices. For $d=2$ the
method is valid for any $a$, and for $d>2$ the method is valid for tensor
valued functions $a$ satisfying an additional constraint. In this way we offer
an efficient method for simulation of sample paths of a non-symmetric
generalized diffusion using MJPs. Each $X_{n}(\cdot)$ can be simulated by
well-known methods.

In the case of classical diffusion determined by an elliptic operator of the
form
$A(\mbox{\boldmath{$x$}})=-\sum_{ij}a_{ij}(\mbox{\boldmath{$x$}})\partial_{i}\partial_{j}+\sum_{i}b_{i}(\mbox{\boldmath{$x$}})\partial_{i}$,
where $a_{ij},b_{i}$ are Hölder continuous on ${{\hbox{\msbm R}}}^{d}$,
approximations by MJPs can be efficiently used to simulate the first exit from
a bounded set of ${{\hbox{\msbm R}}}^{d}$. However, such an approach is one of
several existing possibilities, and the motivation for the construction in
terms of MJPs happens to be of lesser importance, since the process has a
representation in terms of SDE which can be simulated straightforwardly. For
one-dimensional generalized diffusion, there exist representations in terms of
SDE as described by Étoré [Et] and Lejay & Martinez [LM], so that such
representations can be used for simulation. In the case of a process defined
by a differential operator in divergence form on ${{\hbox{\msbm R}}}^{d}$,
$d\geq 2$, there is no such natural representation, so approximations by MJPs
are essentially the only tool available for simulations.

A class of convergent approximations $X_{n}(\cdot)$ in [SZ] is constructed by
using discretizations of the corresponding Dirichlet form $(v,u)\mapsto
a(v,u)$. The functions $\partial_{i}v,\partial_{i}u$ in $a(v,u)$ are
approximated by forward difference operators in local basis which generally
varies. This approach is anticipated in [MW] without any remarks on the
convergence of the constructed MJPs $X_{n}(\cdot)$. In our approach, the MJPs
$X_{n}(\cdot)$ are constructed in terms of generators $A_{n}({\rm gen})$ on
the grids,

$G_{n}\ =\
\Big{\\{}h_{n}\,\sum_{i=1}^{d}\,k_{i}\mbox{\boldmath{$e$}}_{i}\,:\,k_{i}\in{{\hbox{\msbm
Z}}}\Big{\\}},\quad h_{n}=2^{-n},$ (1)

with a fixed basis $\\{\mbox{\boldmath{$e$}}_{i}\\}_{1}^{d}$. The index set of
grid-knots is denoted by $I_{n}$. In order to simplify expressions, we often
write $h$ instead of $h_{n}$. The generators $A_{n}({\rm gen})$ are
constructed explicitly from a general principle which is not directly related
to forward difference operators. Then the discretizations $h^{d}a_{n}(v,u)$ of
the original Dirichlet form are associated to the constructed generators
$A_{n}({\rm gen})$. It turns out that $h^{d}a_{n}(v,u)$ cannot be simply
obtained from discretizations of original form by using forward/backward
difference operators. The obtained class of $A_{n}({\rm gen})$ is not included
among generators constructed in [SZ].

Advantages and drawbacks of the present construction can be briefly described
as follows.

Advantages: The convergence of MJPs is proved for a non-symmetric generalized
diffusion. The generators $A_{n}({\rm gen})$ are explicitly given in terms of
values of functions $a_{ij},b_{i}$. For the case of $d=2$, the matrix entries
of $A_{n}({\rm gen})$ can be easily implemented into a computer code because
rotations of coordinates are avoided. For $d=2$, the construction holds for a
general matrix valued function $a(\cdot)$ on ${{\hbox{\msbm R}}}^{2}$.

Disadvantages: For $d\geq 3$, the proposed construction is not valid for all
$a(\cdot)$ on ${{\hbox{\msbm R}}}^{d}$. The restriction is defined by an
auxiliary matrix valued function $\hat{a}(\cdot)$, with the following matrix
entries:

$\hat{a}_{ii}\ =\ a_{ii},\quad\hat{a}_{ij}\ =\ -|a_{ij}|,\ i\neq j.$ (2)

The here proposed construction of $A_{n}({\rm gen})$ is valid only if
$\hat{a}(\cdot)$ is strictly positive definite on ${{\hbox{\msbm R}}}^{d}$.
For $d\geq 3$, there are simple examples of pairs $a,\hat{a}$, where $a$ is
definite and $\hat{a}$ is indefinite. We can say that (2) is valid if the off-
diagonal entries of $a(\cdot)$ are small in comparison with the diagonal ones.
When $\hat{a}(\mbox{\boldmath{$x$}})$ becomes indefinite, it is necessary to
apply a local rotation of coordinates, that ensures a diminishing of the off-
diagonal entries, thus ensuring the positive definiteness of
$\hat{a}(\mbox{\boldmath{$x$}})$ in the new coordinate system.

Now we can describe basic steps in the construction and the proofs. Let
$U(\cdot)$ be the strongly continuous semigroup in the Banach space
$\dot{C}({{\hbox{\msbm R}}}^{d})$ (continuous functions vanishing at infinity)
which is associated with the diffusion process $X(\cdot)$ [EK]. For each $n$
the space of discretizations of $\dot{C}({{\hbox{\msbm R}}}^{d})$ in terms of
grid-functions on $G_{n}$ is denoted by $\dot{l}_{\infty}(G_{n})$. The
generators $A_{n}({\rm gen})$ determine semigroups $U_{n}(\cdot)$ in
$\dot{l}_{\infty}(G_{n})$. There exist continuous mappings
$\Phi_{n}:\dot{l}_{\infty}(G_{n})\mapsto\dot{C}({{\hbox{\msbm R}}}^{d})$ and
$\Phi_{n}^{-1}:\dot{C}({{\hbox{\msbm R}}}^{d})\mapsto\dot{l}_{\infty}(G_{n})$
with the following properties: $\|\Phi_{n}\|=1,\|\Phi_{n}^{-1}\|\leq 1$ and
$\Phi_{n}^{-1}\Phi_{n}=I$ on $\dot{l}_{\infty}(G_{n})$. For
$f\in\dot{C}({{\hbox{\msbm R}}}^{d})$, the grid-function ${\bf
f}_{n}=\Phi_{n}^{-1}f$ is called the discretization of $f$. Now, for each
$t\geq 0$, the grid-functions ${\bf u}_{n}(t)=\Phi_{n}^{-1}U(t)f$ and
$U_{n}(t){\bf f}_{n}$ can be compared. We need the following relation:

$\lim_{n\to\infty}\>\sup\\{\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{n}(t){\bf
f}_{n}-\Phi_{n}^{-1}U(t)f\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}\,:\,t\in[0,1]\\}\
=\ 0,$ (3)

where
$\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,\cdot\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}$
is the $l_{\infty}(G_{n})$-norm. By Theorem 2.11 of [EK], this relation is a
sufficient condition for the convergence in distribution of MJPs
$X_{n}(\cdot)$ to the diffusion process $X(\cdot)$.

Definitions of various objects are given in Section 2. In Section 3, a class
of generators $A_{n}({\rm gen})$ is constructed explicitly. Grid-functions
${\bf u}_{n}(t)=\exp(A_{n}({\rm gen})t){\bf f}_{n}$ and their images
$u(n,t)=\Phi_{n}{\bf u}_{n}(t)$ in $\dot{C}({{\hbox{\msbm R}}}^{d})$ are
studied in Section 4. The convergence $u(n,t)\mapsto U(t)f$ in
$\dot{C}({{\hbox{\msbm R}}}^{d})$ is proved initially for the case of smooth
coefficients $a_{ij},b_{i}$, and then for general coefficients $a_{ij},b_{i}$.
A numerical example is given in the last section.

##  2 PRELIMINARIES

The Euclidean norm in ${{\hbox{\msbm R}}}^{d}$ is denoted by $|\cdot|$. All
the open subsets of ${{\hbox{\msbm R}}}^{d}$, considered in this work, are
bounded and connected open sets with Lipshitz boundary [Ma, Se]. We call a
subset of this kind a Lipshitz domain, denote it by D, and its boundary by
$\partial D$.

The Banach spaces of functions $C^{(k)}({{\hbox{\msbm
R}}}^{d}),C_{0}^{(k)}({{\hbox{\msbm R}}}^{d})=C_{0}({{\hbox{\msbm
R}}}^{d})\cap C^{(k)}({{\hbox{\msbm R}}}^{d})$ are defined as usual,
$C_{0}({{\hbox{\msbm R}}}^{d})$ being the linear space of continuous functions
with compact support. Their norms are denoted by $\|\cdot\|_{\infty}^{(k)}$.
The closure of functions in $C({{\hbox{\msbm R}}}^{d})$ with compact support
determines the Banach space $\dot{C}({{\hbox{\msbm R}}}^{d})$. The Hölder
space of parameter $k+\alpha,k\in{{\hbox{\msbm N}}}_{0},\alpha\in(0,1)$, is
denoted by $C^{(k+\alpha)}({{\hbox{\msbm R}}}^{d})$ and defined as the
completion of $C_{0}^{(\infty)}({{\hbox{\msbm R}}}^{d})$ in the norm:

$\|u\|^{(k+\alpha)}\ =\
\|u\|_{\infty}^{(k)}\,+\,\sup\Big{\\{}\frac{|\partial^{k}u(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$z$}})-\partial^{k}u(\mbox{\boldmath{$x$}})|}{|\mbox{\boldmath{$z$}}|^{\alpha}}\>:\>\mbox{\boldmath{$x$}}\in{{\hbox{\msbm
R}}}^{d},0<|\mbox{\boldmath{$z$}}|\leq 1\Big{\\}}.$

The $L_{p}$-spaces as well as Sobolev $W_{p}^{1}$-spaces are defined in a
standard way [Ma, Se]. Their norms are denoted by $\|\cdot\|_{p}$ and
$\|\cdot\|_{p,1}$, respectively. For each $p,\>1\leq p\leq\infty$, the norm of
$W_{p}^{1}({{\hbox{\msbm R}}}^{d})$ is defined by
$\|u\|_{p,1}=\big{(}\|u\|_{p}^{2}+\|\nabla u\|_{p}^{2}\big{)}^{1/2}$, where
$\|\,\nabla
u\,\|_{p}\,=\,\big{(}\sum_{j=1}^{d}\,\|\,\partial_{j}u\,\|_{p}^{2}\big{)}^{1/2}$.
We say that a function $f$ on ${{\hbox{\msbm R}}}^{d}$ belongs to a class
$C^{(k+\alpha)}$ on ${{\hbox{\msbm R}}}^{d}$ if $\|f\|^{(k+\alpha)}<\infty$.

In this article we consider a $2^{{\rm nd}}$-order elliptic operator on
${{\hbox{\msbm R}}}^{d}$,

$A(\mbox{\boldmath{$x$}})=-\>\sum_{i,j=1}^{d}\partial_{i}a_{ij}(\mbox{\boldmath{$x$}})\partial_{j}\>+\>\sum_{j=1}^{d}b_{j}(\mbox{\boldmath{$x$}})\,\partial_{j},$
(4)

for which the coefficients fulfill the following:

###### Assumption 2.1

The functions $a_{ij}=a_{ji}$, $b_{i},i,j=1,2,\dots,d$, are measurable on
${{\hbox{\msbm R}}}^{d}$ and have the following properties:

a) There exist positive numbers
$\underline{M},\overline{M},\>0<\underline{M}\leq\overline{M}$, such that the
strict ellipticity is valid:

$\underline{M}\,|\mbox{\boldmath{$z$}}|^{2}\,\leq\,\sum_{i,j=1}^{d}a_{ij}(\mbox{\boldmath{$x$}})z_{i}\bar{z}_{j}\leq\,\overline{M}\,|\mbox{\boldmath{$z$}}|^{2},\quad\mbox{\boldmath{$x$}},\mbox{\boldmath{$z$}}\in{{\hbox{\msbm
R}}}^{d}.$ (5)

b) The functions $b_{i}$ are bounded on ${{\hbox{\msbm R}}}^{d}$.

In our analysis we regularly use the notation
$A_{0}(\mbox{\boldmath{$x$}})=-\sum_{ij=1}^{d}\partial_{i}a_{ij}(\mbox{\boldmath{$x$}})\partial_{j},B(\mbox{\boldmath{$x$}})=\sum_{j=1}^{d}b_{j}(\mbox{\boldmath{$x$}})\partial_{j}$
and
$A(\mbox{\boldmath{$x$}})=A_{0}(\mbox{\boldmath{$x$}})+B(\mbox{\boldmath{$x$}})$.
The following real bilinear form on $W_{q}^{1}(D)\times W_{p}^{1}(D)$,
$1/p+1/q=1$, defined by:

$a(v,u)\ =\
\sum_{i,j=1}^{d}\>\int\>a_{ij}(\mbox{\boldmath{$x$}})\>\partial_{i}v(\mbox{\boldmath{$x$}})\>\partial_{j}u(\mbox{\boldmath{$x$}})\>d\mbox{\boldmath{$x$}}-\sum_{i=1}^{d}\>\int\>b_{i}(\mbox{\boldmath{$x$}})\>v(\mbox{\boldmath{$x$}})\,\partial_{i}u(\mbox{\boldmath{$x$}})\>d\mbox{\boldmath{$x$}}$
(6)

is associated with the differential operator $A(\mbox{\boldmath{$x$}})$.

A basic result towards a proof of (3) is the following theorem [St]:

###### THEOREM 2.1

Let the differential operator $A(\mbox{\boldmath{$x$}})$ be defined by (4) and
Assumption 2.1. Then $-A(\mbox{\boldmath{$x$}})$ has the closure in
$\dot{C}({{\hbox{\msbm R}}}^{d})$ and generates a Feller semigroup in
$\dot{C}({{\hbox{\msbm R}}}^{d})$.

The semigroup of this theorem has the restriction to a closed subspace of
$\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})=\dot{C}({{\hbox{\msbm
R}}}^{d})\cap C^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ which is also a strongly
continuous semigroup.

###### COROLLARY 2.1

Let $A(\mbox{\boldmath{$x$}})$ be defined as in Theorem 2.1. There exists an
$\alpha\in(0,1)$ and $\sigma\geq 0$ that depend only on
$\underline{M},\overline{M}$ and $\|\mbox{\boldmath{$b$}}\|_{\infty}$, such
that the following two assertions are valid:

(i) The operators $U(t)$ in $\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ are
bounded uniformly on segments of $[0,\infty)$.

(ii) There exists a closed subspace $F^{(\alpha)}({{\hbox{\msbm
R}}}^{d})\subseteq\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ such that the
closure of $-A(\mbox{\boldmath{$x$}})$ in $F^{(\alpha)}({{\hbox{\msbm
R}}}^{d})$ generates a strongly continuous semigroup $U(\cdot)$ in
$F^{(\alpha)}({{\hbox{\msbm R}}}^{d})$, with
$\|U(t)\|_{\infty}^{(\alpha)}\leq\exp(\sigma t)$.

Proof:  The Feller semigroup of Theorem 2.1 can be represented as
$t\mapsto\exp(t)V(t)$, where $V(\cdot)$ is the strongly continuous semigroup
generated by the closure of $-(I+A(\mbox{\boldmath{$x$}}))$ in
$\dot{C}({{\hbox{\msbm R}}}^{d})$. Due to Theorem 2.1 the differential
operator $I+A(\mbox{\boldmath{$x$}})$ has the bounded inverse $(I+A)^{-1}$ in
$\dot{C}({{\hbox{\msbm R}}}^{d})$. For any $f\in\dot{C}({{\hbox{\msbm
R}}}^{d})$ the continuous function $(t,\mbox{\boldmath{$x$}})\mapsto
u(t,\mbox{\boldmath{$x$}})$, defined by $u(t)=V(t)f$, is a solution to the
following IVP for PDE

$\partial_{t}u(t,\mbox{\boldmath{$x$}})+(I+A(\mbox{\boldmath{$x$}}))u(t,\mbox{\boldmath{$x$}})\
=\ 0,$ (7)

with the initial condition
$u(0,\mbox{\boldmath{$x$}})=f(\mbox{\boldmath{$x$}})$. In addition,
$\|u(t)\|_{\infty}\leq\|f\|_{\infty}$ for all $t\geq 0$. It turns out that a
solution to (7), for each $t>0$, is an element of
$\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$. For an initial condition
$f\in\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ we have a stronger result, a
solution to (7) is an element of $\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$
uniformly with respect to $t\in[0,1]$. This follows from the following
arguments. Let us consider the balls
$B_{1}(\mbox{\boldmath{$x$}})\subset{{\hbox{\msbm R}}}^{d}$, centered at
$\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{d}$ and having the radius equal
to 1, and let us consider the restrictions $u(t)|B_{1}(\mbox{\boldmath{$x$}})$
for any $\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{d}$ and any $t\in[0,1]$.
According to the basic theorem of Section 10, Chapter 3 of [LSU], there exists
$\alpha\in(0,1)$ depending on
$\underline{M},\overline{M},\|\mbox{\boldmath{$b$}}\|_{\infty}$ and $\beta$
depending on $\underline{M},\overline{M},\|\mbox{\boldmath{$b$}}\|_{\infty}$
and $\|f\|_{\infty}^{(\alpha)}$ such that a solution $u$ to (7) on
$B_{1}(\mbox{\boldmath{$x$}})$ is an element of
$C^{(\alpha)}(B_{1}(\mbox{\boldmath{$x$}}))$, and

$\sup_{t\in[0,1]}\,\|u(t)|B_{1}(\mbox{\boldmath{$x$}})\|^{(\alpha)}\,+\,\sup_{t,s\in[0,1]}\,\sup_{\mbox{\boldmath{\scriptsize$y$}}\in
B_{1}(\mbox{\boldmath{\scriptsize$x$}})}\,\frac{|u(t,\mbox{\boldmath{$y$}})-u(s,\mbox{\boldmath{$y$}})|}{|t-s|^{\alpha}}\
\leq\ \beta.$ (8)

This inequality implies that the operators $V(t),t\in[0,1]$, are linear in
$\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$, with a norm which is uniformly
bounded with respect to $t\in[0,1]$. Hence, the operators $V(t)$ define a
semigroup of bounded operators in $\dot{C}^{(\alpha)}({{\hbox{\msbm
R}}}^{d})$. This proves the assertion (i) of the corollary.

For any $g\in\dot{C}({{\hbox{\msbm R}}}^{d})$ the function $f=(I+A)^{-1}g$
belongs to the domain $\mathfrak{D}(A)$ of the closure of
$A(\mbox{\boldmath{$x$}})$ in $\dot{C}({{\hbox{\msbm R}}}^{d})$, and
$V(t)f=(I+A)^{-1}V(t)g$. The linear space spanned by
$f=(I+A)^{-1}g,\>g\in\dot{C}({{\hbox{\msbm R}}}^{d})$, is denoted by
$E({{\hbox{\msbm R}}}^{d})$. Let us assume now that we can prove
$f=(I+A)^{-1}g\in\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ for any
$g\in\dot{C}({{\hbox{\msbm R}}}^{d})$, i.e.
$\|f\|^{(\alpha)}\leq\beta\|g\|_{\infty}$, where $\beta$ is a number depending
on $\underline{M},\overline{M}$ and $\|\mbox{\boldmath{$b$}}\|_{\infty}$. This
would imply that $E({{\hbox{\msbm R}}}^{d})$ consists of elements in
$\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ and its completion in the
$\|\cdot\|^{(\alpha)}$-norm is a closed space $F^{(\alpha)}({{\hbox{\msbm
R}}}^{d})\subseteq\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$. In addition,
the above assumption would imply the following inequality:

$\|(V(t)-I)f\|^{(\alpha)}\ \leq\ \beta\,\|(V(t)-I)g\|_{\infty},$ (9)

i.e. the continuity of the function $t\mapsto V(t)$ on a dense linear space
$E({{\hbox{\msbm R}}}^{d})\subset F^{(\alpha)}({{\hbox{\msbm R}}}^{d})$. Since
the operators $V(t)$ are bounded in $\dot{C}^{(\alpha)}({{\hbox{\msbm
R}}}^{d})$ uniformly with respect to $t\in[0,1]$, the obtained inequality
would imply that the operators $V(\cdot)$ define a strongly continuous
semigroup in the Banach space $F^{(\alpha)}({{\hbox{\msbm R}}}^{d})$. Then the
semigroup $U(\cdot)$ of (ii) is also strongly continuous on
$F^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ and the existence of $\sigma\geq 0$ in
(ii) follows from the theory of strongly continuous semigroups.

Therefore, in order to prove (ii) it remains to show that $(I+A)^{-1}$ maps
$\dot{C}({{\hbox{\msbm R}}}^{d})$ into $\dot{C}^{(\alpha)}({{\hbox{\msbm
R}}}^{d})$. For this we consider the elliptic problem:

$\big{(}I\,+\,A(\mbox{\boldmath{$x$}})\big{)}\,u(\mbox{\boldmath{$x$}})\ =\
f(\mbox{\boldmath{$x$}}),\quad\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{d},$
(10)

where $f\in\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$. The existence of
solution $u=(I+A)^{-1}f$ in $\dot{C}({{\hbox{\msbm R}}}^{d})$ is granted by
Theorem 2.1, and in addition we have $\|u\|_{\infty}\leq\|f\|_{\infty}$. In
order to prove that $u\in\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$, it is
sufficient to verify that $u|B_{1}(\mbox{\boldmath{$x$}})\in
C^{(\alpha)}(B_{1}(\mbox{\boldmath{$x$}}))$ for any $x$ and a fixed $\alpha$.
To obtain this, we apply the theorem in Section 14, Chapter 3 of [LU]. The
following estimate is valid:

$\|u|B_{1}(\mbox{\boldmath{$x$}})\|^{(\alpha)}\ \leq\ \beta\,\|f\|_{\infty},$

where both, $\alpha\in(0,1)$ and $\beta$, depend on
$\underline{M},\overline{M},\|\mbox{\boldmath{$b$}}\|_{\infty}$.

The parameters $\alpha\in(0,1)$ for the parabolic (7) and elliptic problems
(10) are not necessarily equal. We choose the minimum to obtain (9), and
consequently the assertion (ii). QED

The assertion (i) of Corollary 2.1 is used in the proof of Theorem 4.1, which
establishes the basic inequality (3).

The linear space of grid-functions on $G_{n}$ is denoted by $l(G_{n})$.
Elements of $l(G_{n})$ are also called columns. Columns are denoted by ${\bf
u},{\bf v}$ etc, while their entries are denoted by
$u_{\mbox{\boldmath{\scriptsize$k$}}},v_{\mbox{\boldmath{\scriptsize$l$}}}$
etc. Thus a column ${\bf u}_{n}$ has entries $({\bf
u}_{n})_{\mbox{\boldmath{\scriptsize$k$}}},\mbox{\boldmath{$k$}}\in I_{n}$.
Columns with a finite number of components span a linear space $l_{0}(G_{n})$.
The completion of $l_{0}(G_{n})$ in the
$\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,\cdot\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}$-norm
is denoted by $\dot{l}(G_{n})$. The corresponding $l_{p}$-spaces are denoted
by $l_{p}(G_{n})$ and their norms by
$\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,\cdot\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{p}$.
For $p<\infty$ this norm is defined by
$\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,{\bf
u}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{p}=\big{[}\sum_{\mbox{\boldmath{\scriptsize$k$}}\in
I_{n}}\>|u_{\mbox{\boldmath{\scriptsize$k$}}}|^{p}\big{]}^{1/p}$, and for
$p=\infty$ by $\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,{\bf
u}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}=\sup\\{|u_{\mbox{\boldmath{\scriptsize$k$}}}|:\mbox{\boldmath{$k$}}\in{{\hbox{\msbm
Z}}}^{d}\\}$. The scalar product in $l_{2}(G_{n})$ is denoted by
$\langle\,\cdot|\cdot\,\rangle$ and sometimes by $(\cdot|\cdot)$.

Let $f\in C({{\hbox{\msbm R}}}^{d})$ and the column ${\bf f}_{n}\in l(G_{n})$
be defined by its components $({\bf
f}_{n})_{\mbox{\boldmath{\scriptsize$k$}}}=f(h\mbox{\boldmath{$k$}})$, where
$h=h_{n}=2^{-n}$. Then the mapping $C({{\hbox{\msbm R}}}^{d})\owns
f\mapsto{\bf f}_{n}\in l(G_{n})$ is called the discretization of $f$.

The shift operator
$Z(\mbox{\boldmath{$x$}}),\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{d}$,
acting on functions $f:{{\hbox{\msbm R}}}^{d}\mapsto{{\hbox{\msbm R}}}$, is
defined by
$\big{(}Z(\mbox{\boldmath{$z$}})f\big{)}(\mbox{\boldmath{$x$}})=f(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$z$}})$.
Similarly we define the discretized shift operator. The shift operator
$Z_{n}(r,i)$ by $r$ units in the direction $\mbox{\boldmath{$e$}}_{i}$ is
defined by $\big{(}Z_{n}(r,i){\bf
u}_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}}=({\bf
u}_{n})_{\mbox{\boldmath{\scriptsize$l$}}}$, where
$\mbox{\boldmath{$l$}}=\mbox{\boldmath{$k$}}+r\mbox{\boldmath{$e$}}_{i}$. The
partial derivatives of $u\in C^{(1)}({{\hbox{\msbm R}}}^{d})$ with respect to
a grid step $h$ are discretized by forward/backward finite difference
operators in the usual way,

$\begin{array}[]{c}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(th)u(\mbox{\boldmath{$x$}})\
=\
\frac{1}{th}\big{(}u(\mbox{\boldmath{$x$}}+th\mbox{\boldmath{$e$}}_{i})\>-\>u(\mbox{\boldmath{$x$}})\big{)},\\\
\widehat{\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,}_{i}(th)u(\mbox{\boldmath{$x$}})\
=\
\frac{1}{th}\big{(}u(\mbox{\boldmath{$x$}})\>-\>u(\mbox{\boldmath{$x$}}-th\mbox{\boldmath{$e$}}_{i})\big{)},\end{array}\quad\mbox{\boldmath{$x$}}\in{{\hbox{\msbm
R}}}^{d},\ t>0.$ (11)

Let $r\in{{\hbox{\msbm Z}}}\setminus\\{0\\}$. Discretizations of the functions
$\partial_{i}u$ on $G_{n}$, denoted by $U_{i}(r){\bf u}_{n},V_{i}(r){\bf
u}_{n}$, are defined by:

$\big{(}U_{i}(r)\,{\bf
u}_{n}\big{)}_{\mbox{\boldmath{\scriptsize$m$}}}\>=\>\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}\big{(}rh\big{)}\>u(\mbox{\boldmath{$x$}}_{\mbox{\boldmath{\scriptsize$m$}}}),\quad\big{(}V_{i}(r)\,{\bf
u}_{n}\big{)}_{\mbox{\boldmath{\scriptsize$m$}}}\>=\>\widehat{\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,}_{i}\big{(}rh\big{)}\>u(\mbox{\boldmath{$x$}}_{\mbox{\boldmath{\scriptsize$m$}}}),$

where $\mbox{\boldmath{$x$}}_{m}\in G_{n}$. Then

$\begin{array}[]{l}U_{i}(r)\ =\ (rh)^{-1}(Z_{n}(r,i)\,-\,I),\\\ V_{i}(r)\ =\
(rh)^{-1}\big{(}I-Z_{n}(-r,i)\big{)}\>=\>U_{i}(-r)\>=\>\,-\,U_{i}(r)^{T}.\end{array}$
(12)

Therefore we have $U_{i}(-r)=U_{i}(r)\,Z_{n}(-r,i)=Z_{n}(-r,i)\,U_{i}(r)$, and
similarly for $V_{i}(\cdot)$. We use the abbreviations
$U_{i}=U_{i}(1),V_{i}=V_{i}(1)$.

A matrix $A_{n}$ on $G_{n}$ is called a matrix of positive type if the
diagonal entries of $A_{n}$ are positive, off-diagonal entries are non-
positive and the row sums are non-negative. If $A_{n}(gen)$ is the generator
of a MJP in $G_{n}$, then $-A_{n}(gen)$ is a matrix of positive type.

##  3 CONSTRUCTION OF GENERATORS OF MJPs

To discretize $A(\mbox{\boldmath{$x$}})$ means to associate to
$A(\mbox{\boldmath{$x$}})$ a sequence of matrices $A_{n}$ on
$G_{n},n\in{{\hbox{\msbm N}}}$, with the following properties:

$a(v,u)\ =\ \lim_{n\to\infty}\>h^{d}\langle\,{\bf v}_{n}|A_{n}{\bf
u}_{n}\,\rangle,\quad v,u\in C_{0}^{(1)}({{\hbox{\msbm R}}}^{d}).$

The terminology ”discretizations” of $A(\mbox{\boldmath{$x$}})$ instead of
approximations of $A(\mbox{\boldmath{$x$}})$ appears to be more suitable at
the beginning of the construction, since the convergence analysis is postponed
to Section 4.

We wish to emphasize that discretizations $A_{n}$ of the differential operator
$A_{0}(\mbox{\boldmath{$x$}})$ are derived from a general principle, similar
to the one exploited in [LR2]. This method is not based on finite difference
formulas. Nevertheless, bilinear forms need to be associated to $A_{n}$ so
that $A_{n}$ can be derived from the corresponding variational equalities. The
constructed bilinear forms can be considered as the discretizations of the
original form (6). These forms are basic objects in our proof of convergence
of discretizations.

Discretizations to be considered in this section are possible if certain
conditions on $a_{ij}$ are fulfilled. The required conditions are stronger
than those given in Assumption 2.1. By relaxing them gradually as $n\to\infty$
we obtain discretizations for a general $A(\mbox{\boldmath{$x$}})$ given in
Assumption 2.1.

To each pair $\mbox{\boldmath{$v$}}\in
G_{n},\mbox{\boldmath{$p$}}\in{{\hbox{\msbm N}}}^{d}$, there is associated a
rectangle
$C_{n}(\mbox{\boldmath{$p$}},\mbox{\boldmath{$v$}})=\prod_{i=1}^{d}[v_{i},v_{i}+hp_{i})$
with the ”lower left” vertex $v$ and the edge of size $hp_{i}$ in the $i$-th
coordinate direction. These rectangles define a partition of ${{\hbox{\msbm
R}}}^{d}$. Apart from these rectangles, we will need the closed rectangles,

$S_{n}(\mbox{\boldmath{$p$}},\mbox{\boldmath{$v$}})\ =\
\prod_{i=1}^{d}\>[v_{i}-hp_{i}\,,\,v_{i}+hp_{i}],$ (13)

which are centered at the grid-knots $v$. Evidently,
$S_{n}(\mbox{\boldmath{$p$}},\mbox{\boldmath{$v$}})$ is the union of closures
of those rectangles $C_{n}(\mbox{\boldmath{$p$}},\mbox{\boldmath{$x$}})$ which
contain the grid-knot $v$.

A discretization $A_{n}$ is defined in terms of its matrix entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$,
where $h\mbox{\boldmath{$k$}},h\mbox{\boldmath{$l$}}\in G_{n}$. For a fixed
$\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}\in G_{n}$ the set of all the
grid-knots $\mbox{\boldmath{$y$}}=h\mbox{\boldmath{$l$}}$ such that
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}\neq
0$ is denoted by ${\cal N}(\mbox{\boldmath{$x$}})$ and called the numerical
neighborhood of $A_{n}$ at $\mbox{\boldmath{$x$}}\in G_{n}$. The set ${\cal
N}(\mbox{\boldmath{$x$}})$ contains always a subset consisting of $x$ and $2d$
elements $\mbox{\boldmath{$x$}}\pm h\mbox{\boldmath{$e$}}_{i},i=1,2,\ldots,d$.
Additional elements of ${\cal N}(\mbox{\boldmath{$x$}})$ are possible as we
shall see, depending on the sign of $a_{ij},i\neq j$. In terms of the MJP
$X_{n}(\cdot)$, the set ${\cal N}(\mbox{\boldmath{$x$}})$ consists of the
states of possible jumps from the state $x$. Let us point out that the sets
${\cal N}(\mbox{\boldmath{$x$}})$ vary with $n$, that is for two grids
$G_{n},G_{m},n\neq m$ and $\mbox{\boldmath{$x$}}\in G_{n}\cap G_{m}$, the
corresponding numerical neighborhoods ${\cal N}(\mbox{\boldmath{$x$}})$ are
different.

### General setup

In order to give a comprehensive insight into the proposed construction of the
discretizations, it is convenient to initially consider a differential
operator $A_{0}=-\sum_{ij}a_{ij}\partial_{i}\partial_{j}$ having a constant
diffusion tensor $a=\\{a_{ij}\\}_{11}^{dd}$. In this case, the matrices
$A_{n}$ need to have a property which is called the consistency. Let
$\mbox{\boldmath{$x$}}\mapsto
p(\mbox{\boldmath{$x$}})=p_{0}+\sum_{i}p_{i}x_{i}+\sum_{ij}p_{ij}x_{i}x_{j}$,
$p_{0},p_{i},p_{ij}\in{\hbox{\msbm R}}$, be a second degree polynomial in
arguments $x_{i}$, and let ${\bf p}_{n}$ be its discretization on grid
$G_{n}$. Then the consistency holds if the following identities are valid:

$\delta(A,n,\mbox{\boldmath{$x$}})p(\mbox{\boldmath{$x$}})\ :=\
A(\mbox{\boldmath{$x$}})p(\mbox{\boldmath{$x$}})-\big{(}A_{n}{\bf
p}_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}}\ =\
0,\quad\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}\in G_{n}.$ (14)

These consistency conditions are sufficient for proving the convergence of
$U_{n}(t){\bf f}_{n}$ to $\Phi_{n}U(t)f$ in the Banach space of continuous
functions on ${{\hbox{\msbm R}}}^{d}$ as specified in (3). Actually, the first
step of construction of matrices $A_{n}$ begins with a search for those
matrices $A_{n}$ which are simultaneously of positive type and fulfill the
conditions $\delta(A,n,\mbox{\boldmath{$x$}})p(\mbox{\boldmath{$x$}})=0$ on
$G_{n}$ [LR1].

We request that the discretizations $A_{n}$ have the following properties:

a) The numerical neighborhoods ${\cal N}(\mbox{\boldmath{$x$}})\subset G_{n}$
resemble each other, that is ${\cal
N}(\mbox{\boldmath{$x$}})=\mbox{\boldmath{$x$}}+{\cal N}({\bf 0})$.

b) The numerical neighborhoods ${\cal N}({\bf 0})\subset
G_{n},n\in{{\hbox{\msbm N}}}$, resemble each other, i.e. if
$\mbox{\boldmath{$y$}}=h_{n}\mbox{\boldmath{$l$}}\in{\cal N}({\bf 0})\subset
G_{n}$ then $\mbox{\boldmath{$y$}}^{\prime}=h_{m}\mbox{\boldmath{$l$}}\in{\cal
N}({\bf 0})\subset G_{m}$, $n,m\in{\hbox{\msbm N}},n\neq m$.

c) The matrices $A_{n}$ are symmetric.

d) The matrices $A_{n}$ are consistent discretizations of $A_{0}$.

Now we have the following result.

###### LEMMA 3.1

If $\hat{a}$ is positive definite then there exist matrices $A_{n}$ of
positive type fulfilling the conditions a)-d). The non-trivial entries of
$A_{n}$ are defined in terms of $d$ natural numbers
$r_{1},r_{2},\ldots,r_{d}\in{\hbox{\msbm N}}$ by the following expressions:

$\begin{array}[]{llll}\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\,\mbox{\boldmath{\scriptsize$k$}}}&=&-\,\sum_{h\mbox{\boldmath{\scriptsize$l$}}\in{\cal
N}(h\mbox{\boldmath{\scriptsize$k$}}),\
\mbox{\boldmath{\scriptsize$l$}}\neq\mbox{\boldmath{\scriptsize$k$}}}\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}},&\\\
\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\,\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}&=&-\,\frac{1}{h^{2}}\,\left[a_{ii}-\sum_{m\neq
i}\frac{r_{i}}{r_{m}}|a_{im}|\right],&i=1,2,\ldots,d,\\\
\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\,\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$z$}}(i,j,\pm)}&=&-\,\frac{1}{h^{2}r_{i}r_{j}}\,|a_{ij}|,&i,j=1,2,\ldots,d,\end{array}$
(15)

where in the last line,
$\mbox{\boldmath{$z$}}(i,j,+)=r_{i}\mbox{\boldmath{$e$}}_{i}+r_{j}\mbox{\boldmath{$e$}}_{j}$
corresponds to the case $a_{ij}>0$, and
$\mbox{\boldmath{$z$}}(i,j,-)=r_{i}\mbox{\boldmath{$e$}}_{i}-r_{j}\mbox{\boldmath{$e$}}_{j}$
corresponds to the case $a_{ij}<0$.

Proof:  For a matrix $A_{n}$ defined by (15) the conditions a)- c) are
obviously satisfied. (Recall that ${\cal N}(h\mbox{\boldmath{$k$}})$ consists
of those grid-knots $h\mbox{\boldmath{$l$}}\in{\cal
N}(h\mbox{\boldmath{$k$}})$ for which
$\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\,\mbox{\boldmath{\scriptsize$l$}}}$
are non-trivial). It remains to prove the condition d), and that $A_{n}$ is of
positive type if $\hat{a}$ is a positive definite matrix.

The condition d) is equivalent to the following property. Let
$\mbox{\boldmath{$x$}}\mapsto
p(\mbox{\boldmath{$x$}})=\sum_{ij}p_{ij}x_{i}x_{j}+\sum_{i}p_{i}x_{i}+p_{0}$
be a second degree polynomial. Then d) is valid if and only if

$\sum_{h\mbox{\boldmath{\scriptsize$k$}}\in{\cal N}({\bf
0})}\>\big{(}A_{n}\big{)}_{{\bf
0}\mbox{\boldmath{\scriptsize$k$}}}\,p(h\mbox{\boldmath{$k$}})\ =\
-2\,\sum_{i\leq j}a_{ij}p_{ij}.$

For $A_{n}$ defined by (15) this identity can be easily verified by using the
monomials $p(\mbox{\boldmath{$x$}})=x_{i}x_{j}$ with $i\neq j$ and
$p(\mbox{\boldmath{$x$}})=x_{i}^{2}$, for all $i,j=1,2,\ldots,d$.

In the last step of the proof we show that there exists a sequence
$r_{1},r_{2},\ldots,r_{d}\in{\hbox{\msbm N}}$ for which the brackets in the
second line of (15) are positive, thus ensuring the matrix $A_{n}$ to be of
positive type. For this purpose we assume that the matrix $\hat{a}$ is
positive definite. Let us consider the eigenvalue problem
$\hat{a}\mbox{\boldmath{$w$}}=\lambda\mbox{\boldmath{$w$}}$, and let us also
assume that the matrix $\hat{a}$ is irreducible. For $\mu>0$ sufficiently
large the irreducible matrix $\mu I+\hat{a}$ has the inverse $(\mu
I+\hat{a})^{-1}$ with positive entries, so the Perron-Frobenius theorem can be
applied to $(\mu I+\hat{a})^{-1}$. Thus the eigenvector corresponding to its
maximal eigenvalue is positive. This result can be also formulated in terms of
the problem $\hat{a}\mbox{\boldmath{$w$}}=\lambda_{1}\mbox{\boldmath{$w$}}$
for the minimal eigenvalue $\lambda_{1}$ of $\hat{a}$. We have
$\lambda_{1}\mbox{\boldmath{$w$}}>0$ and consequently
$\hat{a}\mbox{\boldmath{$w$}}>0$. This inequality can be rewritten as
$a_{ii}-\sum_{m\neq i}(q_{i}/q_{m})|a_{im}|>0$ where $q_{i}=w_{i}^{-1}$. One
can find rational approximations $r_{i}/r_{m}$ of $q_{i}/q_{m}$ which also
fulfil the obtained inequalities. If $\hat{a}$ is not irreducible, it can be
rewritten in a block diagonal matrix form with irreducible blocks. The
previous construction can be applied to each block. QED

In the case of $d=2$, either both $a,\hat{a}$ are positive definite or neither
is. For $d=3$, there are symmetric positive definite matrices $a$ for which
$\hat{a}$ are indefinite. For instance, the symmetric matrix $a$ of order
$d=3$ defined by $a_{ii}=1,a_{12}=a_{23}=-1/\sqrt{2}$ has positive eigenvalues
for the case of $a_{13}>0$ and a negative eigenvalue for the case of
$a_{13}<0$.

For $d=2$ and $r_{1}=3,r_{2}=1$, two possible numerical neighborhoods ${\cal\
N}(\mbox{\boldmath{$x$}})$ are illustrated in Figure 1.

###### Figure 1

units ¡1.0pt,1.0pt¿ x from -40 to 300, y from 20 to 100 40 60, 80 60 / 60 40,
60 80 / 0 40, 120 80 / 200 60, 240 60 / 220 40, 220 80 / 160 80, 280 40 /
=-A_n$isageneratorofaMJPin$G_n$.Letusconsidernowfunctions$a_ij$on$ℝ^d$definingadiffusiontensorateach$$x$∈ℝ^d$.Letusassumeinadditionthatthematrix$^a($x$)={^a_ij($x$)}_11^dd$ispositivedefiniteateach$$x$∈ℝ^d$.Onecouldreplacethenumbers$a_ij$of(\ref{exp2.13})bythenumbers$a_ij(h$k$)$.Theresulting$A_n(gen)=-A_n$wouldagainbeageneratorofaMJPin$G_n$.However,thematrices$A_n$thusobtainedarediscretizationsofthedifferentialoperator$A($x$)=-∑_ija_ij($x$)∂_i
∂_j$,aswillbeseeninSection\ref{sec5},andnotof$A($x$)=-∑_ij∂_i a_ij($x$)
∂_j$.Herewepresentamethodofconstructionofdiscretizationsofadifferentialoperatorindivergenceform,$A($x$)=-∑_ij∂_i
a_ij
($x$)∂_j$,resultinginmatricesofpositivetype,withastructuresimilarto(\ref{exp2.13}).Theoperatorindivergenceform,$A($x$)$,naturallycorrespondstothebilinearform$a(v,u)=∑_ij∫a_ij($x$)∂_iv($x$)
∂_ju($x$)d$x$$.Consequently,itsdiscretizations$A_n$willcorrespondtoasequenceofdiscretizedforms$a_n(v,u)$ongrids$G_n.

### Construction for $d=2$

Initially one considers the constant coefficients $a_{ij}$ and constructs the
forms $a_{n}(v,u)=\langle\,{\bf v}_{n}|A_{n}{\bf u}_{n}\,\rangle$ for which
the matrices $A_{n}$ coincide with (15). Then the obtained expression of
$a_{n}(v,u)$ is generalized to the case of non-constant coefficients $a_{ij}$.
The corresponding matrices $A_{n}$ are obtained from the variational method in
the standard way.

For the case of constant coefficients $a_{ij}$, the bilinear form
$h^{2}a_{n}(v,u)$ that discretizes the form
$a(v,u)=\sum_{i,j=1}^{2}a_{ij}\partial_{i}v\partial_{j}u$, will be a second
order polynomial in the quantities
$\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(r_{i}h)v(\mbox{\boldmath{$x$}})$,
$\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{j}(r_{j}h)u(\mbox{\boldmath{$x$}})$
with a certain choice of $r_{i},r_{j}\in{\hbox{\msbm N}}$ and
$\mbox{\boldmath{$x$}}\in G_{n}$. In order to write down the form as simply as
possible, we use the following abbreviations:

$\begin{array}[]{lllll}u_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})&=&\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(r_{i}h)u(h\mbox{\boldmath{$k$}})&=&(r_{i}h)^{-1}\,[u(h\mbox{\boldmath{$k$}}+r_{i}h\mbox{\boldmath{$e$}}_{i})-u(h\mbox{\boldmath{$k$}})],\\\
\widehat{u}_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})&=&\widehat{\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,}_{i}(r_{i}h)u(h\mbox{\boldmath{$k$}})&=&(r_{i}h)^{-1}\,[u(h\mbox{\boldmath{$k$}})-u(h\mbox{\boldmath{$k$}}-r_{i}h\mbox{\boldmath{$e$}}_{i})].\end{array}$

For $a_{12}<0$ the form is defined by:

$\begin{array}[]{ll}\displaystyle a_{n}^{(-)}(v,u)\ =\
\>\sum_{\mbox{\boldmath{\scriptsize$k$}}}&\displaystyle\left(\sum_{i=1}^{2}\>a_{ii}\,v_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)\,u_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)+\>\sum_{i\neq
j}\>a_{ij}v_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})\,u_{\mbox{\boldmath{\scriptsize$k$}}}(j,r_{j})\right.\\\
&+\left.\displaystyle\sum_{i\neq
j}\>a_{ij}\,\frac{r_{i}}{r_{j}}\>\big{[}v_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)\,u_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)-v_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})\,u_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})\big{]}\right).\end{array}$
(16)

If the last summand of (16) were omitted, then the resulting $A_{n}$ would
have neighborhoods ${\cal N}(\mbox{\boldmath{$x$}})$ containing seven grid-
knots $\mbox{\boldmath{$x$}},\mbox{\boldmath{$x$}}\pm
h\mbox{\boldmath{$e$}}_{i},\mbox{\boldmath{$x$}}\pm
h(r_{1}\mbox{\boldmath{$e$}}_{1}-r_{2}\mbox{\boldmath{$e$}}_{2})$ as
illustrated in Figure 1, as well as the following four additional grid-knots,
$\mbox{\boldmath{$x$}}\pm hr_{i}\mbox{\boldmath{$e$}}_{i}$. The second line of
(16) causes cancellation of those matrix entries which would correspond to
four additional grid-knots.

For $a_{12}>0$ the form is obtained from the previous one by changing
$v_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),u_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),i=2$,
into
$\widehat{v}_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),\widehat{u}_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),i=2$,
respectively.

For the case of functions $\mbox{\boldmath{$x$}}\mapsto
a_{ij}(\mbox{\boldmath{$x$}})$ one naturally starts from the just obtained
expression for $a_{n}^{(\pm)}(v,u)$, changing the numbers $a_{ij}$ into the
numbers $a_{ij}(h\mbox{\boldmath{$k$}})$. In fact, the numbers
$a_{ij}(h\mbox{\boldmath{$k$}}+\mbox{\boldmath{$t$}})$, where $t$ are
appropriately selected element of ${{\hbox{\msbm R}}}^{2}$, are acceptable.
The fastest convergence of $h_{n}^{2}a_{n}(v,u)\to a(v,u)$ is a criterion
which helps us to choose the vectors $t$. It turns out that the best choice is
$\mbox{\boldmath{$t$}}(\mbox{\boldmath{$r$}})=(h/2)(r_{1}\mbox{\boldmath{$e$}}_{1}+r_{2}\mbox{\boldmath{$e$}}_{2})$
[TS], i.e. the values for which
$h\mbox{\boldmath{$k$}}+\mbox{\boldmath{$t$}}(\mbox{\boldmath{$r$}})$ are the
mid-point of the rectangle
$C_{n}(\mbox{\boldmath{$r$}},h\mbox{\boldmath{$k$}})$, with the lower left
vertex at $h\mbox{\boldmath{$k$}}$ and the upper right vertex at
$h\mbox{\boldmath{$k$}}+r_{1}\mbox{\boldmath{$e$}}_{1}+r_{2}\mbox{\boldmath{$e$}}_{2}$.
Let us remark that the finite difference operators
$\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(r_{i}h)u(h\mbox{\boldmath{$k$}})$
are defined in terms of function values at the vertices of
$C_{n}(\mbox{\boldmath{$r$}},h\mbox{\boldmath{$k$}})$. Thus, if
$a_{12}(\mbox{\boldmath{$x$}})\leq 0,\mbox{\boldmath{$x$}}\in{{\hbox{\msbm
R}}}^{2}$, we get:

$\begin{array}[]{ll}a_{n}(v,u)\ =\ \sum_{\mbox{\boldmath{\scriptsize$x$}}\in
G_{n}}\>\Big{(}&\displaystyle\sum_{i=1}^{2}\>a_{ii}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}({\bf
1}))\,v_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)\,u_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)\\\
\displaystyle&+\>\sum_{i\neq
j}\>a_{ij}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}(\mbox{\boldmath{$r$}}))\,v_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})\,u_{\mbox{\boldmath{\scriptsize$k$}}}(j,r_{j})\\\
&\displaystyle+\>\sum_{i\neq
j}\>a_{ij}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}(\mbox{\boldmath{$r$}}))\,\frac{r_{i}}{r_{j}}\>\big{[}v_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)\,u_{\mbox{\boldmath{\scriptsize$k$}}}(i,1)-v_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})\,u_{\mbox{\boldmath{\scriptsize$k$}}}(i,r_{i})\big{]}\Big{)},\end{array}$

where $\mbox{\boldmath{$r$}}=(r_{1},r_{2}),{\bf 1}=(1,1)$. By using the
variational method, one obtains the entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}},\mbox{\boldmath{$k$}},\mbox{\boldmath{$l$}}\in
I_{n}$. At the present step of construction it is not necessary to write down
all the entries. In order to describe the influence of the parameters
$r_{1},r_{2}$ on the structure of entries we consider one group of entries:

$\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\,\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$e$}}_{1}}\
=\
-\,\frac{1}{h^{2}}\>\Big{[}a_{11}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}({\bf
1}))\>-\>\frac{r_{1}}{r_{2}}\,\big{|}a_{12}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}(\mbox{\boldmath{$r$}}))\big{|}\Big{]}.$

The matrices $A_{n}$ are of positive type iff the bracket has positive sign
for each $\mbox{\boldmath{$x$}}\in G_{n}$. This is a condition on the
functions $a_{ij}$, which is implied by a particular choice of
$r_{1},r_{2}\in{\hbox{\msbm N}}$.

For $a_{12}\geq 0$ the form is obtained from the constructed one by changing
the following quantities. The finite differences
$v_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),u_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),i=2$,
should be changed into
$\widehat{v}_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),\widehat{u}_{\mbox{\boldmath{\scriptsize$k$}}}(i,\cdot),i=2$,
as in the case of constant coefficients, and
$\mbox{\boldmath{$t$}}(\mbox{\boldmath{$m$}})$ should be changed into
$\mbox{\boldmath{$s$}}(\mbox{\boldmath{$m$}})=(h/2)(m_{1}\mbox{\boldmath{$e$}}_{1}-m_{2}\mbox{\boldmath{$e$}}_{2})$.
With these changes, the obtained $A_{n}$ are symmetric matrices, possibly of
positive type. We say “possibly of positive type” since this depends on the
choice of $r_{1},r_{2}$.

In the general case, the sign of $a_{12}$ is not constant on ${{\hbox{\msbm
R}}}^{2}$. Therefore, we partition ${{\hbox{\msbm R}}}^{2}$ into two subsets
$\\{\mbox{\boldmath{$x$}}\in{{\hbox{\msbm
R}}}^{2}:a_{12}(\mbox{\boldmath{$x$}})\leq 0\\}$ and
$\\{\mbox{\boldmath{$x$}}\in{{\hbox{\msbm
R}}}^{2}:a_{12}(\mbox{\boldmath{$x$}})>0\\}$. Then each of these sets has to
be partitioned further, where each of the partitioned classes is characterized
by a pair $r_{1},r_{2}$, so that the resulting entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}$
have negative values. The construction is carried out for a class of functions
$a_{ij}$ with moderate discontinuities.

###### Assumption 3.1

Let there exist a finite index set $\mathcal{L}$, a partition ${{\hbox{\msbm
R}}}^{2}=\cup_{l\in\mathcal{L}}D_{l}$ and a diffusion tensor
$a=\\{a_{ij}\\}_{11}^{22}$ on ${{\hbox{\msbm R}}}^{2}$ satisfying the strict
ellipticity conditions (5) and the following additional discretization
conditions:

a) There exist $\mbox{\boldmath{$q$}}\in{{\hbox{\msbm N}}}^{2}$,
$n_{0}\in{{\hbox{\msbm N}}}$ and the corresponding $h_{0}=2^{-n_{0}}$ so that
each set $D_{l}$ is a connected union of cubes of the form
$\prod_{i=1}^{d}[x_{i},x_{i}+q_{i}h_{0})$. The matrix-valued function
$\mbox{\boldmath{$x$}}\mapsto a(\mbox{\boldmath{$x$}})$ is continuous on ${\rm
cl}(D_{l})$ and either $a_{12}\geq 0$ or $a_{12}\leq 0$ on $D_{l}$.

b) To each $D_{l}$ there is associated a parameter
$\mbox{\boldmath{$r$}}(l)\in{{\hbox{\msbm N}}}^{2}$, such that the following
inequality is valid:

$\begin{array}[]{ll}\displaystyle\omega(a)\
=&\displaystyle\inf_{n}\>\min_{l\in{\cal
L}}\,\min_{i=1,2}\,\inf\Big{\\{}\,\inf_{\mbox{\boldmath{\scriptsize$z$}}\in
S_{n}(\mbox{\boldmath{\scriptsize$r$}}(l),\mbox{\boldmath{\scriptsize$x$}})\cap
D_{l}}\,a_{ii}(\mbox{\boldmath{$z$}})\\\
&\displaystyle-\>\frac{r_{i}(l)}{r_{m(i)}(l)}\,\sup_{\mbox{\boldmath{\scriptsize$z$}}\in
S_{n}(\mbox{\boldmath{\scriptsize$r$}}(l),\mbox{\boldmath{\scriptsize$x$}})\cap
D_{l}}\,|a_{im(i)}(\mbox{\boldmath{$z$}})|\>:\>\mbox{\boldmath{$x$}}\in
G_{n}\Big{\\}}\ >0.\end{array}$ (17)

where $m(i)=3-i$ and the rectangles
$S_{n}(\mbox{\boldmath{$r$}}(l),\mbox{\boldmath{$x$}})$ are defined by (13).

Condition b) is crucial in our construction of discretizations $A_{n}$ with
the structure of matrices of positive type.

The set $\mathcal{L}$ of Assumption 3.1 is partitioned into the subsets
$\mathcal{L}_{\mp}$, where $l\in\mathcal{L}_{-}$ means that $a_{12}\leq 0$ on
$D_{l}$, and $l\in\mathcal{L}_{+}$ means that $a_{12}\geq 0$ on $D_{l}$. It is
convenient to use a representation
$a_{n}(v,u)=a_{n}^{(-)}(v,u)+a_{n}^{(+)}(v,u)$, where the form
$a_{n}^{(-)}(v,u)$ contains the sums over the grid-knots in
$D_{l},l\in\mathcal{L}_{-}$ and $a_{n}^{(+)}(v,u)$ contains the sums over the
grid-knots in $D_{l},l\in\mathcal{L}_{+}$. Let us define

$\begin{array}[]{l}\mbox{\boldmath{$t$}}^{(\pm+)}(\mbox{\boldmath{$r$}})\ =\
\frac{h}{2}\,\big{(}\pm\,r_{1}(l)\mbox{\boldmath{$e$}}_{1}\,+\,r_{2}(l)\mbox{\boldmath{$e$}}_{2}\big{)}\
\in S_{n}(\mbox{\boldmath{$r$}},{\bf 0}),\\\
\mbox{\boldmath{$t$}}^{(\pm-)}(\mbox{\boldmath{$r$}})\ =\
\frac{h}{2}\,\big{(}\pm\,r_{1}(l)\mbox{\boldmath{$e$}}_{1}\,-\,r_{2}(l)\mbox{\boldmath{$e$}}_{2}\big{)}\
\in S_{n}(\mbox{\boldmath{$r$}},{\bf 0}).\end{array}$ (18)

where $\mbox{\boldmath{$r$}}=(r_{1},r_{2})$. The form
$a_{n}^{(-)}(\cdot,\cdot)$ is defined by:

$\begin{array}[]{ll}a_{n}^{(-)}(v,u)=\
\sum_{l\in\mathcal{L}_{-}}&\displaystyle\Big{\\{}\sum_{\mbox{\boldmath{\scriptsize$x$}}\in
G_{n}(l)}\>\Big{(}\sum_{i=1}^{2}\>a_{ii}(l,\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}^{(++)}({\bf
1}))\,\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(h)v\big{)}(\mbox{\boldmath{$x$}})\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(h)u\big{)}(\mbox{\boldmath{$x$}})\\\
&\displaystyle+\>\sum_{i,j=1,2,i\neq
j}\>a_{ij}(l,\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}^{(++)}(\mbox{\boldmath{$r$}}))\,\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(r_{i}(l)h)v\big{)}(\mbox{\boldmath{$x$}})\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{j}(r_{j}(l)h)u\big{)}(\mbox{\boldmath{$x$}})\\\
&\displaystyle+\>\sum_{i,j=1,2,i\neq
j}\>a_{ij}(l,\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}^{(++)}(\mbox{\boldmath{$r$}}))\,\frac{r_{i}(l)}{r_{j}(l)}\>\Big{[}\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(h)v\big{)}(\mbox{\boldmath{$x$}})\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(h)u\big{)}(\mbox{\boldmath{$x$}})\\\
&\>-\>\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(r_{i}(l)h)v\big{)}(\mbox{\boldmath{$x$}})\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(r_{i}(l)h)u\big{)}(\mbox{\boldmath{$x$}})\Big{]}\
\Big{)}\Big{\\}},\end{array}$ (19)

where, as already noted, $\mbox{\boldmath{$r$}}=(r_{1},r_{2})$, so that ${\bf
1}=(1,1)$. The form $a_{n}^{(+)}(v,u)$ is obtained from $a_{n}^{(-)}(v,u)$
formally by replacing $\mathcal{L}_{-}$ with $\mathcal{L}_{+}$, then
$\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,(mh)f$
with
$\widehat{\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,}(mh)f$
for each $m\in\\{1,r_{2}\\}$ and $f=v,u$, and $\mbox{\boldmath{$t$}}^{(++)}$
with $\mbox{\boldmath{$t$}}^{(+-)}$. Observe that the forms
$a_{n}^{(\mp)}(v,u)$ are the second degree polynomials in the quantities
$\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(qh),\widehat{\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,}_{i}(qh)$
with $q\in\\{1,r_{1},r_{2}\\}$.

We say that $\mbox{\boldmath{$x$}}\in G_{n}$ is an internal grid-knot if
${\cal N}(\mbox{\boldmath{$x$}})\subset G_{n}\cap D_{l}$ for some $l$. All the
other grid-knots are called boundary grid-knots. For an internal grid-knot $x$
the acceptable expressions for
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$
follow directly from the definition of corresponding discrete forms ${\bf
v},{\bf u}\mapsto\langle\,{\bf v}|A_{n}{\bf u}\,\rangle$. For a boundary grid-
knot $\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}$ the obtained expressions
are complex, and the calculated
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$
could break down the structure of matrices of positive type. Therefore, one
seeks simpler procedures for constructing entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$
for boundary grid-knots $\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}$. The
results of such a construction must be matrices $A_{n}$ which determine MJPs
and the convergence (3) of MJPs determined by $A_{n}$ should also be ensured.

In order to write down the entries of $A_{n}$ corresponding to internal grid-
knots, we use (18) and the following abbreviations:

$\begin{array}[]{lll}\mbox{\boldmath{$w$}}^{(\pm)}(l)&=&r_{1}(l)\mbox{\boldmath{$e$}}_{1}\,\pm\,r_{2}(l)\mbox{\boldmath{$e$}}_{2}\
\in I_{n},\\\
a_{ij}^{(\alpha\beta)}(\mbox{\boldmath{$r$}})&=&a_{ij}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}^{(\alpha\beta)}(\mbox{\boldmath{$r$}})),\quad\alpha,\beta\in\\{+,-\\},\\\
\check{a}_{12}^{(-+)}(\mbox{\boldmath{$r$}})&=&a_{12}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}^{(++)}(\mbox{\boldmath{$r$}})-h\mbox{\boldmath{$e$}}_{1}),\\\
\check{a}_{12}^{(+-)}(\mbox{\boldmath{$r$}})&=&a_{12}(\mbox{\boldmath{$x$}}+\mbox{\boldmath{$t$}}^{(++)}(\mbox{\boldmath{$r$}})-h\mbox{\boldmath{$e$}}_{2}),\\\
\check{a}_{ii}^{(++)}(\mbox{\boldmath{$r$}})&=&a_{ii}^{(++)}(\mbox{\boldmath{$r$}}),\quad\check{a}_{ii}^{(--)}(\mbox{\boldmath{$r$}})\
=\ a_{ii}^{(--)}(\mbox{\boldmath{$r$}}).\end{array}$

When we apply variational method to the form $a_{n}^{(-)}$ defined by (19) and
the corresponding form $a_{n}^{(+)}$, we obtain entries of $A_{n}$. Thus, for
an internal grid-knot $h\mbox{\boldmath{$k$}}$ the nontrivial off-diagonal
entries of $A_{n}$ are:

$\begin{array}[]{ll}\displaystyle\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{1}}&\displaystyle=\
-\,\frac{1}{h^{2}}\>\left\\{\begin{array}[]{lll}\displaystyle
a_{11}^{(\pm+)}({\bf
1})\>-\>\frac{r_{1}(l)}{r_{2}(l)}\,\big{|}\check{a}_{12}^{(\pm+)}(\mbox{\boldmath{$r$}})\big{|}&{\rm
for}&a_{12}\leq 0,\ {\rm on}\ D_{l},\\\ \displaystyle a_{11}^{(\pm-)}({\bf
1})\>-\>\frac{r_{1}(l)}{r_{2}(l)}\,\big{|}\check{a}_{12}^{(\pm-)}(\mbox{\boldmath{$r$}})\big{|}&{\rm
for}&a_{12}\geq 0,\ {\rm on}\ D_{l},\end{array}\right.\\\
\displaystyle\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{2}}&\displaystyle=\
-\,\frac{1}{h^{2}}\>\Big{[}a_{22}^{(+\pm)}({\bf
1})\>-\>\frac{r_{2}(l)}{r_{1}(l)}\,\big{|}\check{a}_{12}^{(+\pm)}(\mbox{\boldmath{$r$}})\big{|}\Big{]}.\end{array}$
(20)

The entries corresponding to the grid-knots in the plane spanned by
$\mbox{\boldmath{$e$}}_{1},\mbox{\boldmath{$e$}}_{2}$ have the structure:

$\begin{array}[]{l}\begin{array}[]{l}\displaystyle\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$w$}}^{(-)}(l)}\
=\
-\,\frac{1}{h^{2}r_{1}(l)r_{2}(l)}\>\big{|}a_{12}^{(+-)}(\mbox{\boldmath{$r$}})\big{|}\\\
\displaystyle\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}-\mbox{\boldmath{\scriptsize$w$}}^{(-)}(l)}\
=\
-\,\frac{1}{h^{2}r_{1}(l)r_{2}(l)}\>\big{|}a_{12}^{(-+)}(\mbox{\boldmath{$r$}})\big{|}\end{array}\quad{\rm
for}\quad a_{12}\leq 0\ {\rm on}\ D_{l},\\\
\begin{array}[]{l}\displaystyle\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$w$}}^{(+)}(l)}\
=\
-\,\frac{1}{h^{2}r_{1}(l)r_{2}(l)}\>\big{|}a_{12}^{(++)}(\mbox{\boldmath{$r$}})\big{|}\\\
\displaystyle\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}-\mbox{\boldmath{\scriptsize$w$}}^{(+)}(l)}\
=\
-\,\frac{1}{h^{2}r_{1}(l)r_{2}(l)}\>\big{|}a_{12}^{(--)}(\mbox{\boldmath{$r$}})\big{|}\\\
\end{array}\quad{\rm for}\quad a_{12}\geq 0\ {\rm on}\ D_{l},\end{array}$ (21)

If the quantities
$\check{a}_{ij}^{(\alpha\beta)}(\mbox{\boldmath{$r$}}),\alpha,\beta\in\\{+,-\\}$
in (20) are replaced with $a_{ij}^{(\alpha\beta)}(\mbox{\boldmath{$r$}})$, the
convergence is still preserved. However, the quantities
$a_{ij}^{(\alpha\beta)}(\mbox{\boldmath{$r$}}),\check{a}_{ij}^{(\alpha\beta)}(\mbox{\boldmath{$r$}})$
should not be replaced with $a_{ij}(h\mbox{\boldmath{$k$}})$ since the
resulting
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$
would be discretizations of $-\sum_{ij}a_{ij}\partial_{i}\partial_{j}$ and not
of $-\sum_{ij}\partial_{i}a_{ij}\partial_{j}$. This assertion can be easily
proved for the case of dimension $d=1$ and the diffusion tensor $a$, that is
twice continuously differentiable. We intend to compare the expressions
$A(x)u(x)=-(a(x)u(x)^{\prime})^{\prime}$ and
$A^{cs}(x)u(x)=-a(x)u^{\prime\prime}(x)$ and their discretizations on grids
$G_{n}=\\{hk:k\in{\hbox{\msbm Z}}\\}\subset{\hbox{\msbm R}}$. The
discretizations of $A^{cs}(x)u(x)$ are given by Lemma 3.1:

$\big{(}A^{cs}\big{)}_{kk\pm
1}\,=\,-h^{-2}a(hk),\qquad\big{(}A^{cs}\big{)}_{kk}\,=\,2h^{-2}a(hk),$

and consequently:

$\big{(}(A^{cs})_{n}{\bf u}_{n}\big{)}_{k}\ =\
a(hk)\,\frac{2u(hk)-u(hk+h)-u(hk-h)}{h^{2}}.$

For the discretizations of $A(x)u(x)$ we consider (20). Let us represent
$a(x\pm h/2)$ by its Taylor polynomial of the second degree, $a(x\pm
h/2)=a(x)\pm(h/2)a^{\prime}(x)+(h/2)^{2}a^{\prime\prime}(x)+r(\pm h,x)$, where
the remainder $r(\pm h,x)$ has the property $\lim_{h\to 0}h^{-2}r(\pm h,x)=0$.
Therefore we get:

$\big{(}A\big{)}_{kk\pm
1}\,=\,-\,\Big{(}1+\frac{h^{2}}{4}\frac{a^{\prime\prime}(hk)}{a(hk)}\Big{)}\,\frac{1}{h^{2}}\,a(hk)\,\mp\frac{1}{2h}a^{\prime}(hk)+o(\pm
h,x),$

where $\lim_{h\to 0}o(\pm h,x)=0$. Hence:

$\begin{array}[]{lll}\big{(}A_{n}{\bf
u}_{n}\big{)}_{k}&=&\displaystyle\Big{(}1+\frac{h^{2}}{4}\frac{a^{\prime\prime}(hk)}{a(hk)}\Big{)}\,a(hk)\,\frac{2u(hk)-u(hk+h)-u(hk-h)}{h^{2}}\\\
&+&\displaystyle
a^{\prime}(hk)\frac{u(hk+h)-u(hk-h)}{2h}\,+\,\tilde{o}(h,x).\end{array}$

From the obtained expressions we have:

$\big{(}A_{n}{\bf
u}_{n}\big{)}_{k}\,-\,\Big{(}1+\frac{h^{2}}{4}\frac{a^{\prime\prime}(hk)}{a(hk)}\Big{)}\big{(}(A^{cs})_{n}{\bf
u}_{n}\big{)}_{k}\ =\ \big{(}B_{n}{\bf u}_{n}\big{)}_{k}\,+\,\tilde{o}(h,x),$

where $\big{(}B_{n}{\bf u}_{n}\big{)}_{k}$ are the discretizations of
$-a^{\prime}(x)u^{\prime}(x)$. In other words, we have obtained
discretizations of the expression $A(x)=A^{cs}(x)-a^{\prime}(x)(d/dx)$ as we
should have.

Now we can describe the construction which gives a satisfactory result for
both types of grid-knots, internal and boundary ones. For each
$\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}\in{\rm cls}(D_{l})$ the entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$
are constructed by the rules (20), (21). If ${\cal
N}(\mbox{\boldmath{$x$}})\subset{\rm cls}(D_{l})$ there is nothing more to
adjust. If ${\cal N}(\mbox{\boldmath{$x$}})\cap{\rm
cls}(D_{l})\subseteq/\hskip 5.69054pt{\cal N}(\mbox{\boldmath{$x$}})$ then the
quantities
$a_{12}^{(\alpha\beta)}(\mbox{\boldmath{$r$}}),\check{a}_{12}^{(\alpha\beta)}(\mbox{\boldmath{$r$}})$,
where $\alpha,\beta\in\\{+,-\\}$, should be replaced by zeros in all the
expressions (20), (21). Let us point out that otherwise the entries (20) for
the case of ${\cal N}(\mbox{\boldmath{$x$}})\cap{\rm
cls}(D_{l})\subseteq/\hskip 5.69054pt{\cal N}(\mbox{\boldmath{$x$}})$ could
have a wrong sign, i.e. it could happen that
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}>0$
for some $\mbox{\boldmath{$k$}}\neq\mbox{\boldmath{$l$}}$. In order to avoid
such undesired features, one therefore omits the terms of entries in the
expressions (20), (21), which are proportional to $|a_{12}|$. This adjustment
procedure is equivalent to the assumption that the function $a_{12}$ is zero
in a neighborhood of the boundary $\Gamma=\cup_{l}\partial D_{l}$.

This determines the rules of construction of discretizations of a generalized
diffusion for which the diffusion tensor satisfies Assumption 3.1. The above
described construction of entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$
at boundary grid-knots can be justified as follows.

A numerical neighborhood ${\cal N}(\mbox{\boldmath{$x$}})$, where
$\mbox{\boldmath{$x$}}\in D_{l}\cap G_{n},l\in{\cal L}$ is contained in the
closed rectangle

$S_{n}(\mbox{\boldmath{$r$}}(l),\mbox{\boldmath{$x$}})\ =\
\prod_{i=1}^{d}\>[x_{i}-h_{n}r_{i}(l)\,,\,x_{i}+h_{n}r_{i}(l)].$ (22)

The union of all such rectangles centered at boundary grid-knots
$\mbox{\boldmath{$x$}}\in\Gamma\cap G_{n}$ is denoted by $S_{n}(\Gamma)$. This
is a closed set covering $\Gamma$. Because of $G_{n}\subset G_{n+1}$ we have
$S_{n+1}(\Gamma)\subset S_{n}(\Gamma)$ and the identity
$\Gamma=\cap_{n}S_{n}(\Gamma)$. For each $n\in{\hbox{\msbm N}}$, we
approximate the original diffusion tensor $a$ by the tensor $a^{(n)}$ defined
as follows:

$\begin{array}[]{lllll}a_{ij}^{(n)}(\mbox{\boldmath{$x$}})&=&a_{ij}(\mbox{\boldmath{$x$}})&{\rm
for}&\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{2}\setminus
S_{n}(\Gamma),\quad i,j\in\\{1,2\\},\\\
a_{ii}^{(n)}(\mbox{\boldmath{$x$}})&=&a_{ii}(\mbox{\boldmath{$x$}})&{\rm
for}&\mbox{\boldmath{$x$}}\in S_{n}(\Gamma),\quad i\in\\{1,2\\},\\\
a_{12}^{(n)}(\mbox{\boldmath{$x$}})&=&0&{\rm for}&\mbox{\boldmath{$x$}}\in
S_{n}(\Gamma).\end{array}$ (23)

The defined diffusion tensors $a^{(n)}$ determine differential operators
$A_{0}^{(n)}(\mbox{\boldmath{$x$}})$ which approximate the differential
operator $A_{0}(\mbox{\boldmath{$x$}})$. Each
$A_{0}^{(n)}(\mbox{\boldmath{$x$}})$ determines a generalized diffusion
$X^{(n)}(\cdot)$. The generalized diffusions $X^{(n)}(\cdot)$ converge in
distribution to the original diffusion $X(\cdot)$. Therefore it suffices to
consider the discretizations $A_{n}^{(n)}$ which are defined in terms od
$a_{ij}^{(n)}$ as described above. The convergence in distribution of
diffusions $X^{(n)}(\cdot)$ to $X(\cdot)$ follows from the convergence of the
corresponding semigroups $U^{(n)}(\cdot)$ to the semigroup $U(\cdot)$ in
$\dot{C}({{\hbox{\msbm R}}}^{d})$. To see this, we first apply Corollary 2.1
in order to prove the uniform boundedness $\sup_{t\leq
1}\|U^{(n)}(t)\|^{(\alpha)}\leq\beta$, where $\alpha,\beta$ do not depend on
$n$. Then we use the standard methods of Sobolev spaces in order to prove the
convergence of the semigroups in $L_{2}({{\hbox{\msbm R}}}^{d})$. The two
above properties are combined in order to prove the convergence of semigroups
in $\dot{C}({{\hbox{\msbm R}}}^{d})$ as in the last step of proof of Theorem
4.1.

###  3.1 Construction for $d\geq 3$

The goal of the overall analysis is to find those discretizations $A_{n}$ of
the differential operator $A(\mbox{\boldmath{$x$}})$ which have the structure
of matrices of positive type. Here we describe a general approach, which is
based on reduction to a finite number of two-dimensional problems.

The index set of pairs $I(d)=\\{\\{ij\\}:i<j,\>i,j=1,2,\ldots,d,i\neq j\\}$
has the cardinal number $m(d)=d(d-1)/2$. To each index $\\{kl\\}\in I(d)$ we
associate three coefficients,

$a_{kk}^{\\{kl\\}}=\frac{1}{d-1}\,a_{kk},\quad
a_{ll}^{\\{kl\\}}=\frac{1}{d-1}\,a_{ll},\quad a_{kl}^{\\{kl\\}}=a_{kl},$ (24)

and a two-dimensional bilinear form $a^{\\{kl\\}}(\cdot,\cdot)$,

$a^{\\{kl\\}}(v,u)\ =\sum_{i,j\in\\{k,l\\}}\>\int_{{{\hbox{\msbm
R}}}^{d}}\>a_{ij}^{\\{kl\\}}(\mbox{\boldmath{$x$}})\>\partial_{i}v(\mbox{\boldmath{$x$}})\,\partial_{j}u(\mbox{\boldmath{$x$}})\,d\mbox{\boldmath{$x$}}.$

Clearly, for each pair $v,u\in C_{0}^{(1)}({{\hbox{\msbm R}}}^{d})$ the
following identity is valid:

$a(v,u)\ =\ \sum_{\\{kl\\}\in I(d)}\,a^{\\{kl\\}}(v,u).$ (25)

To each of the forms $a^{\\{kl\\}}(\cdot,\cdot)$ we associate a sequence of
forms $a_{n}^{\\{kl\\}}(\cdot,\cdot)$ and matrices $A_{n}^{\\{kl\\}}$
constructed by using schemes in two dimensions from the previous subsection.
Then the matrix

$A_{n}\ =\ \sum_{\\{kl\\}\in I}\>A_{n}^{\\{kl\\}},$ (26)

is a discretization of $A_{0}(\mbox{\boldmath{$x$}})$. If each
$A_{n}^{\\{kl\\}}$ has the structure of a matrix of positive type, then
$A_{n}$ is also a matrix of positive type. However, $A_{n}$ can have the
structure of a matrix of positive type, even though no $A_{n}^{\\{kl\\}}$ is a
matrix of positive type. This important property, which enables us to
construct matrices $A_{n}$ of positive type for the case $d\geq 3$, can be
proved from the structure of entries
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}$.

First we consider the entry
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$e$}}_{1}}$
defined by (26) in the case that $a_{ij}<0$ for all $i\neq j,\
i,j=1,2,\ldots,d$. In addition, in order to write down expressions as simply
as possible, the index $l\in{\cal L}$ is omitted from the notations. The
contribution from the sum of entries
$\big{(}A_{n}^{\\{kl\\}}\big{)}{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$e$}}_{1}}$
to the entry
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$e$}}_{1}}$
has the following form:

$\frac{1}{d-1}\sum_{s\geq
2}a_{11}\Big{(}\mbox{\boldmath{$x$}}+\frac{h}{2}(\mbox{\boldmath{$e$}}_{1}+\mbox{\boldmath{$e$}}_{s})\Big{)}\
-\ {\rm terms~{}containing~{}~{}}a_{12},a_{13},\ldots,a_{1d}.$

Similarly we can describe the terms containing $a_{ii}$ for any
$i=1,2,\ldots,d$, and any $l\in{\cal L}$. The corresponding sum contributing
to
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}+\mbox{\boldmath{\scriptsize$e$}}_{i}}$
has the following general form:

$\omega_{n}(a_{ii},\mbox{\boldmath{$x$}})\ =\ \frac{1}{d-1}\>\sum_{s=1,s\neq
i}^{d}\>a_{ii}(h\mbox{\boldmath{$k$}}+h\mbox{\boldmath{$m$}}_{ii}(l,s)),$ (27)

where $\mbox{\boldmath{$m$}}_{ii}(l,s)$ are defined by the rules of
construction of (26), (20) and (21). The terms proportional to $a_{is},s\neq
i$, are summed with the just defined
$\omega_{n}(a_{ii},\mbox{\boldmath{$x$}})$ as shown in the next description of
the obtained results.

###### Discretization procedure 3.1

Let there be given a partition ${{\hbox{\msbm R}}}^{d}=\cup_{l}D_{l}$ into a
finite number of subsets $D_{l}$ such that all the functions $a_{ij}$ are
uniformly continuous on each $D_{l}$, and the functions $a_{ij},i\neq
j,i,j=1,2,\ldots,d$, do not change sign on $D_{l}$. Let a parameter
$\mbox{\boldmath{$r$}}(l)=(r_{1}(l),r_{2}(l),\ldots,r_{d}(l))\in{{\hbox{\msbm
N}}}^{d}$ be assigned to each $D_{l}$ and the matrices $A_{n}$ on $G_{n}$ be
constructed by the rules (20), (21) and (26). Then their entries have the
following properties:

  1. 1.

Entries of
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}},\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}$,
$\mbox{\boldmath{$k$}}\in{{\hbox{\msbm Z}}}^{d}$,
$h\mbox{\boldmath{$l$}}\in{\cal N}(\mbox{\boldmath{$x$}})$, are linear
combinations of
$a_{ij}(\mbox{\boldmath{$x$}}_{ij}(n,\mbox{\boldmath{$x$}},l))$ where
$\mbox{\boldmath{$x$}}_{ij}(n,\mbox{\boldmath{$x$}},l)=h\mbox{\boldmath{$k$}}+h\mbox{\boldmath{$m$}}_{ij}(l,s)$,
and where the elements $\mbox{\boldmath{$m$}}_{ij}(l,s)\in{{\hbox{\msbm
R}}}^{d}$ for $i,j,s=1,2,\ldots,d,l\in\mathcal{L}$, do not depend on $n$.

  2. 2.

For each grid-knot $\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}$:
$(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}}=-\sum_{\mbox{\boldmath{\scriptsize$l$}}}(A_{n})_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$.

  3. 3.

For each $\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}\in{\rm cls}(D_{l})$
entries in the coordinate directions $\mbox{\boldmath{$x$}}\pm
h\mbox{\boldmath{$e$}}_{i}$ are defined by:

$\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}\>=\>-\,\frac{1}{h^{2}}\,\Big{[}\omega_{n}(a_{ii},\mbox{\boldmath{$x$}})-\,\sum_{m=1,m\neq
i}^{d}\>\frac{r_{i}(l)}{r_{m}(l)}\,|a_{im}(\mbox{\boldmath{$x$}}_{im}(n,\mbox{\boldmath{$x$}},l))|\Big{]}.$

  4. 4.

For each $l\in\mathcal{L}$ the entries of $A_{n}$ in the plane spanned by
$\mbox{\boldmath{$e$}}_{i},\mbox{\boldmath{$e$}}_{j}$ are defined by using
elements
$\mbox{\boldmath{$z$}}_{ij}(l)=r_{i}(l)\mbox{\boldmath{$e$}}_{i}-r_{j}(l)\mbox{\boldmath{$e$}}_{j}\in{{\hbox{\msbm
Z}}}^{d}$ (if $a_{ij}\leq 0$ on $D_{l}$) or elements
$\mbox{\boldmath{$z$}}_{ij}(l)=r_{i}(l)\mbox{\boldmath{$e$}}_{i}+r_{j}(l)\mbox{\boldmath{$e$}}_{j}\in{{\hbox{\msbm
Z}}}^{d}$, (if $a_{ij}\geq 0$ on $D_{l}$), as follows:

$\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$z$}}_{ij}(l)}\
=\
-\>\frac{1}{h^{2}r_{i}(l)r_{j}(l)}\>|a_{ij}(\mbox{\boldmath{$x$}}_{ij}(n,\mbox{\boldmath{$x$}},l))|.$

An appropriate choice of the parameters $\mbox{\boldmath{$r$}}(l)$ follows
from Theorem 3.1.

Some special features regarding the structure of the sets ${\cal
N}(\mbox{\boldmath{$x$}}),\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}\in
G_{n}$, should be pointed out. If $a_{ij}\neq 0,i\neq j$, then the maximal
number of elements in ${\cal N}(\mbox{\boldmath{$x$}})$ is $1+d+d^{2}$. In
this case the set ${\cal N}(\mbox{\boldmath{$x$}})$ consists of its center,
$2d$-grid-knots in the coordinate directions, and 2 grid-knots in each of
$d(d-1)$ two-dimensional planes. Since there can be at most two grid-knots in
a two-dimensional plane (Property 4. of Discretization procedure 3.1), the
entries of $A_{n}^{(rs)}$ have the following property. Let the pairs
$\mbox{\boldmath{$e$}}_{r},\mbox{\boldmath{$e$}}_{s}$ and
$\mbox{\boldmath{$e$}}_{s},\mbox{\boldmath{$e$}}_{t}$ span two-dimensional
planes and let $A_{n}^{(rs)},A_{n}^{(st)}$ be the corresponding
discretizations which are constructed by using parameters
$\mbox{\boldmath{$r$}}^{(rs)},\mbox{\boldmath{$r$}}^{(st)}$. Then for the
construction defined by Discretization procedure 3.1 the following identity
$(\mbox{\boldmath{$r$}}^{(rs)})_{s}=(\mbox{\boldmath{$r$}}^{(st)})_{s}$ must
be valid.

###  3.2 Lower order differential operators

The discretizations of differential operator
$B(\mbox{\boldmath{$x$}})=\mbox{\boldmath{$b$}}(\mbox{\boldmath{$x$}})\nabla$
are denoted by $B_{n}$. The following general rule should be obeyed. A
positive diagonal entry and a non-positive off diagonal entry is associated to
each $\mbox{\boldmath{$x$}}\in G_{n}$ for which
$\mbox{\boldmath{$b$}}(\mbox{\boldmath{$x$}})\neq{\bf 0}$. Let us define the
sets $\mathcal{K}(i,-)=\\{\mbox{\boldmath{$x$}}\in
G_{n}:b_{i}(\mbox{\boldmath{$x$}})<0\\}$ and analogously
$\mathcal{K}(i,+)=\\{\mbox{\boldmath{$x$}}\in
G_{n}:b_{i}(\mbox{\boldmath{$x$}})>0\\}$. Then the discretizations of $(v|Bu)$
are defined by

$b_{n}(v,u)\ =\
\sum_{i}\>\Big{[}\sum_{\mbox{\boldmath{\scriptsize$x$}}\in\mathcal{K}(i,-)}\>b_{i}(\mbox{\boldmath{$x$}})\,v(\mbox{\boldmath{$x$}})\,\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(h)u\big{)}(\mbox{\boldmath{$x$}})\,+\,\sum_{\mbox{\boldmath{\scriptsize$x$}}\in\mathcal{K}(i,+)}\>b_{i}(\mbox{\boldmath{$x$}})\,v(\mbox{\boldmath{$x$}})\,\big{(}\,\rule[0.0pt]{0.28453pt}{3.69885pt}\rule[0.0pt]{3.69885pt}{0.28453pt}\rule[3.69885pt]{3.69885pt}{0.28453pt}\rule[6.54413pt]{3.69885pt}{0.28453pt}\rule[0.0pt]{0.28453pt}{6.54413pt}\,_{i}(-h)u\big{)}(\mbox{\boldmath{$x$}})\Big{]}.$

These forms have to be summed with the forms (25) in order to get
discretizations of the original form (6). If discretizations $(A_{0})_{n}$ of
$A_{0}(\mbox{\boldmath{$x$}})$ have the structure of matrices of positive
type, then obviously $(A_{0})_{n}+B_{n}$ maintain this structure. The so
defined discretizations of $B$ are usually called upwind schemes.

The constructed forms $a_{n}$ of this section are discretizations of the form
(6). At the present level of analysis the constructed discretizations can be
justified by the limit $a(v,u)=\lim_{n}h^{d}a_{n}(v,u)$, being valid for any
pair $v,u\in C_{0}^{(1)}({{\hbox{\msbm R}}}^{d})$.

###  3.3 Summarized results of the construction

###### THEOREM 3.1

Let there be given a partition ${{\hbox{\msbm R}}}^{d}=\cup_{l}D_{l}$ into a
finite number of connected sets $D_{l}$, each being the union of cubes
$C_{m}(\mbox{\boldmath{$p$}},\mbox{\boldmath{$x$}})$ with some fixed $m$, so
that the functions $a_{ij}$ fulfil (5) and the following additional
conditions:

a) The functions $a_{ij}$ are uniformly continuous on $D_{l}$ and
$a_{ij},i\neq j$, do not change sign on $D_{l}$.

b) For each pair $i,j$ the limit
$\lim_{|\mbox{\boldmath{\scriptsize$x$}}|\to\infty}a_{ij}(\mbox{\boldmath{$x$}})$
has a constant value.

c) The matrix-valued function
$\mbox{\boldmath{$x$}}\mapsto\hat{a}(\mbox{\boldmath{$x$}})$, defined by (2),
is strictly positive definite on ${{\hbox{\msbm R}}}^{d}$, i.e.
$(\mbox{\boldmath{$z$}}|\hat{a}(\mbox{\boldmath{$x$}})\mbox{\boldmath{$z$}})\geq\beta|\mbox{\boldmath{$z$}}|^{2}$
for some $\beta>0$ and all
$\mbox{\boldmath{$z$}},\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{d}$.

Then there exist discretizations $A_{n}$ which are constructed by the rules of
Discretization procedure 3.1, such that $A_{n}$ are matrices of positive type.

Proof:  If for each $D_{l}$ we choose the parameters $r_{i}(l)$ so that the
entries
$\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}$
of item 3. of Discretization procedure 3.1 have all negative values, then the
rules of construction (21) ensure the existence of $A_{n}$ with the structure
of matrices of positive type. It remains to justify the existence of such
parameters $\mbox{\boldmath{$r$}}(l)\in{{\hbox{\msbm N}}}^{d}$ for each
$l\in{\cal L}$. Let us consider a set ${\rm cls}(D_{l})$ and the quantity:

$\begin{array}[]{ll}\displaystyle\omega(a)\
=&\displaystyle\inf_{n}\>\min_{l\in{\cal
L}}\,\min_{i=1,2,\ldots,d}\,\inf\Big{\\{}\,\inf_{\mbox{\boldmath{\scriptsize$z$}}\in
S_{n}(\mbox{\boldmath{\scriptsize$r$}}(l),\mbox{\boldmath{\scriptsize$x$}})\cap
D_{l}}\,a_{ii}(\mbox{\boldmath{$z$}})\\\ &\displaystyle-\>\sum_{m=1,m\neq
i}^{d}\,\frac{r_{i}(l)}{r_{m}(l)}\,\sup_{\mbox{\boldmath{\scriptsize$z$}}\in
S_{n}(\mbox{\boldmath{\scriptsize$r$}}(l),\mbox{\boldmath{\scriptsize$x$}})\cap
D_{l}}\,|a_{im}(\mbox{\boldmath{$z$}})|\>:\>\mbox{\boldmath{$x$}}\in
G_{n}\Big{\\}}.\end{array}$ (28)

If $\omega(a)>0$, the chosen parameters $\mbox{\boldmath{$r$}}(l)$ ensure the
positive value of the brackets in item 3. of Discretization procedure 3.1. If
$\omega(a)\leq 0$, the partition should be refined until the condition
$\omega(a)>0$ is achieved. In accordance with Lemma 3.1, for each
$\mbox{\boldmath{$x$}}\in G_{n}$ there exist
$(r_{1}(\mbox{\boldmath{$x$}}),r_{2}(\mbox{\boldmath{$x$}}),\ldots,r_{d}(\mbox{\boldmath{$x$}}))$
such that $a_{ii}(\mbox{\boldmath{$x$}})-\sum_{m\neq
i}(r_{i}(\mbox{\boldmath{$x$}})/r_{m}(\mbox{\boldmath{$x$}}))|a_{im}(\mbox{\boldmath{$x$}})|>0$.
Due to the uniform continuity of the functions $a_{ij}$ on the sets ${\rm
cls}(D_{l})$, and b), the described procedure results with a desirable result
after a finite number of steps. QED

The uniform continuity of functions $a_{ij}$ on $D_{l}$, and the inequality
$\omega(a)>0$, where $\omega(a)$ is defined by (28), imply another important
property of matrices $A_{n}$ that are constructed by the above described
procedure. There exists $\sigma_{0}>0$, independent of $n$, such that

$\Big{|}\big{(}A_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}\Big{|}\
\geq\
\frac{\sigma_{0}}{h^{2}},\quad\mbox{\boldmath{$x$}}=h\mbox{\boldmath{$k$}}\in
G_{n}.$ (29)

The described construction of matrices $A_{n}$ for which (29) holds is called
admissible method, anticipating that the obtained $A_{n}$ have all the
necessary properties for the convergence of corresponding MJPs to generalized
diffusion. Let us recall $U_{i}=U_{i}(1)$.

###### LEMMA 3.2

Let the matrices $A_{n}$ on $G_{n}$ be discretizations of
$A_{0}=-\sum\partial_{i}a_{ij}\partial_{j}$ by an admissible method. Then the
matrices $A_{n}$ are irreducible. If in addition, the matrices $A_{n}$ are
symmetric, then there exist positive numbers $\underline{M},\overline{M}$ such
that the following inequalities

$\underline{M}\,\sum_{i}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{i}{\bf
u}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{2}^{2}\ \leq\ \langle\,{\bf
u}\,|\,A_{n}{\bf u}\,\rangle\ \leq\
\overline{M}\,\sum_{i}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{i}{\bf
u}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{2}^{2},$

are valid uniformly with respect to $n\in{\hbox{\msbm N}}$.

Proof:  It is easy to check that a symmetric matrix $A_{n}$ of positive type
is positive semidefinite, i.e. $\langle\,{\bf u}|A_{n}{\bf u}\,\rangle\geq 0$.
Let us consider the tensor valued functions $a,b$ where
$a(\mbox{\boldmath{$x$}})=\\{a_{ij}(\mbox{\boldmath{$x$}})\\}_{11}^{dd}$ and
$b(\mbox{\boldmath{$x$}})$ is defined by
$b_{ij}(\mbox{\boldmath{$x$}})=a_{ij}(\mbox{\boldmath{$x$}})-\kappa\delta_{ij}$.
The corresponding auxiliary tensors of (2) are denoted by $\hat{a},\hat{b}$.
Due to the strict positive definitness of the matrices $a,\hat{a}$ one can
choose $\kappa>0$ sufficiently small so that $b,\hat{b}$ are also positive
definite on ${{\hbox{\msbm R}}}^{d}$. Let us define matrices $H_{n}$ on
$G_{n}$ by the following non-trivial entries:

$\big{(}H_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}}\
=\
\frac{2d}{h^{2}},\quad\big{(}H_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$k$}}\pm\mbox{\boldmath{\scriptsize$e$}}_{i}}\
=\ -\frac{1}{h^{2}},\ i=1,2,\ldots,d.$

We have $\langle\,{\bf u}|H_{n}{\bf
u}\,\rangle=\sum_{i=1}^{d}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{i}{\bf
u}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{2}^{2}$. In accordance with
the construction of tensors $b,\hat{b}$ and the inequality (29), the matrices
$B_{n}=A_{n}-\kappa H_{n}$ are also of positive type. Since the symmetric
matrix $B_{n}$ of positive type is necessarily positive semidefinite, i.e.
$\langle\,{\bf u}|B_{n}{\bf u}\,\rangle\geq 0$ for any ${\bf u}\in
l_{0}(G_{n})$, we have:

$\langle\,{\bf u}\,|\,A_{n}\,{\bf u}\,\rangle\ =\ \kappa\langle\,{\bf
u}\,|\,H_{n}\,{\bf u}\,\rangle\>+\>\langle\,{\bf u}\,|\,B_{n}\,{\bf
u}\,\rangle\ \geq\
\kappa\sum_{i=1}^{d}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{i}{\bf
u}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{2}^{2},$

proving the lower bound of the assertion. The upper bound follows from (19)
and (25). The irreducibility follows from the graph theory, since any two
grid-knots $\mbox{\boldmath{$x$}}_{0},\mbox{\boldmath{$y$}}\in G_{n}$ can be
connected by a path of the form
$\\{\mbox{\boldmath{$x$}}_{0},\mbox{\boldmath{$x$}}_{1},\ldots,\mbox{\boldmath{$x$}}_{m},\mbox{\boldmath{$y$}}\\}\subset
G_{n}$ such that ${\cal N}(\mbox{\boldmath{$x$}}_{k-1})\cap{\cal
N}(\mbox{\boldmath{$x$}}_{k})\neq\emptyset$ for $k=1,2,\ldots,m$. QED

For a differential operator
$A(\mbox{\boldmath{$x$}})=A_{0}(\mbox{\boldmath{$x$}})+\mbox{\boldmath{$b$}}(\mbox{\boldmath{$x$}})\nabla$
with non-constant coefficients, the functions
$G_{n}\owns\mbox{\boldmath{$x$}}\mapsto\delta(A,n,\mbox{\boldmath{$x$}})p(\mbox{\boldmath{$x$}})$
of (14) are not identically zero on $G_{n}$. In the proof of convergence in
Section 4 the following weaker result is therefore used:

###### LEMMA 3.3

Let the differential operator $A(\mbox{\boldmath{$x$}})$ satisfy the
conditions of Theorem 3.1 and let the following additional conditions be valid
for some $\alpha\in(0,1)$:

a) Functions $a_{ij}$ belong to the class $C^{(1+\alpha)}$ on ${{\hbox{\msbm
R}}}^{d}$.

b) Functions $b_{i}$ belong to the class $C^{(\alpha)}$ on ${{\hbox{\msbm
R}}}^{d}$.

Then there exists a positive constant $\kappa$, depending on
$\overline{M},\|\mbox{\boldmath{$b$}}\|_{\infty}$, $\|a_{ij}\|^{(1+\alpha)}$
and $\|b_{i}\|^{(\alpha)}$, but not on $n$, such that

$\sup_{\mbox{\boldmath{\scriptsize$x$}}\in
G_{n}}\,|\delta(A,n,\mbox{\boldmath{$x$}})f(\mbox{\boldmath{$x$}})|\ \leq\
\kappa\,h^{\alpha}\,\|f\|^{(2+\alpha)}$

for any $f\in C^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})$.

This result can be easily checked by calculating $A_{n}{\bf f}_{n}$ directly.

Let us mention the following. If the functions $a_{ij},b_{i}$ are uniformly of
the class $C^{(1+\alpha)}$ and $C^{(\alpha)}$ on $D_{l}$, respectively, and
$f\in C^{(2+\alpha)}$ then
$\delta(A,n,\mbox{\boldmath{$x$}})f(\mbox{\boldmath{$x$}})$ converges to zero
on grid-knots $\mbox{\boldmath{$x$}}\in{\rm
int}(D_{l})\cap\big{(}\cup_{n}G_{n}\big{)}$. Otherwise,
$n\mapsto\delta(A,n,\mbox{\boldmath{$x$}})f(\mbox{\boldmath{$x$}})$ is not
bounded as $n$ increases. Nevertheless, the convergence in $W_{2}^{1}$-spaces
[LR2] is ensured as usually.

##  4 CONVERGENCE OF MJPs

The convergence of MJPs to generalized diffusion is analyzed here in terms of
the criterion (3). Therefore, we first need to define explicitly the mappings
$\Phi_{n}:\dot{l}_{\infty}(G_{n})\mapsto\dot{C}({{\hbox{\msbm R}}}^{d})$.

An element (column) ${\bf u}_{n}\in l(G_{n})$ can be associated to a
continuous function on ${{\hbox{\msbm R}}}^{d}$ in various ways. In the
current setting we define a mapping $l(G_{n})\mapsto C({{\hbox{\msbm
R}}}^{d})$ in terms of hat functions. Let $\chi$ be the canonical hat function
on ℝ, centered at the origin and having the support $[-1,1]$:

$\chi(x)\ =\ \left\\{\begin{array}[]{lll}1+x&{\rm for}&x\in[-1,0],\\\ 1-x&{\rm
for}&x\in[0,1],\\\ 0&{\rm for}&|x|>1.\end{array}\right.$

Then $z\mapsto\phi(h,x,z)=\chi(h^{-1}(z-x))$ is the hat function on ℝ,
centered at $x\in{{\hbox{\msbm R}}}$ with support $[x-h,x+h]$. The functions
$\mbox{\boldmath{$z$}}\>\mapsto\>\phi_{\mbox{\boldmath{\scriptsize$k$}}}(\mbox{\boldmath{$z$}})\>=\>\prod_{i=1}^{d}\>\phi(h,x_{i},z_{i}),x_{i}=hk_{i},i=1,2,\ldots,d$,
are the corresponding $d$-dimensional hat functions with support $S_{n}({\bf
1},\mbox{\boldmath{$x$}})=\prod_{i}[x_{i}-h,x_{i}+h]$. The functions
$\phi_{\mbox{\boldmath{\scriptsize$k$}}}(\cdot)\in l(G_{n})$, span a linear
space, denoted by $E(n,{{\hbox{\msbm R}}}^{d})$. Let ${\bf u}_{n}\in l(G_{n})$
have the entries $u_{n\mbox{\boldmath{\scriptsize$k$}}}=({\bf
u}_{n})_{\mbox{\boldmath{\scriptsize$k$}}}$. Then the function
$u(n)=\sum_{\mbox{\boldmath{\scriptsize$k$}}\in
I_{n}}u_{n\mbox{\boldmath{\scriptsize$k$}}}\phi_{\mbox{\boldmath{\scriptsize$k$}}}$
is an element of $E(n,{{\hbox{\msbm R}}}^{d})$ and defines an embedding of
grid-functions into the space of continuous functions. We denote the
corresponding mapping by $\Phi_{n}:l(G_{n})\mapsto E(n,{{\hbox{\msbm
R}}}^{d})$ and write

$u(n)\ =\ \Phi_{n}\,{\bf u}_{n}\ =\
\sum_{\mbox{\boldmath{\scriptsize$k$}}}\>\big{(}{\bf
u}_{n}\big{)}_{\mbox{\boldmath{\scriptsize$k$}}}\>\phi_{\mbox{\boldmath{\scriptsize$k$}}}.$
(30)

The inverse mapping $\Phi_{n}^{-1}:E(n,{{\hbox{\msbm R}}}^{d})\mapsto
l(G_{n})$ is defined by $\Phi_{n}^{-1}\big{(}\sum
u_{\mbox{\boldmath{\scriptsize$k$}}}\phi_{\mbox{\boldmath{\scriptsize$k$}}}\big{)}={\bf
u}$, where the column ${\bf u}$ has the entries
$u_{\mbox{\boldmath{\scriptsize$k$}}}$. It is obvious that the spaces
$l(G_{n})$ and $E(n,{{\hbox{\msbm R}}}^{d})$ are isomorphic with respect to
the pair of mappings $\Phi_{n},\Phi_{n}^{-1}$. Since $h_{n}=2^{-n}$, it is
clear that $E(n,{{\hbox{\msbm R}}}^{d})\subset E(n+1,{{\hbox{\msbm R}}}^{d})$
and that the space of functions $\cup_{n}E(n,{{\hbox{\msbm R}}}^{d})$ is dense
in $L_{p}({{\hbox{\msbm R}}}^{d}),p\in[1,\infty)$, as well as in
$\dot{C}({{\hbox{\msbm R}}}^{d})$. Let us mention that
$\sum_{\mbox{\boldmath{\scriptsize$k$}}}\phi_{\mbox{\boldmath{\scriptsize$k$}}}=1$
on ${{\hbox{\msbm R}}}^{d}$. For two functions $v(n),u(n)\in E(n,{{\hbox{\msbm
R}}}^{d})$ we have
$(v(n)|u(n))=h^{d}\sum_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}s_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}v_{\mbox{\boldmath{\scriptsize$k$}}}u_{\mbox{\boldmath{\scriptsize$k$}}}$
where
$s_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}=\|\phi_{\mbox{\boldmath{\scriptsize$k$}}}\|_{1}^{-1}(\phi_{\mbox{\boldmath{\scriptsize$k$}}}|\phi_{\mbox{\boldmath{\scriptsize$l$}}})$.
The numbers
$s_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}$ do not
depend on $n$ and the following identity is valid:
$\sum_{\mbox{\boldmath{\scriptsize$l$}}}s_{\mbox{\boldmath{\scriptsize$k$}}\mbox{\boldmath{\scriptsize$l$}}}=1$.
Thus we have $\Phi_{n}^{-1}\Phi_{n}=I$ in $l(G_{n})$ and
$\Phi_{n}\Phi_{n}^{-1}=I$ in $E(n,{{\hbox{\msbm R}}}^{d})$. Let $P(n)$ be the
projector onto $E(n,{{\hbox{\msbm R}}}^{d})$ defined by $f\mapsto
P(n)f=\Phi_{n}{\bf f}_{n}$. The mapping $\Phi_{n}^{-1}$ can be extended from
$E(n,{{\hbox{\msbm R}}}^{d})$ to $\dot{C}({{\hbox{\msbm R}}}^{d})$ by defining
$\Phi_{n}^{-1}f=\Phi_{n}^{-1}P(n)f$. Thus we have $\Phi_{n}\Phi_{n}^{-1}=P(n)$
in $\dot{C}({{\hbox{\msbm R}}}^{d})$.

If $F_{n}$ is a matrix on $G_{n}$, then $F(n)=\Phi_{n}F_{n}\Phi_{n}^{-1}$ is a
linear operator in the linear space $E(n,{{\hbox{\msbm R}}}^{d})$. It is easy
to verify that
$\|F(n)\|_{p}=\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,F_{n}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{p}$
for $p=1,\infty$, where the norm $\|\cdot\|_{p}$ is induced by the restriction
of $L_{p}({{\hbox{\msbm R}}}^{d})$ to $E(n,{{\hbox{\msbm R}}}^{d})$. By
applying the interpolation Rietz-Thorin theorem [BL] we get
$\|F(n)\|_{p}\leq\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,F_{n}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{p}$
for $p\in[1,\infty]$.

The objective of the analysis in this section is the comparison of the Feller
semigroup $U(\cdot)$ in $\dot{C}({{\hbox{\msbm R}}}^{d})$ and the matrix
semigroups $U(n,t)=\Phi_{n}\exp(-A_{n}t)\Phi_{n}^{-1}$ in
$\dot{E}(n,{{\hbox{\msbm R}}}^{d})=E(n,{{\hbox{\msbm
R}}}^{d})\cap\dot{C}({{\hbox{\msbm R}}}^{d})$, leading to a proof of (3). For
this purpose we consider the following initial value problems (IVP):

$\begin{array}[]{l}\big{(}\partial_{t}+A(\mbox{\boldmath{$x$}})\big{)}u(t,\mbox{\boldmath{$x$}})\
=\ 0,\quad u(0,\mbox{\boldmath{$x$}})=u_{0}(\mbox{\boldmath{$x$}}),\\\
\dot{{\bf u}}_{n}(t)+A_{n}{\bf u}_{n}(t)={\bf 0},\quad{\bf u}_{n}(0)={\bf
u}_{0n},\quad n\in{\hbox{\msbm N}}.\end{array}$ (31)

where ${\bf u}_{0n}$ are the discretizations of $u_{0}$. By using the standard
methods in the Sobolev space $W_{2}^{1}({{\hbox{\msbm R}}}^{d})$, we can first
prove a result which is weaker than (3). Let $A(\mbox{\boldmath{$x$}})$ be
defined by (4) and Assumption 2.1, and let $A_{n}$ be its discretizations on
$G_{n}$, constructed by the rules of Section 3. We consider $u(t)=U(t)f$ for
$f\in\dot{C}({{\hbox{\msbm R}}}^{d})\cap W_{2}^{1}({{\hbox{\msbm R}}}^{d})$,
the function $f(n)$ defined by (30), and $u(n,t)=U(n,t)f(n)$. As proved in
[LR2], for each specified $f$ we have
$\lim_{n\to\infty}\>\|u(t)-u(n,t)\|_{2,1}=0$, uniformly on segments of
$[0,\infty)$.

The function $u(t)=U(t)f$ for $f\in\dot{C}({{\hbox{\msbm R}}}^{d})\cap
W_{2}^{1}({{\hbox{\msbm R}}}^{d})$ is continuous on ${{\hbox{\msbm R}}}^{d}$,
uniformly on segments of $[0,\infty)$, as shown in [LSU]. The corresponding
functions $u(n,t)$ are also continuous on ${{\hbox{\msbm R}}}^{d}$, uniformly
on segments of $[0,\infty)$, as follows from their structure,
$u(n,t)=\Phi_{n}^{-1}\exp(-A_{n}t){\bf f}_{n}$. Due to the just described
convergence in $W_{2}^{1}({{\hbox{\msbm R}}}^{d})$, the sequence
$\mathfrak{U}=\\{u(n,t):n\in{\hbox{\msbm N}}\\}$ is bounded in
$W_{2}^{1}({{\hbox{\msbm R}}}^{d})$, uniformly on segments of $[0,\infty)$.
Hence, this sequence has a subsequence converging a.e. to $u$ on
${{\hbox{\msbm R}}}^{d}$. However, we need to show the uniform convergence. As
is usual in such problems, the uniform convergence can be proved eventually
for an appropriately selected subsequence of $\mathfrak{U}$. In our approach,
the outline of proof of (3) is as follows. The original differential operator
$A(\mbox{\boldmath{$x$}})$ is approximated by differential operators
$A^{(m)}(\mbox{\boldmath{$x$}})$ with smooth coefficients. The corresponding
semigroups are denoted by $U^{(m)}(\cdot)$ in $\dot{C}({{\hbox{\msbm
R}}}^{d})$ and $U^{(m)}(n,\cdot)$ in $\dot{l}_{\infty}(G_{n})$. The limit (3)
is then proved for each $m\in{\hbox{\msbm N}}$. Finally, by applying the
diagonalization argument to the sequence $\\{u^{(m)}(n,t):m,n\in{\hbox{\msbm
N}}\\}$ we get the desired result. The so outlined steps of the proof are
performed in the next two subsections.

### Convergence for smooth coefficients

Let the differential operator (4) have the coefficients $b_{i}$ that belong to
the class $C^{(\alpha)}$ and the coefficients $a_{ij}$ that belong to the
class $C^{(1+\alpha)}$ so that it can be represented as,

$A(\mbox{\boldmath{$x$}})\ =\
-\,\sum_{ij=1}^{d}\,a_{ij}(\mbox{\boldmath{$x$}})\partial_{i}\partial_{j}\,+\,\sum_{i=1}^{d}\,b_{i}^{\prime}(\mbox{\boldmath{$x$}})\partial_{i},$

where
$b_{i}^{\prime}(\mbox{\boldmath{$x$}})=b_{i}(\mbox{\boldmath{$x$}})-\sum_{j}\partial_{j}a_{ij}(\mbox{\boldmath{$x$}})$.
Hence, $A(\mbox{\boldmath{$x$}})$ can be represented as an elliptic operator
in non-divergence form with coefficients belonging to the class $C^{(\alpha)}$
on ${{\hbox{\msbm R}}}^{d}$. This form of $A(\mbox{\boldmath{$x$}})$ makes it
possible to use results on the existence of a strongly continuous semigroup in
$\dot{C}^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})=C^{(2+\alpha)}({{\hbox{\msbm
R}}}^{d})\cap\dot{C}({{\hbox{\msbm R}}}^{d})$ as developed by Solonnikov [So]
(a detailed exposition of results can be found in [LSU], Sections 13 and 14 of
Chapter 4). Thus we have $\|U(t)\|^{(2+\alpha)}\leq\exp(\sigma t)$ for some
$\sigma\geq 0$.

The matrices $A_{n}$ of Section 3 approximate $A(\mbox{\boldmath{$x$}})$ as
described in Lemma 3.3. We have
$\big{(}\Phi_{n}^{-1}A-A_{n}\Phi_{n}^{-1}\big{)}f(h\mbox{\boldmath{$k$}})=\delta(A,n,h\mbox{\boldmath{$k$}})f(h\mbox{\boldmath{$k$}})$
so that:

$\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,\big{(}\Phi_{n}^{-1}A-A_{n}\Phi_{n}^{-1}\big{)}f\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}\leq\kappa
h^{\alpha}\,\|f\|^{(2+\alpha)}.$ (32)

Now we have a straightforward application of this result on approximations:

###### LEMMA 4.1

Let $t\mapsto U(t)$ be a strongly continuous semigroup in
$\dot{C}^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})$,
$\|U(t)\|^{2+\alpha}\leq\exp(\sigma t)$ with some $\sigma\geq 0$, such that
$u(t)=U(t)u_{0}$ solves the first IVP in (31). Let $t\mapsto U_{n}(t)$,
$\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{n}(t)\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}\leq
1$, be semigroups generated by $-A_{n}$ in $\dot{l}_{\infty}(G_{n})$, such
that ${\bf u}_{n}(t)=U_{n}(t){\bf u}_{0n}$ solve the IVPs for ODE in (31). If
the differential operator $A(\mbox{\boldmath{$x$}})$ on ${{\hbox{\msbm
R}}}^{d}$ fulfills the condition (32), then the following assertion is valid:
For each $T>0$ there exists a positive number $\rho(T)$ such that

$\displaystyle\sup_{t\in[0,T]}\,\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,\Phi_{n}^{-1}\,u(t)\,-\,{\bf
u}_{n}(t)\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}\ \leq\
\rho(T)\|u_{0}\|^{2+\alpha}\,h^{\alpha},$

for all $u_{0}\in\dot{C}^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})$.

Proof:  The function $s\mapsto{\bf
f}(s)=U_{n}(t-s)\Phi_{n}^{-1}U(s)u_{0}\in\dot{l}_{\infty}(G_{n})$, for $0\leq
s\leq t,u_{0}\in\dot{C}^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})$, has a
continuous derivative of the form ${\bf
f}^{\prime}(s)=U_{n}(t-s)\big{(}A_{n}\Phi_{n}^{-1}U(s)-\Phi_{n}^{-1}AU(s)\big{)}u_{0}$.
Therefore the following identity must be valid for each
$u_{0}\in\dot{C}^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})$:

$\Big{(}\Phi_{n}^{-1}U(t)-U_{n}(t)\,\Phi_{n}^{-1}\Big{)}u_{0}\ =\
\int_{0}^{t}\
U_{n}(t-s)\,\big{(}A_{n}\,\Phi_{n}^{-1}-\Phi_{n}^{-1}\,A\big{)}\,U(s)\,u_{0}\,ds.$

The $\dot{l}_{\infty}(G_{n})$-norm of the integrand is first estimated from
above by

$\begin{array}[]{c}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,U_{n}(t-s)\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}\,\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,\big{(}A_{n}\,\Phi_{n}^{-1}\,-\,\Phi_{n}^{-1}\,A\big{)}\,U(s)u_{0}\,\rule[-2.84526pt]{1.70717pt}{11.38109pt}\,_{\infty}\
\leq\ \kappa h^{\alpha}\,\|U(s)u_{0}\|^{(2+\alpha)}\\\ \leq\ \kappa
h^{\alpha}\,\exp(\sigma T)\,\|u_{0}\|^{(2+\alpha)},\end{array}$

where (32) is used. Hence,

$\,\rule[-5.69054pt]{1.70717pt}{19.91692pt}\,\,\Big{(}\Phi_{n}^{-1}U(t)\,-\,U_{n}(t)\Phi_{n}^{-1}\Big{)}u_{0}\,\,\rule[-5.69054pt]{1.70717pt}{19.91692pt}\,_{\raisebox{-1.4pt}{\mbox{$\infty$}}}\
\leq\ T\,\exp(\sigma T)\,\kappa\,h^{\alpha}\,\|u_{0}\|^{2+\alpha},$

implying the assertion. QED

Due to the density of $\dot{C}^{(2+\alpha)}({{\hbox{\msbm R}}}^{d})$ in
$\dot{C}({{\hbox{\msbm R}}}^{d})$ we have the following auxiliary result.

###### COROLLARY 4.1

Let $A(\mbox{\boldmath{$x$}})$ be defined by (4) and Assumption 2.1, and let
it fulfil the conditions of Theorem 3.1. If the coefficients $b_{i}$ belong to
the class $C^{(\alpha)}$ on ${{\hbox{\msbm R}}}^{d}$ and $a_{ij}$ belong to
the class $C^{(1+\alpha)}$ on ${{\hbox{\msbm R}}}^{d}$ then:

(i) The operators $U(n,t)P(n)$ converge strongly in $\dot{C}({{\hbox{\msbm
R}}}^{d})$ to $U(t)$, uniformly on segments of $[0,\infty)$.

(ii) The limit (3) is valid.

The so called classical diffusion, i.e. the process in ${{\hbox{\msbm
R}}}^{d}$ determined by the differential operator
$A(\mbox{\boldmath{$x$}})=-\sum_{ij=1}^{d}a_{ij}(\mbox{\boldmath{$x$}})\partial_{i}\partial_{j}$,
is usually simulated by using its representation in trems of stochastic
differential equations. Corollary 4.1 makes it possible to simulate sample
paths of classical diffusion in terms of MJPs. This alternative approach to
simulation gives better results in an estimation of the statistical moments of
the first exit time from open sets at subsets of the boundary with a rapidly
changing normal,

### Convergence in the general case

Let us define a mollifier
$\mbox{\boldmath{$x$}}\mapsto\vartheta(n,\mbox{\boldmath{$x$}})=h^{-d}\vartheta(h^{-1}\mbox{\boldmath{$x$}})$
in terms of a non-negative function $\vartheta(\cdot)$ of class $C^{(2)}$ on
${{\hbox{\msbm R}}}^{d}$ with the support equal to the unit ball
$B_{1}(0)\subset{{\hbox{\msbm R}}}^{d}$ and
$\int\vartheta(\mbox{\boldmath{$x$}})d\mbox{\boldmath{$x$}}=1$. A smoothing
procedure of coefficients $a_{ij},b_{i}$ is determined by replacing these
coefficients with the sequence of coefficients $a_{ij}^{(m)}=\vartheta(m)\ast
a_{ij},b_{i}^{(m)}=\vartheta(m)\ast b_{i}$, $m\in{\hbox{\msbm N}}$, where
$\ast$ denotes the convolution operator. For each $m\in{\hbox{\msbm N}}$ the
resulting differential operator $A^{(m)}(\mbox{\boldmath{$x$}})$ has the
closure in $\dot{C}({{\hbox{\msbm R}}}^{d})$ and generates the Feller
semigroup $U^{(m)}(\cdot)$. Now we consider the mappings:

$\begin{array}[]{ccc}A(\mbox{\boldmath{$x$}})&\longmapsto&\\{A^{(m)}(\mbox{\boldmath{$x$}})\,:\,m\in{\hbox{\msbm
N}}\\},\\\ \downarrow&&\downarrow\\\
A_{n}&\longmapsto&\\{A_{n}^{(m)}\,:\,m\in{\hbox{\msbm N}}\\},\end{array}$ (33)

and the operators $U^{(m)}(t)$ and $U^{(m)}(n,t)P(n)$ in the Banach space
$\dot{C}({{\hbox{\msbm R}}}^{d})$, where the semigroup $U^{(m)}(\cdot)$ is
generated by the closure of $A^{(m)}(\mbox{\boldmath{$x$}})$, and where
$U^{(m)}(n,t)=\exp(-A_{n}^{(m)}t)$. For the double sequence of operators
$U^{(m)}(n,t)P(n),m,n\in{\hbox{\msbm N}}$, we will prove the following limits:

$U^{(m)}(n,t)P(n)\ \begin{array}[]{c}\raisebox{-4.2679pt}{$n\to\infty$}\\\
\longrightarrow\\\ \\\ \end{array}\ U^{(m)}(t)\
\begin{array}[]{c}\raisebox{-4.2679pt}{$m\to\infty$}\\\ \longrightarrow\\\ \\\
\end{array}\ U(t),$ (34)

uniformly on segments of $[0,\infty)$. By using (34) and applying the
diagonalization argument to the sequence
$\\{U^{(m)}(n,\cdot)P(m):(m,n)\in{{\hbox{\msbm N}}}^{2}\\}$, we get the main
result of this article:

###### THEOREM 4.1

Let the differential operator $A(\mbox{\boldmath{$x$}})$ be defined by (4) and
Assumption 2.1. There exists a sequence of pairs $(m,n(m))\in{{\hbox{\msbm
N}}}^{2}$ such that

$\lim_{m}\,\|U_{n(m)}^{(m)}(t)P(n(m))-U(t)\|_{\infty}\ =\ 0,$

uniformly on segments of $[0,\infty)$. Therefore, the asymptotic (3) is valid
for the sequence $\\{U_{n(m)}^{(m)}(\cdot):m\in{\hbox{\msbm N}}\\}$.

Proof:  Since the coefficients $b_{i}^{(m)}$ belong to the class
$C^{(\alpha)}$ on ${{\hbox{\msbm R}}}^{d}$ and the coefficients $a_{ij}^{(m)}$
belong to the class $C^{(1+\alpha)}$ on ${{\hbox{\msbm R}}}^{d}$, the first
limit in (34) follows from Corollary 4.1.

It remains to prove the second limit in (34), that is, the convergence of the
semigroups $U^{(m)}(\cdot)$ to $U(\cdot)$ in $\dot{C}({{\hbox{\msbm
R}}}^{d})$, uniformly on segments of $[0,\infty)$. The convergence of
$U^{(m)}(\cdot)$ to $U(\cdot)$ in $L_{2}({{\hbox{\msbm R}}}^{d})$ can easily
be obtained (for instance, by using Theorem 6.1, Chapter 1 in [EK]). The
convergence of $U^{(m)}(t)$ to $U(t)$ in $\dot{C}({{\hbox{\msbm R}}}^{d})$,
uniformly on segments of $[0,\infty)$, would follow from such convergence on a
dense subspace in $\dot{C}({{\hbox{\msbm R}}}^{d})$. We choose the subspace
$C_{0}^{(\alpha)}({{\hbox{\msbm R}}}^{d})=C_{0}({{\hbox{\msbm R}}}^{d})\cap
C^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ which is dense in both, the space
$L_{2}({{\hbox{\msbm R}}}^{d})$ and the space $\dot{C}({{\hbox{\msbm
R}}}^{d})$. Due to (i) of Corollary 2.1, the operators $U(t)$ and $U^{(m)}(t)$
are bounded in the space $\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$,
uniformly on segments $K\subset[0,\infty)$, i.e.
$\|U^{(m)}(t)\|^{(\alpha)}\leq\beta(K)$, where $\beta(K)$ does not depend on
$m$. Thus we come to the following conclusion. The operators $U(t)-U^{(m)}(t)$
are continuous mappings from $C_{0}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$ into
$\dot{C}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$,
$\|U(t)-U^{(m)}(t)\|^{(\alpha)}\leq\beta(K)$, and they converge to zero in
$L_{2}({{\hbox{\msbm R}}}^{d})$, uniformly on segments $K\subset[0,\infty)$.

Now we apply the following auxiliary result to
$u_{m}(t)=(U(t)-U^{(m)}(t))v,v\in C_{0}^{(\alpha)}({{\hbox{\msbm R}}}^{d})$.
Let $\mathfrak{U}=\\{u_{n}:n\in{{\hbox{\msbm N}}}\\}$ be a sequence of
continuous functions on $[0,1]\times{{\hbox{\msbm R}}}^{d}$ such that:

a) $u_{n}(t)\in L_{2}({{\hbox{\msbm R}}}^{d})$ for each $t\in[0,1]$, and
$\sup\\{\lim_{n}\|u_{n}(t)\|_{2}:t\in[0,1]\\}=0$.

b) The functions $u_{n}$ are uniformly Hölder continuous in the following
sense. There exist $\alpha\in(0,1)$ and $c_{\alpha}>0$, which do not depend on
$t$ or $x$, such that the restrictions $u_{n}(t)|B_{1}(\mbox{\boldmath{$x$}})$
fulfil the following two conditions:

$u_{n}(t)|B_{1}(\mbox{\boldmath{$x$}})\,\in\,C^{(\alpha)}(B_{1}(\mbox{\boldmath{$x$}})),\qquad\|u_{n}(t)\|_{\infty}^{(\alpha)}\,\leq\,c_{\alpha},$

uniformly with respect to $\mbox{\boldmath{$x$}}\in{{\hbox{\msbm R}}}^{d}$ and
$t\in[0,1]$.

Then $\lim_{n}u_{n}(t)=0$ in $\dot{C}({{\hbox{\msbm R}}}^{d})$, uniformly with
respect to $t\in[0,1]$.

A proof of this auxiliary result is simple. If the assertion were not valid,
there would exist a positive number $\delta$ and a sequence of pairs
$(t_{k},\mbox{\boldmath{$x$}}_{k})\in[0,1]\times{{\hbox{\msbm
R}}}^{d},\lim_{k}|\mbox{\boldmath{$x$}}_{k}|=\infty$, such that
$u_{k}(t_{k},\mbox{\boldmath{$x$}}_{k})\geq\delta$. Due to b) the following
must also be valid: $|u_{k}(t_{k})|\geq\delta/2$ on the ball
$B_{r}(\mbox{\boldmath{$x$}}_{k})$ where $r=(\delta/2c_{\alpha})^{1/\alpha}$.
A consequence of these inequalities would be $\|u_{k}(t_{k})\|_{2}\geq
2^{-1}\delta\sqrt{|B_{r}({\bf 0})|}$, contradicting a). QED

We remind the reader that the approximations $A(\mbox{\boldmath{$x$}})\mapsto
A^{(m)}(\mbox{\boldmath{$x$}})$, constructed in the proof of Theorem 4.1, are
not the only ones that are needed for the proof of (3). The approximations
defined by (23) are also needed in order to get matrices $A_{n}$ of positive
type.

##  5 SIMULATION OF SAMPLE PATHS

Here we demonstrate the efficiency of simulation of sample paths of a
generalized diffusion by using MJPs. We intend to estimate the expectation and
the variance of the first exit time from a Lipshitz domain. The differential
operator
$A(\mbox{\boldmath{$x$}})=-\sum_{ij=1}^{2}\partial_{i}a_{ij}(\mbox{\boldmath{$x$}})\partial_{j}$
on ${{\hbox{\msbm R}}}^{2}$ is defined by its diffusion tensor, being a
piecewise constant tensor-valued function of the form,

$a(\mbox{\boldmath{$x$}})\ =\
\left[\begin{array}[]{cc}\sigma^{2}&\alpha(\mbox{\boldmath{$x$}})\\\
\alpha(\mbox{\boldmath{$x$}})&1\end{array}\right],\quad\alpha(\mbox{\boldmath{$x$}})\>=\>\rho{\hbox{1\kern-1.49994pt\vrule
height=6.88889pt,width=0.3pt\vrule width=0.8pt,height=0.25pt\kern
1.49994pt}}_{D_{0}}(\mbox{\boldmath{$x$}}),\quad\rho^{2}<\sigma^{2},$

where $\sigma^{2}$ is a positive number, $\rho$ is a real number and
$D_{0}=(1/4,3/4)^{2}$. Let $X(\cdot)$ be diffusion determined by
$A(\mbox{\boldmath{$x$}})$ starting from
$\mbox{\boldmath{$x$}}_{0}=(1/2,1/2)$. For a bounded Lipshitz domain $D$ the
expectation and the variance of first exit time from $D$ are given by
expressions [Li]:

${\bf E}[\theta]\ =\ \|u\|_{1},\quad{\bf Var}[\theta]\ =\ 2\|v\|_{1}\,-\,{\bf
E}[\theta]^{2},$ (35)

where $u,v$ are the unique solutions of the following boundary value problems,

$\begin{array}[]{cc}\displaystyle
A(\mbox{\boldmath{$x$}})u(\mbox{\boldmath{$x$}})\ =\
\delta(\mbox{\boldmath{$x$}}-\mbox{\boldmath{$x$}}_{0}),&{\bf x}\in D,\\\
u\,|\,\partial D\ =\ 0,&\end{array}\quad\begin{array}[]{cc}\displaystyle
A(\mbox{\boldmath{$x$}})v(\mbox{\boldmath{$x$}})\ =\
u(\mbox{\boldmath{$x$}}),&{\bf x}\in D,\\\ v\,|\,\partial D\ =\
0,&\end{array}$ (36)

and $\delta(\mbox{\boldmath{$x$}})$ is Dirac $\delta$-function at ${\bf 0}$.
We intend to compute ${\bf E}[\theta],{\bf Var}[\theta]$ by simulations and by
using deterministic methods formulated in terms of Expressions (35), (36). In
order to simulate sample paths of $X(\cdot)$ we shall approximate the
diffusion by a MJP $X_{n}(\cdot)$ and simulate sample paths of $X_{n}(\cdot)$.
In this example we choose $\sigma^{2}=0.1,\rho=0.02$, $D=(0,1)^{2}$ and two
cases of discretizations, $h=1/200$ and $h=1/400$.

The generator of the process $X_{n}(\cdot)$ in $G_{n}$ is denoted by
$A_{n}({\rm gen})$. Since the entries $a_{12}$ are non-trivial on $D_{0}$, we
have to use the construction of Section 3 in order to get $A_{n}({\rm
gen})=-A_{n}$. The parameters $r_{1}=3,r_{2}=1$ of construction are
illustrated in Figure 1. These values of parameters ensure $A_{n}$ to have the
structure of a matrix of positive type.

Two boundary value problems of (36) have unique solutions $u,v\in L_{1}(D)$ as
proved in [BO, LR1]. An efficient numerical method is constructed and the
convergence in $L_{1}(D)$ is proved in [LR2]. This numerical method is based
on the construction of $A_{n}$ which is described in Section 3. Efficiency of
constructed methods is demonstrated by examples in which solutions in closed
forms are compared with numerical solutions.

Results of computation are expressed in terms of relative errors:

$\varepsilon_{exp}\>=\>\frac{<\theta>_{det}-<\theta>_{sim}}{<\theta>_{det}},\quad\varepsilon_{var}\>=\>\frac{\,<<\theta>>_{det}-\,<<\theta>>_{sim}}{\,<<\theta>>_{det}},$

where $<\theta>_{det},<\theta>_{sim}$ are the estimates of ${\bf E}[\theta]$
obtained by using deterministic methods (35) and Monte Carlo simulations,
respectively. Analogously, $\,<<\theta>>_{det},\,<<\theta>>_{sim}$ are the
corresponding quantities for estimates of ${\bf Var}[\theta]$. Some results of
computations are given in the table bellow. The last column contains the
ratios, $r=t_{det}/t_{sim}$, of computational times $t_{det}$ and $t_{sim}$ of
deterministic and Monte Carlo method, respectively. Sample paths are simulated
20000 times.

$h=1/n$ | $\varepsilon_{exp}$ | $\varepsilon_{var}$ | $r$
---|---|---|---
n = 200 |  -0. | 039 | 0. | 018 | 15. | 7
n = 400 | -0. | 044 | 0. | 007 | 20. | 3

Comparison of results obtained by

deterministic and Monte Carlo methods

As expected, the first two statistical moments of the first exit time can be
estimated by Monte Carlo simulations dozen times faster than by using the
deterministic method formulated by (35) and (36).

Acknowledgement. I would like to thank Vlada Limic for numerous discussions
and helpful comments on previous versions of the manuscript.

## References

  * [BO] Bottaro G. & Oppezzi P., Elliptic equations in divergence form, with right hand side a measure, on unbounded domains, Ann. Mat. Pura Appl., 163, 223–246, 1993
  * [BL] Bergh J. & Löfström J., Interpolation Spaces. An Introduction, Springer Verlag, 1976.
  * [EK] Ethier S. N. & Kurtz T. G., Markov Processes, characteristics and Convergence, Wiley, New York, 1986.
  * [Et] Étoré P., On random walk simulation of one-dimensional diffusion process with discontinuous coefficients, Electronic Journal of Probability, 11, 2006, 249-275.
  * [LSU] Ladyzhenskaya O. A., Solonnikov V. A.& Ural’tseva N. N., Linear and Quasilinear Equations of Parabolic type, American mathematical Society, Providence, Rhode island, 1968.
  * [LU] Ladyzhenskaya O.A. & Ural’tseva N.N., Linear and Quasilinear Elliptic Partial Differential Equations, Academic Press, N.Y., 1968.
  * [LM] Lejay A. & Martinez M., A scheme for simulating one-dimensional diffusion process with discontinuous coefficients, The Annales of Applied Probability 16:1, 107-139, 2006.
  * [Li] Limić N., Monte Carlo Simulations of Random variables, Sequences and processes, Element, Zagreb, 2009.
  * [LR1] Limić N. & Rogina M., Explicit stable methods for second order parabolic systems, Math. Commun. Vol 5, 97–115, 2000
  * [LR2] Limić N. & Rogina M., Monotone numerical schemes for a Dirichlet problem for elliptic operators in divergence form, Math. Meth. Appl. Sci. 32, 1129-1155, 2009.
  * [Ma] Maz’ya V.G., Sobolev spaces, Springer, N.Y., 1985
  * [MW] Motzkyn T.S. & Wasov W., On the approximation of linear elliptic differential equations by difference equations with positive coefficients, J. Mathematical Physics 31, 253-259 1953.
  * [Se] Stein F.M., Singular Integrals and Differentiability Properties of Functions, Princeton Univ. Press, New Jersey, 1970
  * [So] Solonnikov V.A., A prioeri estimates for second order equations of parabolic type, Trudy Math. Inst. Steklov, 70, 133-212, 1964; On a boundary value problem for linear general parabolic system of differential equations, Trudy Math. Inst. Steklov, 83, 1965.
  * [St] Stroock D.W., Diffusion semigroups corresponding to uniformly elliptic divergence form operators, Séminaire de Probabilité XXII, Springer Verlag LNMS # 1321, 316-348, 1988.
  * [SZ] Stroock D.W. & W. Zheng, Markov chain approximations to symmetric diffusions, Annal, Section B, 33, 619-649 1997.
  * [TS] Tikhonov AN $ Samarskii AA, On the stability of difference schemes, Dokl. Akad. Nauk SSSR, 149, 529-531 1963.