diff --git a/NumLib/CMakeLists.txt b/NumLib/CMakeLists.txt
index efd3ad4b46bc388a244828cab766c522f4408e83..c1fae72c19152c09009bd6286f485e697b8c4109 100644
--- a/NumLib/CMakeLists.txt
+++ b/NumLib/CMakeLists.txt
@@ -10,7 +10,6 @@ APPEND_SOURCE_FILES(SOURCES Fem/ShapeFunction)
APPEND_SOURCE_FILES(SOURCES TimeStepping)
APPEND_SOURCE_FILES(SOURCES TimeStepping/Algorithms)
APPEND_SOURCE_FILES(SOURCES Function)
-APPEND_SOURCE_FILES(SOURCES Distribution)
APPEND_SOURCE_FILES(SOURCES ODESolver)
APPEND_SOURCE_FILES(SOURCES Extrapolation)
diff --git a/NumLib/Distribution/Distribution.cpp b/NumLib/Distribution/Distribution.cpp
deleted file mode 100644
index 2aef2ec8da4ea9e22badc94e1dbc679952d1f24e..0000000000000000000000000000000000000000
--- a/NumLib/Distribution/Distribution.cpp
+++ /dev/null
@@ -1,36 +0,0 @@
-/**
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#include "Distribution.h"
-
-#include <algorithm>
-
-#include "MeshLib/Mesh.h"
-#include "MeshLib/Node.h"
-#include "NumLib/Function/ISpatialFunction.h"
-
-namespace NumLib
-{
-
-std::vector<double> generateNodeValueDistribution(
- const NumLib::ISpatialFunction &func,
- const MeshLib::Mesh &msh,
- const std::vector<std::size_t> &vec_node_ids)
-{
- // evaluate a given function with nodal coordinates
- std::vector<double> vec_values(vec_node_ids.size());
- std::transform(vec_node_ids.begin(), vec_node_ids.end(), vec_values.begin(),
- [&func, &msh](std::size_t node_id)
- {
- return func(*msh.getNode(node_id));
- });
- return vec_values;
-}
-
-} // namespace NumLib
diff --git a/NumLib/Distribution/Distribution.h b/NumLib/Distribution/Distribution.h
deleted file mode 100644
index 36771244c9704546afc25a323ddc08352735152e..0000000000000000000000000000000000000000
--- a/NumLib/Distribution/Distribution.h
+++ /dev/null
@@ -1,36 +0,0 @@
-/**
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <vector>
-
-namespace MeshLib
-{
-class Mesh;
-}
-
-namespace NumLib
-{
-class ISpatialFunction;
-
-/**
- * Generate distributed node values from a given function
- *
- * @param func a spatial function object
- * @param msh a mesh object
- * @param vec_node_ids a vector of mesh node ids where the function is evaluated
- * @return a vector of nodal values
- */
-std::vector<double> generateNodeValueDistribution(
- const NumLib::ISpatialFunction &func,
- const MeshLib::Mesh &msh,
- const std::vector<std::size_t> &vec_node_ids);
-
-} // namespace NumLib
diff --git a/NumLib/Function/LinearInterpolationAlongPolyline.cpp b/NumLib/Function/LinearInterpolationAlongPolyline.cpp
deleted file mode 100644
index 8eff055a784e4c53aa51107ccc35ab7f7dcee7ab..0000000000000000000000000000000000000000
--- a/NumLib/Function/LinearInterpolationAlongPolyline.cpp
+++ /dev/null
@@ -1,70 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date 2014-03-18
- *
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#include "LinearInterpolationAlongPolyline.h"
-
-#include "GeoLib/Polyline.h"
-#include "MeshLib/Mesh.h"
-#include "MeshGeoToolsLib/MeshNodesAlongPolyline.h"
-
-namespace NumLib
-{
-
-LinearInterpolationAlongPolyline::LinearInterpolationAlongPolyline(
- const GeoLib::Polyline& ply,
- const std::vector<std::size_t>& vec_interpolate_point_ids,
- const std::vector<double>& vec_interpolate_point_values,
- const double search_length, const double default_value)
-: _ply(ply),
- _interpolation(createInterpolation(ply, vec_interpolate_point_ids, vec_interpolate_point_values)),
- _search_length(search_length), _default_value(default_value)
-{}
-
-MathLib::PiecewiseLinearInterpolation LinearInterpolationAlongPolyline::createInterpolation(
- const GeoLib::Polyline& ply,
- const std::vector<std::size_t>& vec_interpolate_point_ids,
- const std::vector<double>& vec_interpolate_point_values)
-{
- std::vector<double> vec_known_dist;
- std::vector<double> vec_known_values;
- vec_known_dist.reserve(vec_interpolate_point_ids.size());
- vec_known_values.reserve(vec_interpolate_point_ids.size());
- for (std::size_t i=0; i<vec_interpolate_point_ids.size(); i++)
- {
- const std::size_t pnt_id = vec_interpolate_point_ids[i];
- if (!ply.isPointIDInPolyline(pnt_id))
- {
- continue;
- }
-
- for (std::size_t j=0; j<ply.getNumberOfPoints(); j++)
- {
- if (pnt_id == ply.getPointID(j))
- {
- vec_known_dist.push_back(ply.getLength(j));
- vec_known_values.push_back(vec_interpolate_point_values[i]);
- break;
- }
- }
- }
-
- return MathLib::PiecewiseLinearInterpolation{std::move(vec_known_dist),
- std::move(vec_known_values)};
-}
-
-double LinearInterpolationAlongPolyline::operator()(const MathLib::Point3d& pnt) const
-{
- const double dist = _ply.getDistanceAlongPolyline(pnt, _search_length);
- return dist>=0 ? _interpolation.getValue(dist) : _default_value;
-}
-
-} // namespace NumLib
diff --git a/NumLib/Function/LinearInterpolationAlongPolyline.h b/NumLib/Function/LinearInterpolationAlongPolyline.h
deleted file mode 100644
index 8e733c589ac30b45ba3d6930499ad99f3df489b4..0000000000000000000000000000000000000000
--- a/NumLib/Function/LinearInterpolationAlongPolyline.h
+++ /dev/null
@@ -1,78 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date 2014-03-18
- *
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <vector>
-
-#include "MathLib/InterpolationAlgorithms/PiecewiseLinearInterpolation.h"
-#include "GeoLib/Point.h"
-#include "ISpatialFunction.h"
-
-namespace GeoLib
-{
-class Polyline;
-}
-
-namespace MeshLib
-{
-class Mesh;
-}
-
-namespace NumLib
-{
-/**
- * Class for linearly interpolating values along a GeoLib::Polyline object
- */
-class LinearInterpolationAlongPolyline : public ISpatialFunction
-{
-public:
- /**
- * Constructor
- * @param ply a polyline object
- * @param vec_interpolate_point_ids a vector of point IDs where values are known
- * @param vec_interpolate_point_values a vector of values at points
- * @param search_length distance threshold to decide whether a point is located on a polyline or not
- * @param default_value default value when a given point is not located on a polyine
- */
- LinearInterpolationAlongPolyline(
- const GeoLib::Polyline& ply,
- const std::vector<std::size_t>& vec_interpolate_point_ids,
- const std::vector<double>& vec_interpolate_point_values,
- const double search_length, const double default_value);
-
- /**
- * interpolate a value at the given point
- * @param pnt a point object
- * @return interpolated value. A default value is returned if the given point
- * is not located on a polyline
- */
- double operator()(const MathLib::Point3d& pnt) const override;
-
-private:
- /// construct an interpolation algorithm
- static MathLib::PiecewiseLinearInterpolation createInterpolation(
- const GeoLib::Polyline& ply,
- const std::vector<std::size_t>& vec_interpolate_point_ids,
- const std::vector<double>& vec_interpolate_point_values);
-
- /// a polyline object
- const GeoLib::Polyline& _ply;
- /// an interpolation algorithm
- const MathLib::PiecewiseLinearInterpolation _interpolation;
- /// search length
- const double _search_length;
- /// default value
- const double _default_value;
-};
-
-} // namespace NumLib
diff --git a/NumLib/Function/LinearInterpolationOnSurface.cpp b/NumLib/Function/LinearInterpolationOnSurface.cpp
deleted file mode 100644
index 67e653cf5fc03f6d4ec8306a3ba9d5b082fb9934..0000000000000000000000000000000000000000
--- a/NumLib/Function/LinearInterpolationOnSurface.cpp
+++ /dev/null
@@ -1,121 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date 2014-03-18
- *
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#include "LinearInterpolationOnSurface.h"
-
-#include <limits>
-#include <array>
-#include <cassert>
-#include <numeric>
-
-#include "GeoLib/Surface.h"
-#include "GeoLib/Triangle.h"
-#include "GeoLib/AnalyticalGeometry.h"
-#include "MathLib/Vector3.h"
-#include "MathLib/GeometricBasics.h"
-#include "MeshLib/Mesh.h"
-#include "MeshLib/Node.h"
-#include "MeshGeoToolsLib/MeshNodesAlongSurface.h"
-
-namespace NumLib
-{
-
-LinearInterpolationOnSurface::LinearInterpolationOnSurface(
- const GeoLib::Surface& sfc,
- const std::vector<std::size_t>& vec_interpolate_point_ids,
- const std::vector<double>& vec_interpolate_point_values,
- const double default_value)
-: _sfc(sfc), _vec_interpolate_point_ids(vec_interpolate_point_ids),
- _vec_interpolate_point_values(vec_interpolate_point_values), _default_value(default_value)
-{
- assert(vec_interpolate_point_ids.size()==vec_interpolate_point_values.size());
-}
-
-double LinearInterpolationOnSurface::operator()(const MathLib::Point3d& pnt) const
-{
- if (!_sfc.isPntInBoundingVolume(pnt, 0))
- {
- return _default_value;
- }
- auto* tri = _sfc.findTriangle(pnt);
- if (tri == nullptr)
- {
- return _default_value;
- }
-
- std::array<double, 3> pnt_values;
- for (unsigned j=0; j<3; j++) {
- auto itr = std::find(_vec_interpolate_point_ids.begin(), _vec_interpolate_point_ids.end(), (*tri)[j]);
- if (itr != _vec_interpolate_point_ids.end()) {
- const std::size_t index = std::distance(_vec_interpolate_point_ids.begin(), itr);
- pnt_values[j] = _vec_interpolate_point_values[index];
- } else {
- pnt_values[j] = _default_value;
- }
- }
- double val = interpolateInTri(*tri, pnt_values.data(), pnt);
- return val;
-}
-
-double LinearInterpolationOnSurface::interpolateInTri(
- const GeoLib::Triangle &tri,
- double const* const vertex_values,
- MathLib::Point3d const& pnt) const
-{
- std::vector<GeoLib::Point> pnts;
- for (unsigned i = 0; i < 3; i++)
- {
- pnts.emplace_back(*tri.getPoint(i));
- }
- pnts.emplace_back(pnt, -1);
- std::vector<GeoLib::Point*> p_pnts = {{&pnts[0], &pnts[1], &pnts[2], &pnts[3]}};
- GeoLib::rotatePointsToXY(p_pnts.begin(), p_pnts.begin()+3, p_pnts.begin(), p_pnts.end());
-
- GeoLib::Point const& v1(pnts[0]);
- GeoLib::Point const& v2(pnts[1]);
- GeoLib::Point const& v3(pnts[2]);
- GeoLib::Point const& v_pnt(pnts[3]);
- const double area = MathLib::calcTriangleArea(v1, v2, v3);
-
- if (area==.0) {
- // take average if all points have the same coordinates
- return std::accumulate(vertex_values, vertex_values+3, 0.0) / 3;
- }
-
- // interpolate as
- // u(x,y) = sum_i[N_i(x,y)*u_i] (i=1,2,3)
- // N_i(x,y) = 1/(2A)*(a_i + b_i*x + c_i*y)
-
- double a[3], b[3], c[3];
- // 1st vertex
- a[0] = 0.5/area*(v2[0]*v3[1]-v3[0]*v2[1]);
- b[0] = 0.5/area*(v2[1]-v3[1]);
- c[0] = 0.5/area*(v3[0]-v2[0]);
- // 2nd vertex
- a[1] = 0.5/area*(v3[0]*v1[1]-v1[0]*v3[1]);
- b[1] = 0.5/area*(v3[1]-v1[1]);
- c[1] = 0.5/area*(v1[0]-v3[0]);
- // 3rd vertex
- a[2] = 0.5/area*(v1[0]*v2[1]-v2[0]*v1[1]);
- b[2] = 0.5/area*(v1[1]-v2[1]);
- c[2] = 0.5/area*(v2[0]-v1[0]);
-
- double val = .0;
- for (unsigned i = 0; i < 3; i++)
- {
- val += (a[i] + b[i] * v_pnt[0] + c[i] * v_pnt[1]) * vertex_values[i];
- }
-
- return val;
-}
-
-} // namespace NumLib
diff --git a/NumLib/Function/LinearInterpolationOnSurface.h b/NumLib/Function/LinearInterpolationOnSurface.h
deleted file mode 100644
index ac020a81393c0dd448fad9607c06805094ec0956..0000000000000000000000000000000000000000
--- a/NumLib/Function/LinearInterpolationOnSurface.h
+++ /dev/null
@@ -1,97 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date 2014-03-18
- *
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <vector>
-
-#include "GeoLib/Point.h"
-#include "ISpatialFunction.h"
-
-namespace GeoLib
-{
-class Surface;
-class Triangle;
-}
-
-namespace MeshLib
-{
-class Mesh;
-}
-
-namespace NumLib
-{
-/**
- * Class for linearly interpolating values along a GeoLib::Surface object
- */
-class LinearInterpolationOnSurface : public ISpatialFunction
-{
-public:
- /**
- * Constructor
- * @param sfc a surface object
- * @param vec_interpolate_point_ids a vector of point IDs where values are known
- * @param vec_interpolate_point_values a vector of values at points
- * @param default_value default value when a given point is not located on a surface
- */
- LinearInterpolationOnSurface(
- const GeoLib::Surface& sfc,
- const std::vector<std::size_t>& vec_interpolate_point_ids,
- const std::vector<double>& vec_interpolate_point_values,
- const double default_value);
-
- /**
- * interpolate a value at the given point
- * @param pnt a point object
- * @return interpolated value. A default value is returned if the given point
- * is not located on a surface
- */
- double operator()(const MathLib::Point3d& pnt) const override;
-
-private:
- /// rotate a triangle to XY plane
- void rotate(std::vector<GeoLib::Point> &_pnts) const;
-
- /**
- * do an interpolation within a triangle. Interpolation is done by shape functions
- * for triangle elements based on areal coordinates as
- * \f[
- * u(x,y) = sum_i[N_i(x,y)*u_i] (i=1,2,3)
- * \f]
- * where \f$N_i\f$ is a shape function for node i, \f$u_i\f$ is a value at node i.
- * The shape function is given by
- * \f[
- * N_i(x,y) = 1/(2A)*(a_i + b_i*x + c_i*y)
- * \f]
- * with an element area \f$A\f$ and geometric parameters \f$a_i\f$, \f$b_i\f$ and \f$c_i\f$.
- * reference: http://kratos-wiki.cimne.upc.edu/index.php/Two-dimensional_Shape_Functions
- *
- * @param tri
- * @param vertex_values
- * @param pnt
- * @return
- */
- double interpolateInTri(const GeoLib::Triangle &tri,
- double const* const vertex_values,
- MathLib::Point3d const& pnt) const;
-
- /// a surface object
- const GeoLib::Surface& _sfc;
- /// a vector of point IDs where values are known
- const std::vector<std::size_t>& _vec_interpolate_point_ids;
- /// a vector of values at points
- const std::vector<double>& _vec_interpolate_point_values;
- /// default value
- const double _default_value;
-};
-
-} // namespace NumLib
diff --git a/NumLib/Function/SpatialFunctionLinear.h b/NumLib/Function/SpatialFunctionLinear.h
deleted file mode 100644
index 16d673d1ee15f0badca2e35f9f0bfac5cf875b02..0000000000000000000000000000000000000000
--- a/NumLib/Function/SpatialFunctionLinear.h
+++ /dev/null
@@ -1,24 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date 2013-08-13
- *
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include "MathLib/LinearFunction.h"
-#include "TemplateSpatialFunction.h"
-
-namespace NumLib
-{
-
-/// Representation of linear functions f: R^3 -> R; f({x, y, z}) = a + bx + cy + dz
-using SpatialFunctionLinear =
- TemplateSpatialFunction<MathLib::LinearFunction<double, 3>>;
-}
diff --git a/NumLib/Function/TemplateSpatialFunction.h b/NumLib/Function/TemplateSpatialFunction.h
deleted file mode 100644
index 6535669311b7a6fe52d15644d5f5102a383f7796..0000000000000000000000000000000000000000
--- a/NumLib/Function/TemplateSpatialFunction.h
+++ /dev/null
@@ -1,51 +0,0 @@
-/**
- * \author Norihiro Watanabe
- * \date 2013-08-13
- *
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#pragma once
-
-#include <utility>
-
-#include "ISpatialFunction.h"
-
-namespace NumLib
-{
-
-/**
- * Template class to relate pure mathematical functions with ISpatialFunction
- *
- * @tparam T_FUNCTION Function object which takes "double const* const" as an argument.
- */
-template <class T_FUNCTION>
-class TemplateSpatialFunction : public ISpatialFunction
-{
-public:
- /**
- * Constructor
- * @param f a function object
- */
- explicit TemplateSpatialFunction(T_FUNCTION f) : _f(std::move(f)) {}
- /**
- * evaluate a function
- * @param pnt a point object
- * @return evaluated value
- */
- double operator()(const MathLib::Point3d& pnt) const override
- {
- return _f(pnt.getCoords());
- }
-
-private:
- /// a function object
- const T_FUNCTION _f;
-};
-
-} // namespace NumLib
diff --git a/Tests/NumLib/TestDistribution.cpp b/Tests/NumLib/TestDistribution.cpp
deleted file mode 100644
index 5afc6b503a37da0a377a37aca5f80703f881454e..0000000000000000000000000000000000000000
--- a/Tests/NumLib/TestDistribution.cpp
+++ /dev/null
@@ -1,235 +0,0 @@
-/**
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#include <vector>
-#include <memory>
-
-#include <gtest/gtest.h>
-
-#include "GeoLib/GEOObjects.h"
-#include "MeshLib/Mesh.h"
-#include "MeshLib/MeshGenerators/MeshGenerator.h"
-#include "MeshGeoToolsLib/HeuristicSearchLength.h"
-#include "MeshGeoToolsLib/MeshNodeSearcher.h"
-#include "MeshGeoToolsLib/MeshNodesAlongPolyline.h"
-#include "MeshGeoToolsLib/MeshNodesAlongSurface.h"
-#include "NumLib/Function/LinearInterpolationAlongPolyline.h"
-#include "NumLib/Function/LinearInterpolationOnSurface.h"
-#include "NumLib/Function/SpatialFunctionLinear.h"
-#include "NumLib/Distribution/Distribution.h"
-
-#include "Tests/TestTools.h"
-
-class NumLibDistributionQuad : public testing::Test
-{
-public:
- NumLibDistributionQuad()
- : _msh(MeshLib::MeshGenerator::generateRegularQuadMesh(
- _geometric_size, _number_of_subdivisions_per_direction)),
- _project_name("test"),
- _mshNodesSearcher(*_msh,
- std::make_unique<MeshGeoToolsLib::SearchLength>(),
- MeshGeoToolsLib::SearchAllNodes::Yes)
- {
- // create geometry
- auto pnts = std::make_unique<std::vector<GeoLib::Point*>>();
- pnts->push_back(new GeoLib::Point(0.0, 0.0, 0.0));
- pnts->push_back(new GeoLib::Point(_geometric_size, 0.0, 0.0));
- pnts->push_back(new GeoLib::Point(_geometric_size, _geometric_size, 0.0));
- pnts->push_back(new GeoLib::Point(0.0, _geometric_size, 0.0));
-
- auto plys = std::make_unique<std::vector<GeoLib::Polyline*>>();
- _ply0 = new GeoLib::Polyline(*pnts);
- _ply0->addPoint(0);
- _ply0->addPoint(1);
- plys->push_back(_ply0);
-
- auto* ply1 = new GeoLib::Polyline(*pnts);
- ply1->addPoint(0);
- ply1->addPoint(1);
- ply1->addPoint(2);
- ply1->addPoint(3);
- ply1->addPoint(0);
- plys->push_back(ply1);
-
- auto sfcs = std::make_unique<std::vector<GeoLib::Surface*>>();
- _sfc1 = new GeoLib::Surface(*pnts);
- _sfc1->addTriangle(0, 1, 2);
- _sfc1->addTriangle(0, 2, 3);
- sfcs->push_back(_sfc1);
-
- _geo_objs.addPointVec(std::move(pnts), _project_name);
- _geo_objs.addPolylineVec(std::move(plys), _project_name);
- _geo_objs.addSurfaceVec(std::move(sfcs), _project_name);
- }
-
-protected:
- const double _geometric_size{10.0};
- const std::size_t _number_of_subdivisions_per_direction{10};
- std::unique_ptr<MeshLib::Mesh> _msh;
- GeoLib::GEOObjects _geo_objs;
- std::string _project_name;
- MeshGeoToolsLib::MeshNodeSearcher _mshNodesSearcher;
- GeoLib::Polyline* _ply0{nullptr};
- GeoLib::Surface* _sfc1;
-};
-
-class NumLibDistributionHex : public testing::Test
-{
-public:
- NumLibDistributionHex()
- : _msh(MeshLib::MeshGenerator::generateRegularHexMesh(
- _geometric_size, _number_of_subdivisions_per_direction)),
- _project_name("test"),
- _mshNodesSearcher(*_msh,
- std::make_unique<MeshGeoToolsLib::SearchLength>(),
- MeshGeoToolsLib::SearchAllNodes::Yes)
- {
- // create geometry
- auto pnts = std::make_unique<std::vector<GeoLib::Point*>>();
- pnts->push_back(new GeoLib::Point(0.0, 0.0, 0.0));
- pnts->push_back(new GeoLib::Point(_geometric_size, 0.0, 0.0));
- pnts->push_back(new GeoLib::Point(_geometric_size, _geometric_size, 0.0));
- pnts->push_back(new GeoLib::Point(0.0, _geometric_size, 0.0));
- pnts->push_back(new GeoLib::Point(0.0, 0.0, _geometric_size));
- pnts->push_back(new GeoLib::Point(_geometric_size, 0.0, _geometric_size));
- pnts->push_back(new GeoLib::Point(_geometric_size, _geometric_size, _geometric_size));
- pnts->push_back(new GeoLib::Point(0.0, _geometric_size, _geometric_size));
-
- auto plys = std::make_unique<std::vector<GeoLib::Polyline*>>();
- _ply0 = new GeoLib::Polyline(*pnts); // vertical polyline
- _ply0->addPoint(0);
- _ply0->addPoint(4);
- plys->push_back(_ply0);
- auto* ply1 = new GeoLib::Polyline(*pnts); // polygon for left surface
- ply1->addPoint(0);
- ply1->addPoint(3);
- ply1->addPoint(7);
- ply1->addPoint(4);
- ply1->addPoint(0);
- plys->push_back(ply1);
-
- auto sfcs = std::make_unique<std::vector<GeoLib::Surface*>>();
- _sfc1 = new GeoLib::Surface(*pnts);
- _sfc1->addTriangle(0, 3, 7);
- _sfc1->addTriangle(0, 7, 4);
- sfcs->push_back(_sfc1);
-
- _geo_objs.addPointVec(std::move(pnts) ,_project_name);
- _geo_objs.addPolylineVec(std::move(plys), _project_name);
- _geo_objs.addSurfaceVec(std::move(sfcs), _project_name);
- }
-
-protected:
- const double _geometric_size{10.0};
- const std::size_t _number_of_subdivisions_per_direction{10};
- std::unique_ptr<MeshLib::Mesh> _msh;
- GeoLib::GEOObjects _geo_objs;
- std::string _project_name;
- MeshGeoToolsLib::MeshNodeSearcher _mshNodesSearcher;
- GeoLib::Polyline* _ply0{nullptr};
- GeoLib::Surface* _sfc1;
-};
-
-TEST_F(NumLibDistributionQuad, Linear)
-{
- // f(x,y,z) = 1 + 2x + 3y + 4z
- std::array<double,4> f_coeff = {{1, 2, 3, 4}};
- MathLib::LinearFunction<double,3> linear_f(f_coeff);
- std::vector<double> expected = {{1, 3, 5, 7, 9, 11, 13, 15, 17, 19, 21}};
-
- NumLib::SpatialFunctionLinear f(linear_f);
- const std::vector<std::size_t>& vec_node_ids = _mshNodesSearcher.getMeshNodesAlongPolyline(*_ply0).getNodeIDs();
- std::vector<double> interpolated_values = NumLib::generateNodeValueDistribution(f, *_msh, vec_node_ids);
- ASSERT_ARRAY_NEAR(expected, interpolated_values, expected.size(), std::numeric_limits<double>::epsilon());
-}
-
-TEST_F(NumLibDistributionHex, Linear)
-{
- // f(x,y,z) = 1 + 2x + 3y + 4z
- std::array<double,4> f_coeff = {{1, 2, 3, 4}};
- MathLib::LinearFunction<double,3> linear_f(f_coeff);
- std::vector<double> expected(static_cast<std::size_t>(std::pow(_number_of_subdivisions_per_direction+1, 2)));
- const double dL = _geometric_size / _number_of_subdivisions_per_direction;
- for (std::size_t i=0; i<expected.size(); i++) {
- double x = 0;
- double y = (i%(_number_of_subdivisions_per_direction+1)) * dL;
- double z = (i/(_number_of_subdivisions_per_direction+1)) * dL;
- expected[i] = f_coeff[0] + f_coeff[1]*x + f_coeff[2]*y + f_coeff[3]*z;
- }
-
- NumLib::SpatialFunctionLinear f(linear_f);
- const std::vector<std::size_t>& vec_node_ids = _mshNodesSearcher.getMeshNodesAlongSurface(*_sfc1).getNodeIDs();
- std::vector<double> interpolated_values = NumLib::generateNodeValueDistribution(f, *_msh, vec_node_ids);
-
- ASSERT_ARRAY_NEAR(expected, interpolated_values, expected.size(), std::numeric_limits<double>::epsilon());
-}
-
-TEST_F(NumLibDistributionQuad, InterpolationPolyline)
-{
- const std::vector<std::size_t> vec_point_ids = {{0, 1}};
- const std::vector<double> vec_point_values = {{0., 100.}};
- std::vector<double> expected = {{0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100}};
-
- NumLib::LinearInterpolationAlongPolyline interpolate(
- *_ply0, vec_point_ids, vec_point_values, std::numeric_limits<double>::epsilon(), std::numeric_limits<double>::max());
- const std::vector<std::size_t>& vec_node_ids = _mshNodesSearcher.getMeshNodesAlongPolyline(*_ply0).getNodeIDs();
- std::vector<double> interpolated_values = NumLib::generateNodeValueDistribution(interpolate, *_msh, vec_node_ids);
-
- ASSERT_ARRAY_NEAR(expected, interpolated_values, expected.size(), std::numeric_limits<double>::epsilon());
-}
-
-TEST_F(NumLibDistributionHex, InterpolationPolyline)
-{
- const std::vector<std::size_t> vec_point_ids = {{0, 4}};
- const std::vector<double> vec_point_values = {{0., 100.}};
- std::vector<double> expected = {{0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100}};
-
- NumLib::LinearInterpolationAlongPolyline interpolate(
- *_ply0, vec_point_ids, vec_point_values, std::numeric_limits<double>::epsilon(), std::numeric_limits<double>::max());
- const std::vector<std::size_t>& vec_node_ids = _mshNodesSearcher.getMeshNodesAlongPolyline(*_ply0).getNodeIDs();
- std::vector<double> interpolated_values = NumLib::generateNodeValueDistribution(interpolate, *_msh, vec_node_ids);
-
- ASSERT_ARRAY_NEAR(expected, interpolated_values, expected.size(), std::numeric_limits<double>::epsilon());
-}
-
-TEST_F(NumLibDistributionQuad, InterpolationSurface)
-{
- const std::vector<std::size_t> vec_point_ids = {{0, 1, 2, 3}};
- const std::vector<double> vec_point_values = {{0., 100., 100., 0.}};
- std::vector<double> expected(_msh->getNumberOfNodes());
- for (std::size_t i=0; i<_msh->getNumberOfNodes(); i++) {
- expected[i] = static_cast<double>((i%(_number_of_subdivisions_per_direction+1)) * 10);
- }
-
- NumLib::LinearInterpolationOnSurface interpolate(*_sfc1, vec_point_ids, vec_point_values, std::numeric_limits<double>::max());
- const std::vector<std::size_t>& vec_node_ids = _mshNodesSearcher.getMeshNodesAlongSurface(*_sfc1).getNodeIDs();
- std::vector<double> interpolated_values = NumLib::generateNodeValueDistribution(interpolate, *_msh, vec_node_ids);
-
- // the machine epsilon for double is too small for this test
- ASSERT_ARRAY_NEAR(expected, interpolated_values, expected.size(), std::numeric_limits<float>::epsilon());
-}
-
-TEST_F(NumLibDistributionHex, InterpolationSurface)
-{
- const std::vector<std::size_t> vec_point_ids = {{0, 3, 7, 4}};
- const std::vector<double> vec_point_values = {{0., 100., 100., 0.}};
- std::vector<double> expected(static_cast<std::size_t>(std::pow(_number_of_subdivisions_per_direction+1, 2)));
- for (std::size_t i=0; i<expected.size(); i++) {
- expected[i] = static_cast<double>((i%(_number_of_subdivisions_per_direction+1)) * 10);
- }
-
- NumLib::LinearInterpolationOnSurface interpolate(*_sfc1, vec_point_ids, vec_point_values, std::numeric_limits<double>::max());
- const std::vector<std::size_t>& vec_node_ids = _mshNodesSearcher.getMeshNodesAlongSurface(*_sfc1).getNodeIDs();
- std::vector<double> interpolated_values = NumLib::generateNodeValueDistribution(interpolate, *_msh, vec_node_ids);
-
- ASSERT_ARRAY_NEAR(expected, interpolated_values, expected.size(), std::numeric_limits<float>::epsilon());
-}
-
-
diff --git a/Tests/NumLib/TestSpatialFunction.cpp b/Tests/NumLib/TestSpatialFunction.cpp
deleted file mode 100644
index 5307bd20373ecc43c3ab994ac0a5b000cdbcd746..0000000000000000000000000000000000000000
--- a/Tests/NumLib/TestSpatialFunction.cpp
+++ /dev/null
@@ -1,117 +0,0 @@
-/**
- * \copyright
- * Copyright (c) 2012-2019, OpenGeoSys Community (http://www.opengeosys.org)
- * Distributed under a Modified BSD License.
- * See accompanying file LICENSE.txt or
- * http://www.opengeosys.org/project/license
- *
- */
-
-#include <vector>
-#include <memory>
-
-#include <gtest/gtest.h>
-
-#include "GeoLib/GEOObjects.h"
-#include "MeshLib/Mesh.h"
-#include "MeshLib/MeshGenerators/MeshGenerator.h"
-#include "NumLib/Function/LinearInterpolationAlongPolyline.h"
-#include "NumLib/Function/LinearInterpolationOnSurface.h"
-#include "NumLib/Function/SpatialFunctionLinear.h"
-
-#include "Tests/TestTools.h"
-
-TEST(NumLib, SpatialFunctionLinear)
-{
- // f(x,y,z) = 1 + 2x + 3y + 4z
- std::array<double,4> f_coeff = {{1, 2, 3, 4}};
- MathLib::LinearFunction<double,3> linear_f(f_coeff);
-
- NumLib::SpatialFunctionLinear f(linear_f);
-
- ASSERT_DOUBLE_EQ(1., f(GeoLib::Point(0,0,0)));
- ASSERT_DOUBLE_EQ(10., f(GeoLib::Point(1,1,1)));
- ASSERT_DOUBLE_EQ(-8, f(GeoLib::Point(-1,-1,-1)));
- for (std::size_t k(0); k < 5; ++k)
- {
- GeoLib::Point pt(0, 0, 0);
- for (unsigned i = 0; i < 3; ++i)
- {
- pt[i] = (double)rand() - (RAND_MAX / 2.0);
- }
- double expected = 1+2*pt[0]+3*pt[1]+4*pt[2];
- ASSERT_DOUBLE_EQ(expected, f(pt));
- }
-}
-
-TEST(NumLib, SpatialFunctionInterpolationPolyline)
-{
- // create geometry
- GeoLib::Point pt1(0.0, 0.0, 0.0);
- GeoLib::Point pt2(10.0, 0.0, 0.0);
- std::vector<GeoLib::Point*> pnts = {&pt1, &pt2};
- GeoLib::Polyline ply0(pnts);
- ply0.addPoint(0);
- ply0.addPoint(1);
-
- // define a function
- const std::vector<std::size_t> vec_point_ids = {{0, 1}};
- const std::vector<double> vec_point_values = {{0., 100.}};
- NumLib::LinearInterpolationAlongPolyline interpolate(
- ply0, vec_point_ids, vec_point_values, std::numeric_limits<double>::epsilon(), std::numeric_limits<double>::max());
-
- // normal
- for (unsigned k = 0; k < 10; k++)
- {
- ASSERT_DOUBLE_EQ(k * 10., interpolate(GeoLib::Point(k, 0, 0)));
- }
-
- // failure
- // x
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(-1,0,0)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(11,0,0)));
- // y
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(0,1,0)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(0,-1,0)));
- // z
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(0,0,1)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(0,0,-1)));
-}
-
-TEST(NumLib, SpatialFunctionInterpolationSurface)
-{
- // create geometry
- GeoLib::Point pt1(0.0, 0.0, 0.0);
- GeoLib::Point pt2(10.0, 0.0, 0.0);
- GeoLib::Point pt3(10.0, 10.0, 0.0);
- GeoLib::Point pt4(0.0, 10.0, 0.0);
- std::vector<GeoLib::Point*> pnts = {&pt1, &pt2, &pt3, &pt4};
- GeoLib::Surface sfc1(pnts);
- sfc1.addTriangle(0, 1, 2);
- sfc1.addTriangle(0, 2, 3);
-
- // define a function
- const std::vector<std::size_t> vec_point_ids = {{0, 1, 2, 3}};
- const std::vector<double> vec_point_values = {{0., 100., 100., 0.}};
- NumLib::LinearInterpolationOnSurface interpolate(
- sfc1, vec_point_ids, vec_point_values,
- std::numeric_limits<double>::max());
-
- // normal
- for (unsigned k=0; k<10; k++) {
- ASSERT_DOUBLE_EQ(k*10., interpolate(GeoLib::Point(k,0,0)));
- ASSERT_DOUBLE_EQ(k*10., interpolate(GeoLib::Point(k,5,0)));
- ASSERT_DOUBLE_EQ(k*10., interpolate(GeoLib::Point(k,10,0)));
- }
-
- // failure
- // x, y
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(-1,-1,0)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(11,-1,0)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(11,11,0)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(-1,11,0)));
- // z
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(0,0,1)));
- ASSERT_DOUBLE_EQ(std::numeric_limits<double>::max(), interpolate(GeoLib::Point(0,0,-1)));
-}
-