/[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.28 - (show annotations) (download)
Sat Feb 2 02:34:49 2008 UTC (17 years, 5 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint59o, checkpoint59n
Changes since 1.27: +26 -1 lines
introduce isopycnal diffusion coefficient control.

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

  ViewVC Help
Powered by ViewVC 1.1.22