Christoph Borgers.

A domain decomposition Laplace solver for internal combustion engine modeling online

. (page 1 of 2)
Online LibraryChristoph BorgersA domain decomposition Laplace solver for internal combustion engine modeling → online text (page 1 of 2)
Font size
QR-code for this ebook


, , j^nj-ytstv.'.Ttv'g-nvi*:




TECHNICAL REPORT



A Domain Decomposition Laplace Solver
For Internal Combustion Engine Modeling

Christoph Borgers t
OlofB. Widlundt



Technical Report #31 5
August 1987



3 0)
rtfO
ft) I—'



n
o

3

cr
c
en
rr
(-■•
O
3

ro

D



> 03

o

Qj M
O OD

3 m nl
0) i-( ol

H- CO sj
3 ' -TJ
cnl

Qi n o

fO 3-
O f^
O ( - •
CO

rr
O



<

O 3

O

i-h CO 13
O H- s*
i-i rr

H-

O

3



I



Ull



o



EW YORK UNIVERSITY




^ Department of Computer Science
Courant Institute of Mathematical Sciences

' .,251 MERCER STREET, NEW YORK, N.Y. 10012



A Domain Decomposition Laplace Solver
For Internal Combustion Engine Modeling

Chrlstoph Borgers t
OlofB. Widlundt



Technical Report #31 5

August 1987



t Supported in part by the Applied Mathematical Sciences Subprogram of the Office
of Energy Research, U.S. Department of Energy, under contract DE-AC03-76SF00098,
at the Lawrence Berkeley Laboratory

t Supported in part by the National Science Foundation under grant NSF-DCR-
8405506 and by the US Department of Energy under contract DE-AC02-76ER0377-V
at the Co u rant Mathematics and Computing Laboratory.



A Domain Decomposition Laplace Solver for Internal Combustion Engine Modeling

Christoph Borgers *
• Department of Mathematics

University of Michigan
Ann Arbor, MI 4 ,9

and

OlofB. Widlund"^

Courant Institute of Mathematical Sciences

New York University

New York, NY 10012

Keywords: Laplace solver, domain decomposition, multigrid methods

AMS subject classification: 65N20, 65F10

Abbreviation of title: Domain Decomposition Laplace Solver



* Supported in part by the Applied Mathematical Sciences Subprogram of the Office of Energy
Research. U.S. Department of Energy, under contract DE-AC03-76SF00098, at the Lawrence Berke-
ley Laboratory.

** Supported in part by the National Science Foundation under grant NSF-DCR-8405506 and by the
U.S. Department of Energy under contract DE-AC02-76ER03077-V at the Courant Mathematics and
Computing Laboratory. ^



Abstract: Discrete elliptic problems can often naturally be partitioned into subproblems correspond-
ing to subregions into which the region has been partitioned or from which it was originally assem-
bled. Fast conjugate gradient methods have been developed to handle the continuity conditions, across
the interfaces, that need to be satisfied by the solution of the entire problem.

In this paper, we focus on a special family of problems that arises in combustion engine model-
ing. We demonstrate that efficient algorithms can conveniently be built using standard software to
solve elliptic problems on simple regions. We give a convergence analysis which is much simpler
than the general theory for domain decomposition algorithms. We also show that the use of multigrid
cycles to solve the subproblems inexactly can be very efficient.



A Domain Decomposition Laplace Solver for Internal Combustion Engine Modeling

1. Introduction

The interest in domain decomposition methods has increased considerably in recent years, to a
large extent because of their promise on parallel computers. There are frequent conferences; cf.
Glowinski and Periaux (1987). Much has been learned on how to design optimal or almost optimal
iterative algorithms to handle the interaction between solvers on the subregions. In this paper, we will
examine the use of one of these methods for a particular problem arising in combustion engine model-
ing. We note, however, that the ideas explored here can be adopted to other classes of problems in
fluid dynamics etc.. By focusing on a quite specific problem, we are able to demonstrate that the
theory sometimes can be simplified, that the algorithms can be speeded up by using special features,
and that standard software can be used to assemble a solver on relatively simple but non-trivial
regions.

We describe a Laplace solver on a domain Q of the shape shown in Fig. 1. This solver will be
used 10 compute the potential flow component in a random vortex method applied to a problem with a
moving piston; see Sethian (1987). Related work is described, e.g., in Sethian (1984), Sethian and
Ghoniem(1986).

In Fig. 1, the large rectangle fi^^' represents a cylinder, and the small attached quadrilaterals,
the union of which we denote by ii*", represent so-called inlets. There are no inlets attached to the
lower side of fi^'. The number of inlets attached to the left, right and upper sides of £i*^\ their posi-
tions, and the angles which they form with ii^' are left variable in our program. Sethian's calcula-
tions are time dependent with the lower side of Q'^' moving. It is a fractional step method, and one of
the steps requires the solution of a Neumann problem of the fonn

-A(|) = Oonn (1.1)

|^ = ^onaa ' (1.2)

an
Here is a velocity potential. The Neumann data g and the position of the lower side of Q*^',

representing the top of the piston, vary with time. The lowest inlet always remains above the lower



-2-

Equations (1.1) and (1.2) are discretized with piecewise linear finite elements on a quasi-
uniform mesh as described in Section 2. The region has re-entrant comers, and this is known to limit
the regularity of the solution, cf. Grisvard (1985). As a consequence, the discretization will yield only
a crude approximation to 0, unless the mesh width is very small. On the other hand, the resulting sys-
tem of linear equations has a relatively simple structure and can be solved very efficiently. We will
focus on this case in this paper, but discuss possible remedies through local refinement of the mesh in
the last section. One of the domain decomposition methods studied by Bj^rstad and Widlund
(1984,1986), the so-called "excellent" or Neumann- Dirichlet algorithm, will be used. A careful
motivation of this choice is given in those papers. For the problem at hand, we use a fast Poisson
solver for the rectangular region and band Cholesky for the relatively narrow inlets. The general con-
vergence analysis of this method requires tools firom the theory of elliptic boundary value problems;
see Bj^rstad and Widlund (1986), Widlund (1987b). However, for the geometry and meshes under
consideration, a much simpler convergence analysis, using elementary mappings and reflections, can
be given. This results in quantative and realistic bounds; see Section 4, where results of numerical
experiments are also presented.

Recently there has been an increased interest in the use of inexact solvers for the problems
defined on the subregions, see Bramble, Pasciak and Schatz (1986) and Widlund (1987a). In the fifth
section, we demonstrate that the exact solver for Q°^ can be replaced by a multigrid cycle. The
number of outer iterations grows modestly, but the cost per step is considerably decreased. We also
give bounds for the rate of convergence of this variant of the algorithm. In the final section, we collect
a number of ideas that could considerably further enhance the performance of the iterative method as
well as the accuracy of the discrelizalior..

The code used in the numerical experiments is available from the authors.
2. Triangulation of the region



-3-

The region is the union of a rectangle ii*^' and N,^u small quadrilaterals called inlets. Q'^' is
covered by a regular, rectangular mesh, with mesh widths h^ in the x-direciion and hy in the y-
direciion. We use logically rectangular meshes on the inlets, see Fig. 2, with common mesh points on
the interface between any mlet and the rectangle.

The triangulaiion of the rectangle is obtained by cutting each cell of the rectangular mesh into
two triangles. We do not need to specify the direction of these cuts, since the same five point formula
results when we use piecewise linear finite elements.

We triangulate an inlet by defining a bilinear mapping y, from a rectangle R, onto the inlet, /,,
as indicated in Fig. 2. On the rectangle R, a regular, rectangular mesh is defined. While the mesh
width in one coordinate direction is determined by the mesh on the large rectangle, the other can be
chosen freely. We triangulate R, by dividing each of its mesh cells into two triangles. This triangula-
tion induces a thangulation of the inlet as shown in Fig. 2.

3. A preconditioned conjugate gradient method

Let Ae/?*"^ be a symmetric, positive semidefiniie matrix, and let beR" be a vector orthogonal
to the kernel of A. Let AeR"*' be symmetric and positive semidefiniie, with ker(A)c ker(i4). There
are several possible preconditioned conjugate gradient algorithms for

Au = b (3.1)

using the preconditioner A, which are mathematically equivalent but algorithmically different The
version given below is particularly useful for domain decomposition methods. The matrix A only
appears in a vector-matrix multiplication with A-A. In our application, this is much cheaper than mul-
tiplying with A. In addition, we need a procedure which gives us the product of the Moore-Penrose
pseudo-inverse A of A with an arbitrary vector.

Preconditioned conjugate gradient algorithm:

Letuf°':-g€/?".



-4 -
Replace ^'"^ by its orthogonal projection onto the range of -4.

^(0).= ^(0)
-(0) -■>■ (0)

g ■=A g^">



-(0) .(0)

d ■.= g



For /=0, 1,2....:

aO' := [,0)^|C.)



^^''((A-A)d^'+d(^':J



(/)



^0*1) :=^0)_aO)[(^_^)j^\^0)]



