Function repousse

scripts/repousse.py:323–381

(mesh, logger)

    return pymesh.form_mesh(vertices, mesh.faces);

def repousse(mesh, logger):
    cell_ids = mesh.get_attribute("cell").ravel().astype(int);
    mesh.add_attribute("edge_length");
    tol = np.amax(mesh.get_attribute("edge_length")) * 0.1;

    bbox_min, bbox_max = mesh.bbox;
    scaling = 2.0 / norm(bbox_max - bbox_min);

    start_time = time();
    num_cells = np.amax(cell_ids)+1;
    results = [];
    for i in range(num_cells):
        to_keep = np.arange(mesh.num_faces, dtype=int)[cell_ids == i];
        if not np.any(to_keep):
            continue;

        cut_mesh = pymesh.submesh(mesh, to_keep, 0);
        pymesh.save_mesh("debug.msh", cut_mesh);
        cut_mesh, __ = pymesh.remove_degenerated_triangles(cut_mesh, 100);
        cut_mesh, __ = pymesh.split_long_edges(cut_mesh, tol);

        dof = cut_mesh.num_vertices;
        assembler = pymesh.Assembler(cut_mesh);
        L = assembler.assemble("laplacian");
        M = assembler.assemble("mass");

        L_rhs = M * np.ones(dof) * -0.5;
        bd_indices = cut_mesh.boundary_vertices;
        n = len(bd_indices);
        C = scipy.sparse.coo_matrix((np.ones(n),
            (np.arange(n, dtype=int), bd_indices)), shape=(n, dof));
        C_rhs = np.zeros(n);

        A = scipy.sparse.bmat([
            [-L, C.T],
            [C, None]
            ]);
        rhs = np.concatenate((L_rhs.ravel(), C_rhs));

        solver = pymesh.SparseSolver.create("SparseLU");
        solver.compute(A);
        x = solver.solve(rhs);

        z = x[:dof].reshape((-1, 1));

        vertices = np.hstack((cut_mesh.vertices, z));
        out_mesh = pymesh.form_mesh(vertices, cut_mesh.faces);
        results.append(out_mesh);

    finish_time = time();
    t = finish_time - start_time;
    logger.info("Repousse running time: {}".format(t));

    mesh = pymesh.merge_meshes(results);

    vertices = mesh.vertices[:,:2];
    mesh_2d = pymesh.form_mesh(vertices, mesh.faces);
    pymesh.save_mesh("out_2d.msh", mesh_2d)