Skip to content
GitLab
Menu
Projects
Groups
Snippets
Loading...
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in / Register
Toggle navigation
Menu
Open sidebar
Development
openfoam
Commits
7f22e310
Commit
7f22e310
authored
Dec 11, 2008
by
henry
Browse files
Added correctPhi to compressibleInterDyMFoam.
parent
a2c058d5
Changes
3
Hide whitespace changes
Inline
Side-by-side
applications/solvers/multiphase/compressibleInterDyMFoam/compressibleInterDyMFoam.C
View file @
7f22e310
...
...
@@ -74,24 +74,26 @@ int main(int argc, char *argv[])
Info
<<
"Time = "
<<
runTime
.
timeName
()
<<
nl
<<
endl
;
scalar
timeBeforeMeshUpdate
=
runTime
.
elapsedCpuTime
();
{
// Store divU from the previous mesh for the correctPhi
volScalarField
divU
=
fvc
::
div
(
phi
);
// Do any mesh changes
mesh
.
update
();
scalar
timeBeforeMeshUpdate
=
runTime
.
elapsedCpuTime
();
if
(
mesh
.
changing
())
{
Info
<<
"Execution time for mesh.update() = "
<<
runTime
.
elapsedCpuTime
()
-
timeBeforeMeshUpdate
<<
" s"
<<
endl
;
// Do any mesh changes
mesh
.
update
();
gh
=
g
&
mesh
.
C
();
ghf
=
g
&
mesh
.
Cf
();
}
if
(
mesh
.
changing
())
{
Info
<<
"Execution time for mesh.update() = "
<<
runTime
.
elapsedCpuTime
()
-
timeBeforeMeshUpdate
<<
" s"
<<
endl
;
}
if
(
mesh
.
changing
()
&&
correctPhi
)
{
//***HGW#include "correctPhi.H"
if
(
mesh
.
changing
()
&&
correctPhi
)
{
#include "correctPhi.H"
}
}
// Make the fluxes relative to the mesh motion
...
...
@@ -102,6 +104,12 @@ int main(int argc, char *argv[])
#include "meshCourantNo.H"
}
if
(
mesh
.
changing
())
{
gh
=
g
&
mesh
.
C
();
ghf
=
g
&
mesh
.
Cf
();
}
turbulence
->
correct
();
// --- Outer-corrector loop
...
...
applications/solvers/multiphase/compressibleInterDyMFoam/correctPhi.H
0 → 100644
View file @
7f22e310
{
#include "continuityErrs.H"
volScalarField
pcorr
(
IOobject
(
"pcorr"
,
runTime
.
timeName
(),
mesh
,
IOobject
::
NO_READ
,
IOobject
::
NO_WRITE
),
mesh
,
dimensionedScalar
(
"pcorr"
,
pd
.
dimensions
(),
0
.
0
),
pcorrTypes
);
dimensionedScalar
rAUf
(
"(1|A(U))"
,
dimTime
/
rho
.
dimensions
(),
1
.
0
);
adjustPhi
(
phi
,
U
,
pcorr
);
for
(
int
nonOrth
=
0
;
nonOrth
<=
nNonOrthCorr
;
nonOrth
++
)
{
fvScalarMatrix
pcorrEqn
(
fvm
::
laplacian
(
rAUf
,
pcorr
)
==
fvc
::
div
(
phi
)
-
divU
);
pcorrEqn
.
solve
();
if
(
nonOrth
==
nNonOrthCorr
)
{
phi
-=
pcorrEqn
.
flux
();
}
}
#include "continuityErrs.H"
}
applications/solvers/multiphase/compressibleInterDyMFoam/createFields.H
View file @
7f22e310
...
...
@@ -150,3 +150,14 @@
(
incompressible
::
turbulenceModel
::
New
(
U
,
phi
,
twoPhaseProperties
)
);
wordList
pcorrTypes
(
pd
.
boundaryField
().
types
());
for
(
label
i
=
0
;
i
<
pd
.
boundaryField
().
size
();
i
++
)
{
if
(
pd
.
boundaryField
()[
i
].
fixesValue
())
{
pcorrTypes
[
i
]
=
fixedValueFvPatchScalarField
::
typeName
;
}
}
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment