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

Contents of /MITgcm/pkg/seaice/seaice_model.F

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


Revision 1.11 - (show annotations) (download)
Fri Oct 24 05:29:36 2003 UTC (20 years, 7 months ago) by edhill
Branch: MAIN
CVS Tags: checkpoint51o_pre, checkpoint51t_post, checkpoint51s_post, checkpoint51q_post, checkpoint51r_post, checkpoint51o_post, checkpoint51p_post
Branch point for: branch-nonh
Changes since 1.10: +2 -1 lines
 o undid all of the cp51 checkin pending some ongoing code cleanups
   and discussion

1 C $Header: /u/u3/gcmpack/MITgcm/pkg/seaice/seaice_model.F,v 1.9 2003/10/09 04:19:20 edhill Exp $
2 C $Name: $
3
4 #include "SEAICE_OPTIONS.h"
5
6 CStartOfInterface
7 SUBROUTINE seaice_model( myTime, myIter, myThid )
8 C /===========================================================\
9 C | SUBROUTINE SEAICE_MODEL |
10 C | o Time stepping of a dynamic/thermodynamic sea ice model. |
11 C | Dynamics solver: Zhang/Hibler, JGR, 102, 8691-8702, 1997 |
12 C | and Zhang/Rothrock, JGR, 105, 3325-3338, 2000 |
13 C | Thermodynamics: Hibler, MWR, 108, 1943-1973, 1980 |
14 C | Rheology: Hibler, JPO, 9, 815- 846, 1979 |
15 C | Snow: Zhang et al. , JPO, 28, 191- 217, 1998 |
16 C | Parallel forward ice model written by Jinlun Zhang PSC/UW|
17 C | & coupled into MITgcm by Dimitris Menemenlis (JPL) 2/2001|
18 C | zhang@apl.washington.edu / menemenlis@jpl.nasa.gov |
19 C |===========================================================|
20 C \===========================================================/
21 IMPLICIT NONE
22
23 C === Global variables ===
24 #include "SIZE.h"
25 #include "EEPARAMS.h"
26 #include "DYNVARS.h"
27 #include "PARAMS.h"
28 #include "FFIELDS.h"
29 #include "SEAICE.h"
30 #include "SEAICE_PARAMS.h"
31 #include "SEAICE_EXTERNAL.h"
32
33 #ifdef SEAICE_EXTERNAL_FORCING
34 # include "SEAICE_FFIELDS.h"
35 #endif
36
37 #ifdef ALLOW_AUTODIFF_TAMC
38 # include "tamc.h"
39 #endif
40
41 C === Routine arguments ===
42 C myTime - Simulation time
43 C myIter - Simulation timestep number
44 C myThid - Thread no. that called this routine.
45 _RL myTime
46 INTEGER myIter
47 INTEGER myThid
48 CEndOfInterface
49
50 #ifdef ALLOW_SEAICE
51
52 C === Local variables ===
53 C i,j,bi,bj - Loop counters
54
55 INTEGER i, j, bi, bj
56
57 #ifdef SEAICE_EXTERNAL_FORCING
58 C-- Atmospheric state and runoff are from
59 C pkg/exf, which does not update edges.
60 _EXCH_XY_R8( uwind, myThid )
61 _EXCH_XY_R8( vwind, myThid )
62 _EXCH_XY_R8( atemp, myThid )
63 _EXCH_XY_R8( aqh, myThid )
64 _EXCH_XY_R8( lwdown, myThid )
65 _EXCH_XY_R8( swdown, mythid )
66 _EXCH_XY_R8( precip, myThid )
67 _EXCH_XY_R8( evap, myThid )
68 _EXCH_XY_R8( runoff, myThid )
69 #else /* SEAICE_EXTERNAL_FORCING */
70 C-- Load atmospheric state and runoff.
71 CALL SEAICE_GET_FORCING ( myTime, myIter, myThid )
72 #endif /* SEAICE_EXTERNAL_FORCING */
73
74 C-- Third level model velocity is used as proxy for geostrophic velocity
75 DO bj=myByLo(myThid),myByHi(myThid)
76 DO bi=myBxLo(myThid),myBxHi(myThid)
77 DO j=0,sNy+1
78 DO i=0,sNx+1
79 GWATX(I,J,bi,bj)=HALF*(uVel(i,j,3,bi,bj)
80 & +uVel(i,j-1,3,bi,bj))
81 GWATY(I,J,bi,bj)=HALF*(vVel(i,j,3,bi,bj)
82 & +vVel(i-1,j,3,bi,bj))
83 #ifdef SEAICE_DEBUG
84 c write(*,'(2i4,2i2,f7.1,7f12.3)')
85 c & ,i,j,bi,bj,UVM(I,J,bi,bj)
86 c & ,GWATX(I,J,bi,bj),GWATY(I,J,bi,bj)
87 c & ,uVel(i+1,j,3,bi,bj),uVel(i+1,j+1,3,bi,bj)
88 c & ,vVel(i,j+1,3,bi,bj),vVel(i+1,j+1,3,bi,bj)
89 #endif
90 ENDDO
91 ENDDO
92 ENDDO
93 ENDDO
94
95 #ifdef ALLOW_AUTODIFF_TAMC
96 CADJ STORE uwind = comlev1, key = ikey_dynamics
97 CADJ STORE vwind = comlev1, key = ikey_dynamics
98 # ifdef SEAICE_ALLOW_DYNAMICS
99 CADJ STORE heff = comlev1, key = ikey_dynamics
100 CADJ STORE area = comlev1, key = ikey_dynamics
101 # endif
102 #endif /* ALLOW_AUTODIFF_TAMC */
103
104 C solve ice momentum equations and calculate ocean surface stress
105 CALL DYNSOLVER ( myTime, myIter, myThid )
106
107 #ifdef ALLOW_AUTODIFF_TAMC
108 # ifdef SEAICE_ALLOW_DYNAMICS
109 CADJ STORE heff = comlev1, key = ikey_dynamics
110 CADJ STORE area = comlev1, key = ikey_dynamics
111 CADJ STORE uice = comlev1, key = ikey_dynamics
112 CADJ STORE vice = comlev1, key = ikey_dynamics
113 # endif
114 #endif /* ALLOW_AUTODIFF_TAMC */
115 C NOW DO ADVECTION
116 CALL ADVECT( UICE, VICE, HEFF, HEFFM, myThid )
117 CALL ADVECT( UICE, VICE, AREA, HEFFM, myThid )
118
119 C NOW DO GROWTH
120 C MUST CALL GROWTH ONLY AFTER CALLING ADVECTION
121 CALL GROWTH( myTime, myIter, myThid)
122
123 C-- Update overlap regions for a bunch of stuff
124 _BARRIER
125 CALL SEAICE_EXCH( HEFF, myThid )
126 CALL SEAICE_EXCH( AREA, myThid )
127 _EXCH_XY_R4(fu , myThid )
128 _EXCH_XY_R4(fv , myThid )
129 _EXCH_XY_R4(EmPmR, myThid )
130 _EXCH_XY_R4(Qnet , myThid )
131 _EXCH_XY_R4(surfaceTendencyTice, myThid )
132 #ifdef SHORTWAVE_HEATING
133 _EXCH_XY_R4(Qsw , myThid )
134 #endif
135 _EXCH_XYZ_R8(theta , myThid )
136
137 C-- Sea ice diagnostics.
138 CALL SEAICE_DO_DIAGS( myTime, myIter, myThid )
139
140 C-- Write sea ice restart files
141 CALL SEAICE_WRITE_PICKUP ( .FALSE.,
142 & myTime+deltaTClock, myIter+1, myThid )
143
144 C---------------------------------------------------
145 C OOH NOOOO we need to move the whole stuff
146 C---------------------------------------------------
147 #ifdef ALLOW_AUTODIFF_TAMC
148 CRG CADJ store UICE,VICE,AREA,HEFF,fu,fv,EmPmR,Qnet,Qsw = comlev1_bibj
149 #endif
150
151 C-- Call sea-ice cost function routine
152 CRG CALL SEAICE_COST( myTime, myIter, myThid )
153
154 #endif /* ALLOW_SEAICE */
155
156 RETURN
157 END

  ViewVC Help
Powered by ViewVC 1.1.22