A Finite-Volume discretization of viscoelastic Saint-Venant equations for FENE-P fluids

0.23 MB·

by  Sébastien Boyaval

A Finite-Volume discretization of viscoelastic
Saint-Venant equations for FENE-P fluids
Se´bastienBoyaval
7
1
0
2

n
a
J
Abstract Saint-Venant equations can be generalized to account for a viscoelastic
7
rheologyinshallowflows.AFinite-Volumediscretizationforthe1DSaint-Venant
1
system generalized to Upper-Convected Maxwell (UCM) fluids was proposed in

] [Bouchut & Boyaval, 2013], which preserved a physically-natural stability prop-
A
erty(i.e.free-energydissipation)ofthefullsystem.Itinvokedarelaxationscheme
N ofSuliciutypeforthenumericalcomputationofapproximatesolutiontoRiemann
. problems.Here,theapproachisextendedtothe1DSaint-Venantsystemgeneralized
h
t tothefinitely-extensiblenonlinearelasticfluidsofPeterlin(FENE-P).Wearecur-
a
rentlynotabletoensureallstabilityconditionsapriori,butnumericalsimulations
m
wentsmoothlyinapracticallyusefulrangeofparameters.
[

1 Key words: Saint-Venant equations, FENE-P viscoelastic fluids, Finite-Volume,
v simpleRiemannsolver,Suliciurelaxationscheme
9
3 MSC(2010):65M08,65N08,35Q30
6
4
0
. 1 Introduction
1
0
7
Saint-Venantequationsstandardlymodelshallowfree-surfacegravityflowsandcan
1
be generalized to account for the viscoelastic rheology of non-Newtonian fluids
:
v [6], Upper-Convected Maxwell (UCM) fluids in particular [5]. Here, we consider
i
X a generalized Saint-Venant (gSV) system for finitely-extensible nonlinear elastic
fluidswithPeterlinclosure(FENE-Pfluids)inCartesiancoordinates
r
a
Se´bastienBoyaval
Laboratoired’hydrauliqueSaint-Venant(EcoledesPontsParisTech–EDFR&D–CEREMA)
Universite´Paris-Est,EDF’lab6quaiWatier78401ChatouCedexFrance,&INRIAParisMATH-
ERIALSe-mail:[email protected]
1
2 Se´bastienBoyaval
∂h+∂ (hu)=0 (1)
t x
∂(hu)+∂ (cid:0)hu2+gh2/2+hN(cid:1)=0 (2)
t x
λ(∂σ +u∂ σ +2(ζ−1)σ ∂ u)=1−σ /(1−(σ +σ )/(cid:96)) (3)
t xx x xx xx x xx zz xx
λ(∂σ +u∂ σ +2(1−ζ)σ ∂ u)=1−σ /(1−(σ +σ )/(cid:96)) (4)
t zz x zz zz x zz zz xx
for1De -translationinvariantflowalonge underauniformgravityfield−ge with
y x z
• meanflowdepthh(t,x)>0(incaseofanon-rugousflatbottom),
• meanflowvelocityu(t,x)(foruniformcrosssections),and
• a normal-stress difference N =G(σ −σ )/(1−(σ +σ )/(cid:96)) given by con-
zz xx zz xx
formation variables σ ,σ >0 constrained by 0<σ +σ <(cid:96), a relaxation
zz xx zz xx
timeλ ≥0andanelasticitymodulusG>0.
Note that (1-2-3-4) formally reduces to the standard viscous Saint-Venant system
withviscosityν ≡2λG≥0when(cid:96)→∞,λ →0andGλ <∞.Moreoverwehave
usedthequitegeneralGordon-Schowalterderivativeswithslipparameterζ ∈[0,1]
2
constrainedbythehyperbolicityofthesystem(1-2-3-4).(Thisfollowsafteraneasy
computationsimilarto[8].)
In this work, we discuss a Finite-Volume method to solve (numerically) the
Cauchyproblemforthenonlinearhyperbolic1Dsystem(1-2-3-4).Standardly,we
need to consider weak solutions (in fact, to (6-7-8-9), see below) plus admissibil-
ityconstraintsthatarephysically-meaningfuldissipationrulesformalizingthether-
modynamics second principle close to an equilibrium [9]. Here, we consider the
inequalityassociatedwiththecompanionconservationlawforthefree-energy
(cid:18)u2 gh G (cid:19)
F =h + − ((cid:96)log(((cid:96)−(σ +σ ))/((cid:96)−2))+log(σ σ ))
xx zz xx zz
2 2 2(1−ζ)
thatis,ondenotingtheimpulsebyP=gh2/2+hN,
Gh (cid:32) (cid:18) σ (cid:19)2 (cid:18) σ (cid:19)2(cid:33)
− σ−1 1− xx +σ−1 1− zz
2(1−ζ)λ xx 1−(σ +σ )/(cid:96) zz 1−(σ +σ )/(cid:96)
zz xx zz xx
=:D≥∂F+∂ (u(F+P)) (5)
t x
wheretheleft-hand-sideisobviouslynon-positiveontheadmissibilitydomain
U(cid:96):={0<h,0<σ ,0<σ ,σ +σ <(cid:96)}.
xx zz xx zz
Notethatwedonotconsiderthevacuumstateh=0asadmissiblehere,see[8].
Saint-Venant/FENE-Pshallowflows 3
2 Finite-VolumediscretizationofFENE-P/Saint-Venant
Piecewise-constantapproximatesolutionstotheCauchyproblemon(t,x)∈[0,T)×
R for the gSV system can be defined by a Finite-Volume (FV) method. With a
viewtopreservingU(cid:96) andthedissipation(5)afterdiscretizationbyaFVmethod,
wechooseq=(h,hu,hσ ,hσ )asdiscretizationvariable.Indeed,thefree-energy
xx zz
functionalFisconvexontheconvexdomainU(cid:96)(cid:51)q(thisfollowsafteraneasycom-
putationfrom[4,Lemma1.3])whileitisnotconvexinthevariable(h,hu,hΠ,hΣ)
whateversmoothinvertiblefunctionsϖ,ς areusedforthereformulationofgSV
∂h+∂ (hu)=0 (6)
t x
(cid:18) gh2 (cid:19)
∂(hu)+∂ hu2+ +hN =0 (7)
t x
2
(cid:32) (cid:33)
h3−2ζϖ(cid:48)(σ h2(1−ζ)) σ
xx xx
∂(hΠ)+∂ (huΠ)= 1− (8)
t x λ 1−σzz+σxx
(cid:96)
(cid:32) (cid:33)
h2ζ−1ς(cid:48)(σ h2(ζ−1)) σ
zz zz
∂(hΣ)+∂ (huΣ)= 1− (9)
t x λ 1−σzz+σxx
(cid:96)
with Π =ϖ(σ h2(1−ζ)), Σ =ς(σ h2(ζ−1)) (computations are similar to [5, Ap-
xx zz
pendix]).Inthesequel,wethereforediscretizeaquasilinearsystemwithsource
∂q+A(q)∂ q=S(q), (10)
t x
whichwerecallisnotambiguoushere(forthosediscontinuoussolutionsbuiltusing
aRiemannsolver,atleast)thankstothedissipationrule(5),see[11,2,8].
2.1 Splitting-in-time
In cell (x ,x ), i∈Z, with volume ∆x =x −x >0 and center
i−1/2 i+1/2 i i+1/2 i−1/2
x =(x +x )/2,weapproximateqsolutionto(10)onR ×R(cid:51)(t,x)by
i i−1/2 i+1/2 ≥0
qn+1≈ 1 (cid:90) xi+1/2q(t,x)dx, i∈Z,t∈(tn,tn+1]
i ∆xi xi−1/2
onatimegrid0=t0<t1<...<tn<tn+1<...<tN =T where∆tn=|tn+1−tn|
willbechosensmallenoughcomparedwith∆x=supi∈Z∆xi<∞toensurestability.
More precisly, having in mind the numerical approximation of a (well-posed)
Cauchyproblemfor(10)onR ×Rwithinitialconditionq(t→0+)=q0∈L∞(R),
≥0
andthereforestartingfromapproximationsq0≈ 1 (cid:82)xi+1/2q0(x)dx,i∈Z,weshall
definethecellvaluesqnintwostepsforeachin=∆1x,i..x.i−,1N/2:
i
(i)anapproximatesolutiontothehomogeneousgSVsystem(i.e.withoutthesource
4 Se´bastienBoyaval
termS)on[tn,tn+1)isfirstcomputedbyanexplicitthree-pointscheme
qn+1/2=qn−∆tn(cid:0)F(qn,qn )−F(qn ,qn)(cid:1), (11)
i i ∆x l i i+1 r i−1 i
i
(ii)anapproximatesolutiontothefullgSVsystemon(tn,tn+1]isnextcomputedas
qn+1=qn+1/2+∆tnS(qn+1). (12)
i i i
Then,theschemeisconsistentwithweaksolutionsof(1–2)equiv.(6–7)
qn+1=qn−∆tn(cid:0)F(qn,qn )−F(qn ,qn)(cid:1)+∆tnS(qn+1) (13)
i i ∆x l i i+1 r i−1 i i
i
providedthetwofirstfluxcomponentsfortheconservativepart(h,hu)ofthevari-
able q (actually solutions to conservation laws) are conservative F =F :=F ,
l,h r,h h
F =F :=F andconsistentF (q,q)=hu| ,F (q,q)=(hu2+gh2/2+hN)|
l,hu r,hu hu h q hu q
asusual,andwiththeconservativeinterpretation(8–9)of(3–4)insofaraswenext
defineF andF usingasimpleapproximateRiemannsolver[10]for(6–7–8–9).
l r
Moreover,withaviewtopreservingU(cid:96)andadiscreteversionof(5)
F(qn+1/2)−F(qn)+∆tn(cid:0)G(qn,qn )−G(qn ,qn))(cid:1)≤0 (14)
i i ∆x i i+1 i−1 i
i
for a numerical free-energy flux function consistent with G(q,q)=u(F+P)| in
q
(5), in the sequel, we shall discuss the relaxation technique introduced by Suliciu
as simple Riemann solver in step (i), because it proved satisfying for other close
systems [3, 4, 5] equipped with an “entropy” convex in the discretization variable
likeF here.Intheend,forthefullscheme(13),aconsistentfree-energydissipation
F(qn+1)−F(qn)+∆tn(cid:0)G(qn,qn )−G(qn ,qn))(cid:1)≤∆tnD(qn+1) (15)
i i ∆x i i+1 i−1 i i
i
thenholdsinsofarhn+1/2=hn+1,un+1/2=un+1andtheconvexityofF imply
i i i i
F(qn+1)−F(qn+1/2)≤∆tnD(qn+1)≤0. (16)
i i i
Proof. Onnotinghn+1/2=hn+1,un+1/2=un+1itsufficestoshowthat
i i i i
(cid:16) (cid:17)
λ σn+1−σn /∆tn =1−σn+1/(1−(σn+1+σn+1)/(cid:96))
xx,i xx,i xx,i zz,i xx,i
(cid:16) (cid:17)
λ σn+1−σn /∆tn =1−σn+1/(1−(σn+1+σn+1)/(cid:96))
zz,i zz,i zz,i zz,i xx,i
imply(16).Now,thisisobvious,onnotingtheconvexityofF| in(σ ,σ )and
h,u xx zz
∇ F| ·S=D
(σxx,σzz) h,u
Saint-Venant/FENE-Pshallowflows 5
since∇ F·(h2(ζ−1)S ,h2(1−ζ)S )=Dbydesign.
(σxxh2(1−ζ),σzzh2(ζ−1)) hσxx hσzz
2.2 SuliciurelaxationoftheRiemannproblemwithoutsource
Foralltimerangest∈[tn,tn+1),n=0...N−1,letusnowdefineateachinterface
x ,i∈Z,betweencellsiandi+1thenumericalfluxfunctionsF andF
i+1 l r
2
(cid:16) (cid:17)
F(q,q )=F (q)−(cid:82)0 R(ξ,q,q )−q dξ,
l l r 0 l −∞ l r l
(cid:16) (cid:17) (17)
F(q,q )=F (q )+(cid:82)∞ R(ξ,q,q )−q dξ.
r l r 0 r 0 l r r
invokinganapproximatesolutionR(cid:0)(x−x )/(t−tn),qn,qn (cid:1)totheRiemann
i+1/2 i i+1
problemfor(10)withinitialconditionqn1 +1 qn att=tn,andanyF .
i x<0 x>0 i+1 0
Inthiswork,weproposeasapproximatesolutionthatgivenbySuliciurelaxation
R(ξ,q,q )=LR(ξ,Q,Q ) (18)
l r l r
i.e. the projection (operator L) onto q ≡ (h,hu,hσ ,hσ ) of the exact solution
xx zz
R(ξ,Q,Q )oftheRiemannproblemforthesystemwithrelaxedpressure
l r
 ∂h+∂ (hu)=0
 ∂t(σxxh2(1−ζ∂)t)(+huu)∂+x(∂σtxx(xhhu22(x1+−ζπ)))==00
