Remeshing advice

Question

Dear all,

I am playing with some domain changes, and I need to remesh my domain. Right now it's not smart, I just create a new mesh. However, the domain changes are really small (per step), so I could easily remesh a portion of it.

Can you give me some advice? I've never done this with Fenics, so I'm pretty newbie.

I will need to generate the domain (a cracked rectangle) using C++, with the equivalent of "PolygonalMeshGenerator.generate".

Thanks!

Answer 1

There is currently no support for remeshing of subdomains so you will need to remesh the whole thing, or: manually extract the portion you want to remesh, then remesh it so that it matches the remaining part of the mesh, then stitch the meshes back together again. That could make a nice utility function.

Comments

So, that would be ok, if I can do it :)

In essence:

1. I need to mark some cells, and that I think I can do: I did it in python for assigning boundary conditions, correct me if I'm wrong, but wouldn't that be the same concept?
2. Then, as far as I understand, I should remove the cells storing them apart, creating a de facto new domain.
3. After that, I must remesh my subdomain with I think some constraints on the boundary, leaving the nodes and edges untouched.
4. Finally, merge the new mesh with the old one.

Can you point me in the right direction? What classes do you think I should be looking for?

Obviously, if the computation is distributed, something will be different, but I don't know anyting about the inner workings of fenics... I'm a newbie!

Thanks!

---

My first comment would be that this is a relatively hard function to code up. For such a function to really shine, it should be implemented in C++, made available to Python via SWIG and also run seamlessly in parallel.

However: I would suggest coding this up in Python: mark the region you want to remesh using a CellFunction, then extract that region as a SubMesh, then take the BoundaryMesh of that, manually create a polygon from it, send that polygon to the mesh generator, get a new mesh, then finally stitch it together with the rest using MeshEditor to create a new Mesh.

---

Thanks @logg,

I'm trying this approach, but I am stuck with the polygon generation. I get an error on the orientation, it seems that the vertices are not oriented CCW: Unable to create polygon. Polygonv vertices must be given in counter clockwise order. This error was encountered inside CSGPrimitives2D.cpp..

I've tried with std::vector and an array (for no reason at all), but the error persists.

Do you have any recommendation?

Thanks!

    // Extract a subdomain
    SubMesh submesh(simple, marker_extractor, 1);

    // Get the boundary
    BoundaryMesh boundary(submesh, "local");

    // Vector of Points from the boundary
    std::vector<Point> vv;

    // Iterate and add to vv each point
    for (CellIterator c(boundary); !c.end(); ++c)
    {
        double v[2];
        //std::vector<double> v;
        c->get_vertex_coordinates(v);

        std::cout << "B idx " << c->index() << " coords " << v[0] << ", " << v[1];

//        for(auto p : v)
//        {
//            std::cout << p << ",";
//        }

        std::cout << std::endl;

        vv.push_back(Point(2, v));
    }

    // New polygon
    mshr::Polygon boundary_polygon(vv);