/[MITgcm]/MITgcm_contrib/ifenty/seaiceAdjointCode/cost_ice_test.F
ViewVC logotype

Annotation of /MITgcm_contrib/ifenty/seaiceAdjointCode/cost_ice_test.F

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


Revision 1.1 - (hide annotations) (download)
Fri Jun 29 18:54:03 2007 UTC (18 years ago) by gforget
Branch: MAIN
CVS Tags: HEAD
Ian Fenty sea-ice growth routines (adjointable, in developement)

1 gforget 1.1 C $Header: /u/gcmpack/MITgcm/pkg/seaice/cost_ice_test.F,v 1.3 2006/12/15 18:02:17 heimbach Exp $
2    
3     #include "SEAICE_OPTIONS.h"
4    
5     subroutine cost_ice_test( mytime, myiter, mythid )
6    
7     c ==================================================================
8     c SUBROUTINE cost_ice_test
9     c ==================================================================
10     c
11     c o Compute sea-ice cost function. The following options can
12     c be selected with data.cost variable cost_ice_flag:
13     c
14     c cost_ice_flag = 1
15     c - compute mean sea-ice volume
16     c costIceStart < mytime < costIceEnd
17     c
18     c cost_ice_flag = 2
19     c - compute mean sea-ice area
20     c costIceStart < mytime < costIceEnd
21     c
22     c cost_ice_flag = 3
23     c - heat content of top level plus latent heat of sea-ice
24     c costIceStart < mytime < costIceEnd
25     c
26     c cost_ice_flag = 4
27     c - heat content of top level
28     c costIceStart < mytime < costIceEnd
29     c
30     c cost_ice_flag = 5
31     c - heat content of top level plus sea-ice plus latent heat of snow
32     c costIceStart < mytime < costIceEnd
33     c
34     c cost_ice_flag = 6
35     c - quadratic cost function measuring difference between pkg/seaice
36     c AREA variable and simulated sea-ice measurements at every time
37     c step.
38     c
39     c ==================================================================
40     c
41     c started: menemenlis@jpl.nasa.gov 26-Feb-2003
42     c
43     c ==================================================================
44     c SUBROUTINE cost_ice_test
45     c ==================================================================
46    
47     implicit none
48    
49     c == global variables ==
50     #ifdef ALLOW_COST_ICE
51     #include "EEPARAMS.h"
52     #include "SIZE.h"
53     #include "GRID.h"
54     #include "PARAMS.h"
55     #include "SEAICE_COST.h"
56     #include "SEAICE.h"
57     #include "DYNVARS.h"
58     #include "cost.h"
59     #endif /* ALLOW_COST_ICE */
60    
61     c == routine arguments ==
62    
63     _RL mytime
64     integer myiter
65     integer mythid
66    
67     #if (defined (ALLOW_SEAICE) && defined (ALLOW_COST_ICE))
68    
69     c == local variables ==
70    
71     c msgBuf - Informational/error message buffer
72     CHARACTER*(MAX_LEN_MBUF) msgBuf
73     integer bi,bj,i,j
74     _RL tempVar
75    
76     c == external functions ==
77    
78     integer ilnblnk
79     external ilnblnk
80    
81     c == end of interface ==
82     print *,'if-ice A1', myTime,endTime,lastinterval
83     print *,'if-ice A2', endTime-lastinterval,deltaTClock
84     if ( myTime .GT. (endTime - lastinterval) ) then
85     tempVar = 1. /
86     & ( ( 1. + min(endTime-startTime,lastinterval) )
87     & / deltaTClock )
88    
89     cph(
90     print *, 'ph-ice B ', myiter, theta(4,4,1,1,1),
91     & area(4,4,1,1,1), heff(4,4,1,1,1)
92     print *,'if-ice B', myTime,endTime-lastinterval,endTime,tempVar
93     cph)
94     if ( cost_ice_flag .eq. 1 ) then
95    
96     c sea-ice volume
97     do bj=myByLo(myThid),myByHi(myThid)
98     do bi=myBxLo(myThid),myBxHi(myThid)
99     do j = 1,sny
100     do i = 1,snx
101     objf_ice(bi,bj) = objf_ice(bi,bj) +
102     & tempVar * rA(i,j,bi,bj) * HEFF(i,j,1,bi,bj)
103     enddo
104     enddo
105     enddo
106     enddo
107    
108     elseif ( cost_ice_flag .eq. 2 ) then
109     c sea-ice area
110     do bj=myByLo(myThid),myByHi(myThid)
111     do bi=myBxLo(myThid),myBxHi(myThid)
112     do j = 1,sny
113     do i = 1,snx
114     objf_ice(bi,bj) = objf_ice(bi,bj) +
115     & tempVar * rA(i,j,bi,bj) * AREA(i,j,1,bi,bj)
116     enddo
117     enddo
118     enddo
119     enddo
120    
121     c heat content of top level:
122     c theta * delZ * (sea water heat capacity = 3996 J/kg/K)
123     c * (density of sea-water = 1026 kg/m^3)
124     c
125     c heat content of sea-ice:
126     c tice * heff * (sea ice heat capacity = 2090 J/kg/K)
127     c * (density of sea-ice = 910 kg/m^3)
128     c
129     c note: to remove mass contribution to heat content,
130     c which is not properly accounted for by volume converving
131     c ocean model, theta and tice are referenced to freezing
132     c temperature of sea-ice, here -1.96 deg C
133     c
134     c latent heat content of sea-ice:
135     c - heff * (latent heat of fusion = 334000 J/kg)
136     c * (density of sea-ice = 910 kg/m^3)
137     c
138     c latent heat content of snow:
139     c - hsnow * (latent heat of fusion = 334000 J/kg)
140     c * (density of snow = 330 kg/m^3)
141    
142     elseif ( cost_ice_flag .eq. 3 ) then
143     c heat content of top level plus latent heat of sea-ice
144     do bj=myByLo(myThid),myByHi(myThid)
145     do bi=myBxLo(myThid),myBxHi(myThid)
146     do j = 1,sny
147     do i = 1,snx
148     objf_ice(bi,bj) = objf_ice(bi,bj) +
149     & tempVar * rA(i,j,bi,bj) * (
150     & (THETA(i,j,1,bi,bj) + 1.96 ) *
151     & drF(1) * 3996 * 1026 -
152     & HEFF(i,j,1,bi,bj) * 334000 * 910 )
153     enddo
154     enddo
155     enddo
156     enddo
157    
158     elseif ( cost_ice_flag .eq. 4 ) then
159     c heat content of top level
160     do bj=myByLo(myThid),myByHi(myThid)
161     do bi=myBxLo(myThid),myBxHi(myThid)
162     do j = 1,sny
163     do i = 1,snx
164     objf_ice(bi,bj) = objf_ice(bi,bj) +
165     & tempVar * rA(i,j,bi,bj) * (
166     & (THETA(i,j,1,bi,bj) + 1.96 ) *
167     & drF(1) * 3996 * 1026 )
168     enddo
169     enddo
170     enddo
171     enddo
172    
173     elseif ( cost_ice_flag .eq. 5 ) then
174     c heat content of top level plus sea-ice plus latent heat of snow
175     do bj=myByLo(myThid),myByHi(myThid)
176     do bi=myBxLo(myThid),myBxHi(myThid)
177     do j = 1,sny
178     do i = 1,snx
179     objf_ice(bi,bj) = objf_ice(bi,bj) +
180     & tempVar * rA(i,j,bi,bj) * (
181     & (THETA(i,j,1,bi,bj) + 1.96 ) *
182     & drF(1) * 3996 * 1026 +
183     & (TICE(i,j,bi,bj) - 273.15 + 1.96 ) *
184     & HEFF(I,J,1,bi,bj) * 2090 * 910 -
185     & HEFF(i,j,1,bi,bj) * 334000 * 910 -
186     & HSNOW(I,J,bi,bj) * 334000 * 330 )
187     enddo
188     enddo
189     enddo
190     enddo
191    
192     elseif ( cost_ice_flag .eq. 6 ) then
193     c Qadratic cost function measuring difference between pkg/seaice
194     c AREA variable and simulated sea-ice measurements at every time
195     c step. For time being no measurements are read-in. It is
196     c assumed that measurements are AREA=0.5 at all times everywhere.
197     do bj=myByLo(myThid),myByHi(myThid)
198     do bi=myBxLo(myThid),myBxHi(myThid)
199     do j = 1,sny
200     do i = 1,snx
201     objf_ice(bi,bj) = objf_ice(bi,bj) +
202     & ( AREA(i,j,1,bi,bj) - 0.5 ) *
203     & ( AREA(i,j,1,bi,bj) - 0.5 )
204     enddo
205     enddo
206     enddo
207     enddo
208    
209     else
210     WRITE(msgBuf,'(A)')
211     & 'COST_ICE: invalid cost_ice_flag'
212     CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
213     & SQUEEZE_RIGHT , myThid )
214     STOP 'ABNORMAL END: S/R COST_ICE'
215     endif
216     endif
217    
218     cph(
219     print *, 'ph-ice C ', myiter, objf_ice(1,1)
220     cph)
221    
222     #endif /* ALLOW_COST_ICE */
223    
224     end

  ViewVC Help
Powered by ViewVC 1.1.22