Ignore:
Timestamp:
2018-06-19T21:54:00+02:00 (2 years ago)
Author:
antsia
Message:

delete iscplhsb, add iscpldiv and make the code readable

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_isf_divg_corr_GO6_package_r9385/NEMOGCM/NEMO/OPA_SRC/DOM/iscplrst.F90

    r9630 r9813  
    2424   USE lbclnk          ! communication 
    2525   USE iscplini        ! ice sheet coupling: initialisation 
    26    USE iscplhsb        ! ice sheet coupling: conservation 
     26   !USE iscplhsb        ! ice sheet coupling: conservation 
     27   USE iscpldiv 
    2728 
    2829   IMPLICIT NONE 
     
    8081      CALL iscpl_rst_interpol( ztmask_b, zumask_b, zvmask_b, zsmask_b, ze3t_b, ze3u_b, ze3v_b, zdepw_b ) 
    8182 
    82       !! compute correction if conservation needed 
    83       IF ( ln_hsb ) THEN 
    84          IF( iscpl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'rst_iscpl : unable to allocate rst_iscpl arrays' ) 
    85          CALL iscpl_cons(ztmask_b, zsmask_b, ze3t_b, htsc_iscpl, hdiv_iscpl, rdt_iscpl) 
    86       END IF 
    87           
    8883      !! print mesh/mask 
    8984      IF( nmsh /= 0 .AND. ln_iscpl )   CALL dom_wri      ! Create a domain file 
    90  
    91       IF ( ln_hsb ) THEN 
    92          cfile='correction' 
    93          cfile = TRIM( cfile ) 
    94          CALL iom_open  ( cfile, inum0, ldwrt = .TRUE., kiolib = jprstlib ) 
    95          CALL iom_rstput( 0, 0, inum0, 'vol_cor', hdiv_iscpl(:,:,:) ) 
    96          CALL iom_rstput( 0, 0, inum0, 'tem_cor', htsc_iscpl(:,:,:,jp_tem) ) 
    97          CALL iom_rstput( 0, 0, inum0, 'sal_cor', htsc_iscpl(:,:,:,jp_sal) ) 
    98          CALL iom_close ( inum0 ) 
    99       END IF 
    10085 
    10186      CALL wrk_dealloc(jpi,jpj,jpk, ztmask_b,zumask_b,zvmask_b )   
     
    150135      !! 
    151136      REAL(wp):: summsk, zsum, zsum1, zarea, zsumn, zsumb, zhdivdiff, zfp_ui, zfm_ui, zfp_vi, zfm_vi 
    152       REAL(wp):: zdz, zdzm1, zdzp1, zvol 
     137      REAL(wp):: zdz, zdzm1, zdzp1, zvol, zuflux_sum, zvflux_sum 
    153138      REAL(wp),DIMENSION(2)::ztrdtsb 
    154139      !! 
     
    156141      REAL(wp), DIMENSION(:,:    ), POINTER :: zbub   , zbvb   , zbun  , zbvn 
    157142      REAL(wp), DIMENSION(:,:    ), POINTER :: zssh0  , zssh1, zhu1, zhv1 
    158       REAL(wp), DIMENSION(:,:    ), POINTER :: zsmask0, zsmask1 
    159       REAL(wp), DIMENSION(:,:,:  ), POINTER :: ztmask0, ztmask1, ztrp, zuflux, zvflux 
     143      REAL(wp), DIMENSION(:,:    ), POINTER :: zsmask0, zsmask1, isfmask 
     144      REAL(wp), DIMENSION(:,:,:  ), POINTER :: ztmask0, ztmask1, ztrp, zuflux, zvflux, zhdiv, zhdiv2 
    160145      REAL(wp), DIMENSION(:,:,:  ), POINTER :: zwmaskn, zwmaskb, ztmp3d 
    161146      REAL(wp), DIMENSION(:,:,:,:), POINTER :: zts0 
     
    165150      CALL wrk_alloc(jpi,jpj,jpk,2, zts0                                   ) 
    166151      CALL wrk_alloc(jpi,jpj,jpk,   ztmask0, ztmask1 , ztrp, ztmp3d        )  
    167       CALL wrk_alloc(jpi,jpj,jpk,   zwmaskn, zwmaskb, zuflux, zvflux                       )  
    168       CALL wrk_alloc(jpi,jpj,       zsmask0, zsmask1                       )  
     152      CALL wrk_alloc(jpi,jpj,jpk,   zwmaskn, zwmaskb, zuflux, zvflux, zhdiv, zhdiv2               )  
     153      CALL wrk_alloc(jpi,jpj,       zsmask0, zsmask1, isfmask                       )  
    169154      CALL wrk_alloc(jpi,jpj,       zdmask , zdsmask, zvcorr, zucorr, zde3t)  
    170155      CALL wrk_alloc(jpi,jpj,       zbub   , zbvb    , zbun , zbvn         )  
     
    277262      END IF 
    278263 
    279 !============================================================================= 
    280 ! compute velocity 
    281  
    282       ub(:,:,:)=un(:,:,:) 
    283       vb(:,:,:)=vn(:,:,:) 
    284  
    285       un(:,:,:)=ub(:,:,:)*umask(:,:,:) 
    286       vn(:,:,:)=vb(:,:,:)*vmask(:,:,:) 
    287      
    288  
    289       !-----initialise the divergence and tracer correction------ 
    290       rhdivdiff(:,:,:) = 0._wp 
    291       rhdivdiff_trac(:,:,:,:) = 0._wp       
    292  
    293       !-----the change of flux in each velocity cell  
    294       zuflux(:,:,:) = (pe3u_b(:,:,:)*pumask_b(:,:,:) - fse3u_n(:,:,:)*umask(:,:,:))*ub(:,:,:) 
    295       zvflux(:,:,:) = (pe3v_b(:,:,:)*pvmask_b(:,:,:) - fse3v_n(:,:,:)*vmask(:,:,:))*vb(:,:,:) 
    296       
    297  
    298       ! compute divergence and tracer correction (whether the horizontal velocity cell is closed or open)     
    299       DO jj = 1,jpj 
    300          DO ji = 1,jpim1   ! jpim1 rather than jpi for the correct lbc_lnk synchronisation 
    301             !--------u-cell is closed-----------   
    302             DO jk  = 1,jpk 
    303                  !------the contribution to the west side of u-cell 
    304                  IF(tmask(ji,jj,jk) == 1) THEN   ! the west cell is still active  
    305                     zvol = e1e2t(ji,jj)*fse3t(ji,jj,jk) 
    306                     rhdivdiff(ji,jj,jk) = rhdivdiff(ji,jj,jk) + zuflux(ji,jj,jk)*e2u(ji,jj)/zvol  
    307                     rhdivdiff_trac(ji,jj,jk,1:2) =  rhdivdiff_trac(ji,jj,jk,1:2) - zuflux(ji,jj,jk)*tsb(ji,jj,jk,1:2)*e2u(ji,jj)/zvol 
    308                  ELSEIF(mikt(ji,jj) > 1) THEN  ! the west cell becomes dry, put the waterflux to the top of the west column 
    309                     jtop=mikt(ji,jj) 
    310                     zvol=e1e2t(ji,jj)*fse3t(ji,jj,jtop) 
    311                     rhdivdiff(ji,jj,jtop) = rhdivdiff(ji,jj,jtop) + zuflux(ji,jj,jk)*e2u(ji,jj)/zvol 
    312                     rhdivdiff_trac(ji,jj,jtop,1:2) = rhdivdiff_trac(ji,jj,jtop,1:2)-zuflux(ji,jj,jk)*tsb(ji,jj,jtop,1:2)*e2u(ji,jj)/zvol 
    313                  !ELSE ! for the closing T-column : nothing to be done  
    314                  ENDIF 
    315  
    316                  !------the contribution to the east side of u-cell 
    317                  IF(tmask(ji+1,jj,jk) == 1) THEN  ! the east cell is still active 
    318                     zvol = e1e2t(ji+1,jj)*fse3t(ji+1,jj,jk) 
    319                     rhdivdiff(ji+1,jj,jk) = rhdivdiff(ji+1,jj,jk) - zuflux(ji,jj,jk)*e2u(ji,jj)/zvol 
    320                     rhdivdiff_trac(ji+1,jj,jk,1:2) = rhdivdiff_trac(ji+1,jj,jk,1:2) + zuflux(ji,jj,jk)*tsb(ji+1,jj,jk,1:2)*e2u(ji,jj)/zvol 
    321                  ELSEIF(mikt(ji+1,jj) > 1) THEN    ! the east cell becomes dry, put the waterflux to the top of the east column 
    322                     jtop = mikt(ji+1,jj) 
    323                     zvol = e1e2t(ji+1,jj)*fse3t(ji+1,jj,jtop) 
    324                     rhdivdiff(ji+1,jj,jtop) = rhdivdiff(ji+1,jj,jtop) - zuflux(ji,jj,jk)*e2u(ji,jj)/zvol 
    325                     rhdivdiff_trac(ji+1,jj,jtop,1:2) = rhdivdiff_trac(ji+1,jj,jtop,1:2) + zuflux(ji,jj,jk)*tsb(ji+1,jj,jtop,1:2)*e2u(ji,jj)/zvol 
    326                  !ELSE  ! for the closing T-column : nothing to be done 
    327                  ENDIF 
    328             ENDDO 
    329          END DO 
    330       END DO 
    331  
    332  
    333       DO jj = 1,jpjm1  ! jpim1 rather than jpi for the correct lbc_lnk synchronisation 
    334          DO ji = 1,jpi    
    335             !-----------v-cell is closed------------- 
    336             DO jk  = 1,jpk 
    337                  !------the contribution to the south side of v-cell 
    338                  IF(tmask(ji,jj,jk) == 1) THEN  ! the south cell is still active 
    339                    zvol = e1e2t(ji,jj)*fse3t(ji,jj,jk) 
    340                    rhdivdiff(ji,jj,jk) = rhdivdiff(ji,jj,jk) + zvflux(ji,jj,jk)*e1v(ji,jj)/zvol 
    341                    rhdivdiff_trac(ji,jj,jk,1:2) =  rhdivdiff_trac(ji,jj,jk,1:2) - zvflux(ji,jj,jk)*tsb(ji,jj,jk,1:2)*e1v(ji,jj)/zvol 
    342                  ELSEIF(mikt(ji,jj) > 1) THEN  !the south cell becomes dry, put the waterflux to the top of the south column 
    343                    jtop = mikt(ji,jj) 
    344                    zvol = e1e2t(ji,jj)*fse3t(ji,jj,jtop) 
    345                    rhdivdiff(ji,jj,jtop) = rhdivdiff(ji,jj,jtop) + zvflux(ji,jj,jk)*e1v(ji,jj)/zvol 
    346                    rhdivdiff_trac(ji,jj,jtop,1:2) = rhdivdiff_trac(ji,jj,jtop,1:2) - zvflux(ji,jj,jk)*tsb(ji,jj,jtop,1:2)*e1v(ji,jj)/zvol 
    347                  ENDIF 
    348  
    349                  !------the contribution to the north side of v-cell 
    350                  IF(tmask(ji,jj+1,jk) == 1) THEN ! the north cell is still active 
    351                     zvol = e1e2t(ji,jj+1)*fse3t(ji,jj+1,jk) 
    352                     rhdivdiff(ji,jj+1,jk) = rhdivdiff(ji,jj+1,jk) - zvflux(ji,jj,jk)*e1v(ji,jj)/zvol 
    353                     rhdivdiff_trac(ji,jj+1,jk,1:2) = rhdivdiff_trac(ji,jj+1,jk,1:2) + zvflux(ji,jj,jk)*tsb(ji,jj+1,jk,1:2)*e1v(ji,jj)/zvol 
    354                  ELSEIF(mikt(ji,jj+1) > 1) THEN  ! the north cell becomes dry, put the waterflux to the top of the north column 
    355                     jtop = mikt(ji,jj+1) 
    356                     zvol = e1e2t(ji,jj+1)*fse3t(ji,jj+1,jtop) 
    357                     rhdivdiff(ji,jj+1,jtop) = rhdivdiff(ji,jj+1,jtop) - zvflux(ji,jj,jk)*e1v(ji,jj)/zvol 
    358                     rhdivdiff_trac(ji,jj+1,jtop,1:2) = rhdivdiff_trac(ji,jj+1,jtop,1:2) + zvflux(ji,jj,jk)*tsb(ji,jj+1,jtop,1:2)*e1v(ji,jj)/zvol 
    359                  ENDIF 
    360             ENDDO 
    361  
    362          END DO 
    363       END DO  
    364        
    365  
    366       CALL lbc_lnk(rhdivdiff,   'T',1._wp) 
    367       CALL lbc_lnk(rhdivdiff_trac(:,:,:,1),   'T',1._wp) 
    368       CALL lbc_lnk(rhdivdiff_trac(:,:,:,2),   'T',1._wp) 
    369  
     264      !----compute and save the divergence difference (before and after remeshing)  
     265      CALL iscpl_div_corr(ptmask_b, pe3t_b, pe3u_b, pe3v_b) 
    370266 
    371267      ! compute temp and salt 
     
    473369      CALL wrk_dealloc(jpi,jpj,jpk,2, zts0                                   ) 
    474370      CALL wrk_dealloc(jpi,jpj,jpk,   ztmask0, ztmask1 , ztrp                )  
    475       CALL wrk_dealloc(jpi,jpj,jpk,   zwmaskn, zwmaskb , ztmp3d, zuflux, zvflux    )  
    476       CALL wrk_dealloc(jpi,jpj,       zsmask0, zsmask1                       )  
     371      CALL wrk_dealloc(jpi,jpj,jpk,   zwmaskn, zwmaskb , ztmp3d, zuflux, zvflux, zhdiv, zhdiv2    )  
     372      CALL wrk_dealloc(jpi,jpj,       zsmask0, zsmask1, isfmask                     )  
    477373      CALL wrk_dealloc(jpi,jpj,       zdmask , zdsmask, zvcorr, zucorr, zde3t)  
    478374      CALL wrk_dealloc(jpi,jpj,       zbub   , zbvb    , zbun  , zbvn        )  
Note: See TracChangeset for help on using the changeset viewer.