zhyang-dev/openfoam
1
0
Fork 0
mirror of https://develop.openfoam.com/Development/openfoam.git synced 2025-11-28 03:28:01 +00:00
Code Issues Actions Packages Projects Releases Wiki Activity

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

Browse Source
This commit is contained in:
andy
2011-07-07 09:15:40 +01:00
parent b2d939833a 073f4960f7
commit ac26124658
44 changed files with 781 additions and 296 deletions
Download Patch File Download Diff File Expand all files Collapse all files

109
src/OSspecific/POSIX/POSIX.C
View File

@ -57,6 +57,14 @@ Description
#include <netinet/in.h>
#ifdef USE_RANDOM
# include <climits>
# if INT_MAX != 2147483647
# error "INT_MAX != 2147483647"
# error "The random number generator may not work!"
# endif
#endif
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
defineTypeNameAndDebug(Foam::POSIX, 0);
@ -68,16 +76,19 @@ pid_t Foam::pid()
return ::getpid();
}
pid_t Foam::ppid()
{
return ::getppid();
}
pid_t Foam::pgid()
{
return ::getpgrp();
}
bool Foam::env(const word& envName)
{
return ::getenv(envName.c_str()) != NULL;
@ -253,10 +264,11 @@ bool Foam::chDir(const fileName& dir)
}
Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
Foam::fileNameList Foam::findEtcFiles(const fileName& name, bool mandatory)
{
//
// search for user files in
fileNameList results;
// Search for user files in
// * ~/.OpenFOAM/VERSION
// * ~/.OpenFOAM
//
@ -266,19 +278,17 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
fileName fullName = searchDir/FOAMversion/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
fullName = searchDir/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
}
//
// search for group (site) files in
// Search for group (site) files in
// * $WM_PROJECT_SITE/VERSION
// * $WM_PROJECT_SITE
//
@ -290,19 +300,18 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
fileName fullName = searchDir/FOAMversion/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
fullName = searchDir/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
}
}
else
{
//
// OR search for group (site) files in
// * $WM_PROJECT_INST_DIR/site/VERSION
// * $WM_PROJECT_INST_DIR/site
@ -313,20 +322,18 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
fileName fullName = searchDir/"site"/FOAMversion/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
fullName = searchDir/"site"/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
}
}
//
// search for other (shipped) files in
// Search for other (shipped) files in
// * $WM_PROJECT_DIR/etc
//
searchDir = getEnv("WM_PROJECT_DIR");
@ -335,24 +342,41 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
fileName fullName = searchDir/"etc"/name;
if (isFile(fullName))
{
return fullName;
results.append(fullName);
}
}
// Not found
// abort if the file is mandatory, otherwise return null
if (mandatory)
if (results.empty())
{
std::cerr
<< "--> FOAM FATAL ERROR in Foam::findEtcFile() :"
" could not find mandatory file\n '"
<< name.c_str() << "'\n\n" << std::endl;
::exit(1);
// Abort if the file is mandatory, otherwise return null
if (mandatory)
{
std::cerr
<< "--> FOAM FATAL ERROR in Foam::findEtcFiles() :"
" could not find mandatory file\n '"
<< name.c_str() << "'\n\n" << std::endl;
::exit(1);
}
}
// Return null-constructed fileName rather than fileName::null
// to avoid cyclic dependencies in the construction of globals
return fileName();
// Return list of matching paths or empty list if none found
return results;
}
Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
{
fileNameList results(findEtcFiles(name, mandatory));
if (results.size())
{
return results[0];
}
else
{
return fileName();
}
}
@ -890,7 +914,6 @@ bool Foam::mvBak(const fileName& src, const std::string& ext)
}
// Remove a file, returning true if successful otherwise false
bool Foam::rm(const fileName& file)
{
@ -1221,4 +1244,34 @@ Foam::fileNameList Foam::dlLoaded()
}
void Foam::osRandomSeed(const label seed)
{
#ifdef USE_RANDOM
srandom((unsigned int)seed);
#else
srand48(seed);
#endif
}
Foam::label Foam::osRandomInteger()
{
#ifdef USE_RANDOM
return random();
#else
return lrand48();
#endif
}
Foam::scalar Foam::osRandomDouble()
{
#ifdef USE_RANDOM
return (scalar)random()/INT_MAX;
#else
return drand48();
#endif
}
// ************************************************************************* //

8
src/OpenFOAM/Make/files
View File

@ -1,6 +1,6 @@
global/global.Cver
global/constants/constants.C
global/constants/dimensionedConstants.C
/* global/constants/constants.C in global.Cver */
/* global/constants/dimensionedConstants.C in global.Cver */
global/argList/argList.C
global/clock/clock.C
@ -136,7 +136,7 @@ $(StringStreams)/StringStreamsPrint.C
Pstreams = $(Streams)/Pstreams
$(Pstreams)/UIPstream.C
$(Pstreams)/IPstream.C
$(Pstreams)/UPstream.C
/* $(Pstreams)/UPstream.C in global.Cver */
$(Pstreams)/UPstreamCommsStruct.C
$(Pstreams)/Pstream.C
$(Pstreams)/UOPstream.C
@ -181,7 +181,7 @@ $(IOobject)/IOobjectReadHeader.C
$(IOobject)/IOobjectWriteHeader.C
regIOobject = db/regIOobject
$(regIOobject)/regIOobject.C
/* $(regIOobject)/regIOobject.C in global.Cver */
$(regIOobject)/regIOobjectRead.C
$(regIOobject)/regIOobjectWrite.C

21
src/OpenFOAM/db/dictionary/functionEntries/codeStream/codeStream.C
View File

