/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see .
\*---------------------------------------------------------------------------*/
#include "MRFZone.H"
#include "fvMesh.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "fvMatrices.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template
void Foam::MRFZone::makeRelativeRhoFlux
(
const RhoFieldType& rho,
surfaceScalarField& phi
) const
{
if (!active_)
{
return;
}
const surfaceVectorField& Cf = mesh_.Cf();
const surfaceVectorField& Sf = mesh_.Sf();
const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
const vectorField& Cfi = Cf;
const vectorField& Sfi = Sf;
scalarField& phii = phi.primitiveFieldRef();
// Internal faces
forAll(internalFaces_, i)
{
label facei = internalFaces_[i];
phii[facei] -= rho[facei]*(Omega ^ (Cfi[facei] - origin_)) & Sfi[facei];
}
makeRelativeRhoFlux(rho.boundaryField(), phi.boundaryFieldRef());
}
template
void Foam::MRFZone::makeRelativeRhoFlux
(
const RhoFieldType& rho,
FieldField& phi
) const
{
if (!active_)
{
return;
}
const surfaceVectorField& Cf = mesh_.Cf();
const surfaceVectorField& Sf = mesh_.Sf();
const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
// Included patches
forAll(includedFaces_, patchi)
{
forAll(includedFaces_[patchi], i)
{
label patchFacei = includedFaces_[patchi][i];
phi[patchi][patchFacei] = 0.0;
}
}
// Excluded patches
forAll(excludedFaces_, patchi)
{
forAll(excludedFaces_[patchi], i)
{
label patchFacei = excludedFaces_[patchi][i];
phi[patchi][patchFacei] -=
rho[patchi][patchFacei]
* (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
& Sf.boundaryField()[patchi][patchFacei];
}
}
}
template
void Foam::MRFZone::makeRelativeRhoFlux
(
const RhoFieldType& rho,
Field& phi,
const label patchi
) const
{
if (!active_)
{
return;
}
const surfaceVectorField& Cf = mesh_.Cf();
const surfaceVectorField& Sf = mesh_.Sf();
const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
// Included patches
forAll(includedFaces_[patchi], i)
{
label patchFacei = includedFaces_[patchi][i];
phi[patchFacei] = 0.0;
}
// Excluded patches
forAll(excludedFaces_[patchi], i)
{
label patchFacei = excludedFaces_[patchi][i];
phi[patchFacei] -=
rho[patchFacei]
* (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
& Sf.boundaryField()[patchi][patchFacei];
}
}
template
void Foam::MRFZone::makeAbsoluteRhoFlux
(
const RhoFieldType& rho,
surfaceScalarField& phi
) const
{
if (!active_)
{
return;
}
const surfaceVectorField& Cf = mesh_.Cf();
const surfaceVectorField& Sf = mesh_.Sf();
const vector Omega = omega_->value(mesh_.time().timeOutputValue())*axis_;
const vectorField& Cfi = Cf;
const vectorField& Sfi = Sf;
scalarField& phii = phi.primitiveFieldRef();
// Internal faces
forAll(internalFaces_, i)
{
label facei = internalFaces_[i];
phii[facei] += rho[facei]*(Omega ^ (Cfi[facei] - origin_)) & Sfi[facei];
}
surfaceScalarField::Boundary& phibf = phi.boundaryFieldRef();
// Included patches
forAll(includedFaces_, patchi)
{
forAll(includedFaces_[patchi], i)
{
label patchFacei = includedFaces_[patchi][i];
phibf[patchi][patchFacei] +=
rho.boundaryField()[patchi][patchFacei]
* (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
& Sf.boundaryField()[patchi][patchFacei];
}
}
// Excluded patches
forAll(excludedFaces_, patchi)
{
forAll(excludedFaces_[patchi], i)
{
label patchFacei = excludedFaces_[patchi][i];
phibf[patchi][patchFacei] +=
rho.boundaryField()[patchi][patchFacei]
* (Omega ^ (Cf.boundaryField()[patchi][patchFacei] - origin_))
& Sf.boundaryField()[patchi][patchFacei];
}
}
}
// ************************************************************************* //