1 !-------------------------------------------------------- 2 ! Initialize the Probe Point Arrays and write the Header 3 ! ------------------------------------------------------- 4 subroutine initProbePoints() 5 !Tests if the probe point file xyzts.dat exists, loads probe 6 !point locations, initializes a number of arrays used by 7 !timedataC, and writes the initial header for the output file. 8 ! 9 ! Rewritten by: Nicholas Mati 2014-04-18 10 ! Revision history: 11 ! 2014-04-18 Code moved from itrdrv to here 12 13 use timedataC 14 include "common.h" 15 include "mpif.h" 16 17 logical :: exVarts 18 19 !Test if xyzts.dat exists and broadcast the result. 20 if(myrank.eq.master) inquire(file='xyzts.dat',exist=exts) 21 if(numpe .gt. 1) then 22 call MPI_BARRIER(MPI_COMM_WORLD, ierr) 23 call MPI_BCAST(exts,1,MPI_INTEGER,master,MPI_COMM_WORLD,ierr) 24 endif 25 26 if(.not. exts) return 27 call readProbePoints 28 29 allocate (statptts(ntspts,2)) 30 allocate (parptts( ntspts,nsd)) 31 allocate (varts( ntspts,ndof)) 32 33 statptts(:,:) = 0 34 parptts(:,:) = zero 35 varts(:,:) = zero 36 ivartsbuff = 0 37 vartsResetbuffer = .false. 38 39 allocate (ivarts( ntspts*ndof)) 40 allocate (ivartsg( ntspts*ndof)) 41 allocate (vartssoln( ntspts*ndof)) 42 allocate (vartssolng(ntspts*ndof)) 43 allocate (vartsbuff( ntspts,ndof,nbuff)) 44 allocate (vartsbuffstep(nbuff)) 45 46 !test if the varts folder exists. If it doesn't create it. 47 if(myrank .eq. master) then 48 inquire(file="./varts/.", exist=exVarts) 49 if(.not. exVarts) then 50 call system("mkdir varts") !Only works on *nix, but we 51 !never really run on Windows 52 !anymore so... 53 endif 54 endif 55 56! initProbePoints = exts 57! end function 58 end subroutine 59 60 61 !------------------------ 62 ! Read Probe Point Input 63 !------------------------ 64 subroutine readProbePoints 65 ! Original Code written by: ?? ????-??-?? 66 ! Rewritten by: Nicholas Mati 2014-04-18 67 ! Revision history: 68 ! 2014-04-18 Rewritten code moved from itrdrv to here. 69 ! 70 !Reads the file xyzts.dat for probe point locations, write 71 !frequency, tolerance, ... The file is expected to have the 72 !form: 73 ! ntspts freq tolpt iterat nbuff 74 ! x1 y1 z1 75 ! x2 y2 z2 76 ! ... 77 ! xN yN zN 78 ! 79 ! ... where ntspts is the number of probe points and freq is the 80 ! number of steps to take before flushing data. If nbuff is 81 ! zero, the number of time steps between restarts, ntout, is 82 ! used. 83 84 use timedataC 85 include "common.h" 86 include "mpif.h" 87 88 if(myrank.eq.master) then 89 open(unit=626,file='xyzts.dat',status='old') 90 read(626,*) ntspts, freq, tolpt, iterat, nbuff 91 endif 92 93 !Broadcase out the header of xyzts.dat. These should probably 94 !be combined into two calls, but this is quick and dirty. 95 if(numpe .gt. 1) then 96 call MPI_BARRIER(MPI_COMM_WORLD, ierr) 97 call MPI_Bcast(ntspts, 1, MPI_INTEGER, master, 98 & MPI_COMM_WORLD, ierr) 99 call MPI_Bcast(freq, 1, MPI_INTEGER, master, 100 & MPI_COMM_WORLD, ierr) 101 call MPI_Bcast(tolpt, 1, MPI_DOUBLE_PRECISION, master, 102 & MPI_COMM_WORLD, ierr) 103 call MPI_Bcast(iterat, 1, MPI_INTEGER, master, 104 & MPI_COMM_WORLD, ierr) 105 call MPI_Bcast(nbuff, 1, MPI_INTEGER, master, 106 & MPI_COMM_WORLD, ierr) 107 endif 108 109 allocate (ptts( ntspts,nsd)) 110 111 !Read probe point coordinates and broadcast to the rest of the 112 !processors 113 if(myrank .eq. master) then 114 do jj=1,ntspts ! read coordinate data where solution desired 115 read(626,*) ptts(jj,1), ptts(jj,2), ptts(jj,3) 116 enddo 117 close(626) 118 endif 119 120 if(numpe .gt. 1) then 121 call MPI_BARRIER(MPI_COMM_WORLD, ierr) 122 call MPI_Bcast(ptts, ntspts*nsd, MPI_DOUBLE_PRECISION, 123 & master, MPI_COMM_WORLD, ierr) 124 endif 125 126 if (nbuff .eq. 0) 127 & nbuff=ntout 128 end subroutine 129 130 131 !----------------------------- 132 ! Write the Header varts file 133 !----------------------------- 134 subroutine TD_writeHeader(fvarts) 135 !Creates the file fvarts and writes the data header. 136 !fvarts: Name The file to create 137 138 use timedataC 139 include "common.h" 140 141 character(len=*) fvarts 142 143 !Open the output varts file and write the header 144 if (myrank .eq. master) then 145 146 !fvarts='varts/varts' 147 !fvarts=trim(fvarts)//trim(cname2(lstep)) 148 !fvarts=trim(fvarts)//'.dat' 149 open(unit=1001, file=fvarts, status='unknown') 150 151 !Write the time step 152 write(1001, *) "Time Step: ", Delt(1) 153 write(1001, *) 154 155 !Write the probe locations to varts.ts.dat so that post 156 !processing tools actually know what point goes where. 157 !From experience, it's difficult to keep this straight. 158 write(1001, *) 159 & "Probe ID x y z" 160 do jj = 1, ntspts 161 write(1001, "(I5, T12, 3(F16.12))") jj, ptts(jj,1:3) 162 enddo 163 write(1001, *) 164 165 !write the output format string. This can't be hard 166 !coded because ntspts is not known in advance. 167 write(vartsIOFrmtStr, '("(I8, ", I4, "(E15.7e2))")') 168 & ndof*ntspts 169 170 !Header to delinieate the probe locations with the data. 171 write(1001, *) "Probe Data:" 172 close(unit=1001) 173 endif ! if(myrank .eq. master) 174 end subroutine 175 176 177 178 !------------------------ 179 ! Accumulate Probe Data 180 !------------------------ 181 subroutine TD_bufferData() 182 183 use timedataC 184 include "common.h" 185 include "mpif.h" 186 187 integer :: icheck, istp, nstp 188 189 if (mod(lstep,freq).eq.0) then 190 if(vartsResetBuffer) then 191 ivartsBuff = 0 192 vartsResetBuffer = .false. 193 endif 194 195 !------------------------ 196 !Merge Data across parts 197 !------------------------ 198 if (numpe > 1) then 199 !load the contents of varts into vartssoln 200 do jj = 1, ntspts 201 vartssoln((jj-1)*ndof+1:jj*ndof)=varts(jj,:) 202 ivarts=zero 203 enddo 204 205 !mark which points have been computed on this processor 206 do k=1,ndof*ntspts 207 if(vartssoln(k).ne.zero) ivarts(k)=1 208 enddo 209 210 !merge the solution 211 call MPI_REDUCE(vartssoln, vartssolng, ndof*ntspts, 212 & MPI_DOUBLE_PRECISION, MPI_SUM, master, 213 & MPI_COMM_WORLD, ierr) 214 215 call MPI_REDUCE(ivarts, ivartsg, ndof*ntspts, 216 & MPI_INTEGER, MPI_SUM, master, 217 & MPI_COMM_WORLD, ierr) 218 219 !if the probe point happened to span multiple partitions, 220 !divide the sum by the number of contributing partitions. 221 if (myrank.eq.master) then 222 do jj = 1, ntspts 223 indxvarts = (jj-1)*ndof 224 do k=1,ndof 225 if(ivartsg(indxvarts+k).ne.0) then ! none of the vartssoln(parts) were non zero 226 varts(jj,k) = 227 & vartssolng(indxvarts+k) / ivartsg(indxvarts+k) 228 endif 229 enddo !loop over states / DoF 230 enddo !loop over probe points 231 endif !only on master 232 endif !only if numpe > 1 233 234 ivartsBuff = ivartsBuff + 1 235 if (myrank.eq.master) then 236 if(ivartsBuff .gt. nbuff) then 237 write(*,*) "WARNING: vartsbuff has overflowed" 238 ivartsBuff = nbuff 239 endif 240 241 vartsBuffStep(ivartsBuff) = lstep 242 do jj = 1, ntspts 243 vartsbuff(jj,1:ndof, ivartsBuff) = varts(jj,1:ndof) 244 enddo 245 endif 246 endif 247 248 varts(:,:) = zero 249 250 end subroutine 251 252 253 !------------ 254 ! Write Data 255 !------------ 256 subroutine TD_writeData(fvarts, forceFlush) 257 !writes the probe point data accumulated durring calls to 258 !TD_bufferData(). Note that actual file IO only occurs when the 259 !buffer is full or when DT_writeData is called with forceFlush 260 !set to true. Also note that TD_writeHeader must be called prior 261 !to calling DT_writeData. 262 use timedataC 263 include "common.h" 264 265 character(len=*) :: fvarts 266 logical :: forceFlush 267! logical, optional :: forceflush 268 logical :: flush 269 integer :: k, ibuf 270 271 if (myrank.eq.master) then 272 273 !if provided, use the default value passed in to determine 274 !wheather to flush the buffer 275! if(present(forceFlush)) then !optional version breaks the 276 flush = forceFlush !compiler on Bluegene? 277! else 278! flush = .false. !set the default value 279! endif 280 281 !make sure incomplete buffers get purged at the end of a run 282 !regardless of the default. 283! if(ivartsBuff .eq. nbuff) flush = .true. 284 if(mod(lstep, nbuff) .eq. 0) flush = .true. 285 if(vartsResetBuffer) flush = .false. !Prevent repeated calls without updating 286 !the buffer from writting multiple times. 287 288 if(flush) then !flush the buffer to disc 289 open(unit=1001, file = fvarts, status = "old", 290 & position = "append", action = "write") 291 do ibuf = 1,ivartsBuff 292 write(1001, vartsIOFrmtStr) 293 & vartsBuffStep(ibuf), !write the time step in the first column. 294 & ((vartsbuff(jj,k,ibuf), k=1, ndof) !loop over the variables that you actually want to output. 295 & , jj=1, ntspts) !loop over probe points 296 enddo 297 298 close(1001) 299 300 vartsResetBuffer = .true. 301! ivartsBuff = 0 !need to reset ivartsBuff because 302! !writeDate can be called consecutively 303 endif !only dump when buffer full 304 endif !only on master 305 306! call flush(1001) 307! call fsync(1001) 308 309 !Code for writting one file per probe point 310! do jj = 1, ntspts !loop through probe points 311! ifile = 1000+jj 312! do ibuf=1,nbuff 313! write(ifile,555) lstep-1 -nbuff+ibuf, 314! & (vartsbuff(jj,k,ibuf) , k=1, ndof) 315!! & vartsbuff(jj,:,ibuf) 316! 317! enddo ! buff empty 318! 319! call flush(ifile) 320! enddo ! jj ntspts 321 322 323! varts(:,:) = zero ! reset the array for next step !MOVED FOR Mach Control 324! 555 format(i6,6(2x,E12.5e2)) 325 326 end subroutine 327 328 329 subroutine TD_finalize() 330 use timedataC 331 332 deallocate(ivarts) 333 deallocate(ivartsg) 334 deallocate(vartssoln) 335 deallocate(vartssolng) 336 deallocate(vartsbuff) 337 deallocate(vartsbuffstep) 338 339 deallocate(ptts) 340 deallocate(varts) 341 end subroutine 342 343 344 !--------------------- 345 ! allocate the arrays 346 !--------------------- 347 subroutine sTD 348 !Allocates the arrays statptts, ptts, parptts, and varts for use 349 !in itrdrv and ?? 350 !Subroutine is Depricated. 351 352 use timedataC 353 include "common.h" 354 355 allocate (statptts(ntspts,2)) 356 allocate (ptts(ntspts,nsd)) 357 allocate (parptts(ntspts,nsd)) 358 allocate (varts(ntspts,ndof)) 359 360 return 361 end 362 363 !------------------- 364 ! delete the arrays 365 !------------------- 366 subroutine dTD 367 !Deallocates ptts and varts 368 use timedataC 369 370 deallocate (ptts) 371 deallocate (varts) 372 373 return 374 end 375