Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev

This commit is contained in:
mattijs
2013-09-02 11:12:34 +01:00
20 changed files with 473 additions and 113 deletions

View File

@ -1,5 +1,4 @@
EXE_INC = \
-I../rhoPorousMRFPimpleFoam \
-I.. \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \

View File

@ -10,5 +10,4 @@ EXE_INC = \
LIB_LIBS = \
-linterfaceProperties \
-lincompressibleTransportModels \
-lcompressibleMultiphaseEulerianInterfacialModels \
-lfiniteVolume

View File

@ -64,6 +64,7 @@ initialPointsMethod/bodyCentredCubic/bodyCentredCubic.C
initialPointsMethod/faceCentredCubic/faceCentredCubic.C
initialPointsMethod/pointFile/pointFile.C
initialPointsMethod/autoDensity/autoDensity.C
initialPointsMethod/rayShooting/rayShooting.C
relaxationModel/relaxationModel/relaxationModel.C
relaxationModel/adaptiveLinear/adaptiveLinear.C

View File

@ -1249,7 +1249,7 @@ void Foam::conformalVoronoiMesh::move()
if
(
(
(vA->internalPoint() && vB->internalPoint())
(vA->internalPoint() || vB->internalPoint())
&& (!vA->referred() || !vB->referred())
// ||
// (

View File

@ -670,6 +670,32 @@ Foam::Field<bool> Foam::conformationSurfaces::wellInOutSide
continue;
}
const searchableSurface& surface(allGeometry_[surfaces_[s]]);
if (!surface.hasVolumeType())
{
pointField sample(1, samplePts[i]);
scalarField nearestDistSqr(1, GREAT);
List<pointIndexHit> info;
surface.findNearest(sample, nearestDistSqr, info);
vector hitDir = info[0].rawPoint() - samplePts[i];
hitDir /= mag(hitDir) + SMALL;
if
(
findSurfaceAnyIntersection
(
samplePts[i],
info[0].rawPoint() - 1e-3*mag(hitDir)*(hitDir)
)
)
{
continue;
}
}
if (surfaceVolumeTests[s][i] == volumeType::OUTSIDE)
{
if
@ -694,44 +720,6 @@ Foam::Field<bool> Foam::conformationSurfaces::wellInOutSide
break;
}
}
// else
// {
// // Surface volume type is unknown
// Info<< "UNKNOWN" << endl;
// // Get nearest face normal
//
// pointField sample(1, samplePts[i]);
// scalarField nearestDistSqr(1, GREAT);
// List<pointIndexHit> info;
// vectorField norms(1);
//
// surface.findNearest(sample, nearestDistSqr, info);
// surface.getNormal(info, norms);
//
// vector fN = norms[0];
// fN /= mag(fN);
//
// vector hitDir = info[0].rawPoint() - samplePts[i];
// hitDir /= mag(hitDir);
//
// if ((fN & hitDir) < 0)
// {
// // Point is OUTSIDE
//
// if
// (
// normalVolumeTypes_[regionI]
// == extendedFeatureEdgeMesh::OUTSIDE
// )
// {
// }
// else
// {
// insidePoint[i] = false;
// break;
// }
// }
// }
}
}

View File

