source: issm/trunk-jpl/src/c/modules/ModelProcessorx/ElementsAndVerticesPartitioning.cpp@ 23504

Last change on this file since 23504 was 23504, checked in by Mathieu Morlighem, 6 years ago

CHG: cosmetics
File size: 5.0 KB
Line 
1/*!\file: ElementsAndVerticesPartitioning.cpp
2 * \brief: partition elements and nodes and vertices
3 */
4
5#ifdef HAVE_CONFIG_H
6 #include <config.h>
7#else
8#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
9#endif
10
11#include <string.h>
12#include "../../classes/classes.h"
13#include "../../shared/shared.h"
14#include "../MeshPartitionx/MeshPartitionx.h"
15#include "../ModelProcessorx/ModelProcessorx.h"
16
17void ElementsAndVerticesPartitioning(bool** pmy_elements, int** pmy_vertices, IoModel* iomodel){
18
19 int numberofelements2d;
20 int numberofvertices2d;
21 int numlayers;
22 int numvertex_pairing;
23
24 /*intermediary: */
25 int *epart = NULL; //element partitioning.
26 int *npart = NULL; //node partitioning.
27 int elements_width; //number of columns in elements (2d->3, 3d->6)
28 int el1,el2;
29 int *elements2d = NULL;
30
31 /*Get my_rank:*/
32 int my_rank = IssmComm::GetRank();
33 int num_procs = IssmComm::GetSize();
34
35 /*First, check that partitioning has not yet been carryed out. Just check whether my_elements pointers is not already assigned a value: */
36 if(*pmy_elements) return;
37
38 /*Number of vertices per elements, needed to correctly retrieve data: */
39 /*Determine parallel partitioning of elements: we use Metis for now. First load the data, then partition*/
40 switch(iomodel->meshelementtype){
41 case TriaEnum:
42 elements_width=3;
43 numberofelements2d = 0;
44 numberofvertices2d = 0;
45 numlayers = 0;
46 break;
47 case TetraEnum:
48 elements_width=4;
49 numberofelements2d = 0;
50 numberofvertices2d = 0;
51 numlayers = 0;
52 break;
53 case PentaEnum:
54 elements_width=6;
55 iomodel->FetchData(&elements2d,NULL,NULL,"md.mesh.elements2d");
56 iomodel->FindConstant(&numberofelements2d,"md.mesh.numberofelements2d");
57 iomodel->FindConstant(&numberofvertices2d,"md.mesh.numberofvertices2d");
58 iomodel->FindConstant(&numlayers,"md.mesh.numberoflayers");
59 break;
60 default:
61 _error_("mesh elements "<< EnumToStringx(iomodel->meshelementtype) <<" not supported yet");
62 }
63
64 /*Partition and free resouces*/
65 MeshPartitionx(&epart,&npart,iomodel->numberofelements,iomodel->numberofvertices,iomodel->elements,numberofelements2d,numberofvertices2d,elements2d,numlayers,elements_width,iomodel->meshelementtype,num_procs);
66 xDelete<int>(elements2d);
67 xDelete<int>(npart);
68
69 /*Deal with rifts, they have to be included into one partition only, not several: */
70 int numrifts;
71 iomodel->FindConstant(&numrifts,"md.rifts.numrifts");
72 if(numrifts){
73 IssmDouble *riftinfo = NULL;
74 iomodel->FetchData(&riftinfo,&numrifts,NULL,"md.rifts.riftstruct");
75 for(int i=0;i<numrifts;i++){
76 const int RIFTINFOSIZE = 12;
77 el1=reCast<int>(*(riftinfo+RIFTINFOSIZE*i+2))-1; //matlab indexing to c indexing
78 el2=reCast<int>(*(riftinfo+RIFTINFOSIZE*i+3))-1; //matlab indexing to c indexing
79 epart[el2]=epart[el1]; //ensures that this pair of elements will be in the same partition, as well as the corresponding vertices;
80 }
81 iomodel->DeleteData(riftinfo,"md.rifts.riftstruct");
82 }
83
84 /*Create my_vertices and my_elements, used by each partition */
85 bool *my_elements = xNewZeroInit<bool>(iomodel->numberofelements);
86 int *my_vertices = xNewZeroInit<int>(iomodel->numberofvertices);
87
88 /*Start figuring out, out of the partition, which elements belong to this cpu: */
89 for(int i=0;i<iomodel->numberofelements;i++){
90
91 /*!All elements have been partitioned above, only deal with elements for this cpu: */
92 if(my_rank==epart[i]){
93 my_elements[i]=true;
94 /*Now that we are here, we can also start building the list of vertices belonging to this cpu partition: we use
95 *the element index to do this. For each element n, we know index[n][0:2] holds the indices (matlab indexing)
96 into the vertices coordinates. If we start plugging 1 into my_vertices for each index[n][i] (i=0:2), then my_vertices
97 will hold which vertices belong to this partition*/
98 for(int j=0;j<elements_width;j++){
99 my_vertices[iomodel->elements[elements_width*i+j]-1]=1;
100 }
101 }
102 }
103
104 /*We might have vertex_pairing in which case, some vertices have to be cloned:
105 * penpair has 2 nodes that are poointing toward 2 vertices.
106 * The 2 vertices must be in the same cpu as the penpair*/
107 int *vertex_pairing = NULL;
108 iomodel->FetchData(&vertex_pairing,&numvertex_pairing,NULL,"md.stressbalance.vertex_pairing");
109 for(int i=0;i<numvertex_pairing;i++){
110 if(my_vertices[vertex_pairing[2*i+0]-1] && !my_vertices[vertex_pairing[2*i+1]-1]){
111 my_vertices[vertex_pairing[2*i+1]-1]=2; //to know that these elements are not on the partition
112 }
113 }
114 xDelete<int>(vertex_pairing);
115 iomodel->FetchData(&vertex_pairing,&numvertex_pairing,NULL,"md.masstransport.vertex_pairing");
116 for(int i=0;i<numvertex_pairing;i++){
117 if(my_vertices[vertex_pairing[2*i+0]-1] && !my_vertices[vertex_pairing[2*i+1]-1]){
118 my_vertices[vertex_pairing[2*i+1]-1]=2; //to know that these elements are not on the partition
119 }
120 }
121 xDelete<int>(vertex_pairing);
122
123 /*cleanup and assign output pointer*/
124 xDelete<int>(epart);
125 *pmy_elements=my_elements;
126 *pmy_vertices=my_vertices;
127}
Note: See TracBrowser for help on using the repository browser.