Skip to content
Snippets Groups Projects
Commit ef3c554c authored by PATELA Julie's avatar PATELA Julie
Browse files

Displace Heat5PointsAlgorithm in its own file

parent e8fda7b8
No related branches found
No related tags found
No related merge requests found
...@@ -2,6 +2,7 @@ ...@@ -2,6 +2,7 @@
add_library(PugsLanguageAlgorithms add_library(PugsLanguageAlgorithms
AcousticSolverAlgorithm.cpp AcousticSolverAlgorithm.cpp
Heat5PointsAlgorithm.cpp
) )
......
#include <language/algorithms/Heat5PointsAlgorithm.hpp>
#include <algebra/CRSMatrix.hpp>
#include <algebra/LocalRectangularMatrix.hpp>
#include <algebra/PCG.hpp>
#include <algebra/SparseMatrixDescriptor.hpp>
#include <algebra/TinyVector.hpp>
#include <algebra/Vector.hpp>
#include <language/utils/InterpolateItemValue.hpp>
#include <mesh/Connectivity.hpp>
#include <mesh/ConnectivityToDiamondDualConnectivityDataMapper.hpp>
#include <mesh/DiamondDualConnectivityManager.hpp>
#include <mesh/DiamondDualMeshManager.hpp>
#include <mesh/Mesh.hpp>
#include <mesh/MeshData.hpp>
#include <mesh/MeshDataManager.hpp>
#include <output/VTKWriter.hpp>
template <size_t Dimension>
Heat5PointsAlgorithm<Dimension>::Heat5PointsAlgorithm(
std::shared_ptr<const IMesh> i_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list,
const FunctionSymbolId& T_id,
const FunctionSymbolId& kappa_id)
{
using ConnectivityType = Connectivity<Dimension>;
using MeshType = Mesh<ConnectivityType>;
using MeshDataType = MeshData<Dimension>;
std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh);
std::cout << "number of bc descr = " << bc_descriptor_list.size() << '\n';
if constexpr (Dimension > 1) {
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*mesh);
CellValue<double> Tj =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(T_id, mesh_data.xj());
NodeValue<double> Tr(mesh->connectivity());
const NodeValue<const TinyVector<Dimension>>& xr = mesh->xr();
const CellValue<const TinyVector<Dimension>>& xj = mesh_data.xj();
const auto& node_to_cell_matrix = mesh->connectivity().nodeToCellMatrix();
CellValuePerNode<double> w_rj{mesh->connectivity()};
for (NodeId i_node = 0; i_node < mesh->numberOfNodes(); ++i_node) {
Vector<double> b{Dimension + 1};
b[0] = 1;
for (size_t i = 1; i < Dimension + 1; i++) {
b[i] = xr[i_node][i - 1];
}
const auto& node_to_cell = node_to_cell_matrix[i_node];
LocalRectangularMatrix<double> A{Dimension + 1, node_to_cell.size()};
for (size_t j = 0; j < node_to_cell.size(); j++) {
A(0, j) = 1;
}
for (size_t i = 1; i < Dimension + 1; i++) {
for (size_t j = 0; j < node_to_cell.size(); j++) {
const CellId J = node_to_cell[j];
A(i, j) = xj[J][i - 1];
}
}
Vector<double> x{node_to_cell.size()};
x = 0;
LocalRectangularMatrix B = transpose(A) * A;
Vector y = transpose(A) * b;
PCG{y, B, B, x, 10, 1e-12};
Tr[i_node] = 0;
for (size_t j = 0; j < node_to_cell.size(); j++) {
Tr[i_node] += x[j] * Tj[node_to_cell[j]];
w_rj(i_node, j) = x[j];
}
}
{
VTKWriter vtk_writer("T_" + std::to_string(Dimension), 0.01);
vtk_writer.write(mesh,
{NamedItemValue{"Tr", Tr}, NamedItemValue{"temperature", Tj},
NamedItemValue{"coords", mesh->xr()}, NamedItemValue{"xj", mesh_data.xj()},
NamedItemValue{"cell_owner", mesh->connectivity().cellOwner()},
NamedItemValue{"node_owner", mesh->connectivity().nodeOwner()}},
0, true); // forces last output
}
{
std::shared_ptr diamond_mesh = DiamondDualMeshManager::instance().getDiamondDualMesh(mesh);
MeshDataType& diamond_mesh_data = MeshDataManager::instance().getMeshData(*diamond_mesh);
std::shared_ptr mapper =
DiamondDualConnectivityManager::instance().getConnectivityToDiamondDualConnectivityDataMapper(
mesh->connectivity());
NodeValue<double> Trd{diamond_mesh->connectivity()};
mapper->toDualNode(Tr, Tj, Trd);
CellValue<double> kappaj =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(kappa_id,
mesh_data.xj());
CellValue<double> dual_kappaj =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(kappa_id,
diamond_mesh_data
.xj());
VTKWriter vtk_writer("D_" + std::to_string(Dimension), 0.01);
vtk_writer.write(diamond_mesh,
{NamedItemValue{"kappaj", dual_kappaj}, NamedItemValue{"coords", diamond_mesh->xr()},
NamedItemValue{"Vj", diamond_mesh_data.Vj()}, NamedItemValue{"xj", diamond_mesh_data.xj()}},
0, true); // forces last output
const CellValue<const double> dual_Vj = diamond_mesh_data.Vj();
const auto& face_to_node_matrix = mesh->connectivity().faceToNodeMatrix();
const FaceValue<const double> mes_l = [=] {
FaceValue<double> compute_mes_l{mesh->connectivity()};
const NodeValue<const TinyVector<Dimension>>& xr = mesh->xr();
if constexpr (Dimension == 1) {
compute_mes_l.fill(1);
} else if constexpr (Dimension == 2) {
parallel_for(
mesh->numberOfFaces(), PUGS_LAMBDA(FaceId face_id) {
const auto& face_to_node = face_to_node_matrix[face_id];
const NodeId node_id1 = face_to_node[0];
const NodeId node_id2 = face_to_node[1];
const TinyVector<Dimension, double> r = xr[node_id1] - xr[node_id2];
compute_mes_l[face_id] = l2Norm(r);
});
} else {
throw NotImplementedError("Not implemented in 3D");
}
return compute_mes_l;
}();
const CellValue<const double> dual_mes_l_j = [=] {
CellValue<double> compute_mes_j{diamond_mesh->connectivity()};
mapper->toDualCell(mes_l, compute_mes_j);
return compute_mes_j;
}();
FaceValue<const double> alpha_l = [&] {
CellValue<double> alpha_j{diamond_mesh->connectivity()};
parallel_for(
diamond_mesh->numberOfCells(), PUGS_LAMBDA(CellId diamond_cell_id) {
alpha_j[diamond_cell_id] =
dual_mes_l_j[diamond_cell_id] * dual_kappaj[diamond_cell_id] / dual_Vj[diamond_cell_id];
});
FaceValue<double> computed_alpha_l{mesh->connectivity()};
mapper->fromDualCell(alpha_j, computed_alpha_l);
VTKWriter vtk_writer("DD_" + std::to_string(Dimension), 0.01);
vtk_writer.write(diamond_mesh,
{NamedItemValue{"alphaj", alpha_j}, NamedItemValue{"xj", diamond_mesh_data.xj()},
NamedItemValue{"Vl", diamond_mesh_data.Vj()}, NamedItemValue{"|l|", dual_mes_l_j}},
0,
true); // forces last output
return computed_alpha_l;
}();
SparseMatrixDescriptor S{mesh->numberOfCells()};
const auto& face_to_cell_matrix = mesh->connectivity().faceToCellMatrix();
for (FaceId face_id = 0; face_id < mesh->numberOfFaces(); ++face_id) {
const auto& primal_face_to_cell = face_to_cell_matrix[face_id];
const double beta_l = 0.5 * alpha_l[face_id] * mes_l[face_id];
for (size_t i_cell = 0; i_cell < primal_face_to_cell.size(); ++i_cell) {
const CellId cell_id1 = primal_face_to_cell[i_cell];
for (size_t j_cell = 0; j_cell < primal_face_to_cell.size(); ++j_cell) {
const CellId cell_id2 = primal_face_to_cell[j_cell];
if (i_cell == j_cell) {
S(cell_id1, cell_id2) -= beta_l;
} else {
S(cell_id1, cell_id2) += beta_l;
}
}
}
}
const CellValue<const double> primal_Vj = mesh_data.Vj();
CRSMatrix A{S};
Vector<double> b{mesh->numberOfCells()};
double pi = 3.14159265;
for (CellId i_cell = 0; i_cell < mesh->numberOfCells(); ++i_cell) {
if constexpr (Dimension == 2) {
b[i_cell] =
-2 * kappaj[i_cell] * pi * pi * sin(pi * xj[i_cell][0]) * sin(pi * xj[i_cell][1]) * primal_Vj[i_cell];
} else if constexpr (Dimension == 3) {
b[i_cell] = -3 * kappaj[i_cell] * pi * pi * sin(pi * xj[i_cell][0]) * sin(pi * xj[i_cell][1]) *
sin(pi * xj[i_cell][2]) * primal_Vj[i_cell];
}
}
Vector<double> T{mesh->numberOfCells()};
T = 0;
PCG{b, A, A, T, 100, 1e-12};
CellValue<double> Temperature{mesh->connectivity()};
parallel_for(
mesh->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { Temperature[cell_id] = T[cell_id]; });
Vector<double> error{mesh->numberOfCells()};
parallel_for(
mesh->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { error[cell_id] = Temperature[cell_id] - Tj[cell_id]; });
std::cout << "||Error||_2 = " << std::sqrt((error, error)) << "\n";
{
VTKWriter vtk_writer("Temperature_" + std::to_string(Dimension), 0.01);
vtk_writer.write(mesh, {NamedItemValue{"T", Temperature}}, 0,
true); // forces last output
}
}
} else {
throw NotImplementedError("not done in 3d");
}
}
template Heat5PointsAlgorithm<1>::Heat5PointsAlgorithm(
std::shared_ptr<const IMesh>,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&,
const FunctionSymbolId&,
const FunctionSymbolId&);
template Heat5PointsAlgorithm<2>::Heat5PointsAlgorithm(
std::shared_ptr<const IMesh>,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&,
const FunctionSymbolId&,
const FunctionSymbolId&);
template Heat5PointsAlgorithm<3>::Heat5PointsAlgorithm(
std::shared_ptr<const IMesh>,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&,
const FunctionSymbolId&,
const FunctionSymbolId&);
#ifndef HEAT_5POINTS_ALGORITHM_HPP
#define HEAT_5POINTS_ALGORITHM_HPP
#include <language/utils/FunctionSymbolId.hpp>
#include <mesh/IMesh.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
#include <memory>
#include <vector>
template <size_t Dimension>
struct Heat5PointsAlgorithm
{
Heat5PointsAlgorithm(std::shared_ptr<const IMesh> i_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list,
const FunctionSymbolId& T_id,
const FunctionSymbolId& kappa_id);
};
#endif // HEAT_5POINTS_ALGORITHM_HPP
...@@ -7,10 +7,10 @@ ...@@ -7,10 +7,10 @@
#include <algebra/TinyVector.hpp> #include <algebra/TinyVector.hpp>
#include <algebra/Vector.hpp> #include <algebra/Vector.hpp>
#include <language/algorithms/AcousticSolverAlgorithm.hpp> #include <language/algorithms/AcousticSolverAlgorithm.hpp>
#include <language/algorithms/Heat5PointsAlgorithm.hpp>
#include <language/utils/BuiltinFunctionEmbedder.hpp> #include <language/utils/BuiltinFunctionEmbedder.hpp>
#include <language/utils/PugsFunctionAdapter.hpp> #include <language/utils/PugsFunctionAdapter.hpp>
#include <language/utils/TypeDescriptor.hpp> #include <language/utils/TypeDescriptor.hpp>
#include <memory>
#include <mesh/ConnectivityToDiamondDualConnectivityDataMapper.hpp> #include <mesh/ConnectivityToDiamondDualConnectivityDataMapper.hpp>
#include <mesh/DiamondDualConnectivityManager.hpp> #include <mesh/DiamondDualConnectivityManager.hpp>
#include <mesh/DiamondDualMeshBuilder.hpp> #include <mesh/DiamondDualMeshBuilder.hpp>
...@@ -25,6 +25,8 @@ ...@@ -25,6 +25,8 @@
#include <scheme/NamedBoundaryDescriptor.hpp> #include <scheme/NamedBoundaryDescriptor.hpp>
#include <scheme/NumberedBoundaryDescriptor.hpp> #include <scheme/NumberedBoundaryDescriptor.hpp>
#include <memory>
#include <scheme/SymmetryBoundaryConditionDescriptor.hpp> #include <scheme/SymmetryBoundaryConditionDescriptor.hpp>
/////////// TEMPORARY /////////// TEMPORARY
...@@ -250,224 +252,6 @@ struct HeatDiamondScheme ...@@ -250,224 +252,6 @@ struct HeatDiamondScheme
} }
}; };
template <size_t Dimension>
struct Laplace5points
{
using ConnectivityType = Connectivity<Dimension>;
using MeshType = Mesh<ConnectivityType>;
using MeshDataType = MeshData<Dimension>;
using UnknownsType = FiniteVolumesEulerUnknowns<MeshType>;
std::shared_ptr<const MeshType> m_mesh;
Laplace5points(std::shared_ptr<const IMesh> i_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list,
const FunctionSymbolId& T_id,
const FunctionSymbolId& kappa_id)
: m_mesh{std::dynamic_pointer_cast<const MeshType>(i_mesh)}
{
std::cout << "number of bc descr = " << bc_descriptor_list.size() << '\n';
if constexpr (Dimension > 1) {
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*m_mesh);
CellValue<double> Tj =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(T_id, mesh_data.xj());
NodeValue<double> Tr(m_mesh->connectivity());
const NodeValue<const TinyVector<Dimension>>& xr = m_mesh->xr();
const CellValue<const TinyVector<Dimension>>& xj = mesh_data.xj();
const auto& node_to_cell_matrix = m_mesh->connectivity().nodeToCellMatrix();
CellValuePerNode<double> w_rj{m_mesh->connectivity()};
for (NodeId i_node = 0; i_node < m_mesh->numberOfNodes(); ++i_node) {
Vector<double> b{Dimension + 1};
b[0] = 1;
for (size_t i = 1; i < Dimension + 1; i++) {
b[i] = xr[i_node][i - 1];
}
const auto& node_to_cell = node_to_cell_matrix[i_node];
LocalRectangularMatrix<double> A{Dimension + 1, node_to_cell.size()};
for (size_t j = 0; j < node_to_cell.size(); j++) {
A(0, j) = 1;
}
for (size_t i = 1; i < Dimension + 1; i++) {
for (size_t j = 0; j < node_to_cell.size(); j++) {
const CellId J = node_to_cell[j];
A(i, j) = xj[J][i - 1];
}
}
Vector<double> x{node_to_cell.size()};
x = 0;
LocalRectangularMatrix B = transpose(A) * A;
Vector y = transpose(A) * b;
PCG{y, B, B, x, 10, 1e-12};
Tr[i_node] = 0;
for (size_t j = 0; j < node_to_cell.size(); j++) {
Tr[i_node] += x[j] * Tj[node_to_cell[j]];
w_rj(i_node, j) = x[j];
}
}
{
VTKWriter vtk_writer("T_" + std::to_string(Dimension), 0.01);
vtk_writer.write(m_mesh,
{NamedItemValue{"Tr", Tr}, NamedItemValue{"temperature", Tj},
NamedItemValue{"coords", m_mesh->xr()}, NamedItemValue{"xj", mesh_data.xj()},
NamedItemValue{"cell_owner", m_mesh->connectivity().cellOwner()},
NamedItemValue{"node_owner", m_mesh->connectivity().nodeOwner()}},
0, true); // forces last output
}
{
std::shared_ptr diamond_mesh = DiamondDualMeshManager::instance().getDiamondDualMesh(m_mesh);
MeshDataType& diamond_mesh_data = MeshDataManager::instance().getMeshData(*diamond_mesh);
std::shared_ptr mapper =
DiamondDualConnectivityManager::instance().getConnectivityToDiamondDualConnectivityDataMapper(
m_mesh->connectivity());
NodeValue<double> Trd{diamond_mesh->connectivity()};
mapper->toDualNode(Tr, Tj, Trd);
CellValue<double> kappaj =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(kappa_id,
mesh_data.xj());
CellValue<double> dual_kappaj =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(kappa_id,
diamond_mesh_data
.xj());
VTKWriter vtk_writer("D_" + std::to_string(Dimension), 0.01);
vtk_writer.write(diamond_mesh,
{NamedItemValue{"kappaj", dual_kappaj}, NamedItemValue{"coords", diamond_mesh->xr()},
NamedItemValue{"Vj", diamond_mesh_data.Vj()}, NamedItemValue{"xj", diamond_mesh_data.xj()}},
0, true); // forces last output
const CellValue<const double> dual_Vj = diamond_mesh_data.Vj();
const auto& face_to_node_matrix = m_mesh->connectivity().faceToNodeMatrix();
const FaceValue<const double> mes_l = [=] {
FaceValue<double> compute_mes_l{m_mesh->connectivity()};
const NodeValue<const TinyVector<Dimension>>& xr = m_mesh->xr();
if constexpr (Dimension == 1) {
compute_mes_l.fill(1);
} else if constexpr (Dimension == 2) {
parallel_for(
m_mesh->numberOfFaces(), PUGS_LAMBDA(FaceId face_id) {
const auto& face_to_node = face_to_node_matrix[face_id];
const NodeId node_id1 = face_to_node[0];
const NodeId node_id2 = face_to_node[1];
const TinyVector<Dimension, double> r = xr[node_id1] - xr[node_id2];
compute_mes_l[face_id] = l2Norm(r);
});
} else {
throw NotImplementedError("Not implemented in 3D");
}
return compute_mes_l;
}();
const CellValue<const double> dual_mes_l_j = [=] {
CellValue<double> compute_mes_j{diamond_mesh->connectivity()};
mapper->toDualCell(mes_l, compute_mes_j);
return compute_mes_j;
}();
FaceValue<const double> alpha_l = [&] {
CellValue<double> alpha_j{diamond_mesh->connectivity()};
parallel_for(
diamond_mesh->numberOfCells(), PUGS_LAMBDA(CellId diamond_cell_id) {
alpha_j[diamond_cell_id] =
dual_mes_l_j[diamond_cell_id] * dual_kappaj[diamond_cell_id] / dual_Vj[diamond_cell_id];
});
FaceValue<double> computed_alpha_l{m_mesh->connectivity()};
mapper->fromDualCell(alpha_j, computed_alpha_l);
VTKWriter vtk_writer("DD_" + std::to_string(Dimension), 0.01);
vtk_writer.write(diamond_mesh,
{NamedItemValue{"alphaj", alpha_j}, NamedItemValue{"xj", diamond_mesh_data.xj()},
NamedItemValue{"Vl", diamond_mesh_data.Vj()}, NamedItemValue{"|l|", dual_mes_l_j}},
0,
true); // forces last output
return computed_alpha_l;
}();
SparseMatrixDescriptor S{m_mesh->numberOfCells()};
const auto& face_to_cell_matrix = m_mesh->connectivity().faceToCellMatrix();
for (FaceId face_id = 0; face_id < m_mesh->numberOfFaces(); ++face_id) {
const auto& primal_face_to_cell = face_to_cell_matrix[face_id];
const double beta_l = 0.5 * alpha_l[face_id] * mes_l[face_id];
for (size_t i_cell = 0; i_cell < primal_face_to_cell.size(); ++i_cell) {
const CellId cell_id1 = primal_face_to_cell[i_cell];
for (size_t j_cell = 0; j_cell < primal_face_to_cell.size(); ++j_cell) {
const CellId cell_id2 = primal_face_to_cell[j_cell];
if (i_cell == j_cell) {
S(cell_id1, cell_id2) -= beta_l;
} else {
S(cell_id1, cell_id2) += beta_l;
}
}
}
}
const CellValue<const double> primal_Vj = mesh_data.Vj();
CRSMatrix A{S};
Vector<double> b{m_mesh->numberOfCells()};
double pi = 3.14159265;
for (CellId i_cell = 0; i_cell < m_mesh->numberOfCells(); ++i_cell) {
if constexpr (Dimension == 2) {
b[i_cell] =
-2 * kappaj[i_cell] * pi * pi * sin(pi * xj[i_cell][0]) * sin(pi * xj[i_cell][1]) * primal_Vj[i_cell];
} else if constexpr (Dimension == 3) {
b[i_cell] = -3 * kappaj[i_cell] * pi * pi * sin(pi * xj[i_cell][0]) * sin(pi * xj[i_cell][1]) *
sin(pi * xj[i_cell][2]) * primal_Vj[i_cell];
}
}
Vector<double> T{m_mesh->numberOfCells()};
T = 0;
PCG{b, A, A, T, 100, 1e-12};
CellValue<double> Temperature{m_mesh->connectivity()};
parallel_for(
m_mesh->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { Temperature[cell_id] = T[cell_id]; });
Vector<double> error{m_mesh->numberOfCells()};
parallel_for(
m_mesh->numberOfCells(),
PUGS_LAMBDA(CellId cell_id) { error[cell_id] = Temperature[cell_id] - Tj[cell_id]; });
std::cout << "||Error||_2 = " << std::sqrt((error, error)) << "\n";
{
VTKWriter vtk_writer("Temperature_" + std::to_string(Dimension), 0.01);
vtk_writer.write(m_mesh, {NamedItemValue{"T", Temperature}}, 0,
true); // forces last output
}
}
} else {
throw NotImplementedError("not done in 3d");
}
}
};
SchemeModule::SchemeModule() SchemeModule::SchemeModule()
{ {
this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IBoundaryDescriptor>>); this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IBoundaryDescriptor>>);
...@@ -655,7 +439,7 @@ SchemeModule::SchemeModule() ...@@ -655,7 +439,7 @@ SchemeModule::SchemeModule()
)); ));
this->_addBuiltinFunction("laplace5points", this->_addBuiltinFunction("heat5points",
std::make_shared<BuiltinFunctionEmbedder< std::make_shared<BuiltinFunctionEmbedder<
void(std::shared_ptr<const IMesh>, void(std::shared_ptr<const IMesh>,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&, const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&,
...@@ -667,15 +451,15 @@ SchemeModule::SchemeModule() ...@@ -667,15 +451,15 @@ SchemeModule::SchemeModule()
const FunctionSymbolId& T_id, const FunctionSymbolId& kappa_id) -> void { const FunctionSymbolId& T_id, const FunctionSymbolId& kappa_id) -> void {
switch (p_mesh->dimension()) { switch (p_mesh->dimension()) {
case 1: { case 1: {
Laplace5points<1>{p_mesh, bc_descriptor_list, T_id, kappa_id}; Heat5PointsAlgorithm<1>{p_mesh, bc_descriptor_list, T_id, kappa_id};
break; break;
} }
case 2: { case 2: {
Laplace5points<2>{p_mesh, bc_descriptor_list, T_id, kappa_id}; Heat5PointsAlgorithm<2>{p_mesh, bc_descriptor_list, T_id, kappa_id};
break; break;
} }
case 3: { case 3: {
Laplace5points<3>{p_mesh, bc_descriptor_list, T_id, kappa_id}; Heat5PointsAlgorithm<3>{p_mesh, bc_descriptor_list, T_id, kappa_id};
break; break;
} }
default: { default: {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment