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

Contents of /MITgcm/pkg/ctrl/ctrl_map_ini.F

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


Revision 1.18 - (show annotations) (download)
Tue Jan 3 17:10:35 2006 UTC (19 years, 6 months ago) by heimbach
Branch: MAIN
Changes since 1.17: +11 -1 lines
o Adding special treatment of seeding control variables for OpenAD
  to facilitate a working version (meant as a temporary solution).

1 C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_map_ini.F,v 1.17 2005/04/07 23:38:43 heimbach Exp $
2
3 #include "CTRL_CPPOPTIONS.h"
4
5 CBOP
6 C !ROUTINE: ctrl_map_ini
7 C !INTERFACE:
8 subroutine ctrl_map_ini( mythid )
9
10 C !DESCRIPTION: \bv
11 c *=================================================================
12 c | SUBROUTINE ctrl_map_ini
13 c | Add the temperature, salinity, and diffusivity parts of the
14 c | control vector to the model state and update the tile halos.
15 c | The control vector is defined in the header file "ctrl.h".
16 c *=================================================================
17 C \ev
18
19 C !USES:
20 implicit none
21
22 c == global variables ==
23 #include "SIZE.h"
24 #include "EEPARAMS.h"
25 #include "PARAMS.h"
26 #include "GRID.h"
27 #include "DYNVARS.h"
28 #include "FFIELDS.h"
29 #include "ctrl.h"
30 #include "ctrl_dummy.h"
31 #include "optim.h"
32 #ifdef ALLOW_PTRACERS
33 # include "PTRACERS_SIZE.h"
34 # include "PTRACERS.h"
35 #endif
36 #ifdef ALLOW_ECCO
37 # include "ecco_cost.h"
38 #endif
39
40 C !INPUT/OUTPUT PARAMETERS:
41 c == routine arguments ==
42 integer mythid
43
44 C !LOCAL VARIABLES:
45 c == local variables ==
46
47 integer bi,bj
48 integer i,j,k
49 integer itlo,ithi
50 integer jtlo,jthi
51 integer jmin,jmax
52 integer imin,imax
53 integer il
54
55 logical equal
56 logical doglobalread
57 logical ladinit
58
59 character*( 80) fnamegeneric
60
61 _RL fac
62 _RL tmptest
63
64 c == external ==
65 integer ilnblnk
66 external ilnblnk
67
68 c == end of interface ==
69 CEOP
70
71 jtlo = mybylo(mythid)
72 jthi = mybyhi(mythid)
73 itlo = mybxlo(mythid)
74 ithi = mybxhi(mythid)
75 jmin = 1
76 jmax = sny
77 imin = 1
78 imax = snx
79
80 doglobalread = .false.
81 ladinit = .false.
82
83 equal = .true.
84
85 if ( equal ) then
86 fac = 1. _d 0
87 else
88 fac = 0. _d 0
89 endif
90
91 #ifdef ALLOW_THETA0_CONTROL
92 c-- Temperature field.
93 il=ilnblnk( xx_theta_file )
94 write(fnamegeneric(1:80),'(2a,i10.10)')
95 & xx_theta_file(1:il),'.',optimcycle
96 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
97 & doglobalread, ladinit, optimcycle,
98 & mythid, xx_theta_dummy )
99
100 do bj = jtlo,jthi
101 do bi = itlo,ithi
102 do k = 1,nr
103 do j = jmin,jmax
104 do i = imin,imax
105 #ifdef ALLOW_ECCO
106 IF (abs(tmpfld3d(i,j,k,bi,bj)).gt.
107 $ 2.0/sqrt(wtheta(k,bi,bj)))
108 $ tmpfld3d(i,j,k,bi,bj)=
109 $ sign(2.0/sqrt(wtheta(k,bi,bj)),tmpfld3d(i,j,k,bi,bj))
110 #endif
111 #ifdef ALLOW_OPENAD
112 theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
113 & fac*xx_theta(i,j,k,bi,bj)
114 #else
115 theta(i,j,k,bi,bj) = theta(i,j,k,bi,bj) +
116 & fac*tmpfld3d(i,j,k,bi,bj)
117 #endif
118 if(theta(i,j,k,bi,bj).lt.-2.0)
119 & theta(i,j,k,bi,bj)= -2.0
120 enddo
121 enddo
122 enddo
123 enddo
124 enddo
125
126 #endif
127
128 #ifdef ALLOW_SALT0_CONTROL
129 c-- Temperature field.
130 il=ilnblnk( xx_salt_file )
131 write(fnamegeneric(1:80),'(2a,i10.10)')
132 & xx_salt_file(1:il),'.',optimcycle
133 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
134 & doglobalread, ladinit, optimcycle,
135 & mythid, xx_salt_dummy )
136
137 do bj = jtlo,jthi
138 do bi = itlo,ithi
139 do k = 1,nr
140 do j = jmin,jmax
141 do i = imin,imax
142 #ifdef ALLOW_ECCO
143 IF (abs(tmpfld3d(i,j,k,bi,bj)).gt.
144 $ 2.0/sqrt(wsalt(k,bi,bj)))
145 $ tmpfld3d(i,j,k,bi,bj)=
146 $ sign(2.0/sqrt(wsalt(k,bi,bj)),tmpfld3d(i,j,k,bi,bj))
147 #endif
148 #ifdef ALLOW_OPENAD
149 salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
150 & fac*xx_salt(i,j,k,bi,bj)
151 #else
152 salt(i,j,k,bi,bj) = salt(i,j,k,bi,bj) +
153 & fac*tmpfld3d(i,j,k,bi,bj)
154 #endif
155
156 enddo
157 enddo
158 enddo
159 enddo
160 enddO
161 #endif
162
163 #ifdef ALLOW_TR10_CONTROL
164 #ifdef ALLOW_PTRACERS
165 c-- Temperature field.
166 il=ilnblnk( xx_tr1_file )
167 write(fnamegeneric(1:80),'(2a,i10.10)')
168 & xx_tr1_file(1:il),'.',optimcycle
169 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
170 & doglobalread, ladinit, optimcycle,
171 & mythid, xx_tr1_dummy )
172
173 do bj = jtlo,jthi
174 do bi = itlo,ithi
175 do k = 1,nr
176 do j = jmin,jmax
177 do i = imin,imax
178 ptracer(i,j,k,bi,bj,1) = ptracer(i,j,k,bi,bj,1) +
179 & fac*tmpfld3d(i,j,k,bi,bj)
180 enddo
181 enddo
182 enddo
183 enddo
184 enddo
185 #endif
186 #endif
187
188 #ifdef ALLOW_SST0_CONTROL
189 c-- sst0.
190 il=ilnblnk( xx_sst_file )
191 write(fnamegeneric(1:80),'(2a,i10.10)')
192 & xx_sst_file(1:il),'.',optimcycle
193 call active_read_xy_loc ( fnamegeneric, tmpfld2d, 1,
194 & doglobalread, ladinit, optimcycle,
195 & mythid, xx_sst_dummy )
196 do bj = jtlo,jthi
197 do bi = itlo,ithi
198 do j = jmin,jmax
199 do i = imin,imax
200 cph sst(i,j,bi,bj) = sst(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
201 theta(i,j,1,bi,bj) = theta(i,j,1,bi,bj)
202 & + tmpfld2d(i,j,bi,bj)
203 enddo
204 enddo
205 enddo
206 enddo
207 #endif
208
209 #ifdef ALLOW_SSS0_CONTROL
210 c-- sss0.
211 il=ilnblnk( xx_sss_file )
212 write(fnamegeneric(1:80),'(2a,i10.10)')
213 & xx_sss_file(1:il),'.',optimcycle
214 call active_read_xy_loc ( fnamegeneric, tmpfld2d, 1,
215 & doglobalread, ladinit, optimcycle,
216 & mythid, xx_sss_dummy )
217 do bj = jtlo,jthi
218 do bi = itlo,ithi
219 do j = jmin,jmax
220 do i = imin,imax
221 cph sss(i,j,bi,bj) = sss(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
222 salt(i,j,1,bi,bj) = salt(i,j,1,bi,bj)
223 & + tmpfld2d(i,j,bi,bj)
224 enddo
225 enddo
226 enddo
227 enddo
228 #endif
229
230 #ifdef ALLOW_DIFFKR_CONTROL
231 c-- diffkr.
232 il=ilnblnk( xx_diffkr_file )
233 write(fnamegeneric(1:80),'(2a,i10.10)')
234 & xx_diffkr_file(1:il),'.',optimcycle
235 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
236 & doglobalread, ladinit, optimcycle,
237 & mythid, xx_diffkr_dummy )
238 do bj = jtlo,jthi
239 do bi = itlo,ithi
240 do k = 1,nr
241 do j = jmin,jmax
242 do i = imin,imax
243 diffkr(i,j,k,bi,bj) = diffkr(i,j,k,bi,bj) +
244 & tmpfld3d(i,j,k,bi,bj)
245 enddo
246 enddo
247 enddo
248 enddo
249 enddo
250 #endif
251
252 #ifdef ALLOW_KAPGM_CONTROL
253 c-- kapgm.
254 il=ilnblnk( xx_kapgm_file )
255 write(fnamegeneric(1:80),'(2a,i10.10)')
256 & xx_kapgm_file(1:il),'.',optimcycle
257 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
258 & doglobalread, ladinit, optimcycle,
259 & mythid, xx_kapgm_dummy )
260 do bj = jtlo,jthi
261 do bi = itlo,ithi
262 do k = 1,nr
263 do j = jmin,jmax
264 do i = imin,imax
265 kapgm(i,j,k,bi,bj) = kapgm(i,j,k,bi,bj) +
266 & tmpfld3d(i,j,k,bi,bj)
267 enddo
268 enddo
269 enddo
270 enddo
271 enddo
272 #endif
273
274 #ifdef ALLOW_EFLUXY0_CONTROL
275 c-- y-component EP-flux field.
276 il=ilnblnk( xx_efluxy_file )
277 write(fnamegeneric(1:80),'(2a,i10.10)')
278 & xx_efluxy_file(1:il),'.',optimcycle
279 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
280 & doglobalread, ladinit, optimcycle,
281 & mythid, xx_efluxy_dummy )
282
283 do bj = jtlo,jthi
284 do bi = itlo,ithi
285 do k = 1,nr
286 do j = jmin,jmax
287 do i = imin,imax
288 EfluxY(i,j,k,bi,bj) = EfluxY(i,j,k,bi,bj)
289 & - fac*tmpfld3d(i,j,k,bi,bj)
290 & *maskS(i,j,k,bi,bj)
291 cph EfluxY(i,j,k,bi,bj) = EfluxY(i,j,k,bi,bj)
292 cph & - rSphere*cosFacU(J,bi,bj)
293 cph & *fac*tmpfld3d(i,j,k,bi,bj)
294 enddo
295 enddo
296 enddo
297 enddo
298 enddo
299 #endif
300
301 #ifdef ALLOW_EFLUXP0_CONTROL
302 c-- p-component EP-flux field.
303 il=ilnblnk( xx_efluxp_file )
304 write(fnamegeneric(1:80),'(2a,i10.10)')
305 & xx_efluxp_file(1:il),'.',optimcycle
306 call active_read_xyz_loc( fnamegeneric, tmpfld3d, 1,
307 & doglobalread, ladinit, optimcycle,
308 & mythid, xx_efluxp_dummy )
309
310 do bj = jtlo,jthi
311 do bi = itlo,ithi
312 do k = 1,nr
313 do j = jmin,jmax
314 do i = imin,imax
315 EfluxP(i,j,k,bi,bj) = EfluxP(i,j,k,bi,bj)
316 & + fCori(i,j,bi,bj)
317 & *fac*tmpfld3d(i,j,k,bi,bj)
318 & *hFacV(i,j,k,bi,bj)
319 cph EfluxP(i,j,k,bi,bj) = EfluxP(i,j,k,bi,bj)
320 cph & + fCori(i,j,bi,bj)
321 cph & *rSphere*cosFacU(J,bi,bj)
322 cph & *fac*tmpfld3d(i,j,k,bi,bj)
323 enddo
324 enddo
325 enddo
326 enddo
327 enddo
328 #endif
329
330 #ifdef ALLOW_BOTTOMDRAG_CONTROL
331 c-- bottom drag
332 il=ilnblnk( xx_bottomdrag_file )
333 write(fnamegeneric(1:80),'(2a,i10.10)')
334 & xx_bottomdrag_file(1:il),'.',optimcycle
335 call active_read_xy_loc ( fnamegeneric, tmpfld2d, 1,
336 & doglobalread, ladinit, optimcycle,
337 & mythid, xx_bottomdrag_dummy )
338 do bj = jtlo,jthi
339 do bi = itlo,ithi
340 do j = jmin,jmax
341 do i = imin,imax
342 bottomdragfld(i,j,bi,bj) = bottomdragfld(i,j,bi,bj)
343 & + tmpfld2d(i,j,bi,bj)
344 enddo
345 enddo
346 enddo
347 enddo
348 #endif
349
350 #ifdef ALLOW_EDTAUX_CONTROL
351 c-- zonal eddy stress : edtaux
352 il=ilnblnk( xx_edtaux_file )
353 write(fnamegeneric(1:80),'(2a,i10.10)')
354 & xx_edtaux_file(1:il),'.',optimcycle
355 call active_read_xyz( fnamegeneric, tmpfld3d, 1,
356 & doglobalread, ladinit, optimcycle,
357 & mythid, xx_edtaux_dummy )
358 do bj = jtlo,jthi
359 do bi = itlo,ithi
360 do k = 1,nr
361 do j = jmin,jmax
362 do i = imin,imax
363 Eddytaux(i,j,k,bi,bj) = Eddytaux(i,j,k,bi,bj) +
364 & tmpfld3d(i,j,k,bi,bj)
365 enddo
366 enddo
367 enddo
368 enddo
369 enddo
370 #endif
371
372 #ifdef ALLOW_EDTAUY_CONTROL
373 c-- meridional eddy stress : edtauy
374 il=ilnblnk( xx_edtauy_file )
375 write(fnamegeneric(1:80),'(2a,i10.10)')
376 & xx_edtauy_file(1:il),'.',optimcycle
377 call active_read_xyz( fnamegeneric, tmpfld3d, 1,
378 & doglobalread, ladinit, optimcycle,
379 & mythid, xx_edtauy_dummy )
380 do bj = jtlo,jthi
381 do bi = itlo,ithi
382 do k = 1,nr
383 do j = jmin,jmax
384 do i = imin,imax
385 Eddytauy(i,j,k,bi,bj) = Eddytauy(i,j,k,bi,bj) +
386 & tmpfld3d(i,j,k,bi,bj)
387 enddo
388 enddo
389 enddo
390 enddo
391 enddo
392 #endif
393
394 #ifdef ALLOW_UVEL0_CONTROL
395 c-- initial zonal velocity
396 il=ilnblnk( xx_uvel_file )
397 write(fnamegeneric(1:80),'(2a,i10.10)')
398 & xx_uvel_file(1:il),'.',optimcycle
399 call active_read_xyz( fnamegeneric, tmpfld3d, 1,
400 & doglobalread, ladinit, optimcycle,
401 & mythid, xx_uvel_dummy )
402 do bj = jtlo,jthi
403 do bi = itlo,ithi
404 do k = 1,nr
405 do j = jmin,jmax
406 do i = imin,imax
407 uVel(i,j,k,bi,bj) = uVel(i,j,k,bi,bj) +
408 & tmpfld3d(i,j,k,bi,bj)
409 enddo
410 enddo
411 enddo
412 enddo
413 enddo
414 #endif
415
416 #ifdef ALLOW_VVEL0_CONTROL
417 c-- initial merid. velocity
418 il=ilnblnk( xx_vvel_file )
419 write(fnamegeneric(1:80),'(2a,i10.10)')
420 & xx_vvel_file(1:il),'.',optimcycle
421 call active_read_xyz( fnamegeneric, tmpfld3d, 1,
422 & doglobalread, ladinit, optimcycle,
423 & mythid, xx_vvel_dummy )
424 do bj = jtlo,jthi
425 do bi = itlo,ithi
426 do k = 1,nr
427 do j = jmin,jmax
428 do i = imin,imax
429 vVel(i,j,k,bi,bj) = vVel(i,j,k,bi,bj) +
430 & tmpfld3d(i,j,k,bi,bj)
431 enddo
432 enddo
433 enddo
434 enddo
435 enddo
436 #endif
437
438 #ifdef ALLOW_ETAN0_CONTROL
439 c-- initial Eta.
440 il=ilnblnk( xx_etan_file )
441 write(fnamegeneric(1:80),'(2a,i10.10)')
442 & xx_etan_file(1:il),'.',optimcycle
443 call active_read_xy_loc ( fnamegeneric, tmpfld2d, 1,
444 & doglobalread, ladinit, optimcycle,
445 & mythid, xx_etan_dummy )
446 do bj = jtlo,jthi
447 do bi = itlo,ithi
448 do j = jmin,jmax
449 do i = imin,imax
450 etaN(i,j,bi,bj) = etaN(i,j,bi,bj) + tmpfld2d(i,j,bi,bj)
451 enddo
452 enddo
453 enddo
454 enddo
455 #endif
456
457 #ifdef ALLOW_RELAXSST_CONTROL
458 c-- SST relaxation coefficient.
459 il=ilnblnk( xx_relaxsst_file )
460 write(fnamegeneric(1:80),'(2a,i10.10)')
461 & xx_relaxsst_file(1:il),'.',optimcycle
462 call active_read_xy_loc ( fnamegeneric, tmpfld2d, 1,
463 & doglobalread, ladinit, optimcycle,
464 & mythid, xx_relaxsst_dummy )
465 do bj = jtlo,jthi
466 do bi = itlo,ithi
467 do j = jmin,jmax
468 do i = imin,imax
469 lambdaThetaClimRelax(i,j,bi,bj) =
470 & lambdaThetaClimRelax(i,j,bi,bj)
471 & + tmpfld2d(i,j,bi,bj)
472 enddo
473 enddo
474 enddo
475 enddo
476 #endif
477
478 #ifdef ALLOW_RELAXSSS_CONTROL
479 c-- SSS relaxation coefficient.
480 il=ilnblnk( xx_relaxsss_file )
481 write(fnamegeneric(1:80),'(2a,i10.10)')
482 & xx_relaxsss_file(1:il),'.',optimcycle
483 call active_read_xy_loc ( fnamegeneric, tmpfld2d, 1,
484 & doglobalread, ladinit, optimcycle,
485 & mythid, xx_relaxsss_dummy )
486 do bj = jtlo,jthi
487 do bi = itlo,ithi
488 do j = jmin,jmax
489 do i = imin,imax
490 lambdaSaltClimRelax(i,j,bi,bj) =
491 & lambdaSaltClimRelax(i,j,bi,bj)
492 & + tmpfld2d(i,j,bi,bj)
493 enddo
494 enddo
495 enddo
496 enddo
497 #endif
498
499 c-- Update the tile edges.
500
501 #if (defined (ALLOW_THETA0_CONTROL) || defined (ALLOW_SST0_CONTROL))
502 _EXCH_XYZ_R8( theta, mythid )
503 #endif
504 #if (defined (ALLOW_SALT0_CONTROL) || defined (ALLOW_SSS0_CONTROL))
505 _EXCH_XYZ_R8( salt, mythid )
506 #endif
507 #ifdef ALLOW_TR10_CONTROL
508 #ifdef ALLOW_PTRACERS
509 _EXCH_XYZ_R8(pTracer(1-Olx,1-Oly,1,1,1,1),myThid)
510 #endif
511 #endif
512 #ifdef ALLOW_DIFFKR_CONTROL
513 _EXCH_XYZ_R8( diffkr, mythid)
514 #endif
515 #ifdef ALLOW_KAPGM_CONTROL
516 _EXCH_XYZ_R8( kapgm, mythid)
517 #endif
518 #ifdef ALLOW_EFLUXY0_CONTROL
519 _EXCH_XYZ_R8( EfluxY, mythid )
520 #endif
521 #ifdef ALLOW_EFLUXP0_CONTROL
522 _EXCH_XYZ_R8( EfluxP, mythid )
523 #endif
524 #ifdef ALLOW_BOTTOMDRAG_CONTROL
525 _EXCH_XY_R8( bottomdragfld, mythid )
526 #endif
527
528 #if (defined (ALLOW_EDTAUX_CONTROL) && defined (ALLOW_EDTAUY_CONTROL))
529 CALL EXCH_UV_XYZ_RS(Eddytaux,Eddytauy,.TRUE.,myThid)
530 #elif (defined (ALLOW_EDTAUX_CONTROL) || defined (ALLOW_EDTAUY_CONTROL))
531 STOP 'ctrl_map_forcing: need BOTH ALLOW_EDTAU[X,Y]_CONTROL'
532 #endif
533
534 #ifdef ALLOW_UVEL0_CONTROL
535 _EXCH_XYZ_R8( uVel, mythid)
536 #endif
537
538 #ifdef ALLOW_VVEL0_CONTROL
539 _EXCH_XYZ_R8( vVel, mythid)
540 #endif
541
542 #ifdef ALLOW_ETAN0_CONTROL
543 _EXCH_XY_R8( etaN, mythid )
544 #endif
545
546 #ifdef ALLOW_RELAXSST_CONTROL
547 _EXCH_XY_R4( lambdaThetaClimRelax, mythid )
548 #endif
549
550 #ifdef ALLOW_RELAXSSS_CONTROL
551 _EXCH_XY_R4( lambdaThetaClimRelax, mythid )
552 #endif
553
554 return
555 end
556

  ViewVC Help
Powered by ViewVC 1.1.22