diff --git a/src/language/modules/SchemeModule.cpp b/src/language/modules/SchemeModule.cpp
index e4b52dc2277b95381213bff69af42749783fb655..08df9ba86ecb63a18bbb9e55fae1adfa544fcb23 100644
--- a/src/language/modules/SchemeModule.cpp
+++ b/src/language/modules/SchemeModule.cpp
@@ -91,6 +91,27 @@ SchemeModule::SchemeModule()
));
+ this->_addBuiltinFunction("integrate",
+ std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
+ const IDiscreteFunction>(std::shared_ptr<const IMesh>,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>&,
+ std::shared_ptr<const IQuadratureDescriptor>,
+ std::shared_ptr<const IDiscreteFunctionDescriptor>,
+ const std::vector<FunctionSymbolId>&)>>(
+ [](std::shared_ptr<const IMesh> mesh,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>& integration_zone_list,
+ std::shared_ptr<const IQuadratureDescriptor> quadrature_descriptor,
+ std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
+ const std::vector<FunctionSymbolId>& function_id_list)
+ -> std::shared_ptr<const IDiscreteFunction> {
+ return DiscreteFunctionVectorIntegrator{mesh, integration_zone_list,
+ quadrature_descriptor,
+ discrete_function_descriptor, function_id_list}
+ .integrate();
+ }
+
+ ));
+
this
->_addBuiltinFunction("integrate",
std::make_shared<BuiltinFunctionEmbedder<
@@ -110,6 +131,20 @@ SchemeModule::SchemeModule()
));
+ this->_addBuiltinFunction(
+ "integrate",
+ std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
+ const IDiscreteFunction>(std::shared_ptr<const IMesh>, const std::vector<std::shared_ptr<const IZoneDescriptor>>&,
+ std::shared_ptr<const IQuadratureDescriptor>, const FunctionSymbolId&)>>(
+ [](std::shared_ptr<const IMesh> mesh,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>& integration_zone_list,
+ std::shared_ptr<const IQuadratureDescriptor> quadrature_descriptor,
+ const FunctionSymbolId& function_id) -> std::shared_ptr<const IDiscreteFunction> {
+ return DiscreteFunctionIntegrator{mesh, integration_zone_list, quadrature_descriptor, function_id}.integrate();
+ }
+
+ ));
+
this->_addBuiltinFunction("integrate",
std::make_shared<BuiltinFunctionEmbedder<
std::shared_ptr<const IDiscreteFunction>(std::shared_ptr<const IMesh>,
@@ -136,6 +171,24 @@ SchemeModule::SchemeModule()
));
+ this->_addBuiltinFunction("interpolate",
+ std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
+ const IDiscreteFunction>(std::shared_ptr<const IMesh>,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>&,
+ std::shared_ptr<const IDiscreteFunctionDescriptor>,
+ const std::vector<FunctionSymbolId>&)>>(
+ [](std::shared_ptr<const IMesh> mesh,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>& interpolation_zone_list,
+ std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
+ const std::vector<FunctionSymbolId>& function_id_list)
+ -> std::shared_ptr<const IDiscreteFunction> {
+ return DiscreteFunctionVectorInterpoler{mesh, interpolation_zone_list,
+ discrete_function_descriptor, function_id_list}
+ .interpolate();
+ }
+
+ ));
+
this->_addBuiltinFunction(
"interpolate",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
@@ -159,24 +212,6 @@ SchemeModule::SchemeModule()
));
- this->_addBuiltinFunction("interpolate",
- std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
- const IDiscreteFunction>(std::shared_ptr<const IMesh>,
- const std::vector<std::shared_ptr<const IZoneDescriptor>>&,
- std::shared_ptr<const IDiscreteFunctionDescriptor>,
- const std::vector<FunctionSymbolId>&)>>(
- [](std::shared_ptr<const IMesh> mesh,
- const std::vector<std::shared_ptr<const IZoneDescriptor>>& interpolation_zone_list,
- std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
- const std::vector<FunctionSymbolId>& function_id_list)
- -> std::shared_ptr<const IDiscreteFunction> {
- return DiscreteFunctionVectorInterpoler{mesh, interpolation_zone_list,
- discrete_function_descriptor, function_id_list}
- .interpolate();
- }
-
- ));
-
this->_addBuiltinFunction("interpolate",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IDiscreteFunction>(std::shared_ptr<const IMesh>,
diff --git a/src/scheme/DiscreteFunctionIntegrator.cpp b/src/scheme/DiscreteFunctionIntegrator.cpp
index 91cf921df2764ad13f4470687fc4cc5f4d6bced1..fe0757cc0ef74e52b887e751fa678e7be6bdb186 100644
--- a/src/scheme/DiscreteFunctionIntegrator.cpp
+++ b/src/scheme/DiscreteFunctionIntegrator.cpp
@@ -1,13 +1,75 @@
#include <scheme/DiscreteFunctionIntegrator.hpp>
#include <language/utils/IntegrateCellValue.hpp>
+#include <mesh/MeshCellZone.hpp>
#include <scheme/DiscreteFunctionP0.hpp>
#include <utils/Exceptions.hpp>
template <size_t Dimension, typename DataType, typename ValueType>
std::shared_ptr<IDiscreteFunction>
-DiscreteFunctionIntegrator::_integrate() const
+DiscreteFunctionIntegrator::_integrateOnZoneList() const
+{
+ static_assert(std::is_convertible_v<DataType, ValueType>);
+ Assert(m_zone_list.size() > 0, "no zone list provided");
+
+ using MeshType = Mesh<Connectivity<Dimension>>;
+
+ std::shared_ptr p_mesh = std::dynamic_pointer_cast<const MeshType>(m_mesh);
+
+ CellValue<bool> is_in_zone{p_mesh->connectivity()};
+ is_in_zone.fill(false);
+
+ size_t number_of_cells = 0;
+ for (const auto& zone : m_zone_list) {
+ auto mesh_cell_zone = getMeshCellZone(*p_mesh, *zone);
+ const auto& cell_list = mesh_cell_zone.cellList();
+ for (size_t i_cell = 0; i_cell < cell_list.size(); ++i_cell) {
+ const CellId cell_id = cell_list[i_cell];
+ if (is_in_zone[cell_id]) {
+ std::ostringstream os;
+ os << "cell " << cell_id << " (number " << p_mesh->connectivity().cellNumber()[cell_id]
+ << ") is present multiple times in zone list";
+ throw NormalError(os.str());
+ }
+ ++number_of_cells;
+ is_in_zone[cell_id] = true;
+ }
+ }
+
+ Array<CellId> zone_cell_list{number_of_cells};
+ {
+ size_t i_cell = 0;
+ for (CellId cell_id = 0; cell_id < p_mesh->numberOfCells(); ++cell_id) {
+ if (is_in_zone[cell_id]) {
+ zone_cell_list[i_cell++] = cell_id;
+ }
+ }
+ }
+
+ Array<DataType> array =
+ IntegrateCellValue<DataType(TinyVector<Dimension>)>::template integrate<MeshType>(m_function_id,
+ *m_quadrature_descriptor, *p_mesh,
+ zone_cell_list);
+
+ CellValue<ValueType> cell_value{p_mesh->connectivity()};
+ if constexpr (is_tiny_vector_v<ValueType> or is_tiny_matrix_v<ValueType>) {
+ cell_value.fill(zero);
+ } else {
+ cell_value.fill(0);
+ }
+
+ parallel_for(
+ zone_cell_list.size(), PUGS_LAMBDA(const size_t i_cell) { cell_value[zone_cell_list[i_cell]] = array[i_cell]; });
+
+ return std::make_shared<DiscreteFunctionP0<Dimension, ValueType>>(p_mesh, cell_value);
+}
+
+template <size_t Dimension, typename DataType, typename ValueType>
+std::shared_ptr<IDiscreteFunction>
+DiscreteFunctionIntegrator::_integrateGlobally() const
{
+ Assert(m_zone_list.size() == 0, "invalid call when zones are defined");
+
using MeshType = Mesh<Connectivity<Dimension>>;
std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(m_mesh);
@@ -33,6 +95,17 @@ DiscreteFunctionIntegrator::_integrate() const
}
}
+template <size_t Dimension, typename DataType, typename ValueType>
+std::shared_ptr<IDiscreteFunction>
+DiscreteFunctionIntegrator::_integrate() const
+{
+ if (m_zone_list.size() == 0) {
+ return this->_integrateGlobally<Dimension, DataType, ValueType>();
+ } else {
+ return this->_integrateOnZoneList<Dimension, DataType, ValueType>();
+ }
+}
+
template <size_t Dimension>
std::shared_ptr<IDiscreteFunction>
DiscreteFunctionIntegrator::_integrate() const
diff --git a/src/scheme/DiscreteFunctionIntegrator.hpp b/src/scheme/DiscreteFunctionIntegrator.hpp
index 1784f891eaa6fd7471b17006b3c9f277c31e51dd..21a461d54561d67f1f745d3257ec76ee514d264a 100644
--- a/src/scheme/DiscreteFunctionIntegrator.hpp
+++ b/src/scheme/DiscreteFunctionIntegrator.hpp
@@ -4,6 +4,7 @@
#include <analysis/IQuadratureDescriptor.hpp>
#include <language/utils/FunctionSymbolId.hpp>
#include <mesh/IMesh.hpp>
+#include <mesh/IZoneDescriptor.hpp>
#include <scheme/IDiscreteFunction.hpp>
#include <memory>
@@ -12,9 +13,16 @@ class DiscreteFunctionIntegrator
{
private:
std::shared_ptr<const IMesh> m_mesh;
+ const std::vector<std::shared_ptr<const IZoneDescriptor>> m_zone_list;
std::shared_ptr<const IQuadratureDescriptor> m_quadrature_descriptor;
const FunctionSymbolId m_function_id;
+ template <size_t Dimension, typename DataType, typename ValueType = DataType>
+ std::shared_ptr<IDiscreteFunction> _integrateOnZoneList() const;
+
+ template <size_t Dimension, typename DataType, typename ValueType = DataType>
+ std::shared_ptr<IDiscreteFunction> _integrateGlobally() const;
+
template <size_t Dimension, typename DataType, typename ValueType = DataType>
std::shared_ptr<IDiscreteFunction> _integrate() const;
@@ -30,6 +38,13 @@ class DiscreteFunctionIntegrator
: m_mesh{mesh}, m_quadrature_descriptor{quadrature_descriptor}, m_function_id{function_id}
{}
+ DiscreteFunctionIntegrator(const std::shared_ptr<const IMesh>& mesh,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>& zone_list,
+ const std::shared_ptr<const IQuadratureDescriptor>& quadrature_descriptor,
+ const FunctionSymbolId& function_id)
+ : m_mesh{mesh}, m_zone_list{zone_list}, m_quadrature_descriptor{quadrature_descriptor}, m_function_id{function_id}
+ {}
+
DiscreteFunctionIntegrator(const DiscreteFunctionIntegrator&) = delete;
DiscreteFunctionIntegrator(DiscreteFunctionIntegrator&&) = delete;
diff --git a/src/scheme/DiscreteFunctionVectorIntegrator.cpp b/src/scheme/DiscreteFunctionVectorIntegrator.cpp
index 6d8937f84aeba6ae96ac5fc6e5ba523e5781f239..105dd0ff9ef05a9af8bf6389ac886f810421f6d1 100644
--- a/src/scheme/DiscreteFunctionVectorIntegrator.cpp
+++ b/src/scheme/DiscreteFunctionVectorIntegrator.cpp
@@ -1,14 +1,76 @@
#include <scheme/DiscreteFunctionVectorIntegrator.hpp>
#include <language/utils/IntegrateCellArray.hpp>
-
+#include <mesh/MeshCellZone.hpp>
#include <scheme/DiscreteFunctionP0Vector.hpp>
#include <utils/Exceptions.hpp>
template <size_t Dimension, typename DataType>
std::shared_ptr<IDiscreteFunction>
-DiscreteFunctionVectorIntegrator::_integrate() const
+DiscreteFunctionVectorIntegrator::_integrateOnZoneList() const
+{
+ Assert(m_zone_list.size() > 0, "no zone list provided");
+
+ std::shared_ptr p_mesh = std::dynamic_pointer_cast<const Mesh<Connectivity<Dimension>>>(m_mesh);
+
+ CellValue<bool> is_in_zone{p_mesh->connectivity()};
+ is_in_zone.fill(false);
+
+ size_t number_of_cells = 0;
+ for (const auto& zone : m_zone_list) {
+ auto mesh_cell_zone = getMeshCellZone(*p_mesh, *zone);
+ const auto& cell_list = mesh_cell_zone.cellList();
+ for (size_t i_cell = 0; i_cell < cell_list.size(); ++i_cell) {
+ const CellId cell_id = cell_list[i_cell];
+ if (is_in_zone[cell_id]) {
+ std::ostringstream os;
+ os << "cell " << cell_id << " (number " << p_mesh->connectivity().cellNumber()[cell_id]
+ << ") is present multiple times in zone list";
+ throw NormalError(os.str());
+ }
+ ++number_of_cells;
+ is_in_zone[cell_id] = true;
+ }
+ }
+
+ Array<CellId> zone_cell_list{number_of_cells};
+ {
+ size_t i_cell = 0;
+ for (CellId cell_id = 0; cell_id < p_mesh->numberOfCells(); ++cell_id) {
+ if (is_in_zone[cell_id]) {
+ zone_cell_list[i_cell++] = cell_id;
+ }
+ }
+ }
+
+ Table<const DataType> table =
+ IntegrateCellArray<DataType(TinyVector<Dimension>)>::template integrate(m_function_id_list,
+ *m_quadrature_descriptor, *p_mesh,
+ zone_cell_list);
+
+ CellArray<DataType> cell_array{p_mesh->connectivity(), m_function_id_list.size()};
+ if constexpr (is_tiny_vector_v<DataType> or is_tiny_matrix_v<DataType>) {
+ cell_array.fill(zero);
+ } else {
+ cell_array.fill(0);
+ }
+
+ parallel_for(
+ zone_cell_list.size(), PUGS_LAMBDA(const size_t i_cell) {
+ for (size_t i = 0; i < table.numberOfColumns(); ++i) {
+ cell_array[zone_cell_list[i_cell]][i] = table[i_cell][i];
+ }
+ });
+
+ return std::make_shared<DiscreteFunctionP0Vector<Dimension, DataType>>(p_mesh, cell_array);
+}
+
+template <size_t Dimension, typename DataType>
+std::shared_ptr<IDiscreteFunction>
+DiscreteFunctionVectorIntegrator::_integrateGlobally() const
{
+ Assert(m_zone_list.size() == 0, "invalid call when zones are defined");
+
std::shared_ptr mesh = std::dynamic_pointer_cast<const Mesh<Connectivity<Dimension>>>(m_mesh);
return std::make_shared<
@@ -17,6 +79,17 @@ DiscreteFunctionVectorIntegrator::_integrate() const
*m_quadrature_descriptor, *mesh));
}
+template <size_t Dimension, typename DataType>
+std::shared_ptr<IDiscreteFunction>
+DiscreteFunctionVectorIntegrator::_integrate() const
+{
+ if (m_zone_list.size() == 0) {
+ return this->_integrateGlobally<Dimension, DataType>();
+ } else {
+ return this->_integrateOnZoneList<Dimension, DataType>();
+ }
+}
+
template <size_t Dimension>
std::shared_ptr<IDiscreteFunction>
DiscreteFunctionVectorIntegrator::_integrate() const
diff --git a/src/scheme/DiscreteFunctionVectorIntegrator.hpp b/src/scheme/DiscreteFunctionVectorIntegrator.hpp
index 01e086e7efd9c2334d07dfdb3cf8dce1f33e49f1..55f5fe3572cbb237e28b9bf86ff7bda19db20a98 100644
--- a/src/scheme/DiscreteFunctionVectorIntegrator.hpp
+++ b/src/scheme/DiscreteFunctionVectorIntegrator.hpp
@@ -4,6 +4,7 @@
#include <analysis/IQuadratureDescriptor.hpp>
#include <language/utils/FunctionSymbolId.hpp>
#include <mesh/IMesh.hpp>
+#include <mesh/IZoneDescriptor.hpp>
#include <scheme/IDiscreteFunction.hpp>
#include <scheme/IDiscreteFunctionDescriptor.hpp>
@@ -14,10 +15,17 @@ class DiscreteFunctionVectorIntegrator
{
private:
std::shared_ptr<const IMesh> m_mesh;
+ const std::vector<std::shared_ptr<const IZoneDescriptor>> m_zone_list;
std::shared_ptr<const IQuadratureDescriptor> m_quadrature_descriptor;
std::shared_ptr<const IDiscreteFunctionDescriptor> m_discrete_function_descriptor;
const std::vector<FunctionSymbolId> m_function_id_list;
+ template <size_t Dimension, typename DataType>
+ std::shared_ptr<IDiscreteFunction> _integrateOnZoneList() const;
+
+ template <size_t Dimension, typename DataType>
+ std::shared_ptr<IDiscreteFunction> _integrateGlobally() const;
+
template <size_t Dimension, typename DataType>
std::shared_ptr<IDiscreteFunction> _integrate() const;
@@ -38,6 +46,19 @@ class DiscreteFunctionVectorIntegrator
m_function_id_list{function_id_list}
{}
+ DiscreteFunctionVectorIntegrator(
+ const std::shared_ptr<const IMesh>& mesh,
+ const std::vector<std::shared_ptr<const IZoneDescriptor>>& zone_list,
+ const std::shared_ptr<const IQuadratureDescriptor>& quadrature_descriptor,
+ const std::shared_ptr<const IDiscreteFunctionDescriptor>& discrete_function_descriptor,
+ const std::vector<FunctionSymbolId>& function_id_list)
+ : m_mesh{mesh},
+ m_zone_list{zone_list},
+ m_quadrature_descriptor(quadrature_descriptor),
+ m_discrete_function_descriptor{discrete_function_descriptor},
+ m_function_id_list{function_id_list}
+ {}
+
DiscreteFunctionVectorIntegrator(const DiscreteFunctionVectorIntegrator&) = delete;
DiscreteFunctionVectorIntegrator(DiscreteFunctionVectorIntegrator&&) = delete;