1 subroutine getvel (y,ilwork, iBC, 2 & nsons, ifath, velbar) 3 4 5 include "common.h" 6 include "mpif.h" 7 include "auxmpi.h" 8 9c 10 dimension nsons(nfath), rinvsons(nfath), 11 & velo(numnp,nflow), velf(nfath,nflow), 12 & velft(nfath,nflow), 13 & velbar(nfath,nflow), 14 & y(numnp,ndof), 15 & ifath(numnp), 16 & ilwork(nlwork), iBC(numnp) 17 18c 19c for now keep the compressible numbering in velbar 20c 21 velo(:,1)=y(:,4) 22 velo(:,2)=y(:,1) 23 velo(:,3)=y(:,2) 24 velo(:,4)=y(:,3) 25 if(nflow.eq.5) velo(:,5)=y(:,5) 26 nsonmax=maxval(nsons) 27 if ((nsonmax.eq.1)) then ! we are doing local clipping -no homog dir 28 velft=velo 29 else 30c 31c zero on processor periodic nodes so that they will not be added twice 32c 33 where(btest(iBC,10).or.btest(iBC,12)) 34 velo(:,1)=zero 35 velo(:,2)=zero 36 velo(:,3)=zero 37 velo(:,4)=zero 38 endwhere 39 if(nflow.eq.5) then 40 where(btest(iBC,10).or.btest(iBC,12)) 41 velo(:,5)=zero 42 endwhere 43 endif 44 if (numpe.gt.1) then 45 46 numtask = ilwork(1) 47 itkbeg = 1 48 49c zero the nodes that are "solved" on the other processors 50 do itask = 1, numtask 51 52 iacc = ilwork (itkbeg + 2) 53 numseg = ilwork (itkbeg + 4) 54 55 if (iacc .eq. 0) then 56 do is = 1,numseg 57 isgbeg = ilwork (itkbeg + 3 + 2*is) 58 lenseg = ilwork (itkbeg + 4 + 2*is) 59 isgend = isgbeg + lenseg - 1 60 velo(isgbeg:isgend,:) = zero 61 velo(isgbeg:isgend,:) = zero 62 enddo 63 endif 64 65 itkbeg = itkbeg + 4 + 2*numseg 66 67 enddo !itask 68 69 endif ! numpe.gt.1 70 71 72 velf = zero 73c 74c accumulate sum of sons to the fathers 75c 76 do i = 1,numnp 77 ifathi=ifath(i) 78 velf(ifathi,1:nflow) = velf(ifathi,1:nflow) 79 & + velo(i,1:nflow) 80 enddo 81 82c 83c Now the true fathers and serrogates combine results and update 84c each other. 85c 86 if(numpe .gt. 1) then 87 call drvAllreduce(velf, velft,nfath*nflow) 88 else 89 velft=velf 90 endif 91c 92c xvelf is the sum of the sons for each father on this processor 93c 94c xvelft is the sum of the sons for each father on all processor combined 95c (the same as if we had not partitioned the mesh for each processor) 96c 97c divide by # of sons to get average father for this step 98c 99 rinvsons = one/nsons ! division is expensive 100 velft(:,1) = velft(:,1) * rinvsons(:) ! / nsons(:) 101 velft(:,2) = velft(:,2) * rinvsons(:) ! / nsons(:) 102 velft(:,3) = velft(:,3) * rinvsons(:) ! / nsons(:) 103 velft(:,4) = velft(:,4) * rinvsons(:) ! / nsons(:) 104 if(nflow.eq.5) velft(:,5) = velft(:,5) * rinvsons(:) 105 endif ! end of homog direction averaging 106 denom=max(one*(lstep),one) 107 if(wtavei.lt.0) then 108 tavef=one/denom 109 else 110 tavef=wtavei 111 endif 112 velbar(:,1)=velft(:,1) 113 velbar(:,2)=velft(:,2) 114 velbar(:,3)=velft(:,3) 115 velbar(:,4)=velft(:,4) 116 117 if(istep.eq.0) then 118 velbar(:,:)=velft(:,:) 119 else 120 velbar(:,:)=tavef*velft(:,:)+(one-tavef)*velbar(:,:) 121 endif 122 123 return 124 end 125 126 127 128 129 130 131