Skip to content
Snippets Groups Projects
Commit abd6eb6f authored by Emmanuel Labourasse's avatar Emmanuel Labourasse
Browse files

add a constructor through vectors to TinyMatrix

parent 26a05415
No related branches found
No related tags found
1 merge request!167Improve fluxing based remapping
...@@ -192,7 +192,43 @@ template <> ...@@ -192,7 +192,43 @@ template <>
FaceValue<double> FaceValue<double>
FluxingAdvectionSolver<3>::_computeAlgebraicFluxingVolume() FluxingAdvectionSolver<3>::_computeAlgebraicFluxingVolume()
{ {
throw NotImplementedError("ViensViensViens"); const auto face_to_node_matrix = m_old_mesh->connectivity().faceToNodeMatrix();
FaceValue<double> algebraic_fluxing_volume(m_new_mesh->connectivity());
auto old_xr = m_old_mesh->xr();
auto new_xr = m_new_mesh->xr();
parallel_for(
m_new_mesh->numberOfFaces(), PUGS_LAMBDA(FaceId face_id) {
const auto& face_to_node = face_to_node_matrix[face_id];
if (face_to_node.size() == 4) {
const Rd& x0 = old_xr[face_to_node[0]];
const Rd& x1 = old_xr[face_to_node[1]];
const Rd& x2 = old_xr[face_to_node[2]];
const Rd& x3 = old_xr[face_to_node[3]];
const Rd& x4 = new_xr[face_to_node[1]];
const Rd& x5 = new_xr[face_to_node[0]];
const Rd& x6 = new_xr[face_to_node[1]];
const Rd& x7 = new_xr[face_to_node[0]];
const Rd& a1 = x6 - x1;
const Rd& a2 = x6 - x3;
const Rd& a3 = x6 - x4;
const Rd& b1 = x7 - x0;
const Rd& b2 = x5 - x0;
const Rd& b3 = x2 - x0;
TinyMatrix<3> M1(a1 + b1, a2, a3);
TinyMatrix<3> M2(b1, a2 + b2, a3);
TinyMatrix<3> M3(a1, b2, a3 + b3);
algebraic_fluxing_volume[face_id] = (det(M1) + det(M2) + det(M3)) / 12;
} else {
throw NotImplementedError("Not implemented for non quad faces");
}
});
return algebraic_fluxing_volume;
} }
template <size_t Dimension> template <size_t Dimension>
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment