/[MITgcm]/MITgcm/optim/optim_readdata.F
ViewVC logotype

Annotation of /MITgcm/optim/optim_readdata.F

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


Revision 1.1.2.3 - (hide annotations) (download)
Fri Dec 6 01:47:35 2002 UTC (21 years, 3 months ago) by heimbach
Branch: ecco-branch
CVS Tags: icebear5, icebear4, icebear3, icebear2, ecco_c44_e26, ecco_c44_e27
Branch point for: icebear
Changes since 1.1.2.2: +25 -27 lines
o lsopt:
   changed BLAS calls from single prec. (SDOT, SNRM2,SSCAL)
   to double prec. (DDOT, DNRM2, DSCAL)
   for compatibility with IBM SP3/SP4
 o optim:
   bringing optim_readdata/optim_writedata formats up-to-date
   with latest ctrl_pack/ctrl_unpack formats.

1 heimbach 1.1.2.1
2     subroutine optim_readdata(
3     I nn,
4     I dfile,
5     I lheaderonly,
6     O ff,
7     O vv
8     & )
9    
10     c ==================================================================
11     c SUBROUTINE optim_readdata
12     c ==================================================================
13     c
14     c o Read the data written by the MITgcmUV state estimation setup and
15     c join them to one vector that is subsequently used by the minimi-
16     c zation algorithm "lsopt". Depending on the specified file name
17     c either the control vector or the gradient vector can be read.
18     c
19     c *dfile* should be the radix of the file: ecco_ctrl or ecco_cost
20     c
21     c started: Christian Eckert eckert@mit.edu 12-Apr-2000
22     c
23     c changed: Patrick Heimbach heimbach@mit.edu 19-Jun-2000
24     c - finished, revised and debugged
25     c
26     c ==================================================================
27     c SUBROUTINE optim_readdata
28     c ==================================================================
29    
30     implicit none
31    
32     c == global variables ==
33    
34     #include "EEPARAMS.h"
35     #include "SIZE.h"
36 heimbach 1.1.2.2 cgg Include ECCO_CPPOPTIONS because the ecco_ctrl,cost files
37     cgg have headers with options for OBCS masks.
38     #include "ECCO_CPPOPTIONS.h"
39 heimbach 1.1.2.1
40     #include "ecco.h"
41     #include "ctrl.h"
42     #include "optim.h"
43     #include "minimization.h"
44    
45     c == routine arguments ==
46    
47     integer nn
48     _RL ff
49    
50     #if defined (DYNAMIC)
51     _RL vv(nn)
52     #elif defined (USE_POINTER) || (MAX_INDEPEND == 0)
53     _RL vv
54     pointer (pvv,vv(1))
55     #else
56     integer nmax
57     parameter( nmax = MAX_INDEPEND )
58     _RL vv(nmax)
59     #endif
60    
61     character*(9) dfile
62     logical lheaderonly
63    
64     c == local variables ==
65    
66     integer bi,bj
67     integer biG,bjG
68     integer i,j
69     integer ii,k
70     integer icvar
71     integer icvrec
72     integer icvcomp
73     integer icvoffset
74     integer nopt
75     integer funit
76    
77     integer cbuffindex
78     _RL cbuff( sNx*nSx*nPx*sNy*nSy*nPy )
79    
80     character*(128) fname
81    
82     integer filei
83     integer filej
84     integer filek
85     integer filenopt
86     integer fileig
87     integer filejg
88     integer filensx
89     integer filensy
90     _RL fileff
91    
92 heimbach 1.1.2.2 cgg(
93     _RL gg
94     integer igg
95     integer iobcs
96     cgg)
97    
98 heimbach 1.1.2.1 c == end of interface ==
99    
100     c-- The reference i/o unit.
101     funit = 20
102    
103     c-- Next optimization cycle.
104     nopt = optimcycle
105    
106     if ( dfile .eq. ctrlname ) then
107     print*
108     print*,' OPTIM_READDATA: Reading control vector'
109     print*,' for optimization cycle: ',nopt
110     print*
111     else if ( dfile .eq. costname ) then
112     print*
113     print*,' OPTIM_READDATA: Reading cost function and'
114     print*,' gradient of cost function'
115     print*,' for optimization cycle: ',nopt
116     print*
117     else
118     print*
119     print*,' OPTIM_READDATA: subroutine called by a false *dfile*'
120     print*,' argument. *dfile* = ',dfile
121     print*
122     stop ' ... stopped in OPTIM_READDATA.'
123     endif
124    
125     c-- Read the data.
126    
127     bjG = 1 + (myygloballo - 1)/sny
128     biG = 1 + (myxgloballo - 1)/snx
129    
130     c-- Generate file name and open the file.
131     write(fname(1:128),'(4a,i4.4)')
132     & dfile,'_',expId(1:10),'.opt', nopt
133     open( funit, file = fname,
134     & status = 'old',
135     & form = 'unformatted',
136     & access = 'sequential' )
137     print*, 'opened file ', fname
138    
139     c-- Read the header.
140     read( funit ) nvartype
141     read( funit ) nvarlength
142     read( funit ) expId
143     read( funit ) filenopt
144     read( funit ) fileff
145     read( funit ) fileiG
146     read( funit ) filejG
147     read( funit ) filensx
148     read( funit ) filensy
149    
150     cph(
151     print *,'ph-opt 1 ', nvartype, nvarlength, filensx, filensy
152     cph)
153    
154 heimbach 1.1.2.3 read( funit ) (nWetcGlobal(k), k=1,nr)
155     read( funit ) (nWetsGlobal(k), k=1,nr)
156     read( funit ) (nWetwGlobal(k), k=1,nr)
157     read( funit ) (nWetvGlobal(k), k=1,nr)
158 heimbach 1.1.2.2
159     cgg( Add OBCS Mask information into the header section for optimization.
160     #ifdef ALLOW_OBCSN_CONTROL
161 heimbach 1.1.2.3 read( funit ) ((nWetobcsnGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
162 heimbach 1.1.2.2 #endif
163     #ifdef ALLOW_OBCSS_CONTROL
164 heimbach 1.1.2.3 read( funit ) ((nWetobcssGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
165 heimbach 1.1.2.2 #endif
166     #ifdef ALLOW_OBCSW_CONTROL
167 heimbach 1.1.2.3 read( funit ) ((nWetobcswGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
168 heimbach 1.1.2.2 #endif
169     #ifdef ALLOW_OBCSE_CONTROL
170 heimbach 1.1.2.3 read( funit ) ((nWetobcseGlo(k,iobcs), k=1,nr),iobcs= 1,nobcs)
171 heimbach 1.1.2.2 #endif
172     cgg)
173 heimbach 1.1.2.1 read( funit ) (ncvarindex(i), i=1,maxcvars)
174     read( funit ) (ncvarrecs(i), i=1,maxcvars)
175     read( funit ) (ncvarxmax(i), i=1,maxcvars)
176     read( funit ) (ncvarymax(i), i=1,maxcvars)
177     read( funit ) (ncvarnrmax(i), i=1,maxcvars)
178     read( funit ) (ncvargrd(i), i=1,maxcvars)
179     read( funit )
180    
181     cph(
182     cph if (lheaderonly) then
183     print *, 'pathei: nvartype ', nvartype
184     print *, 'pathei: nvarlength ', nvarlength
185     print *, 'pathei: expId ', expId
186     print *, 'pathei: filenopt ', filenopt
187     print *, 'pathei: fileff ', fileff
188     print *, 'pathei: fileiG ', fileiG
189     print *, 'pathei: filejG ', filejG
190     print *, 'pathei: filensx ', filensx
191     print *, 'pathei: filensy ', filensy
192    
193 heimbach 1.1.2.3 print *, 'pathei: nWetcGlobal ',
194     & (nWetcGlobal(k), k=1,nr)
195     print *, 'pathei: nWetsGlobal ',
196     & (nWetsGlobal(k), k=1,nr)
197     print *, 'pathei: nWetwGlobal ',
198     & (nWetwGlobal(k), k=1,nr)
199     print *, 'pathei: nWetvGlobal ',
200     & (nWetvGlobal(k), k=1,nr)
201 heimbach 1.1.2.1 print *, 'pathei: ncvarindex ',
202     & (ncvarindex(i), i=1,maxcvars)
203     print *, 'pathei: ncvarrecs ',
204     & (ncvarrecs(i), i=1,maxcvars)
205     print *, 'pathei: ncvarxmax ',
206     & (ncvarxmax(i), i=1,maxcvars)
207     print *, 'pathei: ncvarymax ',
208     & (ncvarymax(i), i=1,maxcvars)
209     print *, 'pathei: ncvarnrmax ',
210     & (ncvarnrmax(i), i=1,maxcvars)
211     print *, 'pathei: ncvargrd ',
212     & (ncvargrd(i), i=1,maxcvars)
213     cph end if
214     cph)
215     c-- Check the header information for consistency.
216    
217     cph if ( filenopt .ne. nopt ) then
218     cph print*
219     cph print*,' READ_HEADER: Input data belong to the wrong'
220     cph print*,' optimization cycle.'
221     cph print*,' optimization cycle = ',nopt
222     cph print*,' input optim cycle = ',filenopt
223     cph print*
224     cph stop ' ... stopped in READ_HEADER.'
225     cph endif
226    
227     if ( (fileiG .ne. biG) .or. (filejG .ne. bjG) ) then
228     print*
229     print*,' READ_HEADER: Tile indices of loop and data '
230     print*,' do not match.'
231     print*,' loop x/y component = ',
232     & biG,bjG
233     print*,' data x/y component = ',
234     & fileiG,filejG
235     print*
236     stop ' ... stopped in READ_HEADER.'
237     endif
238    
239     if ( (filensx .ne. nsx) .or. (filensy .ne. nsy) ) then
240     print*
241     print*,' READ_HEADER: Numbers of tiles do not match.'
242     print*,' parameter x/y no. of tiles = ',
243     & bi,bj
244     print*,' data x/y no. of tiles = ',
245     & filensx,filensy
246     print*
247     stop ' ... stopped in READ_HEADER.'
248     endif
249    
250     ce Add some more checks. ...
251    
252     if (.NOT. lheaderonly) then
253     c-- Read the data.
254     icvoffset = 0
255     do icvar = 1,maxcvars
256     if ( ncvarindex(icvar) .ne. -1 ) then
257     do icvrec = 1,ncvarrecs(icvar)
258     do bj = 1,nsy
259     do bi = 1,nsx
260     read( funit ) ncvarindex(icvar)
261     read( funit ) filej
262     read( funit ) filei
263     do k = 1,ncvarnrmax(icvar)
264     cbuffindex = 0
265     if (ncvargrd(icvar) .eq. 'c') then
266 heimbach 1.1.2.3 cbuffindex = nWetcGlobal(k)
267 heimbach 1.1.2.1 else if (ncvargrd(icvar) .eq. 's') then
268 heimbach 1.1.2.3 cbuffindex = nWetsGlobal(k)
269 heimbach 1.1.2.1 else if (ncvargrd(icvar) .eq. 'w') then
270 heimbach 1.1.2.3 cbuffindex = nWetwGlobal(k)
271     else if (ncvargrd(icvar) .eq. 'v') then
272     cbuffindex = nWetvGlobal(k)
273 heimbach 1.1.2.2 cgg( O.B. points have the grid mask "m".
274     else if (ncvargrd(icvar) .eq. 'm') then
275     cgg From "icvrec", calculate what iobcs must be.
276     gg = (icvrec-1)/nobcs
277     igg = int(gg)
278     iobcs= icvrec - igg*nobcs
279     #ifdef ALLOW_OBCSN_CONTROL
280     if (icvar .eq. 11) then
281 heimbach 1.1.2.3 cbuffindex = nWetobcsnGlo(k,iobcs)
282 heimbach 1.1.2.2 endif
283     #endif
284     #ifdef ALLOW_OBCSS_CONTROL
285     if (icvar .eq. 12) then
286 heimbach 1.1.2.3 cbuffindex = nWetobcssGlo(k,iobcs)
287 heimbach 1.1.2.2 endif
288     #endif
289     #ifdef ALLOW_OBCSW_CONTROL
290     if (icvar .eq. 13) then
291 heimbach 1.1.2.3 cbuffindex = nWetobcswGlo(k,iobcs)
292 heimbach 1.1.2.2 endif
293     #endif
294     #ifdef ALLOW_OBCSE_CONTROL
295     if (icvar .eq. 14) then
296 heimbach 1.1.2.3 cbuffindex = nWetobcseGlo(k,iobcs)
297 heimbach 1.1.2.2 endif
298     #endif
299     cgg)
300 heimbach 1.1.2.1 endif
301     if (cbuffindex .gt. 0) then
302     read( funit ) cbuffindex
303     read( funit ) filek
304     read( funit ) (cbuff(ii), ii=1,cbuffindex)
305     do icvcomp = 1,cbuffindex
306     vv(icvoffset+icvcomp) = cbuff(icvcomp)
307 heimbach 1.1.2.2 cgg( Right now, the changes to the open boundary velocities are not balanced.
308     cgg( The model will crash due to physical reasons.
309     cgg( However, we can optimize with respect to just O.B. T and S if the
310     cgg( next two lines are uncommented.
311     cgg if (iobcs .eq. 3) vv(icvoffset+icvcomp)=0.
312     cgg if (iobcs .eq. 4) vv(icvoffset+icvcomp)=0.
313 heimbach 1.1.2.1 enddo
314     icvoffset = icvoffset + cbuffindex
315     endif
316     enddo
317     enddo
318     enddo
319     enddo
320     endif
321     enddo
322    
323     else
324    
325     c-- Assign the number of control variables.
326     nn = nvarlength
327    
328     endif
329    
330     close( funit )
331    
332     c-- Assign the cost function value in case we read the cost file.
333    
334     if ( dfile .eq. ctrlname ) then
335     ff = 0. d 0
336     else if ( dfile .eq. costname ) then
337     ff = fileff
338     endif
339    
340     return
341     end

  ViewVC Help
Powered by ViewVC 1.1.22