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

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

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

revision 1.17 by cnh, Sun Feb 4 16:46:44 2001 UTC revision 1.20 by adcroft, Thu Sep 27 18:14:52 2001 UTC
# Line 5  C $Name$ Line 5  C $Name$
5    
6  #undef  USE_BACKWARD_COMPATIBLE_GRID  #undef  USE_BACKWARD_COMPATIBLE_GRID
7    
8  CStartOfInterface  CBOP
9    C     !ROUTINE: INI_SPHERICAL_POLAR_GRID
10    C     !INTERFACE:
11        SUBROUTINE INI_SPHERICAL_POLAR_GRID( myThid )        SUBROUTINE INI_SPHERICAL_POLAR_GRID( myThid )
12    C     !DESCRIPTION: \bv
13  C     /==========================================================\  C     /==========================================================\
14  C     | SUBROUTINE INI_SPHERICAL_POLAR_GRID                      |  C     | SUBROUTINE INI_SPHERICAL_POLAR_GRID                      |
15  C     | o Initialise model coordinate system                     |  C     | o Initialise model coordinate system arrays              |
16  C     |==========================================================|  C     |==========================================================|
17  C     | These arrays are used throughout the code in evaluating  |  C     | These arrays are used throughout the code in evaluating  |
18  C     | gradients, integrals and spatial avarages. This routine  |  C     | gradients, integrals and spatial avarages. This routine  |
19  C     | is called separately by each thread and initialise only  |  C     | is called separately by each thread and initialise only  |
20  C     | the region of the domain it is "responsible" for.        |  C     | the region of the domain it is "responsible" for.        |
 C     | Notes:                                                   |  
 C     | Two examples are included. One illustrates the           |  
 C     | initialisation of a cartesian grid. The other shows the  |  
 C     | inialisation of a spherical polar grid. Other orthonormal|  
 C     | grids can be fitted into this design. In this case       |  
 C     | custom metric terms also need adding to account for the  |  
 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 projections 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.                                         |  