@ -197,8 +197,25 @@ Foam::functionEntries::codeStream::getFunction
}
}
// all processes must wait for compile to finish
reduce(create, orOp<bool>());
//- We don't know whether this code was from IOdictionary
// (possibly read on master only) or from e.g. Field so cannot
// decide here.
//// all processes must wait for compile to finish - except if this
//// file is only read on the master
//bool masterOnly =
// (
// regIOobject::fileModificationChecking
// == regIOobject::timeStampMaster
// )
// || (
// regIOobject::fileModificationChecking
// == regIOobject::inotifyMaster
// );
//
//if (!masterOnly)
//{
reduce(create, orOp<bool>());
//}
if (isA<IOdictionary>(topDict(parentDict)))
{

13
src/OpenFOAM/global/debug/debug.C
View File

@ -75,10 +75,15 @@ Foam::dictionary& Foam::debug::controlDict()
{
if (!controlDictPtr_)
{
controlDictPtr_ = new dictionary
(
IFstream(findEtcFile("controlDict", true))()
);
fileNameList controlDictFiles = findEtcFiles("controlDict", true);
controlDictPtr_ = new dictionary();
forAllReverse(controlDictFiles, cdfi)
{
controlDictPtr_->merge
(
dictionary(IFstream(controlDictFiles[cdfi])())
);
}
}
return *controlDictPtr_;

18
src/OpenFOAM/global/global.Cver
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -65,6 +65,22 @@ bool Foam::JobInfo::constructed(false);
#include "debug.C"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Read file modification checking switches
#include "regIOobject.C"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Read parallel communication switches
#include "UPstream.C"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Read constants
#include "constants.C"
#include "dimensionedConstants.C"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Read and set cell models

23
src/OpenFOAM/include/OSspecific.H
View File

@ -93,7 +93,7 @@ fileName cwd();
// else return false
bool chDir(const fileName& dir);
//- Search for a file from user/group/shipped directories.
//- Search for files from user/group/shipped directories.
// The search scheme allows for version-specific and
// version-independent files using the following hierarchy:
// - \b user settings:
@ -108,7 +108,14 @@ bool chDir(const fileName& dir);
// - \b other (shipped) settings:
// - $WM_PROJECT_DIR/etc/
//
// \return The full path name or fileName() if the name cannot be found
// \return The list of full paths of all the matching files or
// an empty list if the name cannot be found.
// Optionally abort if the file cannot be found
fileNameList findEtcFiles(const fileName&, bool mandatory=false);
//- Search for a file using findEtcFiles.
// \return The full path name of the first file found which in the
// search hierarchy or an empty fileName if the name cannot be found.
// Optionally abort if the file cannot be found
fileName findEtcFile(const fileName&, bool mandatory=false);
@ -200,6 +207,18 @@ bool dlSymFound(void* handle, const std::string& symbol);
fileNameList dlLoaded();
// Low level random numbers. Use Random class instead.
//- Seed random number generator.
void osRandomSeed(const label seed);
//- Return random integer (uniform distribution between 0 and 2^31)
label osRandomInteger();
//- Return random double precision (uniform distribution between 0 and 1)
scalar osRandomDouble();
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam

33
src/OpenFOAM/primitives/random/Random/Random.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -24,6 +24,7 @@ License
\*---------------------------------------------------------------------------*/
#include "Random.H"
#include "OSspecific.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -37,11 +38,6 @@ namespace Foam
# error "The random number generator may not work!"
#endif
#ifdef USE_RANDOM
# include <climits>
#else
# include <cstdlib>
#endif
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
@ -57,22 +53,13 @@ Random::Random(const label seed)
Seed = 1;
}
# ifdef USE_RANDOM
srandom((unsigned int)Seed);
# else
srand48(Seed);
# endif
osRandomSeed(Seed);
}
int Random::bit()
{
# ifdef USE_RANDOM
if (random() > INT_MAX/2)
# else
if (lrand48() > INT_MAX/2)
# endif
if (osRandomInteger() > INT_MAX/2)
{
return 1;
}
@ -85,11 +72,7 @@ int Random::bit()
scalar Random::scalar01()
{
# ifdef USE_RANDOM
return (scalar)random()/INT_MAX;
# else
return drand48();
# endif
return osRandomDouble();
}
@ -140,11 +123,7 @@ tensor Random::tensor01()
label Random::integer(const label lower, const label upper)
{
# ifdef USE_RANDOM
return lower + (random() % (upper+1-lower));
# else
return lower + (lrand48() % (upper+1-lower));
# endif
return lower + (osRandomInteger() % (upper+1-lower));
}

10
src/OpenFOAM/primitives/random/cachedRandom/cachedRandom.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2010-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2010-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -24,7 +24,7 @@ License
\*---------------------------------------------------------------------------*/
#include "cachedRandom.H"
#include <cstdlib>
#include "OSspecific.H"
#if INT_MAX != 2147483647
# error "INT_MAX != 2147483647"
@ -37,7 +37,7 @@ Foam::scalar Foam::cachedRandom::scalar01()
{
if (sampleI_ < 0)
{
return drand48();
return osRandomDouble();
}
if (sampleI_ == samples_.size() - 1)
@ -76,7 +76,7 @@ Foam::cachedRandom::cachedRandom(const label seed, const label count)
}
// Initialise samples
srand48(seed_);
osRandomSeed(seed_);
forAll(samples_, i)
{
samples_[i] = drand48();
@ -98,7 +98,7 @@ Foam::cachedRandom::cachedRandom(const cachedRandom& cr, const bool reset)
) << "Copy constructor called, but samples not being cached. "
<< "This may lead to non-repeatable behaviour" << endl;
srand48(seed_);
osRandomSeed(seed_);
}
else if (reset)
{

73
src/finiteVolume/fields/fvPatchFields/derived/fanPressure/fanPressureFvPatchScalarField.C
View File

@ -59,10 +59,7 @@ Foam::fanPressureFvPatchScalarField::fanPressureFvPatchScalarField
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(p, iF),
phiName_("phi"),
rhoName_("rho"),
p0_(p.size(), 0.0),
totalPressureFvPatchScalarField(p, iF),
fanCurve_(),
direction_(ffdOut)
{}
@ -76,10 +73,7 @@ Foam::fanPressureFvPatchScalarField::fanPressureFvPatchScalarField
const fvPatchFieldMapper& mapper
)
:
fixedValueFvPatchScalarField(ptf, p, iF, mapper),
phiName_(ptf.phiName_),
rhoName_(ptf.rhoName_),
p0_(ptf.p0_, mapper),
totalPressureFvPatchScalarField(ptf, p, iF, mapper),
fanCurve_(ptf.fanCurve_),
direction_(ptf.direction_)
{}
@ -92,10 +86,7 @@ Foam::fanPressureFvPatchScalarField::fanPressureFvPatchScalarField
const dictionary& dict
)
:
fixedValueFvPatchScalarField(p, iF),
phiName_(dict.lookupOrDefault<word>("phi", "phi")),
rhoName_(dict.lookupOrDefault<word>("rho", "rho")),
p0_("p0", dict, p.size()),
totalPressureFvPatchScalarField(p, iF),
fanCurve_(dict),
direction_(fanFlowDirectionNames_.read(dict.lookup("direction")))
{
@ -109,10 +100,7 @@ Foam::fanPressureFvPatchScalarField::fanPressureFvPatchScalarField
const fanPressureFvPatchScalarField& pfopsf
)
:
fixedValueFvPatchScalarField(pfopsf),
phiName_(pfopsf.phiName_),
rhoName_(pfopsf.rhoName_),
p0_(pfopsf.p0_),
totalPressureFvPatchScalarField(pfopsf),
fanCurve_(pfopsf.fanCurve_),
direction_(pfopsf.direction_)
{}
@ -124,10 +112,7 @@ Foam::fanPressureFvPatchScalarField::fanPressureFvPatchScalarField
const DimensionedField<scalar, volMesh>& iF
)
:
fixedValueFvPatchScalarField(pfopsf, iF),
phiName_(pfopsf.phiName_),
rhoName_(pfopsf.rhoName_),
p0_(pfopsf.p0_),
totalPressureFvPatchScalarField(pfopsf, iF),
fanCurve_(pfopsf.fanCurve_),
direction_(pfopsf.direction_)
{}
@ -144,7 +129,7 @@ void Foam::fanPressureFvPatchScalarField::updateCoeffs()
// Retrieve flux field
const surfaceScalarField& phi =
db().lookupObject<surfaceScalarField>(phiName_);
db().lookupObject<surfaceScalarField>(phiName());
const fvsPatchField<scalar>& phip =
patch().patchField<surfaceScalarField, scalar>(phi);
@ -161,7 +146,7 @@ void Foam::fanPressureFvPatchScalarField::updateCoeffs()
else if (phi.dimensions() == dimVelocity*dimArea*dimDensity)
{
const scalarField& rhop =
patch().lookupPatchField<volScalarField, scalar>(rhoName_);
patch().lookupPatchField<volScalarField, scalar>(rhoName());
aveFlowRate = dir*gSum(phip/rhop)/gSum(patch().magSf());
}
else
@ -174,51 +159,23 @@ void Foam::fanPressureFvPatchScalarField::updateCoeffs()
<< exit(FatalError);
}
// Normal flow through fan
if (aveFlowRate >= 0.0)
{
// Pressure drop for this flow rate
const scalar pdFan = fanCurve_(aveFlowRate);
// Pressure drop for this flow rate
const scalar pdFan = fanCurve_(max(aveFlowRate, 0.0));
operator==(p0_ - dir*pdFan);
}
// Reverse flow
else
{
// Assume that fan has stalled if flow reversed
// i.e. apply dp for zero flow rate
const scalar pdFan = fanCurve_(0);
// Flow speed across patch
scalarField Up = phip/(patch().magSf());
// Pressure drop associated withback flow = dynamic pressure
scalarField pdBackFlow = 0.5*magSqr(Up);
if (phi.dimensions() == dimVelocity*dimArea*dimDensity)
{
const scalarField& rhop =
patch().lookupPatchField<volScalarField, scalar>(rhoName_);
pdBackFlow /= rhop;
}
operator==(p0_ - dir*(pdBackFlow + pdFan));
}
fixedValueFvPatchScalarField::updateCoeffs();
totalPressureFvPatchScalarField::updateCoeffs
(
p0() - dir*pdFan,
patch().lookupPatchField<volVectorField, vector>(UName())
);
}
void Foam::fanPressureFvPatchScalarField::write(Ostream& os) const
{
fvPatchScalarField::write(os);
os.writeKeyword("phi") << phiName_ << token::END_STATEMENT << nl;
os.writeKeyword("rho") << rhoName_ << token::END_STATEMENT << nl;
totalPressureFvPatchScalarField::write(os);
fanCurve_.write(os);
os.writeKeyword("direction")
<< fanFlowDirectionNames_[direction_] << token::END_STATEMENT << nl;
p0_.writeEntry("p0", os);
writeEntry("value", os);
}

