1 files changed, 353 insertions, 0 deletions
diff --git a/src/functors/analytical/analyticalF.h b/src/functors/analytical/analyticalF.h
new file mode 100644
index 0000000..0abd70f
--- /dev/null
+++ b/src/functors/analytical/analyticalF.h
@@ -0,0 +1,353 @@
+/*  This file is part of the OpenLB library
+ *
+ *  Copyright (C) 2012 Lukas Baron, Tim Dornieden, Mathias J. Krause, Albert Mink
+ *  E-mail contact: info@openlb.net
+ *  The most recent release of OpenLB can be downloaded at
+ *  <http://www.openlb.net/>
+ *
+ *  This program 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 2
+ *  of the License, or (at your option) any later version.
+ *
+ *  This program 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 this program; if not, write to the Free
+ *  Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
+ *  Boston, MA  02110-1301, USA.
+*/
+
+#ifndef ANALYTICAL_F_H
+#define ANALYTICAL_F_H
+
+#include <vector>
+
+#include "analyticalBaseF.h"
+#include "indicator/smoothIndicatorF2D.h"
+#include "indicator/smoothIndicatorF3D.h"
+
+
+/**
+ *  The functor dimensions are given by F: S^m -> T^n  (S=source, T=target)
+ *  and are implemented via GenericF(n,m).
+ *  Don't get confused by the flipped order of source and target.
+ */
+
+namespace olb {
+
+template<typename T, typename S, bool> class SmoothIndicatorSphere3D;
+template<typename T,  typename DESCRIPTOR> class RadiativeUnitConverter;
+
+////////////////////////////////////////////////////////////////////////////////
+////////implementation of several 1d,2d,3d functors (analyticalFXD)/////////////
+////////////////////////////////////////////////////////////////////////////////
+
+
+//////////////////////////////////1D////////////////////////////////////////////
+
+/// AnalyticalConst1D: 1D -> XD, where XD is defined by value.size()
+template <typename T, typename S>
+class AnalyticalConst1D : public AnalyticalF1D<T,S> {
+private:
+  // is constant return value of operator()
+  std::vector<T> _c;
+public:
+  AnalyticalConst1D(T value);
+  AnalyticalConst1D(const std::vector<T>& value);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalLinear1D: 1D -> 1D troughout given points (x0,v0) and (x1,v1)
+//  Punktsteigungsform
+template <typename T, typename S>
+class AnalyticalLinear1D : public AnalyticalF1D<T,S> {
+private:
+  T _a;
+  T _b;
+public:
+  AnalyticalLinear1D(T a, T b);
+  AnalyticalLinear1D(S x0, T v0, S x1, T v1);
+  bool operator() (T output[], const S x[]) override; ///< returns line _a*x + _b
+};
+
+/// AnalyticalRandom1D: 1D -> 1D with random image in (0,1)
+template <typename T, typename S>
+class AnalyticalRandom1D : public AnalyticalF1D<T,S> {
+public:
+  AnalyticalRandom1D();
+  bool operator() (T output[], const S x[]) override;
+};
+
+
+/// represents an inverse parabola profile like it is used in Poiseuille inflow
+/// note: output depends only on first parameter, maps 1D,2D,3D->1D
+template <typename T, typename S>
+class AnalyticalSquare1D : public AnalyticalF1D<T,S> {
+private:
+  S _cp;
+  S _r;
+  T _maxi;
+public:
+  AnalyticalSquare1D(S cp, S r, T maxi);
+  bool operator() (T output[], const S x[]) override;
+};
+
+
+/// SinusStartScale: 1D -> 1D a start curve based on sinus for a continuous transition at 0 and 1
+template <typename T, typename S>
+class SinusStartScale : public AnalyticalF1D<T,S> {
+protected:
+  S _numTimeSteps;
+  T _maxValue;
+  T _pi;
+public:
+  SinusStartScale(int numTimeSteps=1, T maxValue=1);
+  bool operator() (T output[], const S x[]) override;
+};
+
+
+/// PolynomialStartScale: 1D -> 1D a start curve based on a polynomial fifth order for a continuous transition at 0 and 1: maxValue*(6*y^5-15*y^4+10*y^3)
+template <typename T, typename S>
+class PolynomialStartScale : public AnalyticalF1D<T,S> {
+protected:
+  S _numTimeSteps;
+  T _maxValue;
+public:
+  PolynomialStartScale(S numTimeSteps=S(1), T maxValue=T(1));
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// Derivative of a given 1D functor computed with a finite difference
+template <typename T>
+class AnalyticalDiffFD1D : public AnalyticalF1D<T,T> {
+protected:
+  AnalyticalF1D<T,T>& _f;
+  T _eps;
+public:
+  AnalyticalDiffFD1D(AnalyticalF1D<T,T>& f, T eps = 1.e-10);
+  bool operator() (T output[], const T input[]) override;
+};
+
+//////////////////////////////////2D////////////////////////////////////////////
+
+template <typename T, typename S>
+class AnalyticalComposed2D final : public AnalyticalF2D<T,S> {
+private:
+  AnalyticalF2D<T,S>& _f0;
+  AnalyticalF2D<T,S>& _f1;
+public:
+  AnalyticalComposed2D(AnalyticalF2D<T,S>& f0, AnalyticalF2D<T,S>& f1);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalConst2D: 2D -> XD, where XD is defined by value.size()
+template <typename T, typename S>
+class AnalyticalConst2D final : public AnalyticalF2D<T,S> {
+private:
+  // is constant return value of operator()
+  std::vector<T> _c;
+public:
+  AnalyticalConst2D(T value);
+  AnalyticalConst2D(T value0, T value1);
+  AnalyticalConst2D(T value0, T value1, T value2);
+  AnalyticalConst2D(const std::vector<T>& value);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalLinear2D: 2D -> 1D troughout given points (x0,y0,v0), (x1,y1,v1), (x2,y2,v2)
+template <typename T, typename S>
+class AnalyticalLinear2D final : public AnalyticalF2D<T,S> {
+protected:
+  T _a;
+  T _b;
+  T _c;
+public:
+  AnalyticalLinear2D(T a, T b, T c);
+  AnalyticalLinear2D(S x0, S y0, T v0, S x1, S y1, T v1, S x2, S y2, T v2);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalRandom2D: 2D -> 1D with random image in (0,1)
+template <typename T, typename S>
+class AnalyticalRandom2D final : public AnalyticalF2D<T,S> {
+public:
+  AnalyticalRandom2D();
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalRandom2D: 2D -> 1D with maxValue in the center decreasing linearly with the distrance to the center to zero at the radius and zero outside
+template <typename T, typename S>
+class AnalyticalParticleAdsorptionLinear2D final : public AnalyticalF2D<T,S> {
+protected:
+  T _center[2];
+  T _radius;
+  T _maxValue;
+public:
+  AnalyticalParticleAdsorptionLinear2D(T center[], T radius, T maxValue);
+  bool operator() (T output[], const S x[]);
+};
+
+/** Computes resulting velocity of an object from translational and rotational velocity.
+ * \param indicator Class defining the object (needs to be a SmoothIndicatorF2D<T,T,true>)
+ * \param u translational velocity of the object - expected in lattice units
+ * \param omega rotational velocity of the object - expected in lattice units
+ */
+template <typename T, typename S, typename DESCRIPTOR>
+class ParticleU2D : public AnalyticalF2D<T,S> {
+protected:
+  SmoothIndicatorF2D<T,T,true>& _indicator;
+  UnitConverter<T,DESCRIPTOR> const& _converter;
+public:
+  ParticleU2D(SmoothIndicatorF2D<T,T,true>& indicator, UnitConverter<T,DESCRIPTOR> const& converter);
+  bool operator()(T output[], const S input[]);
+};
+
+
+//////////////////////////////////3D////////////////////////////////////////////
+template <typename T, typename S>
+class AnalyticalComposed3D final : public AnalyticalF3D<T,S> {
+private:
+  AnalyticalF3D<T,S>& _f0;
+  AnalyticalF3D<T,S>& _f1;
+  AnalyticalF3D<T,S>& _f2;
+public:
+  AnalyticalComposed3D(AnalyticalF3D<T,S>& f0, AnalyticalF3D<T,S>& f1, AnalyticalF3D<T,S>& f2);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalConst3D: 3D -> XD, where XD is defined by value.size()
+template <typename T, typename S>
+class AnalyticalConst3D final : public AnalyticalF3D<T,S> {
+private:
+  // is constant return value of operator()
+  std::vector<T> _c;
+public:
+  AnalyticalConst3D(T value);
+  AnalyticalConst3D(T value0, T value1);
+  AnalyticalConst3D(T value0, T value1, T value2);
+  AnalyticalConst3D(const std::vector<T>& value);
+  AnalyticalConst3D(const Vector<T,3>& value);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalLinear3D: 3D -> 1D troughout given points (x0,y0,z0,v0), (x1,y1,z1,v1), (x2,y2,z2,v2), (x3,y3,z3,v3)
+template <typename T, typename S>
+class AnalyticalLinear3D final : public AnalyticalF3D<T,S> {
+protected:
+  T _a;
+  T _b;
+  T _c;
+  T _d;
+public:
+  AnalyticalLinear3D(T a, T b, T c, T d);
+  AnalyticalLinear3D(S x0, S y0, S z0, T v0, S x1, S y1, S z1, T v1, S x2, S y2,
+                     S z2, T v2, S x3, S y3, S z3, T v3);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalRandom3D: 3D -> 1D with random image in (0,1)
+template <typename T, typename S>
+class AnalyticalRandom3D final : public AnalyticalF3D<T,S> {
+public:
+  AnalyticalRandom3D();
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// AnalyticalScaled3D: 3D -> Image(AnalyticalF) scales AnalyticalF by _scale
+template <typename T, typename S>
+class AnalyticalScaled3D final : public AnalyticalF3D<T,S> {
+private:
+  AnalyticalF3D<T,S>& _f;
+  T _scale;
+public:
+  AnalyticalScaled3D(AnalyticalF3D<T,S>& f, T scale);
+  bool operator() (T output[], const S x[]) override;
+};
+
+/// see Mink et al. 2016 in Sec.3.1.
+template <typename T, typename S, typename DESCRIPTOR>
+class PLSsolution3D : public AnalyticalF3D<T,S> {
+private:
+  T _physSigmaEff;
+  T _physDiffusionCoefficient;
+public:
+  PLSsolution3D(RadiativeUnitConverter<T,DESCRIPTOR> const& converter);
+  bool operator()(T output[1], const S x[3]) override;
+};
+
+/// light source as a cylinder along z-axis
+template <typename T, typename S, typename DESCRIPTOR>
+class LightSourceCylindrical3D : public AnalyticalF3D<T,S> {
+private:
+  T _physSigmaEff;
+  T _physDiffusionCoefficient;
+  Vector<T,3> _center;
+public:
+  LightSourceCylindrical3D(RadiativeUnitConverter<T,DESCRIPTOR> const& converter, Vector<T,3> center = {T(0), T(0), T(0)});
+  bool operator()(T output[1], const S x[3]) override;
+};
+
+/**
+* \param _position of light source
+* \param _orientation direction of light source (normalized)
+* \param _falloff is power of the cosine
+*/
+template <typename T, typename S>
+class Spotlight : public AnalyticalF3D<T,S> {
+private:
+  Vector<T,3> const _position;
+  Vector<T,3> const _orientation;
+  T const _falloff;
+public:
+  Spotlight(Vector<T,3> position, Vector<T,3> direction, T falloff);
+  bool operator()(T output[1], const S x[3]) override;
+};
+
+/// 8.6.1 Gauss Hill inital values
+template <typename T, typename S>
+class GaussianHill2D : public AnalyticalF2D<T,S> {
+private:
+  T _sigma;
+  Vector<T,2> _x0;
+  T _c0;
+public:
+  GaussianHill2D(T sigma, Vector<T,2> x0, T c0);
+  bool operator()(T output[1], const S x[2]) override;
+};
+
+/// 8.6.1 Gauss Hill time evolution
+template <typename T, typename S>
+class GaussianHillTimeEvolution2D : public AnalyticalF2D<T,S> {
+private:
+  T _sigma02;
+  T _D;
+  T _t;
+  Vector<T,2> _x0;
+  Vector<T,2> _u;
+  T _c0;
+public:
+  GaussianHillTimeEvolution2D(T sigma0, T D, T t, Vector<T,2> x0, Vector<T,2> u, T c0);
+  bool operator()(T output[1], const S x[2]) override;
+};
+
+/** Computes resulting velocity of an object from translational and rotational velocity.
+ * \param indicator Class defining the object (needs to be a SmoothIndicatorF3D)
+ * \param u translational velocity of the object - expected in lattice units
+ * \param omega rotational velocity of the object - expected in lattice units
+ */
+template <typename T, typename S, typename DESCRIPTOR>
+class ParticleU3D : public AnalyticalF3D<T,S> {
+protected:
+  SmoothIndicatorF3D<T, T, true>& _indicator;
+  UnitConverter<T,DESCRIPTOR> const& _converter;
+public:
+  ParticleU3D(SmoothIndicatorF3D<T, T, true>& indicator, UnitConverter<T,DESCRIPTOR> const& converter);
+  bool operator()(T output[], const S input[]);
+};
+
+} // end namespace olb
+#endif