Wildmeshing Toolkit
Loading...
Searching...
No Matches
Swap32EnergyBeforeInvariantDouble.cpp
Go to the documentation of this file.
1#include "Swap32EnergyBeforeInvariantDouble.hpp"
2#include <wmtk/Mesh.hpp>
3#include <wmtk/function/utils/amips.hpp>
4#include <wmtk/utils/orient.hpp>
5
6namespace wmtk {
7Swap32EnergyBeforeInvariantDouble::Swap32EnergyBeforeInvariantDouble(
8 const Mesh& m,
9 const attribute::TypedAttributeHandle<double>& coordinate,
10 double eps)
11 : Invariant(m, true, false, false)
12 , m_coordinate_handle(coordinate)
13 , m_eps(eps)
14{}
15
16bool Swap32EnergyBeforeInvariantDouble::before(const simplex::Simplex& t) const
17{
18 constexpr static PrimitiveType PV = PrimitiveType::Vertex;
19 constexpr static PrimitiveType PE = PrimitiveType::Edge;
20 constexpr static PrimitiveType PF = PrimitiveType::Triangle;
21 constexpr static PrimitiveType PT = PrimitiveType::Tetrahedron;
22
23 auto accessor = mesh().create_const_accessor(m_coordinate_handle);
24
25 // get the coords of the vertices
26 // input edge end points
27 const Tuple v0 = t.tuple();
28 const Tuple v1 = mesh().switch_tuple(v0, PV);
29 // other 3 vertices
30 const Tuple v2 = mesh().switch_tuples(v0, {PE, PV});
31 const Tuple v3 = mesh().switch_tuples(v0, {PF, PE, PV});
32 const Tuple v4 = mesh().switch_tuples(v0, {PT, PF, PE, PV});
33
34
35 std::array<Eigen::Vector3<double>, 5> positions = {
36 {accessor.const_vector_attribute(v0),
37 accessor.const_vector_attribute(v1),
38 accessor.const_vector_attribute(v2),
39 accessor.const_vector_attribute(v3),
40 accessor.const_vector_attribute(v4)}};
41 std::array<Eigen::Vector3d, 5> positions_double = {
42 {accessor.const_vector_attribute(v0),
43 accessor.const_vector_attribute(v1),
44 accessor.const_vector_attribute(v2),
45 accessor.const_vector_attribute(v3),
46 accessor.const_vector_attribute(v4)}};
47
48 std::array<std::array<int, 4>, 3> old_tets = {{{{0, 1, 2, 3}}, {{0, 1, 3, 4}}, {{0, 1, 4, 2}}}};
49 std::array<std::array<int, 4>, 2> new_tets = {{{{2, 3, 4, 0}}, {{2, 4, 3, 1}}}};
50
51
52 double old_energy_max = std::numeric_limits<double>::lowest();
53 double new_energy_max = std::numeric_limits<double>::lowest();
54
55 for (int i = 0; i < 3; ++i) {
56 if (utils::wmtk_orient3d(
57 positions[old_tets[i][0]],
58 positions[old_tets[i][1]],
59 positions[old_tets[i][2]],
60 positions[old_tets[i][3]]) > 0) {
61 auto energy = wmtk::function::utils::Tet_AMIPS_energy({{
62 positions_double[old_tets[i][0]][0],
63 positions_double[old_tets[i][0]][1],
64 positions_double[old_tets[i][0]][2],
65 positions_double[old_tets[i][1]][0],
66 positions_double[old_tets[i][1]][1],
67 positions_double[old_tets[i][1]][2],
68 positions_double[old_tets[i][2]][0],
69 positions_double[old_tets[i][2]][1],
70 positions_double[old_tets[i][2]][2],
71 positions_double[old_tets[i][3]][0],
72 positions_double[old_tets[i][3]][1],
73 positions_double[old_tets[i][3]][2],
74 }});
75
76
77 if (energy > old_energy_max) old_energy_max = energy;
78 } else {
79 auto energy = wmtk::function::utils::Tet_AMIPS_energy({{
80 positions_double[old_tets[i][1]][0],
81 positions_double[old_tets[i][1]][1],
82 positions_double[old_tets[i][1]][2],
83 positions_double[old_tets[i][0]][0],
84 positions_double[old_tets[i][0]][1],
85 positions_double[old_tets[i][0]][2],
86 positions_double[old_tets[i][2]][0],
87 positions_double[old_tets[i][2]][1],
88 positions_double[old_tets[i][2]][2],
89 positions_double[old_tets[i][3]][0],
90 positions_double[old_tets[i][3]][1],
91 positions_double[old_tets[i][3]][2],
92 }});
93
94
95 if (energy > old_energy_max) old_energy_max = energy;
96 }
97 }
98
99 for (int i = 0; i < 2; ++i) {
100 if (utils::wmtk_orient3d(
101 positions[new_tets[i][0]],
102 positions[new_tets[i][1]],
103 positions[new_tets[i][2]],
104 positions[new_tets[i][3]]) > 0) {
105 auto energy = wmtk::function::utils::Tet_AMIPS_energy({{
106 positions_double[new_tets[i][0]][0],
107 positions_double[new_tets[i][0]][1],
108 positions_double[new_tets[i][0]][2],
109 positions_double[new_tets[i][1]][0],
110 positions_double[new_tets[i][1]][1],
111 positions_double[new_tets[i][1]][2],
112 positions_double[new_tets[i][2]][0],
113 positions_double[new_tets[i][2]][1],
114 positions_double[new_tets[i][2]][2],
115 positions_double[new_tets[i][3]][0],
116 positions_double[new_tets[i][3]][1],
117 positions_double[new_tets[i][3]][2],
118 }});
119
120
121 if (energy > new_energy_max) new_energy_max = energy;
122 } else {
123 auto energy = wmtk::function::utils::Tet_AMIPS_energy({{
124 positions_double[new_tets[i][1]][0],
125 positions_double[new_tets[i][1]][1],
126 positions_double[new_tets[i][1]][2],
127 positions_double[new_tets[i][0]][0],
128 positions_double[new_tets[i][0]][1],
129 positions_double[new_tets[i][0]][2],
130 positions_double[new_tets[i][2]][0],
131 positions_double[new_tets[i][2]][1],
132 positions_double[new_tets[i][2]][2],
133 positions_double[new_tets[i][3]][0],
134 positions_double[new_tets[i][3]][1],
135 positions_double[new_tets[i][3]][2],
136 }});
137
138
139 if (energy > new_energy_max) new_energy_max = energy;
140 }
141 }
142
143 return old_energy_max > new_energy_max * m_eps;
144}
145
146} // namespace wmtk
wmtk::Mesh::create_const_accessor
const attribute::Accessor< T, Mesh, D > create_const_accessor(const attribute::MeshAttributeHandle &handle) const
wmtk::Mesh::switch_tuples
Tuple switch_tuples(const Tuple &tuple, const ContainerType &op_sequence) const
Definition Mesh.hpp:953
wmtk::Mesh::switch_tuple
virtual Tuple switch_tuple(const Tuple &tuple, PrimitiveType type) const =0
switch the orientation of the Tuple of the given dimension
wmtk::Swap32EnergyBeforeInvariantDouble::before
bool before(const simplex::Simplex &t) const override
Definition Swap32EnergyBeforeInvariantDouble.cpp:16
wmtk::Swap32EnergyBeforeInvariantDouble::m_coordinate_handle
const attribute::TypedAttributeHandle< double > m_coordinate_handle
Definition Swap32EnergyBeforeInvariantDouble.hpp:19
wmtk::Swap32EnergyBeforeInvariantDouble::Swap32EnergyBeforeInvariantDouble
Swap32EnergyBeforeInvariantDouble(const Mesh &m, const attribute::TypedAttributeHandle< double > &coordinate, double eps=1.0001)
Definition Swap32EnergyBeforeInvariantDouble.cpp:7
wmtk::Tuple
The Tuple is the basic navigation tool in our mesh data structure.
Definition Tuple.hpp:19
wmtk::attribute::TypedAttributeHandle
Handle that represents attributes for some mesh.
Definition TypedAttributeHandle.hpp:28
wmtk::invariants::Invariant::mesh
const Mesh & mesh() const
Definition Invariant.cpp:35
wmtk::simplex::Simplex::tuple
const Tuple & tuple() const
Definition Simplex.hpp:53
PT
constexpr wmtk::PrimitiveType PT
Definition faces_single_dimension.cpp:6
PF
constexpr wmtk::PrimitiveType PF
Definition faces_single_dimension.cpp:5
wmtk::function::utils::Tet_AMIPS_energy
double Tet_AMIPS_energy(const std::array< double, 12 > &T)
Definition amips.cpp:380
wmtk::utils::wmtk_orient3d
int wmtk_orient3d(const Eigen::Ref< const Eigen::Vector3< Rational > > &p0, const Eigen::Ref< const Eigen::Vector3< Rational > > &p1, const Eigen::Ref< const Eigen::Vector3< Rational > > &p2, const Eigen::Ref< const Eigen::Vector3< Rational > > &p3)
Definition orient.cpp:75
wmtk::PE
constexpr PrimitiveType PE
Definition EnvelopeInvariant.cpp:25
wmtk::PV
constexpr PrimitiveType PV
Definition EnvelopeInvariant.cpp:24