31
src/finiteVolume/fields/fvPatchFields/derived/fanPressure/fanPressureFvPatchScalarField.H
View File

@ -25,7 +25,7 @@ Class
Foam::fanPressureFvPatchScalarField
Description
Assigns pressure inlet or outlet condition for a fan.
Assigns pressure inlet or outlet total pressure condition for a fan.
User specifies:
- pressure drop vs volumetric flow rate table (fan curve) file name;
@ -56,8 +56,8 @@ Description
\endverbatim
See Also
Foam::interpolationTable and
Foam::timeVaryingFlowRateInletVelocityFvPatchVectorField
Foam::totalPressureFvPatchScalarField and
Foam::interpolationTable
SourceFiles
fanPressureFvPatchScalarField.C
@ -67,8 +67,7 @@ SourceFiles
#ifndef fanPressureFvPatchScalarField_H
#define fanPressureFvPatchScalarField_H
#include "fvPatchFields.H"
#include "fixedValueFvPatchFields.H"
#include "totalPressureFvPatchScalarField.H"
#include "interpolationTable.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@ -82,21 +81,10 @@ namespace Foam
class fanPressureFvPatchScalarField
:
public fixedValueFvPatchScalarField
public totalPressureFvPatchScalarField
{
// Private data
//- Name of the flux transporting the field
word phiName_;
//- Name of the density field
word rhoName_;
//- Total pressure
scalarField p0_;
//- Tabulated fan curve
interpolationTable<scalar> fanCurve_;
public:
//- Fan flow direction
enum fanFlowDirection
@ -107,6 +95,13 @@ class fanPressureFvPatchScalarField
static const NamedEnum<fanFlowDirection, 2> fanFlowDirectionNames_;
private:
// Private data
//- Tabulated fan curve
interpolationTable<scalar> fanCurve_;
//- Direction of flow through the fan relative to patch
fanFlowDirection direction_;

2
src/finiteVolume/fields/fvPatchFields/derived/rotatingTotalPressure/rotatingTotalPressureFvPatchScalarField.C
View File

@ -112,7 +112,7 @@ void Foam::rotatingTotalPressureFvPatchScalarField::updateCoeffs()
+ rotationVelocity
);
totalPressureFvPatchScalarField::updateCoeffs(Up);
totalPressureFvPatchScalarField::updateCoeffs(p0(), Up);
}

