1 |
C $Header: /u/gcmpack/MITgcm/model/src/ini_vertical_grid.F,v 1.17 2006/11/29 04:39:06 jmc Exp $ |
2 |
C $Name: $ |
3 |
|
4 |
#include "CPP_OPTIONS.h" |
5 |
|
6 |
CBOP |
7 |
C !ROUTINE: INI_VERTICAL_GRID |
8 |
C !INTERFACE: |
9 |
SUBROUTINE INI_VERTICAL_GRID( myThid ) |
10 |
|
11 |
C !DESCRIPTION: \bv |
12 |
C *==========================================================* |
13 |
C | SUBROUTINE INI_VERTICAL_GRID |
14 |
C | o Initialise vertical gridding arrays |
15 |
C *==========================================================* |
16 |
C \ev |
17 |
|
18 |
C !USES: |
19 |
IMPLICIT NONE |
20 |
C === Global variables === |
21 |
#include "SIZE.h" |
22 |
#include "EEPARAMS.h" |
23 |
#include "PARAMS.h" |
24 |
#include "GRID.h" |
25 |
|
26 |
C !INPUT/OUTPUT PARAMETERS: |
27 |
C == Routine arguments == |
28 |
C myThid :: my Thread Id number |
29 |
INTEGER myThid |
30 |
|
31 |
C !LOCAL VARIABLES: |
32 |
C == Local variables == |
33 |
C k :: loop index |
34 |
C msgBuf :: Informational/error meesage buffer |
35 |
INTEGER k |
36 |
_RL tmpRatio, checkRatio1, checkRatio2 |
37 |
CHARACTER*(MAX_LEN_MBUF) msgBuf |
38 |
CEOP |
39 |
|
40 |
_BEGIN_MASTER(myThid) |
41 |
|
42 |
WRITE(msgBuf,'(A,2(A,L5))') 'Enter INI_VERTICAL_GRID:', |
43 |
& ' setInterFDr=', setInterFDr, |
44 |
& ' ; setCenterDr=', setCenterDr |
45 |
CALL PRINT_MESSAGE( msgBuf, standardMessageUnit, |
46 |
& SQUEEZE_RIGHT, myThid ) |
47 |
|
48 |
C-- Set factors required for mixing pressure and meters as vertical coordinate. |
49 |
C rkSign is a "sign" parameter which is used where the orientation of the vertical |
50 |
C coordinate (pressure or meters) relative to the vertical index (k) is important. |
51 |
C rkSign = -1 applies when k and the coordinate are in the opposite sense. |
52 |
C rkSign = 1 applies when k and the coordinate are in the same sense. |
53 |
rkSign = -1. _d 0 |
54 |
gravitySign = -1. _d 0 |
55 |
IF ( usingPCoords ) THEN |
56 |
gravitySign = 1. _d 0 |
57 |
ENDIF |
58 |
|
59 |
IF ( .NOT.(setInterFDr.OR.setCenterDr) ) THEN |
60 |
WRITE(msgBuf,'(A)') |
61 |
& 'S/R INI_VERTICAL_GRID: neither delR nor delRc are defined' |
62 |
CALL PRINT_ERROR( msgBuf, myThid ) |
63 |
WRITE(msgBuf,'(A)') |
64 |
& 'S/R INI_VERTICAL_GRID: Need at least 1 of the 2 (delR,delRc)' |
65 |
CALL PRINT_ERROR( msgBuf, myThid ) |
66 |
STOP 'ABNORMAL END: S/R INI_VERTICAL_GRID' |
67 |
ENDIF |
68 |
|
69 |
C--- Set Level r-thickness (drF) and Center r-distances (drC) |
70 |
|
71 |
IF (setInterFDr) THEN |
72 |
C-- Interface r-distances are defined: |
73 |
DO k=1,Nr |
74 |
drF(k) = delR(k) |
75 |
ENDDO |
76 |
C- Check that all thickness are > 0 : |
77 |
DO k=1,Nr |
78 |
IF (delR(k).LE.0.) THEN |
79 |
WRITE(msgBuf,'(A,I4,A,E16.8)') |
80 |
& 'S/R INI_VERTICAL_GRID: delR(k=',k,' )=',delR(k) |
81 |
CALL PRINT_ERROR( msgBuf, myThid ) |
82 |
WRITE(msgBuf,'(A)') |
83 |
& 'S/R INI_VERTICAL_GRID: Vert. grid spacing MUST BE > 0' |
84 |
CALL PRINT_ERROR( msgBuf, myThid ) |
85 |
STOP 'ABNORMAL END: S/R INI_VERTICAL_GRID' |
86 |
ENDIF |
87 |
ENDDO |
88 |
ELSE |
89 |
C-- Interface r-distances undefined: |
90 |
C assume Interface at middle between 2 Center |
91 |
drF(1) = delRc(1) |
92 |
DO k=2,Nr |
93 |
drF(k-1) = 0.5 _d 0 *delRc(k) + drF(k-1) |
94 |
drF( k ) = 0.5 _d 0 *delRc(k) |
95 |
ENDDO |
96 |
drF(Nr) = delRc(Nr+1) + drF(Nr) |
97 |
ENDIF |
98 |
|
99 |
IF (setCenterDr) THEN |
100 |
C-- Cell Center r-distances are defined: |
101 |
DO k=1,Nr |
102 |
drC(k) = delRc(k) |
103 |
ENDDO |
104 |
C- Check that all thickness are > 0 : |
105 |
DO k=1,Nr+1 |
106 |
IF (delRc(k).LE.0.) THEN |
107 |
WRITE(msgBuf,'(A,I4,A,E16.8)') |
108 |
& 'S/R INI_VERTICAL_GRID: delRc(k=',k,' )=',delRc(k) |
109 |
CALL PRINT_ERROR( msgBuf, myThid ) |
110 |
WRITE(msgBuf,'(A)') |
111 |
& 'S/R INI_VERTICAL_GRID: Vert. grid spacing MUST BE > 0' |
112 |
CALL PRINT_ERROR( msgBuf, myThid ) |
113 |
STOP 'ABNORMAL END: S/R INI_VERTICAL_GRID' |
114 |
ENDIF |
115 |
ENDDO |
116 |
ELSE |
117 |
C-- Cell Center r-distances undefined: |
118 |
C assume Center at middle between 2 Interfaces |
119 |
drC(1) = 0.5 _d 0 *delR(1) |
120 |
DO k=2,Nr |
121 |
drC(k) = 0.5 _d 0 *(delR(k-1)+delR(k)) |
122 |
ENDDO |
123 |
ENDIF |
124 |
|
125 |
C--- Set r-position of interFace (rF) and cell-Center (rC): |
126 |
rF(1) = Ro_SeaLevel |
127 |
DO k=1,Nr |
128 |
rF(k+1) = rF(k) + rkSign*drF(k) |
129 |
ENDDO |
130 |
rC(1) = rF(1) + rkSign*drC(1) |
131 |
DO k=2,Nr |
132 |
rC(k) = rC(k-1) + rkSign*drC(k) |
133 |
ENDDO |
134 |
|
135 |
C--- Check vertical discretization : |
136 |
checkRatio2 = 100. |
137 |
checkRatio1 = 1. _d 0 / checkRatio2 |
138 |
DO k=1,Nr |
139 |
tmpRatio = 0. |
140 |
IF ( (rC(k)-rF(k+1)) .NE. 0. ) |
141 |
& tmpRatio = (rF(k)-rC(k)) / (rC(k)-rF(k+1)) |
142 |
IF ( tmpRatio.LT.checkRatio1 .OR. tmpRatio.GT.checkRatio2 ) THEN |
143 |
c write(0,*) 'drF=',drF |
144 |
c write(0,*) 'drC=',drC |
145 |
c write(0,*) 'rF=',rF |
146 |
c write(0,*) 'rC=',rC |
147 |
WRITE(msgBuf,'(A,I4,A,E16.8)') |
148 |
& 'S/R INI_VERTICAL_GRID: Invalid relative position, level k=', |
149 |
& k, ' :', tmpRatio |
150 |
CALL PRINT_ERROR( msgBuf, myThid ) |
151 |
WRITE(msgBuf,'(A,1PE14.6,A,2E14.6)') |
152 |
& 'S/R INI_VERTICAL_GRID: rC=', rC(k), |
153 |
& ' , rF(k,k+1)=',rF(k),rF(k+1) |
154 |
CALL PRINT_ERROR( msgBuf, myThid ) |
155 |
STOP 'ABNORMAL END: S/R INI_VERTICAL_GRID' |
156 |
ENDIF |
157 |
ENDDO |
158 |
|
159 |
C- Calculate reciprol vertical grid spacing : |
160 |
DO k=1,Nr |
161 |
recip_drC(k) = 1. _d 0/drC(k) |
162 |
recip_drF(k) = 1. _d 0/drF(k) |
163 |
ENDDO |
164 |
|
165 |
C-- Calculate horizontal grid factor for the deep model (<=> dropping the |
166 |
C shallow atmosphere approximation): only function of the vertical index |
167 |
C- first: initialise deep-model grid factor: |
168 |
DO k=1,Nr |
169 |
deepFacC(k) = 1. _d 0 |
170 |
deepFac2C(k)= 1. _d 0 |
171 |
recip_deepFacC(k) = 1. _d 0 |
172 |
recip_deepFac2C(k)= 1. _d 0 |
173 |
ENDDO |
174 |
DO k=1,Nr+1 |
175 |
deepFacF(k) = 1. _d 0 |
176 |
deepFac2F(k)= 1. _d 0 |
177 |
recip_deepFacF(k) = 1. _d 0 |
178 |
recip_deepFac2F(k)= 1. _d 0 |
179 |
ENDDO |
180 |
IF ( deepAtmosphere ) THEN |
181 |
C- set deep-model grid factor: |
182 |
IF ( usingZCoords ) THEN |
183 |
DO k=1,Nr |
184 |
deepFacC(k) = (rSphere+rC(k))*recip_rSphere |
185 |
deepFac2C(k) = deepFacC(k)*deepFacC(k) |
186 |
ENDDO |
187 |
DO k=1,Nr+1 |
188 |
deepFacF(k) = (rSphere+rF(k))*recip_rSphere |
189 |
deepFac2F(k) = deepFacF(k)*deepFacF(k) |
190 |
ENDDO |
191 |
ELSE |
192 |
STOP 'INI_VERTICAL_GRID: setting deepFac is not coded' |
193 |
ENDIF |
194 |
C- set reciprocal of deep-model grid factor: |
195 |
DO k=1,Nr |
196 |
recip_deepFacC(k) = 1. _d 0/deepFacC(k) |
197 |
recip_deepFac2C(k)= 1. _d 0/deepFac2C(k) |
198 |
ENDDO |
199 |
DO k=1,Nr+1 |
200 |
recip_deepFacF(k) = 1. _d 0/deepFacF(k) |
201 |
recip_deepFac2F(k)= 1. _d 0/deepFac2F(k) |
202 |
ENDDO |
203 |
ENDIF |
204 |
|
205 |
_END_MASTER(myThid) |
206 |
_BARRIER |
207 |
|
208 |
RETURN |
209 |
END |