zhyang-dev/OpenFOAM-12
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Function1: Implemented integral evaluations

Browse Source 
Integral evaluations have been implemented for all the ramp function1-s,
as well as the sine and square wave. Bounds handling has also been added
to the integration of table-type functions.

In addition, the sine wave "t0" paramater has been renamed "start" for
consistency with the ramp functions.
This commit is contained in:
Will Bainbridge
2019-10-17 19:11:05 +01:00
parent ef04eed403
commit dcf4d0c505
54 changed files with 958 additions and 768 deletions
Download Patch File Download Diff File Expand all files Collapse all files

139
applications/test/Function1/Test-Function1.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -29,114 +29,61 @@ Description
\*---------------------------------------------------------------------------*/
#include "fvCFD.H"
#include "Function1.H"
#include "IOdictionary.H"
#include "linearInterpolationWeights.H"
#include "splineInterpolationWeights.H"
#include "IFstream.H"
#include "OFstream.H"
#include "ListOps.H"
#include "argList.H"
using namespace Foam;
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
int main(int argc, char *argv[])
{
#include "setRootCase.H"
#include "createTime.H"
#include "createMesh.H"
argList::validArgs.append("dictionary");
argList args(argc, argv);
const word dictName(args[1]);
Info<< "Reading " << dictName << nl << endl;
const dictionary dict = IFstream(dictName)();
{
scalarField samples(4);
samples[0] = 0;
samples[1] = 1;
samples[2] = 2;
samples[3] = 3;
scalarField values(4);
values = 1.0;
// values[0] = 0.0;
// values[1] = 1.0;
Info<< "Constructing the function\n" << endl;
const autoPtr<Function1<scalar>> functionPtr =
Function1<scalar>::New("function", dict);
const Function1<scalar>& function = functionPtr();
linearInterpolationWeights interpolator
// splineInterpolationWeights interpolator
const scalar x0 = readScalar(dict.lookup("x0"));
const scalar x1 = readScalar(dict.lookup("x1"));
const label nX = readLabel(dict.lookup("nX"));
const scalar dx = (x1 - x0)/(nX - 1);
const scalarField xs
(
samples
);
labelList indices;
scalarField weights;
interpolator.integrationWeights(1.1, 1.2, indices, weights);
Pout<< "indices:" << indices << endl;
Pout<< "weights:" << weights << endl;
scalar baseSum = interpolator.weightedSum
(
weights,
UIndirectList<scalar>(values, indices)
);
Pout<< "baseSum=" << baseSum << nl << nl << endl;
// interpolator.integrationWeights(-0.01, 0, indices, weights);
// scalar partialSum = interpolator.weightedSum
// (
// weights,
// UIndirectList<scalar>(values, indices)
// );
// Pout<< "partialSum=" << partialSum << nl << nl << endl;
//
//
// interpolator.integrationWeights(-0.01, 1, indices, weights);
// // Pout<< "samples:" << samples << endl;
// // Pout<< "indices:" << indices << endl;
// // Pout<< "weights:" << weights << endl;
// scalar sum = interpolator.weightedSum
// (
// weights,
// UIndirectList<scalar>(values, indices)
// );
// Pout<< "integrand=" << sum << nl << nl << endl;
return 1;
}
IOdictionary function1Properties
(
IOobject
(
"function1Properties",
runTime.constant(),
mesh,
IOobject::MUST_READ_IF_MODIFIED,
IOobject::NO_WRITE
)
x0 + dx*List<scalar>(identity(nX))
);
autoPtr<Function1<scalar>> function1
(
Function1<scalar>::New
(
"function1",
function1Properties
)
);
Info<< "Calculating values\n" << endl;
const scalarField ys(function.value(xs));
const scalarField integralYs(function.integrate(scalarField(nX, x0), xs));
scalarField trapezoidIntegralYs(nX, 0);
for (label i = 1; i < nX; ++ i)
{
trapezoidIntegralYs[i] =
trapezoidIntegralYs[i-1] + (ys[i] + ys[i-1])/2*dx;
}
scalar x0 = readScalar(function1Properties.lookup("x0"));
scalar x1 = readScalar(function1Properties.lookup("x1"));
OFstream dataFile(dictName + ".dat");
Info<< "Writing to " << dataFile.name() << "\n" << endl;
dataFile<< "# x y integralY trapezoidIntegralY" << endl;
forAll(xs, i)
{
dataFile
<< xs[i] << ' '
<< ys[i] << ' '
<< integralYs[i] << ' '
<< trapezoidIntegralYs[i] << endl;
}
Info<< "Data entry type: " << function1().type() << nl << endl;
Info<< "Inputs" << nl
<< " x0 = " << x0 << nl
<< " x1 = " << x1 << nl
<< endl;
Info<< "Interpolation" << nl
<< " f(x0) = " << function1().value(x0) << nl
<< " f(x1) = " << function1().value(x1) << nl
<< endl;
Info<< "Integration" << nl
<< " int(f(x)) lim(x0->x1) = " << function1().integrate(x0, x1) << nl
<< endl;
Info<< "End\n" << endl;
return 0;
}

25
applications/test/Function1/function1Properties
View File

@ -1,25 +0,0 @@
/*--------------------------------*- C++ -*----------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Version: dev
\\/ M anipulation |
\*---------------------------------------------------------------------------*/
FoamFile
{
version 2.0;
format ascii;
class dictionary;
location "constant";
object function1Properties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
x0 0.5;
x1 1;
function1 table ((0 0)(10 1));
// ************************************************************************* //

8
applications/test/Function1/linearRamp Normal file
View File

@ -0,0 +1,8 @@
function linearRamp;
start 1;
duration 3;
x0 0;
x1 5;
nX 51;

8
applications/test/Function1/quadraticRamp Normal file
View File

@ -0,0 +1,8 @@
function quadraticRamp;
start 1;
duration 3;
x0 0;
x1 5;
nX 51;

8
applications/test/Function1/quarterCosineRamp Normal file
View File

@ -0,0 +1,8 @@
function quarterCosineRamp;
start 1;
duration 3;
x0 0;
x1 5;
nX 51;

8
applications/test/Function1/quarterSineRamp Normal file
View File

@ -0,0 +1,8 @@
function quarterSineRamp;
start 1;
duration 3;
x0 0;
x1 5;
nX 51;

10
applications/test/Function1/reverseQuarterSineRamp Normal file
View File

@ -0,0 +1,10 @@
function reverseRamp;
ramp quarterSineRamp;
start 1;
duration 3;
x0 0;
x1 5;
nX 51;

13
applications/test/Function1/scaleHalfCosineRamp Normal file
View File

@ -0,0 +1,13 @@
function scale;
scale 5;
xScale 1.5;
value halfCosineRamp;
start 1;
duration 3;
x0 0;
x1 5;
nX 51;

10
applications/test/Function1/sine Normal file
View File

@ -0,0 +1,10 @@
function sine;
amplitude 2;
frequency 0.5;
phase -3.14;
level 1;
x0 0;
x1 5;
nX 51;

11
applications/test/Function1/square Normal file
View File

@ -0,0 +1,11 @@
function square;
amplitude 2;
frequency 0.6667;
phase 0;
level 1;
markSpace 0.25;
x0 0;
x1 5;
nX 101;

15
applications/test/Function1/table Normal file
View File

@ -0,0 +1,15 @@
function table
(
(0 0.1)
(0.4 0.31)
(0.8 0.59)
(1.2 0.81)
(1.6 0.95)
(2 1)
);
outOfBounds repeat;
x0 -0.5;
x1 2.5;
nX 50;

9
applications/test/rigidBodyDynamics/Test-rigidBodyDynamics.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2018-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -70,10 +70,11 @@ int main(int argc, char *argv[])
}
}
}
Info<< endl;
// Set up motion fields
const label nDoF = motion.nDoF();
Info << nDoF << " degrees of freedom" << endl;
Info<< nDoF << " degrees of freedom" << nl << endl;
scalarField tau(nDoF, Zero);
Field<spatialVector> fx(motion.nBodies(), Zero);
@ -142,7 +143,7 @@ int main(int argc, char *argv[])
animationFile<< endl;
}
}
Info << endl;
Info<< nl << endl;
// Expand the bound box
{
@ -179,6 +180,8 @@ int main(int argc, char *argv[])
}
animationFile<< endl << " pause " << deltaT << endl << "}" << endl;
Info<< "End" << nl << endl;
return 0;
}

