The finite element formulation of the Lagrangian form of the mass conservation equation is given
by:
図 1 .
d
ρ
d
t
|
x
=
−
(
ρ
/
V
)
d
V
d
t
|
x
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbwvMCKf
MBHbqefqvATv2CG4uz3bIuV1wyUbqedmvETj2BSbqefm0B1jxALjhi
ov2DaebbnrfifHhDYfgasaacH8srps0lbbf9q8WrFfeuY=Hhbbf9v8
qqaqFr0xc9pk0xbba9q8WqFfea0=yr0RYxir=Jbba9q8aq0=yq=He9
q8qqQ8frFve9Fve9Ff0dmeaacaGacmGadaWaaiqacaabaiaafaaake
aadaWcaaqaaiaadsgacqaHbpGCaeaacaWGKbGaamiDaaaadaabbaqa
aiaadIhaaiaawEa7aiabg2da9iabgkHiTmaabmaabaGaeqyWdiNaai
4laiaadAfaaiaawIcacaGLPaaadaWcaaqaaiaadsgacaWGwbaabaGa
amizaiaadshaaaWaaqqaaeaacaWG4baacaGLhWoaaaa@4D14@
When transformed into the ALE formulation it gives:
図 2 .
Applying a Galerkin variation form for the solution of
式 2 :
図 3 .
Where,
ψ
is the Weighting function.
Using a finite volume formulation
Where,
ψ
=1
ρ
= constant density over control volume
V
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn
hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr
4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9
vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x
fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOvaaaa@36D1@
.
Therefore:
図 4 .
∫
V
∂
ρ
∂
t
d
V
+
∫
V
ρ
∂
ν
K
∂
x
K
d
V
=
0
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbwvMCKf
MBHbqefqvATv2CG4uz3bIuV1wyUbqedmvETj2BSbqefm0B1jxALjhi
ov2DaebbnrfifHhDYfgasaacH8srps0lbbf9q8WrFfeuY=Hhbbf9v8
qqaqFr0xc9pk0xbba9q8WqFfea0=yr0RYxir=Jbba9q8aq0=yq=He9
q8qqQ8frFve9Fve9Ff0dmeaacaGacmGadaWaaiqacaabaiaafaaake
aadaWdrbqaamaalaaabaGaeyOaIyRaeqyWdihabaGaeyOaIyRaamiD
aaaacaWGKbGaamOvaiabgUcaRmaapefabaGaeqyWdihaleaacaWGwb
aabeqdcqGHRiI8aaWcbaGaamOvaaqab0Gaey4kIipakmaalaaabaGa
eyOaIyRaeqyVd42aaSbaaSqaaiaadUeaaeqaaaGcbaGaeyOaIyRaam
iEamaaBaaaleaacaWGlbaabeaaaaGccaWGKbGaamOvaiabg2da9iaa icdaaaa@5447@
Using the divergence theorem leads to:
図 5 .
∫
V
∂
ρ
∂
t
d
V
+
∫
S
ρ
(
ν
⋅
j
n
j
)
d
S
=
0
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbwvMCKf
MBHbqefqvATv2CG4uz3bIuV1wyUbqedmvETj2BSbqefm0B1jxALjhi
ov2DaebbnrfifHhDYfgasaacH8srps0lbbf9q8WrFfeuY=Hhbbf9v8
qqaqFr0xc9pk0xbba9q8WqFfea0=yr0RYxir=Jbba9q8aq0=yq=He9
q8qqQ8frFve9Fve9Ff0dmeaacaGacmGadaWaaiqacaabaiaafaaake
aadaWdrbqaamaalaaabaGaeyOaIyRaeqyWdihabaGaeyOaIyRaamiD
aaaacaWGKbGaamOvaiabgUcaRmaapefabaGaeqyWdihaleaacaWGtb
aabeqdcqGHRiI8aaWcbaGaamOvaaqab0Gaey4kIipakmaabmaabaGa
eqyVd42aaSraaSqaaiaadQgaaeqaaOGaeyyXICTaaGPaVlaayIW7ca
WGUbWaaSbaaSqaaiaadQgaaeqaaaGccaGLOaGaayzkaaGaamizaiaa
dofacqGH9aqpcaaIWaaaaa@5889@
Further expansion gives:
図 6 .
d
d
t
∫
v
ρ
d
V
=
∫
s
ρ
(
w
j
−
v
j
)
n
j
d
S
︸
m
a
s
s
f
l
u
x
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbwvMCKf
MBHbqefqvATv2CG4uz3bIuV1wyUbqedmvETj2BSbqefm0B1jxALjhi
ov2DaebbnrfifHhDYfgasaacH8srps0lbbf9q8WrFfeuY=Hhbbf9v8
qqaqFr0xc9pk0xbba9q8WqFfea0=yr0RYxir=Jbba9q8aq0=yq=He9
q8qqQ8frFve9Fve9Ff0dmeaacaGacmGadaWaaiqacaabaiaafaaake
aadaWcaaqaaiaadsgaaeaacaWGKbGaamiDaaaadaWdrbqaaiabeg8a
YjaadsgacaWGwbGaeyypa0Zaa8quaeaadaagaaqaaiabeg8aYnaabm
aabaGaam4DamaaBaaaleaacaWGQbaabeaakiabgkHiTiaadAhadaWg
aaWcbaGaamOAaaqabaaakiaawIcacaGLPaaacaWGUbWaaSbaaSqaai
aadQgaaeqaaOGaamizaiaadofaaSqaaiaad2gacaWGHbGaam4Caiaa
dohacaaMc8UaamOzaiaadYgacaWG1bGaamiEaaGccaGL44paaSqaai
aadohaaeqaniabgUIiYdaaleaacaWG2baabeqdcqGHRiI8aaaa@5E39@
This formula is still valid if density
ρ
is not assumed uniform over volume
V
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn
hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr
4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9
vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x
fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaamOvaaaa@36D1@
.
図 7 . Mass Flux Across a Surface
The density,
ρ
i
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn
hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr
4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9
vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x
fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS
baaSqaaiaadMgaaeqaaaaa@38D0@
, is given computed:
図 8 .
ρ
i
=
1
2
ρ
I
{
1
+
η
s
i
g
n
(
ϕ
i
)
}
+
1
2
ρ
J
{
1
−
η
s
i
g
n
(
ϕ
i
)
}
Where,
0
≤
η
≤
1
is the upwind coefficient given on the input card.
If
η
=0, there is no upwind.
Therefore,
ρ
i
=
ρ
I
+
ρ
J
2
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn
hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr
4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9
vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x
fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS
baaSqaaiaadMgaaeqaaOGaeyypa0ZaaSaaaeaacqaHbpGCdaWgaaWc
baGaamysaaqabaGccqGHRaWkcqaHbpGCdaWgaaWcbaGaamOsaaqaba aakeaacaaIYaaaaaaa@4117@
.
If
η
=1, there is full upwind.
The smaller the upwind factor, the faster the solution; however, the solution is more stable with
a large upwind factor. This upwind coefficient can be tuned with parameter from
/UPWIND keyword (not recommended, this keyword has been obsolete as of
version 2018).
For a free surface:
ρ
J
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn
hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr
4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9
vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x
fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS
baaSqaaiaadMgaaeqaaaaa@38D0@
=
ρ
I
MathType@MTEF@5@5@+=
feaagKart1ev2aaatCvAUfeBSjuyZL2yd9gzLbvyNv2CaerbuLwBLn
hiov2DGi1BTfMBaeXatLxBI9gBaerbd9wDYLwzYbItLDharqqtubsr
4rNCHbGeaGqiVu0Je9sqqrpepC0xbbL8F4rqqrFfpeea0xe9Lq=Jc9
vqaqpepm0xbba9pwe9Q8fs0=yqaqpepae9pg0FirpepeKkFr0xfr=x
fr=xb9adbaqaaeGaciGaaiaabeqaamaabaabaaGcbaGaeqyWdi3aaS
baaSqaaiaadMgaaeqaaaaa@38D0@