Merge branch 'develop' of develop.openfoam.com:Development/OpenFOAM-plus into develop

This commit is contained in:
Andrew Heather
2017-06-26 13:33:19 +01:00
72 changed files with 7072 additions and 53 deletions

View File

@ -86,9 +86,25 @@ if (mesh.changing())
fvm::laplacian(rAUf, pcorr) == fvc::div(phi)
);
//pcorrEqn.setReference(refCellI2, 0.0, true);
scalarList values(nZones, 0.0);
pcorrEqn.setReferences(refCells, values, true);
// Only set reference for cells that are CALCULATED
{
DynamicList<label> validCells(refCells.size());
forAll(refCells, zoneId)
{
if (refCells[zoneId] != -1)
{
validCells.append(refCells[zoneId]);
}
}
pcorrEqn.setReferences
(
validCells,
scalarList(validCells.size(), 0.0),
true
);
}
const dictionary& d = mesh.solver
(
@ -97,6 +113,7 @@ if (mesh.changing())
pimple.finalInnerIter()
)
);
//Bypass virtual layer
mesh.fvMesh::solve(pcorrEqn, d);
if (pimple.finalNonOrthogonalIter())

View File

@ -96,9 +96,25 @@
fvm::laplacian(rAUf, pcorr) == fvc::div(phi)
);
//pcorrEqn.setReference(refCellI2, 0, true);
scalarList values(nZones, 0.0);
pcorrEqn.setReferences(refCells, values, true);
// Only set reference for cells that are CALCULATED
{
DynamicList<label> validCells(refCells.size());
forAll(refCells, zoneId)
{
if (refCells[zoneId] != -1)
{
validCells.append(refCells[zoneId]);
}
}
pcorrEqn.setReferences
(
validCells,
scalarList(validCells.size(), 0.0),
true
);
}
const dictionary& d = mesh.solver
(
@ -117,9 +133,9 @@
}
}
if (runTime.outputTime())
{
volScalarField("contPhiPcorr", fvc::div(phi)).write();
pcorr.write();
}
//if (runTime.outputTime())
//{
// volScalarField("contPhiPcorr", fvc::div(phi)).write();
// pcorr.write();
//}
}

View File

@ -0,0 +1,3 @@
Test-externalCoupler.C
EXE = $(FOAM_USER_APPBIN)/Test-externalCoupler

View File

@ -0,0 +1,7 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/lumpedPointMotion/lnInclude
EXE_LIBS = \
-lfiniteVolume \
-llumpedPointMotion

View File

@ -0,0 +1,97 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 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/>.
Application
Test-externalCoupler
Description
Test of master/slave communication etc.
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "externalCoupler.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Main program:
int main(int argc, char *argv[])
{
argList::noParallel();
argList::addOption("max", "N", "max number of calls (default: 1000)");
argList::addBoolOption("slave", "run as slave");
#include "setRootCase.H"
const label maxCount = args.optionLookupOrDefault<label>("max", 1000);
externalCoupler coupler;
if (args.optionFound("slave"))
{
const word role = "slave";
Info<< "Running as " << role << " max=" << maxCount << endl;
for (label count = 0; count < maxCount; ++count)
{
// Wait for master, but stop if status=done was seen
Info<< role << ": waiting for master" << endl;
if (!coupler.waitForMaster())
{
Info<< role << ": stopping. status=done was detected" << endl;
break;
}
Info<< role << ": switch to master" << endl;
coupler.useMaster();
}
}
else
{
const word role = "master";
Info<< "Running as " << role << " max=" << maxCount << endl;
for (label count = 0; count < maxCount; ++count)
{
// Wait for slave, but stop if status=done was seen
Info<< role << ": waiting for slave" << endl;
if (!coupler.waitForSlave())
{
Info<< role << ": stopping. status=done was detected" << endl;
break;
}
Info<< role << ": switch to slave" << endl;
coupler.useSlave();
}
}
return 0;
}
// ************************************************************************* //

View File

@ -0,0 +1,3 @@
lumpedPointForces.C
EXE = $(FOAM_APPBIN)/lumpedPointForces

View File

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/lumpedPointMotion/lnInclude
EXE_LIBS = \
-lmeshTools \
-lfiniteVolume \
-llumpedPointMotion

View File

@ -0,0 +1,177 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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/>.
Application
lumpedPointForces
Description
Extract force/moment information from existing calculations based
on the segmentation of the pressure integration zones.
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "Time.H"
#include "timeSelector.H"
#include "volFields.H"
#include "IOobjectList.H"
#include "lumpedPointTools.H"
#include "lumpedPointIOMovement.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
template<class GeoFieldType>
autoPtr<GeoFieldType> loadField
(
const volMesh::Mesh& mesh,
const IOobject* io
)
{
if (io && io->headerClassName() == GeoFieldType::typeName)
{
Info<< "Reading " << GeoFieldType::typeName
<< ' ' << io->name() << endl;
return autoPtr<GeoFieldType>
(
new GeoFieldType
(
IOobject
(
io->name(),
io->instance(),
io->local(),
io->db(),
IOobject::MUST_READ,
IOobject::AUTO_WRITE,
io->registerObject()
),
mesh
)
);
}
return autoPtr<GeoFieldType>();
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addNote
(
"Extract force/moment information from existing calculations based "
"on the lumpedPoints pressure zones."
);
argList::addBoolOption
(
"vtk",
"create visualization files of the forces"
);
timeSelector::addOptions(true, false);
argList::noFunctionObjects(); // Never use function objects
#include "addRegionOption.H"
#include "setRootCase.H"
const bool withVTK = args.optionFound("vtk");
#include "createTime.H"
instantList timeDirs = timeSelector::select0(runTime, args);
#include "createNamedMesh.H"
autoPtr<lumpedPointIOMovement> movement = lumpedPointIOMovement::New
(
runTime
);
if (!movement.valid())
{
Info<< "no valid movement given" << endl;
return 1;
}
const labelList patchLst = lumpedPointTools::lumpedPointPatchList(mesh);
if (patchLst.empty())
{
Info<< "no patch list found" << endl;
return 2;
}
movement().setMapping(mesh, patchLst, lumpedPointTools::points0Field(mesh));
List<vector> forces, moments;
forAll(timeDirs, timei)
{
runTime.setTime(timeDirs[timei], timei);
Info<< "Time = " << runTime.timeName() << endl;
if (mesh.readUpdate())
{
Info<< " Read new mesh" << nl;
}
// Search for list of objects for this time
IOobjectList objects(mesh, runTime.timeName());
// Pressure field
autoPtr<volScalarField> pField
= loadField<volScalarField>(mesh, objects.lookup("p"));
// The forces per zone
if (movement().forcesAndMoments(mesh, forces, moments))
{
Info<<"forces per zone: " << forces << endl;
Info<<"moments per zone: " << moments << endl;
if (withVTK && Pstream::master())
{
const word outputName =
Foam::name("forces_%06d.vtp", runTime.timeIndex());
Info<<" " << outputName << endl;
movement().writeForcesAndMomentsVTP
(
outputName,
forces,
moments
);
}
}
}
Info<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -0,0 +1,3 @@
lumpedPointMovement.C
EXE = $(FOAM_APPBIN)/lumpedPointMovement

View File

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/lumpedPointMotion/lnInclude
EXE_LIBS = \
-lmeshTools \
-lfiniteVolume \
-llumpedPointMotion

View File

@ -0,0 +1,263 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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/>.
Application
lumpedPointMovement
Description
Thus utility can be used to produce VTK files to visualize the response
points/rotations and the corresponding movement of the building surfaces.
Uses the tabulated responses from the specified file.
Optionally, it can also be used to a dummy responder for the
externalCoupler logic, which makes it useful as a debugging facility
as well demonstrating how an external application could communicate
with the lumpedPointMovement point-patch boundary condition.
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "Time.H"
#include "timeSelector.H"
#include "OFstream.H"
#include "lumpedPointTools.H"
#include "lumpedPointIOMovement.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addNote
(
"Visualize lumpedPoint movements or provide a slave responder "
"for diagnostic purposes."
);
argList::noParallel();
argList::noFunctionObjects(); // Never use function objects
argList::addOption
(
"max",
"N",
"maximum number of outputs"
);
argList::addOption
(
"span",
"N",
"increment each input by factor N (default: 1)"
);
argList::addOption
(
"scale",
"factor",
"relaxation/scaling factor for movement (default: 1)"
);
argList::addBoolOption
(
"removeLock",
"remove lock-file on termination of slave"
);
argList::addBoolOption
(
"slave",
"invoke as a slave responder for testing"
);
argList::validArgs.insert("responseFile");
#include "setRootCase.H"
const label maxOut =
Foam::max(0, args.optionLookupOrDefault<label>("max", 0));
const label span =
Foam::max(1, args.optionLookupOrDefault<label>("span", 1));
const scalar relax = args.optionLookupOrDefault<scalar>("scale", 1);
const bool slave = args.optionFound("slave");
const bool removeLock = args.optionFound("removeLock");
#include "createTime.H"
autoPtr<lumpedPointIOMovement> movement = lumpedPointIOMovement::New
(
runTime
);
if (!movement.valid())
{
Info<< "no valid movement given" << endl;
return 1;
}
List<lumpedPointStateTuple> responseTable =
lumpedPointTools::lumpedPointStates(args[1]);
Info<< "Using response table with " << responseTable.size()
<< " entries" << endl;
Info << "Increment input by " << span << nl;
if (maxOut)
{
Info<< "Stopping after " << maxOut << " outputs" << endl;
}
if (slave)
{
Info<< "Running as slave responder" << endl;
externalCoupler& coupler = movement().coupler();
label count = 0;
for (label index = 0; index < responseTable.size(); index += span)
{
Info<< args.executable() << ": waiting for master" << endl;
// Wait for master, but stop if status=done was seen
if (!coupler.waitForMaster())
{
Info<< args.executable()
<< ": stopping status=done was detected" << endl;
break;
}
lumpedPointState state = responseTable[index].second();
state.relax(relax, movement().state0());
// Generate input for OpenFOAM
OFstream os(coupler.resolveFile(movement().inputName()));
if
(
movement().inputFormat()
== lumpedPointState::inputFormatType::PLAIN
)
{
state.writePlain(os);
}
else
{
os.writeEntry("time", responseTable[index].first());
state.writeDict(os);
}
Info<< args.executable()
<< ": updated to state " << index
<< " - switch to master"
<< endl;
// Let OpenFOAM know that it can continue
coupler.useMaster();
if (maxOut && ++count >= maxOut)
{
Info<< args.executable()
<<": stopping after " << maxOut << " outputs" << endl;
break;
}
}
if (removeLock)
{
Info<< args.executable() <<": removing lock file" << endl;
coupler.useSlave(); // This removes the lock-file
}
}
else
{
runTime.setTime(instant(0, runTime.constant()), 0);
#include "createNamedPolyMesh.H"
const labelList patchLst = lumpedPointTools::lumpedPointPatchList(mesh);
if (patchLst.empty())
{
Info<< "no patch list found" << endl;
return 2;
}
pointIOField points0 = lumpedPointTools::points0Field(mesh);
movement().setBoundBox(mesh, patchLst, points0);
label index = 0;
// Initial geometry
movement().writeVTP("geom_init.vtp", mesh, patchLst, points0);
forAll(responseTable, i)
{
const bool output = ((i % span) == 0);
lumpedPointState state = responseTable[i].second();
state.relax(relax, movement().state0());
if (output)
{
Info<<"output [" << i << "/"
<< responseTable.size() << "]" << endl;
}
else
{
continue;
}
// State/response = what comes back from FEM
{
const word outputName = Foam::name("state_%06d.vtp", index);
Info<<" " << outputName << endl;
state.writeVTP(outputName, movement().axis());
}
{
const word outputName = Foam::name("geom_%06d.vtp", index);
Info<<" " << outputName << endl;
movement().writeVTP(outputName, state, mesh, patchLst, points0);
}
{
++index;
bool canOutput = !maxOut || (index <= maxOut);
if (!canOutput)
{
Info<<"stopping output after "
<< maxOut << " outputs" << endl;
break;
}
}
}
}
Info<< args.executable() << ": End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -0,0 +1,3 @@
lumpedPointZones.C
EXE = $(FOAM_APPBIN)/lumpedPointZones

View File

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude \
-I$(LIB_SRC)/lumpedPointMotion/lnInclude
EXE_LIBS = \
-lmeshTools \
-lfiniteVolume \
-llumpedPointMotion

View File

@ -0,0 +1,107 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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/>.
Application
lumpedPointZones
Description
Produces a VTK PolyData file \c lumpedPointZones.vtp in which the
segmentation of the pressure integration zones can be visualized
for diagnostic purposes. Does not use external coupling.
\*---------------------------------------------------------------------------*/
#include "argList.H"
#include "Time.H"
#include "timeSelector.H"
#include "lumpedPointTools.H"
#include "lumpedPointIOMovement.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
argList::addNote
(
"Create lumpedPointZones.vtp to verify the segmentation of "
"pressure integration zones used by lumpedPoint BC."
);
argList::noParallel(); // The VTP writer is not yet in parallel
argList::noFunctionObjects(); // Never use function objects
argList::addBoolOption
(
"verbose",
"increased verbosity"
);
#include "addRegionOption.H"
#include "setRootCase.H"
// const bool verbose = args.optionFound("verbose");
#include "createTime.H"
runTime.setTime(instant(0, runTime.constant()), 0);
#include "createNamedPolyMesh.H"
autoPtr<lumpedPointIOMovement> movement = lumpedPointIOMovement::New
(
runTime
);
if (!movement.valid())
{
Info<< "no valid movement found" << endl;
return 1;
}
const labelList patchLst = lumpedPointTools::lumpedPointPatchList(mesh);
if (patchLst.empty())
{
Info<< "no patch list found" << endl;
return 2;
}
pointIOField points0 = lumpedPointTools::points0Field(mesh);
movement().setMapping(mesh, patchLst, points0);
// Initial geometry, but with zone colouring
movement().writeZonesVTP("lumpedPointZones.vtp", mesh, points0);
// Initial positions/rotations
movement().writeStateVTP("initialState.vtp");
Info<< nl
<< "wrote 'lumpedPointZones.vtp'" << nl
<< "wrote 'initialState.vtp'" << nl
<< "End\n" << endl;
return 0;
}
// ************************************************************************* //

View File

@ -42,8 +42,10 @@ options:
-block use blockMesh reader (uses .blockMesh extension)
-case <dir> specify alternative case directory, default is the cwd
-region <name> specify alternative mesh region
-touch only create the file (eg, .blockMesh, .OpenFOAM, etc)
-touchAll create .blockMesh, .OpenFOAM files (and for all regions)
-touch create the file (eg, .blockMesh, .OpenFOAM, .foam, ...)
-touch-all | -touchAll
create .blockMesh, .OpenFOAM files (and for all regions)
-touch-proc same as '-touch' but for each processor
-vtk | -builtin use VTK builtin OpenFOAM reader (uses .foam extension)
-help print the usage
--help paraview help
@ -102,11 +104,16 @@ do
unset plugin
shift
;;
-touchAll)
-touch-all | -touchAll)
optTouch=all
unset plugin
shift
;;
-touch-proc*)
optTouch=processor
unset plugin
shift
;;
--)
shift
break # Stop here, treat balance as paraview options
@ -205,6 +212,16 @@ all)
done
exit 0
;;
processor)
for i in processor*
do
(
cd $i 2> /dev/null && touch "${caseFile%.*}#${i#processor}.$extension"
)
done
echo "Created '$caseFile' for processor directories" 1>&2
exit 0
;;
true)
touch "$caseFile"
echo "Created '$caseFile'" 1>&2

View File

@ -76,6 +76,7 @@ wmake $targetType regionCoupled
functionObjects/Allwmake $targetType $*
wmake $targetType lumpedPointMotion
wmake $targetType sixDoFRigidBodyMotion
wmake $targetType rigidBodyDynamics
wmake $targetType rigidBodyMeshMotion

View File