1
src/OpenFOAM/Make/files
View File

@ -97,6 +97,7 @@ primitives/functions/Function1/quadraticRamp/quadraticRamp.C
primitives/functions/Function1/quarterSineRamp/quarterSineRamp.C
primitives/functions/Function1/quarterCosineRamp/quarterCosineRamp.C
primitives/functions/Function1/halfCosineRamp/halfCosineRamp.C
primitives/functions/Function1/Table/tableBase.C
primitives/subModelBase/subModelBase.C

8
src/OpenFOAM/db/Time/TimeState.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -61,4 +61,10 @@ Foam::scalar Foam::TimeState::timeToUserTime(const scalar t) const
}
Foam::scalar Foam::TimeState::timeToUserTimeRatio() const
{
return 1;
}
// ************************************************************************* //

3
src/OpenFOAM/db/Time/TimeState.H
View File

@ -83,6 +83,9 @@ public:
//- Convert the real-time (s) into user-time (e.g. CA deg)
virtual scalar timeToUserTime(const scalar t) const;
//- Ratio between real-time and user-time
virtual scalar timeToUserTimeRatio() const;
//- Return current time value
inline scalar timeOutputValue() const;

7
src/OpenFOAM/interpolations/interpolationWeights/interpolationWeights/interpolationWeights.H
View File

@ -55,11 +55,14 @@ class objectRegistry;
class interpolationWeights
{
protected:
// Protected Member Data
//- Cached samples
const scalarField& samples_;
public:
//- Runtime type information
@ -124,7 +127,7 @@ public:
scalarField& weights
) const = 0;
//- Helper: weighted sum
//- Weighted sum helper method
template<class ListType1, class ListType2>
static typename outerProduct
<

