femwell icon indicating copy to clipboard operation
femwell copied to clipboard

mesh_from_Dict does not handle MultiLineStrings()

Open duarte-jfs opened this issue 2 months ago • 4 comments

I have noticed an issue when using mesh_from_Dict(), stemming from the following lines:

       # Add overall bounding box
        total = unary_union(list(shapes_dict.values()))
        all_polygons = [total, *list(shapes_dict.values())]
        # Break up all shapes so that plane is tiled with non-overlapping layers, get the maximal number of fragments
        # Equivalent to BooleanFragments
        np.seterr(invalid="ignore")
        listpoly = [a.intersection(b) for a, b in combinations(all_polygons, 2)]
        np.seterr(**initial_settings)
        rings = [
            LineString(list(object.exterior.coords))
            for object in listpoly
            if not (object.geom_type in ["Point", "LineString"] or object.is_empty)
        ]

This is due to the fact that listpoly will have a MultiLineString when intersecting a LinearString with total. The change that could (probably) solve this is:

       # Add overall bounding box
        total = unary_union(list(shapes_dict.values()))
        all_polygons = [total, *list(shapes_dict.values())]
        # Break up all shapes so that plane is tiled with non-overlapping layers, get the maximal number of fragments
        # Equivalent to BooleanFragments
        np.seterr(invalid="ignore")
        listpoly = [a.intersection(b) for a, b in combinations(all_polygons, 2)]
        np.seterr(**initial_settings)
        rings = [
            LineString(list(object.exterior.coords))
            for object in listpoly
            if not (object.geom_type in ["Point", "LineString", "MultiLineString"] or object.is_empty)
        ]

Minimal working example


from collections import OrderedDict

import numpy as np

from shapely.geometry import box, LineString, MultiLineString
from skfem.io.meshio import from_meshio

from femwell.mesh import mesh_from_OrderedDict,mesh_from_Dict

box1 = box(-1,-1,1,1)
box2 = box(-1,5,1,6)

surface = box1.buffer(0.5).exterior
polygons = dict(box1 = box1,
               box2 = box2,
               surface = surface)

resolutions = dict(box1 = {'resolution': 0.2, 'distance':1},
                   box2 =  {'resolution': 0.2, 'distance':1},
                   surface =  {'resolution': 0.2, 'distance':1})

mesh = from_meshio(mesh_from_Dict(polygons, resolutions))

mesh.draw()

duarte-jfs avatar Apr 19 '24 16:04 duarte-jfs