@ -51,9 +51,20 @@ void Foam::functionObjects::timeControl::readControls()
{
timeEnd_ = time_.userTimeToTime(timeEnd_);
}
deltaTCoeff_ = dict_.lookupOrDefault("deltaTCoeff", GREAT);
dict_.readIfPresent("nStepsToStartTimeChange", nStepsToStartTimeChange_);
deltaTCoeff_ = GREAT;
if (dict_.readIfPresent("deltaTCoeff", deltaTCoeff_))
{
nStepsToStartTimeChange_ = labelMax;
}
else
{
nStepsToStartTimeChange_ = 3;
dict_.readIfPresent
(
"nStepsToStartTimeChange",
nStepsToStartTimeChange_
);
}
}
@ -73,7 +84,6 @@ Foam::scalar Foam::functionObjects::timeControl::calcExpansion
)
{
scalar ratio = startRatio;
// This function is used to calculate the 'expansion ratio' or
// 'time step growth factor'. Inputs:
// - ratio: wanted ratio
@ -226,7 +236,6 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
// has reached one that is possible. Returns the new expansionratio.
scalar y = timeToNextMultiple/wantedDT;
label requiredSteps = nSteps;
scalar ratioEstimate = deltaTCoeff_;
@ -240,6 +249,7 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
if (!rampDirectionUp)
{
ratioEstimate = 1/ratioEstimate;
ratioMax = 1/ratioMax;
requiredSteps *= -1;
}
// Provision for fall-back value if we can't satisfy requirements
@ -253,16 +263,18 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
y,
requiredSteps
);
const scalar deltaTLoop = wantedDT/pow(newRatio, requiredSteps-1);
const scalar deltaTLoop =
wantedDT
/ pow(newRatio, mag(requiredSteps)-1);
scalar firstDeltaRatio = deltaTLoop/deltaT0_;
if (debug)
{
// Avoid division by zero for ratio = 1.0
scalar Sn =
deltaTLoop
*(pow(newRatio, requiredSteps)-1)
*(pow(newRatio, mag(requiredSteps))-1)
/(newRatio-1+SMALL);
if (debug)
{
Info<< " nSteps " << requiredSteps
<< " ratioEstimate " << ratioEstimate
<< " a0 " << deltaTLoop
@ -281,7 +293,7 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
)
{
y += 1;
requiredSteps = nSteps;
requiredSteps = mag(nSteps);
if (debug)
{
Info<< "firstDeltaRatio " << firstDeltaRatio << " rampDir"
@ -291,7 +303,7 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
continue;
}
if (firstDeltaRatio > ratioMax)
if (firstDeltaRatio > ratioMax && rampDirectionUp)
{
requiredSteps++;
if (debug)
@ -302,7 +314,22 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
<< endl;
}
}
else if (newRatio > ratioMax)
else if (firstDeltaRatio < ratioMax && !rampDirectionUp)
{
requiredSteps--;
if (debug)
{
Info<< "First ratio " << firstDeltaRatio
<< " exceeds threshold " << ratioMax << nl
<< "Decreasing required steps " << requiredSteps
<< endl;
}
}
else if
(
(newRatio > ratioMax && rampDirectionUp)
|| (newRatio < ratioMax && !rampDirectionUp)
)
{
y += 1;
requiredSteps = nSteps;
@ -319,8 +346,19 @@ void Foam::functionObjects::timeControl::calcDeltaTCoeff
// Reset future series expansion at last step
if (requiredSteps == 1)
{
Sn = deltaT0_*firstDeltaRatio;
seriesDTCoeff_ = GREAT;
}
// Situations where we achieve wantedDT but fail to achieve
// multiple of writeInterval
scalar jumpError =
remainder(Sn + presentTime, wantedDT) / wantedDT;
if (mag(jumpError) > ROOTSMALL)
{
requiredSteps = label(timeToNextWrite/wantedDT);
firstDeltaRatio = timeToNextWrite/requiredSteps/deltaT0_;
}
if (debug)
{
Info<< "All conditions satisfied deltaT0_:" << deltaT0_
@ -365,10 +403,7 @@ Foam::functionObjects::timeControl::timeControl
dict_(dict),
timeStart_(-VGREAT),
timeEnd_(VGREAT),
nStepsToStartTimeChange_
(
dict.lookupOrDefault("nStepsToStartTimeChange", 3)
),
nStepsToStartTimeChange_(labelMax),
executeControl_(t, dict, "execute"),
writeControl_(t, dict, "write"),
foPtr_(functionObject::New(name, t, dict_)),
@ -494,7 +529,14 @@ bool Foam::functionObjects::timeControl::adjustTimeStep()
scalar nSteps = timeToNextWrite/deltaT;
// For tiny deltaT the label can overflow!
if (nSteps < labelMax)
if
(
nSteps < labelMax
&& (
deltaTCoeff_ != GREAT
|| nSteps < nStepsToStartTimeChange_
)
)
{
// nSteps can be < 1 so make sure at least 1
@ -517,6 +559,7 @@ bool Foam::functionObjects::timeControl::adjustTimeStep()
clipThreshold = 1/clipThreshold;
deltaT = max(newDeltaT, clipThreshold*deltaT);
}
const_cast<Time&>(time_).setDeltaT(deltaT, false);
}
}

View File

@ -0,0 +1,13 @@
externalCoupler.C
lumpedPointMovement.C
lumpedPointMovementWriter.C
lumpedPointState.C
lumpedPointStateWriter.C
lumpedPointIOMovement.C
lumpedPointDisplacementPointPatchVectorField.C
lumpedPointTools.C
LIB = $(FOAM_LIBBIN)/liblumpedPointMotion

View File

@ -0,0 +1,10 @@
EXE_INC = \
-I$(LIB_SRC)/finiteVolume/lnInclude \
-I$(LIB_SRC)/dynamicMesh/lnInclude \
-I$(LIB_SRC)/fileFormats/lnInclude \
-I$(LIB_SRC)/meshTools/lnInclude
LIB_LIBS = \
-lfiniteVolume \
-ldynamicMesh \
-lmeshTools

View File

@ -0,0 +1,309 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015-2017 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 "externalCoupler.H"
#include "Pstream.H"
#include "PstreamReduceOps.H"
#include "OSspecific.H"
#include "IFstream.H"
#include "OFstream.H"
#include "Switch.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(externalCoupler, 0);
}
Foam::word Foam::externalCoupler::lockName = "OpenFOAM";
// file-scope
// check file (must exist) for "status=done" content
static bool checkIsDone(const std::string& lck)
{
std::string content;
std::ifstream is(lck.c_str());
is >> content;
return (content.find("done") != std::string::npos);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
const Foam::fileName& Foam::externalCoupler::baseDir() const
{
return commsDir_;
}
Foam::fileName Foam::externalCoupler::lockFile() const
{
return resolveFile(lockName + ".lock");
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::externalCoupler::externalCoupler()
:
runState_(NONE),
commsDir_("$FOAM_CASE/comms"),
waitInterval_(1u),
timeOut_(100u),
slaveFirst_(false),
log(false)
{
commsDir_.expand();
commsDir_.clean();
}
Foam::externalCoupler::externalCoupler(const dictionary& dict)
:
externalCoupler()
{
readDict(dict);
if (Pstream::master())
{
mkDir(baseDir());
}
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::externalCoupler::~externalCoupler()
{
shutdown();
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::externalCoupler::readDict(const dictionary& dict)
{
// Normally cannot change directory or initialization
// if things have already been initialized
dict.lookup("commsDir") >> commsDir_;
commsDir_.expand();
commsDir_.clean();
slaveFirst_ = dict.lookupOrDefault<bool>("initByExternal", false);
waitInterval_ = dict.lookupOrDefault("waitInterval", 1u);
if (!waitInterval_)
{
// Enforce non-zero sleep
waitInterval_ = 1u;
}
timeOut_ = dict.lookupOrDefault("timeOut", 100*waitInterval_);
log = dict.lookupOrDefault<bool>("log", false);
return true;
}
bool Foam::externalCoupler::initialized() const
{
return runState_ != NONE;
}
bool Foam::externalCoupler::slaveFirst() const
{
return slaveFirst_;
}
void Foam::externalCoupler::useMaster(const bool wait) const
{
const bool wasInit = initialized();
runState_ = MASTER;
if (Pstream::master())
{
if (!wasInit)
{
// First time
Foam::mkDir(commsDir_);
}
const fileName lck(lockFile());
// Create lock file - only if it doesn't already exist
if (!Foam::isFile(lck))
{
Log << type() << ": creating lock file" << endl;
OFstream os(lck);
os << "status=openfoam\n";
os.flush();
}
}
if (wait)
{
waitForMaster();
}
}
void Foam::externalCoupler::useSlave(const bool wait) const
{
const bool wasInit = initialized();
runState_ = SLAVE;
if (Pstream::master())
{
if (!wasInit)
{
// First time
Foam::mkDir(commsDir_);
}
Log << type() << ": removing lock file" << endl;
Foam::rm(lockFile());
}
if (wait)
{
waitForSlave();
}
}
bool Foam::externalCoupler::waitForMaster() const
{
if (!initialized())
{
useMaster(); // was not initialized
}
bool isDone = false;
if (Pstream::master())
{
const fileName lck(lockFile());
double prevTime = -1;
double modTime = 0;
// Wait until file disappears (modTime == 0)
// But also check for status=done content in the file
while ((modTime = highResLastModified(lck)) > 0)
{
if (prevTime < modTime)
{
prevTime = modTime;
isDone = checkIsDone(lck);
if (isDone)
{
break;
}
}
sleep(waitInterval_);
}
}
// MPI barrier
Pstream::scatter(isDone);
return !isDone;
}
bool Foam::externalCoupler::waitForSlave() const
{
if (!initialized())
{
useSlave(); // was not initialized
}
bool isDone = false;
if (Pstream::master())
{
const fileName lck(lockFile());
unsigned totalTime = 0;
Log << type() << ": waiting for lock file to appear " << lck << endl;
while (!Foam::isFile(lck))
{
sleep(waitInterval_);
if (timeOut_ && (totalTime += waitInterval_) > timeOut_)
{
FatalErrorInFunction
<< "Wait time exceeded timeout of " << timeOut_
<< " s" << abort(FatalError);
}
Log << type() << ": wait time = " << totalTime << endl;
}
// But check for status=done content in the file
isDone = checkIsDone(lck);
Log << type() << ": found lock file " << lck << endl;
}
// MPI barrier
Pstream::scatter(isDone);
return !isDone;
}
Foam::fileName Foam::externalCoupler::resolveFile
(
const word& file
) const
{
return fileName(baseDir()/file);
}
void Foam::externalCoupler::shutdown() const
{
if (Pstream::master() && runState_ == MASTER && Foam::isDir(commsDir_))
{
const fileName lck(lockFile());
Log << type() << ": lock file status=done" << endl;
OFstream os(lck);
os << "status=done\n";
os.flush();
}
runState_ = DONE; // Avoid re-triggering in destructor
}
// ************************************************************************* //

View File

@ -0,0 +1,207 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015-2017 OpenCFD Ltd.
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify i
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::externalCoupler
Description
Encapsulates the logic for coordinating between OpenFOAM and an
external application.
This class provides a simple interface for explicit coupling with an
external application using plain text files located in the user-specified
communications directory.
These files are to be read/written on the master processor only.
The explicit coupling follows a master/slave model in which OpenFOAM
is the 'master' and the external application is the 'slave'.
The readiness to exchange information in either direction is handled
by a lock file (ie, OpenFOAM.lock).
If the lock file is present, the slave (external application) should wait
for the master (OpenFOAM) to complete.
When the master is finished its tasks and has prepared data for
the slave, the lock file is removed, instructing the external
source to take control of the program execution.
When the slave has completed its tasks, it will reinstate the lock file.
Example of the communication specification:
\verbatim
communication
{
commsDir "${FOAM_CASE}/comms";
waitInterval 1;
timeOut 100;
initByExternal no;
}
\endverbatim
SourceFiles
externalCoupler.C
\*---------------------------------------------------------------------------*/
#ifndef externalCoupler_H
#define externalCoupler_H
#include "fileName.H"
#include "dictionary.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class externalCoupler Declaration
\*---------------------------------------------------------------------------*/
class externalCoupler
{
//- The run state (ie, who is currently in charge)
enum runState
{
NONE, //!< Not initialized
MASTER, //!< The master (OpenFOAM) is in charge
SLAVE, //!< The slave (external program) is in charge
DONE //!< Finished
};
// Private data
//- The current run (and initialization) state
mutable runState runState_;
//- Local path to communications directory
fileName commsDir_;
//- Interval time between checking for return data [s]
unsigned waitInterval_;
//- Timeout [s] while waiting for the external application
unsigned timeOut_;
//- Flag to indicate values are initialized by external application
bool slaveFirst_;
//- Local logging/verbosity flag
bool log;
// Private Member Functions
//- Return the file path to the base communications directory
const fileName& baseDir() const;
//- Return the file path to the lock file
fileName lockFile() const;
//- Disallow default bitwise copy construc
externalCoupler(const externalCoupler&) = delete;
//- Disallow default bitwise assignmen
void operator=(const externalCoupler&) = delete;
public:
//- Runtime type information
TypeName("externalCoupler");
// Static data members
//- Name of the lock file
static word lockName;
// Constructors
//- Construct null using standard defaults
externalCoupler();
//- Construct from dictionary
externalCoupler(const dictionary& dict);
//- Destructor
virtual ~externalCoupler();
// Member Functions
// Access
//- True if state has been initialized
bool initialized() const;
//- External application provides initial values
bool slaveFirst() const;
//- Create lock file to indicate that OpenFOAM is in charge
// Optionally wait for master as well.
void useMaster(const bool wait=false) const;
//- Wait for indication that OpenFOAM has supplied output.
// This is when the lock file disappears, or it exists but with
// "status=done" content.
// \return False if lock file contains "status=done"
bool waitForMaster() const;
//- Remove lock file to indicate that the external program is in charge
// Optionally wait for slave as well.
void useSlave(const bool wait=false) const;
//- Wait for indication that the external program has supplied input.
// This is when the lock file appears.
// \return False if lock file contains "status=done"
bool waitForSlave() const;
//- Return the file path in the communications directory
fileName resolveFile(const word& file) const;
//- Generate status=done in lock (only when run-state = master)
void shutdown() const;
//- Remove files written by OpenFOAM
void removeDirectory() const;
// Edit
//- Read communication settings from dictionary
bool readDict(const dictionary& dict);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,285 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointDisplacementPointPatchVectorField.H"
#include "lumpedPointMovement.H"
#include "lumpedPointIOMovement.H"
#include "addToRunTimeSelectionTable.H"
#include "pointFields.H"
#include "surfaceFields.H"
#include "volFields.H"
#include "Time.H"
#include "polyMesh.H"
#include "displacementMotionSolver.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
makePointPatchTypeField
(
pointPatchVectorField,
lumpedPointDisplacementPointPatchVectorField
);
}
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
Foam::labelList
Foam::lumpedPointDisplacementPointPatchVectorField::patchIds
(
const pointVectorField& pvf
)
{
const pointVectorField::Boundary& bf = pvf.boundaryField();
DynamicList<label> patchLst(bf.size());
forAll(bf, patchi)
{
// All patches of this type
if (isA<patchType>(bf[patchi]))
{
patchLst.append(patchi);
// or patchLst.append(bf[patchi].patch().index());
}
}
return patchLst.shrink();
}
// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
const Foam::pointField&
Foam::lumpedPointDisplacementPointPatchVectorField::points0() const
{
const objectRegistry& obr = this->patch().boundaryMesh().mesh().db();
// Obtain starting locations from the motionSolver
return obr.lookupObject<displacementMotionSolver>
(
"dynamicMeshDict"
).points0();
}
const Foam::lumpedPointMovement&
Foam::lumpedPointDisplacementPointPatchVectorField::movement() const
{
const objectRegistry& obr = this->patch().boundaryMesh().mesh().db();
const lumpedPointIOMovement* ptr =
lumpedPointIOMovement::lookupInRegistry(obr);
if (ptr)
{
return *ptr; // Already exists
}
// create and register with this patch as the owner
autoPtr<lumpedPointIOMovement> obj = lumpedPointIOMovement::New
(
obr,
this->patch().index()
);
return objectRegistry::store(obj);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::lumpedPointDisplacementPointPatchVectorField::
lumpedPointDisplacementPointPatchVectorField
(
const pointPatch& p,
const DimensionedField<vector, pointMesh>& iF
)
:
fixedValuePointPatchField<vector>(p, iF)
{}
Foam::lumpedPointDisplacementPointPatchVectorField::
lumpedPointDisplacementPointPatchVectorField
(
const pointPatch& p,
const DimensionedField<vector, pointMesh>& iF,
const dictionary& dict
)
:
fixedValuePointPatchField<vector>(p, iF, dict)
{}
Foam::lumpedPointDisplacementPointPatchVectorField::
lumpedPointDisplacementPointPatchVectorField
(
const lumpedPointDisplacementPointPatchVectorField& pf,
const pointPatch& p,
const DimensionedField<vector, pointMesh>& iF,
const pointPatchFieldMapper& mapper
)
:
fixedValuePointPatchField<vector>(pf, p, iF, mapper)
{}
Foam::lumpedPointDisplacementPointPatchVectorField::
lumpedPointDisplacementPointPatchVectorField
(
const lumpedPointDisplacementPointPatchVectorField& pf,
const DimensionedField<vector, pointMesh>& iF
)
:
fixedValuePointPatchField<vector>(pf, iF)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::lumpedPointDisplacementPointPatchVectorField::
~lumpedPointDisplacementPointPatchVectorField()
{
// de-register movement if in use and managed by this patch
const lumpedPointIOMovement* ptr = lumpedPointIOMovement::lookupInRegistry
(
this->patch().boundaryMesh().mesh().db()
);
if (ptr && ptr->ownerId() == this->patch().index())
{
movement().coupler().shutdown();
const_cast<lumpedPointIOMovement*>(ptr)->checkOut();
}
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::lumpedPointDisplacementPointPatchVectorField::updateCoeffs()
{
if (this->updated())
{
return;
}
const bool masterPatch = (movement().ownerId() == this->patch().index());
if (masterPatch)
{
if (lumpedPointIOMovement::debug)
{
Pout<<"masterPatch: " << this->patch().index() << endl;
}
const polyMesh& mesh = this->patch().boundaryMesh().mesh().mesh();
// need face 'zones' for calculating forces
// likely need bounding box for the movement
// -> do both now if required
if (!movement().hasMapping())
{
const_cast<lumpedPointMovement&>(movement()).setMapping
(
mesh,
// All patches of this type
patchIds
(
static_cast<const pointVectorField&>
(
this->internalField()
)
),
this->points0()
);
}
if
(
movement().coupler().initialized()
|| !movement().coupler().slaveFirst()
)
{
// Synchronized for all processes
List<vector> forces, moments;
movement().forcesAndMoments(mesh, forces, moments);
if (lumpedPointIOMovement::debug)
{
Pout<<"gatherForces: " << forces << " called from patch "
<< this->patch().index() << endl;
if (Pstream::master())
{
Pout<<"output forces to file: "
<< movement().locations() << " called from patch "
<< this->patch().index() << nl
<<"# " << forces.size() << " force entries" << nl
<<"# fx fy fz" << nl
<<"output forces to file: "
<< forces << " called from patch "
<< this->patch().index() << endl;
}
}
if (Pstream::master())
{
movement().writeData(forces, moments);
// signal external source to execute
movement().coupler().useSlave();
}
}
// Wait for slave to provide data - includes MPI barrier
movement().coupler().waitForSlave();
// Read data passed back from external source - includes MPI barrier
const_cast<lumpedPointMovement&>(movement()).readState();
}
tmp<pointField> tdisp = movement().displacePoints
(
this->points0(),
this->patch().meshPoints()
);
this->operator==(tdisp);
fixedValuePointPatchField<vector>::updateCoeffs();
}
void Foam::lumpedPointDisplacementPointPatchVectorField::write(Ostream& os)
const
{
pointPatchField<vector>::write(os);
writeEntry("value", os);
}
// ************************************************************************* //

View File

@ -0,0 +1,200 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 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::lumpedPointDisplacementPointPatchVectorField
Description
Interpolates pre-specified motion with motion specified as
pointVectorFields.
This is the point-patch responsible for managing the force
integration on a 'lumped-point' basis, waiting for the external
application, reading back the response from the external program
and updating the locations of the associated patch points
accordingly.
The internal patch type name is 'lumpedPointDisplacement'.
\heading Patch usage
Example:
\verbatim
walls
{
type lumpedPointDisplacement;
value uniform (0 0 0);
fieldName wantedDisplacement;
interpolationScheme linear;
}
\endverbatim
This will scan the case for \a wantedDisplacement pointVectorFields
and interpolate those in time (using \c linear interpolation) to
obtain the current displacement.
The advantage of specifying displacement in this way is that it
automatically works through decomposePar.
SourceFiles
lumpedPointDisplacementPointPatchVectorField.C
\*---------------------------------------------------------------------------*/
#ifndef lumpedPointDisplacementPointPatchVectorField_H
#define lumpedPointDisplacementPointPatchVectorField_H
#include "fixedValuePointPatchField.H"
#include "lumpedPointMovement.H"
#include "lumpedPointState.H"
#include "lumpedPointIOMovement.H"
#include "labelList.H"
#include "tmp.H"
#include "pointField.H"
#include "pointFieldsFwd.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class interpolationWeights;
/*---------------------------------------------------------------------------*\
Class lumpedPointDisplacementPointPatchVectorField Declaration
\*---------------------------------------------------------------------------*/
class lumpedPointDisplacementPointPatchVectorField
:
public fixedValuePointPatchField<vector>
{
// Private data
//- Convenience typedefs
typedef lumpedPointDisplacementPointPatchVectorField patchType;
typedef DimensionedField<vector, pointMesh> fieldType;
// Private Member Functions
protected:
//- The starting locations (obtained from the motionSolver).
const pointField& points0() const;
//- The auto-vivifying singleton for movement.
const lumpedPointMovement& movement() const;
public:
//- Runtime type information
TypeName("lumpedPointDisplacement");
// Constructors
//- Construct from patch and internal field
lumpedPointDisplacementPointPatchVectorField
(
const pointPatch& p,
const DimensionedField<vector, pointMesh>& iF
);
//- Construct from patch, internal field and dictionary
lumpedPointDisplacementPointPatchVectorField
(
const pointPatch& p,
const DimensionedField<vector, pointMesh>& iF,
const dictionary& dict
);
//- Construct by mapping given patchField<vector> onto a new patch
lumpedPointDisplacementPointPatchVectorField
(
const lumpedPointDisplacementPointPatchVectorField& pf,
const pointPatch& p,
const DimensionedField<vector, pointMesh>& iF,
const pointPatchFieldMapper& mapper
);
//- Construct and return a clone
virtual autoPtr<pointPatchField<vector>> clone() const
{
return autoPtr<pointPatchField<vector>>
(
new lumpedPointDisplacementPointPatchVectorField
(
*this
)
);
}
//- Construct as copy setting internal field reference
lumpedPointDisplacementPointPatchVectorField
(
const lumpedPointDisplacementPointPatchVectorField& pf,
const DimensionedField<vector, pointMesh>& iF
);
//- Construct and return a clone setting internal field reference
virtual autoPtr<pointPatchField<vector>> clone
(
const DimensionedField<vector, pointMesh>& iF
) const
{
return autoPtr<pointPatchField<vector>>
(
new lumpedPointDisplacementPointPatchVectorField
(
*this,
iF
)
);
}
//- Destructor
virtual ~lumpedPointDisplacementPointPatchVectorField();
// Member functions
//- The ids for all patches of this type
static labelList patchIds(const pointVectorField& pvf);
//- Update the coefficients associated with the patch field
virtual void updateCoeffs();
//- Write
virtual void write(Ostream&) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,142 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointIOMovement.H"
#include "Time.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(lumpedPointIOMovement, 0);
}
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
const Foam::lumpedPointIOMovement*
Foam::lumpedPointIOMovement::lookupInRegistry(const objectRegistry& obr)
{
return obr.lookupObjectPtr<lumpedPointIOMovement>
(
lumpedPointMovement::dictionaryName
);
}
Foam::autoPtr<Foam::lumpedPointIOMovement>
Foam::lumpedPointIOMovement::New
(
const objectRegistry& obr,
label ownerId
)
{
return autoPtr<lumpedPointIOMovement>
(
new lumpedPointIOMovement
(
IOobject
(
lumpedPointMovement::dictionaryName,
obr.time().caseSystem(),
obr,
IOobject::MUST_READ,
IOobject::NO_WRITE,
true // register object
),
ownerId // tag this patch as owner too
)
);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::lumpedPointIOMovement::lumpedPointIOMovement
(
const IOobject& io,
label ownerId
)
:
lumpedPointMovement(),
regIOobject(io)
{
bool ok =
(
readOpt() == IOobject::MUST_READ
|| readOpt() == IOobject::MUST_READ_IF_MODIFIED
);
if (ok)
{
ok = readData(readStream(typeName));
close();
if (ok)
{
this->ownerId(ownerId);
}
}
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::lumpedPointIOMovement::~lumpedPointIOMovement()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
bool Foam::lumpedPointIOMovement::readData(Istream& is)
{
dictionary dict(is);
readDict(dict);
return is.check(FUNCTION_NAME);
}
bool Foam::lumpedPointIOMovement::writeData(Ostream& os) const
{
os << *this;
return os.good();
}
// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
Foam::Ostream& Foam::operator<<(Ostream& os, const lumpedPointIOMovement& obj)
{
obj.writeDict(os);
os.check(FUNCTION_NAME);
return os;
}
// ************************************************************************* //

View File

@ -0,0 +1,133 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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::lumpedPointIOMovement
Description
IO-registered version of lumpedPointMovement.
SourceFiles
lumpedPointMovement.C
\*---------------------------------------------------------------------------*/
#ifndef lumpedPointIOMovement_H
#define lumpedPointIOMovement_H
#include "lumpedPointMovement.H"
#include "regIOobject.H"
#include "className.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward declarations
class lumpedPointIOMovement;
Ostream& operator<<(Ostream& os, const lumpedPointIOMovement& obj);
/*---------------------------------------------------------------------------*\
Class lumpedPointIOMovement Declaration
\*---------------------------------------------------------------------------*/
class lumpedPointIOMovement
:
public lumpedPointMovement,
public regIOobject
{
// Private Member Functions
//- Disallow default bitwise copy construct
lumpedPointIOMovement(const lumpedPointIOMovement&) = delete;
//- Disallow default bitwise assignment
void operator=(const lumpedPointIOMovement&) = delete;
public:
//- Runtime type information
TypeName("lumpedPointMovement");
// Static Member Functions
//- Lookup pointer to object in the object-registry,
// return nullptr if found.
static const lumpedPointIOMovement* lookupInRegistry
(
const objectRegistry& obr
);
//- Create a new object in the registry by reading system dictionary
static autoPtr<lumpedPointIOMovement> New
(
const objectRegistry& obr,
label ownerId = -1
);
// Constructors
//- Construct from IOobject, optionally with some owner information
explicit lumpedPointIOMovement
(
const IOobject& io,
label ownerId = -1
);
//- Destructor
~lumpedPointIOMovement();
// Member Functions
//- readData member function used by regIOobject
bool readData(Istream& is);
//- writeData member function required by regIOobject
bool writeData(Ostream& os) const;
// IOstream Operators
friend Ostream& operator<<
(
Ostream& os,
const lumpedPointIOMovement& obj
);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,790 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointMovement.H"
#include "lumpedPointIOMovement.H"
#include "demandDrivenData.H"
#include "linearInterpolationWeights.H"
#include "IFstream.H"
#include "OFstream.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "PtrMap.H"
#include "faceZoneMesh.H"
#include "PstreamReduceOps.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::Enum<Foam::lumpedPointMovement::outputFormatType>
Foam::lumpedPointMovement::formatNames
{
{ outputFormatType::PLAIN, "plain" },
{ outputFormatType::DICTIONARY, "dictionary" }
};
const Foam::word
Foam::lumpedPointMovement::dictionaryName("lumpedPointMovement");
namespace Foam
{
//! \cond fileScope
//- Write list content with size, bracket, content, bracket one-per-line.
// This makes for consistent for parsing, regardless of the list length.
template <class T>
static void writeList(Ostream& os, const string& header, const UList<T>& lst)
{
// Header string
os << header.c_str() << nl;
// Write size and start delimiter
os << lst.size() << nl
<< token::BEGIN_LIST << nl;
// Write contents
forAll(lst, i)
{
os << lst[i] << nl;
}
// Write end delimiter
os << token::END_LIST << token::END_STATEMENT << nl << nl;
}
//! \endcond
}
// namespace Foam
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
Foam::IOobject Foam::lumpedPointMovement::selectIO
(
const IOobject& io,
const fileName& f
)
{
return
(
f.size()
? IOobject // construct from filePath instead
(
f,
io.db(),
io.readOpt(),
io.writeOpt(),
io.registerObject(),
io.globalObject()
)
: io
);
}
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::lumpedPointMovement::calcThresholds() const
{
thresholdPtr_ = new scalarField(locations_);
scalarField& thresh = *thresholdPtr_;
for (label i=1; i < thresh.size(); ++i)
{
thresh[i-1] =
(
locations_[i-1]
+ (division_ * (locations_[i] - locations_[i-1]))
);
}
}
Foam::label Foam::lumpedPointMovement::threshold(scalar pos) const
{
const scalarField& thresh = this->thresholds();
const label n = thresh.size();
// Clamp above and below
//
// ? may need special treatment to clamp values below the first threshold ?
// ? special treatment when 'axis' is negative ?
for (label i=0; i < n; ++i)
{
if (pos < thresh[i])
{
return i;
}
}
// everything above goes into the last zone too
return n-1;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::lumpedPointMovement::lumpedPointMovement()
:
axis_(0,0,1),
locations_(),
division_(0),
relax_(1),
interpolationScheme_(linearInterpolationWeights::typeName),
ownerId_(-1),
boundBox_(),
centre_(),
autoCentre_(true),
forcesDict_(),
coupler_(),
inputName_("positions.in"),
outputName_("forces.out"),
inputFormat_(lumpedPointState::inputFormatType::DICTIONARY),
outputFormat_(outputFormatType::DICTIONARY),
state0_(),
state_(),
thresholdPtr_(0),
interpolatorPtr_()
{}
Foam::lumpedPointMovement::lumpedPointMovement
(
const dictionary& dict,
label ownerId
)
:
axis_(0,0,1),
locations_(),
division_(0),
relax_(1),
interpolationScheme_
(
dict.lookupOrDefault<word>
(
"interpolationScheme",
linearInterpolationWeights::typeName
)
),
ownerId_(ownerId),
boundBox_(),
centre_(),
autoCentre_(true),
forcesDict_(),
coupler_(),
state0_(),
state_(),
thresholdPtr_(0),
interpolatorPtr_()
{
readDict(dict);
}
// * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * //
Foam::lumpedPointMovement::~lumpedPointMovement()
{
deleteDemandDrivenData(thresholdPtr_);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::lumpedPointMovement::readDict(const dictionary& dict)
{
// assume the worst
deleteDemandDrivenData(thresholdPtr_);
dict.lookup("axis") >> axis_;
division_ = 0;
if (dict.readIfPresent("division", division_))
{
if (division_ < 0)
{
division_ = 0;
}
else if (division_ > 1)
{
division_ = 1;
}
}
dict.readIfPresent("relax", relax_);
dict.lookup("locations") >> locations_;
if (dict.readIfPresent("interpolationScheme", interpolationScheme_))
{
interpolatorPtr_.clear();
}
autoCentre_ = !dict.readIfPresent("centre", centre_);
// Initial state from initial locations, with zero rotation
tmp<vectorField> pts = locations_ * axis_;
state0_ = lumpedPointState(pts);
state_ = state0_;
forcesDict_.merge(dict.subOrEmptyDict("forces"));
const dictionary& commDict = dict.subDict("communication");
coupler_.readDict(commDict);
// TODO: calcFrequency_ = dict.lookupOrDefault("calcFrequency", 1);
commDict.lookup("inputName") >> inputName_;
commDict.lookup("outputName") >> outputName_;
inputFormat_ = lumpedPointState::formatNames.lookup
(
"inputFormat",
commDict
);
outputFormat_ = lumpedPointMovement::formatNames.lookup
(
"outputFormat",
commDict
);
}
void Foam::lumpedPointMovement::setBoundBox
(
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
)
{
boundBox_ = boundBox::invertedBox;
const polyBoundaryMesh& patches = mesh.boundaryMesh();
for (const label patchId : patchIds)
{
boundBox_.add(points0, patches[patchId].meshPoints());
}
boundBox_.reduce();
if (autoCentre_)
{
centre_ = boundBox_.midpoint();
centre_ -= (centre_ & axis_) * axis_;
if (lumpedPointIOMovement::debug)
{
Pout<<"autoCentre on " << centre_ << " boundBox: " << boundBox_ << endl;
}
}
else
{
// User-specified centre
if (lumpedPointIOMovement::debug)
{
Pout<<"centre on " << centre_ << " boundBox: " << boundBox_ << endl;
}
}
}
void Foam::lumpedPointMovement::setMapping
(
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
)
{
setBoundBox(mesh, patchIds, points0);
faceZones_.clear();
// Local storage location (boundary faces only)
const label firstFace = mesh.nInternalFaces();
labelList faceToZoneID(mesh.nFaces() - firstFace, -1);
// Number of faces per zone
labelList nFaces(thresholds().size(), 0);
const polyBoundaryMesh& patches = mesh.boundaryMesh();
for (const label patchId : patchIds)
{
const polyPatch& pp = patches[patchId];
// Classify faces into respective zones according to their centres
label meshFaceI = pp.start();
forAll(pp, i)
{
const face& f = mesh.faces()[meshFaceI];
label zoneId = this->threshold(f.centre(points0));
// localize storage location
faceToZoneID[meshFaceI - firstFace] = zoneId;
++nFaces[zoneId];
++meshFaceI;
}
}
if (lumpedPointIOMovement::debug)
{
Pout<<"faces per zone:" << nFaces << endl;
}
faceZones_.setSize(nFaces.size());
label nZonedFaces = 0;
forAll(nFaces, zonei)
{
label n = nFaces[zonei];
labelList& addr = faceZones_[zonei];
addr.setSize(n);
n = 0;
forAll(faceToZoneID, faceI)
{
if (faceToZoneID[faceI] == zonei)
{
// from local storage to global mesh face
label meshFaceI = faceI + firstFace;
addr[n] = meshFaceI;
++n;
}
}
addr.setSize(n);
nZonedFaces += n;
if (lumpedPointIOMovement::debug)
{
Pout<< "Created faceZone[" << zonei << "] "
<< returnReduce(n, sumOp<label>()) << " faces, "
<< thresholds()[zonei] << " threshold" << nl;
}
}
}
bool Foam::lumpedPointMovement::forcesAndMoments
(
const polyMesh& pmesh,
List<vector>& forces,
List<vector>& moments
) const
{
forces.setSize(faceZones_.size());
moments.setSize(faceZones_.size());
if (faceZones_.empty())
{
WarningInFunction
<< "Attempted to calculate forces without setMapping()"
<< nl << endl;
forces.setSize(thresholds().size(), Zero);
moments.setSize(thresholds().size(), Zero);
return false;
}
// Initialize with zero
forces = Zero;
moments = Zero;
// The mass centres
const pointField& lumpedCentres = state().points();
const polyBoundaryMesh& patches = pmesh.boundaryMesh();
const word pName = forcesDict_.lookupOrDefault<word>("p", "p");
scalar pRef = forcesDict_.lookupOrDefault<scalar>("pRef", 0.0);
scalar rhoRef = forcesDict_.lookupOrDefault<scalar>("rhoRef", 1.0);
// Calculated force per patch - cache
PtrMap<vectorField> forceOnPatches;
const volScalarField* pPtr = pmesh.lookupObjectPtr<volScalarField>(pName);
// fvMesh and has pressure field
if (isA<fvMesh>(pmesh) && pPtr)
{
const fvMesh& mesh = dynamicCast<const fvMesh>(pmesh);
const volScalarField& p = *pPtr;
// face centres (on patches)
const surfaceVectorField::Boundary& patchCf =
mesh.Cf().boundaryField();
// face areas (on patches)
const surfaceVectorField::Boundary& patchSf =
mesh.Sf().boundaryField();
// Pressure (on patches)
const volScalarField::Boundary& patchPress =
p.boundaryField();
// rhoRef if the pressure field is dynamic, i.e. p/rho otherwise 1
rhoRef = (p.dimensions() == dimPressure ? 1.0 : rhoRef);
// Scale pRef by density for incompressible simulations
pRef /= rhoRef;
forAll(faceZones_, zonei)
{
const labelList& zn = faceZones_[zonei];
forAll(zn, localFaceI)
{
const label meshFaceI = zn[localFaceI];
// locate which patch this face is on,
// and its offset within the patch
const label patchI = patches.whichPatch(meshFaceI);
if (patchI == -1)
{
// Should not happen - could also warn
continue;
}
const label patchFaceI = meshFaceI - patches[patchI].start();
if (!forceOnPatches.found(patchI))
{
// Patch faces are +ve outwards,
// so the forces (exerted by fluid on solid)
// already have the correct sign
forceOnPatches.set
(
patchI,
(
rhoRef
* patchSf[patchI] * (patchPress[patchI] - pRef)
).ptr()
);
}
const vectorField& forceOnPatch = *forceOnPatches[patchI];
// Force from the patch-face into sum
forces[zonei] += forceOnPatch[patchFaceI];
// effective torque arm:
// - translated into the lumped-points coordinate system
// prior to determining the distance
const vector lever =
(
patchCf[patchI][patchFaceI] - centre_
- lumpedCentres[zonei]
);
// Moment from the patch-face into sum
moments[zonei] += lever ^ forceOnPatch[patchFaceI];
}
}
}
else
{
Info<<"no pressure" << endl;
}
Pstream::listCombineGather(forces, plusEqOp<vector>());
Pstream::listCombineScatter(forces);
Pstream::listCombineGather(moments, plusEqOp<vector>());
Pstream::listCombineScatter(moments);
return true;
}
// \cond fileScope
template<class T>
static inline T interpolatedValue
(
const Foam::UList<T> input,
const Foam::labelUList& indices,
const Foam::scalarField& weights
)
{
T result = weights[0] * input[indices[0]];
for (Foam::label i=1; i < indices.size(); ++i)
{
result += weights[i] * input[indices[i]];
}
return result;
}
// \endcond
Foam::tmp<Foam::pointField>
Foam::lumpedPointMovement::displacePoints
(
const pointField& points0,
const labelList& pointLabels
) const
{
return displacePoints(state(), points0, pointLabels);
}
Foam::tmp<Foam::pointField>
Foam::lumpedPointMovement::displacePoints
(
const lumpedPointState& state,
const pointField& points0,
const labelList& pointLabels
) const
{
// storage for the interpolation indices/weights
labelList indices;
scalarField weights;
const interpolationWeights& interp = interpolator();
// The mass centres
const pointField& lumpedCentres = state.points();
// local-to-global transformation tensor
const tensorField& localToGlobal = state.rotations();
// could also verify the sizes (state vs original)
tmp<pointField> tdisp(new pointField(pointLabels.size()));
pointField& disp = tdisp.ref();
forAll(pointLabels, ptI)
{
const point& p0 = points0[pointLabels[ptI]];
// Interpolation factor based on the axis component
// of the original point (location)
scalar where = p0 & axis_;
interp.valueWeights(where, indices, weights);
vector origin = interpolatedValue(lumpedCentres, indices, weights);
tensor rotTensor = interpolatedValue(localToGlobal, indices, weights);
if (indices.size() == 1)
{
// out-of-bounds, use corresponding end-point
where = locations_[indices[0]];
}
// local point:
// - translated into the lumped-points coordinate system
// - relative to the interpolation (where) location
const vector local = (p0 - (where * axis_)) - centre_;
// local-to-global rotation and translation
// - translate back out of the lumped-points coordinate system
//
// -> subtract p0 to return displacements
disp[ptI] = (rotTensor & local) + origin + centre_ - p0;
}
return tdisp;
}
void Foam::lumpedPointMovement::writeDict(Ostream& os) const
{
os.writeEntry("axis", axis_);
os.writeEntry("locations", locations_);
os.writeEntry("division", division_);
// os.writeEntry("interpolationScheme", interpolationScheme_);
}
bool Foam::lumpedPointMovement::readState()
{
lumpedPointState prev = state_;
bool status = state_.readData
(
inputFormat_,
coupler().resolveFile(inputName_)
);
state_.relax(relax_, prev);
return status;
}
bool Foam::lumpedPointMovement::writeData
(
const UList<vector>& forces
) const
{
if (!Pstream::master())
{
return false;
}
const fileName output(coupler().resolveFile(outputName_));
OFstream os(output); // ASCII
if (outputFormat_ == outputFormatType::PLAIN)
{
os <<"# output from OpenFOAM" << nl
<<"# N, points, forces" << nl
<< this->size() << nl;
const char* zeroVector = "0 0 0";
forAll(locations_, i)
{
const vector pos = locations_[i] * axis_;
os << pos.x() << ' '
<< pos.y() << ' '
<< pos.z() << ' ';
if (i < forces.size())
{
os << forces[i].x() << ' '
<< forces[i].y() << ' '
<< forces[i].z();
}
else
{
os << zeroVector;
}
os << nl;
}
}
else
{
// Make it easier for external programs to parse
// - exclude the usual OpenFOAM 'FoamFile' header
// - ensure lists have consistent format
os <<"// output from OpenFOAM" << nl << nl;
writeList(os, "points", (locations_*axis_)());
writeList(os, "forces", forces);
}
return true;
}
bool Foam::lumpedPointMovement::writeData
(
const UList<vector>& forces,
const UList<vector>& moments
) const
{
if (!Pstream::master())
{
return false;
}
const fileName output(coupler().resolveFile(outputName_));
OFstream os(output); // ASCII
if (outputFormat_ == outputFormatType::PLAIN)
{
os <<"# output from OpenFOAM" << nl
<<"# N, points, forces, moments" << nl
<< this->size() << nl;
const char* zeroVector = "0 0 0";
forAll(locations_, i)
{
const vector pos = locations_[i] * axis_;
os << pos.x() << ' '
<< pos.y() << ' '
<< pos.z() << ' ';
if (i < forces.size())
{
os << forces[i].x() << ' '
<< forces[i].y() << ' '
<< forces[i].z() << ' ';
}
else
{
os << zeroVector << ' ';
}
if (i < moments.size())
{
os << moments[i].x() << ' '
<< moments[i].y() << ' '
<< moments[i].z();
}
else
{
os << zeroVector;
}
os << nl;
}
}
else
{
// Make it easier for external programs to parse
// - exclude the usual OpenFOAM 'FoamFile' header
// - ensure lists have consistent format
os <<"// output from OpenFOAM" << nl << nl;
writeList(os, "points", (locations_*axis_)());
writeList(os, "forces", forces);
writeList(os, "moments", moments);
}
return true;
}
const Foam::interpolationWeights&
Foam::lumpedPointMovement::interpolator() const
{
if (!interpolatorPtr_.valid())
{
interpolatorPtr_ = interpolationWeights::New
(
interpolationScheme_,
locations_
);
}
return interpolatorPtr_();
}
// ************************************************************************* //

View File

@ -0,0 +1,399 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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::lumpedPointMovement
Description
The movement 'driver' that describes initial point locations, the
segmentation for pressure integration, the current state of the
points/rotations, and forwarding to the communication
coordinator. The 'lumpedPointIOMovement' class is simply a
registered version of the same.
SourceFiles
lumpedPointMovement.C
\*---------------------------------------------------------------------------*/
#ifndef lumpedPointMovement_H
#define lumpedPointMovement_H
#include "dictionary.H"
#include "scalarList.H"
#include "scalarField.H"
#include "pointField.H"
#include "vectorField.H"
#include "tensorField.H"
#include "vector.H"
#include "interpolationWeights.H"
#include "IOobject.H"
#include "tmp.H"
#include "faceZoneMeshFwd.H"
#include "externalCoupler.H"
#include "lumpedPointState.H"
#include "boundBox.H"
#include "Enum.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// Forward declarations
class polyMesh;
/*---------------------------------------------------------------------------*\
Class lumpedPointMovement Declaration
\*---------------------------------------------------------------------------*/
class lumpedPointMovement
{
public:
//- Output format types
enum class outputFormatType
{
PLAIN,
DICTIONARY
};
// Static data
//- Names for the output format types
static const Enum<outputFormatType> formatNames;
private:
// Private data
//- Reference axis for the locations
vector axis_;
//- The locations of the lumped points within the reference axis
// The interpolator needs scalarField, not scalarList.
scalarField locations_;
//- The division when calculating pressure forces (0-1)
scalar division_;
//- The relaxation factor when moving points (default: 1)
scalar relax_;
//- The interpolation type (linear|spline)
word interpolationScheme_;
//- Optional owner information
label ownerId_;
//- The bounding box (if set)
boundBox boundBox_;
//- The offset of patch points to compared to the locations
point centre_;
//- Calculate centre based on the bounding box
bool autoCentre_;
//- Dictionary of controls for force calculation
dictionary forcesDict_;
//- Communication control
externalCoupler coupler_;
//- File io
word inputName_;
word outputName_;
lumpedPointState::inputFormatType inputFormat_;
outputFormatType outputFormat_;
// Demand-driven private data
//- The initial state (positions, rotations)
lumpedPointState state0_;
//- The current state (positions, rotations)
// Eg, as response from external application
lumpedPointState state_;
//- Threshold locations for pressure forces
mutable scalarField* thresholdPtr_;
//- User-specified interpolator
mutable autoPtr<interpolationWeights> interpolatorPtr_;
//- Pressure zones (only used from the master patch)
mutable List<labelList> faceZones_;
// Private Member Functions
//- Calculate threshold locations
void calcThresholds() const;
//- Classify the position to be located in one of the threshold zones
label threshold(scalar pos) const;
//- Disallow default bitwise copy constructor
lumpedPointMovement(const lumpedPointMovement&) = delete;
//- Disallow default bitwise assignment
void operator=(const lumpedPointMovement&) = delete;
public:
// Static data members
//- The standard dictionary name
static const word dictionaryName;
//- Helper: return IOobject with optionally absolute path provided
static IOobject selectIO(const IOobject&, const fileName&);
// Constructors
//- Construct null
lumpedPointMovement();
//- Construct from dictionary, optionally with some owner information
lumpedPointMovement(const dictionary& dict, label ownerId=-1);
//- Destructor
virtual ~lumpedPointMovement();
// Member Functions
//- Update settings from dictionary
void readDict(const dictionary& dict);
//- If no number of lumped points (locations) were specified
inline bool empty() const;
//- The number of lumped points (number of locations)
inline label size() const;
//- The normalized reference axis
inline const vector& axis() const;
//- Read access to the locations
inline const scalarField& locations() const;
//- The division (0-1) when calculating pressure forces
inline scalar division() const;
//- An owner Id, if needed for bookkeeping purposes
inline label ownerId() const;
//- Change the owner id, if needed for bookkeeping purposes
inline void ownerId(label id);
//- Threshold locations for pressure forces
inline const scalarField& thresholds() const;
//- Classify the position to be located in one of the threshold zones
inline label threshold(const point& position) const;
//- Communication control
inline const externalCoupler& coupler() const;
//- Communication control
inline externalCoupler& coupler();
//- The initial state (positions/rotations)
inline const lumpedPointState& state0() const;
//- The current state (positions/rotations)
inline const lumpedPointState& state() const;
//- The relaxation factor when changing states
inline scalar relax() const;
//- The relaxation factor when changing states
inline scalar& relax();
//- The input (state) file name
inline const word& inputName() const;
//- The output (forces) file name
inline const word& outputName() const;
//- The input (state) file format
inline lumpedPointState::inputFormatType inputFormat() const;
//- The output (forces) file format
inline lumpedPointMovement::outputFormatType outputFormat() const;
//- Define the bounding-box required to enclose the specified patches.
// Calculates the centre as required.
//
// \param mesh The volume mesh reference
// \param patchIds The patch ids to be included in the bounding box.
// \param points0 The initial mesh points, prior to movement
void setBoundBox
(
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
);
//- Define the pressure-zones mapping for faces in the specified
// patches.
// The face centres are compared to the threshold positions,
// which are determined by locations along the defined axis.
//
// \param mesh The volume mesh reference
// \param patchIds The patch ids to be included in the mapping
// \param points0 The initial mesh points, prior to movement
void setMapping
(
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
);
//- True if the pressure-zones mapping has already been performed
inline bool hasMapping() const;
//- Return the pressure-zones mapping with the associated
// patch face ids.
inline const List<labelList>& zones() const;
//- The forces and moments acting on each pressure-zone.
// The zones must be previously defined via setMapping.
bool forcesAndMoments
(
const polyMesh& pmesh,
List<vector>& forces,
List<vector>& moments
) const;
//- Displace points according to the current state
tmp<pointField> displacePoints
(
const pointField& points0,
const labelList& pointLabels
) const;
//- Displace points according to specified state
tmp<pointField> displacePoints
(
const lumpedPointState& state,
const pointField& points0,
const labelList& pointLabels
) const;
//- Interpolation weights
const interpolationWeights& interpolator() const;
//- Write axis, locations, division as a dictionary
void writeDict(Ostream& os) const;
//- Write points, forces
bool writeData
(
const UList<vector>& forces
) const;
//- Write points, forces, moments
bool writeData
(
const UList<vector>& forces,
const UList<vector>& moments
) const;
//- Read state from file, applying relaxation as requested
bool readState();
//- Write state as VTK PolyData format.
void writeStateVTP(const fileName& file) const;
//- Write forces on points as VTK PolyData format.
void writeForcesAndMomentsVTP
(
const fileName& file,
const UList<vector>& forces,
const UList<vector>& moments
) const;
//- Write pressure-zones geometry, write as VTK PolyData format.
void writeZonesVTP
(
const fileName& file,
const polyMesh& mesh,
const pointField& points0
) const;
//- Write displaced geometry according to the current state,
// write as VTK PolyData format.
void writeVTP
(
const fileName& file,
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
) const;
//- Write displaced geometry according to the specified state,
// write as VTK PolyData format.
void writeVTP
(
const fileName& file,
const lumpedPointState& state,
const polyMesh& mesh,
const labelUList& patchLst,
const pointField& points0
) const;
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "lumpedPointMovementI.H"
#endif
// ************************************************************************* //

View File

@ -0,0 +1,162 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 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/>.
\*---------------------------------------------------------------------------*/
inline bool Foam::lumpedPointMovement::empty() const
{
return locations_.empty();
}
inline Foam::label Foam::lumpedPointMovement::size() const
{
return locations_.size();
}
inline const Foam::vector& Foam::lumpedPointMovement::axis() const
{
return axis_;
}
inline const Foam::scalarField& Foam::lumpedPointMovement::locations() const
{
return locations_;
}
inline Foam::scalar Foam::lumpedPointMovement::division() const
{
return division_;
}
inline Foam::label Foam::lumpedPointMovement::ownerId() const
{
return ownerId_;
}
inline void Foam::lumpedPointMovement::ownerId(label id)
{
ownerId_ = id;
}
inline const Foam::scalarField& Foam::lumpedPointMovement::thresholds() const
{
if (!thresholdPtr_)
{
calcThresholds();
}
return *thresholdPtr_;
}
inline Foam::label
Foam::lumpedPointMovement::threshold(const point& position) const
{
return threshold(position & axis_);
}
inline const Foam::externalCoupler& Foam::lumpedPointMovement::coupler() const
{
return coupler_;
}
inline Foam::externalCoupler& Foam::lumpedPointMovement::coupler()
{
return coupler_;
}
//- The initial state (positions/rotations)
inline const Foam::lumpedPointState& Foam::lumpedPointMovement::state0() const
{
return state0_;
}
inline const Foam::lumpedPointState& Foam::lumpedPointMovement::state() const
{
return state_;
}
inline Foam::scalar Foam::lumpedPointMovement::relax() const
{
return relax_;
}
inline Foam::scalar& Foam::lumpedPointMovement::relax()
{
return relax_;
}
inline const Foam::word& Foam::lumpedPointMovement::inputName() const
{
return inputName_;
}
inline const Foam::word& Foam::lumpedPointMovement::outputName() const
{
return outputName_;
}
inline Foam::lumpedPointState::inputFormatType
Foam::lumpedPointMovement::inputFormat() const
{
return inputFormat_;
}
Foam::lumpedPointMovement::outputFormatType
Foam::lumpedPointMovement::outputFormat() const
{
return outputFormat_;
}
inline bool Foam::lumpedPointMovement::hasMapping() const
{
return !faceZones_.empty();
}
inline const Foam::List<Foam::labelList>&
Foam::lumpedPointMovement::zones() const
{
return faceZones_;
}
// ************************************************************************* //

View File

@ -0,0 +1,447 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointMovement.H"
#include "polyMesh.H"
#include "OFstream.H"
#include "uindirectPrimitivePatch.H"
#include "foamVtkOutput.H"
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void Foam::lumpedPointMovement::writeStateVTP(const fileName& file) const
{
state().writeVTP(file, axis());
}
void Foam::lumpedPointMovement::writeForcesAndMomentsVTP
(
const fileName& file,
const UList<vector>& forces,
const UList<vector>& moments
) const
{
OFstream fos(file);
std::ostream& os = fos.stdStream();
autoPtr<vtk::formatter> format = vtk::newFormatter
(
os,
vtk::formatType::INLINE_ASCII
);
format().xmlHeader()
.beginVTKFile(vtk::fileTag::POLY_DATA, "0.1");
//
// The 'spine' of lumped mass points
//
const label nPoints = state().points().size();
{
format()
.openTag(vtk::fileTag::PIECE)
.xmlAttr(vtk::fileAttr::NUMBER_OF_POINTS, nPoints)
.xmlAttr(vtk::fileAttr::NUMBER_OF_VERTS, nPoints)
.closeTag();
// 'points'
{
const uint64_t payLoad = (nPoints*3* sizeof(float));
format()
.tag(vtk::fileTag::POINTS)
.openDataArray<float, 3>(vtk::dataArrayAttr::POINTS)
.closeTag();
format().writeSize(payLoad);
vtk::writeList(format(), state().points());
format().flush();
format()
.endDataArray()
.endTag(vtk::fileTag::POINTS);
}
// <Verts>
format().tag(vtk::fileTag::VERTS);
//
// 'connectivity'
//
{
const uint64_t payLoad = (nPoints*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::CONNECTIVITY)
.closeTag();
format().writeSize(payLoad);
for (label i=0; i<nPoints; ++i)
{
format().write(i);
}
format().flush();
format().endDataArray();
}
//
// 'offsets' (connectivity offsets)
// = linear mapping onto points (with 1 offset)
//
{
const uint64_t payLoad = (nPoints*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::OFFSETS)
.closeTag();
format().writeSize(payLoad);
for (label i=0; i<nPoints; ++i)
{
format().write(i+1);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::VERTS);
// </Verts>
}
format().tag(vtk::fileTag::POINT_DATA);
// forces
if (forces.size() == nPoints)
{
const uint64_t payLoad = (nPoints * 3 * sizeof(float));
format().openDataArray<float, 3>("forces")
.closeTag();
format().writeSize(payLoad);
vtk::writeList(format(), forces);
format().flush();
format().endDataArray();
}
// moments
if (moments.size() == nPoints)
{
const uint64_t payLoad = (nPoints * 3 * sizeof(float));
format().openDataArray<float, 3>("moments")
.closeTag();
format().writeSize(payLoad);
vtk::writeList(format(), moments);
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::POINT_DATA);
format().endTag(vtk::fileTag::PIECE);
format().endTag(vtk::fileTag::POLY_DATA);
format().endVTKFile();
}
void Foam::lumpedPointMovement::writeZonesVTP
(
const fileName& file,
const polyMesh& mesh,
const pointField& points0
) const
{
OFstream fos(file);
std::ostream& os = fos.stdStream();
autoPtr<vtk::formatter> format = vtk::newFormatter
(
os,
vtk::formatType::INLINE_ASCII
);
format().xmlHeader()
.beginVTKFile(vtk::fileTag::POLY_DATA, "0.1");
forAll(faceZones_, zoneI)
{
uindirectPrimitivePatch pp
(
UIndirectList<face>(mesh.faces(), faceZones_[zoneI]),
points0
);
format()
.openTag(vtk::fileTag::PIECE)
.xmlAttr(vtk::fileAttr::NUMBER_OF_POINTS, pp.nPoints())
.xmlAttr(vtk::fileAttr::NUMBER_OF_POLYS, pp.size())
.closeTag();
// 'points'
{
const uint64_t payLoad = (pp.nPoints()*3* sizeof(float));
format()
.tag(vtk::fileTag::POINTS)
.openDataArray<float, 3>(vtk::dataArrayAttr::POINTS)
.closeTag();
format().writeSize(payLoad);
vtk::writeList(format(), pp.localPoints());
format().flush();
format()
.endDataArray()
.endTag(vtk::fileTag::POINTS);
}
// <Polys>
format().tag(vtk::fileTag::POLYS);
//
// 'connectivity'
//
{
uint64_t payLoad = 0;
forAll(pp, facei)
{
const face& f = pp[facei];
payLoad += f.size();
}
format().openDataArray<label>(vtk::dataArrayAttr::CONNECTIVITY)
.closeTag();
format().writeSize(payLoad * sizeof(label));
for (const face& f : pp.localFaces())
{
vtk::writeList(format(), f);
}
format().flush();
format().endDataArray();
}
//
// 'offsets' (connectivity offsets)
//
{
const uint64_t payLoad = (pp.size() * sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::OFFSETS)
.closeTag();
format().writeSize(payLoad);
label off = 0;
forAll(pp, facei)
{
const face& f = pp[facei];
off += f.size();
format().write(off);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::POLYS);
format().tag(vtk::fileTag::CELL_DATA);
// zone Id
{
const uint64_t payLoad = (pp.size() * sizeof(label));
format().openDataArray<label>("zoneId")
.closeTag();
format().writeSize(payLoad);
forAll(pp, facei)
{
format().write(zoneI);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::CELL_DATA);
format().endTag(vtk::fileTag::PIECE);
}
format().endTag(vtk::fileTag::POLY_DATA);
format().endVTKFile();
}
void Foam::lumpedPointMovement::writeVTP
(
const fileName& file,
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
) const
{
writeVTP(file, state(), mesh, patchIds, points0);
}
void Foam::lumpedPointMovement::writeVTP
(
const fileName& file,
const lumpedPointState& state,
const polyMesh& mesh,
const labelUList& patchIds,
const pointField& points0
) const
{
const polyBoundaryMesh& boundaryMesh = mesh.boundaryMesh();
OFstream fos(file);
std::ostream& os = fos.stdStream();
autoPtr<vtk::formatter> format = vtk::newFormatter
(
os,
vtk::formatType::INLINE_ASCII
);
format().xmlHeader()
.beginVTKFile(vtk::fileTag::POLY_DATA, "0.1");
for (const label patchId : patchIds)
{
const polyPatch& pp = boundaryMesh[patchId];
format()
.openTag(vtk::fileTag::PIECE)
.xmlAttr(vtk::fileAttr::NUMBER_OF_POINTS, pp.nPoints())
.xmlAttr(vtk::fileAttr::NUMBER_OF_POLYS, pp.size())
.closeTag();
// 'points'
{
const uint64_t payLoad = (pp.nPoints()*3* sizeof(float));
format()
.tag(vtk::fileTag::POINTS)
.openDataArray<float, 3>(vtk::dataArrayAttr::POINTS)
.closeTag();
// Could be more efficient, but not often needed
tmp<pointField> tpts = displacePoints
(
state,
points0,
pp.meshPoints()
) + pointField(points0, pp.meshPoints());
const pointField& pts = tpts();
format().writeSize(payLoad);
vtk::writeList(format(), pts);
format().flush();
format()
.endDataArray()
.endTag(vtk::fileTag::POINTS);
}
// <Polys>
format().tag(vtk::fileTag::POLYS);
//
// 'connectivity'
//
{
uint64_t payLoad = 0;
for (const face& f : pp)
{
payLoad += f.size();
}
format().openDataArray<label>(vtk::dataArrayAttr::CONNECTIVITY)
.closeTag();
format().writeSize(payLoad * sizeof(label));
for (const face& f : pp.localFaces())
{
vtk::writeList(format(), f);
}
format().flush();
format().endDataArray();
}
//
// 'offsets' (connectivity offsets)
//
{
const uint64_t payLoad = (pp.size() * sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::OFFSETS)
.closeTag();
format().writeSize(payLoad);
label off = 0;
forAll(pp, facei)
{
const face& f = pp[facei];
off += f.size();
format().write(off);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::POLYS);
format().endTag(vtk::fileTag::PIECE);
}
format().endTag(vtk::fileTag::POLY_DATA);
format().endVTKFile();
}
// ************************************************************************* //

View File

@ -0,0 +1,362 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointState.H"
#include "demandDrivenData.H"
#include "EulerCoordinateRotation.H"
#include "unitConversion.H"
#include "ISstream.H"
#include "IFstream.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
const Foam::Enum<Foam::lumpedPointState::inputFormatType>
Foam::lumpedPointState::formatNames
{
{ inputFormatType::PLAIN, "plain" },
{ inputFormatType::DICTIONARY, "dictionary" }
};
// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
//! \cond fileScope
static Foam::string getLineNoComment(Foam::ISstream& is)
{
Foam::string line;
do
{
is.getLine(line);
}
while ((line.empty() || line[0] == '#') && is.good());
return line;
}
//! \endcond
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::lumpedPointState::calcRotations() const
{
rotationPtr_ = new tensorField(angles_.size());
forAll(angles_, itemi)
{
rotationPtr_->operator[](itemi) = EulerCoordinateRotation
(
angles_[itemi],
degrees_ // true=degrees, false=radians
).R();
}
}
void Foam::lumpedPointState::readDict(const dictionary& dict)
{
dict.lookup("points") >> points_;
dict.lookup("angles") >> angles_;
degrees_ = dict.lookupOrDefault("degrees", false);
deleteDemandDrivenData(rotationPtr_);
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::lumpedPointState::lumpedPointState()
:
points_(0),
angles_(0),
degrees_(false),
rotationPtr_(nullptr)
{}
Foam::lumpedPointState::lumpedPointState(const lumpedPointState& rhs)
:
points_(rhs.points_),
angles_(rhs.angles_),
degrees_(rhs.degrees_),
rotationPtr_(nullptr)
{}
Foam::lumpedPointState::lumpedPointState
(
const pointField& pts
)
:
points_(pts),
angles_(points_.size(), Zero),
degrees_(false),
rotationPtr_(nullptr)
{}
Foam::lumpedPointState::lumpedPointState
(
tmp<pointField>& pts
)
:
points_(pts),
angles_(points_.size(), Zero),
degrees_(false),
rotationPtr_(nullptr)
{}
Foam::lumpedPointState::lumpedPointState
(
const dictionary& dict
)
:
points_(0),
angles_(0),
degrees_(false),
rotationPtr_(nullptr)
{
readDict(dict);
}
// * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * //
Foam::lumpedPointState::~lumpedPointState()
{
deleteDemandDrivenData(rotationPtr_);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::lumpedPointState::operator=(const lumpedPointState& rhs)
{
points_ = rhs.points_;
angles_ = rhs.angles_;
degrees_ = rhs.degrees_;
deleteDemandDrivenData(rotationPtr_);
}
void Foam::lumpedPointState::relax
(
const scalar alpha,
const lumpedPointState& prev
)
{
points_ = prev.points_ + alpha*(points_ - prev.points_);
scalar convert = 1.0;
if (degrees_ != prev.degrees_)
{
if (prev.degrees_)
{
// Was degrees, now radians
convert = degToRad(1);
}
else
{
// Was radians, now degrees
convert = radToDeg(1);
}
}
angles_ = convert*prev.angles_ + alpha*(angles_ - convert*prev.angles_);
deleteDemandDrivenData(rotationPtr_);
}
bool Foam::lumpedPointState::readPlain(Istream& is)
{
// Assume generic input stream so we can do line-based input
ISstream& iss = dynamic_cast<ISstream&>(is);
string line = getLineNoComment(iss);
label count;
{
IStringStream isstr(line);
isstr >> count;
}
points_.setSize(count);
angles_.setSize(count);
count = 0;
forAll(points_, i)
{
iss.getLine(line);
IStringStream isstr(line);
isstr
>> points_[count].x() >> points_[count].y() >> points_[count].z()
>> angles_[count].x() >> angles_[count].y() >> angles_[count].z();
++count;
}
points_.setSize(count);
angles_.setSize(count);
degrees_ = false;
deleteDemandDrivenData(rotationPtr_);
return count;
}
bool Foam::lumpedPointState::readData(Istream& is)
{
dictionary dict(is);
readDict(dict);
return points_.size();
}
bool Foam::lumpedPointState::writeData(Ostream& os) const
{
writeDict(os);
return true;
}
void Foam::lumpedPointState::writeDict(Ostream& os) const
{
os.writeEntry("points", points_);
os.writeEntry("angles", angles_);
if (degrees_)
{
os.writeKeyword("degrees") << "true;" << nl;
}
}
void Foam::lumpedPointState::writePlain(Ostream& os) const
{
os <<"# input for OpenFOAM\n"
<<"# N, points, angles\n"
<< points_.size() << "\n";
forAll(points_, i)
{
const vector& pt = points_[i];
os << pt.x() << ' '
<< pt.y() << ' '
<< pt.z() << ' ';
if (i < angles_.size())
{
os << angles_[i].x() << ' '
<< angles_[i].y() << ' '
<< angles_[i].z() << '\n';
}
else
{
os << "0 0 0\n";
}
}
}
bool Foam::lumpedPointState::readData
(
const inputFormatType& fmt,
const fileName& file
)
{
bool ok = false;
if (Pstream::master())
{
IFstream is(file);
if (fmt == inputFormatType::PLAIN)
{
ok = this->readPlain(is);
}
else
{
ok = this->readData(is);
}
}
if (Pstream::parRun())
{
// Scatter master data using communication scheme
const List<Pstream::commsStruct>& comms =
(
(Pstream::nProcs() < Pstream::nProcsSimpleSum)
? Pstream::linearCommunication()
: Pstream::treeCommunication()
);
// Get my communication order
const Pstream::commsStruct& myComm = comms[Pstream::myProcNo()];
// Receive from up
if (myComm.above() != -1)
{
IPstream fromAbove
(
UPstream::commsTypes::scheduled,
myComm.above(),
0,
Pstream::msgType(),
Pstream::worldComm
);
fromAbove >> points_ >> angles_ >> degrees_;
}
// Send to downstairs neighbours
forAllReverse(myComm.below(), belowI)
{
OPstream toBelow
(
UPstream::commsTypes::scheduled,
myComm.below()[belowI],
0,
Pstream::msgType(),
Pstream::worldComm
);
toBelow << points_ << angles_ << degrees_;
}
deleteDemandDrivenData(rotationPtr_);
// MPI barrier
Pstream::scatter(ok);
}
return ok;
}
// ************************************************************************* //

View File

@ -0,0 +1,194 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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::lumpedPointState
Description
The 'state' of lumped points corresponds to positions and rotations.
Encapsulates the response from the external application.
Entry place for applying relaxation and sub-stepping etc.
SourceFiles
lumpedPointState.C
\*---------------------------------------------------------------------------*/
#ifndef lumpedPointState_H
#define lumpedPointState_H
#include "dictionary.H"
#include "scalarList.H"
#include "pointField.H"
#include "scalarField.H"
#include "vectorField.H"
#include "tensorField.H"
#include "Enum.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class Istream;
class Ostream;
/*---------------------------------------------------------------------------*\
Class lumpedPointState Declaration
\*---------------------------------------------------------------------------*/
//- Bundling of position/rotation
class lumpedPointState
{
public:
//- Input format types
enum class inputFormatType
{
PLAIN,
DICTIONARY
};
// Static data
//- Names for the input format types
static const Enum<inputFormatType> formatNames;
private:
// Private data
//- Positions of lumped points
pointField points_;
//- Orientation of lumped points (as Euler angles)
vectorField angles_;
//- Euler angles in degrees instead radians
bool degrees_;
//- Tensor rotation of lumped points
mutable tensorField* rotationPtr_;
// Private Member Functions
void calcRotations() const;
void readDict(const dictionary& dict);
public:
// Constructors
//- Construct null
lumpedPointState();
//- Copy constructor
lumpedPointState(const lumpedPointState& rhs);
//- Construct from points with zero-rotation
lumpedPointState(const pointField& pts);
//- Construct from points with zero-rotation
lumpedPointState(tmp<pointField>& pts);
//- Construct from dictionary
lumpedPointState(const dictionary& dict);
//- Destructor
virtual ~lumpedPointState();
// Member Functions
//- Has positions and consistent number of rotations?
inline bool valid() const;
//- If no points were specified
inline bool empty() const;
//- The number of points
inline label size() const;
//- The points corresponding to mass centres
inline const pointField& points() const;
//- The orientation of the points (mass centres)
inline const vectorField& angles() const;
//- The local-to-global transformation for each point
inline const tensorField& rotations() const;
//- Relax the state
// alpha = 1 : no relaxation
// alpha < 1 : relaxation
// alpha = 0 : do nothing
void relax(const scalar alpha, const lumpedPointState& prev);
//- Read input as dictionary content
bool readData(Istream& is);
//- Output as dictionary content
bool writeData(Ostream& os) const;
//- Output as dictionary content
void writeDict(Ostream& os) const;
//- Read input as plain content
bool readPlain(Istream& is);
//- Output as plain content
void writePlain(Ostream& os) const;
//- Read input from file (master process only) using specified format
bool readData(const inputFormatType& fmt, const fileName& file);
//- Output as VTK file for debugging/visualization
// The points are joined as lines, the rotation is visualized
// by planes, write as VTK PolyData format.
void writeVTP(const fileName& outputFile, const vector& axis) const;
// Member Operators
//- Assignment operator
void operator=(const lumpedPointState& rhs);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "lumpedPointStateI.H"
#endif
// ************************************************************************* //

View File

@ -0,0 +1,67 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2017 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/>.
\*---------------------------------------------------------------------------*/
inline bool Foam::lumpedPointState::valid() const
{
return points_.size() && points_.size() == angles_.size();
}
inline bool Foam::lumpedPointState::empty() const
{
return points_.empty();
}
inline Foam::label Foam::lumpedPointState::size() const
{
return points_.size();
}
inline const Foam::pointField& Foam::lumpedPointState::points() const
{
return points_;
}
inline const Foam::vectorField& Foam::lumpedPointState::angles() const
{
return angles_;
}
inline const Foam::tensorField& Foam::lumpedPointState::rotations() const
{
if (!rotationPtr_)
{
calcRotations();
}
return *rotationPtr_;
}
// ************************************************************************* //

View File

@ -0,0 +1,326 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointState.H"
#include "OFstream.H"
#include "axesRotation.H"
#include "foamVtkOutput.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
// file-local
const static Foam::FixedList<Foam::point, 4> standardCorners
{
{-0.1, -0.1, 0},
{+0.1, -0.1, 0},
{+0.1, +0.1, 0},
{-0.1, +0.1, 0}
};
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::lumpedPointState::writeVTP
(
const fileName& outputFile,
const vector& axis
) const
{
// local-to-global transformation tensors
const tensorField& localToGlobal = rotations();
OFstream fos(outputFile);
std::ostream& os = fos.stdStream();
autoPtr<vtk::formatter> format = vtk::newFormatter
(
os,
vtk::formatType::INLINE_ASCII
);
format().xmlHeader()
.beginVTKFile(vtk::fileTag::POLY_DATA, "0.1");
//
// The 'spine' of lumped mass points
//
{
format()
.openTag(vtk::fileTag::PIECE)
.xmlAttr(vtk::fileAttr::NUMBER_OF_POINTS, points_.size())
.xmlAttr(vtk::fileAttr::NUMBER_OF_VERTS, points_.size())
.xmlAttr(vtk::fileAttr::NUMBER_OF_LINES, 1)
.closeTag();
// 'points'
{
const uint64_t payLoad = (points_.size()*3* sizeof(float));
format()
.tag(vtk::fileTag::POINTS)
.openDataArray<float, 3>(vtk::dataArrayAttr::POINTS)
.closeTag();
format().writeSize(payLoad);
vtk::writeList(format(), points_);
format().flush();
format()
.endDataArray()
.endTag(vtk::fileTag::POINTS);
}
// <Verts>
format().tag(vtk::fileTag::VERTS);
//
// 'connectivity'
//
{
const uint64_t payLoad = (points_.size()*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::CONNECTIVITY)
.closeTag();
format().writeSize(payLoad);
forAll(points_, i)
{
format().write(i);
}
format().flush();
format().endDataArray();
}
//
// 'offsets' (connectivity offsets)
// = linear mapping onto points (with 1 offset)
//
{
const uint64_t payLoad = (points_.size()*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::OFFSETS)
.closeTag();
format().writeSize(payLoad);
forAll(points_, i)
{
format().write(i+1);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::VERTS);
// </Verts>
// <Lines>
format().tag(vtk::fileTag::LINES);
//
// 'connectivity'
//
{
const uint64_t payLoad = (points_.size()*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::CONNECTIVITY)
.closeTag();
format().writeSize(payLoad);
forAll(points_, i)
{
format().write(i);
}
format().flush();
format().endDataArray();
}
//
// 'offsets' (connectivity offsets)
// = single line
//
{
const uint64_t payLoad = (1*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::OFFSETS)
.closeTag();
format().writeSize(payLoad);
format().write(points_.size());
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::LINES);
format().endTag(vtk::fileTag::PIECE);
}
// Standard corners in local axis
const axesRotation cornerTransform(axis);
FixedList<point, 4> corners;
forAll(standardCorners, corni)
{
corners[corni] = cornerTransform.transform(standardCorners[corni]);
}
//
// Planes to visualize location/rotation
//
{
const label nPoints = 4*points_.size(); // 4 points per quad
const label nPolys = points_.size();
format()
.openTag(vtk::fileTag::PIECE)
.xmlAttr(vtk::fileAttr::NUMBER_OF_POINTS, nPoints)
.xmlAttr(vtk::fileAttr::NUMBER_OF_POLYS, nPolys)
.closeTag();
// 'points'
{
const uint64_t payLoad = (nPoints*3*sizeof(float));
format()
.tag(vtk::fileTag::POINTS)
.openDataArray<float, 3>(vtk::dataArrayAttr::POINTS)
.closeTag();
format().writeSize(payLoad);
forAll(points_, posI)
{
const point& origin = points_[posI];
const tensor& rotTensor =
(
posI < localToGlobal.size()
? localToGlobal[posI]
: pTraits<tensor>::I
);
for (const point& cornerPt : corners)
{
// Local-to-global rotation and translation
const point pt = (rotTensor & cornerPt) + origin;
vtk::write(format(), pt);
}
}
format().flush();
format()
.endDataArray()
.endTag(vtk::fileTag::POINTS);
}
// <Polys>
format().tag(vtk::fileTag::POLYS);
//
// 'connectivity' - 4 points (ie, quad)
//
{
const uint64_t payLoad = (4*nPolys*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::CONNECTIVITY)
.closeTag();
format().writeSize(payLoad);
for (label i=0; i < 4*nPolys; ++i)
{
format().write(i);
}
format().flush();
format().endDataArray();
}
//
// 'offsets' (connectivity offsets)
// = single quad
//
{
const uint64_t payLoad = (nPolys*sizeof(label));
format().openDataArray<label>(vtk::dataArrayAttr::OFFSETS)
.closeTag();
format().writeSize(payLoad);
for (label i=0; i < nPolys; ++i)
{
const label off = 4 * (i+1);
format().write(off);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::POLYS);
#if 0
format().tag(vtk::fileTag::CELL_DATA);
// zone Id
{
const uint64_t payLoad = (points_.size()*sizeof(label));
format().openDataArray<label>("zoneId")
.closeTag();
format().writeSize(payLoad);
for (label i=0; i < nPolys; ++i)
{
format().write(i);
}
format().flush();
format().endDataArray();
}
format().endTag(vtk::fileTag::CELL_DATA);
#endif
format().endTag(vtk::fileTag::PIECE);
}
// Finally
// could add a 'ghost' level above to visualize extrapolated values
// draw as two triangles to distingush from real levels ...
format().endTag(vtk::fileTag::POLY_DATA);
format().endVTKFile();
}
// ************************************************************************* //

View File

@ -0,0 +1,161 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016-2017 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 "lumpedPointTools.H"
#include "IFstream.H"
#include "IOobjectList.H"
#include "volFields.H"
#include "lumpedPointDisplacementPointPatchVectorField.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
// file-scope
template<class GeoFieldType>
static autoPtr<GeoFieldType> loadPointField
(
const pointMesh::Mesh& mesh,
const IOobject* io
)
{
if (io && io->headerClassName() == GeoFieldType::typeName)
{
Info<< "Reading " << GeoFieldType::typeName
<< ' ' << io->name() << endl;
return autoPtr<GeoFieldType>
(
new GeoFieldType
(
IOobject
(
io->name(),
io->instance(),
io->local(),
io->db(),
IOobject::MUST_READ,
IOobject::AUTO_WRITE,
io->registerObject()
),
mesh
)
);
}
return autoPtr<GeoFieldType>();
}
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Foam::List<Foam::lumpedPointStateTuple>
Foam::lumpedPointTools::lumpedPointStates(Istream& is)
{
dictionary contents(is);
List<dictionary> entries(contents.lookup("response"));
DynamicList<Tuple2<scalar, lumpedPointState>> states(entries.size());
for (const dictionary& dict : entries)
{
states.append
(
lumpedPointStateTuple
(
readScalar(dict.lookup("time")),
lumpedPointState(dict)
)
);
}
return states.shrink();
}
Foam::List<Foam::lumpedPointStateTuple>
Foam::lumpedPointTools::lumpedPointStates(const fileName& file)
{
IFstream is(file);
return lumpedPointStates(is);
}
Foam::pointIOField
Foam::lumpedPointTools::points0Field(const polyMesh& mesh)
{
pointIOField pts
(
IOobject
(
"points",
mesh.time().constant(),
polyMesh::meshSubDir,
mesh,
IOobject::MUST_READ,
IOobject::NO_WRITE,
false // Do not re-register
)
);
return pts;
}
Foam::labelList
Foam::lumpedPointTools::lumpedPointPatchList(const pointVectorField& pvf)
{
return lumpedPointDisplacementPointPatchVectorField::patchIds(pvf);
}
Foam::labelList Foam::lumpedPointTools::lumpedPointPatchList
(
const polyMesh& mesh
)
{
IOobjectList objects0(mesh, "0");
pointMesh pMesh(mesh);
autoPtr<pointVectorField> displacePtr = loadPointField<pointVectorField>
(
pMesh,
objects0.lookup("pointDisplacement")
);
if (!displacePtr.valid())
{
Info<< "no valid pointDisplacement" << endl;
return labelList();
}
return lumpedPointPatchList(displacePtr());
}
// ************************************************************************* //

View File

@ -0,0 +1,86 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2016 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/>.
Namespace
Foam::lumpedPointTools
Description
A collection of utility functions for handling IO related to the
lumped-mass movement.
SourceFiles
lumpedPointTools.C
\*---------------------------------------------------------------------------*/
#ifndef lumpedPointTools_H
#define lumpedPointTools_H
#include "labelList.H"
#include "polyMesh.H"
#include "pointMesh.H"
#include "pointFields.H"
#include "Tuple2.H"
#include "lumpedPointState.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
typedef Tuple2<scalar, lumpedPointState> lumpedPointStateTuple;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace lumpedPointTools
{
//- Load a list of states from an Istream
List<lumpedPointStateTuple> lumpedPointStates(Istream& is);
//- Load a list of states from a file
List<lumpedPointStateTuple> lumpedPointStates(const fileName& file);
//- Return the 0 or constant points field
pointIOField points0Field(const polyMesh& mesh);
//- Return the patch-ids associated a "lumpedPointDisplacement" type
labelList lumpedPointPatchList(const pointVectorField& pvf);
//- Get the "pointDisplacement" at time 0 and use that to determine which
// patches have a "lumpedPointDisplacement" type
labelList lumpedPointPatchList(const polyMesh& mesh);
} // End namespace lumpedPointTools
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View File

@ -0,0 +1,9 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
(cd steady && ./Allclean)
rm -rf movement
rm -rf transient
#------------------------------------------------------------------------------

View File

@ -0,0 +1,92 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# 1) First run steady-state to establish a good initial field.
# 2) Copy the latest state-state results for the transient case,
# but need to copy the pointDisplacement from the 0/ directory
# since it will not have been used for the steady-state case
# 3) Relocate this initial solution to coincide with the first deltaT
# to avoid overwriting the 0/ directory at all.
#
# copyParallelPointDisplacement caseDir timeName
#
# Copy pointDisplacement from caseDir/0/ to caseDir/timeName/
#
copyParallelPointDisplacement()
{
local src=$1
local dstTime=$2
local file=pointDisplacement
[ -d "$src" ] || {
echo "Error: no directory: $src"
return 1
}
# Copy select directories
echo " copy processor '$file' from 0/ -> $dstTime"
if [ -n "$dstTime" ]
then
(
cd $src || exit 1
for proc in processor*
do
[ -d "$proc/0" -a -d "$proc/$dstTime" ] && \
cp $proc/0/$file $proc/$dstTime/$file
done
)
else
echo " no destination time"
fi
# Restart from latestTime
foamDictionary $src/system/controlDict \
-entry startFrom -set latestTime
deltaT=$(foamDictionary $src/system/controlDict -entry deltaT -value)
latestTime=$(foamListTimes -case $src -noZero -latestTime -processor)
# Restart using steady results as first deltaT interval
echo "deltaT=$deltaT latestTime=$latestTime"
if [ -n "$latestTime" -a "$deltaT" != "$latestTime" ]
then
(
cd $src || exit 1
for proc in processor*
do
if [ -d "$proc/$latestTime" -a ! -d "$proc/$deltaT" ]
then
mv $proc/$latestTime $proc/$deltaT
\rm -rf $proc/$deltaT/uniform
fi
done
)
fi
return 0
}
# Do steady-state case
(cd steady && foamRunTutorials)
if ! isTest $@
then
latestTime=$(cd steady && foamListTimes -noZero -latestTime -processor)
# Clone the steady-state case to transient
cloneParallelCase steady transient 0 $latestTime
copyParallelPointDisplacement transient $latestTime
# Do the transient case
\cp files/Allrun.transient transient/Allrun
(cd transient && foamRunTutorials)
fi
#------------------------------------------------------------------------------

View File

@ -0,0 +1,59 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# 1) Run meshing
# 2) Reconstruct
# 3) Test input zones and movement
#
# linkParallelCase srcDir dstDir
#
linkParallelCase()
{
local src=$1
local dst=$2
shift 2
if [ -e "$dst" ]
then
echo "Case already linked: remove case directory $dst prior to linking"
return 1
elif [ ! -d "$src" ]
then
echo "Error: no directory to link: $src"
return 1
fi
echo "Linking $dst parallel case from $src"
mkdir $dst
for i in constant system
do
cp -r $src/$i $dst
done
echo " link processor directories with $# times: $@"
for proc in $(cd $src && \ls -d processor*)
do
( cd $dst && ln -sf ../$src/$proc . )
done
return 0
}
# Do steady-state case
(cd steady && ./Allrun.pre)
if ! isTest $@
then
# Copy/link the steady-state case to movement
linkParallelCase steady movement
# Test movement
\cp files/Allrun.movement movement/Allrun
(cd movement && foamRunTutorials)
fi
#------------------------------------------------------------------------------

View File

@ -0,0 +1,16 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# The 0/ field only
runApplication reconstructPar -withZero -time 0
# Check the location of the pressure zones
# runParallel lumpedPointZones <<- Parallel file writing not yet done
runApplication lumpedPointZones
# Simulated external solver
# Using -scale=1 to see the excessively large movements
runApplication lumpedPointMovement -span 25 -scale 1 ../files/response.txt
#------------------------------------------------------------------------------

View File

@ -0,0 +1,16 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# If OpenFOAM stops prematurely, trigger the external solver to stop
trap '[ -e comms/OpenFOAM.lock ] && echo "status=done" > comms/OpenFOAM.lock' EXIT TERM INT
# Simulated external solver.
# Using -scale since the input movement table is excessively large
runApplication -overwrite \
lumpedPointMovement -scale 0.01 -removeLock -slave ../files/response.txt &
# Run OpenFOAM
runParallel $(getApplication)
#------------------------------------------------------------------------------

View File

@ -0,0 +1,50 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 1 -1 0 0 0 0];
internalField uniform $flowVelocity;
boundaryField
{
#includeEtc "caseDicts/setConstraintTypes"
#include "include/fixedInlet"
outlet
{
type inletOutlet;
inletValue uniform (0 0 0);
value $internalField;
}
// the ground
z_
{
type noSlip;
}
"(?i).*building.*"
{
type noSlip;
}
#include "include/environ"
}
// ************************************************************************* //

View File

@ -0,0 +1,51 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 2 -3 0 0 0 0];
internalField uniform $turbulentEpsilon;
boundaryField
{
#includeEtc "caseDicts/setConstraintTypes"
#include "include/fixedInlet"
outlet
{
type inletOutlet;
inletValue $internalField;
value $internalField;
}
z_
{
type epsilonWallFunction;
value $internalField;
}
"(?i).*building.*"
{
type epsilonWallFunction;
value $internalField;
}
#include "include/environ"
}
// ************************************************************************* //

View File

@ -0,0 +1,24 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
_y
{
type slip;
}
y_
{
type slip;
}
_z
{
type slip;
}
// ************************************************************************* //

View File

@ -0,0 +1,15 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
inlet
{
type fixedValue;
value $internalField;
}
// ************************************************************************* //

View File

@ -0,0 +1,17 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
flowVelocity (100 0 0);
pressure 0;
turbulentKE 37;
turbulentOmega 32;
turbulentEpsilon 30;
#inputMode merge
// ************************************************************************* //

View File

@ -0,0 +1,51 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object k;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 2 -2 0 0 0 0];
internalField uniform $turbulentKE;
boundaryField
{
#includeEtc "caseDicts/setConstraintTypes"
#include "include/fixedInlet"
outlet
{
type inletOutlet;
inletValue $internalField;
value $internalField;
}
z_
{
type kqRWallFunction;
value $internalField;
}
"(?i).*building.*"
{
type kqRWallFunction;
value $internalField;
}
#include "include/environ"
}
// ************************************************************************* //

View File

@ -0,0 +1,44 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object nut;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 2 -1 0 0 0 0];
internalField uniform 0;
boundaryField
{
#includeEtc "caseDicts/setConstraintTypes"
z_
{
type nutkWallFunction;
value uniform 0;
}
"(?i).*building.*"
{
type nutkWallFunction;
value uniform 0;
}
".*"
{
type calculated;
value uniform 0;
}
}
// ************************************************************************* //

View File

@ -0,0 +1,51 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object epsilon;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 0 -1 0 0 0 0];
internalField uniform $turbulentOmega;
boundaryField
{
#includeEtc "caseDicts/setConstraintTypes"
#include "include/fixedInlet"
outlet
{
type inletOutlet;
inletValue $internalField;
value $internalField;
}
z_
{
type omegaWallFunction;
value $internalField;
}
"(?i).*building.*"
{
type omegaWallFunction;
value $internalField;
}
#include "include/environ"
}
// ************************************************************************* //

View File

@ -0,0 +1,51 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "include/initialConditions"
dimensions [0 2 -2 0 0 0 0];
internalField uniform $pressure;
boundaryField
{
#includeEtc "caseDicts/setConstraintTypes"
inlet
{
type zeroGradient;
}
outlet
{
type fixedValue;
value $internalField;
}
z_
{
type zeroGradient;
}
"(?i).*building.*"
{
type zeroGradient;
}
#include "include/environ"
}
// ************************************************************************* //

View File

@ -0,0 +1,73 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class pointVectorField;
object pointDisplacement;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions [0 1 0 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
// Specify directly (without include) for benefit of the paraview reader
processor
{
type processor;
value uniform (0 0 0);
}
inlet
{
type uniformFixedValue;
uniformValue (0 0 0);
}
outlet
{
type uniformFixedValue;
uniformValue (0 0 0);
}
_y
{
type uniformFixedValue;
uniformValue (0 0 0);
}
y_
{
type uniformFixedValue;
uniformValue (0 0 0);
}
_z
{
type uniformFixedValue;
uniformValue (0 0 0);
}
// the ground
z_
{
type slip;
}
"(?i).*building.*"
{
type lumpedPointDisplacement;
value uniform (0 0 0);
}
}
// ************************************************************************* //

View File

@ -0,0 +1,8 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/CleanFunctions
cleanCase
\rm -rf 0
#------------------------------------------------------------------------------

View File

@ -0,0 +1,24 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/RunFunctions
./Allrun.pre
unset parallel
parallel=true
if [ "${parallel:-false}" = false ]
then
# Serial
runApplication simpleFoam
else
# Parallel
runParallel simpleFoam
fi
#------------------------------------------------------------------------------

View File

@ -0,0 +1,45 @@
#!/bin/sh
cd ${0%/*} || exit 1 # Run from this directory
. $WM_PROJECT_DIR/bin/tools/RunFunctions
# Copy building from resources directory
mkdir -p constant/triSurface/
cp $FOAM_TUTORIALS/resources/geometry/building_wtc2.obj constant/triSurface/
# runApplication surfaceFeatureExtract
runApplication blockMesh
\rm -f constant/polyMesh/*Level
unset parallel
parallel=true
# Dummy 0 directory
mkdir 0
if [ "${parallel:-false}" = false ]
then
# Serial
runApplication snappyHexMesh -overwrite
\rm -f constant/polyMesh/refinementHistory*
restore0Dir
runApplication renumberMesh -overwrite
else
# Parallel
runApplication decomposePar -force
runParallel snappyHexMesh -overwrite
\ls -d processor* | xargs -I {} \rm -f ./{}/constant/polyMesh/refinementHistory
restore0Dir -processor
runParallel renumberMesh -overwrite
fi
#------------------------------------------------------------------------------

View File

@ -0,0 +1,31 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object dynamicMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dynamicFvMesh dynamicMotionSolverFvMesh;
motionSolverLibs ("libfvMotionSolvers.so");
solver displacementLaplacian;
displacementLaplacianCoeffs
{
diffusivity inverseDistance ( targetBuilding );
}
// ************************************************************************* //

View File

@ -0,0 +1,21 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
transportModel Newtonian;
nu [0 2 -1 0 0 0 0] 1.5e-05;
// ************************************************************************* //

View File

@ -0,0 +1,29 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object turbulenceProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
simulationType RAS;
RAS
{
// RASModel kOmegaSST;
RASModel kEpsilon;
turbulence on;
printCoeffs on;
}
// ************************************************************************* //

View File

@ -0,0 +1,96 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
convertToMeters 1;
vertices
(
( -1.6065 -1.428 0.0)
( 1.6065 -1.428 0.0)
( 1.6065 1.428 0.0)
( -1.6065 1.428 0.0)
( -1.6065 -1.428 2.021031)
( 1.6065 -1.428 2.021031)
( 1.6065 1.428 2.021031)
( -1.6065 1.428 2.021031)
);
blocks
(
// High resolution
// hex (0 1 2 3 4 5 6 7) (153 136 168) simpleGrading (1 1 5)
// Low resolution
hex (0 1 2 3 4 5 6 7) (77 68 84) simpleGrading (1 1 5)
);
edges
(
);
boundary
(
inlet // -ve X
{
type patch;
faces
(
( 0 4 7 3 )
);
}
outlet // +x X
{
type patch;
faces
(
( 1 2 6 5 )
);
}
y_ // -ve Y
{
type wall;
faces
(
( 0 1 5 4)
);
}
_y // +ve Y
{
type wall;
faces
(
( 3 7 6 2)
);
}
z_ // -ve Z = ground
{
type wall;
faces
(
( 0 3 2 1)
);
}
_z // +ve Z = sky
{
type wall;
faces
(
( 4 5 6 7)
);
}
);
// ************************************************************************* //

View File

@ -0,0 +1,75 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
libs ("liblumpedPointMotion.so");
application pimpleDyMFoam;
startFrom startTime;
startTime 0;
stopAt endTime;
endTime 0.01;
deltaT 1e-4;
writeControl timeStep;
writeInterval 10;
purgeWrite 0;
writeFormat binary;
writePrecision 8;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
adjustTimeStep yes;
// These can be a bit larger when we restart from steady-state
maxCo 0.75;
maxDeltaT 0.01;
// Embed steady-state settings (simpleFoam) without changeDictionary
_simpleFoam
{
endTime 500;
writeInterval 100;
deltaT 1;
adjustTimeStep no;
}
${_${FOAM_EXECUTABLE}};
#remove _simpleFoam
functions
{
}
// ************************************************************************* //

View File

@ -0,0 +1,30 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
note "mesh decomposition control dictionary";
object decomposeParDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
numberOfSubdomains 12;
method scotch;
// method hierarchical;
hierarchicalCoeffs
{
n (3 2 1);
delta 0.001;
order xyz;
}
// ************************************************************************* //

View File

@ -0,0 +1,72 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
ddtSchemes
{
default Euler;
}
gradSchemes
{
default Gauss linear;
grad(U) cellLimited Gauss linear 1;
}
divSchemes
{
default none;
div(phi,U) bounded Gauss linearUpwindV grad(U);
div(phi,k) bounded Gauss upwind;
div(phi,omega) bounded Gauss upwind;
div(phi,epsilon) bounded Gauss upwind;
div((nuEff*dev2(T(grad(U))))) Gauss linear;
}
laplacianSchemes
{
default Gauss linear corrected;
}
interpolationSchemes
{
default linear;
}
snGradSchemes
{
default corrected;
}
wallDist
{
method meshWave;
}
// Embed steady-state settings (simpleFoam) without changeDictionary
_simpleFoam
{
ddtSchemes
{
default steadyState;
}
}
${_${FOAM_EXECUTABLE}};
#remove _simpleFoam
// ************************************************************************* //

View File

@ -0,0 +1,134 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
solvers
{
p
{
solver GAMG;
tolerance 1e-7;
relTol 0.01;
smoother GaussSeidel;
nPreSweeps 0;
nPostSweeps 2;
cacheAgglomeration true;
agglomerator faceAreaPair;
nCellsInCoarsestLevel 10;
mergeLevels 1;
}
"(cellDisplacement)"
{
$p;
tolerance 1e-6;
relTol 0.001;
minIter 1;
}
pFinal
{
$p;
tolerance 1e-6;
relTol 0;
}
cellDisplacementFinal
{
$cellDisplacement;
tolerance 1e-6;
relTol 0;
}
Phi
{
$p;
}
"(U|k|epsilon|omega)"
{
solver smoothSolver;
smoother GaussSeidel;
tolerance 1e-8;
relTol 0.1;
nSweeps 1;
}
"(U|k|epsilon|omega)Final"
{
$U;
relTol 0;
}
}
SIMPLE
{
nNonOrthogonalCorrectors 0;
consistent yes;
pRefCell 0;
pRefValue 0;
}
PIMPLE
{
nOuterCorrectors 2;
nCorrectors 1;
nNonOrthogonalCorrectors 0;
pRefCell 0;
pRefValue 0;
}
potentialFlow
{
nNonOrthogonalCorrectors 10;
}
cache
{
grad(U);
}
// Default relaxation for transient
relaxationFactors
{
equations
{
".*" 1;
}
}
// Embed steady-state settings (simpleFoam) without changeDictionary
_simpleFoam
{
fields
{
p 1.0;
}
equations
{
"(U|k|epsilon)" 0.9;
}
}
relaxationFactors
{
${_${FOAM_EXECUTABLE}};
}
#remove _simpleFoam
// ************************************************************************* //

View File

@ -0,0 +1,75 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "system";
object lumpedPointMovement;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Reference axis for the locations
axis (0 0 1);
// Locations of the lumped points
locations 11(0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5);
// Division for pressure forces (0-1)
division 0.5;
//- If present, the offset of patch points compared to the locations
// Otherwise determined from the bounding box
// centre (0 0 0);
//- The interpolation scheme
interpolationScheme linear;
//- Relaxation/scaling factor when updating positions
relax 1.0;
forces
{
//- The pressure name (default: p)
p p;
//- Reference pressure [Pa] (default: 0)
pRef 0;
//- Reference density for incompressible calculations (default: 1)
rhoRef 1;
}
communication
{
commsDir "comms";
log on;
waitInterval 1;
timeOut 100;
initByExternal false;
// Input file of positions/rotation, written by external application
inputName positions.in;
// Output file of forces, written by OpenFOAM
outputName forces.out;
inputFormat dictionary;
outputFormat dictionary;
debugTable "$FOAM_CASE/output.txt";
}
// ************************************************************************* //

View File

@ -0,0 +1,24 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object meshQualityDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Include defaults parameters from master dictionary
#includeEtc "caseDicts/meshQualityDict"
//- minFaceWeight (0 -> 0.5)
minFaceWeight 0.02;
// ************************************************************************* //

View File

@ -0,0 +1,305 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object snappyHexMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
// Which of the steps to run
castellatedMesh true;
snap true;
addLayers false;
// Geometry. Definition of all surfaces. All surfaces are of class
// searchableSurface.
// Surfaces are used
// - to specify refinement for any mesh cell intersecting it
// - to specify refinement for any mesh cell inside/outside/near
// - to 'snap' the mesh boundary to the surface
geometry
{
building_wtc2.obj
{
type triSurfaceMesh;
name building;
}
};
// Settings for the castellatedMesh generation.
castellatedMeshControls
{
// Refinement parameters
// ~~~~~~~~~~~~~~~~~~~~~
// If local number of cells is >= maxLocalCells on any processor
// switches from from refinement followed by balancing
// (current method) to (weighted) balancing before refinement.
maxLocalCells 100000;
// Overall cell limit (approximately). Refinement will stop immediately
// upon reaching this number so a refinement level might not complete.
// Note that this is the number of cells before removing the part which
// is not 'visible' from the keepPoint. The final number of cells might
// actually be a lot less.
maxGlobalCells 2000000;
// The surface refinement loop might spend lots of iterations refining just a
// few cells. This setting will cause refinement to stop if <= minimumRefine
// are selected for refinement. Note: it will at least do one iteration
// (unless the number of cells to refine is 0)
minRefinementCells 10;
// Allow a certain level of imbalance during refining
// (since balancing is quite expensive)
// Expressed as fraction of perfect balance (= overall number of cells /
// nProcs). 0=balance always.
maxLoadUnbalance 0.10;
// Number of buffer layers between different levels.
// 1 means normal 2:1 refinement restriction, larger means slower
// refinement.
nCellsBetweenLevels 15;
// Explicit feature edge refinement
// ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
// Specifies a level for any cell intersected by its edges.
// This is a featureEdgeMesh, read from constant/triSurface for now.
features
(
);
// Surface based refinement
// ~~~~~~~~~~~~~~~~~~~~~~~~
// Specifies two levels for every surface. The first is the minimum level,
// every cell intersecting a surface gets refined up to the minimum level.
// The second level is the maximum level. Cells that 'see' multiple
// intersections where the intersections make an
// angle > resolveFeatureAngle get refined up to the maximum level.
refinementSurfaces
{
building
{
// Surface-wise min and max refinement level
level (2 2);
// Optional specification of patch type (default is wall). No
// constraint types (cyclic, symmetry) etc. are allowed.
patchInfo
{
type wall;
inGroups (targetBuilding);
}
}
}
// Resolve sharp angles
resolveFeatureAngle 30;
// Region-wise refinement
// ~~~~~~~~~~~~~~~~~~~~~~
// Specifies refinement level for cells in relation to a surface. One of
// three modes
// - distance. 'levels' specifies per distance to the surface the
// wanted refinement level. The distances need to be specified in
// descending order.
// - inside. 'levels' is only one entry and only the level is used. All
// cells inside the surface get refined up to the level. The surface
// needs to be closed for this to be possible.
// - outside. Same but cells outside.
refinementRegions
{
}
// Mesh selection
// ~~~~~~~~~~~~~~
// After refinement patches get added for all refinementSurfaces and
// all cells intersecting the surfaces get put into these patches. The
// section reachable from the locationInMesh is kept.
// NOTE: This point should never be on a face, always inside a cell, even
// after refinement.
locationInMesh (0 0 1.99);
// Whether any faceZones (as specified in the refinementSurfaces)
// are only on the boundary of corresponding cellZones or also allow
// free-standing zone faces. Not used if there are no faceZones.
allowFreeStandingZoneFaces true;
}
// Settings for the snapping.
snapControls
{
//- Number of patch smoothing iterations before finding correspondence
// to surface
nSmoothPatch 3;
//- Relative distance for points to be attracted by surface feature point
// or edge. True distance is this factor times local
// maximum edge length.
tolerance 2.0;
//- Number of mesh displacement relaxation iterations.
nSolveIter 30;
//- Maximum number of snapping relaxation iterations. Should stop
// before upon reaching a correct mesh.
nRelaxIter 5;
// Feature snapping
//- Number of feature edge snapping iterations.
// Leave out altogether to disable.
nFeatureSnapIter 10;
//- Detect (geometric only) features by sampling the surface
// (default=false).
implicitFeatureSnap false;
//- Use castellatedMeshControls::features (default = true)
explicitFeatureSnap true;
//- Detect points on multiple surfaces (only for explicitFeatureSnap)
multiRegionFeatureSnap false;
}
// Settings for the layer addition.
addLayersControls
{
// Are the thickness parameters below relative to the undistorted
// size of the refined cell outside layer (true) or absolute sizes (false).
relativeSizes true;
// Per final patch (so not geometry!) the layer information
layers
{
}
// Expansion factor for layer mesh
expansionRatio 1.0;
// Wanted thickness of final added cell layer. If multiple layers
// is the thickness of the layer furthest away from the wall.
// Relative to undistorted size of cell outside layer.
// See relativeSizes parameter.
finalLayerThickness 0.3;
// Minimum thickness of cell layer. If for any reason layer
// cannot be above minThickness do not add layer.
// Relative to undistorted size of cell outside layer.
minThickness 0.1;
// If points get not extruded do nGrow layers of connected faces that are
// also not grown. This helps convergence of the layer addition process
// close to features.
// Note: changed(corrected) w.r.t 17x! (didn't do anything in 17x)
nGrow 0;
// Advanced settings
// When not to extrude surface. 0 is flat surface, 90 is when two faces
// are perpendicular
featureAngle 60;
// At non-patched sides allow mesh to slip if extrusion direction makes
// angle larger than slipFeatureAngle.
slipFeatureAngle 30;
// Maximum number of snapping relaxation iterations. Should stop
// before upon reaching a correct mesh.
nRelaxIter 3;
// Number of smoothing iterations of surface normals
nSmoothSurfaceNormals 1;
// Number of smoothing iterations of interior mesh movement direction
nSmoothNormals 3;
// Smooth layer thickness over surface patches
nSmoothThickness 10;
// Stop layer growth on highly warped cells
maxFaceThicknessRatio 0.5;
// Reduce layer growth where ratio thickness to medial
// distance is large
maxThicknessToMedialRatio 0.3;
// Angle used to pick up medial axis points
// Note: changed(corrected) w.r.t 17x! 90 degrees corresponds to 130 in 17x.
minMedianAxisAngle 90;
// Create buffer region for new layer terminations
nBufferCellsNoExtrude 0;
// Overall max number of layer addition iterations. The mesher will exit
// if it reaches this number of iterations; possibly with an illegal
// mesh.
nLayerIter 50;
}
// Generic mesh quality settings. At any undoable phase these determine
// where to undo.
meshQualityControls
{
#include "meshQualityDict"
// Advanced
//- Number of error distribution iterations
nSmoothScale 4;
//- Amount to scale back displacement at error points
errorReduction 0.75;
}
// Advanced
// Write flags
writeFlags
(
scalarLevels
layerSets
layerFields // write volScalarField for layer coverage
);
// Merge tolerance. Is fraction of overall bounding box of initial mesh.
// Note: the write tolerance needs to be higher than this.
mergeTolerance 1e-6;
// ************************************************************************* //

View File

@ -0,0 +1,33 @@
/*--------------------------------*- C++ -*----------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: plus |
| \\ / A nd | Web: www.OpenFOAM.com |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object surfaceFeatureExtractDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
building_wtc2.obj
{
// How to obtain raw features (extractFromFile || extractFromSurface)
extractionMethod extractFromSurface;
// Mark edges whose adjacent surface normals are at an angle less
// than includedAngle as features
// - 0 : selects no edges
// - 180: selects all edges
includedAngle 85;
// Write features to obj format for postprocessing
writeObj false;
}
// ************************************************************************* //

View File

@ -32,7 +32,7 @@ boundaryField
wall
{
type nutkRoughWallFunction;
type nutkWallFunction;
Ks uniform 100e-6;
Cs uniform 0.5;
value $internalField;

View File

@ -33,6 +33,7 @@ boundaryField
wall
{
type omegaWallFunction;
blended true;
value $internalField;
}
}

View File

@ -28,7 +28,7 @@ startTime 0;
stopAt endTime;
endTime 2;
endTime 0.6;
deltaT 0.001;

View File

@ -37,7 +37,7 @@ divSchemes
laplacianSchemes
{
default Gauss linear corrected;
default Gauss linear limited corrected 0.333;
}
interpolationSchemes
@ -47,7 +47,7 @@ interpolationSchemes
snGradSchemes
{
default corrected;
default limited corrected 0.333;
}
wallDist

View File

@ -73,10 +73,8 @@ PIMPLE
{
momentumPredictor yes;
nCorrectors 3;
nNonOrthogonalCorrectors 0;
pRefPoint (0.51 0.51 0.51);
pRefValue 0;
nOuterCorrectors 1;
nNonOrthogonalCorrectors 1;
}

View File

@ -0,0 +1,319 @@
# Wavefront OBJ file
# source https://grabcad.com/library/1-world-trade-center-1
#
# Terms of Usage:
# One World Trade Center by CaitlinLe is licensed under the
# Creative Commons - Attribution - Share Alike license.
#
o targetBuilding_wtc2
# points : 101
# faces : 198
# zones : 1
# 0 building (nFaces: 198)
# <points count="101">
v 0.21015 0.243555 7.5e-06
v 0.21015 0.056445 7.5e-06
v 0.38985 0.056445 7.5e-06
v 0.38985 0.243555 7.5e-06
v 0.2109 0.0616875 0.111383
v 0.21522 0.0571875 0.111383
v 0.2109 0.239858 0.111383
v 0.3891 0.0601425 0.111383
v 0.386258 0.0571875 0.111383
v 0.213735 0.242812 0.111383
v 0.3891 0.242812 0.111383
v 0.21015 0.056445 0.111383
v 0.21015 0.243555 0.111383
v 0.38985 0.056445 0.111383
v 0.38985 0.243555 0.111383
v 0.3891 0.242812 0.119212
v 0.2109 0.150773 0.728745
v 0.300742 0.0571875 0.72873
v 0.262478 0.11247 0.742507
v 0.255555 0.120915 0.742507
v 0.270915 0.105547 0.742507
v 0.250403 0.130545 0.742507
v 0.280545 0.100395 0.742507
v 0.290992 0.0972375 0.742507
v 0.247237 0.140992 0.742507
v 0.246173 0.151853 0.742507
v 0.30186 0.0961575 0.742507
v 0.247237 0.162712 0.742507
v 0.31272 0.0972375 0.742507
v 0.250403 0.17316 0.742507
v 0.323168 0.100395 0.742507
v 0.255555 0.18279 0.742507
v 0.332798 0.105547 0.742507
v 0.262478 0.191235 0.742507
v 0.341235 0.11247 0.742507
v 0.348157 0.120915 0.742507
v 0.270915 0.198157 0.742507
v 0.280545 0.20331 0.742507
v 0.35331 0.130545 0.742507
v 0.290992 0.206467 0.742507
v 0.356475 0.140992 0.742507
v 0.35754 0.151853 0.742507
v 0.30186 0.20754 0.742507
v 0.31272 0.206467 0.742507
v 0.356475 0.162712 0.742507
v 0.35331 0.17316 0.742507
v 0.323168 0.20331 0.742507
v 0.348157 0.18279 0.742507
v 0.332798 0.198157 0.742507
v 0.341235 0.191235 0.742507
v 0.212805 0.150773 0.742507
v 0.3891 0.151207 0.742507
v 0.3891 0.152872 0.742507
v 0.300742 0.059175 0.742507
v 0.301163 0.242812 0.742507
v 0.30276 0.242812 0.742507
v 0.3 0.149993 0.764783
v 0.301853 0.150495 0.764783
v 0.298148 0.150495 0.764783
v 0.303218 0.151853 0.764783
v 0.296783 0.151853 0.764783
v 0.30372 0.153713 0.764783
v 0.29628 0.153713 0.764783
v 0.303218 0.155565 0.764783
v 0.296783 0.155565 0.764783
v 0.301853 0.156922 0.764783
v 0.298148 0.156922 0.764783
v 0.3 0.157425 0.764783
v 0.262478 0.11247 0.764783
v 0.255555 0.120915 0.764783
v 0.270915 0.105547 0.764783
v 0.250403 0.130545 0.764783
v 0.280545 0.100395 0.764783
v 0.247237 0.140992 0.764783
v 0.290992 0.0972375 0.764783
v 0.246173 0.151853 0.764783
v 0.30186 0.0961575 0.764783
v 0.247237 0.162712 0.764783
v 0.31272 0.0972375 0.764783
v 0.250403 0.17316 0.764783
v 0.323168 0.100395 0.764783
v 0.255555 0.18279 0.764783
v 0.332798 0.105547 0.764783
v 0.262478 0.191235 0.764783
v 0.341235 0.11247 0.764783
v 0.348157 0.120915 0.764783
v 0.270915 0.198157 0.764783
v 0.35331 0.130545 0.764783
v 0.280545 0.20331 0.764783
v 0.290992 0.206467 0.764783
v 0.356475 0.140992 0.764783
v 0.35754 0.151853 0.764783
v 0.30186 0.20754 0.764783
v 0.31272 0.206467 0.764783
v 0.356475 0.162712 0.764783
v 0.35331 0.17316 0.764783
v 0.323168 0.20331 0.764783
v 0.348157 0.18279 0.764783
v 0.332798 0.198157 0.764783
v 0.341235 0.191235 0.764783
v 0.3 0.153713 0.95784
# </points>
# <faces count="198">
g building
f 1 2 12
f 13 1 12
f 12 2 14
f 3 14 2
f 3 2 4
f 2 1 4
f 4 15 14
f 4 14 3
f 15 4 13
f 1 13 4
f 13 7 10
f 44 94 97
f 44 97 47
f 13 10 15
f 11 15 10
f 44 43 93
f 44 93 94
f 11 10 16
f 47 56 44
f 11 8 15
f 14 15 8
f 14 8 9
f 11 16 8
f 46 96 95
f 46 95 45
f 52 53 42
f 17 7 5
f 95 92 42
f 95 42 45
f 18 9 8
f 98 64 96
f 95 96 62
f 5 12 6
f 9 6 14
f 41 52 42
f 101 64 66
f 92 41 42
f 60 62 101
f 62 64 101
f 12 14 6
f 12 5 13
f 7 13 5
f 65 78 63
f 63 74 61
f 72 70 61
f 67 84 65
f 55 43 56
f 46 45 53
f 45 42 53
f 16 56 53
f 26 76 78
f 26 78 28
f 78 80 30
f 78 30 28
f 16 53 8
f 52 8 53
f 101 66 68
f 28 51 26
f 28 30 51
f 101 58 60
f 48 53 50
f 56 50 53
f 50 100 98
f 50 98 48
f 46 53 48
f 46 48 96
f 98 96 48
f 83 58 81
f 81 58 79
f 85 58 83
f 99 66 100
f 100 64 98
f 50 56 49
f 97 94 66
f 64 100 66
f 62 96 64
f 91 92 62
f 97 66 99
f 99 49 97
f 62 60 91
f 58 85 60
f 88 60 86
f 88 91 60
f 50 49 99
f 50 99 100
f 90 68 93
f 101 68 67
f 56 47 49
f 49 47 97
f 38 40 55
f 77 79 57
f 86 60 85
f 92 95 62
f 57 79 58
f 91 41 92
f 88 39 41
f 88 41 91
f 86 36 39
f 86 39 88
f 85 36 86
f 41 39 52
f 35 52 36
f 39 36 52
f 94 93 68
f 89 67 90
f 87 67 89
f 84 82 65
f 82 80 65
f 78 65 80
f 84 67 87
f 87 37 34
f 87 34 84
f 18 8 52
f 34 32 82
f 34 82 84
f 89 38 37
f 89 37 87
f 57 58 101
f 90 40 38
f 90 38 89
f 19 69 70
f 24 77 75
f 71 73 59
f 57 59 73
f 75 57 73
f 77 57 75
f 73 23 24
f 73 24 75
f 71 21 23
f 71 23 73
f 69 19 21
f 69 21 71
f 23 54 24
f 54 21 19
f 51 25 26
f 51 55 17
f 55 10 7
f 17 55 7
f 56 16 55
f 10 55 16
f 5 6 18
f 54 51 17
f 18 54 17
f 17 5 18
f 18 6 9
f 52 35 54
f 35 33 54
f 52 54 18
f 36 85 35
f 33 31 54
f 83 33 35
f 83 35 85
f 31 29 54
f 81 31 33
f 81 33 83
f 79 29 31
f 79 31 81
f 77 27 29
f 77 29 79
f 77 24 27
f 27 24 54
f 54 29 27
f 72 61 74
f 70 69 61
f 59 61 69
f 69 71 59
f 20 22 51
f 19 70 20
f 25 51 22
f 74 25 22
f 74 22 72
f 74 76 25
f 26 25 76
f 22 20 70
f 22 70 72
f 23 21 54
f 54 19 51
f 20 51 19
f 93 43 40
f 93 40 90
f 80 82 30
f 32 30 82
f 37 38 55
f 34 37 55
f 32 34 51
f 55 51 34
f 30 32 51
f 44 56 43
f 43 55 40
f 101 65 63
f 101 63 61
f 59 101 61
f 101 67 65
f 59 57 101
f 68 90 67
f 78 76 63
f 68 66 94
f 76 74 63
# </faces>