/[MITgcm]/MITgcm/pkg/atm2d/relax_add.F
ViewVC logotype

Contents of /MITgcm/pkg/atm2d/relax_add.F

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


Revision 1.5 - (show annotations) (download)
Fri Jan 20 20:32:57 2012 UTC (12 years, 4 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint64, checkpoint65, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, HEAD
Changes since 1.4: +5 -3 lines
multiply instead of dividing by recip. But should use rhoConst instead of rhoNil.

1 C $Header: /u/gcmpack/MITgcm/pkg/atm2d/relax_add.F,v 1.4 2010/06/16 21:04:22 jscott Exp $
2 C $Name: $
3
4 #include "ctrparam.h"
5 #include "ATM2D_OPTIONS.h"
6
7 C !INTERFACE:
8 SUBROUTINE RELAX_ADD( wght0, wght1,
9 & intime0, intime1, iftime, myIter, myThid)
10 C *==========================================================*
11 C | Adds restoring terms to surface forcing. Note that: |
12 C | - restoring is phased out as restor (or act.) SST <2C |
13 C | - if nsTypeRelax NE 0, salt rest. phased out nr poles |
14 C | - if ntTypeRelax NE 0, temp rest. phased out nr poles |
15 C *==========================================================*
16 IMPLICIT NONE
17
18 #include "ATMSIZE.h"
19 #include "SIZE.h"
20 #include "EEPARAMS.h"
21 #include "PARAMS.h"
22 #include "GRID.h"
23 #include "THSICE_VARS.h"
24 #include "ATM2D_VARS.h"
25
26 c include ocean and seaice vars
27
28 C !INPUT/OUTPUT PARAMETERS:
29 C === Routine arguments ===
30 C wght0, wght1 - weight of first and second month, respectively
31 C intime0,intime1- month id # for first and second months
32 C iftime - true -> prompts a reloading of data from disk
33 C myIter - Ocean iteration number
34 C myThid - Thread no. that called this routine.
35 _RL wght0
36 _RL wght1
37 INTEGER intime0
38 INTEGER intime1
39 LOGICAL iftime
40 INTEGER myIter
41 INTEGER myThid
42
43 C LOCAL VARIABLES:
44 C Save below so that continual file reloads aren't necessary
45 COMMON /OCEANRELAX/
46 & sst0, sst1, sss0, sss1
47
48 _RS sst0(1-Olx:sNx+Olx,1-Oly:sNy+Oly,1,1)
49 _RS sst1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,1,1)
50 _RS sss0(1-Olx:sNx+Olx,1-Oly:sNy+Oly,1,1)
51 _RS sss1(1-Olx:sNx+Olx,1-Oly:sNy+Oly,1,1)
52 _RL lambdaTheta,lambdaSalt
53 _RS nearIce ! constant used to phase out rest near frz point
54 _RL qrelflux, frelflux
55 _RL sstRelax(1:sNx,1:sNy) ! relaxation sst as computed from file
56 _RL sssRelax(1:sNx,1:sNy) ! relaxation sss as computed from file
57 INTEGER i,j
58
59 IF (ifTime) THEN
60
61 C If the above condition is met then we need to read in
62 C data for the period ahead and the period behind current time.
63
64 WRITE(*,*) 'S/R RELAX_ADD: Reading new data'
65 IF ( thetaRelaxFile .NE. ' ' ) THEN
66 CALL READ_REC_XY_RS( thetaRelaxFile,sst0,intime0,
67 & myIter,myThid )
68 CALL READ_REC_XY_RS( thetaRelaxFile,sst1,intime1,
69 & myIter,myThid )
70 ENDIF
71 IF ( saltRelaxFile .NE. ' ' ) THEN
72 CALL READ_REC_XY_RS( saltRelaxFile,sss0,intime0,
73 & myIter,myThid )
74 CALL READ_REC_XY_RS( saltRelaxFile,sss1,intime1,
75 & myIter,myThid )
76 ENDIF
77
78 ENDIF
79
80 IF ((thetaRelaxFile.NE.' ').OR.(saltRelaxFile.NE.' ')) THEN
81
82 C-- Interpolate and add to anomaly
83 DO j=1,sNy
84
85 IF (ntTypeRelax .EQ. 0) THEN
86 lambdaTheta = r_tauThetaRelax
87 ELSE
88 lambdaTheta = r_tauThetaRelax *
89 & max(cos(1.5 _d 0*yC(1,j,1,1)*deg2rad),0. _d 0)
90 ENDIF
91 IF (nsTypeRelax .EQ. 0) THEN
92 lambdaSalt = r_tauSaltRelax
93 ELSE
94 lambdaSalt = r_tauSaltRelax *
95 & max(cos(1.5 _d 0*yC(1,j,1,1)*deg2rad),0. _d 0)
96 ENDIF
97
98 DO i=1,sNx
99
100 IF (maskC(i,j,1,1,1) .EQ. 1.) THEN
101
102 IF (thetaRelaxFile.NE.' ') THEN
103 sstRelax(i,j)= (wght0*sst0(i,j,1,1) + wght1*sst1(i,j,1,1))
104 ELSE !no T restoring; use actual SST to determine if nr freezing
105 sstRelax(i,j)= sstFromOcn(i,j)
106 ENDIF
107
108 IF (saltRelaxFile.NE.' ') THEN
109 sssRelax(i,j)= (wght0*sss0(i,j,1,1) + wght1*sss1(i,j,1,1))
110 ELSE ! no S restoring; this ensures frelflux=0
111 sssRelax(i,j)= sssFromOcn(i,j)
112 ENDIF
113
114
115 C Next lines: linearly phase out SST restoring between 2C and -1C
116 C ONLY if seaice is present
117 IF ((sstRelax(i,j).GT.2. _d 0).OR.
118 & (iceMask(i,j,1,1) .EQ. 0. _d 0)) THEN
119 nearIce=1.0
120 ELSEIF (sstRelax(i,j) .LE. -1. _d 0) THEN
121 nearIce=0.0
122 ELSE
123 nearIce=(sstRelax(i,j)+1.0)/3.0
124 endif
125
126 qrelflux= lambdaTheta*(sstFromOcn(i,j)-sstRelax(i,j))
127 & * (HeatCapacity_Cp*rhoNil*drF(1))*nearIce
128 C- should use rhoConst instead of rhoNil:
129 c & * (HeatCapacity_Cp*rhoConst*drF(1))*nearIce
130
131 C no longer restore on top of ice, but effectively full gridpoint UNDER ice
132 C (unless gridpoint is fully ice covered)
133 IF (iceMask(i,j,1,1) .LT. 0.999 _d 0) THEN
134 qneto_2D(i,j)= qneto_2D(i,j) + qrelflux
135 & / (1. _d 0 - iceMask(i,j,1,1))
136 ENDIF
137 C qneto_2D(i,j)= qneto_2D(i,j) + qrelflux
138 C qneti_2D(i,j)= qneti_2D(i,j) + qrelflux
139
140 frelflux= -lambdaSalt*(sssFromOcn(i,j)-sssRelax(i,j))/
141 & (convertFW2Salt *recip_drF(1))*nearIce
142
143 C or use actual salt instead of convertFW2salt above?
144
145 IF (frelflux .GT. 0. _d 0) THEN
146 evapo_2D(i,j)= evapo_2D(i,j) - frelflux
147 C note most of the time, evapi=0 when iceMask>0 anyway
148 C (i.e., only when relaxing SST >2 but ocn still ice-covered)
149 IF (iceMask(i,j,1,1).GT.0. _d 0)
150 & evapi_2D(i,j)= evapi_2D(i,j) - frelflux
151 ELSE
152 precipo_2D(i,j)= precipo_2D(i,j) + frelflux
153 IF (iceMask(i,j,1,1).GT.0. _d 0)
154 & precipi_2D(i,j)= precipi_2D(i,j) + frelflux
155 ENDIF
156
157 C IF (iceMask(i,j,1,1) .GT. 0. _d 0) THEN
158 C PRINT *,'Frelflux',frelflux,precipi_2D(i,j),atm_precip(j+1)
159 C ENDIF
160
161 C Diagnostics
162 sum_qrel(i,j)= sum_qrel(i,j) + qrelflux*dtatmo
163 sum_frel(i,j)= sum_frel(i,j) + frelflux*dtatmo
164
165 ENDIF
166 ENDDO
167 ENDDO
168 ENDIF
169
170 C PRINT *,'***bottom of relaxadd',wght0,wght1,intime0,intime1
171 C PRINT *,'evapo_2d: ',evapo_2D(JBUGI,JBUGJ)
172 C PRINT *,'precipo_2d: ',precipo_2D(JBUGI,JBUGJ)
173 C PRINT *,'qneto_2d: ',qneto_2D(JBUGI,JBUGJ)
174 C PRINT *,'SStfrom Ocn: ',sstfromocn(JBUGI,JBUGJ)
175 C PRINT *,'SSSfrom Ocn: ',sssfromocn(JBUGI,JBUGJ)
176
177 RETURN
178 END

  ViewVC Help
Powered by ViewVC 1.1.22