From 9c32b96e0173eafadfb085b45d3af3fa1b7caefc Mon Sep 17 00:00:00 2001
From: graham <graham.macpherson@strath.ac.uk>
Date: Fri, 6 Mar 2009 15:53:01 +0000
Subject: [PATCH] Added equipartition internal energy distribution sampling
function and calling on particle intialisation, Maxwellian wall collision and
free stream injection. Interface with WallInteractionModel modified to pass
a scalar& for the internal energy and a label to specify the typeId.
---
.../clouds/Templates/DsmcCloud/DsmcCloud.C | 49 +++++++++++++++++--
.../clouds/Templates/DsmcCloud/DsmcCloud.H | 11 ++++-
.../parcels/Templates/DsmcParcel/DsmcParcel.C | 13 +++--
.../LarsenBorgnakkeVariableHardSphere.C | 4 +-
.../FreeStream/FreeStream.C | 8 +--
.../MaxwellianThermal/MaxwellianThermal.C | 9 +++-
.../MaxwellianThermal/MaxwellianThermal.H | 3 +-
.../SpecularReflection/SpecularReflection.C | 3 +-
.../SpecularReflection/SpecularReflection.H | 3 +-
.../WallInteractionModel.H | 3 +-
10 files changed, 86 insertions(+), 20 deletions(-)
diff --git a/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.C b/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.C
index 0e39c2ff33a..927f2199ec9 100644
--- a/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.C
+++ b/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.C
@@ -172,9 +172,11 @@ void Foam::DsmcCloud<ParcelType>::initialise
cP.mass()
);
- scalar Ei =
- 0.5*cP.internalDegreesOfFreedom()
- *kb*temperature;
+ scalar Ei = equipartitionInternalEnergy
+ (
+ temperature,
+ cP.internalDegreesOfFreedom()
+ );
U += velocity;
@@ -690,6 +692,47 @@ Foam::vector Foam::DsmcCloud<ParcelType>::equipartitionLinearVelocity
}
+template<class ParcelType>
+Foam::scalar Foam::DsmcCloud<ParcelType>::equipartitionInternalEnergy
+(
+ scalar temperature,
+ scalar iDof
+)
+{
+ scalar Ei = 0.0;
+
+ if
+ (
+ iDof < 2.0 + SMALL
+ && iDof > 2.0 - SMALL
+ )
+ {
+ // Special case for iDof = 2, i.e. diatomics;
+ Ei = -log(rndGen_.scalar01())*kb*temperature;
+ }
+ else
+ {
+ scalar a = 0.5*iDof - 1;
+
+ scalar energyRatio;
+
+ scalar P;
+
+ do
+ {
+ energyRatio = 10*rndGen_.scalar01();
+
+ P = pow((energyRatio/a), a)*exp(a - energyRatio);
+
+ } while (P < rndGen_.scalar01());
+
+ Ei = energyRatio*kb*temperature;
+ }
+
+ return Ei;
+}
+
+
template<class ParcelType>
void Foam::DsmcCloud<ParcelType>::dumpParticlePositions() const
{
diff --git a/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.H b/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.H
index 8e8fd5cbd74..b246b049904 100644
--- a/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.H
+++ b/src/lagrangian/dsmc/clouds/Templates/DsmcCloud/DsmcCloud.H
@@ -281,8 +281,17 @@ public:
scalar mass
);
- //- From the Maxwellian distribution:
+ //- Generate a random internal energy, sampled from the
+ // equilibrium distribution (Bird eqn 11.22 and 11.23 and
+ // adapting code from DSMC3.FOR)
+ scalar equipartitionInternalEnergy
+ (
+ scalar temperature,
+ scalar internalDegreesOfFreedom
+ );
+
+ // From the Maxwellian distribution:
//- Average particle speed
inline scalar maxwellianAverageSpeed
(
diff --git a/src/lagrangian/dsmc/parcels/Templates/DsmcParcel/DsmcParcel.C b/src/lagrangian/dsmc/parcels/Templates/DsmcParcel/DsmcParcel.C
index 93058e076e4..f9c13571894 100644
--- a/src/lagrangian/dsmc/parcels/Templates/DsmcParcel/DsmcParcel.C
+++ b/src/lagrangian/dsmc/parcels/Templates/DsmcParcel/DsmcParcel.C
@@ -109,25 +109,28 @@ void Foam::DsmcParcel<ParcelType>::hitWallPatch
{
const constantProperties& constProps(td.cloud().constProps(typeId_));
+ scalar m = constProps.mass();
+
// pre-interaction energy
- scalar preIE = 0.5*constProps.mass()*(U_ & U_) + Ei_;
+ scalar preIE = 0.5*m*(U_ & U_) + Ei_;
// pre-interaction momentum
- vector preIMom = constProps.mass()*U_;
+ vector preIMom = m*U_;
td.cloud().wallInteraction().correct
(
wpp,
this->face(),
U_,
- constProps.mass()
+ Ei_,
+ typeId_
);
// post-interaction energy
- scalar postIE = 0.5*constProps.mass()*(U_ & U_) + Ei_;
+ scalar postIE = 0.5*m*(U_ & U_) + Ei_;
// post-interaction momentum
- vector postIMom = constProps.mass()*U_;
+ vector postIMom = m*U_;
label wppIndex = wpp.index();
diff --git a/src/lagrangian/dsmc/submodels/BinaryCollisionModel/LarsenBorgnakkeVariableHardSphere/LarsenBorgnakkeVariableHardSphere.C b/src/lagrangian/dsmc/submodels/BinaryCollisionModel/LarsenBorgnakkeVariableHardSphere/LarsenBorgnakkeVariableHardSphere.C
index 1cf6fd90dcf..ff75da32b4e 100644
--- a/src/lagrangian/dsmc/submodels/BinaryCollisionModel/LarsenBorgnakkeVariableHardSphere/LarsenBorgnakkeVariableHardSphere.C
+++ b/src/lagrangian/dsmc/submodels/BinaryCollisionModel/LarsenBorgnakkeVariableHardSphere/LarsenBorgnakkeVariableHardSphere.C
@@ -81,9 +81,7 @@ Foam::scalar Foam::LarsenBorgnakkeVariableHardSphere<CloudType>::energyRatio
ChiBMinusOne
);
}
- }
-
- while (P < rndGen.scalar01());
+ } while (P < rndGen.scalar01());
return energyRatio;
}
diff --git a/src/lagrangian/dsmc/submodels/InflowBoundaryModel/FreeStream/FreeStream.C b/src/lagrangian/dsmc/submodels/InflowBoundaryModel/FreeStream/FreeStream.C
index a69f74aa994..1367198efac 100644
--- a/src/lagrangian/dsmc/submodels/InflowBoundaryModel/FreeStream/FreeStream.C
+++ b/src/lagrangian/dsmc/submodels/InflowBoundaryModel/FreeStream/FreeStream.C
@@ -237,9 +237,11 @@ void Foam::FreeStream<CloudType>::inflow()
U += velocity_;
- scalar Ei =
- 0.5*cloud.constProps(typeId).internalDegreesOfFreedom()
- *CloudType::kb*temperature_;
+ scalar Ei = cloud.equipartitionInternalEnergy
+ (
+ temperature_,
+ cloud.constProps(typeId).internalDegreesOfFreedom()
+ );
cloud.addNewParcel
(
diff --git a/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.C b/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.C
index a6ed8e57183..d8973fc5001 100644
--- a/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.C
+++ b/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.C
@@ -54,7 +54,8 @@ void Foam::MaxwellianThermal<CloudType>::correct
const wallPolyPatch& wpp,
const label faceId,
vector& U,
- scalar mass
+ scalar& Ei,
+ label typeId
)
{
label wppIndex = wpp.index();
@@ -102,6 +103,10 @@ void Foam::MaxwellianThermal<CloudType>::correct
scalar T = cloud.T().boundaryField()[wppIndex][wppLocalFace];
+ scalar mass = cloud.constProps(typeId).mass();
+
+ scalar iDof = cloud.constProps(typeId).internalDegreesOfFreedom();
+
U =
sqrt(CloudType::kb*T/mass)
*(
@@ -111,6 +116,8 @@ void Foam::MaxwellianThermal<CloudType>::correct
);
U += cloud.U().boundaryField()[wppIndex][wppLocalFace];
+
+ Ei = cloud.equipartitionInternalEnergy(T, iDof);
}
diff --git a/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.H b/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.H
index 1542d20ff05..a52a066aec2 100644
--- a/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.H
+++ b/src/lagrangian/dsmc/submodels/WallInteractionModel/MaxwellianThermal/MaxwellianThermal.H
@@ -78,7 +78,8 @@ public:
const wallPolyPatch& wpp,
const label faceId,
vector& U,
- scalar mass
+ scalar& Ei,
+ label typeId
);
};
diff --git a/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.C b/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.C
index 38e6621d30b..4de5507bf00 100644
--- a/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.C
+++ b/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.C
@@ -56,7 +56,8 @@ void Foam::SpecularReflection<CloudType>::correct
const wallPolyPatch& wpp,
const label faceId,
vector& U,
- scalar mass
+ scalar& Ei,
+ label typeId
)
{
vector nw = wpp.faceAreas()[wpp.whichFace(faceId)];
diff --git a/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.H b/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.H
index 67d0adf8b07..d2208bd3503 100644
--- a/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.H
+++ b/src/lagrangian/dsmc/submodels/WallInteractionModel/SpecularReflection/SpecularReflection.H
@@ -76,7 +76,8 @@ public:
const wallPolyPatch& wpp,
const label faceId,
vector& U,
- scalar mass
+ scalar& Ei,
+ label typeId
);
};
diff --git a/src/lagrangian/dsmc/submodels/WallInteractionModel/WallInteractionModel/WallInteractionModel.H b/src/lagrangian/dsmc/submodels/WallInteractionModel/WallInteractionModel/WallInteractionModel.H
index d8667201a7a..582bf2b1065 100644
--- a/src/lagrangian/dsmc/submodels/WallInteractionModel/WallInteractionModel/WallInteractionModel.H
+++ b/src/lagrangian/dsmc/submodels/WallInteractionModel/WallInteractionModel/WallInteractionModel.H
@@ -130,7 +130,8 @@ public:
const wallPolyPatch& wpp,
const label faceId,
vector& U,
- scalar mass
+ scalar& Ei,
+ label typeId
) = 0;
};
--
GitLab