(mesh)
| 22 | return parser.parse_args(); |
| 23 | |
| 24 | def compute_distortion_energies_3D(mesh): |
| 25 | if mesh.num_voxels > 0 and mesh.vertex_per_voxel != 4: |
| 26 | raise RuntimeError("Only tet mesh is supported for distortion computation"); |
| 27 | |
| 28 | regular_tet = pymesh.generate_regular_tetrahedron(); |
| 29 | assembler = pymesh.Assembler(regular_tet); |
| 30 | G = assembler.assemble("gradient"); |
| 31 | |
| 32 | vertices = mesh.vertices; |
| 33 | tets = mesh.voxels; |
| 34 | Js = [ G * vertices[tet] for tet in tets ]; |
| 35 | |
| 36 | J_F = np.array([np.trace(np.dot(J.T, J)) for J in Js]); |
| 37 | J_det = np.array([numpy.linalg.det(J) for J in Js]); |
| 38 | invert_J = lambda args: np.full((3,3), np.inf) if args[1] == 0 else numpy.linalg.inv(args[0]) |
| 39 | J_inv = map(invert_J, zip(Js, J_det)); |
| 40 | J_inv_F = np.array([np.trace(np.dot(Ji.T, Ji)) for Ji in J_inv]); |
| 41 | |
| 42 | conformal_amips = np.divide(J_F, np.cbrt(np.square(J_det))); |
| 43 | finite_conformal_amips = np.isfinite(conformal_amips); |
| 44 | symmetric_dirichlet = J_F + J_inv_F; |
| 45 | finite_symmetric_dirichlet = np.isfinite(symmetric_dirichlet); |
| 46 | orientations = pymesh.get_tet_orientations(mesh); |
| 47 | orientations[orientations > 0] = 1; |
| 48 | orientations[orientations < 0] = -1; |
| 49 | |
| 50 | num_degenerate_tets = np.count_nonzero(orientations==0); |
| 51 | num_inverted_tets = np.count_nonzero(orientations<0); |
| 52 | num_nonfinite_amips = np.count_nonzero(np.logical_not(finite_conformal_amips)); |
| 53 | num_nonfinite_dirichlet =\ |
| 54 | np.count_nonzero(np.logical_not(finite_symmetric_dirichlet)); |
| 55 | logger = logging.getLogger("Distorsion"); |
| 56 | if num_degenerate_tets > 0: |
| 57 | logger.warn("degenerate tets: {}".format(num_degenerate_tets)); |
| 58 | if num_inverted_tets > 0: |
| 59 | logger.warn("inverted tets: {}".format(num_inverted_tets)); |
| 60 | if num_nonfinite_amips > 0: |
| 61 | logger.warn("Non-finite conformal AMIPS: {}".format( |
| 62 | num_nonfinite_amips)); |
| 63 | if num_nonfinite_dirichlet > 0: |
| 64 | logger.warn("Non-finite symmetric Dirichlet: {}".format( |
| 65 | num_nonfinite_dirichlet)); |
| 66 | |
| 67 | mesh.add_attribute("conformal_AMIPS"); |
| 68 | mesh.set_attribute("conformal_AMIPS", conformal_amips); |
| 69 | mesh.add_attribute("finite_conformal_AMIPS"); |
| 70 | mesh.set_attribute("finite_conformal_AMIPS", finite_conformal_amips); |
| 71 | mesh.add_attribute("symmetric_Dirichlet"); |
| 72 | mesh.set_attribute("symmetric_Dirichlet", symmetric_dirichlet); |
| 73 | mesh.add_attribute("finite_symmetric_Dirichlet"); |
| 74 | mesh.set_attribute("finite_symmetric_Dirichlet", finite_symmetric_dirichlet); |
| 75 | mesh.add_attribute("orientations"); |
| 76 | mesh.set_attribute("orientations", orientations); |
| 77 | |
| 78 | def compute_distortion_energies_2D(mesh): |
| 79 | assert(mesh.dim == 2); |
no test coverage detected