Changeset 1218 for trunk/NEMO/OPA_SRC/SBC/sbccpl.F90
- Timestamp:
- 2008-10-28T10:12:16+01:00 (16 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/SBC/sbccpl.F90
r1156 r1218 2 2 !!====================================================================== 3 3 !! *** MODULE sbccpl *** 4 !! Ocean forcing: momentum, heat and freshwater coupled formulation5 !!===================================================================== 6 !! History : 9.0 ! 06-07 (R. Redler, N. Keenlyside, W. Park)7 !! Original code split into flxmod & taumod8 !! 9.0 ! 06-07 (G. Madec) surface module4 !! Surface Boundary Condition : momentum, heat and freshwater fluxes in coupled mode 5 !!====================================================================== 6 !! History : 2.0 ! 06-2007 (R. Redler, N. Keenlyside, W. Park) Original code split into flxmod & taumod 7 !! 3.0 ! 02-2008 (G. Madec, C Talandier) surface module 8 !! - ! 08-2008 (S. Masson, E. .... ) generic coupled interface 9 9 !!---------------------------------------------------------------------- 10 #if defined key_ sbc_cpl10 #if defined key_oasis3 || defined key_oasis4 11 11 !!---------------------------------------------------------------------- 12 !! 'key_ sbc_cpl'Coupled Ocean/Atmosphere formulation12 !! 'key_oasis3' or 'key_oasis4' Coupled Ocean/Atmosphere formulation 13 13 !!---------------------------------------------------------------------- 14 !! namsbc_cpl : coupled formulation namlist 15 !! sbc_cpl_init : initialisation of the coupled exchanges 16 !! sbc_cpl_rcv : receive fields from the atmosphere over the ocean (ocean only) 17 !! receive stress from the atmosphere over the ocean (ocean-ice case) 18 !! sbc_cpl_ice_tau : receive stress from the atmosphere over ice 19 !! sbc_cpl_ice_flx : receive fluxes from the atmosphere over ice 20 !! sbc_cpl_snd : send fields to the atmosphere 14 21 !!---------------------------------------------------------------------- 15 !! namsbc_cpl : coupled formulation namlist16 !! sbc_cpl : coupled formulation for the ocean surface boundary condition17 !!----------------------------------------------------------------------18 USE oce ! ocean dynamics and tracers19 22 USE dom_oce ! ocean space and time domain 20 USE phycst ! physical constants 23 USE sbc_oce ! Surface boundary condition: ocean fields 24 USE sbc_ice ! Surface boundary condition: ice fields 25 USE ice_oce ! Shared variables between ice and ocean 26 #if defined key_lim3 27 USE par_ice ! ice parameters 28 #endif 29 USE cpl_oasis3 ! OASIS3 coupling 30 USE geo2ocean ! 31 USE restart ! 32 USE oce , ONLY : tn, un, vn 33 USE phycst, ONLY : rt0 34 USE albedo ! 21 35 USE in_out_manager ! I/O manager 36 USE iom ! NetCDF library 22 37 USE lib_mpp ! distribued memory computing library 23 38 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 24 USE daymod ! calendar 25 26 USE cpl_oasis3 ! OASIS3 coupling (to ECHAM5) 27 USE cpl_oasis4 ! OASIS4 coupling (to ECHAM5) 28 USE geo2ocean, ONLY : repere, repcmo 29 USE ice_2, only : frld ! : leads fraction = 1-a/totalarea 30 31 USE sbc_oce ! Surface boundary condition: ocean fields 32 33 USE iom ! NetCDF library 34 39 USE mod_prism_proto ! OASIS3 prism module: PRISM_* variables... 35 40 IMPLICIT NONE 36 41 PRIVATE 37 42 38 PUBLIC sbc_cpl ! routine called by step.F90 39 40 LOGICAL, PUBLIC :: lk_sbc_cpl = .TRUE. !: coupled formulation flag 41 42 INTEGER , PARAMETER :: jpfld = 5 ! maximum number of files to read 43 INTEGER , PARAMETER :: jp_taux = 1 ! index of wind stress (i-component) file 44 INTEGER , PARAMETER :: jp_tauy = 2 ! index of wind stress (j-component) file 45 INTEGER , PARAMETER :: jp_qtot = 3 ! index of total (non solar+solar) heat file 46 INTEGER , PARAMETER :: jp_qsr = 4 ! index of solar heat file 47 INTEGER , PARAMETER :: jp_emp = 5 ! index of evaporation-precipation file 43 PUBLIC sbc_cpl_rcv ! routine called by sbc_ice_lim(_2).F90 44 PUBLIC sbc_cpl_snd ! routine called by step.F90 45 PUBLIC sbc_cpl_ice_tau ! routine called by sbc_ice_lim(_2).F90 46 PUBLIC sbc_cpl_ice_flx ! routine called by sbc_ice_lim(_2).F90 48 47 49 !!wonsun 50 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: &51 taux, tauy & !: surface stress components in (i,j) referential52 53 54 USE sbc_ice, only : dqns_ice , & ! : derivative of non solar heat flux on sea ice55 qsr_ice , & ! : solar flux over ice56 qns_ice , & ! : total non solar heat flux (Longwave downward radiation) over ice57 tn_ice , & ! : ice surface temperature58 alb_ice , & ! : albedo of ice59 sprecip , & ! : solid (snow) precipitation over water (!) what about ice?60 tprecip , & ! : total precipitation ( or liquid precip minus evaporation in coupled mode)61 calving , & ! : calving62 rrunoff , & ! : monthly runoff (kg/m2/s)63 fr1_i0 , & ! : 1st part of the fraction of sol.rad. which penetrate inside the ice cover64 fr2_i0 ! : 2nd part of the fraction of sol.rad. which penetrate inside the ice cover65 66 USE ice_2, only : hicif , & ! : ice thickness67 frld , & ! : leads fraction = 1-a/totalarea68 hsnif , & ! : snow thickness69 u_ice , v_ice ! : ice velocity70 71 USE sbc_oce, only : sst_m ! : sea surface temperature72 73 REAL(wp), PUBLIC :: & !!! surface fluxes namelist (namflx)74 q0 = 0.e0, & ! net heat flux75 qsr0 = 0.e0, & ! solar heat flux76 emp0 = 0.e0, & ! net freshwater flux77 dqdt0 = -40., & ! coefficient for SST damping (W/m2/K)78 deds0 = 27.7 ! coefficient for SSS damping (mm/day)48 INTEGER, PARAMETER :: jpr_otx1 = 1 ! 3 atmosphere-ocean stress components on grid 1 49 INTEGER, PARAMETER :: jpr_oty1 = 2 ! 50 INTEGER, PARAMETER :: jpr_otz1 = 3 ! 51 INTEGER, PARAMETER :: jpr_otx2 = 4 ! 3 atmosphere-ocean stress components on grid 2 52 INTEGER, PARAMETER :: jpr_oty2 = 5 ! 53 INTEGER, PARAMETER :: jpr_otz2 = 6 ! 54 INTEGER, PARAMETER :: jpr_itx1 = 7 ! 3 atmosphere-ice stress components on grid 1 55 INTEGER, PARAMETER :: jpr_ity1 = 8 ! 56 INTEGER, PARAMETER :: jpr_itz1 = 9 ! 57 INTEGER, PARAMETER :: jpr_itx2 = 10 ! 3 atmosphere-ice stress components on grid 2 58 INTEGER, PARAMETER :: jpr_ity2 = 11 ! 59 INTEGER, PARAMETER :: jpr_itz2 = 12 ! 60 INTEGER, PARAMETER :: jpr_qsroce = 13 ! Qsr above the ocean 61 INTEGER, PARAMETER :: jpr_qsrice = 14 ! Qsr above the ice 62 INTEGER, PARAMETER :: jpr_qsrmix = jpr_qsroce ! Qsr above ocean+ice 63 INTEGER, PARAMETER :: jpr_qnsoce = 15 ! Qns above the ocean 64 INTEGER, PARAMETER :: jpr_qnsice = 16 ! Qns above the ice 65 INTEGER, PARAMETER :: jpr_qnsmix = jpr_qnsoce ! Qns above ocean+ice 66 INTEGER, PARAMETER :: jpr_rain = 17 ! total liquid precipitation (rain) 67 INTEGER, PARAMETER :: jpr_snow = 18 ! solid precipitation over the ocean (snow) 68 INTEGER, PARAMETER :: jpr_tevp = 19 ! total evaporation 69 INTEGER, PARAMETER :: jpr_ievp = 20 ! solid evaporation (sublimation) 70 INTEGER, PARAMETER :: jpr_prsb = 21 ! total precipitation (liquid + solid) 71 INTEGER, PARAMETER :: jpr_semp = 22 ! solid freshwater budget (sublimation - snow) 72 INTEGER, PARAMETER :: jpr_oemp = 23 ! ocean freshwater budget (evap - precip) 73 INTEGER, PARAMETER :: jpr_w10m = 24 ! 74 INTEGER, PARAMETER :: jpr_dqnsdt = 25 ! 75 INTEGER, PARAMETER :: jpr_rnf = 26 ! 76 INTEGER, PARAMETER :: jpr_cal = 27 ! 77 INTEGER, PARAMETER :: jprcv = 27 ! total number of fields recieved 79 78 80 REAL(wp), DIMENSION(jpi,jpj) :: qsr_oce_recv , qsr_ice_recv 81 REAL(wp), DIMENSION(jpi,jpj) :: qns_oce_recv, qns_ice_recv 82 REAL(wp), DIMENSION(jpi,jpj) :: dqns_ice_recv 83 REAL(wp), DIMENSION(jpi,jpj) :: tprecip_recv , precip_recv 84 REAL(wp), DIMENSION(jpi,jpj) :: fr1_i0_recv , fr2_i0_recv 85 REAL(wp), DIMENSION(jpi,jpj) :: rrunoff_recv , calving_recv 86 #if defined key_cpl_ocevel 87 REAL(wp), DIMENSION(jpi,jpj) :: un_weighted, vn_weighted 88 REAL(wp), DIMENSION(jpi,jpj) :: un_send , vn_send 89 #endif 90 REAL(wp), DIMENSION(jpi,jpj) :: zrunriv ! river discharge into ocean 91 REAL(wp), DIMENSION(jpi,jpj) :: zruncot ! continental discharge into ocean 92 93 REAL(wp), DIMENSION(jpi,jpj) :: zpew ! P-E over water 94 REAL(wp), DIMENSION(jpi,jpj) :: zpei ! P-E over ice 95 REAL(wp), DIMENSION(jpi,jpj) :: zpsol ! surface downward snow fall 96 REAL(wp), DIMENSION(jpi,jpj) :: zevice ! surface upward snow flux where sea ice 97 !!wonsun 98 99 !! * Substitutions 100 # include "domzgr_substitute.h90" 79 INTEGER, PARAMETER :: jps_fice = 1 ! ice fraction 80 INTEGER, PARAMETER :: jps_toce = 2 ! ocean temperature 81 INTEGER, PARAMETER :: jps_tice = 3 ! ice temperature 82 INTEGER, PARAMETER :: jps_tmix = 4 ! mixed temperature (ocean+ice) 83 INTEGER, PARAMETER :: jps_albice = 5 ! ice albedo 84 INTEGER, PARAMETER :: jps_albmix = 6 ! mixed albedo 85 INTEGER, PARAMETER :: jps_hice = 7 ! ice thickness 86 INTEGER, PARAMETER :: jps_hsnw = 8 ! snow thickness 87 INTEGER, PARAMETER :: jps_ocx1 = 9 ! ocean current on grid 1 88 INTEGER, PARAMETER :: jps_ocy1 = 10 ! 89 INTEGER, PARAMETER :: jps_ocz1 = 11 ! 90 INTEGER, PARAMETER :: jps_ivx1 = 12 ! ice current on grid 1 91 INTEGER, PARAMETER :: jps_ivy1 = 13 ! 92 INTEGER, PARAMETER :: jps_ivz1 = 14 ! 93 INTEGER, PARAMETER :: jpsnd = 14 ! total number of fields sended 94 95 ! !!** namelist namsbc_cpl ** 96 ! Send to the atmosphere ! 97 CHARACTER(len=100) :: cn_snd_temperature = 'oce only' ! 'oce only' 'weighted oce and ice' or 'mixed oce-ice' 98 CHARACTER(len=100) :: cn_snd_albedo = 'none' ! 'none' 'weighted ice' or 'mixed oce-ice' 99 CHARACTER(len=100) :: cn_snd_thickness = 'none' ! 'none' or 'weighted ice and snow' 100 CHARACTER(len=100) :: cn_snd_crt_nature = 'none' ! 'none' 'oce only' 'weighted oce and ice' or 'mixed oce-ice' 101 CHARACTER(len=100) :: cn_snd_crt_refere = 'spherical' ! 'spherical' or 'cartesian' 102 CHARACTER(len=100) :: cn_snd_crt_orient = 'local grid' ! 'eastward-northward' or 'local grid' 103 CHARACTER(len=100) :: cn_snd_crt_grid = 'T' ! always at 'T' point 104 105 ! Recieved from the atmosphere ! 106 CHARACTER(len=100) :: cn_rcv_tau_nature = 'oce only' ! 'oce only' 'oce and ice' or 'mixed oce-ice' 107 CHARACTER(len=100) :: cn_rcv_tau_refere = 'spherical' ! 'spherical' or 'cartesian' 108 CHARACTER(len=100) :: cn_rcv_tau_orient = 'local grid' ! 'eastward-northward' or 'local grid' 109 CHARACTER(len=100) :: cn_rcv_tau_grid = 'T' ! 'T', 'U,V', 'U,V,I', 'T,I', or 'T,U,V' 110 CHARACTER(len=100) :: cn_rcv_w10m = 'none' ! 'none' or 'coupled' 111 CHARACTER(len=100) :: cn_rcv_dqnsdt = 'none' ! 'none' or 'coupled' 112 CHARACTER(len=100) :: cn_rcv_qsr = 'oce only' ! 'oce only' 'conservative' 'oce and ice' or 'mixed oce-ice' 113 CHARACTER(len=100) :: cn_rcv_qns = 'oce only' ! 'oce only' 'conservative' 'oce and ice' or 'mixed oce-ice' 114 CHARACTER(len=100) :: cn_rcv_emp = 'oce only' ! 'oce only' 'conservative' 'oce and ice' or 'mixed oce-ice' 115 CHARACTER(len=100) :: cn_rcv_rnf = 'coupled' ! 'coupled' 'climato' or 'mixed' 116 CHARACTER(len=100) :: cn_rcv_cal = 'none' ! 'none' or 'coupled' 117 118 !! CHARACTER(len=100), PUBLIC :: cn_rcv_rnf !: ??? ==>> !!gm treat this case in a different maner 119 120 CHARACTER(len=100), DIMENSION(4) :: cn_snd_crt ! array combining cn_snd_crt_* 121 CHARACTER(len=100), DIMENSION(4) :: cn_rcv_tau ! array combining cn_rcv_tau_* 122 123 REAL(wp), DIMENSION(jpi,jpj) :: albedo_oce_mix ! ocean albedo sent to atmosphere (mix clear/overcast sky) 124 125 REAL(wp), DIMENSION(jpi,jpj,jprcv) :: frcv ! all fields recieved from the atmosphere 126 INTEGER , DIMENSION( jprcv) :: nrcvinfo ! OASIS info argument 127 128 !! Substitution 129 # include "vectopt_loop_substitute.h90" 101 130 !!---------------------------------------------------------------------- 102 !! OPA 9.0 , LOCEAN-IPSL (2006)103 !! $Id $131 !! NEMO/OPA 3.0 , LOCEAN-IPSL (2008) 132 !! $Id:$ 104 133 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 105 134 !!---------------------------------------------------------------------- 106 135 107 136 CONTAINS 108 109 SUBROUTINE sbc_cpl( kt ) 137 138 SUBROUTINE sbc_cpl_init( k_ice ) 139 !!---------------------------------------------------------------------- 140 !! *** ROUTINE sbc_cpl_init *** 141 !! 142 !! ** Purpose : Initialisation of send and recieved information from 143 !! the atmospheric component 144 !! 145 !! ** Method : * Read namsbc_cpl namelist 146 !! * define the receive interface 147 !! * define the send interface 148 !! * initialise the OASIS coupler 149 !!---------------------------------------------------------------------- 150 INTEGER, INTENT(in) :: k_ice ! ice management in the sbc (=0/1/2/3) 151 !! 152 INTEGER :: jn ! dummy loop index 153 REAL(wp), DIMENSION(jpi,jpj) :: zacs, zaos ! 2D workspace (clear & overcast sky albedos) 154 !! 155 NAMELIST/namsbc_cpl/ cn_snd_temperature, cn_snd_albedo , cn_snd_thickness, & 156 cn_snd_crt_nature, cn_snd_crt_refere , cn_snd_crt_orient, cn_snd_crt_grid , cn_rcv_w10m, & 157 cn_rcv_tau_nature, cn_rcv_tau_refere , cn_rcv_tau_orient, cn_rcv_tau_grid , & 158 cn_rcv_dqnsdt , cn_rcv_qsr , cn_rcv_qns , cn_rcv_emp , cn_rcv_rnf , cn_rcv_cal 110 159 !!--------------------------------------------------------------------- 111 !! *** ROUTINE sbc_cpl *** 112 !! 113 !! ** Purpose : provide at each time step the surface ocean fluxes 114 !! (momentum, heat, freshwater and runoff) in coupled mode 115 !! 116 !! ** Method : - Recieve from a Atmospheric model via OASIS coupler : 117 !! i-component of the stress taux (N/m2) 118 !! j-component of the stress tauy (N/m2) 119 !! net downward heat flux qtot (watt/m2) 120 !! net downward radiative flux qsr (watt/m2) 121 !! net upward freshwater (evapo - precip) emp (kg/m2/s) 122 !! - send to the Atmospheric model via OASIS coupler : 123 !! 124 !! ** Action : update at each time-step the two components of the 125 !! surface stress in both (i,j) and geographical ref. 126 !! 127 !! 128 !! CAUTION : - never mask the surface stress fields 129 !! 130 !! ** Action : update at each time-step 131 !! - taux & tauy : stress components in (i,j) referential 132 !! - qns : non solar heat flux 133 !! - qsr : solar heat flux 134 !! - emp : evap - precip (volume flux) 135 !! - emps : evap - precip (concentration/dillution) 136 !! 137 !! References : The OASIS User Guide, Version 3.0 and 4.0 138 !!---------------------------------------------------------------------- 139 INTEGER, INTENT(in) :: kt ! ocean time step 140 !! 141 INTEGER :: ji, jj ! dummy loop indices 142 #if defined key_cpl_ocevel 143 INTEGER :: ikchoix 144 #endif 145 INTEGER :: var_id, info 146 INTEGER :: date !???? !!gm bug this is a real !!! 147 REAL(wp) :: zfacflx, zfacwat, zfact 148 149 REAL(wp), DIMENSION(jpi,jpj) :: ztaueuw, ztauevw ! eastward wind stress over water at U and V-points 150 REAL(wp), DIMENSION(jpi,jpj) :: ztaunuw, ztaunvw ! northward wind stress over water at U and V-points 151 REAL(wp), DIMENSION(jpi,jpj) :: ztaueui, ztauevi ! eastward wind stress over ice at U and V-points 152 REAL(wp), DIMENSION(jpi,jpj) :: ztaunui, ztaunvi ! northward wind stress over ice at U and V-points 153 REAL(wp), DIMENSION(jpi,jpj) :: ztaueu , ztauev ! eastward wind stress combined 154 REAL(wp), DIMENSION(jpi,jpj) :: ztaunu , ztaunv ! northward wind stress combined 155 !!--------------------------------------------------------------------- 156 157 date = ( kt - nit000 ) * rdttra(1) ! date of exxhanges 158 ! ! Conversion factor (ocean units are W/m2 and Kg/m2/s] 159 zfacflx = 1.e0 ! no conversion [W/m2] ! W/m2 heat fluxes are send by the Atmosphere 160 zfacwat = 1.e3 ! convert [m/s] to [kg/m2/s] ! m/s freshwater fluxes are send by the atmosphere 161 162 163 ! ! =========================== ! 164 ! ! Send Coupling fields ! 165 ! ! =========================== ! 160 161 ! ================================ ! 162 ! Namelist informations ! 163 ! ================================ ! 164 165 REWIND( numnam ) ! ... read namlist namsbc_cpl 166 READ ( numnam, namsbc_cpl ) 167 168 IF(lwp) THEN ! control print 169 WRITE(numout,*) 170 WRITE(numout,*)'sbc_cpl_init : namsbc_cpl namelist ' 171 WRITE(numout,*)'~~~~~~~~~~~~' 172 WRITE(numout,*)' received fields' 173 WRITE(numout,*)' 10m wind module cn_rcv_w10m = ', cn_rcv_w10m 174 WRITE(numout,*)' surface stress - nature cn_rcv_tau_nature = ', cn_rcv_tau_nature 175 WRITE(numout,*)' - referential cn_rcv_tau_refere = ', cn_rcv_tau_refere 176 WRITE(numout,*)' - orientation cn_rcv_tau_orient = ', cn_rcv_tau_orient 177 WRITE(numout,*)' - mesh cn_rcv_tau_grid = ', cn_rcv_tau_grid 178 WRITE(numout,*)' non-solar heat flux sensitivity cn_rcv_dqnsdt = ', cn_rcv_dqnsdt 179 WRITE(numout,*)' solar heat flux cn_rcv_qsr = ', cn_rcv_qsr 180 WRITE(numout,*)' non-solar heat flux cn_rcv_qns = ', cn_rcv_qns 181 WRITE(numout,*)' freshwater budget cn_rcv_emp = ', cn_rcv_emp 182 WRITE(numout,*)' runoffs cn_rcv_rnf = ', cn_rcv_rnf 183 WRITE(numout,*)' calving cn_rcv_cal = ', cn_rcv_cal 184 WRITE(numout,*)' sent fields' 185 WRITE(numout,*)' surface temperature cn_snd_temperature = ', cn_snd_temperature 186 WRITE(numout,*)' albedo cn_snd_albedo = ', cn_snd_albedo 187 WRITE(numout,*)' ice/snow thickness cn_snd_thickness = ', cn_snd_thickness 188 WRITE(numout,*)' surface current - nature cn_snd_crt_nature = ', cn_snd_crt_nature 189 WRITE(numout,*)' - referential cn_snd_crt_refere = ', cn_snd_crt_refere 190 WRITE(numout,*)' - orientation cn_snd_crt_orient = ', cn_snd_crt_orient 191 WRITE(numout,*)' - mesh cn_snd_crt_grid = ', cn_snd_crt_grid 192 ENDIF 193 194 ! save current & stress in an array and suppress possible blank in the name 195 cn_snd_crt(1) = TRIM( cn_snd_crt_nature ) ; cn_snd_crt(2) = TRIM( cn_snd_crt_refere ) 196 cn_snd_crt(3) = TRIM( cn_snd_crt_orient ) ; cn_snd_crt(4) = TRIM( cn_snd_crt_grid ) 197 cn_rcv_tau(1) = TRIM( cn_rcv_tau_nature ) ; cn_rcv_tau(2) = TRIM( cn_rcv_tau_refere ) 198 cn_rcv_tau(3) = TRIM( cn_rcv_tau_orient ) ; cn_rcv_tau(4) = TRIM( cn_rcv_tau_grid ) 199 200 ! ================================ ! 201 ! Define the receive interface ! 202 ! ================================ ! 203 nrcvinfo(:) = PRISM_NotDef ! needed by nrcvinfo(jpr_otx1) if we do not receive ocea stress 204 205 ! for each field: define the OASIS name (srcv(:)%clname) 206 ! define receive or not from the namelist parameters (srcv(:)%laction) 207 ! define the north fold type of lbc (srcv(:)%nsgn) 208 209 ! default definitions of srcv 210 srcv(:)%laction = .FALSE. ; srcv(:)%clgrid = 'T' ; srcv(:)%nsgn = 1 211 212 ! ! ------------------------- ! 213 ! ! ice and ocean wind stress ! 214 ! ! ------------------------- ! 215 ! ! Name 216 srcv(jpr_otx1)%clname = 'O_OTaux1' ! 1st ocean component on grid ONE (T or U) 217 srcv(jpr_oty1)%clname = 'O_OTauy1' ! 2nd - - - - 218 srcv(jpr_otz1)%clname = 'O_OTauz1' ! 3rd - - - - 219 srcv(jpr_otx2)%clname = 'O_OTaux2' ! 1st ocean component on grid TWO (V) 220 srcv(jpr_oty2)%clname = 'O_OTauy2' ! 2nd - - - - 221 srcv(jpr_otz2)%clname = 'O_OTauz2' ! 3rd - - - - 222 ! 223 srcv(jpr_itx1)%clname = 'O_ITaux1' ! 1st ice component on grid ONE (T, F, I or U) 224 srcv(jpr_ity1)%clname = 'O_ITauy1' ! 2nd - - - - 225 srcv(jpr_itz1)%clname = 'O_ITauz1' ! 3rd - - - - 226 srcv(jpr_itx2)%clname = 'O_ITaux2' ! 1st ice component on grid TWO (V) 227 srcv(jpr_ity2)%clname = 'O_ITauy2' ! 2nd - - - - 228 srcv(jpr_itz2)%clname = 'O_ITauz2' ! 3rd - - - - 166 229 ! 167 !!gm bug ? here send instantaneous SST, not mean over the coupling period.... 168 var_id = send_id(1) ; CALL cpl_prism_send( var_id, date, tn(:,:,1)+rt0, info ) ! ocean surface temperature [K] 169 var_id = send_id(2) ; CALL cpl_prism_send( var_id, date, 1.0-frld , info ) ! fraction of ice-cover 170 #if defined key_cpl_albedo 171 DO jj = 1, jpj 172 DO ji = 1, jpi 173 IF( ( tn_ice(ji,jj) < 50 .OR. tn_ice(ji,jj) > 400 ) .AND. frld(ji,jj) < 1. ) THEN 174 WRITE(numout,*) ' tn_ice, ERROR ', ji, jj, ' = ', tn_ice(ji,jj), & 175 & ' qns_ice_recv=', qns_ice_recv(ji,jj), ' dqns_ice_recv=', dqns_ice_recv(ji,jj) 230 srcv(jpr_otx1:jpr_itz2)%nsgn = -1 ! Vectors: change of sign at north fold 231 232 ! ! Set grid and action 233 SELECT CASE( TRIM( cn_rcv_tau(4) ) ) ! 'T', 'U,V', 'U,V,I', 'U,V,F', 'T,I', 'T,F', or 'T,U,V' 234 CASE( 'T' ) 235 srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point 236 srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 237 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 238 CASE( 'U,V' ) 239 srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point 240 srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point 241 srcv(jpr_itx1:jpr_itz1)%clgrid = 'U' ! ice components given at U-point 242 srcv(jpr_itx2:jpr_itz2)%clgrid = 'V' ! and V-point 243 srcv(jpr_otx1:jpr_itz2)%laction = .TRUE. ! receive oce and ice components on both grid 1 & 2 244 CASE( 'U,V,T' ) 245 srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point 246 srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point 247 srcv(jpr_itx1:jpr_itz1)%clgrid = 'T' ! ice components given at T-point 248 srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 249 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only 250 CASE( 'U,V,I' ) 251 srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point 252 srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point 253 srcv(jpr_itx1:jpr_itz1)%clgrid = 'I' ! ice components given at I-point 254 srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 255 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only 256 CASE( 'U,V,F' ) 257 srcv(jpr_otx1:jpr_otz1)%clgrid = 'U' ! oce components given at U-point 258 srcv(jpr_otx2:jpr_otz2)%clgrid = 'V' ! and V-point 259 srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point 260 srcv(jpr_otx1:jpr_otz2)%laction = .TRUE. ! receive oce components on grid 1 & 2 261 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 only 262 CASE( 'T,I' ) 263 srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point 264 srcv(jpr_itx1:jpr_itz1)%clgrid = 'I' ! ice components given at I-point 265 srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 266 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 267 CASE( 'T,F' ) 268 srcv(jpr_otx1:jpr_itz2)%clgrid = 'T' ! oce and ice components given at T-point 269 srcv(jpr_itx1:jpr_itz1)%clgrid = 'F' ! ice components given at F-point 270 srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 271 srcv(jpr_itx1:jpr_itz1)%laction = .TRUE. ! receive ice components on grid 1 272 CASE( 'T,U,V' ) 273 srcv(jpr_otx1:jpr_otz1)%clgrid = 'T' ! oce components given at T-point 274 srcv(jpr_itx1:jpr_itz1)%clgrid = 'U' ! ice components given at U-point 275 srcv(jpr_itx2:jpr_itz2)%clgrid = 'V' ! and V-point 276 srcv(jpr_otx1:jpr_otz1)%laction = .TRUE. ! receive oce components on grid 1 only 277 srcv(jpr_itx1:jpr_itz2)%laction = .TRUE. ! receive ice components on grid 1 & 2 278 CASE default 279 CALL ctl_stop( 'sbc_cpl_init: wrong definition of cn_rcv_tau(4)' ) 280 END SELECT 281 ! 282 IF( TRIM( cn_rcv_tau(2) ) == 'spherical' ) & ! spherical: 3rd component not received 283 & srcv( (/jpr_otz1, jpr_otz2, jpr_itz1, jpr_itz2/) )%laction = .FALSE. 284 ! 285 IF( TRIM( cn_rcv_tau(1) ) /= 'oce and ice' ) THEN ! 'oce and ice' case ocean stress on ocean mesh used 286 srcv(jpr_itx1:jpr_itz2)%laction = .FALSE. ! ice components not received 287 srcv(jpr_itx1)%clgrid = 'U' ! ocean stress used after its transformation 288 srcv(jpr_ity1)%clgrid = 'V' ! i.e. it is always at U- & V-points for i- & j-comp. resp. 289 ENDIF 290 291 ! ! ------------------------- ! 292 ! ! freshwater budget ! E-P 293 ! ! ------------------------- ! 294 ! we suppose that atmosphere modele do not make the difference between precipiration (liquide or solid) 295 ! over ice of free ocean within the same atmospheric cell.cd 296 srcv(jpr_rain)%clname = 'OTotRain' ! Rain = liquid precipitation 297 srcv(jpr_snow)%clname = 'OTotSnow' ! Snow = solid precipitation 298 srcv(jpr_tevp)%clname = 'OTotEvap' ! total evaporation (over oce + ice sublimation) 299 srcv(jpr_ievp)%clname = 'OIceEvap' ! evaporation over ice = sublimation 300 srcv(jpr_prsb)%clname = 'OPre-Sub' ! liquid precipitation + solid precipitation - sublimation 301 srcv(jpr_semp)%clname = 'OISub-Sn' ! ice solid water budget = sublimation - solid precipitation 302 srcv(jpr_oemp)%clname = 'OOEva-Pr' ! ocean water budget = ocean Evap - ocean precip 303 SELECT CASE( TRIM( cn_rcv_emp ) ) 304 CASE( 'oce only' ) ; srcv( jpr_oemp )%laction = .TRUE. 305 CASE( 'conservative' ) ; srcv( (/jpr_rain, jpr_snow, jpr_ievp, jpr_tevp/) )%laction = .TRUE. 306 CASE( 'oce and ice' ) ; srcv( (/ jpr_prsb, jpr_semp, jpr_oemp/) )%laction = .TRUE. 307 CASE( 'mixed oce-ice' ) ; srcv( (/jpr_rain, jpr_semp, jpr_tevp/) )%laction = .TRUE. 308 CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of cn_rcv_emp' ) 309 END SELECT 310 311 ! ! ------------------------- ! 312 ! ! Runoffs & Calving ! 313 ! ! ------------------------- ! 314 srcv(jpr_rnf )%clname = 'O_Runoff' ; IF( TRIM( cn_rcv_rnf ) == 'coupled' ) srcv(jpr_rnf)%laction = .TRUE. 315 IF( TRIM( cn_rcv_rnf ) == 'climato' ) THEN ; ln_rnf = .TRUE. 316 ELSE ; ln_rnf = .FALSE. 317 ENDIF 318 srcv(jpr_cal )%clname = 'OCalving' ; IF( TRIM( cn_rcv_cal ) == 'coupled' ) srcv(jpr_cal)%laction = .TRUE. 319 320 ! ! ------------------------- ! 321 ! ! non solar radiation ! Qns 322 ! ! ------------------------- ! 323 srcv(jpr_qnsoce)%clname = 'O_QnsOce' 324 srcv(jpr_qnsice)%clname = 'O_QnsIce' 325 srcv(jpr_qnsmix)%clname = 'O_QnsMix' 326 SELECT CASE( TRIM( cn_rcv_qns ) ) 327 CASE( 'oce only' ) ; srcv( jpr_qnsoce )%laction = .TRUE. 328 CASE( 'conservative' ) ; srcv( (/jpr_qnsice, jpr_qnsmix/) )%laction = .TRUE. 329 CASE( 'oce and ice' ) ; srcv( (/jpr_qnsice, jpr_qnsoce/) )%laction = .TRUE. 330 CASE( 'mixed oce-ice' ) ; srcv( jpr_qnsmix )%laction = .TRUE. 331 CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of cn_rcv_qns' ) 332 END SELECT 333 334 ! ! ------------------------- ! 335 ! ! solar radiation ! Qsr 336 ! ! ------------------------- ! 337 srcv(jpr_qsroce)%clname = 'O_QsrOce' 338 srcv(jpr_qsrice)%clname = 'O_QsrIce' 339 srcv(jpr_qsrmix)%clname = 'O_QsrMix' 340 SELECT CASE( TRIM( cn_rcv_qsr ) ) 341 CASE( 'oce only' ) ; srcv( jpr_qsroce )%laction = .TRUE. 342 CASE( 'conservative' ) ; srcv( (/jpr_qsrice, jpr_qsrmix/) )%laction = .TRUE. 343 CASE( 'oce and ice' ) ; srcv( (/jpr_qsrice, jpr_qsroce/) )%laction = .TRUE. 344 CASE( 'mixed oce-ice' ) ; srcv( jpr_qsrmix )%laction = .TRUE. 345 CASE default ; CALL ctl_stop( 'sbc_cpl_init: wrong definition of cn_rcv_qsr' ) 346 END SELECT 347 348 ! ! ------------------------- ! 349 ! ! non solar sensitivity ! d(Qns)/d(T) 350 ! ! ------------------------- ! 351 srcv(jpr_dqnsdt)%clname = 'O_dQnsdT' 352 IF( TRIM( cn_rcv_dqnsdt ) == 'coupled' ) srcv(jpr_dqnsdt)%laction = .TRUE. 353 354 ! ! ------------------------- ! 355 ! ! Ice Qsr penetration ! 356 ! ! ------------------------- ! 357 ! fraction of net shortwave radiation which is not absorbed in the thin surface layer 358 ! and penetrates inside the ice cover ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 359 ! Coupled case: since cloud cover is not received from atmosphere 360 ! ===> defined as constant value -> definition done in sbc_cpl_init 361 fr1_i0(:,:) = 0.18 362 fr2_i0(:,:) = 0.82 363 ! ! ------------------------- ! 364 ! ! 10m wind module ! 365 ! ! ------------------------- ! 366 srcv(jpr_w10m )%clname = 'O_Wind10' ; IF( TRIM(cn_rcv_w10m) == 'coupled' ) srcv(jpr_w10m)%laction = .TRUE. 367 ! ! +++ ---> A brancher et a blinder dans tke si cn_rcv_w10m == 'none' 368 369 370 ! ================================ ! 371 ! Define the send interface ! 372 ! ================================ ! 373 ! for each field: define the OASIS name (srcv(:)%clname) 374 ! define send or not from the namelist parameters (srcv(:)%laction) 375 ! define the north fold type of lbc (srcv(:)%nsgn) 376 377 ! default definitions of nsnd 378 ssnd(:)%laction = .FALSE. ; ssnd(:)%clgrid = 'T' ; ssnd(:)%nsgn = 1 379 380 ! ! ------------------------- ! 381 ! ! Surface temperature ! 382 ! ! ------------------------- ! 383 ssnd(jps_toce)%clname = 'O_SSTSST' 384 ssnd(jps_tice)%clname = 'O_TepIce' 385 ssnd(jps_tmix)%clname = 'O_TepMix' 386 SELECT CASE( TRIM( cn_snd_temperature ) ) 387 CASE( 'oce only' ) ; ssnd( jps_toce )%laction = .TRUE. 388 CASE( 'weighted oce and ice' ) ; ssnd( (/jps_toce, jps_tice/) )%laction = .TRUE. 389 CASE( 'mixed oce-ice' ) ; ssnd( jps_tmix )%laction = .TRUE. 390 END SELECT 391 392 ! ! ------------------------- ! 393 ! ! Albedo ! 394 ! ! ------------------------- ! 395 ssnd(jps_albice)%clname = 'O_AlbIce' 396 ssnd(jps_albmix)%clname = 'O_AlbMix' 397 SELECT CASE( TRIM( cn_snd_albedo ) ) 398 CASE( 'none' ) ! nothing to do 399 CASE( 'weighted ice' ) ; ssnd(jps_albice)%laction = .TRUE. 400 CASE( 'mixed oce-ice' ) ; ssnd(jps_albmix)%laction = .TRUE. 401 CALL albedo_oce( zaos, zacs ) 402 ! Due to lack of information on nebulosity : mean clear/overcast sky 403 albedo_oce_mix(:,:) = ( zacs(:,:) + zaos(:,:) ) * 0.5 404 END SELECT 405 406 ! ! ------------------------- ! 407 ! ! Ice fraction & Thickness ! 408 ! ! ------------------------- ! 409 ssnd(jps_fice)%clname = 'OIceFrac' 410 ssnd(jps_hice)%clname = 'O_IceTck' 411 ssnd(jps_hsnw)%clname = 'O_SnwTck' 412 IF( k_ice /= 0 ) ssnd(jps_fice)%laction = .TRUE. ! if ice treated in the ocean (even in climato case) 413 IF( TRIM( cn_snd_thickness ) == 'weighted ice and snow' ) ssnd( (/jps_hice, jps_hsnw/) )%laction = .TRUE. 414 415 ! ! ------------------------- ! 416 ! ! Surface current ! 417 ! ! ------------------------- ! 418 ! ocean currents ! ice velocities 419 ssnd(jps_ocx1)%clname = 'O_OCurx1' ; ssnd(jps_ivx1)%clname = 'O_IVelx1' 420 ssnd(jps_ocy1)%clname = 'O_OCury1' ; ssnd(jps_ivy1)%clname = 'O_IVely1' 421 ssnd(jps_ocz1)%clname = 'O_OCurz1' ; ssnd(jps_ivz1)%clname = 'O_IVelz1' 422 ! 423 ssnd(jps_ocx1:jps_ivz2)%nsgn = -1 ! vectors: change of the sign at the north fold 424 425 IF( cn_snd_crt(4) /= 'T' ) CALL ctl_stop( 'cn_snd_crt(4) must be equal to T' ) 426 ssnd(jps_ocx1:jps_ivz2)%clgrid = 'T' ! all oce and ice components on the same unique grid 427 ssnd(jps_ocx1:jps_ocz1)%laction = .TRUE. ! oce components on 1 grid 428 ssnd(jps_ivx1:jps_ivz1)%laction = .TRUE. ! ice components on 1 grid 429 430 IF( TRIM( cn_snd_crt(2) ) == 'spherical' ) & ! 3rd component not used 431 & srcv( (/jps_otz1, jps_otz2, jps_itz1, jps_itz2/) )%laction = .FALSE. 432 ! 433 IF( TRIM( cn_snd_crt(1) ) /= 'oce only' .OR. 'oce and ice' ) & ! ice components not used 434 & srcv(jps_itx1:jps_itz2)%laction = FALSE. 435 436 SELECT CASE( TRIM( cn_snd_crt(1) ) ) 437 CASE( 'none' ) ; ssnd(jps_ocx1:jps_ivz2)%laction = .FALSE. 438 CASE( 'oce only' ) ; ssnd(jps_ivx1:jps_ivz2)%laction = .FALSE. 439 CASE( 'weighted oce and ice' ) ! nothing to do 440 CASE( 'mixed oce-ice' ) ; ssnd(jps_ivx1:jps_ivz2)%laction = .FALSE. 441 END SELECT 442 443 ! ================================ ! 444 ! initialisation of the coupler ! 445 ! ================================ ! 446 CALL cpl_prism_define 447 ! 448 END SUBROUTINE sbc_cpl_init 449 450 451 SUBROUTINE sbc_cpl_rcv( kt, k_fsbc, k_ice ) 452 !!---------------------------------------------------------------------- 453 !! *** ROUTINE sbc_cpl_rcv *** 454 !! 455 !! ** Purpose : provide the stress over the ocean and, if no sea-ice, 456 !! provide the ocean heat and freshwater fluxes. 457 !! 458 !! ** Method : - Receive all the atmospheric fields (stored in frcv array). called at each time step. 459 !! OASIS controls if there is something do receive or not. nrcvinfo contains the info 460 !! to know if the field was really received or not 461 !! 462 !! --> If ocean stress was really received: 463 !! 464 !! - transform the received ocean stress vector from the received 465 !! referential and grid into an atmosphere-ocean stress in 466 !! the (i,j) ocean referencial and at the ocean velocity point. 467 !! The received stress are : 468 !! - defined by 3 components (if cartesian coordinate) 469 !! or by 2 components (if spherical) 470 !! - oriented along geographical coordinate (if eastward-northward) 471 !! or along the local grid coordinate (if local grid) 472 !! - given at U- and V-point, resp. if received on 2 grids 473 !! or at T-point if received on 1 grid 474 !! Therefore and if necessary, they are successively 475 !! processed in order to obtain them 476 !! first as 2 components on the sphere 477 !! second as 2 components oriented along the local grid 478 !! third as 2 components on the U,V grid 479 !! 480 !! --> 481 !! 482 !! - In 'ocean only' case, non solar and solar ocean heat fluxes 483 !! and total ocean freshwater fluxes 484 !! 485 !! ** Method : receive all fields from the atmosphere and transform 486 !! them into ocean surface boundary condition fields 487 !! 488 !! ** Action : update utau, vtau ocean stress at U,V grid 489 !! qns , qsr non solar and solar ocean heat fluxes ('ocean only case) 490 !! emp = emps evap. - precip. (- runoffs) (- calving) ('ocean only case) 491 !! wind10m 10m wind speed !!!!gm to be checked 492 !!---------------------------------------------------------------------- 493 INTEGER, INTENT(in) :: kt ! ocean model time step index 494 INTEGER, INTENT(in) :: k_fsbc ! frequency of sbc (-> ice model) computation 495 INTEGER, INTENT(in) :: k_ice ! ice management in the sbc (=0/1/2/3) 496 !! 497 INTEGER :: ji, jj, jn ! dummy loop indices 498 INTEGER :: isec ! number of seconds since nit000 (assuming rdttra did not change since nit000) 499 REAL(wp) :: zcumulneg, zcumulpos ! temporary scalars 500 REAL(wp), DIMENSION(jpi,jpj) :: ztx, ztx ! 2D workspace 501 !!---------------------------------------------------------------------- 502 503 IF( kt == nit000 ) CALL sbc_cpl_init( k_ice ) ! initialisation 504 505 ! ! Receive all the atmos. fields (including ice information) 506 isec = ( kt - nit000 ) * NINT( rdttra(1) ) ! date of exchanges 507 DO jn = 1, jprcv ! received fields sent by the atmosphere 508 IF( srcv(jn)%laction ) CALL cpl_prism_rcv( jn, isec, frcv(:,:,jn), nrcvinfo(jn) ) 509 END DO 510 511 ! ! ========================= ! 512 IF( srcv(jpr_otx1)%laction ) THEN ! ocean stress ! 513 ! ! ========================= ! 514 ! define frcv(:,:,jpr_otx1) and frcv(:,:,jpr_oty1): stress at U/V point along model grid 515 ! => need to be done only when we receive the field 516 IF( nrcvinfo(jpr_otx1) == PRISM_Recvd .OR. nrcvinfo(jpr_otx1) == PRISM_FromRest .OR. & 517 & nrcvinfo(jpr_otx1) == PRISM_RecvOut .OR. nrcvinfo(jpr_otx1) == PRISM_FromRestOut ) THEN 518 ! 519 IF( TRIM( cn_rcv_tau(2) ) == 'cartesian' ) THEN ! 2 components on the sphere 520 ! ! (cartesian to spherical -> 3 to 2 components) 521 ! 522 CALL geo2oce( frcv(:,:,jpr_otx1), frcv(:,:,jpr_oty1), frcv(:,:,jpr_otz1), & 523 & srcv(jpr_otx1)%clgrid, ztx, zty ) 524 frcv(:,:,jpr_otx1) = ztx(:,:) ! overwrite 1st comp. on the 1st grid 525 frcv(:,:,jpr_oty1) = zty(:,:) ! overwrite 2nd comp. on the 1st grid 526 ! 527 IF( srcv(jpr_otx2)%laction ) THEN 528 CALL geo2oce( frcv(:,:,jpr_otx2), frcv(:,:,jpr_oty2), frcv(:,:,jpr_otz2), & 529 & srcv(jpr_otx2)%clgrid, ztx, zty ) 530 frcv(:,:,jpr_otx2) = ztx(:,:) ! overwrite 1st comp. on the 2nd grid 531 frcv(:,:,jpr_oty2) = zty(:,:) ! overwrite 2nd comp. on the 2nd grid 532 ENDIF 533 ! 176 534 ENDIF 177 END DO 178 END DO 179 var_id = send_id(3) ; CALL cpl_prism_send( var_id, date, tn_ice , info ) ! ice surface temperature [K] 180 var_id = send_id(4) ; CALL cpl_prism_send( var_id, date, alb_ice , info ) ! ice albedo [%] 181 #else 182 var_id = send_id(3) ; CALL cpl_prism_send( var_id, date, hicif , info ) ! ice thickness [m] 183 var_id = send_id(4) ; CALL cpl_prism_send( var_id, date, hsnif , info ) ! snow thickness [m] 184 #endif 185 #if defined key_cpl_ocevel 186 !!gm bug??? I have to check the grid point position... 187 !! a priori there is a error here as un, vn are not at the same grid point.... 188 !! there should be a averaged to set u and v at T-point.... with caution for sea-ice velocity at I-point.... 189 un_weighted = un(:,:,1) * frld + u_ice * ( 1. - frld ) 190 vn_weighted = vn(:,:,1) * frld + v_ice * ( 1. - frld ) 191 ikchoix = - 1 ! converte from (i,j) to geographic referential 192 CALL repere( un_weighted, vn_weighted, un_send, vn_send, ikchoix ) 193 !!gm bug : at lbc_lnk is to be added on un_send and vn_send 194 var_id = send_id(5) ; CALL cpl_prism_send( var_id, date, un_send , info ) ! surface current [m/s] 195 var_id = send_id(6) ; CALL cpl_prism_send( var_id, date, vn_send , info ) ! surface current [m/s] 196 #endif 197 198 ! ! =========================== ! 199 ! ! Recieve Momentum fluxes ! 200 ! ! =========================== ! 201 ! 202 ! ... Receive wind stress fields in geographic component over water and ice 203 var_id = recv_id(1) ; CALL cpl_prism_recv( var_id, date, ztaueuw, info ) ! ??? 204 var_id = recv_id(2) ; CALL cpl_prism_recv( var_id, date, ztaunuw, info ) 205 var_id = recv_id(3) ; CALL cpl_prism_recv( var_id, date, ztaueui, info ) 206 var_id = recv_id(4) ; CALL cpl_prism_recv( var_id, date, ztaunui, info ) 207 var_id = recv_id(5) ; CALL cpl_prism_recv( var_id, date, ztauevw, info ) 208 var_id = recv_id(6) ; CALL cpl_prism_recv( var_id, date, ztaunvw, info ) 209 var_id = recv_id(7) ; CALL cpl_prism_recv( var_id, date, ztauevi, info ) 210 var_id = recv_id(8) ; CALL cpl_prism_recv( var_id, date, ztaunvi, info ) 211 ! 212 !!gm bug : keep separate ice and ocean stress ! 213 ! ... combine water / ice stresses 214 ztaueu(:,:) = ztaueuw(:,:) * frld(:,:) + ztaueui(:,:) * ( 1.0 - frld(:,:) ) 215 ztaunu(:,:) = ztaunuw(:,:) * frld(:,:) + ztaunui(:,:) * ( 1.0 - frld(:,:) ) 216 ztauev(:,:) = ztauevw(:,:) * frld(:,:) + ztauevi(:,:) * ( 1.0 - frld(:,:) ) 217 ztaunv(:,:) = ztaunvw(:,:) * frld(:,:) + ztaunvi(:,:) * ( 1.0 - frld(:,:) ) 218 ! 219 ! ... rotate vector components from geographic to (i,j) referential 220 CALL repcmo ( ztaueu, ztaunu, ztauev, ztaunv, utau, vtau, kt ) 221 ! 222 !!gm bug?? not sure but put that for security 223 CALL lbc_lnk( utau , 'U', -1. ) 224 CALL lbc_lnk( vtau , 'V', -1. ) 225 !!gm end bug?? 226 ! 227 ! ! =========================== ! 228 ! ! Recieve heat fluxes ! 229 ! ! =========================== ! 230 ! 231 var_id = recv_id(13) ; CALL cpl_prism_recv( var_id, date, qsr_oce_recv , info ) ! ocean surface net downward shortwave flux 232 var_id = recv_id(14) ; CALL cpl_prism_recv( var_id, date, qns_oce_recv , info ) ! ocean surface downward non-solar heat flux 233 var_id = recv_id(15) ; CALL cpl_prism_recv( var_id, date, qsr_ice_recv , info ) ! ice solar heat flux 234 var_id = recv_id(16) ; CALL cpl_prism_recv( var_id, date, qns_ice_recv , info ) ! ice non-solar heat flux 235 var_id = recv_id(17) ; CALL cpl_prism_recv( var_id, date, dqns_ice_recv, info ) ! ice non-solar heat flux sensitivity 236 237 qsr_oce_recv (:,:) = qsr_oce_recv (:,:) * tmask(:,:,1) * zfacflx 238 qns_oce_recv (:,:) = qns_oce_recv (:,:) * tmask(:,:,1) * zfacflx 239 qsr_ice_recv (:,:) = qsr_ice_recv (:,:) * tmask(:,:,1) * zfacflx 240 qns_ice_recv (:,:) = qns_ice_recv (:,:) * tmask(:,:,1) * zfacflx 241 dqns_ice_recv(:,:) = dqns_ice_recv(:,:) * tmask(:,:,1) * zfacflx 242 243 IF( kt == nit000 ) THEN ! set once for all qsr penetration in sea-ice 244 ! ! Since cloud cover catm not transmitted from atmosphere, it is set to 0. 245 ! ! i.e. constant penetration fractions of 0.18 and 0.82 246 ! fraction of net shortwave radiation which is not absorbed in the thin surface layer and penetrates 247 ! inside the ice cover ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 248 fr1_i0_recv(:,:) = 0.18 249 fr2_i0_recv(:,:) = 0.82 535 ! 536 IF( TRIM( cn_rcv_tau(3) ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid 537 ! ! (geographical to local grid -> rotate the components) 538 CALL rot_rep( frcv(:,:,jpr_otx1), frcv(:,:,jpr_oty1), srcv(jpr_otx1)%clgrid, 'en->i', ztx ) 539 frcv(:,:,jpr_otx1) = ztx(:,:) ! overwrite 1st component on the 1st grid 540 IF( srcv(jpr_otx2)%laction ) THEN 541 CALL rot_rep( frcv(:,:,jpr_otx2), frcv(:,:,jpr_oty2), srcv(jpr_otx2)%clgrid, 'en->j', zty ) 542 ELSE 543 CALL rot_rep( frcv(:,:,jpr_otx1), frcv(:,:,jpr_oty1), srcv(jpr_otx1)%clgrid, 'en->j', zty ) 544 ENDIF 545 frcv(:,:,jpr_oty1) = zty(:,:) ! overwrite 2nd component on the 2nd grid 546 ENDIF 547 ! 548 IF( srcv(jpr_otx1)%clgrid == 'T' ) THEN 549 DO jj = 2, jpjm1 ! T ==> (U,V) 550 DO ji = fs_2, fs_jpim1 ! vector opt. 551 frcv(ji,jj,jpr_otx1) = 0.5 * ( frcv(ji+1,jj ,jpr_otx1) + frcv(ji,jj,jpr_otx1) ) 552 frcv(ji,jj,jpr_oty1) = 0.5 * ( frcv(ji ,jj+1,jpr_oty1) + frcv(ji,jj,jpr_oty1) ) 553 END DO 554 END DO 555 CALL lbc_lnk( frcv(:,:,jpr_otx1), 'U', -1. ) ; CALL lbc_lnk( frcv(:,:,jpr_oty1), 'V', -1. ) 556 ENDIF 557 ENDIF 558 ! ! ========================= ! 559 ELSE ! No dynamical coupling ! 560 ! ! ========================= ! 561 frcv(:,:,jpr_otx1) = 0.e0 ! here simply set to zero 562 frcv(:,:,jpr_oty1) = 0.e0 ! an external read in a file can be added instead 563 ! 250 564 ENDIF 251 !252 ! ! =========================== !253 ! ! Recieve freshwater fluxes !254 ! ! =========================== !255 !256 var_id = recv_id( 9) ; CALL cpl_prism_recv( var_id, date, zpew , info ) ! P-E over water257 var_id = recv_id(10) ; CALL cpl_prism_recv( var_id, date, zpei , info ) ! P-E over ice258 var_id = recv_id(11) ; CALL cpl_prism_recv( var_id, date, zpsol , info ) ! Snow fall over water and ice259 var_id = recv_id(12) ; CALL cpl_prism_recv( var_id, date, zevice, info ) ! Evaporation over ice (sublimination)260 !261 ! ... calculate water flux (P-E over open ocean and ice) and solid precipitation (positive upward)262 tprecip_recv(:,:) = ( zpew (:,:) + zpei (:,:) ) * tmask(:,:,1) * zfacwat263 sprecip_recv(:,:) = ( zpsol(:,:) + zevice(:,:) ) * tmask(:,:,1) * zfacwat264 565 265 ! ... Control print & check 266 IF(ln_ctl) THEN 267 WRITE(numout,*) ' flx:tprecip_recv - Minimum value is ', MINVAL( tprecip_recv ) 268 WRITE(numout,*) ' flx:tprecip_recv - Maximum value is ', MAXVAL( tprecip_recv ) 269 WRITE(numout,*) ' flx:tprecip_recv - Sum value is ', SUM ( tprecip_recv ) 566 ! u(v)tau will be modified by ice model -> need to be reset before each call of the ice/fsbc 567 IF( MOD( kt-1, k_fsbc ) == 0 ) THEN 568 utau(:,:) = frcv(:,:,jpr_otx1) 569 vtau(:,:) = frcv(:,:,jpr_oty1) 270 570 ENDIF 271 !!gm bug in mpp SUM require a mmp_sum call 272 !!gm further more this test is quite expensive ... only needed at the first time-step??? 273 IF( SUM( zpew*e1t*e2t ) /= SUM( zpew*e1t*e2t*tmask(:,:,1) ) ) THEN 274 WRITE(numout,*) ' flx: Forcing values outside Orca mask' 275 WRITE(numout,*) ' flx: Losses in water conservation' 276 WRITE(numout,*) ' flx: Masked ', SUM(zpew*e1t*e2t*tmask(:,:,1)) 277 WRITE(numout,*) ' flx: Unmasked ', SUM(zpew*e1t*e2t) 278 WRITE(numout,*) ' flx: Simulation STOP' 279 CALL FLUSH(numout) 280 STOP 281 END IF 282 ! 283 #if defined key_cpl_discharge 284 ! Runoffs 285 var_id = recv_id(18) ; CALL cpl_prism_recv ( var_id, date, calving_recv, info ) ! ice discharge into ocean 286 var_id = recv_id(19) ; CALL cpl_prism_recv ( var_id, date, zrunriv , info ) ! river discharge into ocean 287 var_id = recv_id(20) ; CALL cpl_prism_recv ( var_id, date, zruncot , info ) ! continental discharge into ocean 288 289 DO jj = 1, jpj 290 DO ji = 1, jpi 291 zfact = zfacwat * tmask(ji,jj,1) 292 calving_recv(ji,jj) = calving_recv(ji,jj) * zfact 293 rrunoff_recv(ji,jj) = ( zrunriv(ji,jj) + zruncot(ji,jj) ) * zfact 294 END DO 295 END DO 296 #else 297 calving_recv(:,:) = 0. 298 rrunoff_recv(:,:) = 0. 299 #endif 300 301 !!gm bug : this is not valid in mpp 302 !!gm and I presum this is not required at all as a lbc_lnk is applied to all the fields at the end of the routine 303 ! Oasis mask shift and update lateral boundary conditions (E. Maisonnave) 304 ! not tested when mpp is used, W. Park 305 !WSPTEST 306 qsr_oce_recv (jpi-1,:) = qsr_oce_recv (1,:) 307 qsr_ice_recv (jpi-1,:) = qsr_ice_recv (1,:) 308 qns_oce_recv (jpi-1,:) = qns_oce_recv (1,:) 309 qns_ice_recv (jpi-1,:) = qns_ice_recv (1,:) 310 dqns_ice_recv(jpi-1,:) = dqns_ice_recv(1,:) 311 tprecip_recv (jpi-1,:) = tprecip_recv (1,:) 312 sprecip_recv (jpi-1,:) = sprecip_recv (1,:) 313 fr1_i0_recv (jpi-1,:) = fr1_i0_recv (1,:) 314 fr2_i0_recv (jpi-1,:) = fr2_i0_recv (1,:) 315 rrunoff_recv (jpi-1,:) = rrunoff_recv (1,:) 316 calving_recv (jpi-1,:) = calving_recv (1,:) 317 !!gm end bug 318 319 qsr (:,:) = qsr_oce_recv (:,:) ! ocean surface boundary condition 320 qns (:,:) = qns_oce_recv (:,:) 321 emp (:,:) = zpew (:,:) 322 emps (:,:) = zpew (:,:) 323 324 qsr_ice (:,:) = qsr_ice_recv (:,:) ! ice forcing fields 325 qns_ice (:,:) = qns_ice_recv (:,:) 326 dqns_ice(:,:) = dqns_ice_recv(:,:) 327 tprecip (:,:) = tprecip_recv (:,:) 328 sprecip (:,:) = sprecip_recv (:,:) 329 fr1_i0 (:,:) = fr1_i0_recv (:,:) 330 fr2_i0 (:,:) = fr2_i0_recv (:,:) 331 332 !WSP rrunoff = rrunoff_recv 333 !WSP calving = calving_recv 334 rrunoff (:,:) = 0.e0 !WSP runoff and calving included in tprecip 335 calving (:,:) = 0.e0 !WSP 336 337 IF(ln_ctl) THEN 338 WRITE(numout,*) 'flx:qsr_oce - Minimum value is ', MINVAL( qsr_oce ) 339 WRITE(numout,*) 'flx:qsr_oce - Maximum value is ', MAXVAL( qsr_oce ) 340 WRITE(numout,*) 'flx:qsr_oce - Sum value is ', SUM ( qsr_oce ) 341 ! 342 WRITE(numout,*) 'flx:tprecip - Minimum value is ', MINVAL( tprecip ) 343 WRITE(numout,*) 'flx:tprecip - Maximum value is ', MAXVAL( tprecip ) 344 WRITE(numout,*) 'flx:tprecip - Sum value is ', SUM ( tprecip ) 571 ! ! ========================= ! 572 IF( k_ice <= 1 ) THEN ! heat & freshwater fluxes ! (Ocean only case) 573 ! ! ========================= ! 574 ! 575 ! ! non solar heat flux over the ocean (qns) 576 IF( srcv(jpr_qnsoce)%laction ) qns(:,:) = frcv(:,:,jpr_qnsoce) 577 IF( srcv(jpr_qnsmix)%laction ) qns(:,:) = frcv(:,:,jpr_qnsmix) 578 ! ! solar flux over the ocean (qsr) 579 IF( srcv(jpr_qsroce)%laction ) qsr(:,:) = frcv(:,:,jpr_qsroce) 580 IF( srcv(jpr_qsrmix)%laction ) qsr(:,:) = frcv(:,:,jpr_qsrmix) 581 ! 582 ! ! total freshwater fluxes over the ocean (emp, emps) 583 SELECT CASE( TRIM( cn_rcv_emp ) ) ! evaporation - precipitation 584 CASE( 'conservative' ) 585 emp(:,:) = frcv(:,:,jpr_tevp) - ( frcv(:,:,jpr_rain) + frcv(:,:,jpr_snow) ) 586 CASE( 'mixed oce-ice' ) 587 emp(:,:) = frcv(:,:,jpr_tevp) - ( frcv(:,:,jpr_rain) + frcv(:,:,jpr_semp) ) 588 CASE( 'ocean only', 'oce and ice' ) 589 emp(:,:) = frcv(:,:,jpr_oemp) 590 END SELECT 591 ! 592 ! ! runoffs and calving (added in emp) 593 IF( srcv(jpr_rnf)%laction ) emp(:,:) = emp(:,:) - frcv(:,:,jpr_rnf) 594 IF( srcv(jpr_cal)%laction ) emp(:,:) = emp(:,:) - ABS( frcv(:,:,jpr_cal) ) 595 ! 596 !!gm : this seems to be internal cooking, not sure to need that in a generic interface 597 !!gm at least should be optional... 598 !! IF( TRIM( cn_rcv_rnf ) == 'coupled' ) THEN ! add to the total freshwater budget 599 !! ! remove negative runoff 600 !! zcumulpos = SUM( MAX( frcv(:,:,jpr_rnf), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) ) 601 !! zcumulneg = SUM( MIN( frcv(:,:,jpr_rnf), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) ) 602 !! IF( lk_mpp ) CALL mpp_sum( zcumulpos ) ! sum over the global domain 603 !! IF( lk_mpp ) CALL mpp_sum( zcumulneg ) 604 !! IF( zcumulpos /= 0. ) THEN ! distribute negative runoff on positive runoff grid points 605 !! zcumulneg = 1.e0 + zcumulneg / zcumulpos 606 !! frcv(:,:,jpr_rnf) = MAX( frcv(:,:,jpr_rnf), 0.e0 ) * zcumulneg 607 !! ENDIF 608 !! ! add runoff to e-p 609 !! emp(:,:) = emp(:,:) - frcv(:,:,jpr_rnf) 610 !! ENDIF 611 !!gm end of internal cooking 612 ! 613 emps(:,:) = emp(:,:) ! concentration/dilution = emp 614 615 ! ! 10 m wind speed 616 IF( srcv(jpr_w10m)%laction ) wind10m(:,:) = frcv(:,:,jpr_w10m) 617 !!gm ---> blinder dans tke si cn_rcv_w10m == 'none' 618 ! 345 619 ENDIF 346 347 CALL lbc_lnk( qsr_oce , 'T', 1. ) 348 CALL lbc_lnk( qsr_ice , 'T', 1. ) 349 CALL lbc_lnk( qns_oce , 'T', 1. ) 350 CALL lbc_lnk( qns_ice , 'T', 1. ) 351 CALL lbc_lnk( tprecip , 'T', 1. ) 352 CALL lbc_lnk( sprecip , 'T', 1. ) 353 CALL lbc_lnk( rrunoff , 'T', 1. ) 354 CALL lbc_lnk( dqns_ice, 'T', 1. ) 355 CALL lbc_lnk( calving , 'T', 1. ) 356 CALL lbc_lnk( fr1_i0 , 'T', 1. ) 357 CALL lbc_lnk( fr2_i0 , 'T', 1. ) 358 359 IF(ln_ctl) THEN 360 WRITE(numout,*) 'flx(af lbc_lnk):qsr_oce - Minimum value is ', MINVAL( qsr_oce ) 361 WRITE(numout,*) 'flx(af lbc_lnk):qsr_oce - Maximum value is ', MAXVAL( qsr_oce ) 362 WRITE(numout,*) 'flx(af lbc_lnk):qsr_oce - Sum value is ', SUM ( qsr_oce ) 363 ! 364 WRITE(numout,*) 'flx(af lbc_lnk):tprecip - Minimum value is ', MINVAL( tprecip ) 365 WRITE(numout,*) 'flx(af lbc_lnk):tprecip - Maximum value is ', MAXVAL( tprecip ) 366 WRITE(numout,*) 'flx(af lbc_lnk):tprecip - Sum value is ', SUM ( tprecip ) 620 ! 621 END SUBROUTINE sbc_cpl_rcv 622 623 624 SUBROUTINE sbc_cpl_ice_tau( p_taui, p_tauj ) 625 !!---------------------------------------------------------------------- 626 !! *** ROUTINE sbc_cpl_ice_tau *** 627 !! 628 !! ** Purpose : provide the stress over sea-ice in coupled mode 629 !! 630 !! ** Method : transform the received stress from the atmosphere into 631 !! an atmosphere-ice stress in the (i,j) ocean referencial 632 !! and at the velocity point of the sea-ice model (cice_grid): 633 !! 'C'-grid : i- (j-) components given at U- (V-) point 634 !! 'B'-grid : both components given at I-point 635 !! 636 !! The received stress are : 637 !! - defined by 3 components (if cartesian coordinate) 638 !! or by 2 components (if spherical) 639 !! - oriented along geographical coordinate (if eastward-northward) 640 !! or along the local grid coordinate (if local grid) 641 !! - given at U- and V-point, resp. if received on 2 grids 642 !! or at a same point (T or I) if received on 1 grid 643 !! Therefore and if necessary, they are successively 644 !! processed in order to obtain them 645 !! first as 2 components on the sphere 646 !! second as 2 components oriented along the local grid 647 !! third as 2 components on the cice_grid point 648 !! 649 !! In 'oce and ice' case, only one vector stress field 650 !! is received. It has already been processed in sbc_cpl_rcv 651 !! so that it is now defined as (i,j) components given at U- 652 !! and V-points, respectively. Therefore, here only the third 653 !! transformation is done and only if the ice-grid is a 'B'-grid. 654 !! 655 !! ** Action : return ptau_i, ptau_j, the stress over the ice at cice_grid point 656 !!---------------------------------------------------------------------- 657 REAL(wp), INTENT(out), DIMENSION(jpi,jpj) :: p_taui ! i- & j-components of atmos-ice stress [N/m2] 658 REAL(wp), INTENT(out), DIMENSION(jpi,jpj) :: p_tauj ! at I-point (B-grid) or U & V-point (C-grid) 659 !! 660 INTEGER :: ji, jj ! dummy loop indices 661 INTEGER :: itx ! index of taux over ice 662 REAL(wp), DIMENSION(jpi,jpj) :: ztx, ztx ! 2D workspace 663 !!---------------------------------------------------------------------- 664 665 IF( srcv(jpr_itx1)%laction ) THEN ; itx = jpr_itx1 666 ELSE ; itx = jpr_otx1 367 667 ENDIF 368 ! 369 END SUBROUTINE sbc_cpl 370 668 669 ! do something only if we just received the stress from atmosphere 670 IF( nrcvinfo(itx) == PRISM_Recvd .OR. nrcvinfo(itx) == PRISM_FromRest .OR. & 671 & nrcvinfo(itx) == PRISM_RecvOut .OR. nrcvinfo(itx) == PRISM_FromRestOut ) THEN 672 673 ! ! ======================= ! 674 IF( srcv(jpr_itx1)%laction ) THEN ! ice stress received ! 675 ! ! ======================= ! 676 ! 677 IF( TRIM( cn_rcv_tau(2) ) == 'cartesian' ) THEN ! 2 components on the sphere 678 ! ! (cartesian to spherical -> 3 to 2 components) 679 CALL geo2oce( frcv(:,:,jpr_itx1), frcv(:,:,jpr_ity1), frcv(:,:,jpr_itz1), & 680 & srcv(jpr_itx1)%clgrid, ztx, zty ) 681 frcv(:,:,jpr_itx1) = ztx(:,:) ! overwrite 1st comp. on the 1st grid 682 frcv(:,:,jpr_itx1) = zty(:,:) ! overwrite 2nd comp. on the 1st grid 683 ! 684 IF( srcv(jpr_itx2)%laction ) THEN 685 CALL geo2oce( frcv(:,:,jpr_itx2), frcv(:,:,jpr_ity2), frcv(:,:,jpr_itz2), & 686 & srcv(jpr_itx2)%clgrid, ztx, zty ) 687 frcv(:,:,jpr_itx2) = ztx(:,:) ! overwrite 1st comp. on the 2nd grid 688 frcv(:,:,jpr_ity2) = zty(:,:) ! overwrite 2nd comp. on the 2nd grid 689 ENDIF 690 ! 691 ENDIF 692 ! 693 IF( TRIM( cn_rcv_tau(3) ) == 'eastward-northward' ) THEN ! 2 components oriented along the local grid 694 ! ! (geographical to local grid -> rotate the components) 695 CALL rot_rep( frcv(:,:,jpr_itx1), frcv(:,:,jpr_ity1), srcv(jpr_itx1)%clgrid, 'en->i', ztx ) 696 frcv(:,:,jpr_itx1) = ztx(:,:) ! overwrite 1st component on the 1st grid 697 IF( srcv(jpr_itx2)%laction ) THEN 698 CALL rot_rep( frcv(:,:,jpr_itx2), frcv(:,:,jpr_ity2), srcv(jpr_itx2)%clgrid, 'en->j', zty ) 699 ELSE 700 CALL rot_rep( frcv(:,:,jpr_itx1), frcv(:,:,jpr_ity1), srcv(jpr_itx1)%clgrid, 'en->j', zty ) 701 ENDIF 702 frcv(:,:,jpr_ity1) = zty(:,:) ! overwrite 2nd component on the 1st grid 703 ENDIF 704 ! ! ======================= ! 705 ELSE ! use ocean stress ! 706 ! ! ======================= ! 707 frcv(:,:,jpr_itx1) = frcv(:,:,jpr_otx1) 708 frcv(:,:,jpr_ity1) = frcv(:,:,jpr_oty1) 709 ! 710 ENDIF 711 712 ! ! ======================= ! 713 ! ! put on ice grid ! 714 ! ! ======================= ! 715 ! 716 ! j+1 j -----V---F 717 ! ice stress on ice velocity point (cice_grid) ! | 718 ! (C-grid ==>(U,V) or B-grid ==> I) j | T U 719 ! | | 720 ! j j-1 -I-------| 721 ! (for I) | | 722 ! i-1 i i 723 ! i i+1 (for I) 724 SELECT CASE ( cice_grid ) 725 ! 726 CASE( 'B' ) ! B-grid ==> I 727 SELECT CASE ( srcv(jpr_itx1)%clgrid ) 728 CASE( 'U' ) 729 DO jj = 2, jpjm1 ! (U,V) ==> I 730 DO ji = fs_2, fs_jpim1 ! vector opt. 731 p_taui(ji,jj) = 0.5 * ( frcv(ji-1,jj ,jpr_itx1) + frcv(ji-1,jj-1,jpr_itx1) ) 732 p_tauj(ji,jj) = 0.5 * ( frcv(ji ,jj-1,jpr_ity1) + frcv(ji-1,jj-1,jpr_ity1) ) 733 END DO 734 END DO 735 CASE( 'F' ) 736 DO jj = 2, jpjm1 ! F ==> I 737 DO ji = fs_2, fs_jpim1 ! vector opt. 738 p_taui(ji,jj) = frcv(ji-1,jj-1,jpr_itx1) 739 p_tauj(ji,jj) = frcv(ji-1,jj-1,jpr_ity1) 740 END DO 741 END DO 742 CASE( 'T' ) 743 DO jj = 2, jpjm1 ! T ==> I 744 DO ji = fs_2, fs_jpim1 ! vector opt. 745 p_taui(ji,jj) = 0.25 * ( frcv(ji,jj ,jpr_itx1) + frcv(ji-1,jj ,jpr_itx1) & 746 & + frcv(ji,jj-1,jpr_itx1) + frcv(ji-1,jj-1,jpr_itx1) ) 747 p_tauj(ji,jj) = 0.25 * ( frcv(ji,jj ,jpr_ity1) + frcv(ji-1,jj ,jpr_ity1) & 748 & + frcv(ji,jj-1,jpr_ity1) + frcv(ji-1,jj-1,jpr_ity1) ) 749 END DO 750 END DO 751 CASE( 'I' ) 752 p_taui(:,:) = frcv(:,:,jpr_itx1) ! I ==> I 753 p_tauj(:,:) = frcv(:,:,jpr_ity1) 754 END SELECT 755 IF( srcv(jpr_itx1)%clgrid /= 'I' ) THEN 756 CALL lbc_lnk( p_taui, 'I', -1. ) ; CALL lbc_lnk( p_tauj, 'I', -1. ) 757 ENDIF 758 ! 759 CASE( 'C' ) ! C-grid ==> U,V 760 SELECT CASE ( srcv(jpr_itx1)%clgrid ) 761 CASE( 'U' ) 762 p_taui(:,:) = frcv(:,:,jpr_itx1) ! (U,V) ==> (U,V) 763 p_tauj(:,:) = frcv(:,:,jpr_ity1) 764 CASE( 'F' ) 765 DO jj = 2, jpjm1 ! F ==> (U,V) 766 DO ji = fs_2, fs_jpim1 ! vector opt. 767 p_taui(ji,jj) = 0.5 * ( frcv(ji,jj,jpr_itx1) + frcv(ji ,jj-1,jpr_itx1) ) 768 p_tauj(ji,jj) = 0.5 * ( frcv(ji,jj,jpr_ity1) + frcv(ji-1,jj ,jpr_ity1) ) 769 END DO 770 END DO 771 CASE( 'T' ) 772 DO jj = 2, jpjm1 ! T ==> (U,V) 773 DO ji = fs_2, fs_jpim1 ! vector opt. 774 p_taui(ji,jj) = 0.5 * ( frcv(ji+1,jj ,jpr_itx1) + frcv(ji,jj,jpr_itx1) ) 775 p_tauj(ji,jj) = 0.5 * ( frcv(ji ,jj+1,jpr_ity1) + frcv(ji,jj,jpr_ity1) ) 776 END DO 777 END DO 778 CASE( 'I' ) 779 DO jj = 2, jpjm1 ! I ==> (U,V) 780 DO ji = fs_2, fs_jpim1 ! vector opt. 781 p_taui(ji,jj) = 0.5 * ( frcv(ji+1,jj+1,jpr_itx1) + frcv(ji+1,jj ,jpr_itx1) ) 782 p_tauj(ji,jj) = 0.5 * ( frcv(ji+1,jj+1,jpr_ity1) + frcv(ji ,jj+1,jpr_ity1) ) 783 END DO 784 END DO 785 END SELECT 786 IF( srcv(jpr_itx1)%clgrid /= 'U' ) THEN 787 CALL lbc_lnk( p_taui, 'U', -1. ) ; CALL lbc_lnk( p_tauj, 'V', -1. ) 788 ENDIF 789 END SELECT 790 791 !!gm Should be useless as sbc_cpl_ice_tau only called at coupled frequency 792 ! The receive stress are transformed such that in all case frcv(:,:,jpr_itx1), frcv(:,:,jpr_ity1) 793 ! become the i-component and j-component of the stress at the right grid point 794 !!gm frcv(:,:,jpr_itx1) = p_taui(:,:) 795 !!gm frcv(:,:,jpr_ity1) = p_tauj(:,:) 796 !!gm 797 ENDIF 798 ! 799 END SUBROUTINE sbc_cpl_ice_tau 800 801 802 SUBROUTINE sbc_cpl_ice_flx( p_frld, palbi , psst , pist, & 803 & pqns_tot, pqns_ice, & 804 & pqsr_tot, pqsr_ice, & 805 & pemp_tot, pemp_ice, psprecip ) 806 !!---------------------------------------------------------------------- 807 !! *** ROUTINE sbc_cpl_ice_flx_rcv *** 808 !! 809 !! ** Purpose : provide the heat and freshwater fluxes of the 810 !! ocean-ice system. 811 !! 812 !! ** Method : transform the fields received from the atmosphere into 813 !! surface heat and fresh water boundary condition for the 814 !! ice-ocean system. The following fields are provided: 815 !! * total non solar, solar and freshwater fluxes (qns_tot, 816 !! qsr_tot and emp_tot) (total means weighted ice-ocean flux) 817 !! NB: emp_tot include runoffs and calving. 818 !! * fluxes over ice (qns_ice, qsr_ice, emp_ice) where 819 !! emp_ice = sublimation - solid precipitation as liquid 820 !! precipitation are re-routed directly to the ocean and 821 !! runoffs and calving directly enter the ocean. 822 !! * solid precipitation (sprecip), used to add to qns_tot 823 !! the heat lost associated to melting solid precipitation 824 !! over the ocean fraction. 825 !! ===>> CAUTION here this changes the net heat flux received from 826 !! the atmosphere 827 !! * 10m wind module (wind10m) 828 !! 829 !! N.B. - fields over sea-ice are passed in argument so that 830 !! the module can be compile without sea-ice. 831 !! - the fluxes have been separated from the stress as 832 !! (a) they are updated at each ice time step compare to 833 !! an update at each coupled time step for the stress, and 834 !! (b) the conservative computation of the fluxes over the 835 !! sea-ice area requires the knowledge of the ice fraction 836 !! after the ice advection and before the ice thermodynamics, 837 !! so that the stress is updated before the ice dynamics 838 !! while the fluxes are updated after it. 839 !! 840 !! ** Action : update at each nf_ice time step: 841 !! pqns_tot, pqsr_tot non-solar and solar total heat fluxes 842 !! pqns_ice, pqsr_ice non-solar and solar heat fluxes over the ice 843 !! pemp_tot total evaporation - precipitation(liquid and solid) (-runoff)(-calving) 844 !! pemp_ice ice sublimation - solid precipitation over the ice 845 !! sprecip solid precipitation over the ocean 846 !! wind10m 10m wind module 847 !!---------------------------------------------------------------------- 848 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) :: p_frld ! lead fraction [0 to 1] 849 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) :: palbi ! ice albedo 850 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) :: psst ! sea surface temperature [Celcius] 851 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) :: pist ! ice surface temperature [Celcius] 852 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: pqns_tot ! total non solar heat flux [W/m2] 853 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: pqns_ice ! ice non solar heat flux [W/m2] 854 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: pqsr_tot ! total solar heat flux [W/m2] 855 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: pqsr_ice ! ice solar heat flux [W/m2] 856 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: pemp_tot ! total freshwater budget [Kg/m2/s] 857 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: pemp_ice ! ice solid freshwater budget [Kg/m2/s] 858 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: psprecip ! solid precipitation [Kg/m2/s] 859 !! 860 INTEGER :: ji, jj ! dummy loop indices 861 INTEGER :: isec, info ! temporary integer 862 REAL(wp):: zcoef, ztsurf ! temporary scalar 863 REAL(wp), DIMENSION(jpi,jpj) :: ztx, ztx ! 2D workspace 864 !!---------------------------------------------------------------------- 865 ! 866 ! ! ========================= ! 867 ! ! freshwater budget ! (emp) 868 ! ! ========================= ! 869 ! 870 ! ! total Precipitations - total Evaporation (emp_tot) 871 ! ! solid precipitation - sublimation (emp_ice) 872 ! ! solid Precipitation (sprecip) 873 SELECT CASE( TRIM( cn_rcv_emp ) ) 874 CASE( 'conservative' ) ! received fields: jpr_rain, jpr_snow, jpr_ievp, jpr_tevp 875 pemp_tot(:,:) = frcv(:,:,jpr_tevp) - frcv(:,:,jpr_rain) - frcv(:,:,jpr_snow) 876 pemp_ice(:,:) = frcv(:,:,jpr_ievp) - frcv(:,:,jpr_snow) 877 psprecip(:,:) = frcv(:,:,jpr_snow) 878 CASE( 'oce and ice' ) ! received fields: jpr_prsb, jpr_semp, jpr_oemp 879 pemp_tot(:,:) = p_frld(:,:) * frcv(:,:,jpr_oemp) + (1.- p_frld(:,:)) * frcv(:,:,jpr_semp) !!sm: rain over ice is missing?? 880 pemp_ice(:,:) = frcv(:,:,jpr_semp) 881 psprecip(:,:) = frcv(:,:,jpr_semp) !!gm here error due to sublimation 882 CASE( 'mixed oce-ice' ) ! received fields: jpr_rain, jpr_semp, jpr_tevp 883 pemp_tot(:,:) = frcv(:,:,jpr_tevp) - frcv(:,:,jpr_rain) + frcv(:,:,jpr_semp) !!gm here sublimation error ??? 884 pemp_ice(:,:) = frcv(:,:,jpr_semp) 885 psprecip(:,:) = frcv(:,:,jpr_semp) !!gm here error due to sublimation 886 END SELECT 887 ! 888 889 ! ! runoffs and calving (put in emp_tot) 890 IF( srcv(jpr_rnf)%laction ) pemp_tot(:,:) = pemp_tot(:,:) - frcv(:,:,jpr_rnf) 891 IF( srcv(jpr_cal)%laction ) pemp_tot(:,:) = pemp_tot(:,:) - ABS( frcv(:,:,jpr_cal) ) 892 ! 893 !!gm : this seems to be internal cooking, not sure to need that in a generic interface 894 !!gm at least should be optional... 895 !! ! remove negative runoff ! sum over the global domain 896 !! zcumulpos = SUM( MAX( frcv(:,:,jpr_rnf), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) ) 897 !! zcumulneg = SUM( MIN( frcv(:,:,jpr_rnf), 0.e0 ) * e1t(:,:) * e2t(:,:) * tmask_i(:,:) ) 898 !! IF( lk_mpp ) CALL mpp_sum( zcumulpos ) 899 !! IF( lk_mpp ) CALL mpp_sum( zcumulneg ) 900 !! IF( zcumulpos /= 0. ) THEN ! distribute negative runoff on positive runoff grid points 901 !! zcumulneg = 1.e0 + zcumulneg / zcumulpos 902 !! frcv(:,:,jpr_rnf) = MAX( frcv(:,:,jpr_rnf), 0.e0 ) * zcumulneg 903 !! ENDIF 904 !! pemp_tot(:,:) = pemp_tot(:,:) - frcv(:,:,jpr_rnf) ! add runoff to e-p 905 !! 906 !!gm end of internal cooking 907 908 909 ! ! ========================= ! 910 SELECT CASE( TRIM( cn_rcv_qns ) ) ! non solar heat fluxes ! (qns) 911 ! ! ========================= ! 912 CASE( 'conservative' ) ! the required fields are directly provided 913 pqns_tot(:,:) = frcv(:,:,jpr_qnsmix) 914 pqns_ice(:,:) = frcv(:,:,jpr_qnsice) 915 CASE( 'oce and ice' ) ! the total flux is computed from ocean and ice fluxes 916 pqns_tot(:,:) = p_frld(:,:) * frcv(:,:,jpr_qnsoce) + ( 1.- p_frld(:,:) ) * frcv(:,:,jpr_qnsice) 917 pqns_ice(:,:) = frcv(:,:,jpr_qnsice) 918 CASE( 'mixed oce-ice' ) ! the ice flux is cumputed from the total flux, the SST and ice informations 919 pqns_tot(:,:) = frcv(:,:,jpr_qnsmix) 920 pqns_ice(:,:) = frcv(:,:,jpr_qnsmix) & 921 & + frcv(:,:,jpr_dqnsdt) * ( pist(:,:) - psst(:,:) ) * ( 1. - p_frld(:,:) ) 922 END SELECT 923 ! ! snow melting heat flux .... 924 ! energy for melting solid precipitation over free ocean 925 zcoef = xlsn / rhosn 926 pqns_tot(:,:) = pqns_tot(:,:) - p_frld(:,:) * psprecip(:,:) * zcoef 927 !!gm 928 !! currently it is taken into account in leads budget but not in the qns_tot, and thus not in 929 !! the flux that enter the ocean.... 930 !! moreover 1 - it is not diagnose anywhere.... 931 !! 2 - it is unclear for me whether this heat lost is taken into account in the atmosphere or not... 932 !! 933 !! similar job should be done for snow and precipitation temperature 934 935 ! ! ========================= ! 936 SELECT CASE( TRIM( cn_rcv_qsr ) ) ! solar heat fluxes ! (qsr) 937 ! ! ========================= ! 938 CASE( 'conservative' ) 939 pqsr_tot(:,:) = frcv(:,:,jpr_qsrmix) 940 pqsr_ice(:,:) = frcv(:,:,jpr_qsrice) 941 CASE( 'oce and ice' ) 942 pqsr_tot(:,:) = p_frld(:,:) * frcv(:,:,jpr_qsroce) + ( 1.- p_frld(:,:) ) * frcv(:,:,jpr_qsrice) 943 pqsr_ice(:,:) = frcv(:,:,jpr_qsrice) 944 CASE( 'mixed oce-ice' ) 945 pqsr_tot(:,:) = frcv(:,:,jpr_qsrmix) 946 !!gm cpl_albedo ???? kezako ????? je pige pas grand chose ici.... 947 pqsr_ice(:,:) = qsr_mix(:,:) * ( 1.- palbi(:,:) ) & 948 & / ( 1.- ( cpl_ocean_albedo(ji,jj) * ( 1.- p_frld(ji,jj) ) & 949 & + palbi (ji,jj) * p_frld(ji,jj) ) ) 950 END SELECT 951 952 953 ! ! ========================= ! 954 ! ! 10 m wind speed ! (wind10m) 955 ! ! ========================= ! 956 ! 957 IF( srcv(jpr_w10m )%laction ) wind10m(:,:) = frcv(:,:,jpr_w10m) 958 !!gm ---> blinder dans tke si cn_rcv_w10m == 'none' 959 ! 960 END SUBROUTINE sbc_cpl_ice_flx_rcv 961 962 963 SUBROUTINE sbc_cpl_snd( kt ) 964 !!---------------------------------------------------------------------- 965 !! *** ROUTINE sbc_cpl_snd *** 966 !! 967 !! ** Purpose : provide the ocean-ice informations to the atmosphere 968 !! 969 !! ** Method : send to the atmosphere through a call to cpl_prism_snd 970 !! all the needed fields (as defined in sbc_cpl_init) 971 !!---------------------------------------------------------------------- 972 INTEGER, INTENT(in) :: kt 973 !! 974 INTEGER :: ji, jj ! dummy loop indices 975 INTEGER :: isec, info ! temporary integer 976 REAL(wp), DIMENSION(jpi,jpj) :: zfr_l ! 1. - fr_i(:,:) 977 REAL(wp), DIMENSION(jpi,jpj) :: ztmp1, ztmp2 978 REAL(wp), DIMENSION(jpi,jpj) :: zotx1 , zoty1 , zotz1, zitx1, zity1, zitz1 979 !!---------------------------------------------------------------------- 980 981 isec = ( kt - nit000 ) * NINT(rdttra(1)) ! date of exchanges 982 983 zfr_l(:,:) = 1.- fr_i(:,:) 984 985 ! ! ------------------------- ! 986 ! ! Surface temperature ! in Kelvin 987 ! ! ------------------------- ! 988 SELECT CASE( cn_snd_temperature) 989 CASE( 'oce only' ) ; ztmp1(:,:) = tn(:,:,1) + rt0 990 CASE( 'weighted oce and ice' ) ; ztmp1(:,:) = ( tn(:,:,1) + rt0 ) * zfr_l(:,:) 991 ztmp2(:,:) = tn_ice(:,:) * fr_i(:,:) 992 CASE( 'mixed oce-ice' ) ; ztmp1(:,:) = ( tn(:,:,1) + rt0 ) * zfr_l(:,:) + tn_ice(:,:) * fr_i(:,:) 993 END SELECT 994 IF( ssnd(jps_toce)%laction ) CALL cpl_prism_snd( jps_toce, isec, ztmp1, info ) 995 IF( ssnd(jps_tice)%laction ) CALL cpl_prism_snd( jps_tice, isec, ztmp2, info ) 996 IF( ssnd(jps_tmix)%laction ) CALL cpl_prism_snd( jps_tmix, isec, ztmp1, info ) 997 ! 998 ! ! ------------------------- ! 999 ! ! Albedo ! 1000 ! ! ------------------------- ! 1001 IF( ssnd(jps_albice)%laction ) THEN ! ice 1002 ztmp(:,:) = alb_ice(:,:) * fr_i(:,:) 1003 CALL cpl_prism_snd( jps_albice, isec, ztmp, info ) 1004 ENDIF 1005 IF( ssnd(jps_albmix)%laction ) THEN ! mixed ice-ocean 1006 ztmp(:,:) = albedo_oce_mix(:,:) * zfr_l(:,:) + alb_ice(:,:) * fr_i(:,:) 1007 CALL cpl_prism_snd( jps_albmix, isec, ztmp, info ) 1008 ENDIF 1009 ! ! ------------------------- ! 1010 ! ! Ice fraction & Thickness ! 1011 ! ! ------------------------- ! 1012 IF( ssnd(jps_fice)%laction ) CALL cpl_prism_snd( jps_fice, isec, fr_i , info ) 1013 IF( ssnd(jps_hice)%laction ) CALL cpl_prism_snd( jps_hice, isec, hice(:,:) * fr_i(:,:), info ) 1014 IF( ssnd(jps_hsnw)%laction ) CALL cpl_prism_snd( jps_hsnw, isec, hsnw(:,:) * fr_i(:,:), info ) 1015 ! 1016 ! ! ------------------------- ! 1017 IF( ssnd(jps_ocx1)%laction ) THEN ! Surface current ! 1018 ! ! ------------------------- ! 1019 SELECT CASE( TRIM( cn_snd_crt(1) ) ) 1020 CASE( 'oce only' ) 1021 DO jj = 2, jpjm1 1022 DO ji = fs_2, fs_jpim1 ! vector opt. 1023 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj ,1) ) 1024 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + un(ji ,jj-1,1) ) 1025 END DO 1026 END DO 1027 CASE( 'weighted oce and ice' ) 1028 IF( cice_grid = 'C' ) THEN ! 'C'-grid ice velocity 1029 DO jj = 2, jpjm1 1030 DO ji = fs_2, fs_jpim1 ! vector opt. 1031 zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(:,:) 1032 zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + un (ji ,jj-1,1) ) * zfr_l(:,:) 1033 zitx1(ji,jj) = 0.5 * ( utaui_ice(ji,jj) + utaui_ice(ji-1,jj ) ) * fr_i(:,:) 1034 zity1(ji,jj) = 0.5 * ( vtaui_ice(ji,jj) + vtaui_ice(ji ,jj-1) ) * fr_i(:,:) 1035 END DO 1036 END DO 1037 ELSE ! 'B'-grid ice velocity 1038 DO jj = 2, jpjm1 1039 DO ji = fs_2, fs_jpim1 ! vector opt. 1040 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj-1,1) ) * zfr_l(:,:) 1041 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + un(ji ,jj-1,1) ) * zfr_l(:,:) 1042 zitx1(ji,jj) = 0.25 * ( utaui_ice(ji+1,jj+1) + utaui_ice(ji,jj+1) & 1043 & + utaui_ice(ji+1,jj ) + utaui_ice(ji,jj ) ) * fr_i(:,:) 1044 zity1(ji,jj) = 0.25 * ( vtaui_ice(ji+1,jj+1) + vtaui_ice(ji,jj+1) & 1045 & + vtaui_ice(ji+1,jj ) + vtaui_ice(ji,jj ) ) * fr_i(:,:) 1046 END DO 1047 END DO 1048 ENDIF 1049 CALL lbc_lnk( zitx1, 'T', -1. ) ; CALL lbc_lnk( zity1, 'T', -1. ) 1050 CASE( 'mixed oce-ice' ) 1051 IF( cice_grid = 'C' ) THEN ! 'C'-grid ice velocity 1052 DO jj = 2, jpjm1 1053 DO ji = fs_2, fs_jpim1 ! vector opt. 1054 zotx1(ji,jj) = 0.5 * ( un (ji,jj,1) + un (ji-1,jj ,1) ) * zfr_l(:,:) 1055 & + 0.5 * ( utaui_ice(ji,jj) + utaui_ice(ji-1,jj ) ) * fr_i(:,:) 1056 zoty1(ji,jj) = 0.5 * ( vn (ji,jj,1) + un (ji ,jj-1,1) ) * zfr_l(:,:) 1057 & + 0.5 * ( vtaui_ice(ji,jj) + vtaui_ice(ji ,jj-1) ) * fr_i(:,:) 1058 END DO 1059 END DO 1060 ELSE ! 'B'-grid ice velocity 1061 DO jj = 2, jpjm1 1062 DO ji = fs_2, fs_jpim1 ! vector opt. 1063 zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj-1,1) ) * zfr_l(:,:) 1064 & + 0.25 * ( utaui_ice(ji+1,jj+1) + utaui_ice(ji,jj+1) & 1065 & + utaui_ice(ji+1,jj ) + utaui_ice(ji,jj ) ) * fr_i(:,:) 1066 zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + un(ji ,jj-1,1) ) * zfr_l(:,:) 1067 & + 0.25 * ( vtaui_ice(ji+1,jj+1) + vtaui_ice(ji,jj+1) & 1068 & + vtaui_ice(ji+1,jj ) + vtaui_ice(ji,jj ) ) * fr_i(:,:) 1069 END DO 1070 END DO 1071 ENDIF 1072 END SELECT 1073 CALL lbc_lnk( zotx1, 'T', -1. ) ; CALL lbc_lnk( zoty1, 'T', -1. ) 1074 ! 1075 ! 1076 IF( TRIM( cn_snd_crt(3) ) == 'eastward-northward' ) THEN ! Rotation of the components 1077 ! ! Ocean component 1078 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 1079 CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 1080 zotx1(:,:) = ztmp1(:,:) ! overwrite the components 1081 zoty1(:,:) = ztmp2(:,:) 1082 IF( ssnd(jps_ivx1)%laction ) THEN ! Ice component 1083 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 ) ! 1st component 1084 CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 ) ! 2nd component 1085 zitx1(:,:) = ztmp1(:,:) ! overwrite the components 1086 zity1(:,:) = ztmp2(:,:) 1087 ENDIF 1088 ENDIF 1089 ! 1090 !!gm Eric : Arnaud, je te laisse coder oce2geo ! 1091 ! spherical coordinates to cartesian -> 2 components to 3 components 1092 IF( TRIM( cn_snd_crt(2) ) == 'cartesian' ) THEN 1093 ztmp1(:,:) = zotx1(:,:) ! ocean currents 1094 ztmp2(:,:) = zoty1(:,:) 1095 CALL oce2geo ( ztmp1, ztmp2, 't', glamt, gphit, zotx1, zoty1, zotz1 ) 1096 ! 1097 IF( ssnd(jps_ivx1)%laction ) THEN ! ice velocities 1098 ztmp1(:,:) = zitx1(:,:) 1099 ztmp1(:,:) = zity1(:,:) 1100 CALL oce2geo ( ztmp1, ztmp2, 't', glamt, gphit, zitx1, zity1, zitz1 ) 1101 ENDIF 1102 ENDIF 1103 ! 1104 IF( ssnd(jps_ocx1)%laction ) CALL cpl_prism_snd( jps_ocx1, isec, zotx1, info ) ! ocean x current 1st grid 1105 IF( ssnd(jps_ocy1)%laction ) CALL cpl_prism_snd( jps_ocy1, isec, zoty1, info ) ! ocean y current 1st grid 1106 IF( ssnd(jps_ocz1)%laction ) CALL cpl_prism_snd( jps_ocz1, isec, zotz1, info ) ! ocean z current 1st grid 1107 ! 1108 IF( ssnd(jps_ivx1)%laction ) CALL cpl_prism_snd( jps_ivx1, isec, zitx1, info ) ! ice x current 1st grid 1109 IF( ssnd(jps_ivy1)%laction ) CALL cpl_prism_snd( jps_ivy1, isec, zity1, info ) ! ice y current 1st grid 1110 IF( ssnd(jps_ivz1)%laction ) CALL cpl_prism_snd( jps_ivz1, isec, zitz1, info ) ! ice z current 1st grid 1111 ! 1112 ENDIF 1113 ! 1114 END SUBROUTINE sbc_cpl_snd 1115 371 1116 #else 372 1117 !!---------------------------------------------------------------------- 373 !! Dummy routine NO sea surface restoring1118 !! Dummy module NO coupling 374 1119 !!---------------------------------------------------------------------- 375 LOGICAL, PUBLIC :: lk_sbc_cpl = .FALSE. !: coupled formulation flag1120 USE par_kind ! kind definition 376 1121 CONTAINS 377 SUBROUTINE sbc_cpl( kt ) ! Dummy routine 378 WRITE(*,*) 'sbc_cpl: you should not have seen that print! error?', kt 379 END SUBROUTINE sbc_cpl 1122 SUBROUTINE sbc_cpl_snd( kt ) 1123 WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?', kt 1124 END SUBROUTINE sbc_cpl_snd 1125 ! 1126 SUBROUTINE sbc_cpl_rcv( kt, k_fsbc, k_ice ) 1127 WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?', kt, k_fsbc, k_ice 1128 END SUBROUTINE sbc_cpl_rcv 1129 ! 1130 SUBROUTINE sbc_cpl_ice_tau( p_taui, p_tauj ) 1131 REAL(wp), INTENT(out), DIMENSION(:,:) :: p_taui ! i- & j-components of atmos-ice stress [N/m2] 1132 REAL(wp), INTENT(out), DIMENSION(:,:) :: p_tauj ! at I-point (B-grid) or U & V-point (C-grid) 1133 p_taui(:,:) = 0. ; p_tauj(:,:) = 0. ! stupid definition to avoid warning message when compiling... 1134 WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?' 1135 END SUBROUTINE sbc_cpl_ice_tau 1136 ! 1137 SUBROUTINE sbc_cpl_ice_flx( p_frld, palbi , psst , pist, & 1138 & pqns_tot, pqns_ice, & 1139 & pqsr_tot, pqsr_ice, & 1140 & pemp_tot, pemp_ice, psprecip ) 1141 REAL(wp), INTENT(in ), DIMENSION(:,:) :: p_frld ! lead fraction [0 to 1] 1142 REAL(wp), INTENT(in ), DIMENSION(:,:) :: palbi ! ice albedo 1143 REAL(wp), INTENT(in ), DIMENSION(:,:) :: psst ! sea surface temperature [Celcius] 1144 REAL(wp), INTENT(in ), DIMENSION(:,:) :: pist ! ice surface temperature [Celcius] 1145 REAL(wp), INTENT( out), DIMENSION(:,:) :: pqns_tot ! total non solar heat flux [W/m2] 1146 REAL(wp), INTENT( out), DIMENSION(:,:) :: pqns_ice ! ice non solar heat flux [W/m2] 1147 REAL(wp), INTENT( out), DIMENSION(:,:) :: pqsr_tot ! total solar heat flux [W/m2] 1148 REAL(wp), INTENT( out), DIMENSION(:,:) :: pqsr_ice ! ice solar heat flux [W/m2] 1149 REAL(wp), INTENT( out), DIMENSION(:,:) :: pemp_tot ! total freshwater budget [Kg/m2/s] 1150 REAL(wp), INTENT( out), DIMENSION(:,:) :: pemp_ice ! ice solid freshwater budget [Kg/m2/s] 1151 REAL(wp), INTENT( out), DIMENSION(:,:) :: psprecip ! solid precipitation [Kg/m2/s] 1152 WRITE(*,*) 'sbc_cpl_snd: You should not have seen this print! error?', p_frld(1,1), palbi(1,1), psst(1,1), pist(1,1) 1153 ! stupid definition to avoid warning message when compiling... 1154 pqns_tot(:,:) = 0. ; pqns_ice(:,:) = 0. 1155 pqsr_tot(:,:) = 0. ; pqsr_ice(:,:) = 0. 1156 pemp_tot(:,:) = 0. ; pemp_ice(:,:) = 0. ; psprecip(:,:) = 0. 1157 END SUBROUTINE sbc_cpl_ice_flx 1158 380 1159 #endif 381 1160
Note: See TracChangeset
for help on using the changeset viewer.