siefer92 |
February 19, 2020 14:58 |
Issue with *Cloud.evolve()
Howdy Yall,
I was able to get the kinematic cloud merged with Sonic foam and achieved some visually agreeable results for the ma 3 forward facing step tutorial. Momentum is coupled and the results look, upon visual inspection, good.
The problem is that I am trying to incorporate heat transfer now and I switched from kinematicCloud -> thermoCloud. From my understanding thermoCloud should just add some functionality and act as a wrapper to kinematicCloud. I made some modifications to the creatFields.H file and coupled the cloud to the continuum energy equation through (+ dustCloud.Sh(e))
My code compiles and I go to run my test case, the simulation proceeds as normal except for the fact that once the particle injection time is reached, the simulation crashes without an error.
Specifically it crashes on the evolution step of a particular time step.
I am attaching my case files in this top post and will make a comment that has my code files
Pardon the language in my log text.
Code:
Evolving dustCloud
Did OpenFOAM Fucking Freeze Again?
Solving 2-D cloud dustCloud
Cloud: dustCloud
Current number of parcels = 0
Current mass in system = 0
Linear momentum = (0 0 0)
|Linear momentum| = 0
Linear kinetic energy = 0
Average particle per parcel = 0
Injector model1:
- parcels added = 0
- mass introduced = 0
Parcel fate: system (number, mass)
- escape = 0, 0
Parcel fate: patch (number, mass) inlet
- escape = 0, 0
- stick = 0, 0
Parcel fate: patch (number, mass) outlet
- escape = 0, 0
- stick = 0, 0
Parcel fate: patch (number, mass) bottom
- escape = 0, 0
- stick = 0, 0
Parcel fate: patch (number, mass) top
- escape = 0, 0
- stick = 0, 0
Parcel fate: patch (number, mass) obstacle
- escape = 0, 0
- stick = 0, 0
Parcel fate: patch (number, mass) defaultFaces
- escape = 0, 0
- stick = 0, 0
Temperature min/max = 0, 0
Not Yet.
ExecutionTime = 10.93 s ClockTime = 11 s
Time = 0.1002
Courant Number mean: 0.0477675 max: 0.0615826
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
PIMPLE: iteration 1
smoothSolver: Solving for Ux, Initial residual = 0.00116844, Final residual = 3.81639e-07, No Iterations 1
smoothSolver: Solving for Uy, Initial residual = 0.00166205, Final residual = 1.54533e-06, No Iterations 1
smoothSolver: Solving for e, Initial residual = 0.00103972, Final residual = 7.14875e-07, No Iterations 1
smoothSolver: Solving for p, Initial residual = 0.000908961, Final residual = 6.32946e-07, No Iterations 1
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 1.6622e-07, global = -4.03005e-08, cumulative = 1.31019e-05
Evolving dustCloud
Did OpenFOAM Fucking Freeze Again?
Solving 2-D cloud dustCloud
-------------------------------------------------------
Primary job terminated normally, but 1 process returned
a non-zero exit code.. Per user-direction, the job has been aborted.
-------------------------------------------------------
--------------------------------------------------------------------------
mpirun detected that one or more processes exited with non-zero status, thus causing
the job to be terminated. The first process to do so was:
Process name: [[18412,1],1]
Exit code: 144
--------------------------------------------------------------------------
It should be noted that this error is just an MPI error and running in serial will have the code crash on "Solving 2-D cloud dustCloud"
I've seen it mentioned on the forums that the lagragian wall rebound treatment can couse hanging, though I am crashing, if the mesh is trash near the wall. Despite not using the rebound model I checked the mesh just to make sure it wasn't the cause using> checkMesh -allGeometry
The output of the check is:
Code:
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: v1912 |
| \\ / A nd | Website: www.openfoam.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Build : _f3950763fe-20191219 OPENFOAM=1912
Arch : "LSB;label=64;scalar=64"
Exec : checkMesh -allGeometry
Date : Feb 19 2020
Time : 12:37:15
Host : AEK997
PID : 27748
I/O : uncollated
Case : /forwardStepLPTT
nProcs : 1
trapFpe: Floating point exception trapping enabled (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster (fileModificationSkew 10)
allowSystemOperations : Allowing user-supplied system call operations
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time
Create mesh for time = 0
Enabling all geometry checks.
Time = 0
Mesh stats
points: 32898
internal points: 0
faces: 64832
internal faces: 31936
cells: 16128
faces per cell: 6
boundary patches: 6
point zones: 0
face zones: 0
cell zones: 0
Overall number of cells of each type:
hexahedra: 16128
prisms: 0
wedges: 0
pyramids: 0
tet wedges: 0
tetrahedra: 0
polyhedra: 0
Checking topology...
Boundary definition OK.
Cell to face addressing OK.
Point usage OK.
Upper triangular ordering OK.
Face vertices OK.
Number of regions: 1 (OK).
Checking patch topology for multiply connected surfaces...
Patch Faces Points Surface topology Bounding box
inlet 80 162 ok (non-closed singly connected) (0 0 -0.05) (0 1 0.05)
outlet 64 130 ok (non-closed singly connected) (3 0.2 -0.05) (3 1 0.05)
bottom 48 98 ok (non-closed singly connected) (0 0 -0.05) (0.6 0 0.05)
top 240 482 ok (non-closed singly connected) (0 1 -0.05) (3 1 0.05)
obstacle 208 418 ok (non-closed singly connected) (0.6 0 -0.05) (3 0.2 0.05)
defaultFaces 32256 32898 ok (non-closed singly connected) (0 0 -0.05) (3 1 0.05)
Checking faceZone topology for multiply connected surfaces...
No faceZones found.
Checking basic cellZone addressing...
No cellZones found.
Checking geometry...
Overall domain bounding box (0 0 -0.05) (3 1 0.05)
Mesh has 2 geometric (non-empty/wedge) directions (1 1 0)
Mesh has 2 solution (non-empty) directions (1 1 0)
All edges aligned with or perpendicular to non-empty directions.
Boundary openness (1.18817e-18 6.76884e-17 -2.78861e-15) OK.
Max cell openness = 1.73472e-16 OK.
Max aspect ratio = 1 OK.
Minimum face area = 0.00015625. Maximum face area = 0.00125. Face area magnitudes OK.
Min volume = 1.5625e-05. Max volume = 1.5625e-05. Total volume = 0.252. Cell volumes OK.
Mesh non-orthogonality Max: 0 average: 0
Non-orthogonality check OK.
Face pyramids OK.
Max skewness = 2.13163e-13 OK.
Coupled point location match (average 0) OK.
Face tets OK.
Min/max edge length = 0.0125 0.1 OK.
All angles in faces OK.
Face flatness (1 = flat, 0 = butterfly) : min = 1 average = 1
All face flatness OK.
Cell determinant (wellposedness) : minimum: 0.125 average: 0.490118
Cell determinant check OK.
Concave cell check OK.
Face interpolation weight : minimum: 0.5 average: 0.5
Face interpolation weight check OK.
Face volume ratio : minimum: 1 average: 1
Face volume ratio check OK.
Mesh OK.
End
here is my constant/dustCloudProperties file
Code:
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 4.x |
| \\ / A nd | Web: www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object kinematicCloudProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
solution
{
active true;
coupled true;
transient yes;
cellValueSourceCorrection off;
maxCo 0.3;
interpolationSchemes
{
rho cell;
U cellPoint;
thermo:mu cell;
T cell;
Cp cell;
kappa cell;
p cell;
G cell;
DUcDt cell;
}
integrationSchemes
{
U Euler;
T analytical;
}
sourceTerms
{
schemes
{
rho semiImplicit 1;
U semiImplicit 1;
Yi semiImplicit 1;
h semiImplicit 1;
}
}
}
constantProperties
{
rho0 8800;
youngsModulus 1.3e5;
poissonsRatio 0.35;
T0 1;
Cp0 4187;
epsilon0 1;
f0 0.5;
}
subModels
{
particleForces
{
sphereDrag;
gravity;
pressureGradient
{
U U;
}
}
injectionModels
{
model1
{
type patchInjection;
patch inlet;
duration 1;
parcelsPerSecond 1390885;
massTotal 40;
parcelBasisType fixed;
flowRateProfile constant 1;
nParticle 1;
SOI 0.1;
U0 (3 0 0);
T0 3;
sizeDistribution
{
type fixedValue;
fixedValueDistribution
{
value 0.00007;
}
}
}
}
dispersionModel none;
patchInteractionModel standardWallInteraction;
standardWallInteractionCoeffs
{
type rebound;
e 0.97;
mu 0.09;
}
surfaceFilmModel none;
stochasticCollisionModel none;
radiation off;
heatTransferModel none;
RanzMarshallCoeffs
{
BirdCorrection true;
}
collisionModel none;
pairCollisionCoeffs
{
maxInteractionDistance 0.00007;
writeReferredParticleCloud no;
pairModel pairSpringSliderDashpot;
pairSpringSliderDashpotCoeffs
{
useEquivalentSize no;
alpha 0.12;
b 1.5;
mu 0.52;
cohesionEnergyDensity 0;
collisionResolutionSteps 12;
};
wallModel wallSpringSliderDashpot;
wallSpringSliderDashpotCoeffs
{
useEquivalentSize no;
collisionResolutionSteps 12;
youngsModulus 1e10;
poissonsRatio 0.23;
alpha 0.12;
b 1.5;
mu 0.43;
cohesionEnergyDensity 0;
};
}
}
cloudFunctions
{
voidFraction1
{
type voidFraction;
}
}
UEQN
Code:
// Solve the Momentum equation
MRF.correctBoundaryVelocity(U);
fvVectorMatrix UEqn
(
fvm::ddt(rho, U) + fvm::div(phi, U)
+ MRF.DDt(rho, U)
+ turbulence->divDevRhoReff(U)
==
fvOptions(rho, U) + dustParcels.SU(U)
);
UEqn.relax();
fvOptions.constrain(UEqn);
if (pimple.momentumPredictor())
{
solve(UEqn == -fvc::grad(p));
fvOptions.correct(U);
K = 0.5*magSqr(U);
}
EEQN
Code:
{
fvScalarMatrix EEqn
(
fvm::ddt(rho, e) + fvm::div(phi, e)
+ fvc::ddt(rho, K) + fvc::div(phi, K)
+ fvc::div(fvc::absolute(phi/fvc::interpolate(rho), U), p, "div(phiv,p)")
- fvm::laplacian(turbulence->alphaEff(), e)
==
fvOptions(rho, e) + dustParcels.Sh(e)
);
EEqn.relax();
fvOptions.constrain(EEqn);
EEqn.solve();
fvOptions.correct(e);
thermo.correct();
}
My Create Fields File
Code:
{
#include "readGravitationalAcceleration.H"
Info<< "Reading thermophysical properties\n" << endl;
autoPtr<psiThermo> pThermo(psiThermo::New(mesh));
psiThermo& thermo = pThermo();
thermo.validate(args.executable(), "h", "e");
SLGThermo slgThermo(mesh, thermo);
volScalarField& p = thermo.p();
// const volScalarField& T = thermo.T();
// const volScalarField& psi = thermo.psi();
multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
fields.add(thermo.he());
volScalarField rho
(
IOobject
(
"rho",
runTime.timeName(),
mesh
),
thermo.rho()
);
Info<< "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
#include "compressibleCreatePhi.H"
mesh.setFluxRequired(p.name());
Info<< "Creating turbulence model\n" << endl;
autoPtr<compressible::turbulenceModel> turbulence
(
compressible::turbulenceModel::New
(
rho,
U,
phi,
thermo
)
);
#include "createK.H"
#include "createMRF.H"
Info<< "Reading transportProperties\n" <<endl;
IOdictionary transportProperties
(
IOobject
(
"transportProperties",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
);
dimensionedScalar rhoInfValue
(
transportProperties.lookup("rhoInf")
);
dimensionedScalar invrhoInf("invrhoInf",(1.0/rhoInfValue));
volScalarField rhoInf
(
IOobject
(
"rho",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
rhoInfValue
);
volScalarField mu
(
IOobject
(
"mu",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
turbulence->nu()*rhoInfValue
);
const word thermoCloudName
(
args.optionLookupOrDefault<word>("CloudName", "dustCloud")
);
Info<< "\nConstructing dust cloud" << endl;
basicThermoCloud dustParcels
(
thermoCloudName,
rhoInf,
U,
g,
slgThermo
);
#include "createFvOptions.H"
}
and the modified solver file
Code:
/*---------------------------------------------------------------------------
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "SLGThermo.H"
#include "psiThermo.H"
#include "turbulentFluidThermoModel.H"
#include "pimpleControl.H"
#include "fvOptions.H"
#include "basicThermoCloud.H"
#include "radiationModel.H"
#include "pressureControl.H"
#include "localEulerDdtScheme.H"
#include "fvcSmooth.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addNote
(
"Transient solver for trans-sonic/supersonic, turbulent flow"
" of a compressible gas."
);
#include "postProcess.H"
#include "addCheckCaseOptions.H"
#include "setRootCaseLists.H"
#include "createTime.H"
#include "createMesh.H"
// #include "readGravitationalAcceleration.H"
#include "createControl.H"
#include "createFields.H"
#include "createFieldRefs.H"
#include "initContinuityErrs.H"turbulence->validate();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Info<< "\nStarting time loop\n" << endl;
while (runTime.loop())
{
Info<< "Time = " << runTime.timeName() << nl << endl;
#include "compressibleCourantNo.H"
#include "rhoEqn.H"
// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#include "UEqn.H"
#include "EEqn.H"
// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}
if (pimple.turbCorr())
{
turbulence->correct();
}
}
rho = thermo.rho();
Info << "\nEvolving " << dustParcels.name() << endl;
Info << "\nDid OpenFOAM Fucking Freeze Again?" << endl;
dustParcels.evolve();
Info << "\nNot Yet." << endl;
runTime.write();
runTime.printExecutionTime(Info);
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //
Given the lack of an error code, the complicated nature of the Lagrangian source code and my, relative, inexperience with openFOAM solver programming I'm quite at a loss here.
I know that this has to caused by user error, but I haven't pinned down a specific reason for why within the lagrangian evolve() function
if you want additional files I'll upload them in the comments if it is requested.
Thanks yall!
|