1 | #include <stdio.h> |
---|
2 | #include <stdlib.h> |
---|
3 | #include <string.h> |
---|
4 | #include <time.h> |
---|
5 | #include "MSTK.h" |
---|
6 | #include "MSTK_private.h" |
---|
7 | #ifdef __cplusplus |
---|
8 | extern "C" { |
---|
9 | #endif |
---|
10 | |
---|
11 | |
---|
12 | |
---|
13 | /* |
---|
14 | This function is a collective call |
---|
15 | It assign vertex PType similar as MESH_AssignGlobalIDs(), |
---|
16 | except that here global IDs are given |
---|
17 | |
---|
18 | It assigns all the elements with a POVERLAP vertex as POVERLAP element |
---|
19 | |
---|
20 | It also assign master partition ID |
---|
21 | |
---|
22 | Author(s): Duo Wang, Rao Garimella |
---|
23 | */ |
---|
24 | |
---|
25 | int MESH_LabelPType_Vertex(Mesh_ptr submesh, int rank, int num, MPI_Comm comm); |
---|
26 | int MESH_LabelPType_Face(Mesh_ptr submesh, int rank, int num, MPI_Comm comm); |
---|
27 | int MESH_LabelPType_Region(Mesh_ptr submesh, int rank, int num, MPI_Comm comm); |
---|
28 | |
---|
29 | int MESH_LabelPType(Mesh_ptr submesh, int rank, int num, MPI_Comm comm) { |
---|
30 | int nf, nr; |
---|
31 | RepType rtype; |
---|
32 | |
---|
33 | /* basic mesh information */ |
---|
34 | rtype = MESH_RepType(submesh); |
---|
35 | nf = MESH_Num_Faces(submesh); |
---|
36 | nr = MESH_Num_Regions(submesh); |
---|
37 | /* first label vertex */ |
---|
38 | MESH_LabelPType_Vertex(submesh, rank, num, comm); |
---|
39 | if (nr) |
---|
40 | MESH_LabelPType_Region(submesh, rank, num, comm); |
---|
41 | else if(nf) |
---|
42 | MESH_LabelPType_Face(submesh, rank, num, comm); |
---|
43 | else { |
---|
44 | MSTK_Report("MESH_LabelPType()","only send volume or surface mesh",MSTK_ERROR); |
---|
45 | exit(-1); |
---|
46 | } |
---|
47 | return 1; |
---|
48 | } |
---|
49 | |
---|
50 | int MESH_LabelPType_Vertex(Mesh_ptr submesh, int rank, int num, MPI_Comm comm) { |
---|
51 | int i, j, k, nv, nbv, ne, nf, mesh_info[10], nevs, nfes, nfv, natt, nset, ncomp, dir; |
---|
52 | int nfe; |
---|
53 | MVertex_ptr mv; |
---|
54 | MEdge_ptr me; |
---|
55 | MFace_ptr mf; |
---|
56 | List_ptr boundary_verts, fverts, mfedges; |
---|
57 | RepType rtype; |
---|
58 | char attname[256], msetname[256]; |
---|
59 | MType mtype; |
---|
60 | MAttType att_type; |
---|
61 | MAttrib_ptr attrib; |
---|
62 | MSet_ptr mset; |
---|
63 | int *list_attr, *list_mset; |
---|
64 | char *list_attr_names, *list_mset_names; |
---|
65 | int *loc; |
---|
66 | int iloc, num_ghost_verts, global_id; |
---|
67 | double coor[3]; |
---|
68 | for (i = 0; i < 10; i++) mesh_info[i] = 0; |
---|
69 | |
---|
70 | /* mesh_info store the mesh reptype, nv, ne, nf, nbv */ |
---|
71 | rtype = MESH_RepType(submesh); |
---|
72 | nv = MESH_Num_Vertices(submesh); |
---|
73 | ne = MESH_Num_Edges(submesh); |
---|
74 | nf = MESH_Num_Faces(submesh); |
---|
75 | |
---|
76 | mesh_info[0] = rtype; |
---|
77 | mesh_info[1] = nv; |
---|
78 | mesh_info[2] = ne; |
---|
79 | mesh_info[3] = nf; |
---|
80 | |
---|
81 | /* calculate number of boundary vertices */ |
---|
82 | nbv = 0; |
---|
83 | boundary_verts = List_New(10); |
---|
84 | for(i = 0; i < nv; i++) { |
---|
85 | mv = MESH_Vertex(submesh,i); |
---|
86 | if (MV_PType(mv)) { |
---|
87 | List_Add(boundary_verts,mv); |
---|
88 | nbv++; |
---|
89 | } |
---|
90 | } |
---|
91 | mesh_info[4] = nbv; |
---|
92 | |
---|
93 | /* sort boundary vertices based on global ID, for binary search */ |
---|
94 | List_Sort(boundary_verts,nbv,sizeof(MVertex_ptr),compareGlobalID); |
---|
95 | |
---|
96 | /* |
---|
97 | gather submeshes information |
---|
98 | right now we only need nv and nbv, and later num_ghost_verts, but we gather all mesh_info |
---|
99 | */ |
---|
100 | int *global_mesh_info = (int *)MSTK_malloc(10*num*sizeof(int)); |
---|
101 | MPI_Allgather(mesh_info,10,MPI_INT,global_mesh_info,10,MPI_INT,comm); |
---|
102 | |
---|
103 | /* get largest number of boundary vertices of all the processors */ |
---|
104 | int max_nbv = 0; |
---|
105 | for(i = 0; i < num; i++) |
---|
106 | if(max_nbv < global_mesh_info[10*i+4]) |
---|
107 | max_nbv = global_mesh_info[10*i+4]; |
---|
108 | |
---|
109 | int *list_boundary_vertex = (int *)MSTK_malloc(max_nbv*sizeof(int)); |
---|
110 | int *recv_list_vertex = (int *)MSTK_malloc(num*max_nbv*sizeof(int)); |
---|
111 | |
---|
112 | /* only global ID are sent */ |
---|
113 | int index_nbv = 0; |
---|
114 | for(i = 0; i < nbv; i++) { |
---|
115 | mv = List_Entry(boundary_verts,i); |
---|
116 | list_boundary_vertex[index_nbv] = MV_GlobalID(mv); |
---|
117 | index_nbv++; |
---|
118 | } |
---|
119 | |
---|
120 | /* gather boundary vertices */ |
---|
121 | MPI_Allgather(list_boundary_vertex,max_nbv,MPI_INT,recv_list_vertex,max_nbv,MPI_INT,comm); |
---|
122 | |
---|
123 | /* indicate if a vertex is overlapped */ |
---|
124 | int *vertex_ov_label = (int *)MSTK_malloc(num*max_nbv*sizeof(int)); |
---|
125 | |
---|
126 | /* |
---|
127 | store the local boundary id on ov processor |
---|
128 | it is used to assign global id of local ghost vertices |
---|
129 | do not need to store master partition id, MV_MasterParID(mv) is assigned |
---|
130 | */ |
---|
131 | for (i = 0; i < num*max_nbv; i++) |
---|
132 | vertex_ov_label[i] = 0; |
---|
133 | num_ghost_verts = 0; |
---|
134 | /* for processor other than 0 */ |
---|
135 | for(i = 0; i < nbv; i++) { |
---|
136 | mv = List_Entry(boundary_verts,i); |
---|
137 | global_id = MV_GlobalID(mv); |
---|
138 | /* check which processor has the same coordinate vertex |
---|
139 | Different from assigning global id, we check all the other processors, from large to small |
---|
140 | by rank, whenever a vertex is found, mark rank and j as neighbors, and the masterparid is still |
---|
141 | the smallest rank processor |
---|
142 | */ |
---|
143 | for(j = num-1; j >= 0; j--) { |
---|
144 | if(j == rank) continue; |
---|
145 | /* since each processor has sorted the boundary vertices, use binary search */ |
---|
146 | loc = (int *)bsearch(&global_id, |
---|
147 | &recv_list_vertex[max_nbv*j], |
---|
148 | global_mesh_info[10*j+4], |
---|
149 | sizeof(int), |
---|
150 | compareINT); |
---|
151 | /* if found the vertex on previous processors */ |
---|
152 | if(loc) { |
---|
153 | /* here the location iloc is relative to the beginning of the jth processor */ |
---|
154 | iloc = (int)(loc - &recv_list_vertex[max_nbv*j]); |
---|
155 | MV_Set_PType(mv,PGHOST); |
---|
156 | MV_Set_MasterParID(mv,j); |
---|
157 | MESH_Flag_Has_Ghosts_From_Prtn(submesh,j,MVERTEX); |
---|
158 | num_ghost_verts++; |
---|
159 | /* label the original vertex as overlapped */ |
---|
160 | vertex_ov_label[max_nbv*j+iloc] |= 1; |
---|
161 | /* if found on processor j, no need to test other processors*/ |
---|
162 | break; |
---|
163 | } |
---|
164 | } |
---|
165 | } |
---|
166 | |
---|
167 | mesh_info[5] = num_ghost_verts; |
---|
168 | |
---|
169 | /* update ghost verts number */ |
---|
170 | MPI_Allgather(mesh_info,10,MPI_INT,global_mesh_info,10,MPI_INT,comm); |
---|
171 | /* since this is a or reduction, we can use MPI_IN_PLACE, send buffer same as recv buffer */ |
---|
172 | MPI_Allreduce(MPI_IN_PLACE,vertex_ov_label,num*max_nbv,MPI_INT,MPI_LOR,comm); |
---|
173 | |
---|
174 | for(i = 0; i < nv; i++) { |
---|
175 | mv = MESH_Vertex(submesh,i); |
---|
176 | if (MV_PType(mv) == PGHOST) |
---|
177 | continue; |
---|
178 | MV_Set_MasterParID(mv,rank); |
---|
179 | } |
---|
180 | |
---|
181 | /* label OVERLAP vertices */ |
---|
182 | for(i = 0; i < nbv; i++) { |
---|
183 | if(vertex_ov_label[rank*max_nbv+i]) { |
---|
184 | mv = List_Entry(boundary_verts,i); |
---|
185 | MV_Set_PType(mv,POVERLAP); |
---|
186 | } |
---|
187 | } |
---|
188 | |
---|
189 | List_Delete(boundary_verts); |
---|
190 | MSTK_free(global_mesh_info); |
---|
191 | MSTK_free(vertex_ov_label); |
---|
192 | |
---|
193 | MSTK_free(list_boundary_vertex); |
---|
194 | MSTK_free(recv_list_vertex); |
---|
195 | return 1; |
---|
196 | } |
---|
197 | |
---|
198 | /* |
---|
199 | Right now assume faces are not overlapped across processors |
---|
200 | Label the face that has OVERLAP or Ghost vertex as OVERLAP |
---|
201 | |
---|
202 | */ |
---|
203 | int MESH_LabelPType_Face(Mesh_ptr submesh, int rank, int num, MPI_Comm comm) { |
---|
204 | int i, idx; |
---|
205 | MVertex_ptr mv; |
---|
206 | MFace_ptr mf; |
---|
207 | List_ptr fverts; |
---|
208 | idx = 0; |
---|
209 | while( (mf = MESH_Next_Face(submesh,&idx)) ) { |
---|
210 | fverts = MF_Vertices(mf,1,0); |
---|
211 | for(i = 0; i < List_Num_Entries(fverts); i++) { |
---|
212 | mv = List_Entry(fverts,i); |
---|
213 | if( MV_PType(mv) == PGHOST || MV_PType(mv) == POVERLAP ) { |
---|
214 | MF_Set_PType(mf,POVERLAP); |
---|
215 | break; |
---|
216 | } |
---|
217 | } |
---|
218 | List_Delete(fverts); |
---|
219 | } |
---|
220 | |
---|
221 | return 1; |
---|
222 | } |
---|
223 | |
---|
224 | /* |
---|
225 | right now assume regions are not overlapped across processors |
---|
226 | Label the region that has OVERLAP vertex as OVERLAP |
---|
227 | */ |
---|
228 | |
---|
229 | |
---|
230 | int MESH_LabelPType_Region(Mesh_ptr submesh, int rank, int num, MPI_Comm comm) { |
---|
231 | int i, idx; |
---|
232 | MRegion_ptr mr; |
---|
233 | MVertex_ptr mv; |
---|
234 | List_ptr rverts; |
---|
235 | idx = 0; |
---|
236 | while( (mr=MESH_Next_Region(submesh,&idx)) ) { |
---|
237 | rverts = MR_Vertices(mr); |
---|
238 | for(i = 0; i < List_Num_Entries(rverts); i++) { |
---|
239 | mv = List_Entry(rverts,i); |
---|
240 | if(MV_PType(mv) == PGHOST || MV_PType(mv) == POVERLAP) { |
---|
241 | MR_Set_PType(mr,POVERLAP); |
---|
242 | break; |
---|
243 | } |
---|
244 | } |
---|
245 | List_Delete(rverts); |
---|
246 | } |
---|
247 | |
---|
248 | return 1; |
---|
249 | } |
---|
250 | |
---|
251 | |
---|
252 | |
---|
253 | #ifdef __cplusplus |
---|
254 | } |
---|
255 | #endif |
---|
256 | |
---|