/[MITgcm]/MITgcm_contrib/rpa_layers/layers/layers_init_fixed.F
ViewVC logotype

Annotation of /MITgcm_contrib/rpa_layers/layers/layers_init_fixed.F

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


Revision 1.1 - (hide annotations) (download)
Tue Sep 15 19:16:53 2009 UTC (15 years, 10 months ago) by rpa
Branch: MAIN
importing layers package

1 rpa 1.1 #include "LAYERS_OPTIONS.h"
2    
3     C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
4    
5     SUBROUTINE LAYERS_INIT_FIXED( myThid )
6    
7     C ===================================================================
8     C Initialize LAYERS variables that are kept fixed during the run.
9     C ===================================================================
10    
11     IMPLICIT NONE
12     #include "EEPARAMS.h"
13     #include "SIZE.h"
14     #include "PARAMS.h"
15     #include "GRID.h"
16     #include "LAYERS_SIZE.h"
17     #include "LAYERS.h"
18    
19     C INPUT/OUTPUT PARAMETERS:
20     C myThid :: my Thread Id number
21     INTEGER myThid
22    
23     C LOCAL VARIABLES:
24     C k :: loop index
25     C kk :: fine grid loop index
26     C Zf :: depth at cell boundaries
27     C Zf :: depth at cell centers
28     C ZZf :: depth at cell boundaries (fine grid)
29     C ZZc :: depth at cell centers (fine grid)
30     C msgBuf :: Informational/error meesage buffer
31     INTEGER k,kk
32     _RL Zf(Nr+1)
33     _RL Zc(Nr)
34     _RL ZZf(FineGridMax+1)
35     _RL ZZc(FineGridMax)
36    
37     CHARACTER*(MAX_LEN_MBUF) msgBuf
38    
39     #ifdef ALLOW_MNC
40     #ifdef LAYERS_MNC
41     IF (layers_MNC) THEN
42     CALL LAYERS_MNC_INIT( myThid )
43     ENDIF
44     #endif /* LAYERS_MNC */
45     #endif /* ALLOW_MNC */
46    
47     C Set up the vertical grid
48    
49     C for now, just use up the entire available array for ZZ
50     NZZ = FineGridMax
51    
52     C Z and ZZ are INCREASING DOWNWARD!!!
53     C Maybe this is dumb but it will work as long as we are consistent
54    
55     C find the depths
56     Zf(1) = 0. _d 0
57     Zc(1) = drC(1)
58     DO k=2,Nr
59     Zf(k) = Zf(k-1) + drF(k-1)
60     Zc(k) = Zc(k-1) + drC(k)
61     ENDDO
62     Zf(Nr+1) = Zf(Nr) + drF(Nr)
63    
64    
65     C calculate dZZ based on depth
66     dZZ = Zf(Nr+1) / NZZ
67    
68    
69     C create ZZ
70     ZZf(1) = 0. _d 0
71     ZZc(1) = 0.5 _d 0 * dZZ
72    
73     DO kk=2,NZZ
74     ZZf(kk) = ZZf(kk-1) + dZZ
75     ZZc(kk) = ZZc(kk-1) + dZZ
76     ENDDO
77     ZZf(NZZ+1) = ZZf(NZZ) + dZZ
78    
79     C create the interpolating mapping arrays
80     k = 1
81     DO kk=1,NZZ
82     C see if ZZc point is less than the top Zc point
83     IF ( ZZc(kk) .LT. Zc(1) ) THEN
84     MapIndex(kk) = 1
85     MapFact(kk) = 1.0 _d 0
86     C see if ZZc point is greater than the bottom Zc point
87     ELSE IF ( (ZZc(kk) .GE. Zc(Nr)) .OR. (k .EQ. Nr) ) THEN
88     MapIndex(kk) = Nr - 1
89     MapFact(kk) = 0.0 _d 0
90     C Otherwise we are somewhere in between and need to do interpolation)
91     ELSE IF ( (ZZc(kk) .GE. Zc(k))
92     & .AND. (ZZc(kk) .LT. Zc(Nr)) ) THEN
93     C Find the proper k value
94     DO WHILE (ZZc(kk) .GE. Zc(k+1))
95     k = k + 1
96     ENDDO
97     C If the loop stopped, that means Zc(k) <= ZZc(kk) < ZZc(k+1)
98     MapIndex(kk) = k
99     MapFact(kk) = 1.0 - (( ZZc(kk) - Zc(k) ) / drC(k+1))
100     ELSE
101     C This means there was a problem
102     WRITE(msgBuf,'(A,I4,A,I4,A,1E14.6,A,2E14.6)')
103     & 'S/R LAYERS_INIT_FIXED: kk=', kk, ' k=', k,
104     & 'ZZc(kk)=', ZZc(kk),' , Zc(k)=',Zc(k)
105     CALL PRINT_ERROR( msgBuf, myThid )
106     STOP 'ABNORMAL END: S/R LAYERS_INIT_FIXED'
107     END IF
108    
109     C See which grid box the point lies in
110     IF (ZZc(kk) < Zf(MapIndex(kk)+1)) THEN
111     CellIndex(kk) = MapIndex(kk)
112     ELSE
113     CellIndex(kk) = MapIndex(kk)+1
114     END IF
115     ENDDO
116    
117     IF (debugLevel .GT. 0) THEN
118     CALL PRINT_MESSAGE( 'LAYERS_INIT_FIXED Debugging:',
119     & standardMessageUnit,SQUEEZE_RIGHT , 1)
120     DO kk=1,NZZ
121     WRITE(msgBuf,'(A,1F6.1,A,I3,A,I3,A,I3,A,1F6.4)')
122     & '// ZZc=', ZZc(kk),
123     & ', MapIndex(',kk,')=',MapIndex(kk),
124     & ', MapFact(',kk,')=',MapFact(kk)
125     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
126     & SQUEEZE_RIGHT , 1)
127     ENDDO
128     END IF
129    
130     RETURN
131     END
132    
133    
134    
135    
136    
137    
138    
139    
140    
141    
142    
143    
144    
145    
146    
147    
148    
149    
150    
151    

  ViewVC Help
Powered by ViewVC 1.1.22