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 6020 for branches/UKMO/icebergs_restart_single_file/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90 – NEMO

Ignore:
Timestamp:
2015-12-08T12:39:53+01:00 (8 years ago)
Author:
timgraham
Message:

Merged with head of trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/icebergs_restart_single_file/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r6019 r6020  
    1717   !!                 ! 2005-11  (V. Garnier) Surface pressure gradient organization 
    1818   !!            3.2  ! 2008-11  (B. Lemaire) creation from old diawri 
     19   !!            3.7  ! 2014-01  (G. Madec) remove eddy induced velocity from no-IOM output 
     20   !!                 !                     change name of output variables in dia_wri_state 
    1921   !!---------------------------------------------------------------------- 
    2022 
     
    2729   USE dynadv, ONLY: ln_dynadv_vec 
    2830   USE zdf_oce         ! ocean vertical physics 
    29    USE ldftra_oce      ! ocean active tracers: lateral physics 
    30    USE ldfdyn_oce      ! ocean dynamics: lateral physics 
    31    USE traldf_iso_grif, ONLY : psix_eiv, psiy_eiv 
    32    USE sol_oce         ! solver variables 
     31   USE ldftra          ! lateral physics: eddy diffusivity coef. 
    3332   USE sbc_oce         ! Surface boundary condition: ocean fields 
    3433   USE sbc_ice         ! Surface boundary condition: ice fields 
     
    4645   USE iom 
    4746   USE ioipsl 
    48    USE dynspg_oce, ONLY: un_adv, vn_adv ! barotropic velocities      
    4947 
    5048#if defined key_lim2 
     
    206204         CALL iom_put( "sbu", z2d )                ! bottom i-current 
    207205      ENDIF 
    208 #if defined key_dynspg_ts 
    209       CALL iom_put(  "ubar", un_adv(:,:)      )    ! barotropic i-current 
    210 #else 
    211       CALL iom_put(  "ubar", un_b(:,:)        )    ! barotropic i-current 
    212 #endif 
     206 
     207      IF ( ln_dynspg_ts ) THEN 
     208         CALL iom_put(  "ubar", un_adv(:,:)      )    ! barotropic i-current 
     209      ELSE 
     210         CALL iom_put(  "ubar", un_b(:,:)        )    ! barotropic i-current 
     211      ENDIF 
    213212       
    214213      CALL iom_put( "voce", vn(:,:,:)         )    ! 3D j-current 
     
    223222         CALL iom_put( "sbv", z2d )                ! bottom j-current 
    224223      ENDIF 
    225 #if defined key_dynspg_ts 
    226       CALL iom_put(  "vbar", vn_adv(:,:)      )    ! barotropic j-current 
    227 #else 
    228       CALL iom_put(  "vbar", vn_b(:,:)        )    ! barotropic j-current 
    229 #endif 
     224 
     225      IF ( ln_dynspg_ts ) THEN 
     226         CALL iom_put(  "vbar", vn_adv(:,:)      )    ! barotropic j-current 
     227      ELSE 
     228         CALL iom_put(  "vbar", vn_b(:,:)        )    ! barotropic j-current 
     229      ENDIF 
    230230 
    231231      CALL iom_put( "woce", wn )                   ! vertical velocity 
    232232      IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN   ! vertical mass transport & its square value 
    233233         ! Caution: in the VVL case, it only correponds to the baroclinic mass transport. 
    234          z2d(:,:) = rau0 * e12t(:,:) 
     234         z2d(:,:) = rau0 * e1e2t(:,:) 
    235235         DO jk = 1, jpk 
    236236            z3d(:,:,jk) = wn(:,:,jk) * z2d(:,:) 
     
    247247         DO jj = 2, jpjm1                                    ! sst gradient 
    248248            DO ji = fs_2, fs_jpim1   ! vector opt. 
    249                zztmp      = tsn(ji,jj,1,jp_tem) 
    250                zztmpx     = ( tsn(ji+1,jj  ,1,jp_tem) - zztmp ) / e1u(ji,jj) + ( zztmp - tsn(ji-1,jj  ,1,jp_tem) ) / e1u(ji-1,jj  ) 
    251                zztmpy     = ( tsn(ji  ,jj+1,1,jp_tem) - zztmp ) / e2v(ji,jj) + ( zztmp - tsn(ji  ,jj-1,1,jp_tem) ) / e2v(ji  ,jj-1) 
     249               zztmp  = tsn(ji,jj,1,jp_tem) 
     250               zztmpx = ( tsn(ji+1,jj,1,jp_tem) - zztmp ) * r1_e1u(ji,jj) + ( zztmp - tsn(ji-1,jj  ,1,jp_tem) ) * r1_e1u(ji-1,jj) 
     251               zztmpy = ( tsn(ji,jj+1,1,jp_tem) - zztmp ) * r1_e2v(ji,jj) + ( zztmp - tsn(ji  ,jj-1,1,jp_tem) ) * r1_e2v(ji,jj-1) 
    252252               z2d(ji,jj) = 0.25 * ( zztmpx * zztmpx + zztmpy * zztmpy )   & 
    253253                  &              * umask(ji,jj,1) * umask(ji-1,jj,1) * vmask(ji,jj,1) * umask(ji,jj-1,1) 
     
    401401      !!      Each nwrite time step, output the instantaneous or mean fields 
    402402      !!---------------------------------------------------------------------- 
    403       !! 
    404       INTEGER, INTENT( in ) ::   kt      ! ocean time-step index 
    405       !! 
     403      INTEGER, INTENT( in ) ::   kt   ! ocean time-step index 
     404      ! 
    406405      LOGICAL ::   ll_print = .FALSE.                        ! =T print and flush numout 
    407406      CHARACTER (len=40) ::   clhstnam, clop, clmx           ! local names 
     
    412411      INTEGER  ::   jn, ierror                               ! local integers 
    413412      REAL(wp) ::   zsto, zout, zmax, zjulian, zdt           ! local scalars 
    414       !! 
     413      ! 
    415414      REAL(wp), POINTER, DIMENSION(:,:)   :: zw2d       ! 2D workspace 
    416415      REAL(wp), POINTER, DIMENSION(:,:,:) :: zw3d       ! 3D workspace 
     
    419418      IF( nn_timing == 1 )   CALL timing_start('dia_wri') 
    420419      ! 
    421       CALL wrk_alloc( jpi , jpj      , zw2d ) 
    422       IF ( ln_traldf_gdia .OR. lk_vvl )  call wrk_alloc( jpi , jpj , jpk  , zw3d ) 
     420                     CALL wrk_alloc( jpi,jpj      , zw2d ) 
     421      IF( lk_vvl )   CALL wrk_alloc( jpi,jpj,jpk  , zw3d ) 
    423422      ! 
    424423      ! Output the initial state and forcings 
     
    438437      zdt = rdt 
    439438      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 
     439      clop = "x"         ! no use of the mask value (require less cpu time and otherwise the model crashes) 
    443440#if defined key_diainstant 
    444441      zsto = nwrite * zdt 
     
    659656          
    660657         clmx ="l_max(only(x))"    ! max index on a period 
    661          CALL histdef( nid_T, "sobowlin", "Bowl Index"                         , "W-point",   &  ! bowl INDEX  
    662             &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clmx, zsto, zout ) 
     658!         CALL histdef( nid_T, "sobowlin", "Bowl Index"                         , "W-point",   &  ! bowl INDEX  
     659!            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clmx, zsto, zout ) 
    663660#if defined key_diahth 
    664661         CALL histdef( nid_T, "sothedep", "Thermocline Depth"                  , "m"      ,   & ! hth 
     
    684681         CALL histdef( nid_U, "vozocrtx", "Zonal Current"                      , "m/s"    ,   &  ! un 
    685682            &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
    686          IF( ln_traldf_gdia ) THEN 
    687             CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv 
    688                  &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
    689          ELSE 
    690 #if defined key_diaeiv 
    691             CALL histdef( nid_U, "vozoeivu", "Zonal EIV Current"                  , "m/s"    ,   &  ! u_eiv 
    692             &          jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 
    693 #endif 
    694          END IF 
    695683         !                                                                                      !!! nid_U : 2D 
    696684         CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! utau 
     
    702690         CALL histdef( nid_V, "vomecrty", "Meridional Current"                 , "m/s"    ,   &  ! vn 
    703691            &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
    704          IF( ln_traldf_gdia ) THEN 
    705             CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv 
    706                  &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
    707          ELSE  
    708 #if defined key_diaeiv 
    709             CALL histdef( nid_V, "vomeeivv", "Meridional EIV Current"             , "m/s"    ,   &  ! v_eiv 
    710             &          jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 
    711 #endif 
    712          END IF 
    713692         !                                                                                      !!! nid_V : 2D 
    714693         CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! vtau 
     
    720699         CALL histdef( nid_W, "vovecrtz", "Vertical Velocity"                  , "m/s"    ,   &  ! wn 
    721700            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
    722          IF( ln_traldf_gdia ) THEN 
    723             CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity"              , "m/s"    ,   &  ! w_eiv 
    724                  &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
    725          ELSE 
    726 #if defined key_diaeiv 
    727             CALL histdef( nid_W, "voveeivw", "Vertical EIV Velocity"              , "m/s"    ,   &  ! w_eiv 
    728                  &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
    729 #endif 
    730          END IF 
    731701         CALL histdef( nid_W, "votkeavt", "Vertical Eddy Diffusivity"          , "m2/s"   ,   &  ! avt 
    732702            &          jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 
     
    739709         ENDIF 
    740710         !                                                                                      !!! nid_W : 2D 
    741 #if defined key_traldf_c2d 
    742          CALL histdef( nid_W, "soleahtw", "lateral eddy diffusivity"           , "m2/s"   ,   &  ! ahtw 
    743             &          jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout ) 
    744 # if defined key_traldf_eiv  
    745             CALL histdef( nid_W, "soleaeiw", "eddy induced vel. coeff. at w-point", "m2/s",   &  ! aeiw 
    746                &       jpi, jpj, nh_W, 1  , 1, 1  , - 99, 32, clop, zsto, zout ) 
    747 # endif 
    748 #endif 
    749  
    750711         CALL histend( nid_W, snc4chunks=snc4set ) 
    751712 
     
    855816 
    856817      CALL histwrite( nid_U, "vozocrtx", it, un            , ndim_U , ndex_U )    ! i-current 
    857       IF( ln_traldf_gdia ) THEN 
    858          IF (.not. ALLOCATED(psix_eiv))THEN 
    859             ALLOCATE( psix_eiv(jpi,jpj,jpk) , psiy_eiv(jpi,jpj,jpk) , STAT=ierr ) 
    860             IF( lk_mpp   )   CALL mpp_sum ( ierr ) 
    861             IF( ierr > 0 )   CALL ctl_stop('STOP', 'diawri: unable to allocate psi{x,y}_eiv') 
    862             psix_eiv(:,:,:) = 0.0_wp 
    863             psiy_eiv(:,:,:) = 0.0_wp 
    864          ENDIF 
    865          DO jk=1,jpkm1 
    866             zw3d(:,:,jk) = (psix_eiv(:,:,jk+1) - psix_eiv(:,:,jk))/fse3u(:,:,jk)  ! u_eiv = -dpsix/dz 
    867          END DO 
    868          zw3d(:,:,jpk) = 0._wp 
    869          CALL histwrite( nid_U, "vozoeivu", it, zw3d, ndim_U , ndex_U )           ! i-eiv current 
    870       ELSE 
    871 #if defined key_diaeiv 
    872          CALL histwrite( nid_U, "vozoeivu", it, u_eiv, ndim_U , ndex_U )          ! i-eiv current 
    873 #endif 
    874       ENDIF 
    875818      CALL histwrite( nid_U, "sozotaux", it, utau          , ndim_hU, ndex_hU )   ! i-wind stress 
    876819 
    877820      CALL histwrite( nid_V, "vomecrty", it, vn            , ndim_V , ndex_V  )   ! j-current 
    878       IF( ln_traldf_gdia ) THEN 
    879          DO jk=1,jpk-1 
    880             zw3d(:,:,jk) = (psiy_eiv(:,:,jk+1) - psiy_eiv(:,:,jk))/fse3v(:,:,jk)  ! v_eiv = -dpsiy/dz 
    881          END DO 
    882          zw3d(:,:,jpk) = 0._wp 
    883          CALL histwrite( nid_V, "vomeeivv", it, zw3d, ndim_V , ndex_V )           ! j-eiv current 
    884       ELSE 
    885 #if defined key_diaeiv 
    886          CALL histwrite( nid_V, "vomeeivv", it, v_eiv, ndim_V , ndex_V )          ! j-eiv current 
    887 #endif 
    888       ENDIF 
    889821      CALL histwrite( nid_V, "sometauy", it, vtau          , ndim_hV, ndex_hV )   ! j-wind stress 
    890822 
    891823      CALL histwrite( nid_W, "vovecrtz", it, wn             , ndim_T, ndex_T )    ! vert. current 
    892       IF( ln_traldf_gdia ) THEN 
    893          DO jk=1,jpk-1 
    894             DO jj = 2, jpjm1 
    895                DO ji = fs_2, fs_jpim1  ! vector opt. 
    896                   zw3d(ji,jj,jk) = (psiy_eiv(ji,jj,jk) - psiy_eiv(ji,jj-1,jk))/e2v(ji,jj) + & 
    897                        &    (psix_eiv(ji,jj,jk) - psix_eiv(ji-1,jj,jk))/e1u(ji,jj) ! w_eiv = dpsiy/dy + dpsiy/dx 
    898                END DO 
    899             END DO 
    900          END DO 
    901          zw3d(:,:,jpk) = 0._wp 
    902          CALL histwrite( nid_W, "voveeivw", it, zw3d          , ndim_T, ndex_T )    ! vert. eiv current 
    903       ELSE 
    904 #   if defined key_diaeiv 
    905          CALL histwrite( nid_W, "voveeivw", it, w_eiv          , ndim_T, ndex_T )    ! vert. eiv current 
    906 #   endif 
    907       ENDIF 
    908824      CALL histwrite( nid_W, "votkeavt", it, avt            , ndim_T, ndex_T )    ! T vert. eddy diff. coef. 
    909825      CALL histwrite( nid_W, "votkeavm", it, avmu           , ndim_T, ndex_T )    ! T vert. eddy visc. coef. 
     
    911827         CALL histwrite( nid_W, "voddmavs", it, fsavs(:,:,:), ndim_T, ndex_T )    ! S vert. eddy diff. coef. 
    912828      ENDIF 
    913 #if defined key_traldf_c2d 
    914       CALL histwrite( nid_W, "soleahtw", it, ahtw          , ndim_hT, ndex_hT )   ! lateral eddy diff. coef. 
    915 # if defined key_traldf_eiv 
    916          CALL histwrite( nid_W, "soleaeiw", it, aeiw       , ndim_hT, ndex_hT )   ! EIV coefficient at w-point 
    917 # endif 
    918 #endif 
    919829 
    920830      ! 3. Close all files 
     
    927837      ENDIF 
    928838      ! 
    929       CALL wrk_dealloc( jpi , jpj      , zw2d ) 
    930       IF ( ln_traldf_gdia .OR. lk_vvl )  call wrk_dealloc( jpi , jpj , jpk  , zw3d ) 
     839                     CALL wrk_dealloc( jpi , jpj        , zw2d ) 
     840      IF( lk_vvl )   CALL wrk_dealloc( jpi , jpj , jpk  , zw3d ) 
    931841      ! 
    932842      IF( nn_timing == 1 )   CALL timing_stop('dia_wri') 
     
    960870      !!---------------------------------------------------------------------- 
    961871      !  
    962 !     IF( nn_timing == 1 )   CALL timing_start('dia_wri_state') ! not sure this works for routines not called in first timestep 
    963  
    964872      ! 0. Initialisation 
    965873      ! ----------------- 
     
    1020928         CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      ,   &   ! t-point depth 
    1021929            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
    1022       END IF 
     930      ENDIF 
    1023931 
    1024932#if defined key_lim2 
     
    1044952      CALL histwrite( id_i, "vomecrty", kt, vn               , jpi*jpj*jpk, idex )    ! now j-velocity 
    1045953      CALL histwrite( id_i, "vovecrtz", kt, wn               , jpi*jpj*jpk, idex )    ! now k-velocity 
    1046       CALL histwrite( id_i, "sowaflup", kt, (emp-rnf )       , jpi*jpj    , idex )    ! freshwater budget 
     954      CALL histwrite( id_i, "sowaflup", kt, emp-rnf          , jpi*jpj    , idex )    ! freshwater budget 
    1047955      CALL histwrite( id_i, "sohefldo", kt, qsr + qns        , jpi*jpj    , idex )    ! total heat flux 
    1048956      CALL histwrite( id_i, "soshfldo", kt, qsr              , jpi*jpj    , idex )    ! solar heat flux 
     
    1062970      ENDIF 
    1063971#endif 
    1064         
    1065 !     IF( nn_timing == 1 )   CALL timing_stop('dia_wri_state') ! not sure this works for routines not called in first timestep 
    1066972      !  
    1067  
    1068973   END SUBROUTINE dia_wri_state 
    1069974   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.