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

Annotation of /MITgcm/model/src/find_alpha.F

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


Revision 1.16 - (hide annotations) (download)
Tue May 1 04:09:36 2007 UTC (17 years, 1 month ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint60, checkpoint61, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint59e, checkpoint59d, checkpoint59g, checkpoint59f, checkpoint59a, checkpoint59c, checkpoint59b, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59h, checkpoint59k, checkpoint59j, checkpoint61b, checkpoint61a
Changes since 1.15: +13 -9 lines
 add myThid to all kpp routines (long overdue)

1 mlosch 1.16 C $Header: /u/gcmpack/MITgcm/model/src/find_alpha.F,v 1.15 2003/10/27 22:32:55 heimbach Exp $
2 cnh 1.5 C $Name: $
3 adcroft 1.1
4     #include "CPP_OPTIONS.h"
5     #define USE_FACTORIZED_POLY
6    
7 cnh 1.5 CBOP
8     C !ROUTINE: FIND_ALPHA
9     C !INTERFACE:
10     SUBROUTINE FIND_ALPHA (
11 mlosch 1.10 I bi, bj, iMin, iMax, jMin, jMax, k, kRef,
12 mlosch 1.16 O alphaloc,
13     I myThid )
14 adcroft 1.1
15 cnh 1.5 C !DESCRIPTION: \bv
16     C *==========================================================*
17     C | o SUBROUTINE FIND_ALPHA
18     C | Calculates [drho(S,T,z) / dT] of a horizontal slice
19     C *==========================================================*
20     C |
21     C | k - is the Theta/Salt level
22     C | kRef - determines pressure reference level
23     C | (not used in 'LINEAR' mode)
24     C |
25     C | alphaloc - drho / dT (kg/m^3/C)
26     C |
27     C *==========================================================*
28     C \ev
29    
30     C !USES:
31     IMPLICIT NONE
32 jmc 1.12 C === Global variables ===
33 adcroft 1.1 #include "SIZE.h"
34     #include "DYNVARS.h"
35     #include "EEPARAMS.h"
36     #include "PARAMS.h"
37 mlosch 1.6 #include "EOS.h"
38     #include "GRID.h"
39 adcroft 1.1
40 cnh 1.5 C !INPUT/OUTPUT PARAMETERS:
41 jmc 1.12 C == Routine arguments ==
42 mlosch 1.16 C k :: Level of Theta/Salt slice
43     C kRef :: Pressure reference level
44     c myThid :: thread number for this instance of the routine
45     INTEGER myThid
46 jmc 1.12 INTEGER bi,bj,iMin,iMax,jMin,jMax
47     INTEGER k
48     INTEGER kRef
49 adcroft 1.1 _RL alphaloc(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
50    
51 cnh 1.5 C !LOCAL VARIABLES:
52 jmc 1.12 C == Local variables ==
53     INTEGER i,j
54 adcroft 1.1 _RL refTemp,refSalt,tP,sP
55 mlosch 1.9 _RL t1, t2, t3, s1, s3o2, p1, p2, sp5, p1t1
56 mlosch 1.6 _RL drhoP0dtheta, drhoP0dthetaFresh, drhoP0dthetaSalt
57     _RL dKdtheta, dKdthetaFresh, dKdthetaSalt, dKdthetaPres
58 jmc 1.12 _RL locPres(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
59     _RL rhoP0 (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
60 mlosch 1.6 _RL bulkMod(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
61 mlosch 1.9 _RL dnum_dtheta, dden_dtheta
62 jmc 1.12 _RL rhoDen (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
63     _RL rhoLoc (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
64 cnh 1.5 CEOP
65 adcroft 1.1
66 mlosch 1.11 #ifdef CHECK_SALINITY_FOR_NEGATIVE_VALUES
67     CALL LOOK_FOR_NEG_SALINITY( bi, bj, iMin, iMax, jMin, jMax, k,
68     & sFld, myThid )
69     #endif
70    
71 jmc 1.12 IF (equationOfState.EQ.'LINEAR') THEN
72 adcroft 1.1
73 jmc 1.12 DO j=jMin,jMax
74     DO i=iMin,iMax
75 adcroft 1.1 alphaloc(i,j) = -rhonil * tAlpha
76 jmc 1.12 ENDDO
77     ENDDO
78 adcroft 1.1
79 jmc 1.12 ELSEIF (equationOfState.EQ.'POLY3') THEN
80 adcroft 1.1
81     refTemp=eosRefT(kRef)
82     refSalt=eosRefS(kRef)
83    
84 jmc 1.12 DO j=jMin,jMax
85     DO i=iMin,iMax
86 adcroft 1.1 tP=theta(i,j,k,bi,bj)-refTemp
87     sP=salt(i,j,k,bi,bj)-refSalt
88     #ifdef USE_FACTORIZED_POLY
89     alphaloc(i,j) =
90     & ( eosC(6,kRef)
91     & *tP*3.
92     & +(eosC(7,kRef)*sP + eosC(3,kRef))*2.
93     & )*tP
94     & +(eosC(8,kRef)*sP + eosC(4,kRef) )*sP + eosC(1,kRef)
95     &
96     #else
97     alphaloc(i,j) =
98     & eosC(1,kRef) +
99     & eosC(3,kRef)*tP*2. +
100     & eosC(4,kRef) *sP +
101     & eosC(6,kRef)*tP*tP*3. +
102     & eosC(7,kRef)*tP*2. *sP +
103     & eosC(8,kRef) *sP*sP
104     #endif
105 jmc 1.12 ENDDO
106     ENDDO
107 adcroft 1.1
108 jmc 1.12 ELSEIF ( equationOfState(1:5).EQ.'JMD95'
109     & .OR. equationOfState.EQ.'UNESCO' ) THEN
110 mlosch 1.6 C nonlinear equation of state in pressure coordinates
111    
112 jmc 1.12 CALL PRESSURE_FOR_EOS(
113     I bi, bj, iMin, iMax, jMin, jMax, kRef,
114     O locPres,
115     I myThid )
116    
117 mlosch 1.6 CALL FIND_RHOP0(
118     I bi, bj, iMin, iMax, jMin, jMax, k,
119     I theta, salt,
120     O rhoP0,
121     I myThid )
122    
123     CALL FIND_BULKMOD(
124 jmc 1.12 I bi, bj, iMin, iMax, jMin, jMax, k,
125     I locPres, theta, salt,
126 mlosch 1.6 O bulkMod,
127     I myThid )
128    
129 jmc 1.12 DO j=jMin,jMax
130     DO i=iMin,iMax
131 mlosch 1.6
132     C abbreviations
133 mlosch 1.9 t1 = theta(i,j,k,bi,bj)
134     t2 = t1*t1
135     t3 = t2*t1
136 mlosch 1.6
137 mlosch 1.9 s1 = salt(i,j,k,bi,bj)
138 jmc 1.14 IF ( s1 .GT. 0. _d 0 ) THEN
139 jmc 1.12 s3o2 = SQRT(s1*s1*s1)
140 jmc 1.14 ELSE
141     s1 = 0. _d 0
142     s3o2 = 0. _d 0
143     ENDIF
144 mlosch 1.6
145 jmc 1.12 p1 = locPres(i,j)*SItoBar
146 mlosch 1.9 p2 = p1*p1
147 mlosch 1.6
148     C d(rho)/d(theta)
149     C of fresh water at p = 0
150     drhoP0dthetaFresh =
151     & eosJMDCFw(2)
152 mlosch 1.9 & + 2.*eosJMDCFw(3)*t1
153 mlosch 1.6 & + 3.*eosJMDCFw(4)*t2
154     & + 4.*eosJMDCFw(5)*t3
155 mlosch 1.9 & + 5.*eosJMDCFw(6)*t3*t1
156 mlosch 1.6 C of salt water at p = 0
157     drhoP0dthetaSalt =
158 mlosch 1.9 & s1*(
159 mlosch 1.6 & eosJMDCSw(2)
160 mlosch 1.9 & + 2.*eosJMDCSw(3)*t1
161 mlosch 1.6 & + 3.*eosJMDCSw(4)*t2
162     & + 4.*eosJMDCSw(5)*t3
163     & )
164     & + s3o2*(
165     & + eosJMDCSw(7)
166 mlosch 1.9 & + 2.*eosJMDCSw(8)*t1
167 mlosch 1.6 & )
168     C d(bulk modulus)/d(theta)
169     C of fresh water at p = 0
170     dKdthetaFresh =
171     & eosJMDCKFw(2)
172 mlosch 1.9 & + 2.*eosJMDCKFw(3)*t1
173 mlosch 1.6 & + 3.*eosJMDCKFw(4)*t2
174     & + 4.*eosJMDCKFw(5)*t3
175     C of sea water at p = 0
176     dKdthetaSalt =
177 mlosch 1.9 & s1*( eosJMDCKSw(2)
178     & + 2.*eosJMDCKSw(3)*t1
179 mlosch 1.6 & + 3.*eosJMDCKSw(4)*t2
180     & )
181     & + s3o2*( eosJMDCKSw(6)
182 mlosch 1.9 & + 2.*eosJMDCKSw(7)*t1
183 mlosch 1.6 & )
184     C of sea water at p
185     dKdthetaPres =
186 mlosch 1.9 & p1*( eosJMDCKP(2)
187     & + 2.*eosJMDCKP(3)*t1
188 mlosch 1.6 & + 3.*eosJMDCKP(4)*t2
189     & )
190 mlosch 1.9 & + p1*s1*( eosJMDCKP(6)
191     & + 2.*eosJMDCKP(7)*t1
192 mlosch 1.6 & )
193     & + p2*( eosJMDCKP(10)
194 mlosch 1.9 & + 2.*eosJMDCKP(11)*t1
195 mlosch 1.6 & )
196 mlosch 1.9 & + p2*s1*( eosJMDCKP(13)
197     & + 2.*eosJMDCKP(14)*t1
198 mlosch 1.6 & )
199    
200     drhoP0dtheta = drhoP0dthetaFresh
201     & + drhoP0dthetaSalt
202     dKdtheta = dKdthetaFresh
203     & + dKdthetaSalt
204     & + dKdthetaPres
205     alphaloc(i,j) =
206     & ( bulkmod(i,j)**2*drhoP0dtheta
207 mlosch 1.9 & - bulkmod(i,j)*p1*drhoP0dtheta
208     & - rhoP0(i,j)*p1*dKdtheta )
209     & /( bulkmod(i,j) - p1 )**2
210 mlosch 1.6
211    
212 jmc 1.12 ENDDO
213     ENDDO
214     ELSEIF ( equationOfState.EQ.'MDJWF' ) THEN
215    
216     CALL PRESSURE_FOR_EOS(
217     I bi, bj, iMin, iMax, jMin, jMax, kRef,
218     O locPres,
219     I myThid )
220    
221 heimbach 1.15 c$$$ CALL FIND_RHO( bi, bj, iMin, iMax, jMin, jMax, k,
222     c$$$ & kRef, theta, salt, rhoLoc, myThid )
223    
224     CALL FIND_RHONUM( bi, bj, iMin, iMax, jMin, jMax, k,
225     & locPres, theta, salt, rhoLoc, myThid )
226 heimbach 1.13
227 jmc 1.12 CALL FIND_RHODEN( bi, bj, iMin, iMax, jMin, jMax, k,
228     & locPres, theta, salt, rhoDen, myThid )
229 mlosch 1.9
230 jmc 1.12 DO j=jMin,jMax
231     DO i=iMin,iMax
232 mlosch 1.9 t1 = theta(i,j,k,bi,bj)
233     t2 = t1*t1
234     s1 = salt(i,j,k,bi,bj)
235 jmc 1.14 IF ( s1 .GT. 0. _d 0 ) THEN
236 jmc 1.12 sp5 = SQRT(s1)
237 jmc 1.14 ELSE
238     s1 = 0. _d 0
239     sp5 = 0. _d 0
240     ENDIF
241 mlosch 1.9
242 jmc 1.12 p1 = locPres(i,j)*SItodBar
243 mlosch 1.9 p1t1 = p1*t1
244    
245     dnum_dtheta = eosMDJWFnum(1)
246     & + t1*(2.*eosMDJWFnum(2) + 3.*eosMDJWFnum(3)*t1)
247     & + eosMDJWFnum(5)*s1
248 jmc 1.14 & + p1t1*(2.*eosMDJWFnum(8) + 2.*eosMDJWFnum(11)*p1)
249 mlosch 1.9
250     dden_dtheta = eosMDJWFden(1)
251     & + t1*(2.*eosMDJWFden(2)
252     & + t1*(3.*eosMDJWFden(3)
253     & + 4.*eosMDJWFden(4)*t1 ) )
254     & + s1*(eosMDJWFden(6)
255     & + t1*(3.*eosMDJWFden(7)*t1
256     & + 2.*eosMDJWFden(9)*sp5 ) )
257     & + p1*p1*(3.*eosMDJWFden(11)*t2 + eosMDJWFden(12)*p1)
258    
259     alphaLoc(i,j) = rhoDen(i,j)*(dnum_dtheta
260 heimbach 1.15 & - (rhoLoc(i,j)*rhoDen(i,j))*dden_dtheta)
261 mlosch 1.9
262 jmc 1.12 ENDDO
263     ENDDO
264 mlosch 1.9
265 jmc 1.12 ELSE
266     WRITE(*,*) 'FIND_ALPHA: equationOfState = ',equationOfState
267     STOP 'FIND_ALPHA: "equationOfState" has illegal value'
268     ENDIF
269 adcroft 1.1
270 jmc 1.12 RETURN
271     END
272 adcroft 1.1
273 jmc 1.12 SUBROUTINE FIND_BETA (
274 mlosch 1.10 I bi, bj, iMin, iMax, jMin, jMax, k, kRef,
275 mlosch 1.16 O betaloc,
276     I myThid )
277 adcroft 1.1 C /==========================================================\
278     C | o SUBROUTINE FIND_BETA |
279     C | Calculates [drho(S,T,z) / dS] of a horizontal slice |
280     C |==========================================================|
281     C | |
282     C | k - is the Theta/Salt level |
283     C | kRef - determines pressure reference level |
284     C | (not used in 'LINEAR' mode) |
285     C | |
286     C | betaloc - drho / dS (kg/m^3/PSU) |
287     C | |
288     C \==========================================================/
289 jmc 1.12 IMPLICIT NONE
290 adcroft 1.1
291 jmc 1.12 C === Global variables ===
292 adcroft 1.1 #include "SIZE.h"
293     #include "DYNVARS.h"
294     #include "EEPARAMS.h"
295     #include "PARAMS.h"
296 mlosch 1.6 #include "EOS.h"
297     #include "GRID.h"
298 adcroft 1.1
299 jmc 1.12 C == Routine arguments ==
300 mlosch 1.16 C k :: Level of Theta/Salt slice
301     C kRef :: Pressure reference level
302     c myThid :: thread number for this instance of the routine
303     INTEGER myThid
304 jmc 1.12 INTEGER bi,bj,iMin,iMax,jMin,jMax
305     INTEGER k
306     INTEGER kRef
307 adcroft 1.1 _RL betaloc(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
308    
309 jmc 1.12 C == Local variables ==
310     INTEGER i,j
311 adcroft 1.1 _RL refTemp,refSalt,tP,sP
312 mlosch 1.9 _RL t1, t2, t3, s1, s3o2, p1, sp5, p1t1
313 mlosch 1.6 _RL drhoP0dS
314     _RL dKdS, dKdSSalt, dKdSPres
315 jmc 1.12 _RL locPres(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
316     _RL rhoP0 (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
317 mlosch 1.6 _RL bulkMod(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
318 mlosch 1.9 _RL dnum_dsalt, dden_dsalt
319 jmc 1.12 _RL rhoDen (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
320     _RL rhoLoc (1-Olx:sNx+Olx,1-Oly:sNy+Oly)
321 mlosch 1.6 CEOP
322    
323 mlosch 1.11 #ifdef CHECK_SALINITY_FOR_NEGATIVE_VALUES
324     CALL LOOK_FOR_NEG_SALINITY( bi, bj, iMin, iMax, jMin, jMax, k,
325     & sFld, myThid )
326     #endif
327    
328 jmc 1.12 IF (equationOfState.EQ.'LINEAR') THEN
329 adcroft 1.1
330 jmc 1.12 DO j=jMin,jMax
331     DO i=iMin,iMax
332 adcroft 1.1 betaloc(i,j) = rhonil * sBeta
333 jmc 1.12 ENDDO
334     ENDDO
335 adcroft 1.1
336 jmc 1.12 ELSEIF (equationOfState.EQ.'POLY3') THEN
337 adcroft 1.1
338     refTemp=eosRefT(kRef)
339     refSalt=eosRefS(kRef)
340    
341 jmc 1.12 DO j=jMin,jMax
342     DO i=iMin,iMax
343 adcroft 1.1 tP=theta(i,j,k,bi,bj)-refTemp
344     sP=salt(i,j,k,bi,bj)-refSalt
345     #ifdef USE_FACTORIZED_POLY
346     betaloc(i,j) =
347     & ( eosC(9,kRef)*sP*3. + eosC(5,kRef)*2. )*sP + eosC(2,kRef)
348     & + ( eosC(7,kRef)*tP
349     & +eosC(8,kRef)*sP*2. + eosC(4,kRef)
350     & )*tP
351     #else
352     betaloc(i,j) =
353     & eosC(2,kRef) +
354     & eosC(4,kRef)*tP +
355     & eosC(5,kRef) *sP*2. +
356     & eosC(7,kRef)*tP*tP +
357     & eosC(8,kRef)*tP *sP*2. +
358     & eosC(9,kRef) *sP*sP*3.
359     #endif
360 jmc 1.12 ENDDO
361     ENDDO
362 adcroft 1.1
363 jmc 1.12 ELSEIF ( equationOfState(1:5).EQ.'JMD95'
364     & .OR. equationOfState.EQ.'UNESCO' ) THEN
365 mlosch 1.6 C nonlinear equation of state in pressure coordinates
366    
367 jmc 1.12 CALL PRESSURE_FOR_EOS(
368     I bi, bj, iMin, iMax, jMin, jMax, kRef,
369     O locPres,
370     I myThid )
371    
372 mlosch 1.6 CALL FIND_RHOP0(
373     I bi, bj, iMin, iMax, jMin, jMax, k,
374     I theta, salt,
375     O rhoP0,
376     I myThid )
377    
378     CALL FIND_BULKMOD(
379 jmc 1.12 I bi, bj, iMin, iMax, jMin, jMax, k,
380     I locPres, theta, salt,
381 mlosch 1.6 O bulkMod,
382     I myThid )
383    
384 jmc 1.12 DO j=jMin,jMax
385     DO i=iMin,iMax
386 mlosch 1.6
387     C abbreviations
388 mlosch 1.9 t1 = theta(i,j,k,bi,bj)
389     t2 = t1*t1
390     t3 = t2*t1
391 mlosch 1.6
392 mlosch 1.9 s1 = salt(i,j,k,bi,bj)
393 jmc 1.14 IF ( s1 .GT. 0. _d 0 ) THEN
394 jmc 1.12 s3o2 = 1.5*SQRT(s1)
395 jmc 1.14 ELSE
396     s1 = 0. _d 0
397     s3o2 = 0. _d 0
398     ENDIF
399 mlosch 1.6
400 jmc 1.12 p1 = locPres(i,j)*SItoBar
401 mlosch 1.6
402     C d(rho)/d(S)
403     C of fresh water at p = 0
404     drhoP0dS = 0. _d 0
405     C of salt water at p = 0
406     drhoP0dS = drhoP0dS
407     & + eosJMDCSw(1)
408 mlosch 1.9 & + eosJMDCSw(2)*t1
409 mlosch 1.6 & + eosJMDCSw(3)*t2
410     & + eosJMDCSw(4)*t3
411 mlosch 1.9 & + eosJMDCSw(5)*t3*t1
412 mlosch 1.6 & + s3o2*(
413     & eosJMDCSw(6)
414 mlosch 1.9 & + eosJMDCSw(7)*t1
415 mlosch 1.6 & + eosJMDCSw(8)*t2
416     & )
417 mlosch 1.9 & + 2*eosJMDCSw(9)*s1
418 mlosch 1.6 C d(bulk modulus)/d(S)
419     C of fresh water at p = 0
420     dKdS = 0. _d 0
421     C of sea water at p = 0
422     dKdSSalt =
423     & eosJMDCKSw(1)
424 mlosch 1.9 & + eosJMDCKSw(2)*t1
425 mlosch 1.6 & + eosJMDCKSw(3)*t2
426     & + eosJMDCKSw(4)*t3
427     & + s3o2*( eosJMDCKSw(5)
428 mlosch 1.9 & + eosJMDCKSw(6)*t1
429 mlosch 1.6 & + eosJMDCKSw(7)*t2
430     & )
431    
432     C of sea water at p
433     dKdSPres =
434 mlosch 1.9 & p1*( eosJMDCKP(5)
435     & + eosJMDCKP(6)*t1
436 mlosch 1.6 & + eosJMDCKP(7)*t2
437     & )
438 mlosch 1.9 & + s3o2*p1*eosJMDCKP(8)
439     & + p1*p1*( eosJMDCKP(12)
440     & + eosJMDCKP(13)*t1
441 mlosch 1.6 & + eosJMDCKP(14)*t2
442     & )
443    
444     dKdS = dKdSSalt + dKdSPres
445    
446     betaloc(i,j) =
447     & ( bulkmod(i,j)**2*drhoP0dS
448 mlosch 1.9 & - bulkmod(i,j)*p1*drhoP0dS
449     & - rhoP0(i,j)*p1*dKdS )
450     & /( bulkmod(i,j) - p1 )**2
451 mlosch 1.6
452    
453 jmc 1.12 ENDDO
454     ENDDO
455     ELSEIF ( equationOfState.EQ.'MDJWF' ) THEN
456    
457     CALL PRESSURE_FOR_EOS(
458     I bi, bj, iMin, iMax, jMin, jMax, kRef,
459     O locPres,
460     I myThid )
461    
462 heimbach 1.15 c$$$ CALL FIND_RHO( bi, bj, iMin, iMax, jMin, jMax, k,
463     c$$$ & kRef, theta, salt, rhoLoc, myThid )
464    
465     CALL FIND_RHONUM( bi, bj, iMin, iMax, jMin, jMax, k,
466     & locPres, theta, salt, rhoLoc, myThid )
467 heimbach 1.13
468 jmc 1.12 CALL FIND_RHODEN( bi, bj, iMin, iMax, jMin, jMax, k,
469     & locPres, theta, salt, rhoDen, myThid )
470 mlosch 1.9
471 jmc 1.12 DO j=jMin,jMax
472     DO i=iMin,iMax
473 mlosch 1.9 t1 = theta(i,j,k,bi,bj)
474     t2 = t1*t1
475     s1 = salt(i,j,k,bi,bj)
476 jmc 1.14 IF ( s1 .GT. 0. _d 0 ) THEN
477 jmc 1.12 sp5 = SQRT(s1)
478 jmc 1.14 ELSE
479     s1 = 0. _d 0
480     sp5 = 0. _d 0
481     ENDIF
482 mlosch 1.9
483 jmc 1.12 p1 = locPres(i,j)*SItodBar
484 mlosch 1.9 p1t1 = p1*t1
485    
486     dnum_dsalt = eosMDJWFnum(4)
487     & + eosMDJWFnum(5)*t1
488     & + 2.*eosMDJWFnum(6)*s1 + eosMDJWFnum(9)*p1
489     dden_dsalt = eosMDJWFden(5)
490     & + t1*( eosMDJWFden(6) + eosMDJWFden(7)*t2 )
491     & + 1.5*sp5*(eosMDJWFden(8) + eosMDJWFden(9)*t2)
492    
493     betaLoc(i,j) = rhoDen(i,j)*( dnum_dsalt
494 heimbach 1.15 & - (rhoLoc(i,j)*rhoDen(i,j))*dden_dsalt )
495 mlosch 1.9
496 jmc 1.12 ENDDO
497     ENDDO
498 mlosch 1.9
499 jmc 1.12 ELSE
500     WRITE(*,*) 'FIND_BETA: equationOfState = ',equationOfState
501     STOP 'FIND_BETA: "equationOfState" has illegal value'
502     ENDIF
503 adcroft 1.1
504 jmc 1.12 RETURN
505     END

  ViewVC Help
Powered by ViewVC 1.1.22