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

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

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


Revision 1.13 - (hide annotations) (download)
Tue Nov 1 04:09:46 2005 UTC (18 years, 6 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint57y_post, checkpoint58, checkpoint57z_post, checkpoint57y_pre, checkpoint57w_post, checkpoint57x_post
Changes since 1.12: +9 -0 lines
Completely restructured the arpack2model interface.
Now (again) only 1-d wetpoint vector is passed to ARPACK.
ctrl_unpack/pack are mimiced by admtlm_dsvd2model/model2dsvd

1 heimbach 1.2
2     #include "CTRL_CPPOPTIONS.h"
3    
4     subroutine ctrl_set_unpack_xyz(
5 heimbach 1.9 & lxxadxx, cunit, ivartype, fname, masktype, weighttype,
6 heimbach 1.2 & weightfld, nwetglobal, mythid)
7    
8     c ==================================================================
9     c SUBROUTINE ctrl_set_unpack_xyz
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 Use a more precise nondimensionalization that depends on (x,y)
15     c Added weighttype to the argument list so that I can geographically
16     c vary the nondimensionalization.
17     c gebbie@mit.edu, 18-Mar-2003
18 heimbach 1.2 c
19     c ==================================================================
20    
21     implicit none
22    
23     c == global variables ==
24    
25     #include "EEPARAMS.h"
26     #include "SIZE.h"
27     #include "PARAMS.h"
28     #include "GRID.h"
29    
30     #include "ctrl.h"
31     #include "optim.h"
32    
33     c == routine arguments ==
34    
35 heimbach 1.9 logical lxxadxx
36 heimbach 1.2 integer cunit
37     integer ivartype
38     character*( 80) fname
39 heimbach 1.11 character*( 9) masktype
40 heimbach 1.3 character*( 80) weighttype
41 heimbach 1.2 _RL weightfld( nr,nsx,nsy )
42     integer nwetglobal(nr)
43     integer mythid
44    
45     c == local variables ==
46    
47     integer bi,bj
48     integer ip,jp
49     integer i,j,k
50     integer ii
51     integer il
52     integer irec
53     integer itlo,ithi
54     integer jtlo,jthi
55     integer jmin,jmax
56     integer imin,imax
57    
58     integer cbuffindex
59    
60     _RL globmsk ( snx,nsx,npx,sny,nsy,npy,nr )
61     _RL globfld3d( snx,nsx,npx,sny,nsy,npy,nr )
62 heimbach 1.3 #ifdef CTRL_UNPACK_PRECISE
63     _RL weightfld3d( snx,nsx,npx,sny,nsy,npy,nr )
64     #endif
65 heimbach 1.9 real*4 cbuff ( snx*nsx*npx*sny*nsy*npy )
66     real*4 globfldtmp2( snx,nsx,npx,sny,nsy,npy )
67     real*4 globfldtmp3( snx,nsx,npx,sny,nsy,npy )
68 heimbach 1.2
69     character*(128) cfile
70 heimbach 1.3 character*(80) weightname
71 heimbach 1.2
72 heimbach 1.9 _RL delZnorm
73     integer reclen, irectrue
74     integer cunit2, cunit3
75     character*(80) cfile2, cfile3
76    
77 heimbach 1.2 c == external ==
78    
79     integer ilnblnk
80     external ilnblnk
81    
82     cc == end of interface ==
83    
84     jtlo = 1
85     jthi = nsy
86     itlo = 1
87     ithi = nsx
88     jmin = 1
89     jmax = sny
90     imin = 1
91     imax = snx
92    
93 heimbach 1.9 #ifdef CTRL_DELZNORM
94     delZnorm = 0.
95     do k = 1, Nr
96     delZnorm = delZnorm + delR(k)/FLOAT(Nr)
97     enddo
98     #endif
99    
100 heimbach 1.2 c Initialise temporary file
101     do k = 1,nr
102     do jp = 1,nPy
103     do bj = jtlo,jthi
104     do j = jmin,jmax
105     do ip = 1,nPx
106     do bi = itlo,ithi
107     do i = imin,imax
108 heimbach 1.9 globfld3d (i,bi,ip,j,bj,jp,k) = 0. _d 0
109     globmsk (i,bi,ip,j,bj,jp,k) = 0. _d 0
110     globfldtmp2(i,bi,ip,j,bj,jp) = 0.
111     globfldtmp3(i,bi,ip,j,bj,jp) = 0.
112 heimbach 1.2 enddo
113     enddo
114     enddo
115     enddo
116     enddo
117     enddo
118     enddo
119    
120     c-- Only the master thread will do I/O.
121     _BEGIN_MASTER( mythid )
122    
123 heimbach 1.9 #ifdef CTRL_DELZNORM
124     do k = 1, nr
125     print *, 'ph-delznorm ', k, delZnorm, delR(k)
126     print *, 'ph-weight ', weightfld(k,1,1)
127     enddo
128     #endif
129    
130     if ( doPackDiag ) then
131     write(cfile2(1:80),'(80a)') ' '
132     write(cfile3(1:80),'(80a)') ' '
133     if ( lxxadxx ) then
134     write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
135     & 'diag_unpack_nondim_ctrl_',
136     & ivartype, '_', optimcycle, '.bin'
137     write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
138     & 'diag_unpack_dimens_ctrl_',
139     & ivartype, '_', optimcycle, '.bin'
140     else
141     write(cfile2(1:80),'(a,I2.2,a,I4.4,a)')
142     & 'diag_unpack_nondim_grad_',
143     & ivartype, '_', optimcycle, '.bin'
144     write(cfile3(1:80),'(a,I2.2,a,I4.4,a)')
145     & 'diag_unpack_dimens_grad_',
146     & ivartype, '_', optimcycle, '.bin'
147     endif
148    
149     reclen = FLOAT(snx*nsx*npx*sny*nsy*npy*4)
150     call mdsfindunit( cunit2, mythid )
151 heimbach 1.10 open( cunit2, file=cfile2, status='unknown',
152 heimbach 1.9 & access='direct', recl=reclen )
153     call mdsfindunit( cunit3, mythid )
154 heimbach 1.10 open( cunit3, file=cfile3, status='unknown',
155 heimbach 1.9 & access='direct', recl=reclen )
156     endif
157    
158 heimbach 1.3 #ifdef CTRL_UNPACK_PRECISE
159     il=ilnblnk( weighttype)
160     write(weightname(1:80),'(80a)') ' '
161     write(weightname(1:80),'(a)') weighttype(1:il)
162    
163     call MDSREADFIELD_3D_GL(
164     & weightname, ctrlprec, 'RL',
165     & Nr, weightfld3d, 1, mythid)
166     #endif
167    
168 heimbach 1.2 call MDSREADFIELD_3D_GL(
169     & masktype, ctrlprec, 'RL',
170     & Nr, globmsk, 1, mythid)
171    
172     do irec = 1, ncvarrecs(ivartype)
173 heimbach 1.13 #ifndef ALLOW_ADMTLM
174 heimbach 1.2 read(cunit) filencvarindex(ivartype)
175     if (filencvarindex(ivartype) .NE. ncvarindex(ivartype))
176     & then
177     print *, 'ctrl_set_unpack_xyz:WARNING: wrong ncvarindex ',
178     & filencvarindex(ivartype), ncvarindex(ivartype)
179     STOP 'in S/R ctrl_unpack'
180     endif
181     read(cunit) filej
182     read(cunit) filei
183 heimbach 1.13 #endif /* ALLOW_ADMTLM */
184 heimbach 1.2 do k = 1, Nr
185 heimbach 1.9 irectrue = (irec-1)*nr + k
186 heimbach 1.10 if ( doZscaleUnpack ) then
187 heimbach 1.12 cph delZnorm = SQRT(delR(1)/delR(k))
188     delZnorm = delR(1)/delR(k)
189 heimbach 1.10 else
190     delZnorm = 1. _d 0
191     endif
192 heimbach 1.2 cbuffindex = nwetglobal(k)
193     if ( cbuffindex .gt. 0 ) then
194 heimbach 1.13 #ifndef ALLOW_ADMTLM
195 heimbach 1.2 read(cunit) filencbuffindex
196     if (filencbuffindex .NE. cbuffindex) then
197     print *, 'WARNING: wrong cbuffindex ',
198     & filencbuffindex, cbuffindex
199     STOP 'in S/R ctrl_unpack'
200     endif
201     read(cunit) filek
202     if (filek .NE. k) then
203     print *, 'WARNING: wrong k ',
204     & filek, k
205     STOP 'in S/R ctrl_unpack'
206     endif
207 heimbach 1.13 #endif /* ALLOW_ADMTLM */
208 heimbach 1.2 read(cunit) (cbuff(ii), ii=1,cbuffindex)
209     endif
210 heimbach 1.13 c
211 heimbach 1.2 cbuffindex = 0
212     do jp = 1,nPy
213     do bj = jtlo,jthi
214     do j = jmin,jmax
215     do ip = 1,nPx
216     do bi = itlo,ithi
217     do i = imin,imax
218     if ( globmsk(i,bi,ip,j,bj,jp,k) .ne. 0. ) then
219     cbuffindex = cbuffindex + 1
220     globfld3d(i,bi,ip,j,bj,jp,k) = cbuff(cbuffindex)
221 heimbach 1.9 cph(
222     globfldtmp2(i,bi,ip,j,bj,jp) = cbuff(cbuffindex)
223     cph)
224 heimbach 1.13 #ifdef ALLOW_ADMTLM
225     nveccount = nveccount + 1
226     phtmpadmtlm(nveccount) = cbuff(cbuffindex)
227     #endif
228 heimbach 1.2 #ifdef ALLOW_NONDIMENSIONAL_CONTROL_IO
229 heimbach 1.9 if ( lxxadxx ) then
230 heimbach 1.10 globfld3d(i,bi,ip,j,bj,jp,k) = delZnorm
231     & * globfld3d(i,bi,ip,j,bj,jp,k)
232 heimbach 1.9 # ifdef CTRL_UNPACK_PRECISE
233 heimbach 1.10 & / sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
234 heimbach 1.9 # else
235 heimbach 1.10 & / sqrt(weightfld(k,bi,bj))
236 heimbach 1.9 # endif
237     else
238 heimbach 1.10 globfld3d(i,bi,ip,j,bj,jp,k) = delZnorm
239     & * globfld3d(i,bi,ip,j,bj,jp,k)
240 heimbach 1.3 # ifdef CTRL_UNPACK_PRECISE
241 heimbach 1.9 & * sqrt(weightfld3d(i,bi,ip,j,bj,jp,k))
242 heimbach 1.3 # else
243 heimbach 1.9 & * sqrt(weightfld(k,bi,bj))
244 heimbach 1.3 # endif
245 heimbach 1.9 endif
246 heimbach 1.3 #endif /* ALLOW_NONDIMENSIONAL_CONTROL_IO */
247 heimbach 1.2 else
248     globfld3d(i,bi,ip,j,bj,jp,k) = 0. _d 0
249     endif
250 heimbach 1.9 cph(
251     globfldtmp3(i,bi,ip,j,bj,jp) =
252     & globfld3d(i,bi,ip,j,bj,jp,k)
253     cph)
254 heimbach 1.2 enddo
255     enddo
256     enddo
257     enddo
258     enddo
259     enddo
260     c
261 heimbach 1.9 if ( doPackDiag ) then
262     write(cunit2,rec=irectrue) globfldtmp2
263     write(cunit3,rec=irectrue) globfldtmp3
264     endif
265     c
266 heimbach 1.2 enddo
267    
268     call MDSWRITEFIELD_3D_GL( fname, ctrlprec, 'RL',
269     & Nr, globfld3d,
270     & irec, optimcycle, mythid)
271    
272     enddo
273    
274 heimbach 1.9 if ( doPackDiag ) then
275     close ( cunit2 )
276     close ( cunit3 )
277     endif
278    
279 heimbach 1.2 _END_MASTER( mythid )
280    
281     return
282     end
283    

  ViewVC Help
Powered by ViewVC 1.1.22