- Timestamp:
- 2020-07-30T13:20:57+02:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/tickets_icb_1900/src/OCE/ICB/icbutl.F90
r13265 r13357 19 19 !!---------------------------------------------------------------------- 20 20 USE par_oce ! ocean parameters 21 USE oce, ONLY: ts, uu, vv 21 22 USE dom_oce ! ocean domain 22 23 USE in_out_manager ! IO parameters … … 34 35 PRIVATE 35 36 37 INTERFACE icb_utl_bilin_h 38 MODULE PROCEDURE icb_utl_bilin_2d_h, icb_utl_bilin_3d_h 39 END INTERFACE 40 36 41 PUBLIC icb_utl_copy ! routine called in icbstp module 42 PUBLIC icb_utl_getkb ! routine called in icbdyn and icbthm modules 37 43 PUBLIC icb_utl_interp ! routine called in icbdyn, icbthm modules 38 44 PUBLIC icb_utl_bilin_h ! routine called in icbdyn module … … 56 62 CONTAINS 57 63 58 SUBROUTINE icb_utl_copy( )64 SUBROUTINE icb_utl_copy( Kmm ) 59 65 !!---------------------------------------------------------------------- 60 66 !! *** ROUTINE icb_utl_copy *** … … 64 70 !! ** Method : - blah blah 65 71 !!---------------------------------------------------------------------- 72 REAL(wp), DIMENSION(0:jpi+1,0:jpj+1) :: ztmp 66 73 #if defined key_si3 67 74 REAL(wp), DIMENSION(jpi,jpj) :: zssh_lead_m ! ocean surface (ssh_m) if ice is not embedded 68 75 ! ! ocean surface in leads if ice is embedded 69 76 #endif 77 INTEGER :: jk ! vertical loop index 78 INTEGER :: Kmm ! ocean time levelindex 79 ! 70 80 ! copy nemo forcing arrays into iceberg versions with extra halo 71 81 ! only necessary for variables not on T points 72 82 ! and ssh which is used to calculate gradients 73 74 uo_e(1:jpi,1:jpj) = ssu_m(:,:) * umask(:,:,1) 75 vo_e(1:jpi,1:jpj) = ssv_m(:,:) * vmask(:,:,1) 76 tt_e(1:jpi,1:jpj) = sst_m(:,:) 77 fr_e(1:jpi,1:jpj) = fr_i (:,:) 78 ua_e(1:jpi,1:jpj) = utau (:,:) * umask(:,:,1) ! maybe mask useless because mask applied in sbcblk 79 va_e(1:jpi,1:jpj) = vtau (:,:) * vmask(:,:,1) ! maybe mask useless because mask applied in sbcblk 80 ! 81 CALL lbc_lnk_icb( 'icbutl', uo_e, 'U', -1._wp, 1, 1 ) 82 CALL lbc_lnk_icb( 'icbutl', vo_e, 'V', -1._wp, 1, 1 ) 83 ! 84 ! surface forcing 85 ! 86 ssu_e(1:jpi,1:jpj) = ssu_m(:,:) * umask(:,:,1) 87 ssv_e(1:jpi,1:jpj) = ssv_m(:,:) * vmask(:,:,1) 88 sst_e(1:jpi,1:jpj) = sst_m(:,:) 89 fr_e (1:jpi,1:jpj) = fr_i (:,:) 90 ua_e (1:jpi,1:jpj) = utau (:,:) * umask(:,:,1) ! maybe mask useless because mask applied in sbcblk 91 va_e (1:jpi,1:jpj) = vtau (:,:) * vmask(:,:,1) ! maybe mask useless because mask applied in sbcblk 92 ! 93 CALL lbc_lnk_icb( 'icbutl', ssu_e, 'U', -1._wp, 1, 1 ) 94 CALL lbc_lnk_icb( 'icbutl', ssv_e, 'V', -1._wp, 1, 1 ) 83 95 CALL lbc_lnk_icb( 'icbutl', ua_e, 'U', -1._wp, 1, 1 ) 84 96 CALL lbc_lnk_icb( 'icbutl', va_e, 'V', -1._wp, 1, 1 ) … … 98 110 #endif 99 111 ! 112 ! (PM) could be improve with a 3d lbclnk gathering both variables 113 ! should be done once extra haloe generalised 114 IF ( ln_M2016 ) THEN 115 DO jk = 1,jpk 116 ! uoce 117 ztmp(1:jpi,1:jpj) = uu(:,:,jk,Kmm) 118 CALL lbc_lnk_icb( 'icbutl', ztmp, 'U', -1._wp, 1, 1 ) 119 uoce_e(:,:,jk) = ztmp(:,:) 120 ! 121 ! voce 122 ztmp(1:jpi,1:jpj) = vv(:,:,jk,Kmm) 123 CALL lbc_lnk_icb( 'icbutl', ztmp, 'V', -1._wp, 1, 1 ) 124 voce_e(:,:,jk) = ztmp(:,:) 125 END DO 126 toce_e(1:jpi,1:jpj,1:jpk) = ts(:,:,:,1,Kmm) 127 e3t_e (1:jpi,1:jpj,1:jpk) = e3t(:,:,:,Kmm) 128 END IF 129 ! 100 130 END SUBROUTINE icb_utl_copy 101 131 102 132 103 SUBROUTINE icb_utl_interp( pi, pj, pe1, puo, pui, pua, pssh_i, & 104 & pe2, pvo, pvi, pva, pssh_j, & 105 & psst, pcn, phi, pff, plon, plat) 133 SUBROUTINE icb_utl_interp( pi, pj, pe1, pssu, pui, pua, pssh_i, & 134 & pe2, pssv, pvi, pva, pssh_j, & 135 & psst, pcn, phi, pff, plon, plat, & 136 & ptoce, puoce, pvoce, pe3t ) 106 137 !!---------------------------------------------------------------------- 107 138 !! *** ROUTINE icb_utl_interp *** … … 122 153 REAL(wp), INTENT(in ) :: pi , pj ! position in (i,j) referential 123 154 REAL(wp), INTENT( out), OPTIONAL :: pe1, pe2 ! i- and j scale factors 124 REAL(wp), INTENT( out), OPTIONAL :: p uo, pvo, pui, pvi, pua, pva ! ocean, ice and wind speeds155 REAL(wp), INTENT( out), OPTIONAL :: pssu, pssv, pui, pvi, pua, pva ! ocean, ice and wind speeds 125 156 REAL(wp), INTENT( out), OPTIONAL :: pssh_i, pssh_j ! ssh i- & j-gradients 126 157 REAL(wp), INTENT( out), OPTIONAL :: psst, pcn, phi, pff ! SST, ice concentration, ice thickness, Coriolis 127 158 REAL(wp), INTENT( out), OPTIONAL :: plat, plon ! position 159 REAL(wp), DIMENSION(jpk), INTENT( out), OPTIONAL :: ptoce, puoce, pvoce, pe3t ! 3D variables 128 160 ! 129 161 REAL(wp), DIMENSION(4) :: zwT , zwU , zwV , zwF ! interpolation weight … … 147 179 IF ( PRESENT(plat) ) plat= icb_utl_bilin_h( rlat_e, iiT, ijT, zwT, .false. ) 148 180 ! 149 IF ( PRESENT(p uo ) ) puo = icb_utl_bilin_h( uo_e, iiU, ijU, zwU , .false. ) ! ocean velocities150 IF ( PRESENT(p vo ) ) pvo = icb_utl_bilin_h( vo_e, iiV, ijV, zwV , .false. ) !151 IF ( PRESENT(psst) ) psst = icb_utl_bilin_h( tt_e, iiT, ijT, zwT * zmskT, .false. ) ! sst152 IF ( PRESENT(pcn ) ) pcn = icb_utl_bilin_h( fr_e , iiT, ijT, zwT * zmskT, .false. ) ! ice concentration153 IF ( PRESENT(pff ) ) pff = icb_utl_bilin_h( ff_e , iiF, ijF, zwF , .false. ) ! Coriolis parameter181 IF ( PRESENT(pssu) ) pssu = icb_utl_bilin_h( ssu_e, iiU, ijU, zwU , .false. ) ! ocean velocities 182 IF ( PRESENT(pssv) ) pssv = icb_utl_bilin_h( ssv_e, iiV, ijV, zwV , .false. ) ! 183 IF ( PRESENT(psst) ) psst = icb_utl_bilin_h( sst_e, iiT, ijT, zwT * zmskT, .false. ) ! sst 184 IF ( PRESENT(pcn ) ) pcn = icb_utl_bilin_h( fr_e , iiT, ijT, zwT * zmskT, .false. ) ! ice concentration 185 IF ( PRESENT(pff ) ) pff = icb_utl_bilin_h( ff_e , iiF, ijF, zwF , .false. ) ! Coriolis parameter 154 186 ! 155 187 IF ( PRESENT(pua) .AND. PRESENT(pva) ) THEN … … 188 220 & icb_utl_bilin_h( ssh_e, iiTm, ijTm, zwTm*zmskTm, .false. ) ) / ( 0.2_wp * pe2 ) 189 221 END IF 222 ! 223 ! 3d interpolation 224 IF ( PRESENT(puoce) .AND. PRESENT(pvoce) ) THEN 225 puoce(:) = icb_utl_bilin_h( uoce_e , iiU, ijU, zwU ) 226 pvoce(:) = icb_utl_bilin_h( voce_e , iiV, ijV, zwV ) 227 END IF 228 IF ( PRESENT(ptoce) ) ptoce(:) = icb_utl_bilin_h( toce_e , iiT, ijT, zwT * zmskT ) 229 IF ( PRESENT(pe3t) ) pe3t(:) = e3t_e(iiT,ijT,:) ! as in Nacho tarball need to be fix once we are able to reproduce Nacho results 190 230 ! 191 231 END SUBROUTINE icb_utl_interp … … 291 331 END SUBROUTINE icb_utl_pos 292 332 293 REAL(wp) FUNCTION icb_utl_bilin_ h( pfld, pii, pij, pw, pllon )333 REAL(wp) FUNCTION icb_utl_bilin_2d_h( pfld, pii, pij, pw, pllon ) 294 334 !!---------------------------------------------------------------------- 295 335 !! *** FUNCTION icb_utl_bilin *** … … 321 361 ! 322 362 ! compute interpolated value 323 icb_utl_bilin_h = ( zdat(1)*pw(1) + zdat(2)*pw(2) + zdat(3)*pw(3) + zdat(4)*pw(4) ) / MAX(1.e-20, pw(1)+pw(2)+pw(3)+pw(4)) 324 ! 325 IF( pllon .AND. icb_utl_bilin_h > 180._wp ) icb_utl_bilin_h = icb_utl_bilin_h - 360._wp 326 ! 327 END FUNCTION icb_utl_bilin_h 363 icb_utl_bilin_2d_h = ( zdat(1)*pw(1) + zdat(2)*pw(2) + zdat(3)*pw(3) + zdat(4)*pw(4) ) / MAX(1.e-20, pw(1)+pw(2)+pw(3)+pw(4)) 364 ! 365 IF( pllon .AND. icb_utl_bilin_2d_h > 180._wp ) icb_utl_bilin_2d_h = icb_utl_bilin_2d_h - 360._wp 366 ! 367 END FUNCTION icb_utl_bilin_2d_h 368 369 FUNCTION icb_utl_bilin_3d_h( pfld, pii, pij, pw ) 370 !!---------------------------------------------------------------------- 371 !! *** FUNCTION icb_utl_bilin *** 372 !! 373 !! ** Purpose : bilinear interpolation at berg location depending on the grid-point type 374 !! this version deals with extra halo points 375 !! 376 !! !!gm CAUTION an optional argument should be added to handle 377 !! the slip/no-slip conditions ==>>> to be done later 378 !! 379 !!---------------------------------------------------------------------- 380 REAL(wp), DIMENSION(0:jpi+1,0:jpj+1, jpk), INTENT(in) :: pfld ! field to be interpolated 381 REAL(wp), DIMENSION(4) , INTENT(in) :: pw ! weight 382 INTEGER , INTENT(in) :: pii, pij ! bottom left corner 383 REAL(wp), DIMENSION(jpk) :: icb_utl_bilin_3d_h 384 ! 385 REAL(wp), DIMENSION(4,jpk) :: zdat ! input data 386 !!---------------------------------------------------------------------- 387 ! 388 ! data 389 zdat(1,:) = pfld(pii ,pij ,:) 390 zdat(2,:) = pfld(pii+1,pij ,:) 391 zdat(3,:) = pfld(pii ,pij+1,:) 392 zdat(4,:) = pfld(pii+1,pij+1,:) 393 ! 394 ! compute interpolated value 395 icb_utl_bilin_3d_h(:) = ( zdat(1,:)*pw(1) + zdat(2,:)*pw(2) + zdat(3,:)*pw(3) + zdat(4,:)*pw(4) ) / MAX(1.e-20, pw(1)+pw(2)+pw(3)+pw(4)) 396 ! 397 END FUNCTION icb_utl_bilin_3d_h 328 398 329 399 REAL(wp) FUNCTION icb_utl_bilin_e( pet, peu, pev, pef, pi, pj ) … … 405 475 END FUNCTION icb_utl_bilin_e 406 476 477 SUBROUTINE icb_utl_getkb( kb, pe3, pD ) 478 !!---------------------------------------------------------------------- 479 !! *** ROUTINE icb_utl_getkb *** 480 !! 481 !! ** Purpose : compute the latest level affected by icb 482 !! 483 !!---------------------------------------------------------------------- 484 INTEGER, INTENT(out) :: kb 485 REAL(wp), DIMENSION(:), INTENT(in) :: pe3 486 REAL(wp), INTENT(in) :: pD 487 !! 488 INTEGER :: jk 489 REAL(wp) :: zdepw 490 !! 491 zdepw = 0.0 492 kb = 1 493 DO WHILE ( zdepw < pD) 494 zdepw = zdepw + pe3(kb) 495 kb = kb + 1 496 END DO 497 kb = kb - 1 498 END SUBROUTINE 407 499 408 500 SUBROUTINE icb_utl_add( bergvals, ptvals ) … … 593 685 WRITE(numicb, 9200) kt, berg%number(1), & 594 686 pt%xi, pt%yj, pt%lon, pt%lat, pt%uvel, pt%vvel, & 595 pt% uo, pt%vo, pt%ua, pt%va, pt%ui, pt%vi687 pt%ssu, pt%ssv, pt%ua, pt%va, pt%ui, pt%vi 596 688 CALL flush( numicb ) 597 689 9200 FORMAT(5x,i5,2x,i10,6(2x,2f10.4)) … … 619 711 WRITE(numicb,'(a," pe=(",i3,")")' ) cd_label, narea 620 712 WRITE(numicb,'(a8,4x,a6,12x,a5,15x,a7,19x,a3,17x,a5,17x,a5,17x,a5)' ) & 621 & 'timestep', 'number', 'xi,yj','lon,lat','u,v',' uo,vo','ua,va','ui,vi'713 & 'timestep', 'number', 'xi,yj','lon,lat','u,v','ssu,ssv','ua,va','ui,vi' 622 714 ENDIF 623 715 DO WHILE( ASSOCIATED(this) )
Note: See TracChangeset
for help on using the changeset viewer.
