/[MITgcm]/MITgcm/model/src/ini_grid.F
ViewVC logotype

Diff of /MITgcm/model/src/ini_grid.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.7.2.2 by jmc, Wed May 16 20:54:49 2001 UTC revision 1.35 by jmc, Sun Feb 17 02:18:16 2013 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3    
4    #include "PACKAGES_CONFIG.h"
5  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
6    
7  CStartOfInterface  CBOP
8    C     !ROUTINE: INI_GRID
9    
10    C     !INTERFACE:
11        SUBROUTINE INI_GRID( myThid )        SUBROUTINE INI_GRID( myThid )
12  C     /==========================================================\  C     !DESCRIPTION:
13  C     | SUBROUTINE INI_GRID                                      |  C     These arrays are used throughout the code in evaluating gradients,
14  C     | o Initialise model grid                                  |  C     integrals and spatial avarages. This routine is called separately
15  C     |==========================================================|  C     by each thread and initializes only the region of the domain it is
16  C     | These arrays are used throughout the code in evaluating  |  C     "responsible" for.
17  C     | gradients, integrals and spatial avarages. This routine  |  
18  C     | is called separately by each thread and initialise only  |  C     !CALLING SEQUENCE:
19  C     | the region of the domain it is "responsible" for.        |  C     INI_GRID
20  C     | Notes:                                                   |  C      |   -- LOAD_GRID_SPACING
21  C     | Two examples are shown in this code. One illustrates the |  C      |   -- INI_VERTICAL_GRID
22  C     | initialisation of a cartesian grid. The other shows the  |  C      |    / INI_CARTESIAN_GRID
23  C     | inialisation of a spherical polar grid. Other orthonormal|  C      |   /  INI_SPHERICAL_POLAR_GRID
24  C     | grids can be fitted into this design. In this case       |  C      |   \  INI_CURVILINEAR_GRID
25  C     | custom metric terms also need adding to account for the  |  C      |    \ INI_CYLINDER_GRID
 C     | projections of velocity vectors onto these grids.        |  
 C     | The structure used here also makes it possible to        |  
 C     | implement less regular grid mappings. In particular      |  
 C     | o Schemes which leave out blocks of the domain that are  |  
 C     |   all land could be supported.                           |  
 C     | o Multi-level schemes such as icosohedral or cubic       |  
 C     |   grid projectedions onto a sphere can also be fitted    |  
 C     |   within the strategy we use.                            |  
 C     |   Both of the above also require modifying the support   |  
 C     |   routines that map computational blocks to simulation   |  
 C     |   domain blocks.                                         |  
 C     \==========================================================/  
       IMPLICIT NONE  
26    
27  C     === Global variables ===  C     !USES:
28          IMPLICIT NONE
29  #include "SIZE.h"  #include "SIZE.h"
30  #include "EEPARAMS.h"  #include "EEPARAMS.h"
31  #include "PARAMS.h"  #include "PARAMS.h"
32  #include "GRID.h"  #include "GRID.h"
33    #ifdef ALLOW_MNC
34    #include "MNC_PARAMS.h"
35    #endif
36    #ifdef ALLOW_MONITOR
37    #include "MONITOR.h"
38    #endif
39    
40  C     == Routine arguments ==  C     !INPUT/OUTPUT PARAMETERS:
41  C     myThid -  Number of this instance of INI_GRID  C     myThid  :: my Thread Id number
42        INTEGER myThid        INTEGER myThid
 CEndOfInterface  