2
src/finiteVolume/fields/fvPatchFields/derived/rotatingTotalPressure/rotatingTotalPressureFvPatchScalarField.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License

20
src/finiteVolume/fields/fvPatchFields/derived/totalPressure/totalPressureFvPatchScalarField.C
View File

@ -153,7 +153,11 @@ void Foam::totalPressureFvPatchScalarField::rmap
}
void Foam::totalPressureFvPatchScalarField::updateCoeffs(const vectorField& Up)
void Foam::totalPressureFvPatchScalarField::updateCoeffs
(
const scalarField& p0p,
const vectorField& Up
)
{
if (updated())
{
@ -165,7 +169,7 @@ void Foam::totalPressureFvPatchScalarField::updateCoeffs(const vectorField& Up)
if (psiName_ == "none" && rhoName_ == "none")
{
operator==(p0_ - 0.5*(1.0 - pos(phip))*magSqr(Up));
operator==(p0p - 0.5*(1.0 - pos(phip))*magSqr(Up));
}
else if (rhoName_ == "none")
{
@ -178,7 +182,7 @@ void Foam::totalPressureFvPatchScalarField::updateCoeffs(const vectorField& Up)
operator==
(
p0_
p0p
/pow
(
(1.0 + 0.5*psip*gM1ByG*(1.0 - pos(phip))*magSqr(Up)),
@ -188,7 +192,7 @@ void Foam::totalPressureFvPatchScalarField::updateCoeffs(const vectorField& Up)
}
else
{
operator==(p0_/(1.0 + 0.5*psip*(1.0 - pos(phip))*magSqr(Up)));
operator==(p0p/(1.0 + 0.5*psip*(1.0 - pos(phip))*magSqr(Up)));
}
}
else if (psiName_ == "none")
@ -196,7 +200,7 @@ void Foam::totalPressureFvPatchScalarField::updateCoeffs(const vectorField& Up)
const fvPatchField<scalar>& rho =
patch().lookupPatchField<volScalarField, scalar>(rhoName_);
operator==(p0_ - 0.5*rho*(1.0 - pos(phip))*magSqr(Up));
operator==(p0p - 0.5*rho*(1.0 - pos(phip))*magSqr(Up));
}
else
{
@ -220,7 +224,11 @@ void Foam::totalPressureFvPatchScalarField::updateCoeffs(const vectorField& Up)
void Foam::totalPressureFvPatchScalarField::updateCoeffs()
{
updateCoeffs(patch().lookupPatchField<volVectorField, vector>(UName_));
updateCoeffs
(
p0(),
patch().lookupPatchField<volVectorField, vector>(UName())
);
}

49
src/finiteVolume/fields/fvPatchFields/derived/totalPressure/totalPressureFvPatchScalarField.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -157,6 +157,45 @@ public:
return UName_;
}
//- Return the name of the flux field
const word& phiName() const
{
return phiName_;
}
//- Return reference to the name of the flux field
// to allow adjustment
word& phiName()
{
return phiName_;
}
//- Return the name of the density field
const word& rhoName() const
{
return rhoName_;
}
//- Return reference to the name of the density field
// to allow adjustment
word& rhoName()
{
return rhoName_;
}
//- Return the name of the compressibility field
const word& psiName() const
{
return psiName_;
}
//- Return reference to the name of the compressibility field
// to allow adjustment
word& psiName()
{
return psiName_;
}
//- Return the heat capacity ratio
scalar gamma() const
{
@ -201,8 +240,12 @@ public:
// Evaluation functions
//- Update the coefficients associated with the patch field
// using the given patch velocity field
virtual void updateCoeffs(const vectorField& Up);
// using the given patch total pressure and velocity fields
virtual void updateCoeffs
(
const scalarField& p0p,
const vectorField& Up
);
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();

3
src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.C
View File

@ -214,7 +214,8 @@ void Foam::PairSpringSliderDashpot<CloudType>::evaluatePair
rHat_AB
*(kN*pow(normalOverlapMag, b_) - etaN*(U_AB & rHat_AB));
// Cohesion force
// Cohesion force, energy density multiplied by the area of
// particle-particle overlap
if (cohesion_)
{
fN_AB +=

42
src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallLocalSpringSliderDashpot/WallLocalSpringSliderDashpot.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2008-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -76,7 +76,8 @@ void Foam::WallLocalSpringSliderDashpot<CloudType>::evaluateWall
typename CloudType::parcelType& p,
const point& site,
const WallSiteData<vector>& data,
scalar pREff
scalar pREff,
bool cohesion
) const
{
// wall patch index
@ -88,14 +89,18 @@ void Foam::WallLocalSpringSliderDashpot<CloudType>::evaluateWall
scalar alpha = alpha_[wPI];
scalar b = b_[wPI];
scalar mu = mu_[wPI];
scalar cohesionEnergyDensity = cohesionEnergyDensity_[wPI];
cohesion = cohesion && cohesion_[wPI];
vector r_PW = p.position() - site;
vector U_PW = p.U() - data.wallData();
scalar normalOverlapMag = max(pREff - mag(r_PW), 0.0);
scalar r_PW_mag = mag(r_PW);
vector rHat_PW = r_PW/(mag(r_PW) + VSMALL);
scalar normalOverlapMag = max(pREff - r_PW_mag, 0.0);
vector rHat_PW = r_PW/(r_PW_mag + VSMALL);
scalar kN = (4.0/3.0)*sqrt(pREff)*Estar;
@ -105,6 +110,16 @@ void Foam::WallLocalSpringSliderDashpot<CloudType>::evaluateWall
rHat_PW
*(kN*pow(normalOverlapMag, b) - etaN*(U_PW & rHat_PW));
// Cohesion force, energy density multiplied by the area of wall/particle
// overlap
if (cohesion)
{
fN_PW +=
-cohesionEnergyDensity
*mathematical::pi*(sqr(pREff) - sqr(r_PW_mag))
*rHat_PW;
}
p.f() += fN_PW;
vector USlip_PW =
@ -168,6 +183,8 @@ Foam::WallLocalSpringSliderDashpot<CloudType>::WallLocalSpringSliderDashpot
alpha_(),
b_(),
mu_(),
cohesionEnergyDensity_(),
cohesion_(),
patchMap_(),
maxEstarIndex_(-1),
collisionResolutionSteps_
@ -212,6 +229,8 @@ Foam::WallLocalSpringSliderDashpot<CloudType>::WallLocalSpringSliderDashpot
alpha_.setSize(nWallPatches);
b_.setSize(nWallPatches);
mu_.setSize(nWallPatches);
cohesionEnergyDensity_.setSize(nWallPatches);
cohesion_.setSize(nWallPatches);
scalar maxEstar = -GREAT;
@ -238,6 +257,13 @@ Foam::WallLocalSpringSliderDashpot<CloudType>::WallLocalSpringSliderDashpot
mu_[wPI] = readScalar(patchCoeffDict.lookup("mu"));
cohesionEnergyDensity_[wPI] = readScalar
(
patchCoeffDict.lookup("cohesionEnergyDensity")
);
cohesion_[wPI] = (mag(cohesionEnergyDensity_[wPI]) > VSMALL);
if (Estar_[wPI] > maxEstar)
{
maxEstarIndex_ = wPI;
@ -325,20 +351,22 @@ void Foam::WallLocalSpringSliderDashpot<CloudType>::evaluateWall
p,
flatSitePoints[siteI],
flatSiteData[siteI],
pREff
pREff,
true
);
}
forAll(sharpSitePoints, siteI)
{
// Treating sharp sites like flat sites
// Treating sharp sites like flat sites, except suppress cohesion
evaluateWall
(
p,
sharpSitePoints[siteI],
sharpSiteData[siteI],
pREff
pREff,
false
);
}
}

9
src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallLocalSpringSliderDashpot/WallLocalSpringSliderDashpot.H
View File

@ -65,6 +65,12 @@ class WallLocalSpringSliderDashpot
//- Coefficient of friction in for tangential sliding
scalarList mu_;
//- Cohesion energy density [J/m^3]
scalarList cohesionEnergyDensity_;
// Switch cohesion on and off
boolList cohesion_;
//- Mapping the patch index to the model data
labelList patchMap_;
@ -115,7 +121,8 @@ class WallLocalSpringSliderDashpot
typename CloudType::parcelType& p,
const point& site,
const WallSiteData<vector>& data,
scalar pREff
scalar pREff,
bool cohesion
) const;

36
src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2008-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2008-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -77,16 +77,19 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
const point& site,
const WallSiteData<vector>& data,
scalar pREff,
scalar kN
scalar kN,
bool cohesion
) const
{
vector r_PW = p.position() - site;
vector U_PW = p.U() - data.wallData();
scalar normalOverlapMag = max(pREff - mag(r_PW), 0.0);
scalar r_PW_mag = mag(r_PW);
vector rHat_PW = r_PW/(mag(r_PW) + VSMALL);
scalar normalOverlapMag = max(pREff - r_PW_mag, 0.0);
vector rHat_PW = r_PW/(r_PW_mag + VSMALL);
scalar etaN = alpha_*sqrt(p.mass()*kN)*pow025(normalOverlapMag);
@ -94,6 +97,16 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
rHat_PW
*(kN*pow(normalOverlapMag, b_) - etaN*(U_PW & rHat_PW));
// Cohesion force, energy density multiplied by the area of wall/particle
// overlap
if (cohesion)
{
fN_PW +=
-cohesionEnergyDensity_
*mathematical::pi*(sqr(pREff) - sqr(r_PW_mag))
*rHat_PW;
}
p.f() += fN_PW;
vector USlip_PW =
@ -157,6 +170,11 @@ Foam::WallSpringSliderDashpot<CloudType>::WallSpringSliderDashpot
alpha_(readScalar(this->coeffDict().lookup("alpha"))),
b_(readScalar(this->coeffDict().lookup("b"))),
mu_(readScalar(this->coeffDict().lookup("mu"))),
cohesionEnergyDensity_
(
readScalar(this->coeffDict().lookup("cohesionEnergyDensity"))
),
cohesion_(false),
collisionResolutionSteps_
(
readScalar
@ -183,6 +201,8 @@ Foam::WallSpringSliderDashpot<CloudType>::WallSpringSliderDashpot
Estar_ = 1/((1 - sqr(pNu))/pE + (1 - sqr(nu))/E);
Gstar_ = 1/(2*((2 + pNu - sqr(pNu))/pE + (2 + nu - sqr(nu))/E));
cohesion_ = (mag(cohesionEnergyDensity_) > VSMALL);
}
@ -266,13 +286,14 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
flatSitePoints[siteI],
flatSiteData[siteI],
pREff,
kN
kN,
cohesion_
);
}
forAll(sharpSitePoints, siteI)
{
// Treating sharp sites like flat sites
// Treating sharp sites like flat sites, except suppress cohesion
evaluateWall
(
@ -280,7 +301,8 @@ void Foam::WallSpringSliderDashpot<CloudType>::evaluateWall
sharpSitePoints[siteI],
sharpSiteData[siteI],
pREff,
kN
kN,
false
);
}
}

