- Timestamp:
- 2019-07-17T15:29:15+02:00 (5 years ago)
- Location:
- branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC
- Files:
-
- 13 edited
-
cpl_oasis3.F90 (modified) (21 diffs)
-
geo2ocean.F90 (modified) (2 diffs)
-
sbc_oce.F90 (modified) (5 diffs)
-
sbcapr.F90 (modified) (5 diffs)
-
sbcblk_core.F90 (modified) (3 diffs)
-
sbcblk_mfs.F90 (modified) (1 diff)
-
sbccpl.F90 (modified) (45 diffs)
-
sbcflx.F90 (modified) (5 diffs)
-
sbcice_cice.F90 (modified) (2 diffs)
-
sbcmod.F90 (modified) (9 diffs)
-
sbcrnf.F90 (modified) (3 diffs)
-
sbcssm.F90 (modified) (9 diffs)
-
sbcwave.F90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/cpl_oasis3.F90
r8058 r11277 24 24 !! cpl_finalize : finalize the coupled mode communication 25 25 !!---------------------------------------------------------------------- 26 #if defined key_oasis3 26 #if defined key_oasis3 || defined key_oasis3mct 27 27 USE mod_oasis ! OASIS3-MCT module 28 28 #endif … … 32 32 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 33 33 34 #if defined key_cpl_rootexchg 35 USE lib_mpp, only : mppsync 36 USE lib_mpp, only : mppscatter,mppgather 37 #endif 38 34 39 IMPLICIT NONE 35 40 PRIVATE … … 41 46 PUBLIC cpl_freq 42 47 PUBLIC cpl_finalize 48 #if defined key_mpp_mpi 49 INCLUDE 'mpif.h' 50 #endif 51 52 INTEGER, PARAMETER :: localRoot = 0 53 LOGICAL :: commRank ! true for ranks doing OASIS communication 54 #if defined key_cpl_rootexchg 55 LOGICAL :: rootexchg =.true. ! logical switch 56 #else 57 LOGICAL :: rootexchg =.false. ! logical switch 58 #endif 43 59 44 60 INTEGER, PUBLIC :: OASIS_Rcv = 1 !: return code if received field … … 46 62 INTEGER :: ncomp_id ! id returned by oasis_init_comp 47 63 INTEGER :: nerror ! return error code 48 #if ! defined key_oasis3 64 #if ! defined key_oasis3 && ! defined key_oasis3mct 49 65 ! OASIS Variables not used. defined only for compilation purpose 50 66 INTEGER :: OASIS_Out = -1 … … 65 81 INTEGER :: nsnd ! total number of fields sent 66 82 INTEGER :: ncplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 67 INTEGER, PUBLIC, PARAMETER :: nmaxfld=5 0! Maximum number of coupling fields83 INTEGER, PUBLIC, PARAMETER :: nmaxfld=55 ! Maximum number of coupling fields 68 84 INTEGER, PUBLIC, PARAMETER :: nmaxcat=5 ! Maximum number of coupling fields 69 85 INTEGER, PUBLIC, PARAMETER :: nmaxcpl=5 ! Maximum number of coupling fields … … 82 98 83 99 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: exfld ! Temporary buffer for receiving 100 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: tbuf ! Temporary buffer for sending / receiving 101 INTEGER, PUBLIC :: localComm 84 102 85 103 !!---------------------------------------------------------------------- … … 120 138 IF ( nerror /= OASIS_Ok ) & 121 139 CALL oasis_abort (ncomp_id, 'cpl_init','Failure in oasis_get_localcomm' ) 140 localComm = kl_comm 122 141 ! 123 142 END SUBROUTINE cpl_init … … 149 168 IF(lwp) WRITE(numout,*) 150 169 170 commRank = .false. 171 IF ( rootexchg ) THEN 172 IF ( nproc == localRoot ) commRank = .true. 173 ELSE 174 commRank = .true. 175 ENDIF 176 151 177 ncplmodel = kcplmodel 152 178 IF( kcplmodel > nmaxcpl ) THEN … … 172 198 ishape(:,2) = (/ 1, nlej-nldj+1 /) 173 199 ! 174 ! ... Allocate memory for data exchange175 !176 ALLOCATE(exfld(nlei-nldi+1, nlej-nldj+1), stat = nerror)177 IF( nerror > 0 ) THEN178 CALL oasis_abort ( ncomp_id, 'cpl_define', 'Failure in allocating exfld') ; RETURN179 ENDIF180 200 ! 181 201 ! ----------------------------------------------------------------- … … 183 203 ! ----------------------------------------------------------------- 184 204 205 IF ( rootexchg ) THEN 206 207 paral(1) = 2 ! box partitioning 208 paral(2) = 0 ! NEMO lower left corner global offset 209 paral(3) = jpiglo ! local extent in i 210 paral(4) = jpjglo ! local extent in j 211 paral(5) = jpiglo ! global extent in x 212 213 ELSE 185 214 paral(1) = 2 ! box partitioning 186 215 paral(2) = jpiglo * (nldj-1+njmpp-1) + (nldi-1+nimpp-1) ! NEMO lower left corner global offset … … 196 225 ENDIF 197 226 198 CALL oasis_def_partition ( id_part, paral, nerror ) 227 ENDIF 228 IF ( commRank ) CALL oasis_def_partition ( id_part, paral, nerror ) 229 230 ! ... Allocate memory for data exchange 231 ! 232 ALLOCATE(exfld(paral(3), paral(4)), stat = nerror) 233 IF( nerror > 0 ) THEN 234 CALL oasis_abort ( ncomp_id, 'cpl_define', 'Failure in allocating exfld') ; RETURN 235 ENDIF 236 IF ( rootexchg ) THEN 237 ! Should possibly use one of the work arrays for tbuf really 238 ALLOCATE(tbuf(jpi, jpj, jpnij), stat = nerror) 239 IF( nerror > 0 ) THEN 240 CALL oasis_abort ( ncomp_id, 'cpl_define', 'Failure in allocating tbuf') ; RETURN 241 ENDIF 242 ENDIF 243 ! 244 IF (commRank ) THEN 199 245 ! 200 246 ! ... Announce send variables. … … 288 334 ENDIF 289 335 END DO 336 ! 337 ENDIF ! commRank=true 290 338 291 339 !------------------------------------------------------------------ … … 293 341 !------------------------------------------------------------------ 294 342 295 CALL oasis_enddef(nerror) 296 IF( nerror /= OASIS_Ok ) CALL oasis_abort ( ncomp_id, 'cpl_define', 'Failure in oasis_enddef') 343 IF ( commRank ) THEN 344 CALL oasis_enddef(nerror) 345 IF( nerror /= OASIS_Ok ) CALL oasis_abort ( ncomp_id, 'cpl_define', 'Failure in oasis_enddef') 346 ENDIF 297 347 ! 298 348 END SUBROUTINE cpl_define … … 311 361 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pdata 312 362 !! 313 INTEGER :: j c,jm! local loop index363 INTEGER :: jn,jc,jm ! local loop index 314 364 !!-------------------------------------------------------------------- 315 365 ! … … 320 370 321 371 IF( ssnd(kid)%nid(jc,jm) /= -1 ) THEN 322 CALL oasis_put ( ssnd(kid)%nid(jc,jm), kstep, pdata(nldi:nlei, nldj:nlej,jc), kinfo ) 372 IF ( rootexchg ) THEN 373 ! 374 ! collect data on the local root process 375 ! 376 CALL mppgather (pdata(:,:,jc),localRoot,tbuf) 377 CALL mppsync 378 379 IF ( nproc == localRoot ) THEN 380 DO jn = 1, jpnij 381 exfld(nimppt(jn)-1+nldit(jn):nimppt(jn)+nleit(jn)-1,njmppt(jn)-1+nldjt(jn):njmppt(jn)+nlejt(jn)-1)= & 382 tbuf(nldit(jn):nleit(jn),nldjt(jn):nlejt(jn),jn) 383 ENDDO 384 ! snd data to OASIS3 385 CALL oasis_put ( ssnd(kid)%nid(jc,jm), kstep, exfld, kinfo ) 386 ENDIF 387 ELSE 388 ! snd data to OASIS3 389 CALL oasis_put ( ssnd(kid)%nid(jc,jm), kstep, pdata(nldi:nlei, nldj:nlej,jc), kinfo ) 390 ENDIF 323 391 324 392 IF ( ln_ctl ) THEN … … 358 426 INTEGER , INTENT( out) :: kinfo ! OASIS3 info argument 359 427 !! 360 INTEGER :: j c,jm! local loop index428 INTEGER :: jn,jc,jm ! local loop index 361 429 LOGICAL :: llaction, llfisrt 362 430 !!-------------------------------------------------------------------- … … 372 440 373 441 IF( srcv(kid)%nid(jc,jm) /= -1 ) THEN 374 375 CALL oasis_get ( srcv(kid)%nid(jc,jm), kstep, exfld, kinfo ) 442 ! 443 ! receive data from OASIS3 444 ! 445 IF ( commRank ) CALL oasis_get ( srcv(kid)%nid(jc,jm), kstep, exfld, kinfo ) 446 IF ( rootexchg ) CALL MPI_BCAST ( kinfo, 1, MPI_INTEGER, localRoot, localComm, nerror ) 376 447 377 448 llaction = kinfo == OASIS_Recvd .OR. kinfo == OASIS_FromRest .OR. & … … 384 455 kinfo = OASIS_Rcv 385 456 IF( llfisrt ) THEN 386 pdata(nldi:nlei,nldj:nlej,jc) = exfld(:,:) * pmask(nldi:nlei,nldj:nlej,jm) 457 IF ( rootexchg ) THEN 458 ! distribute data to processes 459 ! 460 IF ( nproc == localRoot ) THEN 461 DO jn = 1, jpnij 462 tbuf(nldit(jn):nleit(jn),nldjt(jn):nlejt(jn),jn)= & 463 exfld(nimppt(jn)-1+nldit(jn):nimppt(jn)+nleit(jn)-1,njmppt(jn)-1+nldjt(jn):njmppt(jn)+nlejt(jn)-1) 464 ! NOTE: we are missing combining this with pmask (see else below) 465 ENDDO 466 ENDIF 467 CALL mppscatter(tbuf,localRoot,pdata(:,:,jc)) 468 CALL mppsync 469 ELSE 470 pdata(nldi:nlei, nldj:nlej, jc) = exfld(:,:) * pmask(nldi:nlei,nldj:nlej,jm) 471 ENDIF 387 472 llfisrt = .FALSE. 388 473 ELSE … … 462 547 #if defined key_oa3mct_v3 463 548 CALL oasis_get_freqs(id, mop, 1, itmp, info) 464 #else 549 #endif 550 #if defined key_oasis3 465 551 CALL oasis_get_freqs(id, 1, itmp, info) 466 552 #endif 467 553 cpl_freq = itmp(1) 554 #if defined key_oasis3mct 555 cpl_freq = namflddti( id ) 556 #endif 468 557 ENDIF 469 558 ! … … 481 570 ! 482 571 DEALLOCATE( exfld ) 572 IF ( rootexchg ) DEALLOCATE ( tbuf ) 483 573 IF (nstop == 0) THEN 484 574 CALL oasis_terminate( nerror ) … … 489 579 END SUBROUTINE cpl_finalize 490 580 491 #if ! defined key_oasis3 581 #if ! defined key_oasis3 && ! defined key_oasis3mct 492 582 493 583 !!---------------------------------------------------------------------- -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/geo2ocean.F90
r8058 r11277 51 51 52 52 SUBROUTINE repcmo ( pxu1, pyu1, pxv1, pyv1, & 53 px2 , py2 )53 px2 , py2, kchoix ) 54 54 !!---------------------------------------------------------------------- 55 55 !! *** ROUTINE repcmo *** … … 67 67 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: px2 ! i-componante (defined at u-point) 68 68 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: py2 ! j-componante (defined at v-point) 69 !!---------------------------------------------------------------------- 70 71 ! Change from geographic to stretched coordinate 72 ! ---------------------------------------------- 73 CALL rot_rep( pxu1, pyu1, 'U', 'en->i',px2 ) 74 CALL rot_rep( pxv1, pyv1, 'V', 'en->j',py2 ) 69 INTEGER, INTENT( IN ) :: kchoix ! type of transformation 70 ! = 1 change from geographic to model grid. 71 ! =-1 change from model to geographic grid 72 !!---------------------------------------------------------------------- 73 74 SELECT CASE (kchoix) 75 CASE ( 1) 76 ! Change from geographic to stretched coordinate 77 ! ---------------------------------------------- 78 79 CALL rot_rep( pxu1, pyu1, 'U', 'en->i',px2 ) 80 CALL rot_rep( pxv1, pyv1, 'V', 'en->j',py2 ) 81 CASE (-1) 82 ! Change from stretched to geographic coordinate 83 ! ---------------------------------------------- 84 85 CALL rot_rep( pxu1, pyu1, 'U', 'ij->e',px2 ) 86 CALL rot_rep( pxv1, pyv1, 'V', 'ij->n',py2 ) 87 END SELECT 75 88 76 89 END SUBROUTINE repcmo -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90
r8059 r11277 35 35 LOGICAL , PUBLIC :: ln_blk_core !: CORE bulk formulation 36 36 LOGICAL , PUBLIC :: ln_blk_mfs !: MFS bulk formulation 37 #if defined key_oasis3 37 #if defined key_oasis3 || defined key_oasis3mct 38 38 LOGICAL , PUBLIC :: lk_oasis = .TRUE. !: OASIS used 39 39 #else … … 42 42 LOGICAL , PUBLIC :: ln_cpl !: ocean-atmosphere coupled formulation 43 43 LOGICAL , PUBLIC :: ln_mixcpl !: ocean-atmosphere forced-coupled mixed formulation 44 LOGICAL , PUBLIC :: ln_wavcpl !: ocean-wave forced-coupled mixed formulation 45 LOGICAL , PUBLIC :: ll_purecpl !: ocean-atmosphere or ocean-wave pure coupled formulation 44 46 LOGICAL , PUBLIC :: ln_dm2dc !: Daily mean to Diurnal Cycle short wave (qsr) 45 47 LOGICAL , PUBLIC :: ln_rnf !: runoffs / runoff mouths 46 48 LOGICAL , PUBLIC :: ln_ssr !: Sea Surface restoring on SST and/or SSS 47 49 LOGICAL , PUBLIC :: ln_apr_dyn !: Atmospheric pressure forcing used on dynamics (ocean & ice) 50 LOGICAL, PUBLIC :: ln_usernfmask = .false. ! use a runoff mask file to merge data received from several models 51 ! -> file rnfmask.nc with the float variable called rnfmask (jpi,jpj,nn_cplmodel) 48 52 INTEGER , PUBLIC :: nn_ice !: flag for ice in the surface boundary condition (=0/1/2/3) 49 53 INTEGER , PUBLIC :: nn_isf !: flag for isf in the surface boundary condition (=0/1/2/3/4) … … 64 68 LOGICAL , PUBLIC :: ln_wave !: true if some coupling with wave model 65 69 LOGICAL , PUBLIC :: ln_cdgw !: true if neutral drag coefficient from wave model 66 LOGICAL , PUBLIC :: ln_sdw !: true if 3d stokes drift from wave model 70 LOGICAL , PUBLIC :: ln_sdw !: true if 3d Stokes drift from wave model 71 LOGICAL , PUBLIC :: ln_stcor !: true if Stokes-Coriolis term is used 72 LOGICAL , PUBLIC :: ln_tauoc !: true if normalized stress from wave is used 73 LOGICAL , PUBLIC :: ln_tauw !: true if ocean stress components from wave is used 74 LOGICAL , PUBLIC :: ln_phioc !: true if wave energy to ocean is used 75 LOGICAL , PUBLIC :: ln_rough !: true if wave roughness equals significant wave height 76 LOGICAL , PUBLIC :: ln_zdfqiao !: Enhanced wave vertical mixing Qiao(2010) formulation flag 77 INTEGER , PUBLIC :: nn_drag ! type of formula to calculate wind stress from wind components 78 INTEGER , PUBLIC :: nn_wmix ! type of wave breaking mixing 67 79 ! 68 80 LOGICAL , PUBLIC :: ln_icebergs !: Icebergs … … 91 103 INTEGER , PUBLIC, PARAMETER :: jp_iam_sas = 2 !: Multi executable configuration - SAS component 92 104 ! (internal OASIS coupling) 105 !!---------------------------------------------------------------------- 106 !! wind stress definition 107 !!---------------------------------------------------------------------- 108 INTEGER, PUBLIC, PARAMETER :: jp_ukmo = 0 ! UKMO SHELF formulation 109 INTEGER, PUBLIC, PARAMETER :: jp_std = 1 ! standard formulation with forced or coupled drag coefficient 110 INTEGER, PUBLIC, PARAMETER :: jp_const = 2 ! standard formulation with constant drag coefficient 111 INTEGER, PUBLIC, PARAMETER :: jp_mcore = 3 ! momentum calculated from core forcing fields 112 113 !!---------------------------------------------------------------------- 114 !! Wave mixing vertical parameterization 115 !!---------------------------------------------------------------------- 116 INTEGER, PUBLIC, PARAMETER :: jp_craigbanner = 0 ! Craig and Banner formulation (original NEMO formulation - 117 ! direct conversion of mechanical to turbulent energy) 118 INTEGER, PUBLIC, PARAMETER :: jp_janssen = 1 ! Janssen formulation - no assumption on direct energy conversion 93 119 !!---------------------------------------------------------------------- 94 120 !! Ocean Surface Boundary Condition fields … … 128 154 #endif 129 155 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xcplmask !: coupling mask for ln_mixcpl (warning: allocated in sbccpl) 156 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xrnfmask !: coupling mask for ln_usernfmask (warning: allocated in sbcrnf) 130 157 131 158 !!---------------------------------------------------------------------- -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcapr.F90
r8058 r11277 26 26 27 27 ! !!* namsbc_apr namelist (Atmospheric PRessure) * 28 LOGICAL, PUBLIC :: cpl_mslp = .FALSE. ! Is the pressure read from coupling? 28 29 LOGICAL, PUBLIC :: ln_apr_obc !: inverse barometer added to OBC ssh data 29 30 LOGICAL, PUBLIC :: ln_ref_apr !: ref. pressure: global mean Patm (F) or a constant (F) 30 REAL(wp) 31 REAL(wp),PUBLIC :: rn_pref ! reference atmospheric pressure [N/m2] 31 32 32 33 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: ssh_ib ! Inverse barometer now sea surface height [m] … … 34 35 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: apr ! atmospheric pressure at kt [N/m2] 35 36 36 REAL(wp) :: tarea ! whole domain mean masked ocean surface37 REAL(wp) :: r1_grau ! = 1.e0 / (grav * rau0)37 REAL(wp), PUBLIC :: tarea ! whole domain mean masked ocean surface 38 REAL(wp), PUBLIC :: r1_grau ! = 1.e0 / (grav * rau0) 38 39 39 40 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_apr ! structure of input fields (file informations, fields read) … … 85 86 IF(lwm) WRITE ( numond, namsbc_apr ) 86 87 ! 87 ALLOCATE( sf_apr(1), STAT=ierror ) !* allocate and fill sf_sst (forcing structure) with sn_sst 88 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_apr: unable to allocate sf_apr structure' ) 89 ! 90 CALL fld_fill( sf_apr, (/ sn_apr /), cn_dir, 'sbc_apr', 'Atmospheric pressure ', 'namsbc_apr' ) 91 ALLOCATE( sf_apr(1)%fnow(jpi,jpj,1) ) 92 IF( sn_apr%ln_tint ) ALLOCATE( sf_apr(1)%fdta(jpi,jpj,1,2) ) 88 IF( .NOT. cpl_mslp ) THEN 89 ALLOCATE( sf_apr(1), STAT=ierror ) !* allocate and fill sf_sst (forcing structure) with sn_sst 90 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_apr: unable to allocate sf_apr structure' ) 91 ! 92 CALL fld_fill( sf_apr, (/ sn_apr /), cn_dir, 'sbc_apr', 'Atmospheric pressure ', 'namsbc_apr' ) 93 ALLOCATE( sf_apr(1)%fnow(jpi,jpj,1) ) 94 IF( sn_apr%ln_tint ) ALLOCATE( sf_apr(1)%fdta(jpi,jpj,1,2) ) 95 ENDIF 93 96 ALLOCATE( ssh_ib(jpi,jpj) , ssh_ibb(jpi,jpj) ) 94 97 ALLOCATE( apr (jpi,jpj) ) … … 96 99 IF(lwp) THEN !* control print 97 100 WRITE(numout,*) 98 WRITE(numout,*) ' Namelist namsbc_apr : Atmospheric PRessure as extrenal forcing' 101 IF( cpl_mslp ) THEN 102 WRITE(numout,*) ' Namelist namsbc_apr : Atmospheric Pressure as extrenal coupling' 103 ELSE 104 WRITE(numout,*) ' Namelist namsbc_apr : Atmospheric Pressure as extrenal forcing' 105 ENDIF 99 106 WRITE(numout,*) ' ref. pressure: global mean Patm (T) or a constant (F) ln_ref_apr = ', ln_ref_apr 100 107 ENDIF … … 119 126 ENDIF 120 127 121 ! ! ========================== ! 122 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN ! At each sbc time-step ! 123 ! ! ===========+++============ ! 124 ! 125 IF( kt /= nit000 ) ssh_ibb(:,:) = ssh_ib(:,:) !* Swap of ssh_ib fields 126 ! 127 CALL fld_read( kt, nn_fsbc, sf_apr ) !* input Patm provided at kt + nn_fsbc/2 128 ! 129 ! !* update the reference atmospheric pressure (if necessary) 130 IF( ln_ref_apr ) rn_pref = glob_sum( sf_apr(1)%fnow(:,:,1) * e1e2t(:,:) ) / tarea 131 ! 132 ! !* Patm related forcing at kt 133 ssh_ib(:,:) = - ( sf_apr(1)%fnow(:,:,1) - rn_pref ) * r1_grau ! equivalent ssh (inverse barometer) 134 apr (:,:) = sf_apr(1)%fnow(:,:,1) ! atmospheric pressure 135 ! 136 CALL iom_put( "ssh_ib", ssh_ib ) !* output the inverse barometer ssh 137 ENDIF 138 139 ! ! ---------------------------------------- ! 140 IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! 141 ! ! ---------------------------------------- ! 142 ! !* Restart: read in restart file 143 IF( ln_rstart .AND. iom_varid( numror, 'ssh_ibb', ldstop = .FALSE. ) > 0 ) THEN 144 IF(lwp) WRITE(numout,*) 'sbc_apr: ssh_ibb read in the restart file' 145 CALL iom_get( numror, jpdom_autoglo, 'ssh_ibb', ssh_ibb ) ! before inv. barometer ssh 128 IF( .NOT. cpl_mslp ) THEN ! ========================== ! 129 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN ! At each sbc time-step ! 130 ! ! ===========+++============ ! 146 131 ! 147 ELSE !* no restart: set from nit000 values 148 IF(lwp) WRITE(numout,*) 'sbc_apr: ssh_ibb set to nit000 values' 149 ssh_ibb(:,:) = ssh_ib(:,:) 132 IF( kt /= nit000 ) ssh_ibb(:,:) = ssh_ib(:,:) !* Swap of ssh_ib fields 133 ! 134 CALL fld_read( kt, nn_fsbc, sf_apr ) !* input Patm provided at kt + nn_fsbc/2 135 ! 136 ! !* update the reference atmospheric pressure (if necessary) 137 IF( ln_ref_apr ) rn_pref = glob_sum( sf_apr(1)%fnow(:,:,1) * e1e2t(:,:) ) / tarea 138 ! 139 ! !* Patm related forcing at kt 140 ssh_ib(:,:) = - ( sf_apr(1)%fnow(:,:,1) - rn_pref ) * r1_grau ! equivalent ssh (inverse barometer) 141 apr (:,:) = sf_apr(1)%fnow(:,:,1) ! atmospheric pressure 142 ! 143 CALL iom_put( "ssh_ib", ssh_ib ) !* output the inverse barometer ssh 150 144 ENDIF 151 ENDIF 152 ! ! ---------------------------------------- ! 153 IF( lrst_oce ) THEN ! Write in the ocean restart file ! 154 ! ! ---------------------------------------- ! 155 IF(lwp) WRITE(numout,*) 156 IF(lwp) WRITE(numout,*) 'sbc_apr : ssh_ib written in ocean restart file at it= ', kt,' date= ', ndastp 157 IF(lwp) WRITE(numout,*) '~~~~' 158 CALL iom_rstput( kt, nitrst, numrow, 'ssh_ibb' , ssh_ib ) 145 146 ! ! ---------------------------------------- ! 147 IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! 148 ! ! ---------------------------------------- ! 149 ! !* Restart: read in restart file 150 IF( ln_rstart .AND. iom_varid( numror, 'ssh_ibb', ldstop = .FALSE. ) > 0 ) THEN 151 IF(lwp) WRITE(numout,*) 'sbc_apr: ssh_ibb read in the restart file' 152 CALL iom_get( numror, jpdom_autoglo, 'ssh_ibb', ssh_ibb ) ! before inv. barometer ssh 153 ! 154 ELSE !* no restart: set from nit000 values 155 IF(lwp) WRITE(numout,*) 'sbc_apr: ssh_ibb set to nit000 values' 156 ssh_ibb(:,:) = ssh_ib(:,:) 157 ENDIF 158 ENDIF 159 ! ! ---------------------------------------- ! 160 IF( lrst_oce ) THEN ! Write in the ocean restart file ! 161 ! ! ---------------------------------------- ! 162 IF(lwp) WRITE(numout,*) 163 IF(lwp) WRITE(numout,*) 'sbc_apr : ssh_ib written in ocean restart file at it= ', kt,' date= ', ndastp 164 IF(lwp) WRITE(numout,*) '~~~~' 165 CALL iom_rstput( kt, nitrst, numrow, 'ssh_ibb' , ssh_ib ) 166 ENDIF 159 167 ENDIF 160 168 ! -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r8058 r11277 319 319 & Cd, Ch, Ce, zt_zu, zq_zu ) 320 320 321 IF ( ln_cdgw .AND. nn_drag==jp_std ) Cd(:,:) = cdn_wave(:,:) 321 322 ! ... tau module, i and j component 322 323 DO jj = 1, jpj … … 737 738 738 739 !! Neutral coefficients at 10m: 739 IF( ln_cdgw ) THEN ! wave drag case740 IF( ln_cdgw .AND. nn_drag==jp_mcore ) THEN ! wave drag case 740 741 cdn_wave(:,:) = cdn_wave(:,:) + rsmall * ( 1._wp - tmask(:,:,1) ) 741 742 ztmp0 (:,:) = cdn_wave(:,:) … … 783 784 END IF 784 785 785 IF( ln_cdgw ) THEN ! surface wave case786 IF( ln_cdgw .AND. nn_drag==jp_mcore ) THEN ! surface wave case 786 787 sqrt_Cd = vkarmn / ( vkarmn / sqrt_Cd_n10 - zpsi_m_u ) 787 788 Cd = sqrt_Cd * sqrt_Cd -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_mfs.F90
r8058 r11277 24 24 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 25 25 USE prtctl ! Print control 26 USE sbcwave,ONLY : cdn_wave !wave module27 26 28 27 IMPLICIT NONE -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
r8058 r11277 23 23 USE sbcapr 24 24 USE sbcdcy ! surface boundary condition: diurnal cycle 25 USE sbcwave ! surface boundary condition: waves 25 26 USE phycst ! physical constants 26 27 #if defined key_lim3 … … 41 42 USE timing ! Timing 42 43 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 44 #if defined key_oasis3 || defined key_oasis3mct 45 USE mod_oasis ! OASIS3-MCT module 46 #endif 43 47 USE eosbn2 44 48 USE sbcrnf , ONLY : l_rnfcpl … … 105 109 INTEGER, PARAMETER :: jpr_e3t1st = 41 ! first T level thickness 106 110 INTEGER, PARAMETER :: jpr_fraqsr = 42 ! fraction of solar net radiation absorbed in the first ocean level 107 INTEGER, PARAMETER :: jprcv = 42 ! total number of fields received 111 INTEGER, PARAMETER :: jpr_mslp = 43 ! mean sea level pressure 112 INTEGER, PARAMETER :: jpr_hsig = 44 ! Hsig 113 INTEGER, PARAMETER :: jpr_phioc = 45 ! Wave=>ocean energy flux 114 INTEGER, PARAMETER :: jpr_sdrftx = 46 ! Stokes drift on grid 1 115 INTEGER, PARAMETER :: jpr_sdrfty = 47 ! Stokes drift on grid 2 116 INTEGER, PARAMETER :: jpr_wper = 48 ! Mean wave period 117 INTEGER, PARAMETER :: jpr_wnum = 49 ! Mean wavenumber 118 INTEGER, PARAMETER :: jpr_tauoc = 50 ! Stress fraction adsorbed by waves 119 INTEGER, PARAMETER :: jpr_wdrag = 51 ! Neutral surface drag coefficient 120 INTEGER, PARAMETER :: jpr_wfreq = 52 ! Wave peak frequency 121 INTEGER, PARAMETER :: jpr_tauwx = 53 ! x component of the ocean stress from waves 122 INTEGER, PARAMETER :: jpr_tauwy = 54 ! y component of the ocean stress from waves 123 INTEGER, PARAMETER :: jprcv = 54 ! total number of fields received 108 124 109 125 INTEGER, PARAMETER :: jps_fice = 1 ! ice fraction sent to the atmosphere … … 135 151 INTEGER, PARAMETER :: jps_e3t1st = 27 ! first level depth (vvl) 136 152 INTEGER, PARAMETER :: jps_fraqsr = 28 ! fraction of solar net radiation absorbed in the first ocean level 137 INTEGER, PARAMETER :: jpsnd = 28 ! total number of fields sended 153 INTEGER, PARAMETER :: jps_ficet = 29 ! total ice fraction 154 INTEGER, PARAMETER :: jps_ocxw = 30 ! currents on grid 1 155 INTEGER, PARAMETER :: jps_ocyw = 31 ! currents on grid 2 156 INTEGER, PARAMETER :: jps_wlev = 32 ! water level 157 INTEGER, PARAMETER :: jpsnd = 32 ! total number of fields sent 138 158 139 159 ! !!** namelist namsbc_cpl ** … … 149 169 ! Received from the atmosphere ! 150 170 TYPE(FLD_C) :: sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf 151 TYPE(FLD_C) :: sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2 171 TYPE(FLD_C) :: sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp 172 ! Send to waves 173 TYPE(FLD_C) :: sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev 174 ! Received from waves 175 TYPE(FLD_C) :: sn_rcv_hsig,sn_rcv_phioc,sn_rcv_sdrft,sn_rcv_wper, & 176 sn_rcv_wfreq,sn_rcv_wnum,sn_rcv_tauoc,sn_rcv_tauw, & 177 sn_rcv_wdrag 152 178 ! Other namelist parameters ! 153 179 INTEGER :: nn_cplmodel ! Maximum number of models to/from which NEMO is potentialy sending/receiving data … … 161 187 162 188 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: albedo_oce_mix ! ocean albedo sent to atmosphere (mix clear/overcast sky) 163 189 164 190 INTEGER , ALLOCATABLE, SAVE, DIMENSION( :) :: nrcvinfo ! OASIS info argument 165 191 … … 179 205 !! *** FUNCTION sbc_cpl_alloc *** 180 206 !!---------------------------------------------------------------------- 181 INTEGER :: ierr( 3)207 INTEGER :: ierr(4) 182 208 !!---------------------------------------------------------------------- 183 209 ierr(:) = 0 … … 188 214 ALLOCATE( a_i(jpi,jpj,1) , STAT=ierr(2) ) ! used in sbcice_if.F90 (done here as there is no sbc_ice_if_init) 189 215 #endif 190 ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) ) 191 ! 216 ! ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) ) 217 ! Hardwire three models as nn_cplmodel has not been read in from the namelist yet. 218 ALLOCATE( xcplmask(jpi,jpj,0:3) , STAT=ierr(3) ) 219 ! 220 IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(4) ) 221 192 222 sbc_cpl_alloc = MAXVAL( ierr ) 193 223 IF( lk_mpp ) CALL mpp_sum ( sbc_cpl_alloc ) … … 216 246 REAL(wp), POINTER, DIMENSION(:,:) :: zacs, zaos 217 247 !! 218 NAMELIST/namsbc_cpl/ sn_snd_temp, sn_snd_alb , sn_snd_thick, sn_snd_crt , sn_snd_co2, & 219 & sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau , sn_rcv_dqnsdt, sn_rcv_qsr, & 220 & sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal , sn_rcv_iceflx, & 221 & sn_rcv_co2 , nn_cplmodel , ln_usecplmask 248 NAMELIST/namsbc_cpl/ sn_snd_temp , sn_snd_alb , sn_snd_thick , sn_snd_crt , sn_snd_co2, & 249 & sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau , sn_rcv_dqnsdt, sn_rcv_qsr, & 250 & sn_snd_ifrac, sn_snd_crtw , sn_snd_wlev , sn_rcv_hsig , sn_rcv_phioc , & 251 & sn_rcv_sdrft, sn_rcv_wper , sn_rcv_wnum , sn_rcv_wfreq, sn_rcv_tauoc, & 252 & sn_rcv_wdrag, sn_rcv_qns , sn_rcv_emp , sn_rcv_rnf , sn_rcv_cal , & 253 & sn_rcv_iceflx, sn_rcv_co2 , sn_rcv_mslp , sn_rcv_tauw , & 254 & nn_cplmodel, ln_usecplmask, ln_usernfmask 222 255 !!--------------------------------------------------------------------- 223 256 ! … … 260 293 WRITE(numout,*)' sea ice heat fluxes = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 261 294 WRITE(numout,*)' atm co2 = ', TRIM(sn_rcv_co2%cldes ), ' (', TRIM(sn_rcv_co2%clcat ), ')' 295 WRITE(numout,*)' mean sea level pressure = ', TRIM(sn_rcv_mslp%cldes ), ' (', TRIM(sn_rcv_mslp%clcat ), ')' 296 WRITE(numout,*)' significant wave heigth = ', TRIM(sn_rcv_hsig%cldes ), ' (', TRIM(sn_rcv_hsig%clcat ), ')' 297 WRITE(numout,*)' wave to oce energy flux = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')' 298 WRITE(numout,*)' Surface Stokes drift u,v = ', TRIM(sn_rcv_sdrft%cldes ), ' (', TRIM(sn_rcv_sdrft%clcat ), ')' 299 WRITE(numout,*)' Mean wave period = ', TRIM(sn_rcv_wper%cldes ), ' (', TRIM(sn_rcv_wper%clcat ), ')' 300 WRITE(numout,*)' Mean wave number = ', TRIM(sn_rcv_wnum%cldes ), ' (', TRIM(sn_rcv_wnum%clcat ), ')' 301 WRITE(numout,*)' Wave peak frequency = ', TRIM(sn_rcv_wfreq%cldes ), ' (', TRIM(sn_rcv_wfreq%clcat ), ')' 302 WRITE(numout,*)' Stress frac adsorbed by waves = ', TRIM(sn_rcv_tauoc%cldes ), ' (', TRIM(sn_rcv_tauoc%clcat ), ')' 303 WRITE(numout,*)' Stress components by waves = ', TRIM(sn_rcv_tauw%cldes ), ' (', TRIM(sn_rcv_tauw%clcat ), ')' 304 WRITE(numout,*)' Neutral surf drag coefficient = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' 262 305 WRITE(numout,*)' sent fields (multiple ice categories)' 263 306 WRITE(numout,*)' surface temperature = ', TRIM(sn_snd_temp%cldes ), ' (', TRIM(sn_snd_temp%clcat ), ')' 264 307 WRITE(numout,*)' albedo = ', TRIM(sn_snd_alb%cldes ), ' (', TRIM(sn_snd_alb%clcat ), ')' 265 308 WRITE(numout,*)' ice/snow thickness = ', TRIM(sn_snd_thick%cldes ), ' (', TRIM(sn_snd_thick%clcat ), ')' 309 WRITE(numout,*)' total ice fraction = ', TRIM(sn_snd_ifrac%cldes ), ' (', TRIM(sn_snd_ifrac%clcat ), ')' 266 310 WRITE(numout,*)' surface current = ', TRIM(sn_snd_crt%cldes ), ' (', TRIM(sn_snd_crt%clcat ), ')' 267 311 WRITE(numout,*)' - referential = ', sn_snd_crt%clvref … … 269 313 WRITE(numout,*)' - mesh = ', sn_snd_crt%clvgrd 270 314 WRITE(numout,*)' oce co2 flux = ', TRIM(sn_snd_co2%cldes ), ' (', TRIM(sn_snd_co2%clcat ), ')' 315 WRITE(numout,*)' water level = ', TRIM(sn_snd_wlev%cldes ), ' (', TRIM(sn_snd_wlev%clcat ), ')' 316 WRITE(numout,*)' surface current to waves = ', TRIM(sn_snd_crtw%cldes ), ' (', TRIM(sn_snd_crtw%clcat ), ')' 317 WRITE(numout,*)' - referential = ', sn_snd_crtw%clvref 318 WRITE(numout,*)' - orientation = ', sn_snd_crtw%clvor 319 WRITE(numout,*)' - mesh = ', sn_snd_crtw%clvgrd 271 320 WRITE(numout,*)' nn_cplmodel = ', nn_cplmodel 272 321 WRITE(numout,*)' ln_usecplmask = ', ln_usecplmask 322 WRITE(numout,*)' ln_usernfmask = ', ln_usernfmask 273 323 ENDIF 274 324 275 325 ! ! allocate sbccpl arrays 276 IF( sbc_cpl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' )326 !IF( sbc_cpl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' ) 277 327 278 328 ! ================================ ! … … 307 357 ! 308 358 ! Vectors: change of sign at north fold ONLY if on the local grid 359 IF( TRIM( sn_rcv_tau%cldes ) == 'oce only' .OR. TRIM(sn_rcv_tau%cldes ) == 'oce and ice') THEN ! avoid working with the atmospheric fields if they are not coupled 309 360 IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' ) srcv(jpr_otx1:jpr_itz2)%nsgn = -1. 310 361 … … 337 388 srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point 338 389 srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point 339 srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 390 !srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 391 ! Currently needed for HadGEM3 - but shouldn't affect anyone else for the moment 392 srcv(jpr_otx1)%laction = .TRUE. 393 srcv(jpr_oty1)%laction = .TRUE. 394 ! 340 395 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only 341 396 CASE( 'T,I' ) … … 374 429 srcv(jpr_ity1)%clgrid = 'V' ! i.e. it is always at U- & V-points for i- & j-comp. resp. 375 430 ENDIF 431 ENDIF 376 432 377 433 ! ! ------------------------- ! … … 403 459 IF( TRIM( sn_rcv_rnf%cldes ) == 'coupled' ) THEN 404 460 srcv(jpr_rnf)%laction = .TRUE. 405 l_rnfcpl = . TRUE.! -> no need to read runoffs in sbcrnf461 l_rnfcpl = .NOT. ln_usernfmask ! -> no need to read runoffs in sbcrnf 406 462 ln_rnf = nn_components /= jp_iam_sas ! -> force to go through sbcrnf if not sas 407 463 IF(lwp) WRITE(numout,*) … … 470 526 ! ! ------------------------- ! 471 527 srcv(jpr_co2 )%clname = 'O_AtmCO2' ; IF( TRIM(sn_rcv_co2%cldes ) == 'coupled' ) srcv(jpr_co2 )%laction = .TRUE. 528 529 ! ! ------------------------- ! 530 ! ! Mean Sea Level Pressure ! 531 ! ! ------------------------- ! 532 srcv(jpr_mslp)%clname = 'O_MSLP' 533 IF( TRIM(sn_rcv_mslp%cldes ) == 'coupled' ) THEN 534 srcv(jpr_mslp)%laction = .TRUE. 535 cpl_mslp = .TRUE. 536 ENDIF 537 472 538 ! ! ------------------------- ! 473 539 ! ! topmelt and botmelt ! … … 483 549 srcv(jpr_topm:jpr_botm)%laction = .TRUE. 484 550 ENDIF 551 ! ! ------------------------- ! 552 ! ! Wave breaking ! 553 ! ! ------------------------- ! 554 srcv(jpr_hsig)%clname = 'O_Hsigwa' ! significant wave height 555 IF( TRIM(sn_rcv_hsig%cldes ) == 'coupled' ) THEN 556 srcv(jpr_hsig)%laction = .TRUE. 557 cpl_hsig = .TRUE. 558 ENDIF 559 srcv(jpr_phioc)%clname = 'O_PhiOce' ! wave to ocean energy 560 IF( TRIM(sn_rcv_phioc%cldes ) == 'coupled' ) THEN 561 srcv(jpr_phioc)%laction = .TRUE. 562 cpl_phioc = .TRUE. 563 ENDIF 564 srcv(jpr_sdrftx)%clname = 'O_Sdrfx' ! Stokes drift in the u direction 565 srcv(jpr_sdrfty)%clname = 'O_Sdrfy' ! Stokes drift in the v direction 566 IF( TRIM(sn_rcv_sdrft%cldes ) == 'coupled' ) THEN 567 srcv(jpr_sdrftx)%laction = .TRUE. 568 srcv(jpr_sdrfty)%laction = .TRUE. 569 cpl_sdrft = .TRUE. 570 ENDIF 571 srcv(jpr_wper)%clname = 'O_WPer' ! mean wave period 572 IF( TRIM(sn_rcv_wper%cldes ) == 'coupled' ) THEN 573 srcv(jpr_wper)%laction = .TRUE. 574 cpl_wper = .TRUE. 575 ENDIF 576 srcv(jpr_wfreq)%clname = 'O_WFreq' ! wave peak frequency 577 IF( TRIM(sn_rcv_wfreq%cldes ) == 'coupled' ) THEN 578 srcv(jpr_wfreq)%laction = .TRUE. 579 cpl_wfreq = .TRUE. 580 ENDIF 581 srcv(jpr_wnum)%clname = 'O_WNum' ! mean wave number 582 IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' ) THEN 583 srcv(jpr_wnum)%laction = .TRUE. 584 cpl_wnum = .TRUE. 585 ENDIF 586 srcv(jpr_tauoc)%clname = 'O_TauOce' ! stress fraction adsorbed by the wave 587 IF( TRIM(sn_rcv_tauoc%cldes ) == 'coupled' ) THEN 588 srcv(jpr_tauoc)%laction = .TRUE. 589 cpl_tauoc = .TRUE. 590 ENDIF 591 srcv(jpr_tauwx)%clname = 'O_Tauwx' ! ocean stress from wave in the x direction 592 srcv(jpr_tauwy)%clname = 'O_Tauwy' ! ocean stress from wave in the y direction 593 IF( TRIM(sn_rcv_tauw%cldes ) == 'coupled' ) THEN 594 srcv(jpr_tauwx)%laction = .TRUE. 595 srcv(jpr_tauwy)%laction = .TRUE. 596 cpl_tauw = .TRUE. 597 ENDIF 598 srcv(jpr_wdrag)%clname = 'O_WDrag' ! neutral surface drag coefficient 599 IF( TRIM(sn_rcv_wdrag%cldes ) == 'coupled' ) THEN 600 srcv(jpr_wdrag)%laction = .TRUE. 601 cpl_wdrag = .TRUE. 602 ENDIF 603 ! 604 IF( srcv(jpr_tauoc)%laction .AND. srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction ) & 605 CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 606 '(sn_rcv_tauoc=coupled and sn_rcv_tauw=coupled)' ) 607 ! 608 ! 485 609 ! ! ------------------------------- ! 486 610 ! ! OPA-SAS coupling - rcv by opa ! … … 637 761 ! ! ------------------------- ! 638 762 ssnd(jps_fice)%clname = 'OIceFrc' 763 ssnd(jps_ficet)%clname = 'OIceFrcT' 639 764 ssnd(jps_hice)%clname = 'OIceTck' 640 765 ssnd(jps_hsnw)%clname = 'OSnwTck' … … 645 770 ENDIF 646 771 772 IF (TRIM( sn_snd_ifrac%cldes ) == 'coupled') ssnd(jps_ficet)%laction = .TRUE. 773 647 774 SELECT CASE ( TRIM( sn_snd_thick%cldes ) ) 648 775 CASE( 'none' ) ! nothing to do … … 665 792 ssnd(jps_ocy1)%clname = 'O_OCury1' ; ssnd(jps_ivy1)%clname = 'O_IVely1' 666 793 ssnd(jps_ocz1)%clname = 'O_OCurz1' ; ssnd(jps_ivz1)%clname = 'O_IVelz1' 794 ssnd(jps_ocxw)%clname = 'O_OCurxw' 795 ssnd(jps_ocyw)%clname = 'O_OCuryw' 667 796 ! 668 797 ssnd(jps_ocx1:jps_ivz1)%nsgn = -1. ! vectors: change of the sign at the north fold … … 685 814 END SELECT 686 815 816 ssnd(jps_ocxw:jps_ocyw)%nsgn = -1. ! vectors: change of the sign at the north fold 817 818 IF( sn_snd_crtw%clvgrd == 'U,V' ) THEN 819 ssnd(jps_ocxw)%clgrid = 'U' ; ssnd(jps_ocyw)%clgrid = 'V' 820 ELSE IF( sn_snd_crtw%clvgrd /= 'T' ) THEN 821 CALL ctl_stop( 'sn_snd_crtw%clvgrd must be equal to T' ) 822 ENDIF 823 IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) ssnd(jps_ocxw:jps_ocyw)%nsgn = 1. 824 SELECT CASE( TRIM( sn_snd_crtw%cldes ) ) 825 CASE( 'none' ) ; ssnd(jps_ocxw:jps_ocyw)%laction = .FALSE. 826 CASE( 'oce only' ) ; ssnd(jps_ocxw:jps_ocyw)%laction = .TRUE. 827 CASE( 'weighted oce and ice' ) ! nothing to do 828 CASE( 'mixed oce-ice' ) ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE. 829 CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crtw%cldes' ) 830 END SELECT 831 687 832 ! ! ------------------------- ! 688 833 ! ! CO2 flux ! 689 834 ! ! ------------------------- ! 690 835 ssnd(jps_co2)%clname = 'O_CO2FLX' ; IF( TRIM(sn_snd_co2%cldes) == 'coupled' ) ssnd(jps_co2 )%laction = .TRUE. 836 837 ! ! ------------------------- ! 838 ! ! Sea surface height ! 839 ! ! ------------------------- ! 840 ssnd(jps_wlev)%clname = 'O_Wlevel' ; IF( TRIM(sn_snd_wlev%cldes) == 'coupled' ) ssnd(jps_wlev)%laction = .TRUE. 691 841 692 842 ! ! ------------------------------- ! … … 783 933 IF( ln_dm2dc .AND. ln_cpl .AND. ncpl_qsr_freq /= 86400 ) & 784 934 & CALL ctl_stop( 'sbc_cpl_init: diurnal cycle reconstruction (ln_dm2dc) needs daily couping for solar radiation' ) 785 ncpl_qsr_freq = 86400 / ncpl_qsr_freq935 IF( ln_dm2dc .AND. ln_cpl ) ncpl_qsr_freq = 86400 / ncpl_qsr_freq 786 936 787 937 CALL wrk_dealloc( jpi,jpj, zacs, zaos ) … … 837 987 !! emp upward mass flux [evap. - precip. (- runoffs) (- calving)] (ocean only case) 838 988 !!---------------------------------------------------------------------- 989 USE sbcflx , ONLY : ln_shelf_flx 990 USE sbcssm , ONLY : sbc_ssm_cpl 991 USE lib_fortran ! distributed memory computing library 992 839 993 INTEGER, INTENT(in) :: kt ! ocean model time step index 840 994 INTEGER, INTENT(in) :: k_fsbc ! frequency of sbc (-> ice model) computation … … 845 999 INTEGER :: ji, jj, jn ! dummy loop indices 846 1000 INTEGER :: isec ! number of seconds since nit000 (assuming rdttra did not change since nit000) 1001 INTEGER :: ikchoix 847 1002 REAL(wp) :: zcumulneg, zcumulpos ! temporary scalars 848 1003 REAL(wp) :: zcoef ! temporary scalar … … 850 1005 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 851 1006 REAL(wp) :: zzx, zzy ! temporary variables 852 REAL(wp), POINTER, DIMENSION(:,:) :: ztx, zty, z msk, zemp, zqns, zqsr1007 REAL(wp), POINTER, DIMENSION(:,:) :: ztx, zty, ztx2, zty2, zmsk, zemp, zqns, zqsr 853 1008 !!---------------------------------------------------------------------- 854 1009 ! 855 1010 IF( nn_timing == 1 ) CALL timing_start('sbc_cpl_rcv') 856 1011 ! 857 CALL wrk_alloc( jpi,jpj, ztx, zty, z msk, zemp, zqns, zqsr )1012 CALL wrk_alloc( jpi,jpj, ztx, zty, ztx2, zty2, zmsk, zemp, zqns, zqsr ) 858 1013 ! 859 1014 IF( ln_mixcpl ) zmsk(:,:) = 1. - xcplmask(:,:,0) … … 893 1048 IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid 894 1049 ! ! (geographical to local grid -> rotate the components) 895 CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx ) 896 IF( srcv(jpr_otx2)%laction ) THEN 897 CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty ) 898 ELSE 899 CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty ) 1050 IF( srcv(jpr_otx1)%clgrid == 'U' .AND. (.NOT. srcv(jpr_otx2)%laction) ) THEN 1051 ! Temporary code for HadGEM3 - will be removed eventually. 1052 ! Only applies when we have only taux on U grid and tauy on V grid 1053 DO jj=2,jpjm1 1054 DO ji=2,jpim1 1055 ztx(ji,jj)=0.25*vmask(ji,jj,1) & 1056 *(frcv(jpr_otx1)%z3(ji,jj,1)+frcv(jpr_otx1)%z3(ji-1,jj,1) & 1057 +frcv(jpr_otx1)%z3(ji,jj+1,1)+frcv(jpr_otx1)%z3(ji-1,jj+1,1)) 1058 zty(ji,jj)=0.25*umask(ji,jj,1) & 1059 *(frcv(jpr_oty1)%z3(ji,jj,1)+frcv(jpr_oty1)%z3(ji+1,jj,1) & 1060 +frcv(jpr_oty1)%z3(ji,jj-1,1)+frcv(jpr_oty1)%z3(ji+1,jj-1,1)) 1061 ENDDO 1062 ENDDO 1063 1064 ikchoix = 1 1065 CALL repcmo(frcv(jpr_otx1)%z3(:,:,1),zty,ztx,frcv(jpr_oty1)%z3(:,:,1),ztx2,zty2,ikchoix) 1066 CALL lbc_lnk (ztx2,'U', -1. ) 1067 CALL lbc_lnk (zty2,'V', -1. ) 1068 frcv(jpr_otx1)%z3(:,:,1)=ztx2(:,:) 1069 frcv(jpr_oty1)%z3(:,:,1)=zty2(:,:) 1070 ELSE 1071 CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx ) 1072 frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st component on the 1st grid 1073 IF( srcv(jpr_otx2)%laction ) THEN 1074 CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty ) 1075 ELSE 1076 CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty ) 1077 ENDIF 1078 frcv(jpr_oty1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd component on the 2nd grid 900 1079 ENDIF 901 frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:) ! overwrite 1st component on the 1st grid902 frcv(jpr_oty1)%z3(:,:,1) = zty(:,:) ! overwrite 2nd component on the 2nd grid903 1080 ENDIF 904 1081 ! … … 919 1096 ELSE ! No dynamical coupling ! 920 1097 ! ! ========================= ! 1098 ! it is possible that the momentum is calculated from the winds (ln_shelf_flx) and a coupled drag coefficient 1099 IF( srcv(jpr_wdrag)%laction .AND. ln_shelf_flx .AND. ln_cdgw .AND. nn_drag == jp_std ) THEN 1100 DO jj = 1, jpjm1 1101 DO ji = 1, jpim1 1102 ! here utau and vtau should contain the wind components as read from the forcing files, 1103 ! and the wind module should already be properly calculated 1104 frcv(jpr_otx1)%z3(ji,jj,1) = zrhoa * 0.5 * ( frcv(jpr_wdrag)%z3(ji,jj,1) + frcv(jpr_wdrag)%z3(ji+1,jj,1) ) * & 1105 utau(ji,jj) * 0.5 * ( wndm(ji,jj) + wndm(ji+1,jj) ) 1106 frcv(jpr_oty1)%z3(ji,jj,1) = zrhoa * 0.5 * ( frcv(jpr_wdrag)%z3(ji,jj,1) + frcv(jpr_wdrag)%z3(ji,jj+1,1) ) * & 1107 vtau(ji,jj) * 0.5 * ( wndm(ji,jj) + wndm(ji,jj+1) ) 1108 utau(ji,jj) = frcv(jpr_otx1)%z3(ji,jj,1) 1109 vtau(ji,jj) = frcv(jpr_oty1)%z3(ji,jj,1) 1110 END DO 1111 END DO 1112 CALL lbc_lnk_multi( frcv(jpr_otx1)%z3(:,:,1), 'U', -1. , frcv(jpr_oty1)%z3(:,:,1), 'V', -1. , & 1113 utau(:,:), 'U', -1. , vtau(:,:), 'V', -1. ) 1114 llnewtx = .TRUE. 1115 ELSE 921 1116 frcv(jpr_otx1)%z3(:,:,1) = 0.e0 ! here simply set to zero 922 1117 frcv(jpr_oty1)%z3(:,:,1) = 0.e0 ! an external read in a file can be added instead 923 1118 llnewtx = .TRUE. 1119 ENDIF 924 1120 ! 925 1121 ENDIF … … 941 1137 END DO 942 1138 CALL lbc_lnk( frcv(jpr_taum)%z3(:,:,1), 'T', 1. ) 1139 IF( .NOT. srcv(jpr_otx1)%laction .AND. srcv(jpr_wdrag)%laction .AND. & 1140 ln_shelf_flx .AND. ln_cdgw .AND. nn_drag == jp_std ) & 1141 taum(:,:) = frcv(jpr_taum)%z3(:,:,1) 943 1142 llnewtau = .TRUE. 944 1143 ELSE … … 955 1154 ! ! ========================= ! 956 1155 ! ! 10 m wind speed ! (wndm) 1156 ! ! include wave drag coef ! (wndm) 957 1157 ! ! ========================= ! 958 1158 ! … … 965 1165 !CDIR NOVERRCHK 966 1166 DO ji = 1, jpi 1167 IF( ln_shelf_flx ) THEN ! the 10 wind module is properly calculated before if ln_shelf_flx 1168 frcv(jpr_w10m)%z3(ji,jj,1) = wndm(ji,jj) 1169 ELSE 967 1170 frcv(jpr_w10m)%z3(ji,jj,1) = SQRT( frcv(jpr_taum)%z3(ji,jj,1) * zcoef ) 1171 ENDIF 968 1172 END DO 969 1173 END DO … … 975 1179 IF( MOD( kt-1, k_fsbc ) == 0 ) THEN 976 1180 ! 1181 ! if ln_wavcpl, the fields already contain the right information from forcing even if not ln_mixcpl 977 1182 IF( ln_mixcpl ) THEN 978 utau(:,:) = utau(:,:) * xcplmask(:,:,0) + frcv(jpr_otx1)%z3(:,:,1) * zmsk(:,:) 979 vtau(:,:) = vtau(:,:) * xcplmask(:,:,0) + frcv(jpr_oty1)%z3(:,:,1) * zmsk(:,:) 980 taum(:,:) = taum(:,:) * xcplmask(:,:,0) + frcv(jpr_taum)%z3(:,:,1) * zmsk(:,:) 981 wndm(:,:) = wndm(:,:) * xcplmask(:,:,0) + frcv(jpr_w10m)%z3(:,:,1) * zmsk(:,:) 982 ELSE 1183 IF( srcv(jpr_otx1)%laction ) THEN 1184 utau(:,:) = utau(:,:) * xcplmask(:,:,0) + frcv(jpr_otx1)%z3(:,:,1) * zmsk(:,:) 1185 vtau(:,:) = vtau(:,:) * xcplmask(:,:,0) + frcv(jpr_oty1)%z3(:,:,1) * zmsk(:,:) 1186 ENDIF 1187 IF( srcv(jpr_taum)%laction ) & 1188 taum(:,:) = taum(:,:) * xcplmask(:,:,0) + frcv(jpr_taum)%z3(:,:,1) * zmsk(:,:) 1189 IF( srcv(jpr_w10m)%laction ) & 1190 wndm(:,:) = wndm(:,:) * xcplmask(:,:,0) + frcv(jpr_w10m)%z3(:,:,1) * zmsk(:,:) 1191 ELSE IF( ll_purecpl ) THEN 983 1192 utau(:,:) = frcv(jpr_otx1)%z3(:,:,1) 984 1193 vtau(:,:) = frcv(jpr_oty1)%z3(:,:,1) … … 996 1205 IF( srcv(jpr_co2)%laction ) atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1) 997 1206 #endif 1207 ! 1208 ! ! ========================= ! 1209 ! ! Mean Sea Level Pressure ! (taum) 1210 ! ! ========================= ! 1211 ! 1212 IF( srcv(jpr_mslp)%laction ) THEN ! UKMO SHELF effect of atmospheric pressure on SSH 1213 IF( kt /= nit000 ) ssh_ibb(:,:) = ssh_ib(:,:) !* Swap of ssh_ib fields 1214 1215 ! !* update the reference atmospheric pressure (if necessary) 1216 IF( ln_ref_apr ) rn_pref = glob_sum( frcv(jpr_mslp)%z3(:,:,1) * e1e2t(:,:) ) / tarea 1217 1218 ssh_ib(:,:) = - ( frcv(jpr_mslp)%z3(:,:,1) - rn_pref ) * r1_grau ! equivalent ssh (inverse barometer) 1219 apr (:,:) = frcv(jpr_mslp)%z3(:,:,1) !atmospheric pressure 1220 ! 1221 CALL iom_put( "ssh_ib", ssh_ib ) !* output the inverse barometer ssh 1222 1223 ! ! ---------------------------------------- ! 1224 IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! 1225 ! ! ---------------------------------------- ! 1226 !* Restart: read in restart file 1227 IF( ln_rstart .AND. iom_varid( numror, 'ssh_ibb', ldstop = .FALSE. ) > 0 ) THEN 1228 IF(lwp) WRITE(numout,*) 'sbc_cpl: ssh_ibb read in the restart file' 1229 CALL iom_get( numror, jpdom_autoglo, 'ssh_ibb', ssh_ibb ) ! before inv. barometer ssh 1230 ELSE !* no restart: set from nit000 values 1231 IF(lwp) WRITE(numout,*) 'sbc_cpl: ssh_ibb set to nit000 values' 1232 ssh_ibb(:,:) = ssh_ib(:,:) 1233 ENDIF 1234 ENDIF 1235 ! ! ---------------------------------------- ! 1236 IF( lrst_oce ) THEN ! Write in the ocean restart file ! 1237 ! ! ---------------------------------------- ! 1238 IF(lwp) WRITE(numout,*) 1239 IF(lwp) WRITE(numout,*) 'sbc_cpl : ssh_ib written in ocean restart file at it= ', kt,' date= ', ndastp 1240 IF(lwp) WRITE(numout,*) '~~~~' 1241 CALL iom_rstput( kt, nitrst, numrow, 'ssh_ibb' , ssh_ib ) 1242 ENDIF 1243 1244 ! Update mean ssh 1245 CALL sbc_ssm_cpl( kt ) 1246 END IF 1247 ! 1248 IF( ln_sdw ) THEN ! Stokes Drift correction activated 1249 ! ! ========================= ! 1250 ! ! Stokes drift u,v ! 1251 ! ! ========================= ! 1252 IF( srcv(jpr_sdrftx)%laction .AND. srcv(jpr_sdrfty)%laction ) THEN 1253 ut0sd(:,:) = frcv(jpr_sdrftx)%z3(:,:,1) 1254 vt0sd(:,:) = frcv(jpr_sdrfty)%z3(:,:,1) 1255 ENDIF 1256 ! 1257 ! ! ========================= ! 1258 ! ! Wave mean period ! 1259 ! ! ========================= ! 1260 IF( srcv(jpr_wper)%laction ) wmp(:,:) = frcv(jpr_wper)%z3(:,:,1) 1261 ! 1262 ! ! ========================= ! 1263 ! ! Significant wave height ! 1264 ! ! ========================= ! 1265 IF( srcv(jpr_hsig)%laction ) hsw(:,:) = frcv(jpr_hsig)%z3(:,:,1) 1266 ! 1267 ! ! ========================= ! 1268 ! ! Wave peak frequency ! 1269 ! ! ========================= ! 1270 IF( srcv(jpr_wfreq)%laction ) wfreq(:,:) = frcv(jpr_wfreq)%z3(:,:,1) 1271 ! 1272 ! ! ========================= ! 1273 ! ! Vertical mixing Qiao ! 1274 ! ! ========================= ! 1275 IF( srcv(jpr_wnum)%laction .AND. ln_zdfqiao ) wnum(:,:) = frcv(jpr_wnum)%z3(:,:,1) 1276 1277 ! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode 1278 IF( (srcv(jpr_sdrftx)%laction .AND. srcv(jpr_sdrfty)%laction) .OR. srcv(jpr_wper)%laction & 1279 .OR. srcv(jpr_hsig)%laction .OR. srcv(jpr_wfreq)%laction) & 1280 CALL sbc_stokes() 1281 ENDIF 1282 ! ! ========================= ! 1283 ! ! Stress adsorbed by waves ! 1284 ! ! ========================= ! 1285 IF( srcv(jpr_tauoc)%laction .AND. ln_tauoc ) THEN 1286 tauoc_wave(:,:) = frcv(jpr_tauoc)%z3(:,:,1) 1287 ! cap the value of tauoc 1288 WHERE(tauoc_wave < 0.0 ) tauoc_wave = 1.0 1289 WHERE(tauoc_wave > 100.0 ) tauoc_wave = 1.0 1290 ENDIF 1291 ! ! ========================= ! 1292 ! ! Stress component by waves ! 1293 ! ! ========================= ! 1294 IF( srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction .AND. ln_tauw ) THEN 1295 tauw_x(:,:) = frcv(jpr_tauwx)%z3(:,:,1) 1296 tauw_y(:,:) = frcv(jpr_tauwy)%z3(:,:,1) 1297 ! cap the value of tauoc 1298 WHERE(tauw_x < -100.0 ) tauw_x = 0.0 1299 WHERE(tauw_x > 100.0 ) tauw_x = 0.0 1300 WHERE(tauw_y < -100.0 ) tauw_y = 0.0 1301 WHERE(tauw_y > 100.0 ) tauw_y = 0.0 1302 ENDIF 1303 1304 ! ! ========================= ! 1305 ! ! Wave to ocean energy ! 1306 ! ! ========================= ! 1307 IF( srcv(jpr_phioc)%laction .AND. ln_phioc ) THEN 1308 rn_crban(:,:) = 29.0 * frcv(jpr_phioc)%z3(:,:,1) 1309 WHERE( rn_crban < 0.0 ) rn_crban = 0.0 1310 WHERE( rn_crban > 1000.0 ) rn_crban = 1000.0 1311 ENDIF 998 1312 999 1313 ! Fields received by SAS when OASIS coupling … … 1067 1381 CALL ctl_stop( 'sbc_cpl_rcv: wrong definition of sn_rcv_emp%cldes' ) 1068 1382 END SELECT 1069 ELSE 1383 ELSE IF( ll_purecpl ) THEN 1070 1384 zemp(:,:) = 0._wp 1071 1385 ENDIF … … 1075 1389 IF( srcv(jpr_cal)%laction ) zemp(:,:) = zemp(:,:) - frcv(jpr_cal)%z3(:,:,1) 1076 1390 1077 IF( ln_mixcpl ) THEN ; emp(:,:) = emp(:,:) * xcplmask(:,:,0) + zemp(:,:) * zmsk(:,:) 1078 ELSE ; emp(:,:) = zemp(:,:) 1391 IF( ln_mixcpl .AND. ( srcv(jpr_oemp)%laction .OR. srcv(jpr_rain)%laction )) THEN 1392 emp(:,:) = emp(:,:) * xcplmask(:,:,0) + zemp(:,:) * zmsk(:,:) 1393 ELSE IF( ll_purecpl ) THEN ; emp(:,:) = zemp(:,:) 1079 1394 ENDIF 1080 1395 ! … … 1091 1406 ENDIF 1092 1407 ENDIF 1093 IF( ln_mixcpl ) THEN ; qns(:,:) = qns(:,:) * xcplmask(:,:,0) + zqns(:,:) * zmsk(:,:) 1094 ELSE ; qns(:,:) = zqns(:,:) 1408 IF( ln_mixcpl .AND. ( srcv(jpr_qnsoce)%laction .OR. srcv(jpr_qnsmix)%laction )) THEN 1409 qns(:,:) = qns(:,:) * xcplmask(:,:,0) + zqns(:,:) * zmsk(:,:) 1410 ELSE IF( ll_purecpl ) THEN ; qns(:,:) = zqns(:,:) 1095 1411 ENDIF 1096 1412 … … 1101 1417 ENDIF 1102 1418 IF( ln_dm2dc .AND. ln_cpl ) zqsr(:,:) = sbc_dcy( zqsr ) ! modify qsr to include the diurnal cycle 1103 IF( ln_mixcpl ) THEN ; qsr(:,:) = qsr(:,:) * xcplmask(:,:,0) + zqsr(:,:) * zmsk(:,:) 1104 ELSE ; qsr(:,:) = zqsr(:,:) 1419 IF( ln_mixcpl .AND. ( srcv(jpr_qsroce)%laction .OR. srcv(jpr_qsrmix)%laction )) THEN 1420 qsr(:,:) = qsr(:,:) * xcplmask(:,:,0) + zqsr(:,:) * zmsk(:,:) 1421 ELSE IF( ll_purecpl ) THEN ; qsr(:,:) = zqsr(:,:) 1105 1422 ENDIF 1106 1423 ! … … 1113 1430 ENDIF 1114 1431 ! 1115 CALL wrk_dealloc( jpi,jpj, ztx, zty, z msk, zemp, zqns, zqsr )1432 CALL wrk_dealloc( jpi,jpj, ztx, zty, ztx2, zty2, zmsk, zemp, zqns, zqsr ) 1116 1433 ! 1117 1434 IF( nn_timing == 1 ) CALL timing_stop('sbc_cpl_rcv') … … 1708 2025 ! 1709 2026 INTEGER :: ji, jj, jl ! dummy loop indices 2027 INTEGER :: ikchoix 1710 2028 INTEGER :: isec, info ! local integer 1711 2029 REAL(wp) :: zumax, zvmax … … 1736 2054 ! 1737 2055 SELECT CASE( sn_snd_temp%cldes) 1738 CASE( 'oce only' ) ; ztmp1(:,:) = ztmp1(:,:) + rt02056 CASE( 'oce only' ) ; ztmp1(:,:) = (ztmp1(:,:) + rt0) * tmask(:,:,1) 1739 2057 CASE( 'oce and ice' ) ; ztmp1(:,:) = ztmp1(:,:) + rt0 1740 2058 SELECT CASE( sn_snd_temp%clcat ) … … 1765 2083 ztmp1(:,:) = ztmp1(:,:) + tn_ice(:,:,jl) * a_i(:,:,jl) 1766 2084 ENDDO 2085 CASE( 'none' ) ! nothing to do 1767 2086 CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%cldes' ) 1768 2087 END SELECT … … 1889 2208 ! j+1 j -----V---F 1890 2209 ! surface velocity always sent from T point ! | 1891 ! 2210 ! [except for HadGEM3] j | T U 1892 2211 ! | | 1893 2212 ! j j-1 -I-------| … … 1901 2220 SELECT CASE( TRIM( sn_snd_crt%cldes ) ) 1902 2221 CASE( 'oce only' ) ! C-grid ==> T 1903 DO jj = 2, jpjm1 1904 DO ji = fs_2, fs_jpim1 ! vector opt. 1905 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) 1906 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) 1907 END DO 1908 END DO 2222 IF ( TRIM( sn_snd_crt%clvgrd ) == 'T' ) THEN 2223 DO jj = 2, jpjm1 2224 DO ji = fs_2, fs_jpim1 ! vector opt. 2225 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) 2226 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji,jj-1,1) ) 2227 END DO 2228 END DO 2229 ELSE 2230 ! Temporarily Changed for UKV 2231 DO jj = 2, jpjm1 2232 DO ji = 2, jpim1 2233 zotx1(ji,jj) = un(ji,jj,1) 2234 zoty1(ji,jj) = vn(ji,jj,1) 2235 END DO 2236 END DO 2237 ENDIF 1909 2238 CASE( 'weighted oce and ice' ) 1910 2239 SELECT CASE ( cp_ice_msh ) … … 1930 2259 END DO 1931 2260 CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T 1932 DO jj = 2, jpjm1 1933 DO ji = 2, jpim1 ! NO vector opt. 1934 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) 1935 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) 1936 zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) & 1937 & + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj) 1938 zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) & 1939 & + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj) 1940 END DO 1941 END DO 2261 IF ( TRIM( sn_snd_crt%clvgrd ) == 'T' ) THEN 2262 DO jj = 2, jpjm1 2263 DO ji = 2, jpim1 ! NO vector opt. 2264 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj,1) ) * zfr_l(ji,jj) & 2265 & + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) & 2266 & + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj) 2267 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji,jj-1,1) ) * zfr_l(ji,jj) & 2268 & + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) & 2269 & + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj) 2270 END DO 2271 END DO 2272 #if defined key_cice 2273 ELSE 2274 ! Temporarily Changed for HadGEM3 2275 DO jj = 2, jpjm1 2276 DO ji = 2, jpim1 ! NO vector opt. 2277 zotx1(ji,jj) = (1.0-fr_iu(ji,jj)) * un(ji,jj,1) & 2278 & + fr_iu(ji,jj) * 0.5 * ( u_ice(ji,jj-1) + u_ice(ji,jj) ) 2279 zoty1(ji,jj) = (1.0-fr_iv(ji,jj)) * vn(ji,jj,1) & 2280 & + fr_iv(ji,jj) * 0.5 * ( v_ice(ji-1,jj) + v_ice(ji,jj) ) 2281 END DO 2282 END DO 2283 #endif 2284 ENDIF 1942 2285 END SELECT 1943 2286 CALL lbc_lnk( zitx1, 'T', -1. ) ; CALL lbc_lnk( zity1, 'T', -1. ) … … 1984 2327 IF( TRIM( sn_snd_crt%clvor ) == 'eastward-northward' ) THEN ! Rotation of the components 1985 2328 ! ! Ocean component 1986 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 1987 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 1988 zotx1(:,:) = ztmp1(:,:) ! overwrite the components 1989 zoty1(:,:) = ztmp2(:,:) 1990 IF( ssnd(jps_ivx1)%laction ) THEN ! Ice component 1991 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 1992 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 1993 zitx1(:,:) = ztmp1(:,:) ! overwrite the components 1994 zity1(:,:) = ztmp2(:,:) 1995 ENDIF 2329 IF ( TRIM( sn_snd_crt%clvgrd ) == 'T' ) THEN 2330 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 2331 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 2332 zotx1(:,:) = ztmp1(:,:) ! overwrite the components 2333 zoty1(:,:) = ztmp2(:,:) 2334 IF( ssnd(jps_ivx1)%laction ) THEN ! Ice component 2335 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 2336 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 2337 zitx1(:,:) = ztmp1(:,:) ! overwrite the components 2338 zity1(:,:) = ztmp2(:,:) 2339 ENDIF 2340 ELSE 2341 ! Temporary code for HadGEM3 - will be removed eventually. 2342 ! Only applies when we want uvel on U grid and vvel on V grid 2343 ! Rotate U and V onto geographic grid before sending. 2344 2345 DO jj=2,jpjm1 2346 DO ji=2,jpim1 2347 ztmp1(ji,jj)=0.25*vmask(ji,jj,1) & 2348 *(zotx1(ji,jj)+zotx1(ji-1,jj) & 2349 +zotx1(ji,jj+1)+zotx1(ji-1,jj+1)) 2350 ztmp2(ji,jj)=0.25*umask(ji,jj,1) & 2351 *(zoty1(ji,jj)+zoty1(ji+1,jj) & 2352 +zoty1(ji,jj-1)+zoty1(ji+1,jj-1)) 2353 ENDDO 2354 ENDDO 2355 2356 ! Ensure any N fold and wrap columns are updated 2357 CALL lbc_lnk(ztmp1, 'V', -1.0) 2358 CALL lbc_lnk(ztmp2, 'U', -1.0) 2359 2360 ikchoix = -1 2361 CALL repcmo(zotx1,ztmp2,ztmp1,zoty1,zotx1,zoty1,ikchoix) 2362 ENDIF 1996 2363 ENDIF 1997 2364 ! … … 2019 2386 ENDIF 2020 2387 ! 2388 ! ! ------------------------- ! 2389 ! ! Surface current to waves ! 2390 ! ! ------------------------- ! 2391 IF( ssnd(jps_ocxw)%laction .OR. ssnd(jps_ocyw)%laction ) THEN 2392 ! 2393 ! j+1 j -----V---F 2394 ! surface velocity always sent from T point ! | 2395 ! j | T U 2396 ! | | 2397 ! j j-1 -I-------| 2398 ! (for I) | | 2399 ! i-1 i i 2400 ! i i+1 (for I) 2401 SELECT CASE( TRIM( sn_snd_crtw%cldes ) ) 2402 CASE( 'oce only' ) ! C-grid ==> T 2403 DO jj = 2, jpjm1 2404 DO ji = fs_2, fs_jpim1 ! vector opt. 2405 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) 2406 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji , jj-1,1) ) 2407 END DO 2408 END DO 2409 CASE( 'weighted oce and ice' ) 2410 SELECT CASE ( cp_ice_msh ) 2411 CASE( 'C' ) ! Ocean and Ice on C-grid ==> T 2412 DO jj = 2, jpjm1 2413 DO ji = fs_2, fs_jpim1 ! vector opt. 2414 zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(ji,jj) 2415 zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + vn (ji ,jj-1,1) ) * zfr_l(ji,jj) 2416 zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj) 2417 zity1(ji,jj) = 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj) 2418 END DO 2419 END DO 2420 CASE( 'I' ) ! Ocean on C grid, Ice on I-point (B-grid) ==> T 2421 DO jj = 2, jpjm1 2422 DO ji = 2, jpim1 ! NO vector opt. 2423 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) 2424 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) 2425 zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1) & 2426 & + u_ice(ji+1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj) 2427 zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1) & 2428 & + v_ice(ji+1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj) 2429 END DO 2430 END DO 2431 CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T 2432 DO jj = 2, jpjm1 2433 DO ji = 2, jpim1 ! NO vector opt. 2434 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) 2435 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) 2436 zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) & 2437 & + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj) 2438 zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) & 2439 & + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj) 2440 END DO 2441 END DO 2442 END SELECT 2443 CALL lbc_lnk( zitx1, 'T', -1. ) ; CALL lbc_lnk( zity1, 'T', -1. ) 2444 CASE( 'mixed oce-ice' ) 2445 SELECT CASE ( cp_ice_msh ) 2446 CASE( 'C' ) ! Ocean and Ice on C-grid ==> T 2447 DO jj = 2, jpjm1 2448 DO ji = fs_2, fs_jpim1 ! vector opt. 2449 zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) & 2450 & + 0.5 * ( u_ice(ji,jj ) + u_ice(ji-1,jj ) ) * fr_i(ji,jj) 2451 zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) & 2452 & + 0.5 * ( v_ice(ji,jj ) + v_ice(ji ,jj-1 ) ) * fr_i(ji,jj) 2453 END DO 2454 END DO 2455 CASE( 'I' ) ! Ocean on C grid, Ice on I-point (B-grid) ==> T 2456 DO jj = 2, jpjm1 2457 DO ji = 2, jpim1 ! NO vector opt. 2458 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) & 2459 & + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1) & 2460 & + u_ice(ji+1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj) 2461 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) & 2462 & + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1) & 2463 & + v_ice(ji+1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj) 2464 END DO 2465 END DO 2466 CASE( 'F' ) ! Ocean on C grid, Ice on F-point (B-grid) ==> T 2467 DO jj = 2, jpjm1 2468 DO ji = 2, jpim1 ! NO vector opt. 2469 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) * zfr_l(ji,jj) & 2470 & + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1) & 2471 & + u_ice(ji-1,jj ) + u_ice(ji,jj ) ) * fr_i(ji,jj) 2472 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji ,jj-1,1) ) * zfr_l(ji,jj) & 2473 & + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1) & 2474 & + v_ice(ji-1,jj ) + v_ice(ji,jj ) ) * fr_i(ji,jj) 2475 END DO 2476 END DO 2477 END SELECT 2478 END SELECT 2479 CALL lbc_lnk( zotx1, ssnd(jps_ocxw)%clgrid, -1. ) ; CALL lbc_lnk( zoty1, ssnd(jps_ocyw)%clgrid, -1. ) 2480 ! 2481 ! 2482 IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) THEN ! Rotation of the components 2483 ! ! Ocean component 2484 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->e', ztmp1 ) ! 1st component 2485 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 2486 zotx1(:,:) = ztmp1(:,:) ! overwrite the components 2487 zoty1(:,:) = ztmp2(:,:) 2488 IF( ssnd(jps_ivx1)%laction ) THEN ! Ice component 2489 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 2490 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 2491 zitx1(:,:) = ztmp1(:,:) ! overwrite the components 2492 zity1(:,:) = ztmp2(:,:) 2493 ENDIF 2494 ENDIF 2495 ! 2496 ! ! spherical coordinates to cartesian -> 2 components to 3 components 2497 ! IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN 2498 ! ztmp1(:,:) = zotx1(:,:) ! ocean currents 2499 ! ztmp2(:,:) = zoty1(:,:) 2500 ! CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 ) 2501 ! ! 2502 ! IF( ssnd(jps_ivx1)%laction ) THEN ! ice velocities 2503 ! ztmp1(:,:) = zitx1(:,:) 2504 ! ztmp1(:,:) = zity1(:,:) 2505 ! CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 ) 2506 ! ENDIF 2507 ! ENDIF 2508 ! 2509 IF( ssnd(jps_ocxw)%laction ) CALL cpl_snd( jps_ocxw, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info ) ! ocean x current 1st grid 2510 IF( ssnd(jps_ocyw)%laction ) CALL cpl_snd( jps_ocyw, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info ) ! ocean y current 1st grid 2511 ! 2512 ENDIF 2513 ! 2514 IF( ssnd(jps_ficet)%laction ) THEN 2515 CALL cpl_snd( jps_ficet, isec, RESHAPE ( fr_i, (/jpi,jpj,1/) ), info ) 2516 END IF 2517 ! ! ------------------------- ! 2518 ! ! Water levels to waves ! 2519 ! ! ------------------------- ! 2520 IF( ssnd(jps_wlev)%laction ) THEN 2521 IF( ln_apr_dyn ) THEN 2522 IF( kt /= nit000 ) THEN 2523 ztmp1(:,:) = sshb(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 2524 ELSE 2525 ztmp1(:,:) = sshb(:,:) 2526 ENDIF 2527 ELSE 2528 ztmp1(:,:) = sshn(:,:) 2529 ENDIF 2530 CALL cpl_snd( jps_wlev , isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info ) 2531 END IF 2021 2532 ! 2022 2533 ! Fields sent by OPA to SAS when doing OPA<->SAS coupling -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90
r8918 r11277 20 20 USE iom ! IOM library 21 21 USE in_out_manager ! I/O manager 22 USE sbcwave ! wave physics 22 23 USE lib_mpp ! distribued memory computing library 23 24 USE lbclnk ! ocean lateral boundary conditions (or mpp link) … … 87 88 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 88 89 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 90 REAL(wp) :: totwind ! UKMO SHELF: Module of wind speed 89 91 REAL(wp) :: ztx, zty, zmod, zcoef ! temporary variables 90 92 REAL :: cs ! UKMO SHELF: Friction co-efficient at surface … … 115 117 IF(lwm) WRITE ( numond, namsbc_flx ) 116 118 ! 119 IF(lwp) THEN ! Namelist print 120 WRITE(numout,*) 121 WRITE(numout,*) 'sbc_flx : Flux forcing' 122 WRITE(numout,*) '~~~~~~~~~~~' 123 WRITE(numout,*) ' Namelist namsbc_flx : shelf seas configuration (force with winds instead of momentum)' 124 WRITE(numout,*) ' shelf seas configuration ln_shelf_flx = ', ln_shelf_flx 125 WRITE(numout,*) ' relative wind speed ln_rel_wind = ', ln_rel_wind 126 WRITE(numout,*) ' wind multiplication factor rn_wfac = ', rn_wfac 127 ENDIF 117 128 ! ! check: do we plan to use ln_dm2dc with non-daily forcing? 118 129 IF( ln_dm2dc .AND. sn_qsr%nfreqh /= 24 ) & … … 210 221 END DO 211 222 END DO 223 ! ! add modification due to drag coefficient read from wave forcing 224 ! ! this code is inefficient but put here to allow merging with another UKMO branch 225 IF( ln_shelf_flx ) THEN 226 ! calculate first the wind module, as it will be used later 227 DO jj = 2, jpjm1 228 DO ji = fs_2, fs_jpim1 ! vect. opt. 229 ztx = zwnd_i(ji-1,jj ) + zwnd_i(ji,jj) 230 zty = zwnd_j(ji ,jj-1) + zwnd_j(ji,jj) 231 wndm(ji,jj) = 0.5 * SQRT( ztx * ztx + zty * zty ) 232 END DO 233 END DO 234 CALL lbc_lnk( wndm(:,:), 'T', 1. ) 235 236 IF( ln_cdgw .AND. nn_drag == jp_std ) THEN 237 IF( cpl_wdrag ) THEN 238 ! reset utau and vtau to the wind components: the momentum will 239 ! be calculated from the coupled value of the drag coefficient 240 DO jj = 1, jpj 241 DO ji = 1, jpi 242 utau(ji,jj) = zwnd_i(ji,jj) 243 vtau(ji,jj) = zwnd_j(ji,jj) 244 END DO 245 END DO 246 ELSE 247 DO jj = 1, jpjm1 248 DO ji = 1, jpim1 249 utau(ji,jj) = zrhoa * 0.5 * ( cdn_wave(ji,jj) + cdn_wave(ji+1,jj) ) * zwnd_i(ji,jj) * & 250 0.5 * ( wndm(ji,jj) + wndm(ji+1,jj) ) 251 vtau(ji,jj) = zrhoa * 0.5 * ( cdn_wave(ji,jj) + cdn_wave(ji,jj+1) ) * zwnd_j(ji,jj) * & 252 0.5 * ( wndm(ji,jj) + wndm(ji,jj+1) ) 253 END DO 254 END DO 255 CALL lbc_lnk_multi( utau(:,:), 'U', -1., vtau(:,:), 'V', -1. ) 256 ENDIF 257 ELSE IF( nn_drag == jp_const ) THEN 258 DO jj = 1, jpjm1 259 DO ji = 1, jpim1 260 utau(ji,jj) = zrhoa * zcdrag * zwnd_i(ji,jj) * 0.5 * ( wndm(ji,jj) + wndm(ji+1,jj) ) 261 vtau(ji,jj) = zrhoa * zcdrag * zwnd_j(ji,jj) * 0.5 * ( wndm(ji,jj) + wndm(ji,jj+1) ) 262 END DO 263 END DO 264 CALL lbc_lnk_multi( utau(:,:), 'U', -1., vtau(:,:), 'V', -1. ) 265 ENDIF 266 ENDIF 212 267 ! ! add to qns the heat due to e-p 213 268 qns(:,:) = qns(:,:) - emp(:,:) * sst_m(:,:) * rcp ! mass flux is at SST … … 230 285 zmod = 0.5 * SQRT( ztx * ztx + zty * zty ) 231 286 taum(ji,jj) = zmod 287 IF ( .NOT. (ln_shelf_flx .AND. ln_cpl)) THEN 232 288 wndm(ji,jj) = SQRT( zmod * zcoef ) 289 ENDIF 233 290 END DO 234 291 END DO -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90
r8058 r11277 284 284 CALL wrk_dealloc( jpi,jpj, ztmp1, ztmp2 ) 285 285 ! 286 ! In coupled mode get extra fields from CICE for passing back to atmosphere 287 IF ( ksbc == jp_purecpl ) CALL cice_sbc_hadgam(nit000) 288 ! 286 289 IF( nn_timing == 1 ) CALL timing_stop('cice_sbc_init') 287 290 ! … … 708 711 IF( nn_timing == 1 ) CALL timing_start('cice_sbc_hadgam') 709 712 ! 710 IF( kt == nit000 ) THEN711 IF(lwp) WRITE(numout,*)'cice_sbc_hadgam'712 IF( sbc_cpl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' )713 ENDIF714 715 713 ! ! =========================== ! 716 714 ! ! Prepare Coupling fields ! -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
r8058 r11277 88 88 NAMELIST/namsbc/ nn_fsbc , ln_ana , ln_flx, ln_blk_clio, ln_blk_core, ln_mixcpl, & 89 89 & ln_blk_mfs, ln_apr_dyn, nn_ice, nn_ice_embd, ln_dm2dc , ln_rnf , & 90 & ln_ssr , nn_isf , nn_fwb, ln_ cdgw , ln_wave , ln_sdw , &91 & nn_lsm , nn_limflx , nn_components, ln_cpl90 & ln_ssr , nn_isf , nn_fwb, ln_wave, nn_lsm , nn_limflx, & 91 nn_components, ln_cpl , ln_wavcpl, nn_drag 92 92 INTEGER :: ios 93 93 INTEGER :: ierr, ierr0, ierr1, ierr2, ierr3, jpm 94 LOGICAL :: ll_purecpl95 94 !!---------------------------------------------------------------------- 96 95 … … 131 130 WRITE(numout,*) ' MFS bulk formulation ln_blk_mfs = ', ln_blk_mfs 132 131 WRITE(numout,*) ' ocean-atmosphere coupled formulation ln_cpl = ', ln_cpl 133 WRITE(numout,*) ' forced-coupled mixed formulation ln_mixcpl = ', ln_mixcpl 132 WRITE(numout,*) ' forced-coupled atm mixed formulation ln_mixcpl = ', ln_mixcpl 133 WRITE(numout,*) ' forced-coupled wav mixed formulation ln_wavcpl = ', ln_wavcpl 134 WRITE(numout,*) ' wave physics ln_wave = ', ln_wave 134 135 WRITE(numout,*) ' OASIS coupling (with atm or sas) lk_oasis = ', lk_oasis 135 136 WRITE(numout,*) ' components of your executable nn_components = ', nn_components 136 137 WRITE(numout,*) ' Multicategory heat flux formulation (LIM3) nn_limflx = ', nn_limflx 138 WRITE(numout,*) ' momentum formulation nn_drag = ', nn_drag 137 139 WRITE(numout,*) ' Misc. options of sbc : ' 138 140 WRITE(numout,*) ' Patm gradient added in ocean & ice Eqs. ln_apr_dyn = ', ln_apr_dyn … … 164 166 IF ( nn_components == jp_iam_opa .AND. ln_cpl ) & 165 167 & CALL ctl_stop( 'STOP', 'sbc_init : OPA-SAS coupled via OASIS, but ln_cpl = T in OPA' ) 166 IF ( nn_components == jp_iam_opa .AND. ln_mixcpl ) &167 & CALL ctl_stop( 'STOP', 'sbc_init : OPA-SAS coupled via OASIS, but ln_mixcpl = T in OPA' )168 IF ( nn_components == jp_iam_opa .AND. ( ln_mixcpl .OR. ln_wavcpl) ) & 169 & CALL ctl_stop( 'STOP', 'sbc_init : OPA-SAS coupled via OASIS, but ln_mixcpl or ln_wavcpl = T in OPA' ) 168 170 IF ( ln_cpl .AND. .NOT. lk_oasis ) & 169 171 & CALL ctl_stop( 'STOP', 'sbc_init : OASIS-coupled atmosphere model, but key_oasis3 disabled' ) 170 IF( ln_mixcpl .AND. .NOT. lk_oasis ) &171 & CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl ) requires the cpp key key_oasis3' )172 IF( ln_mixcpl .AND. .NOT. ln_cpl ) &173 & CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl ) requires ln_cpl = T' )174 IF( ln_mixcpl .AND. nn_components /= jp_iam_nemo ) &175 & CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl ) is not yet working with sas-opa coupling via oasis' )172 IF( ( ln_mixcpl .OR. ln_wavcpl ) .AND. .NOT. lk_oasis ) & 173 & CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl or ln_wavcpl) requires the cpp key key_oasis3' ) 174 IF( ( ln_mixcpl .OR. ln_wavcpl ) .AND. .NOT. ln_cpl ) & 175 & CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl or ln_wavcpl) requires ln_cpl = T' ) 176 IF( ( ln_mixcpl .OR. ln_wavcpl ) .AND. nn_components /= jp_iam_nemo ) & 177 & CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl or ln_wavcpl) is not yet working with sas-opa coupling via oasis' ) 176 178 177 179 ! ! allocate sbc arrays … … 214 216 IF( ln_dm2dc .AND. .NOT.( ln_flx .OR. ln_blk_core ) .AND. nn_components /= jp_iam_opa ) & 215 217 & CALL ctl_stop( 'diurnal cycle into qsr field from daily values requires a flux or core-bulk formulation' ) 216 217 IF ( ln_wave ) THEN218 !Activated wave module but neither drag nor stokes drift activated219 IF ( .NOT.(ln_cdgw .OR. ln_sdw) ) THEN220 CALL ctl_warn( 'Ask for wave coupling but nor drag coefficient (ln_cdgw=F) neither stokes drift activated (ln_sdw=F)' )221 !drag coefficient read from wave model definable only with mfs bulk formulae and core222 ELSEIF (ln_cdgw .AND. .NOT.(ln_blk_mfs .OR. ln_blk_core) ) THEN223 CALL ctl_stop( 'drag coefficient read from wave model definable only with mfs bulk formulae and core')224 ENDIF225 ELSE226 IF ( ln_cdgw .OR. ln_sdw ) &227 & CALL ctl_stop('Not Activated Wave Module (ln_wave=F) but &228 & asked coupling with drag coefficient (ln_cdgw =T) or Stokes drift (ln_sdw=T) ')229 ENDIF230 218 ! ! Choice of the Surface Boudary Condition (set nsbc) 231 ll_purecpl = ln_cpl .AND. .NOT. ln_mixcpl 219 ll_purecpl = ln_cpl .AND. .NOT. ln_mixcpl .AND. .NOT. ln_wavcpl 232 220 ! 233 221 icpt = 0 … … 261 249 IF( nsbc == jp_mfs ) WRITE(numout,*) ' MFS Bulk formulation' 262 250 IF( nsbc == jp_none ) WRITE(numout,*) ' OPA coupled to SAS via oasis' 263 IF( ln_mixcpl ) WRITE(numout,*) ' + forced-coupled mixed formulation' 251 IF( ln_mixcpl ) WRITE(numout,*) ' + forced-coupled mixed atm formulation' 252 IF( ln_wavcpl ) WRITE(numout,*) ' + forced-coupled mixed wav formulation' 264 253 IF( nn_components/= jp_iam_nemo ) & 265 254 & WRITE(numout,*) ' + OASIS coupled SAS' 266 255 ENDIF 267 256 ! 268 IF( lk_oasis ) CALL sbc_cpl_init (nn_ice) ! OASIS initialisation. must be done before: (1) first time step 269 ! ! (2) the use of nn_fsbc 270 257 IF( lk_oasis ) THEN 258 IF( sbc_cpl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' ) 259 CALL sbc_cpl_init (nn_ice) ! OASIS initialisation. must be done before: (1) first time step 260 ! (2) the use of nn_fsbc 261 ENDIF 271 262 ! nn_fsbc initialization if OPA-SAS coupling via OASIS 272 263 ! sas model time step has to be declared in OASIS (mandatory) -> nn_fsbc has to be modified accordingly … … 305 296 306 297 IF( nn_ice == 4 ) CALL cice_sbc_init( nsbc ) ! CICE initialisation 298 ! 299 IF( ln_wave ) CALL sbc_wave_init ! surface wave initialisation 300 ! 307 301 308 302 END SUBROUTINE sbc_init … … 325 319 !! - updte the ice fraction : fr_i 326 320 !!---------------------------------------------------------------------- 321 USE bdydta, ONLY: bdy_dta 322 ! 327 323 INTEGER, INTENT(in) :: kt ! ocean time step 328 324 !!--------------------------------------------------------------------- … … 345 341 ! ! ---------------------------------------- ! 346 342 ! 347 IF ( .NOT. lk_bdy ) then 348 IF( ln_apr_dyn ) CALL sbc_apr( kt ) ! atmospheric pressure provided at kt+0.5*nn_fsbc 349 ENDIF 343 344 IF( ln_apr_dyn ) CALL sbc_apr( kt ) ! atmospheric pressure provided at kt+0.5*nn_fsbc 350 345 ! (caution called before sbc_ssm) 351 346 ! … … 382 377 END SELECT 383 378 384 IF( ln_mixcpl ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice ) ! forced-coupled mixed formulation after forcing 385 379 IF( ln_mixcpl .OR. ln_wavcpl ) CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice ) ! forced-coupled mixed formulation after forcing 380 381 IF( ln_wave .AND. (ln_tauoc .OR. ln_tauw) ) CALL sbc_stress( ) ! Wave stress update 382 IF( lk_bdy ) CALL bdy_dta ( kt, time_offset=+1 ) ! update dynamic & tracer data at open boundaries 383 ! (caution called after sbc_ssm[_cpl] and before ice) 386 384 387 385 ! !== Misc. Options ==! -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90
r8058 r11277 82 82 ALLOCATE( rnfmsk(jpi,jpj) , rnfmsk_z(jpk) , & 83 83 & h_rnf (jpi,jpj) , nk_rnf (jpi,jpj) , & 84 & rnf_tsc_b(jpi,jpj,jpts) , rnf_tsc (jpi,jpj,jpts) , STAT=sbc_rnf_alloc ) 84 & rnf_tsc_b(jpi,jpj,jpts) , rnf_tsc (jpi,jpj,jpts) , & 85 & xrnfmask(jpi,jpj,1) , STAT=sbc_rnf_alloc ) 85 86 ! 86 87 IF( lk_mpp ) CALL mpp_sum ( sbc_rnf_alloc ) … … 128 129 IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN 129 130 ! 130 IF( .NOT. l_rnfcpl ) rnf(:,:) = rn_rfact * ( sf_rnf(1)%fnow(:,:,1) ) ! updated runoff value at time step kt 131 IF( .NOT. l_rnfcpl ) & 132 rnf(:,:) = rnf(:,:) * (1. - xrnfmask(:,:,1)) + rn_rfact * sf_rnf(1)%fnow(:,:,1) * xrnfmask(:,:,1) ! updated runoff value at time step kt 131 133 ! 132 134 ! ! set temperature & salinity content of runoffs … … 442 444 ENDIF 443 445 ! 446 xrnfmask(:,:,:) = 1. ! default value for points using river forcing 447 IF (ln_usernfmask) THEN 448 IF(lwp) WRITE(numout,*) 449 IF(lwp) WRITE(numout,*) ' runoff mask read in a file' 450 CALL iom_open( 'rnfmask', inum ) 451 CALL iom_get( inum, jpdom_data, 'rnfmask', xrnfmask(:,:,1), 1) 452 CALL iom_close( inum ) 453 ENDIF 454 ! 444 455 rnf(:,:) = 0._wp ! runoff initialisation 445 456 rnf_tsc(:,:,:) = 0._wp ! runoffs temperature & salinty contents initilisation -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssm.F90
r8058 r11277 26 26 27 27 PUBLIC sbc_ssm ! routine called by step.F90 28 PUBLIC sbc_ssm_cpl ! routine called by sbccpl.F90 28 29 PUBLIC sbc_ssm_init ! routine called by sbcmod.F90 29 30 … … 77 78 sss_m(:,:) = zts(:,:,jp_sal) 78 79 ! ! removed inverse barometer ssh when Patm forcing is used (for sea-ice dynamics) 79 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 80 ELSE ; ssh_m(:,:) = sshn(:,:) 80 IF( .NOT. cpl_mslp ) THEN 81 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 82 ELSE ; ssh_m(:,:) = sshn(:,:) 83 ENDIF 81 84 ENDIF 82 85 ! … … 99 102 sss_m(:,:) = zcoef * zts(:,:,jp_sal) 100 103 ! ! removed inverse barometer ssh when Patm forcing is used (for sea-ice dynamics) 101 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = zcoef * ( sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) ) 102 ELSE ; ssh_m(:,:) = zcoef * sshn(:,:) 104 IF( .NOT. cpl_mslp ) THEN 105 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = zcoef * ( sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) ) 106 ELSE ; ssh_m(:,:) = zcoef * sshn(:,:) 107 ENDIF 103 108 ENDIF 104 109 ! … … 113 118 sst_m(:,:) = 0.e0 114 119 sss_m(:,:) = 0.e0 115 ssh_m(:,:) = 0.e0120 IF( .NOT. cpl_mslp ) ssh_m(:,:) = 0.e0 116 121 IF( lk_vvl ) e3t_m(:,:) = 0.e0 117 122 frq_m(:,:) = 0.e0 … … 127 132 sss_m(:,:) = sss_m(:,:) + zts(:,:,jp_sal) 128 133 ! ! removed inverse barometer ssh when Patm forcing is used (for sea-ice dynamics) 129 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = ssh_m(:,:) + sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 130 ELSE ; ssh_m(:,:) = ssh_m(:,:) + sshn(:,:) 134 IF( .NOT. cpl_mslp ) THEN 135 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = ssh_m(:,:) + sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 136 ELSE ; ssh_m(:,:) = ssh_m(:,:) + sshn(:,:) 137 ENDIF 131 138 ENDIF 132 139 ! … … 143 150 ssu_m(:,:) = ssu_m(:,:) * zcoef ! mean suface current [m/s] 144 151 ssv_m(:,:) = ssv_m(:,:) * zcoef ! 145 ssh_m(:,:) = ssh_m(:,:) * zcoef ! mean SSH [m]152 IF( .NOT. cpl_mslp ) ssh_m(:,:) = ssh_m(:,:) * zcoef ! mean SSH [m] 146 153 IF( lk_vvl ) e3t_m(:,:) = fse3t_m(:,:) * zcoef ! mean vertical scale factor [m] 147 154 frq_m(:,:) = frq_m(:,:) * zcoef ! mean fraction of solar net radiation absorbed in the 1st T level [-] … … 161 168 CALL iom_rstput( kt, nitrst, numrow, 'sst_m' , sst_m ) 162 169 CALL iom_rstput( kt, nitrst, numrow, 'sss_m' , sss_m ) 163 CALL iom_rstput( kt, nitrst, numrow, 'ssh_m' , ssh_m )170 IF( .NOT. cpl_mslp ) CALL iom_rstput( kt, nitrst, numrow, 'ssh_m' , ssh_m ) 164 171 IF( lk_vvl ) CALL iom_rstput( kt, nitrst, numrow, 'e3t_m' , e3t_m ) 165 172 CALL iom_rstput( kt, nitrst, numrow, 'frq_m' , frq_m ) … … 174 181 CALL iom_put( 'sst_m', sst_m ) 175 182 CALL iom_put( 'sss_m', sss_m ) 176 CALL iom_put( 'ssh_m', ssh_m )183 IF( .NOT. cpl_mslp ) CALL iom_put( 'ssh_m', ssh_m ) 177 184 IF( lk_vvl ) CALL iom_put( 'e3t_m', e3t_m ) 178 185 CALL iom_put( 'frq_m', frq_m ) … … 180 187 ! 181 188 END SUBROUTINE sbc_ssm 189 190 SUBROUTINE sbc_ssm_cpl( kt ) 191 !!--------------------------------------------------------------------- 192 !! *** ROUTINE sbc_ssm_cpl *** 193 !! 194 !! ** Purpose : provide ocean surface variable to sea-surface boundary 195 !! condition computation when pressure is read from coupling 196 !! 197 !! ** Method : The inverse barometer ssh (i.e. ssh associated with Patm) 198 !! is added to ssh_m when ln_apr_dyn = T. Required for sea-ice dynamics. 199 !!--------------------------------------------------------------------- 200 INTEGER, INTENT(in) :: kt ! ocean time step 201 ! 202 REAL(wp) :: zcoef ! local scalar 203 !!--------------------------------------------------------------------- 204 ! 205 IF( nn_fsbc == 1 ) THEN ! Instantaneous surface fields ! 206 ! ! ---------------------------------------- ! 207 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 208 ELSE ; ssh_m(:,:) = sshn(:,:) 209 ENDIF 210 ELSE 211 ! ! ----------------------------------------------- ! 212 IF( kt == nit000 .AND. .NOT. l_ssm_mean ) THEN ! Initialisation: 1st time-step, no input means ! 213 ! ! ----------------------------------------------- ! 214 IF(lwp) WRITE(numout,*) 215 IF(lwp) WRITE(numout,*) '~~~~~~~ mean ssh field initialised to instantaneous values' 216 zcoef = REAL( nn_fsbc - 1, wp ) 217 zcoef = REAL( nn_fsbc - 1, wp ) 218 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = zcoef * ( sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) ) 219 ELSE ; ssh_m(:,:) = zcoef * sshn(:,:) 220 ENDIF 221 ! ! ---------------------------------------- ! 222 ELSEIF( MOD( kt - 2 , nn_fsbc ) == 0 ) THEN ! Initialisation: New mean computation ! 223 ! ! ---------------------------------------- ! 224 ssh_m(:,:) = 0.e0 225 ENDIF 226 227 IF( ln_apr_dyn ) THEN ; ssh_m(:,:) = ssh_m(:,:) + sshn(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) ) 228 ELSE ; ssh_m(:,:) = ssh_m(:,:) + sshn(:,:) 229 ENDIF 230 ! ! ---------------------------------------- ! 231 IF( MOD( kt - 1 , nn_fsbc ) == 0 ) THEN ! Mean value at each nn_fsbc time-step ! 232 ! ! ---------------------------------------- ! 233 zcoef = 1. / REAL( nn_fsbc, wp ) 234 ssh_m(:,:) = ssh_m(:,:) * zcoef ! mean SSH [m] 235 ENDIF 236 ! ! ---------------------------------------- ! 237 IF( lrst_oce ) THEN ! Write in the ocean restart file ! 238 ! ! ---------------------------------------- ! 239 IF(lwp) WRITE(numout,*) 240 IF(lwp) WRITE(numout,*) 'sbc_ssm_cpl : ssh mean field written in ocean restart file ', & 241 & 'at it= ', kt,' date= ', ndastp 242 IF(lwp) WRITE(numout,*) '~~~~~~~' 243 CALL iom_rstput( kt, nitrst, numrow, 'ssh_m' , ssh_m ) 244 ENDIF 245 ENDIF 246 ! 247 IF( MOD( kt - 1 , nn_fsbc ) == 0 ) THEN ! Mean value at each nn_fsbc time-step ! 248 CALL iom_put( 'ssh_m', ssh_m ) 249 ENDIF 250 ! 251 END SUBROUTINE sbc_ssm_cpl 182 252 183 253 SUBROUTINE sbc_ssm_init -
branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
r8058 r11277 4 4 !! Wave module 5 5 !!====================================================================== 6 !! History : 3.3.1 ! 2011-09 (Adani M) Original code: Drag Coefficient 7 !! : 3.4 ! 2012-10 (Adani M) Stokes Drift 6 !! History : 3.3 ! 2011-09 (Adani M) Original code: Drag Coefficient 7 !! : 3.4 ! 2012-10 (Adani M) Stokes Drift 8 !! 3.6 ! 2014-09 (Clementi E, Oddo P)New Stokes Drift Computation 9 !! - ! 2016-12 (G. Madec, E. Clementi) update Stoke drift computation 10 !! + add sbc_wave_ini routine 8 11 !!---------------------------------------------------------------------- 9 USE iom ! I/O manager library 10 USE in_out_manager ! I/O manager 11 USE lib_mpp ! distribued memory computing library 12 USE fldread ! read input fields 13 USE oce 14 USE sbc_oce ! Surface boundary condition: ocean fields 15 USE domvvl 16 17 12 18 13 !!---------------------------------------------------------------------- 19 !! sbc_wave : read drag coefficient from wave model in netcdf files 14 !! sbc_stokes : calculate 3D Stokes-drift velocities 15 !! sbc_wave : wave data from wave model in netcdf files 16 !! sbc_wave_init : initialisation fo surface waves 20 17 !!---------------------------------------------------------------------- 18 USE oce ! ocean variables 19 USE sbc_oce ! Surface boundary condition: ocean fields 20 USE bdy_oce ! open boundary condition variables 21 USE domvvl ! domain: variable volume layers 22 ! 23 USE iom ! I/O manager library 24 USE in_out_manager ! I/O manager 25 USE lib_mpp ! distribued memory computing library 26 USE fldread ! read input fields 27 USE wrk_nemo ! 28 USE phycst ! physical constants 21 29 22 30 IMPLICIT NONE 23 31 PRIVATE 24 32 25 PUBLIC sbc_wave ! routine called in sbc_blk_core or sbc_blk_mfs 33 PUBLIC sbc_stokes ! routine called in sbccpl 34 PUBLIC sbc_stress ! routine called in sbcmod 35 PUBLIC sbc_wave ! routine called in sbcmod 36 PUBLIC sbc_wave_init ! routine called in sbcmod 26 37 27 INTEGER , PARAMETER :: jpfld = 3 ! maximum number of files to read for srokes drift 28 INTEGER , PARAMETER :: jp_usd = 1 ! index of stokes drift (i-component) (m/s) at T-point 29 INTEGER , PARAMETER :: jp_vsd = 2 ! index of stokes drift (j-component) (m/s) at T-point 30 INTEGER , PARAMETER :: jp_wn = 3 ! index of wave number (1/m) at T-point 31 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 32 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 33 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:) :: cdn_wave 34 REAL(wp),ALLOCATABLE,DIMENSION (:,:) :: usd2d,vsd2d,uwavenum,vwavenum 35 REAL(wp),PUBLIC,ALLOCATABLE,DIMENSION (:,:,:) :: usd3d,vsd3d,wsd3d 36 37 !! * Substitutions 38 # include "domzgr_substitute.h90" 38 ! Variables checking if the wave parameters are coupled (if not, they are read from file) 39 LOGICAL, PUBLIC :: cpl_hsig = .FALSE. 40 LOGICAL, PUBLIC :: cpl_phioc = .FALSE. 41 LOGICAL, PUBLIC :: cpl_sdrft = .FALSE. 42 LOGICAL, PUBLIC :: cpl_wper = .FALSE. 43 LOGICAL, PUBLIC :: cpl_wfreq = .FALSE. 44 LOGICAL, PUBLIC :: cpl_wnum = .FALSE. 45 LOGICAL, PUBLIC :: cpl_tauoc = .FALSE. 46 LOGICAL, PUBLIC :: cpl_tauw = .FALSE. 47 LOGICAL, PUBLIC :: cpl_wdrag = .FALSE. 48 49 INTEGER :: nn_sdrift ! type of parameterization to calculate vertical Stokes drift 50 INTEGER, PARAMETER :: jp_breivik = 0 ! Breivik 2015: v_z=v_0*[exp(2*k*z)/(1-8*k*z)] 51 INTEGER, PARAMETER :: jp_phillips = 1 ! Phillips: v_z=v_o*[exp(2*k*z)-beta*sqrt(-2*k*pi*z)*erfc(sqrt(-2*k*z))] 52 INTEGER, PARAMETER :: jp_peakph = 2 ! Phillips using the peak wave number read from wave model instead of the inverse depth scale 53 54 INTEGER :: jpfld ! number of files to read for stokes drift 55 INTEGER :: jp_usd ! index of stokes drift (i-component) (m/s) at T-point 56 INTEGER :: jp_vsd ! index of stokes drift (j-component) (m/s) at T-point 57 INTEGER :: jp_hsw ! index of significant wave hight (m) at T-point 58 INTEGER :: jp_wmp ! index of mean wave period (s) at T-point 59 INTEGER :: jp_wfr ! index of wave peak frequency (s^-1) at T-point 60 61 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_cd ! structure of input fields (file informations, fields read) Drag Coefficient 62 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_sd ! structure of input fields (file informations, fields read) Stokes Drift 63 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_wn ! structure of input fields (file informations, fields read) wave number for Qiao 64 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 65 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tauw ! structure of input fields (file informations, fields read) ocean stress components from wave model 66 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_phioc ! structure of input fields (file informations, fields read) wave to ocean energy 67 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: cdn_wave !: 68 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: hsw, wmp, wnum !: 69 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: wfreq !: 70 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: rn_crban !: Craig and Banner constant for surface breaking waves mixing 71 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauoc_wave !: 72 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauw_x !: 73 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tauw_y !: 74 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: tsd2d !: 75 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: div_sd !: barotropic stokes drift divergence 76 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) :: ut0sd, vt0sd !: surface Stokes drift velocities at t-point 77 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:) :: usd , vsd , wsd !: Stokes drift velocities at u-, v- & w-points, resp. 78 79 # include "vectopt_loop_substitute.h90" 39 80 !!---------------------------------------------------------------------- 40 81 !! NEMO/OPA 4.0 , NEMO Consortium (2011) … … 44 85 CONTAINS 45 86 87 SUBROUTINE sbc_stokes( ) 88 !!--------------------------------------------------------------------- 89 !! *** ROUTINE sbc_stokes *** 90 !! 91 !! ** Purpose : compute the 3d Stokes Drift according to Breivik et al., 92 !! 2014 (DOI: 10.1175/JPO-D-14-0020.1) 93 !! 94 !! ** Method : - Calculate Stokes transport speed 95 !! - Calculate horizontal divergence 96 !! - Integrate the horizontal divergenze from the bottom 97 !! ** action 98 !!--------------------------------------------------------------------- 99 INTEGER :: jj, ji, jk ! dummy loop argument 100 INTEGER :: ik ! local integer 101 REAL(wp) :: ztransp, zfac, ztemp 102 REAL(wp) :: zdep_u, zdep_v, zkh_u, zkh_v, zda_u, zda_v 103 REAL(wp), DIMENSION(:,:) , POINTER :: zk_t, zk_u, zk_v, zu0_sd, zv0_sd ! 2D workspace 104 REAL(wp), DIMENSION(:,:,:), POINTER :: ze3divh ! 3D workspace 105 106 !!--------------------------------------------------------------------- 107 ! 108 109 CALL wrk_alloc( jpi,jpj,jpk, ze3divh ) 110 CALL wrk_alloc( jpi,jpj, zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 111 ! 112 ! select parameterization for the calculation of vertical Stokes drift 113 ! exp. wave number at t-point 114 IF( nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips ) THEN ! (Eq. (19) in Breivick et al. (2014) ) 115 zfac = 2.0_wp * rpi / 16.0_wp 116 DO jj = 1, jpj 117 DO ji = 1, jpi 118 ! Stokes drift velocity estimated from Hs and Tmean 119 ztransp = zfac * hsw(ji,jj)*hsw(ji,jj) / MAX( wmp(ji,jj), 0.0000001_wp ) 120 ! Stokes surface speed 121 tsd2d(ji,jj) = SQRT( ut0sd(ji,jj)*ut0sd(ji,jj) + vt0sd(ji,jj)*vt0sd(ji,jj)) 122 ! Wavenumber scale 123 zk_t(ji,jj) = ABS( tsd2d(ji,jj) ) / MAX( ABS( 5.97_wp*ztransp ), 0.0000001_wp ) 124 END DO 125 END DO 126 DO jj = 1, jpjm1 ! exp. wave number & Stokes drift velocity at u- & v-points 127 DO ji = 1, jpim1 128 zk_u(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji+1,jj) ) 129 zk_v(ji,jj) = 0.5_wp * ( zk_t(ji,jj) + zk_t(ji,jj+1) ) 130 ! 131 zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 132 zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 133 END DO 134 END DO 135 ELSE IF( nn_sdrift==jp_peakph ) THEN ! peak wave number calculated from the peak frequency received by the wave model 136 DO jj = 1, jpjm1 137 DO ji = 1, jpim1 138 zk_u(ji,jj) = 0.5_wp * ( wfreq(ji,jj)*wfreq(ji,jj) + wfreq(ji+1,jj)*wfreq(ji+1,jj) ) / grav 139 zk_v(ji,jj) = 0.5_wp * ( wfreq(ji,jj)*wfreq(ji,jj) + wfreq(ji,jj+1)*wfreq(ji,jj+1) ) / grav 140 ! 141 zu0_sd(ji,jj) = 0.5_wp * ( ut0sd(ji,jj) + ut0sd(ji+1,jj) ) 142 zv0_sd(ji,jj) = 0.5_wp * ( vt0sd(ji,jj) + vt0sd(ji,jj+1) ) 143 END DO 144 END DO 145 ENDIF 146 ! 147 ! !== horizontal Stokes Drift 3D velocity ==! 148 IF( nn_sdrift==jp_breivik ) THEN 149 DO jk = 1, jpkm1 150 DO jj = 2, jpjm1 151 DO ji = 2, jpim1 152 zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) 153 zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) 154 ! 155 zkh_u = zk_u(ji,jj) * zdep_u ! k * depth 156 zkh_v = zk_v(ji,jj) * zdep_v 157 ! ! Depth attenuation 158 zda_u = EXP( -2.0_wp*zkh_u ) / ( 1.0_wp + 8.0_wp*zkh_u ) 159 zda_v = EXP( -2.0_wp*zkh_v ) / ( 1.0_wp + 8.0_wp*zkh_v ) 160 ! 161 usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 162 vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 163 END DO 164 END DO 165 END DO 166 ELSE IF( nn_sdrift==jp_phillips .OR. nn_sdrift==jp_peakph ) THEN 167 DO jk = 1, jpkm1 168 DO jj = 2, jpjm1 169 DO ji = 2, jpim1 170 zdep_u = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji+1,jj,jk) ) 171 zdep_v = 0.5_wp * ( gdept_n(ji,jj,jk) + gdept_n(ji,jj+1,jk) ) 172 ! 173 zkh_u = zk_u(ji,jj) * zdep_u ! k * depth 174 zkh_v = zk_v(ji,jj) * zdep_v 175 ! ! Depth attenuation 176 zda_u = EXP( -2.0_wp*zkh_u ) - SQRT(2.0_wp*rpi*zkh_u) * ERFC(SQRT(2.0_wp*zkh_u)) 177 zda_v = EXP( -2.0_wp*zkh_v ) - SQRT(2.0_wp*rpi*zkh_v) * ERFC(SQRT(2.0_wp*zkh_v)) 178 ! 179 usd(ji,jj,jk) = zda_u * zu0_sd(ji,jj) * umask(ji,jj,jk) 180 vsd(ji,jj,jk) = zda_v * zv0_sd(ji,jj) * vmask(ji,jj,jk) 181 END DO 182 END DO 183 END DO 184 ENDIF 185 186 CALL lbc_lnk( usd(:,:,:), 'U', vsd(:,:,:), 'V', -1. ) 187 ! 188 ! !== vertical Stokes Drift 3D velocity ==! 189 ! 190 DO jk = 1, jpkm1 ! Horizontal e3*divergence 191 DO jj = 2, jpj 192 DO ji = fs_2, jpi 193 ze3divh(ji,jj,jk) = ( e2u(ji ,jj) * e3u_n(ji ,jj,jk) * usd(ji, jj,jk) & 194 & - e2u(ji-1,jj) * e3u_n(ji-1,jj,jk) * usd(ji-1,jj,jk) & 195 & + e1v(ji,jj ) * e3v_n(ji,jj ,jk) * vsd(ji,jj ,jk) & 196 & - e1v(ji,jj-1) * e3v_n(ji,jj-1,jk) * vsd(ji,jj-1,jk) ) * r1_e12t(ji,jj) 197 END DO 198 END DO 199 END DO 200 ! 201 IF( .NOT. AGRIF_Root() ) THEN 202 IF( nbondi == 1 .OR. nbondi == 2 ) ze3divh(nlci-1, : ,:) = 0._wp ! east 203 IF( nbondi == -1 .OR. nbondi == 2 ) ze3divh( 2 , : ,:) = 0._wp ! west 204 IF( nbondj == 1 .OR. nbondj == 2 ) ze3divh( : ,nlcj-1,:) = 0._wp ! north 205 IF( nbondj == -1 .OR. nbondj == 2 ) ze3divh( : , 2 ,:) = 0._wp ! south 206 ENDIF 207 ! 208 CALL lbc_lnk( ze3divh, 'T', 1. ) 209 ! 210 IF( .NOT. lk_vvl ) THEN ; ik = 1 ! none zero velocity through the sea surface 211 ELSE ; ik = 2 ! w=0 at the surface (set one for all in sbc_wave_init) 212 ENDIF 213 DO jk = jpkm1, ik, -1 ! integrate from the bottom the hor. divergence (NB: at k=jpk w is always zero) 214 wsd(:,:,jk) = wsd(:,:,jk+1) - ze3divh(:,:,jk) 215 END DO 216 #if defined key_bdy 217 IF( lk_bdy ) THEN 218 DO jk = 1, jpkm1 219 wsd(:,:,jk) = wsd(:,:,jk) * bdytmask(:,:) 220 END DO 221 ENDIF 222 #endif 223 ! !== Horizontal divergence of barotropic Stokes transport ==! 224 div_sd(:,:) = 0._wp 225 DO jk = 1, jpkm1 ! 226 div_sd(:,:) = div_sd(:,:) + ze3divh(:,:,jk) 227 END DO 228 ! 229 CALL iom_put( "ustokes", usd ) 230 CALL iom_put( "vstokes", vsd ) 231 CALL iom_put( "wstokes", wsd ) 232 ! 233 CALL wrk_dealloc( jpi,jpj,jpk, ze3divh ) 234 CALL wrk_dealloc( jpi,jpj, zk_t, zk_u, zk_v, zu0_sd, zv0_sd ) 235 ! 236 END SUBROUTINE sbc_stokes 237 238 239 SUBROUTINE sbc_stress( ) 240 !!--------------------------------------------------------------------- 241 !! *** ROUTINE sbc_stress *** 242 !! 243 !! ** Purpose : Updates the ocean momentum modified by waves 244 !! 245 !! ** Method : - Calculate u,v components of stress depending on stress 246 !! model 247 !! - Calculate the stress module 248 !! - The wind module is not modified by waves 249 !! ** action 250 !!--------------------------------------------------------------------- 251 INTEGER :: jj, ji ! dummy loop argument 252 ! 253 IF( ln_tauoc ) THEN 254 utau(:,:) = utau(:,:)*tauoc_wave(:,:) 255 vtau(:,:) = vtau(:,:)*tauoc_wave(:,:) 256 taum(:,:) = taum(:,:)*tauoc_wave(:,:) 257 ENDIF 258 ! 259 IF( ln_tauw ) THEN 260 DO jj = 1, jpjm1 261 DO ji = 1, jpim1 262 ! Stress components at u- & v-points 263 utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) ) 264 vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) ) 265 ! 266 ! Stress module at t points 267 taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) ) 268 END DO 269 END DO 270 CALL lbc_lnk_multi( utau(:,:), 'U', -1. , vtau(:,:), 'V', -1. , taum(:,: ), 'T', -1. ) 271 ENDIF 272 ! 273 END SUBROUTINE sbc_stress 274 275 46 276 SUBROUTINE sbc_wave( kt ) 47 277 !!--------------------------------------------------------------------- 48 !! *** ROUTINE sbc_ apr***49 !! 50 !! ** Purpose : read drag coefficientfrom wave model in netcdf files.278 !! *** ROUTINE sbc_wave *** 279 !! 280 !! ** Purpose : read wave parameters from wave model in netcdf files. 51 281 !! 52 282 !! ** Method : - Read namelist namsbc_wave 53 283 !! - Read Cd_n10 fields in netcdf files 54 284 !! - Read stokes drift 2d in netcdf files 55 !! - Read wave number in netcdf files 56 !! - Compute 3d stokes drift using monochromatic 57 !! ** action : 58 !! 59 !!--------------------------------------------------------------------- 60 USE oce, ONLY : un,vn,hdivn,rotn 61 USE divcur 62 USE wrk_nemo 63 #if defined key_bdy 64 USE bdy_oce, ONLY : bdytmask 65 #endif 66 INTEGER, INTENT( in ) :: kt ! ocean time step 67 INTEGER :: ierror ! return error code 68 INTEGER :: ifpr, jj,ji,jk 69 INTEGER :: ios ! Local integer output status for namelist read 70 REAL(wp),DIMENSION(:,:,:),POINTER :: udummy,vdummy,hdivdummy,rotdummy 71 REAL :: z2dt,z1_2dt 72 TYPE(FLD_N), DIMENSION(jpfld) :: slf_i ! array of namelist informations on the fields to read 285 !! - Read wave number in netcdf files 286 !! - Compute 3d stokes drift using Breivik et al.,2014 287 !! formulation 288 !! ** action 289 !!--------------------------------------------------------------------- 290 INTEGER, INTENT(in ) :: kt ! ocean time step 291 !!--------------------------------------------------------------------- 292 ! 293 IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN !== Neutral drag coefficient ==! 294 CALL fld_read( kt, nn_fsbc, sf_cd ) ! read from external forcing 295 cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 296 ! check that the drag coefficient contains proper information even if 297 ! the masks do not match - the momentum stress is not masked! 298 WHERE( cdn_wave < 0.0 ) cdn_wave = 1.5e-3 299 WHERE( cdn_wave > 1.0 ) cdn_wave = 1.5e-3 300 ENDIF 301 302 IF( ln_tauoc .AND. .NOT. cpl_tauoc ) THEN !== Wave induced stress ==! 303 CALL fld_read( kt, nn_fsbc, sf_tauoc ) ! read wave norm stress from external forcing 304 tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 305 WHERE( tauoc_wave < 0.0 ) tauoc_wave = 1.0 306 WHERE( tauoc_wave > 100.0 ) tauoc_wave = 1.0 307 ENDIF 308 309 IF( ln_tauw .AND. .NOT. cpl_tauw ) THEN !== Wave induced stress ==! 310 CALL fld_read( kt, nn_fsbc, sf_tauw ) ! read ocean stress components from external forcing (T grid) 311 tauw_x(:,:) = sf_tauw(1)%fnow(:,:,1) 312 WHERE( tauw_x < -100.0 ) tauw_x = 0.0 313 WHERE( tauw_x > 100.0 ) tauw_x = 0.0 314 315 tauw_y(:,:) = sf_tauw(2)%fnow(:,:,1) 316 WHERE( tauw_y < -100.0 ) tauw_y = 0.0 317 WHERE( tauw_y > 100.0 ) tauw_y = 0.0 318 ENDIF 319 320 IF( ln_phioc .AND. .NOT. cpl_phioc ) THEN !== Wave to ocean energy ==! 321 CALL fld_read( kt, nn_fsbc, sf_phioc ) ! read wave to ocean energy from external forcing 322 rn_crban(:,:) = 29.0 * sf_phioc(1)%fnow(:,:,1) ! ! Alfa is phioc*sqrt(rau0/zrhoa) : rau0=water density, zhroa= air density 323 WHERE( rn_crban > 1.e8 ) rn_crban = 0.0 !remove first mask mistmatch points, then cap values in case of low friction velocity 324 WHERE( rn_crban < 0.0 ) rn_crban = 0.0 325 WHERE( rn_crban > 1000.0 ) rn_crban = 1000.0 326 ENDIF 327 328 IF( ln_sdw .OR. ln_rough ) THEN !== Computation of the 3d Stokes Drift ==! 329 ! 330 IF( jpfld > 0 ) THEN ! Read from file only if the field is not coupled 331 CALL fld_read( kt, nn_fsbc, sf_sd ) ! read wave parameters from external forcing 332 IF( jp_hsw > 0 ) THEN 333 hsw (:,:) = sf_sd(jp_hsw)%fnow(:,:,1) ! significant wave height 334 WHERE( hsw > 100.0 ) hsw = 0.0 335 WHERE( hsw < 0.0 ) hsw = 0.0 336 ENDIF 337 IF( jp_wmp > 0 ) THEN 338 wmp (:,:) = sf_sd(jp_wmp)%fnow(:,:,1) ! wave mean period 339 WHERE( wmp > 100.0 ) wmp = 0.0 340 WHERE( wmp < 0.0 ) wmp = 0.0 341 ENDIF 342 IF( jp_wfr > 0 ) THEN 343 wfreq(:,:) = sf_sd(jp_wfr)%fnow(:,:,1) ! Peak wave frequency 344 WHERE( wfreq < 0.0 ) wfreq = 0.0 345 WHERE( wfreq > 100.0 ) wfreq = 0.0 346 ENDIF 347 IF( jp_usd > 0 ) THEN 348 ut0sd(:,:) = sf_sd(jp_usd)%fnow(:,:,1) ! 2D zonal Stokes Drift at T point 349 WHERE( ut0sd < -100.0 ) ut0sd = 1.0 350 WHERE( ut0sd > 100.0 ) ut0sd = 1.0 351 ENDIF 352 IF( jp_vsd > 0 ) THEN 353 vt0sd(:,:) = sf_sd(jp_vsd)%fnow(:,:,1) ! 2D meridional Stokes Drift at T point 354 WHERE( vt0sd < -100.0 ) vt0sd = 1.0 355 WHERE( vt0sd > 100.0 ) vt0sd = 1.0 356 ENDIF 357 ENDIF 358 ENDIF 359 ! 360 IF( ln_sdw ) THEN 361 ! Read also wave number if needed, so that it is available in coupling routines 362 IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 363 CALL fld_read( kt, nn_fsbc, sf_wn ) ! read wave parameters from external forcing 364 wnum(:,:) = sf_wn(1)%fnow(:,:,1) 365 ENDIF 366 367 ! !== Computation of the 3d Stokes Drift ==! 368 ! 369 IF( ((nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) .AND. & 370 jp_hsw>0 .AND. jp_wmp>0 .AND. jp_usd>0 .AND. jp_vsd>0) .OR. & 371 (nn_sdrift==jp_peakph .AND. jp_wfr>0 .AND. jp_usd>0 .AND. jp_vsd>0) ) & 372 CALL sbc_stokes() ! Calculate only if required fields are read 373 ! ! In coupled wave model-NEMO case the call is done after coupling 374 ! 375 ENDIF 376 ! 377 END SUBROUTINE sbc_wave 378 379 380 SUBROUTINE sbc_wave_init 381 !!--------------------------------------------------------------------- 382 !! *** ROUTINE sbc_wave_init *** 383 !! 384 !! ** Purpose : read wave parameters from wave model in netcdf files. 385 !! 386 !! ** Method : - Read namelist namsbc_wave 387 !! - Read Cd_n10 fields in netcdf files 388 !! - Read stokes drift 2d in netcdf files 389 !! - Read wave number in netcdf files 390 !! - Compute 3d stokes drift using Breivik et al.,2014 391 !! formulation 392 !! ** action 393 !!--------------------------------------------------------------------- 394 INTEGER :: ierror, ios ! local integer 395 INTEGER :: ifpr 396 !! 73 397 CHARACTER(len=100) :: cn_dir ! Root directory for location of drag coefficient files 74 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, sn_wn ! informations about the fields to be read 75 !!--------------------------------------------------------------------- 76 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_wn 77 !!--------------------------------------------------------------------- 78 79 !!---------------------------------------------------------------------- 80 ! 81 ! 82 ! ! -------------------- ! 83 IF( kt == nit000 ) THEN ! First call kt=nit000 ! 84 ! ! -------------------- ! 85 REWIND( numnam_ref ) ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 86 READ ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 87 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 88 89 REWIND( numnam_cfg ) ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 90 READ ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 91 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 92 IF(lwm) WRITE ( numond, namsbc_wave ) 93 ! 94 95 IF ( ln_cdgw ) THEN 398 TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) :: slf_i, slf_j ! array of namelist informations on the fields to read 399 TYPE(FLD_N) :: sn_cdg, sn_usd, sn_vsd, sn_phioc, & 400 & sn_hsw, sn_wmp, sn_wfr, sn_wnum , & 401 & sn_tauoc, sn_tauwx, sn_tauwy ! information about the fields to be read 402 ! 403 NAMELIST/namsbc_wave/ sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, sn_wnum, sn_phioc, & 404 sn_tauoc, sn_tauwx, sn_tauwy, & 405 ln_cdgw, ln_sdw, ln_stcor, ln_phioc, ln_tauoc, ln_tauw, ln_zdfqiao, ln_rough, & 406 nn_sdrift, nn_wmix 407 !!--------------------------------------------------------------------- 408 ! 409 REWIND( numnam_ref ) ! Namelist namsbc_wave in reference namelist : File for drag coeff. from wave model 410 READ ( numnam_ref, namsbc_wave, IOSTAT = ios, ERR = 901) 411 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in reference namelist', lwp ) 412 413 REWIND( numnam_cfg ) ! Namelist namsbc_wave in configuration namelist : File for drag coeff. from wave model 414 READ ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 ) 415 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp ) 416 IF(lwm) WRITE ( numond, namsbc_wave ) 417 ! 418 IF(lwp) THEN !* Parameter print 419 WRITE(numout,*) 420 WRITE(numout,*) 'sbc_wave_init: wave physics' 421 WRITE(numout,*) '~~~~~~~~' 422 WRITE(numout,*) ' Namelist namsbc_wave : set wave physics parameters' 423 WRITE(numout,*) ' Stokes drift corr. to vert. velocity ln_sdw = ', ln_sdw 424 WRITE(numout,*) ' vertical parametrization nn_sdrift = ', nn_sdrift 425 WRITE(numout,*) ' Stokes coriolis term ln_stcor = ', ln_stcor 426 WRITE(numout,*) ' wave modified ocean stress ln_tauoc = ', ln_tauoc 427 WRITE(numout,*) ' wave modified ocean stress components ln_tauw = ', ln_tauw 428 WRITE(numout,*) ' wave to ocean energy ln_phioc = ', ln_phioc 429 WRITE(numout,*) ' vertical mixing parametrization nn_wmix = ', nn_wmix 430 WRITE(numout,*) ' neutral drag coefficient ln_cdgw = ', ln_cdgw 431 WRITE(numout,*) ' wave roughness length modification ln_rough = ', ln_rough 432 WRITE(numout,*) ' Qiao vertical mixing formulation ln_zdfqiao = ', ln_zdfqiao 433 ENDIF 434 435 IF ( ln_wave ) THEN 436 ! Activated wave physics but no wave physics components activated 437 IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_tauw .OR. ln_stcor .OR. ln_phioc & 438 .OR. ln_rough .OR. ln_zdfqiao) ) THEN 439 CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_tauw=F, ln_stcor=F ', & 440 'ln_phioc=F, ln_rough=F, ln_zdfqiao=F' ) 441 ELSE 442 IF (ln_stcor .AND. .NOT. ln_sdw) & 443 CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T') 444 IF ( ln_cdgw .AND. .NOT.(nn_drag==jp_ukmo .OR. nn_drag==jp_std .OR. nn_drag==jp_const .OR. nn_drag==jp_mcore) ) & 445 CALL ctl_stop( 'The chosen nn_drag for momentum calculation must be 0, 1, 2, or 3') 446 IF ( ln_cdgw .AND. ln_blk_core .AND. nn_drag==0 ) & 447 CALL ctl_stop( 'The chosen nn_drag for momentum calculation in core forcing must be 1, 2, or 3') 448 IF ( ln_cdgw .AND. ln_flx .AND. nn_drag==3 ) & 449 CALL ctl_stop( 'The chosen nn_drag for momentum calculation in direct forcing must be 0, 1, or 2') 450 IF( ln_phioc .AND. .NOT.(nn_wmix==jp_craigbanner .OR. nn_wmix==jp_janssen) ) & 451 CALL ctl_stop( 'The chosen nn_wmix for wave vertical mixing must be 0, or 1' ) 452 IF( ln_sdw .AND. .NOT.(nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips .OR. nn_sdrift==jp_peakph) ) & 453 CALL ctl_stop( 'The chosen nn_sdrift for Stokes drift vertical velocity must be 0, 1, or 2' ) 454 IF( ln_zdfqiao .AND. .NOT.ln_sdw ) & 455 CALL ctl_stop( 'Qiao vertical mixing can not be used without Stokes drift (ln_sdw)' ) 456 IF( ln_tauoc .AND. ln_tauw ) & 457 CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 458 '(ln_tauoc=.true. and ln_tauw=.true.)' ) 459 IF( ln_tauoc ) & 460 CALL ctl_warn( 'You are subtracting the wave stress to the ocean (ln_tauoc=.true.)' ) 461 IF( ln_tauw ) & 462 CALL ctl_warn( 'The wave modified ocean stress components are used (ln_tauw=.true.) ', & 463 'This will override any other specification of the ocean stress' ) 464 ENDIF 465 ELSE 466 IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_tauw .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao ) & 467 & CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ', & 468 & 'with drag coefficient (ln_cdgw =T) ' , & 469 & 'or Stokes Drift (ln_sdw=T) ' , & 470 & 'or Stokes-Coriolis term (ln_stcor=T)', & 471 & 'or ocean stress modification due to waves (ln_tauoc=T) ', & 472 & 'or ocean stress components from waves (ln_tauw=T) ', & 473 & 'or wave to ocean energy modification (ln_phioc=T) ', & 474 & 'or wave surface roughness (ln_rough=T) ', & 475 & 'or Qiao vertical mixing formulation (ln_zdfqiao=T) ' ) 476 ENDIF 477 ! 478 IF( ln_cdgw ) THEN 479 IF( .NOT. cpl_wdrag ) THEN 96 480 ALLOCATE( sf_cd(1), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 97 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave : unable to allocate sf_wavestructure' )481 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_cd structure' ) 98 482 ! 99 483 ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1) ) 100 484 IF( sn_cdg%ln_tint ) ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 101 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 102 ALLOCATE( cdn_wave(jpi,jpj) ) 103 cdn_wave(:,:) = 0.0 104 ENDIF 105 IF ( ln_sdw ) THEN 106 slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; slf_i(jp_wn) = sn_wn 107 ALLOCATE( sf_sd(3), STAT=ierror ) !* allocate and fill sf_wave with sn_cdg 108 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 485 CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 486 ENDIF 487 ALLOCATE( cdn_wave(jpi,jpj) ) 488 ENDIF 489 490 IF( ln_tauoc ) THEN 491 IF( .NOT. cpl_tauoc ) THEN 492 ALLOCATE( sf_tauoc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_tauoc 493 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauoc structure' ) 109 494 ! 110 DO ifpr= 1, jpfld 111 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 112 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 113 END DO 114 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 115 ALLOCATE( usd2d(jpi,jpj),vsd2d(jpi,jpj),uwavenum(jpi,jpj),vwavenum(jpi,jpj) ) 116 ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 117 usd2d(:,:) = 0.0 ; vsd2d(:,:) = 0.0 ; uwavenum(:,:) = 0.0 ; vwavenum(:,:) = 0.0 118 usd3d(:,:,:) = 0.0 ;vsd3d(:,:,:) = 0.0 ; wsd3d(:,:,:) = 0.0 119 ENDIF 120 ENDIF 495 ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1) ) 496 IF( sn_tauoc%ln_tint ) ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 497 CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 498 ENDIF 499 ALLOCATE( tauoc_wave(jpi,jpj) ) 500 ENDIF 501 502 IF( ln_tauw ) THEN 503 IF( .NOT. cpl_tauw ) THEN 504 ALLOCATE( sf_tauw(2), STAT=ierror ) !* allocate and fill sf_wave with sn_tauwx/y 505 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauw structure' ) 506 ! 507 ALLOCATE( slf_j(2) ) 508 slf_j(1) = sn_tauwx 509 slf_j(2) = sn_tauwy 510 ALLOCATE( sf_tauw(1)%fnow(jpi,jpj,1) ) 511 ALLOCATE( sf_tauw(2)%fnow(jpi,jpj,1) ) 512 IF( slf_j(1)%ln_tint ) ALLOCATE( sf_tauw(1)%fdta(jpi,jpj,1,2) ) 513 IF( slf_j(2)%ln_tint ) ALLOCATE( sf_tauw(2)%fdta(jpi,jpj,1,2) ) 514 CALL fld_fill( sf_tauw, (/ slf_j /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 515 ENDIF 516 ALLOCATE( tauw_x(jpi,jpj) ) 517 ALLOCATE( tauw_y(jpi,jpj) ) 518 ENDIF 519 520 IF( ln_phioc ) THEN 521 IF( .NOT. cpl_phioc ) THEN 522 ALLOCATE( sf_phioc(1), STAT=ierror ) !* allocate and fill sf_wave with sn_phioc 523 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_phioc structure' ) 524 ! 525 ALLOCATE( sf_phioc(1)%fnow(jpi,jpj,1) ) 526 IF( sn_phioc%ln_tint ) ALLOCATE( sf_phioc(1)%fdta(jpi,jpj,1,2) ) 527 CALL fld_fill( sf_phioc, (/ sn_phioc /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 528 ENDIF 529 ALLOCATE( rn_crban(jpi,jpj) ) 530 ENDIF 531 532 ! Find out how many fields have to be read from file if not coupled 533 jpfld=0 534 jp_usd=0 ; jp_vsd=0 ; jp_hsw=0 ; jp_wmp=0 ; jp_wfr=0 535 IF( ln_sdw ) THEN 536 IF( .NOT. cpl_sdrft ) THEN 537 jpfld = jpfld + 1 538 jp_usd = jpfld 539 jpfld = jpfld + 1 540 jp_vsd = jpfld 541 ENDIF 542 IF( .NOT. cpl_hsig .AND. (nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) ) THEN 543 jpfld = jpfld + 1 544 jp_hsw = jpfld 545 ENDIF 546 IF( .NOT. cpl_wper .AND. (nn_sdrift==jp_breivik .OR. nn_sdrift==jp_phillips) ) THEN 547 jpfld = jpfld + 1 548 jp_wmp = jpfld 549 ENDIF 550 IF( .NOT. cpl_wfreq .AND. nn_sdrift==jp_peakph ) THEN 551 jpfld = jpfld + 1 552 jp_wfr = jpfld 553 ENDIF 554 ENDIF 555 556 IF( ln_rough .AND. .NOT. cpl_hsig .AND. jp_hsw==0 ) THEN 557 jpfld = jpfld + 1 558 jp_hsw = jpfld 559 ENDIF 560 561 ! Read from file only the non-coupled fields 562 IF( jpfld > 0 ) THEN 563 ALLOCATE( slf_i(jpfld) ) 564 IF( jp_usd > 0 ) slf_i(jp_usd) = sn_usd 565 IF( jp_vsd > 0 ) slf_i(jp_vsd) = sn_vsd 566 IF( jp_hsw > 0 ) slf_i(jp_hsw) = sn_hsw 567 IF( jp_wmp > 0 ) slf_i(jp_wmp) = sn_wmp 568 IF( jp_wfr > 0 ) slf_i(jp_wfr) = sn_wfr 569 ALLOCATE( sf_sd(jpfld), STAT=ierror ) !* allocate and fill sf_sd with stokes drift 570 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_sd structure' ) 121 571 ! 572 DO ifpr= 1, jpfld 573 ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 574 IF( slf_i(ifpr)%ln_tint ) ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 575 END DO 122 576 ! 123 IF ( ln_cdgw ) THEN 124 CALL fld_read( kt, nn_fsbc, sf_cd ) !* read drag coefficient from external forcing 125 cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 126 ENDIF 127 IF ( ln_sdw ) THEN 128 CALL fld_read( kt, nn_fsbc, sf_sd ) !* read drag coefficient from external forcing 129 130 ! Interpolate wavenumber, stokes drift into the grid_V and grid_V 131 !------------------------------------------------- 132 133 DO jj = 1, jpjm1 134 DO ji = 1, jpim1 135 uwavenum(ji,jj)=0.5 * ( 2. - umask(ji,jj,1) ) * ( sf_sd(3)%fnow(ji,jj,1) * tmask(ji,jj,1) & 136 & + sf_sd(3)%fnow(ji+1,jj,1) * tmask(ji+1,jj,1) ) 137 138 vwavenum(ji,jj)=0.5 * ( 2. - vmask(ji,jj,1) ) * ( sf_sd(3)%fnow(ji,jj,1) * tmask(ji,jj,1) & 139 & + sf_sd(3)%fnow(ji,jj+1,1) * tmask(ji,jj+1,1) ) 140 141 usd2d(ji,jj) = 0.5 * ( 2. - umask(ji,jj,1) ) * ( sf_sd(1)%fnow(ji,jj,1) * tmask(ji,jj,1) & 142 & + sf_sd(1)%fnow(ji+1,jj,1) * tmask(ji+1,jj,1) ) 143 144 vsd2d(ji,jj) = 0.5 * ( 2. - vmask(ji,jj,1) ) * ( sf_sd(2)%fnow(ji,jj,1) * tmask(ji,jj,1) & 145 & + sf_sd(2)%fnow(ji,jj+1,1) * tmask(ji,jj+1,1) ) 146 END DO 147 END DO 148 149 !Computation of the 3d Stokes Drift 150 DO jk = 1, jpk 151 DO jj = 1, jpj-1 152 DO ji = 1, jpi-1 153 usd3d(ji,jj,jk) = usd2d(ji,jj)*exp(2.0*uwavenum(ji,jj)*(-MIN( gdept_0(ji,jj,jk) , gdept_0(ji+1,jj ,jk)))) 154 vsd3d(ji,jj,jk) = vsd2d(ji,jj)*exp(2.0*vwavenum(ji,jj)*(-MIN( gdept_0(ji,jj,jk) , gdept_0(ji ,jj+1,jk)))) 155 END DO 156 END DO 157 usd3d(jpi,:,jk) = usd2d(jpi,:)*exp( 2.0*uwavenum(jpi,:)*(-gdept_0(jpi,:,jk)) ) 158 vsd3d(:,jpj,jk) = vsd2d(:,jpj)*exp( 2.0*vwavenum(:,jpj)*(-gdept_0(:,jpj,jk)) ) 159 END DO 160 161 CALL wrk_alloc( jpi,jpj,jpk,udummy,vdummy,hdivdummy,rotdummy) 162 163 udummy(:,:,:)=un(:,:,:) 164 vdummy(:,:,:)=vn(:,:,:) 165 hdivdummy(:,:,:)=hdivn(:,:,:) 166 rotdummy(:,:,:)=rotn(:,:,:) 167 un(:,:,:)=usd3d(:,:,:) 168 vn(:,:,:)=vsd3d(:,:,:) 169 CALL div_cur(kt) 170 ! !------------------------------! 171 ! ! Now Vertical Velocity ! 172 ! !------------------------------! 173 z2dt = 2._wp * rdt ! set time step size (Euler/Leapfrog) 174 175 z1_2dt = 1.e0 / z2dt 176 DO jk = jpkm1, 1, -1 ! integrate from the bottom the hor. divergence 177 ! - ML - need 3 lines here because replacement of fse3t by its expression yields too long lines otherwise 178 wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - fse3t_n(:,:,jk) * hdivn(:,:,jk) & 179 & - ( fse3t_a(:,:,jk) - fse3t_b(:,:,jk) ) & 180 & * tmask(:,:,jk) * z1_2dt 181 #if defined key_bdy 182 wsd3d(:,:,jk) = wsd3d(:,:,jk) * bdytmask(:,:) 183 #endif 184 END DO 185 hdivn(:,:,:)=hdivdummy(:,:,:) 186 rotn(:,:,:)=rotdummy(:,:,:) 187 vn(:,:,:)=vdummy(:,:,:) 188 un(:,:,:)=udummy(:,:,:) 189 CALL wrk_dealloc( jpi,jpj,jpk,udummy,vdummy,hdivdummy,rotdummy) 190 ENDIF 191 END SUBROUTINE sbc_wave 192 577 CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' ) 578 ENDIF 579 580 IF( ln_sdw ) THEN 581 ALLOCATE( usd (jpi,jpj,jpk), vsd (jpi,jpj,jpk), wsd(jpi,jpj,jpk) ) 582 ALLOCATE( wmp (jpi,jpj) ) 583 ALLOCATE( wfreq (jpi,jpj) ) 584 ALLOCATE( ut0sd(jpi,jpj) , vt0sd(jpi,jpj) ) 585 ALLOCATE( div_sd(jpi,jpj) ) 586 ALLOCATE( tsd2d (jpi,jpj) ) 587 usd(:,:,:) = 0._wp 588 vsd(:,:,:) = 0._wp 589 wsd(:,:,:) = 0._wp 590 ! Wave number needed only if ln_zdfqiao=T 591 IF( ln_zdfqiao .AND. .NOT.cpl_wnum ) THEN 592 ALLOCATE( sf_wn(1), STAT=ierror ) !* allocate and fill sf_wave with sn_wnum 593 IF( ierror > 0 ) CALL ctl_stop( 'STOP', 'sbc_wave_init: unable toallocate sf_wn structure' ) 594 ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1) ) 595 IF( sn_wnum%ln_tint ) ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 596 CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'read wave input', 'namsbc_wave' ) 597 ENDIF 598 ALLOCATE( wnum(jpi,jpj) ) 599 ENDIF 600 601 IF( ln_sdw .OR. ln_rough ) THEN 602 ALLOCATE( hsw (jpi,jpj) ) 603 ENDIF 604 ! 605 END SUBROUTINE sbc_wave_init 606 193 607 !!====================================================================== 194 608 END MODULE sbcwave
Note: See TracChangeset
for help on using the changeset viewer.