43    
44  C     == Local variables ==  C     !FUNCTIONS:
45  C     msgBuf - Used for informational I/O.        LOGICAL  MASTER_CPU_IO
46          EXTERNAL MASTER_CPU_IO
47    
48    C     !LOCAL VARIABLES:
49    C     bi, bj  :: tile indices
50    C     i, j    :: Loop counters
51    C     msgBuf  :: Informational/error message buffer
52          INTEGER bi, bj
53          INTEGER i, j
54        CHARACTER*(MAX_LEN_MBUF) msgBuf        CHARACTER*(MAX_LEN_MBUF) msgBuf
55    CEOP
56    
57    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
58    
59    C--   Load grid spacing (vector) from files
60          CALL LOAD_GRID_SPACING( myThid )
61    
62  C--   Set up vertical grid and coordinate system  C--   Set up vertical grid and coordinate system
63        CALL INI_VERTICAL_GRID( myThid )        CALL INI_VERTICAL_GRID( myThid )
64    
65    C--   Initialise (everywhere) all horizontal grid array to null value
66    C     Note: some arrays are not defined in some parts of the halo
67    C     region. We set them to zero here for safety. If they are ever
68    C     referred to, especially in the denominator then it is a mistake!
69          DO bj = myByLo(myThid), myByHi(myThid)
70           DO bi = myBxLo(myThid), myBxHi(myThid)
71            DO j=1-OLy,sNy+OLy
72             DO i=1-OLx,sNx+OLx
73              xC(i,j,bi,bj)  = 0.
74              yC(i,j,bi,bj)  = 0.
75              xG(i,j,bi,bj)  = 0.
76              yG(i,j,bi,bj)  = 0.
77              dxC(i,j,bi,bj) = 0.
78              dyC(i,j,bi,bj) = 0.
79              dxG(i,j,bi,bj) = 0.
80              dyG(i,j,bi,bj) = 0.
81              dxF(i,j,bi,bj) = 0.
82              dyF(i,j,bi,bj) = 0.
83              dxV(i,j,bi,bj) = 0.
84              dyU(i,j,bi,bj) = 0.
85              rA(i,j,bi,bj)  = 0.
86              rAz(i,j,bi,bj) = 0.
87              rAw(i,j,bi,bj) = 0.
88              rAs(i,j,bi,bj) = 0.
89              recip_dxG(i,j,bi,bj) = 0.
90              recip_dyG(i,j,bi,bj) = 0.
91              recip_dxC(i,j,bi,bj) = 0.
92              recip_dyC(i,j,bi,bj) = 0.
93              recip_dxF(i,j,bi,bj) = 0.
94              recip_dyF(i,j,bi,bj) = 0.
95              recip_dxV(i,j,bi,bj) = 0.
96              recip_dyU(i,j,bi,bj) = 0.
97              recip_rA (i,j,bi,bj) = 0.
98              recip_rAs(i,j,bi,bj) = 0.
99              recip_rAw(i,j,bi,bj) = 0.
100              recip_rAz(i,j,bi,bj) = 0.
101              tanPhiAtU(i,j,bi,bj) = 0.
102              tanPhiAtV(i,j,bi,bj) = 0.
103              angleCosC(i,j,bi,bj) = 1.
104              angleSinC(i,j,bi,bj) = 0.
105              u2zonDir(i,j,bi,bj)  = 1.
106              v2zonDir(i,j,bi,bj)  = 0.
107             ENDDO
108             cosFacU(j,bi,bj)   = 1.
109             cosFacV(j,bi,bj)   = 1.
110             sqCosFacU(j,bi,bj) = 1.
111             sqCosFacV(j,bi,bj) = 1.
112            ENDDO
113           ENDDO
114          ENDDO
115    
116    C     Two examples are shown in this code. One illustrates the
117    C     initialization of a cartesian grid. The other shows the
118    C     inialization of a spherical polar grid. Other orthonormal grids
119    C     can be fitted into this design. In this case custom metric terms
120    C     also need adding to account for the projections of velocity
121    C     vectors onto these grids.  The structure used here also makes it
122    C     possible to implement less regular grid mappings. In particular:
123    C      o Schemes which leave out blocks of the domain that are
124    C        all land could be supported.
125    C      o Multi-level schemes such as icosohedral or cubic
126    C        grid projectedions onto a sphere can also be fitted
127    C       within the strategy we use.
128    C        Both of the above also require modifying the support
129    C        routines that map computational blocks to simulation
130    C        domain blocks.
131    
132  C--   Set up horizontal grid and coordinate system  C--   Set up horizontal grid and coordinate system
133        IF ( usingCartesianGrid ) THEN        IF ( usingCartesianGrid ) THEN
134         CALL INI_CARTESIAN_GRID( myThid )          CALL INI_CARTESIAN_GRID( myThid )
135        ELSEIF ( usingSphericalPolarGrid ) THEN        ELSEIF ( usingSphericalPolarGrid ) THEN
136         CALL INI_SPHERICAL_POLAR_GRID( myThid )          CALL INI_SPHERICAL_POLAR_GRID( myThid )
137        ELSEIF ( usingCurvilinearGrid ) THEN        ELSEIF ( usingCurvilinearGrid ) THEN
138         CALL INI_CURVILINEAR_GRID( myThid )          CALL INI_CURVILINEAR_GRID( myThid )
139          ELSEIF ( usingCylindricalGrid ) THEN
140            CALL INI_CYLINDER_GRID( myThid )
141        ELSE        ELSE
142         _BEGIN_MASTER(myThid)          _BEGIN_MASTER(myThid)
143          WRITE(msgBuf,'(A)')          WRITE(msgBuf,'(2A)') 'S/R INI_GRID: ',
144       &   'S/R INI_GRID: No grid coordinate system has been selected'       &       'No grid coordinate system has been selected'
145           CALL PRINT_ERROR( msgBuf , myThid)          CALL PRINT_ERROR( msgBuf , myThid)
146           STOP 'ABNORMAL END: S/R INI_GRID'          CALL ALL_PROC_DIE( 0 )
147         _END_MASTER(myThid)          STOP 'ABNORMAL END: S/R INI_GRID'
148            _END_MASTER(myThid)
149          ENDIF
150    
151    C--   Calculate reciprocals grid lengths (formerly part of INI_MASKS_ETC)
152          DO bj = myByLo(myThid), myByHi(myThid)
153           DO bi = myBxLo(myThid), myBxHi(myThid)
154            DO j=1-OLy,sNy+OLy
155             DO i=1-OLx,sNx+OLx
156              IF ( dxG(i,j,bi,bj) .NE. 0. )
157         &    recip_dxG(i,j,bi,bj) = 1. _d 0/dxG(i,j,bi,bj)
158              IF ( dyG(i,j,bi,bj) .NE. 0. )
159         &    recip_dyG(i,j,bi,bj) = 1. _d 0/dyG(i,j,bi,bj)
160              IF ( dxC(i,j,bi,bj) .NE. 0. )
161         &    recip_dxC(i,j,bi,bj) = 1. _d 0/dxC(i,j,bi,bj)
162              IF ( dyC(i,j,bi,bj) .NE. 0. )
163         &    recip_dyC(i,j,bi,bj) = 1. _d 0/dyC(i,j,bi,bj)
164              IF ( dxF(i,j,bi,bj) .NE. 0. )
165         &    recip_dxF(i,j,bi,bj) = 1. _d 0/dxF(i,j,bi,bj)
166              IF ( dyF(i,j,bi,bj) .NE. 0. )
167         &    recip_dyF(i,j,bi,bj) = 1. _d 0/dyF(i,j,bi,bj)
168              IF ( dxV(i,j,bi,bj) .NE. 0. )
169         &    recip_dxV(i,j,bi,bj) = 1. _d 0/dxV(i,j,bi,bj)
170              IF ( dyU(i,j,bi,bj) .NE. 0. )
171         &    recip_dyU(i,j,bi,bj) = 1. _d 0/dyU(i,j,bi,bj)
172              IF ( rA (i,j,bi,bj) .NE. 0. )
173         &    recip_rA (i,j,bi,bj) = 1. _d 0/rA (i,j,bi,bj)
174              IF ( rAs(i,j,bi,bj) .NE. 0. )
175         &    recip_rAs(i,j,bi,bj) = 1. _d 0/rAs(i,j,bi,bj)
176              IF ( rAw(i,j,bi,bj) .NE. 0. )
177         &    recip_rAw(i,j,bi,bj) = 1. _d 0/rAw(i,j,bi,bj)
178              IF ( rAz(i,j,bi,bj) .NE. 0. )
179         &    recip_rAz(i,j,bi,bj) = 1. _d 0/rAz(i,j,bi,bj)
180             ENDDO
181            ENDDO
182           ENDDO
183          ENDDO
184    
185    C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
186    
187    #ifdef ALLOW_MONITOR
188          IF ( MASTER_CPU_IO(myThid) ) THEN
189    C--   only the master thread is allowed to switch On/Off mon_write_stdout
190    C     & mon_write_mnc (since it is the only thread that uses those flags):
191    
192            IF (monitor_stdio) THEN
193              mon_write_stdout = .TRUE.
194            ELSE
195              mon_write_stdout = .FALSE.
196            ENDIF
197            mon_write_mnc = .FALSE.
198    #ifdef ALLOW_MNC
199            IF (useMNC .AND. monitor_mnc) THEN
200              DO i = 1,MAX_LEN_MBUF
201                mon_fname(i:i) = ' '
202              ENDDO
203              mon_fname(1:12) = 'monitor_grid'
204              CALL MNC_CW_SET_UDIM(mon_fname, 1, myThid)
205              mon_write_mnc = .TRUE.
206            ENDIF
207    #endif /*  ALLOW_MNC  */
208    
209          ENDIF
210    
211    C     Print out statistics of each horizontal grid array (helps when debugging)
212          CALL MON_PRINTSTATS_RS(1,xC,'XC',myThid)
213          CALL MON_PRINTSTATS_RS(1,xG,'XG',myThid)
214          CALL MON_PRINTSTATS_RS(1,dxC,'DXC',myThid)
215          CALL MON_PRINTSTATS_RS(1,dxF,'DXF',myThid)
216          CALL MON_PRINTSTATS_RS(1,dxG,'DXG',myThid)
217          CALL MON_PRINTSTATS_RS(1,dxV,'DXV',myThid)
218          CALL MON_PRINTSTATS_RS(1,yC,'YC',myThid)
219          CALL MON_PRINTSTATS_RS(1,yG,'YG',myThid)
220          CALL MON_PRINTSTATS_RS(1,dyC,'DYC',myThid)
221          CALL MON_PRINTSTATS_RS(1,dyF,'DYF',myThid)
222          CALL MON_PRINTSTATS_RS(1,dyG,'DYG',myThid)
223          CALL MON_PRINTSTATS_RS(1,dyU,'DYU',myThid)
224          CALL MON_PRINTSTATS_RS(1,rA,'RA',myThid)
225          CALL MON_PRINTSTATS_RS(1,rAw,'RAW',myThid)
226          CALL MON_PRINTSTATS_RS(1,rAs,'RAS',myThid)
227          CALL MON_PRINTSTATS_RS(1,rAz,'RAZ',myThid)
228          CALL MON_PRINTSTATS_RS(1,angleCosC,'AngleCS',myThid)
229          CALL MON_PRINTSTATS_RS(1,angleSinC,'AngleSN',myThid)
230    
231          IF ( MASTER_CPU_IO(myThid) ) THEN
232            mon_write_stdout = .FALSE.
233            mon_write_mnc    = .FALSE.
234        ENDIF        ENDIF
235    #endif /* ALLOW_MONITOR */
236    
237  C--   Write certain grid data to files (useful for creating netCDF  C--   Everyone else must wait for the grid to be set
238  C     and general post-analysis)        _BARRIER
       CALL WRITE_FLD_XY_RS( 'XC',' ',XC,0,myThid)  
       CALL WRITE_FLD_XY_RS( 'YC',' ',YC,0,myThid)  
       CALL WRITE_FLD_XY_RS( 'XG',' ',XG,0,myThid)  
       CALL WRITE_FLD_XY_RS( 'YG',' ',YG,0,myThid)  
       CALL WRITE_FLD_XY_RS( 'AC',' ',rA,0,myThid)  
239    
240        RETURN        RETURN
241        END        END

Legend:
Removed from v.1.7.2.2  
changed lines
  Added in v.1.35

  ViewVC Help
Powered by ViewVC 1.1.22