9
src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/WallModel/WallSpringSliderDashpot/WallSpringSliderDashpot.H
View File

@ -65,6 +65,12 @@ class WallSpringSliderDashpot
//- Coefficient of friction in for tangential sliding
scalar mu_;
//- Cohesion energy density [J/m^3]
scalar cohesionEnergyDensity_;
// Switch cohesion on and off
bool cohesion_;
//- The number of steps over which to resolve the minimum
// harmonic approximation of the collision period
scalar collisionResolutionSteps_;
@ -110,7 +116,8 @@ class WallSpringSliderDashpot
const point& site,
const WallSiteData<vector>& data,
scalar pREff,
scalar kN
scalar kN,
bool cohesion
) const;

3
src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoRefineDriver.C
View File

@ -744,7 +744,8 @@ void Foam::autoRefineDriver::doRefine
100 // maxIter
);
// Introduce baffles at surface intersections
// Introduce baffles at surface intersections. Remove sections unreachable
// from keepPoint.
baffleAndSplitMesh(refineParams, prepareForSnapping, motionDict);
// Mesh is at its finest. Do optional zoning.

1
src/meshTools/Make/files
View File

@ -62,6 +62,7 @@ $(searchableSurface)/searchableSurfaces.C
$(searchableSurface)/searchableSurfacesQueries.C
$(searchableSurface)/searchableSurfaceWithGaps.C
$(searchableSurface)/triSurfaceMesh.C
$(searchableSurface)/closedTriSurfaceMesh.C
topoSets = sets/topoSets
$(topoSets)/cellSet.C