@ -0,0 +1,218 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013 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 <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "rayShooting.H"
#include "addToRunTimeSelectionTable.H"
#include "triSurfaceMesh.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(rayShooting, 0);
addToRunTimeSelectionTable(initialPointsMethod, rayShooting, dictionary);
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
rayShooting::rayShooting
(
const dictionary& initialPointsDict,
const Time& runTime,
Random& rndGen,
const conformationSurfaces& geometryToConformTo,
const cellShapeControl& cellShapeControls,
const autoPtr<backgroundMeshDecomposition>& decomposition
)
:
initialPointsMethod
(
typeName,
initialPointsDict,
runTime,
rndGen,
geometryToConformTo,
cellShapeControls,
decomposition
),
randomiseInitialGrid_(detailsDict().lookup("randomiseInitialGrid")),
randomPerturbationCoeff_
(
readScalar(detailsDict().lookup("randomPerturbationCoeff"))
)
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
List<Vb::Point> rayShooting::initialPoints() const
{
// Loop over surface faces
const searchableSurfaces& surfaces = geometryToConformTo().geometry();
const labelList& surfacesToConformTo = geometryToConformTo().surfaces();
const scalar maxRayLength = surfaces.bounds().mag();
// Initialise points list
label initialPointsSize = 0;
forAll(surfaces, surfI)
{
initialPointsSize += surfaces[surfI].size();
}
DynamicList<Vb::Point> initialPoints(initialPointsSize);
forAll(surfacesToConformTo, surfI)
{
const searchableSurface& s = surfaces[surfacesToConformTo[surfI]];
tmp<pointField> faceCentresTmp(s.coordinates());
const pointField& faceCentres = faceCentresTmp();
Info<< " Shoot rays from " << s.name() << nl
<< " nRays = " << faceCentres.size() << endl;
forAll(faceCentres, fcI)
{
const Foam::point& fC = faceCentres[fcI];
if
(
Pstream::parRun()
&& !decomposition().positionOnThisProcessor(fC)
)
{
continue;
}
const scalar pert =
randomPerturbationCoeff_
*cellShapeControls().cellSize(fC);
pointIndexHit surfHitStart;
label hitSurfaceStart;
// Face centres should be on the surface so search distance can be
// small
geometryToConformTo().findSurfaceNearest
(
fC,
sqr(pert),
surfHitStart,
hitSurfaceStart
);
vectorField normStart(1, vector::min);
geometryToConformTo().getNormal
(
hitSurfaceStart,
List<pointIndexHit>(1, surfHitStart),
normStart
);
pointIndexHit surfHitEnd;
label hitSurfaceEnd;
geometryToConformTo().findSurfaceNearestIntersection
(
fC - normStart[0]*pert,
fC - normStart[0]*maxRayLength,
surfHitEnd,
hitSurfaceEnd
);
if (surfHitEnd.hit())
{
vectorField normEnd(1, vector::min);
geometryToConformTo().getNormal
(
hitSurfaceEnd,
List<pointIndexHit>(1, surfHitEnd),
normEnd
);
if ((normStart[0] & normEnd[0]) < 0)
{
line<point, point> l(fC, surfHitEnd.hitPoint());
if (Pstream::parRun())
{
// Clip the line in parallel
pointIndexHit procIntersection =
decomposition().findLine
(
l.start(),
l.end()
);
if (procIntersection.hit())
{
l =
line<point, point>
(
l.start(),
procIntersection.hitPoint()
);
}
}
Foam::point midPoint(l.centre());
const scalar minDistFromSurfaceSqr =
minimumSurfaceDistanceCoeffSqr_
*sqr(cellShapeControls().cellSize(midPoint));
if (randomiseInitialGrid_)
{
midPoint.x() += pert*(rndGen().scalar01() - 0.5);
midPoint.y() += pert*(rndGen().scalar01() - 0.5);
midPoint.z() += pert*(rndGen().scalar01() - 0.5);
}
if
(
magSqr(midPoint - l.start()) > minDistFromSurfaceSqr
&& magSqr(midPoint - l.end()) > minDistFromSurfaceSqr
)
{
initialPoints.append(toPoint(midPoint));
}
}
}
}
}
return initialPoints.shrink();
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// ************************************************************************* //

View File

@ -0,0 +1,103 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2013 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 <http://www.gnu.org/licenses/>.
Class
Foam::rayShooting
Description
SourceFiles
rayShooting.C
\*---------------------------------------------------------------------------*/
#ifndef rayShooting_H
#define rayShooting_H
#include "initialPointsMethod.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class rayShooting Declaration
\*---------------------------------------------------------------------------*/
class rayShooting
:
public initialPointsMethod
{
private:
// Private data
//- Should the initial positions be randomised
Switch randomiseInitialGrid_;
//- Randomise the initial positions by fraction of the initialCellSize_
scalar randomPerturbationCoeff_;
public:
//- Runtime type information
TypeName("rayShooting");
// Constructors
//- Construct from components
rayShooting
(
const dictionary& initialPointsDict,
const Time& runTime,
Random& rndGen,
const conformationSurfaces& geometryToConformTo,
const cellShapeControl& cellShapeControls,
const autoPtr<backgroundMeshDecomposition>& decomposition
);
//- Destructor
virtual ~rayShooting()
{}
// Member Functions
//- Return the initial points for the conformalVoronoiMesh
virtual List<Vb::Point> initialPoints() const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -61,6 +61,9 @@ void reactingOneDim::readReactingOneDimControls()
coeffs().lookup("gasHSource") >> gasHSource_;
coeffs().lookup("QrHSource") >> QrHSource_;
useChemistrySolvers_ =
coeffs().lookupOrDefault<bool>("useChemistrySolvers", true);
}
@ -321,6 +324,8 @@ void reactingOneDim::solveEnergy()
(
fvm::ddt(rho_, h_)
- fvm::laplacian(alpha, h_)
+ fvc::laplacian(alpha, h_)
- fvc::laplacian(kappa(), T())
==
chemistrySh_
- fvm::Sp(solidChemistry_->RRg(), h_)
@ -462,7 +467,8 @@ reactingOneDim::reactingOneDim(const word& modelType, const fvMesh& mesh)
totalGasMassFlux_(0.0),
totalHeatRR_(dimensionedScalar("zero", dimEnergy/dimTime, 0.0)),
gasHSource_(false),
QrHSource_(false)
QrHSource_(false),
useChemistrySolvers_(true)
{
if (active_)
{
@ -560,7 +566,8 @@ reactingOneDim::reactingOneDim
totalGasMassFlux_(0.0),
totalHeatRR_(dimensionedScalar("zero", dimEnergy/dimTime, 0.0)),
gasHSource_(false),
QrHSource_(false)
QrHSource_(false),
useChemistrySolvers_(true)
{
if (active_)
{
@ -681,11 +688,18 @@ void reactingOneDim::evolveRegion()
{
Info<< "\nEvolving pyrolysis in region: " << regionMesh().name() << endl;
solidChemistry_->solve
(
time().value() - time().deltaTValue(),
time().deltaTValue()
);
if (useChemistrySolvers_)
{
solidChemistry_->solve
(
time().value() - time().deltaTValue(),
time().deltaTValue()
);
}
else
{
solidChemistry_->calculate();
}
solveContinuity();

View File

@ -154,6 +154,9 @@ protected:
//- Add in depth radiation source term
bool QrHSource_;
//- Use chemistry solvers (ode or sequential)
bool useChemistrySolvers_;
// Protected member functions

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -140,6 +140,12 @@ public:
const label i
) const = 0;
//- Return access to chemical source terms [kg/m3/s]
virtual DimensionedField<scalar, volMesh>& RR
(
const label i
) = 0;
// Chemistry solution

View File

@ -211,6 +211,12 @@ public:
const label i
) const;
//- Return non const access to chemical source terms [kg/m3/s]
virtual DimensionedField<scalar, volMesh>& RR
(
const label i
);
//- Solve the reaction system for the given start time and time
// step and return the characteristic time
virtual scalar solve(const scalar t0, const scalar deltaT);

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -78,5 +78,15 @@ Foam::chemistryModel<CompType, ThermoType>::RR
return RR_[i];
}
template<class CompType, class ThermoType>
Foam::DimensionedField<Foam::scalar, Foam::volMesh>&
Foam::chemistryModel<CompType, ThermoType>::RR
(
const label i
)
{
return RR_[i];
}
// ************************************************************************* //

View File

@ -53,4 +53,35 @@ Foam::basicSolidChemistryModel::~basicSolidChemistryModel()
{}
const Foam::DimensionedField<Foam::scalar, Foam::volMesh>&
Foam::basicSolidChemistryModel::RR(const label i) const
{
notImplemented
(
"const Foam::DimensionedField<Foam::scalar, Foam::volMesh>&"
"basicSolidChemistryModel::RR(const label)"
);
return (DimensionedField<scalar, volMesh>::null());
}
Foam::DimensionedField<Foam::scalar, Foam::volMesh>&
Foam::basicSolidChemistryModel::RR(const label i)
{
notImplemented
(
"Foam::DimensionedField<Foam::scalar, Foam::volMesh>&"
"basicSolidChemistryModel::RR(const label)"
);
return dynamic_cast<DimensionedField<scalar, volMesh>&>
(
const_cast<DimensionedField<scalar, volMesh>& >
(
DimensionedField<scalar, volMesh>::null()
)
);
}
// ************************************************************************* //

View File

@ -149,6 +149,15 @@ public:
//- Calculates the reaction rates
virtual void calculate() = 0;
//- Return const access to the total source terms
virtual const DimensionedField<scalar, volMesh>& RR
(
const label i
) const;
//- Return non-const access to the total source terms
virtual DimensionedField<scalar, volMesh>& RR(const label i);
};

View File

@ -535,14 +535,21 @@ updateConcsInReactionI
c[si] = max(0.0, c[si]);
}
scalar sr = 0.0;
forAll(R.rhs(), s)
{
label si = R.rhs()[s].index;
const scalar rhoR = this->solidThermo_[si].rho(p, T);
const scalar sr = rhoR/rhoL;
sr = rhoR/rhoL;
c[si] += dt*sr*omeg;
c[si] = max(0.0, c[si]);
}
forAll(R.grhs(), g)
{
label gi = R.grhs()[g].index;
c[gi + this->nSolids_] += (1.0 - sr)*omeg*dt;
}
}
@ -561,24 +568,11 @@ updateRRInReactionI
simpleMatrix<scalar>& RR
) const
{
const Reaction<SolidThermo>& R = this->reactions_[index];
scalar rhoL = 0.0;
forAll(R.lhs(), s)
{
label si = R.lhs()[s].index;
rhoL = this->solidThermo_[si].rho(p, T);
RR[si][rRef] -= pr*corr;
RR[si][lRef] += pf*corr;
}
forAll(R.rhs(), s)
{
label si = R.rhs()[s].index;
const scalar rhoR = this->solidThermo_[si].rho(p, T);
const scalar sr = rhoR/rhoL;
RR[si][lRef] -= sr*pf*corr;
RR[si][rRef] += sr*pr*corr;
}
notImplemented
(
"void Foam::pyrolysisChemistryModel<CompType, SolidThermo,GasThermo>::"
"updateRRInReactionI"
);
}
@ -666,7 +660,6 @@ Foam::pyrolysisChemistryModel<CompType, SolidThermo, GasThermo>::solve
scalar newCp = 0.0;
scalar newhi = 0.0;
scalar invRho = 0.0;
scalarList dcdt = (c - c0)/dt;
for (label i=0; i<this->nSolids_; i++)
@ -675,7 +668,6 @@ Foam::pyrolysisChemistryModel<CompType, SolidThermo, GasThermo>::solve
scalar Yi = c[i]/cTot;
newCp += Yi*this->solidThermo_[i].Cp(pi, Ti);
newhi -= dYi*this->solidThermo_[i].Hc();
invRho += Yi/this->solidThermo_[i].rho(pi, Ti);
}
scalar dTi = (newhi/newCp)*dt;

