The SIMPLE algorithm

The SIMPLE algorithm for the veloity-pressure oupling has been used in this work. It

wasdeveloped by Patankar and Spalding [66℄ and has sine then been rened by a number

of authors.

Variantsrejoiing in the names of SIMPLEC [94℄, [79℄, SIMPLEX [79℄, SIMPLEN [92℄,

SIMPLER [66℄ and PISO [39℄ aim to improve the ouplingof the momentum and pressure

equations via minor modiationsof the SIMPLE algorithm.

Theshemeisapreditor-orretormethod,withaninitialestimatefortheveloityeld

from the Navier-Stokes equations being orreted with the ontinuity equation to fore the

onservation of mass. The predition and orretion operationsare enlosed in an iterative

loopwhihonverges togivea solutionthat satises all the equationsin the system.

The initialsheme by Patankar and Spaldingwasfor a staggered Cartesian mesh, with

the veloity values being loated at the faes of nite volume ells, and the pressure,

tem-perature and other salar variables being loated at the ell enters. Rhie and Chow [80℄

extended the method to use olloated grids, where the veloities and the other variables

are all loated at the ell enters, and this has been further developed by Peri¢ and other

authors ([71℄, [28℄). Suh a grid allows an easier onversion tonon-Cartesian meshes. Here

it willbe desribed forCartesian meshes.

The rst stage of the alulation proess is the resolution of the disretized versions

of the momentum equations 2.2 using the urrent estimate of the pressure eld, and using

a ell fae mass ux that is interpolated from the urrent estimate of the veloity eld

(this interpolation proedure is disussed in more detail in Setion 2.6.4). The momentum

equations 2.2 are in the same general form as the generi transport equation. For a given

mass uxeld and pressure eld they an bedisretized intoequationsof the form:

a _P u ^∗ _P + X

n = nb

a _n u ^∗ _n = S _x + A _x dp dx a P v _P ^∗ + X

n = nb

a n v _n ^∗ = S y + A y

dp

dy

^(2.50)

a P w ^∗ _P + X

n = nb

a n w ^∗ _n = S z + A z

dp dz

where

u ^∗ _p

^,

v ^∗ _p

^and

w _p ^∗

^{are the new estimates of veloity in the x, y and z. axis,and nb refers}

to the neighboring ells.

A _x , A _y

^and

A _z

^{represent the areas of the ell faes normal to the}

x, y

^and

z

^{axis. The pressure gradient an be found by} interpolating the pressure eld at the ell faes using a linear interpolation, and then approximating the gradient aross the

ell with aentered dierene as:

dp

dx ≈ p e − p w

∆x , dp

dy ≈ p n − p s

∆y ,

^(2.51)

dp

dz ≈ p t − p b

∆z ,

for a regularmesh where

∆x, ∆y

^and

∆z

^{are the ell}dimensions.

After the alulationof the veloity eld estimates, the ellfaeveloities an be

inter-polatedfromtheir values atthe ellenters, and the ellfaemass uxes anbealulated.

Forthe easternfaeof aellthe veloity normaltothe faeis

u e

^{,whilstthe faehasanarea}

A e

^{. The mass ow aross the faeis:}

m e = ρA e u e

^(2.52)

In generalthis interpolatedveloity eld willnot bemass onserving (i.e.,itwillnot havea

disrete divergene of zero), and so willnot satisfy the disretized version of the ontinuity

equation:

m e − m w + m n − m s + m t − m b = 0,

^(2.53)

whih an be writtenfor the fae veloities ona Cartesian meshesas:

ρA e u e − ρA w u w + ρA n v n − ρA s v s + ρA t w t − ρA b w b = 0,

^(2.54)

Wethereforewishtoalulateaorreted veloityeld

u ^∗∗

^,

v ^∗∗

^and

w ^∗∗

^{that ismass}

onserv-ing, together with a orresponding pressure eld

p ^∗∗

^{. We do so by adding a veloity and}

pressure orretion tothe originalestimation of the veloity and pressure elds:

u ^∗∗ = u ^∗ + u ^′ ,

v ^∗∗ = v ^∗ + v ^′ ,

^(2.55)

w ^∗∗ = w ^∗ + w ^′ , p ^∗∗ = p + p ^′ ,

where a dash

′

signiesthe orretion eld.

The expressions inequation 2.55 are substituted into the

u

^{equation inequation 2.50}

a _P (u ^∗ _P + u ^′ _P ) + X

n = nb

a _n (u ^∗ _n + u ^′ _n ) = S _x + A _x d d x

(p + p ^′ ),

^(2.56)

and the sum of the neighboring veloity terms approximated by:

X

n=nb

a _n (u ^∗ _n + u ^′ _n ) ≈ A _x dp ^′

dx

^(2.57)

whih should be validas

p ^′ → 0

^and

u ^′ → 0

^{. Subtrating the momentum equationgives an}

expression relatingthe orretion pressure and veloity eld toeah other,

a _P u ^′ _P ≈ A _x dp ^′

dx ,

^(2.58)

or

u ^′ _P = A x

a P

dp ^′ dx , v _P ^′ = A y

a P

dp ^′

dy ,

^(2.59)

w _P ^′ = A z

a _P dp ^′

dz ,

where similar approximations have been made for the

x, y

^and

z

^{omponents of veloity.}

By interpolating the expressions in equation 2.59 to the faes of the ell, the orreted ell

fae veloities normalto the faeare given by:

u ^′ _e = A _xe a P e

p ^′ _E − p ^′ _P δx e

,

v _n ^′ = A yn

a P n

p ^′ _N − p ^′ _P δy n

,

^(2.60)

w ^′ _t = A zt

a P t

p ^′ _T − p ^′ _P δz n

,

with the

a P

^{terms being} approximated atthe faes by a linear interpolation

a P e = a _{P P} + a _{P E}

2 ,

a P n = a _{P P} + a _{P N}

2 ,

^(2.61)

a _{P t} = a _{P P} + a _{P T}

2 ,

In this interpolation

a P P

^{is the}

a P

^{termin the equation for the ell}

P

^{, whilst}

a P E

^{is the}

a P

^{term in the equation for ell}

E

^. Substituting the equations in 2.55 into the disretized ontinuity equation 2.54 yields

A e (u ^∗ _e − u ^′ _e ) − A w (u ^∗ _w + u ^′ _w ) + A n (u ^∗ _n + u ^′ _n )

−A s (u ^∗ _s + u ^′ _s ) + A t (u ^∗ _t + u ^′ _t ) − A b (u ^∗ _b + u ^′ _b ) = 0

(2.62)

Usingthe expressions for

u ^′

^{fromequation2.60and fatorizingyieldsfollowingequation}

for the pressure orretion:

b P p ^′ _P + b E p ^′ _E + b W p ^′ _W + b N p ^′ _N + b S p ^′ _S + b T p ^′ _T + b B p ^′ _B = c

^(2.63)

where

b E = _a ^{A 2 e}

Pe , b W = _a ^{A 2 w}

Pw , b N = _a ^{A 2 n}

Pn , b S = _a ^{A 2 s}

Ps , b T = _a ^{A 2 t}

Pt , b B = _a ^{A 2 b}

Pb

b P = −(b E + b W + b N + b S + b T + b B ) c = _ρ ¹ (m w − m e + m n − m s + m t − m b )

This an be solved for the pressure orretion

p ^′

^{, whih is then used to update the ell}

enter and ell fae veloities using equations 2.60 and 2.55, the resulting fae veloities

satisfying the ontinuity equation. The pressure eld is updated using equation 2.55 and

then theproess isrepeated,withthe newveloityandpressure eldbeingusedtoalulate

the

u ^∗

^veloities.

TheSIMPLEalgorithmissummarizedinFigure2.1forthesolutionofathermallydriven

three dimensionalow. Toobtain theveloity,pressure and temperatureelds, anestimate

oftheveloityeld

u ^∗ , v ^∗

^and

w ^∗

^{isalulatedthroughequation2.50. Theellfaeveloities}

are theninterpolatedfromthe veloityeld andtheellfaemassuxes alulated

m ^∗

^{. The}

pressure orretion equation 2.63 is solved for

p ^′

^{, and the veloity, mass ux, and pressure}

elds are updated using equation 2.55. The resultingmass onserving veloity eld is then

used to solve any transport equations for auxiliary salar elds, suh as temperature and

speies onentration.

•

^{set initialelds for}

u ^∗ , v ^∗ w ^∗ , p

^and

T

•

interpolate tond ell faemass uxes F. repeat

solve Equation 2.50 for

u ^∗ , v ^∗ , w ^∗

interpolateto nd ellfae mass uxes

m ^∗

alulatepressure orretion

p ^′

^{from equation2.63}

update u, v,w, p and m using equation 2.55.

alulate

T

^{and any other salarelds.}

hek for onvergene. Ifonverged, halt.

Figure 2.1.: The SIMPLE veloity-pressure ouplingalgorithm

The divergene ofthe

m ^∗

^{massuxeld(alulatedasthesoureterminequation2.63)}

is usually used as a onvergene riteria. The system at this stage satises the momentum

equations (they havingjustbeen solved)whilst adivergene free

m ^∗

^{mass uxeldsignies}

that thesolutionalsosatisestheontinuityequation. After theupdateof the

m

^eldswith

the pressure orretion

p ^′

^{, the divergene should be equal to zero, a non-zero divergene}

signifying an inorretly solved pressure orretion equation. However, the orretions to

the veloity elds means that the momentum equations may no longer be satised, and so

the algorithm repeats. To aid the onvergene of the method, the veloity, pressure and

salar eld updates an be underrelaxed using some relaxation parameter. There are two

obvious methods of under-relaxation, either by relaxing the update:

u = u ^∗ + α u u ^′ , v = v ^∗ + α v v ^′ , w = w ^∗ + α w w ^′ ,

p = p ^∗ + α p p ^′ ,

(2.64)

where

α u

^,

α v

^,

α w

^and

α p

^{are the relaxation parameters for the veloity and pressure elds}

respetively,orby relaxingthe diagonaloeient ofthe linearsystem foreahvariable,for

example the u veloity equation 2.50 being modied to:

a P

α _u u ^∗ _P + X

n = nb

a n u n = S x + A x

dp

dx

^(2.65)

The pressure equation relaxation should always be of the form of equation 2.64, and

no relaxation should be used in the update of the fae veloities and mass uxes. Unlike

the transport equation, the pressure equation is always diagonally dominant and does not

require the under-relaxationof the diagonal. By avoidingthe relaxation ofthe fae veloity

updates,amassonservinguxeldisalulated,whihensurestheonservationofenthalpy,

momentum and other properties. It an also be noted that if a relaxation of the form of

equation 2.65 is used for the momentum equations, then a non-relaxed value of

a P

^should

beused in the pressure orretionand veloityinterpolation operations.

TomodelatransientproblemtheiterativeproessoutlinedinFigure2.1isarriedoutat

every timestep, using the veloity,pressure and salar elds fromthe previoustime step as

the initialguess for the values atthe new time step. This an be quitetime onsuming and

more eient time stepping proedures, whih will not be disussed here, are used instead

[28℄.

Forsteady stateproblems, theSIMPLEouplingshemean beonsideredasa

pseudo-transient proess, with an impliit alulation of the momentum equations being orreted

via an expliit pressure orretion proess, eah iteration of the sheme orresponding to

a pseudo time step. When modeling a transient ow, eah time step omprises a number

of pseudo-transient time steps to obtain the onverged solution for the physially real time

step.

Im Dokument Numerical Simulation of Reacting Flows under Laminar Conditions with Detailed Chemistry (Seite 30-35)