1 static char help[] = "The main goal of this code is to retrieve the original element numbers as found in the " 2 "initial partitions (sInitialPartition)... but after the call to DMPlexDistribute"; 3 4 #include <petsc.h> 5 6 /* Coordinates of a 2x5 rectangular mesh of quads : */ 7 PetscReal sCoords2x5Mesh[18][2] = { 8 {0.00000000000000000e+00, 0.00000000000000000e+00}, 9 {2.00000000000000000e+00, 0.00000000000000000e+00}, 10 {0.00000000000000000e+00, 1.00000000000000000e+00}, 11 {2.00000000000000000e+00, 1.00000000000000000e+00}, 12 {9.99999999997387978e-01, 0.00000000000000000e+00}, 13 {9.99999999997387978e-01, 1.00000000000000000e+00}, 14 {0.00000000000000000e+00, 2.00000000000000011e-01}, 15 {0.00000000000000000e+00, 4.00000000000000022e-01}, 16 {0.00000000000000000e+00, 5.99999999999999978e-01}, 17 {0.00000000000000000e+00, 8.00000000000000044e-01}, 18 {2.00000000000000000e+00, 2.00000000000000011e-01}, 19 {2.00000000000000000e+00, 4.00000000000000022e-01}, 20 {2.00000000000000000e+00, 5.99999999999999978e-01}, 21 {2.00000000000000000e+00, 8.00000000000000044e-01}, 22 {9.99999999997387756e-01, 2.00000000000000011e-01}, 23 {9.99999999997387978e-01, 4.00000000000000022e-01}, 24 {9.99999999997387978e-01, 6.00000000000000089e-01}, 25 {9.99999999997388089e-01, 8.00000000000000044e-01} 26 }; 27 28 /* Connectivity of a 2x5 rectangular mesh of quads : */ 29 const PetscInt sConnectivity2x5Mesh[10][4] = { 30 {0, 4, 14, 6 }, 31 {6, 14, 15, 7 }, 32 {7, 15, 16, 8 }, 33 {8, 16, 17, 9 }, 34 {9, 17, 5, 2 }, 35 {4, 1, 10, 14}, 36 {14, 10, 11, 15}, 37 {15, 11, 12, 16}, 38 {16, 12, 13, 17}, 39 {17, 13, 3, 5 } 40 }; 41 42 /* Partitions of a 2x5 rectangular mesh of quads : */ 43 const PetscInt sInitialPartition2x5Mesh[2][5] = { 44 {0, 2, 4, 6, 8}, 45 {1, 3, 5, 7, 9} 46 }; 47 48 const PetscInt sNLoclCells2x5Mesh = 5; 49 const PetscInt sNGlobVerts2x5Mesh = 18; 50 51 /* Prisms mesh */ 52 PetscReal sCoordsPrismsMesh[125][3] = { 53 {2.24250931694056355e-01, 0.00000000000000000e+00, 0.00000000000000000e+00}, 54 {2.20660660151932697e-01, 2.87419338850266937e-01, 0.00000000000000000e+00}, 55 {0.00000000000000000e+00, 0.00000000000000000e+00, 2.70243537720639027e-01}, 56 {2.32445727460992402e-01, 0.00000000000000000e+00, 2.60591845015572310e-01}, 57 {2.41619971105419079e-01, 2.69894910706158231e-01, 2.42844781736072490e-01}, 58 {0.00000000000000000e+00, 2.46523339883120779e-01, 2.69072907562752262e-01}, 59 {0.00000000000000000e+00, 0.00000000000000000e+00, 0.00000000000000000e+00}, 60 {1.00000000000000000e+00, 2.75433417601945563e-01, 0.00000000000000000e+00}, 61 {1.00000000000000000e+00, 0.00000000000000000e+00, 2.33748605950385602e-01}, 62 {7.32445727460992457e-01, 0.00000000000000000e+00, 2.42344379130445597e-01}, 63 {1.00000000000000000e+00, 2.78258478013028610e-01, 2.57379172987105553e-01}, 64 {1.00000000000000000e+00, 0.00000000000000000e+00, 0.00000000000000000e+00}, 65 {7.49586880891153995e-01, 1.00000000000000000e+00, 0.00000000000000000e+00}, 66 {1.00000000000000000e+00, 1.00000000000000000e+00, 2.51949651266657582e-01}, 67 {7.41619971105419107e-01, 7.69894910706158120e-01, 2.33697838509081768e-01}, 68 {1.00000000000000000e+00, 7.78258478013028610e-01, 2.66479695645241543e-01}, 69 {7.55042653233710115e-01, 1.00000000000000000e+00, 2.58019637386860512e-01}, 70 {1.00000000000000000e+00, 1.00000000000000000e+00, 0.00000000000000000e+00}, 71 {0.00000000000000000e+00, 7.59235710423095789e-01, 0.00000000000000000e+00}, 72 {0.00000000000000000e+00, 1.00000000000000000e+00, 2.17232187874490473e-01}, 73 {0.00000000000000000e+00, 7.46523339883120807e-01, 2.42567232639677999e-01}, 74 {2.55042653233710115e-01, 1.00000000000000000e+00, 2.40660905690776916e-01}, 75 {0.00000000000000000e+00, 1.00000000000000000e+00, 0.00000000000000000e+00}, 76 {2.38934376044866809e-01, 0.00000000000000000e+00, 1.00000000000000000e+00}, 77 {2.18954188589218168e-01, 2.26916038449581692e-01, 1.00000000000000000e+00}, 78 {2.39787449636397643e-01, 0.00000000000000000e+00, 7.60591845015572310e-01}, 79 {2.40766735324061815e-01, 2.39643260505815608e-01, 7.42844781736072490e-01}, 80 {0.00000000000000000e+00, 2.57448248192627016e-01, 7.69072907562752262e-01}, 81 {0.00000000000000000e+00, 0.00000000000000000e+00, 1.00000000000000000e+00}, 82 {1.00000000000000000e+00, 2.38666970143638080e-01, 1.00000000000000000e+00}, 83 {7.39787449636397643e-01, 0.00000000000000000e+00, 7.42344379130445597e-01}, 84 {1.00000000000000000e+00, 2.59875254283874868e-01, 7.57379172987105553e-01}, 85 {1.00000000000000000e+00, 0.00000000000000000e+00, 1.00000000000000000e+00}, 86 {7.76318984007844159e-01, 1.00000000000000000e+00, 1.00000000000000000e+00}, 87 {7.40766735324061787e-01, 7.39643260505815636e-01, 7.33697838509081768e-01}, 88 {1.00000000000000000e+00, 7.59875254283874924e-01, 7.66479695645241543e-01}, 89 {7.68408704792055142e-01, 1.00000000000000000e+00, 7.58019637386860512e-01}, 90 {1.00000000000000000e+00, 1.00000000000000000e+00, 1.00000000000000000e+00}, 91 {0.00000000000000000e+00, 7.81085527042108207e-01, 1.00000000000000000e+00}, 92 {0.00000000000000000e+00, 7.57448248192627016e-01, 7.42567232639678054e-01}, 93 {2.68408704792055197e-01, 1.00000000000000000e+00, 7.40660905690776916e-01}, 94 {0.00000000000000000e+00, 1.00000000000000000e+00, 1.00000000000000000e+00}, 95 {7.24250931694056410e-01, 0.00000000000000000e+00, 0.00000000000000000e+00}, 96 {7.24250931694056410e-01, 2.75433417601945563e-01, 0.00000000000000000e+00}, 97 {4.44911591845989052e-01, 2.87419338850266937e-01, 0.00000000000000000e+00}, 98 {4.64891454921984804e-01, 0.00000000000000000e+00, 2.50940152310505593e-01}, 99 {4.74065698566411453e-01, 2.69894910706158231e-01, 2.33193089031005774e-01}, 100 {4.48501863388112709e-01, 0.00000000000000000e+00, 0.00000000000000000e+00}, 101 {7.20660660151932753e-01, 7.87419338850266937e-01, 0.00000000000000000e+00}, 102 {7.20660660151932753e-01, 5.62852756452212555e-01, 0.00000000000000000e+00}, 103 {2.20660660151932697e-01, 5.46655049273362614e-01, 0.00000000000000000e+00}, 104 {4.83239942210838158e-01, 5.39789821412316462e-01, 2.15446025751505982e-01}, 105 {7.41619971105419107e-01, 5.48153388719186951e-01, 2.48227882887665785e-01}, 106 {2.41619971105419079e-01, 5.16418250589278927e-01, 2.41674151578185781e-01}, 107 {4.41321320303865394e-01, 5.74838677700533873e-01, 0.00000000000000000e+00}, 108 {1.00000000000000000e+00, 7.75433417601945507e-01, 0.00000000000000000e+00}, 109 {1.00000000000000000e+00, 5.56516956026057219e-01, 2.81009740023825560e-01}, 110 {7.32445727460992457e-01, 2.78258478013028610e-01, 2.65974946167165549e-01}, 111 {1.00000000000000000e+00, 5.50866835203891125e-01, 0.00000000000000000e+00}, 112 {0.00000000000000000e+00, 2.59235710423095733e-01, 0.00000000000000000e+00}, 113 {0.00000000000000000e+00, 4.93046679766241558e-01, 2.67902277404865552e-01}, 114 {2.55042653233710115e-01, 7.46523339883120807e-01, 2.65995950455964469e-01}, 115 {0.00000000000000000e+00, 5.18471420846191466e-01, 0.00000000000000000e+00}, 116 {2.49586880891154023e-01, 1.00000000000000000e+00, 0.00000000000000000e+00}, 117 {2.49586880891154023e-01, 7.59235710423095789e-01, 0.00000000000000000e+00}, 118 {4.70247541043086748e-01, 7.87419338850266937e-01, 0.00000000000000000e+00}, 119 {5.10085306467420230e-01, 1.00000000000000000e+00, 2.64089623507063387e-01}, 120 {4.96662624339129222e-01, 7.69894910706158231e-01, 2.39767824629284698e-01}, 121 {4.99173761782308045e-01, 1.00000000000000000e+00, 0.00000000000000000e+00}, 122 {0.00000000000000000e+00, 0.00000000000000000e+00, 7.70243537720639027e-01}, 123 {2.40640523227928449e-01, 0.00000000000000000e+00, 5.21183690031144620e-01}, 124 {2.62579282058905461e-01, 2.52370482562049525e-01, 4.85689563472144981e-01}, 125 {0.00000000000000000e+00, 2.33810969343145825e-01, 5.38145815125504523e-01}, 126 {0.00000000000000000e+00, 0.00000000000000000e+00, 5.40487075441278053e-01}, 127 {1.00000000000000000e+00, 0.00000000000000000e+00, 7.33748605950385602e-01}, 128 {7.40640523227928504e-01, 0.00000000000000000e+00, 4.84688758260891195e-01}, 129 {1.00000000000000000e+00, 2.81083538424111656e-01, 5.14758345974211107e-01}, 130 {1.00000000000000000e+00, 0.00000000000000000e+00, 4.67497211900771203e-01}, 131 {7.38934376044866781e-01, 0.00000000000000000e+00, 1.00000000000000000e+00}, 132 {4.79574899272795285e-01, 0.00000000000000000e+00, 7.50940152310505593e-01}, 133 {4.77868752089733617e-01, 0.00000000000000000e+00, 1.00000000000000000e+00}, 134 {1.00000000000000000e+00, 1.00000000000000000e+00, 7.51949651266657582e-01}, 135 {7.62579282058905461e-01, 7.52370482562049525e-01, 4.67395677018163536e-01}, 136 {1.00000000000000000e+00, 7.81083538424111712e-01, 5.32959391290483087e-01}, 137 {7.60498425576266124e-01, 1.00000000000000000e+00, 5.16039274773721024e-01}, 138 {1.00000000000000000e+00, 1.00000000000000000e+00, 5.03899302533315163e-01}, 139 {7.18954188589218113e-01, 7.26916038449581636e-01, 1.00000000000000000e+00}, 140 {4.81533470648123629e-01, 4.79286521011631217e-01, 7.15446025751505954e-01}, 141 {4.57888564634085005e-01, 2.26916038449581692e-01, 1.00000000000000000e+00}, 142 {4.95273172597062383e-01, 7.26916038449581636e-01, 1.00000000000000000e+00}, 143 {4.37908377178436337e-01, 4.53832076899163384e-01, 1.00000000000000000e+00}, 144 {1.00000000000000000e+00, 7.38666970143638135e-01, 1.00000000000000000e+00}, 145 {1.00000000000000000e+00, 5.19750508567749736e-01, 7.81009740023825616e-01}, 146 {7.38934376044866781e-01, 2.38666970143638080e-01, 1.00000000000000000e+00}, 147 {7.18954188589218113e-01, 4.65583008593219771e-01, 1.00000000000000000e+00}, 148 {1.00000000000000000e+00, 4.77333940287276159e-01, 1.00000000000000000e+00}, 149 {0.00000000000000000e+00, 1.00000000000000000e+00, 7.17232187874490501e-01}, 150 {0.00000000000000000e+00, 7.33810969343145825e-01, 4.85134465279355998e-01}, 151 {2.60498425576266179e-01, 1.00000000000000000e+00, 4.81321811381553832e-01}, 152 {0.00000000000000000e+00, 1.00000000000000000e+00, 4.34464375748980947e-01}, 153 {0.00000000000000000e+00, 2.81085527042108152e-01, 1.00000000000000000e+00}, 154 {0.00000000000000000e+00, 5.14896496385254032e-01, 7.67902277404865607e-01}, 155 {2.76318984007844215e-01, 7.81085527042108207e-01, 1.00000000000000000e+00}, 156 {2.18954188589218168e-01, 5.08001565491689844e-01, 1.00000000000000000e+00}, 157 {0.00000000000000000e+00, 5.62171054084216304e-01, 1.00000000000000000e+00}, 158 {2.76318984007844215e-01, 1.00000000000000000e+00, 1.00000000000000000e+00}, 159 {5.36817409584110394e-01, 1.00000000000000000e+00, 7.64089623507063331e-01}, 160 {5.52637968015688430e-01, 1.00000000000000000e+00, 1.00000000000000000e+00}, 161 {5.03219805286833965e-01, 2.52370482562049525e-01, 4.66386178062011547e-01}, 162 {4.80554184960459430e-01, 2.39643260505815608e-01, 7.33193089031005774e-01}, 163 {4.81281046455856898e-01, 0.00000000000000000e+00, 5.01880304621011186e-01}, 164 {7.62579282058905461e-01, 5.33454020986161126e-01, 4.96455765775331570e-01}, 165 {2.62579282058905461e-01, 4.86181451905195350e-01, 4.83348303156371562e-01}, 166 {7.40766735324061787e-01, 4.99518514789690449e-01, 7.48227882887665841e-01}, 167 {2.40766735324061815e-01, 4.97091508698442541e-01, 7.41674151578185725e-01}, 168 {5.25158564117810922e-01, 5.04740965124099050e-01, 4.30892051503011964e-01}, 169 {7.40640523227928504e-01, 2.81083538424111656e-01, 5.31949892334331098e-01}, 170 {7.39787449636397643e-01, 2.59875254283874868e-01, 7.65974946167165549e-01}, 171 {1.00000000000000000e+00, 5.62167076848223313e-01, 5.62019480047651121e-01}, 172 {2.60498425576266179e-01, 7.33810969343145825e-01, 5.31991900911928939e-01}, 173 {2.68408704792055197e-01, 7.57448248192627016e-01, 7.65995950455964469e-01}, 174 {0.00000000000000000e+00, 4.67621938686291649e-01, 5.35804554809731104e-01}, 175 {5.23077707635171585e-01, 7.52370482562049525e-01, 4.79535649258569396e-01}, 176 {5.09175440116116929e-01, 7.39643260505815636e-01, 7.39767824629284698e-01}, 177 {5.20996851152532359e-01, 1.00000000000000000e+00, 5.28179247014126774e-01} 178 }; 179 180 const PetscInt sConnectivityPrismsMesh[128][6] = { 181 /* rank 0 */ 182 {11, 7, 42, 8, 10, 9 }, 183 {47, 42, 43, 45, 9, 57 }, 184 {8, 10, 9, 77, 76, 75 }, 185 {45, 9, 57, 110, 75, 116}, 186 {17, 48, 55, 13, 14, 15 }, 187 {58, 55, 49, 56, 15, 52 }, 188 {13, 14, 15, 85, 82, 83 }, 189 {56, 15, 52, 118, 83, 111}, 190 {6, 0, 1, 2, 3, 4 }, 191 {54, 1, 44, 51, 4, 46 }, 192 {2, 3, 4, 73, 70, 71 }, 193 {51, 4, 46, 115, 71, 108}, 194 {58, 49, 43, 56, 52, 57 }, 195 {47, 43, 44, 45, 57, 46 }, 196 {56, 52, 57, 118, 111, 116}, 197 {45, 57, 46, 110, 116, 108}, 198 {77, 76, 75, 74, 31, 30 }, 199 {110, 75, 116, 79, 30, 117}, 200 {74, 31, 30, 32, 29, 78 }, 201 {79, 30, 117, 80, 78, 93 }, 202 {85, 82, 83, 81, 34, 35 }, 203 {118, 83, 111, 92, 35, 113}, 204 {81, 34, 35, 37, 86, 91 }, 205 {92, 35, 113, 95, 91, 94 }, 206 {73, 70, 71, 69, 25, 26 }, 207 {115, 71, 108, 87, 26, 109}, 208 {69, 25, 26, 28, 23, 24 }, 209 {87, 26, 109, 90, 24, 88 }, 210 {118, 111, 116, 92, 113, 117}, 211 {110, 116, 108, 79, 117, 109}, 212 {92, 113, 117, 95, 94, 93 }, 213 {79, 117, 109, 80, 93, 88 }, 214 {22, 18, 63, 19, 20, 21 }, 215 {68, 63, 64, 66, 21, 61 }, 216 {19, 20, 21, 99, 97, 98 }, 217 {66, 21, 61, 124, 98, 119}, 218 {6, 1, 59, 2, 4, 5 }, 219 {62, 59, 50, 60, 5, 53 }, 220 {2, 4, 5, 73, 71, 72 }, 221 {60, 5, 53, 121, 72, 112}, 222 {17, 12, 48, 13, 16, 14 }, 223 {54, 48, 65, 51, 14, 67 }, 224 {13, 16, 14, 85, 84, 82 }, 225 {51, 14, 67, 115, 82, 122}, 226 {62, 50, 64, 60, 53, 61 }, 227 {68, 64, 65, 66, 61, 67 }, 228 {60, 53, 61, 121, 112, 119}, 229 {66, 61, 67, 124, 119, 122}, 230 {99, 97, 98, 96, 39, 40 }, 231 {124, 98, 119, 106, 40, 120}, 232 {96, 39, 40, 41, 38, 105}, 233 {106, 40, 120, 107, 105, 102}, 234 {73, 71, 72, 69, 26, 27 }, 235 {121, 72, 112, 101, 27, 114}, 236 {69, 26, 27, 28, 24, 100}, 237 {101, 27, 114, 104, 100, 103}, 238 {85, 84, 82, 81, 36, 34 }, 239 {115, 82, 122, 87, 34, 123}, 240 {81, 36, 34, 37, 33, 86 }, 241 {87, 34, 123, 90, 86, 89 }, 242 {121, 112, 119, 101, 114, 120}, 243 {124, 119, 122, 106, 120, 123}, 244 {101, 114, 120, 104, 103, 102}, 245 {106, 120, 123, 107, 102, 89 }, 246 /* rank 1 */ 247 {58, 43, 7, 56, 57, 10 }, 248 {7, 43, 42, 10, 57, 9 }, 249 {56, 57, 10, 118, 116, 76 }, 250 {10, 57, 9, 76, 116, 75 }, 251 {54, 49, 48, 51, 52, 14 }, 252 {48, 49, 55, 14, 52, 15 }, 253 {51, 52, 14, 115, 111, 82 }, 254 {14, 52, 15, 82, 111, 83 }, 255 {47, 44, 0, 45, 46, 3 }, 256 {0, 44, 1, 3, 46, 4 }, 257 {45, 46, 3, 110, 108, 70 }, 258 {3, 46, 4, 70, 108, 71 }, 259 {54, 44, 49, 51, 46, 52 }, 260 {49, 44, 43, 52, 46, 57 }, 261 {51, 46, 52, 115, 108, 111}, 262 {52, 46, 57, 111, 108, 116}, 263 {118, 116, 76, 92, 117, 31 }, 264 {76, 116, 75, 31, 117, 30 }, 265 {92, 117, 31, 95, 93, 29 }, 266 {31, 117, 30, 29, 93, 78 }, 267 {115, 111, 82, 87, 113, 34 }, 268 {82, 111, 83, 34, 113, 35 }, 269 {87, 113, 34, 90, 94, 86 }, 270 {34, 113, 35, 86, 94, 91 }, 271 {110, 108, 70, 79, 109, 25 }, 272 {70, 108, 71, 25, 109, 26 }, 273 {79, 109, 25, 80, 88, 23 }, 274 {25, 109, 26, 23, 88, 24 }, 275 {115, 108, 111, 87, 109, 113}, 276 {111, 108, 116, 113, 109, 117}, 277 {87, 109, 113, 90, 88, 94 }, 278 {113, 109, 117, 94, 88, 93 }, 279 {62, 64, 18, 60, 61, 20 }, 280 {18, 64, 63, 20, 61, 21 }, 281 {60, 61, 20, 121, 119, 97 }, 282 {20, 61, 21, 97, 119, 98 }, 283 {54, 50, 1, 51, 53, 4 }, 284 {1, 50, 59, 4, 53, 5 }, 285 {51, 53, 4, 115, 112, 71 }, 286 {4, 53, 5, 71, 112, 72 }, 287 {68, 65, 12, 66, 67, 16 }, 288 {12, 65, 48, 16, 67, 14 }, 289 {66, 67, 16, 124, 122, 84 }, 290 {16, 67, 14, 84, 122, 82 }, 291 {54, 65, 50, 51, 67, 53 }, 292 {50, 65, 64, 53, 67, 61 }, 293 {51, 67, 53, 115, 122, 112}, 294 {53, 67, 61, 112, 122, 119}, 295 {121, 119, 97, 101, 120, 39 }, 296 {97, 119, 98, 39, 120, 40 }, 297 {101, 120, 39, 104, 102, 38 }, 298 {39, 120, 40, 38, 102, 105}, 299 {115, 112, 71, 87, 114, 26 }, 300 {71, 112, 72, 26, 114, 27 }, 301 {87, 114, 26, 90, 103, 24 }, 302 {26, 114, 27, 24, 103, 100}, 303 {124, 122, 84, 106, 123, 36 }, 304 {84, 122, 82, 36, 123, 34 }, 305 {106, 123, 36, 107, 89, 33 }, 306 {36, 123, 34, 33, 89, 86 }, 307 {115, 122, 112, 87, 123, 114}, 308 {112, 122, 119, 114, 123, 120}, 309 {87, 123, 114, 90, 89, 103}, 310 {114, 123, 120, 103, 89, 102} 311 }; 312 313 /* Partitions of prisms mesh : */ 314 const PetscInt sInitialPartitionPrismsMesh[2][64] = { 315 {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63 }, 316 {64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 317 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127} 318 }; 319 320 const PetscInt sNLoclCellsPrismsMesh = 64; 321 const PetscInt sNGlobVertsPrismsMesh = 125; 322 323 int main(int argc, char **argv) 324 { 325 PetscInt Nc = 0; 326 const PetscInt *InitPartForRank[2]; 327 DM dm, idm, ddm; 328 PetscSF sfVert, sfMig, sfPart; 329 PetscPartitioner part; 330 PetscSection s; 331 PetscInt *cells, c; 332 PetscMPIInt size, rank; 333 PetscBool box = PETSC_FALSE, field = PETSC_FALSE, quadsmesh = PETSC_FALSE, prismsmesh = PETSC_FALSE; 334 335 PetscFunctionBeginUser; 336 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 337 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 338 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 339 PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_SUP, "This is a 2 processors example only"); 340 PetscCall(PetscOptionsGetBool(NULL, NULL, "-field", &field, NULL)); 341 PetscCall(PetscOptionsGetBool(NULL, NULL, "-box", &box, NULL)); 342 PetscCall(PetscOptionsGetBool(NULL, NULL, "-quadsmesh", &quadsmesh, NULL)); 343 PetscCall(PetscOptionsGetBool(NULL, NULL, "-prismsmesh", &prismsmesh, NULL)); 344 PetscCheck(1 == (box ? 1 : 0) + (quadsmesh ? 1 : 0) + (prismsmesh ? 1 : 0), PETSC_COMM_WORLD, PETSC_ERR_SUP, "Specify one and only one of -box, -quadsmesh or -prismsmesh"); 345 346 PetscCall(DMPlexCreate(PETSC_COMM_WORLD, &dm)); 347 if (box) { 348 PetscCall(DMSetType(dm, DMPLEX)); 349 PetscCall(DMSetFromOptions(dm)); 350 } else { 351 if (quadsmesh) { 352 Nc = sNLoclCells2x5Mesh; //Same on each rank for this example... 353 PetscInt Nv = sNGlobVerts2x5Mesh; 354 InitPartForRank[0] = &sInitialPartition2x5Mesh[0][0]; 355 InitPartForRank[1] = &sInitialPartition2x5Mesh[1][0]; 356 const PetscInt(*Conn)[4] = sConnectivity2x5Mesh; 357 358 const PetscInt Ncor = 4; 359 const PetscInt dim = 2; 360 361 PetscCall(PetscMalloc1(Nc * Ncor, &cells)); 362 for (c = 0; c < Nc; ++c) { 363 PetscInt cell = (InitPartForRank[rank])[c], cor; 364 365 for (cor = 0; cor < Ncor; ++cor) cells[c * Ncor + cor] = Conn[cell][cor]; 366 } 367 PetscCall(DMSetDimension(dm, dim)); 368 PetscCall(DMPlexBuildFromCellListParallel(dm, Nc, PETSC_DECIDE, Nv, Ncor, cells, &sfVert, NULL)); 369 } else if (prismsmesh) { 370 Nc = sNLoclCellsPrismsMesh; //Same on each rank for this example... 371 PetscInt Nv = sNGlobVertsPrismsMesh; 372 InitPartForRank[0] = &sInitialPartitionPrismsMesh[0][0]; 373 InitPartForRank[1] = &sInitialPartitionPrismsMesh[1][0]; 374 const PetscInt(*Conn)[6] = sConnectivityPrismsMesh; 375 376 const PetscInt Ncor = 6; 377 const PetscInt dim = 3; 378 379 PetscCall(PetscMalloc1(Nc * Ncor, &cells)); 380 for (c = 0; c < Nc; ++c) { 381 PetscInt cell = (InitPartForRank[rank])[c], cor; 382 383 for (cor = 0; cor < Ncor; ++cor) cells[c * Ncor + cor] = Conn[cell][cor]; 384 } 385 PetscCall(DMSetDimension(dm, dim)); 386 PetscCall(DMPlexBuildFromCellListParallel(dm, Nc, PETSC_DECIDE, Nv, Ncor, cells, &sfVert, NULL)); 387 } 388 PetscCall(PetscSFDestroy(&sfVert)); 389 PetscCall(PetscFree(cells)); 390 PetscCall(DMPlexInterpolate(dm, &idm)); 391 PetscCall(DMDestroy(&dm)); 392 dm = idm; 393 } 394 PetscCall(DMSetUseNatural(dm, PETSC_TRUE)); 395 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 396 397 if (field) { 398 const PetscInt Nf = 1; 399 const PetscInt numBC = 0; 400 const PetscInt numComp[1] = {1}; 401 PetscInt numDof[4] = {0, 0, 0, 0}; 402 PetscInt dim; 403 404 PetscCall(DMGetDimension(dm, &dim)); 405 numDof[dim] = 1; 406 407 PetscCall(DMSetNumFields(dm, Nf)); 408 PetscCall(DMPlexCreateSection(dm, NULL, numComp, numDof, numBC, NULL, NULL, NULL, NULL, &s)); 409 PetscCall(DMSetLocalSection(dm, s)); 410 PetscCall(PetscSectionView(s, PETSC_VIEWER_STDOUT_WORLD)); 411 PetscCall(PetscSectionDestroy(&s)); 412 } 413 414 PetscCall(DMPlexGetPartitioner(dm, &part)); 415 PetscCall(PetscPartitionerSetFromOptions(part)); 416 417 PetscCall(DMPlexDistribute(dm, 0, &sfMig, &ddm)); 418 PetscCall(PetscSFView(sfMig, PETSC_VIEWER_STDOUT_WORLD)); 419 PetscCall(PetscSFCreateInverseSF(sfMig, &sfPart)); 420 PetscCall(PetscObjectSetName((PetscObject)sfPart, "Inverse Migration SF")); 421 PetscCall(PetscSFView(sfPart, PETSC_VIEWER_STDOUT_WORLD)); 422 423 Vec lGlobalVec, lNatVec; 424 PetscScalar *lNatVecArray; 425 426 { 427 PetscSection s; 428 429 PetscCall(DMGetGlobalSection(dm, &s)); 430 PetscCall(PetscSectionView(s, PETSC_VIEWER_STDOUT_WORLD)); 431 } 432 PetscCall(DMGetGlobalVector(dm, &lNatVec)); 433 PetscCall(PetscObjectSetName((PetscObject)lNatVec, "Natural Vector (initial partition)")); 434 435 //Copying the initial partition into the "natural" vector: 436 PetscCall(VecZeroEntries(lNatVec)); 437 PetscCall(VecGetArray(lNatVec, &lNatVecArray)); 438 for (c = 0; c < Nc; ++c) lNatVecArray[c] = (InitPartForRank[rank])[c]; 439 PetscCall(VecRestoreArray(lNatVec, &lNatVecArray)); 440 441 PetscCall(DMGetGlobalVector(ddm, &lGlobalVec)); 442 PetscCall(PetscObjectSetName((PetscObject)lGlobalVec, "Global Vector (reordered element numbers in the petsc distributed order)")); 443 PetscCall(VecZeroEntries(lGlobalVec)); 444 445 // The call to DMPlexNaturalToGlobalBegin/End does not produce our expected result... 446 // In lGlobalVec, we expect to have: 447 /* 448 * Process [0] 449 * 2. 450 * 4. 451 * 8. 452 * 3. 453 * 9. 454 * Process [1] 455 * 1. 456 * 5. 457 * 7. 458 * 0. 459 * 6. 460 * 461 * but we obtained: 462 * 463 * Process [0] 464 * 2. 465 * 4. 466 * 8. 467 * 0. 468 * 0. 469 * Process [1] 470 * 0. 471 * 0. 472 * 0. 473 * 0. 474 * 0. 475 */ 476 477 { 478 PetscSF nsf; 479 480 PetscCall(DMPlexGetGlobalToNaturalSF(ddm, &nsf)); 481 PetscCall(PetscSFView(nsf, NULL)); 482 } 483 PetscCall(DMPlexNaturalToGlobalBegin(ddm, lNatVec, lGlobalVec)); 484 PetscCall(DMPlexNaturalToGlobalEnd(ddm, lNatVec, lGlobalVec)); 485 486 PetscCall(VecView(lNatVec, PETSC_VIEWER_STDOUT_WORLD)); 487 PetscCall(VecView(lGlobalVec, PETSC_VIEWER_STDOUT_WORLD)); 488 489 PetscCall(DMRestoreGlobalVector(dm, &lNatVec)); 490 PetscCall(DMRestoreGlobalVector(ddm, &lGlobalVec)); 491 492 const PetscBool lUseCone = PETSC_FALSE; 493 const PetscBool lUseClosure = PETSC_TRUE; 494 PetscCall(DMSetBasicAdjacency(ddm, lUseCone, lUseClosure)); 495 const PetscInt lNbCellsInOverlap = 1; 496 PetscSF lSFMigrationOvl; 497 DM ddm_with_overlap; 498 499 PetscCall(DMPlexDistributeOverlap(ddm, lNbCellsInOverlap, &lSFMigrationOvl, &ddm_with_overlap)); 500 501 IS lISCellWithOvl = 0; 502 /* This is the buggy call with prisms since commit 5ae96e2b862 */ 503 PetscCall(DMPlexCreateCellNumbering(ddm_with_overlap, PETSC_TRUE, &lISCellWithOvl)); 504 /* Here, we can see the elements in the overlap within the IS: they are the ones with negative indices */ 505 PetscCall(ISView(lISCellWithOvl, PETSC_VIEWER_STDOUT_WORLD)); 506 PetscCall(ISDestroy(&lISCellWithOvl)); 507 508 PetscCall(PetscSFDestroy(&lSFMigrationOvl)); 509 PetscCall(DMDestroy(&ddm_with_overlap)); 510 PetscCall(PetscSFDestroy(&sfMig)); 511 PetscCall(PetscSFDestroy(&sfPart)); 512 PetscCall(DMDestroy(&dm)); 513 PetscCall(DMDestroy(&ddm)); 514 PetscCall(PetscFinalize()); 515 return 0; 516 } 517 518 /*TEST 519 520 testset: 521 args: -field -petscpartitioner_type simple 522 nsize: 2 523 524 test: 525 suffix: 0 526 args: -quadsmesh 527 output_file: output/ex47_0.out 528 529 test: 530 suffix: 1 531 args: -box -dm_plex_simplex 0 -dm_plex_box_faces 2,5 -dm_distribute 532 output_file: output/ex47_1.out 533 534 test: 535 suffix: 2 536 args: -prismsmesh 537 output_file: output/ex47_2.out 538 539 TEST*/ 540