Enhancements for replacing grids functionality

[jhabriel]

The replacing grids and interfaces functionality has two important restrictions:

1. Some of the methods assume that the new grids are Gmsh-generated, which is no longer true in the case of grids created via AMR routines.
2. Replacement is limited to mixed-dimensional grids whose ambient dimension is 2.

The first restriction can be tackled with careful geometric sanity checks, whereas the second restriction is more challenging, although it seems that the current implementation can be extended with some effort to the 3d case.

Another point to consider is whether it makes sense to centralise this functionality into a new class that manages the grid replacement process. Currently, the functionality depends on several methods/functions that are located in different parts. These are: grids.md_grid.replace_subdomains and interfaces(), grids.mortar_grid.update_mortar, grids.mortar_grid.update_primary, grids.mortar_grid.update_secondary, and grids.match_grids. Not sure whether the update methods are used in other parts or only for replacement.

PD: Since I will have to address the two aforementioned restrictions for my AMR project, it could be a good opportunity to incorporate these changes into PorePy.

[keileg]

Sounds very interesting, @jhabriel. Let's stay in touch.

[keileg]

Notes after a discussion today:

Updates of mortar grids

MortarGrid.update_secondary() must be updated to accommodate incremental modifications of grid. Specifically, lines 444-9,

        self._secondary_to_mortar_avg = sps.bmat(matrix_avg, format="csc")
        self._secondary_to_mortar_int = sps.bmat(matrix_int, format="csc")

must be updated to include the existing mapping. A template can be found in the similar method update_primary(), lines 528-9. It is important that the existing mapping is inserted at the right place (either left or right multiplication, doing this wrong will most likely lead to an array mismatch).

Matching of 2d grids with embedded fractures

In match_grids.py, line 423-5 faces along a line fracture are found by combining geometric information (faces_by_hit) and tags (g_new.tags["fracture_faces"]), where the latter is implicitly assumed available. This will not necessarily be the case for grids not created by the standard mesh construction workflow in PorePy, for instance grids modified by AMR. There are two natural approaches:

1. Use a distance criterion (point to polygon for all points of a face) to assign fracture face tags also for AMR grids.
2. Relax the test in match_grids to consider face tags only if these are available, and then delete the item from the set of tags in AMR grids.

Of these, only option 1 will leave a grid which is fully compatible with the wider PorePy functionality.

[IvarStefansson]

Would there, at least in theory, be a third option: Make the tag mandatory?

[keileg]

See also #1118