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

Contents of /MITgcm/model/src/find_rho.F

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


Revision 1.8 - (show annotations) (download)
Fri Jun 12 19:33:33 1998 UTC (25 years, 10 months ago) by cnh
Branch: MAIN
CVS Tags: checkpoint11, checkpoint10, checkpoint13, checkpoint15, checkpoint14, checkpoint17, checkpoint7, checkpoint9, checkpoint8, checkpoint12, checkpoint16, branch-point-rdot
Branch point for: checkpoint7-4degree-ref, branch-rdot
Changes since 1.7: +2 -2 lines
Chages to make default setup correct for 4 degreee global comparison

1 C $Header: /u/gcmpack/models/MITgcmUV/model/src/find_rho.F,v 1.7 1998/06/09 16:48:02 cnh Exp $
2
3 #include "CPP_EEOPTIONS.h"
4 #define USE_FACTORIZED_POLY
5
6 ! ==============================================================================
7 subroutine FIND_RHO(
8 I bi, bj, iMin, iMax, jMin, jMax, k, kRef, eqn,
9 O rholoc,
10 I myThid )
11 C /==========================================================\
12 C | o SUBROUTINE FIND_RHO |
13 C | Calculates [rho(S,T,z)-Rhonil] of a slice |
14 C |==========================================================|
15 C | |
16 C | k - is the Theta/Salt level |
17 C | kRef - determines pressure reference level |
18 C | (not used in 'LINEAR' mode) |
19 C | eqn - determines the eqn. of state: 'LINEAR' or 'POLY3' |
20 C | |
21 C \==========================================================/
22 implicit none
23 ! Common
24 #include "SIZE.h"
25 #include "DYNVARS.h"
26 #include "EEPARAMS.h"
27 #include "PARAMS.h"
28 ! Arguments
29 integer bi,bj,iMin,iMax,jMin,jMax
30 integer k ! Level of Theta/Salt slice
31 integer kRef ! Pressure reference level
32 character*(*) eqn
33 _RL rholoc(1-Olx:sNx+Olx,1-Oly:sNy+Oly)
34 integer myThid
35 ! Local
36 integer i,j
37 _RL refTemp,refSalt,sigRef,tP,sP,deltaSig
38 ! ------------------------------------------------------------------------------
39
40 if (eqn.eq.'LINEAR') then
41
42 C ***NOTE***
43 C In the linear EOS, to make the static stability calculation meaningful
44 C we alway calculate the perturbation with respect to the surface level.
45 C **********
46 refTemp=tRef(kRef)
47 refSalt=sRef(kRef)
48
49 do j=jMin,jMax
50 do i=iMin,iMax
51 rholoc(i,j)=rhonil*(
52 & sBeta*( salt(i,j,k,bi,bj)-refSalt)
53 & -tAlpha*(theta(i,j,k,bi,bj)-refTemp) )
54 enddo
55 enddo
56
57 elseif (eqn.eq.'POLY3') then
58
59 refTemp=eosRefT(kRef)
60 refSalt=eosRefS(kRef)
61 sigRef=eosSig0(kRef) + (1000.-Rhonil)
62
63 do j=jMin,jMax
64 do i=iMin,iMax
65 tP=theta(i,j,k,bi,bj)-refTemp
66 sP=salt(i,j,k,bi,bj)-refSalt
67 #ifdef USE_FACTORIZED_POLY
68 deltaSig=
69 & (( eosC(9,kRef)*sP + eosC(5,kRef) )*sP + eosC(2,kRef) )*sP
70 & + ( ( eosC(6,kRef)
71 & *tP
72 & +eosC(7,kRef)*sP + eosC(3,kRef)
73 & )*tP
74 & +(eosC(8,kRef)*sP + eosC(4,kRef) )*sP + eosC(1,kRef)
75 & )*tP
76 #else
77 deltaSig=
78 & eosC(1,kRef)*tP
79 & +eosC(2,kRef) *sP
80 & +eosC(3,kRef)*tP*tP
81 & +eosC(4,kRef)*tP *sP
82 & +eosC(5,kRef) *sP*sP
83 & +eosC(6,kRef)*tP*tP*tP
84 & +eosC(7,kRef)*tP*tP *sP
85 & +eosC(8,kRef)*tP *sP*sP
86 & +eosC(9,kRef) *sP*sP*sP
87 #endif
88 rholoc(i,j)=sigRef+deltaSig
89 enddo
90 enddo
91
92 else
93 write(0,*) 'FIND_RHO: eqn = ',eqn
94 stop 'FIND_RHO: argument "eqn" has illegal value'
95 endif
96
97 ! ------------------------------------------------------------------------------
98 return
99 end
100 ! ==============================================================================

  ViewVC Help
Powered by ViewVC 1.1.22