pkg/regrid/regrid_init_varia.F

aa076db465 Ed H*0001 #include "REGRID_OPTIONS.h"
                0002
                0003 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0004 CBOP 0
                0005 C     !ROUTINE: REGRID_INIT_VARIA
                0006
                0007 C     !INTERFACE:
                0008       SUBROUTINE REGRID_INIT_VARIA( myThid )
                0009
                0010 C     !DESCRIPTION:
                0011 C     Initialize REGRID variables
                0012
                0013 C     !USES:
                0014       IMPLICIT NONE
                0015 #include "SIZE.h"
                0016 #include "EEPARAMS.h"
                0017 #include "PARAMS.h"
                0018 #include "GRID.h"
                0019 #include "REGRID_SIZE.h"
                0020 #include "REGRID.h"
                0021 #ifdef ALLOW_EXCH2
f9f661930b Jean*0022 #include "W2_EXCH2_SIZE.h"
aa076db465 Ed H*0023 #include "W2_EXCH2_TOPOLOGY.h"
                0024 #endif
                0025       INTEGER  ILNBLNK
                0026       EXTERNAL ILNBLNK
                0027
                0028 C     !INPUT/OUTPUT PARAMETERS:
                0029 C     myThid ::  my Thread Id number
                0030       INTEGER myThid
                0031 CEOP
                0032
                0033 C     !LOCAL VARIABLES:
                0034       INTEGER i,k, iface, uniq_tnum, bi,bj
                0035       INTEGER irx, isrc, idst, nFx,nFy, init_nlpts,nlpts
                0036       INTEGER iUnit, errIO, nnb
                0037       INTEGER iminx,iminy, imaxx,imaxy
                0038       _RL wt
                0039       CHARACTER*(MAX_LEN_FNAM) fname
                0040       CHARACTER*(MAX_LEN_MBUF) msgbuf
                0041       LOGICAL  exst
                0042
                0043 C     Regrid files contain information on a per-face basis.  This is
                0044 C     convenient in two respects: (1) the domain can be re-tiled without
                0045 C     changing any of the files [since the ordering with respect to
                0046 C     tiles is performed here in the model] and (2) when faces are
                0047 C     removed or added only the corresponding per-face files will need
                0048 C     to be removed or added [and all the other per-face files remain
                0049 C     unchanged provided the face numbers do not change].
                0050 C
c424ee7cc7 Jean*0051 C     The convention is: "points cycle most quickly in X and then Y"
aa076db465 Ed H*0052 C
                0053 C        +-------------------+
                0054 C        |  Face             |
                0055 C        |                   |
                0056 C        |       +-----+     |
                0057 C      Y |       |Tile |     |
                0058 C        |       +-----+     |
                0059 C        |                   |
                0060 C        |123...             |
                0061 C        +-------------------+
                0062 C                X
                0063
                0064       _BEGIN_MASTER( myThid )
                0065
c424ee7cc7 Jean*0066       WRITE(msgBuf,'(a)')
aa076db465 Ed H*0067      &     '// ======================================================='
                0068       CALL PRINT_MESSAGE(msgBuf,standardMessageUnit,
                0069      &     SQUEEZE_RIGHT,myThid)
c424ee7cc7 Jean*0070       WRITE(msgBuf,'(a)')
aa076db465 Ed H*0071      &     '// Begin reading the per-face REGRID information'
                0072       CALL PRINT_MESSAGE(msgBuf,standardMessageUnit,
                0073      &     SQUEEZE_RIGHT,myThid)
c424ee7cc7 Jean*0074       WRITE(msgBuf,'(a)')
aa076db465 Ed H*0075      &     '// ======================================================='
                0076       CALL PRINT_MESSAGE(msgBuf,standardMessageUnit,
                0077      &     SQUEEZE_RIGHT,myThid)
                0078
                0079       nlpts = 0
                0080
                0081       CALL MDSFINDUNIT(iUnit, myThid)
                0082
                0083       DO bj = myByLo(myThid), myByHi(myThid)
                0084         DO bi = myBxLo(myThid), myBxHi(myThid)
c424ee7cc7 Jean*0085
aa076db465 Ed H*0086 #ifdef ALLOW_EXCH2
                0087 C         EXCH2 domains
