19 #include <petsc/finclude/petsc.h>
24 integer(kind=intType),
intent(in) :: level, sps
27 integer(kind=intType) :: i, j, k, l, ii, jj, kk, iBeg, iEnd, jBeg
28 integer(kind=intType) :: jEnd, kBeg, kEnd, nn, mm, iDim, symOnFace(6)
29 integer(kind=intType) :: count, countLocal, symOnFaceLocal(6), procDims(3)
30 integer(kind=intType) :: cellDims(3), nNodes, nCells, ierr, globalInd
31 integer(kind=intType) :: nSubI, nSubJ, nSeed, iSeed, nChanged, loopIter
32 integer(kind=intType) :: iSize, jSize, kSize, nChangedLocal, stackPointer
33 integer(kind=intType),
dimension(:),
allocatable :: indices
34 integer(kind=intType),
dimension(:, :),
allocatable :: lSizes
35 integer(kind=intType),
dimension(:, :),
allocatable :: stack, floodSeeds
37 real(kind=realtype),
dimension(:),
pointer :: cartpointer
38 real(kind=realtype),
dimension(:, :, :),
pointer :: xx
39 real(kind=realtype),
dimension(:, :, :),
pointer :: arrvals, changed
40 real(kind=realtype),
dimension(:),
allocatable :: values
41 real(kind=realtype),
dimension(3) :: xminlocal, xmaxlocal, xmin, xmax, sss
42 real(kind=realtype),
dimension(3) :: pt1, pt2, pt3, pt4, newpt, v1, v2
43 real(kind=realtype) :: arealocal, area, areaavg, err1, err2, h
44 real(kind=realtype) :: cooravg, scalesize, length, u, v
48 vec cartvecglobal, cartveclocal, blankvec, blankveclocal, changedveclocal, changedvecglobal
50 vecscatter blankscatter, blankscatterlocal
83 do j = jbeg + 1, jend + 1
84 do i = ibeg + 1, iend + 1
86 xminlocal(idim) = min(xminlocal(idim), xx(i, j, idim))
87 xmaxlocal(idim) = max(xmaxlocal(idim), xx(i, j, idim))
95 v1 = xx(i, j, :) - xx(i - 1, j - 1, :)
96 v2 = xx(i - 1, j, :) - xx(i, j - 1, :)
98 sss(1) = (v1(2) * v2(3) - v1(3) * v2(2))
99 sss(2) = (v1(3) * v2(1) - v1(1) * v2(3))
100 sss(3) = (v1(1) * v2(2) - v1(2) * v2(1))
102 arealocal = arealocal +
half * (sss(1)**2 + sss(2)**2 + sss(3)**2)
103 countlocal = countlocal + 1
110 call mpi_allreduce(xminlocal, xmin, 3, adflow_real, mpi_min,
adflow_comm_world, ierr)
111 call echk(ierr, __file__, __line__)
113 call mpi_allreduce(xmaxlocal, xmax, 3, adflow_real, mpi_max,
adflow_comm_world, ierr)
114 call echk(ierr, __file__, __line__)
116 call mpi_allreduce(arealocal, area, 1, adflow_real, mpi_sum,
adflow_comm_world, ierr)
117 call echk(ierr, __file__, __line__)
119 call mpi_allreduce(countlocal, count, 1, adflow_integer, mpi_sum,
adflow_comm_world, ierr)
120 call echk(ierr, __file__, __line__)
123 areaavg = area / count
142 scalesize =
mynorm2(xmax - xmin)
173 idim = maxloc(abs(
bcdata(mm)%symNorm), 1)
176 cooravg = sum(xx(ibeg + 1:iend + 1, jbeg + 1:jend + 1, idim)) &
177 / ((iend - ibeg + 1) * (jend - jbeg + 1))
180 err1 = abs(cooravg - xmin(idim)) / scalesize
181 err2 = abs(cooravg - xmax(idim)) / scalesize
183 if (err1 < 1e-8)
then
184 symonfacelocal(2 * idim - 1) = 1
187 if (err2 < 1e-8)
then
188 symonfacelocal(2 * idim) = 1
197 call mpi_allreduce(symonfacelocal, symonface, 6, adflow_integer, mpi_max,
adflow_comm_world, ierr)
198 call echk(ierr, __file__, __line__)
204 if (symonface(2 * idim - 1) == 0)
then
205 xmin(idim) = xmin(idim) - 2 * h
208 if (symonface(2 * idim) == 0)
then
209 xmax(idim) = xmax(idim) + 2 * h
219 celldims = floor((xmax - xmin) / h) + 1
225 if (celldims(1) >= celldims(2) .and. celldims(2) >= celldims(3))
then
226 procdims = (/procdims(3), procdims(2), procdims(1)/)
227 else if (celldims(1) >= celldims(3) .and. celldims(3) >= celldims(2))
then
228 procdims = (/procdims(3), procdims(1), procdims(2)/)
230 else if (celldims(2) >= celldims(1) .and. celldims(1) >= celldims(3))
then
231 procdims = (/procdims(2), procdims(3), procdims(1)/)
232 else if (celldims(2) >= celldims(3) .and. celldims(3) >= celldims(1))
then
233 procdims = (/procdims(1), procdims(3), procdims(2)/)
235 else if (celldims(3) >= celldims(1) .and. celldims(1) >= celldims(2))
then
236 procdims = (/procdims(2), procdims(1), procdims(3)/)
237 else if (celldims(3) >= celldims(2) .and. celldims(2) >= celldims(1))
then
238 procdims = (/procdims(1), procdims(2), procdims(3)/)
245 allocate (lsizes(maxval(procdims), 3))
247 do j = 1, procdims(idim)
249 ii = celldims(idim) / procdims(idim)
250 if (mod(celldims(idim), procdims(idim)) > j - 1)
then
257 print *,
'xmin:', xmin
258 print *,
'xmax:', xmax
259 print *,
'dims:', celldims
261 print *,
'procDims:', procdims
262 print *,
'I sizes:', lsizes(1:procdims(1), 1)
263 print *,
'J sizes:', lsizes(1:procdims(2), 2)
264 print *,
'K sizes:', lsizes(1:procdims(3), 3)
268 dm_boundary_ghosted, dmda_stencil_star, celldims(1), celldims(2), &
269 celldims(3), procdims(1), procdims(2), procdims(3), 1, 1, &
270 lsizes(1:procdims(1), 1), lsizes(1:procdims(2), 2), lsizes(1:procdims(3), 3), &
272 call echk(ierr, __file__, __line__)
287 allocate (indices(i))
320 ii = int((xx(i + 1, j + 1, 1) - xmin(1)) / h)
321 jj = int((xx(i + 1, j + 1, 2) - xmin(2)) / h)
322 kk = int((xx(i + 1, j + 1, 3) - xmin(3)) / h)
324 globalind = celldims(1) * celldims(2) * kk + celldims(1) * jj + ii
326 indices(count) = globalind
332 do i = ibeg + 1, iend
335 pt1 = xx(i, j + 1, :)
336 pt2 = xx(i + 1, j + 1, :)
339 nsubi = int(length / h)
343 newpt = pt1 + dble(k) / (nsubi + 1) * (pt2 - pt1)
344 ii = int((newpt(1) - xmin(1)) / h)
345 jj = int((newpt(2) - xmin(2)) / h)
346 kk = int((newpt(3) - xmin(3)) / h)
348 globalind = celldims(1) * celldims(2) * kk + celldims(1) * jj + ii
350 indices(count) = globalind
356 do j = jbeg + 1, jend
360 pt1 = xx(i + 1, j, :)
361 pt2 = xx(i + 1, j + 1, :)
364 nsubj = int(length / h)
367 newpt = pt1 + dble(k) / (nsubj + 1) * (pt2 - pt1)
368 ii = int((newpt(1) - xmin(1)) / h)
369 jj = int((newpt(2) - xmin(2)) / h)
370 kk = int((newpt(3) - xmin(3)) / h)
372 globalind = celldims(1) * celldims(2) * kk + celldims(1) * jj + ii
374 indices(count) = globalind
380 do j = jbeg + 1, jend
381 do i = ibeg + 1, iend
385 pt2 = xx(i + 1, j, :)
387 pt3 = xx(i + 1, j + 1, :)
388 pt4 = xx(i, j + 1, :)
392 nsubi = int(length / h)
395 nsubi = max(nsubi, int(length / h))
399 nsubj = int(length / h)
402 nsubj = max(nsubj, int(length / h))
406 u = dble(k) / (nsubi + 1)
407 v = dble(l) / (nsubj + 1)
409 newpt = (
one - u) * (
one - v) * pt1 + u * (
one - v) * pt2 + &
410 u * v * pt3 + (
one - u) * v * pt4
412 ii = int((newpt(1) - xmin(1)) / h)
413 jj = int((newpt(2) - xmin(2)) / h)
414 kk = int((newpt(3) - xmin(3)) / h)
416 globalind = celldims(1) * celldims(2) * kk + celldims(1) * jj + ii
418 indices(count) = globalind
427 call dmcreateglobalvector(cartarray, cartvecglobal, ierr)
428 call echk(ierr, __file__, __line__)
431 call vecset(cartvecglobal,
one, ierr)
432 call echk(ierr, __file__, __line__)
434 call dmdagetao(cartarray, cartordering, ierr)
435 call echk(ierr, __file__, __line__)
438 call aoapplicationtopetsc(cartordering,
size(indices), indices, ierr)
439 call echk(ierr, __file__, __line__)
442 allocate (values(count))
445 call vecsetvalues(cartvecglobal, count, indices, values, insert_values, ierr)
446 call echk(ierr, __file__, __line__)
447 deallocate (values, indices)
450 call vecassemblybegin(cartvecglobal, ierr)
451 call echk(ierr, __file__, __line__)
453 call vecassemblyend(cartvecglobal, ierr)
454 call echk(ierr, __file__, __line__)
457 call dmgetlocalvector(cartarray, cartveclocal, ierr)
458 call echk(ierr, __file__, __line__)
460 call dmgetlocalvector(cartarray, changedveclocal, ierr)
461 call echk(ierr, __file__, __line__)
463 call vecset(changedveclocal,
zero, ierr)
464 call echk(ierr, __file__, __line__)
466 call dmgetglobalvector(cartarray, changedvecglobal, ierr)
467 call echk(ierr, __file__, __line__)
470 call dmglobaltolocalbegin(cartarray, cartvecglobal, insert_values, cartveclocal, ierr)
471 call echk(ierr, __file__, __line__)
473 call dmglobaltolocalend(cartarray, cartvecglobal, insert_values, cartveclocal, ierr)
474 call echk(ierr, __file__, __line__)
477 call dmdavecgetarrayf90(cartarray, cartveclocal, arrvals, ierr)
478 call echk(ierr, __file__, __line__)
480 ibeg = lbound(arrvals, 1)
481 jbeg = lbound(arrvals, 2)
482 kbeg = lbound(arrvals, 3)
484 iend = ubound(arrvals, 1)
485 jend = ubound(arrvals, 2)
486 kend = ubound(arrvals, 3)
488 isize = iend - ibeg + 1
489 jsize = jend - jbeg + 1
490 ksize = kend - kbeg + 1
492 call dmdavecrestorearrayf90(cartarray, cartveclocal, arrvals, ierr)
493 call echk(ierr, __file__, __line__)
496 allocate (stack(3, 6 * isize * jsize * ksize + 1))
497 allocate (floodseeds(3, 2 * isize * jsize + 2 * isize * ksize + 2 * jsize * ksize))
501 call dmdavecgetarrayf90(cartarray, cartveclocal, arrvals, ierr)
502 call echk(ierr, __file__, __line__)
504 call dmdavecgetarrayf90(cartarray, changedveclocal, changed, ierr)
505 call echk(ierr, __file__, __line__)
516 if (loopiter == 1)
then
536 do k = kbeg + 1, kend - 1
537 do j = jbeg + 1, jend - 1
538 if (int(changed(ibeg, j, k)) == 1)
then
541 if (int(changed(iend, j, k)) == 1)
then
548 do k = kbeg + 1, kend - 1
549 do i = ibeg + 1, iend - 1
550 if (int(changed(i, jbeg, k)) == 1)
then
553 if (int(changed(i, jend, k)) == 1)
then
560 do j = jbeg + 1, jend - 1
561 do i = ibeg + 1, iend - 1
562 if (int(changed(i, j, kbeg)) == 1)
then
565 if (int(changed(i, j, kend)) == 1)
then
576 stack(:, 1) = floodseeds(:, iseed)
583 do while (stackpointer > 0)
586 i = stack(1, stackpointer)
587 j = stack(2, stackpointer)
588 k = stack(3, stackpointer)
589 stackpointer = stackpointer - 1
591 if (int(arrvals(i, j, k)) == 1)
then
594 changed(i, j, k) =
one
600 nchangedlocal = nchangedlocal + 1
604 arrvals(i, j, k) =
zero
609 if (i - 1 >= ibeg)
then
610 stackpointer = stackpointer + 1
611 stack(:, stackpointer) = (/i - 1, j, k/)
614 if (i + 1 <= iend)
then
615 stackpointer = stackpointer + 1
616 stack(:, stackpointer) = (/i + 1, j, k/)
619 if (j - 1 >= jbeg)
then
620 stackpointer = stackpointer + 1
621 stack(:, stackpointer) = (/i, j - 1, k/)
624 if (j + 1 <= jend)
then
625 stackpointer = stackpointer + 1
626 stack(:, stackpointer) = (/i, j + 1, k/)
629 if (k - 1 >= kbeg)
then
630 stackpointer = stackpointer + 1
631 stack(:, stackpointer) = (/i, j, k - 1/)
634 if (k + 1 <= kend)
then
635 stackpointer = stackpointer + 1
636 stack(:, stackpointer) = (/i, j, k + 1/)
642 call dmdavecrestorearrayf90(cartarray, cartveclocal, arrvals, ierr)
643 call echk(ierr, __file__, __line__)
645 call dmdavecrestorearrayf90(cartarray, changedveclocal, changed, ierr)
646 call echk(ierr, __file__, __line__)
649 call dmlocaltoglobalbegin(cartarray, changedveclocal, insert_values, changedvecglobal, ierr)
650 call echk(ierr, __file__, __line__)
652 call dmlocaltoglobalend(cartarray, changedveclocal, insert_values, changedvecglobal, ierr)
653 call echk(ierr, __file__, __line__)
655 call dmglobaltolocalbegin(cartarray, changedvecglobal, insert_values, changedveclocal, ierr)
656 call echk(ierr, __file__, __line__)
658 call dmglobaltolocalend(cartarray, changedvecglobal, insert_values, changedveclocal, ierr)
659 call echk(ierr, __file__, __line__)
662 call mpi_allreduce(nchangedlocal, nchanged, 1, adflow_integer, mpi_sum, &
664 call echk(ierr, __file__, __line__)
666 print *,
'Cart Flood Iteration:', loopiter,
'Blanked ', nchanged,
'Interior Cells.'
669 if (nchanged == 0)
then
670 exit parallelsyncloop
673 loopiter = loopiter + 1
674 end do parallelsyncloop
676 deallocate (stack, floodseeds)
681 call dmdavecgetarrayf90(cartarray, cartveclocal, arrvals, ierr)
682 call echk(ierr, __file__, __line__)
684 do k = kbeg + 1, kend - 1
685 do j = jbeg + 1, jend - 1
686 do i = ibeg + 1, iend - 1
687 if (int(arrvals(i, j, k)) == 1)
then
688 arrvals(i, j, k) = -
three
694 call dmdavecrestorearrayf90(cartarray, cartveclocal, arrvals, ierr)
695 call echk(ierr, __file__, __line__)
698 call dmrestorelocalvector(cartarray, cartveclocal, ierr)
699 call echk(ierr, __file__, __line__)
701 call dmrestorelocalvector(cartarray, changedveclocal, ierr)
702 call echk(ierr, __file__, __line__)
704 call dmrestoreglobalvector(cartarray, changedvecglobal, ierr)
705 call echk(ierr, __file__, __line__)
708 call dmlocaltoglobalbegin(cartarray, cartveclocal, insert_values, cartvecglobal, ierr)
709 call echk(ierr, __file__, __line__)
711 call dmlocaltoglobalend(cartarray, cartveclocal, insert_values, cartvecglobal, ierr)
712 call echk(ierr, __file__, __line__)
718 i = celldims(1) * celldims(2) * celldims(3)
720 call echk(ierr, __file__, __line__)
723 call vecgetownershiprange(blankvec, i, j, ierr)
724 call echk(ierr, __file__, __line__)
726 call echk(ierr, __file__, __line__)
728 call isduplicate(is2, is1, ierr)
729 call aoapplicationtopetscis(cartordering, is1, ierr)
730 call echk(ierr, __file__, __line__)
732 call vecscattercreate(cartvecglobal, is1, blankvec, is2, blankscatter, ierr)
733 call echk(ierr, __file__, __line__)
735 call vecscatterbegin(blankscatter, cartvecglobal, blankvec, &
736 insert_values, scatter_forward, ierr)
737 call echk(ierr, __file__, __line__)
739 call vecscatterend(blankscatter, cartvecglobal, blankvec, &
740 insert_values, scatter_forward, ierr)
741 call echk(ierr, __file__, __line__)
745 call vecscatterdestroy(blankscatter, ierr)
746 call echk(ierr, __file__, __line__)
748 call isdestroy(is1, ierr)
749 call echk(ierr, __file__, __line__)
751 call isdestroy(is2, ierr)
752 call echk(ierr, __file__, __line__)
754 call vecdestroy(cartvecglobal, ierr)
755 call echk(ierr, __file__, __line__)
757 call dmdestroy(cartarray, ierr)
758 call echk(ierr, __file__, __line__)
766 call vecscattercreatetoall(blankvec, blankscatterlocal, blankveclocal, ierr)
767 call echk(ierr, __file__, __line__)
769 call vecscatterbegin(blankscatterlocal, blankvec, blankveclocal, insert_values, scatter_forward, ierr)
770 call echk(ierr, __file__, __line__)
772 call vecscatterend(blankscatterlocal, blankvec, blankveclocal, insert_values, scatter_forward, ierr)
773 call echk(ierr, __file__, __line__)
775 call vecscatterdestroy(blankscatterlocal, ierr)
776 call echk(ierr, __file__, __line__)
778 call vecgetarrayf90(blankveclocal, cartpointer, ierr)
779 call echk(ierr, __file__, __line__)
789 x(i - 1, j - 1, k - 1, idim) + &
790 x(i, j - 1, k - 1, idim) + &
791 x(i - 1, j, k - 1, idim) + &
792 x(i, j, k - 1, idim) + &
793 x(i - 1, j - 1, k, idim) + &
794 x(i, j - 1, k, idim) + &
795 x(i - 1, j, k, idim) + &
800 ii = int((pt1(1) - xmin(1)) / h)
801 jj = int((pt1(2) - xmin(2)) / h)
802 kk = int((pt1(3) - xmin(3)) / h)
805 ii = min(max(0, ii), celldims(1) - 1)
806 jj = min(max(0, jj), celldims(2) - 1)
807 kk = min(max(0, kk), celldims(3) - 1)
810 globalind = celldims(1) * celldims(2) * kk + celldims(1) * jj + ii
816 if (length > 1.73205080 * h)
then
817 if (int(cartpointer(globalind + 1)) == -3)
then
823 if (int(cartpointer(globalind + 1)) == -3)
then
832 call vecrestorearrayf90(blankveclocal, cartpointer, ierr)
833 call echk(ierr, __file__, __line__)
836 call vecdestroy(blankvec, ierr)
837 call echk(ierr, __file__, __line__)
839 call vecdestroy(blankveclocal, ierr)
840 call echk(ierr, __file__, __line__)
843 domainloop:
do nn = 1, ndom
844 flowdoms(nn, level, sps)%intCommVars(1)%var => &
845 flowdoms(nn, level, sps)%iblank(:, :, :)
856 integer(kind=intType),
intent(in) :: i, j, k
858 floodseeds(:, nseed) = (/i, j, k/)
866 integer(kind=intType),
intent(in) :: i, j, k
869 if (i >= 0 .and. i <= iend - 1 .and. j >= 0 .and. j <= jend - 1 .and. k >= 0 .and. k <= kend - 1)
then
879 #include <petsc/finclude/petsc.h>
884 integer(kind=intType),
intent(in),
dimension(3) :: cellDims
885 real(kind=realtype),
intent(in),
dimension(3) :: xmin
886 real(kind=realtype),
intent(in) :: h
889 integer(kind=intType) :: ierr
892 character(len=32) :: coorNames(3)
893 integer(kind=intType) :: base, zoneID, coordID, cg, zone, iField, iSol, i, j, k, iDim
894 real(kind=realtype),
dimension(:, :, :, :),
allocatable :: xtmp
895 real(kind=realtype),
dimension(:),
pointer :: cartpointer
899 vecscatter blankscatterlocal
901 coornames(1) =
"CoordinateX"
902 coornames(2) =
"CoordinateY"
903 coornames(3) =
"CoordinateZ"
905 call vecscattercreatetozero(blankvec, blankscatterlocal, blankveclocal, ierr)
906 call echk(ierr, __file__, __line__)
908 call vecscatterbegin(blankscatterlocal, blankvec, blankveclocal, &
909 insert_values, scatter_forward, ierr)
910 call echk(ierr, __file__, __line__)
912 call vecscatterend(blankscatterlocal, blankvec, blankveclocal, &
913 insert_values, scatter_forward, ierr)
914 call echk(ierr, __file__, __line__)
918 call cg_open_f(
"cartblock.cgns", mode_write, cg, ierr)
920 call cg_base_write_f(cg,
"Base#1", 3, 3, base, ierr)
922 call cg_zone_write_f(cg, base,
"cartblock", int((/celldims(1) + 1, celldims(2) + 1, celldims(3) + 1, &
923 celldims(1), celldims(2), celldims(3), 0, 0, 0/), &
924 cgsize_t), structured, zoneid, ierr)
926 allocate (xtmp(celldims(1) + 1, celldims(2) + 1, celldims(3) + 1, 3))
928 do k = 1, celldims(3) + 1
929 do j = 1, celldims(2) + 1
930 do i = 1, celldims(1) + 1
931 xtmp(i, j, k, 1) = xmin(1) + (i - 1) * h
936 do k = 1, celldims(3) + 1
937 do j = 1, celldims(2) + 1
938 do i = 1, celldims(1) + 1
939 xtmp(i, j, k, 2) = xmin(2) + (j - 1) * h
944 do k = 1, celldims(3) + 1
945 do j = 1, celldims(2) + 1
946 do i = 1, celldims(1) + 1
947 xtmp(i, j, k, 3) = xmin(3) + (k - 1) * h
953 call cg_coord_write_f(cg, base, zoneid, realdouble, coornames(idim), &
954 xtmp(:, :, :, idim), coordid, ierr)
958 call cg_sol_write_f(cg, base, zoneid,
"flowSolution", cellcenter, isol, ierr)
959 call vecgetarrayf90(blankveclocal, cartpointer, ierr)
960 call echk(ierr, __file__, __line__)
962 call cg_field_write_f(cg, base, zoneid, isol, realdouble,
"iBlank", &
963 cartpointer, ifield, ierr)
965 call vecrestorearrayf90(blankveclocal, cartpointer, ierr)
966 call echk(ierr, __file__, __line__)
967 call cg_close_f(cg, ierr)
970 call vecdestroy(blankveclocal, ierr)
971 call echk(ierr, __file__, __line__)
973 call vecscatterdestroy(blankscatterlocal, ierr)
974 call echk(ierr, __file__, __line__)
983 integer(kind=intType),
intent(in) :: N
984 integer(kind=intType),
intent(out) :: s(3)
987 integer(kind=intType) :: a, b, a1, b1, a2, b2, s1(3), s2(3)
1004 if (s1(1) > s2(1))
then
1014 integer(kind=intType) :: N, f1, f2, i, j, s
1017 s = int(sqrt(dble(n)))
1019 if (mod(n, j) == 0)
then
subroutine addseed(i, j, k)
logical function onblock(i, j, k)
integer(kind=inttype), dimension(:, :, :), pointer iblank
integer(kind=inttype), dimension(:), pointer bcfaceid
integer(kind=inttype) nbocos
integer(kind=inttype), dimension(:), pointer bctype
real(kind=realtype), dimension(:, :, :, :), pointer x
real(kind=realtype), dimension(:, :, :, :), pointer xseed
subroutine triplefactor(N, s)
subroutine largefactor(N, f1, f2)
subroutine writecartmesh(blankVec, cellDims, xMin, h)
subroutine createcartmesh(level, sps)
type(internalcommtype), dimension(:), allocatable internalcell_2nd
type(commtype), dimension(:), allocatable commpatterncell_2nd
integer adflow_comm_world
integer(kind=inttype), parameter eulerwall
real(kind=realtype), parameter zero
real(kind=realtype), parameter three
integer(kind=inttype), parameter imax
integer(kind=inttype), parameter kmin
integer(kind=inttype), parameter jmax
real(kind=realtype), parameter eighth
integer(kind=inttype), parameter nswalladiabatic
integer(kind=inttype), parameter symm
real(kind=realtype), parameter one
real(kind=realtype), parameter half
integer(kind=inttype), parameter imin
integer(kind=inttype), parameter nswallisothermal
integer(kind=inttype), parameter ibcgroupwalls
real(kind=realtype), parameter large
integer(kind=inttype), parameter kmax
integer(kind=inttype), parameter jmin
subroutine whalo1to1intgeneric(nVar, level, sps, commPattern, internal)
subroutine getwallsize(famList, nNodes, nCells, dualMesh)
subroutine qsortintegers(arr, nn)
type(bcgrouptype), dimension(nfamexchange) bcfamgroups
real(kind=realtype) function mynorm2(x)
subroutine echk(errorcode, file, line)
subroutine setpointers(nn, mm, ll)