New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 7731 for branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC – NEMO

Ignore:
Timestamp:
2017-02-23T14:23:32+01:00 (7 years ago)
Author:
dford
Message:

Merge in revisions 6625:7726 of dev_r5518_v3.4_asm_nemovar_community, so this branch will be identical to revison 7726 of dev_r5518_v3.6_asm_nemovar_community.

Location:
branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC
Files:
48 edited
2 copied

Legend:

Unmodified
Added
Removed

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ASM/asmbkg.F90

    r7730 r7731  
    5050   USE ice 
    5151#endif 
     52   USE asminc, ONLY: ln_avgbkg 
    5253   IMPLICIT NONE 
    5354   PRIVATE 
    5455    
    5556   PUBLIC   asm_bkg_wri   !: Write out the background state 
     57 
     58  !! * variables for calculating time means 
     59   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   tn_tavg  , sn_tavg   
     60   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   un_tavg  , vn_tavg 
     61   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   avt_tavg 
     62#if defined key_zdfgls || key_zdftke 
     63   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:,:) ::   en_tavg 
     64#endif 
     65   REAL(wp),SAVE, ALLOCATABLE,   DIMENSION(:,:)   ::   sshn_tavg 
     66   REAL(wp),SAVE :: numtimes_tavg     ! No of times to average over 
    5667 
    5768   !!---------------------------------------------------------------------- 
    … …  
    8192      INTEGER :: inum          ! File unit number 
    8293      REAL(wp) :: zdate        ! Date 
     94      INTEGER :: ierror 
    8395      !!----------------------------------------------------------------------- 
    8496 
    85       !                                !------------------------------------------- 
    86       IF( kt == nitbkg_r ) THEN        ! Write out background at time step nitbkg_r 
    87          !                             !-----------------------------------======== 
     97      ! If creating an averaged assim bkg, initialise on first timestep 
     98      IF ( ln_avgbkg .AND. kt == ( nn_it000 - 1) ) THEN 
     99         ! Allocate memory  
     100         ALLOCATE( tn_tavg(jpi,jpj,jpk), STAT=ierror ) 
     101         IF( ierror > 0 ) THEN 
     102            CALL ctl_stop( 'asm_wri_bkg: unable to allocate tn_tavg' )   ;   RETURN 
     103         ENDIF 
     104         tn_tavg(:,:,:)=0 
     105         ALLOCATE( sn_tavg(jpi,jpj,jpk), STAT=ierror ) 
     106         IF( ierror > 0 ) THEN 
     107            CALL ctl_stop( 'asm_wri_bkg: unable to allocate sn_tavg' )   ;   RETURN 
     108         ENDIF 
     109         sn_tavg(:,:,:)=0 
     110         ALLOCATE( un_tavg(jpi,jpj,jpk), STAT=ierror ) 
     111         IF( ierror > 0 ) THEN 
     112            CALL ctl_stop( 'asm_wri_bkg: unable to allocate un_tavg' )   ;   RETURN 
     113         ENDIF 
     114         un_tavg(:,:,:)=0 
     115         ALLOCATE( vn_tavg(jpi,jpj,jpk), STAT=ierror ) 
     116         IF( ierror > 0 ) THEN 
     117            CALL ctl_stop( 'asm_wri_bkg: unable to allocate vn_tavg' )   ;   RETURN 
     118         ENDIF 
     119         vn_tavg(:,:,:)=0 
     120         ALLOCATE( sshn_tavg(jpi,jpj), STAT=ierror ) 
     121         IF( ierror > 0 ) THEN 
     122            CALL ctl_stop( 'asm_wri_bkg: unable to allocate sshn_tavg' )   ;   RETURN 
     123         ENDIF 
     124         sshn_tavg(:,:)=0 
     125#if defined key_zdftke 
     126         ALLOCATE( en_tavg(jpi,jpj,jpk), STAT=ierror ) 
     127         IF( ierror > 0 ) THEN 
     128            CALL ctl_stop( 'asm_wri_bkg: unable to allocate en_tavg' )   ;   RETURN 
     129         ENDIF 
     130         en_tavg(:,:,:)=0 
     131#endif 
     132         ALLOCATE( avt_tavg(jpi,jpj,jpk), STAT=ierror ) 
     133         IF( ierror > 0 ) THEN 
     134            CALL ctl_stop( 'asm_wri_bkg: unable to allocate avt_tavg' )   ;   RETURN 
     135         ENDIF 
     136         avt_tavg(:,:,:)=0 
     137          
     138         numtimes_tavg = REAL ( nitavgbkg_r -  nn_it000 + 1 ) 
     139      ENDIF    
     140 
     141      ! If creating an averaged assim bkg, sum the contribution every timestep 
     142      IF ( ln_avgbkg ) THEN  
     143         IF (lwp) THEN 
     144              WRITE(numout,*) 'asm_wri_bkg : Summing assim bkg fields at timestep ',kt 
     145              WRITE(numout,*) '~~~~~~~~~~~~ ' 
     146         ENDIF 
     147 
     148         tn_tavg(:,:,:)        = tn_tavg(:,:,:) + tsn(:,:,:,jp_tem) / numtimes_tavg 
     149         sn_tavg(:,:,:)        = sn_tavg(:,:,:) + tsn(:,:,:,jp_sal) / numtimes_tavg 
     150         sshn_tavg(:,:)        = sshn_tavg(:,:) + sshn (:,:) / numtimes_tavg 
     151         un_tavg(:,:,:)        = un_tavg(:,:,:) + un(:,:,:) / numtimes_tavg 
     152         vn_tavg(:,:,:)        = vn_tavg(:,:,:) + vn(:,:,:) / numtimes_tavg 
     153         avt_tavg(:,:,:)        = avt_tavg(:,:,:) + avt(:,:,:) / numtimes_tavg 
     154#if defined key_zdftke 
     155         en_tavg(:,:,:)       = en_tavg(:,:,:) + en(:,:,:) / numtimes_tavg 
     156#endif 
     157      ENDIF 
     158      
     159 
     160      ! Write out background at time step nitbkg_r or nitavgbkg_r 
     161      IF ( ( .NOT. ln_avgbkg .AND. (kt == nitbkg_r) ) .OR. & 
     162      &          ( ln_avgbkg .AND. (kt == nitavgbkg_r) ) ) THEN 
    88163         ! 
    89164         WRITE(cl_asmbkg, FMT='(A,".nc")' ) TRIM( c_asmbkg ) 
    … …  
    97172            CALL iom_open( c_asmbkg, inum, ldwrt = .TRUE., kiolib = jprstlib) 
    98173            ! 
    99             IF( nitbkg_r == nit000 - 1 ) THEN      ! Treat special case when nitbkg = 0 
    100                zdate = REAL( ndastp ) 
    101 #if defined key_zdftke 
    102                ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
    103                IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...' 
    104                CALL tke_rst( nit000, 'READ' )               ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
    105  
    106 #endif 
     174            ! 
     175            ! Write the information 
     176            IF ( ln_avgbkg ) THEN 
     177               IF( nitavgbkg_r == nit000 - 1 ) THEN      ! Treat special case when nitavgbkg = 0 
     178                  zdate = REAL( ndastp ) 
     179#if defined key_zdftke 
     180                  ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
     181                  IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...' 
     182                  CALL tke_rst( nit000, 'READ' )               ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
     183 
     184#endif 
     185               ELSE 
     186                  zdate = REAL( ndastp ) 
     187               ENDIF 
     188               CALL iom_rstput( kt, nitavgbkg_r, inum, 'rdastp' , zdate   ) 
     189               CALL iom_rstput( kt, nitavgbkg_r, inum, 'un'     , un_tavg ) 
     190               CALL iom_rstput( kt, nitavgbkg_r, inum, 'vn'     , vn_tavg ) 
     191               CALL iom_rstput( kt, nitavgbkg_r, inum, 'tn'     , tn_tavg ) 
     192               CALL iom_rstput( kt, nitavgbkg_r, inum, 'sn'     , sn_tavg ) 
     193               CALL iom_rstput( kt, nitavgbkg_r, inum, 'sshn'   , sshn_tavg) 
     194#if defined key_zdftke 
     195               CALL iom_rstput( kt, nitavgbkg_r, inum, 'en'     , en_tavg ) 
     196#endif 
     197               CALL iom_rstput( kt, nitavgbkg_r, inum, 'avt'    , avt_tavg) 
     198               ! 
    107199            ELSE 
    108                zdate = REAL( ndastp ) 
     200               IF( nitbkg_r == nit000 - 1 ) THEN      ! Treat special case when nitbkg = 0 
     201                  zdate = REAL( ndastp ) 
     202#if defined key_zdftke 
     203                  ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
     204                  IF(lwp) WRITE(numout,*) ' Reading TKE (en) from restart...' 
     205                  CALL tke_rst( nit000, 'READ' )               ! lk_zdftke=T :   Read turbulent kinetic energy ( en ) 
     206 
     207#endif 
     208               ELSE 
     209                  zdate = REAL( ndastp ) 
     210               ENDIF 
     211               CALL iom_rstput( kt, nitbkg_r, inum, 'rdastp' , zdate   ) 
     212               CALL iom_rstput( kt, nitbkg_r, inum, 'un'     , un                ) 
     213               CALL iom_rstput( kt, nitbkg_r, inum, 'vn'     , vn                ) 
     214               CALL iom_rstput( kt, nitbkg_r, inum, 'tn'     , tsn(:,:,:,jp_tem) ) 
     215               CALL iom_rstput( kt, nitbkg_r, inum, 'sn'     , tsn(:,:,:,jp_sal) ) 
     216               CALL iom_rstput( kt, nitbkg_r, inum, 'sshn'   , sshn              ) 
     217#if defined key_zdftke 
     218               CALL iom_rstput( kt, nitbkg_r, inum, 'en'     , en                ) 
     219#endif 
     220               CALL iom_rstput( kt, nitbkg_r, inum, 'avt'    , avt               ) 
     221               ! 
    109222            ENDIF 
    110             ! 
    111             !                                      ! Write the information 
    112             CALL iom_rstput( kt, nitbkg_r, inum, 'rdastp' , zdate             ) 
    113             CALL iom_rstput( kt, nitbkg_r, inum, 'un'     , un                ) 
    114             CALL iom_rstput( kt, nitbkg_r, inum, 'vn'     , vn                ) 
    115             CALL iom_rstput( kt, nitbkg_r, inum, 'tn'     , tsn(:,:,:,jp_tem) ) 
    116             CALL iom_rstput( kt, nitbkg_r, inum, 'sn'     , tsn(:,:,:,jp_sal) ) 
    117             CALL iom_rstput( kt, nitbkg_r, inum, 'sshn'   , sshn              ) 
    118 #if defined key_zdftke 
    119             CALL iom_rstput( kt, nitbkg_r, inum, 'en'     , en                ) 
    120 #endif 
    121             CALL iom_rstput( kt, nitbkg_r, inum, 'gcx'    , gcx               ) 
    122             ! 
     223             
    123224            CALL iom_close( inum ) 
     225 
    124226         ENDIF 
    125227         ! 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90

    r7730 r7731  
    2222   !!   ssh_asm_inc    : Apply the SSH increment 
    2323   !!   seaice_asm_inc : Apply the seaice increment 
     24   !!   logchl_asm_inc : Apply the logchl increment 
    2425   !!---------------------------------------------------------------------- 
    2526   USE wrk_nemo         ! Memory Allocation 
    … …  
    4041#endif 
    4142   USE sbc_oce          ! Surface boundary condition variables. 
     43   USE zdfmxl, ONLY :  &   
     44   &  hmld_tref,       &    
     45#if defined key_karaml 
     46   &  hmld_kara,       & 
     47   &  ln_kara,         & 
     48#endif    
     49   &  hmld,            &  
     50   &  hmlp,            & 
     51   &  hmlpt 
     52#if defined key_bdy  
     53   USE bdy_oce, ONLY: bdytmask   
     54#endif   
     55#if defined key_top 
     56   USE trc, ONLY: & 
     57      & trn,      & 
     58      & trb 
     59   USE par_trc, ONLY: & 
     60      & jptra 
     61#endif 
     62#if defined key_fabm 
     63   USE asmlogchlbal_ersem, ONLY: & 
     64      & asm_logchl_bal_ersem 
     65   USE par_fabm 
     66#elif defined key_medusa && defined key_foam_medusa 
     67   USE asmlogchlbal_medusa, ONLY: & 
     68      & asm_logchl_bal_medusa 
     69   USE par_medusa 
     70#elif defined key_hadocc 
     71   USE asmlogchlbal_hadocc, ONLY: & 
     72      & asm_logchl_bal_hadocc 
     73   USE par_hadocc 
     74#endif 
    4275 
    4376   IMPLICIT NONE 
    … …  
    5083   PUBLIC   ssh_asm_inc    !: Apply the SSH increment 
    5184   PUBLIC   seaice_asm_inc !: Apply the seaice increment 
     85   PUBLIC   logchl_asm_inc !: Apply the logchl increment 
    5286 
    5387#if defined key_asminc 
    … …  
    5791#endif 
    5892   LOGICAL, PUBLIC :: ln_bkgwri = .FALSE.      !: No output of the background state fields 
     93   LOGICAL, PUBLIC :: ln_balwri = .FALSE.      !: No output of the assimilation balancing increments 
     94   LOGICAL, PUBLIC :: ln_avgbkg = .FALSE.      !: No output of the mean background state fields 
    5995   LOGICAL, PUBLIC :: ln_asmiau = .FALSE.      !: No applying forcing with an assimilation increment 
    6096   LOGICAL, PUBLIC :: ln_asmdin = .FALSE.      !: No direct initialization 
    … …  
    6298   LOGICAL, PUBLIC :: ln_dyninc = .FALSE.      !: No dynamics (u and v) assimilation increments 
    6399   LOGICAL, PUBLIC :: ln_sshinc = .FALSE.      !: No sea surface height assimilation increment 
    64    LOGICAL, PUBLIC :: ln_seaiceinc             !: No sea ice concentration increment 
     100   LOGICAL, PUBLIC :: ln_seaiceinc = .FALSE.   !: No sea ice concentration increment 
     101   LOGICAL, PUBLIC :: ln_logchltotinc = .FALSE. !: No total log10(chlorophyll) increment 
     102   LOGICAL, PUBLIC :: ln_logchlpftinc = .FALSE. !: No PFT   log10(chlorophyll) increment 
    65103   LOGICAL, PUBLIC :: ln_salfix = .FALSE.      !: Apply minimum salinity check 
    66104   LOGICAL, PUBLIC :: ln_temnofreeze = .FALSE. !: Don't allow the temperature to drop below freezing 
    … …  
    80118   INTEGER , PUBLIC ::   nitiaustr   !: Time step of the start of the IAU interval  
    81119   INTEGER , PUBLIC ::   nitiaufin   !: Time step of the end of the IAU interval 
     120   INTEGER , PUBLIC ::   nitavgbkg   !: Number of timesteps to average assim bkg [0,nitavgbkg] 
    82121   !  
    83122   INTEGER , PUBLIC ::   niaufn      !: Type of IAU weighing function: = 0   Constant weighting 
    … …  
    87126   REAL(wp), DIMENSION(:,:), ALLOCATABLE ::   ssh_bkg, ssh_bkginc   ! Background sea surface height and its increment 
    88127   REAL(wp), DIMENSION(:,:), ALLOCATABLE ::   seaice_bkginc         ! Increment to the background sea ice conc 
     128 
     129   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE :: logchl_bkginc !: Increment to background logchl 
     130#if defined key_top 
     131   REAL(wp), PUBLIC, DIMENSION(:,:,:,:), ALLOCATABLE :: logchl_balinc  !: Increment to BGC variables from logchl assim 
     132#endif 
     133 
     134   INTEGER :: mld_choice        = 4   !: choice of mld criteria to use for physics assimilation 
     135                                      !: 1) hmld      - Turbocline/mixing depth                           [W points] 
     136                                      !: 2) hmlp      - Density criterion (0.01 kg/m^3 change from 10m)   [W points] 
     137                                      !: 3) hmld_kara - Kara MLD                                          [Interpolated] 
     138                                      !: 4) hmld_tref - Temperature criterion (0.2 K change from surface) [T points] 
     139 
     140   INTEGER :: mld_choice_bgc    = 4   !: choice of mld criteria to use for physics assimilation 
     141                                      !: 1) hmld      - Turbocline/mixing depth                           [W points] 
     142                                      !: 2) hmlp      - Density criterion (0.01 kg/m^3 change from 10m)   [W points] 
     143                                      !: 3) hmld_kara - Kara MLD                                          [Interpolated] 
     144                                      !: 4) hmld_tref - Temperature criterion (0.2 K change from surface) [T points] 
     145                                      !: 5) hmlpt     - Density criterion (0.01 kg/m^3 change from 10m)   [T points] 
     146 
     147   INTEGER :: nn_asmpfts        = 0   !: number of logchl PFTs assimilated 
    89148 
    90149   !! * Substitutions 
    … …  
    119178      INTEGER :: iitiaustr_date  ! Date YYYYMMDD of IAU interval start time step 
    120179      INTEGER :: iitiaufin_date  ! Date YYYYMMDD of IAU interval final time step 
     180      INTEGER :: isurfstat       ! Local integer for status of reading surft variable 
     181      INTEGER :: iitavgbkg_date  ! Date YYYYMMDD of end of assim bkg averaging period 
    121182      ! 
    122183      REAL(wp) :: znorm        ! Normalization factor for IAU weights 
    … …  
    127188      REAL(wp) :: zdate_inc    ! Time axis in increments file 
    128189      ! 
     190      REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: &  
     191          &       t_bkginc_2d  ! file for reading in 2D   
     192      !                        ! temperature increments  
     193      REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: &  
     194          &       z_mld     ! Mixed layer depth  
     195           
    129196      REAL(wp), POINTER, DIMENSION(:,:) ::   hdiv   ! 2D workspace 
    130       !! 
    131       NAMELIST/nam_asminc/ ln_bkgwri,                                      & 
     197      ! 
     198      LOGICAL :: lk_surft      ! Logical: T => Increments file contains surft variable  
     199                               !               so only apply surft increments. 
     200      ! 
     201      CHARACTER(LEN=2) :: cl_pftstr 
     202      !! 
     203      NAMELIST/nam_asminc/ ln_bkgwri, ln_balwri, ln_avgbkg,                & 
    132204         &                 ln_trainc, ln_dyninc, ln_sshinc,                & 
     205         &                 ln_logchltotinc, ln_logchlpftinc,               & 
    133206         &                 ln_asmdin, ln_asmiau,                           & 
    134207         &                 nitbkg, nitdin, nitiaustr, nitiaufin, niaufn,   & 
    135          &                 ln_salfix, salfixmin, nn_divdmp 
     208         &                 ln_salfix, salfixmin, nn_divdmp, nitavgbkg, mld_choice, & 
     209         &                 mld_choice_bgc 
    136210      !!---------------------------------------------------------------------- 
    137211 
    … …  
    139213      ! Read Namelist nam_asminc : assimilation increment interface 
    140214      !----------------------------------------------------------------------- 
     215 
     216      ! Set default values 
     217      ln_bkgwri = .FALSE. 
     218      ln_balwri = .FALSE. 
     219      ln_avgbkg = .FALSE. 
     220      ln_trainc = .FALSE. 
     221      ln_dyninc = .FALSE. 
     222      ln_sshinc = .FALSE. 
    141223      ln_seaiceinc = .FALSE. 
     224      ln_logchltotinc = .FALSE. 
     225      ln_logchlpftinc = .FALSE. 
     226      ln_asmdin = .FALSE. 
     227      ln_asmiau = .TRUE. 
     228      ln_salfix = .FALSE. 
    142229      ln_temnofreeze = .FALSE. 
     230      salfixmin = -9999 
     231      nitbkg    = 0 
     232      nitdin    = 0       
     233      nitiaustr = 1 
     234      nitiaufin = 150 
     235      niaufn    = 0 
     236      nitavgbkg = 1 
     237#if defined key_fabm 
     238      nn_asmpfts = 4 
     239#elif defined key_medusa && defined key_foam_medusa 
     240      nn_asmpfts = 2 
     241#elif defined key_hadocc 
     242      nn_asmpfts = 1 
     243#else 
     244      nn_asmpfts = 0 
     245#endif 
    143246 
    144247      REWIND( numnam_ref )              ! Namelist nam_asminc in reference namelist : Assimilation increment 
    … …  
    150253902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_asminc in configuration namelist', lwp ) 
    151254      IF(lwm) WRITE ( numond, nam_asminc ) 
     255 
     256      IF ( ln_logchltotinc ) THEN 
     257         nn_asmpfts = 1 
     258      ELSE IF ( .NOT.( ln_logchlpftinc ) ) THEN 
     259         nn_asmpfts = 0 
     260      ENDIF 
    152261 
    153262      ! Control print 
    … …  
    158267         WRITE(numout,*) '   Namelist namasm : set assimilation increment parameters' 
    159268         WRITE(numout,*) '      Logical switch for writing out background state          ln_bkgwri = ', ln_bkgwri 
     269         WRITE(numout,*) '      Logical switch for writing out balancing increments      ln_balwri = ', ln_balwri 
     270         WRITE(numout,*) '      Logical switch for writing mean background state         ln_avgbkg = ', ln_avgbkg 
    160271         WRITE(numout,*) '      Logical switch for applying tracer increments            ln_trainc = ', ln_trainc 
    161272         WRITE(numout,*) '      Logical switch for applying velocity increments          ln_dyninc = ', ln_dyninc 
    … …  
    163274         WRITE(numout,*) '      Logical switch for Direct Initialization (DI)            ln_asmdin = ', ln_asmdin 
    164275         WRITE(numout,*) '      Logical switch for applying sea ice increments        ln_seaiceinc = ', ln_seaiceinc 
     276         WRITE(numout,*) '      Logical switch for applying total logchl incs      ln_logchltotinc = ', ln_logchltotinc 
     277         WRITE(numout,*) '      Logical switch for applying PFT   logchl incs      ln_logchlpftinc = ', ln_logchlpftinc 
     278         WRITE(numout,*) '      Number of logchl PFTs assimilated                       nn_asmpfts = ', nn_asmpfts 
    165279         WRITE(numout,*) '      Logical switch for Incremental Analysis Updating (IAU)   ln_asmiau = ', ln_asmiau 
    166280         WRITE(numout,*) '      Timestep of background in [0,nitend-nit000-1]            nitbkg    = ', nitbkg 
    … …  
    168282         WRITE(numout,*) '      Timestep of start of IAU interval in [0,nitend-nit000-1] nitiaustr = ', nitiaustr 
    169283         WRITE(numout,*) '      Timestep of end of IAU interval in [0,nitend-nit000-1]   nitiaufin = ', nitiaufin 
     284         WRITE(numout,*) '      Number of timesteps to average assim bkg [0,nitavgbkg]   nitavgbkg = ', nitavgbkg 
    170285         WRITE(numout,*) '      Type of IAU weighting function                           niaufn    = ', niaufn 
    171286         WRITE(numout,*) '      Logical switch for ensuring that the sa > salfixmin      ln_salfix = ', ln_salfix 
    172287         WRITE(numout,*) '      Minimum salinity after applying the increments           salfixmin = ', salfixmin 
     288         WRITE(numout,*) '      Choice of MLD for physics assimilation                  mld_choice = ', mld_choice 
     289         WRITE(numout,*) '      Choice of MLD for BGC assimilation                  mld_choice_bgc = ', mld_choice_bgc 
    173290      ENDIF 
    174291 
    … …  
    177294      nitiaustr_r = nitiaustr + nit000 - 1  ! Start of IAU interval referenced to nit000 
    178295      nitiaufin_r = nitiaufin + nit000 - 1  ! End of IAU interval referenced to nit000 
     296      nitavgbkg_r = nitavgbkg + nit000 - 1  ! Averaging period referenced to nit000 
    179297 
    180298      iiauper = nitiaufin_r - nitiaustr_r + 1  ! IAU interval length 
    … …  
    186304      CALL calc_date( nit000, nitiaustr_r, ndate0, iitiaustr_date )     ! IAU start time referenced to ndate0 
    187305      CALL calc_date( nit000, nitiaufin_r, ndate0, iitiaufin_date )     ! IAU end time referenced to ndate0 
     306      CALL calc_date( nit000, nitavgbkg_r, ndate0, iitavgbkg_date )     ! End of assim bkg averaging period referenced to ndate0 
    188307      ! 
    189308      IF(lwp) THEN 
    … …  
    197316         WRITE(numout,*) '       nitiaustr_r = ', nitiaustr_r 
    198317         WRITE(numout,*) '       nitiaufin_r = ', nitiaufin_r 
     318         WRITE(numout,*) '       nitavgbkg_r = ', nitavgbkg_r 
    199319         WRITE(numout,*) 
    200320         WRITE(numout,*) '   Dates referenced to current cycle:' 
    … …  
    206326         WRITE(numout,*) '       iitiaustr_date = ', iitiaustr_date 
    207327         WRITE(numout,*) '       iitiaufin_date = ', iitiaufin_date 
     328         WRITE(numout,*) '       iitavgbkg_date = ', iitavgbkg_date 
    208329      ENDIF 
    209330 
    … …  
    218339 
    219340      IF (      ( ( .NOT. ln_asmdin ).AND.( .NOT. ln_asmiau ) ) & 
    220            .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ) .OR. ( ln_seaiceinc) )) & 
    221          & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc and ln_seaiceinc is set to .true.', & 
     341         & .AND.( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ).OR. & 
     342         &        ( ln_logchltotinc ).OR.( ln_logchlpftinc ) )) & 
     343         & CALL ctl_stop( ' One or more of ln_trainc, ln_dyninc, ln_sshinc, ln_seaiceinc,', & 
     344         &                ' ln_logchltotinc, and ln_logchlpftinc is set to .true.', & 
    222345         &                ' but ln_asmdin and ln_asmiau are both set to .false. :', & 
    223346         &                ' Inconsistent options') 
    … …  
    228351 
    229352      IF ( ( .NOT. ln_trainc ).AND.( .NOT. ln_dyninc ).AND.( .NOT. ln_sshinc ).AND.( .NOT. ln_seaiceinc ) & 
    230          &                     )  & 
    231          & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc and ln_seaiceinc are set to .false. :', & 
     353         & .AND.( .NOT. ln_logchltotinc ).AND.( .NOT. ln_logchlpftinc ) )  & 
     354         & CALL ctl_warn( ' ln_trainc, ln_dyninc, ln_sshinc, ln_seaiceinc,', & 
     355         &                ' ln_logchltotinc, and ln_logchlpftinc are set to .false. :', & 
    232356         &                ' The assimilation increments are not applied') 
    233357 
    … …  
    249373         &                ' Background time step for Direct Initialization is outside', & 
    250374         &                ' the cycle interval') 
     375 
     376      IF ( nitavgbkg_r > nitend ) & 
     377         & CALL ctl_stop( ' nitavgbkg_r :',  & 
     378         &                ' Assim bkg averaging period is outside', & 
     379         &                ' the cycle interval') 
     380 
     381      IF ( ( ln_logchltotinc ).AND.( ln_logchlpftinc ) ) THEN 
     382         CALL ctl_stop( ' ln_logchltotinc and ln_logchlpftinc both set:', & 
     383            &           ' These options are not compatible') 
     384      ENDIF 
     385 
     386      IF ( ( ln_balwri ).AND.( .NOT. ( ( ln_logchltotinc ).OR.( ln_logchlpftinc ) ) ) ) THEN 
     387         CALL ctl_warn( ' Balancing increments are only calculated for logchl', & 
     388            &           ' Not assimilating logchl, so ln_balwri will be set to .false.') 
     389         ln_balwri = .FALSE. 
     390      ENDIF 
    251391 
    252392      IF ( nstop > 0 ) RETURN       ! if there are any errors then go no further 
    … …  
    327467      !-------------------------------------------------------------------- 
    328468 
    329       ALLOCATE( t_bkginc(jpi,jpj,jpk) ) 
    330       ALLOCATE( s_bkginc(jpi,jpj,jpk) ) 
    331       ALLOCATE( u_bkginc(jpi,jpj,jpk) ) 
    332       ALLOCATE( v_bkginc(jpi,jpj,jpk) ) 
    333       ALLOCATE( ssh_bkginc(jpi,jpj)   ) 
    334       ALLOCATE( seaice_bkginc(jpi,jpj)) 
     469      IF ( ln_trainc ) THEN 
     470         ALLOCATE( t_bkginc(jpi,jpj,jpk) ) 
     471         ALLOCATE( s_bkginc(jpi,jpj,jpk) ) 
     472         t_bkginc(:,:,:) = 0.0 
     473         s_bkginc(:,:,:) = 0.0 
     474      ENDIF 
     475      IF ( ln_dyninc ) THEN  
     476         ALLOCATE( u_bkginc(jpi,jpj,jpk) ) 
     477         ALLOCATE( v_bkginc(jpi,jpj,jpk) ) 
     478         u_bkginc(:,:,:) = 0.0 
     479         v_bkginc(:,:,:) = 0.0 
     480      ENDIF 
     481      IF ( ln_sshinc ) THEN 
     482         ALLOCATE( ssh_bkginc(jpi,jpj)   ) 
     483         ssh_bkginc(:,:) = 0.0 
     484      ENDIF 
     485      IF ( ln_seaiceinc ) THEN  
     486         ALLOCATE( seaice_bkginc(jpi,jpj)) 
     487         seaice_bkginc(:,:) = 0.0 
     488      ENDIF 
    335489#if defined key_asminc 
    336490      ALLOCATE( ssh_iau(jpi,jpj)      ) 
    337 #endif 
    338       t_bkginc(:,:,:) = 0.0 
    339       s_bkginc(:,:,:) = 0.0 
    340       u_bkginc(:,:,:) = 0.0 
    341       v_bkginc(:,:,:) = 0.0 
    342       ssh_bkginc(:,:) = 0.0 
    343       seaice_bkginc(:,:) = 0.0 
    344 #if defined key_asminc 
    345491      ssh_iau(:,:)    = 0.0 
    346492#endif 
    347       IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ) ) THEN 
     493      IF ( ( ln_logchltotinc ).OR.( ln_logchlpftinc ) ) THEN 
     494         ALLOCATE( logchl_bkginc(jpi,jpj,nn_asmpfts)) 
     495         logchl_bkginc(:,:,:) = 0.0 
     496#if defined key_top 
     497         ALLOCATE( logchl_balinc(jpi,jpj,jpk,jptra) ) 
     498         logchl_balinc(:,:,:,:) = 0.0 
     499#endif 
     500      ENDIF 
     501      IF ( ( ln_trainc ).OR.( ln_dyninc ).OR.( ln_sshinc ).OR.( ln_seaiceinc ) & 
     502         &  .OR.( ln_logchltotinc ).OR.( ln_logchlpftinc ) ) THEN 
    348503 
    349504         !-------------------------------------------------------------------- 
    … …  
    378533 
    379534         IF ( ln_trainc ) THEN    
    380             CALL iom_get( inum, jpdom_autoglo, 'bckint', t_bkginc, 1 ) 
    381             CALL iom_get( inum, jpdom_autoglo, 'bckins', s_bkginc, 1 ) 
    382             ! Apply the masks 
    383             t_bkginc(:,:,:) = t_bkginc(:,:,:) * tmask(:,:,:) 
    384             s_bkginc(:,:,:) = s_bkginc(:,:,:) * tmask(:,:,:) 
    385             ! Set missing increments to 0.0 rather than 1e+20 
    386             ! to allow for differences in masks 
    387             WHERE( ABS( t_bkginc(:,:,:) ) > 1.0e+10 ) t_bkginc(:,:,:) = 0.0 
    388             WHERE( ABS( s_bkginc(:,:,:) ) > 1.0e+10 ) s_bkginc(:,:,:) = 0.0 
     535             
     536            !Test if the increments file contains the surft variable. 
     537            isurfstat = iom_varid( inum, 'bckinsurft', ldstop = .FALSE. ) 
     538            IF ( isurfstat == -1 ) THEN 
     539               lk_surft = .FALSE. 
     540            ELSE 
     541               lk_surft = .TRUE. 
     542               CALL ctl_warn( ' Applying 2D temperature increment to bottom of ML: ', & 
     543            &                 ' bckinsurft found in increments file.' ) 
     544            ENDIF              
     545 
     546            IF (lk_surft) THEN  
     547                 
     548               ALLOCATE(z_mld(jpi,jpj))  
     549               SELECT CASE(mld_choice)  
     550               CASE(1)  
     551                  z_mld = hmld  
     552               CASE(2)  
     553                  z_mld = hmlp  
     554               CASE(3)  
     555#if defined key_karaml 
     556                  IF ( ln_kara ) THEN 
     557                     z_mld = hmld_kara 
     558                  ELSE 
     559                     CALL ctl_stop("Kara mixed layer not calculated as ln_kara=.false.") 
     560                  ENDIF 
     561#else 
     562                  CALL ctl_stop("Kara mixed layer not defined in current version of NEMO")  ! JW: Safety feature, should be removed 
     563                                                                                            ! once the Kara mixed layer is available  
     564#endif 
     565               CASE(4)  
     566                  z_mld = hmld_tref  
     567               END SELECT        
     568                       
     569               ALLOCATE( t_bkginc_2d(jpi,jpj) )  
     570               CALL iom_get( inum, jpdom_autoglo, 'bckinsurft', t_bkginc_2d, 1)  
     571#if defined key_bdy                 
     572               DO jk = 1,jpkm1  
     573                  WHERE( z_mld(:,:) > fsdepw(:,:,jk) )  
     574                     t_bkginc(:,:,jk) = t_bkginc_2d(:,:) * 0.5 * & 
     575                     &       ( 1 + cos( (fsdept(:,:,jk)/z_mld(:,:) ) * rpi ) ) 
     576                      
     577                     t_bkginc(:,:,jk) = t_bkginc(:,:,jk) * bdytmask(:,:)  
     578                  ELSEWHERE  
     579                     t_bkginc(:,:,jk) = 0.  
     580                  ENDWHERE  
     581               ENDDO  
     582#else  
     583               t_bkginc(:,:,:) = 0.  
     584#endif                 
     585               s_bkginc(:,:,:) = 0.  
     586                 
     587               DEALLOCATE(z_mld, t_bkginc_2d)  
     588             
     589            ELSE  
     590                
     591               CALL iom_get( inum, jpdom_autoglo, 'bckint', t_bkginc, 1 ) 
     592               CALL iom_get( inum, jpdom_autoglo, 'bckins', s_bkginc, 1 ) 
     593               ! Apply the masks 
     594               t_bkginc(:,:,:) = t_bkginc(:,:,:) * tmask(:,:,:) 
     595               s_bkginc(:,:,:) = s_bkginc(:,:,:) * tmask(:,:,:) 
     596               ! Set missing increments to 0.0 rather than 1e+20 
     597               ! to allow for differences in masks 
     598               WHERE( ABS( t_bkginc(:,:,:) ) > 1.0e+10 ) t_bkginc(:,:,:) = 0.0 
     599               WHERE( ABS( s_bkginc(:,:,:) ) > 1.0e+10 ) s_bkginc(:,:,:) = 0.0 
     600          
     601            ENDIF 
     602          
    389603         ENDIF 
    390604 
    … …  
    417631            ! to allow for differences in masks 
    418632            WHERE( ABS( seaice_bkginc(:,:) ) > 1.0e+10 ) seaice_bkginc(:,:) = 0.0 
     633         ENDIF 
     634 
     635         IF ( ln_logchltotinc ) THEN 
     636            CALL iom_get( inum, jpdom_autoglo, 'bckinlogchl', logchl_bkginc(:,:,1), 1 ) 
     637            ! Apply the masks 
     638            logchl_bkginc(:,:,1) = logchl_bkginc(:,:,1) * tmask(:,:,1) 
     639            ! Set missing increments to 0.0 rather than 1e+20 
     640            ! to allow for differences in masks 
     641            WHERE( ABS( logchl_bkginc(:,:,1) ) > 1.0e+10 ) logchl_bkginc(:,:,1) = 0.0 
     642         ENDIF 
     643 
     644         IF ( ln_logchlpftinc ) THEN 
     645            DO jt = 1, nn_asmpfts 
     646               WRITE(cl_pftstr,'(I2.2)') jt 
     647               CALL iom_get( inum, jpdom_autoglo, 'bckinlogchl'//cl_pftstr, logchl_bkginc(:,:,jt), 1 ) 
     648               ! Apply the masks 
     649               logchl_bkginc(:,:,jt) = logchl_bkginc(:,:,jt) * tmask(:,:,1) 
     650               ! Set missing increments to 0.0 rather than 1e+20 
     651               ! to allow for differences in masks 
     652               WHERE( ABS( logchl_bkginc(:,:,jt) ) > 1.0e+10 ) logchl_bkginc(:,:,jt) = 0.0 
     653            END DO 
    419654         ENDIF 
    420655 
    … …  
    11461381 
    11471382   END SUBROUTINE seaice_asm_inc 
     1383 
     1384   SUBROUTINE logchl_asm_inc( kt ) 
     1385      !!---------------------------------------------------------------------- 
     1386      !!                    ***  ROUTINE logchl_asm_inc  *** 
     1387      !!           
     1388      !! ** Purpose : Apply the chlorophyll assimilation increments. 
     1389      !! 
     1390      !! ** Method  : Calculate increments to state variables using nitrogen 
     1391      !!              balancing. 
     1392      !!              Direct initialization or Incremental Analysis Updating. 
     1393      !! 
     1394      !! ** Action  :  
     1395      !!---------------------------------------------------------------------- 
     1396      INTEGER, INTENT(IN) :: kt   ! Current time step 
     1397      ! 
     1398      INTEGER  :: jk              ! Loop counter 
     1399      INTEGER  :: it              ! Index 
     1400      REAL(wp) :: zincwgt         ! IAU weight for current time step 
     1401      REAL(wp) :: zincper         ! IAU interval in seconds 
     1402      !!---------------------------------------------------------------------- 
     1403 
     1404      IF ( kt <= nit000 ) THEN 
     1405 
     1406         zincper = (nitiaufin_r - nitiaustr_r + 1) * rdt 
     1407 
     1408#if defined key_fabm 
     1409         CALL asm_logchl_bal_ersem( ln_logchlpftinc, nn_asmpfts, mld_choice_bgc, & 
     1410            &                       logchl_bkginc, logchl_balinc ) 
     1411#elif defined key_medusa && defined key_foam_medusa 
     1412         !CALL asm_logchl_bal_medusa() 
     1413         CALL ctl_stop( 'Attempting to assimilate logchl into MEDUSA, ', & 
     1414            &           'but not fully implemented yet' ) 
     1415#elif defined key_hadocc 
     1416         !CALL asm_logchl_bal_hadocc() 
     1417         CALL ctl_stop( 'Attempting to assimilate logchl into HadOCC, ', & 
     1418            &           'but not fully implemented yet' ) 
     1419#else 
     1420         CALL ctl_stop( 'Attempting to assimilate logchl, ', & 
     1421            &           'but not defined a biogeochemical model' ) 
     1422#endif 
     1423 
     1424      ENDIF 
     1425 
     1426      IF ( ln_asmiau ) THEN 
     1427 
     1428         !-------------------------------------------------------------------- 
     1429         ! Incremental Analysis Updating 
     1430         !-------------------------------------------------------------------- 
     1431 
     1432         IF ( ( kt >= nitiaustr_r ).AND.( kt <= nitiaufin_r ) ) THEN 
     1433 
     1434            it = kt - nit000 + 1 
     1435            zincwgt = wgtiau(it)      ! IAU weight for the current time step 
     1436            ! note this is not a tendency so should not be divided by rdt 
     1437 
     1438            IF(lwp) THEN 
     1439               WRITE(numout,*)  
     1440               WRITE(numout,*) 'logchl_asm_inc : logchl IAU at time step = ', & 
     1441                  &  kt,' with IAU weight = ', wgtiau(it) 
     1442               WRITE(numout,*) '~~~~~~~~~~~~' 
     1443            ENDIF 
     1444 
     1445            ! Update the biogeochemical variables 
     1446            ! Add directly to trn and trb, rather than to tra, as not a tendency 
     1447#if defined key_fabm 
     1448            DO jk = 1, jpkm1 
     1449               trn(:,:,jk,jp_fabm0:jp_fabm1) = trn(:,:,jk,jp_fabm0:jp_fabm1) + & 
     1450                  &                            logchl_balinc(:,:,jk,jp_fabm0:jp_fabm1) * zincwgt 
     1451               trb(:,:,jk,jp_fabm0:jp_fabm1) = trb(:,:,jk,jp_fabm0:jp_fabm1) + & 
     1452                  &                            logchl_balinc(:,:,jk,jp_fabm0:jp_fabm1) * zincwgt 
     1453            END DO 
     1454#elif defined key_medusa && defined key_foam_medusa 
     1455            DO jk = 1, jpkm1 
     1456               trn(:,:,jk,jp_msa0:jp_msa1) = trn(:,:,jk,jp_msa0:jp_msa1) + & 
     1457                  &                          logchl_balinc(:,:,jk,jp_msa0:jp_msa1) * zincwgt 
     1458               trb(:,:,jk,jp_msa0:jp_msa1) = trb(:,:,jk,jp_msa0:jp_msa1) + & 
     1459                  &                          logchl_balinc(:,:,jk,jp_msa0:jp_msa1) * zincwgt 
     1460            END DO 
     1461#elif defined key_hadocc 
     1462            DO jk = 1, jpkm1 
     1463               trn(:,:,jk,jp_had0:jp_had1) = trn(:,:,jk,jp_had0:jp_had1) + & 
     1464                  &                          logchl_balinc(:,:,jk,jp_had0:jp_had1) * zincwgt 
     1465               trb(:,:,jk,jp_had0:jp_had1) = trb(:,:,jk,jp_had0:jp_had1) + & 
     1466                  &                          logchl_balinc(:,:,jk,jp_had0:jp_had1) * zincwgt 
     1467            END DO 
     1468#endif 
     1469            
     1470            ! Do not deallocate arrays - needed by asm_bal_wri 
     1471            ! which is called at end of model run 
     1472 
     1473         ENDIF 
     1474 
     1475      ELSEIF ( ln_asmdin ) THEN  
     1476 
     1477         !-------------------------------------------------------------------- 
     1478         ! Direct Initialization 
     1479         !-------------------------------------------------------------------- 
     1480          
     1481         IF ( kt == nitdin_r ) THEN 
     1482 
     1483            neuler = 0                    ! Force Euler forward step 
     1484 
     1485#if defined key_fabm 
     1486            ! Initialize the now fields with the background + increment 
     1487            ! Background currently is what the model is initialised with 
     1488            CALL ctl_warn( ' Doing direct initialisation of ERSEM with chlorophyll assimilation', & 
     1489               &           ' Background state is taken from model rather than background file' ) 
     1490            trn(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) + & 
     1491               &                           logchl_balinc(:,:,:,jp_fabm0:jp_fabm1) 
     1492            trb(:,:,:,jp_fabm0:jp_fabm1) = trn(:,:,:,jp_fabm0:jp_fabm1) 
     1493#elif defined key_medusa && defined key_foam_medusa 
     1494            ! Initialize the now fields with the background + increment 
     1495            ! Background currently is what the model is initialised with 
     1496            CALL ctl_warn( ' Doing direct initialisation of MEDUSA with chlorophyll assimilation', & 
     1497               &           ' Background state is taken from model rather than background file' ) 
     1498            trn(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) + & 
     1499               &                         logchl_balinc(:,:,:,jp_msa0:jp_msa1) 
     1500            trb(:,:,:,jp_msa0:jp_msa1) = trn(:,:,:,jp_msa0:jp_msa1) 
     1501#elif defined key_hadocc 
     1502            ! Initialize the now fields with the background + increment 
     1503            ! Background currently is what the model is initialised with 
     1504            CALL ctl_warn( ' Doing direct initialisation of HadOCC with chlorophyll assimilation', & 
     1505               &           ' Background state is taken from model rather than background file' ) 
     1506            trn(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) + & 
     1507               &                         logchl_balinc(:,:,:,jp_had0:jp_had1) 
     1508            trb(:,:,:,jp_had0:jp_had1) = trn(:,:,:,jp_had0:jp_had1) 
     1509#endif 
     1510  
     1511            ! Do not deallocate arrays - needed by asm_bal_wri 
     1512            ! which is called at end of model run 
     1513         ENDIF 
     1514         ! 
     1515      ENDIF 
     1516      ! 
     1517   END SUBROUTINE logchl_asm_inc 
    11481518    
    11491519   !!====================================================================== 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ASM/asmpar.F90

    r7730 r7731  
    2020      & c_asmtrj = 'assim_trj',                  & !: Filename for storing the  
    2121                                                   !: reference trajectory 
    22       & c_asminc = 'assim_background_increments'  !: Filename for storing the  
     22      & c_asminc = 'assim_background_increments', & !: Filename for storing the  
    2323                                                   !: increments to the background 
    2424                                                   !: state 
     25      & c_asmbal = 'assim.balincs'                 !: Filename for storing the  
     26                                                   !: balancing increments calculated 
     27                                                   !: for biogeochemistry 
    2528 
    2629   INTEGER, PUBLIC :: nitbkg_r      !: Background time step referenced to nit000 
    … …  
    2932   INTEGER, PUBLIC :: nitiaufin_r   !: IAU final time step referenced to nit000 
    3033   INTEGER, PUBLIC :: nittrjfrq     !: Frequency of trajectory output for 4D-VAR 
     34   INTEGER, PUBLIC :: nitavgbkg_r   !: Averaging period for assim bkg referenced to nit000 
    3135 
    3236   !!---------------------------------------------------------------------- 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90

    r7730 r7731  
    430430      CHARACTER(len=100)                     ::   cn_dir        ! Root directory for location of data files 
    431431      CHARACTER(len=100), DIMENSION(nb_bdy)  ::   cn_dir_array  ! Root directory for location of data files 
     432      CHARACTER(len = 256)::   clname                           ! temporary file name 
    432433      LOGICAL                                ::   ln_full_vel   ! =T => full velocities in 3D boundary data 
    433434                                                                ! =F => baroclinic velocities in 3D boundary data 
    … …  
    669670            ! sea ice 
    670671            IF( nn_ice_lim_dta(ib_bdy) .eq. 1 ) THEN 
    671  
    672                ! Test for types of ice input (lim2 or lim3)  
    673                CALL iom_open ( bn_a_i%clname, inum ) 
    674                id1 = iom_varid ( inum, bn_a_i%clvar, kdimsz=zdimsz, kndims=zndims, ldstop = .FALSE. ) 
     672               ! Test for types of ice input (lim2 or lim3) 
     673               ! Build file name to find dimensions  
     674               clname=TRIM(bn_a_i%clname) 
     675               IF( .NOT. bn_a_i%ln_clim ) THEN    
     676                                                  WRITE(clname, '(a,"_y",i4.4)' ) TRIM( bn_a_i%clname ), nyear    ! add year 
     677                  IF( bn_a_i%cltype /= 'yearly' ) WRITE(clname, '(a,"m" ,i2.2)' ) TRIM( clname        ), nmonth   ! add month 
     678               ELSE 
     679                  IF( bn_a_i%cltype /= 'yearly' ) WRITE(clname, '(a,"_m",i2.2)' ) TRIM( bn_a_i%clname ), nmonth   ! add month 
     680               ENDIF 
     681               IF( bn_a_i%cltype == 'daily' .OR. bn_a_i%cltype(1:4) == 'week' ) & 
     682               &                                  WRITE(clname, '(a,"d" ,i2.2)' ) TRIM( clname        ), nday     ! add day 
     683               ! 
     684               CALL iom_open  ( clname, inum ) 
     685               id1 = iom_varid( inum, bn_a_i%clvar, kdimsz=zdimsz, kndims=zndims, ldstop = .FALSE. ) 
    675686               CALL iom_close ( inum ) 
    676                !CALL fld_clopn ( bn_a_i, nyear, nmonth, nday, ldstop=.TRUE. ) 
    677                !CALL iom_open ( bn_a_i%clname, inum ) 
    678                !id1 = iom_varid ( bn_a_i%num, bn_a_i%clvar, kdimsz=zdimsz, kndims=zndims, ldstop = .FALSE. ) 
     687 
    679688                IF ( zndims == 4 ) THEN 
    680689                 ll_bdylim3 = .TRUE.   ! lim3 input 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/BDY/bdydyn2d.F90

    r7730 r7731  
    4949      !!---------------------------------------------------------------------- 
    5050      INTEGER,                      INTENT(in) ::   kt   ! Main time step counter 
    51       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d  
    52       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pub2d, pvb2d 
    53       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: phur, phvr 
    54       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pssh 
     51      REAL(wp), DIMENSION(:,:), INTENT(inout) :: pua2d, pva2d  
     52      REAL(wp), DIMENSION(:,:), INTENT(in   ) :: pub2d, pvb2d 
     53      REAL(wp), DIMENSION(:,:), INTENT(in   ) :: phur, phvr 
     54      REAL(wp), DIMENSION(:,:), INTENT(in   ) :: pssh 
    5555      !! 
    5656      INTEGER                                  ::   ib_bdy ! Loop counter 
    … …  
    9292      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data 
    9393      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index 
    94       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d  
     94      REAL(wp), DIMENSION(:,:), INTENT(inout) :: pua2d, pva2d  
    9595      !! 
    9696      INTEGER  ::   jb, jk         ! dummy loop indices 
    … …  
    147147      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data 
    148148      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index 
    149       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d 
    150       REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh, phur, phvr  
     149      REAL(wp), DIMENSION(:,:), INTENT(inout) :: pua2d, pva2d 
     150      REAL(wp), DIMENSION(:,:), INTENT(in) ::   pssh, phur, phvr  
    151151 
    152152      INTEGER  ::   jb, igrd                         ! dummy loop indices 
    … …  
    237237      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data 
    238238      INTEGER,                      INTENT(in) ::   ib_bdy  ! number of current open boundary set 
    239       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d 
    240       REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pub2d, pvb2d  
     239      REAL(wp), DIMENSION(:,:),    INTENT(inout) :: pua2d, pva2d 
     240      REAL(wp), DIMENSION(:,:),    INTENT(in) :: pub2d, pvb2d  
    241241      LOGICAL,                      INTENT(in) ::   ll_npo  ! flag for NPO version 
    242242 
    … …  
    271271      !! 
    272272      !!---------------------------------------------------------------------- 
    273       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) ::   zssh ! Sea level 
     273      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zssh ! Sea level 
    274274      !! 
    275275      INTEGER  ::   ib_bdy, ib, igrd                        ! local integers 
    276       INTEGER  ::   ii, ij, zcoef, zcoef1, zcoef2, ip, jp   !   "       " 
     276      INTEGER  ::   ii, ij, zcoef, ip, jp   !   "       " 
    277277 
    278278      igrd = 1                       ! Everything is at T-points here 
    … …  
    283283            ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
    284284            ! Set gradient direction: 
    285             zcoef1 = bdytmask(ii-1,ij  ) +  bdytmask(ii+1,ij  ) 
    286             zcoef2 = bdytmask(ii  ,ij-1) +  bdytmask(ii  ,ij+1) 
    287             IF ( zcoef1+zcoef2 == 0 ) THEN 
    288                ! corner 
    289 !               zcoef = tmask(ii-1,ij,1) + tmask(ii+1,ij,1) +  tmask(ii,ij-1,1) +  tmask(ii,ij+1,1) 
    290 !               zssh(ii,ij) = zssh(ii-1,ij  ) * tmask(ii-1,ij  ,1) + & 
    291 !                 &           zssh(ii+1,ij  ) * tmask(ii+1,ij  ,1) + & 
    292 !                 &           zssh(ii  ,ij-1) * tmask(ii  ,ij-1,1) + & 
    293 !                 &           zssh(ii  ,ij+1) * tmask(ii  ,ij+1,1) 
    294                zcoef = bdytmask(ii-1,ij) + bdytmask(ii+1,ij) +  bdytmask(ii,ij-1) +  bdytmask(ii,ij+1) 
    295                zssh(ii,ij) = zssh(ii-1,ij  ) * bdytmask(ii-1,ij  ) + & 
    296                  &           zssh(ii+1,ij  ) * bdytmask(ii+1,ij  ) + & 
    297                  &           zssh(ii  ,ij-1) * bdytmask(ii  ,ij-1) + & 
    298                  &           zssh(ii  ,ij+1) * bdytmask(ii  ,ij+1) 
    299                zssh(ii,ij) = ( zssh(ii,ij) / MAX( 1, zcoef) ) * tmask(ii,ij,1) 
     285            zcoef = bdytmask(ii-1,ij) + bdytmask(ii+1,ij) +  bdytmask(ii,ij-1) +  bdytmask(ii,ij+1) 
     286            IF ( zcoef == 0 ) THEN 
     287               zssh(ii,ij) = 0._wp 
    300288            ELSE 
    301289               ip = bdytmask(ii+1,ij  ) - bdytmask(ii-1,ij  ) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90

    r7730 r7731  
    107107      REAL(wp) ::   zwgt, zwgt1        ! local scalar 
    108108      REAL(wp) ::   ztmelts, zdh 
     109#if  defined key_lim2 && ! defined key_lim2_vp && defined key_agrif 
     110     USE ice_2, vt_s => hsnm 
     111     USE ice_2, vt_i => hicm 
     112#endif 
    109113 
    110114      !!------------------------------------------------------------------------------ 
    … …  
    115119      ! 
    116120#if defined key_lim2 
    117       DO jb = 1, idx%nblen(jgrd) 
     121      DO jb = 1, idx%nblenrim(jgrd) 
    118122         ji    = idx%nbi(jb,jgrd) 
    119123         jj    = idx%nbj(jb,jgrd) 
    … …  
    135139 
    136140      DO jl = 1, jpl 
    137          DO jb = 1, idx%nblen(jgrd) 
     141         DO jb = 1, idx%nblenrim(jgrd) 
    138142            ji    = idx%nbi(jb,jgrd) 
    139143            jj    = idx%nbj(jb,jgrd) 
    … …  
    171175 
    172176      DO jl = 1, jpl 
    173          DO jb = 1, idx%nblen(jgrd) 
     177         DO jb = 1, idx%nblenrim(jgrd) 
    174178            ji    = idx%nbi(jb,jgrd) 
    175179            jj    = idx%nbj(jb,jgrd) 
    … …  
    324328                
    325329               jgrd = 2      ! u velocity 
    326                DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 
     330               DO jb = 1, idx_bdy(ib_bdy)%nblenrim(jgrd) 
    327331                  ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd) 
    328332                  jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd) 
    … …  
    353357                
    354358               jgrd = 3      ! v velocity 
    355                DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 
     359               DO jb = 1, idx_bdy(ib_bdy)%nblenrim(jgrd) 
    356360                  ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd) 
    357361                  jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90

    r7730 r7731  
    7676      INTEGER  ::   ib_bdy, ii, ij, ik, igrd, ib, ir, iseg ! dummy loop indices 
    7777      INTEGER  ::   icount, icountr, ibr_max, ilen1, ibm1  ! local integers 
    78       INTEGER  ::   iw, ie, is, in, inum, id_dummy         !   -       - 
     78      INTEGER  ::   iwe, ies, iso, ino, inum, id_dummy         !   -       - 
    7979      INTEGER  ::   igrd_start, igrd_end, jpbdta           !   -       - 
    8080      INTEGER  ::   jpbdtau, jpbdtas                       !   -       - 
    … …  
    777777!      is = mjg(1) + 1            ! if monotasking and no zoom, is=2 
    778778!      in = mjg(1) + nlcj-1 - 1   ! if monotasking and no zoom, in=jpjm1       
    779       iw = mig(1) - jpizoom + 2         ! if monotasking and no zoom, iw=2 
    780       ie = mig(1) + nlci - jpizoom - 1  ! if monotasking and no zoom, ie=jpim1 
    781       is = mjg(1) - jpjzoom + 2         ! if monotasking and no zoom, is=2 
    782       in = mjg(1) + nlcj - jpjzoom - 1  ! if monotasking and no zoom, in=jpjm1 
     779      iwe = mig(1) - jpizoom + 2         ! if monotasking and no zoom, iw=2 
     780      ies = mig(1) + nlci - jpizoom - 1  ! if monotasking and no zoom, ie=jpim1 
     781      iso = mjg(1) - jpjzoom + 2         ! if monotasking and no zoom, is=2 
     782      ino = mjg(1) + nlcj - jpjzoom - 1  ! if monotasking and no zoom, in=jpjm1 
    783783 
    784784      ALLOCATE( nbondi_bdy(nb_bdy)) 
    … …  
    853853               ENDIF 
    854854               ! check if point is in local domain 
    855                IF(  nbidta(ib,igrd,ib_bdy) >= iw .AND. nbidta(ib,igrd,ib_bdy) <= ie .AND.   & 
    856                   & nbjdta(ib,igrd,ib_bdy) >= is .AND. nbjdta(ib,igrd,ib_bdy) <= in       ) THEN 
     855               IF(  nbidta(ib,igrd,ib_bdy) >= iwe .AND. nbidta(ib,igrd,ib_bdy) <= ies .AND.   & 
     856                  & nbjdta(ib,igrd,ib_bdy) >= iso .AND. nbjdta(ib,igrd,ib_bdy) <= ino      ) THEN 
    857857                  ! 
    858858                  icount = icount  + 1 
    … …  
    890890         com_south_b = 0 
    891891         com_north_b = 0 
     892 
    892893         DO igrd = 1, jpbgrd 
    893894            icount  = 0 
    … …  
    896897               DO ib = 1, nblendta(igrd,ib_bdy) 
    897898                  ! check if point is in local domain and equals ir 
    898                   IF(  nbidta(ib,igrd,ib_bdy) >= iw .AND. nbidta(ib,igrd,ib_bdy) <= ie .AND.   & 
    899                      & nbjdta(ib,igrd,ib_bdy) >= is .AND. nbjdta(ib,igrd,ib_bdy) <= in .AND.   & 
     899                  IF(  nbidta(ib,igrd,ib_bdy) >= iwe .AND. nbidta(ib,igrd,ib_bdy) <= ies .AND.   & 
     900                     & nbjdta(ib,igrd,ib_bdy) >= iso .AND. nbjdta(ib,igrd,ib_bdy) <= ino .AND.   & 
    900901                     & nbrdta(ib,igrd,ib_bdy) == ir  ) THEN 
    901902                     ! 
    … …  
    15941595            ELSE 
    15951596               ! This is a corner 
    1596                WRITE(numout,*) 'Found a South-West corner at (i,j): ', jpiwob(ib), jpjwdt(ib) 
     1597               IF(lwp) WRITE(numout,*) 'Found a South-West corner at (i,j): ', jpiwob(ib), jpjwdt(ib) 
    15971598               CALL bdy_ctl_corn(npckgw(ib), icornw(ib,1)) 
    15981599               itest=itest+1 
    … …  
    16081609            ELSE 
    16091610               ! This is a corner 
    1610                WRITE(numout,*) 'Found a North-West corner at (i,j): ', jpiwob(ib), jpjwft(ib) 
     1611               IF(lwp) WRITE(numout,*) 'Found a North-West corner at (i,j): ', jpiwob(ib), jpjwft(ib) 
    16111612               CALL bdy_ctl_corn(npckgw(ib), icornw(ib,2)) 
    16121613               itest=itest+1 
    … …  
    16381639            ELSE 
    16391640               ! This is a corner 
    1640                WRITE(numout,*) 'Found a South-East corner at (i,j): ', jpieob(ib)+1, jpjedt(ib) 
     1641               IF(lwp) WRITE(numout,*) 'Found a South-East corner at (i,j): ', jpieob(ib)+1, jpjedt(ib) 
    16411642               CALL bdy_ctl_corn(npckge(ib), icorne(ib,1)) 
    16421643               itest=itest+1 
    … …  
    16521653            ELSE 
    16531654               ! This is a corner 
    1654                WRITE(numout,*) 'Found a North-East corner at (i,j): ', jpieob(ib)+1, jpjeft(ib) 
     1655               IF(lwp) WRITE(numout,*) 'Found a North-East corner at (i,j): ', jpieob(ib)+1, jpjeft(ib) 
    16551656               CALL bdy_ctl_corn(npckge(ib), icorne(ib,2)) 
    16561657               itest=itest+1 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/BDY/bdytides.F90

    r7730 r7731  
    416416      ! Absolute time from model initialization:    
    417417      IF( PRESENT(kit) ) THEN   
    418          z_arg = ( kt + (kit+0.5_wp*(time_add-1)) / REAL(nn_baro,wp) ) * rdt 
     418         z_arg = ( kt + (kit+time_add-1) / REAL(nn_baro,wp) ) * rdt 
    419419      ELSE                               
    420420         z_arg = ( kt + time_add ) * rdt 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90

    r7730 r7731  
    9191      ! Calculate the cumulate surface Flux z_cflxemp (m3/s) over all the domain 
    9292      ! ----------------------------------------------------------------------- 
    93       z_cflxemp = SUM ( ( emp(:,:)-rnf(:,:)+rdivisf*fwfisf(:,:) ) * bdytmask(:,:) * e1t(:,:) * e2t(:,:) ) / rau0 
     93      z_cflxemp = SUM ( ( emp(:,:)-rnf(:,:)+fwfisf(:,:) ) * bdytmask(:,:) * e1t(:,:) * e2t(:,:) ) / rau0 
    9494      IF( lk_mpp )   CALL mpp_sum( z_cflxemp )     ! sum over the global domain 
    9595 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DIA/diaharm.F90

    r7730 r7731  
    196196                  DO ji = 1,jpi 
    197197                     ! Elevation 
    198                      ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp*sshn(ji,jj)           *tmask_i(ji,jj)         
    199 #if defined key_dynspg_ts 
    200                      ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp*un_b(ji,jj)*hur(ji,jj)*umask_i(ji,jj) 
    201                      ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp*vn_b(ji,jj)*hvr(ji,jj)*vmask_i(ji,jj) 
    202 #endif 
     198                     ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp*sshn(ji,jj)*tmask_i(ji,jj)         
     199                     ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp*un_b(ji,jj)*umask_i(ji,jj) 
     200                     ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp*vn_b(ji,jj)*vmask_i(ji,jj) 
    203201                  END DO 
    204202               END DO 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90

    r7730 r7731  
    9393      ! 1 - Trends due to forcing ! 
    9494      ! ------------------------- ! 
    95       z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) ) * surf(:,:) ) ! volume fluxes 
     95      z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) * surf(:,:) ) ! volume fluxes 
    9696      z_frc_trd_t =           glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) )                               ! heat fluxes 
    9797      z_frc_trd_s =           glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) )                               ! salt fluxes 
    … …  
    101101      ! Add ice shelf heat & salt input 
    102102      IF( nn_isf .GE. 1 )  THEN 
    103           z_frc_trd_t = z_frc_trd_t & 
    104               &   + glob_sum( ( risf_tsc(:,:,jp_tem) - rdivisf * fwfisf(:,:) * (-1.9) * r1_rau0 ) * surf(:,:) ) 
    105           z_frc_trd_s = z_frc_trd_s + (1.0_wp - rdivisf) * glob_sum( risf_tsc(:,:,jp_sal) * surf(:,:) ) 
     103          z_frc_trd_t = z_frc_trd_t + glob_sum( risf_tsc(:,:,jp_tem) * surf(:,:) ) 
     104          z_frc_trd_s = z_frc_trd_s + glob_sum( risf_tsc(:,:,jp_sal) * surf(:,:) ) 
    106105      ENDIF 
    107106 
    … …  
    200199!      ENDIF 
    201200!!gm end 
    202  
    203201 
    204202      IF( lk_vvl ) THEN 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r7730 r7731  
    438438      zdt = rdt 
    439439      IF( nacc == 1 ) zdt = rdtmin 
    440       IF( ln_mskland )   THEN   ;   clop = "only(x)"   ! put 1.e+20 on land (very expensive!!) 
    441       ELSE                      ;   clop = "x"         ! no use of the mask value (require less cpu time) 
    442       ENDIF 
     440      clop = "x"         ! no use of the mask value (require less cpu time, and otherwise the model crashes) 
    443441#if defined key_diainstant 
    444442      zsto = nwrite * zdt 
    … …  
    10201018         CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      ,   &   ! t-point depth 
    10211019            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
     1020         CALL histdef( id_i, "vovvle3t", "T point thickness"         , "m"      ,   &   ! t-point depth 
     1021            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
    10221022      END IF 
    10231023 
    … …  
    10501050      CALL histwrite( id_i, "sozotaux", kt, utau             , jpi*jpj    , idex )    ! i-wind stress 
    10511051      CALL histwrite( id_i, "sometauy", kt, vtau             , jpi*jpj    , idex )    ! j-wind stress 
     1052      IF( lk_vvl ) THEN 
     1053         CALL histwrite( id_i, "vovvldep", kt, fsdept_n(:,:,:), jpi*jpj*jpk, idex )!  T-cell depth        
     1054         CALL histwrite( id_i, "vovvle3t", kt, fse3t_n (:,:,:), jpi*jpj*jpk, idex )!  T-cell thickness   
     1055      END IF 
    10521056 
    10531057      ! 3. Close the file 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DOM/daymod.F90

    r7730 r7731  
    7373      !!---------------------------------------------------------------------- 
    7474      ! 
     75      ! max number of seconds between each restart 
     76      IF( REAL( nitend - nit000 + 1 ) * rdt > REAL( HUGE( nsec1jan000 ) ) ) THEN 
     77         CALL ctl_stop( 'The number of seconds between each restart exceeds the integer 4 max value: 2^31-1. ',   & 
     78            &           'You must do a restart at higher frequency (or remove this stop and recompile the code in I8)' ) 
     79      ENDIF 
    7580      ! all calendar staff is based on the fact that MOD( rday, rdttra(1) ) == 0 
    7681      IF( MOD( rday     , rdttra(1) ) /= 0. )   CALL ctl_stop( 'the time step must devide the number of second of in a day' ) 
    … …  
    238243               nday_year = 1 
    239244               nsec_year = ndt05 
    240                IF( nsec1jan000 >= 2 * (2**30 - nsecd * nyear_len(1) / 2 ) ) THEN   ! test integer 4 max value 
    241                   CALL ctl_stop( 'The number of seconds between Jan. 1st 00h of nit000 year and Jan. 1st 00h ',   & 
    242                      &           'of the current year is exceeding the INTEGER 4 max VALUE: 2^31-1 -> 68.09 years in seconds', & 
    243                      & 'You must do a restart at higher frequency (or remove this STOP and recompile everything in I8)' ) 
    244                ENDIF 
    245245               nsec1jan000 = nsec1jan000 + nsecd * nyear_len(1) 
    246246               IF( nleapy == 1 )   CALL day_mth 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DOM/domhgr.F90

    r7730 r7731  
    169169            ! 
    170170            ii0 = 282           ;   ii1 = 283        ! Gibraltar Strait (e2u = 20 km) 
    171             ij0 = 201 + isrow   ;   ij1 = 241 - isrow   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  20.e3 
     171            ij0 = 241 - isrow   ;   ij1 = 241 - isrow   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  20.e3 
    172172            IF(lwp) WRITE(numout,*) 
    173173            IF(lwp) WRITE(numout,*) '             orca_r1: Gibraltar : e2u reduced to 20 km' 
    174174 
    175175            ii0 = 314           ;   ii1 = 315        ! Bhosporus Strait (e2u = 10 km) 
    176             ij0 = 208 + isrow   ;   ij1 = 248 - isrow   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
     176            ij0 = 248 - isrow   ;   ij1 = 248 - isrow   ;   e2u( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  10.e3 
    177177            IF(lwp) WRITE(numout,*) 
    178178            IF(lwp) WRITE(numout,*) '             orca_r1: Bhosporus : e2u reduced to 10 km' 
    179179 
    180180            ii0 =  44           ;   ii1 =  44        ! Lombok Strait (e1v = 13 km) 
    181             ij0 = 124 + isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  13.e3 
     181            ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  13.e3 
    182182            IF(lwp) WRITE(numout,*) 
    183183            IF(lwp) WRITE(numout,*) '             orca_r1: Lombok : e1v reduced to 10 km' 
    184184 
    185185            ii0 =  48           ;   ii1 =  48        ! Sumba Strait (e1v = 8 km) [closed from bathy_11 on] 
    186             ij0 = 124 + isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  8.e3 
     186            ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) =  8.e3 
    187187            IF(lwp) WRITE(numout,*) 
    188188            IF(lwp) WRITE(numout,*) '             orca_r1: Sumba : e1v reduced to 8 km' 
    189189 
    190190            ii0 =  53           ;   ii1 =  53        ! Ombai Strait (e1v = 13 km) 
    191             ij0 = 124 + isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3 
     191            ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 13.e3 
    192192            IF(lwp) WRITE(numout,*) 
    193193            IF(lwp) WRITE(numout,*) '             orca_r1: Ombai : e1v reduced to 13 km' 
    194194 
    195195            ii0 =  56           ;   ii1 =  56        ! Timor Passage (e1v = 20 km) 
    196             ij0 = 124 + isrow   ;   ij1 = 145 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3 
     196            ij0 = 164 - isrow   ;   ij1 = 145 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 20.e3 
    197197            IF(lwp) WRITE(numout,*) 
    198198            IF(lwp) WRITE(numout,*) '             orca_r1: Timor Passage : e1v reduced to 20 km' 
    199199 
    200200            ii0 =  55           ;   ii1 =  55        ! West Halmahera Strait (e1v = 30 km) 
    201             ij0 = 141 + isrow   ;   ij1 = 182 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 30.e3 
     201            ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 30.e3 
    202202            IF(lwp) WRITE(numout,*) 
    203203            IF(lwp) WRITE(numout,*) '             orca_r1: W Halmahera : e1v reduced to 30 km' 
    204204 
    205205            ii0 =  58           ;   ii1 =  58        ! East Halmahera Strait (e1v = 50 km) 
    206             ij0 = 141 + isrow   ;   ij1 = 182 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 50.e3 
     206            ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   e1v( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 50.e3 
    207207            IF(lwp) WRITE(numout,*) 
    208208            IF(lwp) WRITE(numout,*) '             orca_r1: E Halmahera : e1v reduced to 50 km' 
    … …  
    544544         IF ( cp_cfg == 'eel' .AND. jp_cfg == 6 ) THEN    ! for EEL6 configuration only 
    545545            IF( .NOT. Agrif_Root() ) THEN 
    546               zphi0 = ppgphi0 - FLOAT( Agrif_Parent(jpjglo)/2)*Agrif_Parent(ppe2_m) / (ra * rad) 
     546              zphi0 = ppgphi0 - FLOAT( Agrif_Parent(jpjglo)/2)*Agrif_Parent(ppe2_m)   &  
     547                    &           / (ra * rad) 
    547548            ENDIF 
    548549         ENDIF 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90

    r7730 r7731  
    413413         IF(lwp) WRITE(numout,*) '      Gibraltar ' 
    414414         ii0 = 282           ;   ii1 = 283        ! Gibraltar Strait  
    415          ij0 = 201 + isrow   ;   ij1 = 241 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
     415         ij0 = 241 - isrow   ;   ij1 = 241 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    416416 
    417417         IF(lwp) WRITE(numout,*) '      Bhosporus ' 
    418418         ii0 = 314           ;   ii1 = 315        ! Bhosporus Strait  
    419          ij0 = 208 + isrow   ;   ij1 = 248 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
     419         ij0 = 248 - isrow   ;   ij1 = 248 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    420420 
    421421         IF(lwp) WRITE(numout,*) '      Makassar (Top) ' 
    422422         ii0 =  48           ;   ii1 =  48        ! Makassar Strait (Top)  
    423          ij0 = 149 + isrow   ;   ij1 = 190 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
     423         ij0 = 189 - isrow   ;   ij1 = 190 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
    424424 
    425425         IF(lwp) WRITE(numout,*) '      Lombok ' 
    426426         ii0 =  44           ;   ii1 =  44        ! Lombok Strait  
    427          ij0 = 124 + isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
     427         ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    428428 
    429429         IF(lwp) WRITE(numout,*) '      Ombai ' 
    430430         ii0 =  53           ;   ii1 =  53        ! Ombai Strait  
    431          ij0 = 124 + isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
     431         ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    432432 
    433433         IF(lwp) WRITE(numout,*) '      Timor Passage ' 
    434434         ii0 =  56           ;   ii1 =  56        ! Timor Passage  
    435          ij0 = 124 + isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
     435         ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    436436 
    437437         IF(lwp) WRITE(numout,*) '      West Halmahera ' 
    438438         ii0 =  58           ;   ii1 =  58        ! West Halmahera Strait  
    439          ij0 = 141 + isrow   ;   ij1 = 182 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
     439         ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
    440440 
    441441         IF(lwp) WRITE(numout,*) '      East Halmahera ' 
    442442         ii0 =  55           ;   ii1 =  55        ! East Halmahera Strait  
    443          ij0 = 141 + isrow   ;   ij1 = 182 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
     443         ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
    444444         ! 
    445445      ENDIF 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90

    r7730 r7731  
    215215         CALL iom_rstput( 0, 0, inum4, 'gdept_1d' , gdept_1d )  !    ! stretched system 
    216216         CALL iom_rstput( 0, 0, inum4, 'gdepw_1d' , gdepw_1d ) 
     217         CALL iom_rstput( 0, 0, inum4, 'gdept_0', gdept_0, ktype = jp_r4 )      
     218         CALL iom_rstput( 0, 0, inum4, 'gdepw_0', gdepw_0, ktype = jp_r4 )      
    217219      ENDIF 
    218220       

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r7730 r7731  
    219219         &  ppsur == pp_to_be_computed           ) THEN 
    220220         ! 
     221#if defined key_agrif 
     222         za1  = (  ppdzmin - pphmax / FLOAT(jpkdta-1)  )                                                   & 
     223            & / ( TANH((1-ppkth)/ppacr) - ppacr/FLOAT(jpkdta-1) * (  LOG( COSH( (jpkdta - ppkth) / ppacr) )& 
     224            &                                                      - LOG( COSH( ( 1  - ppkth) / ppacr) )  )  ) 
     225#else 
    221226         za1  = (  ppdzmin - pphmax / FLOAT(jpkm1)  )                                                      & 
    222227            & / ( TANH((1-ppkth)/ppacr) - ppacr/FLOAT(jpk-1) * (  LOG( COSH( (jpk - ppkth) / ppacr) )      & 
    223228            &                                                   - LOG( COSH( ( 1  - ppkth) / ppacr) )  )  ) 
     229#endif 
    224230         za0  = ppdzmin - za1 *              TANH( (1-ppkth) / ppacr ) 
    225231         zsur =   - za0 - za1 * ppacr * LOG( COSH( (1-ppkth) / ppacr )  ) 
    … …  
    236242              WRITE(numout,*) '            Uniform grid with ',jpk-1,' layers' 
    237243              WRITE(numout,*) '            Total depth    :', zhmax 
     244#if defined key_agrif 
     245              WRITE(numout,*) '            Layer thickness:', zhmax/(jpkdta-1) 
     246#else 
    238247              WRITE(numout,*) '            Layer thickness:', zhmax/(jpk-1) 
     248#endif 
    239249         ELSE 
    240250            IF( ppa1 == 0._wp .AND. ppa0 == 0._wp .AND. ppsur == 0._wp ) THEN 
    … …  
    260270      ! Reference z-coordinate (depth - scale factor at T- and W-points) 
    261271      ! ====================== 
    262       IF( ppkth == 0._wp ) THEN            !  uniform vertical grid        
     272      IF( ppkth == 0._wp ) THEN            !  uniform vertical grid  
     273#if defined key_agrif 
     274         za1 = zhmax / FLOAT(jpkdta-1)  
     275#else 
    263276         za1 = zhmax / FLOAT(jpk-1)  
     277#endif 
    264278         DO jk = 1, jpk 
    265279            zw = FLOAT( jk ) 
    … …  
    18701884             iim1 = MAX( ji-1, 1 ) 
    18711885             ijm1 = MAX( jj-1, 1 ) 
    1872              IF( (bathy(iip1,jj) + bathy(iim1,jj) + bathy(ji,ijp1) + bathy(ji,ijm1) +              & 
    1873         &         bathy(iip1,ijp1) + bathy(iim1,ijm1) + bathy(iip1,ijp1) + bathy(iim1,ijm1)) > 0._wp ) THEN 
    1874                zenv(ji,jj) = rn_sbot_min 
     1886             IF( ( + bathy(iim1,ijp1) + bathy(ji,ijp1) + bathy(iip1,ijp1)  & 
     1887                &  + bathy(iim1,jj  )                  + bathy(iip1,jj  )  & 
     1888                &  + bathy(iim1,ijm1) + bathy(ji,ijm1) + bathy(iip1,ijp1)  ) > 0._wp ) THEN 
     1889                zenv(ji,jj) = rn_sbot_min 
    18751890             ENDIF 
    18761891           ENDIF 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DYN/divcur.F90

    r7730 r7731  
    9797      IF( nn_timing == 1 )  CALL timing_start('div_cur') 
    9898      ! 
    99       CALL wrk_alloc( jpi  , jpj+2, zwu               ) 
    100       CALL wrk_alloc( jpi+4, jpj  , zwv, kistart = -1 ) 
     99      CALL wrk_alloc( jpi  , jpj+2, zwu  ) 
     100      CALL wrk_alloc( jpi+2, jpj  , zwv ) 
    101101      ! 
    102102      IF( kt == nit000 ) THEN 
    … …  
    236236      CALL lbc_lnk( hdivn, 'T', 1. )   ;   CALL lbc_lnk( rotn , 'F', 1. )    ! lateral boundary cond. (no sign change) 
    237237      ! 
    238       CALL wrk_dealloc( jpi  , jpj+2, zwu               ) 
    239       CALL wrk_dealloc( jpi+4, jpj  , zwv, kistart = -1 ) 
     238      CALL wrk_dealloc( jpi  , jpj+2, zwu ) 
     239      CALL wrk_dealloc( jpi+2, jpj  , zwv ) 
    240240      ! 
    241241      IF( nn_timing == 1 )  CALL timing_stop('div_cur') 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90

    r7730 r7731  
    266266               ! Add volume filter correction: compatibility with tracer advection scheme 
    267267               ! => time filter + conservation correction (only at the first level) 
    268                fse3t_b(:,:,1) = fse3t_b(:,:,1) - atfp * rdt * r1_rau0 * ( emp_b(:,:) - emp(:,:) & 
    269                               &                                          -rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1) 
     268               IF ( nn_isf == 0) THEN   ! if no ice shelf melting 
     269                  fse3t_b(:,:,1) = fse3t_b(:,:,1) - atfp * rdt * r1_rau0 * ( emp_b(:,:) - emp(:,:) & 
     270                                 &                                          -rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1) 
     271               ELSE                     ! if ice shelf melting 
     272                  DO jj = 1,jpj 
     273                     DO ji = 1,jpi 
     274                        jk = mikt(ji,jj) 
     275                        fse3t_b(ji,jj,jk) = fse3t_b(ji,jj,jk) - atfp * rdt * r1_rau0                       & 
     276                                          &                          * ( (emp_b(ji,jj)    - emp(ji,jj)   ) & 
     277                                          &                            - (rnf_b(ji,jj)    - rnf(ji,jj)   ) & 
     278                                          &                            + (fwfisf_b(ji,jj) - fwfisf(ji,jj)) ) * tmask(ji,jj,jk) 
     279                     END DO 
     280                  END DO 
     281               END IF 
    270282            ENDIF 
    271283            ! 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r7730 r7731  
    4545   USE agrif_opa_interp ! agrif 
    4646#endif 
    47 #if defined key_asminc    
    48    USE asminc          ! Assimilation increment 
    49 #endif 
    5047 
    5148   IMPLICIT NONE 
    … …  
    187184      ! 
    188185                                                       ! time offset in steps for bdy data update 
    189       IF (.NOT.ln_bt_fw) THEN ; noffset=-2*nn_baro ; ELSE ;  noffset = 0 ; ENDIF 
     186      IF (.NOT.ln_bt_fw) THEN ; noffset=-nn_baro ; ELSE ;  noffset = 0 ; ENDIF 
    190187      ! 
    191188      IF( kt == nit000 ) THEN                !* initialisation 
    … …  
    454451      !                                         ! Surface net water flux and rivers 
    455452      IF (ln_bt_fw) THEN 
    456          zssh_frc(:,:) = zraur * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) ) 
     453         zssh_frc(:,:) = zraur * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) 
    457454      ELSE 
    458455         zssh_frc(:,:) = zraur * z1_2 * (  emp(:,:) + emp_b(:,:) - rnf(:,:) - rnf_b(:,:)   & 
    459                 &                        + rdivisf * ( fwfisf(:,:) + fwfisf_b(:,:) )       ) 
    460       ENDIF 
    461 #if defined key_asminc 
    462       !                                         ! Include the IAU weighted SSH increment 
    463       IF( lk_asminc .AND. ln_sshinc .AND. ln_asmiau ) THEN 
    464          zssh_frc(:,:) = zssh_frc(:,:) - ssh_iau(:,:) 
    465       ENDIF 
    466 #endif 
    467       !                                   !* Fill boundary data arrays with AGRIF 
    468       !                                   ! ------------------------------------- 
     456                &                        + fwfisf(:,:) + fwfisf_b(:,:)                     ) 
     457      ENDIF 
     458      !                                   !* Fill boundary data arrays for AGRIF 
     459      !                                   ! ------------------------------------ 
    469460#if defined key_agrif 
    470461         IF( .NOT.Agrif_Root() ) CALL agrif_dta_ts( kt ) 
    … …  
    523514         ! Update only tidal forcing at open boundaries 
    524515#if defined key_tide 
    525          IF ( lk_bdy .AND. lk_tide )      CALL bdy_dta_tides( kt, kit=jn, time_offset=(noffset+1) ) 
    526          IF ( ln_tide_pot .AND. lk_tide ) CALL upd_tide( kt, kit=jn, koffset=noffset ) 
     516         IF ( lk_bdy .AND. lk_tide ) CALL bdy_dta_tides( kt, kit=jn, time_offset=(noffset+1) ) 
     517         IF ( ln_tide_pot .AND. lk_tide ) CALL upd_tide( kt, kit=jn, time_offset=noffset ) 
    527518#endif 
    528519         ! 
    … …  
    900891#if defined key_agrif 
    901892      ! Save time integrated fluxes during child grid integration 
    902       ! (used to update coarse grid transports) 
    903       ! Useless with 2nd order momentum schemes 
     893      ! (used to update coarse grid transports at next time step) 
    904894      ! 
    905895      IF ( (.NOT.Agrif_Root()).AND.(ln_bt_fw) ) THEN 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

    r7730 r7731  
    3131   USE bdydyn2d        ! bdy_ssh routine 
    3232#if defined key_agrif 
    33    USE agrif_opa_update 
    3433   USE agrif_opa_interp 
    3534#endif 
    … …  
    7574      INTEGER, INTENT(in) ::   kt                      ! time step 
    7675      !  
    77       INTEGER             ::   jk                      ! dummy loop indice 
     76      INTEGER             ::   ji, jj, jk              ! dummy loop indices 
    7877      REAL(wp)            ::   z2dt, z1_rau0           ! local scalars 
    7978      !!---------------------------------------------------------------------- 
    … …  
    9594      z2dt = 2._wp * rdt                              ! set time step size (Euler/Leapfrog) 
    9695      IF( neuler == 0 .AND. kt == nit000 )   z2dt = rdt 
     96 
     97#if defined key_asminc 
     98      !                                                ! Include the IAU weighted SSH increment 
     99      IF( lk_asminc .AND. ln_sshinc .AND. ln_asmiau ) THEN 
     100         CALL ssh_asm_inc( kt ) 
     101#if defined key_vvl 
     102! Don't directly adjust ssh but change hdivn at all levels instead 
     103! In trasbc also add in the heat and salt content associated with these changes at each level   
     104        DO jk = 1, jpkm1                                  
     105                 hdivn(:,:,jk) = hdivn(:,:,jk) - ( ssh_iau(:,:) / ( ht_0(:,:) + 1.0 - ssmask(:,:) ) ) * ( e3t_0(:,:,jk) / fse3t_n(:,:,jk) ) * tmask(:,:,jk)  
     106        END DO 
     107        CALL lbc_lnk( hdivn, 'T', 1. ) ! Not sure that's necessary 
     108      ENDIF 
     109#endif 
     110#endif 
    97111 
    98112      !                                           !------------------------------! 
    … …  
    124138#endif 
    125139 
    126 #if defined key_asminc 
    127       !                                                ! Include the IAU weighted SSH increment 
    128       IF( lk_asminc .AND. ln_sshinc .AND. ln_asmiau ) THEN 
    129          CALL ssh_asm_inc( kt ) 
    130          ssha(:,:) = ssha(:,:) + z2dt * ssh_iau(:,:) 
    131       ENDIF 
    132 #endif 
    133  
    134140      !                                           !------------------------------! 
    135141      !                                           !           outputs            ! 
    … …  
    268274      ELSE                                         !** Leap-Frog time-stepping: Asselin filter + swap 
    269275         sshb(:,:) = sshn(:,:) + atfp * ( sshb(:,:) - 2 * sshn(:,:) + ssha(:,:) )     ! before <-- now filtered 
    270          IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( emp_b(:,:) - emp(:,:) - rnf_b(:,:) + rnf(:,:) ) * ssmask(:,:) 
     276         IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( emp_b(:,:)    - emp(:,:)    & 
     277                                &                                 - rnf_b(:,:)    + rnf(:,:)    & 
     278                                &                                 + fwfisf_b(:,:) - fwfisf(:,:) ) * ssmask(:,:) 
    271279         sshn(:,:) = ssha(:,:)                           ! now <-- after 
    272280      ENDIF 
    273       ! 
    274       ! Update velocity at AGRIF zoom boundaries 
    275 #if defined key_agrif 
    276       IF ( .NOT.Agrif_Root() ) CALL Agrif_Update_Dyn( kt ) 
    277 #endif 
    278281      ! 
    279282      IF(ln_ctl)   CALL prt_ctl( tab2d_1=sshb, clinfo1=' sshb  - : ', mask1=tmask, ovlap=1 ) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/IOM/iom.F90

    r7730 r7731  
    9494      CHARACTER(len=*), INTENT(in)  :: cdname 
    9595#if defined key_iomput 
    96       TYPE(xios_time)   :: dtime    = xios_time(0, 0, 0, 0, 0, 0) 
    97       CHARACTER(len=19) :: cldate  
    98       CHARACTER(len=10) :: clname 
    99       INTEGER           ::   ji 
     96#if ! defined key_xios2 
     97      TYPE(xios_time)     :: dtime    = xios_time(0, 0, 0, 0, 0, 0) 
     98      CHARACTER(len=19)   :: cldate  
     99#else 
     100      TYPE(xios_duration) :: dtime    = xios_duration(0, 0, 0, 0, 0, 0) 
     101      TYPE(xios_date)     :: start_date 
     102#endif 
     103      CHARACTER(len=10)   :: clname 
     104      INTEGER             :: ji 
    100105      ! 
    101106      REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: z_bnds 
    102107      !!---------------------------------------------------------------------- 
    103  
     108#if ! defined key_xios2 
    104109      ALLOCATE( z_bnds(jpk,2) ) 
     110#else 
     111      ALLOCATE( z_bnds(2,jpk) ) 
     112#endif 
    105113 
    106114      clname = cdname 
    … …  
    110118 
    111119      ! calendar parameters 
     120#if ! defined key_xios2 
    112121      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL 
    113122      CASE ( 1)   ;   CALL xios_set_context_attr(TRIM(clname), calendar_type= "Gregorian") 
    … …  
    117126      WRITE(cldate,"(i4.4,'-',i2.2,'-',i2.2,' 00:00:00')") nyear,nmonth,nday  
    118127      CALL xios_set_context_attr(TRIM(clname), start_date=cldate ) 
    119  
     128#else 
     129      ! Calendar type is now defined in xml file  
     130      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL 
     131      CASE ( 1)   ; CALL xios_define_calendar( TYPE = "Gregorian", time_origin = xios_date(1900,01,01,00,00,00), & 
     132          &                                    start_date = xios_date(nyear,nmonth,nday,0,0,0) ) 
     133      CASE ( 0)   ; CALL xios_define_calendar( TYPE = "NoLeap"   , time_origin = xios_date(1900,01,01,00,00,00), & 
     134          &                                    start_date = xios_date(nyear,nmonth,nday,0,0,0) ) 
     135      CASE (30)   ; CALL xios_define_calendar( TYPE = "D360"     , time_origin = xios_date(1900,01,01,00,00,00), & 
     136          &                                    start_date = xios_date(nyear,nmonth,nday,0,0,0) ) 
     137      END SELECT 
     138#endif 
    120139      ! horizontal grid definition 
     140 
     141#if ! defined key_xios2 
    121142      CALL set_scalar 
     143#endif 
    122144 
    123145      IF( TRIM(cdname) == TRIM(cxios_context) ) THEN   
    … …  
    170192 
    171193      ! Add vertical grid bounds 
     194#if ! defined key_xios2 
    172195      z_bnds(:      ,1) = gdepw_1d(:) 
    173196      z_bnds(1:jpkm1,2) = gdepw_1d(2:jpk) 
    174197      z_bnds(jpk:   ,2) = gdepw_1d(jpk) + e3t_1d(jpk) 
     198#else 
     199      z_bnds(1      ,:) = gdepw_1d(:) 
     200      z_bnds(2,1:jpkm1) = gdepw_1d(2:jpk) 
     201      z_bnds(2,jpk:   ) = gdepw_1d(jpk) + e3t_1d(jpk) 
     202#endif 
     203 
    175204      CALL iom_set_axis_attr( "deptht", bounds=z_bnds ) 
    176205      CALL iom_set_axis_attr( "depthu", bounds=z_bnds ) 
    177206      CALL iom_set_axis_attr( "depthv", bounds=z_bnds ) 
    178       z_bnds(:    ,2) = gdept_1d(:) 
    179       z_bnds(2:jpk,1) = gdept_1d(1:jpkm1) 
    180       z_bnds(1    ,1) = gdept_1d(1) - e3w_1d(1) 
     207 
     208#if ! defined key_xios2 
     209      z_bnds(:    ,2)  = gdept_1d(:) 
     210      z_bnds(2:jpk,1)  = gdept_1d(1:jpkm1) 
     211      z_bnds(1    ,1)  = gdept_1d(1) - e3w_1d(1) 
     212#else 
     213      z_bnds(2,:    )  = gdept_1d(:) 
     214      z_bnds(1,2:jpk)  = gdept_1d(1:jpkm1) 
     215      z_bnds(1,1    )  = gdept_1d(1) - e3w_1d(1) 
     216#endif 
    181217      CALL iom_set_axis_attr( "depthw", bounds=z_bnds ) 
     218 
    182219 
    183220# if defined key_floats 
    … …  
    11581195      LOGICAL,  DIMENSION(:,:) , OPTIONAL, INTENT(in) ::   mask 
    11591196 
     1197#if ! defined key_xios2 
    11601198      IF ( xios_is_valid_domain     (cdid) ) THEN 
    11611199         CALL xios_set_domain_attr     ( cdid, ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, jbegin=jbegin, ni=ni, nj=nj,   & 
    … …  
    11641202            &    lonvalue=lonvalue, latvalue=latvalue, mask=mask, nvertex=nvertex, bounds_lon=bounds_lon,                  & 
    11651203            &    bounds_lat=bounds_lat, area=area ) 
    1166       ENDIF 
    1167  
     1204     ENDIF 
    11681205      IF ( xios_is_valid_domaingroup(cdid) ) THEN 
    11691206         CALL xios_set_domaingroup_attr( cdid, ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, jbegin=jbegin, ni=ni, nj=nj,   & 
    … …  
    11731210            &    bounds_lat=bounds_lat, area=area ) 
    11741211      ENDIF 
     1212 
     1213#else 
     1214      IF ( xios_is_valid_domain     (cdid) ) THEN 
     1215         CALL xios_set_domain_attr     ( cdid, ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, jbegin=jbegin, ni=ni, nj=nj,   & 
     1216            &    data_dim=data_dim, data_ibegin=data_ibegin, data_ni=data_ni, data_jbegin=data_jbegin, data_nj=data_nj ,   & 
     1217            &    lonvalue_1D=lonvalue, latvalue_1D=latvalue, mask_2D=mask, nvertex=nvertex, bounds_lon_1D=bounds_lon,      & 
     1218            &    bounds_lat_1D=bounds_lat, area=area, type='curvilinear') 
     1219     ENDIF 
     1220      IF ( xios_is_valid_domaingroup(cdid) ) THEN 
     1221         CALL xios_set_domaingroup_attr( cdid, ni_glo=ni_glo, nj_glo=nj_glo, ibegin=ibegin, jbegin=jbegin, ni=ni, nj=nj,   & 
     1222            &    data_dim=data_dim, data_ibegin=data_ibegin, data_ni=data_ni, data_jbegin=data_jbegin, data_nj=data_nj ,   & 
     1223            &    lonvalue_1D=lonvalue, latvalue_1D=latvalue, mask_2D=mask, nvertex=nvertex, bounds_lon_1D=bounds_lon,      & 
     1224            &    bounds_lat_1D=bounds_lat, area=area, type='curvilinear' ) 
     1225      ENDIF 
     1226#endif 
    11751227      CALL xios_solve_inheritance() 
    11761228 
    11771229   END SUBROUTINE iom_set_domain_attr 
     1230 
     1231#if defined key_xios2 
     1232  SUBROUTINE iom_set_zoom_domain_attr( cdid, ibegin, jbegin, ni, nj) 
     1233     CHARACTER(LEN=*)                   , INTENT(in) ::   cdid 
     1234     INTEGER                  , OPTIONAL, INTENT(in) ::   ibegin, jbegin, ni, nj 
     1235 
     1236     IF ( xios_is_valid_domain     (cdid) ) THEN 
     1237         CALL xios_set_zoom_domain_attr     ( cdid, ibegin=ibegin, jbegin=jbegin, ni=ni,    & 
     1238           &   nj=nj) 
     1239    ENDIF 
     1240  END SUBROUTINE iom_set_zoom_domain_attr 
     1241#endif 
    11781242 
    11791243 
    … …  
    11831247      REAL(wp), DIMENSION(:,:), OPTIONAL, INTENT(in) ::   bounds 
    11841248      IF ( PRESENT(paxis) ) THEN 
     1249#if ! defined key_xios2 
    11851250         IF ( xios_is_valid_axis     (cdid) )   CALL xios_set_axis_attr     ( cdid, size=SIZE(paxis), value=paxis ) 
    11861251         IF ( xios_is_valid_axisgroup(cdid) )   CALL xios_set_axisgroup_attr( cdid, size=SIZE(paxis), value=paxis ) 
     1252#else 
     1253         IF ( xios_is_valid_axis     (cdid) )   CALL xios_set_axis_attr     ( cdid, n_glo=SIZE(paxis), value=paxis ) 
     1254         IF ( xios_is_valid_axisgroup(cdid) )   CALL xios_set_axisgroup_attr( cdid, n_glo=SIZE(paxis), value=paxis ) 
     1255#endif 
    11871256      ENDIF 
    11881257      IF ( xios_is_valid_axis     (cdid) )   CALL xios_set_axis_attr     ( cdid, bounds=bounds ) 
    … …  
    11911260   END SUBROUTINE iom_set_axis_attr 
    11921261 
    1193  
    11941262   SUBROUTINE iom_set_field_attr( cdid, freq_op, freq_offset ) 
    11951263      CHARACTER(LEN=*)          , INTENT(in) ::   cdid 
    1196       CHARACTER(LEN=*),OPTIONAL , INTENT(in) ::   freq_op 
    1197       CHARACTER(LEN=*),OPTIONAL , INTENT(in) ::   freq_offset 
    1198       IF ( xios_is_valid_field     (cdid) )   CALL xios_set_field_attr     ( cdid, freq_op=freq_op, freq_offset=freq_offset ) 
    1199       IF ( xios_is_valid_fieldgroup(cdid) )   CALL xios_set_fieldgroup_attr( cdid, freq_op=freq_op, freq_offset=freq_offset ) 
     1264#if ! defined key_xios2 
     1265      CHARACTER(LEN=*)   ,OPTIONAL , INTENT(in) ::   freq_op 
     1266      CHARACTER(LEN=*)   ,OPTIONAL , INTENT(in) ::   freq_offset 
     1267#else 
     1268      TYPE(xios_duration),OPTIONAL , INTENT(in) ::   freq_op 
     1269      TYPE(xios_duration),OPTIONAL , INTENT(in) ::   freq_offset 
     1270#endif 
     1271      IF ( xios_is_valid_field     (cdid) )   CALL xios_set_field_attr       & 
     1272    &     ( cdid, freq_op=freq_op, freq_offset=freq_offset ) 
     1273      IF ( xios_is_valid_fieldgroup(cdid) )   CALL xios_set_fieldgroup_attr  & 
     1274    &                    ( cdid, freq_op=freq_op, freq_offset=freq_offset ) 
    12001275      CALL xios_solve_inheritance() 
    12011276   END SUBROUTINE iom_set_field_attr 
    1202  
    12031277 
    12041278   SUBROUTINE iom_set_file_attr( cdid, name, name_suffix ) 
    … …  
    12131287   SUBROUTINE iom_get_file_attr( cdid, name, name_suffix, output_freq ) 
    12141288      CHARACTER(LEN=*)          , INTENT(in ) ::   cdid 
    1215       CHARACTER(LEN=*),OPTIONAL , INTENT(out) ::   name, name_suffix, output_freq 
     1289      CHARACTER(LEN=*),OPTIONAL , INTENT(out) ::   name, name_suffix 
     1290#if ! defined key_xios2 
     1291      CHARACTER(LEN=*),OPTIONAL , INTENT(out) ::    output_freq 
     1292#else 
     1293      TYPE(xios_duration)   ,OPTIONAL , INTENT(out) :: output_freq 
     1294#endif   
    12161295      LOGICAL                                 ::   llexist1,llexist2,llexist3 
    12171296      !--------------------------------------------------------------------- 
    12181297      IF( PRESENT( name        ) )   name = ''          ! default values 
    12191298      IF( PRESENT( name_suffix ) )   name_suffix = '' 
     1299#if ! defined key_xios2 
    12201300      IF( PRESENT( output_freq ) )   output_freq = '' 
     1301#else 
     1302      IF( PRESENT( output_freq ) )   output_freq = xios_duration(0,0,0,0,0,0) 
     1303#endif 
    12211304      IF ( xios_is_valid_file     (cdid) ) THEN 
    12221305         CALL xios_solve_inheritance() 
    … …  
    12391322      CHARACTER(LEN=*)                   , INTENT(in) ::   cdid 
    12401323      LOGICAL, DIMENSION(:,:,:), OPTIONAL, INTENT(in) ::   mask 
     1324#if ! defined key_xios2 
    12411325      IF ( xios_is_valid_grid     (cdid) )   CALL xios_set_grid_attr     ( cdid, mask=mask ) 
    12421326      IF ( xios_is_valid_gridgroup(cdid) )   CALL xios_set_gridgroup_attr( cdid, mask=mask ) 
     1327#else 
     1328      IF ( xios_is_valid_grid     (cdid) )   CALL xios_set_grid_attr     ( cdid, mask3=mask ) 
     1329      IF ( xios_is_valid_gridgroup(cdid) )   CALL xios_set_gridgroup_attr( cdid, mask3=mask ) 
     1330#endif 
    12431331      CALL xios_solve_inheritance() 
    12441332   END SUBROUTINE iom_set_grid_attr 
    … …  
    12821370      ni=nlei-nldi+1 ; nj=nlej-nldj+1 
    12831371 
    1284       CALL iom_set_domain_attr("grid_"//cdgrd, ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-1, jbegin=njmpp+nldj-1, ni=ni, nj=nj) 
     1372#if ! defined key_xios2 
     1373     CALL iom_set_domain_attr("grid_"//cdgrd, ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-1, jbegin=njmpp+nldj-1, ni=ni, nj=nj) 
     1374#else 
     1375     CALL iom_set_domain_attr("grid_"//cdgrd, ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-2, jbegin=njmpp+nldj-2, ni=ni, nj=nj) 
     1376#endif      
    12851377      CALL iom_set_domain_attr("grid_"//cdgrd, data_dim=2, data_ibegin = 1-nldi, data_ni = jpi, data_jbegin = 1-nldj, data_nj = jpj) 
    12861378      CALL iom_set_domain_attr("grid_"//cdgrd, lonvalue = RESHAPE(plon(nldi:nlei, nldj:nlej),(/ ni*nj /)),   & 
    … …  
    14301522      ALLOCATE( zlon(ni*nj) )       ;       zlon(:) = 0. 
    14311523 
     1524      CALL dom_ngb( 180., 90., ix, iy, 'T' ) !  i-line that passes near the North Pole : Reference latitude (used in plots) 
     1525#if ! defined key_xios2 
    14321526      CALL iom_set_domain_attr("gznl", ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-1, jbegin=njmpp+nldj-1, ni=ni, nj=nj) 
    14331527      CALL iom_set_domain_attr("gznl", data_dim=2, data_ibegin = 1-nldi, data_ni = jpi, data_jbegin = 1-nldj, data_nj = jpj) 
    … …  
    14351529         &                             latvalue = RESHAPE(plat(nldi:nlei, nldj:nlej),(/ ni*nj /)))   
    14361530      ! 
    1437       CALL dom_ngb( 180., 90., ix, iy, 'T' ) !  i-line that passes near the North Pole : Reference latitude (used in plots) 
    14381531      CALL iom_set_domain_attr ('ptr', zoom_ibegin=ix, zoom_nj=jpjglo) 
     1532#else 
     1533! Pas teste : attention aux indices ! 
     1534      CALL iom_set_domain_attr("ptr", ni_glo=jpiglo, nj_glo=jpjglo, ibegin=nimpp+nldi-2, jbegin=njmpp+nldj-2, ni=ni, nj=nj) 
     1535      CALL iom_set_domain_attr("ptr", data_dim=2, data_ibegin = 1-nldi, data_ni = jpi, data_jbegin = 1-nldj, data_nj = jpj) 
     1536      CALL iom_set_domain_attr("ptr", lonvalue = zlon,   & 
     1537         &                             latvalue = RESHAPE(plat(nldi:nlei, nldj:nlej),(/ ni*nj /)))   
     1538       CALL iom_set_zoom_domain_attr ('ptr', ibegin=ix, nj=jpjglo) 
     1539#endif 
     1540 
    14391541      CALL iom_update_file_name('ptr') 
    14401542      ! 
    … …  
    14551557      zz=REAL(narea,wp) 
    14561558      CALL iom_set_domain_attr('scalarpoint', lonvalue=zz, latvalue=zz) 
    1457  
     1559       
    14581560   END SUBROUTINE set_scalar 
    14591561 
    … …  
    14791581      REAL(wp)        ,DIMENSION( 3) ::   zlonpira                 ! longitudes of pirata moorings 
    14801582      REAL(wp)        ,DIMENSION( 9) ::   zlatpira                 ! latitudes  of pirata moorings 
     1583#if  defined key_xios2 
     1584      TYPE(xios_duration)            ::   f_op, f_of 
     1585#endif 
     1586  
    14811587      !!---------------------------------------------------------------------- 
    14821588      !  
    14831589      ! frequency of the call of iom_put (attribut: freq_op) 
    1484       WRITE(cl1,'(i1)')        1   ;   CALL iom_set_field_attr('field_definition', freq_op = cl1//'ts', freq_offset='0ts') 
    1485       WRITE(cl1,'(i1)')  nn_fsbc   ;   CALL iom_set_field_attr('SBC'             , freq_op = cl1//'ts', freq_offset='0ts') 
    1486       WRITE(cl1,'(i1)')  nn_fsbc   ;   CALL iom_set_field_attr('SBC_scalar'      , freq_op = cl1//'ts', freq_offset='0ts') 
    1487       WRITE(cl1,'(i1)') nn_dttrc   ;   CALL iom_set_field_attr('ptrc_T'          , freq_op = cl1//'ts', freq_offset='0ts') 
    1488       WRITE(cl1,'(i1)') nn_dttrc   ;   CALL iom_set_field_attr('diad_T'          , freq_op = cl1//'ts', freq_offset='0ts') 
     1590#if ! defined key_xios2 
     1591      WRITE(cl1,'(i1)')        1   ;   CALL iom_set_field_attr('field_definition', freq_op=cl1//'ts', freq_offset='0ts') 
     1592      WRITE(cl1,'(i1)')  nn_fsbc   ;   CALL iom_set_field_attr('SBC'             , freq_op=cl1//'ts', freq_offset='0ts') 
     1593      WRITE(cl1,'(i1)')  nn_fsbc   ;   CALL iom_set_field_attr('SBC_scalar'      , freq_op=cl1//'ts', freq_offset='0ts') 
     1594      WRITE(cl1,'(i1)') nn_dttrc   ;   CALL iom_set_field_attr('ptrc_T'          , freq_op=cl1//'ts', freq_offset='0ts') 
     1595      WRITE(cl1,'(i1)') nn_dttrc   ;   CALL iom_set_field_attr('diad_T'          , freq_op=cl1//'ts', freq_offset='0ts') 
     1596#else 
     1597      f_op%timestep = 1        ;  f_of%timestep = 0  ; CALL iom_set_field_attr('field_definition', freq_op=f_op, freq_offset=f_of) 
     1598      f_op%timestep = nn_fsbc  ;  f_of%timestep = 0  ; CALL iom_set_field_attr('SBC'             , freq_op=f_op, freq_offset=f_of) 
     1599      f_op%timestep = nn_fsbc  ;  f_of%timestep = 0  ; CALL iom_set_field_attr('SBC_scalar'      , freq_op=f_op, freq_offset=f_of) 
     1600      f_op%timestep = nn_dttrc ;  f_of%timestep = 0  ; CALL iom_set_field_attr('ptrc_T'          , freq_op=f_op, freq_offset=f_of) 
     1601      f_op%timestep = nn_dttrc ;  f_of%timestep = 0  ; CALL iom_set_field_attr('diad_T'          , freq_op=f_op, freq_offset=f_of) 
     1602#endif 
    14891603        
    14901604      ! output file names (attribut: name) 
    … …  
    15081622         ! Equatorial section (attributs: jbegin, ni, name_suffix) 
    15091623         CALL dom_ngb( 0., 0., ix, iy, cl1 ) 
     1624#if ! defined key_xios2 
    15101625         CALL iom_set_domain_attr ('Eq'//cl1, zoom_jbegin=iy, zoom_ni=jpiglo) 
     1626#else 
     1627         CALL iom_set_zoom_domain_attr ('Eq'//cl1, jbegin=iy-1, ni=jpiglo) 
     1628#endif 
    15111629         CALL iom_get_file_attr   ('Eq'//cl1, name_suffix = clsuff             ) 
    15121630         CALL iom_set_file_attr   ('Eq'//cl1, name_suffix = TRIM(clsuff)//'_Eq') 
    … …  
    15881706               ENDIF 
    15891707               clname = TRIM(ADJUSTL(clat))//TRIM(ADJUSTL(clon)) 
     1708#if ! defined key_xios2 
    15901709               CALL iom_set_domain_attr (TRIM(clname)//cl1, zoom_ibegin= ix, zoom_jbegin= iy) 
     1710#else 
     1711               CALL iom_set_zoom_domain_attr  (TRIM(clname)//cl1, ibegin= ix-1, jbegin= iy-1) 
     1712#endif 
    15911713               CALL iom_get_file_attr   (TRIM(clname)//cl1, name_suffix = clsuff                         ) 
    15921714               CALL iom_set_file_attr   (TRIM(clname)//cl1, name_suffix = TRIM(clsuff)//'_'//TRIM(clname)) 
    … …  
    16171739      REAL(wp)           ::   zsec 
    16181740      LOGICAL            ::   llexist 
    1619       !!---------------------------------------------------------------------- 
     1741#if  defined key_xios2 
     1742      TYPE(xios_duration)   ::   output_freq  
     1743#endif       
     1744      !!---------------------------------------------------------------------- 
     1745 
    16201746 
    16211747      DO jn = 1,2 
    1622  
     1748#if ! defined key_xios2 
    16231749         IF( jn == 1 )   CALL iom_get_file_attr( cdid, name        = clname, output_freq = clfreq ) 
     1750#else 
     1751         output_freq = xios_duration(0,0,0,0,0,0) 
     1752         IF( jn == 1 )   CALL iom_get_file_attr( cdid, name        = clname, output_freq = output_freq ) 
     1753#endif 
    16241754         IF( jn == 2 )   CALL iom_get_file_attr( cdid, name_suffix = clname ) 
    16251755 
    … …  
    16321762            END DO 
    16331763 
     1764#if ! defined key_xios2 
    16341765            idx = INDEX(clname,'@freq@') + INDEX(clname,'@FREQ@') 
    16351766            DO WHILE ( idx /= 0 )  
    … …  
    16441775               idx = INDEX(clname,'@freq@') + INDEX(clname,'@FREQ@') 
    16451776            END DO 
    1646  
     1777#else 
     1778            idx = INDEX(clname,'@freq@') + INDEX(clname,'@FREQ@') 
     1779            DO WHILE ( idx /= 0 )  
     1780              IF ( output_freq%hour /= 0 ) THEN 
     1781                  WRITE(clfreq,'(I19,A1)')INT(output_freq%hour),'h'  
     1782                  itrlen = LEN_TRIM(ADJUSTL(clfreq)) 
     1783              ELSE IF ( output_freq%day /= 0 ) THEN 
     1784                  WRITE(clfreq,'(I19,A1)')INT(output_freq%day),'d'  
     1785                  itrlen = LEN_TRIM(ADJUSTL(clfreq)) 
     1786              ELSE IF ( output_freq%month /= 0 ) THEN    
     1787                  WRITE(clfreq,'(I19,A1)')INT(output_freq%month),'m'  
     1788                  itrlen = LEN_TRIM(ADJUSTL(clfreq)) 
     1789              ELSE IF ( output_freq%year /= 0 ) THEN    
     1790                  WRITE(clfreq,'(I19,A1)')INT(output_freq%year),'y'  
     1791                  itrlen = LEN_TRIM(ADJUSTL(clfreq)) 
     1792              ELSE 
     1793                  CALL ctl_stop('error in the name of file id '//TRIM(cdid),   & 
     1794                     & ' attribute output_freq is undefined -> cannot replace @freq@ in '//TRIM(clname) ) 
     1795              ENDIF 
     1796              clname = clname(1:idx-1)//TRIM(ADJUSTL(clfreq))//clname(idx+6:LEN_TRIM(clname)) 
     1797              idx = INDEX(clname,'@freq@') + INDEX(clname,'@FREQ@') 
     1798            END DO 
     1799#endif 
    16471800            idx = INDEX(clname,'@startdate@') + INDEX(clname,'@STARTDATE@') 
    16481801            DO WHILE ( idx /= 0 )  
    … …  
    16731826            END DO 
    16741827 
     1828            IF( jn == 1 .AND. TRIM(Agrif_CFixed()) /= '0' )   clname = TRIM(Agrif_CFixed())//"_"//TRIM(clname) 
    16751829            IF( jn == 1 )   CALL iom_set_file_attr( cdid, name        = clname ) 
    16761830            IF( jn == 2 )   CALL iom_set_file_attr( cdid, name_suffix = clname ) 
    … …  
    17201874      ENDIF 
    17211875       
     1876!$AGRIF_DO_NOT_TREAT       
     1877! Should be fixed in the conv 
    17221878      IF( llfull ) THEN  
    17231879         clfmt = TRIM(clfmt)//",'_',i2.2,':',i2.2,':',i2.2" 
    … …  
    17301886         WRITE(iom_sdate, '('//TRIM(clfmt)//')') iyear, imonth, iday                          ! date of the end of run 
    17311887      ENDIF 
     1888!$AGRIF_END_DO_NOT_TREAT       
    17321889 
    17331890   END FUNCTION iom_sdate 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/LBC/lib_mpp.F90

    r7730 r7731  
    298298      ENDIF 
    299299 
     300#if defined key_agrif 
     301      IF (Agrif_Root()) THEN 
     302         CALL Agrif_MPI_Init(mpi_comm_opa) 
     303      ELSE 
     304         CALL Agrif_MPI_set_grid_comm(mpi_comm_opa) 
     305      ENDIF 
     306#endif 
     307 
    300308      CALL mpi_comm_rank( mpi_comm_opa, mpprank, ierr ) 
    301309      CALL mpi_comm_size( mpi_comm_opa, mppsize, ierr ) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90

    r7730 r7731  
    3232   PUBLIC   fld_map    ! routine called by tides_init 
    3333   PUBLIC   fld_read, fld_fill   ! called by sbc... modules 
     34   PUBLIC   fld_clopn 
    3435 
    3536   TYPE, PUBLIC ::   FLD_N      !: Namelist field informations 
    … …  
    815816         imonth = kmonth 
    816817         iday = kday 
     818         IF ( sdjf%cltype(1:4) == 'week' ) THEN             ! find the day of the beginning of the week 
     819            isec_week = ksec_week( sdjf%cltype(6:8) )- (86400 * 8 )   
     820            llprevmth  = isec_week > nsec_month             ! longer time since beginning of the week than the month 
     821            llprevyr   = llprevmth .AND. nmonth == 1 
     822            iyear  = nyear  - COUNT((/llprevyr /)) 
     823            imonth = nmonth - COUNT((/llprevmth/)) + 12 * COUNT((/llprevyr /)) 
     824            iday   = nday   + nmonth_len(nmonth-1) * COUNT((/llprevmth/)) - isec_week / NINT(rday) 
     825         ENDIF 
    817826      ELSE                                                  ! use current day values 
    818827         IF ( sdjf%cltype(1:4) == 'week' ) THEN             ! find the day of the beginning of the week 
    … …  
    12811290      CHARACTER(LEN=*)          , INTENT(in   ) ::   lsmfile ! land sea mask file name 
    12821291      !!  
    1283       REAL(wp),DIMENSION(:,:,:),ALLOCATABLE     ::   ztmp_fly_dta,zfieldo                  ! temporary array of values on input grid 
     1292      REAL(wp),DIMENSION(:,:,:),ALLOCATABLE     ::   ztmp_fly_dta                          ! temporary array of values on input grid 
    12841293      INTEGER, DIMENSION(3)                     ::   rec1,recn                             ! temporary arrays for start and length 
    12851294      INTEGER, DIMENSION(3)                     ::   rec1_lsm,recn_lsm                     ! temporary arrays for start and length in case of seaoverland 
    … …  
    13471356 
    13481357 
    1349          itmpi=SIZE(ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:),1) 
    1350          itmpj=SIZE(ztmp_fly_dta(jpi1_lsm:jpi2_lsm,jpj1_lsm:jpj2_lsm,:),2) 
     1358         itmpi=jpi2_lsm-jpi1_lsm+1 
     1359         itmpj=jpj2_lsm-jpj1_lsm+1 
    13511360         itmpz=kk 
    13521361         ALLOCATE(ztmp_fly_dta(itmpi,itmpj,itmpz)) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r7730 r7731  
    403403         CALL iom_put( "qsr_oce" ,   qsr  )                 ! output downward solar heat over the ocean 
    404404         CALL iom_put( "qt_oce"  ,   qns+qsr )              ! output total downward heat over the ocean 
     405         tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac   ! output total precipitation [kg/m2/s] 
     406         sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac   ! output solid precipitation [kg/m2/s] 
     407         CALL iom_put( 'snowpre', sprecip * 86400. )        ! Snow 
     408         CALL iom_put( 'precip' , tprecip * 86400. )        ! Total precipitation 
    405409      ENDIF 
    406410      ! 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r7730 r7731  
    10291029         ssu_m(:,:) = frcv(jpr_ocx1)%z3(:,:,1) 
    10301030         ub (:,:,1) = ssu_m(:,:)                             ! will be used in sbcice_lim in the call of lim_sbc_tau 
     1031         un (:,:,1) = ssu_m(:,:)                             ! will be used in sbc_cpl_snd if atmosphere coupling 
    10311032         CALL iom_put( 'ssu_m', ssu_m ) 
    10321033      ENDIF 
    … …  
    10341035         ssv_m(:,:) = frcv(jpr_ocy1)%z3(:,:,1) 
    10351036         vb (:,:,1) = ssv_m(:,:)                             ! will be used in sbcice_lim in the call of lim_sbc_tau 
     1037         vn (:,:,1) = ssv_m(:,:)                             ! will be used in sbc_cpl_snd if atmosphere coupling 
    10361038         CALL iom_put( 'ssv_m', ssv_m ) 
    10371039      ENDIF 
    … …  
    17431745                     ztmp3(:,:,1) = SUM( tn_ice * a_i, dim=3 ) / SUM( a_i, dim=3 ) 
    17441746                  ELSEWHERE 
    1745                      ztmp3(:,:,1) = rt0 ! TODO: Is freezing point a good default? (Maybe SST is better?) 
     1747                     ztmp3(:,:,1) = rt0 
    17461748                  END WHERE 
    17471749               CASE default   ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%clcat' ) 
    … …  
    17741776      !                                                      ! ------------------------- ! 
    17751777      IF( ssnd(jps_albice)%laction ) THEN                         ! ice  
    1776          SELECT CASE( sn_snd_alb%cldes ) 
    1777          CASE( 'ice'          )   ; ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) 
    1778          CASE( 'weighted ice' )   ; ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl) 
    1779          CASE default             ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_alb%cldes' ) 
     1778          SELECT CASE( sn_snd_alb%cldes ) 
     1779          CASE( 'ice' ) 
     1780             SELECT CASE( sn_snd_alb%clcat ) 
     1781             CASE( 'yes' )    
     1782                ztmp3(:,:,1:jpl) = alb_ice(:,:,1:jpl) 
     1783             CASE( 'no' ) 
     1784                WHERE( SUM( a_i, dim=3 ) /= 0. ) 
     1785                   ztmp1(:,:) = SUM( alb_ice (:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 ) / SUM( a_i(:,:,1:jpl), dim=3 ) 
     1786                ELSEWHERE 
     1787                   ztmp1(:,:) = albedo_oce_mix(:,:) 
     1788                END WHERE 
     1789             CASE default   ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_alb%clcat' ) 
     1790             END SELECT 
     1791          CASE( 'weighted ice' )   ; 
     1792             SELECT CASE( sn_snd_alb%clcat ) 
     1793             CASE( 'yes' )    
     1794                ztmp3(:,:,1:jpl) =  alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl) 
     1795             CASE( 'no' ) 
     1796                WHERE( fr_i (:,:) > 0. ) 
     1797                   ztmp1(:,:) = SUM (  alb_ice(:,:,1:jpl) * a_i(:,:,1:jpl), dim=3 ) 
     1798                ELSEWHERE 
     1799                   ztmp1(:,:) = 0. 
     1800                END WHERE 
     1801             CASE default   ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_ice%clcat' ) 
     1802             END SELECT 
     1803          CASE default      ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_alb%cldes' ) 
    17801804         END SELECT 
    1781          CALL cpl_snd( jps_albice, isec, ztmp3, info ) 
    1782       ENDIF 
     1805 
     1806         SELECT CASE( sn_snd_alb%clcat ) 
     1807            CASE( 'yes' )    
     1808               CALL cpl_snd( jps_albice, isec, ztmp3, info )      !-> MV this has never been checked in coupled mode 
     1809            CASE( 'no'  )    
     1810               CALL cpl_snd( jps_albice, isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )  
     1811         END SELECT 
     1812      ENDIF 
     1813 
    17831814      IF( ssnd(jps_albmix)%laction ) THEN                         ! mixed ice-ocean 
    17841815         ztmp1(:,:) = albedo_oce_mix(:,:) * zfr_l(:,:) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90

    r7730 r7731  
    108108         ! 
    109109         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN 
    110             z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) - snwice_fmass(:,:) ) ) / area   ! sum over the global domain 
     110            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) - snwice_fmass(:,:) ) ) / area   ! sum over the global domain 
    111111            zcoef = z_fwf * rcp 
    112112            emp(:,:) = emp(:,:) - z_fwf              * tmask(:,:,1) 
    … …  
    162162            zsurf_pos = glob_sum( e1e2t(:,:)*ztmsk_pos(:,:) ) 
    163163            !                                                  ! fwf global mean (excluding ocean to ice/snow exchanges)  
    164             z_fwf     = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) - snwice_fmass(:,:) ) ) / area 
     164            z_fwf     = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) - snwice_fmass(:,:) ) ) / area 
    165165            !             
    166166            IF( z_fwf < 0._wp ) THEN         ! spread out over >0 erp area to increase evaporation 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90

    r7730 r7731  
    5353   REAL(wp)   , PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:)     ::  risfLeff               !:effective length (Leff) BG03 nn_isf==2 
    5454   REAL(wp)   , PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:)     ::  ttbl, stbl, utbl, vtbl !:top boundary layer variable at T point 
    55 #if defined key_agrif 
    56    ! AGRIF can not handle these arrays as integers. The reason is a mystery but problems avoided by declaring them as reals 
    57    REAL(wp),    PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:)     ::  misfkt, misfkb         !:Level of ice shelf base 
    58                                                                                           !: (first wet level and last level include in the tbl) 
    59 #else 
    6055   INTEGER,    PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:)     ::  misfkt, misfkb         !:Level of ice shelf base 
    61 #endif 
    6256 
    6357 
    … …  
    9286    REAL(wp)                     ::   rmin 
    9387    REAL(wp)                     ::   zhk 
    94     CHARACTER(len=256)           ::   cfisf, cvarzisf, cvarhisf   ! name for isf file 
     88    REAL(wp)                     ::   zt_frz, zpress 
     89    CHARACTER(len=256)           ::   cfisf , cvarzisf, cvarhisf   ! name for isf file 
    9590    CHARACTER(LEN=256)           :: cnameis                     ! name of iceshelf file 
    9691    CHARACTER (LEN=32)           :: cvarLeff                    ! variable name for efficient Length scale 
    … …  
    176171              DO jj = 1, jpj 
    177172                  jk = 2 
    178                   DO WHILE ( jk .LE. mbkt(ji,jj) .AND. fsdepw(ji,jj,jk) < rzisf_tbl(ji,jj) ) ;  jk = jk + 1 ;  END DO 
     173                  DO WHILE ( jk .LE. mbkt(ji,jj) .AND. gdepw_0(ji,jj,jk) < rzisf_tbl(ji,jj) ) ;  jk = jk + 1 ;  END DO 
    179174                  misfkt(ji,jj) = jk-1 
    180175               END DO 
    … …  
    194189         END IF 
    195190          
     191         ! save initial top boundary layer thickness          
    196192         rhisf_tbl_0(:,:) = rhisf_tbl(:,:) 
     193 
     194      END IF 
     195 
     196      !                                            ! ---------------------------------------- ! 
     197      IF( kt /= nit000 ) THEN                      !          Swap of forcing fields          ! 
     198         !                                         ! ---------------------------------------- ! 
     199         fwfisf_b  (:,:  ) = fwfisf  (:,:  )               ! Swap the ocean forcing fields except at nit000 
     200         risf_tsc_b(:,:,:) = risf_tsc(:,:,:)               ! where before fields are set at the end of the routine 
     201         ! 
     202      ENDIF 
     203 
     204      IF( MOD( kt-1, nn_fsbc) == 0 ) THEN 
    197205 
    198206         ! compute bottom level of isf tbl and thickness of tbl below the ice shelf 
    … …  
    205213 
    206214               ! determine the deepest level influenced by the boundary layer 
    207                ! test on tmask useless ????? 
    208215               DO jk = ikt, mbkt(ji,jj) 
    209216                  IF ( (SUM(fse3t_n(ji,jj,ikt:jk-1)) .LT. rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) ) ikb = jk 
    … …  
    217224            END DO 
    218225         END DO 
    219           
    220       END IF 
    221  
    222       !                                            ! ---------------------------------------- ! 
    223       IF( kt /= nit000 ) THEN                      !          Swap of forcing fields          ! 
    224          !                                         ! ---------------------------------------- ! 
    225          fwfisf_b  (:,:  ) = fwfisf  (:,:  )               ! Swap the ocean forcing fields except at nit000 
    226          risf_tsc_b(:,:,:) = risf_tsc(:,:,:)               ! where before fields are set at the end of the routine 
    227          ! 
    228       ENDIF 
    229  
    230       IF( MOD( kt-1, nn_fsbc) == 0 ) THEN 
    231  
    232226 
    233227         ! compute salf and heat flux 
    … …  
    270264         END IF 
    271265         ! compute tsc due to isf 
    272          ! WARNING water add at temp = 0C, correction term is added in trasbc, maybe better here but need a 3D variable). 
    273          risf_tsc(:,:,jp_tem) = qisf(:,:) * r1_rau0_rcp ! 
     266         ! WARNING water add at temp = 0C, correction term is added, maybe better here but need a 3D variable). 
     267!         zpress = grav*rau0*fsdept(ji,jj,jk)*1.e-04 
     268         zt_frz = -1.9 !eos_fzp( tsn(ji,jj,jk,jp_sal), zpress ) 
     269         risf_tsc(:,:,jp_tem) = qisf(:,:) * r1_rau0_rcp - rdivisf * fwfisf(:,:) * zt_frz * r1_rau0 ! 
    274270          
    275271         ! salt effect already take into account in vertical advection 
    276272         risf_tsc(:,:,jp_sal) = (1.0_wp-rdivisf) * fwfisf(:,:) * stbl(:,:) * r1_rau0 
    277            
     273 
     274         ! output 
     275         IF( iom_use('qisf'  ) )   CALL iom_put('qisf'  , qisf) 
     276         IF( iom_use('fwfisf') )   CALL iom_put('fwfisf', fwfisf * stbl(:,:) / soce ) 
     277 
     278         ! if apply only on the trend and not as a volume flux (rdivisf = 0), fwfisf have to be set to 0 now 
     279         fwfisf(:,:) = rdivisf * fwfisf(:,:)          
     280  
    278281         ! lbclnk 
    279282         CALL lbc_lnk(risf_tsc(:,:,jp_tem),'T',1.) 
    … …  
    295298         ENDIF 
    296299         !  
    297          ! output 
    298          CALL iom_put('qisf'  , qisf) 
    299          IF( iom_use('fwfisf') )   CALL iom_put('fwfisf', fwfisf * stbl(:,:) / soce ) 
    300300      END IF 
    301301   
    … …  
    472472 
    473473                     nit = nit + 1 
    474                      IF (nit .GE. 100) THEN 
    475                         !WRITE(numout,*) "sbcisf : too many iteration ... ", zhtflx, zhtflx_b,zgammat, rn_gammat0, rn_tfri2, nn_gammablk, ji,jj 
    476                         !WRITE(numout,*) "sbcisf : too many iteration ... ", (zhtflx - zhtflx_b)/zhtflx 
    477                         CALL ctl_stop( 'STOP', 'sbc_isf_hol99 : too many iteration ...' ) 
    478                      END IF 
     474                     IF (nit .GE. 100) CALL ctl_stop( 'STOP', 'sbc_isf_hol99 : too many iteration ...' ) 
     475 
    479476! save gammat and compute zhtflx_b 
    480477                     zgammat2d(ji,jj)=zgammat 
    … …  
    794791               ! test on tmask useless ????? 
    795792               DO jk = ikt, mbkt(ji,jj) 
    796 !                  IF ( (SUM(fse3t(ji,jj,ikt:jk-1)) .LT. rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) ) ikb = jk 
     793                  IF ( (SUM(fse3t(ji,jj,ikt:jk-1)) .LT. rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) ) ikb = jk 
    797794               END DO 
    798795               rhisf_tbl(ji,jj) = MIN(rhisf_tbl(ji,jj), SUM(fse3t(ji,jj,ikt:ikb)))  ! limit the tbl to water thickness. 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r7730 r7731  
    179179 
    180180      !                          ! Checks: 
    181       IF( nn_isf .EQ. 0 ) THEN                      ! no specific treatment in vicinity of ice shelf  
     181      IF( nn_isf .EQ. 0 ) THEN                      ! variable initialisation if no ice shelf  
    182182         IF( sbc_isf_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_init : unable to allocate sbc_isf arrays' ) 
    183          fwfisf  (:,:) = 0.0_wp 
    184          fwfisf_b(:,:) = 0.0_wp 
     183         fwfisf  (:,:)   = 0.0_wp ; fwfisf_b  (:,:)   = 0.0_wp 
     184         risf_tsc(:,:,:) = 0.0_wp ; risf_tsc_b(:,:,:) = 0.0_wp 
     185         rdivisf       = 0.0_wp 
    185186      END IF 
    186187      IF( nn_ice == 0 .AND. nn_components /= jp_iam_opa )   fr_i(:,:) = 0.e0 ! no ice in the domain, ice fraction is always zero 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90

    r7730 r7731  
    126126      IF(   ln_rnf_sal   )   CALL fld_read ( kt, nn_fsbc, sf_s_rnf )    ! idem for runoffs salinity    if required 
    127127      ! 
    128       ! Runoff reduction only associated to the ORCA2_LIM configuration 
    129       ! when reading the NetCDF file runoff_1m_nomask.nc 
    130       IF( cp_cfg == 'orca' .AND. jp_cfg == 2 .AND. .NOT. l_rnfcpl )   THEN 
    131          WHERE( 40._wp < gphit(:,:) .AND. gphit(:,:) < 65._wp ) 
    132             sf_rnf(1)%fnow(:,:,1) = 0.85 * sf_rnf(1)%fnow(:,:,1) 
    133          END WHERE 
    134       ENDIF 
    135       ! 
    136128      IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN 
    137129         ! 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/SBC/updtide.F90

    r7730 r7731  
    3131CONTAINS 
    3232 
    33    SUBROUTINE upd_tide( kt, kit, kbaro, koffset ) 
     33   SUBROUTINE upd_tide( kt, kit, time_offset ) 
    3434      !!---------------------------------------------------------------------- 
    3535      !!                 ***  ROUTINE upd_tide  *** 
    … …  
    4242      !!----------------------------------------------------------------------       
    4343      INTEGER, INTENT(in)           ::   kt      ! ocean time-step index 
    44       INTEGER, INTENT(in), OPTIONAL ::   kit     ! external mode sub-time-step index (lk_dynspg_ts=T only) 
    45       INTEGER, INTENT(in), OPTIONAL ::   kbaro   ! number of sub-time-step           (lk_dynspg_ts=T only) 
    46       INTEGER, INTENT(in), OPTIONAL ::   koffset ! time offset in number  
    47                                                  ! of sub-time-steps                 (lk_dynspg_ts=T only) 
     44      INTEGER, INTENT(in), OPTIONAL ::   kit     ! external mode sub-time-step index (lk_dynspg_ts=T) 
     45      INTEGER, INTENT(in), OPTIONAL ::   time_offset ! time offset in number  
     46                                                     ! of internal steps             (lk_dynspg_ts=F) 
     47                                                     ! of external steps             (lk_dynspg_ts=T) 
    4848      ! 
    4949      INTEGER  ::   joffset      ! local integer 
    … …  
    5757      ! 
    5858      joffset = 0 
    59       IF( PRESENT( koffset ) )   joffset = koffset 
     59      IF( PRESENT( time_offset ) )   joffset = time_offset 
    6060      ! 
    61       IF( PRESENT( kit ) .AND. PRESENT( kbaro ) )   THEN 
    62          zt = zt + ( kit + 0.5_wp * ( joffset - 1 ) ) * rdt / REAL( kbaro, wp ) 
     61      IF( PRESENT( kit ) )   THEN 
     62         zt = zt + ( kit +  joffset - 1 ) * rdt / REAL( nn_baro, wp ) 
    6363      ELSE 
    6464         zt = zt + joffset * rdt 
    … …  
    7474      IF( ln_tide_ramp ) THEN         ! linear increase if asked 
    7575         zt = ( kt - nit000 ) * rdt 
    76          IF( PRESENT( kit ) .AND. PRESENT( kbaro ) )   zt = zt + kit * rdt / REAL( kbaro, wp ) 
     76         IF( PRESENT( kit ) )   zt = zt + ( kit + joffset -1) * rdt / REAL( nn_baro, wp ) 
    7777         zramp = MIN(  MAX( zt / (rdttideramp*rday) , 0._wp ) , 1._wp  ) 
    7878         pot_astro(:,:) = zramp * pot_astro(:,:) 
    … …  
    8686  !!---------------------------------------------------------------------- 
    8787CONTAINS 
    88   SUBROUTINE upd_tide( kt, kit, kbaro, koffset )          ! Empty routine 
     88  SUBROUTINE upd_tide( kt, kit, time_offset )  ! Empty routine 
    8989    INTEGER, INTENT(in)           ::   kt      !  integer  arg, dummy routine 
    9090    INTEGER, INTENT(in), OPTIONAL ::   kit     !  optional arg, dummy routine 
    91     INTEGER, INTENT(in), OPTIONAL ::   kbaro   !  optional arg, dummy routine 
    92     INTEGER, INTENT(in), OPTIONAL ::   koffset !  optional arg, dummy routine 
     91    INTEGER, INTENT(in), OPTIONAL ::   time_offset !  optional arg, dummy routine 
    9392    WRITE(*,*) 'upd_tide: You should not have seen this print! error?', kt 
    9493  END SUBROUTINE upd_tide 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/STO/stopar.F90

    r7730 r7731  
    849849 
    850850 
    851    REAL(wp) FUNCTION sto_par_flt_fac( kpasses ) 
     851   FUNCTION sto_par_flt_fac( kpasses ) 
    852852      !!---------------------------------------------------------------------- 
    853853      !!                  ***  FUNCTION sto_par_flt_fac  *** 
    … …  
    858858      !!---------------------------------------------------------------------- 
    859859      INTEGER, INTENT(in) :: kpasses 
     860      REAL(wp) :: sto_par_flt_fac 
    860861      !! 
    861862      INTEGER :: jpasses, ji, jj, jflti, jfltj 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_eiv.F90

    r7730 r7731  
    212212      CHARACTER(len=3) ::   cdtype 
    213213      REAL, DIMENSION(:,:,:) ::   pun, pvn, pwn 
    214       WRITE(*,*) 'tra_adv_eiv: You should not have seen this print! error?', kt, cdtype, pun(1,1,1), pvn(1,1,1), pwn(1,1,1) 
     214      WRITE(*,*) 'tra_adv_eiv: You should not have seen this print! error?', & 
     215          &  kt, cdtype, pun(1,1,1), pvn(1,1,1), pwn(1,1,1) 
    215216   END SUBROUTINE tra_adv_eiv 
    216217#endif 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_tvd.F90

    r7730 r7731  
    326326      CALL wrk_alloc( jpi, jpj, zwx_sav, zwy_sav ) 
    327327      CALL wrk_alloc( jpi, jpj, jpk, zwi, zwz , zhdiv, zwz_sav, zwzts ) 
    328       CALL wrk_alloc( jpi, jpj, jpk, 3, ztrs ) 
     328      CALL wrk_alloc( jpi, jpj, jpk, kjpt+1, ztrs ) 
    329329      ! 
    330330      IF( kt == kit000 )  THEN 
    … …  
    564564      ! 
    565565                   CALL wrk_dealloc( jpi, jpj, jpk, zwi, zwz, zhdiv, zwz_sav, zwzts ) 
    566                    CALL wrk_dealloc( jpi, jpj, jpk, 3, ztrs ) 
     566                   CALL wrk_dealloc( jpi, jpj, jpk, kjpt+1, ztrs ) 
    567567                   CALL wrk_dealloc( jpi, jpj, zwx_sav, zwy_sav ) 
    568568      IF( l_trd )  CALL wrk_dealloc( jpi, jpj, jpk, ztrdx, ztrdy, ztrdz ) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90

    r7730 r7731  
    2828   USE sbc_oce         ! surface boundary condition: ocean 
    2929   USE sbcrnf          ! river runoffs 
     30   USE sbcisf          ! ice shelf melting/freezing 
    3031   USE zdf_oce         ! ocean vertical mixing 
    3132   USE domvvl          ! variable volume 
    … …  
    4647   USE timing          ! Timing 
    4748#if defined key_agrif 
    48    USE agrif_opa_update 
    4949   USE agrif_opa_interp 
    5050#endif 
    … …  
    110110      ! Update after tracer on domain lateral boundaries 
    111111      !  
     112#if defined key_agrif 
     113      CALL Agrif_tra                     ! AGRIF zoom boundaries 
     114#endif 
     115      ! 
    112116      CALL lbc_lnk( tsa(:,:,:,jp_tem), 'T', 1._wp )      ! local domain boundaries  (T-point, unchanged sign) 
    113117      CALL lbc_lnk( tsa(:,:,:,jp_sal), 'T', 1._wp ) 
    … …  
    115119#if defined key_bdy  
    116120      IF( lk_bdy )   CALL bdy_tra( kt )  ! BDY open boundaries 
    117 #endif 
    118 #if defined key_agrif 
    119       CALL Agrif_tra                     ! AGRIF zoom boundaries 
    120121#endif 
    121122  
    … …  
    148149         ELSE                 ;   CALL tra_nxt_fix( kt, nit000,         'TRA', tsb, tsn, tsa, jpts )  ! fixed    volume level  
    149150         ENDIF 
    150       ENDIF  
    151       ! 
    152 #if defined key_agrif 
    153       ! Update tracer at AGRIF zoom boundaries 
    154       IF( .NOT.Agrif_Root() )    CALL Agrif_Update_Tra( kt )      ! children only 
    155 #endif       
    156       ! 
    157       ! trends computation 
     151      ENDIF      
     152      ! 
     153     ! trends computation 
    158154      IF( l_trdtra ) THEN      ! trend of the Asselin filter (tb filtered - tb)/dt      
    159155         DO jk = 1, jpkm1 
    … …  
    279275 
    280276      !!      
    281       LOGICAL  ::   ll_tra_hpg, ll_traqsr, ll_rnf   ! local logical 
     277      LOGICAL  ::   ll_tra_hpg, ll_traqsr, ll_rnf, ll_isf   ! local logical 
    282278      INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices 
    283279      REAL(wp) ::   zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar 
    … …  
    295291         ll_traqsr  = ln_traqsr        ! active  tracers case  and  solar penetration 
    296292         ll_rnf     = ln_rnf           ! active  tracers case  and  river runoffs 
     293         IF (nn_isf .GE. 1) THEN  
     294            ll_isf = .TRUE.            ! active  tracers case  and  ice shelf melting/freezing 
     295         ELSE 
     296            ll_isf = .FALSE. 
     297         END IF 
    297298      ELSE                           
    298299         ll_tra_hpg = .FALSE.          ! passive tracers case or NO semi-implicit hpg 
    299300         ll_traqsr  = .FALSE.          ! active  tracers case and NO solar penetration 
    300301         ll_rnf     = .FALSE.          ! passive tracers or NO river runoffs 
     302         ll_isf     = .FALSE.          ! passive tracers or NO ice shelf melting/freezing 
    301303      ENDIF 
    302304      ! 
    … …  
    321323                  ztc_f  = ztc_n  + atfp * ztc_d 
    322324                  ! 
    323                   IF( jk == 1 ) THEN           ! first level  
    324                      ze3t_f = ze3t_f - zfact2 * ( emp_b(ji,jj) - emp(ji,jj) + rnf(ji,jj) - rnf_b(ji,jj) ) 
     325                  IF( jk == mikt(ji,jj) ) THEN           ! first level  
     326                     ze3t_f = ze3t_f - zfact2 * ( (emp_b(ji,jj)    - emp(ji,jj)   )  & 
     327                            &                   - (rnf_b(ji,jj)    - rnf(ji,jj)   )  & 
     328                            &                   + (fwfisf_b(ji,jj) - fwfisf(ji,jj))  ) 
    325329                     ztc_f  = ztc_f  - zfact1 * ( psbc_tc(ji,jj,jn) - psbc_tc_b(ji,jj,jn) ) 
    326330                  ENDIF 
    327331 
    328                   IF( ll_traqsr .AND. jn == jp_tem .AND. jk <= nksr )   &     ! solar penetration (temperature only) 
     332                  ! solar penetration (temperature only) 
     333                  IF( ll_traqsr .AND. jn == jp_tem .AND. jk <= nksr )                            &  
    329334                     &     ztc_f  = ztc_f  - zfact1 * ( qsr_hc(ji,jj,jk) - qsr_hc_b(ji,jj,jk) )  
    330335 
    331                   IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) )   &            ! river runoffs 
     336                  ! river runoff 
     337                  IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) )                                          & 
    332338                     &     ztc_f  = ztc_f  - zfact1 * ( rnf_tsc(ji,jj,jn) - rnf_tsc_b(ji,jj,jn) ) &  
    333339                     &                              * fse3t_n(ji,jj,jk) / h_rnf(ji,jj) 
     340 
     341                  ! ice shelf 
     342                  IF( ll_isf ) THEN 
     343                     ! level fully include in the Losch_2008 ice shelf boundary layer 
     344                     IF ( jk >= misfkt(ji,jj) .AND. jk < misfkb(ji,jj) )                          & 
     345                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) )  & 
     346                               &                 * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj) 
     347                     ! level partially include in Losch_2008 ice shelf boundary layer  
     348                     IF ( jk == misfkb(ji,jj) )                                                   & 
     349                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) )  & 
     350                               &                 * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj) * ralpha(ji,jj) 
     351                  END IF 
    334352 
    335353                  ze3t_f = 1.e0 / ze3t_f 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90

    r7730 r7731  
    3333   USE timing          ! Timing 
    3434   USE eosbn2 
     35#if defined key_asminc    
     36   USE asminc          ! Assimilation increment 
     37#endif 
    3538 
    3639   IMPLICIT NONE 
    … …  
    120123      REAL(wp) ::   zfact, z1_e3t, zdep 
    121124      REAL(wp) ::   zalpha, zhk 
    122       REAL(wp) ::  zt_frz, zpress 
    123125      REAL(wp), POINTER, DIMENSION(:,:,:) ::  ztrdt, ztrds 
    124126      !!---------------------------------------------------------------------- 
    … …  
    232234               DO jk = ikt, ikb - 1 
    233235               ! compute tfreez for the temperature correction (we add water at freezing temperature) 
    234 !                  zpress = grav*rau0*fsdept(ji,jj,jk)*1.e-04 
    235                   zt_frz = -1.9 !eos_fzp( tsn(ji,jj,jk,jp_sal), zpress ) 
    236236               ! compute trend 
    237237                  tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)                                          & 
    238                      &           + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)          & 
    239                      &               - rdivisf * (fwfisf(ji,jj) + fwfisf_b(ji,jj)) * zt_frz * r1_rau0) & 
    240                      &           * r1_hisf_tbl(ji,jj) 
     238                     &           + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)) * r1_hisf_tbl(ji,jj) 
    241239                  tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)                                          & 
    242240                     &           + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj) 
    … …  
    245243               ! level partially include in ice shelf boundary layer  
    246244               ! compute tfreez for the temperature correction (we add water at freezing temperature) 
    247 !               zpress = grav*rau0*fsdept(ji,jj,ikb)*1.e-04 
    248                zt_frz = -1.9 !eos_fzp( tsn(ji,jj,ikb,jp_sal), zpress ) 
    249245               ! compute trend 
    250246               tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem)                                           & 
    251                   &              + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)          & 
    252                   &                  - rdivisf * (fwfisf(ji,jj) + fwfisf_b(ji,jj)) * zt_frz * r1_rau0) &  
    253                   &              * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) 
     247                  &              + zfact * (risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) 
    254248               tsa(ji,jj,ikb,jp_sal) = tsa(ji,jj,ikb,jp_sal)                                           & 
    255249                  &              + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal)) * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)  
    … …  
    287281         END DO   
    288282      ENDIF 
     283 
     284#if defined key_asminc 
     285! WARNING: THIS MAY WELL NOT BE REQUIRED - WE DON'T WANT TO CHANGE T&S BUT THIS MAY COMPENSATE ANOTHER TERM... 
     286! Rate of change in e3t for each level is ssh_iau*e3t_0/ht_0 
     287! Contribution to tsa should be rate of change in level / per m of ocean? (hence the division by fse3t_n) 
     288      IF( ln_sshinc ) THEN         ! input of heat and salt due to assimilation 
     289         DO jj = 2, jpj  
     290            DO ji = fs_2, fs_jpim1 
     291               zdep = ssh_iau(ji,jj) / ( ht_0(ji,jj) + 1.0 - ssmask(ji, jj) ) 
     292               DO jk = 1, jpkm1 
     293                  tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)   & 
     294                                        &            + tsn(ji,jj,jk,jp_tem) * zdep * ( e3t_0(ji,jj,jk) / fse3t_n(ji,jj,jk) ) 
     295                  tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)   & 
     296                                        &            + tsn(ji,jj,jk,jp_sal) * zdep * ( e3t_0(ji,jj,jk) / fse3t_n(ji,jj,jk) ) 
     297               END DO 
     298            END DO   
     299         END DO   
     300      ENDIF 
     301#endif 
    289302  
    290303      IF( l_trdtra )   THEN                      ! send trends for further diagnostics 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRD/trdken.F90

    r7730 r7731  
    117117      ! 
    118118      SELECT CASE( ktrd ) 
    119          CASE( jpdyn_hpg )   ;   CALL iom_put( "ketrd_hpg", zke )    ! hydrostatic pressure gradient 
    120          CASE( jpdyn_spg )   ;   CALL iom_put( "ketrd_spg", zke )    ! surface pressure gradient 
    121          CASE( jpdyn_spgexp );   CALL iom_put( "ketrd_spgexp", zke ) ! surface pressure gradient (explicit) 
    122          CASE( jpdyn_spgflt );   CALL iom_put( "ketrd_spgflt", zke ) ! surface pressure gradient (filter) 
    123          CASE( jpdyn_pvo )   ;   CALL iom_put( "ketrd_pvo", zke )    ! planetary vorticity 
    124          CASE( jpdyn_rvo )   ;   CALL iom_put( "ketrd_rvo", zke )    ! relative  vorticity     (or metric term) 
    125          CASE( jpdyn_keg )   ;   CALL iom_put( "ketrd_keg", zke )    ! Kinetic Energy gradient (or had) 
    126          CASE( jpdyn_zad )   ;   CALL iom_put( "ketrd_zad", zke )    ! vertical   advection 
    127          CASE( jpdyn_ldf )   ;   CALL iom_put( "ketrd_ldf", zke )    ! lateral diffusion 
    128          CASE( jpdyn_zdf )   ;   CALL iom_put( "ketrd_zdf", zke )    ! vertical diffusion  
     119        CASE( jpdyn_hpg )   ;   CALL iom_put( "ketrd_hpg", zke )    ! hydrostatic pressure gradient 
     120        CASE( jpdyn_spg )   ;   CALL iom_put( "ketrd_spg", zke )    ! surface pressure gradient 
     121        CASE( jpdyn_spgexp );   CALL iom_put( "ketrd_spgexp", zke ) ! surface pressure gradient (explicit) 
     122        CASE( jpdyn_spgflt );   CALL iom_put( "ketrd_spgflt", zke ) ! surface pressure gradient (filter) 
     123        CASE( jpdyn_pvo )   ;   CALL iom_put( "ketrd_pvo", zke )    ! planetary vorticity 
     124        CASE( jpdyn_rvo )   ;   CALL iom_put( "ketrd_rvo", zke )    ! relative  vorticity     (or metric term) 
     125        CASE( jpdyn_keg )   ;   CALL iom_put( "ketrd_keg", zke )    ! Kinetic Energy gradient (or had) 
     126        CASE( jpdyn_zad )   ;   CALL iom_put( "ketrd_zad", zke )    ! vertical   advection 
     127        CASE( jpdyn_ldf )   ;   CALL iom_put( "ketrd_ldf", zke )    ! lateral diffusion 
     128        CASE( jpdyn_zdf )   ;   CALL iom_put( "ketrd_zdf", zke )    ! vertical diffusion  
    129129                                 !                                   ! wind stress trends 
    130                                  CALL wrk_alloc( jpi, jpj, z2dx, z2dy, zke2d ) 
    131                            z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1u(:,:) * e2u(:,:) * umask(:,:,1) 
    132                            z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1v(:,:) * e2v(:,:) * vmask(:,:,1) 
    133                            zke2d(1,:) = 0._wp   ;   zke2d(:,1) = 0._wp 
    134                            DO jj = 2, jpj 
    135                               DO ji = 2, jpi 
    136                                  zke2d(ji,jj) = 0.5_wp * (   z2dx(ji,jj) + z2dx(ji-1,jj)   & 
    137                                  &                         + z2dy(ji,jj) + z2dy(ji,jj-1)   ) * r1_bt(ji,jj,1) 
    138                               END DO 
    139                            END DO 
    140                                  CALL iom_put( "ketrd_tau", zke2d ) 
    141                                  CALL wrk_dealloc( jpi, jpj     , z2dx, z2dy, zke2d ) 
    142          CASE( jpdyn_bfr )   ;   CALL iom_put( "ketrd_bfr", zke )    ! bottom friction (explicit case)  
     130                                CALL wrk_alloc( jpi, jpj, z2dx, z2dy, zke2d ) 
     131                     z2dx(:,:) = un(:,:,1) * ( utau_b(:,:) + utau(:,:) ) * e1u(:,:) * e2u(:,:) * umask(:,:,1) 
     132                     z2dy(:,:) = vn(:,:,1) * ( vtau_b(:,:) + vtau(:,:) ) * e1v(:,:) * e2v(:,:) * vmask(:,:,1) 
     133                     zke2d(1,:) = 0._wp   ;   zke2d(:,1) = 0._wp 
     134                     DO jj = 2, jpj 
     135                         DO ji = 2, jpi 
     136                           zke2d(ji,jj) = 0.5_wp * (   z2dx(ji,jj) + z2dx(ji-1,jj)   & 
     137                            &                         + z2dy(ji,jj) + z2dy(ji,jj-1)   ) * r1_bt(ji,jj,1) 
     138                         END DO 
     139                     END DO 
     140                                CALL iom_put( "ketrd_tau", zke2d ) 
     141                                CALL wrk_dealloc( jpi, jpj     , z2dx, z2dy, zke2d ) 
     142        CASE( jpdyn_bfr )   ;   CALL iom_put( "ketrd_bfr", zke )    ! bottom friction (explicit case)  
    143143!!gm TO BE DONE properly 
    144144!!gm only valid if ln_bfrimp=F otherwise the bottom stress as to be recomputed at the end of the computation.... 
    … …  
    162162!         ENDIF 
    163163!!gm end 
    164          CASE( jpdyn_atf )   ;   CALL iom_put( "ketrd_atf", zke )    ! asselin filter trends  
     164        CASE( jpdyn_atf )   ;   CALL iom_put( "ketrd_atf", zke )    ! asselin filter trends  
    165165!! a faire !!!!  idee changer dynnxt pour avoir un appel a jpdyn_bfr avant le swap !!! 
    166166!! reflechir a une possible sauvegarde du "vrai" un,vn pour le calcul de atf.... 
    … …  
    184184!                              CALL iom_put( "ketrd_bfri", zke2d ) 
    185185!         ENDIF 
    186          CASE( jpdyn_ken )   ;   ! kinetic energy 
    187                            ! called in dynnxt.F90 before asselin time filter 
    188                            ! with putrd=ua and pvtrd=va 
    189                            zke(:,:,:) = 0.5_wp * zke(:,:,:) 
    190                            CALL iom_put( "KE", zke ) 
    191                            ! 
    192                            CALL ken_p2k( kt , zke ) 
    193                            CALL iom_put( "ketrd_convP2K", zke )     ! conversion -rau*g*w 
     186        CASE( jpdyn_ken )   ;   ! kinetic energy 
     187                    ! called in dynnxt.F90 before asselin time filter 
     188                    ! with putrd=ua and pvtrd=va 
     189                    zke(:,:,:) = 0.5_wp * zke(:,:,:) 
     190                    CALL iom_put( "KE", zke ) 
     191                    ! 
     192                    CALL ken_p2k( kt , zke ) 
     193                      CALL iom_put( "ketrd_convP2K", zke )     ! conversion -rau*g*w 
    194194         ! 
    195195      END SELECT 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRD/trdmxl.F90

    r7730 r7731  
    165165 
    166166 
    167       SELECT CASE( ktrd ) 
    168       CASE( jptra_npc  )               ! non-penetrative convection: regrouped with zdf 
     167      SELECT CASE( ktrd ) 
     168      CASE( jptra_npc  )               ! non-penetrative convection: regrouped with zdf 
    169169!!gm : to be completed !  
    170 !        IF( .... 
     170!         IF( .... 
    171171!!gm end 
    172       CASE( jptra_zdfp )               ! iso-neutral diffusion: "pure" vertical diffusion 
    173          !                                   ! regroup iso-neutral diffusion in one term 
     172      CASE( jptra_zdfp )               ! iso-neutral diffusion: "pure" vertical diffusion 
     173         !                                   ! regroup iso-neutral diffusion in one term 
    174174         tmltrd(:,:,jpmxl_ldf) = tmltrd(:,:,jpmxl_ldf) + ( tmltrd(:,:,jpmxl_zdf) - tmltrd(:,:,jpmxl_zdfp) ) 
    175175         smltrd(:,:,jpmxl_ldf) = smltrd(:,:,jpmxl_ldf) + ( smltrd(:,:,jpmxl_zdf) - smltrd(:,:,jpmxl_zdfp) ) 
    … …  
    811811 
    812812 
    813       nkstp     = nit000 - 1              ! current time step indicator initialization 
     813      nkstp     = nit000 - 1              ! current time step indicator initialization 
    814814 
    815815 
    … …  
    851851      IF( nn_ctls == 1 ) THEN 
    852852         CALL ctl_opn( inum, 'ctlsurf_idx', 'OLD', 'UNFORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 
    853          READ ( inum ) nbol 
     853         READ ( inum, * ) nbol 
    854854         CLOSE( inum ) 
    855855      END IF 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRD/trdmxl_oce.F90

    r7730 r7731  
    1515 
    1616   !                                                !* mixed layer trend indices 
    17    INTEGER, PUBLIC, PARAMETER ::   jpltrd = 11      !: number of mixed-layer trends arrays 
     17   INTEGER, PUBLIC, PARAMETER ::   jpltrd = 12      !: number of mixed-layer trends arrays 
    1818   INTEGER, PUBLIC            ::   jpktrd           !: max level for mixed-layer trends diag. 
    1919   ! 
    … …  
    2828   INTEGER, PUBLIC, PARAMETER ::   jpmxl_for =  9   !: forcing  
    2929   INTEGER, PUBLIC, PARAMETER ::   jpmxl_dmp = 10   !: internal restoring trend 
    30    INTEGER, PUBLIC, PARAMETER ::   jpmxl_zdfp = 11   !: asselin trend (**MUST BE THE LAST ONE**) 
    31    INTEGER, PUBLIC, PARAMETER ::   jpmxl_atf  = 12   !: asselin trend (**MUST BE THE LAST ONE**) 
     30   INTEGER, PUBLIC, PARAMETER ::   jpmxl_zdfp = 11  !: iso-neutral diffusion:"pure" vertical diffusion 
     31   INTEGER, PUBLIC, PARAMETER ::   jpmxl_atf  = 12  !: asselin trend (**MUST BE THE LAST ONE**) 
    3232   !                                                            !!* Namelist namtrd_mxl:  trend diagnostics in the mixed layer * 
    3333   INTEGER           , PUBLIC ::   nn_ctls  = 0                  !: control surface type for trends vertical integration 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/TRD/trdpen.F90

    r7730 r7731  
    9999                                   CALL wrk_alloc( jpi, jpj, z2d ) 
    100100                                   z2d(:,:) = wn(:,:,1) * ( & 
    101                                       &   - ( rab_n(:,:,1,jp_tem) + rab_pe(:,:,1,jp_tem) ) * tsn(:,:,1,jp_tem)    & 
    102                                       &   + ( rab_n(:,:,1,jp_sal) + rab_pe(:,:,1,jp_sal) ) * tsn(:,:,1,jp_sal)    & 
    103                                       &                  ) / fse3t(:,:,1) 
     101                                       &   - ( rab_n(:,:,1,jp_tem) + rab_pe(:,:,1,jp_tem) ) * tsn(:,:,1,jp_tem)    & 
     102                                       &   + ( rab_n(:,:,1,jp_sal) + rab_pe(:,:,1,jp_sal) ) * tsn(:,:,1,jp_sal)    & 
     103                                       &             ) / fse3t(:,:,1) 
    104104                                   CALL iom_put( "petrd_sad" , z2d ) 
    105105                                   CALL wrk_dealloc( jpi, jpj, z2d ) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90

    r7730 r7731  
    4343   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   avmu , avmv    !: vertical viscosity coef at uw- & vw-pts       [m2/s] 
    4444   REAL(wp), PUBLIC, SAVE, ALLOCATABLE, DIMENSION(:,:,:) ::   avm  , avt     !: vertical viscosity & diffusivity coef at w-pt [m2/s] 
     45   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k , avm_k  ! not enhanced Kz 
     46   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k, avmv_k ! not enhanced Kz 
     47   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en              !: now turbulent kinetic energy   [m2/s2] 
    4548  
    4649   !!---------------------------------------------------------------------- 
    … …  
    6063         &     tfrua(jpi, jpj), tfrva(jpi, jpj)              ,      & 
    6164         &     avmu(jpi,jpj,jpk), avm(jpi,jpj,jpk)           ,      & 
    62          &     avmv(jpi,jpj,jpk), avt(jpi,jpj,jpk)           , STAT = zdf_oce_alloc ) 
     65         &     avmv  (jpi,jpj,jpk), avt   (jpi,jpj,jpk)      ,      & 
     66         &     avt_k (jpi,jpj,jpk), avm_k (jpi,jpj,jpk)      ,      &  
     67         &     avmu_k(jpi,jpj,jpk), avmv_k(jpi,jpj,jpk)      ,      & 
     68         &     en    (jpi,jpj,jpk), STAT = zdf_oce_alloc ) 
    6369         ! 
    6470      IF( zdf_oce_alloc /= 0 )   CALL ctl_warn('zdf_oce_alloc: failed to allocate arrays') 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90

    r7730 r7731  
    4242   LOGICAL , PUBLIC, PARAMETER ::   lk_zdfgls = .TRUE.   !: TKE vertical mixing flag 
    4343   ! 
    44    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en      !: now turbulent kinetic energy 
    4544   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   mxln    !: now mixing length 
    4645   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   zwall   !: wall function 
    47    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k   ! not enhanced Kz 
    48    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avm_k   ! not enhanced Kz 
    49    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k  ! not enhanced Kz 
    50    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmv_k  ! not enhanced Kz 
    5146   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ustars2 !: Squared surface velocity scale at T-points 
    5247   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   ustarb2 !: Squared bottom  velocity scale at T-points 
    … …  
    120115      !!                ***  FUNCTION zdf_gls_alloc  *** 
    121116      !!---------------------------------------------------------------------- 
    122       ALLOCATE( en(jpi,jpj,jpk),  mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) ,     & 
    123          &      avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk),                    & 
    124          &      avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk),                    & 
    125          &      ustars2(jpi,jpj), ustarb2(jpi,jpj)                      , STAT= zdf_gls_alloc ) 
     117      ALLOCATE( mxln(jpi,jpj,jpk), zwall(jpi,jpj,jpk) ,     & 
     118         &      ustars2(jpi,jpj) , ustarb2(jpi,jpj)   , STAT= zdf_gls_alloc ) 
    126119         ! 
    127120      IF( lk_mpp             )   CALL mpp_sum ( zdf_gls_alloc ) 
    … …  
    329322      !  
    330323      ! One level below 
    331       en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+fsdepw(:,:,2))/zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 
     324      en(:,:,2) = rc02r * ustars2(:,:) * (1._wp + rsbc_tke1 * ((zhsro(:,:)+fsdepw(:,:,2)) & 
     325          &            / zhsro(:,:) )**(1.5_wp*ra_sf))**(2._wp/3._wp) 
    332326      en(:,:,2) = MAX(en(:,:,2), rn_emin ) 
    333327      z_elem_a(:,:,2) = 0._wp  
    … …  
    350344      z_elem_a(:,:,2) = 0._wp 
    351345      zkar(:,:)       = (rl_sf + (vkarmn-rl_sf)*(1.-exp(-rtrans*fsdept(:,:,1)/zhsro(:,:)) )) 
    352       zflxs(:,:)      = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 
     346      zflxs(:,:)      = rsbc_tke2 * ustars2(:,:)**1.5_wp * zkar(:,:) & 
     347           &                      * ((zhsro(:,:)+fsdept(:,:,1))/zhsro(:,:) )**(1.5_wp*ra_sf) 
    353348 
    354349      en(:,:,2) = en(:,:,2) + zflxs(:,:)/fse3w(:,:,2) 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfmxl.F90

    r7730 r7731  
    2222   USE timing          ! Timing 
    2323   USE trc_oce, ONLY : lk_offline ! offline flag 
     24   USE lbclnk  
     25   USE eosbn2          ! Equation of state 
    2426 
    2527   IMPLICIT NONE 
    … …  
    2729 
    2830   PUBLIC   zdf_mxl       ! called by step.F90 
     31   PUBLIC   zdf_mxl_tref  ! called by asminc.F90 
     32   PUBLIC   zdf_mxl_alloc ! Used in zdf_tke_init 
    2933 
    3034   INTEGER , PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   nmln    !: number of level in the mixed layer (used by TOP) 
    … …  
    3236   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlp    !: mixed layer depth  (rho=rho0+zdcrit) [m] 
    3337   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmlpt   !: mixed layer depth at t-points        [m] 
     38   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmld_tref  !: mixed layer depth at t-points - temperature criterion [m] 
     39   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hmld_kara  !: mixed layer depth of Kara et al   [m] 
    3440 
    3541   REAL(wp), PUBLIC ::   rho_c = 0.01_wp    !: density criterion for mixed layer depth 
    3642   REAL(wp)         ::   avt_c = 5.e-4_wp   ! Kz criterion for the turbocline depth 
     43 
     44   ! Namelist variables for  namzdf_karaml  
     45  
     46   LOGICAL   :: ln_kara            ! Logical switch for calculating kara mixed 
     47                                     ! layer 
     48   LOGICAL   :: ln_kara_write25h   ! Logical switch for 25 hour outputs 
     49   INTEGER   :: jpmld_type         ! mixed layer type              
     50   REAL(wp)  :: ppz_ref            ! depth of initial T_ref  
     51   REAL(wp)  :: ppdT_crit          ! Critical temp diff   
     52   REAL(wp)  :: ppiso_frac         ! Fraction of ppdT_crit used  
     53    
     54   !Used for 25h mean 
     55   LOGICAL, PRIVATE :: kara_25h_init = .TRUE.   !Logical used to initalise 25h  
     56                                                !outputs. Necissary, because we need to 
     57                                                !initalise the kara_25h on the zeroth 
     58                                                !timestep (i.e in the nemogcm_init call) 
     59   REAL, PRIVATE, ALLOCATABLE, DIMENSION(:,:) :: hmld_kara_25h 
    3760 
    3861   !! * Substitutions 
    … …  
    5174      zdf_mxl_alloc = 0      ! set to zero if no array to be allocated 
    5275      IF( .NOT. ALLOCATED( nmln ) ) THEN 
    53          ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), STAT= zdf_mxl_alloc ) 
     76         ALLOCATE( nmln(jpi,jpj), hmld(jpi,jpj), hmlp(jpi,jpj), hmlpt(jpi,jpj), & 
     77         &                           hmld_tref(jpi,jpj), STAT= zdf_mxl_alloc ) 
    5478         ! 
    5579         IF( lk_mpp             )   CALL mpp_sum ( zdf_mxl_alloc ) 
    … …  
    122146      END DO 
    123147      ! depth of the mixing and mixed layers 
     148 
     149      CALL zdf_mxl_kara( kt ) ! kara MLD 
     150 
    124151      DO jj = 1, jpj 
    125152         DO ji = 1, jpi 
    … …  
    127154            iikn = nmln(ji,jj) 
    128155            imkt = mikt(ji,jj) 
    129             hmld (ji,jj) = ( fsdepw(ji,jj,iiki  ) - fsdepw(ji,jj,imkt )            )  * ssmask(ji,jj)    ! Turbocline depth  
    130             hmlp (ji,jj) = ( fsdepw(ji,jj,iikn  ) - fsdepw(ji,jj,MAX( imkt,nla10 ) ) ) * ssmask(ji,jj)    ! Mixed layer depth 
    131             hmlpt(ji,jj) = ( fsdept(ji,jj,iikn-1) - fsdepw(ji,jj,imkt )            )  * ssmask(ji,jj)    ! depth of the last T-point inside the mixed layer 
     156            hmld (ji,jj) = ( fsdepw(ji,jj,iiki  ) - fsdepw(ji,jj,imkt ) ) * ssmask(ji,jj)    ! Turbocline depth  
     157            hmlp (ji,jj) = ( fsdepw(ji,jj,iikn  ) - fsdepw(ji,jj,imkt ) ) * ssmask(ji,jj)    ! Mixed layer depth 
     158            hmlpt(ji,jj) = ( fsdept(ji,jj,iikn-1) - fsdepw(ji,jj,imkt ) ) * ssmask(ji,jj)    ! depth of the last T-point inside the mixed layer 
    132159         END DO 
    133160      END DO 
    … …  
    145172   END SUBROUTINE zdf_mxl 
    146173 
     174 
     175   SUBROUTINE zdf_mxl_tref() 
     176      !!---------------------------------------------------------------------- 
     177      !!                  ***  ROUTINE zdf_mxl_tref  *** 
     178      !!                    
     179      !! ** Purpose :   Compute the mixed layer depth with temperature criteria. 
     180      !! 
     181      !! ** Method  :   The temperature-defined mixed layer depth is required 
     182      !!                   when assimilating SST in a 2D analysis.  
     183      !! 
     184      !! ** Action  :   hmld_tref 
     185      !!---------------------------------------------------------------------- 
     186      ! 
     187      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     188      REAL(wp) ::   t_ref               ! Reference temperature   
     189      REAL(wp) ::   temp_c = 0.2        ! temperature criterion for mixed layer depth   
     190      !!---------------------------------------------------------------------- 
     191      ! 
     192      ! Initialise array 
     193      IF( zdf_mxl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'zdf_mxl_tref : unable to allocate arrays' ) 
     194       
     195      !For the AMM model assimiation uses a temperature based mixed layer depth   
     196      !This is defined here   
     197      DO jj = 1, jpj   
     198         DO ji = 1, jpi   
     199           hmld_tref(ji,jj)=fsdept(ji,jj,1  )    
     200           IF(ssmask(ji,jj) > 0.)THEN   
     201             t_ref=tsn(ji,jj,1,jp_tem)  
     202             DO jk=2,jpk   
     203               IF(ssmask(ji,jj)==0.)THEN   
     204                  hmld_tref(ji,jj)=fsdept(ji,jj,jk )   
     205                  EXIT   
     206               ELSEIF( ABS(tsn(ji,jj,jk,jp_tem)-t_ref) < temp_c)THEN   
     207                  hmld_tref(ji,jj)=fsdept(ji,jj,jk )   
     208               ELSE   
     209                  EXIT   
     210               ENDIF   
     211             ENDDO   
     212           ENDIF   
     213         ENDDO   
     214      ENDDO 
     215 
     216   END SUBROUTINE zdf_mxl_tref 
     217 
     218   SUBROUTINE zdf_mxl_kara( kt )  
     219      !!----------------------------------------------------------------------------------  
     220      !!                    ***  ROUTINE zdf_mxl_kara  ***  
     221      !                                                                         
     222      !   Calculate mixed layer depth according to the definition of           
     223      !   Kara et al, 2000, JGR, 105, 16803.   
     224      !  
     225      !   If mld_type=1 the mixed layer depth is calculated as the depth at which the   
     226      !   density has increased by an amount equivalent to a temperature difference of   
     227      !   0.8C at the surface.  
     228      !  
     229      !   For other values of mld_type the mixed layer is calculated as the depth at   
     230      !   which the temperature differs by 0.8C from the surface temperature.   
     231      !                                                                         
     232      !   Original version: David Acreman                                       
     233      !  
     234      !!----------------------------------------------------------------------------------- 
     235      
     236      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index  
     237  
     238      NAMELIST/namzdf_karaml/ ln_kara,jpmld_type, ppz_ref, ppdT_crit, & 
     239      &                       ppiso_frac, ln_kara_write25h  
     240  
     241      ! Local variables                                                         
     242      REAL, DIMENSION(jpi,jpk) :: ppzdep      ! depth for use in calculating d(rho)  
     243      REAL(wp), DIMENSION(jpi,jpj,jpts) :: ztsn_2d  !Local version of tsn  
     244      LOGICAL :: ll_found(jpi,jpj)              ! Is T_b to be found by interpolation ?  
     245      LOGICAL :: ll_belowml(jpi,jpj,jpk)        ! Flag points below mixed layer when ll_found=F  
     246      INTEGER :: ji, jj, jk                     ! loop counter  
     247      INTEGER :: ik_ref(jpi,jpj)                ! index of reference level  
     248      INTEGER :: ik_iso(jpi,jpj)                ! index of last uniform temp level  
     249      REAL    :: zT(jpi,jpj,jpk)                ! Temperature or denisty  
     250      REAL    :: zT_ref(jpi,jpj)                ! reference temperature  
     251      REAL    :: zT_b                           ! base temperature  
     252      REAL    :: zdTdz(jpi,jpj,jpk-2)           ! gradient of zT  
     253      REAL    :: zmoddT(jpi,jpj,jpk-2)          ! Absolute temperature difference  
     254      REAL    :: zdz                            ! depth difference  
     255      REAL    :: zdT                            ! temperature difference  
     256      REAL    :: zdelta_T(jpi,jpj)              ! difference critereon  
     257      REAL    :: zRHO1(jpi,jpj), zRHO2(jpi,jpj) ! Densities 
     258      INTEGER, SAVE :: idt, i_steps 
     259      INTEGER, SAVE :: i_cnt_25h     ! Counter for 25 hour means 
     260      INTEGER :: ios                 ! Local integer output status for namelist read 
     261 
     262      
     263  
     264      !!-------------------------------------------------------------------------------------  
     265  
     266      IF( kt == nit000 ) THEN  
     267         ! Default values  
     268         ln_kara      = .FALSE. 
     269         ln_kara_write25h = .FALSE.  
     270         jpmld_type   = 1      
     271         ppz_ref      = 10.0  
     272         ppdT_crit    = 0.2   
     273         ppiso_frac   = 0.1    
     274         ! Read namelist 
     275         REWIND ( numnam_ref ) 
     276         READ   ( numnam_ref, namzdf_karaml, IOSTAT=ios, ERR= 901 ) 
     277901      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_karaml in reference namelist', lwp ) 
     278 
     279         REWIND( numnam_cfg )              ! Namelist nam_diadiaop in configuration namelist  3D hourly diagnostics 
     280         READ  ( numnam_cfg,  namzdf_karaml, IOSTAT = ios, ERR = 902 ) 
     281902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namzdf_karaml in configuration namelist', lwp ) 
     282         IF(lwm) WRITE ( numond, namzdf_karaml ) 
     283  
     284 
     285         WRITE(numout,*) '===== Kara mixed layer ====='  
     286         WRITE(numout,*) 'ln_kara = ',    ln_kara 
     287         WRITE(numout,*) 'jpmld_type = ', jpmld_type  
     288         WRITE(numout,*) 'ppz_ref = ',    ppz_ref  
     289         WRITE(numout,*) 'ppdT_crit = ',  ppdT_crit  
     290         WRITE(numout,*) 'ppiso_frac = ', ppiso_frac 
     291         WRITE(numout,*) 'ln_kara_write25h = ', ln_kara_write25h 
     292         WRITE(numout,*) '============================'  
     293       
     294         IF ( .NOT.ln_kara ) THEN 
     295            WRITE(numout,*) "ln_kara not set; Kara mixed layer not calculated"  
     296         ELSE 
     297            IF (.NOT. ALLOCATED(hmld_kara) ) ALLOCATE( hmld_kara(jpi,jpj) ) 
     298            IF ( ln_kara_write25h .AND. kara_25h_init ) THEN 
     299               i_cnt_25h = 0 
     300               IF (.NOT. ALLOCATED(hmld_kara_25h) ) & 
     301               &   ALLOCATE( hmld_kara_25h(jpi,jpj) ) 
     302               hmld_kara_25h = 0._wp 
     303               IF( nacc == 1 ) THEN 
     304                  idt = INT(rdtmin) 
     305               ELSE 
     306                  idt = INT(rdt) 
     307               ENDIF 
     308               IF( MOD( 3600,idt) == 0 ) THEN  
     309                  i_steps = 3600 / idt   
     310               ELSE  
     311                  CALL ctl_stop('STOP', & 
     312                  & 'zdf_mxl_kara: timestep must give MOD(3600,rdt)'// & 
     313                  & ' = 0 otherwise no hourly values are possible')  
     314               ENDIF   
     315            ENDIF 
     316         ENDIF 
     317      ENDIF 
     318       
     319      IF ( ln_kara ) THEN 
     320          
     321         !set critical ln_kara 
     322         ztsn_2d = tsn(:,:,1,:) 
     323         ztsn_2d(:,:,jp_tem) = ztsn_2d(:,:,jp_tem) + ppdT_crit 
     324      
     325         ! Set the mixed layer depth criterion at each grid point  
     326         ppzdep = 0._wp 
     327         IF( jpmld_type == 1 ) THEN                                          
     328            CALL eos ( tsn(:,:,1,:), & 
     329            &          ppzdep(:,:), zRHO1(:,:) )  
     330            CALL eos ( ztsn_2d(:,:,:), & 
     331            &           ppzdep(:,:), zRHO2(:,:) )  
     332            zdelta_T(:,:) = abs( zRHO1(:,:) - zRHO2(:,:) ) * rau0  
     333            ! RHO from eos (2d version) doesn't calculate north or east halo:  
     334            CALL lbc_lnk( zdelta_T, 'T', 1. )  
     335            zT(:,:,:) = rhop(:,:,:)  
     336         ELSE  
     337            zdelta_T(:,:) = ppdT_crit                       
     338            zT(:,:,:) = tsn(:,:,:,jp_tem)                            
     339         ENDIF  
     340      
     341         ! Calculate the gradient of zT and absolute difference for use later  
     342         DO jk = 1 ,jpk-2  
     343            zdTdz(:,:,jk)  =    ( zT(:,:,jk+1) - zT(:,:,jk) ) / fse3w(:,:,jk+1)  
     344            zmoddT(:,:,jk) = abs( zT(:,:,jk+1) - zT(:,:,jk) )  
     345         END DO  
     346      
     347         ! Find density/temperature at the reference level (Kara et al use 10m).           
     348         ! ik_ref is the index of the box centre immediately above or at the reference level  
     349         ! Find ppz_ref in the array of model level depths and find the ref     
     350         ! density/temperature by linear interpolation.                                    
     351         ik_ref = -1 
     352         DO jk = jpkm1, 2, -1  
     353            WHERE( fsdept(:,:,jk) > ppz_ref )  
     354               ik_ref(:,:) = jk - 1  
     355               zT_ref(:,:) = zT(:,:,jk-1) + & 
     356               &             zdTdz(:,:,jk-1) * ( ppz_ref - fsdept(:,:,jk-1) )  
     357            ENDWHERE  
     358         END DO 
     359         IF ( ANY( ik_ref  < 0 ) .OR. ANY( ik_ref  > jpkm1 ) ) THEN 
     360            CALL ctl_stop( "STOP", & 
     361            & "zdf_mxl_kara: unable to find reference level for kara ML" )  
     362         ELSE 
     363            ! If the first grid box centre is below the reference level then use the  
     364            ! top model level to get zT_ref  
     365            WHERE( fsdept(:,:,1) > ppz_ref )   
     366               zT_ref = zT(:,:,1)  
     367               ik_ref = 1  
     368            ENDWHERE  
     369      
     370            ! Search for a uniform density/temperature region where adjacent levels           
     371            ! differ by less than ppiso_frac * deltaT.                                       
     372            ! ik_iso is the index of the last level in the uniform layer   
     373            ! ll_found indicates whether the mixed layer depth can be found by interpolation  
     374            ik_iso(:,:)   = ik_ref(:,:)  
     375            ll_found(:,:) = .false.  
     376            DO jj = 1, nlcj  
     377               DO ji = 1, nlci  
     378                 !CDIR NOVECTOR  
     379                  DO jk = ik_ref(ji,jj), mbathy(ji,jj)-1  
     380                     IF( zmoddT(ji,jj,jk) > ( ppiso_frac * zdelta_T(ji,jj) ) ) THEN  
     381                        ik_iso(ji,jj)   = jk  
     382                        ll_found(ji,jj) = ( zmoddT(ji,jj,jk) > zdelta_T(ji,jj) )  
     383                        EXIT  
     384                     ENDIF  
     385                  END DO  
     386               END DO  
     387            END DO  
     388      
     389            ! Use linear interpolation to find depth of mixed layer base where possible  
     390            hmld_kara(:,:) = ppz_ref  
     391            DO jj = 1, jpj  
     392               DO ji = 1, jpi  
     393                  IF( ll_found(ji,jj) .and. tmask(ji,jj,1) == 1.0 ) THEN  
     394                     zdz =  abs( zdelta_T(ji,jj) / zdTdz(ji,jj,ik_iso(ji,jj)) )  
     395                     hmld_kara(ji,jj) = fsdept(ji,jj,ik_iso(ji,jj)) + zdz  
     396                  ENDIF  
     397               END DO  
     398            END DO  
     399      
     400            ! If ll_found = .false. then calculate MLD using difference of zdelta_T     
     401            ! from the reference density/temperature  
     402      
     403            ! Prevent this section from working on land points  
     404            WHERE( tmask(:,:,1) /= 1.0 )  
     405               ll_found = .true.  
     406            ENDWHERE  
     407      
     408            DO jk = 1, jpk  
     409               ll_belowml(:,:,jk) = abs( zT(:,:,jk) - zT_ref(:,:) ) >= & 
     410               & zdelta_T(:,:) 
     411            END DO  
     412      
     413            ! Set default value where interpolation cannot be used (ll_found=false)   
     414            DO jj = 1, jpj  
     415               DO ji = 1, jpi  
     416                  IF( .NOT. ll_found(ji,jj) )  & 
     417                  &   hmld_kara(ji,jj) = fsdept(ji,jj,mbathy(ji,jj))  
     418               END DO  
     419            END DO  
     420      
     421            DO jj = 1, jpj  
     422               DO ji = 1, jpi  
     423                  !CDIR NOVECTOR  
     424                  DO jk = ik_ref(ji,jj)+1, mbathy(ji,jj)  
     425                     IF( ll_found(ji,jj) ) EXIT  
     426                     IF( ll_belowml(ji,jj,jk) ) THEN                 
     427                        zT_b = zT_ref(ji,jj) + zdelta_T(ji,jj) * & 
     428                        &      SIGN(1.0, zdTdz(ji,jj,jk-1) )  
     429                        zdT  = zT_b - zT(ji,jj,jk-1)                                       
     430                        zdz  = zdT / zdTdz(ji,jj,jk-1)                                        
     431                        hmld_kara(ji,jj) = fsdept(ji,jj,jk-1) + zdz  
     432                        EXIT                                                    
     433                     ENDIF  
     434                  END DO  
     435               END DO  
     436            END DO  
     437      
     438            hmld_kara(:,:) = hmld_kara(:,:) * tmask(:,:,1)  
     439  
     440            IF(  ln_kara_write25h  ) THEN 
     441               !Double IF required as i_steps not defined if ln_kara_write25h = 
     442               ! FALSE 
     443               IF ( ( MOD( kt, i_steps ) == 0 ) .OR.  kara_25h_init ) THEN 
     444                  hmld_kara_25h = hmld_kara_25h + hmld_kara 
     445                  i_cnt_25h = i_cnt_25h + 1 
     446                  IF ( kara_25h_init ) kara_25h_init = .FALSE. 
     447               ENDIF 
     448            ENDIF 
     449  
     450#if defined key_iomput  
     451            IF( kt /= nit000 ) THEN  
     452               CALL iom_put( "mldkara"  , hmld_kara )    
     453               IF( ( MOD( i_cnt_25h, 25) == 0 ) .AND.  ln_kara_write25h ) & 
     454                  CALL iom_put( "kara25h"  , ( hmld_kara_25h / 25._wp ) ) 
     455            ENDIF 
     456#endif 
     457  
     458         ENDIF 
     459      ENDIF 
     460        
     461   END SUBROUTINE zdf_mxl_kara  
     462 
    147463   !!====================================================================== 
    148464END MODULE zdfmxl 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90

    r7730 r7731  
    5353   USE timing         ! Timing 
    5454   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
     55#if defined key_agrif 
     56   USE agrif_opa_interp 
     57   USE agrif_opa_update 
     58#endif 
     59 
     60 
    5561 
    5662   IMPLICIT NONE 
    … …  
    8591   REAL(wp) ::   rhftau_scl = 1.0_wp       ! scale factor applied to HF part of taum  (nn_etau=3) 
    8692 
    87    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   en             !: now turbulent kinetic energy   [m2/s2] 
    8893   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   htau           ! depth of tke penetration (nn_htau) 
    8994   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   dissl          ! now mixing lenght of dissipation 
    90    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avt_k , avm_k  ! not enhanced Kz 
    91    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   avmu_k, avmv_k ! not enhanced Kz 
    9295#if defined key_c1d 
    9396   !                                                                        !!** 1D cfg only  **   ('key_c1d') 
    … …  
    115118         &      e_pdl(jpi,jpj,jpk) , e_ric(jpi,jpj,jpk) ,                          & 
    116119#endif 
    117          &      en    (jpi,jpj,jpk) , htau  (jpi,jpj)    , dissl(jpi,jpj,jpk) ,     &  
    118          &      avt_k (jpi,jpj,jpk) , avm_k (jpi,jpj,jpk),                          & 
    119          &      avmu_k(jpi,jpj,jpk) , avmv_k(jpi,jpj,jpk), STAT= zdf_tke_alloc      ) 
     120         &      htau  (jpi,jpj)    , dissl(jpi,jpj,jpk) , STAT= zdf_tke_alloc      ) 
    120121         ! 
    121122      IF( lk_mpp             )   CALL mpp_sum ( zdf_tke_alloc ) 
    … …  
    189190      avmv_k(:,:,:) = avmv(:,:,:)  
    190191      ! 
     192#if defined key_agrif 
     193      ! Update child grid f => parent grid  
     194      IF( .NOT.Agrif_Root() )   CALL Agrif_Update_Tke( kt )      ! children only 
     195#endif       
     196     !  
    191197   END SUBROUTINE zdf_tke 
    192198 
    … …  
    317323                  zwlc = zind * rn_lc * zus * SIN( rpi * fsdepw(ji,jj,jk) / zhlc(ji,jj) ) 
    318324                  !                                           ! TKE Langmuir circulation source term 
    319                   en(ji,jj,jk) = en(ji,jj,jk) + rdt * ( zwlc * zwlc * zwlc ) / zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
     325                  en(ji,jj,jk) = en(ji,jj,jk) + rdt * ( 1._wp - fr_i(ji,jj) ) * ( zwlc * zwlc * zwlc ) /   & 
     326                     &   zhlc(ji,jj) * wmask(ji,jj,jk) * tmask(ji,jj,1) 
    320327               END DO 
    321328            END DO 
    … …  
    710717      !!---------------------------------------------------------------------- 
    711718      INTEGER ::   ji, jj, jk   ! dummy loop indices 
    712       INTEGER ::   ios 
     719      INTEGER ::   ios, ierr 
    713720      !! 
    714721      NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin  ,   & 
    … …  
    768775      ENDIF 
    769776       
    770       IF( nn_etau == 2  )   CALL zdf_mxl( nit000 )      ! Initialization of nmln  
     777      IF( nn_etau == 2  ) THEN 
     778          ierr = zdf_mxl_alloc() 
     779          nmln(:,:) = nlb10           ! Initialization of nmln 
     780      ENDIF 
    771781 
    772782      !                               !* depth of penetration of surface tke 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r7730 r7731  
    6161   USE asminc          ! assimilation increments      
    6262   USE asmbkg          ! writing out state trajectory 
     63   USE asmbal          ! writing out assimilation balancing increments 
    6364   USE diaptr          ! poleward transports           (dia_ptr_init routine) 
    6465   USE diadct          ! sections transports           (dia_dct_init routine) 
    … …  
    158159                IF( ln_dyninc ) CALL dyn_asm_inc( nit000 - 1 )    ! Dynamics 
    159160                IF( ln_sshinc ) CALL ssh_asm_inc( nit000 - 1 )    ! SSH 
     161                IF( ln_logchltotinc .OR. ln_logchlpftinc ) CALL logchl_asm_inc( nit000 - 1 ) 
    160162             ENDIF 
    161163          ENDIF 
    162164 
     165#if defined key_agrif 
     166          CALL Agrif_Regrid() 
     167#endif 
     168 
    163169         DO WHILE ( istp <= nitend .AND. nstop == 0 ) 
    164170#if defined key_agrif 
    165             CALL Agrif_Step( stp )           ! AGRIF: time stepping 
     171            CALL stp                         ! AGRIF: time stepping 
    166172#else 
    167173            CALL stp( istp )                 ! standard time stepping 
    … …  
    173179 
    174180      IF( lk_diaobs   )   CALL dia_obs_wri 
     181      ! 
     182      IF( ( lk_asminc ).AND.( ln_balwri ) ) CALL asm_bal_wri( nitend )  ! Output balancing increments 
    175183      ! 
    176184      IF( ln_icebergs )   CALL icb_end( nitend ) 
    … …  
    187195      ! 
    188196#if defined key_agrif 
    189       CALL Agrif_ParentGrid_To_ChildGrid() 
    190       IF( lk_diaobs ) CALL dia_obs_wri 
    191       IF( nn_timing == 1 )   CALL timing_finalize 
    192       CALL Agrif_ChildGrid_To_ParentGrid() 
     197      IF( .NOT. Agrif_Root() ) THEN 
     198         CALL Agrif_ParentGrid_To_ChildGrid() 
     199         IF( lk_diaobs ) CALL dia_obs_wri 
     200         IF( nn_timing == 1 )   CALL timing_finalize 
     201         CALL Agrif_ChildGrid_To_ParentGrid() 
     202      ENDIF 
    193203#endif 
    194204      IF( nn_timing == 1 )   CALL timing_finalize 
    … …  
    334344         jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj   ! second dim. 
    335345#endif 
    336       ENDIF 
     346      ENDIF          
    337347         jpk = jpkdta                                             ! third dim 
     348#if defined key_agrif 
     349         ! simple trick to use same vertical grid as parent 
     350         ! but different number of levels:  
     351         ! Save maximum number of levels in jpkdta, then define all vertical grids 
     352         ! with this number. 
     353         ! Suppress once vertical online interpolation is ok 
     354         IF(.NOT.Agrif_Root()) jpkdta = Agrif_Parent(jpkdta) 
     355#endif 
    338356         jpim1 = jpi-1                                            ! inner domain indices 
    339357         jpjm1 = jpj-1                                            !   "           " 
    … …  
    459477 
    460478      !                                     ! Assimilation increments 
    461       IF( lk_asminc     )   CALL asm_inc_init   ! Initialize assimilation increments 
     479      IF( lk_asminc ) THEN  
     480#if defined key_shelf  
     481         CALL  zdf_mxl_tref()     ! Initialization of hmld_tref 
     482#endif  
     483         CALL asm_inc_init     ! Initialize assimilation increments  
     484      ENDIF 
     485            
    462486      IF(lwp) WRITE(numout,*) 'Euler time step switch is ', neuler 
    463487      ! 
    … …  
    710734      INTEGER :: ifac, jl, inu 
    711735      INTEGER, PARAMETER :: ntest = 14 
    712       INTEGER :: ilfax(ntest) 
    713       ! 
    714       ! lfax contains the set of allowed factors. 
    715       data (ilfax(jl),jl=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256,  & 
    716          &                            128,   64,   32,   16,    8,   4,   2  / 
    717       !!---------------------------------------------------------------------- 
     736      INTEGER, DIMENSION(ntest) :: ilfax 
     737      ! 
     738      ! ilfax contains the set of allowed factors. 
     739      ilfax(:) = (/(2**jl,jl=ntest,1,-1)/) 
     740      !!---------------------------------------------------------------------- 
     741      ! ilfax contains the set of allowed factors. 
     742      ilfax(:) = (/(2**jl,jl=ntest,1,-1)/) 
    718743 
    719744      ! Clear the error flag and initialise output vars 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/step.F90

    r7730 r7731  
    5050 
    5151#if defined key_agrif 
    52    SUBROUTINE stp( ) 
     52   RECURSIVE SUBROUTINE stp( ) 
    5353      INTEGER             ::   kstp   ! ocean time-step index 
    5454#else 
    … …  
    7979#if defined key_agrif 
    8080      kstp = nit000 + Agrif_Nb_Step() 
    81 !      IF ( Agrif_Root() .and. lwp) Write(*,*) '---' 
    82 !      IF (lwp) Write(*,*) 'Grid Number',Agrif_Fixed(),' time step ',kstp 
     81      IF ( lk_agrif_debug ) THEN 
     82         IF ( Agrif_Root() .and. lwp) Write(*,*) '---' 
     83         IF (lwp) Write(*,*) 'Grid Number',Agrif_Fixed(),' time step ',kstp, 'int tstep',Agrif_NbStepint() 
     84      ENDIF 
     85 
    8386      IF ( kstp == (nit000 + 1) ) lk_agrif_fstep = .FALSE. 
     87 
    8488# if defined key_iomput 
    8589      IF( Agrif_Nbstepint() == 0 )   CALL iom_swap( cxios_context ) 
    … …  
    110114      ! Update stochastic parameters and random T/S fluctuations 
    111115      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    112                         CALL sto_par( kstp )          ! Stochastic parameters 
     116       IF( ln_sto_eos ) CALL sto_par( kstp )          ! Stochastic parameters 
     117       IF( ln_sto_eos ) CALL sto_pts( tsn  )          ! Random T/S fluctuations 
    113118 
    114119      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    … …  
    152157      ! 
    153158      IF( lk_ldfslp ) THEN                            ! slope of lateral mixing 
    154          IF(ln_sto_eos ) CALL sto_pts( tsn )          ! Random T/S fluctuations 
    155159                         CALL eos( tsb, rhd, gdept_0(:,:,:) )               ! before in situ density 
    156160         IF( ln_zps .AND. .NOT. ln_isfcav)                               & 
    … …  
    188192          ! Note that the computation of vertical velocity above, hence "after" sea level 
    189193          ! is necessary to compute momentum advection for the rhs of barotropic loop: 
    190             IF(ln_sto_eos ) CALL sto_pts( tsn )                             ! Random T/S fluctuations 
    191194                            CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) ) ! now in situ density for hpg computation 
    192195            IF( ln_zps .AND. .NOT. ln_isfcav)                               & 
    … …  
    200203                                  ua(:,:,:) = 0.e0             ! set dynamics trends to zero 
    201204                                  va(:,:,:) = 0.e0 
    202           IF(  ln_asmiau .AND. & 
     205          IF(  lk_asminc .AND. ln_asmiau .AND. & 
    203206             & ln_dyninc       )  CALL dyn_asm_inc  ( kstp )   ! apply dynamics assimilation increment 
    204207          IF( ln_neptsimp )       CALL dyn_nept_cor ( kstp )   ! subtract Neptune velocities (simplified) 
    … …  
    239242      ! Passive Tracer Model 
    240243      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     244      IF( lk_asminc .AND. ln_asmiau .AND. ( ln_logchltotinc .OR. ln_logchlpftinc ) ) & 
     245         &               CALL logchl_asm_inc( kstp )  ! logchl assimilation 
    241246                         CALL trc_stp( kstp )         ! time-stepping 
    242247#endif 
    … …  
    248253                             tsa(:,:,:,:) = 0.e0            ! set tracer trends to zero 
    249254 
    250       IF(  ln_asmiau .AND. & 
     255      IF(  lk_asminc .AND. ln_asmiau .AND. & 
    251256         & ln_trainc     )   CALL tra_asm_inc( kstp )       ! apply tracer assimilation increment 
    252257                             CALL tra_sbc    ( kstp )       ! surface boundary condition 
    … …  
    270275         IF( ln_zdfnpc   )   CALL tra_npc( kstp )                ! update after fields by non-penetrative convection 
    271276                             CALL tra_nxt( kstp )                ! tracer fields at next time step 
    272             IF( ln_sto_eos ) CALL sto_pts( tsn )                 ! Random T/S fluctuations 
    273277                             CALL eos    ( tsa, rhd, rhop, fsdept_n(:,:,:) )  ! Time-filtered in situ density for hpg computation 
    274278            IF( ln_zps .AND. .NOT. ln_isfcav)                                & 
    … …  
    281285      ELSE                                                  ! centered hpg  (eos then time stepping) 
    282286         IF ( .NOT. lk_dynspg_ts ) THEN                     ! eos already called in time-split case 
    283             IF( ln_sto_eos ) CALL sto_pts( tsn )    ! Random T/S fluctuations 
    284287                             CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
    285288         IF( ln_zps .AND. .NOT. ln_isfcav)                                   & 
    … …  
    298301      ! Dynamics                                    (tsa used as workspace) 
    299302      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
     303 
     304      IF( ln_bkgwri )        CALL asm_bkg_wri( kstp )     ! output background fields 
     305 
    300306      IF( lk_dynspg_ts   )  THEN 
    301307                                                             ! revert to previously computed momentum tendencies 
    … …  
    314320                               va(:,:,:) = 0.e0 
    315321 
    316         IF(  ln_asmiau .AND. & 
     322        IF(  lk_asminc .AND. ln_asmiau .AND. & 
    317323           & ln_dyninc      )  CALL dyn_asm_inc( kstp )     ! apply dynamics assimilation increment 
    318         IF( ln_bkgwri )        CALL asm_bkg_wri( kstp )     ! output background fields 
    319324        IF( ln_neptsimp )      CALL dyn_nept_cor( kstp )    ! subtract Neptune velocities (simplified) 
    320325        IF( lk_bdy          )  CALL bdy_dyn3d_dmp(kstp )    ! bdy damping trends 
    … …  
    335340                               CALL ssh_swp( kstp )         ! swap of sea surface height 
    336341      IF( lk_vvl           )   CALL dom_vvl_sf_swp( kstp )  ! swap of vertical scale factors 
    337  
     342      ! 
     343      IF( lrst_oce         )   CALL rst_write( kstp )       ! write output ocean restart file 
     344 
     345#if defined key_agrif 
     346      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
     347      ! AGRIF 
     348      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<       
     349                               CALL Agrif_Integrate_ChildGrids( stp )   
     350 
     351      IF ( Agrif_NbStepint().EQ.0 ) THEN 
     352                               CALL Agrif_Update_Tra()      ! Update active tracers 
     353                               CALL Agrif_Update_Dyn()      ! Update momentum 
     354      ENDIF 
     355#endif 
    338356      IF( ln_diahsb        )   CALL dia_hsb( kstp )         ! - ML - global conservation diagnostics 
    339357      IF( lk_diaobs  )         CALL dia_obs( kstp )         ! obs-minus-model (assimilation) diagnostics (call after dynamics update) 
    340358 
    341359      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    342       ! Control and restarts 
     360      ! Control 
    343361      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
    344362                               CALL stp_ctl( kstp, indic ) 
    … …  
    352370         IF( lwm.AND.numoni /= -1 ) CALL FLUSH    ( numoni )     ! flush output namelist ice 
    353371      ENDIF 
    354       IF( lrst_oce         )   CALL rst_write    ( kstp )   ! write output ocean restart file 
    355372 
    356373      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
    … …  
    367384      ! 
    368385      IF( nn_timing == 1 .AND.  kstp == nit000  )   CALL timing_reset 
     386      !      
    369387      ! 
    370388   END SUBROUTINE stp 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/step_oce.F90

    r7730 r7731  
    112112#if defined key_agrif 
    113113   USE agrif_opa_sponge ! Momemtum and tracers sponges 
     114   USE agrif_opa_update ! Update (2-way nesting) 
    114115#endif 
    115116#if defined key_top 

  • branches/UKMO/dev_r5518_v3.6_asm_nemovar_community/NEMOGCM/NEMO/OPA_SRC/stpctl.F90

    r7730 r7731  
    1717   USE dom_oce         ! ocean space and time domain variables  
    1818   USE sol_oce         ! ocean space and time domain variables  
     19   USE sbc_oce         ! surface boundary conditions variables 
    1920   USE in_out_manager  ! I/O manager 
    2021   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    … …  
    2223   USE dynspg_oce      ! pressure gradient schemes  
    2324   USE c1d             ! 1D vertical configuration 
     25 
    2426 
    2527   IMPLICIT NONE 
    … …  
    5254      INTEGER, INTENT( inout ) ::   kindic  ! indicator of solver convergence 
    5355      !! 
     56      CHARACTER(len = 32) ::        clfname ! time stepping output file name 
    5457      INTEGER  ::   ji, jj, jk              ! dummy loop indices 
    5558      INTEGER  ::   ii, ij, ik              ! temporary integers 
    … …  
    6366         WRITE(numout,*) 'stp_ctl : time-stepping control' 
    6467         WRITE(numout,*) '~~~~~~~' 
    65          ! open time.step file 
    66          CALL ctl_opn( numstp, 'time.step', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
     68         ! open time.step file with special treatment for SAS 
     69         IF ( nn_components == jp_iam_sas ) THEN 
     70            clfname = 'time.step.sas' 
     71         ELSE 
     72            clfname = 'time.step' 
     73         ENDIF 
     74         CALL ctl_opn( numstp, TRIM(clfname), 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
    6775      ENDIF 
    6876 
Note: See TracChangeset for help on using the changeset viewer.