/[MITgcm]/MITgcm/pkg/ctrl/ctrl_set_unpack_xz.F
ViewVC logotype

Annotation of /MITgcm/pkg/ctrl/ctrl_set_unpack_xz.F

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


Revision 1.9 - (hide annotations) (download)
Tue Nov 16 05:42:12 2004 UTC (19 years, 6 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint57m_post, checkpoint57g_pre, checkpoint57s_post, checkpoint57b_post, checkpoint57g_post, checkpoint56b_post, checkpoint57r_post, checkpoint57d_post, checkpoint57i_post, checkpoint57, checkpoint56, checkpoint57n_post, checkpoint57l_post, checkpoint57t_post, checkpoint57v_post, checkpoint57f_post, checkpoint57a_post, checkpoint57h_pre, checkpoint57h_post, checkpoint57c_post, checkpoint57c_pre, checkpoint57e_post, checkpoint57p_post, checkpint57u_post, checkpoint57q_post, eckpoint57e_pre, checkpoint56a_post, checkpoint57h_done, checkpoint57j_post, checkpoint57f_pre, checkpoint56c_post, checkpoint57a_pre, checkpoint57o_post, checkpoint57k_post
Changes since 1.8: +0 -11 lines
More on dsvd vs. MITgcm interfacing
o handling of g_, ad, via admtlm_vector (mds...vector)
o use ctrl_pack/unpack for admtlm_vector I/O
o use optimcycle for dsvd iteration
o make sure norm is w.r.t. derived quantities

1 heimbach 1.2
2     #include "CTRL_CPPOPTIONS.h"
3    
4     subroutine ctrl_set_unpack_xz(
5 heimbach 1.3 & cunit, ivartype, fname, masktype, weighttype,
6 heimbach 1.2 & weightfld, nwetglobal, mythid)
7    
8     c ==================================================================
9     c SUBROUTINE ctrl_set_unpack_xz
10     c ==================================================================
11     c
12 heimbach 1.3 c o Unpack the control vector such that land points are filled in.
13     c
14     c o Open boundary packing added :
15     c gebbie@mit.edu, 18-Mar-2003
16     c
17     c changed: heimbach@mit.edu 17-Jun-2003
18     c merged Armin's changes to replace write of
19     c nr * globfld2d by 1 * globfld3d
20     c (ad hoc fix to speed up global I/O)
21 heimbach 1.2 c
22     c ==================================================================
23    
24     implicit none
25    
26     c == global variables ==
27    
28     #include "EEPARAMS.h"
29     #include "SIZE.h"
30     #include "PARAMS.h"
31     #include "GRID.h"
32    
33     #include "ctrl.h"
34     #include "optim.h"
35    
36     c == routine arguments ==
37    
38     integer cunit
39     integer ivartype
40     character*( 80) fname
41     character* (9) masktype
42 heimbach 1.3 character*( 80) weighttype
43 heimbach 1.2 _RL weightfld( nr,nobcs )
44     integer nwetglobal(nr,nobcs)
45     integer mythid
46    
47     c == local variables ==
48    
49     integer bi,bj
50     integer ip,jp
51     integer i,j,k
52 heimbach 1.4 integer ii,jj,kk
53 heimbach 1.2 integer il
54 heimbach 1.3 integer irec,iobcs,nrec_nl
55 heimbach 1.2 integer itlo,ithi
56     integer jtlo,jthi
57     integer jmin,jmax
58     integer imin,imax
59    
60     integer cbuffindex
61    
62 heimbach 1.7 real*4 cbuff ( snx*nsx*npx*nsy*npy )
63 heimbach 1.2 _RL globfldxz( snx,nsx,npx,nsy,npy,nr )
64 heimbach 1.3 _RL globfld3d( snx,nsx,npx,sny,nsy,npy,nr )
65     _RL globmskxz( snx,nsx,npx,nsy,npy,nr,nobcs )
66     #ifdef CTRL_UNPACK_PRECISE
67     _RL weightfldxz( snx,nsx,npx,nsy,npy,nr,nobcs )
68     #endif
69 heimbach 1.2
70     cgg(
71     integer igg
72     _RL gg
73 heimbach 1.3 character*(80) weightname
74 heimbach 1.2 cgg)
75    
76     c == external ==
77    
78     integer ilnblnk
79     external ilnblnk
80    
81     cc == end of interface ==
82    
83     jtlo = 1
84     jthi = nsy
85     itlo = 1
86     ithi = nsx
87     jmin = 1
88     jmax = sny
89     imin = 1
90     imax = snx
91    
92     c Initialise temporary file
93     do k = 1,nr
94     do jp = 1,nPy
95     do bj = jtlo,jthi
96     do ip = 1,nPx
97     do bi = itlo,ithi
98     do i = imin,imax
99     globfldxz(i,bi,ip,bj,jp,k) = 0. _d 0
100 heimbach 1.3 do iobcs=1,nobcs
101     globmskxz(i,bi,ip,bj,jp,k,iobcs) = 0. _d 0
102     enddo
103     enddo
104     enddo
105     enddo
106     enddo
107     enddo
108     enddo
109     c Initialise temporary file
110     do k = 1,nr
111     do jp = 1,nPy
112     do bj = jtlo,jthi
113     do j = jmin,jmax
114     do ip = 1,nPx
115     do bi = itlo,ithi
116     do i = imin,imax
117     globfld3d(i,bi,ip,j,bj,jp,k) = 0. _d 0
118     enddo
119 heimbach 1.2 enddo
120     enddo
121     enddo
122     enddo
123     enddo
124     enddo
125    
126     c-- Only the master thread will do I/O.
127     _BEGIN_MASTER( mythid )
128    
129 heimbach 1.3 do iobcs=1,nobcs
130     call MDSREADFIELD_XZ_GL(
131     & masktype, ctrlprec, 'RL',
132     & Nr, globmskxz(1,1,1,1,1,1,iobcs), iobcs,mythid)
133     #ifdef CTRL_UNPACK_PRECISE
134     il=ilnblnk( weighttype)
135     write(weightname(1:80),'(80a)') ' '
136     write(weightname(1:80),'(a)') weighttype(1:il)
137     call MDSREADFIELD_XZ_GL(
138     & weightname, ctrlprec, 'RL',
139     & Nr, weightfldxz(1,1,1,1,1,1,iobcs), iobcs, mythid)
140     CGG One special exception: barotropic velocity should be nondimensionalized
141     cgg differently. Probably introduce new variable.
142     if (iobcs .eq. 3 .or. iobcs .eq. 4) then
143     k = 1
144     do jp = 1,nPy
145     do bj = jtlo,jthi
146     do ip = 1,nPx
147     do bi = itlo,ithi
148     do i = imin,imax
149     weightfldxz(i,bi,ip,bj,jp,k,iobcs) = wbaro
150     enddo
151     enddo
152     enddo
153     enddo
154     enddo
155     endif
156     #endif /* CTRL_UNPACK_PRECISE */
157     enddo
158    
159     nrec_nl=int(ncvarrecs(ivartype)/sny)
160     do irec = 1, nrec_nl
161 heimbach 1.2 cgg do iobcs = 1, nobcs
162 heimbach 1.3 cgg And now back-calculate what iobcs should be.
163     do j=1,sny
164     iobcs= mod((irec-1)*sny+j-1,nobcs)+1
165    
166     read(cunit) filencvarindex(ivartype)
167     if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
168     & then
169     print *, 'ctrl-set_unpack:xz:WARNING: wrong ncvarindex ',
170     & filencvarindex(ivartype), ncvarindex(ivartype)
171     STOP 'in S/R ctrl_unpack'
172     endif
173     read(cunit) filej
174     read(cunit) filei
175     do k = 1, Nr
176     cbuffindex = nwetglobal(k,iobcs)
177     if ( cbuffindex .gt. 0 ) then
178     read(cunit) filencbuffindex
179     if (filencbuffindex .NE. cbuffindex) then
180     print *, 'WARNING: wrong cbuffindex ',
181     & filencbuffindex, cbuffindex
182     STOP 'in S/R ctrl_unpack'
183     endif
184     read(cunit) filek
185     if (filek .NE. k) then
186     print *, 'WARNING: wrong k ',
187     & filek, k
188     STOP 'in S/R ctrl_unpack'
189     endif
190     read(cunit) (cbuff(ii), ii=1,cbuffindex)
191     endif
192     cbuffindex = 0
193 heimbach 1.7 jj=mod((j-1)*nr+k-1,sny)+1
194     kk=int((j-1)*nr+k-1)/sny+1
195 heimbach 1.3 do jp = 1,nPy
196     do bj = jtlo,jthi
197     do ip = 1,nPx
198     do bi = itlo,ithi
199     do i = imin,imax
200     if ( globmskxz(i,bi,ip,bj,jp,k,iobcs) .ne. 0. ) then
201     cbuffindex = cbuffindex + 1
202 heimbach 1.4 globfld3d(i,bi,ip,jj,bj,jp,kk) =
203     & cbuff(cbuffindex)
204 heimbach 1.3 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
205 heimbach 1.4 globfld3d(i,bi,ip,jj,bj,jp,kk) =
206     & globfld3d(i,bi,ip,jj,bj,jp,kk)/
207 heimbach 1.3 # ifdef CTRL_UNPACK_PRECISE
208     & sqrt(weightfldxz(i,bi,ip,bj,jp,k,iobcs))
209     # else
210     & sqrt(weightfld(k,iobcs))
211     # endif
212     #endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
213     else
214 heimbach 1.4 globfld3d(i,bi,ip,jj,bj,jp,kk) = 0. _d 0
215 heimbach 1.3 endif
216     enddo
217     enddo
218     enddo
219     enddo
220     enddo
221     c
222     c -- end of k loop --
223     enddo
224     c -- end of j loop --
225     enddo
226    
227     call MDSWRITEFIELD_3D_GL( fname, ctrlprec, 'RL',
228     & Nr, globfld3d, irec,
229     & optimcycle, mythid)
230    
231     c -- end of iobcs loop -- This loop removed. 3-28-02.
232     cgg enddo
233     c -- end of irec loop --
234     enddo
235 heimbach 1.2
236 heimbach 1.3 do irec = nrec_nl*sny+1, ncvarrecs(ivartype)
237     cgg do iobcs = 1, nobcs
238     cgg And now back-calculate what iobcs should be.
239     iobcs= mod(irec-1,nobcs)+1
240 heimbach 1.2
241     read(cunit) filencvarindex(ivartype)
242     if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
243     & then
244     print *, 'ctrl-set_unpack:xz:WARNING: wrong ncvarindex ',
245     & filencvarindex(ivartype), ncvarindex(ivartype)
246     STOP 'in S/R ctrl_unpack'
247     endif
248     read(cunit) filej
249     read(cunit) filei
250     do k = 1, Nr
251     cbuffindex = nwetglobal(k,iobcs)
252     if ( cbuffindex .gt. 0 ) then
253     read(cunit) filencbuffindex
254     if (filencbuffindex .NE. cbuffindex) then
255     print *, 'WARNING: wrong cbuffindex ',
256     & filencbuffindex, cbuffindex
257     STOP 'in S/R ctrl_unpack'
258     endif
259     read(cunit) filek
260     if (filek .NE. k) then
261     print *, 'WARNING: wrong k ',
262     & filek, k
263     STOP 'in S/R ctrl_unpack'
264     endif
265     read(cunit) (cbuff(ii), ii=1,cbuffindex)
266     endif
267     cbuffindex = 0
268     do jp = 1,nPy
269     do bj = jtlo,jthi
270     do ip = 1,nPx
271     do bi = itlo,ithi
272     do i = imin,imax
273 heimbach 1.3 if ( globmskxz(i,bi,ip,bj,jp,k,iobcs) .ne. 0. ) then
274 heimbach 1.2 cbuffindex = cbuffindex + 1
275     globfldxz(i,bi,ip,bj,jp,k) = cbuff(cbuffindex)
276     #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
277     globfldxz(i,bi,ip,bj,jp,k) =
278     & globfldxz(i,bi,ip,bj,jp,k)/
279 heimbach 1.3 # ifdef CTRL_UNPACK_PRECISE
280     & sqrt(weightfldxz(i,bi,ip,bj,jp,k,iobcs))
281     # else
282 heimbach 1.2 & sqrt(weightfld(k,iobcs))
283 heimbach 1.3 # endif
284     #endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
285 heimbach 1.2 else
286     globfldxz(i,bi,ip,bj,jp,k) = 0. _d 0
287     endif
288     enddo
289     enddo
290     enddo
291     enddo
292     enddo
293     c
294 heimbach 1.3 c -- end of k loop --
295 heimbach 1.2 enddo
296    
297     call MDSWRITEFIELD_XZ_GL( fname, ctrlprec, 'RL',
298     & Nr, globfldxz, irec,
299     & optimcycle, mythid)
300    
301     c -- end of iobcs loop -- This loop removed. 3-28-02.
302     cgg enddo
303     c -- end of irec loop --
304     enddo
305    
306     _END_MASTER( mythid )
307    
308     return
309     end
310    
311    
312    
313    
314    

  ViewVC Help
Powered by ViewVC 1.1.22