Changeset 5168


Ignore:
Timestamp:
2015-03-25T10:01:32+01:00 (6 years ago)
Author:
pabouttier
Message:

Fixed several bugs described in Ticket #1360

Location:
branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/tstool_tam.F90

    r3611 r5168  
    1414      & e1u,                 & 
    1515      & e2u,                 & 
    16 #if defined key_zco 
    1716      & e3t_0,               & 
    18 #else 
    1917      & e3u,                 & 
    20 #endif 
    2118      & umask,               & 
    2219      & mig,                 & 
  • branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPA_SRC/ASM/asminc.F90

    r3784 r5168  
    99   !!                 ! 2007-04  (A. Weaver)  Merge with OPAVAR/NEMOVAR 
    1010   !!   NEMO     3.3  ! 2010-05  (D. Lea)  Update to work with NEMO v3.2 
    11    !!             -   ! 2010-05  (D. Lea)  add calc_month_len routine based on day_init  
     11   !!             -   ! 2010-05  (D. Lea)  add calc_month_len routine based on day_init 
    1212   !!            3.4  ! 2012-10  (A. Weaver and K. Mogensen) Fix for direct initialization 
    1313   !!---------------------------------------------------------------------- 
     
    4343   IMPLICIT NONE 
    4444   PRIVATE 
    45     
     45 
     46   PUBLIC asm_inc_rea_nam  !: Read namelist for ASM 
    4647   PUBLIC   asm_inc_init   !: Initialize the increment arrays and IAU weights 
    4748   PUBLIC   calc_date      !: Compute the calendar date YYYYMMDD on a given step 
     
    7071   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE ::   u_bkg   , v_bkg      !: Background u- & v- velocity components 
    7172   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE ::   t_bkginc, s_bkginc   !: Increment to the background T & S 
    72    REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE ::   u_bkginc, v_bkginc   !: Increment to the u- & v-components  
     73   REAL(wp), PUBLIC, DIMENSION(:,:,:), ALLOCATABLE ::   u_bkginc, v_bkginc   !: Increment to the u- & v-components 
    7374   REAL(wp), PUBLIC, DIMENSION(:)    , ALLOCATABLE ::   wgtiau               !: IAU weights for each time step 
    7475#if defined key_asminc 
     
    7879   INTEGER , PUBLIC ::   nitbkg      !: Time step of the background state used in the Jb term 
    7980   INTEGER , PUBLIC ::   nitdin      !: Time step of the background state for direct initialization 
    80    INTEGER , PUBLIC ::   nitiaustr   !: Time step of the start of the IAU interval  
     81   INTEGER , PUBLIC ::   nitiaustr   !: Time step of the start of the IAU interval 
    8182   INTEGER , PUBLIC ::   nitiaufin   !: Time step of the end of the IAU interval 
    82    !  
     83   ! 
    8384   INTEGER , PUBLIC ::   niaufn      !: Type of IAU weighing function: = 0   Constant weighting 
    84    !                                 !: = 1   Linear hat-like, centred in middle of IAU interval  
     85   !                                 !: = 1   Linear hat-like, centred in middle of IAU interval 
    8586   REAL(wp), PUBLIC ::   salfixmin   !: Ensure that the salinity is larger than this  value if (ln_salfix) 
    8687 
     
    100101CONTAINS 
    101102 
    102    SUBROUTINE asm_inc_init 
    103       !!---------------------------------------------------------------------- 
    104       !!                    ***  ROUTINE asm_inc_init  *** 
    105       !!           
    106       !! ** Purpose : Initialize the assimilation increment and IAU weights. 
    107       !! 
    108       !! ** Method  : Initialize the assimilation increment and IAU weights. 
    109       !! 
    110       !! ** Action  :  
    111       !!---------------------------------------------------------------------- 
    112       !! 
    113       !! 
    114       INTEGER :: ji,jj,jk 
    115       INTEGER :: jt 
    116       INTEGER :: imid 
    117       INTEGER :: inum 
    118       INTEGER :: iiauper         ! Number of time steps in the IAU period 
    119       INTEGER :: icycper         ! Number of time steps in the cycle 
    120       INTEGER :: iitend_date     ! Date YYYYMMDD of final time step 
    121       INTEGER :: iitbkg_date     ! Date YYYYMMDD of background time step for Jb term 
    122       INTEGER :: iitdin_date     ! Date YYYYMMDD of background time step for DI 
    123       INTEGER :: iitiaustr_date  ! Date YYYYMMDD of IAU interval start time step 
    124       INTEGER :: iitiaufin_date  ! Date YYYYMMDD of IAU interval final time step 
    125  
    126       REAL(wp) :: znorm        ! Normalization factor for IAU weights 
    127       REAL(wp) :: ztotwgt      ! Value of time-integrated IAU weights  
    128                                ! (should be equal to one) 
    129       REAL(wp) :: z_inc_dateb  ! Start date of interval on which increment is valid 
    130       REAL(wp) :: z_inc_datef  ! End date of interval on which increment is valid 
    131       REAL(wp) :: zdate_bkg    ! Date in background state file for DI 
    132       REAL(wp) :: zdate_inc    ! Time axis in increments file 
    133  
    134       REAL(wp), POINTER, DIMENSION(:,:) :: hdiv 
    135       !! 
     103   SUBROUTINE asm_inc_rea_nam 
    136104      NAMELIST/nam_asminc/ ln_bkgwri,                                      & 
    137105         &                 ln_trainc, ln_dyninc, ln_sshinc,                & 
     
    158126      salfixmin = -9999 
    159127      nitbkg    = 0 
    160       nitdin    = 0       
     128      nitdin    = 0 
    161129      nitiaustr = 1 
    162130      nitiaufin = 150      ! = 10 days with ORCA2 
     
    187155         WRITE(numout,*) '      Minimum salinity after applying the increments           salfixmin = ', salfixmin 
    188156      ENDIF 
     157   END SUBROUTINE 
     158 
     159   SUBROUTINE asm_inc_init 
     160      !!---------------------------------------------------------------------- 
     161      !!                    ***  ROUTINE asm_inc_init  *** 
     162      !! 
     163      !! ** Purpose : Initialize the assimilation increment and IAU weights. 
     164      !! 
     165      !! ** Method  : Initialize the assimilation increment and IAU weights. 
     166      !! 
     167      !! ** Action  : 
     168      !!---------------------------------------------------------------------- 
     169      !! 
     170      !! 
     171      INTEGER :: ji,jj,jk 
     172      INTEGER :: jt 
     173      INTEGER :: imid 
     174      INTEGER :: inum 
     175      INTEGER :: iiauper         ! Number of time steps in the IAU period 
     176      INTEGER :: icycper         ! Number of time steps in the cycle 
     177      INTEGER :: iitend_date     ! Date YYYYMMDD of final time step 
     178      INTEGER :: iitbkg_date     ! Date YYYYMMDD of background time step for Jb term 
     179      INTEGER :: iitdin_date     ! Date YYYYMMDD of background time step for DI 
     180      INTEGER :: iitiaustr_date  ! Date YYYYMMDD of IAU interval start time step 
     181      INTEGER :: iitiaufin_date  ! Date YYYYMMDD of IAU interval final time step 
     182 
     183      REAL(wp) :: znorm        ! Normalization factor for IAU weights 
     184      REAL(wp) :: ztotwgt      ! Value of time-integrated IAU weights 
     185                               ! (should be equal to one) 
     186      REAL(wp) :: z_inc_dateb  ! Start date of interval on which increment is valid 
     187      REAL(wp) :: z_inc_datef  ! End date of interval on which increment is valid 
     188      REAL(wp) :: zdate_bkg    ! Date in background state file for DI 
     189      REAL(wp) :: zdate_inc    ! Time axis in increments file 
     190 
     191      REAL(wp), POINTER, DIMENSION(:,:) :: hdiv 
     192      !! 
    189193 
    190194      nitbkg_r    = nitbkg    + nit000 - 1  ! Background time referenced to nit000 
     
    293297 
    294298            !--------------------------------------------------------- 
    295             ! Constant IAU forcing  
     299            ! Constant IAU forcing 
    296300            !--------------------------------------------------------- 
    297301 
     
    303307 
    304308            !--------------------------------------------------------- 
    305             ! Linear hat-like, centred in middle of IAU interval  
     309            ! Linear hat-like, centred in middle of IAU interval 
    306310            !--------------------------------------------------------- 
    307311 
     
    309313            znorm = 0.0 
    310314            IF ( MOD( iiauper, 2 ) == 0 ) THEN  ! Even number of time steps in IAU interval 
    311                imid = iiauper / 2  
     315               imid = iiauper / 2 
    312316               DO jt = 1, imid 
    313317                  znorm = znorm + REAL( jt ) 
     
    315319               znorm = 2.0 * znorm 
    316320            ELSE                               ! Odd number of time steps in IAU interval 
    317                imid = ( iiauper + 1 ) / 2         
     321               imid = ( iiauper + 1 ) / 2 
    318322               DO jt = 1, imid - 1 
    319323                  znorm = znorm + REAL( jt ) 
     
    342346             DO jt = 1, icycper 
    343347                ztotwgt = ztotwgt + wgtiau(jt) 
    344                 WRITE(numout,*) '         ', jt, '       ', wgtiau(jt)  
    345              END DO    
     348                WRITE(numout,*) '         ', jt, '       ', wgtiau(jt) 
     349             END DO 
    346350             WRITE(numout,*) '         ===================================' 
    347351             WRITE(numout,*) '         Time-integrated weight = ', ztotwgt 
    348352             WRITE(numout,*) '         ===================================' 
    349353          ENDIF 
    350           
     354 
    351355      ENDIF 
    352356 
     
    381385         CALL iom_open( c_asminc, inum ) 
    382386 
    383          CALL iom_get( inum, 'time', zdate_inc )  
     387         CALL iom_get( inum, 'time', zdate_inc ) 
    384388 
    385389         CALL iom_get( inum, 'z_inc_dateb', z_inc_dateb ) 
     
    389393 
    390394         IF(lwp) THEN 
    391             WRITE(numout,*)  
     395            WRITE(numout,*) 
    392396            WRITE(numout,*) 'asm_inc_init : Assimilation increments valid ', & 
    393397               &            ' between dates ', NINT( z_inc_dateb ),' and ',  & 
     
    405409            &                ' not agree with Direct Initialization time' ) 
    406410 
    407          IF ( ln_trainc ) THEN    
     411         IF ( ln_trainc ) THEN 
    408412            CALL iom_get( inum, jpdom_autoglo, 'bckint', t_bkginc, 1 ) 
    409413            CALL iom_get( inum, jpdom_autoglo, 'bckins', s_bkginc, 1 ) 
     
    417421         ENDIF 
    418422 
    419          IF ( ln_dyninc ) THEN    
    420             CALL iom_get( inum, jpdom_autoglo, 'bckinu', u_bkginc, 1 )               
    421             CALL iom_get( inum, jpdom_autoglo, 'bckinv', v_bkginc, 1 )               
     423         IF ( ln_dyninc ) THEN 
     424            CALL iom_get( inum, jpdom_autoglo, 'bckinu', u_bkginc, 1 ) 
     425            CALL iom_get( inum, jpdom_autoglo, 'bckinv', v_bkginc, 1 ) 
    422426            ! Apply the masks 
    423427            u_bkginc(:,:,:) = u_bkginc(:,:,:) * umask(:,:,:) 
     
    428432            WHERE( ABS( v_bkginc(:,:,:) ) > 1.0e+10 ) v_bkginc(:,:,:) = 0.0 
    429433         ENDIF 
    430          
     434 
    431435         IF ( ln_sshinc ) THEN 
    432436            CALL iom_get( inum, jpdom_autoglo, 'bckineta', ssh_bkginc, 1 ) 
     
    448452 
    449453         CALL iom_close( inum ) 
    450   
     454 
    451455      ENDIF 
    452456 
     
    457461      IF ( ln_dyninc .AND. nn_divdmp > 0 ) THEN 
    458462 
    459          CALL wrk_alloc(jpi,jpj,hdiv)  
     463         CALL wrk_alloc(jpi,jpj,hdiv) 
    460464 
    461465         DO  jt = 1, nn_divdmp 
     
    482486                     u_bkginc(ji,jj,jk) = u_bkginc(ji,jj,jk) + 0.2_wp * ( e1t(ji+1,jj)*e2t(ji+1,jj) * hdiv(ji+1,jj)   & 
    483487                                                                        - e1t(ji  ,jj)*e2t(ji  ,jj) * hdiv(ji  ,jj) ) & 
    484                                                                       / e1u(ji,jj) * umask(ji,jj,jk)  
     488                                                                      / e1u(ji,jj) * umask(ji,jj,jk) 
    485489                     v_bkginc(ji,jj,jk) = v_bkginc(ji,jj,jk) + 0.2_wp * ( e1t(ji,jj+1)*e2t(ji,jj+1) * hdiv(ji,jj+1)   & 
    486490                                                                        - e1t(ji,jj  )*e2t(ji,jj  ) * hdiv(ji,jj  ) ) & 
    487                                                                       / e2v(ji,jj) * vmask(ji,jj,jk)  
     491                                                                      / e2v(ji,jj) * vmask(ji,jj,jk) 
    488492                  END DO 
    489493               END DO 
     
    493497         END DO 
    494498 
    495          CALL wrk_dealloc(jpi,jpj,hdiv)  
     499         CALL wrk_dealloc(jpi,jpj,hdiv) 
    496500 
    497501      ENDIF 
     
    523527         CALL iom_open( c_asmdin, inum ) 
    524528 
    525          CALL iom_get( inum, 'rdastp', zdate_bkg )  
    526          
     529         CALL iom_get( inum, 'rdastp', zdate_bkg ) 
     530 
    527531         IF(lwp) THEN 
    528             WRITE(numout,*)  
     532            WRITE(numout,*) 
    529533            WRITE(numout,*) 'asm_inc_init : Assimilation background state valid at : ', & 
    530534               &  NINT( zdate_bkg ) 
     
    536540            &                ' not agree with Direct Initialization time' ) 
    537541 
    538          IF ( ln_trainc ) THEN    
     542         IF ( ln_trainc ) THEN 
    539543            CALL iom_get( inum, jpdom_autoglo, 'tn', t_bkg ) 
    540544            CALL iom_get( inum, jpdom_autoglo, 'sn', s_bkg ) 
     
    543547         ENDIF 
    544548 
    545          IF ( ln_dyninc ) THEN    
     549         IF ( ln_dyninc ) THEN 
    546550            CALL iom_get( inum, jpdom_autoglo, 'un', u_bkg ) 
    547551            CALL iom_get( inum, jpdom_autoglo, 'vn', v_bkg ) 
     
    549553            v_bkg(:,:,:) = v_bkg(:,:,:) * vmask(:,:,:) 
    550554         ENDIF 
    551          
     555 
    552556         IF ( ln_sshinc ) THEN 
    553557            CALL iom_get( inum, jpdom_autoglo, 'sshn', ssh_bkg ) 
     
    565569      !!---------------------------------------------------------------------- 
    566570      !!                    ***  ROUTINE calc_date  *** 
    567       !!           
     571      !! 
    568572      !! ** Purpose : Compute the calendar date YYYYMMDD at a given time step. 
    569573      !! 
    570574      !! ** Method  : Compute the calendar date YYYYMMDD at a given time step. 
    571575      !! 
    572       !! ** Action  :  
     576      !! ** Action  : 
    573577      !!---------------------------------------------------------------------- 
    574578      INTEGER, INTENT(IN) :: kit000  ! Initial time step 
     
    595599      iyea0 =   kdate0 / 10000 
    596600      imon0 = ( kdate0 - ( iyea0 * 10000 ) ) / 100 
    597       iday0 =   kdate0 - ( iyea0 * 10000 ) - ( imon0 * 100 )  
     601      iday0 =   kdate0 - ( iyea0 * 10000 ) - ( imon0 * 100 ) 
    598602 
    599603      ! Check that kt >= kit000 - 1 
     
    610614      ! Compute the number of days from the initial date 
    611615      idaystp = INT( REAL( kt - kit000 ) * rdt / 86400. ) 
    612     
     616 
    613617      iday    = iday0 
    614618      imon    = imon0 
     
    628632            iyea = iyea + 1 
    629633            CALL calc_month_len( iyea, imonth_len )  ! update month lengths 
    630          ENDIF                  
     634         ENDIF 
    631635         idaycnt = idaycnt + 1 
    632636      END DO 
     
    640644      !!---------------------------------------------------------------------- 
    641645      !!                    ***  ROUTINE calc_month_len  *** 
    642       !!           
     646      !! 
    643647      !! ** Purpose : Compute the number of days in a months given a year. 
    644648      !! 
    645       !! ** Method  :  
     649      !! ** Method  : 
    646650      !!---------------------------------------------------------------------- 
    647651      INTEGER, DIMENSION(12) ::   imonth_len    !: length in days of the months of the current year 
     
    650654      ! 
    651655      ! length of the month of the current year (from nleapy, read in namelist) 
    652       IF ( nleapy < 2 ) THEN  
     656      IF ( nleapy < 2 ) THEN 
    653657         imonth_len(:) = (/ 31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31 /) 
    654658         IF ( nleapy == 1 ) THEN   ! we are using calendar with leap years 
     
    667671      !!---------------------------------------------------------------------- 
    668672      !!                    ***  ROUTINE tra_asm_inc  *** 
    669       !!           
     673      !! 
    670674      !! ** Purpose : Apply the tracer (T and S) assimilation increments 
    671675      !! 
    672676      !! ** Method  : Direct initialization or Incremental Analysis Updating 
    673677      !! 
    674       !! ** Action  :  
     678      !! ** Action  : 
    675679      !!---------------------------------------------------------------------- 
    676680      INTEGER, INTENT(IN) :: kt               ! Current time step 
     
    682686      !!---------------------------------------------------------------------- 
    683687 
    684       ! freezing point calculation taken from oc_fz_pt (but calculated for all depths)  
    685       ! used to prevent the applied increments taking the temperature below the local freezing point  
    686  
    687 #if defined key_cice  
     688      ! freezing point calculation taken from oc_fz_pt (but calculated for all depths) 
     689      ! used to prevent the applied increments taking the temperature below the local freezing point 
     690 
     691#if defined key_cice 
    688692        fzptnz(:,:,:) = -1.8_wp 
    689 #else  
     693#else 
    690694        DO jk = 1, jpk 
    691695           DO jj = 1, jpj 
    692696              DO ji = 1, jpk 
    693                  fzptnz (ji,jj,jk) = ( -0.0575_wp + 1.710523e-3_wp * SQRT( tsn(ji,jj,jk,jp_sal) )                   &  
    694                                                   - 2.154996e-4_wp *       tsn(ji,jj,jk,jp_sal)   ) * tsn(ji,jj,jk,jp_sal)  &  
    695                                                   - 7.53e-4_wp * fsdepw(ji,jj,jk)       ! (pressure in dbar)  
     697                 fzptnz (ji,jj,jk) = ( -0.0575_wp + 1.710523e-3_wp * SQRT( tsn(ji,jj,jk,jp_sal) )                   & 
     698                                                  - 2.154996e-4_wp *       tsn(ji,jj,jk,jp_sal)   ) * tsn(ji,jj,jk,jp_sal)  & 
     699                                                  - 7.53e-4_wp * fsdepw(ji,jj,jk)       ! (pressure in dbar) 
    696700              END DO 
    697701           END DO 
    698702        END DO 
    699 #endif  
     703#endif 
    700704 
    701705      IF ( ln_asmiau ) THEN 
     
    711715 
    712716            IF(lwp) THEN 
    713                WRITE(numout,*)  
     717               WRITE(numout,*) 
    714718               WRITE(numout,*) 'tra_asm_inc : Tracer IAU at time step = ', & 
    715719                  &  kt,' with IAU weight = ', wgtiau(it) 
     
    722726                  ! Do not apply negative increments if the temperature will fall below freezing 
    723727                  WHERE(t_bkginc(:,:,jk) > 0.0_wp .OR. & 
    724                      &   tsn(:,:,jk,jp_tem) + tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * wgtiau(it) > fzptnz(:,:,jk) )  
    725                      tsa(:,:,jk,jp_tem) = tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * zincwgt   
     728                     &   t_bkg(:,:,jk) + tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * wgtiau(it) > fzptnz(:,:,jk) ) 
     729                     tsa(:,:,jk,jp_tem) = tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * zincwgt 
    726730                  END WHERE 
    727731               ELSE 
    728                   tsa(:,:,jk,jp_tem) = tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * zincwgt   
     732                  tsa(:,:,jk,jp_tem) = tsa(:,:,jk,jp_tem) + t_bkginc(:,:,jk) * zincwgt 
    729733               ENDIF 
    730734               IF (ln_salfix) THEN 
     
    732736                  ! minimum value salfixmin 
    733737                  WHERE(s_bkginc(:,:,jk) > 0.0_wp .OR. & 
    734                      &   tsn(:,:,jk,jp_sal) + tsa(:,:,jk,jp_sal) + s_bkginc(:,:,jk) * wgtiau(it) > salfixmin )  
     738                     &  s_bkg(:,:,jk) + tsa(:,:,jk,jp_sal) + s_bkginc(:,:,jk) * wgtiau(it) > salfixmin ) 
    735739                     tsa(:,:,jk,jp_sal) = tsa(:,:,jk,jp_sal) + s_bkginc(:,:,jk) * zincwgt 
    736740                  END WHERE 
     
    753757         ! Direct Initialization 
    754758         !-------------------------------------------------------------------- 
    755              
     759 
    756760         IF ( kt == nitdin_r ) THEN 
    757761 
     
    762766               ! Do not apply negative increments if the temperature will fall below freezing 
    763767               WHERE(t_bkginc(:,:,:) > 0.0_wp .OR. & 
    764                   &   tsn(:,:,:,jp_tem) + t_bkginc(:,:,:) > fzptnz(:,:,:) )  
    765                   tsn(:,:,:,jp_tem) = t_bkg(:,:,:) + t_bkginc(:,:,:)    
     768                  &   tsn(:,:,:,jp_tem) + t_bkginc(:,:,:) > fzptnz(:,:,:) ) 
     769                  tsn(:,:,:,jp_tem) = t_bkg(:,:,:) + t_bkginc(:,:,:) 
    766770               END WHERE 
    767771            ELSE 
    768                tsn(:,:,:,jp_tem) = t_bkg(:,:,:) + t_bkginc(:,:,:)    
     772               tsn(:,:,:,jp_tem) = t_bkg(:,:,:) + t_bkginc(:,:,:) 
    769773            ENDIF 
    770774            IF (ln_salfix) THEN 
     
    772776               ! minimum value salfixmin 
    773777               WHERE(s_bkginc(:,:,:) > 0.0_wp .OR. & 
    774                   &   tsn(:,:,:,jp_sal) + s_bkginc(:,:,:) > salfixmin )  
    775                   tsn(:,:,:,jp_sal) = s_bkg(:,:,:) + s_bkginc(:,:,:)    
     778                  &   tsn(:,:,:,jp_sal) + s_bkginc(:,:,:) > salfixmin ) 
     779                  tsn(:,:,:,jp_sal) = s_bkg(:,:,:) + s_bkginc(:,:,:) 
    776780               END WHERE 
    777781            ELSE 
    778                tsn(:,:,:,jp_sal) = s_bkg(:,:,:) + s_bkginc(:,:,:)    
     782               tsn(:,:,:,jp_sal) = s_bkg(:,:,:) + s_bkginc(:,:,:) 
    779783            ENDIF 
    780784 
     
    782786 
    783787            CALL eos( tsb, rhd, rhop )                ! Before potential and in situ densities 
    784           
     788 
    785789            IF( ln_zps .AND. .NOT. lk_c1d ) & 
    786790               &  CALL zps_hde( nit000, jpts, tsb, &  ! Partial steps: before horizontal derivative 
     
    797801            DEALLOCATE( s_bkg    ) 
    798802         ENDIF 
    799          !   
     803         ! 
    800804      ENDIF 
    801805      ! Perhaps the following call should be in step 
     
    808812      !!---------------------------------------------------------------------- 
    809813      !!                    ***  ROUTINE dyn_asm_inc  *** 
    810       !!           
     814      !! 
    811815      !! ** Purpose : Apply the dynamics (u and v) assimilation increments. 
    812816      !! 
    813817      !! ** Method  : Direct initialization or Incremental Analysis Updating. 
    814818      !! 
    815       !! ** Action  :  
     819      !! ** Action  : 
    816820      !!---------------------------------------------------------------------- 
    817821      INTEGER, INTENT(IN) :: kt   ! Current time step 
     
    834838 
    835839            IF(lwp) THEN 
    836                WRITE(numout,*)  
     840               WRITE(numout,*) 
    837841               WRITE(numout,*) 'dyn_asm_inc : Dynamics IAU at time step = ', & 
    838842                  &  kt,' with IAU weight = ', wgtiau(it) 
     
    845849               va(:,:,jk) = va(:,:,jk) + v_bkginc(:,:,jk) * zincwgt 
    846850            END DO 
    847             
     851 
    848852            IF ( kt == nitiaufin_r ) THEN 
    849853               DEALLOCATE( u_bkginc ) 
     
    853857         ENDIF 
    854858 
    855       ELSEIF ( ln_asmdin ) THEN  
     859      ELSEIF ( ln_asmdin ) THEN 
    856860 
    857861         !-------------------------------------------------------------------- 
    858862         ! Direct Initialization 
    859863         !-------------------------------------------------------------------- 
    860           
     864 
    861865         IF ( kt == nitdin_r ) THEN 
    862866 
     
    865869            ! Initialize the now fields with the background + increment 
    866870            un(:,:,:) = u_bkg(:,:,:) + u_bkginc(:,:,:) 
    867             vn(:,:,:) = v_bkg(:,:,:) + v_bkginc(:,:,:)   
     871            vn(:,:,:) = v_bkg(:,:,:) + v_bkginc(:,:,:) 
    868872 
    869873            ub(:,:,:) = un(:,:,:)         ! Update before fields 
    870874            vb(:,:,:) = vn(:,:,:) 
    871   
     875 
    872876            DEALLOCATE( u_bkg    ) 
    873877            DEALLOCATE( v_bkg    ) 
     
    884888      !!---------------------------------------------------------------------- 
    885889      !!                    ***  ROUTINE ssh_asm_inc  *** 
    886       !!           
     890      !! 
    887891      !! ** Purpose : Apply the sea surface height assimilation increment. 
    888892      !! 
    889893      !! ** Method  : Direct initialization or Incremental Analysis Updating. 
    890894      !! 
    891       !! ** Action  :  
     895      !! ** Action  : 
    892896      !!---------------------------------------------------------------------- 
    893897      INTEGER, INTENT(IN) :: kt   ! Current time step 
     
    910914 
    911915            IF(lwp) THEN 
    912                WRITE(numout,*)  
     916               WRITE(numout,*) 
    913917               WRITE(numout,*) 'ssh_asm_inc : SSH IAU at time step = ', & 
    914918                  &  kt,' with IAU weight = ', wgtiau(it) 
     
    938942 
    939943            ! Initialize the now fields the background + increment 
    940             sshn(:,:) = ssh_bkg(:,:) + ssh_bkginc(:,:)   
     944            sshn(:,:) = ssh_bkg(:,:) + ssh_bkginc(:,:) 
    941945 
    942946            ! Update before fields 
    943             sshb(:,:) = sshn(:,:)          
     947            sshb(:,:) = sshn(:,:) 
    944948 
    945949            IF( lk_vvl ) THEN 
     
    961965      !!---------------------------------------------------------------------- 
    962966      !!                    ***  ROUTINE seaice_asm_inc  *** 
    963       !!           
     967      !! 
    964968      !! ** Purpose : Apply the sea ice assimilation increment. 
    965969      !! 
    966970      !! ** Method  : Direct initialization or Incremental Analysis Updating. 
    967971      !! 
    968       !! ** Action  :  
     972      !! ** Action  : 
    969973      !! 
    970974      !! History : 
     
    9981002 
    9991003            it = kt - nit000 + 1 
    1000             zincwgt = wgtiau(it)      ! IAU weight for the current time step  
     1004            zincwgt = wgtiau(it)      ! IAU weight for the current time step 
    10011005            ! note this is not a tendency so should not be divided by rdt (as with the tracer and other increments) 
    10021006 
    10031007            IF(lwp) THEN 
    1004                WRITE(numout,*)  
     1008               WRITE(numout,*) 
    10051009               WRITE(numout,*) 'seaice_asm_inc : sea ice conc IAU at time step = ', & 
    10061010                  &  kt,' with IAU weight = ', wgtiau(it) 
     
    10211025            ! Nudge sea ice depth to bring it up to a required minimum depth 
    10221026 
    1023             WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin )  
    1024                zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt     
     1027            WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin ) 
     1028               zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt 
    10251029            ELSEWHERE 
    10261030               zhicifinc(:,:) = 0.0_wp 
     
    10291033! nudge ice depth 
    10301034            hicif(:,:)=hicif(:,:) + zhicifinc(:,:) 
    1031             phicif(:,:)=phicif(:,:) + zhicifinc(:,:)        
     1035            phicif(:,:)=phicif(:,:) + zhicifinc(:,:) 
    10321036 
    10331037! seaice salinity balancing (to add) 
     
    10711075            zofrld(:,:)=frld(:,:) 
    10721076            zohicif(:,:)=hicif(:,:) 
    1073   
     1077 
    10741078            ! Initialize the now fields the background + increment 
    10751079 
    10761080            frld(:,:) = MIN( MAX( frld(:,:) - seaice_bkginc(:,:), 0.0_wp), 1.0_wp) 
    1077             pfrld(:,:) = frld(:,:)  
     1081            pfrld(:,:) = frld(:,:) 
    10781082            fr_i(:,:) = 1.0_wp - frld(:,:)        ! adjust ice fraction 
    10791083 
     
    10821086            ! Nudge sea ice depth to bring it up to a required minimum depth 
    10831087 
    1084             WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin )  
    1085                zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt     
     1088            WHERE( zseaicendg(:,:) > 0.0_wp .AND. hicif(:,:) < zhicifmin ) 
     1089               zhicifinc(:,:) = (zhicifmin - hicif(:,:)) * zincwgt 
    10861090            ELSEWHERE 
    10871091               zhicifinc(:,:) = 0.0_wp 
     
    10901094! nudge ice depth 
    10911095            hicif(:,:)=hicif(:,:) + zhicifinc(:,:) 
    1092             phicif(:,:)=phicif(:,:)        
     1096            phicif(:,:)=phicif(:,:) 
    10931097 
    10941098! seaice salinity balancing (to add) 
    1095    
    1096 #endif 
    1097   
     1099 
     1100#endif 
     1101 
    10981102#if defined key_cice 
    10991103 
     
    11101114#if defined key_cice 
    11111115 
    1112 ! Zero ice increment tendency into CICE  
     1116! Zero ice increment tendency into CICE 
    11131117            ndaice_da(:,:) = 0.0_wp 
    11141118 
    11151119#endif 
    1116           
     1120 
    11171121         ENDIF 
    11181122 
    11191123!#if defined key_lim2 || defined key_cice 
    11201124! 
    1121 !            IF (ln_seaicebal ) THEN        
     1125!            IF (ln_seaicebal ) THEN 
    11221126!             !! balancing salinity increments 
    11231127!             !! simple case from limflx.F90 (doesn't include a mass flux) 
     
    11311135! 
    11321136!             DO jj = 1, jpj 
    1133 !               DO ji = 1, jpi  
     1137!               DO ji = 1, jpi 
    11341138!           ! calculate change in ice and snow mass per unit area 
    11351139!           ! positive values imply adding salt to the ocean (results from ice formation) 
     
    11421146! 
    11431147!           ! prevent small mld 
    1144 !           ! less than 10m can cause salinity instability  
     1148!           ! less than 10m can cause salinity instability 
    11451149!                 IF (mld < 10) mld=10 
    11461150! 
    1147 !           ! set to bottom of a level  
     1151!           ! set to bottom of a level 
    11481152!                 DO jk = jpk-1, 2, -1 
    1149 !                   IF ((mld > gdepw(ji,jj,jk)) .and. (mld < gdepw(ji,jj,jk+1))) THEN  
     1153!                   IF ((mld > gdepw(ji,jj,jk)) .and. (mld < gdepw(ji,jj,jk+1))) THEN 
    11501154!                     mld=gdepw(ji,jj,jk+1) 
    11511155!                     jkmax=jk 
     
    11541158! 
    11551159!            ! avoid applying salinity balancing in shallow water or on land 
    1156 !            !  
     1160!            ! 
    11571161! 
    11581162!            ! dsal_ocn (psu kg m^-2) / (kg m^-3 * m) 
     
    11651169! 
    11661170!           ! put increments in for levels in the mixed layer 
    1167 !           ! but prevent salinity below a threshold value  
    1168 ! 
    1169 !                   DO jk = 1, jkmax               
    1170 ! 
    1171 !                     IF (dsal_ocn > 0.0_wp .or. sb(ji,jj,jk)+dsal_ocn > sal_thresh) THEN  
     1171!           ! but prevent salinity below a threshold value 
     1172! 
     1173!                   DO jk = 1, jkmax 
     1174! 
     1175!                     IF (dsal_ocn > 0.0_wp .or. sb(ji,jj,jk)+dsal_ocn > sal_thresh) THEN 
    11721176!                           sb(ji,jj,jk) = sb(ji,jj,jk) + dsal_ocn 
    11731177!                           sn(ji,jj,jk) = sn(ji,jj,jk) + dsal_ocn 
     
    11801184!      ! 
    11811185!      !! Adjust fsalt. A +ve fsalt means adding salt to ocean 
    1182 !      !!           fsalt(ji,jj) =  fsalt(ji,jj) + zpmess     ! adjust fsalt   
    1183 !      !!                
    1184 !      !!           emps(ji,jj) = emps(ji,jj) + zpmess        ! or adjust emps (see icestp1d)  
     1186!      !!           fsalt(ji,jj) =  fsalt(ji,jj) + zpmess     ! adjust fsalt 
     1187!      !! 
     1188!      !!           emps(ji,jj) = emps(ji,jj) + zpmess        ! or adjust emps (see icestp1d) 
    11851189!      !!                                                     ! E-P (kg m-2 s-2) 
    11861190!      !            emp(ji,jj) = emp(ji,jj) + zpmess          ! E-P (kg m-2 s-2) 
  • branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90

    r3711 r5168  
    1010   !!            8.0  !  2001-09  (M. Levy, M. Ben Jelloul)  istate_uvg 
    1111   !!   NEMO     1.0  !  2003-08  (G. Madec, C. Talandier)  F90: Free form, modules + EEL R5 
    12    !!             -   !  2004-05  (A. Koch-Larrouy)  istate_gyre  
     12   !!             -   !  2004-05  (A. Koch-Larrouy)  istate_gyre 
    1313   !!            2.0  !  2006-07  (S. Masson)  distributed restart using iom 
    1414   !!            3.3  !  2010-10  (C. Ethe) merge TRC-TRA 
    15    !!            3.4  !  2011-04  (G. Madec) Merge of dtatem and dtasal & suppression of tb,tn/sb,sn  
     15   !!            3.4  !  2011-04  (G. Madec) Merge of dtatem and dtasal & suppression of tb,tn/sb,sn 
    1616   !!---------------------------------------------------------------------- 
    1717 
     
    2424   !!   istate_uvg    : initial velocity in geostropic balance 
    2525   !!---------------------------------------------------------------------- 
    26    USE oce             ! ocean dynamics and active tracers  
    27    USE dom_oce         ! ocean space and time domain  
     26   USE oce             ! ocean dynamics and active tracers 
     27   USE dom_oce         ! ocean space and time domain 
    2828   USE daymod          ! calendar 
    2929   USE eosbn2          ! eq. of state, Brunt Vaisala frequency (eos     routine) 
     
    4646   USE wrk_nemo        ! Memory allocation 
    4747   USE timing          ! Timing 
     48   USE asminc 
    4849 
    4950   IMPLICIT NONE 
     
    6566      !!---------------------------------------------------------------------- 
    6667      !!                   ***  ROUTINE istate_init  *** 
    67       !!  
     68      !! 
    6869      !! ** Purpose :   Initialization of the dynamics and tracer fields. 
    6970      !!---------------------------------------------------------------------- 
     
    8384      rhd  (:,:,:  ) = 0.e0 
    8485      rhop (:,:,:  ) = 0.e0 
    85       rn2  (:,:,:  ) = 0.e0  
    86       tsa  (:,:,:,:) = 0.e0     
     86      rn2  (:,:,:  ) = 0.e0 
     87      tsa  (:,:,:,:) = 0.e0 
    8788 
    8889      IF( ln_rstart ) THEN                    ! Restart from a file 
     
    100101         CALL day_init                           ! model calendar (using both namelist and restart infos) 
    101102         !                                       ! Initialization of ocean to zero 
    102          !   before fields      !       now fields      
     103         !   before fields      !       now fields 
    103104         sshb (:,:)   = 0._wp   ;   sshn (:,:)   = 0._wp 
    104105         ub   (:,:,:) = 0._wp   ;   un   (:,:,:) = 0._wp 
    105          vb   (:,:,:) = 0._wp   ;   vn   (:,:,:) = 0._wp   
     106         vb   (:,:,:) = 0._wp   ;   vn   (:,:,:) = 0._wp 
    106107         rotb (:,:,:) = 0._wp   ;   rotn (:,:,:) = 0._wp 
    107108         hdivb(:,:,:) = 0._wp   ;   hdivn(:,:,:) = 0._wp 
     
    109110         IF( cp_cfg == 'eel' ) THEN 
    110111            CALL istate_eel                      ! EEL   configuration : start from pre-defined U,V T-S fields 
    111          ELSEIF( cp_cfg == 'gyre' ) THEN          
     112         ELSEIF( cp_cfg == 'gyre' ) THEN 
    112113            CALL istate_gyre                     ! GYRE  configuration : start from pre-defined T-S fields 
    113114         ELSEIF( ln_tsd_init      ) THEN         ! Initial T-S fields read in files 
     
    124125            &                                       rhd, gru , grv  )   ! of t,s,rd at ocean bottom 
    125126#endif 
    126          !    
     127         ! 
    127128         ! - ML - sshn could be modified by istate_eel, so that initialization of fse3t_b is done here 
    128129         IF( lk_vvl ) THEN 
     
    133134         !                                       ! define e3u_b, e3v_b from e3t_b initialized in domzgr 
    134135         CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 
    135          !  
     136         ! 
    136137      ENDIF 
    137138      ! 
     
    146147      ENDIF 
    147148      ! 
     149      IF ( lk_asminc .AND. ln_asmdin ) THEN 
     150         neuler = 0            ! restart with an euler from the corrected background 
     151      ENDIF 
     152 
    148153      IF( nn_timing == 1 )  CALL timing_stop('istate_init') 
    149154      ! 
     
    153158      !!--------------------------------------------------------------------- 
    154159      !!                  ***  ROUTINE istate_t_s  *** 
    155       !!    
    156       !! ** Purpose :   Intialization of the temperature field with an  
     160      !! 
     161      !! ** Purpose :   Intialization of the temperature field with an 
    157162      !!      analytical profile or a file (i.e. in EEL configuration) 
    158163      !! 
     
    184189      !!---------------------------------------------------------------------- 
    185190      !!                   ***  ROUTINE istate_eel  *** 
    186       !!  
     191      !! 
    187192      !! ** Purpose :   Initialization of the dynamics and tracers for EEL R5 
    188193      !!      configuration (channel with or without a topographic bump) 
     
    195200      USE divcur     ! hor. divergence & rel. vorticity      (div_cur routine) 
    196201      USE iom 
    197   
     202 
    198203      INTEGER  ::   inum              ! temporary logical unit 
    199204      INTEGER  ::   ji, jj, jk        ! dummy loop indices 
     
    207212      !!---------------------------------------------------------------------- 
    208213 
    209       SELECT CASE ( jp_cfg )  
     214      SELECT CASE ( jp_cfg ) 
    210215         !                                              ! ==================== 
    211216         CASE ( 5 )                                     ! EEL R5 configuration 
     
    244249            ! set the dynamics: U,V, hdiv, rot (and ssh if necessary) 
    245250            ! ---------------- 
    246             ! Start EEL5 configuration with barotropic geostrophic velocities  
     251            ! Start EEL5 configuration with barotropic geostrophic velocities 
    247252            ! according the sshb and sshn SSH imposed. 
    248253            ! we assume a uniform grid (hence the use of e1t(1,1) for delta_y) 
    249             ! we use the Coriolis frequency at mid-channel.    
     254            ! we use the Coriolis frequency at mid-channel. 
    250255            ub(:,:,:) = zueel * umask(:,:,:) 
    251256            un(:,:,:) = ub(:,:,:) 
     
    254259            ! 
    255260            DO jj = 1, jpjglo 
    256                zssh(:,jj) = - (FLOAT(jj)- FLOAT(jpjglo-1)/2.)*zueel*e1t(1,1)*zfcor/grav  
     261               zssh(:,jj) = - (FLOAT(jj)- FLOAT(jpjglo-1)/2.)*zueel*e1t(1,1)*zfcor/grav 
    257262            END DO 
    258263            ! 
     
    273278            sshn(:,:) = sshb(:,:)                   ! set now ssh to the before value 
    274279            ! 
    275             IF( nn_rstssh /= 0 ) THEN   
    276                nn_rstssh = 0                        ! hand-made initilization of ssh  
     280            IF( nn_rstssh /= 0 ) THEN 
     281               nn_rstssh = 0                        ! hand-made initilization of ssh 
    277282               CALL ctl_warn( 'istate_eel: force nn_rstssh = 0' ) 
    278283            ENDIF 
     
    326331      !!---------------------------------------------------------------------- 
    327332      !!                   ***  ROUTINE istate_gyre  *** 
    328       !!  
     333      !! 
    329334      !! ** Purpose :   Initialization of the dynamics and tracers for GYRE 
    330335      !!      configuration (double gyre with rotated domain) 
     
    351356                       &           * (-TANH( (500-fsdept(ji,jj,jk)) / 150 ) + 1) / 2               & 
    352357                       &       + (      15. * ( 1. - TANH( (fsdept(ji,jj,jk)-50.) / 1500.) )       & 
    353                        &                - 1.4 * TANH((fsdept(ji,jj,jk)-100.) / 100.)               &     
    354                        &                + 7.  * (1500. - fsdept(ji,jj,jk)) / 1500.             )   &  
     358                       &                - 1.4 * TANH((fsdept(ji,jj,jk)-100.) / 100.)               & 
     359                       &                + 7.  * (1500. - fsdept(ji,jj,jk)) / 1500.             )   & 
    355360                       &           * (-TANH( (fsdept(ji,jj,jk) - 500) / 150) + 1) / 2 
    356361                  tsn(ji,jj,jk,jp_tem) = tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) 
     
    363368                     &                + 0.2  * TANH( (fsdept(ji,jj,jk) - 35.  ) / 100. )       & 
    364369                     &                + 0.2  * TANH( (fsdept(ji,jj,jk) - 1000.) / 5000.)    )  & 
    365                      &              * (-TANH((fsdept(ji,jj,jk) - 500) / 150) + 1) / 2  
     370                     &              * (-TANH((fsdept(ji,jj,jk) - 500) / 150) + 1) / 2 
    366371                  tsn(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk) 
    367372                  tsb(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) 
     
    379384         ! ---------------------- 
    380385         CALL iom_open ( 'data_tem', inum ) 
    381          CALL iom_get ( inum, jpdom_data, 'votemper', tsn(:,:,:,jp_tem) )  
     386         CALL iom_get ( inum, jpdom_data, 'votemper', tsn(:,:,:,jp_tem) ) 
    382387         CALL iom_close( inum ) 
    383388 
    384          tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) * tmask(:,:,:)  
     389         tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) * tmask(:,:,:) 
    385390         tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) 
    386391 
     
    388393         ! ------------------- 
    389394         CALL iom_open ( 'data_sal', inum ) 
    390          CALL iom_get ( inum, jpdom_data, 'vosaline', tsn(:,:,:,jp_sal) )  
     395         CALL iom_get ( inum, jpdom_data, 'vosaline', tsn(:,:,:,jp_sal) ) 
    391396         CALL iom_close( inum ) 
    392397 
    393          tsn(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) * tmask(:,:,:)  
     398         tsn(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) * tmask(:,:,:) 
    394399         tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 
    395400 
     
    412417      !! ** Purpose :   Compute the geostrophic velocities from (tn,sn) fields 
    413418      !! 
    414       !! ** Method  :   Using the hydrostatic hypothesis the now hydrostatic  
     419      !! ** Method  :   Using the hydrostatic hypothesis the now hydrostatic 
    415420      !!      pressure is computed by integrating the in-situ density from the 
    416421      !!      surface to the bottom. 
     
    429434      CALL wrk_alloc( jpi, jpj, jpk, zprn) 
    430435      ! 
    431       IF(lwp) WRITE(numout,*)  
     436      IF(lwp) WRITE(numout,*) 
    432437      IF(lwp) WRITE(numout,*) 'istate_uvg : Start from Geostrophy' 
    433438      IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
     
    437442 
    438443      zalfg = 0.5 * grav * rau0 
    439        
     444 
    440445      zprn(:,:,1) = zalfg * fse3w(:,:,1) * ( 1 + rhd(:,:,1) )       ! Surface value 
    441446 
     
    443448         zprn(:,:,jk) = zprn(:,:,jk-1)   & 
    444449            &         + zalfg * fse3w(:,:,jk) * ( 2. + rhd(:,:,jk) + rhd(:,:,jk-1) ) 
    445       END DO   
     450      END DO 
    446451 
    447452      ! Compute geostrophic balance 
     
    485490      CALL lbc_lnk( un, 'U', -1. ) 
    486491      CALL lbc_lnk( vn, 'V', -1. ) 
    487        
     492 
    488493      ub(:,:,:) = un(:,:,:) 
    489494      vb(:,:,:) = vn(:,:,:) 
    490        
     495 
    491496      ! WARNING !!!!! 
    492497      ! after initializing u and v, we need to calculate the initial streamfunction bsf. 
     
    516521      un(:,:,:) = ub(:,:,:) 
    517522      vn(:,:,:) = vb(:,:,:) 
    518         
     523 
    519524      ! Compute the divergence and curl 
    520525 
  • branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

    r4600 r5168  
    7777#endif 
    7878   USE tamtrj          ! Output trajectory, needed for TAM 
     79   USE asminc 
    7980 
    8081   IMPLICIT NONE 
     
    110111      !!              Madec, 2008, internal report, IPSL. 
    111112      !!---------------------------------------------------------------------- 
    112       INTEGER ::   istp       ! time step index 
     113      INTEGER ::   istp, jk      ! time step index 
    113114      !!---------------------------------------------------------------------- 
    114115      ! 
     
    155156            CALL Agrif_Step( stp )           ! AGRIF: time stepping 
    156157#else 
     158            !-------------------------------------------------------------! 
     159            ! This trick ensures a minimum consistency between mixing coef. 
     160            ! and other variables 
     161            IF( lk_zdftke .AND. lk_asminc .AND. ln_asmdin .AND. ( istp == nit000) ) THEN 
     162            CALL tke_avn                        ! recompute avt, avm, 
     163                                                ! avmu, avmv and dissl (approximation) 
     164               DO jk = nit000, nit000 + 2 
     165                  CALL zdf_tke( jk )   ; 
     166               END DO 
     167 
     168            END IF 
     169            !-------------------------------------------------------------! 
     170 
    157171            CALL stp( istp )                 ! standard time stepping 
    158172#endif 
     
    322336                            CALL dyn_nept_init  ! simplified form of Neptune effect 
    323337                            CALL flush(numout) 
    324  
     338                            CALL asm_inc_rea_nam ! to fix properly neuler value 
     339                                                 ! in istate_init 
    325340                            CALL  istate_init   ! ocean initial state (Dynamics and tracers) 
    326341 
Note: See TracChangeset for help on using the changeset viewer.