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

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

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


Revision 1.4 - (hide annotations) (download)
Tue Oct 9 00:10:13 2007 UTC (16 years, 7 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint59i, checkpoint59j
Changes since 1.3: +4 -3 lines
add missing cvs $Header:$ or $Name:$

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

  ViewVC Help
Powered by ViewVC 1.1.22