View File

@ -137,8 +137,9 @@ public:
const bool updateC0 = false
) const;
//- Return the reaction rate for reaction r and the reference
// species and charateristic times
//- Return the reaction rate for reaction r
// NOTE: Currently does not calculate reference specie and
// characteristic times (pf, cf,l Ref, etc.)
virtual scalar omega
(
const Reaction<SolidThermo>& r,
@ -153,8 +154,10 @@ public:
label& rRef
) const;
//- Return the reaction rate for iReaction and the reference
// species and charateristic times
//- Return the reaction rate for iReaction
// NOTE: Currently does not calculate reference specie and
// characteristic times (pf, cf,l Ref, etc.)
virtual scalar omegaI
(
label iReaction,
@ -169,6 +172,7 @@ public:
label& rRef
) const;
//- Calculates the reaction rates
virtual void calculate();
@ -186,6 +190,7 @@ public:
//- Update matrix RR for reaction i. Used by EulerImplicit
// (Not implemented)
virtual void updateRRInReactionI
(
const label i,

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2013-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -65,6 +65,9 @@ class solidChemistryModel
{
// Private Member Functions
//- Disallow copy constructor
solidChemistryModel(const solidChemistryModel&);
//- Disallow default bitwise assignment
void operator=(const solidChemistryModel&);
@ -151,6 +154,7 @@ public:
label& rRef
) const = 0;
//- Return the reaction rate for iReaction and the reference
// species and charateristic times
virtual scalar omegaI
@ -167,6 +171,10 @@ public:
label& rRef
) const = 0;
//- Calculates the reaction rates
virtual void calculate() = 0;
//- Update concentrations in reaction i given dt and reaction rate
// omega used by sequential solver
virtual void updateConcsInReactionI
@ -194,11 +202,8 @@ public:
simpleMatrix<scalar>& RR
) const = 0;
//- Calculates the reaction rates
virtual void calculate() = 0;
// Chemistry model functions
// Solid Chemistry model functions
//- Return const access to the chemical source terms for solids
inline const DimensionedField<scalar, volMesh>& RRs
@ -209,13 +214,6 @@ public:
//- Return total solid source term
inline tmp<DimensionedField<scalar, volMesh> > RRs() const;
//- Return const access to the total source terms
inline const DimensionedField<scalar, volMesh>& RR
(
const label i
) const;
//- Solve the reaction system for the given start time and time
// step and return the characteristic time
virtual scalar solve(const scalar t0, const scalar deltaT) = 0;

