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 7795 for branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC – NEMO

Ignore:
Timestamp:
2017-03-15T09:04:30+01:00 (7 years ago)
Author:
cbricaud
Message:

code cleaning and correct bug for wn computing in vvl case

Location:
branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/CRS/crs.F90

    r7398 r7795  
    154154      ! Physical and dynamical ocean fields for output or passing to TOP, time-mean fields 
    155155      REAL(wp), DIMENSION(:,:,:,:), ALLOCATABLE      :: tsb_crs,tsn_crs,tsa_crs,rab_crs_n 
    156       REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: un_crs, vn_crs, wn_crs, rke_crs 
     156      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: un_crs, vn_crs, wn_crs 
    157157      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: ub_crs, vb_crs 
    158158      REAL(wp), DIMENSION(:,:,:)  , ALLOCATABLE      :: hdivb_crs , hdivn_crs     
     
    296296 
    297297      ALLOCATE( ub_crs(jpi_crs,jpj_crs,jpk) , vb_crs(jpi_crs,jpj_crs,jpk) , & 
    298          &      un_crs(jpi_crs,jpj_crs,jpk) , vn_crs(jpi_crs,jpj_crs,jpk)    ,  wn_crs(jpi_crs,jpj_crs,jpk) , & 
     298         &      un_crs(jpi_crs,jpj_crs,jpk) , vn_crs(jpi_crs,jpj_crs,jpk)    ,  wn_crs(jpi_crs,jpj_crs,jpk) ,& 
    299299         &      hdivb_crs(jpi_crs,jpj_crs,jpk) , hdivn_crs(jpi_crs,jpj_crs,jpk) , & 
    300          &      rke_crs(jpi_crs,jpj_crs,jpk), rhop_crs(jpi_crs,jpj_crs,jpk)  , & 
     300         &      rhop_crs(jpi_crs,jpj_crs,jpk)  , & 
    301301         &      rb2_crs(jpi_crs,jpj_crs,jpk) ,rn2_crs(jpi_crs,jpj_crs,jpk) , & 
    302302         &      rhd_crs(jpi_crs,jpj_crs,jpk)   , rab_crs_n(jpi_crs,jpj_crs,jpk,jpts) , & 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/CRS/crsfld.F90

    r7520 r7795  
    141141         emp_b_crs(:,:)        = emp_crs(:,:) 
    142142         rnf_b_crs(:,:)        = rnf_crs(:,:) 
    143          hdivb_crs(:,:,:)      = hdivn_crs(:,:,:) 
    144143      ELSE 
    145144         emp_b_crs(:,:    ) = 0._wp 
    146145         rnf_b_crs(:,:    ) = 0._wp 
    147          hdivb_crs(:,:,:  ) = 0._wp 
    148146      ENDIF 
    149147 
     
    157155      CALL crs_dom_ope( zs, 'VOL', 'T', tmask, zs_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) 
    158156      tsb_crs(:,:,:,jp_sal) = zs_crs(:,:,:) 
    159  
    160       !  U-velocity 
    161       CALL crs_dom_ope( ub, 'SUM', 'U', umask, ub_crs, p_e12=e2u, p_e3=zfse3u, p_surf_crs=e2e3u_msk, psgn=-1.0 ) 
    162  
    163       !  V-velocity 
    164       CALL crs_dom_ope( vb, 'SUM', 'V', vmask, vb_crs, p_e12=e1v, p_e3=zfse3v, p_surf_crs=e1e3v_msk, psgn=-1.0 ) 
    165157 
    166158      ! n2 
     
    266258      CALL iom_put( "voce"  , vn_crs )   ! i-current  
    267259 
    268       !n2 
    269       CALL crs_dom_ope( rn2 , 'VOL', 'W', tmask, rn2_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) 
     260!      !n2 
     261!      CALL crs_dom_ope( rn2 , 'VOL', 'W', tmask, rn2_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) 
    270262      
    271263      !  Horizontal divergence ( following OPA_SRC/DYN/divcur.F90 )  
     
    303295      ! 
    304296      CALL iom_put( "avt", avt_crs )   !  Kz 
    305       
     297 
    306298      !2D fields 
    307       CALL crs_dom_ope( utau , 'SUM', 'U', umask, utau_crs , p_e12=e2u  , p_surf_crs=e2u_crs  , psgn=1.0 ) 
    308       CALL crs_dom_ope( vtau , 'SUM', 'V', vmask, vtau_crs , p_e12=e1v  , p_surf_crs=e1v_crs  , psgn=1.0 ) 
    309       CALL crs_dom_ope( wndm , 'SUM', 'T', tmask, wndm_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) 
    310299      CALL crs_dom_ope( rnf  , 'MAX', 'T', tmask, rnf_crs                                     , psgn=1.0 ) 
    311300      CALL crs_dom_ope( h_rnf, 'MAX', 'T', tmask, h_rnf_crs                                   , psgn=1.0 ) 
     
    319308      CALL crs_dom_ope( emp   ,'SUM', 'T', tmask, emp_crs   , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) 
    320309      CALL crs_dom_ope( fmmflx,'SUM', 'T', tmask, fmmflx_crs, p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) 
    321       CALL crs_dom_ope( sfx   ,'SUM', 'T', tmask, sfx_crs   , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) 
    322310 
    323311      CALL crs_dom_ope( fr_i  ,'SUM', 'T', tmask, fr_i_crs  , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) 
     
    385373      DO jk = jpkm1, 1, -1 
    386374         wn_crs(:,:,jk) = wn_crs(:,:,jk+1)*e1e2w_msk(:,:,jk+1) - (  hdivn_crs(:,:,jk)                                   & 
    387                &                          + z1_2dt * e1e2w_crs(:,:,jk) * ( e3t_a_crs(:,:,jk) - e3t_b_crs(:,:,jk) ) ) * tmask_crs(:,:,jk) 
     375               &                          + z1_2dt * e1e2w_msk(:,:,jk) * ( e3t_a_crs(:,:,jk) - e3t_b_crs(:,:,jk) ) ) * tmask_crs(:,:,jk) 
    388376         WHERE( e1e2w_msk(:,:,jk) .NE. 0._wp )  wn_crs(:,:,jk) =  wn_crs(:,:,jk) /e1e2w_msk(:,:,jk) 
    389377 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_tvd_crs.F90

    r6772 r7795  
    144144         ! Surface value 
    145145         IF( lk_vvl ) THEN   ;   zwz(:,:, 1 ) = 0.e0                         ! volume variable 
    146          ELSE                ;   zwz(:,:, 1 ) = pwn(:,:,1) !cbr * ptb(:,:,1,jn)   ! linear free surface  
     146         ELSE                ;   zwz(:,:, 1 ) = pwn(:,:,1) * ptb(:,:,1,jn)   ! linear free surface  
    147147         ENDIF 
    148148         ! Interior value 
  • branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/TRA/trazdf_imp_crs.F90

    r6772 r7795  
    149149 
    150150#endif 
     151 
    151152            DO jk = 1, jpkm1 
    152153               DO jj = 2, jpj_crs-1 
     
    154155 
    155156#if defined key_vvl 
    156                      ze3ta =  ( 1. - r_vvl ) +        r_vvl   * fse3t_a_crs(ji,jj,jk)   ! after scale factor at T-point 
    157                      ze3tn =         r_vvl   + ( 1. - r_vvl ) * fse3t_n_crs(ji,jj,jk)   ! now   scale factor at T-point 
     157                     ze3ta =  ( 1. - r_vvl ) +        r_vvl   * e3t_a_crs(ji,jj,jk)   ! after scale factor at T-point 
     158                     ze3tn =         r_vvl   + ( 1. - r_vvl ) * e3t_n_crs(ji,jj,jk)   ! now   scale factor at T-point 
    158159#else 
    159160                     ze3ta =  ( 1. - r_vvl ) +        r_vvl   * e3t_0_crs(ji,jj,jk)   ! after scale factor at T-point 
    160161                     ze3tn =         r_vvl   + ( 1. - r_vvl ) * e3t_0_crs(ji,jj,jk)   ! now   scale factor at T-point 
    161162#endif 
    162                      !cbr zwi(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk  ) / ( ze3tn * e3w_1d(jk  ) )  !cc 
    163                      !cbr zws(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk+1) / ( ze3tn * e3w_1d(jk+1) )  !cc 
     163 
    164164                     zwi(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk  ) / ( ze3tn * fse3w_max_crs(ji,jj,jk) ) 
    165165                     zws(ji,jj,jk) = - p2dt(jk) * zwt(ji,jj,jk+1) / ( ze3tn * fse3w_max_crs(ji,jj,jk+1) ) 
     
    209209            DO ji = 2, jpi_crs-1 
    210210#if defined key_vvl 
    211                ze3tb = ( 1. - r_vvl ) + r_vvl * fse3t_b_crs(ji,jj,1) 
    212                ze3tn = ( 1. - r_vvl ) + r_vvl * fse3t_n_crs(ji,jj,1) 
     211               ze3tb = ( 1. - r_vvl ) + r_vvl * e3t_b_crs(ji,jj,1) 
     212               ze3tn = ( 1. - r_vvl ) + r_vvl * e3t_n_crs(ji,jj,1) 
    213213#else 
    214214               ze3tb = ( 1. - r_vvl ) + r_vvl * e3t_0_crs(ji,jj,1) 
     
    223223               DO ji = 2, jpi_crs-1 
    224224#if defined key_vvl 
    225                   ze3tb = ( 1. - r_vvl ) + r_vvl * fse3t_b_crs(ji,jj,jk) 
    226                   ze3tn = ( 1. - r_vvl ) + r_vvl * fse3t_n_crs(ji,jj,jk) 
     225                  ze3tb = ( 1. - r_vvl ) + r_vvl * e3t_b_crs(ji,jj,jk) 
     226                  ze3tn = ( 1. - r_vvl ) + r_vvl * e3t_n_crs(ji,jj,jk) 
    227227#else 
    228228                  ze3tb = ( 1. - r_vvl ) + r_vvl * e3t_0_crs(ji,jj,jk) 
Note: See TracChangeset for help on using the changeset viewer.