An analysis of the modified exclusion method in certain problems

U.S.S.R. Comput.Maths.Math.Phys., Printed in Great Britain

Vo1.25,No.2,pp.101-X6,1985

0341-5553/85 $10.00+0.00 Pergomon Journals Ltd.

ANANALYSISOF THE MODIFIED EXCLUSION METHOD IN CERTAIN PROBLEMS*

V.V. EEPYAEB A quantitative estimate is given of the error of a modified exclusion method with an increased stability to rounding errors when solving symmetric grid equations generated by a second-order scalar elliptic operator.

In /l/, where the application of the modified exclusion method to solving certain grid equations was described, the basis of the algorithm's increased stability to rounding errors had a qualitative character. In this paper a quantitative estimate of the algorithm's error to a first approximation of perturbation theory when solving sysnsetric equations generated The results obtained here can henceforth by a second-order scalar elliptic operator is given. be used as the basis of the increased stability of the modified method when solving a wider class of problems, as shown by numerical experiments /2/. 1 Certain properties of second-order grid equations. Before estimating the error, we shall carry out a preliminary analysis. 1. Suppose the coefficients of the positive-definite matrix A of order n satisfy following conditions:

the

i-n

a,,>& Let us consider

O,,SO,

i,j=i,2

,...,

n,

j+i;

I&= &>O,

i=i,

2,.. .) rz.

(1.1)

I-1

the influence

of a special kind of error on the solution

of the equation

AU--F. We introduce

matrix

A’ of order n, for which aulnu’Sp.

(1.2) the following relations,hold:

a,p-0.

i,j==l.Z,...,n.

aG

(1.3) d,/d,‘
We formulate

and

the following

i=l,

2..

. , n.

statement.

If two symmetric positive-definite matrices A and A’ satisfy Theorem. (1.3), the following expression holds for the arbitrary vector U : zGUU’AU/

Proof. in the form

It is easy to prove

Ed,& i-1

(1.1)

U’A’V<‘B.

that the quadratic

U’AU-

conditions

(1.4) form for the matrix A can be represented

‘c’r,a,,(~i-,,)~. i-l I--r+‘

A similar representation is also possible for U’A’U. u'_-1U and WA’U are sums of non-negative quantities. from (1.5) and (1.3).

(1.5)

It follows from (1.5) and (1.1) that The correctness of Theorem 1 follows

cz
2.

