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

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

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


Revision 1.5 - (show annotations) (download)
Sat Jul 13 02:47:32 2002 UTC (21 years, 10 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint46n_post, checkpoint46l_post, checkpoint46g_pre, checkpoint46f_post, checkpoint46b_post, checkpoint46l_pre, checkpoint47a_post, checkpoint46d_pre, checkpoint46j_pre, checkpoint46a_post, checkpoint46j_post, checkpoint46k_post, checkpoint46e_pre, checkpoint46b_pre, checkpoint46c_pre, checkpoint46, checkpoint46h_pre, checkpoint46m_post, checkpoint46a_pre, checkpoint46g_post, checkpoint46i_post, checkpoint46c_post, checkpoint46e_post, checkpoint47, checkpoint46h_post, checkpoint46d_post
Changes since 1.4: +827 -124 lines
Merging new ctrl package from release1_p5:
o new ctrl package
  - adopted from ECCO environment to enable optimization
  - added Eliassen Palm fluxes to controls

1 C $Header: /u/gcmpack/MITgcm/pkg/ctrl/ctrl_init.F,v 1.3.4.2 2002/07/12 15:43:54 heimbach Exp $
2
3 #include "CTRL_CPPOPTIONS.h"
4
5
6 subroutine ctrl_init( mythid )
7
8 c ==================================================================
9 c SUBROUTINE ctrl_init
10 c ==================================================================
11 c
12 c o Set parts of the vector of control variables and initialize the
13 c rest to zero.
14 c
15 c The vector of control variables is initialized here. The
16 c temperature and salinity contributions are read from file.
17 c Subsequently, the latter are dimensionalized and the tile
18 c edges are updated.
19 c
20 c started: Christian Eckert eckert@mit.edu 30-Jun-1999
21 c
22 c changed: Christian Eckert eckert@mit.edu 23-Feb-2000
23 c - Restructured the code in order to create a package
24 c for the MITgcmUV.
25 c
26 c Patrick Heimbach heimbach@mit.edu 30-May-2000
27 c - diffsec was falsely declared.
28 c
29 c Patrick Heimbach heimbach@mit.edu 06-Jun-2000
30 c - Transferred some filename declarations
31 c from ctrl_pack/ctrl_unpack to here
32 c - Transferred mask-per-tile to here
33 c - computation of control vector length here
34 c
35 c Patrick Heimbach heimbach@mit.edu 16-Jun-2000
36 c - Added call to ctrl_pack
37 c - Alternatively: transfer writing of scale files to
38 c ctrl_unpack
39 c
40 c ==================================================================
41 c SUBROUTINE ctrl_init
42 c ==================================================================
43
44 implicit none
45
46 c == global variables ==
47
48 #include "EEPARAMS.h"
49 #include "SIZE.h"
50 #include "PARAMS.h"
51 #include "GRID.h"
52 #include "ctrl.h"
53
54 #ifdef ALLOW_CALENDAR
55 #include "cal.h"
56 #endif
57 #ifdef ALLOW_OBCS_CONTROL
58 # include "OBCS.h"
59 #endif
60 #ifdef ALLOW_ECCO_OPTIMIZATION
61 #include "optim.h"
62 #endif
63
64 c == routine arguments ==
65
66 integer mythid
67
68 c == local variables ==
69
70 integer bi,bj
71 integer i,j,k
72 integer itlo,ithi
73 integer jtlo,jthi
74 integer jmin,jmax
75 integer imin,imax
76 integer ntmp
77 integer ivarindex
78
79 integer iobcs
80 integer il
81 integer errio
82 integer startrec
83 integer endrec
84 integer difftime(4)
85 _RL diffsecs
86 _RL dummy
87
88 character*(80) ymaskobcs
89 character*(max_len_prec) record
90 character*(max_len_mbuf) msgbuf
91
92 integer nwetc3d
93
94 c == external ==
95
96 integer ilnblnk
97 external ilnblnk
98
99 c == end of interface ==
100
101 c-- Read the namelist input.
102 namelist /ctrl_nml/
103 & xx_theta_file,
104 & xx_salt_file,
105 & xx_hflux_file,
106 & xx_hfluxstartdate1, xx_hfluxstartdate2, xx_hfluxperiod,
107 & xx_sflux_file,
108 & xx_sfluxstartdate1, xx_sfluxstartdate2, xx_sfluxperiod,
109 & xx_tauu_file,
110 & xx_tauustartdate1, xx_tauustartdate2, xx_tauuperiod,
111 & xx_tauv_file,
112 & xx_tauvstartdate1, xx_tauvstartdate2, xx_tauvperiod,
113 & xx_atemp_file,
114 & xx_atempstartdate1, xx_atempstartdate2, xx_atempperiod,
115 & xx_aqh_file,
116 & xx_aqhstartdate1, xx_aqhstartdate2, xx_aqhperiod,
117 & xx_uwind_file,
118 & xx_uwindstartdate1, xx_uwindstartdate2, xx_uwindperiod,
119 & xx_vwind_file,
120 & xx_vwindstartdate1, xx_vwindstartdate2, xx_vwindperiod,
121 & xx_obcsn_file,
122 & xx_obcsnstartdate1, xx_obcsnstartdate2, xx_obcsnperiod,
123 & xx_obcss_file,
124 & xx_obcssstartdate1, xx_obcssstartdate2, xx_obcssperiod,
125 & xx_obcsw_file,
126 & xx_obcswstartdate1, xx_obcswstartdate2, xx_obcswperiod,
127 & xx_obcse_file,
128 & xx_obcsestartdate1, xx_obcsestartdate2, xx_obcseperiod,
129 & xx_diffkr_file,
130 & xx_kapgm_file,
131 & xx_tr1_file,
132 & xx_efluxy_file,
133 & xx_efluxp_file
134
135 namelist /ctrl_packnames/
136 & yadmark, yctrlid,
137 & ctrlname, costname, scalname, maskname, metaname
138
139 jtlo = mybylo(mythid)
140 jthi = mybyhi(mythid)
141 itlo = mybxlo(mythid)
142 ithi = mybxhi(mythid)
143 jmin = 1-oly
144 jmax = sny+oly
145 imin = 1-olx
146 imax = snx+olx
147
148 _BEGIN_MASTER( myThid )
149
150 c-- Set default values.
151 xx_theta_file = ' '
152 xx_salt_file = ' '
153 xx_hfluxstartdate1 = 0
154 xx_hfluxstartdate2 = 0
155 xx_hfluxperiod = 0. _d 0
156 xx_hflux_file = ' '
157 xx_sfluxstartdate1 = 0
158 xx_sfluxstartdate2 = 0
159 xx_sfluxperiod = 0. _d 0
160 xx_sflux_file = ' '
161 xx_tauustartdate1 = 0
162 xx_tauustartdate2 = 0
163 xx_tauuperiod = 0. _d 0
164 xx_tauu_file = ' '
165 xx_tauvstartdate1 = 0
166 xx_tauvstartdate2 = 0
167 xx_tauvperiod = 0. _d 0
168 xx_tauv_file = ' '
169 xx_atempstartdate1 = 0
170 xx_atempstartdate2 = 0
171 xx_atempperiod = 0. _d 0
172 xx_atemp_file = ' '
173 xx_aqhstartdate1 = 0
174 xx_aqhstartdate2 = 0
175 xx_aqhperiod = 0. _d 0
176 xx_aqh_file = ' '
177 xx_uwindstartdate1 = 0
178 xx_uwindstartdate2 = 0
179 xx_uwindperiod = 0. _d 0
180 xx_uwind_file = ' '
181 xx_vwindstartdate1 = 0
182 xx_vwindstartdate2 = 0
183 xx_vwindperiod = 0. _d 0
184 xx_vwind_file = ' '
185 xx_obcsnstartdate1 = 0
186 xx_obcsnstartdate2 = 0
187 xx_obcsnperiod = 0. _d 0
188 xx_obcsn_file = ' '
189 xx_obcssstartdate1 = 0
190 xx_obcssstartdate2 = 0
191 xx_obcssperiod = 0. _d 0
192 xx_obcss_file = ' '
193 xx_obcswstartdate1 = 0
194 xx_obcswstartdate2 = 0
195 xx_obcswperiod = 0. _d 0
196 xx_obcsw_file = ' '
197 xx_obcsestartdate1 = 0
198 xx_obcsestartdate2 = 0
199 xx_obcseperiod = 0. _d 0
200 xx_obcse_file = ' '
201 xx_diffkr_file = ' '
202 xx_kapgm_file = ' '
203 xx_tr1_file = ' '
204 xx_efluxy_file = ' '
205 xx_efluxp_file = ' '
206 yadmark = 'ad'
207 yctrlid = 'MIT_CE_000'
208 ctrlname = ' '
209 costname = ' '
210 scalname = ' '
211 maskname = ' '
212 metaname = ' '
213
214 c-- Check versions.
215
216 open(unit=scrunit1,status='scratch')
217
218 c-- Next, read the ecco data file.
219 open(unit = modeldataunit,file = 'data.ctrl',
220 & status = 'old', iostat = errio)
221 if ( errio .lt. 0 ) then
222 stop ' stopped in ctrl_init'
223 endif
224
225 do while ( .true. )
226 read(modeldataunit, fmt='(a)', end=1001) record
227 il = max(ilnblnk(record),1)
228 if ( record(1:1) .ne. commentcharacter )
229 & write(unit=scrunit1, fmt='(a)') record(:il)
230 enddo
231 1001 continue
232 close( modeldataunit )
233
234 rewind( scrunit1 )
235 read(unit = scrunit1, nml = ctrl_nml)
236 read(unit = scrunit1, nml = ctrl_packnames)
237 close( scrunit1 )
238
239 #ifdef ALLOW_CALENDAR
240
241 c-- Get the complete dates of the control variables.
242 #if (defined (ALLOW_HFLUX_CONTROL))
243 c-- The heat flux contribution.
244 call cal_FullDate( xx_hfluxstartdate1, xx_hfluxstartdate2,
245 & xx_hfluxstartdate , mythid )
246 #elif (defined (ALLOW_ATEMP_CONTROL))
247 c-- Atmos. temperature contribution.
248 call cal_FullDate( xx_atempstartdate1, xx_atempstartdate2,
249 & xx_atempstartdate , mythid )
250 #endif
251
252 #if (defined (ALLOW_SFLUX_CONTROL))
253 c-- The salt flux contribution.
254 call cal_FullDate( xx_sfluxstartdate1, xx_sfluxstartdate2,
255 & xx_sfluxstartdate , mythid )
256 #elif (defined (ALLOW_AQH_CONTROL))
257 c-- Atmospheric humidity contribution.
258 call cal_FullDate( xx_aqhstartdate1, xx_aqhstartdate2,
259 & xx_aqhstartdate , mythid )
260 #endif
261
262 #if (defined (ALLOW_USTRESS_CONTROL))
263 c-- The zonal wind stress contribution.
264 call cal_FullDate( xx_tauustartdate1, xx_tauustartdate2,
265 & xx_tauustartdate, mythid )
266 #elif (defined (ALLOW_UWIND_CONTROL))
267 c-- Zonal wind speed contribution.
268 call cal_FullDate( xx_uwindstartdate1, xx_uwindstartdate2,
269 & xx_uwindstartdate , mythid )
270 #endif
271
272 #if (defined (ALLOW_VSTRESS_CONTROL))
273 c-- The merid. wind stress contribution.
274 call cal_FullDate( xx_tauvstartdate1, xx_tauvstartdate2,
275 & xx_tauvstartdate, mythid )
276 #elif (defined (ALLOW_VWIND_CONTROL))
277 c-- Merid. wind speed contribution.
278 call cal_FullDate( xx_vwindstartdate1, xx_vwindstartdate2,
279 & xx_vwindstartdate , mythid )
280 #endif
281
282 #ifdef ALLOW_OBCS_CONTROL
283 call cal_FullDate( xx_obcsnstartdate1, xx_obcsnstartdate2,
284 & xx_obcsnstartdate, mythid )
285 call cal_FullDate( xx_obcssstartdate1, xx_obcssstartdate2,
286 & xx_obcssstartdate, mythid )
287 call cal_FullDate( xx_obcswstartdate1, xx_obcswstartdate2,
288 & xx_obcswstartdate, mythid )
289 call cal_FullDate( xx_obcsestartdate1, xx_obcsestartdate2,
290 & xx_obcsestartdate, mythid )
291 #endif
292
293 #endif /* ALLOW_CALENDAR */
294
295 c-- Set default values.
296 do ivarindex = 1,maxcvars
297 ncvarindex(ivarindex) = -1
298 ncvarrecs(ivarindex) = 0
299 ncvarxmax(ivarindex) = 0
300 ncvarymax(ivarindex) = 0
301 ncvarnrmax(ivarindex) = 0
302 ncvargrd(ivarindex) = '?'
303 enddo
304
305 write(msgbuf,'(a)') ' '
306 call print_message( msgbuf, standardmessageunit,
307 & SQUEEZE_RIGHT , mythid)
308 write(msgbuf,'(a)')
309 & ' ctrl_init: Initializing temperature and salinity'
310 call print_message( msgbuf, standardmessageunit,
311 & SQUEEZE_RIGHT , mythid)
312 write(msgbuf,'(a)')
313 & ' part of the control vector.'
314 call print_message( msgbuf, standardmessageunit,
315 & SQUEEZE_RIGHT , mythid)
316 write(msgbuf,'(a,a)')
317 & ' The initial surface fluxes are set',
318 & ' to zero.'
319 call print_message( msgbuf, standardmessageunit,
320 & SQUEEZE_RIGHT , mythid)
321 write(msgbuf,'(a)') ' '
322 call print_message( msgbuf, standardmessageunit,
323 & SQUEEZE_RIGHT , mythid)
324 _END_MASTER( mythid )
325
326 _BARRIER
327
328 c-- =====================
329 c-- Initial state fields.
330 c-- =====================
331
332 cph(
333 cph index 7-10 reserved for atmos. state,
334 cph index 11-14 reserved for open boundaries,
335 cph index 15-16 reserved for mixing coeff.
336 cph index 17 reserved for passive tracer TR1
337 cph index 18,19 reserved for sst, sss
338 cph index 20 for hFacC
339 cph index 21,22 for efluxy, efluxp
340 cph)
341
342 #ifdef ALLOW_THETA0_CONTROL
343 c-- Initial state temperature contribution.
344
345 _BEGIN_MASTER( mythid )
346 ivarindex = 1
347 ncvarindex(ivarindex) = 101
348 ncvarrecs(ivarindex) = 1
349 ncvarxmax(ivarindex) = snx
350 ncvarymax(ivarindex) = sny
351 ncvarnrmax(ivarindex) = nr
352 ncvargrd(ivarindex) = 'c'
353 _END_MASTER( mythid )
354
355 #endif /* ALLOW_THETA0_CONTROL */
356
357 #ifdef ALLOW_SALT0_CONTROL
358 c-- Initial state salinity contribution.
359
360 _BEGIN_MASTER( mythid )
361 ivarindex = 2
362 ncvarindex(ivarindex) = 102
363 ncvarrecs(ivarindex) = 1
364 ncvarxmax(ivarindex) = snx
365 ncvarymax(ivarindex) = sny
366 ncvarnrmax(ivarindex) = nr
367 ncvargrd(ivarindex) = 'c'
368 _END_MASTER( mythid )
369
370 #endif /* ALLOW_SALT0_CONTROL */
371
372 c-- ===========================
373 c-- Surface flux contributions.
374 c-- ===========================
375
376 #if (defined (ALLOW_HFLUX_CONTROL))
377 c-- Heat flux.
378
379 _BEGIN_MASTER( mythid )
380 #ifdef ALLOW_CALENDAR
381 call cal_TimePassed( xx_hfluxstartdate, modelstartdate,
382 & difftime, mythid )
383 call cal_ToSeconds ( difftime, diffsecs, mythid )
384 startrec = int((modelstart - diffsecs)/
385 & xx_hfluxperiod) + 1
386 endrec = int((modelend - diffsecs - modelstep)/
387 & xx_hfluxperiod) + 2
388 #else
389 startrec = 1
390 endrec = 1
391 #endif
392 ivarindex = 3
393 ncvarindex(ivarindex) = 103
394 ncvarrecs(ivarindex) = endrec - startrec + 1
395 ncvarrecstart(ivarindex) = startrec
396 ncvarrecsend(ivarindex) = endrec
397 ncvarxmax(ivarindex) = snx
398 ncvarymax(ivarindex) = sny
399 ncvarnrmax(ivarindex) = 1
400 ncvargrd(ivarindex) = 'c'
401 _END_MASTER( mythid )
402
403 #elif (defined (ALLOW_ATEMP_CONTROL))
404 c-- Atmos. temperature
405
406 _BEGIN_MASTER( mythid )
407 #ifdef ALLOW_CALENDAR
408 call cal_TimePassed( xx_atempstartdate, modelstartdate,
409 & difftime, mythid )
410 call cal_ToSeconds ( difftime, diffsecs, mythid )
411 startrec = int((modelstart - diffsecs)/
412 & xx_atempperiod) + 1
413 endrec = int((modelend - diffsecs - modelstep)/
414 & xx_atempperiod) + 2
415 #else
416 startrec = 1
417 endrec = 1
418 #endif
419 ivarindex = 7
420 ncvarindex(ivarindex) = 107
421 ncvarrecs(ivarindex) = endrec - startrec + 1
422 ncvarrecstart(ivarindex) = startrec
423 ncvarrecsend(ivarindex) = endrec
424 ncvarxmax(ivarindex) = snx
425 ncvarymax(ivarindex) = sny
426 ncvarnrmax(ivarindex) = 1
427 ncvargrd(ivarindex) = 'c'
428 _END_MASTER( mythid )
429
430 #elif (defined (ALLOW_HFLUX0_CONTROL))
431 c-- initial forcing only
432 _BEGIN_MASTER( mythid )
433 ncvarindex(3) = 103
434 ncvarrecs(3) = 1
435 ncvarxmax(3) = snx
436 ncvarymax(3) = sny
437 ncvarnrmax(3) = 1
438 ncvargrd(3) = 'c'
439 _END_MASTER( mythid )
440
441 #endif /* ALLOW_HFLUX_CONTROL */
442
443 #if (defined (ALLOW_SFLUX_CONTROL))
444 c-- Salt flux.
445
446 _BEGIN_MASTER( mythid )
447 #ifdef ALLOW_CALENDAR
448 call cal_TimePassed( xx_sfluxstartdate, modelstartdate,
449 & difftime, mythid )
450 call cal_ToSeconds ( difftime, diffsecs, mythid )
451 startrec = int((modelstart - diffsecs)/
452 & xx_sfluxperiod) + 1
453 endrec = int((modelend - diffsecs - modelstep)/
454 & xx_sfluxperiod) + 2
455 #else
456 startrec = 1
457 endrec = 1
458 #endif
459 ivarindex = 4
460 ncvarindex(ivarindex) = 104
461 ncvarrecs(ivarindex) = endrec - startrec + 1
462 ncvarrecstart(ivarindex) = startrec
463 ncvarrecsend(ivarindex) = endrec
464 ncvarxmax(ivarindex) = snx
465 ncvarymax(ivarindex) = sny
466 ncvarnrmax(ivarindex) = 1
467 ncvargrd(ivarindex) = 'c'
468 _END_MASTER( mythid )
469
470 #elif (defined (ALLOW_AQH_CONTROL))
471 c-- Atmos. humidity
472
473 _BEGIN_MASTER( mythid )
474 #ifdef ALLOW_CALENDAR
475 call cal_TimePassed( xx_aqhstartdate, modelstartdate,
476 & difftime, mythid )
477 call cal_ToSeconds ( difftime, diffsecs, mythid )
478 startrec = int((modelstart - diffsecs)/
479 & xx_aqhperiod) + 1
480 endrec = int((modelend - diffsecs - modelstep)/
481 & xx_aqhperiod) + 2
482 #else
483 startrec = 1
484 endrec = 1
485 #endif
486 ivarindex = 8
487 ncvarindex(ivarindex) = 108
488 ncvarrecs(ivarindex) = endrec - startrec + 1
489 ncvarrecstart(ivarindex) = startrec
490 ncvarrecsend(ivarindex) = endrec
491 ncvarxmax(ivarindex) = snx
492 ncvarymax(ivarindex) = sny
493 ncvarnrmax(ivarindex) = 1
494 ncvargrd(ivarindex) = 'c'
495 _END_MASTER( mythid )
496
497 #elif (defined (ALLOW_SFLUX0_CONTROL))
498 c-- initial forcing only
499 _BEGIN_MASTER( mythid )
500 ncvarindex(4) = 104
501 ncvarrecs(4) = 1
502 ncvarxmax(4) = snx
503 ncvarymax(4) = sny
504 ncvarnrmax(4) = 1
505 ncvargrd(4) = 'c'
506 _END_MASTER( mythid )
507
508 #endif /* ALLOW_SFLUX_CONTROL */
509
510 #if (defined (ALLOW_USTRESS_CONTROL))
511 c-- Zonal wind stress.
512
513 _BEGIN_MASTER( mythid )
514 #ifdef ALLOW_CALENDAR
515 call cal_TimePassed( xx_tauustartdate, modelstartdate,
516 & difftime, mythid )
517 call cal_ToSeconds ( difftime, diffsecs, mythid )
518 startrec = int((modelstart - diffsecs)/
519 & xx_tauuperiod) + 1
520 endrec = int((modelend - diffsecs - modelstep)/
521 & xx_tauuperiod) + 2
522 #else
523 startrec = 1
524 endrec = 1
525 #endif
526 ivarindex = 5
527 ncvarindex(ivarindex) = 105
528 ncvarrecs(ivarindex) = endrec - startrec + 1
529 ncvarrecstart(ivarindex) = startrec
530 ncvarrecsend(ivarindex) = endrec
531 ncvarxmax(ivarindex) = snx
532 ncvarymax(ivarindex) = sny
533 ncvarnrmax(ivarindex) = 1
534 ncvargrd(ivarindex) = 'w'
535 _END_MASTER( mythid )
536
537 #elif (defined (ALLOW_UWIND_CONTROL))
538 c-- Zonal wind speed.
539
540 _BEGIN_MASTER( mythid )
541 #ifdef ALLOW_CALENDAR
542 call cal_TimePassed( xx_uwindstartdate, modelstartdate,
543 & difftime, mythid )
544 call cal_ToSeconds ( difftime, diffsecs, mythid )
545 startrec = int((modelstart - diffsecs)/
546 & xx_uwindperiod) + 1
547 endrec = int((modelend - diffsecs - modelstep)/
548 & xx_uwindperiod) + 2
549 #else
550 startrec = 1
551 endrec = 1
552 #endif
553 ivarindex = 9
554 ncvarindex(ivarindex) = 109
555 ncvarrecs(ivarindex) = endrec - startrec + 1
556 ncvarrecstart(ivarindex) = startrec
557 ncvarrecsend(ivarindex) = endrec
558 ncvarxmax(ivarindex) = snx
559 ncvarymax(ivarindex) = sny
560 ncvarnrmax(ivarindex) = 1
561 ncvargrd(ivarindex) = 'w'
562 _END_MASTER( mythid )
563
564 #elif (defined (ALLOW_TAUU0_CONTROL))
565 c-- initial forcing only
566 _BEGIN_MASTER( mythid )
567 ncvarindex(5) = 105
568 ncvarrecs(5) = 1
569 ncvarxmax(5) = snx
570 ncvarymax(5) = sny
571 ncvarnrmax(5) = 1
572 ncvargrd(5) = 'w'
573 _END_MASTER( mythid )
574
575 #endif /* ALLOW_USTRESS_CONTROL */
576
577 #if (defined (ALLOW_VSTRESS_CONTROL))
578 c-- Meridional wind stress.
579
580 _BEGIN_MASTER( mythid )
581 #ifdef ALLOW_CALENDAR
582 call cal_TimePassed( xx_tauvstartdate, modelstartdate,
583 & difftime, mythid )
584 call cal_ToSeconds ( difftime, diffsecs, mythid )
585 startrec = int((modelstart - diffsecs)/
586 & xx_tauvperiod) + 1
587 endrec = int((modelend - diffsecs - modelstep)/
588 & xx_tauvperiod) + 2
589 #else
590 startrec = 1
591 endrec = 1
592 #endif
593 ivarindex = 6
594 ncvarindex(ivarindex) = 106
595 ncvarrecs(ivarindex) = endrec - startrec + 1
596 ncvarrecstart(ivarindex) = startrec
597 ncvarrecsend(ivarindex) = endrec
598 ncvarxmax(ivarindex) = snx
599 ncvarymax(ivarindex) = sny
600 ncvarnrmax(ivarindex) = 1
601 ncvargrd(ivarindex) = 's'
602 _END_MASTER( mythid )
603
604 #elif (defined (ALLOW_VWIND_CONTROL))
605 c-- Meridional wind speed.
606
607 _BEGIN_MASTER( mythid )
608 #ifdef ALLOW_CALENDAR
609 call cal_TimePassed( xx_vwindstartdate, modelstartdate,
610 & difftime, mythid )
611 call cal_ToSeconds ( difftime, diffsecs, mythid )
612 startrec = int((modelstart - diffsecs)/
613 & xx_vwindperiod) + 1
614 endrec = int((modelend - diffsecs - modelstep)/
615 & xx_vwindperiod) + 2
616 #else
617 startrec = 1
618 endrec = 1
619 #endif
620 ivarindex = 10
621 ncvarindex(ivarindex) = 110
622 ncvarrecs(ivarindex) = endrec - startrec + 1
623 ncvarrecstart(ivarindex) = startrec
624 ncvarrecsend(ivarindex) = endrec
625 ncvarxmax(ivarindex) = snx
626 ncvarymax(ivarindex) = sny
627 ncvarnrmax(ivarindex) = 1
628 ncvargrd(ivarindex) = 's'
629 _END_MASTER( mythid )
630
631 #elif (defined (ALLOW_TAUV0_CONTROL))
632 c-- initial forcing only
633 _BEGIN_MASTER( mythid )
634 ncvarindex(6) = 106
635 ncvarrecs(6) = 1
636 ncvarxmax(6) = snx
637 ncvarymax(6) = sny
638 ncvarnrmax(6) = 1
639 ncvargrd(6) = 's'
640 _END_MASTER( mythid )
641
642 #endif /* ALLOW_VSTRESS_CONTROL */
643
644 #ifdef ALLOW_OBCSN_CONTROL
645 c-- Northern obc.
646
647 _BEGIN_MASTER( mythid )
648 #ifdef ALLOW_CALENDAR
649 call cal_TimePassed( xx_obcsnstartdate, modelstartdate,
650 & difftime, mythid )
651 call cal_ToSeconds ( difftime, diffsecs, mythid )
652 startrec = int((modelstart - diffsecs)/
653 & xx_obcsnperiod) + 1
654 endrec = int((modelend - diffsecs)/
655 & xx_obcsnperiod) + 2
656 #else
657 startrec = 1
658 endrec = 1
659 #endif
660 ivarindex = 11
661 ncvarindex(ivarindex) = 111
662 ncvarrecs(ivarindex) = (endrec - startrec + 1)*nobcs
663 ncvarrecstart(ivarindex) = (startrec - 1)*nobcs + 1
664 ncvarrecsend(ivarindex) = endrec*nobcs
665 ncvarxmax(ivarindex) = snx
666 ncvarymax(ivarindex) = 1
667 ncvarnrmax(ivarindex) = nr
668 ncvargrd(ivarindex) = 'm'
669 _END_MASTER( mythid )
670
671 #endif /* ALLOW_OBCSN_CONTROL */
672
673 #ifdef ALLOW_OBCSS_CONTROL
674 c-- Southern obc.
675
676 _BEGIN_MASTER( mythid )
677 #ifdef ALLOW_CALENDAR
678 call cal_TimePassed( xx_obcssstartdate, modelstartdate,
679 & difftime, mythid )
680 call cal_ToSeconds ( difftime, diffsecs, mythid )
681 startrec = int((modelstart - diffsecs)/
682 & xx_obcssperiod) + 1
683 endrec = int((modelend - diffsecs)/
684 & xx_obcssperiod) + 2
685 #else
686 startrec = 1
687 endrec = 1
688 #endif
689 ivarindex = 12
690 ncvarindex(ivarindex) = 112
691 ncvarrecs(ivarindex) = (endrec - startrec + 1)*nobcs
692 ncvarrecstart(ivarindex) = (startrec - 1)*nobcs + 1
693 ncvarrecsend(ivarindex) = endrec*nobcs
694 ncvarxmax(ivarindex) = snx
695 ncvarymax(ivarindex) = 1
696 ncvarnrmax(ivarindex) = nr
697 ncvargrd(ivarindex) = 'm'
698 _END_MASTER( mythid )
699
700 #endif /* ALLOW_OBCSS_CONTROL */
701
702 #ifdef ALLOW_OBCSW_CONTROL
703 c-- Western obc.
704
705 _BEGIN_MASTER( mythid )
706 #ifdef ALLOW_CALENDAR
707 call cal_TimePassed( xx_obcswstartdate, modelstartdate,
708 & difftime, mythid )
709 call cal_ToSeconds ( difftime, diffsecs, mythid )
710 startrec = int((modelstart - diffsecs)/
711 & xx_obcswperiod) + 1
712 endrec = int((modelend - diffsecs)/
713 & xx_obcswperiod) + 2
714 #else
715 startrec = 1
716 endrec = 1
717 #endif
718 ivarindex = 13
719 ncvarindex(ivarindex) = 113
720 ncvarrecs(ivarindex) = (endrec - startrec + 1)*nobcs
721 ncvarrecstart(ivarindex) = (startrec - 1)*nobcs + 1
722 ncvarrecsend(ivarindex) = endrec*nobcs
723 ncvarxmax(ivarindex) = 1
724 ncvarymax(ivarindex) = sny
725 ncvarnrmax(ivarindex) = nr
726 ncvargrd(ivarindex) = 'm'
727 _END_MASTER( mythid )
728
729 #endif /* ALLOW_OBCSW_CONTROL */
730
731 #ifdef ALLOW_OBCSE_CONTROL
732 c-- Eastern obc.
733
734 _BEGIN_MASTER( mythid )
735 #ifdef ALLOW_CALENDAR
736 call cal_TimePassed( xx_obcsestartdate, modelstartdate,
737 & difftime, mythid )
738 call cal_ToSeconds ( difftime, diffsecs, mythid )
739 startrec = int((modelstart - diffsecs)/
740 & xx_obcseperiod) + 1
741 endrec = int((modelend - diffsecs)/
742 & xx_obcseperiod) + 2
743 #else
744 startrec = 1
745 endrec = 1
746 #endif
747 ivarindex = 14
748 ncvarindex(ivarindex) = 114
749 ncvarrecs(ivarindex) = (endrec - startrec + 1)*nobcs
750 ncvarrecstart(ivarindex) = (startrec - 1)*nobcs + 1
751 ncvarrecsend(ivarindex) = endrec*nobcs
752 ncvarxmax(ivarindex) = 1
753 ncvarymax(ivarindex) = sny
754 ncvarnrmax(ivarindex) = nr
755 ncvargrd(ivarindex) = 'm'
756 _END_MASTER( mythid )
757
758 #endif /* ALLOW_OBCSE_CONTROL */
759
760 #ifdef ALLOW_DIFFKR_CONTROL
761 _BEGIN_MASTER( mythid )
762 ivarindex = 15
763 ncvarindex(ivarindex) = 115
764 ncvarrecs (ivarindex) = 1
765 ncvarxmax (ivarindex) = snx
766 ncvarymax (ivarindex) = sny
767 ncvarnrmax(ivarindex) = nr
768 ncvargrd (ivarindex) = 'c'
769 _END_MASTER( mythid )
770 #endif /* ALLOW_DIFFKR_CONTROL */
771
772 #ifdef ALLOW_KAPGM_CONTROL
773 _BEGIN_MASTER( mythid )
774 ivarindex = 16
775 ncvarindex(ivarindex) = 116
776 ncvarrecs (ivarindex) = 1
777 ncvarxmax (ivarindex) = snx
778 ncvarymax (ivarindex) = sny
779 ncvarnrmax(ivarindex) = nr
780 ncvargrd (ivarindex) = 'c'
781 _END_MASTER( mythid )
782 #endif /* ALLOW_KAPGM_CONTROL */
783
784 #ifdef ALLOW_TR10_CONTROL
785 _BEGIN_MASTER( mythid )
786 ivarindex = 17
787 ncvarindex(ivarindex) = 117
788 ncvarrecs (ivarindex) = 1
789 ncvarxmax (ivarindex) = snx
790 ncvarymax (ivarindex) = sny
791 ncvarnrmax(ivarindex) = nr
792 ncvargrd (ivarindex) = 'c'
793 _END_MASTER( mythid )
794 #endif /* ALLOW_TR10_CONTROL */
795
796 #ifdef ALLOW_EFLUXY0_CONTROL
797 _BEGIN_MASTER( mythid )
798 ivarindex = 21
799 ncvarindex(ivarindex) = 121
800 ncvarrecs(ivarindex) = 1
801 ncvarxmax(ivarindex) = snx
802 ncvarymax(ivarindex) = sny
803 ncvarnrmax(ivarindex) = nr
804 ncvargrd(ivarindex) = 's'
805 _END_MASTER( mythid )
806 #endif /* ALLOW_EFLUXY0_CONTROL */
807
808 #ifdef ALLOW_EFLUXP0_CONTROL
809 _BEGIN_MASTER( mythid )
810 ivarindex = 22
811 ncvarindex(ivarindex) = 122
812 ncvarrecs(ivarindex) = 1
813 ncvarxmax(ivarindex) = snx
814 ncvarymax(ivarindex) = sny
815 ncvarnrmax(ivarindex) = nr
816 ncvargrd(ivarindex) = 'v'
817 _END_MASTER( mythid )
818 #endif /* ALLOW_EFLUXP0_CONTROL */
819
820 c-- Determine the number of wet points in each tile:
821 c-- maskc, masks, and maskw.
822
823 c-- Set loop ranges.
824 jmin = 1
825 jmax = sny
826 imin = 1
827 imax = snx
828
829 c-- Initialise the counters.
830 do bj = jtlo,jthi
831 do bi = itlo,ithi
832 do k = 1,nr
833 nwetctile(bi,bj,k) = 0
834 nwetstile(bi,bj,k) = 0
835 nwetwtile(bi,bj,k) = 0
836 nwetvtile(bi,bj,k) = 0
837 enddo
838 enddo
839 enddo
840
841 #ifdef ALLOW_OBCS_CONTROL
842 c-- Initialise obcs counters.
843 do bj = jtlo,jthi
844 do bi = itlo,ithi
845 do k = 1,nr
846 do iobcs = 1,nobcs
847 #ifdef ALLOW_OBCSN_CONTROL
848 nwetobcsn(bi,bj,k,iobcs) = 0
849 #endif
850 #ifdef ALLOW_OBCSS_CONTROL
851 nwetobcss(bi,bj,k,iobcs) = 0
852 #endif
853 #ifdef ALLOW_OBCSW_CONTROL
854 nwetobcsw(bi,bj,k,iobcs) = 0
855 #endif
856 #ifdef ALLOW_OBCSE_CONTROL
857 nwetobcse(bi,bj,k,iobcs) = 0
858 #endif
859 enddo
860 enddo
861 enddo
862 enddo
863 #endif
864
865 c-- Count wet points on each tile.
866 do bj = jtlo,jthi
867 do bi = itlo,ithi
868 do k = 1,nr
869 do j = jmin,jmax
870 do i = imin,imax
871 c-- Center mask.
872 if (hFacC(i,j,k,bi,bj) .ne. 0.) then
873 nwetctile(bi,bj,k) = nwetctile(bi,bj,k) + 1
874 endif
875 c-- South mask.
876 if (maskS(i,j,k,bi,bj) .eq. 1.) then
877 nwetstile(bi,bj,k) = nwetstile(bi,bj,k) + 1
878 endif
879 c-- West mask.
880 if (maskW(i,j,k,bi,bj) .eq. 1.) then
881 nwetwtile(bi,bj,k) = nwetwtile(bi,bj,k) + 1
882 endif
883 #if (defined (ALLOW_EFLUXP0_CONTROL))
884 c-- Vertical mask.
885 if (hFacV(i,j,k,bi,bj) .ne. 0.) then
886 nwetvtile(bi,bj,k) = nwetvtile(bi,bj,k) + 1
887 endif
888 #endif
889 enddo
890 enddo
891 enddo
892 enddo
893 enddo
894
895 #ifdef ALLOW_OBCSN_CONTROL
896 c-- Count wet points at Northern boundary.
897 c-- mask conventions are adopted from obcs_apply_ts, obcs_apply_uv
898 ymaskobcs = 'maskobcsn'
899 call ctrl_mask_set_xz(
900 & 0, OB_Jn, nwetobcsn, ymaskobcs, mythid
901 & )
902 #endif
903
904 #ifdef ALLOW_OBCSS_CONTROL
905 c-- Count wet points at Northern boundary.
906 c-- mask conventions are adopted from obcs_apply_ts, obcs_apply_uv
907 ymaskobcs = 'maskobcss'
908 call ctrl_mask_set_xz(
909 & 1, OB_Js, nwetobcss, ymaskobcs, mythid
910 & )
911 #endif
912
913 #ifdef ALLOW_OBCSW_CONTROL
914 c-- Count wet points at Northern boundary.
915 c-- mask conventions are adopted from obcs_apply_ts, obcs_apply_uv
916 ymaskobcs = 'maskobcsw'
917 call ctrl_mask_set_yz(
918 & 1, OB_Iw, nwetobcsw, ymaskobcs, mythid
919 & )
920 #endif
921
922 #ifdef ALLOW_OBCSE_CONTROL
923 c-- Count wet points at Northern boundary.
924 c-- mask conventions are adopted from obcs_apply_ts, obcs_apply_uv
925 ymaskobcs = 'maskobcse'
926 call ctrl_mask_set_yz(
927 & 0, OB_Ie, nwetobcse, ymaskobcs, mythid
928 & )
929 #endif
930
931 _BEGIN_MASTER( mythid )
932 c-- Determine the total number of control variables.
933 nvartype = 0
934 nvarlength = 0
935 do i = 1,maxcvars
936 c
937 if ( ncvarindex(i) .ne. -1 ) then
938 nvartype = nvartype + 1
939 do bj = jtlo,jthi
940 do bi = itlo,ithi
941 do k = 1,ncvarnrmax(i)
942 if ( ncvargrd(i) .eq. 'c' ) then
943 nvarlength = nvarlength +
944 & ncvarrecs(i)*nwetctile(bi,bj,k)
945 else if ( ncvargrd(i) .eq. 's' ) then
946 nvarlength = nvarlength +
947 & ncvarrecs(i)*nwetstile(bi,bj,k)
948 else if ( ncvargrd(i) .eq. 'w' ) then
949 nvarlength = nvarlength +
950 & ncvarrecs(i)*nwetwtile(bi,bj,k)
951 else if ( ncvargrd(i) .eq. 'v' ) then
952 nvarlength = nvarlength +
953 & ncvarrecs(i)*nwetvtile(bi,bj,k)
954 else if ( ncvargrd(i) .eq. 'm' ) then
955 #ifdef ALLOW_OBCS_CONTROL
956 do iobcs = 1, nobcs
957 if ( i .eq. 11 ) then
958 #ifdef ALLOW_OBCSN_CONTROL
959 nvarlength = nvarlength +
960 & (ncvarrecs(i)/nobcs)
961 & *nwetobcsn(bi,bj,k,iobcs)
962 #endif
963 else if ( i .eq. 12 ) then
964 #ifdef ALLOW_OBCSS_CONTROL
965 nvarlength = nvarlength +
966 & (ncvarrecs(i)/nobcs)
967 & *nwetobcss(bi,bj,k,iobcs)
968 #endif
969 else if ( i .eq. 13 ) then
970 #ifdef ALLOW_OBCSW_CONTROL
971 nvarlength = nvarlength +
972 & (ncvarrecs(i)/nobcs)
973 & *nwetobcsw(bi,bj,k,iobcs)
974 #endif
975 else if ( i .eq. 14 ) then
976 #ifdef ALLOW_OBCSE_CONTROL
977 nvarlength = nvarlength +
978 & (ncvarrecs(i)/nobcs)
979 & *nwetobcse(bi,bj,k,iobcs)
980 #endif
981 end if
982 enddo
983 #endif
984 else
985 print*,'ctrl_init: invalid grid location'
986 print*,' control variable = ',ncvarindex(i)
987 print*,' grid location = ',ncvargrd(i)
988 stop ' ... stopped in ctrl_init'
989 endif
990 enddo
991 enddo
992 enddo
993 endif
994 enddo
995
996 cph(
997 print *, 'ph-wet 1: nvarlength = ', nvarlength
998 print *, 'ph-wet 2: surface wet C = ', nwetctile(1,1,1)
999 print *, 'ph-wet 3: surface wet W = ', nwetwtile(1,1,1)
1000 print *, 'ph-wet 4: surface wet S = ', nwetstile(1,1,1)
1001 print *, 'ph-wet 4a:surface wet V = ', nwetvtile(1,1,1)
1002 nwetc3d = 0
1003 do k = 1, Nr
1004 nwetc3d = nwetc3d + nwetctile(1,1,k)
1005 end do
1006 print *, 'ph-wet 5: 3D wet points = ', nwetc3d
1007 do i = 1, maxcvars
1008 print *, 'ph-wet 6: no recs for i = ', i, ncvarrecs(i)
1009 end do
1010 print *, 'ph-wet 7: ',
1011 & 2*nwetc3d +
1012 & ncvarrecs(3)*nwetctile(1,1,1) +
1013 & ncvarrecs(4)*nwetctile(1,1,1) +
1014 & ncvarrecs(5)*nwetwtile(1,1,1) +
1015 & ncvarrecs(6)*nwetstile(1,1,1)
1016 print *, 'ph-wet 8: ',
1017 & 2*nwetc3d +
1018 & ncvarrecs(7)*nwetctile(1,1,1) +
1019 & ncvarrecs(8)*nwetctile(1,1,1) +
1020 & ncvarrecs(9)*nwetwtile(1,1,1) +
1021 & ncvarrecs(10)*nwetstile(1,1,1)
1022 #ifdef ALLOW_OBCSN_CONTROL
1023 print *, 'ph-wet 9: surface wet obcsn = '
1024 & , nwetobcsn(1,1,1,1), nwetobcsn(1,1,1,2)
1025 & , nwetobcsn(1,1,1,3), nwetobcsn(1,1,1,4)
1026 #endif
1027 #ifdef ALLOW_OBCSS_CONTROL
1028 print *, 'ph-wet 10: surface wet obcss = '
1029 & , nwetobcss(1,1,1,1), nwetobcss(1,1,1,2)
1030 & , nwetobcss(1,1,1,3), nwetobcss(1,1,1,4)
1031 #endif
1032 #ifdef ALLOW_OBCSW_CONTROL
1033 print *, 'ph-wet 11: surface wet obcsw = '
1034 & , nwetobcsw(1,1,1,1), nwetobcsw(1,1,1,2)
1035 & , nwetobcsw(1,1,1,3), nwetobcsw(1,1,1,4)
1036 #endif
1037 #ifdef ALLOW_OBCSE_CONTROL
1038 print *, 'ph-wet 12: surface wet obcse = '
1039 & , nwetobcse(1,1,1,1), nwetobcse(1,1,1,2)
1040 & , nwetobcse(1,1,1,3), nwetobcse(1,1,1,4)
1041 #endif
1042 cph)
1043
1044 CALL GLOBAL_SUM_INT( nvarlength, myThid )
1045
1046 print *, 'ph-wet 13: global nvarlength vor k=', k, nvarlength
1047
1048 c
1049 c Summation of wet point counters
1050 c
1051 do k = 1, nr
1052
1053 ntmp=0
1054 do bj=1,nSy
1055 do bi=1,nSx
1056 ntmp=ntmp+nWetcTile(bi,bj,k)
1057 enddo
1058 enddo
1059 CALL GLOBAL_SUM_INT( ntmp, myThid )
1060 nWetcGlobal(k)=ntmp
1061
1062 print *, 'ph-wet 14a: global nWet... vor k=', k, ntmp
1063
1064 ntmp=0
1065 do bj=1,nSy
1066 do bi=1,nSx
1067 ntmp=ntmp+nWetsTile(bi,bj,k)
1068 enddo
1069 enddo
1070 CALL GLOBAL_SUM_INT( ntmp, myThid )
1071 nWetsGlobal(k)=ntmp
1072
1073 print *, 'ph-wet 14b: global nWet... vor k=', k, ntmp
1074
1075 ntmp=0
1076 do bj=1,nSy
1077 do bi=1,nSx
1078 ntmp=ntmp+nWetwTile(bi,bj,k)
1079 enddo
1080 enddo
1081 CALL GLOBAL_SUM_INT( ntmp, myThid )
1082 nWetwGlobal(k)=ntmp
1083
1084 print *, 'ph-wet 14c: global nWet... vor k=', k, ntmp
1085
1086 ntmp=0
1087 do bj=1,nSy
1088 do bi=1,nSx
1089 ntmp=ntmp+nWetvTile(bi,bj,k)
1090 enddo
1091 enddo
1092 CALL GLOBAL_SUM_INT( ntmp, myThid )
1093 nWetvGlobal(k)=ntmp
1094
1095 print *, 'ph-wet 14d: global nWet... vor k=', k, ntmp
1096
1097 #ifdef ALLOW_OBCSN_CONTROL
1098 do iobcs = 1, nobcs
1099 ntmp=0
1100 do bj=1,nSy
1101 do bi=1,nSx
1102 ntmp=ntmp+nwetobcsn(bi,bj,k,iobcs)
1103 enddo
1104 enddo
1105 CALL GLOBAL_SUM_INT( ntmp, myThid )
1106 nwetobcsnglo(k,iobcs)=ntmp
1107 enddo
1108 #endif
1109 #ifdef ALLOW_OBCSS_CONTROL
1110 do iobcs = 1, nobcs
1111 ntmp=0
1112 do bj=1,nSy
1113 do bi=1,nSx
1114 ntmp=ntmp+nwetobcss(bi,bj,k,iobcs)
1115 enddo
1116 enddo
1117 CALL GLOBAL_SUM_INT( ntmp, myThid )
1118 nwetobcssglo(k,iobcs)=ntmp
1119 enddo
1120 #endif
1121 #ifdef ALLOW_OBCSW_CONTROL
1122 do iobcs = 1, nobcs
1123 ntmp=0
1124 do bj=1,nSy
1125 do bi=1,nSx
1126 ntmp=ntmp+nwetobcsw(bi,bj,k,iobcs)
1127 enddo
1128 enddo
1129 CALL GLOBAL_SUM_INT( ntmp, myThid )
1130 nwetobcswglo(k,iobcs)=ntmp
1131 enddo
1132 #endif
1133 #ifdef ALLOW_OBCSE_CONTROL
1134 do iobcs = 1, nobcs
1135 ntmp=0
1136 do bj=1,nSy
1137 do bi=1,nSx
1138 ntmp=ntmp+nwetobcse(bi,bj,k,iobcs)
1139 enddo
1140 enddo
1141 CALL GLOBAL_SUM_INT( ntmp, myThid )
1142 nwetobcseglo(k,iobcs)=ntmp
1143 enddo
1144 #endif
1145
1146 enddo
1147
1148 print*, 'ctrl_init: no. of control variables: ', nvartype
1149 print*, 'ctrl_init: control vector length: ', nvarlength
1150 _END_MASTER( mythid )
1151
1152 c write masks and weights to files to be read by a master process
1153 c
1154 call active_write_xyz( 'hFacC', hFacC, 1, 0, mythid, dummy)
1155 call active_write_xyz( 'maskW', maskW, 1, 0, mythid, dummy)
1156 call active_write_xyz( 'maskS', maskS, 1, 0, mythid, dummy)
1157 #if (defined (ALLOW_EFLUXP0_CONTROL))
1158 call active_write_xyz( 'hFacV', hFacV, 1, 0, mythid, dummy)
1159 #endif
1160
1161 c-- Summarize the control vector's setup.
1162 _BEGIN_MASTER( mythid )
1163 cph call ctrl_Summary( mythid )
1164 _END_MASTER( mythid )
1165
1166 _BARRIER
1167
1168 return
1169 end
1170

  ViewVC Help
Powered by ViewVC 1.1.22