1 | include(dom.inc)
   2 | 
   3 |       PROGRAM dom2ascii
   4 | 
   5 | !       ================================================================!
   6 | !                                                                       !
   7 | !       dom2ascii.F : Converts PRISSMA *.out files into ascii files     !
   8 | !                     and interpolates the data over a line             !
   9 | !                                                                       !
  10 | !       author      : J. AMAYA (mars 2007)                              !
  11 | !                                                                       !
  12 | !       ================================================================!
  13 | 
  14 |         IMPLICIT NONE
  15 | 
  16 |         include 'dom_constants.h'
  17 | 
  18 |         DOM_INT       :: npoints, i, j, l, typeunit, readwall
  19 |         DOM_INT       :: ntcells, ntnodes, ntfaces, ndirs, nbdyfaces
  20 |         DOM_INT       :: nmaxfaces, nmax
  21 |         DOM_INT       :: ios, ic, nf, iw, jw
  22 | 
  23 |         DOM_REAL      :: x1, y1, z1, x2, y2, z2
  24 |         DOM_REAL      :: xd, yd, zd
  25 |         DOM_REAL      :: segmentsize, d,unity
  26 | 
  27 |         DOM_REAL, ALLOCATABLE, DIMENSION(:)   :: x,  y,  z
  28 |         DOM_REAL, ALLOCATABLE, DIMENSION(:)   :: xn, yn, zn
  29 |         DOM_REAL, ALLOCATABLE, DIMENSION(:,:) :: xw, yw, zw
  30 |         DOM_REAL, ALLOCATABLE, DIMENSION(:)   :: xfc,yfc,zfc
  31 |         DOM_REAL, ALLOCATABLE, DIMENSION(:)   :: Sr,  H,  G,  Qw
  32 |         DOM_REAL, ALLOCATABLE, DIMENSION(:)   :: Srn, Hn, Gn, Qwn 
  33 |         DOM_REAL, ALLOCATABLE, DIMENSION(:,:) :: Qr
  34 |         DOM_REAL, ALLOCATABLE, DIMENSION(:,:) :: Qrn
  35 | 
  36 |         CHARACTER*80  :: infilespath, outfilespath
  37 |         CHARACTER*80  :: GFile, SrFile, HFile, QwFile, QrFile
  38 |         CHARACTER*80  :: clfacefile, gdatafile,nodessfile,cfacesfile
  39 |         CHARACTER*80 :: gfile2,srfile2,hfile2,qwfile2,qrfile2
  40 | 
  41 | !       ------------------!
  42 | !       Read choices file !
  43 | !       ------------------!
  44 | 
  45 |         OPEN (FILE_CHCS , FILE='dom2ascii.choices', FORM='FORMATTED')
  46 |         READ (FILE_CHCS,*) infilespath
  47 |         READ (FILE_CHCS,*) outfilespath
  48 |         READ (FILE_CHCS,*) x1, y1, z1
  49 |         READ (FILE_CHCS,*) x2, y2, z2
  50 |         READ (FILE_CHCS,*) npoints
  51 |         READ (FILE_CHCS,*) typeunit
  52 |         READ (FILE_CHCS,*) readwall
  53 |         CLOSE(FILE_CHCS)
  54 | 
  55 |         IF (typeunit.eq.1) THEN
  56 |           unity=1.
  57 |         ELSE IF (typeunit.eq.2) THEN
  58 |           unity=1000.
  59 |         ENDIF
  60 |   
  61 |         IF (npoints.lt.2) THEN
  62 |           WRITE(*,*) " Fatal error: the minimum number of points is 2"
  63 |           STOP
  64 |         ENDIF
  65 | 
  66 | !       -----------------!
  67 | !       Allocate vectors !
  68 | !       -----------------!
  69 | 
  70 |         ALLOCATE(x(npoints))
  71 |         ALLOCATE(y(npoints))
  72 |         ALLOCATE(z(npoints))
  73 | 
  74 | !       ---------------------------!
  75 | !       Calculate geometrical data !
  76 | !       ---------------------------!
  77 | 
  78 |         d  = sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)+(z1-z2)*(z1-z2))
  79 | 
  80 |         segmentsize = d / real(npoints - 1)
  81 |         xd          = (x2-x1) / d
  82 |         yd          = (y2-y1) / d
  83 |         zd          = (z2-z1) / d
  84 | 
  85 |         DO i=1, npoints
  86 |           x(i) = x1 + (i-1) * segmentsize * xd
  87 |           y(i) = y1 + (i-1) * segmentsize * yd
  88 |           z(i) = z1 + (i-1) * segmentsize * zd
  89 |         ENDDO
  90 | 
  91 | !       --------------------!
  92 | !       Open PRISSMA files !
  93 | !       --------------------!
  94 | 
  95 |         GFile  = trim(outfilespath)//'/G.out'
  96 |         SrFile = trim(outfilespath)//'/Sr.out'
  97 |         HFile  = trim(outfilespath)//'/H.out'
  98 |         QwFile = trim(outfilespath)//'/Qw.out'
  99 |         QrFile = trim(outfilespath)//'/Qr.out'
 100 | 
 101 |         gdatafile  = trim(infilespath)//'/Global.in'
 102 |         clfacefile = trim(infilespath)//'/CLFaces.in'
 103 |         nodessfile = trim(infilespath)//'/Nodelist.in'
 104 |         cfacesfile = trim(infilespath)//'/Centerfaces.in'
 105 | 
 106 |         OPEN(FILE_G ,FILE=Gfile ,FORM='UNFORMATTED')
 107 |         OPEN(FILE_Sr,FILE=SrFile,FORM='UNFORMATTED')
 108 |         OPEN(FILE_H ,FILE=HFile ,FORM='UNFORMATTED')
 109 |         OPEN(FILE_Qw,FILE=QwFile,FORM='UNFORMATTED')
 110 |         OPEN(FILE_Qr,FILE=QrFile,FORM='UNFORMATTED')
 111 | 
 112 |         OPEN(FILE_CLFAC, FILE=clfacefile, FORM='UNFORMATTED')
 113 |         OPEN(FILE_GDATA, FILE=gdatafile , FORM='UNFORMATTED')
 114 |         OPEN(FILE_NODES, FILE=nodessfile, FORM='UNFORMATTED')
 115 |         OPEN(FILE_CFACE, FILE=cfacesfile, FORM='UNFORMATTED')
 116 | 
 117 | !       ----------------!
 118 | !       Get global data !
 119 | !       ----------------!
 120 | 
 121 |         READ(FILE_GDATA) ntcells, ntnodes, ntfaces,nmaxfaces,ndirs
 122 |         READ(FILE_CLFAC) nbdyfaces
 123 | 
 124 |         WRITE(*,*) "Dirs : ", ndirs
 125 |         WRITE(*,*) "Cells: ", ntcells
 126 |         WRITE(*,*) "Nodes: ", ntnodes
 127 |         WRITE(*,*) "Faces: ", ntfaces
 128 |         WRITE(*,*) "Boundary faces: ", nbdyfaces
 129 |         WRITE(*,*)
 130 | 
 131 |         CLOSE(FILE_GDATA)
 132 |         CLOSE(FILE_CLFAC)
 133 | 
 134 | !       -----------------!
 135 | !       Allocate vectors !
 136 | !       -----------------!
 137 | 
 138 |         ALLOCATE(xn (ntnodes))
 139 |         ALLOCATE(yn (ntnodes))
 140 |         ALLOCATE(zn (ntnodes))
 141 | 
 142 |         ALLOCATE(xw (nmaxfaces,ntcells))
 143 |         ALLOCATE(yw (nmaxfaces,ntcells))
 144 |         ALLOCATE(zw (nmaxfaces,ntcells))
 145 | 
 146 |         ALLOCATE(xfc(nbdyfaces))
 147 |         ALLOCATE(yfc(nbdyfaces))
 148 |         ALLOCATE(zfc(nbdyfaces))
 149 | 
 150 |         ALLOCATE(Sr (npoints))
 151 |         ALLOCATE(H  (ntcells))
 152 |         ALLOCATE(G  (npoints))
 153 |         ALLOCATE(Qw (ntcells))
 154 |         ALLOCATE(Srn(ntnodes))
 155 |         ALLOCATE(Hn (ntcells))
 156 |         ALLOCATE(Gn (ntnodes))
 157 |         ALLOCATE(Qwn(ntcells))
 158 | 
 159 |         ALLOCATE(Qr (3,npoints))
 160 |         ALLOCATE(Qrn(3,ntnodes))
 161 | 
 162 | !       -------------------!
 163 | !       Initialise vectors !
 164 | !       -------------------!
 165 | 
 166 |         H   = 0.
 167 |         G   = 0.
 168 |         Qw  = 0.
 169 |         Sr  = 0.
 170 |         Qr  = 0.
 171 |         Hn  = 0.
 172 |         Gn  = 0.
 173 |         Qwn = 0.
 174 |         Srn = 0.
 175 |         Qrn = 0.
 176 | 
 177 | !       ---------------------!
 178 | !       Read geometical data !
 179 | !       ---------------------!
 180 | 
 181 |         xw = 0
 182 |         yw = 0
 183 |         zw = 0
 184 | 
 185 |         DO i = 1, ntcells
 186 |           READ(FILE_CFACE) ic,nf,(xw(j,ic),yw(j,ic),zw(j,ic),j=1,nf)
 187 |         ENDDO
 188 | 
 189 |         DO i = 1, ntnodes
 190 |           READ(FILE_NODES) ic,xn(ic),yn(ic),zn(ic)
 191 |         ENDDO
 192 | 
 193 | 
 194 |         CLOSE(FILE_NODES)
 195 |         CLOSE(FILE_CFACE)
 196 | 
 197 | !       ---------------!
 198 | !       Read wall data !
 199 | !       ---------------!
 200 | 
 201 |         IF (readwall.eq.1) THEN
 202 | 
 203 |           DO i=1,nbdyfaces
 204 | 
 205 |             READ(FILE_Qw) iw,jw,Qwn(i)
 206 |             READ(FILE_H ) iw,jw,Hn(i)
 207 | 
 208 |             xfc(i) = xw(iw,jw)
 209 |             yfc(i) = yw(iw,jw)
 210 |             zfc(i) = zw(iw,jw)
 211 | 
 212 |           ENDDO
 213 | 
 214 |           CLOSE(FILE_Qw)
 215 |           CLOSE(FILE_H )
 216 | 
 217 | !         ----------------------!
 218 | !         Interpolate wall data !
 219 | !         ----------------------!
 220 | 
 221 |           IF(npoints.ge.ntcells) THEN
 222 |             nmax = npoints
 223 |           ELSE
 224 |             nmax = ntcells
 225 |           ENDIF            
 226 | 
 227 |           CALL INTERP3D(npoints,  x,  y,  z,  Qw,                       &
 228 |      &                  nbdyfaces,xfc,yfc,zfc,Qwn,nmax)
 229 | 
 230 |           CALL INTERP3D(npoints,  x,  y,  z,  H,                        &
 231 |      &                  nbdyfaces,xfc,yfc,zfc,Hn,nmax)
 232 | 
 233 |         ENDIF
 234 | 
 235 |         DEALLOCATE(xfc,yfc,zfc,xw,yw,zw)
 236 | 
 237 | !       ---------------!
 238 | !       Read cell data !
 239 | !       ---------------!
 240 | 
 241 |         DO i=1,ntnodes
 242 | 
 243 |           READ(FILE_G ) ic,Gn(ic)
 244 |           READ(FILE_Sr) ic,Srn(ic)
 245 | 
 246 |         ENDDO
 247 | 
 248 |         CLOSE(FILE_G )
 249 |         CLOSE(FILE_Sr)
 250 |         CLOSE(FILE_Qr)
 251 | 
 252 | !       ----------------------!
 253 | !       Interpolate cell data !
 254 | !       ----------------------!
 255 | 
 256 |         IF(npoints.ge.ntnodes) THEN
 257 |           nmax = npoints
 258 |         ELSE
 259 |           nmax = ntnodes
 260 |         ENDIF
 261 | 
 262 |         CALL INTERP3D(npoints,x,y,z,G,ntnodes,xn,yn,zn,Gn,nmax)
 263 |         CALL INTERP3D(npoints,x,y,z,Sr,ntnodes,xn,yn,zn,Srn,nmax)
 264 |         CALL INTERP3D(npoints,x,y,z,Qr(1,:),ntnodes,xn,yn,zn,Qrn(1,:),  &
 265 |      &                nmax)
 266 |         CALL INTERP3D(npoints,x,y,z,Qr(2,:),ntnodes,xn,yn,zn,Qrn(2,:),  &
 267 |      &                nmax)
 268 |         CALL INTERP3D(npoints,x,y,z,Qr(3,:),ntnodes,xn,yn,zn,Qrn(3,:),  &
 269 |      &                nmax)
 270 | 
 271 | !       ------------------!
 272 | !       Write ascii files !
 273 | !       ------------------!
 274 | 
 275 |         gfile2     = trim(outfilespath)//'/G.dat'
 276 |         srfile2    = trim(outfilespath)//'/Sr.dat'
 277 |         hfile2     = trim(outfilespath)//'/H.dat'
 278 |         qwfile2    = trim(outfilespath)//'/Qw.dat'
 279 |         qrfile2    = trim(outfilespath)//'/Qr.dat'
 280 | 
 281 |         OPEN(FILE_G2 ,FILE=gfile2 ,FORM='formatted')
 282 |         OPEN(FILE_Sr2,FILE=srfile2,FORM='formatted')
 283 |         OPEN(FILE_H2 ,FILE=hfile2 ,FORM='formatted')
 284 |         OPEN(FILE_Qw2,FILE=qwfile2,FORM='formatted')
 285 |         OPEN(FILE_Qr2,FILE=qrfile2,FORM='formatted')
 286 | 
 287 |         DO i=1,npoints
 288 | 
 289 |           WRITE(FILE_H2 ,'(3f16.8,1f18.6)') x(i),y(i),z(i),H(i)/unity
 290 |           WRITE(FILE_Qw2,'(3f16.8,1f18.6)') x(i),y(i),z(i),Qw(i)/unity
 291 |           WRITE(FILE_G2 ,'(3f16.8,1f18.6)') x(i),y(i),z(i),G(i)/unity
 292 |           WRITE(FILE_Sr2,'(3f16.8,1f18.6)') x(i),y(i),z(i),Sr(i)/unity
 293 |           WRITE(FILE_Qr2,'(3f16.8,3f18.6)') x(i),y(i),z(i),             &
 294 |      &                                     (Qr(l,i)/unity,l=1,3)
 295 |         ENDDO
 296 | 
 297 |         CLOSE(FILE_H2 )
 298 |         CLOSE(FILE_G2 )
 299 |         CLOSE(FILE_Sr2)
 300 |         CLOSE(FILE_Qr2)
 301 |         CLOSE(FILE_Qw2)
 302 | 
 303 |       END PROGRAM dom2ascii