10 using namespace utils;
16 std::vector<json> flatten_ids(
const json &p_j_boundary_tmp)
20 std::vector<json> j_boundary;
22 for (
size_t i = 0; i < j_boundary_tmp.size(); ++i)
24 const auto &
tmp = j_boundary_tmp[i];
28 j_boundary.push_back(tmp);
31 if (!
tmp.contains(
"id"))
34 if (tmp[
"id"].is_array())
36 for (
size_t j = 0; j <
tmp[
"id"].size(); ++j)
39 newj[
"id"] =
tmp[
"id"][j].get<
int>();
40 j_boundary.push_back(newj);
44 j_boundary.push_back(tmp);
55 double z = pts.size() == 2 ? 0 : pts(2);
75 assert(pts.size() == 2 || pts.size() == 3);
76 double x = pts(0),
y = pts(1),
z = pts.size() == 3 ? pts(2) : 0.0;
94 for (
int i = 0; i < 3; ++i)
99 for (
int i = 0; i < 3; ++i)
103 v.set_unit_type(units.
length());
106 v.set_unit_type(units.
force());
115 v.second.set_unit_type(units.
pressure());
118 for (
int i = 0; i < 3; ++i)
119 v.second[i].set_unit_type(units.
length());
122 for (
int i = 0; i < 3; ++i)
123 v.second[i].set_unit_type(units.
velocity());
126 for (
int i = 0; i < 3; ++i)
130 v.second.set_unit_type(units.
length());
133 v.second.set_unit_type(units.
force());
139 val.resize(pts.rows(), pts.cols());
147 const bool planar = pts.cols() == 2;
148 for (
int i = 0; i < pts.rows(); ++i)
150 for (
int j = 0; j < pts.cols(); ++j)
152 double x = pts(i, 0),
y = pts(i, 1),
z = planar ? 0 : pts(i, 2);
184 val = Eigen::MatrixXd::Zero(pts.rows(), mesh.
dimension());
186 for (
long i = 0; i < pts.rows(); ++i)
191 for (
int d = 0; d <
val.cols(); ++d)
203 for (
int d = 0; d <
val.cols(); ++d)
217 val = Eigen::MatrixXd::Zero(pts.rows(), mesh.
dimension());
219 for (
long i = 0; i < pts.rows(); ++i)
227 for (
int d = 0; d <
val.cols(); ++d)
229 val(i, d) =
forces_[b].eval(pts.row(i), d, t);
240 for (
int d = 0; d <
val.cols(); ++d)
252 val = Eigen::MatrixXd::Zero(pts.rows(), 1);
254 for (
long i = 0; i < pts.rows(); ++i)
279 const bool planar = pts.cols() == 2;
280 val.resize(pts.rows(), pts.cols());
282 for (
int i = 0; i < pts.rows(); ++i)
284 for (
int j = 0; j < pts.cols(); ++j)
286 double x = pts(i, 0),
y = pts(i, 1),
z = planar ? 0 : pts(i, 2);
294 const int size = pts.cols();
295 val.resize(pts.rows(), pts.cols() * size);
299 const bool planar = size == 2;
300 for (
int i = 0; i < pts.rows(); ++i)
302 for (
int j = 0; j < pts.cols() * size; ++j)
304 double x = pts(i, 0),
y = pts(i, 1),
z = planar ? 0 : pts(i, 2);
315 for (
size_t k = 0; k <
val.size(); ++k)
321 if (!isx || !isy || !isz)
338 throw "Invalid boundary id";
344 for (
size_t k = 0; k <
val.size(); ++k)
349 if (!isx || !isy || !isz)
357 for (
size_t k = 0; k <
val.size(); ++k)
360 forces_.back().interpolation.push_back(interp);
376 throw "Invalid boundary id";
379 forces_[index].interpolation.clear();
380 forces_[index].interpolation.push_back(interp);
381 for (
size_t k = 0; k <
val.size(); ++k)
406 throw "Invalid boundary id";
412 void GenericTensorProblem::add_dirichlet_boundary(
const int id,
const std::function<Eigen::MatrixXd(
double x,
double y,
double z,
double t)> &func,
const bool isx,
const bool isy,
const bool isz,
const std::shared_ptr<Interpolation> &interp)
422 if (!isx || !isy || !isz)
439 throw "Invalid boundary id";
448 if (!isx || !isy || !isz)
456 forces_.back().interpolation.push_back(interp);
457 for (
size_t k = 0; k <
forces_.back().value.size(); ++k)
458 forces_.back().value[k].init(func, k);
474 throw "Invalid boundary id";
477 forces_[index].interpolation.clear();
478 forces_[index].interpolation.push_back(interp);
479 for (
size_t k = 0; k <
forces_.back().value.size(); ++k)
504 throw "Invalid boundary id";
514 throw "Val must be an array";
518 if (!interpolation.empty())
521 for (
size_t k = 0; k <
val.size(); ++k)
526 if (!isx || !isy || !isz)
533 throw "Val must be an array";
538 if (!interpolation.empty())
541 for (
size_t k = 0; k <
val.size(); ++k)
550 if (interpolation.empty())
551 pressures_.back().interpolation = std::make_shared<NoInterpolation>();
560 throw "Val must be an array";
572 throw "Invalid boundary id";
576 if (!interpolation.empty())
579 for (
size_t k = 0; k <
val.size(); ++k)
584 if (!isx || !isy || !isz)
591 throw "Val must be an array";
604 throw "Invalid boundary id";
606 forces_[index].interpolation.clear();
607 if (!interpolation.empty())
610 for (
size_t k = 0; k <
val.size(); ++k)
627 throw "Invalid boundary id";
630 if (interpolation.empty())
631 pressures_[index].interpolation = std::make_shared<NoInterpolation>();
654 for (
int d = 0; d <
val.cols(); ++d)
656 val(d) = tmp.eval(pt, d, t);
665 for (
int d = 0; d <
val.cols(); ++d)
667 val(d) = it->second.eval(pt, d, t);
675 for (
int i = 0; i < n_dirichlet.rows(); ++i)
677 if (n_dirichlet(i, 0) == node_id)
679 for (
int d = 0; d <
val.cols(); ++d)
681 val(d) = n_dirichlet(i, d + 1);
705 for (
int i = 0; i < n_dirichlet.rows(); ++i)
707 if (n_dirichlet(i, 0) == n_id)
735 return it->second.dirichlet_dimension(dim);
740 for (
int i = 0; i < n_dirichlet.rows(); ++i)
742 if (n_dirichlet(i, 0) == n_id)
744 return !std::isnan(n_dirichlet(i, dim + 1));
757 for (
int n = 0; n < n_dirichlet.rows(); ++n)
759 const int node_id = in_node_to_node[n_dirichlet(n, 0)];
760 n_dirichlet(n, 0) = node_id;
774 auto rr = params[
"rhs"];
777 for (
size_t k = 0; k < rr.size(); ++k)
782 logger().warn(
"Invalid problem rhs: should be an array.");
789 auto ex = params[
"reference"][
"solution"];
793 for (
size_t k = 0; k < ex.size(); ++k)
804 auto ex = params[
"reference"][
"gradient"];
808 for (
size_t k = 0; k < ex.size(); ++k)
821 std::vector<json> j_boundary = flatten_ids(params[
"dirichlet_boundary"]);
828 if (j_boundary[i - offset].is_string())
830 const std::string path =
resolve_path(j_boundary[i - offset], params[
"root_path"]);
831 if (!std::filesystem::is_regular_file(path))
843 if (j_boundary[i - offset][
"id"] ==
"all")
857 auto ff = j_boundary[i - offset][
"value"];
860 for (
size_t k = 0; k < ff.size(); ++k)
863 if (j_boundary[i - offset].contains(
"time_reference") && j_boundary[i - offset][
"time_reference"].size() > 0)
864 displacements_[i].value[k].set_t(j_boundary[i - offset][
"time_reference"]);
878 if (j_boundary[i - offset].contains(
"dimension"))
881 auto &tmp = j_boundary[i - offset][
"dimension"];
882 assert(tmp.is_array());
883 for (
size_t k = 0; k < tmp.size(); ++k)
890 if (j_boundary[i - offset][
"interpolation"].is_array())
892 for (
int ii = 0; ii < j_boundary[i - offset][
"interpolation"].size(); ++ii)
907 auto j_boundary_tmp = params[
"neumann_boundary"];
908 std::vector<json> j_boundary = flatten_ids(j_boundary_tmp);
911 forces_.resize(offset + j_boundary.size());
917 auto ff = j_boundary[i - offset][
"value"];
918 assert(ff.is_array());
920 for (
size_t k = 0; k < ff.size(); ++k)
923 if (j_boundary[i - offset][
"interpolation"].is_array())
925 for (
int ii = 0; ii < j_boundary[i - offset][
"interpolation"].size(); ++ii)
939 auto j_boundary_tmp = params[
"normal_aligned_neumann_boundary"];
940 std::vector<json> j_boundary = flatten_ids(j_boundary_tmp);
949 auto ff = j_boundary[i - offset][
"value"];
952 if (j_boundary[i - offset].contains(
"interpolation"))
964 auto j_boundary_tmp = params[
"pressure_boundary"];
965 std::vector<json> j_boundary = flatten_ids(j_boundary_tmp);
968 pressures_.resize(offset + j_boundary.size());
974 auto ff = j_boundary[i - offset][
"value"];
977 if (j_boundary[i - offset].contains(
"interpolation"))
980 pressures_[i].interpolation = std::make_shared<NoInterpolation>();
988 auto j_boundary_tmp = params[
"pressure_cavity"];
989 std::vector<json> j_boundary = flatten_ids(j_boundary_tmp);
995 int boundary_id = j_boundary[i - offset][
"id"];
1002 auto ff = j_boundary[i - offset][
"value"];
1010 auto rr = params[
"solution"];
1012 assert(rr.is_array());
1014 for (
size_t k = 0; k < rr.size(); ++k)
1017 const auto v = rr[k][
"value"];
1018 for (
size_t d = 0; d < v.size(); ++d)
1025 auto rr = params[
"velocity"];
1027 assert(rr.is_array());
1029 for (
size_t k = 0; k < rr.size(); ++k)
1032 const auto v = rr[k][
"value"];
1033 for (
size_t d = 0; d < v.size(); ++d)
1040 auto rr = params[
"acceleration"];
1042 assert(rr.is_array());
1044 for (
size_t k = 0; k < rr.size(); ++k)
1047 const auto v = rr[k][
"value"];
1048 for (
size_t d = 0; d < v.size(); ++d)
1056 val.resize(pts.rows(), pts.cols());
1063 const bool planar = pts.cols() == 2;
1064 for (
int i = 0; i < pts.rows(); ++i)
1078 val.row(i).setZero();
1082 for (
int j = 0; j < pts.cols(); ++j)
1089 val.resize(pts.rows(), pts.cols());
1096 const bool planar = pts.cols() == 2;
1097 for (
int i = 0; i < pts.rows(); ++i)
1111 val.row(i).setZero();
1115 for (
int j = 0; j < pts.cols(); ++j)
1122 val.resize(pts.rows(), pts.cols());
1129 const bool planar = pts.cols() == 2;
1130 for (
int i = 0; i < pts.rows(); ++i)
1144 val.row(i).setZero();
1148 for (
int j = 0; j < pts.cols(); ++j)
1174 for (
int i = 0; i <
rhs_.size(); ++i)
1176 for (
int i = 0; i <
exact_.size(); ++i)
1184 :
Problem(name), is_all_(false)
1193 v.set_unit_type(
"");
1196 v.set_unit_type(
"");
1199 v.second.set_unit_type(
"");
1204 for (
int i = 0; i < 3; ++i)
1210 val.resize(pts.rows(), 1);
1216 const bool planar = pts.cols() == 2;
1217 for (
int i = 0; i < pts.rows(); ++i)
1219 double x = pts(i, 0),
y = pts(i, 1),
z = planar ? 0 : pts(i, 2);
1226 val = Eigen::MatrixXd::Zero(pts.rows(), 1);
1228 for (
long i = 0; i < pts.rows(); ++i)
1252 val = Eigen::MatrixXd::Zero(pts.rows(), 1);
1254 for (
long i = 0; i < pts.rows(); ++i)
1262 double x = pts(i, 0),
y = pts(i, 1),
z = pts.cols() == 2 ? 0 : pts(i, 2);
1272 val.resize(pts.rows(), 1);
1279 const bool planar = pts.cols() == 2;
1280 for (
int i = 0; i < pts.rows(); ++i)
1305 const bool planar = pts.cols() == 2;
1306 val.resize(pts.rows(), 1);
1308 for (
int i = 0; i < pts.rows(); ++i)
1310 double x = pts(i, 0),
y = pts(i, 1),
z = pts.cols() == 2 ? 0 : pts(i, 2);
1317 val.resize(pts.rows(), pts.cols());
1321 const bool planar = pts.cols() == 2;
1322 for (
int i = 0; i < pts.rows(); ++i)
1324 for (
int j = 0; j < pts.cols(); ++j)
1326 double x = pts(i, 0),
y = pts(i, 1),
z = pts.cols() == 2 ? 0 : pts(i, 2);
1333 val = Eigen::MatrixXd::Zero(1, 1);
1341 val(0) = tmp.eval(pt, t);
1348 val(0) = it->second.eval(pt, t);
1354 for (
int i = 0; i < n_dirichlet.rows(); ++i)
1356 if (n_dirichlet(i, 0) == node_id)
1358 val(0) = n_dirichlet(i, 1);
1380 for (
int i = 0; i < n_dirichlet.rows(); ++i)
1382 if (n_dirichlet(i, 0) == n_id)
1409 for (
int n = 0; n < n_dirichlet.rows(); ++n)
1411 const int node_id = in_node_to_node[n_dirichlet(n, 0)];
1412 n_dirichlet(n, 0) = node_id;
1431 has_exact_ = !params[
"reference"][
"solution"].empty();
1432 exact_.
init(params[
"reference"][
"solution"]);
1437 auto ex = params[
"reference"][
"gradient"];
1441 for (
size_t k = 0; k < ex.size(); ++k)
1454 std::vector<json> j_boundary = flatten_ids(params[
"dirichlet_boundary"]);
1457 dirichlet_.resize(offset + j_boundary.size());
1461 if (j_boundary[i - offset].is_string())
1463 const std::string path =
resolve_path(j_boundary[i - offset], params[
"root_path"]);
1464 if (!std::filesystem::is_regular_file(path))
1467 Eigen::MatrixXd tmp;
1474 int current_id = -1;
1476 if (j_boundary[i - offset][
"id"] ==
"all")
1491 auto ff = j_boundary[i - offset][
"value"];
1495 if (j_boundary[i - offset][
"interpolation"].is_array())
1497 if (j_boundary[i - offset][
"interpolation"].size() == 0)
1498 dirichlet_[i].interpolation = std::make_shared<NoInterpolation>();
1499 else if (j_boundary[i - offset][
"interpolation"].size() == 1)
1515 auto j_boundary_tmp = params[
"neumann_boundary"];
1516 std::vector<json> j_boundary = flatten_ids(j_boundary_tmp);
1519 neumann_.resize(offset + j_boundary.size());
1525 auto ff = j_boundary[i - offset][
"value"];
1528 if (j_boundary[i - offset][
"interpolation"].is_array())
1530 if (j_boundary[i - offset][
"interpolation"].size() == 0)
1531 neumann_[i].interpolation = std::make_shared<NoInterpolation>();
1532 else if (j_boundary[i - offset][
"interpolation"].size() == 1)
1544 auto rr = params[
"solution"];
1546 assert(rr.is_array());
1548 for (
size_t k = 0; k < rr.size(); ++k)
1577 throw "Invalid boundary id";
1589 neumann_.back().interpolation = interp;
1605 throw "Invalid boundary id";
1609 neumann_[index].interpolation = interp;
1633 throw "Invalid boundary id";
1644 neumann_.back().interpolation = interp;
1660 throw "Invalid boundary id";
1663 neumann_[index].interpolation = interp;
1672 dirichlet_.back().interpolation = std::make_shared<NoInterpolation>();
1683 neumann_.back().interpolation = std::make_shared<NoInterpolation>();
1701 throw "Invalid boundary id";
1706 dirichlet_[index].interpolation = std::make_shared<NoInterpolation>();
1724 throw "Invalid boundary id";
1729 neumann_[index].interpolation = std::make_shared<NoInterpolation>();
std::string force() const
std::string pressure() const
std::string velocity() const
std::string acceleration() const
const std::string & length() const
virtual bool is_fluid() const
bool is_nodal_neumann_boundary(const int n_id, const int tag) override
void dirichlet_bc(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &uv, const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
bool has_nodal_dirichlet() override
void initial_solution(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &pts, Eigen::MatrixXd &val) const override
std::vector< std::pair< int, utils::ExpressionValue > > initial_solution_
void set_parameters(const json ¶ms) override
void rhs(const assembler::Assembler &assembler, const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
std::vector< ScalarBCValue > neumann_
void update_nodes(const Eigen::VectorXi &in_node_to_node) override
bool is_rhs_zero() const override
bool is_nodal_dirichlet_boundary(const int n_id, const int tag) override
void dirichlet_nodal_value(const mesh::Mesh &mesh, const int node_id, const RowVectorNd &pt, const double t, Eigen::MatrixXd &val) const override
utils::ExpressionValue exact_
std::array< utils::ExpressionValue, 3 > exact_grad_
void add_dirichlet_boundary(const int id, const double val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
bool has_nodal_neumann() override
void exact_grad(const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
void add_neumann_boundary(const int id, const double val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
std::map< int, ScalarBCValue > nodal_dirichlet_
void update_dirichlet_boundary(const int id, const double val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
void set_units(const assembler::Assembler &assembler, const Units &units) override
GenericScalarProblem(const std::string &name)
utils::ExpressionValue rhs_
std::map< int, ScalarBCValue > nodal_neumann_
void neumann_nodal_value(const mesh::Mesh &mesh, const int node_id, const RowVectorNd &pt, const Eigen::MatrixXd &normal, const double t, Eigen::MatrixXd &val) const override
bool has_exact_sol() const override
void neumann_bc(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &uv, const Eigen::MatrixXd &pts, const Eigen::MatrixXd &normals, const double t, Eigen::MatrixXd &val) const override
std::vector< ScalarBCValue > dirichlet_
void exact(const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
void update_neumann_boundary(const int id, const double val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
std::vector< Eigen::MatrixXd > nodal_dirichlet_mat_
std::map< int, TensorBCValue > nodal_neumann_
std::unordered_map< int, ScalarBCValue > cavity_pressures_
void neumann_bc(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &uv, const Eigen::MatrixXd &pts, const Eigen::MatrixXd &normals, const double t, Eigen::MatrixXd &val) const override
std::vector< std::pair< int, std::array< utils::ExpressionValue, 3 > > > initial_position_
void add_dirichlet_boundary(const int id, const Eigen::RowVector3d &val, const bool isx, const bool isy, const bool isz, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
void set_parameters(const json ¶ms) override
bool all_dimensions_dirichlet() const override
void update_nodes(const Eigen::VectorXi &in_node_to_node) override
GenericTensorProblem(const std::string &name)
void initial_velocity(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &pts, Eigen::MatrixXd &val) const override
bool all_dimensions_dirichlet_
std::vector< std::pair< int, std::array< utils::ExpressionValue, 3 > > > initial_acceleration_
void set_rhs(double x, double y, double z)
void update_pressure_boundary(const int id, const double val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
std::vector< std::pair< int, std::array< utils::ExpressionValue, 3 > > > initial_velocity_
std::vector< TensorBCValue > displacements_
std::vector< ScalarBCValue > normal_aligned_forces_
bool is_nodal_dirichlet_boundary(const int n_id, const int tag) override
bool is_rhs_zero() const override
bool has_nodal_neumann() override
void update_dirichlet_boundary(const int id, const Eigen::RowVector3d &val, const bool isx, const bool isy, const bool isz, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
void neumann_nodal_value(const mesh::Mesh &mesh, const int node_id, const RowVectorNd &pt, const Eigen::MatrixXd &normal, const double t, Eigen::MatrixXd &val) const override
void initial_acceleration(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &pts, Eigen::MatrixXd &val) const override
void update_neumann_boundary(const int id, const Eigen::RowVector3d &val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
void add_pressure_boundary(const int id, const double val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
bool has_nodal_dirichlet() override
void pressure_bc(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &uv, const Eigen::MatrixXd &pts, const Eigen::MatrixXd &normals, const double t, Eigen::MatrixXd &val) const override
bool is_dimension_dirichet(const int tag, const int dim) const override
void exact(const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
std::map< int, TensorBCValue > nodal_dirichlet_
std::array< utils::ExpressionValue, 3 > rhs_
std::array< utils::ExpressionValue, 9 > exact_grad_
void add_neumann_boundary(const int id, const Eigen::RowVector3d &val, const std::shared_ptr< utils::Interpolation > &interp=std::make_shared< utils::NoInterpolation >())
std::vector< TensorBCValue > forces_
double pressure_cavity_bc(const int boundary_id, const double t) const override
void dirichlet_nodal_value(const mesh::Mesh &mesh, const int node_id, const RowVectorNd &pt, const double t, Eigen::MatrixXd &val) const override
void dirichlet_bc(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &uv, const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
std::vector< ScalarBCValue > pressures_
bool is_nodal_neumann_boundary(const int n_id, const int tag) override
bool has_exact_sol() const override
void exact_grad(const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
bool is_nodal_dimension_dirichlet(const int n_id, const int tag, const int dim) const override
void set_units(const assembler::Assembler &assembler, const Units &units) override
void initial_solution(const mesh::Mesh &mesh, const Eigen::MatrixXi &global_ids, const Eigen::MatrixXd &pts, Eigen::MatrixXd &val) const override
void rhs(const assembler::Assembler &assembler, const Eigen::MatrixXd &pts, const double t, Eigen::MatrixXd &val) const override
std::array< utils::ExpressionValue, 3 > exact_
std::vector< Eigen::MatrixXd > nodal_dirichlet_mat_
std::vector< int > normal_aligned_neumann_boundary_ids_
std::vector< int > pressure_boundary_ids_
std::vector< int > boundary_ids_
virtual void init(const mesh::Mesh &mesh)
std::vector< int > neumann_boundary_ids_
std::vector< int > pressure_cavity_ids_
Abstract mesh class to capture 2d/3d conforming and non-conforming meshes.
virtual int get_body_id(const int primitive) const
Get the volume selection of an element (cell in 3d, face in 2d)
virtual int get_boundary_id(const int primitive) const
Get the boundary selection of an element (face in 3d, edge in 2d)
int dimension() const
utily for dimension
virtual int get_node_id(const int node_id) const
Get the boundary selection of a node.
void set_unit_type(const std::string &unit_type)
void init(const json &vals)
static std::shared_ptr< Interpolation > build(const json ¶ms)
bool read_matrix(const std::string &path, Eigen::Matrix< T, Eigen::Dynamic, Eigen::Dynamic > &mat)
Reads a matrix from a file. Determines the file format based on the path's extension.
std::string resolve_path(const std::string &path, const std::string &input_file_path, const bool only_if_exists=false)
std::vector< T > json_as_array(const json &j)
Return the value of a json object as an array.
bool is_param_valid(const json ¶ms, const std::string &key)
Determine if a key exists and is non-null in a json object.
spdlog::logger & logger()
Retrieves the current logger.
Eigen::Matrix< double, 1, Eigen::Dynamic, Eigen::RowMajor, 1, 3 > RowVectorNd
void log_and_throw_error(const std::string &msg)
std::shared_ptr< utils::Interpolation > interpolation
double eval(const RowVectorNd &pts, const double t) const
utils::ExpressionValue value
std::vector< std::shared_ptr< utils::Interpolation > > interpolation
std::array< utils::ExpressionValue, 3 > value
double eval(const RowVectorNd &pts, const int dim, const double t, const int el_id=-1) const