∂(σ h2(ζ−1))+u∂ (σ h2(ζ−1))=0 (19)
t zz x zz
∂t(cid:0)h(u2/2+eˆ)(cid:1)+∂∂tx((cid:0)hhπu)(+u2∂/x2(h+ueπˆ)++uucπ2(cid:1))==00
∂c+u∂ c=0
t x
andinitialconditiongivenby(o=l,r)
(cid:16) (cid:17)
Q = h ,(hu) ,h1−2ζ(hσ ) ,h2ζ−3(hσ ) ,h P(q ),(hu)2/2h +e(q ),c
o o o o xx o o zz o o o o o o o
(20)
wherec (q,q )arechosensoastoensurestability,thatisthedissipationrule(14)
o l r
here(seebelow).Notethat(19)isahyperbolicsystemwhichfullydecomposesinto
linearly degenerate eigenfields, so R has an analytic expression (see formulas in
[4,5]).Notealso:theRiemannsolverRisconsistentundertheCFLcondition
∆tn≤ 1inf 1 max(cid:0)un−c(qn,qn )/hn,un+c (qn,qn )/hn (cid:1). (21)
2i∈Z∆xi i l i i+1 i i r i i+1 i+1
Itremainstospecifyachoiceoffunctionsc,c preservingU(cid:96)andensuring(14).
l r
6 Se´bastienBoyaval
Althoughitisnotclearwhetherourconstructionallowsonetoapproximateso-
lutionsonanytimerangest ∈[0,T),sincetheseries∑ ∆tn maybeboundeduni-
n
formlyforallspace-gridchoice(sup |un|maygrowunboundedlyasn→∞),speci-
i i
fyingsuchc,c fullydefinesacomputablescheme.Inparticular,(15)thenimplies
l r
that(12)atstep(ii)alwayshasatleastonesolutionqn+1∈U(cid:96) forany∆tn fixedat
i
step(i).(ThiscanbeshownusingBrouwerfixed-pointtheoremlikein[1].)
Note however a difficulty here for FENE-P fluids with c,c . Suliciu relaxation
l r
approach(19)wasretainedatstep(i)becausethesolveroftenallowsonetopreserve
invariantdomainslikeU(cid:96)andadissipationrule(14)throughwell-chosenc,c ,see
l r
e.g. [3, 4, 5]. Indeed, on noting the exact Riemann solution to (19), to get (14) on
choosingG(ql,qr)=u(cid:16)h(cid:0)u22 +eˆ(cid:1)+π(cid:17)|R(0,ql,qr),itisenoughthat∀ql,qr∈U(cid:96)
q :=LR(ξ,Q,Q )∈U(cid:96)andh2∂ | P(q )≤c2, ∀ξ ∈R (22)
ξ l r ξ h h2−2ζσxx,h2ζ−2σzz ξ ξ
usingc =c(q,q )ifξ<=u∗andc =c (q,q )ifξ>u∗withu∗:=clul+πl+crur−πr.
ξ l l r ξ r l r cl+cr
Onecaneasilyproposec,c satisfyingthefirstconditionin(22),i.e.
l r
(cid:18) (cid:19)
1 1 c (u −u)+π −π
r r l l r
= 1+ >0 (23)
h∗ h (c/h)(c +c )
l l l l l r
(cid:18) (cid:19)
1 1 c(u −u)+π −π
l r l r l
= 1+ >0 (24)
h∗ h (c /h )(c +c )
r r r r l r
asusualforSaint-Venantsystems,plustheadmissibilityconditions(o=l/r)
(h∗)2(1−ζ)(h )2(ζ−1)σ +(h∗)2(ζ−1)(h )2(1−ζ)σ <(cid:96) (25)
o o zz,o o o xx,o
foranyσ ,σ >0satisfyingσ +σ <(cid:96)(FENE-Pfluids,seebelow).But
zz,o xx,o zz,o xx,o
(cid:113)
thesecondconditionisusuallytreatedforφ :h→h ∂ | Pmono-
o h ho2−2ζσxx,o,h2oζ−2σzz,o
tone. Unfortunately, a lengthy (but easy) computation shows that the latter is not
monotonehere,sothestandardmethodtochoosec,c aprioridoesnotapply.
l r
2.3 Choiceofrelaxationparameter
(cid:112)
Let us treat the first part of (22) as usual and define c = max(h ∂ P(q ) :=
o o h o
h a ,c˜ ),o=l/rsuchthatthefunctionsc˜ (q,q )ensure(23–24)and(25).
o o o o l r
First,letusinspect(23–24)classicallyfollowing[7,section3.3].DenotingaY =
l l
(u −u ) + (πr−πl)+ ≥0, a Y =(u −u ) + (πl−πr)+ ≥0 so 1 ≥ 1−hoaoYo/co,
l r + hlal+hrar r r l r + hlal+hrar h∗o ho
it then holds (h∗)−1 ≥(h )−1y >0 with y :=1− Yo ∈(αo−1,1] provided
o o o o 1+αoYo αo
one chooses c˜ >0 such that c ≥h a (1+α Y ) for α >1, which yields h∗ ∈
o o o o o o o o
(0,h /y ]thus(23–24)inparticular.
o o
Ontheotherhand,letusnowinspect(25),whichrewriteswithh∗>0
o
Saint-Venant/FENE-Pshallowflows 7
(cid:16) (cid:112) (cid:112) (cid:17)
w A +w−1B <1⇔2A w ∈ 1− 1−4A B ,1+ 1−4A B ⊂R (26)
o o o o o o o o o o >0
withw =(h∗/h )2(1−ζ),A =σ /(cid:96),B =σ /(cid:96)positivesuchthatA +B <1
o o o o zz,o o xx,o o o
(hence A B ≤A (1−A )≤ 1) and 2(1−ζ)∈[1,2]. The upper-bound in (26) is
o o o o 4
1 1
satisfiedwithαo=(w+o)2(1−ζ)/((w+o)2(1−ζ) −1)>1,onnoting
(w+o)2(11−ζ) :=(cid:16)(1+(cid:112)1−4AoBo)/(2Ao)(cid:17)2(11−ζ) ≥ αoα−o1 ≥1/yo≥h∗o/ho. (27)
√
Itremainstoensurethelowerboundin(26).Obviously,w−:= 1− 1−4AoBo <1so
o 2Ao
oneonlyneedstoinspectthecaseh∗≤h .Now,withaW =(u −u) +(πl−πr)+ ≥
o o l l r l + hlal+hrar
0,arWr=(ur−ul)++h(πlarl−+πhlr)a+r ≥0,ifco≥hoaoWo((w−o)−2(11−ζ)−1)−1thenholds
(w−o)2(11−ζ) ≤(1+aohoWo/co)−1≤h∗o/ho.
Intheend,weclaimthefollowingchoices
(cid:18) (cid:18) (cid:19) (cid:18) (cid:19)(cid:19)
(π −π) (π −π )
r l + l r +
c =h max a +α (u −u ) + ,β (u −u) + (28)
l l l l l r + l r l +
ha +h a ha +h a
l l r r l l r r
(cid:18) (cid:18) (cid:19) (cid:18) (cid:19)(cid:19)
(π −π ) (π −π)
l r + r l +
c =h max a +α (u −u ) + ,β (u −u) + (29)
r r r r l r + r r l +
ha +h a ha +h a
l l r r l l r r
(cid:112)
satisfysimultaneously(23–24)and(25)inacompatiblewaywitha = ∂ P(q ),
o h o
1 1 1 1
αo√=max(2,(w+o)2(1−ζ)/((√w+o)2(1−ζ)−1)),βo=(w−o)2(1−ζ)/(1−(w−o)2(1−ζ))),w−o =
(cid:96)− (cid:96)−4σzz,oσxx,o, w+ = (cid:96)+ (cid:96)−4σzz,oσxx,o, for o=l/r. Moreover, note that we have
2σzz,o o 2σzz,o
chosenα suchthatallsubcharacteristicconditions(22)aresatisfiedinthe(cid:96)→∞
o
limit, hence also the free-energy dissipation (15). Indeed, φ is monotone in the
o
(cid:96)→∞limitandonecanthenapplythestandardmethodtochoosec,c [5].
l r
3 Numericalillustrations
Wenumericallyapproximateont∈[0,.1]thesolutiontoaRiemannproblemwith
(cid:40)
(h,u,σ ,σ )=(1,0,1,1) x<.5
l l xx,l zz,l
(h ,u ,σ ,σ )=(.1,0,1,1) x>.5
r r xx,r zz,r
as initial condition when g=10, ζ =0, G=.1, λ =.1. In Fig. 1, we show the
initialconditionandtheresultatt=.1when∆x=2−8 for(cid:96)=10,100,1000.Note
theinfluenceoftheparameter(cid:96)onthestretchσ +σ .Oncomputingnumerically
xx zz
8 Se´bastienBoyaval
Fig.1 Top:h(left)andu(right),bottom:σxxandσzz.
thefree-energydissipationwiththechoiceofrelaxationparameterabove,wehave
neverobservedthewrongsign,whilethetime-stepdidnotgotozero.
References
1. Barrett,J.W.,Boyaval,S.:Existenceandapproximationofa(regularized)Oldroyd-Bmodel.
M3AS21(9),1783–1837(2011).
2. Berthon, C., Coquel, F., LeFloch, P.G.: Why many theories of shock waves are necessary:
kineticrelationsfornon-conservativesystems.ProceedingsoftheRoyalSocietyofEdinburgh:
SectionAMathematics142,1–37(2012)
3. Bouchut,F.:EntropysatisfyingfluxvectorsplittingsandkineticBGKmodels. Numerische
Mathematik94,623–672(2003).
4. Bouchut,F.:Nonlinearstabilityoffinitevolumemethodsforhyperbolicconservationlaws
andwell-balancedschemesforsources. FrontiersinMathematics.Birkha¨userVerlag,Basel
(2004)
5. Bouchut,F.,Boyaval,S.:Anewmodelforshallowviscoelasticfluids. M3AS23(08),1479–
1526(2013).
6. Bouchut,F.,Boyaval,S.:Unifiedderivationofthin-layerreducedmodelsforshallowfree-
surfacegravityflowsofviscousfluids. EuropeanJournalofMechanics-B/Fluids55,Part1,
116–131(2016).
7. Bouchut,F.,Klingenberg,C.,Waagan,K.:AmultiwaveapproximateRiemannsolverforideal
MHDbasedonrelaxationII:numericalimplementationwith3and5waves. Numer.Math.
115(4),647–679(2010).
8. Boyaval,S.:Johnson-Segalman–Saint-Venantequationsforviscoelasticshallowflowsinthe
elasticlimit. ArXive-prints(2016)
9. Dafermos,C.M.:Hyperbolicconservationlawsincontinuumphysics,vol.GM325.Springer-
Verlag,Berlin(2000)
10. Harten,A.,Lax,P.D.,vanLeer,B.:OnupstreamdifferencingandGodunov-typeschemesfor
hyperbolicconservationlaws. SIAMRev.25(1),35–61(1983).
11. LeFloch, P.G.: Hyperbolic systems of conservation laws. Lectures in Mathematics ETH
Zu¨rich.Birkha¨userVerlag,Basel(2002).Thetheoryofclassicalandnonclassicalshockwaves