Replace ^ '^'^ by its orthogonal projection onto the range of A.



pO).=



^(^••■"^oO+l)



i^^'^^^



^0*1) :=gO-i)+|30)^



0)



-0*1 -0*1) om^O)

The vectors u^' converge to a solution of (3.1).

The orthogonal projections onto the range of A have no effect in exact arithmetic, but these
steps are sometimes necessary to obtain convergence in floating point arithmetic.

4. The domain decomposition algorithm

In this section, we describe the Neumann-Dirichlet domain decomposition method as given in
Bj^rstad and Widlund (1984,1986). Because of the simple structure of our problem, we can give a
simplified convergence analysis. We also present numerical results.



4.1 Notation and preliminary remarks

The finite element approximation of eqs. (1.1), (1.2) can be written as a system of linear equa-
tions of the form



-5



Kx = b,



(4.1)



with K, X and b partitioned as follows:



K



Ku K,l




f -




'^1


a: 22 ^^23


, X =


i2


. ^ =


^2


Ki-i K2i AT 33




£3




^3



(4.2)



Here the subscript 1 refers to the nodes of the inlets not belonging lo the interfaces between the inlets
and the rectangle, the subscript 2 to the nodes in the rectangle not belonging to the interfaces, and the
subscript 3 the nodes on the interfaces.



The entries of the stiffness matrix K are integrals of the form



JV4)Vydi.



(4.3)



where on Q'-^^ The quadratic form on the right-hand
side of (4.23) is the Dirichlet integral of


1

Online LibraryChristoph BorgersA domain decomposition Laplace solver for internal combustion engine modeling → online text (page 1 of 2)