5 Analysis of the Operator-Splitting Method for the Model Problem
In this section the properties of the operator-splitting method of the
previous section are considered. In particular the convergence of
the Jacobi-like method on the model problem is analysed.
Returning to the general operator-splitting method of section 4.1,
since the operator B of the previous section is an approximate
inverse of K then it is helpful to write
BK = P - Â
(14)
where Â:X ® Z
can be regarded as the residual operator. It follows from
(12) that
P f(n) = Â f(n-1) + B g
(15)
and that
P f(n) = ( Á + Â + Â2+ ... + Â(n-1)) Bg
(16)
where Á:Z ® Z is the identity operator.
The series (16) is a Neumann series. The iterative method
would always converge if ||Â|| < 1,
in which case (15) would be a contraction mapping
and K f=g would have a unique solution.
If the function g has been obtained through applying the operator
K to a function f then substituting K f for g in (16) gives
P f(n) = ( Á + Â + Â2+ ... + Ân-1) BK f.
and hence that
P f(n) = ( Á + Â + Â2+ ... + Ân-1) (P-Â)f
= (P - Ân) f.
It follows that
the error d(n) in the approximation to f at in iterate n is
given by
P d(n) = P ( f - f(n) ) = Ân f
(17)
and the error e(n) in the corresponding approximation to g is
Returning to the simple diffusion equation of section 3, recall that
the eigenvalues are kN = e-N2 p2 T and hence they
cluster and tend toward zero as N ® ¥. It follows
that the corresponding rN=1-kN lie in [0,1] and tend toward unity
as N ® ¥. Thus, as a result of this, Â is
not a strict contraction mapping.
Let us now consider the the domain of the diffusion equation to
be (x,t) Î [0,1] ×[0,T]
with u(x,0)=f(x), u(x,T)=g(x) and u(0,t)=u(1,t)=0.
Given any particular eigenfunction
fN = sin(N px),
it follows from (18) that the error in iterate n
is
eN(n) = gN - gN(n) = kN rNn fN = e-N2 p2 T (1-e-N2 p2 T)n sin(N px)
(20)
where X=Y=C¥[0,1] and P is the identity operator.
Since rN ® ¥ as N ® ¥ then
(19) is a contraction when N is finite.
In this case, the Jacobi-like method converges
to the original f and the error function for iteration n
is given by (19).
The eigenvalues of the model problem converge rapidly to zero.
If l is the eigenvalue corresponding to
N=1 then the complete set of eigenvalues form the
sequence l, l4, l9, l16,...,
and these correspond to the eigenfunctions
sin(px),
sin(2px),
sin(3px),
sin(4px),... ,
Hence if the first eigenvalue is 0.9 then the the simple
diffusion equation (6)
has the eigenvalues 0.9, 0.6561, 0.387420,
0.185302, 0.071790, 0.022528, 0.005726, 0.001179, 0.000197,
0.000026, ...
to six decimal places, if the first eigenvalue is 0.5 then the
eigenvalues are
0.5, 0.0625, 0.001953, 0.000015,....
Since rN = 1- kN then the sequence rN, N=1,2,...¥
converges quickly to unity with N.
It follows that the error term converges much more slowly to
zero for the higher eigenfunctions. For example if the
eigenvalues are
0.5, 0.0625, 0.001953, 0.000015,...
then the rn are 0.5, 0.9373, 0.998047, 0.999985,.... In
order to return a relative error of less than 1% the number of
iterations required is 6, 71, 2355, and 307009 for the
first four eigenfunctions. Test problems demonstrating the
convergence rates of the method on the eigenfunctions of the
model problem are given in the next section.
The analysis of the previous section suggests that progress
towards the inverse solution of the first
eigenfunction of K is relatively rapid,
with the rate of convergence declining as we pass through the
eigenfunctions. The slower rate of convergence of the
higher eigenfunctions is a positive property of the method. Potentially,
the iteration could be terminated before the higher eigenfunctions
(ultimately the noise) affect the solution. For example the iteration
can be terminated once
iterate f(n) has the property that
||K f(n) - g|| < e
where e represents the uncertainty in g as a result of
measurement or numerical error.
It is in this respect that the method is similar to the truncated
eigenfunction expansion method described in section 2.3. However,
in this method the involvement of the higher eigenfunctions is
naturally delayed by the method, rather than being deliberately removed.