Changeset 11277 for branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
- Timestamp:
- 2019-07-17T15:29:15+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
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
Note: See TracChangeset
for help on using the changeset viewer.