/[MITgcm]/MITgcm/verification/internal_wave/code/obcs_calc.F
ViewVC logotype

Annotation of /MITgcm/verification/internal_wave/code/obcs_calc.F

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


Revision 1.7 - (hide annotations) (download)
Mon Feb 28 16:31:36 2011 UTC (13 years, 1 month ago) by jmc
Branch: MAIN
CVS Tags: checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62x
Changes since 1.6: +35 -22 lines
add a call to S/R OBCS_BALANCE_FLOW.

1 jmc 1.7 C $Header: /u/gcmpack/MITgcm/verification/internal_wave/code/obcs_calc.F,v 1.6 2009/12/15 17:03:29 jahn Exp $
2 adcroft 1.3 C $Name: $
3 adcroft 1.2
4     #include "OBCS_OPTIONS.h"
5    
6 jahn 1.6 SUBROUTINE OBCS_CALC( futureTime, futureIter,
7 jmc 1.7 & uVel, vVel, wVel, theta, salt,
8 adcroft 1.2 & myThid )
9 jmc 1.7 C *==========================================================*
10     C | SUBROUTINE OBCS_CALC
11     C | o Calculate future boundary data at open boundaries
12     C | at time = futureTime
13     C *==========================================================*
14 adcroft 1.2 IMPLICIT NONE
15    
16     C === Global variables ===
17     #include "SIZE.h"
18     #include "EEPARAMS.h"
19     #include "PARAMS.h"
20 jmc 1.7 #include "GRID.h"
21     #include "OBCS.h"
22 mlosch 1.5 #include "EOS.h"
23 adcroft 1.2
24     C == Routine arguments ==
25 jmc 1.4 INTEGER futureIter
26 adcroft 1.2 _RL futureTime
27     _RL uVel (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
28     _RL vVel (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
29     _RL wVel (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
30     _RL theta(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
31     _RL salt (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
32     INTEGER myThid
33    
34     #ifdef ALLOW_OBCS
35    
36     C == Local variables ==
37 jahn 1.6 INTEGER bi, bj
38 adcroft 1.2 INTEGER I, J ,K
39     _RL obTimeScale,Uinflow,rampTime2
40     _RL vertStructWst(Nr)
41     _RL mz,strat,kx
42     _RL tmpsum
43    
44 jmc 1.7 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
45    
46     #ifdef ALLOW_DEBUG
47     IF (debugMode) CALL DEBUG_ENTER('OBCS_CALC',myThid)
48     #endif
49    
50 adcroft 1.2 C Vertical mode number
51 adcroft 1.3 mz=1.0 _d 0
52 adcroft 1.2 C Stratification
53     strat = 1.0 _d -6 / (gravity*tAlpha)
54    
55     C Create a vertical structure function with zero mean
56     tmpsum=0.
57     do K=1,Nr
58     vertStructWst(K)=cos(mz*PI* (rC(K)/rF(Nr+1)) )
59     tmpsum=tmpsum+vertStructWst(K)*drF(K)
60     enddo
61     tmpsum=tmpsum/rF(Nr+1)
62     do K=1,Nr
63     vertStructWst(K)=vertStructWst(K)-tmpsum
64     enddo
65     c
66 adcroft 1.3 obTimeScale = 44567.0 _d 0
67     kx=mz*2. _d 0*pi/400.0 _d 0
68     & *sqrt((2.0 _d 0*pi*2.0 _d 0*pi/(obTimeScale*obTimeScale)
69 adcroft 1.2 & - f0*f0)/(1.0 _d -6
70 adcroft 1.3 & - 2.0 _d 0*pi*2.0 _d 0*pi/(obTimeScale*obTimeScale)))
71     Uinflow = 0.024 _d 0
72     C *NOTE* I have commented out the ramp function below
73     C just to speed things up. You will probably want to use it
74     C for smoother looking solutions.
75     rampTime2 = 4. _d 0*44567.0 _d 0
76 adcroft 1.2
77 jahn 1.6 DO bj=myByLo(myThid),myByHi(myThid)
78     DO bi=myBxLo(myThid),myBxHi(myThid)
79 adcroft 1.2
80     C Eastern OB
81     IF (useOrlanskiEast) THEN
82     CALL ORLANSKI_EAST(
83 jmc 1.7 & bi, bj, futureTime,
84     & uVel, vVel, wVel, theta, salt,
85 adcroft 1.2 & myThid )
86     ELSE
87     DO K=1,Nr
88     DO J=1-Oly,sNy+Oly
89     OBEu(J,K,bi,bj)=0.
90     OBEv(J,K,bi,bj)=0.
91     OBEt(J,K,bi,bj)=tRef(K)
92     OBEs(J,K,bi,bj)=sRef(K)
93     #ifdef ALLOW_NONHYDROSTATIC
94     OBEw(J,K,bi,bj)=0.
95     #endif
96     ENDDO
97     ENDDO
98     ENDIF
99    
100     C Western OB
101     IF (useOrlanskiWest) THEN
102     CALL ORLANSKI_WEST(
103 jmc 1.7 & bi, bj, futureTime,
104     & uVel, vVel, wVel, theta, salt,
105 adcroft 1.2 & myThid )
106     ELSE
107     DO K=1,Nr
108     DO J=1-Oly,sNy+Oly
109 adcroft 1.3 OBWu(J,K,bi,bj)=0. _d 0
110 adcroft 1.2 & +Uinflow
111     & *vertStructWst(K)
112 adcroft 1.3 & *sin(2. _d 0*PI*futureTime/obTimeScale)
113     c & *(exp(futureTime/rampTime2)
114     c & - exp(-futureTime/rampTime2))
115     c & /(exp(futureTime/rampTime2)
116     c & + exp(-futureTime/rampTime2))
117     & *cos(kx*(3. _d 0-2. _d 0-0.5 _d 0)*delX(1))
118     OBWv(J,K,bi,bj)=0. _d 0
119 adcroft 1.2 & +Uinflow
120 adcroft 1.3 & *f0/(2.0 _d 0*PI/obTimeScale)
121 adcroft 1.2 & *vertStructWst(K)
122 adcroft 1.3 & *cos(2. _d 0*PI*futureTime/obTimeScale )
123 adcroft 1.2 & * (exp(futureTime/rampTime2)
124     & - exp(-futureTime/rampTime2))
125     & /(exp(futureTime/rampTime2)
126     & + exp(-futureTime/rampTime2))
127     OBWt(J,K,bi,bj)=tRef(K)
128 adcroft 1.3 & + Uinflow*sin(mz*PI*(float(k)-0.5 _d 0)/float(Nr))
129     & * sin(2.0 _d 0*PI*futureTime/obTimeScale)
130 adcroft 1.2 & *sqrt(strat/(tAlpha*gravity))
131 adcroft 1.3 & *sqrt(2.0 _d 0*PI/obTimeScale*2.0*PI/obTimeScale - f0*f0)
132     & /(2.0 _d 0*PI/obTimeScale)
133     c & * (exp(futureTime/rampTime2)
134     c & - exp(-futureTime/rampTime2))
135     c & /(exp(futureTime/rampTime2)
136     c & + exp(-futureTime/rampTime2))
137 adcroft 1.2 #ifdef ALLOW_NONHYDROSTATIC
138     OBWw(J,K,bi,bj)=-Uinflow
139 adcroft 1.3 & *sqrt(2.0 _d 0*PI/obTimeScale*2.0 _d 0*PI/obTimeScale - f0*f0)
140     & /sqrt(strat*strat -
141     & 2.0 _d 0*PI/obTimeScale*2.0 _d 0*PI/obTimeScale)
142     & *sin(mz*PI*(float(k)-0.5 _d 0)/float(Nr))
143     & *cos(2. _d 0*PI*futureTime/obTimeScale)
144     c & *(exp(futureTime/rampTime2)
145     c & - exp(-futureTime/rampTime2))
146     c & /(exp(futureTime/rampTime2)
147     c & + exp(-futureTime/rampTime2))
148 adcroft 1.2 #endif
149     ENDDO
150     ENDDO
151     ENDIF
152    
153     C Northern OB, template for forcing
154     IF (useOrlanskiNorth) THEN
155     CALL ORLANSKI_NORTH(
156 jmc 1.7 & bi, bj, futureTime,
157     & uVel, vVel, wVel, theta, salt,
158 adcroft 1.2 & myThid )
159     ELSE
160     DO K=1,Nr
161     DO I=1-Olx,sNx+Olx
162     OBNv(I,K,bi,bj)=0.
163     OBNu(I,K,bi,bj)=0.
164     OBNt(I,K,bi,bj)=tRef(K)
165     OBNs(I,K,bi,bj)=sRef(K)
166     #ifdef ALLOW_NONHYDROSTATIC
167     OBNw(I,K,bi,bj)=0.
168     #endif
169     ENDDO
170     ENDDO
171     ENDIF
172    
173     C Southern OB, template for forcing
174 jmc 1.7 IF (useOrlanskiSouth) THEN
175 adcroft 1.2 CALL ORLANSKI_SOUTH(
176 jmc 1.7 & bi, bj, futureTime,
177     & uVel, vVel, wVel, theta, salt,
178 adcroft 1.2 & myThid )
179     ELSE
180     DO K=1,Nr
181     DO I=1-Olx,sNx+Olx
182     OBSu(I,K,bi,bj)=0.
183     OBSv(I,K,bi,bj)=0.
184     OBSt(I,K,bi,bj)=tRef(K)
185     OBSs(I,K,bi,bj)=sRef(K)
186     #ifdef ALLOW_NONHYDROSTATIC
187     OBSw(I,K,bi,bj)=0.
188     #endif
189     ENDDO
190     ENDDO
191     ENDIF
192    
193 jahn 1.6 C-- end bi,bj loops.
194     ENDDO
195 jmc 1.7 ENDDO
196    
197     #ifdef ALLOW_OBCS_BALANCE
198     IF ( useOBCSbalance ) THEN
199     CALL OBCS_BALANCE_FLOW( futureTime, futureIter, myThid )
200     ENDIF
201     #endif /* ALLOW_OBCS_BALANCE */
202 jahn 1.6
203 jmc 1.7 #ifdef ALLOW_DEBUG
204     IF (debugMode) CALL DEBUG_LEAVE('OBCS_CALC',myThid)
205     #endif
206 adcroft 1.2 #endif /* ALLOW_OBCS */
207 jmc 1.7
208 adcroft 1.2 RETURN
209     END

  ViewVC Help
Powered by ViewVC 1.1.22