Commit 1332dad6 authored by mattijs's avatar mattijs
Browse files

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

parents 1d021ecd 4ced40cb
......@@ -102,7 +102,8 @@ Foam::lagrangianFieldDecomposer::lagrangianFieldDecomposer
(
positions_,
ppi.position(),
procCelli
procCelli,
false
)
);
}
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
\\/ 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/>.
Typedef
Foam::floatScalar
Description
single floating point number identical to float
SourceFiles
floatScalar.C
\*---------------------------------------------------------------------------*/
#ifndef floatScalar_H
#define floatScalar_H
#include "label.H"
#include "word.H"
#include <cmath>
#ifdef ibm
float lgamma(float);
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Define floatScalar as a float
namespace Foam
{
typedef float floatScalar;
}
/*
// Define floatScalar as a float
namespace Foam
{
typedef float floatScalar;
// Largest and smallest floatScalar values allowed in certain
// parts of the code (6 is the number of significant figures in an
// IEEE single precision number. See limits.h or float.h)
static const floatScalar GREAT = 1.0e+6;
static const floatScalar VGREAT = 1.0e+37;
static const floatScalar SMALL = 1.0e-6;
static const floatScalar VSMALL = 1.0e-37;
}
*/
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "pTraits.H"
#include "products.H"
#include "direction.H"
namespace Foam
{
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// template specialisation for pTraits<floatScalar>
template<>
class pTraits<floatScalar>
{
floatScalar p_;
public:
//- Component type
typedef floatScalar cmptType;
// Member constants
enum
{
dim = 3, // Dimensionality of space
rank = 0, // Rank od floatScalar is 0
nComponents = 1 // Number of components in floatScalar is 1
};
// Static data members
static const char* const typeName;
static const char* componentNames[];
static const floatScalar zero;
static const floatScalar one;
// Constructors
//- Construct from Istream
pTraits(Istream& is);
// Member Functions
operator floatScalar() const
{
return p_;
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
//template<class Cmpt>
//class typeOfRank<Cmpt, 0>
//{
//public:
//
// typedef Cmpt type;
//};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Return a string representation of a floatScalar
word name(const floatScalar s)
{
return name(scalar(s));
}
#define MAXMINPOW(retType, type1, type2) \
\
MAXMIN(retType, type1, type2) \
\
inline float pow(const type1 s, const type2 e) \
{ \
return ::pow(float(s), float(e)); \
}
//MAXMINPOW(float, float, float)
//MAXMINPOW(float, float, int)
//MAXMINPOW(float, int, float)
//MAXMINPOW(float, float, long)
//MAXMINPOW(float, long, float)
//MAXMINPOW(float, float, int)
//MAXMINPOW(float, int, float)
//MAXMINPOW(float, float, long)
//MAXMINPOW(float, long, float)
#undef MAXMINPOW
inline floatScalar mag(const floatScalar s)
{
return ::fabs(s);
}
inline floatScalar sign(const floatScalar s)
{
return (s >= 0)? 1: -1;
}
inline floatScalar pos(const floatScalar s)
{
return (s >= 0)? 1: 0;
}
inline floatScalar neg(const floatScalar s)
{
return (s < 0)? 1: 0;
}
inline floatScalar limit(const floatScalar s1, const floatScalar s2)
{
return (mag(s1) < mag(s2))? s1: 0.0;
}
inline floatScalar magSqr(const floatScalar s)
{
return s*s;
}
inline floatScalar sqr(const floatScalar s)
{
return s*s;
}
inline floatScalar pow3(const floatScalar s)
{
return s*s*s;
}
inline floatScalar pow4(const floatScalar s)
{
return sqr(sqr(s));
}
inline floatScalar cmptAv(const floatScalar s)
{
return s;
}
inline floatScalar cmptMag(const floatScalar s)
{
return mag(s);
}
inline floatScalar scale(const floatScalar s, const floatScalar d)
{
return s*d;
}
#define transFunc(func) \
inline floatScalar func(const floatScalar s) \
{ \
return ::func(s); \
}
// Standard C++ transcendental functions
transFunc(sqrt)
transFunc(exp)
transFunc(log)
transFunc(log10)
transFunc(sin)
transFunc(cos)
transFunc(tan)
transFunc(asin)
transFunc(acos)
transFunc(atan)
transFunc(sinh)
transFunc(cosh)
transFunc(tanh)
transFunc(asinh)
transFunc(acosh)
transFunc(atanh)
// Standard ANSI-C (but not in <cmath>) transcendental functions
transFunc(erf)
transFunc(erfc)
transFunc(lgamma)
transFunc(j0)
transFunc(j1)
transFunc(y0)
transFunc(y1)
#undef transFunc
// Stabilisation around zero for division
inline floatScalar stabilise(const floatScalar s, const floatScalar small)
{
if (s >= 0)
{
return s + small;
}
else
{
return s - small;
}
}
// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
//floatScalar readScalar(Istream& is);
Istream& operator>>(Istream&, floatScalar&);
Ostream& operator<<(Ostream&, const floatScalar);
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
......@@ -162,7 +162,8 @@ Foam::Particle<ParticleType>::Particle
(
const Cloud<ParticleType>& cloud,
const vector& position,
const label cellI
const label cellI,
bool doCellFacePt
)
:
cloud_(cloud),
......@@ -175,7 +176,10 @@ Foam::Particle<ParticleType>::Particle
origProc_(Pstream::myProcNo()),
origId_(cloud_.getNewParticleID())
{
initCellFacePt();
if (doCellFacePt)
{
initCellFacePt();
}
}
......
......@@ -346,7 +346,8 @@ public:
(
const Cloud<ParticleType>&,
const vector& position,
const label cellI
const label cellI,
bool doCellFacePt = true
);
//- Construct from Istream
......
......@@ -68,15 +68,17 @@ public:
Particle<passiveParticle>(c, position, cellI, tetFaceI, tetPtI)
{}
//- Construct from components, with searching for tetFace and tetPt
//- Construct from components, with searching for tetFace and
// tetPt unless disabled by doCellFacePt = false.
passiveParticle
(
const Cloud<passiveParticle>& c,
const vector& position,
const label cellI
const label cellI,
bool doCellFacePt = true
)
:
Particle<passiveParticle>(c, position, cellI)
Particle<passiveParticle>(c, position, cellI, doCellFacePt)
{}
//- Construct from Istream
......
......@@ -57,6 +57,7 @@ Foam::coalParcel::coalParcel
const label typeId,
const scalar nParticle0,
const scalar d0,
const scalar dTarget0,
const vector& U0,
const vector& f0,
const vector& pi0,
......@@ -78,6 +79,7 @@ Foam::coalParcel::coalParcel
typeId,
nParticle0,
d0,
dTarget0,
U0,
f0,
pi0,
......
......@@ -79,6 +79,7 @@ public:
const label typeId,
const scalar nParticle0,
const scalar d0,
const scalar dTarget0,
const vector& U0,
const vector& f0,
const vector& pi0,
......
......@@ -43,6 +43,52 @@ void Foam::KinematicCloud<ParcelType>::preEvolve()
{
this->dispersion().cacheFields(true);
forces_.cacheFields(true, interpolationSchemes_);
updateCellOccupancy();
}
template<class ParcelType>
void Foam::KinematicCloud<ParcelType>::buildCellOccupancy()
{
if (cellOccupancyPtr_.empty())
{
cellOccupancyPtr_.reset
(
new List<DynamicList<ParcelType*> >(mesh_.nCells())
);
}
else if (cellOccupancyPtr_().size() != mesh_.nCells())
{
// If the size of the mesh has changed, reset the
// cellOccupancy size
cellOccupancyPtr_().setSize(mesh_.nCells());
}
List<DynamicList<ParcelType*> >& cellOccupancy = cellOccupancyPtr_();
forAll(cellOccupancy, cO)
{
cellOccupancy[cO].clear();
}
forAllIter(typename KinematicCloud<ParcelType>, *this, iter)
{
cellOccupancy[iter().cell()].append(&iter());
}
}
template<class ParcelType>
void Foam::KinematicCloud<ParcelType>::updateCellOccupancy()
{
// Only build the cellOccupancy if the pointer is set, i.e. it has
// been requested before.
if (cellOccupancyPtr_.valid())
{
buildCellOccupancy();
}
}
......@@ -80,8 +126,19 @@ void Foam::KinematicCloud<ParcelType>::evolveCloud()
g_.value()
);
label preInjectionSize = this->size();
this->surfaceFilm().inject(td);
// Update the cellOccupancy if the size of the cloud has changed
// during the injection.
if (preInjectionSize != this->size())
{
updateCellOccupancy();
preInjectionSize = this->size();
}
this->injection().inject(td);
if (coupled_)
......@@ -89,6 +146,18 @@ void Foam::KinematicCloud<ParcelType>::evolveCloud()
resetSourceTerms();
}
// Assume that motion will update the cellOccupancy as necessary
// before it is required.
motion(td);
}
template<class ParcelType>
void Foam::KinematicCloud<ParcelType>::motion
(
typename ParcelType::trackData& td
)
{
// Sympletic leapfrog integration of particle forces:
// + apply half deltaV with stored force
// + move positions with new velocity
......@@ -136,6 +205,8 @@ void Foam::KinematicCloud<ParcelType>::moveCollide
// td.part() = ParcelType::trackData::tpRotationalTrack;
// Cloud<ParcelType>::move(td);
updateCellOccupancy();
this->collision().collide();
td.part() = ParcelType::trackData::tpVelocityHalfStep;
......@@ -194,6 +265,7 @@ Foam::KinematicCloud<ParcelType>::KinematicCloud
particleProperties_.lookup("cellValueSourceCorrection")
),
rndGen_(label(0)),
cellOccupancyPtr_(),
rho_(rho),
U_(U),
mu_(mu),
......
......@@ -136,6 +136,9 @@ protected:
//- Random number generator - used by some injection routines
Random rndGen_;
//- Cell occupancy information for each parcel, (demand driven)
autoPtr<List<DynamicList<ParcelType*> > > cellOccupancyPtr_;
// References to the carrier gas fields
......@@ -209,10 +212,20 @@ protected:
//- Pre-evolve
void preEvolve();
//- Build the cellOccupancy
void buildCellOccupancy();
//- Update (i.e. build) the cellOccupancy if it has
// already been used
void updateCellOccupancy();
//- Evolve the cloud
void evolveCloud();
//- Move-collide
//- Particle motion
void motion(typename ParcelType::trackData& td);
//- Move-collide particles
void moveCollide(typename ParcelType::trackData& td);
//- Post-evolve
......@@ -284,6 +297,12 @@ public:
//- Return refernce to the random object
inline Random& rndGen();
//- Return the cell occupancy information for each
// parcel, non-const access, the caller is
// responsible for updating it for its own purposes
// if particles are removed or created.
inline List<DynamicList<ParcelType*> >& cellOccupancy();
// References to the carrier gas fields
......
......@@ -302,6 +302,19 @@ inline Foam::Random& Foam::KinematicCloud<ParcelType>::rndGen()
}
template<class ParcelType>
inline Foam::List<Foam::DynamicList<ParcelType*> >&
Foam::KinematicCloud<ParcelType>::cellOccupancy()
{
if (cellOccupancyPtr_.empty())
{
buildCellOccupancy();
}
return cellOccupancyPtr_();
}
template<class ParcelType>
inline Foam::DimensionedField<Foam::vector, Foam::volMesh>&
Foam::KinematicCloud<ParcelType>::UTrans()
......
......@@ -120,8 +120,19 @@ void Foam::ReactingCloud<ParcelType>::evolveCloud()
this->g().value()
);
label preInjectionSize = this->size();
this->surfaceFilm().inject(td);
// Update the cellOccupancy if the size of the cloud has changed
// during the injection.
if (preInjectionSize != this->size())
{
this->updateCellOccupancy();
preInjectionSize = this->size();
}
this->injection().inject(td);
if (this->coupled())
......@@ -129,7 +140,19 @@ void Foam::ReactingCloud<ParcelType>::evolveCloud()
resetSourceTerms();
}
Cloud<ParcelType>::move(td);
// Assume that motion will update the cellOccupancy as necessary
// before it is required.
motion(td);
}
template<class ParcelType>
void Foam::ReactingCloud<ParcelType>::motion
(
typename ParcelType::trackData& td
)
{
ThermoCloud<ParcelType>::motion(td);
}
......