1 2 subroutine D2wall (dwall, x) 3c 4c---------------------------------------------------------------------- 5c 6c This routine calculates the distance to the nearest wall for all 7c field points. 8c 9c 10c Zdenek Johan, Winter 1991. (Fortran 90) 11c---------------------------------------------------------------------- 12c 13 use pointer_data 14 include "common.h" 15 include "mpif.h" 16c 17 dimension x(numnp,nsd), dwall(numnp), 18 & nwall(numpe), idisp(numpe) 19 20c 21c Allocatable arrays 22c 23 allocatable xwi(:,:,:), xw(:,:,:) 24 25c 26c.... Count the welts (wall-elements) 27c 28 nwalli=0 29 do iblk = 1, nelblb ! loop over boundary elt blocks 30 npro = lcblkb(1,iblk+1) - lcblkb(1,iblk) 31 do j = 1, npro 32 if(miBCB(iblk)%p(j,5).eq.2) nwalli=nwalli+1 33 enddo 34 enddo 35c 36c.... Create wallnode-coord list for welts on processor 37c 38 if (nwalli.eq.0) nwalli=1 ! patch for mpi's lack of imagination 39 allocate (xwi(nwalli,nenb+1,nsd)) 40 xwi = 1.0d32 41 xwi(:,nenb+1,:)=zero 42 nwalli = 0 43 do iblk = 1, nelblb ! loop over boundary elt blocks 44c 45 iel = lcblkb(1,iblk) 46 nenbl = lcblkb(6,iblk) ! no. of vertices per bdry. face 47 npro = lcblkb(1,iblk+1) - iel 48c 49 do j = 1, npro ! loop over belts in this blk 50 if(miBCB(iblk)%p(j,5).eq.2) then 51 nwalli=nwalli+1 52c assemble local coord list 53 do node = 1, nenbl 54 xwi(nwalli,node,1:3)=x(mienb(iblk)%p(j,node),:) 55 enddo 56 do node = 1, nenbl 57 xwi(nwalli,nenb+1,:)=xwi(nwalli,nenb+1,:) 58 & +xwi(nwalli,node,:) 59 enddo 60 xwi(nwalli,nenb+1,:)=xwi(nwalli,nenb+1,:)/nenbl 61c 62 endif 63 enddo ! loop over belts in this blk 64c 65 enddo ! loop over boundary elt blocks 66c 67 if (nwalli.eq.0) xwi=1.0e32 ! fix for mpi's lack of imagination 68 if (nwalli.eq.0) nwalli=1 ! patch for mpi's lack of imagination 69c 70c Pool "number of welts" info from all processors 71c 72 if (numpe.gt.1) then 73 call MPI_ALLGATHER(nwalli,1,MPI_INTEGER,nwall,1, 74 & MPI_INTEGER,MPI_COMM_WORLD,ierr) 75 nwallt=0 76 do j=1,numpe 77 nwallt=nwallt+nwall(j) 78 enddo 79 else 80 nwall=nwalli 81 nwallt=nwalli 82 endif 83c 84c Make all-processor wallnode-coord collage 85c 86 allocate (xw(nwallt,nenb+1,nsd)) 87 if (numpe.gt.1) then 88 idisp(:)=0 89 do j=2,numpe 90 idisp(j)=idisp(j-1)+nwall(j-1) 91 enddo 92 do j=1,nenb+1 93 do k=1,nsd 94 call MPI_ALLGATHERV(xwi(:,j,k),nwalli, 95 & MPI_DOUBLE_PRECISION,xw(:,j,k),nwall,idisp, 96 & MPI_DOUBLE_PRECISION,MPI_COMM_WORLD,ierr) 97 enddo 98 enddo 99 else 100 xw=xwi 101 endif 102 103c 104c For each point, calculate distance to centroid; 105c save the distance in this node's position of dwall if it's 106c shorter than currently stored distance 107c 108 dwall=1.0e32 109 do i=1,numnp 110 do j=1, nwallt 111 do k=1,nenb+1 112 distance = ( x(i,1) - xw(j,k,1) )**2 113 & +( x(i,2) - xw(j,k,2) )**2 114 & +( x(i,3) - xw(j,k,3) )**2 115 if ( dwall(i).gt.distance ) dwall(i) = distance 116 enddo 117 enddo 118 enddo 119 dwall=sqrt(dwall) 120c 121 deallocate(xwi) 122 deallocate(xw) 123c 124 return 125 end 126 127 128 129