Measurement of surface velocity. Created momentum and rhoM field inside the

cloud - duplicated memory with solver for the moment.
This commit is contained in:
graham
2009-07-17 15:40:11 +01:00
parent 9321f7e1e5
commit 103cd27ef5
7 changed files with 225 additions and 131 deletions

View File

@ -490,6 +490,23 @@ void Foam::DsmcCloud<ParcelType>::resetSurfaceDataFields()
{
fDBF[p] = vector::zero;
}
volVectorField::GeometricBoundaryField& momentumBF
(
momentum_.boundaryField()
);
forAll(momentumBF, p)
{
momentumBF[p] = vector::zero;
}
volScalarField::GeometricBoundaryField& rhoMBF(rhoM_.boundaryField());
forAll(rhoMBF, p)
{
rhoMBF[p] = VSMALL;
}
}
@ -591,6 +608,37 @@ Foam::DsmcCloud<ParcelType>::DsmcCloud
vector::zero
)
),
momentum_
(
IOobject
(
this->name() + "momentum_",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedVector
(
"zero",
dimensionSet(1, -2, -1, 0, 0),
vector::zero
)
),
rhoM_
(
IOobject
(
this->name() + "rhoM_",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedScalar("zero", dimensionSet(1, -3, 0, 0, 0), 0.0)
),
constProps_(),
rndGen_(label(149382906) + 7183*Pstream::myProcNo()),
T_(T),
@ -703,6 +751,37 @@ Foam::DsmcCloud<ParcelType>::DsmcCloud
vector::zero
)
),
momentum_
(
IOobject
(
this->name() + "momentum_",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedVector
(
"zero",
dimensionSet(1, -2, -1, 0, 0),
vector::zero
)
),
rhoM_
(
IOobject
(
this->name() + "rhoM_",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
mesh_,
dimensionedScalar("zero", dimensionSet(1, -3, 0, 0, 0), 0.0)
),
constProps_(),
rndGen_(label(971501) + 1526*Pstream::myProcNo()),
T_

View File

@ -110,6 +110,12 @@ class DsmcCloud
//- Force density at surface field
volVectorField fD_;
//- Momentum density field
volVectorField momentum_;
//- Mass density field
volScalarField rhoM_;
//- Parcel constant properties - one for each type
List<typename ParcelType::constantProperties> constProps_;
@ -254,19 +260,19 @@ public:
inline scalar cachedDeltaT() const;
// References to the surface data collection fields
// References to the boundary fields for surface data collection
//- Return heat flux at surface field
inline const volScalarField& q() const;
//- Return non-const heat flux boundary field reference
inline volScalarField::GeometricBoundaryField& qBF();
//- Return non-const heat flux at surface field
inline volScalarField& q();
//- Return non-const force density at boundary field reference
inline volVectorField::GeometricBoundaryField& fDBF();
//- Return force density at surface field
inline const volVectorField& fD() const;
//- Return non-const momentum density boundary field reference
inline volVectorField::GeometricBoundaryField& momentumBF();
//- Return non-const force density at surface field
inline volVectorField& fD();
//- Return non-const mass density boundary field reference
inline volScalarField::GeometricBoundaryField& rhoMBF();
// References to the macroscopic fields
@ -391,17 +397,20 @@ public:
// Fields
//- Return heat flux at surface field
inline const volScalarField& q() const;
//- Return force density at surface field
inline const volVectorField& fD() const;
//- Return the real particle number density field
inline const tmp<volScalarField> rhoN() const;
//- Return the particle mass density field
inline const tmp<volScalarField> rhoM() const;
//- Return the field of number of DSMC particles
inline const tmp<volScalarField> dsmcRhoN() const;
inline const volScalarField& rhoM();
//- Return the momentum density field
inline const tmp<volVectorField> momentum() const;
inline const volVectorField& momentum();
//- Return the total linear kinetic energy (translational and
// thermal density field
@ -413,6 +422,9 @@ public:
//- Return the average internal degrees of freedom field
inline const tmp<volScalarField> iDof() const;
//- Return the field of number of DSMC particles
inline const tmp<volScalarField> dsmcRhoN() const;
// Cloud evolution functions

View File

@ -135,30 +135,34 @@ inline Foam::scalar Foam::DsmcCloud<ParcelType>::cachedDeltaT() const
template<class ParcelType>
inline const Foam::volScalarField& Foam::DsmcCloud<ParcelType>::q() const
inline Foam::volScalarField::GeometricBoundaryField&
Foam::DsmcCloud<ParcelType>::qBF()
{
return q_;
return q_.boundaryField();
}
template<class ParcelType>
inline Foam::volScalarField& Foam::DsmcCloud<ParcelType>::q()
inline Foam::volVectorField::GeometricBoundaryField&
Foam::DsmcCloud<ParcelType>::fDBF()
{
return q_;
return fD_.boundaryField();
}
template<class ParcelType>
inline const Foam::volVectorField& Foam::DsmcCloud<ParcelType>::fD() const
inline Foam::volVectorField::GeometricBoundaryField&
Foam::DsmcCloud<ParcelType>::momentumBF()
{
return fD_;
return momentum_.boundaryField();
}
template<class ParcelType>
inline Foam::volVectorField& Foam::DsmcCloud<ParcelType>::fD()
inline Foam::volScalarField::GeometricBoundaryField&
Foam::DsmcCloud<ParcelType>::rhoMBF()
{
return fD_;
return rhoM_.boundaryField();
}
@ -373,6 +377,20 @@ Foam::DsmcCloud<ParcelType>::maxwellianMostProbableSpeed
}
template<class ParcelType>
inline const Foam::volScalarField& Foam::DsmcCloud<ParcelType>::q() const
{
return q_;
}
template<class ParcelType>
inline const Foam::volVectorField& Foam::DsmcCloud<ParcelType>::fD() const
{
return fD_;
}
template<class ParcelType>
inline const Foam::tmp<Foam::volScalarField>
Foam::DsmcCloud<ParcelType>::rhoN() const
@ -411,116 +429,44 @@ Foam::DsmcCloud<ParcelType>::rhoN() const
template<class ParcelType>
inline const Foam::tmp<Foam::volScalarField>
Foam::DsmcCloud<ParcelType>::rhoM() const
inline const Foam::volScalarField& Foam::DsmcCloud<ParcelType>::rhoM()
{
tmp<volScalarField> trhoM
(
new volScalarField
(
IOobject
(
this->name() + "rhoM",
this->db().time().timeName(),
this->db(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh_,
dimensionedScalar("zero", dimensionSet(1, -3, 0, 0, 0), VSMALL)
)
);
scalarField& rM = rhoM_.internalField();
rM = VSMALL;
scalarField& rhoM = trhoM().internalField();
forAllConstIter(typename DsmcCloud<ParcelType>, *this, iter)
{
const ParcelType& p = iter();
const label cellI = p.cell();
rhoM[cellI] += constProps(p.typeId()).mass();
rM[cellI] += constProps(p.typeId()).mass();
}
rhoM *= nParticle_/mesh().cellVolumes();
rM *= nParticle_/mesh().cellVolumes();
return trhoM;
return rhoM_;
}
template<class ParcelType>
inline const Foam::tmp<Foam::volScalarField>
Foam::DsmcCloud<ParcelType>::dsmcRhoN() const
inline const Foam::volVectorField& Foam::DsmcCloud<ParcelType>::momentum()
{
tmp<volScalarField> tdsmcRhoN
(
new volScalarField
(
IOobject
(
this->name() + "dsmcRhoN",
this->db().time().timeName(),
this->db(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh_,
dimensionedScalar("zero", dimensionSet(0, -3, 0, 0, 0), 0.0)
)
);
vectorField& mom = momentum_.internalField();
mom = vector::zero;
scalarField& dsmcRhoN = tdsmcRhoN().internalField();
forAllConstIter(typename DsmcCloud<ParcelType>, *this, iter)
{
const ParcelType& p = iter();
const label cellI = p.cell();
dsmcRhoN[cellI]++;
mom[cellI] += constProps(p.typeId()).mass()*p.U();
}
return tdsmcRhoN;
}
mom *= nParticle_/mesh().cellVolumes();
template<class ParcelType>
inline const Foam::tmp<Foam::volVectorField>
Foam::DsmcCloud<ParcelType>::momentum() const
{
tmp<volVectorField> tmomentum
(
new volVectorField
(
IOobject
(
this->name() + "momentum",
this->db().time().timeName(),
this->db(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh_,
dimensionedVector
(
"zero",
dimensionSet(1, -2, -1, 0, 0),
vector::zero
)
)
);
vectorField& momentum = tmomentum().internalField();
forAllConstIter(typename DsmcCloud<ParcelType>, *this, iter)
{
const ParcelType& p = iter();
const label cellI = p.cell();
momentum[cellI] += constProps(p.typeId()).mass()*p.U();
}
momentum *= nParticle_/mesh().cellVolumes();
return tmomentum;
return momentum_;
}
@ -636,6 +582,41 @@ Foam::DsmcCloud<ParcelType>::iDof() const
}
template<class ParcelType>
inline const Foam::tmp<Foam::volScalarField>
Foam::DsmcCloud<ParcelType>::dsmcRhoN() const
{
tmp<volScalarField> tdsmcRhoN
(
new volScalarField
(
IOobject
(
this->name() + "dsmcRhoN",
this->db().time().timeName(),
this->db(),
IOobject::NO_READ,
IOobject::NO_WRITE,
false
),
mesh_,
dimensionedScalar("zero", dimensionSet(0, -3, 0, 0, 0), 0.0)
)
);
scalarField& dsmcRhoN = tdsmcRhoN().internalField();
forAllConstIter(typename DsmcCloud<ParcelType>, *this, iter)
{
const ParcelType& p = iter();
const label cellI = p.cell();
dsmcRhoN[cellI]++;
}
return tdsmcRhoN;
}
template<class ParcelType>
inline void Foam::DsmcCloud<ParcelType>::clear()
{

View File

@ -114,10 +114,30 @@ void Foam::DsmcParcel<ParcelType>::hitWallPatch
TrackData& td
)
{
label wppIndex = wpp.index();
label wppLocalFace = wpp.whichFace(this->face());
const scalar fA = mag(wpp.faceAreas()[wppLocalFace]);
const scalar deltaT = td.cloud().cachedDeltaT();
const constantProperties& constProps(td.cloud().constProps(typeId_));
scalar m = constProps.mass();
vector nw = wpp.faceAreas()[wppLocalFace];
nw /= mag(nw);
scalar U_dot_nw = U_ & nw;
vector Ut = U_ - U_dot_nw*nw;
td.cloud().momentumBF()[wppIndex][wppLocalFace] +=
m*Ut/max(mag(U_dot_nw)*fA, VSMALL);
td.cloud().rhoMBF()[wppIndex][wppLocalFace] +=
m/max(mag(U_dot_nw)*fA, VSMALL);
// pre-interaction energy
scalar preIE = 0.5*m*(U_ & U_) + Ei_;
@ -133,27 +153,29 @@ void Foam::DsmcParcel<ParcelType>::hitWallPatch
typeId_
);
U_dot_nw = U_ & nw;
Ut = U_ - U_dot_nw*nw;
td.cloud().momentumBF()[wppIndex][wppLocalFace] +=
m*Ut/max(mag(U_dot_nw)*fA, VSMALL);
td.cloud().rhoMBF()[wppIndex][wppLocalFace] +=
m/max(mag(U_dot_nw)*fA, VSMALL);
// post-interaction energy
scalar postIE = 0.5*m*(U_ & U_) + Ei_;
// post-interaction momentum
vector postIMom = m*U_;
label wppIndex = wpp.index();
label wppLocalFace = wpp.whichFace(this->face());
const scalar fA = mag(wpp.faceAreas()[wppLocalFace]);
const scalar deltaT = td.cloud().cachedDeltaT();
scalar deltaQ = td.cloud().nParticle()*(preIE - postIE)/(deltaT*fA);
vector deltaFD = td.cloud().nParticle()*(preIMom - postIMom)/(deltaT*fA);
td.cloud().q().boundaryField()[wppIndex][wppLocalFace] += deltaQ;
td.cloud().qBF()[wppIndex][wppLocalFace] += deltaQ;
td.cloud().fDBF()[wppIndex][wppLocalFace] += deltaFD;
td.cloud().fD().boundaryField()[wppIndex][wppLocalFace] += deltaFD;
}

View File

@ -68,10 +68,10 @@ void Foam::MaxwellianThermal<CloudType>::correct
nw /= mag(nw);
// Normal velocity magnitude
scalar magUn = U & nw;
scalar U_dot_nw = U & nw;
// Wall tangential velocity (flow direction)
vector Ut = U - magUn*nw;
vector Ut = U - U_dot_nw*nw;
CloudType& cloud(this->owner());
@ -90,9 +90,9 @@ void Foam::MaxwellianThermal<CloudType>::correct
U.z()*(0.8 + 0.2*rndGen.scalar01())
);
magUn = U & nw;
U_dot_nw = U & nw;
Ut = U - magUn*nw;
Ut = U - U_dot_nw*nw;
}
// Wall tangential unit vector

View File

@ -69,7 +69,7 @@ void Foam::MixedDiffuseSpecular<CloudType>::correct
nw /= mag(nw);
// Normal velocity magnitude
scalar magUn = U & nw;
scalar U_dot_nw = U & nw;
CloudType& cloud(this->owner());
@ -80,7 +80,7 @@ void Foam::MixedDiffuseSpecular<CloudType>::correct
// Diffuse reflection
// Wall tangential velocity (flow direction)
vector Ut = U - magUn*nw;
vector Ut = U - U_dot_nw*nw;
while (mag(Ut) < SMALL)
{
@ -95,9 +95,9 @@ void Foam::MixedDiffuseSpecular<CloudType>::correct
U.z()*(0.8 + 0.2*rndGen.scalar01())
);
magUn = U & nw;
U_dot_nw = U & nw;
Ut = U - magUn*nw;
Ut = U - U_dot_nw*nw;
}
// Wall tangential unit vector
@ -128,9 +128,9 @@ void Foam::MixedDiffuseSpecular<CloudType>::correct
{
// Specular reflection
if (magUn > 0.0)
if (U_dot_nw > 0.0)
{
U -= 2.0*magUn*nw;
U -= 2.0*U_dot_nw*nw;
}
}

View File

@ -63,11 +63,11 @@ void Foam::SpecularReflection<CloudType>::correct
vector nw = wpp.faceAreas()[wpp.whichFace(faceId)];
nw /= mag(nw);
scalar magUn = U & nw;
scalar U_dot_nw = U & nw;
if (magUn > 0.0)
if (U_dot_nw > 0.0)
{
U -= 2.0*magUn*nw;
U -= 2.0*U_dot_nw*nw;
}
}