1 subroutine rstatic (res, y, Dy) 2c 3c---------------------------------------------------------------------- 4c 5c This subroutine calculates the statistics of the residual. 6c 7c input: 8c res (nshg,nflow) : preconditioned residual 9c 10c output: 11c The time step, cpu-time and entropy-norm of the residual 12c are printed in the file HISTOR.DAT. 13c 14c 15c Zdenek Johan, Winter 1991. (Fortran 90) 16c---------------------------------------------------------------------- 17c 18 include "common.h" 19 include "mpif.h" 20 include "auxmpi.h" 21c 22 dimension res(nshg,nflow), mproc(1) 23!SCATTER dimension rvec(numpe) 24 dimension rtmp(nshg), nrsmax(1) 25!SCATTER dimension irecvcount(numpe), resvec(numpe) 26 27 real*8 y(nshg,ndof), Dy(nshg,4) 28c integer tmrc 29c 30c$$$ ttim(68) = ttim(68) - tmr() 31c 32c.... compute max delta y 33c 34 rdy1 = zero 35 rdy2 = zero 36 rdy4 = zero 37 rdy5 = zero 38 call sumgatN( abs(gami*Delt(itseq) 39 & * Dy(1:numnp,1:3)),3,rdy1, numnp) 40 call sumgatN( abs( y(1:numnp,1:3)),3,rdy2,numnp) 41 call sumgatN( abs(gami*alfi*Delt(itseq) 42 & * Dy(1:numnp,4)),1,rdy4,numnp) 43 call sumgatN( abs( y(1:numnp,4)), 1,rdy5,numnp) 44 rmaxdyU = rdy1/rdy2 45 rmaxdyP = rdy4/rdy5 46 47c 48c..... Signal to quit if delta is very small. look in itrdrv.f for the 49c completion of the hack. 50c 51 if( rmaxdyU .lt. dtol(1) .and. rmaxdyP .lt. dtol(2)) then 52 istop = 1000 53 endif 54 55 if (numpe == 1) nshgt=nshg ! global = this processor 56c 57c 58c.... -----------------------> Convergence <------------------------- 59c 60c.... compute the maximum residual and the corresponding node number 61c 62 rtmp = zero 63 if (impl(itseq) .ge. 9) then 64 do i = 1, nflow 65 rtmp = rtmp + res(:,i)**2 ! only add continuity and momentum 66 enddo 67 endif 68 69 call sumgat (rtmp, 1, resnrm) 70 71 resmaxl = maxval(rtmp) 72 73 irecvcount = 1 74 resvec = resmaxl 75 if (numpe > 1) then 76 if(impistat.eq.1) then 77 iAllR = iAllR+1 78 elseif(impistat.eq.2) then 79 iAllRScal = iAllRScal+1 80 endif 81 if(impistat2.eq.1) call MPI_BARRIER (MPI_COMM_WORLD, ierr) 82 if(impistat.gt.0) rmpitmr = TMRC() 83 call MPI_ALLREDUCE (resvec, resmax, irecvcount, 84 & MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, ierr) 85 if(impistat.eq.1) then 86 rAllR = rAllR+TMRC()-rmpitmr 87 elseif(impistat.eq.2) then 88 rAllRScal = rAllRScal+TMRC()-rmpitmr 89 endif 90c call MPI_REDUCE_SCATTER (resvec, resmax, irecvcount, 91c & MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, ierr) 92 if (resmax .eq. resvec ) then 93c if (resmax .eq. resvec(1) ) then 94 mproc(1) = myrank 95 nrsmax = maxloc(rtmp) 96 else 97 mproc(1) = -1 98 nrsmax(1) = -1 99 endif 100 resvec = nrsmax(1) 101 if(impistat.eq.1) then 102 iAllR = iAllR+1 103 elseif(impistat.eq.2) then 104 iAllRScal = iAllRScal+1 105 endif 106 if(impistat2.eq.1) call MPI_BARRIER (MPI_COMM_WORLD, ierr) 107 if(impistat.gt.0) rmpitmr = TMRC() 108 call MPI_ALLREDUCE (resvec, rvec, irecvcount, 109 & MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, ierr) 110 if(impistat.eq.1) then 111 rAllR = rAllR+TMRC()-rmpitmr 112 elseif(impistat.eq.2) then 113 rAllRScal = rAllRScal+TMRC()-rmpitmr 114 endif 115c call MPI_REDUCE_SCATTER (resvec, rvec, irecvcount, 116c & MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, ierr) 117 nrsmax = rvec 118c nrsmax = rvec(1) 119 resvec = mproc(1) 120 if(impistat.eq.1) then 121 iAllR = iAllR+1 122 elseif(impistat.eq.2) then 123 iAllRScal = iAllRScal+1 124 endif 125 if(impistat2.eq.1) call MPI_BARRIER (MPI_COMM_WORLD, ierr) 126 if(impistat.gt.0) rmpitmr = TMRC() 127 call MPI_ALLREDUCE (resvec, rvec, irecvcount, 128 & MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, ierr) 129 if(impistat.eq.1) then 130 rAllR = rAllR+TMRC()-rmpitmr 131 elseif(impistat.eq.2) then 132 rAllRScal = rAllRScal+TMRC()-rmpitmr 133 endif 134c call MPI_REDUCE_SCATTER (resvec, rvec, irecvcount, 135c & MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, ierr) 136 mproc(1) = rvec 137c mproc = rvec(1) 138 else 139 resmax = resmaxl 140 nrsmax = maxloc(rtmp) 141 mproc(1) = 0 142 endif 143c 144c.... correct the residuals 145c 146 if (loctim(itseq) .eq. 0) then 147 resnrm = resnrm 148 resmax = resmax 149 else 150 resnrm = resnrm 151 resmax = resmax 152 endif 153c 154c.... approximate the number of entries 155c 156 totres = resnrm / float(nshgt) 157 totres = sqrt(totres) 158 resmax = sqrt(resmax) 159 if (resfrt(1) .eq. zero) resfrt(1) = totres 160 jtotrs = int ( 10.d0 * log10 ( totres / resfrt(1) ) ) 161 jresmx = int ( 10.d0 * log10 ( resmax / totres ) ) 162c 163c.... get the CPU-time 164c 165CAD cputme = (second(0) - ttim(100)) 166 rsec=TMRC() 167 cputme = (rsec - ttim(100)) 168c 169c.... output the result 170c 171c if (numpe > 1) call MPI_BARRIER (MPI_COMM_WORLD, ierr) 172 173 if (myrank .eq. master) then 174c 175c.... results of continuity and momentum 176c 177 178 print 2000, lstep+1, cputme, totres, jtotrs, rmaxdyU, 179 & rmaxdyP,nrsmax, 180 & mproc(1)+1, jresmx, int(statsflow(4)), 181 & int(statsflow(1)) 182 write (ihist,2000) lstep+1, cputme, totres, jtotrs, 183 & rmaxdyU, rmaxdyP, nrsmax, 184 & mproc(1)+1,jresmx,int(statsflow(4)), 185 & int(statsflow(1)) 186 187 call flush(ihist) 188 endif 189c if(numpe>1) call MPI_BARRIER (MPI_COMM_WORLD,ierr) 190 191c$$$ ttim(68) = ttim(68) + tmr() 192 193c 194c.... return 195c 196 return 197c 198 1000 format(1p,i6,5e13.5) 199 2000 format(1p,i6,e10.3,e10.3,2x,'(',i4,')',2x,e10.3,2x,e10.3, 200 & 2x,'<',i6,'-',i5,'|', 201 & i4,'>', ' [', i4,' -',i4,']') 202 3000 format(1p,i6,e10.3,e10.3,3x,'(',i4,')',3x,'<',i6,'-',i5,'|', 203 & i4,'>', ' [', i4,' -',i4,' -',i4,']') 204 205c 206 end 207 208 209 subroutine rstaticSclr (res, y, Dy, icomp) 210c 211c---------------------------------------------------------------------- 212c 213c This subroutine calculates the statistics of the residual 214c 215c---------------------------------------------------------------------- 216c 217 include "common.h" 218 include "mpif.h" 219 include "auxmpi.h" 220c 221 dimension res(nshg) 222 dimension rtmp(nshg) 223 real*8 y(nshg,ndof), Dy(nshg), nrm 224c integer tmrc 225c 226c.... compute max delta y 227c 228 rdy1 = zero 229 rdy2 = zero 230c 231c.... normalize turbulence with molecular viscosity 232c 233 if ( (icomp .eq. 6).and. (iRANS.eq.-1) ) then 234 nrm = datmat(1,2,1) 235 else 236 nrm = zero 237 endif 238 call sumgat( abs(gami*Delt(itseq)*Dy(:)),1,rdy1) 239 call sumgat( abs( y(:,icomp)),1,rdy2) 240 rmaxdyT = rdy1/(rdy2+nrm) 241c 242c.... compute the maximum residual and the corresponding node number 243c 244 rtmp = zero 245 rtmp = rtmp + res**2 ! add temperature also 246 call sumgat (rtmp, 1, resnrm) 247 248 if (numpe == 1) nshgt=nshg ! global = this processor 249 250 totres = resnrm / float(nshgt) 251 totres = sqrt(totres) 252 253c if (mod(impl(1),100)/10 .eq. 0) then !not solving flow 254 if (myrank .eq. master) then 255c 256c.... get the CPU-time 257c 258 rsec=TMRC() 259 cputme = (rsec - ttim(100)) 260 261 print 802, lstep+1, cputme, totres, rmaxdyT, 262 & int(statssclr(1)) 263 write (ihist,802) lstep+1, cputme, totres, 264 & rmaxdyT,int(statssclr(1)) 265 266 call flush(ihist) 267 endif 268c else 269c if (myrank .eq. master) then 270c print 803, totres, rmaxdyT, int(statssclr(1)) 271c write(ihist,803) totres, rmaxdyT, int(statssclr(1)) 272c endif 273c endif 274 275 return 276 277 802 format(1p,i6,e10.3,e10.3,10X,e10.3,31X,'[',i6,']') 278 803 format(1p,16x,e10.3,10x,e10.3,31X,'[',i10,']') 279 280 end 281 282 283 284