1 subroutine proces 2c 3c---------------------------------------------------------------------- 4c 5c This subroutine generates the problem data and calls the solution 6c driver. 7c 8c 9c Zdenek Johan, Winter 1991. (Fortran 90) 10c---------------------------------------------------------------------- 11c 12 use readarrays ! used to access x, iper, ilwork 13 use turbsa ! used to access d2wall 14 use dtnmod 15 use periodicity 16 use pvsQbi 17 include "common.h" 18 include "mpif.h" 19c 20c arrays in the following 2 lines are now dimensioned in readnblk 21c dimension x(numnp,nsd) 22c dimension iper(nshg), ilwork(nlwork) 23c 24 dimension y(nshg,ndof), 25 & iBC(nshg), 26 & BC(nshg,ndofBC), 27 & ac(nshg,ndof) 28c 29c.... shape function declarations 30c 31 dimension shp(MAXTOP,maxsh,MAXQPT), 32 & shgl(MAXTOP,nsd,maxsh,MAXQPT), 33 & shpb(MAXTOP,maxsh,MAXQPT), 34 & shglb(MAXTOP,nsd,maxsh,MAXQPT) 35c 36c stuff for dynamic model s.w.avg and wall model 37c 38c dimension ifath(numnp), velbar(nfath,nflow), 39c & nsons(nfath) 40 41 dimension velbar(nfath,nflow) 42c 43c stuff to interpolate profiles at inlet 44c 45 real*8 bcinterp(100,ndof+1),interp_mask(ndof) 46 logical exlog 47 48c if ((irscale .ge. 0) .and. (myrank.eq.master)) then 49c call setSPEBC(numnp, point2nfath, nsonmax) 50c endif 51c 52c.... generate the geometry and boundary conditions data 53c 54 call gendat (y, ac, point2x, 55 & iBC, BC, 56 & point2iper, point2ilwork, shp, 57 & shgl, shpb, shglb, 58 & point2ifath, velbar, point2nsons ) 59 call setper(nshg) 60 call perprep(iBC,point2iper,nshg) 61 if (iLES/10 .eq. 1) then 62 call keeplhsG ! Allocating space for the mass (Gram) matrix to be 63 ! used for projection filtering and reconstruction 64 ! of the strain rate tensor. 65 66 call setrls ! Allocating space for the resolved Leonard stresses 67c See bardmc.f 68 endif 69c 70c.... time averaged statistics 71c 72 if (ioform .eq. 2) then 73 call initStats(point2x, iBC, point2iper, point2ilwork) 74 endif 75c 76c.... RANS turbulence model 77c 78 if (iRANS .lt. 0) then 79 call initTurb( point2x ) 80 endif 81c 82c.... p vs. Q boundary 83c 84 call initNABI( point2x, shpb ) 85c 86c.... check for execution mode 87c 88 if (iexec .eq. 0) then 89 lstep = 0 90 call restar ('out ', y ,ac) 91 return 92 endif 93c 94c.... initialize AutoSponge 95c 96 if(matflg(5,1).ge.4) then ! cool case (sponge) 97 call initSponge( y,point2x) 98 endif 99c 100c 101c.... adjust BC's to interpolate from file 102c 103 104 inquire(file="inlet.dat",exist=exlog) 105 if(exlog) then 106 open (unit=654,file="inlet.dat",status="old") 107 read(654,*) ninterp,ibcmatch,(interp_mask(j),j=1,ndof) 108 do i=1,ninterp 109 read(654,*) (bcinterp(i,j),j=1,ndof+1) ! distance to wall+ 110 ! ndof but note order of BC's 111 ! p,t,u,v,w,scalars. Also note that must be in 112 ! increasing distance from the wall order. 113 114 enddo 115 do i=1,nshg ! only correct for linears at this time 116 if(mod(iBC(i),1024).eq.ibcmatch) then 117 iupper=0 118 do j=2,ninterp 119 if(bcinterp(j,1).gt.d2wall(i)) then !bound found 120 xi=(d2wall(i)-bcinterp(j-1,1))/ 121 & (bcinterp(j,1)-bcinterp(j-1,1)) 122 iupper=j 123 exit 124 endif 125 enddo 126 if(iupper.eq.0) then 127 write(*,*) 'failure in finterp, ynint=',d2wall(i) 128 stop 129 endif 130 if(interp_mask(1).ne.zero) then 131 BC(i,1)=(xi*bcinterp(iupper,2) 132 & +(one-xi)*bcinterp(iupper-1,2))*interp_mask(1) 133 endif 134 if(interp_mask(2).ne.zero) then 135 BC(i,2)=(xi*bcinterp(iupper,3) 136 & +(one-xi)*bcinterp(iupper-1,3))*interp_mask(2) 137 endif 138 if(interp_mask(3).ne.zero) then 139 BC(i,3)=(xi*bcinterp(iupper,4) 140 & +(one-xi)*bcinterp(iupper-1,4))*interp_mask(3) 141 endif 142 if(interp_masK(4).ne.zero) then 143 BC(i,4)=(xi*bcinterp(iupper,5) 144 & +(one-xi)*bcinterp(iupper-1,5))*interp_mask(4) 145 endif 146 if(interp_mask(5).ne.zero) then 147 BC(i,5)=(xi*bcinterp(iupper,6) 148 & +(one-xi)*bcinterp(iupper-1,6))*interp_mask(5) 149 endif 150 if(interp_mask(6).ne.zero) then 151 BC(i,7)=(xi*bcinterp(iupper,7) 152 & +(one-xi)*bcinterp(iupper-1,7))*interp_mask(6) 153 endif 154 if(interp_mask(7).ne.zero) then 155 BC(i,8)=(xi*bcinterp(iupper,8) 156 & +(one-xi)*bcinterp(iupper-1,8))*interp_mask(7) 157 endif 158 endif 159 enddo 160 endif 161c$$$$$$$$$$$$$$$$$$$$ 162 163c 164c 165c.... call the semi-discrete predictor multi-corrector iterative driver 166c 167 call itrdrv (y, ac, 168 & uold, point2x, 169 & iBC, BC, 170 & point2iper, point2ilwork, shp, 171 & shgl, shpb, shglb, 172 & point2ifath, velbar, point2nsons ) 173c 174c.... return 175c 176c 177c.... stop CPU-timer 178c 179CAD call timer ('End ') 180c 181c.... close echo file 182c 183 184!MR CHANGE 185 if (numpe > 1) call MPI_BARRIER(MPI_COMM_WORLD, ierr) 186 if(myrank.eq.0) then 187 write(*,*) 'process - before closing iecho' 188 endif 189!MR CHANGE END 190 191 close (iecho) 192 193!MR CHANGE 194 if (numpe > 1) call MPI_BARRIER(MPI_COMM_WORLD, ierr) 195 if(myrank.eq.0) then 196 write(*,*) 'process - after closing iecho' 197 endif 198!MR CHANGE END 199 200 201c 202c.... end of the program 203c 204CAD write(6,*) 'Life: ', second(0) - ttim(100) 205 deallocate(point2iper) 206 if(numpe.gt.1) then 207 call Dctypes(point2ilwork(1)) 208 deallocate(point2ilwork) 209 endif 210 deallocate(point2x) 211 deallocate(point2nsons) 212 deallocate(point2ifath) 213 deallocate(uold) 214 deallocate(wnrm) 215 deallocate(otwn) 216 call finalizeDtN 217 call clearper 218 call finalizeNABI 219 220 return 221 end 222 223 224