74
src/meshTools/searchableSurface/closedTriSurfaceMesh.C Normal file
View File

@ -0,0 +1,74 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2011 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/>.
\*---------------------------------------------------------------------------*/
#include "closedTriSurfaceMesh.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(closedTriSurfaceMesh, 0);
addToRunTimeSelectionTable(searchableSurface, closedTriSurfaceMesh, dict);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::closedTriSurfaceMesh::closedTriSurfaceMesh
(
const IOobject& io,
const triSurface& s
)
:
triSurfaceMesh(io, s)
{}
Foam::closedTriSurfaceMesh::closedTriSurfaceMesh(const IOobject& io)
:
triSurfaceMesh(io)
{}
Foam::closedTriSurfaceMesh::closedTriSurfaceMesh
(
const IOobject& io,
const dictionary& dict
)
:
triSurfaceMesh(io, dict)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::closedTriSurfaceMesh::~closedTriSurfaceMesh()
{}
// ************************************************************************* //

111
src/meshTools/searchableSurface/closedTriSurfaceMesh.H Normal file
View File

@ -0,0 +1,111 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2011 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/>.
Class
Foam::closedTriSurfaceMesh
Description
A triSurfaceMesh where it is forced to check volumeTypes, used for surfaces
that are topologically non-manifold (small holes or multiple parts) but are
geometrically essentially closed
SourceFiles
closedTriSurfaceMesh.C
\*---------------------------------------------------------------------------*/
#ifndef closedTriSurfaceMesh_H
#define closedTriSurfaceMesh_H
#include "triSurfaceMesh.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class closedTriSurfaceMesh Declaration
\*---------------------------------------------------------------------------*/
class closedTriSurfaceMesh
:
public triSurfaceMesh
{
private:
//- Disallow default bitwise copy construct
closedTriSurfaceMesh(const closedTriSurfaceMesh&);
//- Disallow default bitwise assignment
void operator=(const closedTriSurfaceMesh&);
public:
//- Runtime type information
TypeName("closedTriSurfaceMesh");
// Constructors
//- Construct from triSurface
closedTriSurfaceMesh(const IOobject&, const triSurface&);
//- Construct read.
closedTriSurfaceMesh(const IOobject& io);
//- Construct from IO and dictionary (used by searchableSurface).
// Dictionary may contain a 'scale' entry (eg, 0.001: mm -> m)
closedTriSurfaceMesh
(
const IOobject& io,
const dictionary& dict
);
// Destructor
virtual ~closedTriSurfaceMesh();
// Member Functions
//- Whether supports volume type, forcing to true to force getVolumeType
// queries for this type
virtual bool hasVolumeType() const
{
return true;
}
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

84
src/meshTools/searchableSurface/triSurfaceMesh.C
View File

@ -331,6 +331,7 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io, const triSurface& s)
),
triSurface(s),
tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
minQuality_(-1),
maxTreeDepth_(10),
surfaceClosed_(-1)
{}
@ -375,6 +376,7 @@ Foam::triSurfaceMesh::triSurfaceMesh(const IOobject& io)
)
),
tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
minQuality_(-1),
maxTreeDepth_(10),
surfaceClosed_(-1)
{}
@ -422,6 +424,7 @@ Foam::triSurfaceMesh::triSurfaceMesh
)
),
tolerance_(indexedOctree<treeDataTriSurface>::perturbTol()),
minQuality_(-1),
maxTreeDepth_(10),
surfaceClosed_(-1)
{
@ -444,6 +447,13 @@ Foam::triSurfaceMesh::triSurfaceMesh
<< tolerance_ << endl;
}
// Have optional minimum quality for normal calculation
if (dict.readIfPresent("minQuality", minQuality_) && minQuality_ > 0)
{
Info<< searchableSurface::name()
<< " : ignoring triangles with quality < "
<< minQuality_ << " for normals calculation." << endl;
}
// Have optional non-standard tree-depth to limit storage.
if (dict.readIfPresent("maxTreeDepth", maxTreeDepth_) && maxTreeDepth_ > 0)
@ -804,22 +814,70 @@ void Foam::triSurfaceMesh::getNormal
{
normal.setSize(info.size());
forAll(info, i)
if (minQuality_ >= 0)
{
if (info[i].hit())
{
label faceI = info[i].index();
//- Cached:
//normal[i] = faceNormals()[faceI];
// Make sure we don't use triangles with low quality since
// normal is not reliable.
//- Uncached
normal[i] = triSurface::operator[](faceI).normal(points());
normal[i] /= mag(normal[i]) + VSMALL;
}
else
const triSurface& s = static_cast<const triSurface&>(*this);
const labelListList& faceFaces = s.faceFaces();
forAll(info, i)
{
// Set to what?
normal[i] = vector::zero;
if (info[i].hit())
{
label faceI = info[i].index();
scalar qual = s[faceI].tri(points()).quality();
if (qual < minQuality_)
{
// Search neighbouring triangles
const labelList& fFaces = faceFaces[faceI];
forAll(fFaces, j)
{
label nbrI = fFaces[j];
scalar nbrQual = s[nbrI].tri(points()).quality();
if (nbrQual > qual)
{
qual = nbrQual;
normal[i] = s[nbrI].normal(points());
}
}
}
else
{
normal[i] = s[faceI].normal(points());
}
normal[i] /= mag(normal[i]);
}
else
{
// Set to what?
normal[i] = vector::zero;
}
}
}
else
{
forAll(info, i)
{
if (info[i].hit())
{
label faceI = info[i].index();
//- Cached:
//normal[i] = faceNormals()[faceI];
//- Uncached
normal[i] = triSurface::operator[](faceI).normal(points());
normal[i] /= mag(normal[i]) + VSMALL;
}
else
{
// Set to what?
normal[i] = vector::zero;
}
}
}
}

