A S able Hyb id Po en ial–SPH Technique o En o ce he Fluid
Incomp essibili y
JUAN J. PEREA
Uni e si y o Se ille
School o Compu e Enginee ing
A da. Reina Me cedes, s/n, 41012 Se illa
SPAIN
[email p o ec ed]
JUAN M. CORDERO
Uni e si y o Se ille
School o Compu e Enginee ing
A da. Reina Me cedes, s/n, 41012 Se illa
SPAIN
[email p o ec ed]
Abs ac : The SPH me hod has ex ensi ely used in luid low simula ion. Th ough SPH he luid is modelled by
a pa icles sys em whose mu ual in e ac ion is weighed by a unc ion, named ke nel unc ion, whose limi s de ine
he neighbou ing o each pa icle. In spi e o he high capabili ies o SPH o simula e complex en i onmen s, i
shows sho comings specially i he luid is subjec o high changes in he p essu e, he eloci y and he densi y
as occu in phenomena such as shock– ube, blas –wa e o in he bounda y and discon inui ies whe e he numbe
o neighbou pa icles is ela i ely low. In his case, he p essu e g adien is inaccu a e. As consequence, he
simula ion is ins able wi h an e a ic beha iou o pa icles. To a oid his p oblem, we p opose a hyb id echnique.
This one consis s in o mula ing he p essu e g adien om a po en ial de ined on each pa icles pai . Thus, he
p essu e g adien is immune o he low numbe o neighbou pa icles. Also, ou p oposal allows en o cing he
luid incomp essibili y. To show he imp o emen s ob ained we will ca y ou a se o simula ions.
Key–Wo ds: luid simula ion, SPH, s abili y, accu acy, po en ial–based.
1 In oduc ion
The luid low simula ion plays an essen ial ole in
se e al disciplines bo h he scien i ic and echnical.
F om an analy ic poin o iew, i is o mula ed by
Fluid Mechanics ha se he pa ial di e en ial equa-
ions desc ibing he luid low. Howe e , i s com-
plexi y only allows ob aining analy ic esul s o a
es ic ed se o cases wi h high spa ial symme y.
In mos cases, he nume ical me hods mus be used.
These ones a e s udied by Compu a ional Fluid Dy-
namic (CFD). Th ough CFD, he pa ial di e en ial
equa ions, a e ans o med in o a se o algeb aic
equa ions, whose solu ion allow simula ing he luid
low. One o hese me hods is Smoo hed Pa icle Hy-
d odynamics (SPH) which is widely used o i s adap -
abili y and e sa ili y [21, 32].
SPH ope a es on luid disc e iza ion h ough pa -
icles ha in e ac among hemsel es. Thus, i is based
on he Lag angian o mula ion. Desc ip i ely, each
pa icle depic s a disc e e po ion o he con inuum
luid and mo es wi h he low. In his me hodology,
he dynamic quan i ies a e i on each pa icle and i s
mu ual in e ac ion is weighed by a scala unc ion de-
penden on dis ance. This unc ion, known as ke -
nel unc ion, mus sa is y a se o analy ical p ope -
ies such as i mus be mono onically dec easing, i
mus be no malized and i s domain, named suppo ed
domain, should be limi ed. Each pa icle has go an
associa ed ke nel unc ion, whose suppo ed domain
limi s he numbe o neighbou pa icles.
The SPH me hod has highligh ea u es such as
he conse a ion o mass [13, 19, 21], i s adap abili y
and e sa ili y o simula e complex en i onmen s [19]
and he simpli ica ion o he o mula ion o sol e he
luid dynamic equa ions [9, 13]. Ne e heless, despi e
hese ad an ages, he simula ion by SPH can show
s abili y p oblems ha a ec he ealism o simula-
ion [3, 6, 28]. The e a e se e al ins abili ies ypes:
he local mixing ins abili y (LMI) [25], he hea ing
ins abili y [25], bu he mos ema kable is ela ed o
p essu e g adien [20, 24, 25]. The inaccu acy o p es-
su e g adien induces a ensile ins abili y [14, 18, 21]
o pai ing o pa icles [2, 3, 24].
The e is a consensus ha ela es he ins abil-
i y o p essu e g adien and he numbe o neigh-
bou ing pa icles, specially when his numbe is el-
a i ely low o he pa icles dis ibu ion is no smoo h
[3, 10, 21, 24, 28]. Acco ding he numbe o neigh-
bou ing pa icles dec eases, he inaccu acy o p es-
su e o ce is inc eased. Consequen ly, he pa icles
a e sca e ed and show e a ic beha iou . Usually, i
has been ega ded ha inc easing he adius o he
suppo ed domain he inaccu acy is educed. Ne -
e heless, he inc ease in he suppo ed domain e-
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
50
Volume 13, 2018
duces he esolu ion and, as consequence, he accu-
acy [3, 21, 32]. Fu he mo e, he e a e egions, such
as he bounda ies o luid o i s discon inui ies, whe e
he numbe o neighbou pa icles canno be inc eased
[28]. O he me hodologies ha e been de eloped o
a oid he p oblems associa ed o he low numbe o
neighbou ing. These echniques ei he e o mula e
he dynamic equa ions modelled by SPH [21, 25] o
posi ioning con ol poin s o educe he e ec s o he
ex e nal loads on he p essu e g adien [4, 33]. Ne e -
heless, hese p oposals show sho comings since can
iola e he conse a ion o momen um [3, 16], pa ic-
ula ly in bounda ies and discon inui y egions. Also,
hey a e only e ec i e unde speci ic condi ions, i
hese condi ions a e no sa is ied, an e oneous o e -
damping is induced [3].
To a oid hese sho comings ela ed o he p es-
su e o ce and he incomp essibili y en o cemen ,
we p opose a hyb id echnique po en ial–SPH. Suc-
cinc ly, ou p oposal is based on a o mula ion o a
conse a i e po en ial ha allow us o de i e he p es-
su e o ce. The cons ain s o ob ain his o ce has
he aim he imp o ing he s abili y and accu acy spe-
cially in he discon inui ies and bounda y egions. On
he o he hand, he o mula ed p essu e o ce allow us
o en o ce incomp essibili y wi hou he need o use
a p edic o –co ec o p ocess wha imp o es he e i-
ciency. To show he imp o emen s ha ou p oposal
allows ob aining, we will ca y ou a se o es usually
used o quan i y he s abili y and accu acy.
2 SPH Me hod
The SPH me hod, ini ially de eloped by Gingold and
Monaghan [5] and Lucy [11], was o mula ed om
he in eg al o mula ion o he del a Di ac unc ion
[13].
( ) = Z ( 0)δ( − 0)d 0,(1)
whe e depic s any scala unc ion and δis he del a
Di ac unc ion.
Acco ding o he in e pola ion heo y, he unc-
ion δcan be eplaced by a unc ion wi h spa ial ex-
ension. Thus, he in eg al exp ession (1) can be e-
o mula ed by he below equa ion.
( ) = ZΩ
( 0)W( − 0, h)d 0+O(h2),(2)
whe e Wis named he ke nel unc ion, Ω e e s o he
olume o he de ini ion domain o Wand his he
adius o he suppo ed domain.
Rega ding he ke nel unc ion, his one mus sa -
is y, a leas , he ollowing cons ain s:
1. I mus be no malized and end o δwhen h ends
o ze o, i.e:
lim
h→0W( − 0, h) = δ( − 0); ZV
WdV 0= 1.
(3)
2. I mus be posi i e and dec ease con inuously
wi h ( − 0).
3. I mus be symme ic wi h espec o ( − 0).
Quan i a i ely, he modelling o a con inuous
medium using a pa icles sys em en ails ans o ming
he in eg al, equa ion (2), by a sum whe e he mass o
each pa icle depic s a olume elemen whose mass is
ρdV , i.e:
h ( )i=ZΩ
( 0)
ρ( 0)W( − 0, h)ρ( 0)d 0
≈X
j∈N(i)
mj
j
ρj
W( j− i, h),(4)
whe e hi e e s o he app oxima ed alue o unc ion
since he second o de e m has no been conside ed,
N(i)depic s he neighbou pa icles jo each pa icle
iand ρjis he mass densi y alloca ed o he pa icle j.
F om he symme y p ope y (i em 3), he g adi-
en o any scala magni ude ( )can be simpli ied as:
h∇ ( )i=∂
∂ Z ( 0)
ρ( 0)W( − 0, h)ρ( 0)d 0
≈X
j∈N(i)
mj
j
ρj∇W( j− i, h).(5)
Gene alizing he equa ion (5), i is possible o o -
mula e a alid equa ion o any di e en ial o de . This
gene alized equa ion is:
h∇l ( i)i=X
j∈N(i)
mj
j
ρj∇lW( j− i, h),(6)
whe e l e e s o he di e en ial o de equa ion.
I is no ed ha he equa ion (6) does no ega d he
second o de e m [21]. Also, i iola es he conse a-
ion o angula momen um [13]. As consequence, al-
hough he equa ion (6) sa is ies he SPH p emises, he
accu acy o he ob ained esul s can be comp omised.
To p e en his issue, his equa ion mus be modi ied
o be applied o luid low simula ion, acco ding is de-
sc ibed by [13]. Ne e heless, his symme ical o -
mula ion can inc ease he ins abili y [2, 7, 21].
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
51
Volume 13, 2018
3 Dynamic Equa ions
The dynamic o an incomp essible and in iscid luid,
whose low e ol es adiaba ically wi h high Reynolds
numbe is desc ibed by he below PDEs sys em:
∇· = 0,(7)
D
D =−1
ρ(∇p) + F, (8)
whe e D/D =∂/∂ + ·∇ e e s o he con ec i e
de i a i e, ρis he mass densi y, is he eloci y, pis
he p essu e and Fdepic s he body o ce ec o pe
uni olume.
To close his sys em o equa ions, i is necessa y
o ega d a ela ionship be ween he p essu e and he
mass densi y. This ela ionship is se by he equa ion
o s a e (EOS). Al hough he e a e se e al o mula-
ions o he EOS, in SPH is usually used Tai ’s equa-
ion exp essed as:
p= (γ−1)ρ, (9)
whe e γis he adiaba ic index.
The modelling o he equa ions (7) and (8) by
SPH shows p oblems o s abili y and accu acy. To
sol e hese p oblems, se e al p oposals ha e been de-
eloped [15, 21, 23, 30]. The mos ou s anding a e he
esea ch de eloped by P ice [21], appoin ed as MPM,
and he echnique desc ibed by Inu suka [7] known as
GPSH.
Rega ding he a ia ional app oach, P ice [21]
models he equa ion (8) as:
d i
d =−X
j∈N(i)
mj(χi∇Wij(hi) + χj∇Wij(hj)) +
X
j∈N(i)
αmi
ρij
s ij ·ˆ
Rij∇˜
Wij +F(10)
whe e he subsc ip iand j e e s o pa icle iand j
espec i ely, χ= (Ω−1p/ρ2),Ωis a co ec ion e m
ela ed o he smoo hing leng h, ˆ
Rij is he uni ec o
de ined by he line joining he pa icles iand j,˜
Wij
is he a e aged ke nel and s=ci+cj− ij ·ˆ
Rij.
In P ice [21] is shown a de ailed explana ion o he
de i a ion o his equa ion.
On he o he hand, he echnique GSPH is seeking
o es o e consis ency o he disc e e densi y es ima e
by he ke nel unc ion con olu ions using a Riemann
sol e . Acco ding o he GSPH echnique, he equa-
ion (8) is o mula ed as:
d i
d =−X
j∈N(i)
mjp∗
ij hV2
ij(hi)∇Wij(hi√2)+
V2
ij(hj)∇Wij(hj√2)i
≈ − X
j∈N(i)
mjp∗
ij 1
ρ2
i∇˜
Wij
1
ρ2
j∇˜
Wij!+F(11)
whe e p∗is he in e media e p essu e a ising om he
solu ion o a Riemman sol e pa icula ized o he pa -
icles pai i, j and Vij a e he in e media e eloci ies
ela ed o he solu ion o Riemann sol e pa icula -
ized o wo pa icles [24]. The pape s [7, 8] show a
de ailed explana ion o he de i a ion o he equa ion
(11).
4 P oposed Model
In his sec ion, i will be desc ibed he ea u es o
ou p oposal o imp o e he s abili y and accu acy o
SPH simula ions. Fi s ly, he undamen al hypo he-
sis o he model will be se . F om his hypo hesis
will be o mula ed ou p oposal o ob ain he p es-
su e o ce. Nex , i will be desc ibed he analy ic ea-
u es ha should be sa is ied by he p essu e o ce o
gua an ee s abili y. These ea u es a e se om he
mos ou s anding s udies in s able simula ions scope
by SPH. F om hese cons ain s and he undamen al
hypo hesis, i will be o mula ed he p essu e o ce.
I s o mula ion will allow en o cing he incomp ess-
ibili y. Once he incomp essibili y is en o ced, i will
be explained he p ocess o calcula e he densi y us-
ing bo h he SPH o mula ion and he p essu e o ce.
This joined o mula ion allows ul illing he incom-
p essibili y cons ain .
Analy ically, when ex e nal o ces ac ing on a
luid, he dis ances be ween adjacen egions change.
Then, in o he luid appea s a o ce ha ends o e-
s o e he balance condi ions [34]. This in e nal o ce
is quan i ied by he p essu e g adien . Ex apola ing
his c i e ion o he modelled luid by a pa icle sys-
em, we can o mula e ou undamen al hypo hesis,
namely: he p essu e o ce is due o change o he dis-
ance be ween pa icles.
Taking in o accoun his p emise, we a oid he
p oblems o s abili y ela ed o modelling o he p es-
su e g adien by SPH, especially when he numbe
o neighbo ing pa icles is low [3, 22]. Fu he mo e,
we can en o ce he incomp essibili y cons ain on he
p essu e o ce and he densi y. On he o he hand, as
he o ce only depends on dis ance, we can o mula e
i om a conse a i e po en ial [34]. Ne e heless,
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
52
Volume 13, 2018
any conse a i e po en ial is no alid o ob ain s able
luid simula ions. To se a sui able o mula ion mus
be aking in o accoun he cons ain s ha gua an ee
he s abili y. These cons ain s i will be desc ibed be-
low.
4.1 P essu e Fo ce
Se e al esea ch ha e been ca ied ou o se he ela-
ionship be ween s abili y–accu acy and he p essu e
g adien . The mos highligh ing s udies asse ha he
p essu e o ce:
1. I should ensu e s abili y o a wide ange o
he numbe o neighbou pa icles. This ea u e
can be deduced om he s udies de eloped by
[3, 12, 21, 28]. These esea ch asse ha i he
numbe o pa icles is oo low, he simula ions
show an e a ic sca e ing o pa icles, mainly in
he bounda y o luid. On he o he hand, i he
numbe o pa icles is high, a clus e ing o pa i-
cles can appea .
2. I should be conca e in espec o he abscissa
axis. This cons ain is deduced om he e-
sea ch de eloped by [1, 3, 24], whe e is de-
sc ibed ha a ke nel unc ion whose g adien ,
ha appea in p essu e o ce, shows an in lec ion
poin and conca i y, imp o e he s abili y.
3. I should show an asymme ic p o ile wi h ega d
o he minimum alue o he o ce. Thus, wo
equi emen s a e sa is ied: i s ly, he epulsi e
o ce is inc eased when he dis ance is dec eased,
as i can be deduced om [9] and las ly, he in-
e ac ion can be weighed by he dis ance, as is
sugges ed by [3, 21]. A ha monic o ce can only
sa is y he i s condi ion bu i would no ensu e
he las one [21].
4. I should allow us o en o ce he incomp essibil-
i y [12, 28], namely, i s o mula ion should limi
he minimum dis ance among pa icles. Thus,
he o mula ed o ce should show an asymp o ic
beha io wi h dis ance. I s aim would be o a oid
ha he dis ance be ween pa icles will be incom-
pa ible wi h incomp essibili y cons ain .
F om hese ea u es, we o mula e he conse a-
i e po en ial, om which, we will ob ain he p essu e
o ce. The p oposed po en ial is:
Vi( i) = 21
24πhiX
j∈N(i) 2
(Rij −0.9ε)2−
7π2
2ζ(Rij −(0.85ε)) + (3π2σ)!,(12)
whe e Rij = j− i,εis named as asymp o ic ac o ,
ha allows us o con ol he incomp essibili y a e, ζ
is appoin ed as dep h ac o , by his magni ude, we
can a oid he pa icles clus e ing and σis named as
ail coe icien , by means o his coe icien , i is pos-
sible o con ol he po en ial a he bounda y o he
suppo ed domain. Empi ically, we ha e e alua ed he
alue ange, o hese pa ame e s, ha o e he bes e-
sul s, hese alues a e: 0.3h≤ε≤h,0.5≤ζ≤2.0
and 0.8≤σ≤2.25. The Figu e 1 shows he p o ile
o he o mula ed po en ial.
Figu e 1: G aph o he p oposed po en ial. The pa-
ame e alues a e h= 3,ε= 0.6,ζ= 1.25 and
σ= 1.0.
Analy ically he po en ial, equa ion (12), only de-
pends on posi ion. As consequence, he o ce ha is
de i ed can be w i en as ~
Fpi=−∇V( i), i.e:
~
Fpi=−5.25
6πhiX
j∈N(i) 3.5π2
ζ(Rij −(0.85ε))2−
4
(Rij −0.9ε)3!(h−Rij)~
Rij
Rij
,(13)
F om he Figu e 1, we can desc ibe he ea u es o
ou p oposed po en ial. Fi s ly, i shows a conca i y in
espec o he abscissa, as i was es ablished in i em 2.
Secondly, i shows an asymme ic p o ile wi h espec
o minimum o o ce, acco ding o he i em 3. Las ly,
as a consequence o asymp o ic limi , o low alues
o , we can con ol he incomp essibili y, as i was
equi ed in i em 4. I is no ed, as consequence o he
p essu e o ce is de i ed om his po en ial, his ea-
u es will be inhe i ed. On he o he hand, analyzing
he equa ion (13) we can say ha he o ce is sa u a ed
wi h he dis ance [34]. Fu he mo e, his sa u a ion is
weighed by he dis ance. Thus, he p essu e o ce sa -
is ies he cons ain sugges ed by [3, 13, 17, 32]. Also,
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
53
Volume 13, 2018
he o ce o mula ed in equa ion (13) sa is ies he con-
se a ion o momen um.
Then, once he p essu e o ce has been ob ained,
we a e going o desc ibe ou p oposal o con ol he
incomp essibili y, aking ad an age o he asymp o ic
limi shown by he o mula ed p essu e o ce.
Due o he ac ha he p essu e o ce is ob ained
om a po en ial, ou echnique will he eina e be e-
e ed o as P–SPH.
4.2 The En o cemen o he Incomp essibil-
i y
Quali a i ely, he incomp essibili y is ela ed o he
maximum o ce ha can be bo ne by he luid, wi hou
changing i s local olume [34]. Acco ding o his ac ,
he incomp essibili y in a luid modelled by a pa icle
sys em can be explained as he minimum dis ance ha
he pa icles can app oach. In his con ex , he o ce
calcula ed by he equa ion (13), allows us o model
he incomp essibili y by he asymp o ic alue when R
is low. Ne e heless, i he incomp essibili y coe i-
cien is lowe han 2%, we mus in oduce a con ol
p ocess o educe he nume ical e o s accumula ed in
each ime s ep.
Ou p oposal o con ol he incomp essibili y
when i s coe icien lowe han 2% can be di ided in o
wo s ages. In he i s s age, he p essu e o ce is
calcula ed by he equa ion (13). In he las s age, we
check he pa icles whose p essu e o ce iola es he
incomp essibili y cons ain . I he p essu e o ce o
any pa icle iola es he incomp essibili y cons ain ,
hen he o ce excess is quan i ied and bo h he dy-
namic and he posi ion o hese pa icles a e modi ied
o sa is y he incomp essibili y cons ain .
To implemen his p ocess, i s ly we calcula e he
maximum comp ession o ce, namely, he alue om
which he incomp essibili y cons ain is iola ed. To
his end, we pa icula ize he equa ion (13) wi h a dis-
ance, whose alue is in line wi h he es densi y luid.
Acco ding o his desc ip ion, he ob ained equa ion
is:
~
FC
j|i=−5.25
6πh 3.5π2
ζ(ξ−(0.85ε))2−4
(ξ−0.9ε)3!~
R
R,
(14)
whe e ~
FC
j|iis he maximum o comp ession o ce ha
he pa icle jcan acqui e in espec o pa icle iand ξ
is he minimum dis ance ha he pa icle jcan ap-
p oxima e o pa icle iwi hou iola es he incom-
p essibili y. The equa ion ha sa is ies his c i e ia is:
ξ=3
s3m
4π(1 −η)ρ0
,(15)
whe e ηis he comp essibili y coe icien and ρ0is he
es densi y.
I |~
FC
j|i|>|~
Fpj|i|, hen, we ha e o con ol he
incomp essibili y. Thus, we change he o ce o he
pa icle jand i s posi ion, o he maximum ole a ed
alues, i.e.:
~
Fpj|i=~
FC
pj|iy ~ 0
j=ξ~
R
R,(16)
whe e ~
Fpj|iis he p essu e o ce o pa icle jin espec
o pa icle iand ~ 0
jis i s new posi ion.This change only
is alid i he di e ence be ween |FC
pj|i|and |FC
pj|i|is
lowe han 5%. Al hough his limi a ion seems es ic-
i e, as he p essu e o ce con ol he incomp essibil-
i y, i allows us ha he di e ence be ween wo o ces
is lowe han his alue.
5 Resul s
In his sec ion, we will ca y ou a se o es s o show
he imp o emen s ha o e ou p oposal. To ha
end, we will simula e a luid whose dynamic is de-
sc ibed by he equa ions (7) and (8). To sol e his sys-
em o equa ions, i will be ega ded he echniques
MPM [21], GSPH [7] and P–SPH. The aim is o com-
pa e he ob ained esul s om each echnique o show
he imp o emen ha can be ob ained by he p o-
posed echnique. I will be implemen ed h ee es s
commonly used o quan i y he s abili y and accu acy
[24, 26, 27], hese a e: se ling o a andom pa icles
dis ibu ion,Sod’s shock– ube and blas –wa e. In
each es , we will ca y ou an analysis o he ange
o he numbe o neighbou pa icles whe e each ech-
nique is s able.
We ha e ega ded he mos o he alues o sim-
ula ion pa ame e s sugges ed by [24]. Each es has
di e en alues ha be speci ied in each simula ion.
Ne e heless, in all hem, he alue o he adiaba ic
pa ame e γ, equa ion (9), will be γ= 1,46. Rega d-
ing ou p oposed p essu e o ce, equa ion (13), he se-
lec ed alues o εand ζa e ε= 0.45 and ζ= 1.0.
Likewise, all simula ions will be implemen ed using
a Wendland C4as ke nel unc ion. We ha e selec ed
his unc ion due o i s analy ic ea u es ha a ou he
s abili y and accu acy, acco ding is highligh ed by [3].
5.1 Se ling o a Random Pa icles Dis ibu-
ion
Se e al esea ch show a consensous which ega ds he
se ling o a andom pa icle dis ibu ion as a sui able
es o quan i y he s abili y in SPH [3, 21, 29]. De-
sc ip i ely, in his es , a pa icles sys em e ol e om
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
54
Volume 13, 2018
P-SPH
Figu e 2: E olu ion om a andom pa icles dis ibu ion o a se led pa icles dis ibu ion.
a andom dis ibu ion o a s a e ha shows a le el o
egula i y. The achie ed egula i y allows quan i ying
he s abili y deg ee o any SPH simula ion [21, 24].
To implemen his es , we ha e used 980 pa icles
whose mass is m= 0.00213. The simula ion limi s
a e [50 ×30] and he ime s ep used ∆ = 1.5 10−4s
and h= 0.023. F om his alues, we ha e ob ained
he simula ion shown in Figu e 2.
To analyse he ange o he numbe o neighbou
pa icle ha allow ob aining s able and accu acy sim-
ula ion, we ha e changed om he lowes numbe o
highes . In his es , he used alues a e:
MPM GSPH P–SPH
Ni120–375 95–375 15–375
Ni180 165 25
Table 1: Range o he numbe o neighbou pa icles.
Ni e e s o ange ha allows s abili y and Niis he
alue ega ded o ob ain he simula ions shown in Fig-
u e 2.
F om he Figu e 2, we can say ha he bes esul
is ob ained om P–SPH. This one allows ob aining
he esul s o he highes egula i y. To quan i y his
la ice deg ee, we ha e o mula ed an o de pa ame e
ς ha sa is ies:
ς=Ps
k=1(|dl−dk
ij|)
s,(17)
whe e s= (m−1)(n−1), and mand na e he pa -
icles numbe in each space di ec ion, dlis he la ice
leng h, dk
ij is he dis ance be ween each pai o neigh-
bou pa icles. F om he o mula ion o he pa ame e
ς, i can be said ha he lowe alue, he mo e la ice
deg ee is achie ed and, consequen ly, he mo e s abil-
i y is ob ained.
The ob ained alues o his pa ame e a e shown
in Table 2.
MPM GSPH P–SPH
ς21.7 16.4 4.7
Table 2: Ob ained alues o he o de pa ame e om
Figu e 2.
5.2 Sod’s Tube–Shock Tes
Sod’s shock– ube es [27] allow quan i ying he accu-
acy and s abili y o SPH [10, 14, 26]. Desc ip i ely,
he shock– ube es is buil by a ube whe e a non–
po ous memb ane is in oduced. This memb ane de-
ines wo sec ions ha a e illed wi h gas. In a sec ion,
ha is named as high sec ion, he gas is a ound 10
imes dense and mo e con ined han in he o he sec-
ion which is named low sec ion. F om hese ini ial
cons ain s, he memb ane is ins an ly aken away, a
= 0. A his momen , due o di e ences o p essu e
and densi y be ween high–low sec ions, a shock wa e
is gene a ed. This one e ol es in o he low sec ion.
A he same ime, a a e ac ion wa e appea s in o he
high sec ion. Be ween each wa e, a con ac discon i-
nui y is a isen. Quan i e i ely, om he ob ained ac-
cu acy in he simula ion o hese h ee egions, i is
possible o quan i y he sui abili y o he used ech-
nique.
To ca y ou his es we ha e used a domain
[-0.5,0.5]. The low sec ion is limi ed be ween
[-0.5,0.0)and he high sec ion is limi ed by
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
55
Volume 13, 2018
(a) (b) (c)
(d) (e) ( )
Figu e 3: Compa a i e g aphs o he exac solu ion and nume ical ob ained alues om MPM, GSPH and P-SPH
o Sod’s shock– ube es .
(0.0,0.5]. The emaining used pa eme e s a e he
sugges ed by [24].
F om hese ini ial alues, we ha e ob ained he e-
sul s shown in Figu e 3 a ime =0.15. I is ega ded
he densi y and he eloci y o e alua e he ob ain ac-
cu acy, acco ding o he sugges ion gi en by [24, 26].
The Figu e 3 shows a compa ison be ween he
exac solu ion, solid line, and he nume ical solu ion
shows as do s.
In each g aphs, i is possible o iden i y h ee e-
gion: he con ac discon inui y, he a e ac ion wa e
and he shock wa e. The con ac discon inui y appea s
a ound xcd = 1.004, he a e ac ion wa e appea om
x<xcd and he shock wa e om x>xcd. Al hough
he posi ion o he con ac discon inui y is he same in
all g aphs, bo h he spa ial wid h o he con ac dis-
con inui y and i s wa ing, shows signi ican disc ep-
ancies.
To quan i y he ob ained accu acy, we ha e se-
lec ed he nume ic alues o each pa o he low,
i.e., a e ac ion wa e, con ac discon inui y and shock
wa e, speci ically he wid h ha a e appoin ed as: W
he wid h o he a e ac ion ail, Wcd he wid h o
con ac discon inui y and Wsw he wid h o he shock
wa e. The ob ained alues a e shown in Tables 3 and
4.
5.3 Blas –Wa e Tes
The blas –wa e es is a mo e ex eme e sion o
Sod’s shock– ube es and is o mula ed o a e y
MPM GSPH P–SPH
20.9651 0.9708 0.9891
W 0.105 0.073 0.0205
Wcd 0.078 0.061 0.026
Wsw 0.128 0.104 0.051
Table 3: The mos ele an alues associa ed wi h he
densi y g aphs o Sod’s shock– ube es , Figs. 3a, 3b
and 3c.
MPM GSPH P–SPH
20.9651 0.9688 0.9891
W 0.126 0.099 0.039
Wcd 0.101 0.075 0.043
Wsw 0.089 0.065 0.037
Table 4: The mos ele an alues associa ed wi h he
eloci y g aphs o Sod’s shock– ube es , Figs. 3d, 3e
and 3 .
high Mach numbe , a ound M= 200. I was de-
eloped by Woodwa d and Colella [31] o highligh
he beha iou o he con ac discon inui y, mainly o
he eloci y, a high densi y and p essu e alues. In
blas –wa e es , he low e olu ion shows he same
ea u es ha Sod’s shock– ube es , namely, a con ac
discon inui y will appea inse ed be ween he shock
wa e and he a e ac ion wa e.
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
56
Volume 13, 2018
(a) (b) (c)
(d) (e) ( )
Figu e 4: Compa a i e g aphs o he exac solu ion and nume ical ob ained alues om MPM, GSPH and P-SPH
o blas –wa e es .
As a esul o he simila i ies wi h Sod’s shock–
ube es , we ha e only changed he ini ial alues
o some simula ion pa ame e s. The new alues
o mass densi y, p essu e and he mal ene gy a e:
(ρh, ph, uh) = (1.0,1 103,0.0) and (ρl, pl, ul) =
(1.0,0.01,0.0), acco ding o he sugges ions gi en by
[26, 31]. The new ime s ep is ∆ = 8 10−5. Us-
ing hese alues, we ha e ob ained he esul s ha a e
shown in Figu e 4 a ime = 7.5 10−4, whe e he
exac solu ions and nume ic solu ions a e depic ed by
solid line and do s, espec i ely.
Simila ly o Sod’s shock– ube es , he h ee e-
gions ha a ise when he memb ane is aken away
can be loca ed in he g aphs o Figu e 4: he con ac
discon inui y posi ion xcd = 0.1081, he a e ac ion
wa e x<xcd and he shock wa e x>xcd. Likewise,
in each g aph o 4 he wid h o con ac discon inu-
i y changes o echnique ob aining he lowes wid h
when ou p oposal is implemen ed.
Analogously as he p e ious sec ion 5.2 o shown
he ob ained accu acy, we ha e selec ed he mos ele-
an alues associa ed o he g aphs o 4. These alues
a e shown in Tables 5 and 6.
6 Conclusions
In his pape , we ha e e iewed some o he mos
ele an ea u es ha in luence on he accu acy and
s abili y o low simula ion by SPH. F om hese ea-
MPM GSPH P–SPH
20.9618 0.9791 0.9903
W 0.094 0.081 0.046
Wcd 0.035 0.018 0.009
Wsw 0.094 0.079 0.038
Table 5: The mos ele an alues associa ed wi h he
densi y g aphs o blas –wa e shock es , Figs. 4a, 4b
and 4c.
MPM GSPH P–SPH
20.9621 0.9788 0.9901
W 0.0945 0.0813 0.0585
Wcd 0.0826 0.0582 0.0348
Wsw 0.087 0.026 0.009
Table 6: The mos ele an alues associa ed wi h he
eloci y g aphs o blas –wa e es , Figs. 4d, 4e and
4 .
u es, we ha e p oposed a hyb id po en ial–SPH based
me hod. To show he imp o emen s ob ained by he
p oposed echnique, we ha e implemen ed a se o
es s using an incomp essible and in iscid luid. In
each implemen ed es we ha e compa ed ou p o-
posal wi h he MPM and GSPH echniques which a e
o mula ed o a oid he ins abili ies ela ed o he p es-
su e g adien . Pa icula ly, analysing he es s’ esul s,
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
57
Volume 13, 2018
we can conclude:
1. The P-SPH echnique allows us o ob ain a
high s abili y. This conclusion can be deduced
om he ob ained esul s in se ling es , Figu e
2, whe e he ob ained esul s by P-SPH show
g ea e la ice deg ee han he ob ained ones by
he o he wo echniques. Quan i a i ely, his im-
p o emen is shown by he esul s o he pa am-
e e ς, Table 2, whe e he ob ained alues by P-
SPH a e he lowes ones.
2. The P-SPH allow us o ob ain highly accu a e e-
sul s in simula ions whe e he shock wa e and
con ac discon inui y appea . This a i ma ion is
based on he ob ained esul s in Sod’s shock– ube
es . The ob ained g aphs o he densi y and e-
loci y, Figu e 3, a e mo e accu a ely simula ed
by P-SPH han by MPM o GSPH. This accu-
acy is mo e appa en in he a e ac ion ail, he
con ac discon inui y and he shock wa e. To e-
in o ce his ac , we ha e selec ed he alues as-
socia ed o hese egions ha a e shown in Tables
3 and 4. All ob ained alues e idence ha he
bes accu acy is o e ed by P-SPH.
3. The imp o emen s ob ained in Sod’s shock– ube
es is also accomplished in a mo e ex eme e -
sion o his es . This conclusion can be ob ained
om he esul s shown in he blas –wa e es .
Quali a i ely, he ob ained esul s a e in ag ee-
men wi h he exac solu ions. Ne e heless, he
mos accu a e esul s a e ob ained by he P-SPH,
as can be seen in he g aphs o Figu e 4. Quan i-
a i ely, he imp o emen is shown in he alues
o Tables 5 and 6.
4. The P-SPH allow us o ob ain s able simula ions
o a wide ange o he numbe o neighbou pa -
icles. Also, he s abili y is achie ed wi h he
lowes numbe , as consequence. The P–SPH no
only a ou he s abili y bu also he e iciency,
om he compu a ional cos poin o iew. This
conclusion can be deduced om Table 1.
F om hese esul s, we can conclude ha he P-
SPH imp o es bo h he s abili y and he accu acy o
luid low simula ions by SPH
Acknowledgemen s: This esea ch has been sup-
po ed by he Pololas p ojec (TIN2016-76953-C3-2-
R) o he Spanish Minis y o Economy and Compe -
i i eness.
Re e ences:
[1] R. M. Cabez´
on, D. Ga c´
ıa-Senz, and A. Re-
lano. A one–pa ame e amily o in e pola -
ing ke nels o smoo hed pa icle hyd odynam-
ics s udies. Jou nal o Compu a ional Physics,
227(19):8523–8540, Oc 2000.
[2] S. J. Cummins and R. M. An imp o ed sph
me hod o modeling liquid sloshing dynamics.
Jou nal o Compu a ional Physics, 152:584–
607, Jul 1999.
[3] W. Dehnen and H. Aly. Imp o ing con e gence
in smoo hed pa icle hyd odynamics simula ions
wi hou pai ing ins abili y. Mon hly No ices o
he Royal As onomical Socie y, 425(2):1068–
1082, Ap 2012.
[4] C. T. Dyka, P. W. Randles, and R. P. Ingel. S ess
poin s o ension ins abili y in sph. In e na-
ional Jou nal o Nume ical Me hods in Engi-
nee ing, 40(13):2325–2341, Jul 1997.
[5] R. Gingold and J. Monaghan. Smoo hed pa -
icle hid odynamics: Theo y and applica ion o
non–sphe ical s a s. Royal As onomical Soci-
e y, 181:375–389, 1977.
[6] M. Ihmsen, J. O hmann, B. Solen hale ,
A. Kolb, and M. Teschne . Sph luids in com-
pu e g aphics. In Eu og aphics 2014 - S a e o
he A Repo s. The Eu og aphics Associa ion,
2014.
[7] S. Inu suka. Re o mula ion o smoo hed pa i-
cle hyd odynamics wi h iemann sol e. Jou nal
o Compu a ional Physics, 179(1):238–267, Jun
2012.
[8] K. Iwasaki and S. Inu suka. Smoo hed pa icle
magne ohyd odynamics wi h a iemann sol e
and he me hod o cha ac e is ics. Mon hly
No ices o he Royal As onomical Socie y,
418(3):1668–1688, Dec 2011.
[9] G. Liu and M. Liu. Smoo hed Pa icle Hyd o-
dynamics:A Mesh ee Pa icle Me hod. Wo ld
Scien i ic, 2003.
[10] M. Liu, G. Liu, and K. Lam. Cons uc ing smoo-
hing unc ions in smoo hed pa icle hyd ody-
namics wi h applica ions. Jou nal o Compu-
a ional and Applied Ma hema ics, 155(2):263–
284, Jun 2003.
[11] L. B. Lucy. A nume ical app oach o he es ing
o he ission hypo hesis. As onomical Jou nal,
82(12):1013–1024, Dec 1977.
WSEAS TRANSACTIONS on FLUID MECHANICS
Juan J. Pe ea, Juan M. Co de o
E-ISSN: 2224-347X
58
Volume 13, 2018