- Email: [email protected]
The numerical solution of elliptic differential equations
U.S.S.R. Comput.Maths.Math.Phys., Printed in Great Britian
~01.24,No.5,pp.
0041-5553/84 $lO.OO+O.OO 0 1985 Pergamon Press Ltd.
158-X5,1984
THE NUMERICAL SOLUTION OF ELLIPTIC DIFFERENTIAL EQUATIONS*
N.YU. BAKAEV
A method of using difference schemes to solve initial and boundary value problems in the case of a linear differential equation in partial derivaWhen tives of the elliptic type in a cylindrical region is considered. this method is used, the solution of the starting (initial value or boundary value) problem of the elliptic type reduces to the solution of The a Cauchy problem for a parabolic equation of the same dimensions. transition to a parabolic problem is realized in such a manner that the application of special methods of computational economization (alternating direction type methods) is possible when difference methods are employed. 1. Introduct ion. At the present time an approach based on a difference approximation to the initial problem and on the subsequent numerical solution of the resulting system of algebraic equations is the generally accepted approach to the solution of elliptic boundary value problems At the same time, iterative methods, based on an analogy using high-speed digital computers. with the method of regulation /l/, are the most general methods for solving a system of Actually, they imply a transition from an difference equations for an elliptic problem. elliptic problem to a parabolic problem at the price of an increase in the dimensions of the In the case of multidimensional problems, in spite of the techniques which problem by unity. have been specially developed for computational economization(alternating-direction, factorization, and overall approximation methods /2, 3/), the procudure for solving a system of difnet is even laborious for modern digital computers. On the ferenceequatibnsfor a fine mesh other hand, by using specific characteristics of the scheme (three-point equations, or fivepoint equations), it is possible to use the matrix pivotal condensation method, which is wellknown to be extremely laborious and, for which, this laboriousness can only be reduced for a Finally, the application of Fourier and reduction methods prenarrow range of problems. supposes both a special structure for the difference scheme and a requirement concerning the number of points in it. Moreover, the laboriousness of the reduction method is proportional to the effort required in inverting the matrices which arise diring the course of the calculation. Another approach is also possible: the solution of the elliptic problem in a cylindrical region is represented in the form of a sum of the direct and inverse problems for a certain parabolic equation and, moreover, the dimensions of the variable space remains the same. However, if the local operator figuring in the initial structure also possesses a simple structure as a differential operator, then the operator in the corresponding parabolic problem will not be local. It will have an exceedingly complex structure and, what is most significant, it will not lend itself to separation into one-dimensional operators. It will be shown below that the above-mentioned difficulties can be overcome by a special organization of the computational process where the initial problem will be a parabolic problem which has one dimension less than the problem which is solved by the method of regulation (i.e. iterative methods).
2. Bounded solutions of an elliptic equation. Let us consider anellipticpartial differential equation We shall write this equation in operator form. dWdt*=Au,
in a functional
t=(O, -),
Banach
space E. (2.1)
where u is the required function which belongs to the space E for each t, t is an explicitly selected variable, and A is an operator which acts on the remaining variables, which are not The operator A also takes account of the boundary conditions apertaining explicitly shown. Let us formulate the following problem: it is required to to the corresponding variables_ find the solution of (2.1) which satisfies the initial condition u(O)=ao and the condition
of the boundedness
(2.2)
of the norm of the solution
Iln(t)lJ,6+f=censt as
t-cm.
12.3)
Let the operator A satisfy the following conditions. 1) A is a linear closed operator with a dense domain of definition D(A). 2) the spectrum of the operator A belongs to the set S=(h; larghlz&/2), and the operator
lZh.vychisl.Mat.
mat.Fiz.,24.10,1538-1547,1984
158
159
(-A) generates the subgroup {T(--A, &)), .+O, of class C.. 3) in any closed sector which belongs to the complement of S, the estimate llR& -4)llz< M/IX] on the resolvent of the operator A is satisfied, where the constant M may depend on the choice of sector. is not an eigenvector of the operator A. 4) the point h=O A wide class of differential operators which satisfy conditions (l)-_(4) can be constructed using the results obtained in /4--6/. Above all we note that, if A satisfies conditions (l)-_(4), problem (2.1)-(2.3) is equivalent to the Cauchy problem for a parabolic equation /7, 8/. duidf=-A’;%,
t=(O, m),
u (0) =uo,
(2.4)
where the operator A’” is a fractional power of the operator A (see /9/l. the problem reduces to finding the value of the operator of the subgroup element no: u(t)=T(-A’“,
In other words, T(-A”,t) on the
t)u,.
(2.5)
Reasons have already been cited in the introduction as to why the direct application of difference methods to problem (2.4) gives rise to difficulties. To avoid these difficulties we use the following relationship /lo/ which was employed in /7/ in a conceptually similar situation: T(-A’h , t) u-
where
Vu&,
I;&s)T(-A,s)uds 0
ft,~~(s)-t(2x")-Ls-Kexp (-t*/&), s>O.
Hence,
C-0,
(2.6)
by virtue of 12.5) and (2.6)
.(t)-I;,&)B(a)ds, 0 where
g(s)
is the solution
of the Cauchy dglds==-Ag,
(2.7)
problem ~(0,
m),
(2.81
g (0) =J&
with the operator A, which possesses "good" properties: simplicity of structure, a local nature, and the possibility of being separated into one-dimensional operators. Economic difference schemes for solving multidimensional problems /2, 3/ may be used to solve problem (2.8). Let us now consider the final stage of the computational process which involves evaluation of the integral (2.7). If the functions g(s) in problem (2.8) are successively calculated on the layers s,,i=O, I,...,N (sp=O, s,,+), then it is possible to supplement the definition of the tabular function z(s) gt,i-o, 1,..., N to include the continuous function and, for example, to write
E (s-s,-,)] aw-r,[n-t-t
(2.9)
Q,(s),
1-S I
ii(t)-= jft,%(S)Bbws,
(2.10)
0
where g, is the computed value of g(s) on the i-th layer, Q,(s) is the indicator of the i-1, the interval is [so,&]) and Z(t) is the approximate value of interval (a,+ a,] (when the required function By calculating (2.10), we obtain n(t) r calculated as indicated above. (2.11) I
*a
(-;)
r,S[
g,-,+~(s-s,_,)]s-“exp
ds=
(2.11)
f-L .,.I
~~{~~(-~,~)-r(-~,~)]+ i-l 4 7 where
(
gc-I- gEsM)[
r(+$r
(&f-)])
r(a,y) is the incomplete gamma-function /II/. Let an error 6 arise in the solution of problem
,
(2.8). 1.e.
llg(s*)-gtllr~~, i=O, I,...,N. Then, by assuming that the function beyond the limits of this error
g(s)
constructed
Ilg(s)-B(s)Iw~~
usrng formuia
:2.:! does not gc
160 we obtain
i.e.
u(t) has the same error as g(s). It is obvious that (2.11) can be reduced
to the form
(2.12) where m,(t) are the weights with which the values g, are summed. Hence, relationship (2.12) may be taken as a basis for the approximate solution of (2.1)-(2.3). Let us evaluate the amount of effort required to carry out the computational process We shall estimate this, as is usal, using the number of multiusing the proposed scheme. plicative operations. It will consist of the number of operations at the stage when the parabolic equation (2.8) is solved and the number of operations during the computation using O(MN) operations are required to calculate the weights formula (2.12). m,(t)S where M is the number of mesh points with respect to the variable t. If the weights m,(t) have already been calculated then WNL operations are required to compute the sum (2.12), where L is the number of mesh points with respect to the variables on which the operator A acts. If we where K is the overall number of mesh points, then we obtain the expression put K=ML, KN+O(.&{&J) for the amount of computational effort required at the second stage. The amount of computational effort required at the first stage, where the parabolic problem is solved, is well-known. We emphasize that the dimensionality of this problem is reduced by one compared with the method of regulation which can be used if one starts out from the direct solution of the boundary value problem of Eq.(2.1). To calculate the solution on any one layer of the net Moreover, we note the following. with respect to t, the number of operations during the second stage is equal to LN+O (MN) and is independent of the site at which this layer is situated. If it is unnecessary to determine the solution of the problem at all mesh points but only at a number, J, of them, then the computational effort required at the second stage is equal to IN+O(MN). Subject to economization methods being used in the calculation of the parabolic problems, the amount of computational effort at the first stage is equal to O(LN) (see /2/j.
3.
Boundary
value
Let us now consider
problem
for an elliptic
a boundary
value problem t=(O,
d=uldP=Au, in the Banach space E, where us assume that the operator
b),
equation.
of the form
the operator A satisfiesconditions (l)-(4). Additionally, (-A) generates a subgroup of the negative type, i.e.
(-A’“)
/9/, the operator
IjT(-A’i,
also
possesses
t)IlaGfc
the analogous
where
t) [i,+n,]++A’“.
j’(-A”‘,
property (3.2)
o,>o.
exp (--NJ),
It has been shown in /9/ that, in this case, the solution in the form
u(t)=;
let
o>o.
(JT(-A, t)IlsGMexp (-at), Consequently
(3.1)
u(b) =ub.
u (0) =uo,
b-t)
of the problem
[&.+&I>
can be represented
06&b,
&=-A-"'u~', &-A-"QJ~‘. From the boundary
conditions, &+T(-A”,
from where
we obtain b) [ub+&] =uo,
the equation i&,+T(-A’“,
for determining
Giland iia:
b) [uo+&] =ub,
the relationships &=[Z+T(-A’“,
2b)] [I-T(A’“,
2T(-A'",b) [I-T(-A’“,
26) ]-‘ub,
Es=-23'(-A'", 6)[ I-T(-A”‘, [Z+T(-A’“,
26)]-‘uo-
26) ][I--T(-A”‘,
2b) ]-‘uo+ Sb)]-‘ub
follow. M, sxp(--2o,b)
161
u(t)=(Z’(-A’“, T(-A”‘, (-I’(-A”‘, 2,-A’“,
We now use relationship
2b)+T(-A’“.
t)[Z+T(-A”‘,
b-t)(T(-A”‘, i)[T(-.,t’“,
46)+...]-
b)+T(-A’“,
3b)+...])uo+
b)+T(-A’“,
36)+...]+
b--t)[Z+T(-A”,
Zb)+T(-A’“,
4b)+...l)o.+
(2.61, on the basis of which _ _
one may write
u(t)=j B&,s)g,(s)ds+j &(&s)g*(s)ds, 0 .
(3.3)
where
E,(t,S)==[fr.‘,:(S)+!~+~a.‘,,~S)t/,+~~.~/,(S)+...l[h_,.‘iiiS)+f‘b--l.‘/:(J‘) +-..I,
Bb(t, S)=-_[l,;b.",(S~+!l+lb.',,(S)+...1+ [h_l.~~~(S)ijJ1-,.'/,(S)+ . ..I. and
g.(s)
and
g*(s)
are
the solutions $=-Ag,,,
of the Cauchy
go(o)=uo
problems
dge
and z=-&,
g,(O)=u,
(3.4)
respectively. If problems (3.4) are calculated over the layers of the mesh s,,i=O, 1,.._, N, then, by using the procedure employed in Sect.2, it is possible to write the approximate solution of the boundary value problem (3.1) is the form
ii(t)=
r: I-0
[rn~(t)g,“+ml”(t)g,“l,
(3.5)
where g,' and g: are the solutions of (3.4) calculated on the layers si and m:(t) and ml”(t) are the corresponding weights used in the summation which can be expressed in terms of an incomplete gamma-function (they are not written out on account of their large size). Hence the procedure for solving the boundary value problem (3.11 reduces to the following. At the first stage, the parabolic problems (3.4) are calculated with the same operator A, and, at the second stage, the coefficients mj"(t) are calculated and the summation m<'(t) and is carried out using formula (3.5). We emphasize that the number of operations for the calculation of l&"(l) and mP(t) is equal to O(MN) where M and N have their previous meaning and all the reasoning concerning the laboriousness of the calculation which has been formulated in Sect.2 is transferred to the case of boundary value problem (3.11. Finally, in calculating Bo(& s) and &(l, s), the infinite series should be replaced by finite sums. It is an elementary matter to estimate the number of terms retained in a series since the explicit form of the function Another method is to use the f0 (s) is known. inverse Fourier transform of a function of the form {I-exp [-2b(io)‘.‘])-‘(esp
[ --t(io)‘“]-exp
[ -(26--t)
(io)‘:]).
(in the case of
&(t,s)). A more general boundary value problem may also be considered contain terms with corresponding derivatives /9/.
when
the boundary
conditions
4. The Cauchy problem for an elliptic equation. Let us investigate problem
the structure dWdt’=Au,
of a numerical t=(U, b),
scheme
for the solution
u(O)=uo, uI (0) L-u,,
of the Cauchy (4.1)
which we shall also consider in a functional Banach space E. We assume that the operator A satisfies conditions (l)-_(4) and the operator (-A) generates the semigroup of operators of negative type. {T(--4, 0) The fundamental characteristic of the Cauchy problem (4.1) is that it is ill-posed and it must be solved numerically on the basis of a family of operations irhich is Tikhonovregularized. The method from /13/ is employed below since it enables us to use the approach developed in the preceeding paragraphs of the present paper. The family of solutions of problem (4.1) which is regularized 1; /13/ has the form
u”(t) ==+-A’L,t)
“D.O
?.,,-Y
(~.--A-“~u,)++~
‘G. 2 j
162
c*(+)”
T(-A’“,
nq) (&+A-‘“u,),
PO,
li-”
where q is the regularization followinq manner:
parameter.
Let us now regroup
the latter expression
in the
u,(t)=~rT(-A.h.t)+exp(~)~ $($)“x n-0
I
T(-A”,nq)
kxp
($)
u,+
c
-A-V(-Aya,t)+
p$q+,
‘A-‘T(-A”.nq)] u,,
q>o.
n-0
The representation
is known /14/ for the operator A-" . BY computing the composition of the integrals the representation for the operator A-'T(-A",t):
(2.6) and
(4.3). it is possible
to obtain
_ A-“T(-As,
t)u=-
j r,:h(s)T(-A,s)u&
(4.4)
0 where ~,;~(s)=(xs)-'"exp(-t*/4s).By using expressions (2.6) and ized solution of problem (4.11 in the following manner:
_ u&)=
(4.4),
we represent
the regular-
m
jC&s)g,(s)ds+jc&,s)gmds,
(4.5)
0
0 where C,tt,s)=~[,,,,h(s)+exp($)C7;i-(~)
-
(-1)”
t
”
fnl~~h(s)]v
-rc,~b~+exp(~~~~ct, ‘rn.,lh(s)],
C,(f,s)--$[
“-0
and
g.(s) and
g,(s) are the solutions
of the Cauchy problems
g,(O)=u, and -$=-Ag,,
%=-Ago,
g4 (0) = UI
(4.6)
respectively. As was done above, by starting out from the fact that the calculations on problem (4.6) are made on the layers s,,i=O, I,...,N,the approximate regularized solution of the Cauchy problem (4.1) can be written in the form
and where g," and g,' are the solutions of (4.6) calculated on the layers sit and m.‘(q, t) m,l(n, f) are the weights used in the summation which are calculated using the functions fl.~,(S) C,(t,s) and C,(t,s) can be calculated either and rr.'h(S). As was done above, the functions by retaining a finite number of terms of the series or by using a Fourier transform. Instead of (4.2), a method may be used follows from the results presented in /15/. u.(t)=+
&I)T.1T
T(-A”,
t) [ uo-A-“%,]+
(-A’+)
[u,+A-“u,],
n=i,
2,.
h-0 In this case a finite sum occurs instead of an infinite series. All of the remarks concerning the laboriousness and efficiency of the computational process which were made above in the case of problems associated with an elliptic equation are valid in the case of the numerical scheme for solving the Cauchy problem (4.1).
Some concluding remarks. 1. In considering one of the problems being studied above, let a variable T be explicrtly
5.
selected from among the variables on which the operator A acts. E possess the property that, when space Ea of the space u=Ea ,
Additionally,leta certainsubthe element T(-A. S)U. which
163 T, satisfies the condition ]I(-A. s)~](~)=O when is realized as a function of the variable S>T, where 0 is to be understood as a function of the remaining variables (the variables apart from c) of the operator A which is equal to zero everywhere. In this case, if the initial conditions and the boundary conditions in problems (2-l)-_(2.3), (3.1), or (4.1) belong to the subspace E. then, instead of an infinite upper limit in formulae (2.7), (3.3), and (4.5), a limit equal to r may be set which also leads to economy in the computational operations. where A0 1sa certain As an example, let us consider an operator A of the form A-a/h+Ao Let us denote the elliptic operator, the variables of which are not explicitly indicated. set of variables of the operator A0 by X. Let E be the Banach space of the functions ~(7.2) which are strongly continuous with respect to r on [--, -1 and satisfy some conditions applying to the set of variables x which guarantee the Banach natrue of the space. We shall not actually state in detail which space E and which operator Ao are chosen but we shall assume that the necessary conditions for the applicability of all the arguments are satisfied. Since Let us select the subspace E. in the following manner: u-=Eo, if Y(T,2)-O when ~40. the operator A0 only acts on the variables x, it is obvious that andthat this function is equal to zero when S>T. Hence, the procedure of cutting off the upper limit in integrals (2.7), (3.3), and (4.5) is correct if the corresponding initial and boundary conditions belong to the subspace Ea which, in this case, is interpreted as a subclass of casual functions. 2. Everywhere above we have studied homogeneous equations. However, the results ob tained in the present paper can be extended to an equation of the form d*u/dt'-Au+f(t). f(l) is a function with values As an example, let us consider a c(t,r) for the boundary value problem case of equation (5.1) may be written
where
(5.1)
in the Banach space E. boundary value problem. By using Green's function /9/, the solution of the boundary value problem in the in the following manner:
By applying (2.6)and (4.4) as was done above, the operator c(t,r) may be expressed in terms of the semigroup T(-A,#). In the general case on such a pathway, however, it is difficulttoderive relationships which enable one to obtain economic computational formulae. Nevertheless, if it is assumed that f(r) is represented in the form f(t)--g(t)h,whereg(r) is a numerical function and h is a certain element of the space E (i.e. h is a function of the remaining variables), then we may write
in which the function p(',I) is represented by a certain integral which contains the functions and g(t). In this case the amount of effort required to carry out the computatrons 1:,a(~), rr.'/*(~) is the same as that for the problems considered above. 3. By applying operator calculus in an analogous manner, more general second-order equations may be treated:
A and B are and the corresponding problems in the case of these equations, if the operators functions of the one and the same operator. 4. It has been shown above that, by applying operator calculus, it is possible :o transfer from the solution of an initial value or boundary value problem for an elliptic differential equation to the solution of an initial value problem for a parabolic equation of the same dimensions (if we interpret the dimensions as the sum of the dimensions with respect co the explicitly selected and unslected variables). The construction which is pro_=osed here does, in fact, correspond to a certain way of constructing Green's function for the initial problem. This method simply involves representing Green's function with the help of operator calculus in terms of the solution of the parabolic problem. At the same time, we avoid having to calculate the trial coefficients (this stage of the pivotal condensation method is equivalent to the solution of a certain initial value problem which is, moreover, non-linear) or the analogous situation in the invariant embedding method of Bellman and Kallabi /16/ and we did not encounter any need for greater dimensions (as in the regulation method) in solving the initial value problem. All these difficulties are absorbed by the problem of calculating the numerical coefficients ft.,,, (8).&(t, I),Mt. a),co(t.I) and c,(t.s). We note that other methods which of Green's functions for second-order equaticns ir. cylindrical are known for the construction domains may turn out to be inapplicable. For example, the method described ir. /17,' f-r a boundary value problem and associated with the solution of two initial value problecs for d second-order equation is not always suitable since the Cauchy problem for a second-cr,d.er equation is, generally speaking, ill-posed and associated with numerical insril::;. The possibility of applying the scheme presented in the present paper rt~ursr~:;-- snouid
164
also be mentioned. 5. It is not easy to carry out a comparative analysis of the laboriousness of the proposed method and methods which already exist if only because the laboriousnessofaspecific calculation of a modification of it may depend in a different manner on the number of mesh points along the different coordinate directions (such examples have been cited in /l/l. In the present paper a net is constructed in an intermediate variable s which does not have analogoues in other methods and makes any comparison difficult It is possible from general If we consider the composition of a difference positions, however, to note the following. method with direct methods for solving the system of difference equations, then the matrix pivotal condensation method is exceedingly laborious in any case and, furthermore, requires Of the direct methods, the Fourier and that certain stability conditions are satisfied. reduction methods seem to be the most preferable. However, the latter methods are only applicable under rather strict conditions so that, for example, their use must be abandoned on changing to a net which is non-uniform with respect to t. On the other hand, the composition of a difference method with an iterative technique for solving the system of difference under quite general conditions and, inparticular equations (the regulation method) can be used when comparing non-uniform nets. However, this approach is equivalent to solving an evolutionary problem, the dimensions of which are one greater then that of the initial problem. Moreover, apart from the dependence on the methods which are used to solve the system of difference equations, the difference scheme itself must be compared on a sufficiently fine mesh in order that the initial continuous problem should be approximated fairly well (to say nothing of possible additional constraints in the case of conditionally stable schemes). The scheme presented in the present paper not only admits of non-uniformity of the mesh with respect to t but is also the most efficient (in the sense of Laboriousness) when this for a conventional difference scheme mesh is sparse, i.e. when the approximation conditions may even be exceedingly poorly satisfied. Moreover, at the stage of the calculation using formulae (2.7). (3.31, and (4.5) it is not obligatory to calculate the solution for all values of the variables which are not explicitly indicated if this is not dictated by the conditions In certain cases this may substantially reduce the laboriousness of the of the problem. calculations. It is very easy to construct a flow version of the method developed here, i.e. a scheme For this purpose it suffices to substitute the derivatives of the for calculating duldf. functions f,.Ob(J), ~~(t.~), B&, ~1.Cl@, 11and c,(t,S) with respect to t in formulae (2.7), (3.31, and (4.5) instead of the functions themselves; the derivatives are calculated in analytical form. We note that the approach developed in this paper enables one to avoid both the solution of a boundary value problem in the case of the formulation described in Sect.2 and the introduction of the corresponding error into the solution of the problem which is associated with the need to set the second boundary condition at a finite point (instead of at the point t--l in the numerical realization. REFERENCES (Metody 1. SAMARSKII A.A. and NIKOLAEV E.S., I4ethods of solving finite difference equations reshenie setochnykh uravneniil, Nauka, Moscow, 1978. 2. SAMARSKII A.A., Introduction to the theory of difference schemes (Vvedenie v teoriyu raznostnykh skhem), Nauka, Moscow, 1971. (Metody vychislitel'noi matematiki), 3. MARCHUK G.I., Methods of computational mathematics Nauka, Moscow, 1977. 4. SOLOMYAK M-Z., The analytical nature of the semigroup generated by an elliptic operator in L, spaces, Dokl. Akad. Nauk SSSR, 127, 1, 37-39, 1959. 5. SOLOMYAK M.Z., Estimate of the norm of the resolvent of an elliptic operator in L, spaces, Uspekhi mat. Nauk, 15, 6(96), 141-148, 1960. 6. KRASNOSEL'SKII 14-A. et al., Integral operators in spaces of summed functions (Integral'nye operproty v prostranstvakh summiryemykh funktsiil, Nauka, MOSCO;~, 1966. pulse in a nomogeneous medium with 7. BAKAEV N.YU., The propagation of an elelctromagnetic conductivity, Zh. vychisl. Mat mat. Fiz., 22, 3, 671-677, 1982. dif8. BAKAEV N.YU., The problem of the continuation of the solutions of a seccnd-order at the All-Union Institute for Scientific ferentialequationin Banach space, Deposited and Technical Information (VINITI), No.1753-80 DEP. 1980. (Annotated: Sib. mat. Zh., 22, 1, 227, 19811. 9. KREIN S.G., Linear differential equations in Banach space (Linei?ye differentsial'nye uravnenie v banchovom prostranstve), Nauka, Moscow, 1967. Berlin (1968) (Russian translation, 10. YOSHIDA K., Functional analysis, Springer-Verlag, Mir, MOSCOW, 19671. Functions (Russian translation, 11. ABRAMOVTIZ M. and STEGUN A., Handbook of Mathematical 1977). Nauka, Moscow, 12. KANTOROVICH L.V. and AKILOV G.P., Functional analysis (Funksticnal'nyi analiz), Nauka, Moscow, 1977. 13. BAKAEV N.YU. and TARASOV R-P., One method of solving the Cauchy problem for an elliptic equation, Sib. mat. Zh. 20, 6, 1198-1205, 1979. analiz), Nauka, lloscow, 14. KREIN S.G. (editor), Functional analysis (SMB) (Funksionai'nyi 1972.
165 trajectories of differential equations 15. BAKAEV N.YU., The problem of the continuationof in Banach space and its application to the equations of mechanics and electrodynamics with separation of the steady (vibrational) and decaying regimes, in: Problems in the modern theory of periodic motions (Problemy sovremennoi teorii periodicheskikh dvizhenii), No.5, 63-75, 1981, Udmursk. Gos. Univ., Izhevsk. equations, Academic 16. BELLMAN R. and ANGEL E., Dynamic programming and partial differential Press, 1972 (Russian translation, Mir, Moscow, 1974). 17. BAKHVALOV N.S.. Numerical methods (Chislennye metodyl, Vol.1, Nauka, Moscow, 1975.
Translated
U.S.S.R. comput.Maths.Math.Phys., Printed in Great Britain
Vol 24,~0.5,pp.165-170,1984
by E.L.S
0041-5553/84 $10.00+0.00 al985 Pergamon Press Ltd.
DIFFRACTION OF PLANE WAVES AT A HORIZONTAL HALF-PLANE IN A COMPRESSIBLE STRATIFIED FLUID* A.K. SHATOV The problem of the scattering of steady plane waves at the edge of a horizontal half-plane submerged in an exponentially stratified compressible fluid is considered. The wave propagation is assumed to be adiabatic but not necessarily isentropic. 1. The propagation of plane waves in a compressible stratified fluid. The problem of wave propagation is discussed with reference to a Cartesian coordinate system (.r,, s,,z) in which the vector of the acceleration due to gravity g is oriented in the negative direction of the &-axis. Small adiabatic motions are described by the system of equations
(la)
(lb) (see /l/l. Here v=(u,, VI, US) is the velocity vector of the fluid particles, p is the disturbance of the density caused by the fluid motion, p is the prssure perturbation, c is the adiabatic velocity of sound, and p&(z) and pO(z) are the stationary distributions of the pressure and density, connnected by the relation The first equation in (la) VPo=P&' describes the law of variation of the momentum of the fluid particles; the second is the equation of continuity, a consequence of the equation of state. Henceforth we discuss the problem assuming that pl=e-**, f,>O (for simplicity, we assume po(O)=I), and a 8, g and c are constants in the adopted approximation of the equation of state. We consider the function \y using the formula Y==pe". By eliminating p, p and v, system (11 can be reduced to the following equation for Y:
(2) where A, is the Laplace operator in the variables (z,,zt) and N is the Weisel-Brent frequency defined by the relation p ---gpd(a)/pe(z)-ga/cZ=2~g-g*/c’ (assuming N*>O, i.e. the fluid The components stratification is stable: this occurs under practical conditions; see /l/J. of the velocity Y, pressure p and density p are connected with function \I!by the expressions (3a)
(3bl where
p-p-w/g while IMlCg. Consider the solution of Eq.(2) with a time-dependence e-'Of. For the amplitude function Y(z,,zr,a,t) (we retain for it the symbolY) we have --$yr+(i-$)A.Y+($-@‘)Y-0.
of the
(41
~01.24,No.5,pp.
0041-5553/84 $lO.OO+O.OO 0 1985 Pergamon Press Ltd.
158-X5,1984
THE NUMERICAL SOLUTION OF ELLIPTIC DIFFERENTIAL EQUATIONS*
N.YU. BAKAEV
A method of using difference schemes to solve initial and boundary value problems in the case of a linear differential equation in partial derivaWhen tives of the elliptic type in a cylindrical region is considered. this method is used, the solution of the starting (initial value or boundary value) problem of the elliptic type reduces to the solution of The a Cauchy problem for a parabolic equation of the same dimensions. transition to a parabolic problem is realized in such a manner that the application of special methods of computational economization (alternating direction type methods) is possible when difference methods are employed. 1. Introduct ion. At the present time an approach based on a difference approximation to the initial problem and on the subsequent numerical solution of the resulting system of algebraic equations is the generally accepted approach to the solution of elliptic boundary value problems At the same time, iterative methods, based on an analogy using high-speed digital computers. with the method of regulation /l/, are the most general methods for solving a system of Actually, they imply a transition from an difference equations for an elliptic problem. elliptic problem to a parabolic problem at the price of an increase in the dimensions of the In the case of multidimensional problems, in spite of the techniques which problem by unity. have been specially developed for computational economization(alternating-direction, factorization, and overall approximation methods /2, 3/), the procudure for solving a system of difnet is even laborious for modern digital computers. On the ferenceequatibnsfor a fine mesh other hand, by using specific characteristics of the scheme (three-point equations, or fivepoint equations), it is possible to use the matrix pivotal condensation method, which is wellknown to be extremely laborious and, for which, this laboriousness can only be reduced for a Finally, the application of Fourier and reduction methods prenarrow range of problems. supposes both a special structure for the difference scheme and a requirement concerning the number of points in it. Moreover, the laboriousness of the reduction method is proportional to the effort required in inverting the matrices which arise diring the course of the calculation. Another approach is also possible: the solution of the elliptic problem in a cylindrical region is represented in the form of a sum of the direct and inverse problems for a certain parabolic equation and, moreover, the dimensions of the variable space remains the same. However, if the local operator figuring in the initial structure also possesses a simple structure as a differential operator, then the operator in the corresponding parabolic problem will not be local. It will have an exceedingly complex structure and, what is most significant, it will not lend itself to separation into one-dimensional operators. It will be shown below that the above-mentioned difficulties can be overcome by a special organization of the computational process where the initial problem will be a parabolic problem which has one dimension less than the problem which is solved by the method of regulation (i.e. iterative methods).
2. Bounded solutions of an elliptic equation. Let us consider anellipticpartial differential equation We shall write this equation in operator form. dWdt*=Au,
in a functional
t=(O, -),
Banach
space E. (2.1)
where u is the required function which belongs to the space E for each t, t is an explicitly selected variable, and A is an operator which acts on the remaining variables, which are not The operator A also takes account of the boundary conditions apertaining explicitly shown. Let us formulate the following problem: it is required to to the corresponding variables_ find the solution of (2.1) which satisfies the initial condition u(O)=ao and the condition
of the boundedness
(2.2)
of the norm of the solution
Iln(t)lJ,6+f=censt as
t-cm.
12.3)
Let the operator A satisfy the following conditions. 1) A is a linear closed operator with a dense domain of definition D(A). 2) the spectrum of the operator A belongs to the set S=(h; larghlz&/2), and the operator
lZh.vychisl.Mat.
mat.Fiz.,24.10,1538-1547,1984
158
159
(-A) generates the subgroup {T(--A, &)), .+O, of class C.. 3) in any closed sector which belongs to the complement of S, the estimate llR& -4)llz< M/IX] on the resolvent of the operator A is satisfied, where the constant M may depend on the choice of sector. is not an eigenvector of the operator A. 4) the point h=O A wide class of differential operators which satisfy conditions (l)-_(4) can be constructed using the results obtained in /4--6/. Above all we note that, if A satisfies conditions (l)-_(4), problem (2.1)-(2.3) is equivalent to the Cauchy problem for a parabolic equation /7, 8/. duidf=-A’;%,
t=(O, m),
u (0) =uo,
(2.4)
where the operator A’” is a fractional power of the operator A (see /9/l. the problem reduces to finding the value of the operator of the subgroup element no: u(t)=T(-A’“,
In other words, T(-A”,t) on the
t)u,.
(2.5)
Reasons have already been cited in the introduction as to why the direct application of difference methods to problem (2.4) gives rise to difficulties. To avoid these difficulties we use the following relationship /lo/ which was employed in /7/ in a conceptually similar situation: T(-A’h , t) u-
where
Vu&,
I;&s)T(-A,s)uds 0
ft,~~(s)-t(2x")-Ls-Kexp (-t*/&), s>O.
Hence,
C-0,
(2.6)
by virtue of 12.5) and (2.6)
.(t)-I;,&)B(a)ds, 0 where
g(s)
is the solution
of the Cauchy dglds==-Ag,
(2.7)
problem ~(0,
m),
(2.81
g (0) =J&
with the operator A, which possesses "good" properties: simplicity of structure, a local nature, and the possibility of being separated into one-dimensional operators. Economic difference schemes for solving multidimensional problems /2, 3/ may be used to solve problem (2.8). Let us now consider the final stage of the computational process which involves evaluation of the integral (2.7). If the functions g(s) in problem (2.8) are successively calculated on the layers s,,i=O, I,...,N (sp=O, s,,+), then it is possible to supplement the definition of the tabular function z(s) gt,i-o, 1,..., N to include the continuous function and, for example, to write
E (s-s,-,)] aw-r,[n-t-t
(2.9)
Q,(s),
1-S I
ii(t)-= jft,%(S)Bbws,
(2.10)
0
where g, is the computed value of g(s) on the i-th layer, Q,(s) is the indicator of the i-1, the interval is [so,&]) and Z(t) is the approximate value of interval (a,+ a,] (when the required function By calculating (2.10), we obtain n(t) r calculated as indicated above. (2.11) I
*a
(-;)
r,S[
g,-,+~(s-s,_,)]s-“exp
ds=
(2.11)
f-L .,.I
~~{~~(-~,~)-r(-~,~)]+ i-l 4 7 where
(
gc-I- gEsM)[
r(+$r
(&f-)])
r(a,y) is the incomplete gamma-function /II/. Let an error 6 arise in the solution of problem
,
(2.8). 1.e.
llg(s*)-gtllr~~, i=O, I,...,N. Then, by assuming that the function beyond the limits of this error
g(s)
constructed
Ilg(s)-B(s)Iw~~
usrng formuia
:2.:! does not gc
160 we obtain
i.e.
u(t) has the same error as g(s). It is obvious that (2.11) can be reduced
to the form
(2.12) where m,(t) are the weights with which the values g, are summed. Hence, relationship (2.12) may be taken as a basis for the approximate solution of (2.1)-(2.3). Let us evaluate the amount of effort required to carry out the computational process We shall estimate this, as is usal, using the number of multiusing the proposed scheme. plicative operations. It will consist of the number of operations at the stage when the parabolic equation (2.8) is solved and the number of operations during the computation using O(MN) operations are required to calculate the weights formula (2.12). m,(t)S where M is the number of mesh points with respect to the variable t. If the weights m,(t) have already been calculated then WNL operations are required to compute the sum (2.12), where L is the number of mesh points with respect to the variables on which the operator A acts. If we where K is the overall number of mesh points, then we obtain the expression put K=ML, KN+O(.&{&J) for the amount of computational effort required at the second stage. The amount of computational effort required at the first stage, where the parabolic problem is solved, is well-known. We emphasize that the dimensionality of this problem is reduced by one compared with the method of regulation which can be used if one starts out from the direct solution of the boundary value problem of Eq.(2.1). To calculate the solution on any one layer of the net Moreover, we note the following. with respect to t, the number of operations during the second stage is equal to LN+O (MN) and is independent of the site at which this layer is situated. If it is unnecessary to determine the solution of the problem at all mesh points but only at a number, J, of them, then the computational effort required at the second stage is equal to IN+O(MN). Subject to economization methods being used in the calculation of the parabolic problems, the amount of computational effort at the first stage is equal to O(LN) (see /2/j.
3.
Boundary
value
Let us now consider
problem
for an elliptic
a boundary
value problem t=(O,
d=uldP=Au, in the Banach space E, where us assume that the operator
b),
equation.
of the form
the operator A satisfiesconditions (l)-(4). Additionally, (-A) generates a subgroup of the negative type, i.e.
(-A’“)
/9/, the operator
IjT(-A’i,
also
possesses
t)IlaGfc
the analogous
where
t) [i,+n,]++A’“.
j’(-A”‘,
property (3.2)
o,>o.
exp (--NJ),
It has been shown in /9/ that, in this case, the solution in the form
u(t)=;
let
o>o.
(JT(-A, t)IlsGMexp (-at), Consequently
(3.1)
u(b) =ub.
u (0) =uo,
b-t)
of the problem
[&.+&I>
can be represented
06&b,
&=-A-"'u~', &-A-"QJ~‘. From the boundary
conditions, &+T(-A”,
from where
we obtain b) [ub+&] =uo,
the equation i&,+T(-A’“,
for determining
Giland iia:
b) [uo+&] =ub,
the relationships &=[Z+T(-A’“,
2b)] [I-T(A’“,
2T(-A'",b) [I-T(-A’“,
26) ]-‘ub,
Es=-23'(-A'", 6)[ I-T(-A”‘, [Z+T(-A’“,
26)]-‘uo-
26) ][I--T(-A”‘,
2b) ]-‘uo+ Sb)]-‘ub
follow. M, sxp(--2o,b)
161
u(t)=(Z’(-A’“, T(-A”‘, (-I’(-A”‘, 2,-A’“,
We now use relationship
2b)+T(-A’“.
t)[Z+T(-A”‘,
b-t)(T(-A”‘, i)[T(-.,t’“,
46)+...]-
b)+T(-A’“,
3b)+...])uo+
b)+T(-A’“,
36)+...]+
b--t)[Z+T(-A”,
Zb)+T(-A’“,
4b)+...l)o.+
(2.61, on the basis of which _ _
one may write
u(t)=j B&,s)g,(s)ds+j &(&s)g*(s)ds, 0 .
(3.3)
where
E,(t,S)==[fr.‘,:(S)+!~+~a.‘,,~S)t/,+~~.~/,(S)+...l[h_,.‘iiiS)+f‘b--l.‘/:(J‘) +-..I,
Bb(t, S)=-_[l,;b.",(S~+!l+lb.',,(S)+...1+ [h_l.~~~(S)ijJ1-,.'/,(S)+ . ..I. and
g.(s)
and
g*(s)
are
the solutions $=-Ag,,,
of the Cauchy
go(o)=uo
problems
dge
and z=-&,
g,(O)=u,
(3.4)
respectively. If problems (3.4) are calculated over the layers of the mesh s,,i=O, 1,.._, N, then, by using the procedure employed in Sect.2, it is possible to write the approximate solution of the boundary value problem (3.1) is the form
ii(t)=
r: I-0
[rn~(t)g,“+ml”(t)g,“l,
(3.5)
where g,' and g: are the solutions of (3.4) calculated on the layers si and m:(t) and ml”(t) are the corresponding weights used in the summation which can be expressed in terms of an incomplete gamma-function (they are not written out on account of their large size). Hence the procedure for solving the boundary value problem (3.11 reduces to the following. At the first stage, the parabolic problems (3.4) are calculated with the same operator A, and, at the second stage, the coefficients mj"(t) are calculated and the summation m<'(t) and is carried out using formula (3.5). We emphasize that the number of operations for the calculation of l&"(l) and mP(t) is equal to O(MN) where M and N have their previous meaning and all the reasoning concerning the laboriousness of the calculation which has been formulated in Sect.2 is transferred to the case of boundary value problem (3.11. Finally, in calculating Bo(& s) and &(l, s), the infinite series should be replaced by finite sums. It is an elementary matter to estimate the number of terms retained in a series since the explicit form of the function Another method is to use the f0 (s) is known. inverse Fourier transform of a function of the form {I-exp [-2b(io)‘.‘])-‘(esp
[ --t(io)‘“]-exp
[ -(26--t)
(io)‘:]).
(in the case of
&(t,s)). A more general boundary value problem may also be considered contain terms with corresponding derivatives /9/.
when
the boundary
conditions
4. The Cauchy problem for an elliptic equation. Let us investigate problem
the structure dWdt’=Au,
of a numerical t=(U, b),
scheme
for the solution
u(O)=uo, uI (0) L-u,,
of the Cauchy (4.1)
which we shall also consider in a functional Banach space E. We assume that the operator A satisfies conditions (l)-_(4) and the operator (-A) generates the semigroup of operators of negative type. {T(--4, 0) The fundamental characteristic of the Cauchy problem (4.1) is that it is ill-posed and it must be solved numerically on the basis of a family of operations irhich is Tikhonovregularized. The method from /13/ is employed below since it enables us to use the approach developed in the preceeding paragraphs of the present paper. The family of solutions of problem (4.1) which is regularized 1; /13/ has the form
u”(t) ==+-A’L,t)
“D.O
?.,,-Y
(~.--A-“~u,)++~
‘G. 2 j
162
c*(+)”
T(-A’“,
nq) (&+A-‘“u,),
PO,
li-”
where q is the regularization followinq manner:
parameter.
Let us now regroup
the latter expression
in the
u,(t)=~rT(-A.h.t)+exp(~)~ $($)“x n-0
I
T(-A”,nq)
kxp
($)
u,+
c
-A-V(-Aya,t)+
p$q+,
‘A-‘T(-A”.nq)] u,,
q>o.
n-0
The representation
is known /14/ for the operator A-" . BY computing the composition of the integrals the representation for the operator A-'T(-A",t):
(2.6) and
(4.3). it is possible
to obtain
_ A-“T(-As,
t)u=-
j r,:h(s)T(-A,s)u&
(4.4)
0 where ~,;~(s)=(xs)-'"exp(-t*/4s).By using expressions (2.6) and ized solution of problem (4.11 in the following manner:
_ u&)=
(4.4),
we represent
the regular-
m
jC&s)g,(s)ds+jc&,s)gmds,
(4.5)
0
0 where C,tt,s)=~[,,,,h(s)+exp($)C7;i-(~)
-
(-1)”
t
”
fnl~~h(s)]v
-rc,~b~+exp(~~~~ct, ‘rn.,lh(s)],
C,(f,s)--$[
“-0
and
g.(s) and
g,(s) are the solutions
of the Cauchy problems
g,(O)=u, and -$=-Ag,,
%=-Ago,
g4 (0) = UI
(4.6)
respectively. As was done above, by starting out from the fact that the calculations on problem (4.6) are made on the layers s,,i=O, I,...,N,the approximate regularized solution of the Cauchy problem (4.1) can be written in the form
and where g," and g,' are the solutions of (4.6) calculated on the layers sit and m.‘(q, t) m,l(n, f) are the weights used in the summation which are calculated using the functions fl.~,(S) C,(t,s) and C,(t,s) can be calculated either and rr.'h(S). As was done above, the functions by retaining a finite number of terms of the series or by using a Fourier transform. Instead of (4.2), a method may be used follows from the results presented in /15/. u.(t)=+
&I)T.1T
T(-A”,
t) [ uo-A-“%,]+
(-A’+)
[u,+A-“u,],
n=i,
2,.
h-0 In this case a finite sum occurs instead of an infinite series. All of the remarks concerning the laboriousness and efficiency of the computational process which were made above in the case of problems associated with an elliptic equation are valid in the case of the numerical scheme for solving the Cauchy problem (4.1).
Some concluding remarks. 1. In considering one of the problems being studied above, let a variable T be explicrtly
5.
selected from among the variables on which the operator A acts. E possess the property that, when space Ea of the space u=Ea ,
Additionally,leta certainsubthe element T(-A. S)U. which
163 T, satisfies the condition ]I(-A. s)~](~)=O when is realized as a function of the variable S>T, where 0 is to be understood as a function of the remaining variables (the variables apart from c) of the operator A which is equal to zero everywhere. In this case, if the initial conditions and the boundary conditions in problems (2-l)-_(2.3), (3.1), or (4.1) belong to the subspace E. then, instead of an infinite upper limit in formulae (2.7), (3.3), and (4.5), a limit equal to r may be set which also leads to economy in the computational operations. where A0 1sa certain As an example, let us consider an operator A of the form A-a/h+Ao Let us denote the elliptic operator, the variables of which are not explicitly indicated. set of variables of the operator A0 by X. Let E be the Banach space of the functions ~(7.2) which are strongly continuous with respect to r on [--, -1 and satisfy some conditions applying to the set of variables x which guarantee the Banach natrue of the space. We shall not actually state in detail which space E and which operator Ao are chosen but we shall assume that the necessary conditions for the applicability of all the arguments are satisfied. Since Let us select the subspace E. in the following manner: u-=Eo, if Y(T,2)-O when ~40. the operator A0 only acts on the variables x, it is obvious that andthat this function is equal to zero when S>T. Hence, the procedure of cutting off the upper limit in integrals (2.7), (3.3), and (4.5) is correct if the corresponding initial and boundary conditions belong to the subspace Ea which, in this case, is interpreted as a subclass of casual functions. 2. Everywhere above we have studied homogeneous equations. However, the results ob tained in the present paper can be extended to an equation of the form d*u/dt'-Au+f(t). f(l) is a function with values As an example, let us consider a c(t,r) for the boundary value problem case of equation (5.1) may be written
where
(5.1)
in the Banach space E. boundary value problem. By using Green's function /9/, the solution of the boundary value problem in the in the following manner:
By applying (2.6)and (4.4) as was done above, the operator c(t,r) may be expressed in terms of the semigroup T(-A,#). In the general case on such a pathway, however, it is difficulttoderive relationships which enable one to obtain economic computational formulae. Nevertheless, if it is assumed that f(r) is represented in the form f(t)--g(t)h,whereg(r) is a numerical function and h is a certain element of the space E (i.e. h is a function of the remaining variables), then we may write
in which the function p(',I) is represented by a certain integral which contains the functions and g(t). In this case the amount of effort required to carry out the computatrons 1:,a(~), rr.'/*(~) is the same as that for the problems considered above. 3. By applying operator calculus in an analogous manner, more general second-order equations may be treated:
A and B are and the corresponding problems in the case of these equations, if the operators functions of the one and the same operator. 4. It has been shown above that, by applying operator calculus, it is possible :o transfer from the solution of an initial value or boundary value problem for an elliptic differential equation to the solution of an initial value problem for a parabolic equation of the same dimensions (if we interpret the dimensions as the sum of the dimensions with respect co the explicitly selected and unslected variables). The construction which is pro_=osed here does, in fact, correspond to a certain way of constructing Green's function for the initial problem. This method simply involves representing Green's function with the help of operator calculus in terms of the solution of the parabolic problem. At the same time, we avoid having to calculate the trial coefficients (this stage of the pivotal condensation method is equivalent to the solution of a certain initial value problem which is, moreover, non-linear) or the analogous situation in the invariant embedding method of Bellman and Kallabi /16/ and we did not encounter any need for greater dimensions (as in the regulation method) in solving the initial value problem. All these difficulties are absorbed by the problem of calculating the numerical coefficients ft.,,, (8).&(t, I),Mt. a),co(t.I) and c,(t.s). We note that other methods which of Green's functions for second-order equaticns ir. cylindrical are known for the construction domains may turn out to be inapplicable. For example, the method described ir. /17,' f-r a boundary value problem and associated with the solution of two initial value problecs for d second-order equation is not always suitable since the Cauchy problem for a second-cr,d.er equation is, generally speaking, ill-posed and associated with numerical insril::;. The possibility of applying the scheme presented in the present paper rt~ursr~:;-- snouid
164
also be mentioned. 5. It is not easy to carry out a comparative analysis of the laboriousness of the proposed method and methods which already exist if only because the laboriousnessofaspecific calculation of a modification of it may depend in a different manner on the number of mesh points along the different coordinate directions (such examples have been cited in /l/l. In the present paper a net is constructed in an intermediate variable s which does not have analogoues in other methods and makes any comparison difficult It is possible from general If we consider the composition of a difference positions, however, to note the following. method with direct methods for solving the system of difference equations, then the matrix pivotal condensation method is exceedingly laborious in any case and, furthermore, requires Of the direct methods, the Fourier and that certain stability conditions are satisfied. reduction methods seem to be the most preferable. However, the latter methods are only applicable under rather strict conditions so that, for example, their use must be abandoned on changing to a net which is non-uniform with respect to t. On the other hand, the composition of a difference method with an iterative technique for solving the system of difference under quite general conditions and, inparticular equations (the regulation method) can be used when comparing non-uniform nets. However, this approach is equivalent to solving an evolutionary problem, the dimensions of which are one greater then that of the initial problem. Moreover, apart from the dependence on the methods which are used to solve the system of difference equations, the difference scheme itself must be compared on a sufficiently fine mesh in order that the initial continuous problem should be approximated fairly well (to say nothing of possible additional constraints in the case of conditionally stable schemes). The scheme presented in the present paper not only admits of non-uniformity of the mesh with respect to t but is also the most efficient (in the sense of Laboriousness) when this for a conventional difference scheme mesh is sparse, i.e. when the approximation conditions may even be exceedingly poorly satisfied. Moreover, at the stage of the calculation using formulae (2.7). (3.31, and (4.5) it is not obligatory to calculate the solution for all values of the variables which are not explicitly indicated if this is not dictated by the conditions In certain cases this may substantially reduce the laboriousness of the of the problem. calculations. It is very easy to construct a flow version of the method developed here, i.e. a scheme For this purpose it suffices to substitute the derivatives of the for calculating duldf. functions f,.Ob(J), ~~(t.~), B&, ~1.Cl@, 11and c,(t,S) with respect to t in formulae (2.7), (3.31, and (4.5) instead of the functions themselves; the derivatives are calculated in analytical form. We note that the approach developed in this paper enables one to avoid both the solution of a boundary value problem in the case of the formulation described in Sect.2 and the introduction of the corresponding error into the solution of the problem which is associated with the need to set the second boundary condition at a finite point (instead of at the point t--l in the numerical realization. REFERENCES (Metody 1. SAMARSKII A.A. and NIKOLAEV E.S., I4ethods of solving finite difference equations reshenie setochnykh uravneniil, Nauka, Moscow, 1978. 2. SAMARSKII A.A., Introduction to the theory of difference schemes (Vvedenie v teoriyu raznostnykh skhem), Nauka, Moscow, 1971. (Metody vychislitel'noi matematiki), 3. MARCHUK G.I., Methods of computational mathematics Nauka, Moscow, 1977. 4. SOLOMYAK M-Z., The analytical nature of the semigroup generated by an elliptic operator in L, spaces, Dokl. Akad. Nauk SSSR, 127, 1, 37-39, 1959. 5. SOLOMYAK M.Z., Estimate of the norm of the resolvent of an elliptic operator in L, spaces, Uspekhi mat. Nauk, 15, 6(96), 141-148, 1960. 6. KRASNOSEL'SKII 14-A. et al., Integral operators in spaces of summed functions (Integral'nye operproty v prostranstvakh summiryemykh funktsiil, Nauka, MOSCO;~, 1966. pulse in a nomogeneous medium with 7. BAKAEV N.YU., The propagation of an elelctromagnetic conductivity, Zh. vychisl. Mat mat. Fiz., 22, 3, 671-677, 1982. dif8. BAKAEV N.YU., The problem of the continuation of the solutions of a seccnd-order at the All-Union Institute for Scientific ferentialequationin Banach space, Deposited and Technical Information (VINITI), No.1753-80 DEP. 1980. (Annotated: Sib. mat. Zh., 22, 1, 227, 19811. 9. KREIN S.G., Linear differential equations in Banach space (Linei?ye differentsial'nye uravnenie v banchovom prostranstve), Nauka, Moscow, 1967. Berlin (1968) (Russian translation, 10. YOSHIDA K., Functional analysis, Springer-Verlag, Mir, MOSCOW, 19671. Functions (Russian translation, 11. ABRAMOVTIZ M. and STEGUN A., Handbook of Mathematical 1977). Nauka, Moscow, 12. KANTOROVICH L.V. and AKILOV G.P., Functional analysis (Funksticnal'nyi analiz), Nauka, Moscow, 1977. 13. BAKAEV N.YU. and TARASOV R-P., One method of solving the Cauchy problem for an elliptic equation, Sib. mat. Zh. 20, 6, 1198-1205, 1979. analiz), Nauka, lloscow, 14. KREIN S.G. (editor), Functional analysis (SMB) (Funksionai'nyi 1972.
165 trajectories of differential equations 15. BAKAEV N.YU., The problem of the continuationof in Banach space and its application to the equations of mechanics and electrodynamics with separation of the steady (vibrational) and decaying regimes, in: Problems in the modern theory of periodic motions (Problemy sovremennoi teorii periodicheskikh dvizhenii), No.5, 63-75, 1981, Udmursk. Gos. Univ., Izhevsk. equations, Academic 16. BELLMAN R. and ANGEL E., Dynamic programming and partial differential Press, 1972 (Russian translation, Mir, Moscow, 1974). 17. BAKHVALOV N.S.. Numerical methods (Chislennye metodyl, Vol.1, Nauka, Moscow, 1975.
Translated
U.S.S.R. comput.Maths.Math.Phys., Printed in Great Britain
Vol 24,~0.5,pp.165-170,1984
by E.L.S
0041-5553/84 $10.00+0.00 al985 Pergamon Press Ltd.
DIFFRACTION OF PLANE WAVES AT A HORIZONTAL HALF-PLANE IN A COMPRESSIBLE STRATIFIED FLUID* A.K. SHATOV The problem of the scattering of steady plane waves at the edge of a horizontal half-plane submerged in an exponentially stratified compressible fluid is considered. The wave propagation is assumed to be adiabatic but not necessarily isentropic. 1. The propagation of plane waves in a compressible stratified fluid. The problem of wave propagation is discussed with reference to a Cartesian coordinate system (.r,, s,,z) in which the vector of the acceleration due to gravity g is oriented in the negative direction of the &-axis. Small adiabatic motions are described by the system of equations
(la)
(lb) (see /l/l. Here v=(u,, VI, US) is the velocity vector of the fluid particles, p is the disturbance of the density caused by the fluid motion, p is the prssure perturbation, c is the adiabatic velocity of sound, and p&(z) and pO(z) are the stationary distributions of the pressure and density, connnected by the relation The first equation in (la) VPo=P&' describes the law of variation of the momentum of the fluid particles; the second is the equation of continuity, a consequence of the equation of state. Henceforth we discuss the problem assuming that pl=e-**, f,>O (for simplicity, we assume po(O)=I), and a 8, g and c are constants in the adopted approximation of the equation of state. We consider the function \y using the formula Y==pe". By eliminating p, p and v, system (11 can be reduced to the following equation for Y:
(2) where A, is the Laplace operator in the variables (z,,zt) and N is the Weisel-Brent frequency defined by the relation p ---gpd(a)/pe(z)-ga/cZ=2~g-g*/c’ (assuming N*>O, i.e. the fluid The components stratification is stable: this occurs under practical conditions; see /l/J. of the velocity Y, pressure p and density p are connected with function \I!by the expressions (3a)
(3bl where
p-p-w/g while IMlCg. Consider the solution of Eq.(2) with a time-dependence e-'Of. For the amplitude function Y(z,,zr,a,t) (we retain for it the symbolY) we have --$yr+(i-$)A.Y+($-@‘)Y-0.
of the
(41