231
src/OpenFOAM/interpolations/interpolationWeights/linearInterpolationWeights/linearInterpolationWeights.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2012-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2012-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -44,32 +44,6 @@ addToRunTimeSelectionTable
);
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
Foam::Pair<Foam::scalar> linearInterpolationWeights::integrationWeights
(
const label i,
const scalar t
) const
{
// t is in range samples_[i] .. samples_[i+1]
scalar s = (t-samples_[i])/(samples_[i+1]-samples_[i]);
if (s < -small || s > 1+small)
{
FatalErrorInFunction
<< "Value " << t << " outside range " << samples_[i]
<< " .. " << samples_[i+1]
<< exit(FatalError);
}
scalar d = samples_[i+1]-t;
return Pair<scalar>(d*0.5*(1-s), d*0.5*(1+s));
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
linearInterpolationWeights::linearInterpolationWeights
@ -91,157 +65,174 @@ bool linearInterpolationWeights::valueWeights
scalarField& weights
) const
{
bool indexChanged = false;
// Check if current timeIndex is still valid
// Check if current index is still valid
bool changed = false;
if
(
index_ >= 0
&& index_ < samples_.size()
&& (
samples_[index_] <= t
&& (index_ == samples_.size()-1 || t <= samples_[index_+1])
(
index_ == -1
&& t <= samples_.first()
)
|| (
index_ >= 0
&& index_ < samples_.size() - 1
&& t >= samples_[index_]
&& t <= samples_[index_ + 1]
)
|| (
index_ == samples_.size() - 1
&& t >= samples_.last()
)
)
{
// index_ still at correct slot
// The index is still in the correct slot
}
else
{
// search for correct index
// The index is no longer in the correct slot, so search for a new one
index_ = findLower(samples_, t);
indexChanged = true;
changed = true;
}
// Calculate the number of indices and weights and resize the result
const label n = index_ == -1 || index_ == samples_.size() - 1 ? 1 : 2;
indices.resize(n);
weights.resize(n);
// Compute the value
if (index_ == -1)
{
// Use first element only
indices.setSize(1);
weights.setSize(1);
// Use the first value
indices[0] = 0;
weights[0] = 1.0;
weights[0] = 1;
}
else if (index_ == samples_.size()-1)
else if (index_ == samples_.size() - 1)
{
// Use last element only
indices.setSize(1);
weights.setSize(1);
indices[0] = samples_.size()-1;
weights[0] = 1.0;
// Use the last value
indices[0] = samples_.size() - 1;
weights[0] = 1;
}
else
{
// Interpolate
indices.setSize(2);
weights.setSize(2);
// Interpolate within the interval
const scalar f =
(t - samples_[index_])/(samples_[index_ + 1] - samples_[index_]);
indices[0] = index_;
indices[1] = index_+1;
scalar t0 = samples_[indices[0]];
scalar t1 = samples_[indices[1]];
scalar deltaT = t1-t0;
weights[0] = (t1-t)/deltaT;
weights[1] = 1.0-weights[0];
weights[0] = 1 - f;
indices[1] = index_ + 1;
weights[1] = f;
}
return indexChanged;
return changed;
}
bool linearInterpolationWeights::integrationWeights
(
const scalar t1,
const scalar t2,
scalar t1,
scalar t2,
labelList& indices,
scalarField& weights
) const
{
if (t2 < t1-vSmall)
// If the arguments are in descending order, then swap them and set the
// weights' sign negative
label sign = +1;
if (t1 > t2)
{
FatalErrorInFunction
<< "Integration should be in positive direction."
<< " t1:" << t1 << " t2:" << t2
<< exit(FatalError);
Swap(t1, t2);
sign = -1;
}
// Currently no fancy logic on cached index like in value
//- Search for lower indices
// Note: currently there is no caching of this search like in valueWeights
const label i1 = findLower(samples_, t1);
const label i2 = findLower(samples_, t2);
const label iClip1 = min(max(i1, 0), samples_.size() - 2);
const label iClip2 = min(max(i2, 0), samples_.size() - 2);
const label n = max(i2 - i1 + (i1 == iClip1) + (i2 == iClip2), 1);
//- Find lower or equal index
label i1 = findLower(samples_, t1, 0, lessEqOp<scalar>());
//- Find lower index
label i2 = findLower(samples_, t2);
// For now just fail if any outside table
if (i1 == -1 || i2 == samples_.size()-1)
// Check if anything changed
bool changed = false;
if (indices.size() != n)
{
FatalErrorInFunction
<< "Integrating outside table " << samples_[0] << ".."
<< samples_.last() << " not implemented."
<< " t1:" << t1 << " t2:" << t2 << exit(FatalError);
}
label nIndices = i2-i1+2;
// Determine if indices already correct
bool anyChanged = false;
if (nIndices != indices.size())
{
anyChanged = true;
changed = true;
}
else
{
// Closer check
label index = i1;
forAll(indices, i)
forAll(indices, indexi)
{
if (indices[i] != index)
if (indices[indexi] == indexi + max(i1, 0))
{
anyChanged = true;
changed = true;
break;
}
index++;
}
}
indices.setSize(nIndices);
weights.setSize(nIndices);
weights = 0.0;
// Resize the result arrays
indices.resize(n);
indices = -1;
weights.resize(n);
weights = 0;
// Sum from i1+1 to i2+1
for (label i = i1+1; i <= i2; i++)
// Add out of bounds interval below the table
if (i1 == -1)
{
scalar d = samples_[i+1]-samples_[i];
indices[i-i1] = i;
weights[i-i1] += 0.5*d;
indices[i+1-i1] = i+1;
weights[i+1-i1] += 0.5*d;
indices[0] = 0;
weights[0] += sign*(samples_[0] - t1);
}
if (i2 == -1)
{
indices[0] = 0;
weights[0] -= sign*(samples_[0] - t2);
}
// Add from i1 to t1
// Add partial interval from t1 to i1 + 1
if (i1 == iClip1)
{
Pair<scalar> i1Tot1 = integrationWeights(i1, t1);
const scalar f = (t1 - samples_[i1])/(samples_[i1 + 1] - samples_[i1]);
const scalar d = samples_[i1 + 1] - t1;
indices[0] = i1;
weights[0] += i1Tot1.first();
indices[1] = i1+1;
weights[1] += i1Tot1.second();
weights[0] += sign*(1 - f)*d/2;
indices[1] = i1 + 1;
weights[1] += sign*(1 + f)*d/2;
}
// Subtract from t2 to i2+1
// Sum whole intervals from i1 + 1 to i2
if (i1 != i2) for (label i = i1 + 1; i <= iClip2; i ++)
{
Pair<scalar> wghts = integrationWeights(i2, t2);
indices[i2-i1] = i2;
weights[i2-i1] += -wghts.first();
indices[i2-i1+1] = i2+1;
weights[i2-i1+1] += -wghts.second();
const scalar d = samples_[i + 1] - samples_[i];
indices[i - iClip1] = i;
weights[i - iClip1] += sign*d/2;
indices[i - iClip1 + 1] = i + 1;
weights[i - iClip1 + 1] += sign*d/2;
}
return anyChanged;
// Subtract partial interval from t2 to i2 + 1
if (i2 == iClip2)
{
const scalar f = (t2 - samples_[i2])/(samples_[i2 + 1] - samples_[i2]);
const scalar d = samples_[i2 + 1] - t2;
indices[n - 2] = i2;
weights[n - 2] -= sign*(1 - f)*d/2;
indices[n - 1] = i2 + 1;
weights[n - 1] -= sign*(1 + f)*d/2;
}
// Add out of bounds interval above the table
if (i1 == samples_.size() - 1)
{
indices[n - 1] = samples_.size() - 1;
weights[n - 1] -= sign*(t1 - samples_.last());
}
if (i2 == samples_.size() - 1)
{
indices[n - 1] = samples_.size() - 1;
weights[n - 1] += sign*(t2 - samples_.last());
}
return changed;
}

14
src/OpenFOAM/interpolations/interpolationWeights/linearInterpolationWeights/linearInterpolationWeights.H
View File

@ -54,21 +54,13 @@ class linearInterpolationWeights
//- Cached index in samples from previous invocation
mutable label index_;
// Private Member Functions
//- Get weights of i and i+1 to calculate integration from t to
// samples_[i+1]
Pair<scalar> integrationWeights
(
const label i,
const scalar t
) const;
public:
//- Runtime type information
TypeName("linear");
// Constructors
//- Construct from components
@ -98,8 +90,8 @@ public:
// from samples. Returns true if indices changed.
virtual bool integrationWeights
(
const scalar t1,
const scalar t2,
scalar t1,
scalar t2,
labelList& indices,
scalarField& weights
) const;

2
src/OpenFOAM/interpolations/interpolationWeights/splineInterpolationWeights/splineInterpolationWeights.H
View File

@ -55,11 +55,13 @@ class splineInterpolationWeights
//- Cached index in samples from previous invocation
mutable label index_;
public:
//- Runtime type information
TypeName("spline");
// Constructors
//- Construct from components. By default make sure samples are

5
src/OpenFOAM/primitives/functions/Function1/Function1/Function1.C
View File

@ -59,11 +59,6 @@ const Foam::word& Foam::Function1<Type>::name() const
}
template<class Type>
void Foam::Function1<Type>::convertTimeBase(const Time&)
{}
template<class Type>
Type Foam::Function1<Type>::value(const scalar x) const
{

6
src/OpenFOAM/primitives/functions/Function1/Function1/Function1.H
View File

@ -125,12 +125,6 @@ public:
const word& name() const;
// Manipulation
//- Convert time
virtual void convertTimeBase(const Time& t);
// Evaluation
//- Return value as a function of (scalar) independent variable

17
src/OpenFOAM/primitives/functions/Function1/PolynomialEntry/PolynomialEntry.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2011-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -128,21 +128,6 @@ Foam::Function1Types::Polynomial<Type>::~Polynomial()
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
void Foam::Function1Types::Polynomial<Type>::convertTimeBase(const Time& t)
{
forAll(coeffs_, i)
{
Type value = coeffs_[i].first();
for (direction cmpt = 0; cmpt < pTraits<Type>::nComponents; cmpt++)
{
setComponent(coeffs_[i].first(), cmpt) =
t.userTimeToTime(component(value, cmpt));
}
}
}
template<class Type>
Type Foam::Function1Types::Polynomial<Type>::value(const scalar x) const
{

6
src/OpenFOAM/primitives/functions/Function1/PolynomialEntry/PolynomialEntry.H
View File

@ -101,12 +101,6 @@ public:
// Member Functions
// Manipulation
//- Convert time
virtual void convertTimeBase(const Time& t);
// Evaluation
//- Return Polynomial value

25
src/OpenFOAM/primitives/functions/Function1/Scale/Scale.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -25,20 +25,26 @@ License
#include "Scale.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
void Foam::Function1Types::Scale<Type>::read(const dictionary& coeffs)
{
xScale_ = coeffs.found("xScale")
? Function1<scalar>::New("xScale", coeffs)
: autoPtr<Function1<scalar>>(new Constant<scalar>("xScale", 1));
scale_ = Function1<scalar>::New("scale", coeffs);
xScale_ =
coeffs.found("xScale")
? Function1<scalar>::New("xScale", coeffs)
: autoPtr<Function1<scalar>>(new Constant<scalar>("xScale", 1));
value_ = Function1<Type>::New("value", coeffs);
integrable_ =
isA<Constant<scalar>>(xScale_())
&& isA<Constant<scalar>>(scale_());
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
Foam::Function1Types::Scale<Type>::Scale
(
@ -56,9 +62,10 @@ template<class Type>
Foam::Function1Types::Scale<Type>::Scale(const Scale<Type>& se)
:
Function1<Type>(se),
xScale_(se.xScale_, false),
scale_(se.scale_, false),
value_(se.value_, false)
xScale_(se.xScale_, false),
value_(se.value_, false),
integrable_(se.integrable_)
{}
@ -78,8 +85,8 @@ void Foam::Function1Types::Scale<Type>::writeData(Ostream& os) const
os << token::END_STATEMENT << nl;
os << indent << word(this->name() + "Coeffs") << nl;
os << indent << token::BEGIN_BLOCK << incrIndent << nl;
xScale_->writeData(os);
scale_->writeData(os);
xScale_->writeData(os);
value_->writeData(os);
os << decrIndent << indent << token::END_BLOCK << endl;
}

18
src/OpenFOAM/primitives/functions/Function1/Scale/Scale.H
View File

@ -139,15 +139,18 @@ class Scale
{
// Private Data
//- Argument scaling function
autoPtr<Function1<scalar>> xScale_;
//- Scalar scaling function
autoPtr<Function1<scalar>> scale_;
//- Argument scaling function
autoPtr<Function1<scalar>> xScale_;
//- Value function
autoPtr<Function1<Type>> value_;
//- Is this function integrable?
bool integrable_;
// Private Member Functions
@ -180,12 +183,15 @@ public:
// Member Functions
//- Return value for time t
virtual inline Type value(const scalar t) const;
//- Write in dictionary format
virtual void writeData(Ostream& os) const;
//- Return value
virtual inline Type value(const scalar x) const;
//- Integrate between two values
virtual inline Type integrate(const scalar x1, const scalar x2) const;
// Member Operators

28
src/OpenFOAM/primitives/functions/Function1/Scale/ScaleI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -28,10 +28,30 @@ License
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
inline Type Foam::Function1Types::Scale<Type>::value(const scalar t) const
inline Type Foam::Function1Types::Scale<Type>::value(const scalar x) const
{
const scalar st = xScale_->value(t)*t;
return scale_->value(st)*value_->value(st);
const scalar sx = xScale_->value(x)*x;
return scale_->value(sx)*value_->value(sx);
}
template<class Type>
inline Type Foam::Function1Types::Scale<Type>::integrate
(
const scalar x1,
const scalar x2
) const
{
if (!integrable_)
{
FatalErrorInFunction
<< "Integration is not defined for " << type() << " functions "
<< "unless all scaling is constant"
<< exit(FatalError);
}
const scalar sx = xScale_->value(NaN);
return scale_->value(NaN)/sx*value_->integrate(sx*x1, sx*x2);
}

28
src/OpenFOAM/primitives/functions/Function1/Sine/Sine.C
View File

@ -25,19 +25,24 @@ License
#include "Sine.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
void Foam::Function1Types::Sine<Type>::read(const dictionary& coeffs)
{
t0_ = coeffs.lookupOrDefault<scalar>("t0", 0);
amplitude_ = Function1<scalar>::New("amplitude", coeffs);
frequency_ = Function1<scalar>::New("frequency", coeffs);
scale_ = Function1<Type>::New("scale", coeffs);
amplitude_ = Function1<Type>::New("amplitude", coeffs);
frequency_ = readScalar(coeffs.lookup("frequency"));
start_ = coeffs.lookupOrDefault<scalar>("start", 0);
level_ = Function1<Type>::New("level", coeffs);
integrable_ =
isA<Constant<Type>>(amplitude_())
&& isA<Constant<Type>>(level_());
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
Foam::Function1Types::Sine<Type>::Sine
(
@ -55,11 +60,11 @@ template<class Type>
Foam::Function1Types::Sine<Type>::Sine(const Sine<Type>& se)
:
Function1<Type>(se),
t0_(se.t0_),
amplitude_(se.amplitude_, false),
frequency_(se.frequency_, false),
scale_(se.scale_, false),
level_(se.level_, false)
frequency_(se.frequency_),
start_(se.start_),
level_(se.level_, false),
integrable_(se.integrable_)
{}
@ -79,10 +84,9 @@ void Foam::Function1Types::Sine<Type>::writeData(Ostream& os) const
os << token::END_STATEMENT << nl;
os << indent << word(this->name() + "Coeffs") << nl;
os << indent << token::BEGIN_BLOCK << incrIndent << nl;
writeEntry(os, "t0", t0_);
amplitude_->writeData(os);
frequency_->writeData(os);
scale_->writeData(os);
writeEntry(os, "frequency", frequency_);
writeEntry(os, "start", start_);
level_->writeData(os);
os << decrIndent << indent << token::END_BLOCK << endl;
}

49
src/OpenFOAM/primitives/functions/Function1/Sine/Sine.H
View File

@ -28,19 +28,17 @@ Description
Templated sine function with support for an offset level.
\f[
a sin(2 \pi f (t - t_0)) s + l
a sin(2 \pi f (x - x_0)) + l
\f]
where
\vartable
symbol | Description | Data type
a | Amplitude | Function1<scalar>
f | Frequency [1/s] | Function1<scalar>
s | Type scale factor | Function1<Type>
l | Type offset level | Function1<Type>
t_0 | Start time [s] | scalar
t | Time [s] | scalar
Symbol | Description | Data type | Default
a | Amplitude | Function1<Type> |
f | Frequency | scalar |
x_0 | Start | scalar | 0
l | Offset level | Function1<Type> |
\endvartable
Example for a scalar:
@ -48,9 +46,9 @@ Description
<entryName> sine;
<entryName>Coeffs
{
frequency 10;
amplitude 0.1;
scale 2e-6;
frequency 10;
start 0;
level 2e-6;
}
\endverbatim
@ -60,9 +58,9 @@ Description
<entryName> sine;
<entryName>Coeffs
{
frequency 10;
amplitude 1;
scale (1 0.1 0);
frequency 10;
start 0;
level (10 1 0);
}
\endverbatim
@ -95,21 +93,21 @@ class Sine
{
// Private Data
//- Start-time for the sin function
scalar t0_;
//- Amplitude of the sine function
autoPtr<Function1<Type>> amplitude_;
//- Scalar amplitude of the sin function
autoPtr<Function1<scalar>> amplitude_;
//- Frequency of the sine function
scalar frequency_;
//- Frequency of the sin function
autoPtr<Function1<scalar>> frequency_;
//- Argument offset
scalar start_;
//- Scaling factor of the sin function
autoPtr<Function1<Type>> scale_;
//- Level to which the sin function is added
//- Level to which the sine function is added
autoPtr<Function1<Type>> level_;
//- Is this function integrable?
bool integrable_;
// Private Member Functions
@ -142,8 +140,11 @@ public:
// Member Functions
//- Return value for time t
virtual inline Type value(const scalar t) const;
//- Return value
virtual inline Type value(const scalar x) const;
//- Integrate between two values
virtual inline Type integrate(const scalar x1, const scalar x2) const;
//- Write in dictionary format
virtual void writeData(Ostream& os) const;

36
src/OpenFOAM/primitives/functions/Function1/Sine/SineI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -29,13 +29,37 @@ License
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
inline Type Foam::Function1Types::Sine<Type>::value(const scalar t) const
inline Type Foam::Function1Types::Sine<Type>::value(const scalar x) const
{
const scalar phi = constant::mathematical::twoPi*frequency_*(x - start_);
return amplitude_->value(x)*sin(phi) + level_->value(x);
}
template<class Type>
inline Type Foam::Function1Types::Sine<Type>::integrate
(
const scalar x1,
const scalar x2
) const
{
using namespace constant::mathematical;
if (!integrable_)
{
FatalErrorInFunction
<< "Integration is not defined for " << type() << " functions "
<< "unless the amplitude is constant"
<< exit(FatalError);
}
const scalar phi1 = constant::mathematical::twoPi*frequency_*(x1 - start_);
const scalar phi2 = constant::mathematical::twoPi*frequency_*(x2 - start_);
return
amplitude_->value(t)
*sin(constant::mathematical::twoPi*frequency_->value(t)*(t - t0_))
*scale_->value(t)
+ level_->value(t);
- amplitude_->value(NaN)*(cos(phi2) - cos(phi1))/(twoPi*frequency_)
+ level_->integrate(x1, x2);
}

34
src/OpenFOAM/primitives/functions/Function1/Square/Square.C
View File

@ -25,20 +25,25 @@ License
#include "Square.H"
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
template<class Type>
void Foam::Function1Types::Square<Type>::read(const dictionary& coeffs)
{
t0_ = coeffs.lookupOrDefault<scalar>("t0", 0);
markSpace_ = coeffs.lookupOrDefault<scalar>("markSpace", 1);
amplitude_ = Function1<scalar>::New("amplitude", coeffs);
frequency_ = Function1<scalar>::New("frequency", coeffs);
scale_ = Function1<Type>::New("scale", coeffs);
amplitude_ = Function1<Type>::New("amplitude", coeffs);
frequency_ = readScalar(coeffs.lookup("frequency"));
start_ = coeffs.lookupOrDefault<scalar>("start", 0);
level_ = Function1<Type>::New("level", coeffs);
markSpace_ = coeffs.lookupOrDefault<scalar>("markSpace", 1);
integrable_ =
isA<Constant<Type>>(amplitude_())
&& isA<Constant<Type>>(level_());
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
Foam::Function1Types::Square<Type>::Square
(
@ -56,12 +61,12 @@ template<class Type>
Foam::Function1Types::Square<Type>::Square(const Square<Type>& se)
:
Function1<Type>(se),
t0_(se.t0_),
markSpace_(se.markSpace_),
amplitude_(se.amplitude_, false),
frequency_(se.frequency_, false),
scale_(se.scale_, false),
level_(se.level_, false)
frequency_(se.frequency_),
start_(se.start_),
level_(se.level_, false),
markSpace_(se.markSpace_),
integrable_(se.integrable_)
{}
@ -81,12 +86,11 @@ void Foam::Function1Types::Square<Type>::writeData(Ostream& os) const
os << token::END_STATEMENT << nl;
os << indent << word(this->name() + "Coeffs") << nl;
os << indent << token::BEGIN_BLOCK << incrIndent << nl;
writeEntry(os, "t0", t0_);
writeEntry(os, "markSpace", markSpace_);
amplitude_->writeData(os);
frequency_->writeData(os);
scale_->writeData(os);
writeEntry(os, "frequency", frequency_);
writeEntry(os, "start", start_);
level_->writeData(os);
writeEntry(os, "markSpace", markSpace_);
os << decrIndent << indent << token::END_BLOCK << endl;
}

59
src/OpenFOAM/primitives/functions/Function1/Square/Square.H
View File

@ -28,23 +28,21 @@ Description
Templated square-wave function with support for an offset level.
\f[
a square(f (t - t_0)) s + l
a square(f (x - x_0) ) + l
\f]
where
\f$ square(t) \f$ is the square-wave function in range \f$ [-1, 1] \f$
with a mark/space ratio of \f$ r \f$
\f$square(x)\f$ is the square-wave function in range \f$[-1, 1]\f$ with a
mark/space ratio of \f$r\f$
\vartable
symbol | Description | Data type | Default
a | Amplitude | Function1<scalar> |
f | Frequency [1/s] | Function1<scalar> |
s | Type scale factor | Function1<Type> |
l | Type offset level | Function1<Type> |
t_0 | Start time [s] | scalar | 0
r | mark/space ratio | scalar | 1
t | Time [s] | scalar
Symbol | Description | Data type | Default
a | Amplitude | Function1<Type> |
f | Frequency | scalar |
x_0 | Phase | scalar | 0
l | Offset level | Function1<Type> |
r | Mark/space ratio | scalar | 1
\endvartable
Example for a scalar:
@ -52,10 +50,11 @@ Description
<entryName> square;
<entryName>Coeffs
{
frequency 10;
amplitude 0.1;
scale 2e-6;
frequency 10;
start 0;
level 2e-6;
markSpace 0.5;
}
\endverbatim
@ -64,10 +63,11 @@ Description
<entryName> square;
<entryName>Coeffs
{
frequency 10;
amplitude 1;
scale (1 0.1 0);
frequency 10;
start 0;
level (10 1 0);
markSpace 0.5;
}
\endverbatim
@ -99,24 +99,26 @@ class Square
{
// Private Data
//- Start-time for the square function
scalar t0_;
//- Mark/space ratio of the square function
scalar markSpace_;
//- Scalar amplitude of the square function
autoPtr<Function1<scalar>> amplitude_;
autoPtr<Function1<Type>> amplitude_;
//- Frequency of the square function
autoPtr<Function1<scalar>> frequency_;
scalar frequency_;
//- Scaling factor of the square function
autoPtr<Function1<Type>> scale_;
//- Argument offset
scalar start_;
//- Level to which the square function is added
autoPtr<Function1<Type>> level_;
//- Mark/space ratio of the square function; the fraction of one period
// in which the function value is at the higher of the two possible
// values.
scalar markSpace_;
//- Is this function integrable?
bool integrable_;
// Private Member Functions
@ -149,8 +151,11 @@ public:
// Member Functions
//- Return value for time t
virtual inline Type value(const scalar t) const;
//- Return value
virtual inline Type value(const scalar x) const;
//- Integrate between two values
virtual inline Type integrate(const scalar x1, const scalar x2) const;
//- Write in dictionary format
virtual void writeData(Ostream& os) const;

59
src/OpenFOAM/primitives/functions/Function1/Square/SquareI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -24,29 +24,56 @@ License
\*---------------------------------------------------------------------------*/
#include "Square.H"
#include "mathematicalConstants.H"
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
inline Type Foam::Function1Types::Square<Type>::value(const scalar t) const
inline Type Foam::Function1Types::Square<Type>::value(const scalar x) const
{
// Number of waves including fractions
scalar waves = frequency_->value(t)*(t - t0_);
const scalar markFraction = markSpace_/(1 + markSpace_);
// Number of complete waves
scalar nWaves;
// Fraction of last incomplete wave
scalar waveFrac = std::modf(waves, &nWaves);
// Mark fraction of a wave
scalar markFrac = markSpace_/(1.0 + markSpace_);
const scalar phi = frequency_*(x - start_);
const scalar fraction = phi - floor(phi);
return
amplitude_->value(t)
*(waveFrac < markFrac ? 1 : -1)
*scale_->value(t)
+ level_->value(t);
amplitude_->value(x)*(fraction < markFraction ? 1 : -1)
+ level_->value(x);
}
template<class Type>
inline Type Foam::Function1Types::Square<Type>::integrate
(
const scalar x1,
const scalar x2
) const
{
if (!integrable_)
{
FatalErrorInFunction
<< "Integration is not defined for " << type() << " functions "
<< "unless the amplitude is constant"
<< exit(FatalError);
}
auto integrate0 = [&](const scalar x)
{
const scalar markFraction = markSpace_/(1 + markSpace_);
const scalar phi = frequency_*(x - start_);
const scalar fraction = phi - floor(phi);
return
(2*amplitude_->value(x)/frequency_)
*(
markFraction*floor(phi + 1 - markFraction)
+ fraction*(fraction < markFraction ? 1 : 0)
- frequency_/2*x
);
};
return integrate0(x2) - integrate0(x1) + level_->integrate(x1, x2);
}

351
src/OpenFOAM/primitives/functions/Function1/Table/TableBase.C
View File

@ -53,138 +53,6 @@ Foam::Function1Types::TableBase<Type>::interpolator() const
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
Foam::Function1Types::TableBase<Type>::TableBase
(
const word& name,
const dictionary& dict
)
:
Function1<Type>(name),
name_(name),
boundsHandling_
(
wordToBoundsHandling
(
dict.lookupOrDefault<word>("outOfBounds", "clamp")
)
),
interpolationScheme_
(
dict.lookupOrDefault<word>("interpolationScheme", "linear")
),
table_()
{}
template<class Type>
Foam::Function1Types::TableBase<Type>::TableBase(const TableBase<Type>& tbl)
:
Function1<Type>(tbl),
name_(tbl.name_),
boundsHandling_(tbl.boundsHandling_),
interpolationScheme_(tbl.interpolationScheme_),
table_(tbl.table_),
tableSamplesPtr_(tbl.tableSamplesPtr_),
interpolatorPtr_(tbl.interpolatorPtr_)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
template<class Type>
Foam::Function1Types::TableBase<Type>::~TableBase()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
Foam::word Foam::Function1Types::TableBase<Type>::boundsHandlingToWord
(
const boundsHandling& bound
) const
{
word enumName("warn");
switch (bound)
{
case ERROR:
{
enumName = "error";
break;
}
case WARN:
{
enumName = "warn";
break;
}
case CLAMP:
{
enumName = "clamp";
break;
}
case REPEAT:
{
enumName = "repeat";
break;
}
}
return enumName;
}
template<class Type>
typename Foam::Function1Types::TableBase<Type>::boundsHandling
Foam::Function1Types::TableBase<Type>::wordToBoundsHandling
(
const word& bound
) const
{
if (bound == "error")
{
return ERROR;
}
else if (bound == "warn")
{
return WARN;
}
else if (bound == "clamp")
{
return CLAMP;
}
else if (bound == "repeat")
{
return REPEAT;
}
else
{
WarningInFunction
<< "bad outOfBounds specifier " << bound << " using 'warn'"
<< endl;
return WARN;
}
}
template<class Type>
typename Foam::Function1Types::TableBase<Type>::boundsHandling
Foam::Function1Types::TableBase<Type>::outOfBounds
(
const boundsHandling& bound
)
{
boundsHandling prev = boundsHandling_;
boundsHandling_ = bound;
return prev;
}
template<class Type>
void Foam::Function1Types::TableBase<Type>::check() const
{
@ -215,146 +83,114 @@ void Foam::Function1Types::TableBase<Type>::check() const
template<class Type>
bool Foam::Function1Types::TableBase<Type>::checkMinBounds
(
const scalar x,
scalar& xDash
) const
Foam::scalar Foam::Function1Types::TableBase<Type>::bound(const scalar x) const
{
if (x < table_[0].first())
const bool under = x < table_.first().first();
const bool over = x > table_.last().first();
auto errorMessage = [&]()
{
return "value (" + name(x) + ") " + (under ? "under" : "over") + "flow";
};
if (under || over)
{
switch (boundsHandling_)
{
case ERROR:
case tableBase::boundsHandling::error:
{
FatalErrorInFunction
<< "value (" << x << ") underflow"
<< exit(FatalError);
<< errorMessage() << nl << exit(FatalError);
break;
}
case WARN:
case tableBase::boundsHandling::warn:
{
WarningInFunction
<< "value (" << x << ") underflow" << nl
<< endl;
// fall-through to 'CLAMP'
[[fallthrough]];
}
case CLAMP:
{
xDash = table_[0].first();
return true;
<< errorMessage() << nl << endl;
break;
}
case REPEAT:
case tableBase::boundsHandling::clamp:
{
// adjust x to >= minX
scalar span = table_.last().first() - table_[0].first();
xDash = fmod(x - table_[0].first(), span) + table_[0].first();
break;
}
case tableBase::boundsHandling::repeat:
{
const scalar t0 = table_.first().first();
const scalar t1 = table_.last().first();
const scalar dt = t1 - t0;
const label n = floor((x - t0)/dt);
return x - n*dt;
}
}
}
else
{
xDash = x;
}
return false;
return x;
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class Type>
bool Foam::Function1Types::TableBase<Type>::checkMaxBounds
Foam::Function1Types::TableBase<Type>::TableBase
(
const scalar x,
scalar& xDash
) const
{
if (x > table_.last().first())
{
switch (boundsHandling_)
{
case ERROR:
{
FatalErrorInFunction
<< "value (" << x << ") overflow"
<< exit(FatalError);
break;
}
case WARN:
{
WarningInFunction
<< "value (" << x << ") overflow" << nl
<< endl;
// fall-through to 'CLAMP'
[[fallthrough]];
}
case CLAMP:
{
xDash = table_.last().first();
return true;
break;
}
case REPEAT:
{
// adjust x to >= minX
scalar span = table_.last().first() - table_[0].first();
xDash = fmod(x - table_[0].first(), span) + table_[0].first();
break;
}
}
}
else
{
xDash = x;
}
return false;
}
const word& name,
const dictionary& dict
)
:
tableBase(),
Function1<Type>(name),
name_(name),
boundsHandling_
(
dict.found("outOfBounds")
? tableBase::boundsHandlingNames_.read(dict.lookup("outOfBounds"))
: tableBase::boundsHandling::warn
),
interpolationScheme_
(
dict.lookupOrDefault<word>("interpolationScheme", "linear")
),
table_()
{}
template<class Type>
void Foam::Function1Types::TableBase<Type>::convertTimeBase(const Time& t)
{
forAll(table_, i)
{
scalar value = table_[i].first();
table_[i].first() = t.userTimeToTime(value);
}
Foam::Function1Types::TableBase<Type>::TableBase(const TableBase<Type>& tbl)
:
tableBase(),
Function1<Type>(tbl),
name_(tbl.name_),
boundsHandling_(tbl.boundsHandling_),
interpolationScheme_(tbl.interpolationScheme_),
table_(tbl.table_),
tableSamplesPtr_(tbl.tableSamplesPtr_),
interpolatorPtr_(tbl.interpolatorPtr_)
{}
tableSamplesPtr_.clear();
interpolatorPtr_.clear();
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
template<class Type>
Foam::Function1Types::TableBase<Type>::~TableBase()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class Type>
Type Foam::Function1Types::TableBase<Type>::value(const scalar x) const
{
scalar xDash = x;
const scalar bx = bound(x);
if (checkMinBounds(x, xDash))
Type y = Zero;
interpolator().valueWeights(bx, indices_, weights_);
forAll(indices_, i)
{
return table_[0].second();
y += weights_[i]*table_[indices_[i]].second();
}
if (checkMaxBounds(xDash, xDash))
{
return table_.last().second();
}
// Use interpolator
interpolator().valueWeights(xDash, currentIndices_, currentWeights_);
Type t = currentWeights_[0]*table_[currentIndices_[0]].second();
for (label i = 1; i < currentIndices_.size(); i++)
{
t += currentWeights_[i]*table_[currentIndices_[i]].second();
}
return t;
return y;
}
@ -365,16 +201,38 @@ Type Foam::Function1Types::TableBase<Type>::integrate
const scalar x2
) const
{
// Use interpolator
interpolator().integrationWeights(x1, x2, currentIndices_, currentWeights_);
const scalar bx1 = bound(x1), bx2 = bound(x2);
Type sum = currentWeights_[0]*table_[currentIndices_[0]].second();
for (label i = 1; i < currentIndices_.size(); i++)
Type sumY = Zero;
interpolator().integrationWeights(bx1, bx2, indices_, weights_);
forAll(indices_, i)
{
sum += currentWeights_[i]*table_[currentIndices_[i]].second();
sumY += weights_[i]*table_[indices_[i]].second();
}
return sum;
if (boundsHandling_ == tableBase::boundsHandling::repeat)
{
const scalar t0 = table_.first().first();
const scalar t1 = table_.last().first();
const scalar dt = t1 - t0;
const label n = floor((x2 - t0)/dt) - floor((x1 - t0)/dt);
if (n != 0)
{
Type sumY01 = Zero;
interpolator().integrationWeights(t0, t1, indices_, weights_);
forAll(indices_, i)
{
sumY01 += weights_[i]*table_[indices_[i]].second();
}
sumY += n*sumY01;
}
}
return sumY;
}
@ -411,9 +269,14 @@ Foam::tmp<Foam::Field<Type>> Foam::Function1Types::TableBase<Type>::y() const
template<class Type>
void Foam::Function1Types::TableBase<Type>::writeEntries(Ostream& os) const
{
if (boundsHandling_ != CLAMP)
if (boundsHandling_ != tableBase::boundsHandling::clamp)
{
writeEntry(os, "outOfBounds", boundsHandlingToWord(boundsHandling_));
writeEntry
(
os,
"outOfBounds",
tableBase::boundsHandlingNames_[boundsHandling_]
);
}
if (interpolationScheme_ != "linear")
{

54
src/OpenFOAM/primitives/functions/Function1/Table/TableBase.H
View File

@ -35,6 +35,7 @@ SourceFiles
#ifndef TableBase_H
#define TableBase_H
#include "tableBase.H"
#include "Function1.H"
#include "Tuple2.H"
@ -55,22 +56,9 @@ namespace Function1Types
template<class Type>
class TableBase
:
public tableBase,
public Function1<Type>
{
public:
// Public data types
//- Enumeration for handling out-of-bound values
enum boundsHandling
{
ERROR, //!< Exit with a FatalError
WARN, //!< Issue warning and clamp value (default)
CLAMP, //!< Clamp value to the start/end value
REPEAT //!< Treat as a repeating list
};
protected:
// Protected data
@ -79,7 +67,7 @@ protected:
const word name_;
//- Enumeration for handling out-of-bound values
const boundsHandling boundsHandling_;
const tableBase::boundsHandling boundsHandling_;
//- Interpolation type
const word interpolationScheme_;
@ -93,10 +81,11 @@ protected:
//- Interpolator method
mutable autoPtr<interpolationWeights> interpolatorPtr_;
//- Cached indices and weights
mutable labelList currentIndices_;
//- Cached indices
mutable labelList indices_;
mutable scalarField currentWeights_;
//- Cached weights
mutable scalarField weights_;
// Protected Member Functions
@ -104,6 +93,14 @@ protected:
//- Return (demand driven) interpolator
const interpolationWeights& interpolator() const;
//- Check the table for size and consistency
void check() const;
//- Bound the argument to the table. Errors or warns, or shifts the
// value if the table repeats. Does not clamp to the ends of the table
// as the interpolator already performs that function.
scalar bound(const scalar x) const;
private:
@ -128,27 +125,6 @@ public:
// Member Functions
//- Return the out-of-bounds handling as a word
word boundsHandlingToWord(const boundsHandling& bound) const;
//- Return the out-of-bounds handling as an enumeration
boundsHandling wordToBoundsHandling(const word& bound) const;
//- Set the out-of-bounds handling from enum, return previous setting
boundsHandling outOfBounds(const boundsHandling& bound);
//- Check the table for size and consistency
void check() const;
//- Check minimum table bounds
bool checkMinBounds(const scalar x, scalar& xDash) const;
//- Check maximum table bounds
bool checkMaxBounds(const scalar x, scalar& xDash) const;
//- Convert time
virtual void convertTimeBase(const Time& t);
//- Return Table value
virtual Type value(const scalar x) const;

54
src/OpenFOAM/primitives/functions/Function1/Table/tableBase.C Normal file
View File

@ -0,0 +1,54 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
\*---------------------------------------------------------------------------*/
#include "tableBase.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
namespace Foam
{
template<>
const char*
NamedEnum<Function1Types::tableBase::boundsHandling, 4>::names[] =
{"error", "warn", "clamp", "repeat"};
}
const Foam::NamedEnum<Foam::Function1Types::tableBase::boundsHandling, 4>
Foam::Function1Types::tableBase::boundsHandlingNames_;
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::Function1Types::tableBase::tableBase()
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::Function1Types::tableBase::~tableBase()
{}
// ************************************************************************* //

87
src/OpenFOAM/primitives/functions/Function1/Table/tableBase.H Normal file
View File

@ -0,0 +1,87 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::tableBase
SourceFiles
tableBase.C
\*---------------------------------------------------------------------------*/
#ifndef tableBase_H
#define tableBase_H
#include "NamedEnum.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace Function1Types
{
/*---------------------------------------------------------------------------*\
Class tableBase Declaration
\*---------------------------------------------------------------------------*/
class tableBase
{
public:
// Public data types
//- Enumeration for handling out-of-bound values
enum class boundsHandling
{
error, // Exit with a FatalError
warn, // Issue warning and clamp value (default)
clamp, // Clamp value to the start/end value
repeat // Treat as a repeating list
};
//- Enumeration names for handling out-of-bound values
static const NamedEnum<boundsHandling, 4> boundsHandlingNames_;
// Constructors
//- Construct null
tableBase();
//- Destructor
virtual ~tableBase();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // namespace Function1Types
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

7
src/OpenFOAM/primitives/functions/Function1/halfCosineRamp/halfCosineRamp.H
View File

@ -82,6 +82,13 @@ public:
//- Return value for time t
virtual inline scalar value(const scalar t) const;
//- Return the integral between times t1 and t2
virtual inline scalar integrate
(
const scalar t1,
const scalar t2
) const;
// Member Operators

18
src/OpenFOAM/primitives/functions/Function1/halfCosineRamp/halfCosineRampI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -37,4 +37,20 @@ inline Foam::scalar Foam::Function1Types::halfCosineRamp::value
}
inline Foam::scalar Foam::Function1Types::halfCosineRamp::integrate
(
const scalar t1,
const scalar t2
) const
{
using namespace constant::mathematical;
const scalar l1 = linearRamp(t1), l2 = linearRamp(t2);
return
0.5*(l2 - l1 - (sin(pi*l2) - sin(pi*l1))/pi)/dLinearRampDt()
+ (t2 - min(t2, end()));
}
// ************************************************************************* //

7
src/OpenFOAM/primitives/functions/Function1/linearRamp/linearRamp.H
View File

@ -82,6 +82,13 @@ public:
//- Return value for time t
virtual inline scalar value(const scalar t) const;
//- Return the integral between times t1 and t2
virtual inline scalar integrate
(
const scalar t1,
const scalar t2
) const;
// Member Operators

16
src/OpenFOAM/primitives/functions/Function1/linearRamp/linearRampI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -36,4 +36,18 @@ inline Foam::scalar Foam::Function1Types::linearRamp::value
}
inline Foam::scalar Foam::Function1Types::linearRamp::integrate
(
const scalar t1,
const scalar t2
) const
{
const scalar l1 = ramp::linearRamp(t1), l2 = ramp::linearRamp(t2);
return
(sqr(l2) - sqr(l1))/2/dLinearRampDt()
+ (t2 - min(t2, end()));
}
// ************************************************************************* //

7
src/OpenFOAM/primitives/functions/Function1/quadraticRamp/quadraticRamp.H
View File

@ -82,6 +82,13 @@ public:
//- Return value for time t
virtual inline scalar value(const scalar t) const;
//- Return the integral between times t1 and t2
virtual inline scalar integrate
(
const scalar t1,
const scalar t2
) const;
// Member Operators

14
src/OpenFOAM/primitives/functions/Function1/quadraticRamp/quadraticRampI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -36,4 +36,16 @@ inline Foam::scalar Foam::Function1Types::quadraticRamp::value
}
inline Foam::scalar Foam::Function1Types::quadraticRamp::integrate
(
const scalar t1,
const scalar t2
) const
{
return
(pow3(linearRamp(t2)) - pow3(linearRamp(t1)))/3/dLinearRampDt()
+ (t2 - min(t2, end()));
}
// ************************************************************************* //

7
src/OpenFOAM/primitives/functions/Function1/quarterCosineRamp/quarterCosineRamp.H
View File

@ -82,6 +82,13 @@ public:
//- Return value for time t
virtual inline scalar value(const scalar t) const;
//- Return the integral between times t1 and t2
virtual inline scalar integrate
(
const scalar t1,
const scalar t2
) const;
// Member Operators

20
src/OpenFOAM/primitives/functions/Function1/quarterCosineRamp/quarterCosineRampI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -33,7 +33,23 @@ inline Foam::scalar Foam::Function1Types::quarterCosineRamp::value
const scalar t
) const
{
return 1 - cos(0.5*constant::mathematical::pi*linearRamp(t));
return 1 - cos(constant::mathematical::piByTwo*linearRamp(t));
}
inline Foam::scalar Foam::Function1Types::quarterCosineRamp::integrate
(
const scalar t1,
const scalar t2
) const
{
using namespace constant::mathematical;
const scalar l1 = linearRamp(t1), l2 = linearRamp(t2);
return
(l2 - l1 - (sin(piByTwo*l2) - sin(piByTwo*l1))/piByTwo)/dLinearRampDt()
+ (t2 - min(t2, end()));
}

7
src/OpenFOAM/primitives/functions/Function1/quarterSineRamp/quarterSineRamp.H
View File

@ -82,6 +82,13 @@ public:
//- Return value for time t
virtual inline scalar value(const scalar t) const;
//- Return the integral between times t1 and t2
virtual inline scalar integrate
(
const scalar t1,
const scalar t2
) const;
// Member Operators

20
src/OpenFOAM/primitives/functions/Function1/quarterSineRamp/quarterSineRampI.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -33,7 +33,23 @@ inline Foam::scalar Foam::Function1Types::quarterSineRamp::value
const scalar t
) const
{
return sin(0.5*constant::mathematical::pi*linearRamp(t));
return sin(constant::mathematical::piByTwo*linearRamp(t));
}
inline Foam::scalar Foam::Function1Types::quarterSineRamp::integrate
(
const scalar t1,
const scalar t2
) const
{
using namespace constant::mathematical;
const scalar l1 = linearRamp(t1), l2 = linearRamp(t2);
return
- (cos(piByTwo*l2) - cos(piByTwo*l1))/piByTwo/dLinearRampDt()
+ (t2 - min(t2, end()));
}

17
src/OpenFOAM/primitives/functions/Function1/ramp/ramp.H
View File

@ -93,13 +93,26 @@ protected:
//- Duration of the ramp function
scalar duration_;
//- Simple linear ramp function
// which form the basis of many more complex ramp functions
//- Return the end-time of the ramp function
inline scalar end() const
{
return start_ + duration_;
}
//- Simple linear ramp function which form the basis of many more
// complex ramp functions
inline scalar linearRamp(const scalar t) const
{
return max(min((t - start_)/duration_, 1), 0);
}
//- The gradient of the linear ramp within the bounds of the ramp
// duration
inline scalar dLinearRampDt() const
{
return 1/duration_;
}
private:

7
src/OpenFOAM/primitives/functions/Function1/reverseRamp/reverseRamp.H
View File

@ -115,6 +115,13 @@ public:
//- Return value for time t
virtual inline scalar value(const scalar t) const;
//- Return the integral between times t1 and t2
virtual inline scalar integrate
(
const scalar t1,
const scalar t2
) const;
//- Write in dictionary format
virtual void writeData(Ostream& os) const;

10
src/OpenFOAM/primitives/functions/Function1/reverseRamp/reverseRampI.H
View File

@ -36,4 +36,14 @@ inline Foam::scalar Foam::Function1Types::reverseRamp::value
}
inline Foam::scalar Foam::Function1Types::reverseRamp::integrate
(
const scalar t1,
const scalar t2
) const
{
return (t2 - t1) - ramp_->integrate(t1, t2);
}
// ************************************************************************* //

85
src/OpenFOAM/primitives/functions/TimeFunction1/TimeFunction1.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2012-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2012-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -30,43 +30,40 @@ License
template<class Type>
Foam::TimeFunction1<Type>::TimeFunction1
(
const Time& t,
const word& name,
const Time& time,
const word& entryName,
const dictionary& dict
)
:
time_(t),
name_(name),
entry_(Function1<Type>::New(name, dict))
{
entry_->convertTimeBase(t);
}
template<class Type>
Foam::TimeFunction1<Type>::TimeFunction1(const Time& t, const word& name)
:
time_(t),
name_(name),
entry_(nullptr)
time_(time),
name_(entryName),
function_(Function1<Type>::New(entryName, dict))
{}
template<class Type>
Foam::TimeFunction1<Type>::TimeFunction1
(
const TimeFunction1<Type>& tde
const Time& time,
const word& entryName
)
:
time_(tde.time_),
name_(tde.name_),
entry_()
{
if (tde.entry_.valid())
{
entry_.reset(tde.entry_->clone().ptr());
}
}
time_(time),
name_(entryName),
function_(nullptr)
{}
template<class Type>
Foam::TimeFunction1<Type>::TimeFunction1
(
const TimeFunction1<Type>& tf
)
:
time_(tf.time_),
name_(tf.name_),
function_(tf.function_, false)
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
@ -81,30 +78,14 @@ Foam::TimeFunction1<Type>::~TimeFunction1()
template<class Type>
void Foam::TimeFunction1<Type>::reset(const dictionary& dict)
{
entry_.reset
(
Function1<Type>::New
(
name_,
dict
).ptr()
);
entry_->convertTimeBase(time_);
}
template<class Type>
const Foam::word& Foam::TimeFunction1<Type>::name() const
{
return entry_->name();
function_.reset(Function1<Type>::New(name_, dict).ptr());
}
template<class Type>
Type Foam::TimeFunction1<Type>::value(const scalar x) const
{
return entry_->value(x);
return function_->value(time_.userTimeToTime(x));
}
@ -115,7 +96,13 @@ Type Foam::TimeFunction1<Type>::integrate
const scalar x2
) const
{
return entry_->integrate(x1, x2);
return
time_.timeToUserTimeRatio()
*function_->integrate
(
time_.userTimeToTime(x1),
time_.userTimeToTime(x2)
);
}
@ -125,17 +112,17 @@ template<class Type>
Foam::Ostream& Foam::operator<<
(
Ostream& os,
const TimeFunction1<Type>& de
const TimeFunction1<Type>& tf
)
{
return de.entry_->operator<<(os, de);
return os << tf.function_();
}
template<class Type>
void Foam::TimeFunction1<Type>::writeData(Ostream& os) const
{
entry_->writeData(os);
function_->writeData(os);
}

20
src/OpenFOAM/primitives/functions/TimeFunction1/TimeFunction1.H
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2012-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2012-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -60,7 +60,6 @@ Ostream& operator<<
template<class Type>
class TimeFunction1
{
protected:
// Protected data
@ -71,26 +70,26 @@ protected:
//- Name of the data entry
const word name_;
//- The underlying Function1
autoPtr<Function1<Type>> entry_;
//- The function
autoPtr<Function1<Type>> function_;
public:
// Constructor
//- Construct from entry name
//- Construct from time, entry name and dictionary
TimeFunction1
(
const Time& t,
const word& name,
const Time& time,
const word& entryName,
const dictionary& dict
);
//- Construct null from entry name
//- Construct null from time and entry name
TimeFunction1
(
const Time& t,
const Time& time,
const word& entryName
);
@ -109,9 +108,6 @@ public:
//- Reset entry by re-reading from dictionary
void reset(const dictionary& dict);
//- Return the name of the entry
const word& name() const;
// Evaluation

8
src/engine/engineTime/crankConRod/crankConRod.C
View File

@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration | Website: https://openfoam.org
\\ / A nd | Copyright (C) 2017-2018 OpenFOAM Foundation
\\ / A nd | Copyright (C) 2017-2019 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@ -201,4 +201,10 @@ Foam::scalar Foam::crankConRod::timeToUserTime(const scalar t) const
}
Foam::scalar Foam::crankConRod::timeToUserTimeRatio() const
{
return timeToDeg(1);
}
// ************************************************************************* //

3
src/engine/engineTime/crankConRod/crankConRod.H
View File

@ -183,6 +183,9 @@ public:
//- Convert the real-time (s) into user-time (CA deg)
virtual scalar timeToUserTime(const scalar t) const;
//- Ratio between real-time and user-time
virtual scalar timeToUserTimeRatio() const;
//- Read the control dictionary and set the write controls etc.
virtual void readDict();

7
tutorials/multiphase/potentialFreeSurfaceFoam/oscillatingBox/0/U
View File

@ -49,10 +49,9 @@ boundaryField
value sine;
valueCoeffs
{
t0 0;
amplitude constant 0.025;
frequency constant 1;
scale constant (0 1 0);
amplitude constant (0 0.025 0);
frequency 1;
start 0;
level constant (0 0 0);
}
}