@@ -910,132 +910,128 @@ def _call_gmsh(pmesh, phases, res):
910910
911911 # ---------------------------------------------------------------------- #
912912 # CREATE GEOMETRY
913- # ---------------------------------------------------------------------- #
914- geom = pg .built_in .Geometry ()
915-
916- # Add points
917- pts = [geom .add_point (_pt3d (pt ), res ) for pt in pmesh .points ]
918- n_dim = len (pmesh .points [0 ])
919-
920- # Add edges to geometry
921- for edge in edge_keys :
922- line = geom .add_line (* [pts [kp ] for kp in edge ])
923- edge_lines .append (line )
913+ # ----------------------------------------------------------------------
914+ with pg .geo .Geometry () as geom :
915+ # Add points
916+ pts = [geom .add_point (_pt3d (pt ), res ) for pt in pmesh .points ]
917+ n_dim = len (pmesh .points [0 ])
918+
919+ # Add edges to geometry
920+ for edge in edge_keys :
921+ line = geom .add_line (* [pts [kp ] for kp in edge ])
922+ edge_lines .append (line )
923+
924+ if n_dim == 2 :
925+ lbl = 'facet-{}' .format (edges_info [edge ]['facets' ][0 ])
926+ if facet_check (edges_info [edge ]['neighbors' ], pmesh , phases ):
927+ geom .add_physical (edge_lines [- 1 ], lbl )
924928
925929 if n_dim == 2 :
926- lbl = 'facet-{}' .format (edges_info [edge ]['facets' ][0 ])
927- if facet_check (edges_info [edge ]['neighbors' ], pmesh , phases ):
928- geom .add_physical (edge_lines [- 1 ], lbl )
929-
930- if n_dim == 2 :
931- # Add surfaces to geometry
932- loops = []
933- surfs = []
934- seed_facets = {}
935- for i , r in enumerate (pmesh .regions ):
936- s = pmesh .seed_numbers [i ]
937- seed_facets .setdefault (s , set ()).symmetric_difference_update (r )
938- for i in seed_facets :
939- region = list (seed_facets [i ])
940- sorted_pairs = _sort_facets ([pmesh .facets [f ] for f in region ])
941- loop = []
942- for facet in sorted_pairs :
943- key = tuple (sorted (facet ))
944- if facet [0 ] == key [0 ]:
945- sgn = 1
946- else :
947- sgn = - 1
930+ # Add surfaces to geometry
931+ loops = []
932+ surfs = []
933+ seed_facets = {}
934+ for i , r in enumerate (pmesh .regions ):
935+ s = pmesh .seed_numbers [i ]
936+ seed_facets .setdefault (s , set ()).symmetric_difference_update (r )
937+ for i in seed_facets :
938+ region = list (seed_facets [i ])
939+ sorted_pairs = _sort_facets ([pmesh .facets [f ] for f in region ])
940+ loop = []
941+ for facet in sorted_pairs :
942+ key = tuple (sorted (facet ))
943+ if facet [0 ] == key [0 ]:
944+ sgn = 1
945+ else :
946+ sgn = - 1
948947
949- n = edges_info [key ]['ind' ]
950- line = edge_lines [n ]
951- if sgn > 0 :
952- loop .append (line )
953- else :
954- loop .append (- line )
948+ n = edges_info [key ]['ind' ]
949+ line = edge_lines [n ]
950+ if sgn > 0 :
951+ loop .append (line )
952+ else :
953+ loop .append (- line )
955954
956- loops .append (geom .add_line_loop (loop ))
957- surfs .append (geom .add_plane_surface (loops [- 1 ]))
958- p_num = pmesh .phase_numbers [i ]
959- mat_type = phases [p_num ].get ('material_type' , 'solid' )
960- if mat_type not in _misc .kw_void :
961- geom .add_physical (surfs [- 1 ], 'seed-' + str (i ))
955+ loops .append (geom .add_curve_loop (loop ))
956+ surfs .append (geom .add_plane_surface (loops [- 1 ]))
957+ p_num = pmesh .phase_numbers [i ]
958+ mat_type = phases [p_num ].get ('material_type' , 'solid' )
959+ if mat_type not in _misc .kw_void :
960+ geom .add_physical (surfs [- 1 ], 'seed-' + str (i ))
962961
963- elif n_dim == 3 :
964- # Add surfaces to geometry
965- loops = []
966- surfs = []
967- seed_surfs = {}
968- seed_phases = dict (zip (pmesh .seed_numbers , pmesh .phase_numbers ))
969- for i in facets_info :
970- info = facets_info [i ]
971- facet_seeds = info ['seeds' ]
972- to_add = len (facet_seeds ) < 2 or facet_seeds [0 ] != facet_seeds [1 ]
973- if not to_add :
974- surfs .append ('' )
975- continue
976-
977- loop = []
978- for n , sgn in zip (info ['edge_numbers' ], info ['edge_signs' ]):
979- line = edge_lines [n ]
980- if sgn > 0 :
981- loop .append (line )
982- else :
983- loop .append (- line )
984- loops .append (geom .add_line_loop (loop ))
985- surfs .append (geom .add_plane_surface (loops [- 1 ]))
986- if facet_check (info ['neighbors' ], pmesh , phases ):
987- geom .add_physical (surfs [- 1 ], 'facet-' + str (i ))
988- for seed_num in facet_seeds :
989- if seed_num not in seed_surfs :
990- seed_surfs [seed_num ] = []
991- seed_surfs [seed_num ].append (surfs [- 1 ])
992-
993- # Add volumes to geometry
994- surf_loops = []
995- volumes = []
996- for seed_num in seed_surfs :
997- surf_loop = seed_surfs [seed_num ]
998- surf_loops .append (geom .add_surface_loop (surf_loop ))
999- volumes .append (geom .add_volume (surf_loops [- 1 ]))
1000-
1001- p_num = seed_phases [seed_num ]
1002- mat_type = phases [p_num ].get ('material_type' , 'solid' )
1003- if mat_type not in _misc .kw_void :
1004- geom .add_physical (volumes [- 1 ], 'seed-' + str (seed_num ))
1005- else :
1006- raise ValueError ('Points cannot have dimension ' + str (n_dim ) + '.' )
962+ elif n_dim == 3 :
963+ # Add surfaces to geometry
964+ loops = []
965+ surfs = []
966+ seed_surfs = {}
967+ seed_phases = dict (zip (pmesh .seed_numbers , pmesh .phase_numbers ))
968+ for i in facets_info :
969+ info = facets_info [i ]
970+ facet_seeds = info ['seeds' ]
971+ to_add = len (facet_seeds ) < 2 or facet_seeds [0 ] != facet_seeds [1 ]
972+ if not to_add :
973+ surfs .append ('' )
974+ continue
975+
976+ loop = []
977+ for n , sgn in zip (info ['edge_numbers' ], info ['edge_signs' ]):
978+ line = edge_lines [n ]
979+ if sgn > 0 :
980+ loop .append (line )
981+ else :
982+ loop .append (- line )
983+ loops .append (geom .add_curve_loop (loop ))
984+ surfs .append (geom .add_plane_surface (loops [- 1 ]))
985+ if facet_check (info ['neighbors' ], pmesh , phases ):
986+ geom .add_physical (surfs [- 1 ], 'facet-' + str (i ))
987+ for seed_num in facet_seeds :
988+ if seed_num not in seed_surfs :
989+ seed_surfs [seed_num ] = []
990+ seed_surfs [seed_num ].append (surfs [- 1 ])
991+
992+ # Add volumes to geometry
993+ surf_loops = []
994+ volumes = []
995+ for seed_num in seed_surfs :
996+ surf_loop = seed_surfs [seed_num ]
997+ surf_loops .append (geom .add_surface_loop (surf_loop ))
998+ volumes .append (geom .add_volume (surf_loops [- 1 ]))
999+
1000+ p_num = seed_phases [seed_num ]
1001+ mat_type = phases [p_num ].get ('material_type' , 'solid' )
1002+ if mat_type not in _misc .kw_void :
1003+ geom .add_physical (volumes [- 1 ], 'seed-' + str (seed_num ))
1004+ else :
1005+ raise ValueError ('Points cannot have dimension ' + str (n_dim ) + '.' )
10071006
1008- # ---------------------------------------------------------------------- #
1009- # CALL GMSH
1010- # ---------------------------------------------------------------------- #
1011- mesh = pg .generate_mesh (geom )
1007+ mesh = geom .generate_mesh ()
10121008
10131009 # ---------------------------------------------------------------------- #
10141010 # CREATE MICROSTRUCTPY.MESHING.TRIMESH
10151011 # ---------------------------------------------------------------------- #
1016- f_key = {2 : 'line' , 3 : 'triangle' }[n_dim ]
1017- e_key = {2 : 'triangle' , 3 : 'tetra' }[n_dim ]
1012+ f_ind = {2 : 0 , 3 : 1 }[n_dim ]
1013+ e_ind = {2 : 1 , 3 : 2 }[n_dim ]
10181014
10191015 pts = np .array (mesh .points )[:, :n_dim ]
1020- facets = mesh .cells_dict [ f_key ]
1016+ facets = mesh .cells [ f_ind ]. data
10211017
10221018 # Sort Element Keypoints for Positive Volume
1023- tets = [e [_sort_element ([mesh .points [k ] for k in e ])]
1024- for e in mesh .cells_dict [ e_key ]]
1025-
1026- facet_phsy2att = {}
1027- seed_phys2att = {}
1028- for key in mesh . field_data :
1029- phys = mesh .field_data [ key ][ 0 ]
1030- if key . startswith ( 'facet-' ):
1031- facet_phsy2att [ phys ] = int ( key .split ('-' )[ - 1 ] )
1032- elif key . startswith ( 'seed-' ):
1033- seed_phys2att [ phys ] = int ( key . split ( '-' )[ - 1 ])
1034-
1035- atts = mesh . cell_data_dict [ 'gmsh:physical' ]
1036- tet_atts = [ seed_phys2att [ k ] for k in atts [ e_key ]]
1037- tet_atts = [ amorph_seeds . get ( sd , sd ) for sd in tet_atts ]
1038- facet_atts = [ facet_phsy2att [ k ] for k in atts [ f_key ]]
1019+ tets = np . array ( [e [_sort_element ([mesh .points [k ] for k in e ])]
1020+ for e in mesh .cells [ e_ind ]. data ])
1021+
1022+ tet_atts = np . array ([ - 1 for tet in tets ])
1023+ facet_atts = np . array ([ - 1 for f in facets ])
1024+
1025+ n_facets = len ( mesh .cells [ f_ind ]. data )
1026+ for key , elem_sets in mesh . cell_sets . items ( ):
1027+ set_kind , set_num_str = key .split ('-' )
1028+ att = int ( set_num_str )
1029+ if set_kind == 'seed' :
1030+ elem_set = elem_sets [ e_ind ] - n_facets
1031+ tet_atts [ elem_set ] = amorph_seeds . get ( att , att )
1032+ elif set_kind == 'facet' :
1033+ elem_set = elem_sets [ f_ind ]
1034+ facet_atts [ elem_set ] = att
10391035
10401036 tri_args = (pts , tets , tet_atts , facets , facet_atts )
10411037 return tri_args
0 commit comments