/[MITgcm]/MITgcm_contrib/SOSE/code_ad/cost_driftw.F
ViewVC logotype

Annotation of /MITgcm_contrib/SOSE/code_ad/cost_driftw.F

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


Revision 1.1 - (hide annotations) (download)
Fri Apr 23 19:55:11 2010 UTC (15 years, 3 months ago) by mmazloff
Branch: MAIN
CVS Tags: HEAD
original files

1 mmazloff 1.1 C $Header: /u/gcmpack/MITgcm/pkg/ecco/cost_driftw.F,v 1.8 2009/04/28 18:13:28 jmc Exp $
2     C $Name: $
3    
4     #include "COST_CPPOPTIONS.h"
5    
6    
7     subroutine cost_driftw(
8     I myiter,
9     I mytime,
10     I mythid
11     & )
12    
13     c ==================================================================
14     c SUBROUTINE cost_Driftw
15     c ==================================================================
16     c
17     c o Evaluate cost function contribution of the w difference
18     c between the first and the last year.
19     c
20     c started: from cost_drift
21     c
22     c Armin Koehl akoehl@ucsd.edu 26-Feb-2002
23     c
24     c ==================================================================
25     c SUBROUTINE cost_Driftw
26     c ==================================================================
27    
28     implicit none
29    
30     c == global variables ==
31    
32     #include "EEPARAMS.h"
33     #include "SIZE.h"
34     #include "GRID.h"
35     #include "DYNVARS.h"
36    
37     #include "cal.h"
38     #include "ecco_cost.h"
39     #include "ctrl.h"
40     #include "ctrl_dummy.h"
41     #include "optim.h"
42     CMM(
43     #ifdef ALLOW_OBCS
44     # include "OBCS.h"
45     #endif
46     CMM)
47    
48     c == routine arguments ==
49    
50     integer myiter
51     _RL mytime
52     integer mythid
53    
54     c == local variables ==
55    
56     _RS one_rs
57     parameter( one_rs = 1. )
58    
59     integer bi,bj
60     integer i,j,k
61     integer itlo,ithi
62     integer jtlo,jthi
63     integer jmin,jmax
64     integer imin,imax
65     integer irec
66     integer ilw
67     integer nf, nl, nfmin
68     integer minrec
69    
70     _RL fctilew
71     _RL fcthread_wdrift
72    
73     character*(80) fnamew
74    
75     logical doglobalread
76     logical ladinit
77    
78     character*(MAX_LEN_MBUF) msgbuf
79    
80     c == external functions ==
81    
82     integer ilnblnk
83     external ilnblnk
84    
85     c == end of interface ==
86    
87     jtlo = mybylo(mythid)
88     jthi = mybyhi(mythid)
89     itlo = mybxlo(mythid)
90     ithi = mybxhi(mythid)
91     jmin = 1
92     jmax = sny
93     imin = 1
94     imax = snx
95    
96     c-- Read tiled data.
97     doglobalread = .false.
98     ladinit = .false.
99    
100     #ifdef ALLOW_DRIFTW_COST_CONTRIBUTION
101    
102     if (optimcycle .ge. 0) then
103     ilw = ilnblnk( wbarfile )
104     write(fnamew(1:80),'(2a,i10.10)')
105     & wbarfile(1:ilw),'.',optimcycle
106     endif
107    
108     fcthread_wdrift = 0. _d 0
109    
110     do bj = jtlo,jthi
111     do bi = itlo,ithi
112     do k = 1,nr
113     do j = jmin,jmax
114     do i = imin,imax
115     wfmean(i,j,k,bi,bj) = 0.0
116     wlmean(i,j,k,bi,bj) = 0.0
117     enddo
118     enddo
119     enddo
120     enddo
121     enddo
122    
123     nf = 0
124     nl = 0
125     c-- Number of full years
126     CMM nfmin = MAX(INT(FLOAT(nmonsrec)/12.),1)
127     c-- Prevent code from crashing if integrated for less than a year
128     CMM minrec = MIN(nmonsrec,12)
129    
130     CMM only do k = 25 = 1348.5m do k = 1,nr
131    
132     k = 25
133    
134     do irec = 1,nmonsrec
135     c-- Read time averages and the monthly mean data.
136     call active_read_xyz( fnamew, wbar, irec,
137     & doglobalread, ladinit,
138     & optimcycle, mythid,
139     & xx_wbar_mean_dummy )
140    
141     nl = nl + 1
142    
143     do bj = jtlo,jthi
144     do bi = itlo,ithi
145     CMM do k = 1,nr
146     do j = jmin,jmax
147     do i = imin,imax
148     wlmean(i,j,k,bi,bj) = wlmean(i,j,k,bi,bj) +
149     & wbar(i,j,k,bi,bj)/nmonsrec
150     enddo
151     enddo
152     CMM enddo
153     enddo
154     enddo
155    
156     enddo
157     CMM(
158     do bj = jtlo,jthi
159     do bi = itlo,ithi
160    
161     c-- Loop over the model layers
162     fctilew = 0. _d 0
163    
164     c-- Compute model misfit and cost function term for
165     c the vertical velovity field. The error is 1e-4 m/s.
166     CMM see if on open boundary north or south
167     CMM assumed using obcs west south north
168     do i = imin,imax
169     if (OB_Jn(I,bi,bj).GT.0) then
170     do j = OB_Jn(I,bi,bj)-10,OB_Jn(I,bi,bj)-1
171     if (_hFacC(i,j,32,bi,bj) .GT. 0. ) then
172     fctilew = fctilew + wlmean(i,j,k,bi,bj)
173     num_wdrift(bi,bj) = num_wdrift(bi,bj) + 1. _d 0
174     endif
175     enddo
176     endif
177     enddo
178    
179     CMM bi bj
180     CMM enddo
181     CMM enddo
182     CMM)
183    
184     if (mult_ageos.gt.0.0) then
185     print*,'CMM:cost w: fctilew, mult_ageo = ', fctilew, mult_ageos
186     print*,'CMM:cost w num terms = ', num_wdrift
187     fctilew = fctilew / mult_ageos / mult_ageos
188     endif
189    
190     c-- End of loop over layers.
191    
192     fcthread_wdrift = fcthread_wdrift + fctilew
193     objf_wdrift(bi,bj) = objf_wdrift(bi,bj) + fctilew
194    
195     #ifdef ECCO_VERBOSE
196     c-- Print cost function for each tile in each thread.
197     write(msgbuf,'(a)') ' '
198     call print_message( msgbuf, standardmessageunit,
199     & SQUEEZE_RIGHT , mythid)
200     write(msgbuf,'(a,i8.8,1x,i3.3,1x,i3.3)')
201     & ' cost_Driftw: irec,bi,bj = ',irec,bi,bj
202     call print_message( msgbuf, standardmessageunit,
203     & SQUEEZE_RIGHT , mythid)
204     write(msgbuf,'(a,d22.15)')
205     & ' cost function (wvel) = ',
206     & fcthread_wdrift
207     call print_message( msgbuf, standardmessageunit,
208     & SQUEEZE_RIGHT , mythid)
209     write(msgbuf,'(a)') ' '
210     call print_message( msgbuf, standardmessageunit,
211     & SQUEEZE_RIGHT , mythid)
212     #endif
213    
214     enddo
215     enddo
216    
217     #ifdef ECCO_VERBOSE
218     c-- Print cost function for all tiles.
219     _GLOBAL_SUM_RL( fcthread_wdrift , myThid )
220     write(msgbuf,'(a)') ' '
221     call print_message( msgbuf, standardmessageunit,
222     & SQUEEZE_RIGHT , mythid)
223     write(msgbuf,'(a,i8.8)')
224     & ' cost_Driftw: irec = ',irec
225     call print_message( msgbuf, standardmessageunit,
226     & SQUEEZE_RIGHT , mythid)
227     write(msgbuf,'(a,a,d22.15)')
228     & ' cost function value',
229     & ' (wvel) = ',fcthread_wdrift
230     call print_message( msgbuf, standardmessageunit,
231     & SQUEEZE_RIGHT , mythid)
232     write(msgbuf,'(a)') ' '
233     call print_message( msgbuf, standardmessageunit,
234     & SQUEEZE_RIGHT , mythid)
235     #endif
236    
237     #endif
238    
239     return
240     end
241    

  ViewVC Help
Powered by ViewVC 1.1.22