21  C     | Under the spherical polar grid mode primitive distances  |  C     | Under the spherical polar grid mode primitive distances  |
22  C     | in X and Y are in degrees. Distance in Z are in m or Pa  |  C     | in X and Y are in degrees. Distance in Z are in m or Pa  |
23  C     | depending on the vertical gridding mode.                 |  C     | depending on the vertical gridding mode.                 |
24  C     \==========================================================/  C     \==========================================================/
25        IMPLICIT NONE  C     \ev
26    
27    C     !USES:
28          IMPLICIT NONE
29  C     === Global variables ===  C     === Global variables ===
30  #include "SIZE.h"  #include "SIZE.h"
31  #include "EEPARAMS.h"  #include "EEPARAMS.h"
32  #include "PARAMS.h"  #include "PARAMS.h"
33  #include "GRID.h"  #include "GRID.h"
34    
35    C     !INPUT/OUTPUT PARAMETERS:
36  C     == Routine arguments ==  C     == Routine arguments ==
37  C     myThid -  Number of this instance of INI_CARTESIAN_GRID  C     myThid -  Number of this instance of INI_CARTESIAN_GRID
38        INTEGER myThid        INTEGER myThid
39  CEndOfInterface  CEndOfInterface
40    
41    C     !LOCAL VARIABLES:
42  C     == Local variables ==  C     == Local variables ==
43  C     xG, yG - Global coordinate location.  C     xG, yG - Global coordinate location.
44  C     xBase  - South-west corner location for process.  C     xBase  - South-west corner location for process.
# Line 72  C     I,J,K Line 62  C     I,J,K
62        _RL lat, dlat, dlon, xG0, yG0        _RL lat, dlat, dlon, xG0, yG0
63    
64    
65  C "Long" real for temporary coordinate calculation  C     "Long" real for temporary coordinate calculation
66  C  NOTICE the extended range of indices!!  C      NOTICE the extended range of indices!!
67        _RL xGloc(1-Olx:sNx+Olx+1,1-Oly:sNy+Oly+1)        _RL xGloc(1-Olx:sNx+Olx+1,1-Oly:sNy+Oly+1)
68        _RL yGloc(1-Olx:sNx+Olx+1,1-Oly:sNy+Oly+1)        _RL yGloc(1-Olx:sNx+Olx+1,1-Oly:sNy+Oly+1)
69    
70  C These functions return the "global" index with valid values beyond  C     The functions iGl, jGl return the "global" index with valid values beyond
71  C halo regions  C     halo regions
72    C     cnh wrote:
73    C     >    I dont understand why we would ever have to multiply the
74    C     >    overlap by the total domain size e.g
75    C     >    OLx*Nx, OLy*Ny.
76    C     >    Can anybody explain? Lines are in ini_spherical_polar_grid.F.
77    C     >    Surprised the code works if its wrong, so I am puzzled.        
78    C     jmc replied:
79    C     Yes, I can explain this since I put this modification to work
80    C     with small domain (where Oly > Ny, as for instance, zonal-average
81    C     case):
82    C     This has no effect on the acuracy of the evaluation of iGl(I,bi)
83    C     and jGl(J,bj) since we take mod(a+OLx*Nx,Nx) and mod(b+OLy*Ny,Ny).
84    C     But in case a or b is negative, then the FORTRAN function "mod"
85    C     does not return the matematical value of the "modulus" function,
86    C     and this is not good for your purpose.
87    C     This is why I add +OLx*Nx and +OLy*Ny to be sure that the 1rst
88    C     argument of the mod function is positive.
89        INTEGER iGl,jGl        INTEGER iGl,jGl
90        iGl(I,bi) = 1+mod(myXGlobalLo-1+(bi-1)*sNx+I+Olx*Nx-1,Nx)        iGl(I,bi) = 1+mod(myXGlobalLo-1+(bi-1)*sNx+I+Olx*Nx-1,Nx)
91        jGl(J,bj) = 1+mod(myYGlobalLo-1+(bj-1)*sNy+J+Oly*Ny-1,Ny)        jGl(J,bj) = 1+mod(myYGlobalLo-1+(bj-1)*sNy+J+Oly*Ny-1,Ny)
92    CEOP
93    
94    
95  C     For each tile ...  C     For each tile ...
96        DO bj = myByLo(myThid), myByHi(myThid)        DO bj = myByLo(myThid), myByHi(myThid)
# Line 176  C         by formula Line 185  C         by formula
185            dlon = delX( iGl(I,bi) )            dlon = delX( iGl(I,bi) )
186            dlat = delY( jGl(J,bj) )            dlat = delY( jGl(J,bj) )
187            dXG(I,J,bi,bj) = rSphere*COS(deg2rad*lat)*dlon*deg2rad            dXG(I,J,bi,bj) = rSphere*COS(deg2rad*lat)*dlon*deg2rad
188              if (dXG(I,J,bi,bj).LT.1.) dXG(I,J,bi,bj)=0.
189            dYG(I,J,bi,bj) = rSphere*dlat*deg2rad            dYG(I,J,bi,bj) = rSphere*dlat*deg2rad
190           ENDDO           ENDDO
191          ENDDO          ENDDO
# Line 309  C--     Calculate trigonometric terms Line 319  C--     Calculate trigonometric terms
319           ENDDO           ENDDO
320          ENDDO          ENDDO
321    
322    C--     Cosine(lat) scaling
323            DO J=1-OLy,sNy+OLy
324             jG = myYGlobalLo + (bj-1)*sNy + J-1
325             jG = min(max(1,jG),Ny)
326             IF (cosPower.NE.0.) THEN
327              cosFacU(J,bi,bj)=COS(yC(1,J,bi,bj)*deg2rad)
328         &                    **cosPower
329              cosFacV(J,bi,bj)=COS((yC(1,J,bi,bj)-0.5*delY(jG))*deg2rad)
330         &                    **cosPower
331              sqcosFacU(J,bi,bj)=sqrt(cosFacU(J,bi,bj))
332              sqcosFacV(J,bi,bj)=sqrt(cosFacV(J,bi,bj))
333             ELSE
334              cosFacU(J,bi,bj)=1.
335              cosFacV(J,bi,bj)=1.
336              sqcosFacU(J,bi,bj)=1.
337              sqcosFacV(J,bi,bj)=1.
338             ENDIF
339            ENDDO
340    
341         ENDDO ! bi         ENDDO ! bi
342        ENDDO ! bj        ENDDO ! bj
343    

Legend:
Removed from v.1.17  
changed lines
  Added in v.1.20

  ViewVC Help
Powered by ViewVC 1.1.22