C $Header: /home/ubuntu/mnt/e9_copy/MITgcm/pkg/ctrl/ctrl_init.F,v 1.41 2011/07/19 12:44:11 mlosch Exp $ C $Name: $ #include "CTRL_CPPOPTIONS.h" subroutine ctrl_init( mythid ) c ================================================================== c SUBROUTINE ctrl_init c ================================================================== c c o Set parts of the vector of control variables and initialize the c rest to zero. c c The vector of control variables is initialized here. The c temperature and salinity contributions are read from file. c Subsequently, the latter are dimensionalized and the tile c edges are updated. c c started: Christian Eckert eckert@mit.edu 30-Jun-1999 c c changed: Christian Eckert eckert@mit.edu 23-Feb-2000 c - Restructured the code in order to create a package c for the MITgcmUV. c c Patrick Heimbach heimbach@mit.edu 30-May-2000 c - diffsec was falsely declared. c c Patrick Heimbach heimbach@mit.edu 06-Jun-2000 c - Transferred some filename declarations c from ctrl_pack/ctrl_unpack to here c - Transferred mask-per-tile to here c - computation of control vector length here c c Patrick Heimbach heimbach@mit.edu 16-Jun-2000 c - Added call to ctrl_pack c - Alternatively: transfer writing of scale files to c ctrl_unpack c c Dimitris Menemenlis menemenlis@mit.edu 7-Mar-2003 c - To be consistent with usage in ctrl_getrec.F, c startrec and endrec need to be referenced to c model time = 0, not to startTime. c Also "- modelstep" -> "+ modelstep/2": c old: startrec = int((modelstart - diffsecs)/ c old: & xx_???period) + 1 c old: endrec = int((modelend - diffsecs - modelstep)/ c old: & xx_???period) + 2 c new: startrec = int((modelstart + startTime - diffsecs)/ c new: & xx_???period) + 1 c new: endrec = int((modelend + startTime - diffsecs + modelstep/2)/ c new: & xx_???period) + 2 c c heimbach@mit.edu totally restructured 28-Oct-2003 c c ================================================================== c SUBROUTINE ctrl_init c ================================================================== implicit none c == global variables == #include "EEPARAMS.h" #include "SIZE.h" #include "PARAMS.h" #include "GRID.h" #include "ctrl.h" #include "optim.h" #ifdef ALLOW_CAL # include "cal.h" #endif #ifdef ALLOW_DIC_CONTROL # include "DIC_CTRL.h" #endif c == routine arguments == integer mythid c == local variables == integer bi,bj integer i,j integer itlo,ithi integer jtlo,jthi integer jmin,jmax integer imin,imax integer ivar integer startrec integer endrec integer diffrec #ifdef ALLOW_OBCS_CONTROL_MODES INTEGER k,length_of_rec,dUnit INTEGER MDS_RECLEN EXTERNAL MDS_RECLEN #endif c == external == integer ilnblnk external ilnblnk c == end of interface == jtlo = mybylo(mythid) jthi = mybyhi(mythid) itlo = mybxlo(mythid) ithi = mybxhi(mythid) jmin = 1-oly jmax = sny+oly imin = 1-olx imax = snx+olx c-- Set default values. do ivar = 1,maxcvars ncvarindex(ivar) = -1 ncvarrecs(ivar) = 0 ncvarxmax(ivar) = 0 ncvarymax(ivar) = 0 ncvarnrmax(ivar) = 0 ncvargrd(ivar) = '?' enddo _BARRIER c-- ===================== c-- Initial state fields. c-- ===================== cph( cph index 7-10 reserved for atmos. state, cph index 11-14 reserved for open boundaries, cph index 15-16 reserved for mixing coeff. cph index 17 reserved for passive tracer TR1 cph index 18,19 reserved for sst, sss cph index 20 for hFacC cph index 21-22 for efluxy, efluxp cph index 23 for bottom drag cph index 24 cph index 25-26 for edtaux, edtauy cph index 27-29 for uvel0, vvel0, etan0 cph index 30-31 for generic 2d, 3d field cph index 32 reserved for precip (atmos. state) cph index 33 reserved for swflux (atmos. state) cph index 34 reserved for swdown (atmos. state) cph 35 lwflux cph 36 lwdown cph 37 evap cph 38 snowprecip cph 39 apressure cph 40 runoff cph 41 seaice SIAREA cph 42 seaice SIHEFF cph 43 seaice SIHSNOW cph 44 gmredi kapredi cph 45 shelfice shifwflx cph) c---------------------------------------------------------------------- c-- #ifdef ALLOW_THETA0_CONTROL c-- Initial state temperature contribution. call ctrl_init_ctrlvar ( & xx_theta_file, 1, 101, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_THETA0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SALT0_CONTROL c-- Initial state salinity contribution. call ctrl_init_ctrlvar ( & xx_salt_file, 2, 102, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_SALT0_CONTROL */ c-- =========================== c-- Surface flux contributions. c-- =========================== c---------------------------------------------------------------------- c-- #if (defined (ALLOW_HFLUX_CONTROL)) c-- Heat flux. call ctrl_init_rec ( I xx_hfluxstartdate1, xx_hfluxstartdate2, xx_hfluxperiod, 1, O xx_hfluxstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_hflux_file, 3, 103, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_ATEMP_CONTROL)) c-- Atmos. temperature call ctrl_init_rec ( I xx_atempstartdate1, xx_atempstartdate2, xx_atempperiod, 1, O xx_atempstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_atemp_file, 7, 107, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_HFLUX0_CONTROL)) c-- initial forcing only call ctrl_init_ctrlvar ( & xx_hflux_file, 3, 103, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_HFLUX_CONTROL */ c---------------------------------------------------------------------- c-- #if (defined (ALLOW_SFLUX_CONTROL)) c-- Salt flux. call ctrl_init_rec ( I xx_sfluxstartdate1, xx_sfluxstartdate2, xx_sfluxperiod, 1, O xx_sfluxstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_sflux_file, 4, 104, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_AQH_CONTROL)) c-- Atmos. humidity call ctrl_init_rec ( I xx_aqhstartdate1, xx_aqhstartdate2, xx_aqhperiod, 1, O xx_aqhstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_aqh_file, 8, 108, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_SFLUX0_CONTROL)) c-- initial forcing only call ctrl_init_ctrlvar ( & xx_sflux_file, 4, 104, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_SFLUX_CONTROL */ c---------------------------------------------------------------------- c-- #if (defined (ALLOW_USTRESS_CONTROL)) c-- Zonal wind stress. call ctrl_init_rec ( I xx_tauustartdate1, xx_tauustartdate2, xx_tauuperiod, 1, O xx_tauustartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_tauu_file, 5, 105, diffrec, startrec, endrec, #ifndef ALLOW_ROTATE_UV_CONTROLS & snx, sny, 1, 'w', 'xy', mythid ) #else & snx, sny, 1, 'c', 'xy', mythid ) #endif #elif (defined (ALLOW_UWIND_CONTROL)) c-- Zonal wind speed. call ctrl_init_rec ( I xx_uwindstartdate1, xx_uwindstartdate2, xx_uwindperiod, 1, O xx_uwindstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_uwind_file, 9, 109, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_TAUU0_CONTROL)) c-- initial forcing only call ctrl_init_ctrlvar ( & xx_tauu_file, 5, 105, 1, 1, 1, & snx, sny, 1, 'w', 'xy', mythid ) #endif /* ALLOW_USTRESS_CONTROL */ c---------------------------------------------------------------------- c-- #if (defined (ALLOW_VSTRESS_CONTROL)) c-- Meridional wind stress. call ctrl_init_rec ( I xx_tauvstartdate1, xx_tauvstartdate2, xx_tauvperiod, 1, O xx_tauvstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_tauv_file, 6, 106, diffrec, startrec, endrec, #ifndef ALLOW_ROTATE_UV_CONTROLS & snx, sny, 1, 's', 'xy', mythid ) #else & snx, sny, 1, 'c', 'xy', mythid ) #endif #elif (defined (ALLOW_VWIND_CONTROL)) c-- Meridional wind speed. call ctrl_init_rec ( I xx_vwindstartdate1, xx_vwindstartdate2, xx_vwindperiod, 1, O xx_vwindstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_vwind_file, 10, 110, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_TAUV0_CONTROL)) c-- initial forcing only call ctrl_init_ctrlvar ( & xx_tauv_file, 6, 106, 1, 1, 1, & snx, sny, 1, 's', 'xy', mythid ) #endif /* ALLOW_VSTRESS_CONTROL */ c-- =========================== c-- Open boundary contributions. c-- =========================== c---------------------------------------------------------------------- c-- #ifdef ALLOW_OBCSN_CONTROL c-- Northern obc. call ctrl_init_rec ( I xx_obcsnstartdate1, xx_obcsnstartdate2, xx_obcsnperiod, 4, O xx_obcsnstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_obcsn_file, 11, 111, diffrec, startrec, endrec, & snx, 1, nr, 'm', 'xz', mythid ) #endif /* ALLOW_OBCSN_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_OBCSS_CONTROL c-- Southern obc. call ctrl_init_rec ( I xx_obcssstartdate1, xx_obcssstartdate2, xx_obcssperiod, 4, O xx_obcssstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_obcss_file, 12, 112, diffrec, startrec, endrec, & snx, 1, nr, 'm', 'xz', mythid ) #endif /* ALLOW_OBCSS_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_OBCSW_CONTROL c-- Western obc. call ctrl_init_rec ( I xx_obcswstartdate1, xx_obcswstartdate2, xx_obcswperiod, 4, O xx_obcswstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_obcsw_file, 13, 113, diffrec, startrec, endrec, & 1, sny, nr, 'm', 'yz', mythid ) #endif /* ALLOW_OBCSW_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_OBCSE_CONTROL c-- Eastern obc. call ctrl_init_rec ( I xx_obcsestartdate1, xx_obcsestartdate2, xx_obcseperiod, 4, O xx_obcsestartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_obcse_file, 14, 114, diffrec, startrec, endrec, & 1, sny, nr, 'm', 'yz', mythid ) #endif /* ALLOW_OBCSE_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_OBCS_CONTROL_MODES cih Get matrices for reconstruction from barotropic-barclinic modes CMM To use modes now hardcoded with ECCO_CPPOPTION. Would be good to have c run-time option and also define filename=baro_invmodes.bin CALL MDSFINDUNIT( dUnit, myThid ) length_of_rec = MDS_RECLEN( 64, NR*NR, myThid ) open(dUnit, file='baro_invmodes.bin', status='old', & access='direct', recl=length_of_rec ) do j = 1,Nr read(dUnit,rec=j) ((modesv(k,i,j), k=1,Nr), i=1,Nr) end do CLOSE( dUnit ) CMM double precision modesv is size [NR,NR,NR] c dim one is z-space c dim two is mode space c dim three is the total depth for which this set of modes applies c so for example modesv(:,2,nr) will be the second mode c in z-space for the full model depth c The modes are to be orthogonal when weighted by dz. c i.e. if f_i(z) = mode i, sum_j(f_i(z_j)*f_j(z_j)*dz_j = delta_ij c first mode should also be constant in depth...barotropic c For a matlab code example how to construct the orthonormal modes, c which are ideally the solution of planetary vertical mode equation c using model mean dRho/dz, see c MITgcm/verification/obcs_ctrl/input/gendata.m c This code is compatible with partial cells #endif c---------------------------------------------------------------------- c-- #ifdef ALLOW_DIFFKR_CONTROL call ctrl_init_ctrlvar ( & xx_diffkr_file, 15, 115, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_DIFFKR_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_KAPGM_CONTROL call ctrl_init_ctrlvar ( & xx_kapgm_file, 16, 116, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_KAPGM_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_TR10_CONTROL call ctrl_init_ctrlvar ( & xx_tr1_file, 17, 117, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_TR10_CONTROL */ c---------------------------------------------------------------------- c-- #if (defined (ALLOW_SST_CONTROL)) call ctrl_init_rec ( I xx_sststartdate1, xx_sststartdate2, xx_sstperiod, 1, O xx_sststartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_sst_file, 18, 118, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_SST0_CONTROL)) call ctrl_init_ctrlvar ( & xx_sst_file, 18, 118, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_SST_CONTROL */ c---------------------------------------------------------------------- c-- #if (defined (ALLOW_SSS_CONTROL)) call ctrl_init_rec ( I xx_sssstartdate1, xx_sssstartdate2, xx_sssperiod, 1, O xx_sssstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_sss_file, 19, 119, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #elif (defined (ALLOW_SSS0_CONTROL)) call ctrl_init_ctrlvar ( & xx_sss_file, 19, 119, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_SSS0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_DEPTH_CONTROL call ctrl_init_ctrlvar ( & xx_depth_file, 20, 120, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_DEPTH_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_EFLUXY0_CONTROL call ctrl_init_ctrlvar ( & xx_efluxy_file, 21, 121, 1, 1, 1, & snx, sny, nr, 's', '3d', mythid ) #endif /* ALLOW_EFLUXY0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_EFLUXP0_CONTROL call ctrl_init_ctrlvar ( & xx_efluxp_file, 22, 122, 1, 1, 1, & snx, sny, nr, 'v', '3d', mythid ) #endif /* ALLOW_EFLUXP0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_BOTTOMDRAG_CONTROL call ctrl_init_ctrlvar ( & xx_bottomdrag_file, 23, 123, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_BOTTOMDRAG_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_HFLUXM_CONTROL call ctrl_init_ctrlvar ( & xx_hfluxm_file, 24, 124, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_HFLUXM_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_EDDYPSI_CONTROL call ctrl_init_ctrlvar ( & xx_edtaux_file, 25, 125, 1, 1, 1, & snx, sny, nr, 'w', '3d', mythid ) call ctrl_init_ctrlvar ( & xx_edtauy_file, 26, 126, 1, 1, 1, & snx, sny, nr, 's', '3d', mythid ) #endif /* ALLOW_EDDYPSI_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_UVEL0_CONTROL call ctrl_init_ctrlvar ( & xx_uvel_file, 27, 127, 1, 1, 1, & snx, sny, nr, 'w', '3d', mythid ) #endif /* ALLOW_UVEL0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_VVEL0_CONTROL call ctrl_init_ctrlvar ( & xx_vvel_file, 28, 128, 1, 1, 1, & snx, sny, nr, 's', '3d', mythid ) #endif /* ALLOW_VVEL0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_ETAN0_CONTROL call ctrl_init_ctrlvar ( & xx_etan_file, 29, 129, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_VVEL0_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_GEN2D_CONTROL call ctrl_init_ctrlvar ( & xx_gen2d_file, 30, 130, 1, 1, 1, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_GEN2D_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_GEN3D_CONTROL call ctrl_init_ctrlvar ( & xx_gen3d_file, 31, 131, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_GEN3D_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_PRECIP_CONTROL c-- Atmos. precipitation call ctrl_init_rec ( I xx_precipstartdate1, xx_precipstartdate2, xx_precipperiod,1, O xx_precipstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_precip_file, 32, 132, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_PRECIP_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SWFLUX_CONTROL c-- Atmos. swflux call ctrl_init_rec ( I xx_swfluxstartdate1, xx_swfluxstartdate2, xx_swfluxperiod, 1, O xx_swfluxstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_swflux_file, 33, 133, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_SWFLUX_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SWDOWN_CONTROL c-- Atmos. swdown call ctrl_init_rec ( I xx_swdownstartdate1, xx_swdownstartdate2, xx_swdownperiod, 1, O xx_swdownstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_swdown_file, 34, 134, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_SWDOWN_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_LWFLUX_CONTROL c-- Atmos. lwflux call ctrl_init_rec ( I xx_lwfluxstartdate1, xx_lwfluxstartdate2, xx_lwfluxperiod, 1, O xx_lwfluxstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_lwflux_file, 35, 135, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_LWFLUX_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_LWDOWN_CONTROL c-- Atmos. lwdown call ctrl_init_rec ( I xx_lwdownstartdate1, xx_lwdownstartdate2, xx_lwdownperiod, 1, O xx_lwdownstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_lwdown_file, 36, 136, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_LWDOWN_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_EVAP_CONTROL c-- Atmos. evap call ctrl_init_rec ( I xx_evapstartdate1, xx_evapstartdate2, xx_evapperiod, 1, O xx_evapstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_evap_file, 37, 137, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_EVAP_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SNOWPRECIP_CONTROL c-- Atmos. snowprecip call ctrl_init_rec ( I xx_snowprecipstartdate1, xx_snowprecipstartdate2, I xx_snowprecipperiod, 1, O xx_snowprecipstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_snowprecip_file, 38, 138, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_SNOWPRECIP_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_APRESSURE_CONTROL c-- Atmos. apressure call ctrl_init_rec ( I xx_apressurestartdate1, xx_apressurestartdate2, I xx_apressureperiod, 1, O xx_apressurestartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_apressure_file, 39, 139, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_APRESSURE_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_RUNOFF_CONTROL c-- Atmos. runoff call ctrl_init_rec ( I xx_runoffstartdate1, xx_runoffstartdate2, xx_runoffperiod, 1, O xx_runoffstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_runoff_file, 40, 140, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_RUNOFF_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SIAREA_CONTROL C-- so far there are no xx_siareastartdate1, etc., so we need to fudge it. CML call ctrl_init_rec ( CML I xx_siareastartdate1, xx_siareastartdate2, xx_siareaperiod, 1, CML O xx_siareastartdate, diffrec, startrec, endrec, CML I mythid ) startrec = 1 endrec = 1 diffrec = endrec - startrec + 1 call ctrl_init_ctrlvar ( & xx_siarea_file, 41, 141, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_siarea_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SIHEFF_CONTROL C-- so far there are no xx_siheffstartdate1, etc., so we need to fudge it. CML call ctrl_init_rec ( CML I xx_siheffstartdate1, xx_siheffstartdate2, xx_siheffperiod, 1, CML O xx_siheffstartdate, diffrec, startrec, endrec, CML I mythid ) startrec = 1 endrec = 1 diffrec = endrec - startrec + 1 call ctrl_init_ctrlvar ( & xx_siheff_file, 42, 142, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_siheff_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_SIHSNOW_CONTROL C-- so far there are no xx_sihsnowstartdate1, etc., so we need to fudge it. CML call ctrl_init_rec ( CML I xx_sihsnowstartdate1, xx_sihsnowstartdate2, xx_sihsnowperiod, 1, CML O xx_sihsnowstartdate, diffrec, startrec, endrec, CML I mythid ) startrec = 1 endrec = 1 diffrec = endrec - startrec + 1 call ctrl_init_ctrlvar ( & xx_sihsnow_file, 43, 143, diffrec, startrec, endrec, & snx, sny, 1, 'c', 'xy', mythid ) #endif /* ALLOW_sihsnow_CONTROL */ c---------------------------------------------------------------------- c-- #ifdef ALLOW_KAPREDI_CONTROL call ctrl_init_ctrlvar ( & xx_kapredi_file, 44, 144, 1, 1, 1, & snx, sny, nr, 'c', '3d', mythid ) #endif /* ALLOW_KAPREDI_CONTROL */ c---------------------------------------------------------------------- c---------------------------------------------------------------------- #ifdef ALLOW_SHIFWFLX_CONTROL c-- freshwater flux underneath ice-shelves call ctrl_init_rec ( I xx_shifwflxstartdate1, xx_shifwflxstartdate2, I xx_shifwflxperiod, 1, O xx_shifwflxstartdate, diffrec, startrec, endrec, I mythid ) call ctrl_init_ctrlvar ( & xx_shifwflx_file, 45, 145, diffrec, startrec, endrec, & snx, sny, 1, 'i', 'xy', mythid ) #endif /* ALLOW_SHIFWFLX_CONTROL */ c---------------------------------------------------------------------- c---------------------------------------------------------------------- call ctrl_init_wet( mythid ) c---------------------------------------------------------------------- c---------------------------------------------------------------------- #ifdef ALLOW_DIC_CONTROL do i = 1, dic_n_control xx_dic(i) = 0. _d 0 enddo #endif c---------------------------------------------------------------------- c---------------------------------------------------------------------- do bj = jtlo,jthi do bi = itlo,ithi do j = jmin,jmax do i = imin,imax wareaunit (i,j,bi,bj) = 1.0 #ifndef ALLOW_ECCO whflux (i,j,bi,bj) = maskC(i,j,1,bi,bj) wsflux (i,j,bi,bj) = maskC(i,j,1,bi,bj) wtauu (i,j,bi,bj) = maskW(i,j,1,bi,bj) wtauv (i,j,bi,bj) = maskS(i,j,1,bi,bj) watemp (i,j,bi,bj) = maskC(i,j,1,bi,bj) waqh (i,j,bi,bj) = maskC(i,j,1,bi,bj) wprecip (i,j,bi,bj) = maskC(i,j,1,bi,bj) wswflux (i,j,bi,bj) = maskC(i,j,1,bi,bj) wswdown (i,j,bi,bj) = maskC(i,j,1,bi,bj) wuwind (i,j,bi,bj) = maskC(i,j,1,bi,bj) wvwind (i,j,bi,bj) = maskC(i,j,1,bi,bj) wlwflux (i,j,bi,bj) = maskC(i,j,1,bi,bj) wlwdown (i,j,bi,bj) = maskC(i,j,1,bi,bj) wevap (i,j,bi,bj) = maskC(i,j,1,bi,bj) wsnowprecip(i,j,bi,bj) = maskC(i,j,1,bi,bj) wapressure(i,j,bi,bj) = maskC(i,j,1,bi,bj) wrunoff (i,j,bi,bj) = maskC(i,j,1,bi,bj) wsst (i,j,bi,bj) = maskC(i,j,1,bi,bj) wsss (i,j,bi,bj) = maskC(i,j,1,bi,bj) #endif enddo enddo enddo enddo return end subroutine ctrl_init_rec( I fldstartdate1, fldstartdate2, fldperiod, nfac, O fldstartdate, diffrec, startrec, endrec, I mythid ) c ================================================================== c SUBROUTINE ctrl_init_rec c ================================================================== c c helper routine to compute the first and last record of a c time dependent control variable c c Martin.Losch@awi.de, 2011-Mar-15 c c ================================================================== c SUBROUTINE ctrl_init_rec c ================================================================== implicit none c == global variables == #include "SIZE.h" #include "EEPARAMS.h" #include "PARAMS.h" #ifdef ALLOW_CAL # include "cal.h" #endif c == input variables == c fldstartdate1/2 : start time (date/time) of fld c fldperod : sampling interval of fld c nfac : factor for the case that fld is an obcs variable c in this case nfac = 4, otherwise nfac = 1 c mythid : thread ID of this instance integer fldstartdate1 integer fldstartdate2 _RL fldperiod integer nfac integer mythid c == output variables == c fldstartdate : full date from fldstartdate1 and 2 c startrec : first record of ctrl variable c startrec : last record of ctrl variable c diffrec : difference between first and last record of ctrl variable integer fldstartdate(4) integer startrec integer endrec integer diffrec c == local variables == integer i #ifdef ALLOW_CAL integer difftime(4) _RL diffsecs #endif /* ALLOW_CAL */ c initialise some output do i = 1,4 fldstartdate(i) = 0 end do startrec = 0 endrec = 0 diffrec = 0 # ifdef ALLOW_CAL call cal_FullDate( fldstartdate1, fldstartdate2, & fldstartdate , mythid ) call cal_TimePassed( fldstartdate, modelstartdate, & difftime, mythid ) call cal_ToSeconds ( difftime, diffsecs, mythid ) if ( fldperiod .EQ. -12. ) then startrec = 1 endrec = 12*nfac elseif ( fldperiod .EQ. 0. ) then startrec = 1 endrec = 1*nfac else startrec = int((modelstart + startTime - diffsecs)/ & fldperiod) + 1 endrec = int((modelend + startTime - diffsecs + modelstep/2)/ & fldperiod) + 2 if ( nfac .ne. 1 ) then c This is the case of obcs. endrec = endrec*nfac + 1 startrec = (startrec - 1)*nfac + 1 endif endif # else /* ndef ALLOW_CAL */ if ( fldperiod .EQ. 0. ) then startrec = 1 endrec = 1*nfac else startrec = 1 endrec = (int((endTime - startTime)/fldperiod) + 1)*nfac endif #endif /* ALLOW_CAL */ diffrec = endrec - startrec + 1 return end