View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2011-2013 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2013-2013 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -95,16 +95,4 @@ Foam::solidChemistryModel<CompType, SolidThermo>::RRs() const
}
template<class CompType, class SolidThermo>
inline const Foam::DimensionedField<Foam::scalar, Foam::volMesh>&
Foam::solidChemistryModel<CompType, SolidThermo>::RR
(
const label i
) const
{
notImplemented("solidChemistryModel::RR(const label)");
return (DimensionedField<scalar, volMesh>::null());
}
// ************************************************************************* //

View File

@ -52,9 +52,8 @@ namespace Foam
defineTemplateTypeNameAndDebugWithName \
( \
SS##Schem##Comp##SThermo##GThermo, \
(#SS"<" + word(Schem::typeName_()) \
+ "<"#Comp"," + SThermo::typeName() \
+ "," + GThermo::typeName() + ">>").c_str(), \
(#SS"<"#Schem"<"#Comp"," + SThermo::typeName() + "," \
+ GThermo::typeName() + ">>").c_str(), \
0 \
); \
\
@ -77,14 +76,6 @@ namespace Foam
GThermo \
); \
\
makeSolidChemistrySolverType \
( \
EulerImplicit, \
SolidChem, \
Comp, \
SThermo, \
GThermo \
); \
\
makeSolidChemistrySolverType \
( \
@ -104,7 +95,6 @@ namespace Foam
GThermo \
);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

View File

@ -30,7 +30,7 @@ Description
\*---------------------------------------------------------------------------*/
#ifndef makeSolidThermo_H
#define makesolidThermo_H
#define makeSolidThermo_H
#include "addToRunTimeSelectionTable.H"