MODULE crsfld !!====================================================================== !! *** MODULE crsdfld *** !! Ocean coarsening : coarse ocean fields !!===================================================================== !! 2012-07 (J. Simeon, C. Calone, G. Madec, C. Ethe) !!---------------------------------------------------------------------- !!---------------------------------------------------------------------- !! crs_fld : create the standard output files for coarse grid and prep !! other variables needed to be passed to TOP !!---------------------------------------------------------------------- USE oce ! ocean dynamics and tracers USE dom_oce ! ocean space and time domain USE ldftra_oce ! ocean active tracers: lateral physics USE sbc_oce ! Surface boundary condition: ocean fields USE sbcrnf USE zdf_oce ! vertical physics: ocean fields USE zdfddm ! vertical physics: double diffusion USe zdfmxl USE lbclnk ! ocean lateral boundary conditions (or mpp link) USE in_out_manager ! I/O manager USE timing ! preformance summary USE wrk_nemo ! working array USE crs USE crsdom USE domvvl USE crslbclnk USE iom USE zdfmxl_crs USE eosbn2 USE zdftke USE ieee_arithmetic IMPLICIT NONE PRIVATE PUBLIC crs_fld ! routines called by step.F90 !! * Substitutions # include "zdfddm_substitute.h90" # include "domzgr_substitute.h90" # include "vectopt_loop_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/OPA 3.3 , NEMO Consortium (2010) !! $Id $ !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS SUBROUTINE crs_fld( kt ) !!--------------------------------------------------------------------- !! *** ROUTINE crs_fld *** !! !! ** Purpose : Basic output of coarsened dynamics and tracer fields !! NETCDF format is used by default !! 1. Accumulate in time the dimensionally-weighted fields !! 2. At time of output, rescale [1] by dimension and time !! to yield the spatial and temporal average. !! !! ** Method : !!---------------------------------------------------------------------- !! INTEGER, INTENT( in ) :: kt ! ocean time-step index !! INTEGER :: ji, jj, jk ! dummy loop indices REAL(wp) :: z2dcrsu, z2dcrsv REAL(wp) :: z1_2dt REAL(wp), POINTER, DIMENSION(:,:,:) :: zfse3t, zfse3u, zfse3v, zfse3w ! 3D workspace for e3 REAL(wp), POINTER, DIMENSION(:,:,:) :: zt, zs , ztmp REAL(wp), POINTER, DIMENSION(:,:) :: z2d,z2d_crs REAL(wp), POINTER, DIMENSION(:,:,:) :: zt_crs, zs_crs !!---------------------------------------------------------------------- IF( nn_timing == 1 ) CALL timing_start('crs_fld') ! Initialize arrays CALL wrk_alloc( jpi, jpj, jpk, zfse3t, zfse3w ) CALL wrk_alloc( jpi, jpj, jpk, zfse3u, zfse3v ) CALL wrk_alloc( jpi, jpj, jpk, zt, zs , ztmp ) CALL wrk_alloc( jpi, jpj, z2d ) ! CALL wrk_alloc( jpi_crs, jpj_crs, jpk, zt_crs, zs_crs ) CALL wrk_alloc( jpi_crs, jpj_crs, z2d_crs ) CALL iom_swap( "nemo_crs" ) ! swap on the coarse grid !--------------------------------------------------------------------------------------------------- !scale factors: before and now !--------------------------------------------------------------------------------------------------- #if defined key_vvl IF( kt /= nit000 )THEN zfse3t(:,:,:) = e3t_b(:,:,:) zfse3u(:,:,:) = e3u_b(:,:,:) zfse3v(:,:,:) = e3v_b(:,:,:) zfse3w(:,:,:) = e3w_b(:,:,:) CALL crs_dom_e3( e1t, e2t, zfse3t, p_sfc_3d_crs=e1e2w_crs, cd_type='T', p_mask=tmask, p_e3_crs=e3t_b_crs, p_e3_max_crs=zs_crs) CALL crs_dom_e3( e1t, e2t, zfse3w, p_sfc_3d_crs=e1e2w_crs, cd_type='W', p_mask=tmask, p_e3_crs=e3w_b_crs, p_e3_max_crs=zs_crs) CALL crs_dom_e3( e1u, e2u, zfse3u, p_sfc_2d_crs=e2u_crs , cd_type='U', p_mask=umask, p_e3_crs=e3u_b_crs, p_e3_max_crs=zs_crs) CALL crs_dom_e3( e1v, e2v, zfse3v, p_sfc_2d_crs=e1v_crs , cd_type='V', p_mask=vmask, p_e3_crs=e3v_b_crs, p_e3_max_crs=zs_crs) DO jk = 1, jpk DO ji = 1, jpi_crs DO jj = 1, jpj_crs IF( e3t_b_crs(ji,jj,jk) == 0._wp ) e3t_b_crs(ji,jj,jk) = e3t_1d(jk) IF( e3w_b_crs(ji,jj,jk) == 0._wp ) e3w_b_crs(ji,jj,jk) = e3w_1d(jk) IF( e3u_b_crs(ji,jj,jk) == 0._wp ) e3u_b_crs(ji,jj,jk) = e3t_1d(jk) IF( e3v_b_crs(ji,jj,jk) == 0._wp ) e3v_b_crs(ji,jj,jk) = e3t_1d(jk) ENDDO ENDDO ENDDO e3t_n_crs(:,:,:) = e3t_a_crs(:,:,:) e3u_n_crs(:,:,:) = e3u_a_crs(:,:,:) e3v_n_crs(:,:,:) = e3v_a_crs(:,:,:) e3w_n_crs(:,:,:) = e3w_a_crs(:,:,:) ENDIF #endif !--------------------------------------------------------------------------------------------------- !variables domaine au temps before : swap !--------------------------------------------------------------------------------------------------- #if defined key_vvl zfse3t(:,:,:) = e3t_b(:,:,:) zfse3u(:,:,:) = e3u_b(:,:,:) zfse3v(:,:,:) = e3v_b(:,:,:) zfse3w(:,:,:) = e3w_b(:,:,:) #else zfse3t(:,:,:) = e3t_0(:,:,:) zfse3u(:,:,:) = e3u_0(:,:,:) zfse3v(:,:,:) = e3v_0(:,:,:) zfse3w(:,:,:) = e3w_0(:,:,:) #endif IF( kt /= nit000 )THEN emp_b_crs(:,:) = emp_crs(:,:) rnf_b_crs(:,:) = rnf_crs(:,:) hdivb_crs(:,:,:) = hdivn_crs(:,:,:) ELSE emp_b_crs(:,: ) = 0._wp rnf_b_crs(:,: ) = 0._wp hdivb_crs(:,:,: ) = 0._wp ENDIF ! Temperature zt(:,:,:) = tsb(:,:,:,jp_tem) ; zt_crs(:,:,:) = 0._wp CALL crs_dom_ope( zt, 'VOL', 'T', tmask, zt_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) tsb_crs(:,:,:,jp_tem) = zt_crs(:,:,:) ! Salinity zs(:,:,:) = tsb(:,:,:,jp_sal) ; zs_crs(:,:,:) = 0._wp CALL crs_dom_ope( zs, 'VOL', 'T', tmask, zs_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) tsb_crs(:,:,:,jp_sal) = zs_crs(:,:,:) ! U-velocity CALL crs_dom_ope( ub, 'SUM', 'U', umask, ub_crs, p_e12=e2u, p_e3=zfse3u, p_surf_crs=e2e3u_msk, psgn=-1.0 ) ! V-velocity CALL crs_dom_ope( vb, 'SUM', 'V', vmask, vb_crs, p_e12=e1v, p_e3=zfse3v, p_surf_crs=e1e3v_msk, psgn=-1.0 ) ! n2 CALL crs_dom_ope( rn2b, 'VOL', 'W', tmask, rb2_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) !ssh zfse3t(:,:,:) = 1._wp CALL crs_dom_ope( sshb , 'VOL', 'T', tmask, sshb_crs , p_e12=e1e2t, p_e3=zfse3t , psgn=1.0 ) !--------------------------------------------------------------------------------------------------- !variables at now time : !--------------------------------------------------------------------------------------------------- #if defined key_vvl zfse3t(:,:,:) = e3t_n(:,:,:) zfse3u(:,:,:) = e3u_n(:,:,:) zfse3v(:,:,:) = e3v_n(:,:,:) zfse3w(:,:,:) = e3w_n(:,:,:) CALL iom_put("e3t",e3t_n_crs) CALL iom_put("e3u",e3u_n_crs) CALL iom_put("e3v",e3v_n_crs) CALL iom_put("e3w",e3w_n_crs) #else zfse3t(:,:,:) = e3t_0(:,:,:) zfse3u(:,:,:) = e3u_0(:,:,:) zfse3v(:,:,:) = e3v_0(:,:,:) zfse3w(:,:,:) = e3w_0(:,:,:) #endif #if defined key_vvl ! surfaces CALL crs_dom_sfc( umask, 'U', e2e3u_crs, e2e3u_msk, p_e2=e2u, p_e3=zfse3u ) CALL crs_dom_sfc( vmask, 'V', e1e3v_crs, e1e3v_msk, p_e1=e1v, p_e3=zfse3v ) !cbr CALL iom_put("e2e3u_crs",e2e3u_crs) !CALL iom_put("e2e3u_msk",e2e3u_msk) !CALL iom_put("e1e3v_crs",e1e3v_crs) !CALL iom_put("e1e3v_msk",e1e3v_msk) ! vertical scale factors CALL crs_dom_e3( e1t, e2t, zfse3t, p_sfc_3d_crs=e1e2w_crs, cd_type='T', p_mask=tmask, p_e3_crs=zs_crs, p_e3_max_crs=e3t_max_n_crs) CALL crs_dom_e3( e1t, e2t, zfse3w, p_sfc_3d_crs=e1e2w_crs, cd_type='W', p_mask=tmask, p_e3_crs=zs_crs, p_e3_max_crs=e3w_max_n_crs) CALL crs_dom_e3( e1u, e2u, zfse3u, p_sfc_2d_crs=e2u_crs , cd_type='U', p_mask=umask, p_e3_crs=zs_crs, p_e3_max_crs=e3u_max_n_crs) CALL crs_dom_e3( e1v, e2v, zfse3v, p_sfc_2d_crs=e1v_crs , cd_type='V', p_mask=vmask, p_e3_crs=zs_crs, p_e3_max_crs=e3v_max_n_crs) DO jk = 1, jpk DO ji = 1, jpi_crs DO jj = 1, jpj_crs IF( e3t_max_n_crs(ji,jj,jk) == 0._wp ) e3t_max_n_crs(ji,jj,jk) = e3t_1d(jk) IF( e3w_max_n_crs(ji,jj,jk) == 0._wp ) e3w_max_n_crs(ji,jj,jk) = e3w_1d(jk) IF( e3u_max_n_crs(ji,jj,jk) == 0._wp ) e3u_max_n_crs(ji,jj,jk) = e3t_1d(jk) IF( e3v_max_n_crs(ji,jj,jk) == 0._wp ) e3v_max_n_crs(ji,jj,jk) = e3t_1d(jk) ENDDO ENDDO ENDDO ! depth CALL crs_dom_ope( gdepw_n, 'MAX', 'T', tmask, gdept_n_crs, p_e3=zfse3t, psgn=1.0 ) CALL crs_dom_ope( gdepw_n, 'MAX', 'W', tmask, gdepw_n_crs, p_e3=zfse3w, psgn=1.0 ) DO jk = 1, jpk DO ji = 1, jpi_crs DO jj = 1, jpj_crs IF( gdept_n_crs(ji,jj,jk) .LE. 0._wp ) gdept_n_crs(ji,jj,jk) = gdept_1d(jk) IF( gdepw_n_crs(ji,jj,jk) .LE. 0._wp ) gdepw_n_crs(ji,jj,jk) = gdepw_1d(jk) ENDDO ENDDO ENDDO ! volume and facvol CALL crs_dom_facvol( tmask, 'T', e1t, e2t, zfse3t, ocean_volume_crs_t, facvol_t ) !cbr CALL iom_put("cvol_crs_t",ocean_volume_crs_t) ! bt_crs(:,:,:) = ocean_volume_crs_t(:,:,:) * facvol_t(:,:,:)*tmask_crs(:,:,:) ! r1_bt_crs(:,:,:) = 0._wp WHERE( bt_crs /= 0._wp ) r1_bt_crs(:,:,:) = 1._wp / bt_crs(:,:,:) CALL crs_dom_facvol( tmask, 'W', e1t, e2t, zfse3w, ocean_volume_crs_w, facvol_w ) #endif ! Temperature zt(:,:,:) = tsn(:,:,:,jp_tem) ; zt_crs(:,:,:) = 0._wp CALL crs_dom_ope( zt, 'VOL', 'T', tmask, zt_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) tsn_crs(:,:,:,jp_tem) = zt_crs(:,:,:) CALL iom_put( "toce", tsn_crs(:,:,:,jp_tem) ) ! temp CALL iom_put( "sst" , tsn_crs(:,:,1,jp_tem) ) ! sst ! Salinity zs(:,:,:) = tsn(:,:,:,jp_sal) ; zs_crs(:,:,:) = 0._wp CALL crs_dom_ope( zs, 'VOL', 'T', tmask, zs_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) tsn_crs(:,:,:,jp_sal) = zs_crs(:,:,:) CALL iom_put( "soce" , tsn_crs(:,:,:,jp_sal) ) ! sal CALL iom_put( "sss" , tsn_crs(:,:,1,jp_sal) ) ! sss ! U-velocity CALL crs_dom_ope( un, 'SUM', 'U', umask, un_crs, p_e12=e2u, p_e3=zfse3u, p_surf_crs=e2e3u_msk, psgn=-1.0 ) CALL iom_put( "uoce" , un_crs ) ! i-current ! V-velocity CALL crs_dom_ope( vn, 'SUM', 'V', vmask, vn_crs, p_e12=e1v, p_e3=zfse3v, p_surf_crs=e1e3v_msk, psgn=-1.0 ) CALL iom_put( "voce" , vn_crs ) ! i-current !n2 CALL crs_dom_ope( rn2 , 'VOL', 'W', tmask, rn2_crs, p_e12=e1e2t, p_e3=zfse3t, psgn=1.0 ) ! Horizontal divergence ( following OPA_SRC/DYN/divcur.F90 ) hdivn_crs(:,:,:)=0._wp DO jk = 1, jpkm1 DO jj = 2,jpj_crs DO ji = 2,jpi_crs z2dcrsu = ( un_crs(ji ,jj ,jk) * e2e3u_msk(ji ,jj ,jk) ) & & - ( un_crs(ji-1,jj ,jk) * e2e3u_msk(ji-1,jj ,jk) ) z2dcrsv = ( vn_crs(ji ,jj ,jk) * e1e3v_msk(ji ,jj ,jk) ) & & - ( vn_crs(ji ,jj-1,jk) * e1e3v_msk(ji ,jj-1,jk) ) hdivn_crs(ji,jj,jk) = ( z2dcrsu + z2dcrsv ) ENDDO ENDDO ENDDO CALL crs_lbc_lnk( hdivn_crs, 'T', 1.0 ) ! CALL iom_put( "hdiv", hdivn_crs ) ! avt, avs SELECT CASE ( nn_crs_kz ) CASE ( 0 ) CALL crs_dom_ope( avt, 'VOL', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=zfse3w, psgn=1.0 ) CASE ( 1 ) CALL crs_dom_ope( avt, 'MAX', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=zfse3w, psgn=1.0 ) CASE ( 2 ) CALL crs_dom_ope( avt, 'MIN', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=zfse3w, psgn=1.0 ) CASE ( 3 ) CALL crs_dom_ope( avt, 'LOGVOL', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=zfse3w, p_mask_crs=tmask_crs, psgn=1.0 ) CASE ( 4 ) CALL crs_dom_ope( avt, 'MED', 'W', tmask, avt_crs, p_e12=e1e2t, p_e3=zfse3w, psgn=1.0 ) END SELECT ! CALL iom_put( "avt", avt_crs ) ! Kz !2D fields CALL crs_dom_ope( utau , 'SUM', 'U', umask, utau_crs , p_e12=e2u , p_surf_crs=e2u_crs , psgn=1.0 ) CALL crs_dom_ope( vtau , 'SUM', 'V', vmask, vtau_crs , p_e12=e1v , p_surf_crs=e1v_crs , psgn=1.0 ) CALL crs_dom_ope( wndm , 'SUM', 'T', tmask, wndm_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) CALL crs_dom_ope( rnf , 'MAX', 'T', tmask, rnf_crs , psgn=1.0 ) CALL crs_dom_ope( h_rnf, 'MAX', 'T', tmask, h_rnf_crs , psgn=1.0 ) z2d=REAL(nk_rnf,wp) CALL crs_dom_ope( z2d , 'MAX', 'T', tmask, z2d_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) nk_rnf_crs=INT(z2d_crs) CALL crs_dom_ope( qsr , 'SUM', 'T', tmask, qsr_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) CALL crs_dom_ope( emp ,'SUM', 'T', tmask, emp_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) CALL crs_dom_ope( fmmflx,'SUM', 'T', tmask, fmmflx_crs, p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) CALL crs_dom_ope( sfx ,'SUM', 'T', tmask, sfx_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) CALL crs_dom_ope( fr_i ,'SUM', 'T', tmask, fr_i_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) fr_i_crs=MAX( 0._wp, MIN( fr_i_crs , 1._wp ) ) z2d=REAL(nmln,wp) CALL crs_dom_ope( z2d , 'MAX', 'T', tmask, z2d_crs , p_e12=e1e2t, p_surf_crs=e1e2t_crs, psgn=1.0 ) nmln_crs=INT(z2d_crs) nmln_crs=MAX(nlb10,nmln_crs) CALL iom_put( "utau" , utau_crs ) ! i-tau output CALL iom_put( "vtau" , vtau_crs ) ! j-tau output CALL iom_put( "wspd" , wndm_crs ) ! wind speed output CALL iom_put( "runoffs" , rnf_crs ) ! runoff output CALL iom_put( "qsr" , qsr_crs ) ! qsr output CALL iom_put( "empmr" , emp_crs ) ! water flux output CALL iom_put( "saltflx" , sfx_crs ) ! salt flux output CALL iom_put( "ice_cover", fr_i_crs ) ! ice cover output zfse3t(:,:,:) = 1._wp CALL crs_dom_ope( sshn , 'VOL', 'T', tmask, sshn_crs , p_e12=e1e2t, p_e3=zfse3t , psgn=1.0 ) CALL iom_put( "ssh" , sshn_crs ) ! ssh output #if defined key_vvl !--------------------------------------------------------------------------------------------------- !variables au temps after !--------------------------------------------------------------------------------------------------- !ssha zfse3t(:,:,:) = 1._wp zt(:,:,:) = tmask(:,:,:) CALL crs_dom_ope( ssha , 'VOL', 'T', zt, ssha_crs , p_e12=e1e2t, p_e3=zfse3t , psgn=1.0 ) !vertical scale factors zfse3t(:,:,:) = e3t_a(:,:,:) zfse3u(:,:,:) = e3u_a(:,:,:) zfse3v(:,:,:) = e3v_a(:,:,:) CALL dom_vvl_interpol( zfse3t(:,:,:), zfse3w(:,:,:), 'W' ) CALL crs_dom_e3( e1t, e2t, zfse3t, p_sfc_3d_crs=e1e2w_crs, cd_type='T', p_mask=tmask, p_e3_crs=e3t_a_crs, p_e3_max_crs=zs_crs) CALL crs_dom_e3( e1t, e2t, zfse3w, p_sfc_3d_crs=e1e2w_crs, cd_type='W', p_mask=tmask, p_e3_crs=e3w_a_crs, p_e3_max_crs=zs_crs) CALL crs_dom_e3( e1u, e2u, zfse3u, p_sfc_2d_crs=e2u_crs , cd_type='U', p_mask=umask, p_e3_crs=e3u_a_crs, p_e3_max_crs=zs_crs) CALL crs_dom_e3( e1v, e2v, zfse3v, p_sfc_2d_crs=e1v_crs , cd_type='V', p_mask=vmask, p_e3_crs=e3v_a_crs, p_e3_max_crs=zs_crs) DO jk = 1, jpk DO ji = 1, jpi_crs DO jj = 1, jpj_crs IF( e3t_a_crs(ji,jj,jk) == 0._wp ) e3t_a_crs(ji,jj,jk) = e3t_1d(jk) IF( e3w_a_crs(ji,jj,jk) == 0._wp ) e3w_a_crs(ji,jj,jk) = e3w_1d(jk) IF( e3u_a_crs(ji,jj,jk) == 0._wp ) e3u_a_crs(ji,jj,jk) = e3t_1d(jk) IF( e3v_a_crs(ji,jj,jk) == 0._wp ) e3v_a_crs(ji,jj,jk) = e3t_1d(jk) ENDDO ENDDO ENDDO #endif #if defined key_vvl z1_2dt = 1._wp / ( 2. * rdt ) ! set time step size (Euler/Leapfrog) IF( neuler == 0 .AND. kt == nit000 ) z1_2dt = 1._wp / rdt wn_crs(:,:,jpk) = 0._wp DO jk = jpkm1, 1, -1 wn_crs(:,:,jk) = wn_crs(:,:,jk+1)*e1e2w_msk(:,:,jk+1) - ( hdivn_crs(:,:,jk) & & + z1_2dt * e1e2w_crs(:,:,jk) * ( e3t_a_crs(:,:,jk) - e3t_b_crs(:,:,jk) ) ) * tmask_crs(:,:,jk) WHERE( e1e2w_msk(:,:,jk) .NE. 0._wp ) wn_crs(:,:,jk) = wn_crs(:,:,jk) /e1e2w_msk(:,:,jk) ENDDO #else IF( ln_crs_wn ) THEN CALL crs_dom_ope( wn, 'SUM', 'W', tmask, wn_crs, p_e12=e1e2t, p_surf_crs=e1e2w_msk, psgn=1.0 ) ELSE wn_crs(:,:,jpk) = 0._wp DO jk = jpkm1, 1, -1 wn_crs(:,:,jk) = e1e2w_msk(:,:,jk+1)*wn_crs(:,:,jk+1) - hdivn_crs(:,:,jk) WHERE( e1e2w_msk(:,:,jk) .NE. 0._wp ) wn_crs(:,:,jk) = wn_crs(:,:,jk) /e1e2w_msk(:,:,jk) ENDDO ENDIF #endif CALL iom_put( "woce", wn_crs ) ! vertical velocity !--------------------------------------------------------------------------------------------------- ! free memory CALL wrk_dealloc( jpi, jpj, jpk, zfse3t, zfse3w ) CALL wrk_dealloc( jpi, jpj, jpk, zfse3u, zfse3v ) CALL wrk_dealloc( jpi, jpj, jpk, zt, zs, ztmp ) CALL wrk_dealloc( jpi, jpj, z2d ) CALL wrk_dealloc( jpi_crs, jpj_crs, jpk, zt_crs, zs_crs ) CALL wrk_dealloc( jpi_crs, jpj_crs, z2d_crs ) ! CALL iom_swap( "nemo" ) ! return back on high-resolution grid ! IF( nn_timing == 1 ) CALL timing_stop('crs_fld') ! END SUBROUTINE crs_fld !!====================================================================== END MODULE crsfld