/[MITgcm]/MITgcm/pkg/seaice/seaice_tracer_phys.F
ViewVC logotype

Annotation of /MITgcm/pkg/seaice/seaice_tracer_phys.F

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


Revision 1.1 - (hide annotations) (download)
Tue Jun 7 03:58:23 2011 UTC (14 years, 1 month ago) by gforget
Branch: MAIN
- introducing ALLOW_SITRACER and seaice_tracer_phys.F to handle generic seaice tracer.
  For now it covers, and was tested for, salinity and age (work in progress).
- introducing siEps (1e-5, parameter, defined in SEAICE_PARAMS.h).

1 gforget 1.1 C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_tracer_growth.F,v 1.1 2011/06/02 16:21:08 gforget Exp $
2     C $Name: $
3    
4     #include "SEAICE_OPTIONS.h"
5    
6     CStartOfInterface
7     SUBROUTINE SEAICE_TRACER_PHYS( myTime, myIter, myThid )
8     C /=======================================================\
9     C | SUBROUTINE seaice_tracer_phys |
10     C | o Time step SItr/SItrEFF as a result of |
11     C | seaice thermodynamics and specific tracer physics |
12     C \=======================================================/
13     IMPLICIT NONE
14    
15     C === Global variables ===
16     #include "SIZE.h"
17     #include "EEPARAMS.h"
18     #include "FFIELDS.h"
19     #include "DYNVARS.h"
20     #include "SEAICE_SIZE.h"
21     #include "SEAICE.h"
22     #include "SEAICE_PARAMS.h"
23     #include "SEAICE_TRACER.h"
24    
25     C === Routine arguments ===
26     C INPUT:
27     C myTime :: Simulation time
28     C myIter :: Simulation timestep number
29     C myThid :: Thread no. that called this routine.
30     C OUTPUT:
31     _RL myTime
32     INTEGER myIter, myThid
33     CEndOfInterface
34    
35     C === Local variables ===
36     #ifdef ALLOW_SITRACER
37    
38     INTEGER iTr, jTh, I, J, bi, bj, ks
39     _RL SItrFromOcean (1:sNx,1:sNy)
40     _RL SItrFromSnow (1:sNx,1:sNy)
41     _RL HEFFprev, HEFFpost, growFact, meltPart, tmpscal1
42    
43     #ifdef ALLOW_SITRACER_DIAG
44     _RL DIAGarray (1:sNx,1:sNy,Nr)
45     #endif
46    
47     cgf for now I do not fully account for ocean-ice fluxes of tracer
48     cgf -> I just prescribe it consistent with age tracer
49     cgf eventually I will need to handle them as function params
50    
51     ks=1
52    
53     DO bj=myByLo(myThid),myByHi(myThid)
54     DO bi=myBxLo(myThid),myBxHi(myThid)
55     DO iTr=1,SItrMaxNum
56    
57     c 0) set ice-ocean and ice-snow exchange values
58     c =============================================
59     DO J=1,sNy
60     DO I=1,sNx
61     SItrFromOcean(i,j)=0. _d 0
62     SItrFromSnow(i,j)=0. _d 0
63     ENDDO
64     ENDDO
65     if (SItrName(iTr).EQ.'age') then
66     c age tracer: not passed from ocean or snow
67     DO J=1,sNy
68     DO I=1,sNx
69     SItrFromOcean(i,j)=0. _d 0
70     ENDDO
71     ENDDO
72     elseif (SItrName(iTr).EQ.'salinity') then
73     c salinity tracer:
74     DO J=1,sNy
75     DO I=1,sNx
76     SItrFromOcean(i,j)=SIsal0
77     #ifdef SEAICE_VARIABLE_SALINITY
78     if (SIsalFRAC.GT.0.)
79     & SItrFromOcean(i,j)=SIsalFRAC*salt(I,j,ks,bi,bj)
80     #endif
81     ENDDO
82     ENDDO
83     elseif (SItrName(iTr).EQ.'one') then
84     c "ice concentration" tracer that should remain .EQ.1.
85     DO J=1,sNy
86     DO I=1,sNx
87     SItrFromOcean(i,j)=1. _d 0
88     SItrFromSnow(i,j)=1. _d 0
89     ENDDO
90     ENDDO
91     endif
92     c 1) seaice thermodynamics processes
93     c ==================================
94     DO J=1,sNy
95     DO I=1,sNx
96     HEFFprev=SItrHEFF(i,j,bi,bj,1)
97     #ifdef ALLOW_SITRACER_DIAG
98     diagArray(I,J,5+(iTr-1)*5) =
99     & HEFFprev*SItracer(i,j,bi,bj,iTr) + SItrBucket(i,j,bi,bj,iTr)
100     #endif
101     c apply the sequence of thermodynamics increments to actual traceur
102     c (see seaice_growth.F)
103     c (jTh=1 tendency due to ice-ocean interaction)
104     c (jTh=2 tendency due to the atmosphere, over ice covered part)
105     c (jTh=3 tendency due to the atmosphere, over open water part)
106     c (jTh=4 tendency due to flooding)
107     DO jTh=1,3
108     HEFFprev=SItrHEFF(i,j,bi,bj,jTh)
109     HEFFpost=SItrHEFF(i,j,bi,bj,jTh+1)
110     c compute ratio in [0. 1.] range for either growth or melt
111     growFact=1. _d 0
112     meltPart=0. _d 0
113     if (HEFFpost.GT.HEFFprev) growFact=HEFFprev/HEFFpost
114     if (HEFFpost.LT.HEFFprev) meltPart=HEFFprev-HEFFpost
115     c update SItr accordingly
116     SItracer(i,j,bi,bj,iTr)=SItracer(i,j,bi,bj,iTr)*growFact
117     & +SItrFromOcean(i,j)*(1. _d 0 - growFact)
118     SItrBucket(i,j,bi,bj,iTr)=SItrBucket(i,j,bi,bj,iTr)
119     & -HEFFpost*SItrFromOcean(i,j)*(1. _d 0 - growFact)
120     SItrBucket(i,j,bi,bj,iTr)=SItrBucket(i,j,bi,bj,iTr)
121     & +meltPart*SItracer(i,j,bi,bj,iTr)
122     ENDDO
123     c apply flooding term
124     growFact=1. _d 0
125     HEFFprev=SItrHEFF(i,j,bi,bj,4)
126     HEFFpost=SItrHEFF(i,j,bi,bj,5)
127     if (HEFFpost.GT.HEFFprev) growFact=HEFFprev/HEFFpost
128     SItracer(i,j,bi,bj,iTr)=SItracer(i,j,bi,bj,iTr)*growFact
129     & +SItrFromSnow(i,j) *(1. _d 0 - growFact)
130     if ((SItrName(iTr).EQ.'age').OR.(SItrName(iTr).EQ.'one')) then
131     c Including this term in the bucket only makes sense for diagnostics purposes
132     c and should not be done for tracers that are exchanged with the ocean (we would
133     c need another bucket for ice-snow exchange, and snow tracers to start with).
134     SItrBucket(i,j,bi,bj,iTr)=SItrBucket(i,j,bi,bj,iTr)
135     & -HEFFpost*SItrFromSnow(i,j)*(1. _d 0 - growFact)
136     endif
137     #ifdef ALLOW_SITRACER_DIAG
138     diagArray(I,J,5+(iTr-1)*5) = HEFFpost*SItracer(i,j,bi,bj,iTr)
139     & +SItrBucket(i,j,bi,bj,iTr)-diagArray(I,J,5+(iTr-1)*5)
140     #endif
141     ENDDO
142     ENDDO
143     c 2) very ice tracer processes
144     c ============================
145     if (SItrName(iTr).EQ.'age') then
146     c age tracer: grow old as time passes by
147     DO J=1,sNy
148     DO I=1,sNx
149     if (SItrHEFF(i,j,bi,bj,5).GT.0. _d 0) then
150     SItracer(i,j,bi,bj,iTr)=
151     & SItracer(i,j,bi,bj,iTr)+SEAICE_deltaTtherm
152     else
153     SItracer(i,j,bi,bj,iTr)=0. _d 0
154     endif
155     ENDDO
156     ENDDO
157     elseif (SItrName(iTr).EQ.'salinity') then
158     c salinity tracer: no specific process
159     elseif (SItrName(iTr).EQ.'one') then
160     c "ice concentration" tracer: no specific process
161     endif
162     c 3) ice-ocean tracer exchange
163     c =============================
164     if (SItrName(iTr).EQ.'age') then
165     c age tracer: not passed to ocean
166     elseif (SItrName(iTr).EQ.'salinity') then
167     c salinity tracer: salt flux
168     c DO J=1,sNy
169     c DO I=1,sNx
170     c saltFlux(I,J,bi,bj) = - SItrBucket(i,j,bi,bj,iTr)
171     c & *HEFFM(I,J,bi,bj)/SEAICE_deltaTtherm*SEAICE_rhoIce
172     c note: at this point of the time step, that is the correct sign
173     c saltPlumeFlux(I,J,bi,bj) = ...
174     c ENDDO
175     c ENDDO
176     elseif (SItrName(iTr).EQ.'one') then
177     c "ice concentration" tracer: not passed to ocean
178     endif
179     DO J=1,sNy
180     DO I=1,sNx
181     #ifdef ALLOW_SITRACER_DIAG
182     diagArray(I,J,4+(iTr-1)*5) = - SItrBucket(i,j,bi,bj,iTr)
183     & *HEFFM(I,J,bi,bj)/SEAICE_deltaTtherm*SEAICE_rhoIce
184     #endif
185     c empty bucket
186     SItrBucket(i,j,bi,bj,iTr)=0. _d 0
187     ENDDO
188     ENDDO
189     c 4) diagnostics
190     c ==============
191     #ifdef ALLOW_SITRACER_DIAG
192     DO J=1,sNy
193     DO I=1,sNx
194     HEFFpost=SItrHEFF(i,j,bi,bj,5)
195     DIAGarray(I,J,1+(iTr-1)*5) = SItracer(i,j,bi,bj,iTr)
196     DIAGarray(I,J,2+(iTr-1)*5) = SItracer(i,j,bi,bj,iTr)*HEFFpost
197     c diagArray(:,:,3) is the term of comparison for diagArray(:,:,2)
198     if (SItrName(iTr).EQ.'age') then
199     DIAGarray(I,J,3+(iTr-1)*5) = IceAgeTr(i,j,bi,bj,2)
200     elseif (SItrName(iTr).EQ.'salinity') then
201     DIAGarray(I,J,3+(iTr-1)*5) = HSALT(i,j,bi,bj)/SEAICE_rhoIce
202     elseif (SItrName(iTr).EQ.'one') then
203     DIAGarray(I,J,3+(iTr-1)*5) = HEFFpost
204     endif
205     c diagArray(:,:,4) allows check of conservation : del(SItrBucket)+del(SItr*HEFF)=0. over do_phys
206     c diagArray(:,:,5) is the tracer flux from the ocean (<0 incr. ocean tracer)
207     ENDDO
208     ENDDO
209     #endif
210     ENDDO
211     #ifdef ALLOW_SITRACER_DIAG
212     CALL DIAGNOSTICS_FILL(DIAGarray,'UDIAG1 ',0,Nr,3,bi,bj,myThid)
213     #endif
214     ENDDO
215     ENDDO
216    
217     #endif /* ALLOW_SITRACER */
218    
219     RETURN
220     END

  ViewVC Help
Powered by ViewVC 1.1.22