1 /******************************************************************************
2 * Copyright 1998-2019 Lawrence Livermore National Security, LLC and other
3 * HYPRE Project Developers. See the top-level COPYRIGHT file for details.
4 *
5 * SPDX-License-Identifier: (Apache-2.0 OR MIT)
6 ******************************************************************************/
7
8 /******************************************************************************
9 * OpenMP Problems
10 *
11 * Need to fix the way these variables are set and incremented in loops:
12 * cnt
13 *
14 ******************************************************************************/
15
16 #include "_hypre_sstruct_ls.h"
17
18 /*--------------------------------------------------------------------------
19 * Finds the physical boundary boxes for all levels. Since the coarse grid's
20 * boundary may not be on the physical bdry, we need to compare the coarse
21 * grid to the finest level boundary boxes. All boxes of the coarse grids
22 * must be checked, not just the bounding box.
23 * Algo:
24 * 1) obtain boundary boxes for the finest grid
25 * i) mark the fboxes that have boundary elements.
26 * 2) loop over coarse levels
27 * i) for a cbox that maps to a fbox that has boundary layers
28 * a) refine the cbox
29 * b) intersect with the cell boundary layers of the fbox
30 * c) coarsen the intersection
31 * ii) determine the var boxes
32 * iii) mark the coarse box
33 *
34 * Concerns: Checking an individual pgrid may give artificial physical
35 * boundaries. Need to check if any other pgrid is adjacent to it.
36 * We omit this case and assume only one part for now.
37 *--------------------------------------------------------------------------*/
38
39 HYPRE_Int
hypre_Maxwell_PhysBdy(hypre_SStructGrid ** grid_l,HYPRE_Int num_levels,hypre_Index rfactors,HYPRE_Int *** BdryRanksl_ptr,HYPRE_Int ** BdryRanksCntsl_ptr)40 hypre_Maxwell_PhysBdy( hypre_SStructGrid **grid_l,
41 HYPRE_Int num_levels,
42 hypre_Index rfactors,
43 HYPRE_Int ***BdryRanksl_ptr,
44 HYPRE_Int **BdryRanksCntsl_ptr )
45 {
46
47 MPI_Comm comm= (grid_l[0]-> comm);
48
49 HYPRE_Int **BdryRanks_l;
50 HYPRE_Int *BdryRanksCnts_l;
51
52 HYPRE_Int *npts;
53 HYPRE_BigInt *ranks, *upper_rank, *lower_rank;
54 hypre_BoxManEntry *boxman_entry;
55
56 hypre_SStructGrid *grid;
57 hypre_SStructPGrid *pgrid;
58 hypre_StructGrid *cell_fgrid, *cell_cgrid, *sgrid;
59
60 hypre_BoxArrayArray ****bdry;
61 hypre_BoxArrayArray *fbdry;
62 hypre_BoxArrayArray *cbdry;
63
64 hypre_BoxArray *box_array;
65 hypre_BoxArray *fboxes, *cboxes;
66
67 hypre_Box *fbox, *cbox;
68 hypre_Box *box, *contract_fbox, rbox;
69 hypre_Box intersect;
70
71 HYPRE_Int **cbox_mapping, **fbox_mapping;
72 HYPRE_Int **boxes_with_bdry;
73
74 HYPRE_Int ndim, nvars;
75 HYPRE_Int nboxes, nfboxes;
76 HYPRE_Int boxi;
77
78 hypre_Index zero_shift, upper_shift, lower_shift;
79 hypre_Index loop_size, start, index, lindex;
80
81 HYPRE_Int i, j, k, l, m, n, p;
82 HYPRE_Int d;
83 HYPRE_Int cnt;
84
85 HYPRE_Int part= 0; /* NOTE, ASSUMING ONE PART */
86 HYPRE_Int matrix_type= HYPRE_PARCSR;
87 HYPRE_Int myproc;
88
89 HYPRE_Int ierr= 0;
90
91 hypre_MPI_Comm_rank(comm, &myproc);
92
93 ndim= hypre_SStructGridNDim(grid_l[0]);
94 hypre_SetIndex3(zero_shift, 0, 0, 0);
95
96 hypre_BoxInit(&intersect, ndim);
97
98 /* bounding global ranks of this processor & allocate boundary box markers. */
99 upper_rank= hypre_CTAlloc(HYPRE_BigInt, num_levels, HYPRE_MEMORY_HOST);
100 lower_rank= hypre_CTAlloc(HYPRE_BigInt, num_levels, HYPRE_MEMORY_HOST);
101
102 boxes_with_bdry= hypre_TAlloc(HYPRE_Int *, num_levels, HYPRE_MEMORY_HOST);
103 for (i= 0; i< num_levels; i++)
104 {
105 grid = grid_l[i];
106 lower_rank[i]= hypre_SStructGridStartRank(grid);
107
108 /* note we are assuming only one part */
109 pgrid= hypre_SStructGridPGrid(grid, part);
110 nvars= hypre_SStructPGridNVars(pgrid);
111 sgrid= hypre_SStructPGridSGrid(pgrid, nvars-1);
112 box_array= hypre_StructGridBoxes(sgrid);
113 box = hypre_BoxArrayBox(box_array, hypre_BoxArraySize(box_array)-1);
114
115 hypre_SStructGridBoxProcFindBoxManEntry(grid, part, nvars-1,
116 hypre_BoxArraySize(box_array)-1, myproc, &boxman_entry);
117 hypre_SStructBoxManEntryGetGlobalCSRank(boxman_entry, hypre_BoxIMax(box),
118 &upper_rank[i]);
119
120 sgrid= hypre_SStructPGridCellSGrid(pgrid);
121 box_array= hypre_StructGridBoxes(sgrid);
122 boxes_with_bdry[i]= hypre_CTAlloc(HYPRE_Int, hypre_BoxArraySize(box_array), HYPRE_MEMORY_HOST);
123 }
124
125 /*-----------------------------------------------------------------------------
126 * construct box_number mapping between levels, and offset strides because of
127 * projection coarsening. Note: from the way the coarse boxes are created and
128 * numbered, to determine the coarse box that matches the fbox, we need to
129 * only check the tail end of the list of cboxes. In fact, given fbox_i,
130 * if it's coarsened extents do not interesect with the first coarse box of the
131 * tail end, then this fbox vanishes in the coarsening.
132 * c/fbox_mapping gives the fine/coarse box mapping between two consecutive levels
133 * of the multilevel hierarchy.
134 *-----------------------------------------------------------------------------*/
135 if (num_levels > 1)
136 {
137 cbox_mapping= hypre_CTAlloc(HYPRE_Int *, num_levels, HYPRE_MEMORY_HOST);
138 fbox_mapping= hypre_CTAlloc(HYPRE_Int *, num_levels, HYPRE_MEMORY_HOST);
139 }
140 for (i= 0; i< (num_levels-1); i++)
141 {
142 grid = grid_l[i];
143 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
144 cell_fgrid= hypre_SStructPGridCellSGrid(pgrid);
145 fboxes= hypre_StructGridBoxes(cell_fgrid);
146 nfboxes= hypre_BoxArraySize(hypre_StructGridBoxes(cell_fgrid));
147 fbox_mapping[i]= hypre_CTAlloc(HYPRE_Int, nfboxes, HYPRE_MEMORY_HOST);
148
149 grid = grid_l[i+1];
150 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
151 cell_cgrid= hypre_SStructPGridCellSGrid(pgrid);
152 cboxes= hypre_StructGridBoxes(cell_cgrid);
153 nboxes= hypre_BoxArraySize(hypre_StructGridBoxes(cell_cgrid));
154
155 cbox_mapping[i+1]= hypre_CTAlloc(HYPRE_Int, nboxes, HYPRE_MEMORY_HOST);
156
157 /* assuming if i1 > i2 and (box j1) is coarsened from (box i1)
158 and (box j2) from (box i2), then j1 > j2. */
159 k= 0;
160 hypre_ForBoxI(j, fboxes)
161 {
162 fbox= hypre_BoxArrayBox(fboxes, j);
163 hypre_CopyBox(fbox, &rbox);
164 hypre_ProjectBox(&rbox, zero_shift, rfactors);
165 hypre_StructMapFineToCoarse(hypre_BoxIMin(&rbox), zero_shift,
166 rfactors, hypre_BoxIMin(&rbox));
167 hypre_StructMapFineToCoarse(hypre_BoxIMax(&rbox), zero_shift,
168 rfactors, hypre_BoxIMax(&rbox));
169
170 /* since the ordering of the cboxes was determined by the fbox
171 ordering, we only have to check if the first cbox in the
172 list intersects with rbox. If not, this fbox vanished in the
173 coarsening. */
174 cbox= hypre_BoxArrayBox(cboxes, k);
175 hypre_IntersectBoxes(&rbox, cbox, &rbox);
176 if (hypre_BoxVolume(&rbox))
177 {
178 cbox_mapping[i+1][k]= j;
179 fbox_mapping[i][j]= k;
180 k++;
181 } /* if (hypre_BoxVolume(&rbox)) */
182 } /* hypre_ForBoxI(j, fboxes) */
183 } /* for (i= 0; i< (num_levels-1); i++) */
184
185 bdry= hypre_TAlloc(hypre_BoxArrayArray ***, num_levels, HYPRE_MEMORY_HOST);
186 npts= hypre_CTAlloc(HYPRE_Int, num_levels, HYPRE_MEMORY_HOST);
187
188 /* finest level boundary determination */
189 grid = grid_l[0];
190 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
191 nvars= hypre_SStructPGridNVars(pgrid);
192 cell_fgrid= hypre_SStructPGridCellSGrid(pgrid);
193 nboxes= hypre_BoxArraySize(hypre_StructGridBoxes(cell_fgrid));
194
195 hypre_Maxwell_PNedelec_Bdy(cell_fgrid, pgrid, &bdry[0]);
196 for (i= 0; i< nboxes; i++)
197 {
198 if (bdry[0][i]) /* boundary layers on box[i] */
199 {
200 for (j= 0; j< nvars; j++)
201 {
202 fbdry= bdry[0][i][j+1]; /*(j+1) since j= 0 stores cell-centred boxes*/
203 hypre_ForBoxArrayI(k, fbdry)
204 {
205 box_array= hypre_BoxArrayArrayBoxArray(fbdry, k);
206 hypre_ForBoxI(p, box_array)
207 {
208 box= hypre_BoxArrayBox(box_array, p);
209 npts[0]+= hypre_BoxVolume(box);
210 }
211 }
212 } /* for (j= 0; j< nvars; j++) */
213
214 boxes_with_bdry[0][i]= 1; /* mark this box as containing boundary layers */
215 } /* if (bdry[0][i]) */
216 }
217 nfboxes= nboxes;
218
219 /* coarser levels */
220 for (i= 1; i< num_levels; i++)
221 {
222 grid = grid_l[i-1];
223 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
224 cell_fgrid= hypre_SStructPGridCellSGrid(pgrid);
225 fboxes= hypre_StructGridBoxes(cell_fgrid);
226
227 grid = grid_l[i];
228 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
229 cell_cgrid= hypre_SStructPGridCellSGrid(pgrid);
230 nvars= hypre_SStructPGridNVars(pgrid);
231 cboxes= hypre_StructGridBoxes(cell_cgrid);
232 nboxes= hypre_BoxArraySize(hypre_StructGridBoxes(cell_cgrid));
233
234 bdry[i]= hypre_TAlloc(hypre_BoxArrayArray **, nboxes, HYPRE_MEMORY_HOST);
235 p= 2*(ndim-1);
236 for (j= 0; j< nboxes; j++)
237 {
238 bdry[i][j]= hypre_TAlloc(hypre_BoxArrayArray *, nvars+1, HYPRE_MEMORY_HOST);
239
240 /* cell grid boxarrayarray */
241 bdry[i][j][0]= hypre_BoxArrayArrayCreate(2*ndim, ndim);
242
243 /* var grid boxarrayarrays */
244 for (k= 0; k< nvars; k++)
245 {
246 bdry[i][j][k+1]= hypre_BoxArrayArrayCreate(p, ndim);
247 }
248 }
249
250 /* check if there are boundary points from the previous level */
251 for (j= 0; j< nfboxes; j++)
252 {
253 /* see if the j box of level (i-1) has any boundary layers */
254 if (boxes_with_bdry[i-1][j])
255 {
256 boxi= fbox_mapping[i-1][j];
257 cbox= hypre_BoxArrayBox(cboxes, boxi);
258 fbox= hypre_BoxArrayBox(fboxes, j);
259
260 /* contract the fbox so that divisible in rfactor */
261 contract_fbox= hypre_BoxContraction(fbox, cell_fgrid, rfactors);
262
263 /* refine the cbox. Expand the refined cbox so that the complete
264 chunk of the fine box that coarsened to it is included. This
265 requires some offsets */
266 hypre_ClearIndex(upper_shift);
267 hypre_ClearIndex(lower_shift);
268 for (k= 0; k< ndim; k++)
269 {
270 m= hypre_BoxIMin(contract_fbox)[k];
271 p= m%rfactors[k];
272
273 if (p > 0 && m > 0)
274 {
275 upper_shift[k]= p-1;
276 lower_shift[k]= p-rfactors[k];
277 }
278 else
279 {
280 upper_shift[k]= rfactors[k]-p-1;
281 lower_shift[k]=-p;
282 }
283 }
284 hypre_BoxDestroy(contract_fbox);
285
286 hypre_CopyBox(cbox, &rbox);
287 hypre_StructMapCoarseToFine(hypre_BoxIMin(&rbox), zero_shift,
288 rfactors, hypre_BoxIMin(&rbox));
289 hypre_StructMapCoarseToFine(hypre_BoxIMax(&rbox), zero_shift,
290 rfactors, hypre_BoxIMax(&rbox));
291
292 hypre_AddIndexes(lower_shift, hypre_BoxIMin(&rbox), 3,
293 hypre_BoxIMin(&rbox));
294 hypre_AddIndexes(upper_shift, hypre_BoxIMax(&rbox), 3,
295 hypre_BoxIMax(&rbox));
296
297 /* Determine, if any, boundary layers for this rbox. Since the
298 boundaries of the coarser levels may not be physical, we cannot
299 use hypre_BoxBoundaryDG. But accomplished through intersecting
300 with the finer level boundary boxes. */
301 fbdry= bdry[i-1][j][0]; /* cell-centred boundary layers of level (i-1) */
302 cbdry= bdry[i][boxi][0]; /* cell-centred boundary layers of level i */
303
304 /* fbdry is the cell-centred box_arrayarray. Contains an array of (2*ndim)
305 boxarrays, one for each direction. */
306 cnt= 0;
307 hypre_ForBoxArrayI(l, fbdry)
308 {
309 /* determine which boundary side we are doing. Depending on the
310 boundary, when we coarsen the refined boundary layer, the
311 extents may need to be changed,
312 e.g., index[lower,j,k]= index[upper,j,k]. */
313 switch(l)
314 {
315 case 0: /* lower x direction, x_upper= x_lower */
316 {
317 n= 1; /* n flags whether upper or lower to be replaced */
318 d= 0; /* x component */
319 break;
320 }
321 case 1: /* upper x direction, x_lower= x_upper */
322 {
323 n= 0; /* n flags whether upper or lower to be replaced */
324 d= 0; /* x component */
325 break;
326 }
327 case 2: /* lower y direction, y_upper= y_lower */
328 {
329 n= 1; /* n flags whether upper or lower to be replaced */
330 d= 1; /* y component */
331 break;
332 }
333 case 3: /* upper y direction, y_lower= y_upper */
334 {
335 n= 0; /* n flags whether upper or lower to be replaced */
336 d= 1; /* y component */
337 break;
338 }
339 case 4: /* lower z direction, z_lower= z_upper */
340 {
341 n= 1; /* n flags whether upper or lower to be replaced */
342 d= 2; /* z component */
343 break;
344 }
345 case 5: /* upper z direction, z_upper= z_lower */
346 {
347 n= 0; /* n flags whether upper or lower to be replaced */
348 d= 2; /* z component */
349 break;
350 }
351 }
352
353 box_array= hypre_BoxArrayArrayBoxArray(fbdry, l);
354 hypre_ForBoxI(p, box_array)
355 {
356 hypre_IntersectBoxes(hypre_BoxArrayBox(box_array, p), &rbox,
357 &intersect);
358 if (hypre_BoxVolume(&intersect))
359 {
360 /* coarsen the refined boundary box and append it to
361 boxarray hypre_BoxArrayArrayBoxArray(cbdry, l) */
362 hypre_ProjectBox(&intersect, zero_shift, rfactors);
363 hypre_StructMapFineToCoarse(hypre_BoxIMin(&intersect),
364 zero_shift, rfactors, hypre_BoxIMin(&intersect));
365 hypre_StructMapFineToCoarse(hypre_BoxIMax(&intersect),
366 zero_shift, rfactors, hypre_BoxIMax(&intersect));
367
368 /* the coarsened intersect box may be incorrect because
369 of the box projecting formulas. */
370 if (n) /* replace upper by lower */
371 {
372 hypre_BoxIMax(&intersect)[d]= hypre_BoxIMin(&intersect)[d];
373 }
374 else /* replace lower by upper */
375 {
376 hypre_BoxIMin(&intersect)[d]= hypre_BoxIMax(&intersect)[d];
377 }
378
379 hypre_AppendBox(&intersect,
380 hypre_BoxArrayArrayBoxArray(cbdry, l));
381 cnt++; /* counter to signal boundary layers for cbox boxi */
382 } /* if (hypre_BoxVolume(&intersect)) */
383 } /* hypre_ForBoxI(p, box_array) */
384 } /* hypre_ForBoxArrayI(l, fbdry) */
385
386 /* All the boundary box_arrayarrays have been checked for coarse boxi.
387 Now get the variable boundary layers if any, count the number of
388 boundary points, and appropriately mark boxi. */
389 if (cnt)
390 {
391 hypre_Maxwell_VarBdy(pgrid, bdry[i][boxi]);
392
393 for (p= 0; p< nvars; p++)
394 {
395 cbdry= bdry[i][boxi][p+1];
396 hypre_ForBoxArrayI(l, cbdry)
397 {
398 box_array= hypre_BoxArrayArrayBoxArray(cbdry, l);
399 hypre_ForBoxI(m, box_array)
400 {
401 cbox= hypre_BoxArrayBox(box_array, m);
402 npts[i]+= hypre_BoxVolume(cbox);
403 }
404 }
405 }
406
407 boxes_with_bdry[i][boxi]= 1; /* mark as containing boundary */
408 }
409
410 } /* if (boxes_with_bdry[i-1][j]) */
411 } /* for (j= 0; j< nfboxes; j++) */
412
413 nfboxes= nboxes;
414 } /* for (i= 1; i< num_levels; i++) */
415
416 /* de-allocate objects that are not needed anymore */
417 for (i= 0; i< (num_levels-1); i++)
418 {
419 if (fbox_mapping[i])
420 {
421 hypre_TFree(fbox_mapping[i], HYPRE_MEMORY_HOST);
422 }
423 if (cbox_mapping[i+1])
424 {
425 hypre_TFree(cbox_mapping[i+1], HYPRE_MEMORY_HOST);
426 }
427
428 grid = grid_l[i+1];
429 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
430 cell_cgrid= hypre_SStructPGridCellSGrid(pgrid);
431 cboxes= hypre_StructGridBoxes(cell_cgrid);
432 nboxes= hypre_BoxArraySize(hypre_StructGridBoxes(cell_cgrid));
433 }
434 if (num_levels > 1)
435 {
436 hypre_TFree(fbox_mapping, HYPRE_MEMORY_HOST);
437 hypre_TFree(cbox_mapping, HYPRE_MEMORY_HOST);
438 }
439
440 /* find the ranks for the boundary points */
441 BdryRanks_l = hypre_TAlloc(HYPRE_Int *, num_levels, HYPRE_MEMORY_HOST);
442 BdryRanksCnts_l= hypre_TAlloc(HYPRE_Int , num_levels, HYPRE_MEMORY_HOST);
443
444 /* loop over levels and extract boundary ranks. Only extract unique
445 ranks */
446 for (i= 0; i< num_levels; i++)
447 {
448 grid= grid_l[i];
449 pgrid= hypre_SStructGridPGrid(grid, 0); /* assuming one part */
450 cell_cgrid= hypre_SStructPGridCellSGrid(pgrid);
451 nvars= hypre_SStructPGridNVars(pgrid);
452 cboxes= hypre_StructGridBoxes(cell_cgrid);
453 nboxes= hypre_BoxArraySize(hypre_StructGridBoxes(cell_cgrid));
454
455 ranks= hypre_TAlloc(HYPRE_BigInt, npts[i], HYPRE_MEMORY_HOST);
456 cnt= 0;
457 for (j= 0; j< nboxes; j++)
458 {
459 if (boxes_with_bdry[i][j])
460 {
461 for (k= 0; k< nvars; k++)
462 {
463 fbdry= bdry[i][j][k+1];
464
465 hypre_ForBoxArrayI(m, fbdry)
466 {
467 box_array= hypre_BoxArrayArrayBoxArray(fbdry, m);
468 hypre_ForBoxI(p, box_array)
469 {
470 box= hypre_BoxArrayBox(box_array, p);
471 hypre_BoxGetSize(box, loop_size);
472 hypre_CopyIndex(hypre_BoxIMin(box), start);
473
474 hypre_SerialBoxLoop0Begin(ndim, loop_size);
475 {
476 zypre_BoxLoopGetIndex(lindex);
477 hypre_SetIndex3(index, lindex[0], lindex[1], lindex[2]);
478 hypre_AddIndexes(index, start, 3, index);
479
480 hypre_SStructGridFindBoxManEntry(grid, part, index,
481 k, &boxman_entry);
482 hypre_SStructBoxManEntryGetGlobalRank(boxman_entry, index,
483 &ranks[cnt], matrix_type);
484 cnt++;
485
486 }
487 hypre_SerialBoxLoop0End();
488 } /* hypre_ForBoxI(p, box_array) */
489 } /* hypre_ForBoxArrayI(m, fbdry) */
490
491 } /* for (k= 0; k< nvars; k++) */
492 } /* if (boxes_with_bdry[i][j]) */
493
494 for (k= 0; k< nvars; k++)
495 {
496 hypre_BoxArrayArrayDestroy(bdry[i][j][k+1]);
497 }
498 hypre_BoxArrayArrayDestroy(bdry[i][j][0]);
499 hypre_TFree(bdry[i][j], HYPRE_MEMORY_HOST);
500
501 } /* for (j= 0; j< nboxes; j++) */
502 hypre_TFree(bdry[i], HYPRE_MEMORY_HOST);
503
504 /* mark all ranks that are outside this processor to -1 */
505 for (j= 0; j< cnt; j++)
506 {
507 if ( (ranks[j] < lower_rank[i]) || (ranks[j] > upper_rank[i]) )
508 {
509 ranks[j]= -1;
510 }
511 }
512
513 /* sort the ranks & extract the unique ones */
514 if (cnt) /* recall that some may not have bdry pts */
515 {
516 hypre_BigQsort0(ranks, 0, cnt-1);
517
518 k= 0;
519 if (ranks[0] < 0) /* remove the off-processor markers */
520 {
521 for (j= 1; j< cnt; j++)
522 {
523 if (ranks[j] > -1)
524 {
525 k= j;
526 break;
527 }
528 }
529 }
530
531 l= 1;
532 for (j= k+1; j< cnt; j++)
533 {
534 if (ranks[j] != ranks[j-1])
535 {
536 l++;
537 }
538 }
539 BdryRanks_l[i]= hypre_TAlloc(HYPRE_Int, l, HYPRE_MEMORY_HOST);
540 BdryRanksCnts_l[i]= l;
541
542 l= 0;
543 BdryRanks_l[i][l]= ranks[k]-lower_rank[i];
544 for (j= k+1; j< cnt; j++)
545 {
546 if (ranks[j] != ranks[j-1])
547 {
548 l++;
549 BdryRanks_l[i][l]= ranks[j]-lower_rank[i]; /* store local ranks */
550 }
551 }
552 }
553
554 else /* set BdryRanks_l[i] to be null */
555 {
556 BdryRanks_l[i]= NULL;
557 BdryRanksCnts_l[i]= 0;
558 }
559
560 hypre_TFree(ranks, HYPRE_MEMORY_HOST);
561 hypre_TFree(boxes_with_bdry[i], HYPRE_MEMORY_HOST);
562
563 } /* for (i= 0; i< num_levels; i++) */
564
565 hypre_TFree(boxes_with_bdry, HYPRE_MEMORY_HOST);
566 hypre_TFree(lower_rank, HYPRE_MEMORY_HOST);
567 hypre_TFree(upper_rank, HYPRE_MEMORY_HOST);
568
569 hypre_TFree(bdry, HYPRE_MEMORY_HOST);
570 hypre_TFree(npts, HYPRE_MEMORY_HOST);
571
572 *BdryRanksl_ptr = BdryRanks_l;
573 *BdryRanksCntsl_ptr= BdryRanksCnts_l;
574
575 return ierr;
576 }
577
578 /*-----------------------------------------------------------------------------
579 * Determine the variable boundary layers using the cell-centred boundary
580 * layers. The cell-centred boundary layers are located in bdry[0], a
581 * hypre_BoxArrayArray of size 2*ndim, one array for the upper side and one
582 * for the lower side, for each direction.
583 *-----------------------------------------------------------------------------*/
584 HYPRE_Int
hypre_Maxwell_VarBdy(hypre_SStructPGrid * pgrid,hypre_BoxArrayArray ** bdry)585 hypre_Maxwell_VarBdy( hypre_SStructPGrid *pgrid,
586 hypre_BoxArrayArray **bdry )
587 {
588 HYPRE_Int ierr = 0;
589 HYPRE_Int nvars= hypre_SStructPGridNVars(pgrid);
590
591 hypre_BoxArrayArray *cell_bdry= bdry[0];
592 hypre_BoxArray *box_array, *box_array2;
593 hypre_Box *bdy_box, *shifted_box;
594
595 HYPRE_SStructVariable *vartypes = hypre_SStructPGridVarTypes(pgrid);
596 hypre_Index varoffset, ishift, jshift, kshift;
597 hypre_Index lower, upper;
598
599 HYPRE_Int ndim = hypre_SStructPGridNDim(pgrid);
600 HYPRE_Int i, k, t;
601
602 hypre_SetIndex3(ishift, 1, 0, 0);
603 hypre_SetIndex3(jshift, 0, 1, 0);
604 hypre_SetIndex3(kshift, 0, 0, 1);
605
606 shifted_box= hypre_BoxCreate(ndim);
607 for (i= 0; i< nvars; i++)
608 {
609 t= vartypes[i];
610 hypre_SStructVariableGetOffset(vartypes[i], ndim, varoffset);
611 switch(t)
612 {
613 case 2: /* xface, boundary i= lower, upper */
614 {
615 /* boundary i= lower */
616 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 0);
617 if (hypre_BoxArraySize(box_array))
618 {
619 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 0);
620 hypre_ForBoxI(k, box_array)
621 {
622 bdy_box= hypre_BoxArrayBox(box_array, k);
623
624 /* bdry boxes */
625 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
626 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
627 hypre_SubtractIndexes(lower, varoffset, 3, lower);
628 hypre_SubtractIndexes(upper, varoffset, 3, upper);
629
630 hypre_BoxSetExtents(shifted_box, lower, upper);
631 hypre_AppendBox(shifted_box, box_array2);
632 }
633 }
634
635 /* boundary i= upper */
636 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 1);
637 if (hypre_BoxArraySize(box_array))
638 {
639 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 1);
640 hypre_ForBoxI(k, box_array)
641 {
642 bdy_box= hypre_BoxArrayBox(box_array, k);
643
644 /* bdry boxes */
645 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
646 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
647
648 hypre_BoxSetExtents(shifted_box, lower, upper);
649 hypre_AppendBox(shifted_box, box_array2);
650 }
651 }
652 break;
653 }
654
655 case 3: /* yface, boundary j= lower, upper */
656 {
657 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 2);
658 if (hypre_BoxArraySize(box_array))
659 {
660 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 0);
661 hypre_ForBoxI(k, box_array)
662 {
663 bdy_box= hypre_BoxArrayBox(box_array, k);
664
665 /* bdry boxes */
666 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
667 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
668 hypre_SubtractIndexes(lower, varoffset, 3, lower);
669 hypre_SubtractIndexes(upper, varoffset, 3, upper);
670
671 hypre_BoxSetExtents(shifted_box, lower, upper);
672 hypre_AppendBox(shifted_box, box_array2);
673 }
674 }
675
676 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 3);
677 if (hypre_BoxArraySize(box_array))
678 {
679 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 1);
680 hypre_ForBoxI(k, box_array)
681 {
682 bdy_box= hypre_BoxArrayBox(box_array, k);
683
684 /* bdry boxes */
685 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
686 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
687
688 hypre_BoxSetExtents(shifted_box, lower, upper);
689 hypre_AppendBox(shifted_box, box_array2);
690 }
691 }
692 break;
693 }
694
695 case 5: /* xedge, boundary z_faces & y_faces */
696 {
697 /* boundary k= lower zface*/
698 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 4);
699 if (hypre_BoxArraySize(box_array))
700 {
701 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 0);
702 hypre_ForBoxI(k, box_array)
703 {
704 bdy_box= hypre_BoxArrayBox(box_array, k);
705
706 /* bdry boxes */
707 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
708 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
709 hypre_SubtractIndexes(lower, varoffset, 3, lower);
710 hypre_SubtractIndexes(upper, kshift, 3, upper);
711
712 hypre_BoxSetExtents(shifted_box, lower, upper);
713 hypre_AppendBox(shifted_box, box_array2);
714 }
715 }
716
717 /* boundary k= upper zface*/
718 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 5);
719 if (hypre_BoxArraySize(box_array))
720 {
721 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 1);
722 hypre_ForBoxI(k, box_array)
723 {
724 bdy_box= hypre_BoxArrayBox(box_array, k);
725
726 /* bdry boxes */
727 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
728 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
729 hypre_SubtractIndexes(lower, jshift, 3, lower);
730
731 hypre_BoxSetExtents(shifted_box, lower, upper);
732 hypre_AppendBox(shifted_box, box_array2);
733 }
734 }
735
736 /* boundary j= lower yface*/
737 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 2);
738 if (hypre_BoxArraySize(box_array))
739 {
740 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 2);
741 hypre_ForBoxI(k, box_array)
742 {
743 bdy_box= hypre_BoxArrayBox(box_array, k);
744
745 /* bdry boxes */
746 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
747 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
748 hypre_SubtractIndexes(lower, varoffset, 3, lower);
749 hypre_SubtractIndexes(upper, jshift, 3, upper);
750
751 hypre_BoxSetExtents(shifted_box, lower, upper);
752 hypre_AppendBox(shifted_box, box_array2);
753 }
754 }
755
756 /* boundary j= upper yface*/
757 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 3);
758 if (hypre_BoxArraySize(box_array))
759 {
760 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 3);
761 hypre_ForBoxI(k, box_array)
762 {
763 bdy_box= hypre_BoxArrayBox(box_array, k);
764
765 /* bdry boxes */
766 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
767 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
768 hypre_SubtractIndexes(lower, kshift, 3, lower);
769
770 hypre_BoxSetExtents(shifted_box, lower, upper);
771 hypre_AppendBox(shifted_box, box_array2);
772 }
773 }
774 break;
775 }
776
777 case 6: /* yedge, boundary z_faces & x_faces */
778 {
779 /* boundary k= lower zface*/
780 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 4);
781 if (hypre_BoxArraySize(box_array))
782 {
783 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 0);
784 hypre_ForBoxI(k, box_array)
785 {
786 bdy_box= hypre_BoxArrayBox(box_array, k);
787
788 /* bdry boxes */
789 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
790 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
791 hypre_SubtractIndexes(lower, varoffset, 3, lower);
792 hypre_SubtractIndexes(upper, kshift, 3, upper);
793
794 hypre_BoxSetExtents(shifted_box, lower, upper);
795 hypre_AppendBox(shifted_box, box_array2);
796 }
797 }
798
799 /* boundary k= upper zface*/
800 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 5);
801 if (hypre_BoxArraySize(box_array))
802 {
803 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 1);
804 hypre_ForBoxI(k, box_array)
805 {
806 bdy_box= hypre_BoxArrayBox(box_array, k);
807
808 /* bdry boxes */
809 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
810 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
811 hypre_SubtractIndexes(lower, ishift, 3, lower);
812
813 hypre_BoxSetExtents(shifted_box, lower, upper);
814 hypre_AppendBox(shifted_box, box_array2);
815 }
816 }
817
818 /* boundary i= lower xface*/
819 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 0);
820 if (hypre_BoxArraySize(box_array))
821 {
822 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 2);
823 hypre_ForBoxI(k, box_array)
824 {
825 bdy_box= hypre_BoxArrayBox(box_array, k);
826
827 /* bdry boxes */
828 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
829 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
830 hypre_SubtractIndexes(lower, varoffset, 3, lower);
831 hypre_SubtractIndexes(upper, ishift, 3, upper);
832
833 hypre_BoxSetExtents(shifted_box, lower, upper);
834 hypre_AppendBox(shifted_box, box_array2);
835 }
836 }
837
838 /* boundary i= upper xface*/
839 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 1);
840 if (hypre_BoxArraySize(box_array))
841 {
842 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 3);
843 hypre_ForBoxI(k, box_array)
844 {
845 bdy_box= hypre_BoxArrayBox(box_array, k);
846
847 /* bdry boxes */
848 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
849 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
850 hypre_SubtractIndexes(lower, kshift, 3, lower);
851
852 hypre_BoxSetExtents(shifted_box, lower, upper);
853 hypre_AppendBox(shifted_box, box_array2);
854 }
855 }
856 break;
857 }
858
859 case 7: /* zedge, boundary y_faces & x_faces */
860 {
861 /* boundary j= lower yface*/
862 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 2);
863 if (hypre_BoxArraySize(box_array))
864 {
865 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 0);
866 hypre_ForBoxI(k, box_array)
867 {
868 bdy_box= hypre_BoxArrayBox(box_array, k);
869
870 /* bdry boxes */
871 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
872 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
873 hypre_SubtractIndexes(lower, varoffset, 3, lower);
874 hypre_SubtractIndexes(upper, jshift, 3, upper);
875
876 hypre_BoxSetExtents(shifted_box, lower, upper);
877 hypre_AppendBox(shifted_box, box_array2);
878 }
879 }
880
881 /* boundary j= upper yface*/
882 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 3);
883 if (hypre_BoxArraySize(box_array))
884 {
885 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 1);
886 hypre_ForBoxI(k, box_array)
887 {
888 bdy_box= hypre_BoxArrayBox(box_array, k);
889
890 /* bdry boxes */
891 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
892 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
893 hypre_SubtractIndexes(lower, ishift, 3, lower);
894
895 hypre_BoxSetExtents(shifted_box, lower, upper);
896 hypre_AppendBox(shifted_box, box_array2);
897 }
898 }
899
900 /* boundary i= lower xface*/
901 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 0);
902 if (hypre_BoxArraySize(box_array))
903 {
904 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 2);
905 hypre_ForBoxI(k, box_array)
906 {
907 bdy_box= hypre_BoxArrayBox(box_array, k);
908
909 /* bdry boxes */
910 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
911 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
912 hypre_SubtractIndexes(lower, varoffset, 3, lower);
913 hypre_SubtractIndexes(upper, ishift, 3, upper);
914
915 hypre_BoxSetExtents(shifted_box, lower, upper);
916 hypre_AppendBox(shifted_box, box_array2);
917 }
918 }
919
920 /* boundary i= upper xface*/
921 box_array= hypre_BoxArrayArrayBoxArray(cell_bdry, 1);
922 if (hypre_BoxArraySize(box_array))
923 {
924 box_array2= hypre_BoxArrayArrayBoxArray(bdry[i+1], 3);
925 hypre_ForBoxI(k, box_array)
926 {
927 bdy_box= hypre_BoxArrayBox(box_array, k);
928
929 /* bdry boxes */
930 hypre_CopyIndex(hypre_BoxIMin(bdy_box), lower);
931 hypre_CopyIndex(hypre_BoxIMax(bdy_box), upper);
932 hypre_SubtractIndexes(lower, jshift, 3, lower);
933
934 hypre_BoxSetExtents(shifted_box, lower, upper);
935 hypre_AppendBox(shifted_box, box_array2);
936 }
937 }
938 break;
939 }
940
941 } /* switch(t) */
942 } /* for (i= 0; i< nvars; i++) */
943
944 hypre_BoxDestroy(shifted_box);
945
946 return ierr;
947 }
948
949