7
src/meshTools/searchableSurface/triSurfaceMesh.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2004-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -31,6 +31,7 @@ Description
- scale : scaling factor.
- tolerance : relative tolerance for doing intersections
(see triangle::intersection)
- minQuality: discard triangles with low quality when getting normal
SourceFiles
triSurfaceMesh.C
@ -70,6 +71,10 @@ private:
//- Optional tolerance to use in searches
scalar tolerance_;
//- Optional min triangle quality. Triangles below this get
// ignored for normal calculation
scalar minQuality_;
//- Optional max tree depth of octree
label maxTreeDepth_;

17
src/postProcessing/functionObjects/IO/partialWrite/partialWrite.C
View File

@ -27,6 +27,7 @@ License
#include "dictionary.H"
#include "Time.H"
#include "IOobjectList.H"
#include "polyMesh.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
@ -123,14 +124,26 @@ void Foam::partialWrite::write()
else
{
// Delete all but marked objects
fileName dbDir;
if (isA<polyMesh>(obr_))
{
dbDir = dynamic_cast<const polyMesh&>(obr_).dbDir();
}
IOobjectList objects(obr_, obr_.time().timeName());
forAllConstIter(HashPtrTable<IOobject>, objects, iter)
{
if (!objectNames_.found(iter()->name()))
{
const fileName f = obr_.time().timePath()/iter()->name();
//Pout<< " rm " << f << endl;
const fileName f =
obr_.time().timePath()
/dbDir
/iter()->name();
if (debug)
{
Pout<< " rm " << f << endl;
}
rm(f);
}
}

2
src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C
View File

@ -203,7 +203,7 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs()
// calculate the forces on the motion object from this data, then
// update the positions
motion_.updatePosition(t.deltaTValue());
motion_.updatePosition(t.deltaTValue(), t.deltaT0Value());
dictionary forcesDict;

92
src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -174,6 +174,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion()
initialQ_(I),
momentOfInertia_(diagTensor::one*VSMALL),
mass_(VSMALL),
cDamp_(0.0),
aLim_(VGREAT),
report_(false)
{}
@ -190,6 +192,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion
const point& initialCentreOfMass,
const tensor& initialQ,
const diagTensor& momentOfInertia,
scalar cDamp,
scalar aLim,
bool report
)
:
@ -211,6 +215,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion
initialQ_(initialQ),
momentOfInertia_(momentOfInertia),
mass_(mass),
cDamp_(cDamp),
aLim_(aLim),
report_(report)
{}
@ -233,6 +239,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion(const dictionary& dict)
),
momentOfInertia_(dict.lookup("momentOfInertia")),
mass_(readScalar(dict.lookup("mass"))),
cDamp_(dict.lookupOrDefault<scalar>("accelerationDampingCoeff", 0.0)),
aLim_(dict.lookupOrDefault<scalar>("accelerationLimit", VGREAT)),
report_(dict.lookupOrDefault<Switch>("report", false))
{
addRestraints(dict);
@ -246,17 +254,19 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion
const sixDoFRigidBodyMotion& sDoFRBM
)
:
motionState_(sDoFRBM.motionState()),
restraints_(sDoFRBM.restraints()),
restraintNames_(sDoFRBM.restraintNames()),
constraints_(sDoFRBM.constraints()),
constraintNames_(sDoFRBM.constraintNames()),
maxConstraintIterations_(sDoFRBM.maxConstraintIterations()),
initialCentreOfMass_(sDoFRBM.initialCentreOfMass()),
initialQ_(sDoFRBM.initialQ()),
momentOfInertia_(sDoFRBM.momentOfInertia()),
mass_(sDoFRBM.mass()),
report_(sDoFRBM.report())
motionState_(sDoFRBM.motionState_),
restraints_(sDoFRBM.restraints_),
restraintNames_(sDoFRBM.restraintNames_),
constraints_(sDoFRBM.constraints_),
constraintNames_(sDoFRBM.constraintNames_),
maxConstraintIterations_(sDoFRBM.maxConstraintIterations_),
initialCentreOfMass_(sDoFRBM.initialCentreOfMass_),
initialQ_(sDoFRBM.initialQ_),
momentOfInertia_(sDoFRBM.momentOfInertia_),
mass_(sDoFRBM.mass_),
cDamp_(sDoFRBM.cDamp_),
aLim_(sDoFRBM.aLim_),
report_(sDoFRBM.report_)
{}
@ -358,7 +368,8 @@ void Foam::sixDoFRigidBodyMotion::addConstraints
void Foam::sixDoFRigidBodyMotion::updatePosition
(
scalar deltaT
scalar deltaT,
scalar deltaT0
)
{
// First leapfrog velocity adjust and motion part, required before
@ -366,9 +377,32 @@ void Foam::sixDoFRigidBodyMotion::updatePosition
if (Pstream::master())
{
v() += 0.5*deltaT*a();
vector aClip = a();
scalar aMag = mag(aClip);
pi() += 0.5*deltaT*tau();
if (aMag > SMALL)
{
aClip /= aMag;
}
if (aMag > aLim_)
{
WarningIn
(
"void Foam::sixDoFRigidBodyMotion::updatePosition"
"("
"scalar deltaT, "
"scalar deltaT0"
")"
)
<< "Limited acceleration " << a()
<< " to " << aClip*aLim_
<< endl;
}
v() += 0.5*(1 - cDamp_)*deltaT0*aClip*min(aMag, aLim_);
pi() += 0.5*(1 - cDamp_)*deltaT0*tau();
// Leapfrog move part
centreOfMass() += deltaT*v();
@ -401,9 +435,33 @@ void Foam::sixDoFRigidBodyMotion::updateForce
applyConstraints(deltaT);
v() += 0.5*deltaT*a();
vector aClip = a();
scalar aMag = mag(aClip);
pi() += 0.5*deltaT*tau();
if (aMag > SMALL)
{
aClip /= aMag;
}
if (aMag > aLim_)
{
WarningIn
(
"void Foam::sixDoFRigidBodyMotion::updateForce"
"("
"const vector& fGlobal, "
"const vector& tauGlobal, "
"scalar deltaT"
")"
)
<< "Limited acceleration " << a()
<< " to " << aClip*aLim_
<< endl;
}
v() += 0.5*(1 - cDamp_)*deltaT*aClip*min(aMag, aLim_);
pi() += 0.5*(1 - cDamp_)*deltaT*tau();
if(report_)
{

19
src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -116,6 +116,15 @@ class sixDoFRigidBodyMotion
//- Mass of the body
scalar mass_;
//- Acceleration damping coefficient. Modify applied acceleration:
// v1 = v0 + a*dt - cDamp*a*dt
// = v0 + dt*f*(1 - cDamp)/m
// Increases effective mass by 1/(1 - cDamp).
scalar cDamp_;
//- Acceleration magnitude limit - clips large accelerations
scalar aLim_;
//- Switch to turn reporting of motion data on and off
Switch report_;
@ -235,6 +244,8 @@ public:
const point& initialCentreOfMass,
const tensor& initialQ,
const diagTensor& momentOfInertia,
scalar cDamp = 0.0,
scalar aLim = VGREAT,
bool report = false
);
@ -260,10 +271,12 @@ public:
void addConstraints(const dictionary& dict);
//- First leapfrog velocity adjust and motion part, required
// before force calculation
// before force calculation. Takes old timestep for variable
// timestep cases.
void updatePosition
(
scalar deltaT
scalar deltaT,
scalar deltaT0
);
//- Second leapfrog velocity adjust part, required after motion and

6
src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionIO.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2009-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -40,6 +40,10 @@ void Foam::sixDoFRigidBodyMotion::write(Ostream& os) const
<< momentOfInertia_ << token::END_STATEMENT << nl;
os.writeKeyword("mass")
<< mass_ << token::END_STATEMENT << nl;
os.writeKeyword("accelerationDampingCoeff")
<< cDamp_ << token::END_STATEMENT << nl;
os.writeKeyword("accelerationLimit")
<< aLim_ << token::END_STATEMENT << nl;
os.writeKeyword("report")
<< report_ << token::END_STATEMENT << nl;

9
src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/linearAxialAngularSpring/linearAxialAngularSpring.C
View File

@ -83,15 +83,20 @@ Foam::sixDoFRigidBodyMotionRestraints::linearAxialAngularSpring::restrain
) const
{
vector refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 1, 0);
vector oldDir = refQ_ & refDir;
vector newDir = motion.orientation() & refDir;
if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95)
{
// Directions getting close to the axis, change reference
refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 0, 1);
oldDir = refQ_ & refDir;
newDir = motion.orientation() & refDir;
vector oldDir = refQ_ & refDir;
vector newDir = motion.orientation() & refDir;
}
// Removing any axis component from oldDir and newDir and normalising

6
src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C
View File

@ -87,6 +87,7 @@ Foam::sixDoFRigidBodyMotionRestraints::tabulatedAxialAngularSpring::restrain
vector refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 1, 0);
vector oldDir = refQ_ & refDir;
vector newDir = motion.orientation() & refDir;
if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95)
@ -95,8 +96,9 @@ Foam::sixDoFRigidBodyMotionRestraints::tabulatedAxialAngularSpring::restrain
refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 0, 1);
oldDir = refQ_ & refDir;
newDir = motion.orientation() & refDir;
vector oldDir = refQ_ & refDir;
vector newDir = motion.orientation() & refDir;
}
// Removing any axis component from oldDir and newDir and normalising

4
src/postProcessing/functionObjects/forces/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2004-2010 OpenCFD Ltd.
\\ / A nd | Copyright (C) 2009-2011 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -147,7 +147,7 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs()
const polyMesh& mesh = this->dimensionedInternalField().mesh()();
const Time& t = mesh.time();
motion_.updatePosition(t.deltaTValue());
motion_.updatePosition(t.deltaTValue(), t.deltaT0Value());
vector gravity = vector::zero;