The estimate (V,'-V,)'~(U,'-C',)stl'r:,'dC*

(1.8)

for the solution (I,and u,'of the equations AU=F. d'u-F where vj-mas(la-11/a, }fi-il is). The correctness of corollary 2 follows from (1.4) and from the known estimate in /4/. 2. Suppose the positive-definite matrix A, which satisfies conditions (1.11, is expanded in products of the lower triangular matrix N and the upper triangular matrix A”-‘:

holds

where

R, R,,R,

Suppose

are

the errors

matrices

A+R=N’(A”-I)‘-(N+R,) (A’-‘+&), of the errors, such that Eq.tl.2)

(1.7) is solved

N’Y=F, (A”-‘)‘U=Y. of the coefficients of the matrix A"-' satisfy

1(a;-‘)‘/a;-‘-iJCe,

lZh.vychisl.Mat.mat.Fiz.,25,4,483-491,1985

101

o,:-‘zo,

in two stages: (1.8)

the condition

(1.9)

102 and the matrix N' is determined by the matrix A”-’ as in the symmetric scheme of Gauss's We shall consider the influence of errors of the type (1.9) on the SOhtiOn Of Eq. method. (1.2). (A"-') by ((al-')')-'", i=i, 2 ,..., n. we By multiplying the i-th row of the matrix reduce II (1.7) to the form /5/ A+R=(L+R,)'(LtR,). A-EL, where &is

the error matrix.

According

to (1.91, the relation

II,;i&,-i1
will hold for the coefficients

(MO) It follows

from

(1.8) that

(L’)‘Y--P,

We shall use the general

formula

(i.if)

L’U- Y. /5, p.180/ for the estimate

(1.12) 11811, H-_((L')').'P-(L')-'F

(1.13) The spherical

norm

l/j1 will be applied

Taking into account (see 3., para.1)

below.

the connection M

between

~ k" j&II

IlLI (

IIVI

A and L, we shall reduce i_k" !$

(1.13) to the form

) -I,

(1.14)

where k is the conditionality number of the matrix A: k-_IJAIIIJA-‘11. Further arguments will be produced in the first approximation for the condition

of perturbation

k’hllRsII/ IIW 1.

In this approximation

(1.14) transforms

follows from (1.12) that L’U’=Y+H we denote by Ati, equals the sum of

(1.16)

such that the error of the solution A,U and AZ,

AJ-(L') -‘I’-L-‘Y, A,U== (L’)-‘H. AIIestimate similar to (1.16) holds for IlA,Ull llAtUllwe shall represent

of (1.21, which

(1.17)

(1.H)

llAJII / Il~l,ll~k”‘ll~4I IILII. For the estimate

(1.15)

to the relation

IIHII1IIYll~k”W~il i IPII.

It

theory

(1.17) in the form

A,U-L-‘H+(

(L’)-‘-L-‘)H.

(l.iQ)

With the notation Z-L-iH, from

(1.19) it follows

z,-((L’)_‘-L-‘)H

(1.20)

that

(~+ll~~ll I llzll) 441 (i+k”llHsII I IlLll).

llbW=ll~+W4l4I+IlZ,ll-ll4I To estimate

llrll we use the conditionality

We can show

(see 3., para.l),

number

(1.21)

of Eq.(1.2)

a-ll-WIII~lI 1 ilW that ~~C~oil>i~YII~lA-‘ll”,a-“.. It follows

from

(1.15),

(1.16),

(1.20)-(1.22)

&II< IIU.lll

(1.22)

that

lb-‘II -.ll~ll < ar’,,k’,,p&II IIA-‘Il’k”* IIYII Il.4 ’

lIA4l / jlUolJ~a’“k’“ijR,(I / [IL/l. The estimate

(1.18) and

(1.23) lead to the relation

x, IlW We reduce

(1.23)

lIW+Il~z~< IlUall

+

(i+a’,,)k”_j!$. 11~11

(I.24

(1.24) to the form,

(IbUll /liU#Se (14-z"') k’“, where E corresponds to (1.9) and (l.lO), and S is a coefficient importance will be given to the value S in Sect.2 as a function error of matrix L.

(1.25) of proportionality. Specific of the peculiarities of the

2. Estimating the error. We shall estimate the error of the modified method when solving Eq.(1.2) in the first approximation of the perturbations theory. In this approximation for estimating the general error one can sum the estimates of the errs of the solution produced by different operations. Rounding errors for operations which can be achieved with double the accuracy are not taken into consideration.

103 1. Let us consider the direct run of matrix run, using the same notation as in /l/. The following are determined:

A.

cJ,‘=(l,..., i), I-th step: hold at the

The following operations

We shall cite the algorithm

of the direct

U=ACl.. I-I

’ -d I+‘ ' al+,.r+l-

(2.1)

r, a,:~, j-t+*

m>l+l, A;" =A,'+(-a,:,.,/a,:,.r+~)A:,l,

(2.2)

d:+‘=d,‘+(-a:+,.Ja:+,,,+,)d:+I,

(2.3

m>1-l-t

The error of the solution, produced by rounding errors with the direct run, equals the In its turn, sum of the errors introduced into the solution at each step of the direct path. the error introduced into the solution in one step of the direct path is determined by summing According the errors introduced into the solution by each of the operations (2.1)-(2.3). to (2.2) and (2.3) three arithmetic operations are required and there are no subtractions, such that /5, 6/ (a:'

l+’

n

(dm

)“=“:‘(i*E,),

) =d,

I+,

to calculate

a$

and

ICI , j>m

d,.

(l*e,),

(2.4)

c-0, i if the computer has a t-digit mantissa when representing the where le,1~E,-3.2-‘, number in a form with a floating point. In (2.4) the approximate values of the coefficients of the matrices and vectors - as distinct from the accurate values - are denoted by primes. We shall continue with this notation below. In (2.1) an arithmetic sum is accumulated. We assume, as in /5, 6/, that the accumulation of algebraic sums must be carried out with double the accuracy with subsequent rounding, therefore (aI:1.1+1)'=a1:,.,+1(IfE*),1 E,l
of (2.2) and

I

(a:;‘

(2.3), it follows

from

(2.5) that

) ‘f’ -i

I

lx'+.' 111,

GE,, $'+O,

m>lfi,

cd++’) “’ I d'+lm

1 GE,, d?+O. I

Since the errors of the coefficients, described in (2.4) and (2.6), were caused by different operations, they can be estimated in one step, on the whole, by the expressions

I-where

(2.7)

=‘+I m

In>

E,-E,+E,--4.2-l. We shall use the splitting

of the matrix

A' into submatrices

and use the notation

(2.8) It follows

from the theorem

and from

(2.7) that

(1+E,)-'~(CT'(G'+')'(I)/ {U'G'+'U}=&1+/& (2.9) PG'+'U#O. Otherwise IY(G'+')'U=O. Thus, in one step of the direct run for matrix A, the rounding errors reduce to the following: the error, estimated by expression (2.9), is introduced into the matrix @"'-'-' and the coefficient a,+,,,+,remains to be calculated in the matrix (A"-')'with an error of E, in accordance with (2.5). 2. Consider the problem of the total effect of the errors introduced into the matrices W"_" I I=0 , 1,..., n-2 to solve Eq.cl.2). The products

if

A=QV. where Q and V are lower and upper form

A-

II

triangular V’

yr.*-r

0

VW!

(2.10)

matrices

II II -

Q”!,‘,:’

respectively,

I”

can be represented

in the

(2.11)

v',Vnvr are upper triangular submatrices. where Q,@-' are lower triangular submatrices and Since the expansion (2.10) is unique, if the values of the diagonal element of matrix V (see matrices the representation /5/) are to be determined, and for symmetric positive-definite A=L'L is possible, it then follows from (2.11) that (2.1’)

104 and A=M'+G',

l-1, 2,...,n-l,

(2.13)

G' is defined in (2.8). The matrix Ml is non-negative since it is represented in (2.12) by the product of a certain matrix and the transposed matrix, and G’ is non-negative by virtue of the fact that the matrix w"-' is positive-definite (see /5/J. It follows from (2.13) and (2.9) that

where

(I+&)-‘d(U’A”U)

A”=M:+‘+

/ (U’AU)
Since during the direct run the n-i matrix w is calculated, the perturbations (2.14) used here reduces to the relation where

(G’fl) 1,

(2.14)

the first approximation

of (2.15)

j(U'A'U)/(li'AU)-lI<(n-I)&. A’ is an initial matrix with the total error. In accordance with corollary 2, it follows from (2.15) that (U,-&‘)‘A

As can be shown /7/, it follows

(2.1(i)

(n-l)%Ts’U~‘AL’o.

(Co’-U,,<

from (2.16) that (2.17)

i!U,-U,'// /Ilv.ii9(n-I)EI(RIA)~.,

where R is the Rayleigh ratio ~=(U,rAU,)/(~U,~~'and I is the minimum eigenvalue of matrix A. 3. The error of the solution, calculated by rounding-off the diagonal elements of the matrix d"-'- which must be considered independently of para. - can be estimated using (1.25) where, obviously, E3 from (2.5) must be taken as E and s must equal one. Summing (1.25) from (2.17), we obtain a complete estimate of the error of the solution from the rounding errors for transformations of the matrix in the first approximation of perturbation theory: IlLI,-U,"JI / j/U,IIG (n-l)E,(R/ii)‘“+(I+*‘-)li”-E,= (2.18) (/t(n-l)(Rli.)“‘+(l+a’~)k”‘)2-‘. 4. We shall estimate the error of the solution from the rounding errors when completing the uirect and reverse run with the vector of the right-hand side. When solving triangular sets of equations it is recommended that the operation of accumulating scalar products with algorithm double the accuracy /5, 6/ be used. Within the framework of the modified-method /l/, one can implement this recommendation with a direct run for the vector of the right-hand A"-' side, having previously constructed the matrix B by separating rows of the matrix into diagonal terms, with the following calculations j-i--L (2.19)

+ z b,iK' . r:-' ==f*" i-1 The error of the solution, caused by rounding ji"-',as shown in /6/, corresponds perturbation of the diagonal element of the matrix N in the expansion A=NA”-‘: which

is equivalent

to the perturbation

In,i'/n,,-lj12-1, of the diagonal element

of the matrix

to the

L’: (2.20)

~z,,‘/l,,-lI~~2-‘. When calculating the elements of matrix B the equivalent of the matrix L' satisfies the relation

I (Z/)‘/ 1,;-11<2-‘, which,

together with

perturbation

of the elements

1,'#0, i+j,

(2.20), reduces to the estimate

(see 3, para.2) (2.21)

II(L’) ‘-LT~I/(ILTila2”*k’“2-’

where

h=max(i-j+i),

1.1

Lj’PO.

The estimate of the error for the reverse run with accumulation is described by expression (2.20), from which it follows that

of scalar products

(2.22)

[IL”-LII / IILil
/6/

of l., par-a.2 separately for the forward and reverse run, and on the (2.21) and (2.22) for II&il//jfJil, we arrive at the relations IIA,crl!;llU,ll~lr"'2-L, Il~,UlllllU,ll~2”‘a”‘k’“h’“2-‘, /IA,U+A,UIl/(!U[/G (lf2"'a"'h'") k”Z-‘.

Summing the errors obtained of the solution of Eq.Cl.2)

(2.23)

by (2.18) and (2.23), we obtain a final estimate of the error to the first approximation of perturbation theory:

((U.'-U,lJ /IlUolj=Z(4(n-l) (R/A)'"+ (2+a”+2’“h’“a”)k”)Z-‘. In the special

case, if a-1

and also, consequently,

(IUu'-U,il /~~U,jj~(4(n-l)+(3+2"h")k"')2-' If the diagonal elements, and also the vector of the right-hand from (2.24). calculated and stored with double the accuracy, the expression for estimating can easily be proved, takes the form IjU/-Uo/i

(2.24)

R/iL=i, we have

/IiU~;~~d(rL-l) (R/h)"'2-'.

side, are the error,

as

105

3. Supplement.

1. Let us prove expression

(1.22). Given AUo=I:,

where A is a positive-definite

matrix

A=L’L,

and L is an upper

diagonal

, a=

L'Y=/*',LC,=)

matrix;

Here the spherical norm is used, and (3.1) holds according to /5, p.65/. o,',(I,'>. . -3o,‘W, and Following /5/, we shall denote the eigenvalues of the matrix A by the corresponding orthogonal set of eigenveetors by u,,...,U, and we shall determine the vectors T,: T,=LU, / IILU,ll. (3.2) It is known /5/ that Il~ll== max

~lL~li~llXll=Il~~~lI~II~~ll-o~.

(3.3)

.z+o It follows from (3.3) that (T,. T,)-O,i=2, then, as is easy to prove,

3,...,n. Indeed,

if j is obtained such that (T,,I',)#@ such that the following relation holds:

7 can be selected

Considering, then, the influence of the operators A and L on the vectors in the subspace etc., i.e. the vectors orthogonal to (I,, we conclude that (T,, T,)=O,i-3,4,...,n 1‘,, i=i, 2,..., n form an orthogonal set. Taking (3.2) into account it follows that LU,=o,T’, Using

the representation

U,=Z,p,U,,

L’T=oiU,.

we arrive

(3.4)

at the relations

(3.5) If the expression

for a, is substituted

into

IlUdl From

(3.6) it follows

(l-22), we can reduce

(3.6)

that ll~*a’ll~I!‘~Il

Substituting

(1.22) to the form

/ lY[l~ilUJ” i IiFll”.

the expressions

for the norms

from

(3.7)

w. (3.5) in

(3.7) and using

(c,b,+...+~.b.)~~(c,'+...+c.')(f~,'+...+b,'),

r,=W*ll,,

the Cauchy

inequality

b,=,%,

we can satisfy ourselves of the correctness of (3.7) and, therefore, of (1.22). 2. Suppose matrix A, satisfying conditions (1.11, is expanded in products of the triangular matrices A=L’L and the errors of the coefficients of matrix L satisfy condition (1.10). We shall show that I( L’-Li!,/jiLillg2”‘h”‘E. i.e.,h

1s an indicator of the tape structure From (1.10) it follows that

l,,+O,

h= max (j-i+l), G of matrix

(3.8)

L and, therefore,

of matrix A. (3.9)

(11L'-Lll,llL'--&) i (IlLllrllL/l~)~E2 by virtue of the relations

(see /5/j

IlLIll -

max z

’

I

IU,

ll~ll~=m~s

x

,

lI,jl.

Since /l/

the following holds: (3.to)

Ir;lla.~2ilr,li*. The inequality connecting IlLI!, and (lLl/t is further necessary. Suppose implemented when .i=e.As is easy to prove, using the well-known inequality ((C,*+...+C.*)/n)'",the following relation holds:

~~~II~x.ll~/ll~.il~~c where .Y. is a vector whose single non-zero it follows from (3.11) that

11,.1-IILII,.

component

is the e-th.

IILll,a~“Wil~. The following

relation

m;ls,P,I1,,/is I(c,+ ...+C.)lnlG

13.11) Since

llLlllS!,LX.& / nix&, (3.12)

/5/ is also valid: IJL’-L(i,l~iIL’-LII,IiL’-LIJ,.

(3.13)

106 Substituting correct.

the estimate

(3.101,

(3.12) and

(3.13) into

(3.9),

we can show that (3.8) is

REFERENCES 1.

2.

3. 4. 5. 6. 7.

REPYAKH V.V., Use of the modified exclusion method to solve certain grid equations. Zh. vychisl. Mat. mat. Fir., Vo1.22, No.3, pp.634-645, 1982. P.EPYAKB V.V., Increasing the stability of Gauss's method in problems of structural mechanics. In: Computers in the analysis and design of buildings. Issue V. Kiev: Zonal'nii NII tipovovo i eksperim. proektirovaniya grazhdanskikh zdanii, pp.126-131, 1976. MIKHLIN S.G., Course in mathematical physics. Moscow, Nauka, 1968. MIKHLIN S.G. Numerical realization of variational methods. Moscow, Nauka, 1966. WILKINSON J.H., The algebraic eigenvalue problem. Moscow, Nauka, 1970. VOEVODIN v.v., The calculational basis of linear algebra. Moscow, Nauka, 1977. REPYAKH V.V., The influence of changes in the stiffness of elements of the change in the deformed state of an elastic hinge-bar system in a position of equilibrium. In: Computers in the analysis and design of buildings. Kiev: Budivel'nik,pp.l54-165, 1972.

Translated

U.S.S.R. Comput.Maths.Math.Phys.,25,2,pp.106-111,1985 Printed in Great Britain

by H.Z.

0041-5553/85 $lO.OO+O.CXl Pergamon Journals Ltd.

THE REGULARIZATIONOF CERTAIN METHODS OF MINIMIZING OF HIGH ORDER WHEN THE INITIAL DATA ARE INACCURATE* F.P. VASIL'EV Regularized versions of Steffenson's methods and other methods of high order for the problem of minimization when the function to be minimized and the functions specifying the type of equations are known, together with the errors, are proposed. Conditions for the regularization parameter and the penalty coefficient to agree with the parameters of the method and with the errors guaranteeing convergence of the methods with respect to the norm to the solution with minimum norm are obtained. Let us consider

the minimization

problem 1(~)+inE,

C=(r:Elf:g,(u)
u=U,

(1)

1=1, 2,...,,n;gz(u)=O, i=m+l,...,s),

are functions which are defined where J(lO, gr(a),...,gs(lc) on the real Hilbert space H. We shall assume that

We shall allow for the limitations

P(u)=

(2) using the penalty

[max(g.(u); ,-I

Oil’+

Frechet

(2) differentiable

1(1z)=I.)ZQ.

U.=(U=U:

I.=X/(u)>--m.

and doubly

(3)

function

k,(u) I p7

p>2.

i-l"+,

U, VEH by (u, v), the norm in H by We shall denote the scalar product of the elements and the space of the linear bounded operators mapping H into H by P(H+IZ). (~rc~~=
/'(u. u)=I~(u), similarly, Tikhonov's

suppose function

u.ue.1;

P'(rc. 1.)is a separated difference for problem

for

P'(u).

We shall introduce

A.N.

(11, (2):

T,(u) =I( 1L)+‘4,2( a) +c+l12,

u=H,

k=O, i,....

Then

TI'((L)=J'(a)+AIP'(u)+2a*u, T,'(u, v)=I'(tr, v)+A,P'(u, v)+BadI, where E is a unit operator on H. Problem (11, (2) is, generally

speaking,

an ill-posed

u,v=H, problem,

and to solve it one must