@@ -107,7 +107,18 @@ int MMG2D_anatri(MMG5_pMesh mesh,MMG5_pSol met,int8_t typchk) {
107107 return 1 ;
108108}
109109
110- // Minimum angle function
110+ /**
111+ * \param mesh pointer to the mesh structure.
112+ * \param k triangle index.
113+ * \return the minimal angle of a triangle or 0 if the triangle size
114+ * is not conform to the tolerance sizes.
115+ *
116+ * Compute the minimal angle of triangle k. If the velocity of each
117+ * node is given, the extrapolated minimal angle of the triangle
118+ * is also calculated. If the edge lengths are too small of too
119+ * large, 0 is returned as we want this triangle to be remeshed.
120+ *
121+ */
111122double MMG2D_minangle (MMG5_pMesh mesh , int k ) {
112123
113124 MMG5_pTria pt = & mesh -> tria [k ];
@@ -129,11 +140,13 @@ double MMG2D_minangle(MMG5_pMesh mesh, int k) {
129140 pp2 = & mesh -> point [p2 ];
130141 pp3 = & mesh -> point [p3 ];
131142
132- double square_root_area = sqrt (0.5 * fabs ((pp2 -> c [0 ]- pp1 -> c [0 ])* (pp3 -> c [1 ]- pp1 -> c [1 ]) - (pp2 -> c [1 ]- pp1 -> c [1 ])* (pp3 -> c [0 ]- pp1 -> c [0 ])));
143+ double square_root_area = sqrt (0.5 * fabs ((pp2 -> c [0 ]- pp1 -> c [0 ])* (pp3 -> c [1 ]- pp1 -> c [1 ]) -
144+ (pp2 -> c [1 ]- pp1 -> c [1 ])* (pp3 -> c [0 ]- pp1 -> c [0 ])));
133145 double mean_velocity = 0. ;
134146
135147 // Velocity is only considered for already existing points
136- if (mesh -> velocity && p1 <= mesh -> mark && p2 <= mesh -> mark && p3 <= mesh -> mark && pp1 -> tmp == -1 && pp2 -> tmp == -1 && pp3 -> tmp == -1 )
148+ if (mesh -> velocity && p1 <= mesh -> mark && p2 <= mesh -> mark && p3 <= mesh -> mark
149+ && pp1 -> tmp == -1 && pp2 -> tmp == -1 && pp3 -> tmp == -1 )
137150 mean_velocity = ( sqrt (pow (mesh -> velocity [p1 - 1 ],2 ) + pow (mesh -> velocity [p1 - 1 + mesh -> mark ],2 )) +
138151 sqrt (pow (mesh -> velocity [p2 - 1 ],2 ) + pow (mesh -> velocity [p2 - 1 + mesh -> mark ],2 )) +
139152 sqrt (pow (mesh -> velocity [p3 - 1 ],2 ) + pow (mesh -> velocity [p3 - 1 + mesh -> mark ],2 )) )/3. ;
@@ -149,6 +162,7 @@ double MMG2D_minangle(MMG5_pMesh mesh, int k) {
149162 y3 = pp3 -> c [1 ];
150163 }
151164 else {
165+ // Extrapolated triangle location based on a fictitious time step
152166 dt = 0.25 * square_root_area / mean_velocity ;
153167 x1 = pp1 -> c [0 ] + mesh -> velocity [p1 - 1 ]* dt ;
154168 x2 = pp2 -> c [0 ] + mesh -> velocity [p2 - 1 ]* dt ;
@@ -158,14 +172,17 @@ double MMG2D_minangle(MMG5_pMesh mesh, int k) {
158172 y3 = pp3 -> c [1 ] + mesh -> velocity [p3 - 1 + mesh -> mark ]* dt ;
159173 }
160174
175+ // Compute the length of the triangle edges
161176 length [0 ] = pow (x1 - x2 ,2. ) + pow (y1 - y2 ,2. );
162177 length [1 ] = pow (x2 - x3 ,2. ) + pow (y2 - y3 ,2. );
163178 length [2 ] = pow (x3 - x1 ,2. ) + pow (y3 - y1 ,2. );
164179
180+ // Compute the length of the fictitious triangle edges
165181 length [3 ] = pow (pp1 -> c [0 ] - pp2 -> c [0 ],2. ) + pow (pp1 -> c [1 ] - pp2 -> c [1 ],2. );
166182 length [4 ] = pow (pp2 -> c [0 ] - pp3 -> c [0 ],2. ) + pow (pp2 -> c [1 ] - pp3 -> c [1 ],2. );
167183 length [5 ] = pow (pp3 -> c [0 ] - pp1 -> c [0 ],2. ) + pow (pp3 -> c [1 ] - pp1 -> c [1 ],2. );
168184
185+ // Check if the edge lengths satisfy the extended tolerances
169186 double min = length [0 ];
170187 double max = length [0 ];
171188 int min_index = 0 ;
@@ -198,10 +215,12 @@ double MMG2D_minangle(MMG5_pMesh mesh, int k) {
198215 if (sqrt (min ) < mesh -> info .hmin /factor_min ) return 0. ;
199216 if (sqrt (max ) > factor_max * mesh -> info .hmax ) return 0. ;
200217
201- // If the triangle is on a boundary, the tolerance is multiplied by two
218+ // If the triangle is on a boundary, the angle tolerance is divided by two
219+ // to extend remeshing close to boundaries where the mesh can be often squeezed
202220 if (mesh -> info .bdy_adaptation ) {
203221 if ((pp1 -> tag & MG_BDY ) || (pp2 -> tag & MG_BDY ) || (pp3 -> tag & MG_BDY )) coef = c1 ;
204222 else {
223+ // The neighbours of boundary triangles also have a reduced tolerance
205224 for (int i = 0 ; i < 3 ; i ++ ) {
206225 int k1 = adja [i ]/3 ;
207226
@@ -222,6 +241,7 @@ double MMG2D_minangle(MMG5_pMesh mesh, int k) {
222241 }
223242 }
224243
244+ // Minimal angle
225245 double cosA2 = ( length [(min_index + 1 )%3 + 3 ] + length [(min_index + 2 )%3 + 3 ] - length [min_index + 3 ] ) /
226246 (2 * pow (length [(min_index + 1 )%3 + 3 ],0.5 )* pow (length [(min_index + 2 )%3 + 3 ],0.5 ));
227247
@@ -872,8 +892,11 @@ MMG5_int MMG2D_movtri(MMG5_pMesh mesh,MMG5_pSol met,int maxit,int8_t improve) {
872892 ier = MMG2D_movedgpt (mesh ,met ,ilist ,list ,improve );
873893 if ( ier ) ns ++ ;
874894 }
895+ else if (mesh -> info .isotropic_pt_relocation ) {
896+ ier = MMG2D_movintpt (mesh ,met ,ilist ,list ,improve );
897+ }
875898 else {
876- if ( met -> size == 3 && met -> m && ! mesh -> info . isotropic )
899+ if ( met -> size == 3 && met -> m )
877900 ier = MMG2D_movintpt_ani (mesh ,met ,ilist ,list ,improve );
878901 else
879902 ier = MMG2D_movintpt (mesh ,met ,ilist ,list ,improve );
0 commit comments