Bug fix: laminar forces

This commit is contained in:
Hrvoje Jasak 2015-05-13 16:30:45 +01:00
parent 9772e43abf
commit 8628a03fe0

View file

@ -154,12 +154,6 @@ Foam::immersedBoundaryForces::calcForcesMoment() const
const fvMesh& mesh = U.mesh(); const fvMesh& mesh = U.mesh();
// Stress from turbulence model. Note: this does not account
// for wall functions on the Immersed Boundary: use wallShearStress
// from Immersed Boundary U wall patch instead.
// HJ, 11/Aug/2014
// volSymmTensorField stress = devRhoReff();
// Scale pRef by density for incompressible simulations // Scale pRef by density for incompressible simulations
scalar pRef = pRef_/rho(p); scalar pRef = pRef_/rho(p);
@ -170,9 +164,9 @@ Foam::immersedBoundaryForces::calcForcesMoment() const
// Check and cast into immersed boundary type // Check and cast into immersed boundary type
if if
( (
isA<immersedBoundaryFvPatchVectorField> isA<immersedBoundaryFvPatch>
( (
U.boundaryField()[patchI] mesh.boundary()[patchI]
) )
) )
{ {
@ -184,25 +178,6 @@ Foam::immersedBoundaryForces::calcForcesMoment() const
mesh.boundary()[patchI] mesh.boundary()[patchI]
); );
// Get immersed boundary velocity. Used to access wall
// shear stress
const immersedBoundaryVelocityWallFunctionFvPatchVectorField&
UPatch =
refCast
<
const
immersedBoundaryVelocityWallFunctionFvPatchVectorField
>
(
U.boundaryField()[patchI]
);
// const immersedBoundaryFvPatchSymmTensorField stressPatch =
// refCast<const immersedBoundaryFvPatchSymmTensorField>
// (
// stress.boundaryField()[patchI]
// );
const immersedBoundaryFvPatchScalarField pPatch = const immersedBoundaryFvPatchScalarField pPatch =
refCast<const immersedBoundaryFvPatchScalarField> refCast<const immersedBoundaryFvPatchScalarField>
( (
@ -223,19 +198,82 @@ Foam::immersedBoundaryForces::calcForcesMoment() const
fm.first().first() += rho(p)*sum(pf); fm.first().first() += rho(p)*sum(pf);
fm.second().first() += rho(p)*sum(Md ^ pf); fm.second().first() += rho(p)*sum(Md ^ pf);
// Old treatment: if
// Shear force is obtained from velocity wall functions (
// and integrated on triangular faces isA<immersedBoundaryVelocityWallFunctionFvPatchVectorField>
vectorField vf = (
sA*ibPatch.toTriFaces(UPatch.wallShearStress()); U.boundaryField()[patchI]
)
)
{
// Turbulent wall functions
// New treatment: normal force calculated on triangles // Get immersed boundary velocity. Used to access wall
// Damir Rigler, 30/Apr/2014 // shear stress
// vectorField vfOld = (Sfb & stressPatch.triValue()); const
// Info<< "Force difference " << sumMag(vf - vfOld) << endl; immersedBoundaryVelocityWallFunctionFvPatchVectorField&
UPatch =
refCast
<
const
immersedBoundaryVelocityWallFunctionFvPatchVectorField
>
(
U.boundaryField()[patchI]
);
fm.first().second() += sum(vf); // Shear force is obtained from velocity wall functions
fm.second().second() += sum(Md ^ vf); // and integrated on triangular faces
vectorField vf =
sA*ibPatch.toTriFaces(UPatch.wallShearStress());
fm.first().second() += sum(vf);
fm.second().second() += sum(Md ^ vf);
}
else if
(
isA<immersedBoundaryFvPatchVectorField>
(
U.boundaryField()[patchI]
)
)
{
// Laminar flow
// Get immersed boundary velocity
const immersedBoundaryFvPatchVectorField& UPatch =
refCast<const immersedBoundaryFvPatchVectorField>
(
U.boundaryField()[patchI]
);
// Look up the viscosity
if (mesh.foundObject<dictionary>("transportProperties"))
{
const dictionary& transportProperties =
mesh.lookupObject<dictionary>
(
"transportProperties"
);
dimensionedScalar nu(transportProperties.lookup("nu"));
vectorField vf =
sA*rho(p)*nu.value()*UPatch.triGrad();
fm.first().second() += sum(vf);
fm.second().second() += sum(Md ^ vf);
}
else
{
InfoIn
(
"immersedBoundaryForces::forcesMoments"
"immersedBoundaryForces::calcForcesMoment() const"
) << "Laminar flow, but cannot find nu. Skipping"
<< endl;
}
}
} }
} }
} }