View source for Large-particle method

A method for computing compressible flows of a continuous medium [[#References|[1]]]. The method is based on splitting the original differential equations in accordance with the physical processes they represent (see [[#References|[3]]]). It can be used to solve systems of evolution equations. The procedure may establish the existence of steady-state solutions. The large-particle method is a development of Harlow's method of  "particle-in-cell" .

The method is widely used to investigate aerodynamic flows, diffraction problems, transonic flows, interaction phenomena of radiation and matter, etc.

The difference scheme of the method may be illustrated by considering the motion of an ideal compressible gas (the continuity, the momentum and the energy equation):

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575701.png" /></td> <td valign="top" style="width:5%;text-align:right;">(1)</td></tr></table>

Here <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575702.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575703.png" /> denotes time, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575704.png" /> the density, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575705.png" /> the velocity, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575706.png" /> the total specific energy, and <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575707.png" /> the pressure. To close the system (1) one uses the equation of state

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575708.png" /></td> </tr></table>

where

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l0575709.png" /></td> </tr></table>

is the specific internal energy.

The solution process for an evolution system (1) is divided into three chronological stages, each of which consists of three steps: an Eulerian, a Lagrangian and a final one. One first considers the variation of the internal state of a subsystem — a  "large particle"  (Eulerian step), and then — the motion of this subsystem with the interior state of the subsystem left unchanged (Lagrangian and final steps).

<img style="border:1px solid;" src="https://www.encyclopediaofmath.org/legacyimages/common_img/l057570a.gif" />

Figure: l057570a

Eulerian step. The domain of integration is covered by a stationary (Eulerian) difference grid of arbitrary form (to abbreviate the exposition, a rectangular grid in a two-dimensional (planar) domain is considered, see Fig.). At this step of the calculation only quantities relating to the cell as a whole are varied, while the fluid is assumed to be momentarily constrained. For this reason, convective terms such as <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757010.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757011.png" />, corresponding to translational effects, are dropped from the equations (1). In the remaining equations of (1), <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757012.png" /> is brought forward outside the differentiation symbol, and the equations (1) are solved for the derivatives of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757013.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757014.png" /> and <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757015.png" /> with respect to the time:

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757016.png" /></td> <td valign="top" style="width:5%;text-align:right;">(2)</td></tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757017.png" /></td> </tr></table>

An elementary finite-difference approximation (central differences) yields the following expressions:

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757018.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757019.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757020.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757021.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757022.png" /></td> </tr></table>

Quantities with fractional subscripts relate to cell boundaries, e.g.

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757023.png" /></td> </tr></table>

<img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757024.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757025.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757026.png" /> are intermediate values of the flow parameters, obtained on the assumption that the field <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757027.png" /> is  "frozen"  in the layer <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757028.png" />. Although the scheme of the Eulerian step in this form is unstable, the scheme as a whole is stable if the next steps are suitably formulated. The Eulerian step may be made stable, for example, by introducing elements of the [[Integral-relation method|integral-relation method]]. When this is done the integrands are approximated in the direction parallel to the axes of the body (see Fig.), i.e. as in scheme 1 of the integral-relation method: The original system of equations is taken to be in integral form, and the following integrals are approximated:

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757029.png" /></td> </tr></table>

Lagrangian step. Now one calculates the mass flow through the cell boundaries at time <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757030.png" />. It is assumed throughout that the mass of the large particle is in motion only owing to the velocity component normal to the boundary. Thus,

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757031.png" /></td> </tr></table>

etc. The angular brackets define the parameters <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757032.png" /> and <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757033.png" /> at the cell boundaries. The choice of these quantities is of paramount importance, since it strongly influences the stability and accuracy of the computation. Various difference representations are possible for <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757034.png" />: of different orders of accuracy, with or without allowance for the flow direction, central differences, ZIP-approximations, etc. The momentum (energy) flows are equal to the products of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757035.png" /> and the corresponding value of the velocity (total specific energy). Approximations have been worked out not only for mass flows, but also for momentum and energy flows.

Final step. At this step one finds the final fields of the Eulerian flow parameters at the time <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757036.png" />. The equations of this step are conservation laws for the mass <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757037.png" />, the momentum <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757038.png" /> and the total energy <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757039.png" />, relating to a given cell (large particle) in finite-difference form <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757040.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757041.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757042.png" />. The final values of the flow parameters <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757043.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757044.png" /> in the next time layer are computed according to the following formulas (the direction of the flow is from left to right and upwards):

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757045.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757046.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757047.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757048.png" /></td> </tr></table>

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757049.png" /></td> </tr></table>

The conservativity and total divergency of the difference scheme (divergent-conservative scheme) are guaranteed by the equation for the total energy <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757050.png" />. At the final step (when a discrete model of the medium is being used) one should carry out an additional calculation of the density; this smoothens out fluctuations and increases the accuracy of the computations. Combining the different representations of the steps, one obtains a series of difference schemes; this series, which constitutes the large-particle method, may be used for a broad range of numerical experiments.

The large-particle method may be interpreted from various points of view: the splitting method, the mixed Euler–Lagrange method, computation in local Lagrangian coordinates (Eulerian step) with scaling on the previous grid (Lagrangian and final steps), the difference notation for conservation laws for a fluid element (large particle), and the Eulerian difference scheme.

The boundary conditions are formulated by introducing series of fictitious cells (so that every computation point becomes an interior point and the same algorithm is maintained for all cells). One layer is sufficient for the first-order approximation scheme, two layers are sufficient for the second order, etc. For example, consider the problem of calculating the flow around axially-symmetric or plane bodies with generators of arbitrary form (see Fig.). The above difference formulas are valid for interior cells of the field, surrounded on all sides by fluid, and for cells bordering on the solid body, the contour of which coincides with the cell boundaries.

There are two possible approaches to the computation of flow around bodies using finite-difference methods: A computation in the coordinates <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757051.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757052.png" />; or an introduction of fractional cells (see [[#References|[2]]]). In the first case one encounters difficulties in treating bodies with fractures and concavities. The second approach is free from these disadvantages.

In the case of fractional cells, the boundary conditions on the body are formulated, as in the case of whole cells, by introducing fictitious cells. Inside the body one has a layer of fictitious cells bordering on the fractional cells. To determine the parameters of the gas in these fictitious cells, one constructs the normal from the centre of each fictitious cell <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757053.png" /> to the contour of the body; if <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757054.png" /> is the point of intersection of the normal and the contour, let <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757055.png" /> be the tangent at <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757056.png" /> (see Fig.). Now let <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757057.png" /> be the cell in the field of the flow, symmetric to <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757058.png" /> with respect to <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757059.png" />. The gas-dynamical parameters <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757060.png" /> in <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757061.png" /> are determined by  "weighting" , <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757062.png" /> (<img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757063.png" />), where the summation is over those cells <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757064.png" /> such that part of the area <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757065.png" /> lies in the cell <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757066.png" />. In order to formulate conditions stating that the body is impermeable, one more parameter has to be determined for each fictitious cell — <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757067.png" />, the inclination of the radius vector that cuts the body contour at <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757068.png" />. If one is using adhesion conditions (both velocity components change sign across the surface of the body), this additional determination of <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757069.png" /> is not necessary. The parameters of the gas in the fictitious cell <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757070.png" /> then are

<table class="eq" style="width:100%;"> <tr><td valign="top" style="width:94%;text-align:center;"><img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757071.png" /></td> </tr></table>

The boundary conditions for a body whose contour coincides with the cell boundaries constitute a particular case of the boundary conditions presented here.

For each fractional cell (see Fig.) one needs to know five geometrical characteristics: <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757072.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757073.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757074.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757075.png" />, and <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757076.png" />, where <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757077.png" /> is the ratio of the volume of the fractional cell to the volume <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757078.png" /> of the whole cell, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757079.png" /> is the fraction of the area of the side <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757080.png" /> open for the flow, etc.

The existence of a solid boundary within a cell implies two special features: It shifts the mass centre from the geometrical centre of the cell, bringing it closer to the boundary; and it reduces the actual size of the cell. With regard to both whole and fractional cells, all flow parameters refer to the centre of mass. It is between the centres of mass that one interpolates the gas-dynamical functions. In whole cells, the centre of mass is either identical to the geometrical centre of the cell (the plane Cartesian coordinate system) or near it (the cylindrical coordinate system). In actual computations, the difference, even for the row of cells along the axis, is at most <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757081.png" />. This does not substantially distort the results of the computations. If fractional cells are introduced in an appropriate way, the displacement of the centre of mass relative to the geometrical centre does not increase this error. The reduction in the effective size of the cells is a more serious matter. To avoid violation of the stability condition <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757082.png" />, where <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757083.png" /> is <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757084.png" /> or <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757085.png" />, those cells with <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757086.png" /> should be combined with adjacent whole cells within the flow, and the resulting combinations should be computed using the formulas for fractional cells. In this case, the geometrical dimensions of an enlarged cell will not be less than those of a whole cell: <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757087.png" />, <img align="absmiddle" border="0" src="https://www.encyclopediaofmath.org/legacyimages/l/l057/l057570/l05757088.png" />, so that the question of the stability of fractional cells is eliminated.

In the planar case, the geometrical characteristics of fractional cells may be determined by direct measurement. In the axially-symmetric case one needs an additional computation, incorporating the distance of each fractional cell from the axis of symmetry. The difference formulas for fractional cells are obtained by a slight modification of the difference formulas for whole cells.

The difference schemes of the large-particle method have been investigated (in regard to approximations, viscosity effects, stability) with the help of differential approximations (see [[#References|[4]]]). The method has also been generalized to the spatial (three-dimensional) case.

====References====
<table><TR><TD valign="top">[1]</TD> <TD valign="top">  O.M. Belotserkovskii,   Yu.M. Davydov,   "The method of large particles in gas dynamics. Numerical experiments" , Moscow  (1982)  (In Russian)</TD></TR><TR><TD valign="top">[2]</TD> <TD valign="top">  Yu.M. Davydov,   "Calculation by the  "coarse particle"  method of the flow past a body of arbitrary shape"  ''USSR Comp. Math. Math. Phys.'' , '''11''' :  4  (1972)  pp. 304–312  ''Zh. Vychisl. Mat. Mat. Fiz.'' , '''11''' :  4  (1972)  pp. 1056–1063</TD></TR><TR><TD valign="top">[3]</TD> <TD valign="top">  G.I. Marchuk,   "Methods of numerical mathematics" , Springer  (1975)  (Translated from Russian)</TD></TR><TR><TD valign="top">[4]</TD> <TD valign="top">  Yu.M. Davydov,   "Differential approximations and representations of difference schemes" , Moscow  (1981)  (In Russian)</TD></TR></table>