c424ee7cc7 Jean*0088           uniq_tnum = W2_myTileList(bi,bj)
aa076db465 Ed H*0089           iface = exch2_myFace(uniq_tnum)
                0090           nFx = exch2_mydnx(uniq_tnum)
                0091           nFy = exch2_mydny(uniq_tnum)
                0092           iminx = exch2_tbasex(uniq_tnum) + 1
                0093           imaxx = iminx + exch2_tnx(uniq_tnum) - 1
                0094           iminy = exch2_tbasey(uniq_tnum) + 1
                0095           imaxy = iminy + exch2_tny(uniq_tnum) - 1
                0096 #else
                0097 C         Global tile number for simple single-face "EXCH1" domains
                0098           iG = bi + (myXGlobalLo-1)/sNx
                0099           jG = bj + (myYGlobalLo-1)/sNy
                0100           uniq_tnum = (jG - 1)*(nPx*nSx) + iG
                0101           iface = 1
                0102           nFx = nSx * sNx
                0103           nFy = sSy * sNy
                0104           iminx = myXGlobalLo
                0105           imaxx = myXGlobalLo + sNx - 1
                0106           iminy = myYGlobalLo
                0107           imaxy = myYGlobalLo + sNy - 1
                0108 #endif
                0109
                0110 C         WRITE(*,*) 'iminx, imaxx, nFx, nFy = ',
                0111 C         &         iminx, imaxx, nFx, nFy
                0112
                0113 C         Read through all the weights files for this tile (face) and
                0114 C         locate the points that belong to this tile
                0115           DO i = 1,regrid_ngrids
                0116
                0117             IF (i .EQ. 1) THEN
                0118               nlpts = 0
                0119             ELSE
                0120               nlpts = REGRID_iend(i,bi,bj)
                0121             ENDIF
                0122             init_nlpts = nlpts
                0123
                0124             DO k = 1,MAX_LEN_FNAM
                0125               fname(k:k) = ' '
                0126             ENDDO
                0127             nnb = ILNBLNK(REGRID_fbname_in(i))
c424ee7cc7 Jean*0128             write(fname,'(a,i3.3,a)')
aa076db465 Ed H*0129      &           REGRID_fbname_in(i)(1:nnb),iface,'.regrid.ascii'
                0130             nnb = ILNBLNK(fname)
                0131             INQUIRE( FILE=fname, EXIST=exst )
                0132             IF (.NOT. exst) THEN
                0133               WRITE(msgBuf,'(A)')  'S/R REGRID_INIT_VARIA()'
                0134               CALL PRINT_ERROR( msgBuf , 1)
c424ee7cc7 Jean*0135               WRITE(msgBuf,'(3A)')  ' File "',
aa076db465 Ed H*0136      &             fname(1:nnb), '" does not exist'
                0137               CALL PRINT_ERROR( msgBuf , 1)
                0138               CLOSE(iUnit)
                0139               STOP ' stopped in REGRID_INIT_VARIA()'
                0140             ENDIF
                0141
                0142             open(unit=iUnit, file=fname, status='old', iostat=errIO)
                0143
                0144             IF (errIO .LT. 0) THEN
                0145               WRITE(msgBuf,'(A)')  'S/R REGRID_INIT_VARIA()'
                0146               CALL PRINT_ERROR( msgBuf , 1)
c424ee7cc7 Jean*0147               WRITE(msgBuf,'(3A)')  'Unable to open file="',
aa076db465 Ed H*0148      &             fname(1:nnb), '"'
                0149               CALL PRINT_ERROR( msgBuf , 1)
                0150               CLOSE(iUnit)
                0151               STOP ' stopped in REGRID_INIT_VARIA()'
                0152             ELSE
                0153               WRITE(msgBuf,'(3a)') 'Reading file "', fname(1:nnb),'"'
                0154               call PRINT_MESSAGE(msgBuf,standardMessageUnit,
                0155      &             SQUEEZE_RIGHT,myThid)
                0156             ENDIF
                0157
                0158             DO WHILE ( .TRUE. )
                0159 C             READ(iUnit,fmt='(2(I10,1X),1P1E23.13E3)',iostat=errIO)
c424ee7cc7 Jean*0160               READ(iUnit,fmt='(2(1X,I10),1X,E28.22)',iostat=errIO)
aa076db465 Ed H*0161      &             isrc, idst, wt
                0162               IF ( errIO .NE. 0 ) THEN
                0163                 GOTO 100
                0164               ENDIF
                0165               irx = MOD(isrc,nFx)
                0166               IF (irx .EQ. 0)  irx = nFx
                0167               IF ((iminx .LE. irx) .AND. (irx .LE. imaxx)) THEN
                0168                 nlpts = nlpts + 1
                0169                 REGRID_i_loc(nlpts,bi,bj) = irx
                0170                 REGRID_j_loc(nlpts,bi,bj) = isrc/nFx + 1
                0171                 REGRID_i_out(nlpts,bi,bj) = idst
                0172                 REGRID_amat(nlpts,bi,bj)  = wt
                0173               ENDIF
                0174
                0175             ENDDO
                0176  100        CONTINUE
                0177             close(iUnit)
                0178             WRITE(msgBuf,'(a,i10)') '  num weights read = ',
                0179      &           (nlpts - init_nlpts)
                0180             call PRINT_MESSAGE(msgBuf,standardMessageUnit,
                0181      &           SQUEEZE_RIGHT,myThid)
                0182
                0183             REGRID_ibeg(i,bi,bj) = init_nlpts + 1
                0184             REGRID_iend(i,bi,bj) = nlpts
                0185           ENDDO
c424ee7cc7 Jean*0186
aa076db465 Ed H*0187         ENDDO
                0188       ENDDO
                0189
                0190       WRITE(msgBuf,'(a)') ' '
                0191       CALL PRINT_MESSAGE(msgBuf,standardMessageUnit,
                0192      &     SQUEEZE_RIGHT,myThid)
                0193
                0194       _END_MASTER( myThid )
                0195
                0196       RETURN
                0197       END