/[MITgcm]/MITgcm/pkg/layers/layers_init_fixed.F
ViewVC logotype

Annotation of /MITgcm/pkg/layers/layers_init_fixed.F

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


Revision 1.2 - (hide annotations) (download)
Fri Sep 18 16:13:33 2009 UTC (14 years, 8 months ago) by jmc
Branch: MAIN
Changes since 1.1: +2 -22 lines
fix syntax (was not passing with pgi)

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

  ViewVC Help
Powered by ViewVC 1.1.22