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

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

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

Revision 1.14 - (hide annotations) (download)
Wed Jun 24 08:23:38 2009 UTC (16 years, 4 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61x
Changes since 1.13: +29 -31 lines 
third step of cleaning up the 3-time levels of UICE,VICE,HEFF,AREA:
extra argument ksize no longer required
1 mlosch 1.14 C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_calc_strainrates.F,v 1.13 2009/06/03 08:48:07 mlosch Exp $
2 mlosch 1.1 C $Name: $
3    
4     #include "SEAICE_OPTIONS.h"
5    
6     CStartOfInterface
7 jmc 1.8 SUBROUTINE SEAICE_CALC_STRAINRATES(
8 mlosch 1.1 I uFld, vFld,
9 mlosch 1.12 O e11Loc, e22Loc, e12Loc,
10 mlosch 1.14 I iStep, myTime, myIter, myThid )
11 mlosch 1.1 C /==========================================================\
12     C | SUBROUTINE SEAICE_CALC_STRAINRATES |
13     C | o compute strain rates from ice velocities |
14     C |==========================================================|
15     C | written by Martin Losch, Apr 2007 |
16     C \==========================================================/
17     IMPLICIT NONE
18    
19     C === Global variables ===
20     #include "SIZE.h"
21     #include "EEPARAMS.h"
22     #include "PARAMS.h"
23     #include "GRID.h"
24     #include "SEAICE_PARAMS.h"
25 mlosch 1.11 #include "SEAICE.h"
26 mlosch 1.1
27     #ifdef ALLOW_AUTODIFF_TAMC
28     # include "tamc.h"
29     #endif
30    
31     C === Routine arguments ===
32 jmc 1.8 C iStep :: Sub-time-step number
33     C myTime :: Simulation time
34     C myIter :: Simulation timestep number
35     C myThid :: My Thread Id. number
36     INTEGER iStep
37     _RL myTime
38     INTEGER myIter
39 mlosch 1.1 INTEGER myThid
40     C ice velocities
41 mlosch 1.14 _RL uFld (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
42     _RL vFld (1-Olx:sNx+Olx,1-Oly:sNy+Oly,nSx,nSy)
43 mlosch 1.1 C strain rate tensor
44 mlosch 1.12 _RL e11Loc (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
45     _RL e22Loc (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
46     _RL e12Loc (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
47 mlosch 1.1 CEndOfInterface
48    
49     #ifdef SEAICE_CGRID
50     #ifdef SEAICE_ALLOW_DYNAMICS
51     C === Local variables ===
52     C i,j,bi,bj - Loop counters
53     INTEGER i, j, bi, bj
54 mlosch 1.2 C hFacU, hFacV - determine the no-slip boundary condition
55     INTEGER k
56 mlosch 1.11 _RS hFacU, hFacV, noSlipFac
57 mlosch 1.2
58 mlosch 1.4 k = 1
59 mlosch 1.11 noSlipFac = 0. _d 0
60     IF ( SEAICE_no_slip ) noSlipFac = 1. _d 0
61 mlosch 1.1 C
62 mlosch 1.11 #ifndef SEAICE_OLD_AND_BAD_DISCRETIZATION
63     DO bj=myByLo(myThid),myByHi(myThid)
64     DO bi=myBxLo(myThid),myBxHi(myThid)
65     DO j=1-Oly,sNy+Oly-1
66     DO i=1-Olx,sNx+Olx-1
67     C evaluate strain rates
68 mlosch 1.12 e11Loc(I,J,bi,bj) = _recip_dxF(I,J,bi,bj) *
69 mlosch 1.14 & (uFld(I+1,J,bi,bj)-uFld(I,J,bi,bj))
70 mlosch 1.11 & +HALF*
71 mlosch 1.14 & (vFld(I,J,bi,bj)+vFld(I,J+1,bi,bj))
72 mlosch 1.11 & * k2AtC(I,J,bi,bj)
73 mlosch 1.12 e22Loc(I,J,bi,bj) = _recip_dyF(I,J,bi,bj) *
74 mlosch 1.14 & (vFld(I,J+1,bi,bj)-vFld(I,J,bi,bj))
75 mlosch 1.11 & +HALF*
76 mlosch 1.14 & (uFld(I,J,bi,bj)+uFld(I+1,J,bi,bj))
77 mlosch 1.11 & * k1AtC(I,J,bi,bj)
78     ENDDO
79     ENDDO
80     DO j=1-Oly+1,sNy+Oly
81     DO i=1-Olx+1,sNx+Olx
82     hFacU = _maskW(i,j,k,bi,bj) - _maskW(i,j-1,k,bi,bj)
83     hFacV = _maskS(i,j,k,bi,bj) - _maskS(i-1,j,k,bi,bj)
84 mlosch 1.12 e12Loc(I,J,bi,bj) = HALF*(
85 mlosch 1.14 & ( uFld(I,J,bi,bj) - uFld(I ,J-1,bi,bj) )
86 mlosch 1.11 & * _recip_dyU(I,J,bi,bj)
87 mlosch 1.14 & + ( vFld(I,J,bi,bj) - vFld(I-1,J ,bi,bj) )
88 mlosch 1.11 & * _recip_dxV(I,J,bi,bj)
89     & - k1AtZ(I,J,bi,bj)
90 mlosch 1.14 & * 0.5 _d 0 * (vFld(I,J,bi,bj)+vFld(I-1,J ,bi,bj))
91 mlosch 1.11 & - k2AtZ(I,J,bi,bj)
92 mlosch 1.14 & * 0.5 _d 0 * (uFld(I,J,bi,bj)+uFld(I ,J-1,bi,bj))
93 mlosch 1.11 & )
94     & *maskC(I ,J ,k,bi,bj)*maskC(I-1,J ,k,bi,bj)
95     & *maskC(I ,J-1,k,bi,bj)*maskC(I-1,J-1,k,bi,bj)
96     & + 2.0 _d 0 * noSlipFac * (
97 mlosch 1.14 & ( uFld(I,J,bi,bj) + uFld(I ,J-1,bi,bj) )
98 mlosch 1.11 & * _recip_dyU(I,J,bi,bj) * hFacU
99 mlosch 1.14 & + ( vFld(I,J,bi,bj) + vFld(I-1,J ,bi,bj) )
100 mlosch 1.11 & * _recip_dxV(I,J,bi,bj) * hFacV
101     & )
102     C no slip at the boundary implies u(j)+u(j-1)=0 and v(i)+v(i-1)=0
103     C accross the boundary; this is already accomplished by masking so
104     C that the following lines are not necessary
105     c$$$ & - hFacV * k1AtZ(I,J,bi,bj)
106 mlosch 1.14 c$$$ & * 0.5 _d 0 * (vFld(I,J,bi,bj)+vFld(I-1,J ,bi,bj))
107 mlosch 1.11 c$$$ & - hFacU * k2AtZ(I,J,bi,bj)
108 mlosch 1.14 c$$$ & * 0.5 _d 0 * (uFld(I,J,bi,bj)+uFld(I ,J-1,bi,bj))
109 mlosch 1.11 ENDDO
110     ENDDO
111    
112     c$$$ ENDIF
113     ENDDO
114     ENDDO
115     #else
116     C this the old and incomplete discretization, here I also erroneously
117     C used finite-volumes to discretize the strain rates
118 mlosch 1.1 DO bj=myByLo(myThid),myByHi(myThid)
119     DO bi=myBxLo(myThid),myBxHi(myThid)
120 mlosch 1.5 DO j=1-Oly,sNy+Oly-1
121     DO i=1-Olx,sNx+Olx-1
122     C evaluate strain rates
123 mlosch 1.12 e11Loc(I,J,bi,bj) = _recip_dxF(I,J,bi,bj) *
124 mlosch 1.14 & (uFld(I+1,J,bi,bj)-uFld(I,J,bi,bj))
125 mlosch 1.1 & -HALF*
126 mlosch 1.14 & (vFld(I,J,bi,bj)+vFld(I,J+1,bi,bj))
127 mlosch 1.1 & * _tanPhiAtU(I,J,bi,bj)*recip_rSphere
128 mlosch 1.12 e22Loc(I,J,bi,bj) = _recip_dyF(I,J,bi,bj) *
129 mlosch 1.14 & (vFld(I,J+1,bi,bj)-vFld(I,J,bi,bj))
130 mlosch 1.1 C one metric term is missing
131 mlosch 1.5 ENDDO
132     ENDDO
133     DO j=1-Oly+1,sNy+Oly
134     DO i=1-Olx+1,sNx+Olx
135 mlosch 1.12 e12Loc(I,J,bi,bj) = HALF*(
136 mlosch 1.14 & (uFld(I ,J ,bi,bj) * _dxC(I ,J ,bi,bj)
137     & -uFld(I ,J-1,bi,bj) * _dxC(I ,J-1,bi,bj)
138     & +vFld(I ,J ,bi,bj) * _dyC(I ,J ,bi,bj)
139     & -vFld(I-1,J ,bi,bj) * _dyC(I-1,J ,bi,bj))
140 mlosch 1.1 & * recip_rAz(I,J,bi,bj)
141     & +
142 mlosch 1.14 & 0.25 _d 0 * (uFld(I,J,bi,bj)+uFld(I ,J-1,bi,bj))
143 mlosch 1.1 & * ( _tanPhiAtU(I,J,bi,bj) + _tanPhiAtU(I,J-1,bi,bj) )
144     & *recip_rSphere
145     & )
146 mlosch 1.4 & *maskC(I ,J ,k,bi,bj)*maskC(I-1,J ,k,bi,bj)
147     & *maskC(I ,J-1,k,bi,bj)*maskC(I-1,J-1,k,bi,bj)
148 mlosch 1.1 C one metric term is missing
149     ENDDO
150     ENDDO
151 mlosch 1.2 IF ( SEAICE_no_slip ) THEN
152 mlosch 1.12 C no slip boundary conditions apply only to e12Loc
153 mlosch 1.5 DO j=1-Oly+1,sNy+Oly
154     DO i=1-Olx+1,sNx+Olx
155 mlosch 1.2 hFacU = _maskW(i,j,k,bi,bj) - _maskW(i,j-1,k,bi,bj)
156     hFacV = _maskS(i,j,k,bi,bj) - _maskS(i-1,j,k,bi,bj)
157    
158 mlosch 1.12 e12Loc(I,J,bi,bj) = e12Loc(I,J,bi,bj)
159 mlosch 1.10 & + recip_rAz(i,j,bi,bj) * 2. _d 0 *
160 mlosch 1.14 & ( hFacU * ( _dxC(i,j-1,bi,bj)*uFld(i,j ,bi,bj)
161     & + _dxC(i,j, bi,bj)*uFld(i,j-1,bi,bj) )
162     & + hFacV * ( _dyC(i-1,j,bi,bj)*vFld(i ,j,bi,bj)
163     & + _dyC(i, j,bi,bj)*vFld(i-1,j,bi,bj) ) )
164 jmc 1.8 & - hFacU
165 mlosch 1.14 & * 0.25 _d 0 * (uFld(I,J,bi,bj)+uFld(I ,J-1,bi,bj))
166 mlosch 1.2 & * ( _tanPhiAtU(I,J,bi,bj) + _tanPhiAtU(I,J-1,bi,bj) )
167     & *recip_rSphere
168     C one metric term is missing
169     ENDDO
170     ENDDO
171    
172     ENDIF
173 mlosch 1.1 ENDDO
174     ENDDO
175 mlosch 1.11 #endif /* SEAICE_OLD_AND_BAD_DISCRETIZATION */
176 mlosch 1.1 #endif /* SEAICE_ALLOW_DYNAMICS */
177     #endif /* SEAICE_CGRID */
178     RETURN
179     END
  ViewVC Help
Powered by ViewVC 1.1.22