Commit 6e8f0dbe authored by sergio's avatar sergio Committed by Andrew Heather
Browse files

INT: org integration

1) rPolynomial Eq of State
2) externalForce and softWall in rigidBodyDynamics

INT: Several minor bug fixes plus
parent 5c74b709
......@@ -69,7 +69,7 @@ int main(int argc, char *argv[])
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
fileName vtkPath(runTime.path()/"VTK");
const fileName vtkPath(runTime.rootPath()/runTime.globalCaseName()/"VTK");
mkDir(vtkPath);
Info<< "Scanning times to determine track data for cloud " << cloudName
......
......@@ -137,7 +137,7 @@ int main(int argc, char *argv[])
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
fileName vtkPath(runTime.path()/"VTK");
const fileName vtkPath(runTime.rootPath()/runTime.globalCaseName()/"VTK");
mkDir(vtkPath);
forAll(timeDirs, timeI)
......@@ -145,7 +145,7 @@ int main(int argc, char *argv[])
runTime.setTime(timeDirs[timeI], timeI);
Info<< "Time = " << runTime.timeName() << endl;
fileName vtkTimePath(runTime.path()/"VTK"/runTime.timeName());
const fileName vtkTimePath(vtkPath/runTime.timeName());
mkDir(vtkTimePath);
Info<< " Reading particle positions" << endl;
......
......@@ -42,12 +42,13 @@ namespace Foam
void Foam::porosityModel::adjustNegativeResistance(dimensionedVector& resist)
{
scalar maxCmpt = max(0, cmptMax(resist.value()));
scalar maxCmpt = cmptMax(resist.value());
if (maxCmpt < 0)
{
FatalErrorInFunction
<< "Negative resistances are invalid, resistance = " << resist
<< "Cannot have all resistances set to negative, resistance = "
<< resist
<< exit(FatalError);
}
else
......
......@@ -202,6 +202,13 @@ void Foam::MPPICCloud<CloudType>::motion
if (dampingModel_->active())
{
if (this->mesh().moving())
{
FatalErrorInFunction
<< "MPPIC damping modelling does not support moving meshes."
<< exit(FatalError);
}
// update averages
td.updateAverages(cloud);
......@@ -226,6 +233,13 @@ void Foam::MPPICCloud<CloudType>::motion
if (packingModel_->active())
{
if (this->mesh().moving())
{
FatalErrorInFunction
<< "MPPIC packing modelling does not support moving meshes."
<< exit(FatalError);
}
// same procedure as for damping
td.updateAverages(cloud);
packingModel_->cacheFields(true);
......
......@@ -109,7 +109,7 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::segregated::K() const
L.primitiveFieldRef() = cbrt(mesh.V());
L.correctBoundaryConditions();
volScalarField I
const volScalarField I
(
alpha1
/max
......@@ -118,7 +118,7 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::segregated::K() const
pair_.phase1().residualAlpha() + pair_.phase2().residualAlpha()
)
);
volScalarField magGradI
const volScalarField magGradI
(
max
(
......@@ -127,28 +127,36 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::segregated::K() const
)
);
volScalarField muI
const volScalarField muI
(
rho1*nu1*rho2*nu2
/(rho1*nu1 + rho2*nu2)
);
volScalarField muAlphaI
const volScalarField limitedAlpha1
(
max(alpha1, pair_.phase1().residualAlpha())
);
const volScalarField limitedAlpha2
(
max(alpha2, pair_.phase2().residualAlpha())
);
const volScalarField muAlphaI
(
alpha1*rho1*nu1*alpha2*rho2*nu2
/(
max(alpha1, pair_.phase1().residualAlpha())*rho1*nu1
+ max(alpha2, pair_.phase2().residualAlpha())*rho2*nu2
)
/(limitedAlpha1*rho1*nu1 + limitedAlpha2*rho2*nu2)
);
volScalarField ReI
const volScalarField ReI
(
pair_.rho()
*pair_.magUr()
/(magGradI*muI)
/(magGradI*limitedAlpha1*limitedAlpha2*muI)
);
volScalarField lambda(m_*ReI + n_*muAlphaI/muI);
const volScalarField lambda(m_*ReI + n_*muAlphaI/muI);
return lambda*sqr(magGradI)*muI;
}
......
......@@ -33,6 +33,7 @@ License
#include "specie.H"
#include "perfectGas.H"
#include "rPolynomial.H"
#include "perfectFluid.H"
#include "rhoConst.H"
......@@ -84,6 +85,19 @@ constTransport
>
> constRefFluidHThermoPhysics;
typedef
constTransport
<
species::thermo
<
hRefConstThermo
<
rPolynomial<specie>
>,
sensibleEnthalpy
>
> constRefrPolFluidHThermoPhysics;
typedef
constTransport
<
......@@ -110,6 +124,19 @@ constTransport
>
> constRefFluidEThermoPhysics;
typedef
constTransport
<
species::thermo
<
eRefConstThermo
<
rPolynomial<specie>
>,
sensibleInternalEnergy
>
> constRefrPolFluidEThermoPhysics;
typedef
constTransport
<
......@@ -151,6 +178,18 @@ makeThermos
specie
);
makeThermos
(
rhoThermo,
heRhoThermo,
pureMixture,
constTransport,
sensibleEnthalpy,
hRefConstThermo,
rPolynomial,
specie
);
makeThermos
(
rhoThermo,
......@@ -190,6 +229,18 @@ makeThermos
specie
);
makeThermos
(
rhoThermo,
heRhoThermo,
pureMixture,
constTransport,
sensibleInternalEnergy,
eRefConstThermo,
rPolynomial,
specie
);
makeThermos
(
rhoThermo,
......@@ -235,6 +286,15 @@ makeThermoPhysicsReactionThermos
constRefFluidEThermoPhysics
);
makeThermoPhysicsReactionThermos
(
rhoThermo,
rhoReactionThermo,
heRhoThermo,
multiComponentMixture,
constRefrPolFluidEThermoPhysics
);
makeThermoPhysicsReactionThermos
(
rhoThermo,
......@@ -265,6 +325,15 @@ makeThermoPhysicsReactionThermos
constRefFluidHThermoPhysics
);
makeThermoPhysicsReactionThermos
(
rhoThermo,
rhoReactionThermo,
heRhoThermo,
multiComponentMixture,
constRefrPolFluidHThermoPhysics
);
makeThermoPhysicsReactionThermos
(
rhoThermo,
......
......@@ -33,6 +33,7 @@ restraints/linearDamper/linearDamper.C
restraints/linearAxialAngularSpring/linearAxialAngularSpring.C
restraints/sphericalAngularDamper/sphericalAngularDamper.C
restraints/prescribedRotation/prescribedRotation.C
restraints/externalForce/externalForce.C
restraints/softWall/softWall.C
rigidBodyModel/rigidBodyModel.C
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2019 OpenFOAM Foundation
-------------------------------------------------------------------------------
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 "externalForce.H"
#include "rigidBodyModel.H"
#include "rigidBodyModelState.H"
#include "OneConstant.H"
#include "addToRunTimeSelectionTable.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
namespace RBD
{
namespace restraints
{
defineTypeNameAndDebug(externalForce, 0);
addToRunTimeSelectionTable
(
restraint,
externalForce,
dictionary
);
}
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::RBD::restraints::externalForce::externalForce
(
const word& name,
const dictionary& dict,
const rigidBodyModel& model
)
:
restraint(name, dict, model),
externalForce_(nullptr)
{
read(dict);
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::RBD::restraints::externalForce::~externalForce()
{}
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void Foam::RBD::restraints::externalForce::restrain
(
scalarField& tau,
Field<spatialVector>& fx,
const rigidBodyModelState& state
) const
{
const vector force = externalForce_().value(state.t());
const vector moment(location_ ^ force);
if (model_.debug)
{
Info<< " location " << location_
<< " force " << force
<< " moment " << moment
<< endl;
}
// Accumulate the force for the restrained body
fx[bodyIndex_] += spatialVector(moment, force);
}
bool Foam::RBD::restraints::externalForce::read
(
const dictionary& dict
)
{
restraint::read(dict);
coeffs_.lookup("location") >> location_;
externalForce_ = Function1<vector>::New("force", coeffs_);
return true;
}
void Foam::RBD::restraints::externalForce::write
(
Ostream& os
) const
{
restraint::write(os);
os.writeEntry("location", location_);
externalForce_().writeData(os);
}
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | www.openfoam.com
\\/ M anipulation |
-------------------------------------------------------------------------------
Copyright (C) 2019 OpenFOAM Foundation
-------------------------------------------------------------------------------
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::RBD::restraints::externalForce
Group
grpRigidBodyDynamicsRestraints
Description
Time-dependent external force restraint using Function1.
Usage
Example applying a constant force to the floatingObject:
\verbatim
restraints
{
force
{
type externalForce;
body floatingObject;
location (0 0 0);
force (100 0 0);
}
}
\endverbatim
SourceFiles
externalForce.C
\*---------------------------------------------------------------------------*/
#ifndef RBD_restraints_externalForce_H
#define RBD_restraints_externalForce_H
#include "rigidBodyRestraint.H"
#include "Function1.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace RBD
{
namespace restraints
{
/*---------------------------------------------------------------------------*\
Class externalForce Declaration
\*---------------------------------------------------------------------------*/
class externalForce
:
public restraint
{
// Private data
//- External force (N)
autoPtr<Function1<vector>> externalForce_;
//- Point of application of the force
vector location_;
public:
//- Runtime type information
TypeName("externalForce");
// Constructors
//- Construct from components
externalForce
(
const word& name,
const dictionary& dict,
const rigidBodyModel& model
);
//- Construct and return a clone
virtual autoPtr<restraint> clone() const
{
return autoPtr<restraint>
(
new externalForce(*this)
);
}
//- Destructor
virtual ~externalForce();
// Member Functions
//- Accumulate the retraint internal joint forces into the tau field and
// external forces into the fx field
virtual void restrain
(
scalarField& tau,
Field<spatialVector>& fx,
const rigidBodyModelState& state
) const;
//- Update properties from given dictionary
virtual bool read(const dictionary& dict);
//- Write
virtual void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace restraints
} // End namespace RBD
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
......@@ -77,7 +77,8 @@ Foam::RBD::restraints::linearAxialAngularSpring::~linearAxialAngularSpring()
void Foam::RBD::restraints::linearAxialAngularSpring::restrain
(
scalarField& tau,
Field<spatialVector>& fx
Field<spatialVector>& fx,
const rigidBodyModelState& state
) const
{
vector refDir = rotationTensor(vector(1, 0, 0), axis_) & vector(0, 1, 0);
......
......@@ -111,7 +111,8 @@ public:
virtual void restrain
(
scalarField& tau,
Field<spatialVector>& fx
Field<spatialVector>& fx,
const rigidBodyModelState& state
) const;
//- Update properties from given dictionary
......
......@@ -76,7 +76,8 @@ Foam::RBD::restraints::linearDamper::~linearDamper()
void Foam::RBD::restraints::linearDamper::restrain
(
scalarField& tau,
Field<spatialVector>& fx
Field<spatialVector>& fx,
const rigidBodyModelState& state
) const
{
vector force = -coeff_*model_.v(model_.master(bodyID_)).l();
......
......@@ -102,7 +102,8 @@ public:
virtual void restrain
(
scalarField& tau,
Field<spatialVector>& fx
Field<spatialVector>& fx,
const rigidBodyModelState& state
) const;
//- Update properties from given dictionary
......
......@@ -76,7 +76,8 @@ Foam::RBD::restraints::linearSpring::~linearSpring()
void Foam::RBD::restraints::linearSpring::restrain
(
scalarField& tau,
Field<spatialVector>& fx