Changeset 2370
- Timestamp:
- 2010-11-10T08:48:54+01:00 (14 years ago)
- Location:
- branches/nemo_v3_3_beta/NEMOGCM
- Files:
-
- 29 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/nemo_v3_3_beta/NEMOGCM/ARCH/arch-ifort_osx.fcm
r2281 r2370 15 15 16 16 17 %NCDF_INC -I/usr/local/ pub/netcdf/4.0.1-ifort/include18 %NCDF_LIB -L /usr/local/ pub/netcdf/4.0.1-ifort/lib -lnetcdf17 %NCDF_INC -I/usr/local/include 18 %NCDF_LIB -L /usr/local/lib -lnetcdf 19 19 %FC mpif90 20 20 %FCFLAGS -r8 -O3 -traceback -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/GYRE/EXP00/namelist
r2364 r2370 131 131 ! =1 use observed ice-cover , 132 132 ! =2 ice-model used ("key_lim3" or "key_lim2) 133 nn_ico_cpl = 0 ! ice-ocean coupling : =0 each nn_fsbc134 ! =1 stresses recomputed each ocean time step ("key_lim3" only)135 ! =2 combination of 0 and 1 cases ("key_lim3" only)136 133 ln_dm2dc = .false. ! daily mean to diurnal cycle short wave (qsr) 137 134 ln_rnf = .false. ! runoffs (T => fill namsbc_rnf) -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/GYRE_LOBSTER/EXP00/namelist
r2364 r2370 131 131 ! =1 use observed ice-cover , 132 132 ! =2 ice-model used ("key_lim3" or "key_lim2) 133 nn_ico_cpl = 0 ! ice-ocean coupling : =0 each nn_fsbc134 ! =1 stresses recomputed each ocean time step ("key_lim3" only)135 ! =2 combination of 0 and 1 cases ("key_lim3" only)136 133 ln_dm2dc = .false. ! daily mean to diurnal cycle short wave (qsr) 137 134 ln_rnf = .false. ! runoffs (T => fill namsbc_rnf) -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/1_namelist
r2367 r2370 131 131 ! =1 use observed ice-cover , 132 132 ! =2 ice-model used ("key_lim3" or "key_lim2) 133 nn_ico_cpl = 0 ! ice-ocean coupling : =0 each nn_fsbc134 ! =1 stresses recomputed each ocean time step ("key_lim3" only)135 ! =2 combination of 0 and 1 cases ("key_lim3" only)136 133 ln_dm2dc = .false. ! daily mean to diurnal cycle short wave (qsr) 137 134 ln_rnf = .false. ! runoffs (T => fill namsbc_rnf) -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist
r2367 r2370 131 131 ! =1 use observed ice-cover , 132 132 ! =2 ice-model used ("key_lim3" or "key_lim2) 133 nn_ico_cpl = 0 ! ice-ocean coupling : =0 each nn_fsbc134 ! =1 stresses recomputed each ocean time step ("key_lim3" only)135 ! =2 combination of 0 and 1 cases ("key_lim3" only)136 133 ln_dm2dc = .false. ! daily mean to diurnal cycle short wave (qsr) 137 134 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/cpp_ORCA2_LIM.fcm
r2369 r2370 1 bld::tool::fppkeys key_trabbl key_vectopt_loop key_orca_r2 key_lim2 key_dynspg_flt key_diaeiv key_ldfslp key_traldf_c2d key_traldf_eiv key_dynldf_c3d key_dtatem key_dtasal key_tradmp key_zdftke key_zdfddm key_iomput key_nproci=1 key_nprocj=11 bld::tool::fppkeys key_trabbl key_vectopt_loop key_orca_r2 key_lim2 key_dynspg_flt key_diaeiv key_ldfslp key_traldf_c2d key_traldf_eiv key_dynldf_c3d key_dtatem key_dtasal key_tradmp key_zdftke key_zdfddm key_iomput key_nproci=1 key_nprocj=1 -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM_PISCES/EXP00/namelist
r2367 r2370 131 131 ! =1 use observed ice-cover , 132 132 ! =2 ice-model used ("key_lim3" or "key_lim2) 133 nn_ico_cpl = 0 ! ice-ocean coupling : =0 each nn_fsbc134 ! =1 stresses recomputed each ocean time step ("key_lim3" only)135 ! =2 combination of 0 and 1 cases ("key_lim3" only)136 133 ln_dm2dc = .false. ! daily mean to diurnal cycle short wave (qsr) 137 134 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist
r2364 r2370 131 131 ! =1 use observed ice-cover , 132 132 ! =2 ice-model used ("key_lim3" or "key_lim2) 133 nn_ico_cpl = 0 ! ice-ocean coupling : =0 each nn_fsbc134 ! =1 stresses recomputed each ocean time step ("key_lim3" only)135 ! =2 combination of 0 and 1 cases ("key_lim3" only)136 133 ln_dm2dc = .false. ! daily mean to diurnal cycle short wave (qsr) 137 134 ln_rnf = .true. ! runoffs (T => fill namsbc_rnf) -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/ice_2.F90
r2319 r2370 4 4 !! Sea Ice physics: diagnostics variables of ice defined in memory 5 5 !!===================================================================== 6 !! History : 2.0 ! 03-08 (C. Ethe) F90: Free form and module6 !! History : 2.0 ! 2003-08 (C. Ethe) F90: Free form and module 7 7 !! 3.3 ! 2009-05 (G.Garric) addition of the lim2_evp cas 8 8 !!---------------------------------------------------------------------- … … 11 11 !! 'key_lim2' : LIM 2.0 sea-ice model 12 12 !!---------------------------------------------------------------------- 13 !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010)14 !! $Id$15 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)16 !!----------------------------------------------------------------------17 !! * Modules used18 13 USE par_ice_2 ! LIM sea-ice parameters 19 14 … … 61 56 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: ust2s !: friction velocity 62 57 63 !!* diagnostic quantities58 !!* Ice Rheology 64 59 # if defined key_lim2_vp 65 60 ! !!* VP rheology * 61 LOGICAL , PUBLIC :: lk_lim2_vp = .TRUE. !: Visco-Plactic reology flag 62 ! 66 63 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: hsnm , hicm !: mean snow and ice thicknesses 67 CHARACTER(len=1), PUBLIC :: cl_grid = 'B' !: type of grid used in ice dynamics, 'C' or 'B'68 64 ! 69 65 # else 70 66 ! !!* EVP rheology * 71 CHARACTER(len=1), PUBLIC :: cl_grid = 'C' !: type of grid used in ice dynamics, 'C' or 'B' 67 LOGICAL , PUBLIC:: lk_lim2_vp = .FALSE. !: Visco-Plactic reology flag 68 ! 72 69 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: stress1_i !: first stress tensor element 73 70 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: stress2_i !: second stress tensor element … … 76 73 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: divu_i !: Divergence of the velocity field [s-1] -> limrhg.F90 77 74 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: shear_i !: Shear of the velocity field [s-1] -> limrhg.F90 75 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: at_i !: ice fraction 78 76 ! 79 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: at_i !:80 77 REAL(wp), PUBLIC, DIMENSION(:,:) , POINTER :: vt_s ,vt_i !: mean snow and ice thicknesses 81 78 REAL(wp), PUBLIC, DIMENSION(jpi,jpj), TARGET :: hsnm , hicm !: target vt_s,vt_i pointers 82 !83 79 #endif 84 80 … … 131 127 #endif 132 128 129 !!---------------------------------------------------------------------- 130 !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 131 !! $Id$ 132 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 133 133 !!====================================================================== 134 134 END MODULE ice_2 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limdyn_2.F90
r2332 r2370 4 4 !! Sea-Ice dynamics : 5 5 !!====================================================================== 6 !! History : 1.0 ! 01-04 (LIM) Original code 7 !! 2.0 ! 02-08 (C. Ethe, G. Madec) F90, mpp 8 !! 2.0 ! 03-08 (C. Ethe) add lim_dyn_init 9 !! 2.0 ! 06-07 (G. Madec) Surface module 6 !! History : 1.0 ! 2001-04 (LIM) Original code 7 !! 2.0 ! 2002-08 (C. Ethe, G. Madec) F90, mpp 8 !! 2.0 ! 2003-08 (C. Ethe) add lim_dyn_init 9 !! 2.0 ! 2006-07 (G. Madec) Surface module 10 !! 3.3 ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 10 11 !!--------------------------------------------------------------------- 11 12 #if defined key_lim2 … … 16 17 !! lim_dyn_init_2 : initialization and namelist read 17 18 !!---------------------------------------------------------------------- 18 USE dom_oce ! ocean space and time domain 19 USE sbc_oce ! 20 USE phycst ! 21 USE ice_2 ! 22 USE dom_ice_2 ! 23 USE limistate_2 ! 24 #if defined key_lim2_vp 25 USE limrhg_2 ! ice rheology 26 #else 27 USE limrhg ! LIM : EVP ice rheology 28 #endif 29 USE lbclnk ! 30 USE lib_mpp ! 31 USE in_out_manager ! I/O manager 32 USE prtctl ! Print control 19 USE dom_oce ! ocean space and time domain 20 USE sbc_oce ! ocean surface boundary condition 21 USE phycst ! physical constant 22 USE ice_2 ! LIM-2: ice variables 23 USE sbc_ice ! Surface boundary condition: sea-ice fields 24 USE dom_ice_2 ! LIM-2: ice domain 25 USE limistate_2 ! LIM-2: initial state 26 USE limrhg_2 ! LIM-2: VP ice rheology 27 USE limrhg ! LIM : EVP ice rheology 28 USE lbclnk ! lateral boundary condition - MPP link 29 USE lib_mpp ! MPP library 30 USE in_out_manager ! I/O manager 31 USE prtctl ! Print control 33 32 34 33 IMPLICIT NONE 35 34 PRIVATE 36 35 37 PUBLIC lim_dyn_2 ! routine called by sbc_ice_lim 38 39 !! * Module variables 40 REAL(wp) :: rone = 1.e0 ! constant value 41 36 PUBLIC lim_dyn_2 ! routine called by sbc_ice_lim 37 38 !! * Substitutions 42 39 # include "vectopt_loop_substitute.h90" 43 40 !!---------------------------------------------------------------------- 44 41 !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 45 42 !! $Id$ 46 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 47 !!---------------------------------------------------------------------- 48 43 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 44 !!---------------------------------------------------------------------- 49 45 CONTAINS 50 46 … … 90 86 i_jpj = jpj 91 87 IF(ln_ctl) CALL prt_ctl_info( 'lim_dyn : i_j1 = ', ivar1=i_j1, clinfo2=' ij_jpj = ', ivar2=i_jpj ) 92 #if defined key_lim2_vp 93 CALL lim_rhg_2( i_j1, i_jpj ) ! VP rheology 94 #else 95 CALL lim_rhg ( i_j1, i_jpj ) ! EVP rheology 96 #endif 88 IF( lk_lim2_vp ) THEN ; CALL lim_rhg_2( i_j1, i_jpj ) ! VP rheology 89 ELSE ; CALL lim_rhg ( i_j1, i_jpj ) ! EVP rheology 90 ENDIF 97 91 ! 98 92 ELSE ! optimization of the computational area … … 112 106 i_j1 = i_j1 + 1 113 107 END DO 114 #if defined key_lim2_vp 115 i_j1 = MAX( 1, i_j1-1 ) 116 IF(ln_ctl) WRITE(numout,*) 'lim_dyn : NH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 117 ! 118 CALL lim_rhg_2( i_j1, i_jpj ) 119 ! 120 #else 121 i_j1 = MAX( 1, i_j1-2 ) 108 IF( lk_lim2_vp ) THEN ! VP rheology 109 i_j1 = MAX( 1, i_j1-1 ) 110 CALL lim_rhg_2( i_j1, i_jpj ) 111 ELSE ! EVP rheology 112 i_j1 = MAX( 1, i_j1-2 ) 113 CALL lim_rhg( i_j1, i_jpj ) 114 ENDIF 122 115 IF(ln_ctl) WRITE(numout,*) 'lim_dyn : NH i_j1 = ', i_j1, 'ij_jpj = ', i_jpj 123 CALL lim_rhg( i_j1, i_jpj ) 124 ! 125 #endif 116 ! 126 117 ! Southern hemisphere 127 118 i_j1 = 1 … … 130 121 i_jpj = i_jpj - 1 131 122 END DO 132 #if defined key_lim2_vp 133 i_jpj = MIN( jpj, i_jpj+2 ) 134 IF(ln_ctl) WRITE(numout,*) 'lim_dyn : SH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 135 ! 136 CALL lim_rhg_2( i_j1, i_jpj ) 137 ! 138 #else 139 i_jpj = MIN( jpj, i_jpj+1 ) 123 IF( lk_lim2_vp ) THEN ! VP rheology 124 i_jpj = MIN( jpj, i_jpj+2 ) 125 CALL lim_rhg_2( i_j1, i_jpj ) 126 ELSE ! EVP rheology 127 i_jpj = MIN( jpj, i_jpj+1 ) 128 CALL lim_rhg( i_j1, i_jpj ) 129 ENDIF 140 130 IF(ln_ctl) WRITE(numout,*) 'lim_dyn : SH i_j1 = ', i_j1, 'ij_jpj = ', i_jpj 141 CALL lim_rhg( i_j1, i_jpj ) 142 ! 143 #endif 144 131 ! 145 132 ELSE ! local domain extends over one hemisphere only 146 133 ! ! Rheology is computed only over the ice cover … … 156 143 i_jpj = i_jpj - 1 157 144 END DO 158 i_jpj = MIN( jpj, i_jpj+2) 159 145 i_jpj = MIN( jpj, i_jpj+2 ) 146 ! 147 IF( lk_lim2_vp ) THEN ! VP rheology 148 i_jpj = MIN( jpj, i_jpj+2 ) 149 CALL lim_rhg_2( i_j1, i_jpj ) ! VP rheology 150 ELSE ! EVP rheology 151 i_j1 = MAX( 1 , i_j1-2 ) 152 i_jpj = MIN( jpj, i_jpj+1 ) 153 CALL lim_rhg ( i_j1, i_jpj ) ! EVP rheology 154 ENDIF 160 155 IF(ln_ctl) WRITE(numout,*) 'lim_dyn : one hemisphere: i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 161 !162 #if defined key_lim2_vp163 i_jpj = MIN( jpj, i_jpj+2)164 CALL lim_rhg_2( i_j1, i_jpj ) ! VP rheology165 #else166 i_j1 = MAX( 1, i_j1-2 )167 i_jpj = MIN( jpj, i_jpj+1)168 CALL lim_rhg ( i_j1, i_jpj ) ! EVP rheology169 #endif170 156 ! 171 157 ENDIF … … 177 163 ! computation of friction velocity 178 164 ! -------------------------------- 179 SELECT CASE( cl_grid ) 180 CASE( 'C' ) ! C-grid ice dynamics (EVP) 181 ! ice-ocean & ice velocity at ocean velocity points 182 zu_io(:,:) = u_ice(:,:) - ssu_m(:,:) 165 SELECT CASE( cp_ice_msh ) ! ice-ocean relative velocity at u- & v-pts 166 CASE( 'C' ) ! EVP : C-grid ice dynamics 167 zu_io(:,:) = u_ice(:,:) - ssu_m(:,:) ! ice-ocean & ice velocity at ocean velocity points 183 168 zv_io(:,:) = v_ice(:,:) - ssv_m(:,:) 184 ! 185 CASE( 'B' ) ! B-grid ice dynamics (VP) 186 ! ice-ocean velocity at U & V-points (u_ice v_ice at I-point ; ssu_m, ssv_m at U- & V-points) 187 DO jj = 1, jpjm1 188 DO ji = 1, jpim1 ! NO vector opt. 189 zu_io(ji,jj) = 0.5 * ( u_ice(ji+1,jj+1) + u_ice(ji+1,jj ) ) - ssu_m(ji,jj) 190 zv_io(ji,jj) = 0.5 * ( v_ice(ji+1,jj+1) + v_ice(ji ,jj+1) ) - ssv_m(ji,jj) 169 CASE( 'I' ) ! VP : B-grid ice dynamics (I-point) 170 DO jj = 1, jpjm1 ! u_ice v_ice at I-point ; ssu_m, ssv_m at U- & V-points 171 DO ji = 1, jpim1 ! NO vector opt. ! 172 zu_io(ji,jj) = 0.5_wp * ( u_ice(ji+1,jj+1) + u_ice(ji+1,jj ) ) - ssu_m(ji,jj) 173 zv_io(ji,jj) = 0.5_wp * ( v_ice(ji+1,jj+1) + v_ice(ji ,jj+1) ) - ssv_m(ji,jj) 191 174 END DO 192 175 END DO … … 194 177 195 178 ! frictional velocity at T-point 179 zcoef = 0.5_wp * cw 196 180 DO jj = 2, jpjm1 197 181 DO ji = 2, jpim1 ! NO vector opt. because of zu_io 198 ust2s(ji,jj) = 0.5 * cw & 199 & * ( zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj) & 200 & + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1) ) * tms(ji,jj) 182 ust2s(ji,jj) = zcoef * ( zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj) & 183 & + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1) ) * tms(ji,jj) 201 184 END DO 202 185 END DO … … 207 190 DO jj = 2, jpjm1 208 191 DO ji = fs_2, fs_jpim1 ! vector opt. 209 ust2s(ji,jj) = zcoef * tms(ji,jj) *SQRT( utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj) &210 & + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1))192 ust2s(ji,jj) = zcoef * SQRT( utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj) & 193 & + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) * tms(ji,jj) 211 194 END DO 212 195 END DO … … 217 200 ! 218 201 IF(ln_ctl) CALL prt_ctl(tab2d_1=ust2s , clinfo1=' lim_dyn : ust2s :') 219 202 ! 220 203 END SUBROUTINE lim_dyn_2 221 204 … … 265 248 WRITE(numout,*) ' coefficient for the solution of int. stresses alphaevp = ', alphaevp 266 249 ENDIF 250 ! 251 IF( angvg /= 0._wp .AND. .NOT.lk_lim2_vp ) THEN 252 CALL ctl_warn( 'lim_dyn_init_2: turning angle for oceanic stress not properly coded for EVP ', & 253 & '(see limsbc_2 module). We force angvg = 0._wp' ) 254 angvg = 0._wp 255 ENDIF 267 256 268 257 ! Initialization -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limmsh_2.F90
r2319 r2370 4 4 !! LIM 2.0 ice model : definition of the ice mesh parameters 5 5 !!====================================================================== 6 !! History : - ! 2001-04 (LIM) original code 7 !! 1.0 ! 2002-08 (C. Ethe, G. Madec) F90, module 8 !! 3.3 ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 9 !!---------------------------------------------------------------------- 6 10 #if defined key_lim2 7 11 !!---------------------------------------------------------------------- … … 10 14 !! lim_msh_2 : definition of the ice mesh 11 15 !!---------------------------------------------------------------------- 12 !! * Modules used13 16 USE phycst 14 17 USE dom_oce … … 20 23 PRIVATE 21 24 22 !! * Accessibility23 25 PUBLIC lim_msh_2 ! routine called by ice_ini_2.F90 24 26 … … 26 28 !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 27 29 !! $Id$ 28 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 29 !!---------------------------------------------------------------------- 30 30 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 31 !!---------------------------------------------------------------------- 31 32 CONTAINS 32 33 … … 43 44 !! 44 45 !! ** Refer. : Deleersnijder et al. Ocean Modelling 100, 7-10 45 !!46 !! ** History :47 !! original : 01-04 (LIM)48 !! addition : 02-08 (C. Ethe, G. Madec)49 46 !!--------------------------------------------------------------------- 50 !! * Local variables51 47 INTEGER :: ji, jj ! dummy loop indices 52 48 REAL(wp) :: zusden ! local scalars 53 54 49 #if defined key_lim2_vp 55 50 REAL(wp) :: zusden2 ! local scalars … … 267 262 END DO 268 263 END DO 269 270 !--lateral boundary conditions 271 CALL lbc_lnk( tmu(:,:), 'I', 1. ) 264 CALL lbc_lnk( tmu(:,:), 'I', 1. ) !--lateral boundary conditions 272 265 #else 273 266 ! EVP rheology : ice velocity point are U- & V-points ; ice vorticity … … 278 271 WHERE( fmask(:,:,1) == 1.e0 ) tmf(:,:) = 1.e0 279 272 #endif 280 273 ! 281 274 ! unmasked and masked area of T-grid cell 282 275 area(:,:) = e1t(:,:) * e2t(:,:) 283 276 ! 284 277 END SUBROUTINE lim_msh_2 285 278 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limrhg_2.F90
r2319 r2370 7 7 !! 1.0 ! 1994-12 (H. Goosse) 8 8 !! 2.0 ! 2003-12 (C. Ethe, G. Madec) F90, mpp 9 !! " "! 2006-08 (G. Madec) surface module, ice-stress at I-point10 !! " "! 2009-09 (G. Madec) Huge verctor optimisation9 !! - ! 2006-08 (G. Madec) surface module, ice-stress at I-point 10 !! - ! 2009-09 (G. Madec) Huge verctor optimisation 11 11 !! 3.3 ! 2009-05 (G.Garric, C. Bricaud) addition of the lim2_evp case 12 12 !!---------------------------------------------------------------------- 13 13 #if defined key_lim2 && defined key_lim2_vp 14 14 !!---------------------------------------------------------------------- 15 !! 'key_lim2' and LIM 2.0sea-ice model15 !! 'key_lim2' AND LIM-2 sea-ice model 16 16 !! 'key_lim2_vp' VP ice rheology 17 !!----------------------------------------------------------------------18 17 !!---------------------------------------------------------------------- 19 18 !! lim_rhg_2 : computes ice velocities … … 36 35 PUBLIC lim_rhg_2 ! routine called by lim_dyn 37 36 38 REAL(wp) :: rzero = 0. e0! constant value: zero39 REAL(wp) :: rone = 1. e0! and one37 REAL(wp) :: rzero = 0._wp ! constant value: zero 38 REAL(wp) :: rone = 1._wp ! and one 40 39 41 40 !! * Substitutions … … 44 43 !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 45 44 !! $Id$ 46 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 47 !!---------------------------------------------------------------------- 48 45 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 46 !!---------------------------------------------------------------------- 49 47 CONTAINS 50 48 … … 67 65 INTEGER :: iter, jter ! temporary integers 68 66 CHARACTER (len=50) :: charout 69 REAL(wp) :: ze11 , ze12 , ze22 , ze21 ! temporaryscalars70 REAL(wp) :: zt11 , zt12 , zt21 , zt22 ! " "71 REAL(wp) :: zvis11, zvis21, zvis12, zvis22 ! " "72 REAL(wp) :: zgphsx, ztagnx, zgsshx, zunw, zur, zusw ! " "73 REAL(wp) :: zgphsy, ztagny, zgsshy, zvnw, zvr ! " "67 REAL(wp) :: ze11 , ze12 , ze22 , ze21 ! local scalars 68 REAL(wp) :: zt11 , zt12 , zt21 , zt22 ! - - 69 REAL(wp) :: zvis11, zvis21, zvis12, zvis22 ! - - 70 REAL(wp) :: zgphsx, ztagnx, zgsshx, zunw, zur, zusw ! - - 71 REAL(wp) :: zgphsy, ztagny, zgsshy, zvnw, zvr ! - - 74 72 REAL(wp) :: zresm, za, zac, zmod 75 73 REAL(wp) :: zmpzas, zstms, zindu, zusdtp, zmassdt, zcorlal … … 91 89 REAL(wp), DIMENSION(jpi,0:jpj+1) :: zzfrld, zztms 92 90 REAL(wp), DIMENSION(jpi,0:jpj+1) :: zi1, zi2, zmasst, zpresh 93 94 91 !!------------------------------------------------------------------- 95 96 !!bug97 !! u_oce(:,:) = 0.e098 !! v_oce(:,:) = 0.e099 !! write(*,*) 'rhg min, max u & v', maxval(u_oce), minval(u_oce), maxval(v_oce), minval(v_oce)100 !!bug101 92 102 93 ! Store initial velocities 103 94 ! ---------------- 104 zztms(:,0 ) = 0.e0 ; zzfrld(:,0 ) = 0.e0 105 zztms(:,jpj+1) = 0.e0 ; zzfrld(:,jpj+1) = 0.e0 106 zu0(:,0 ) = 0.e0 ; zv0(:,0 ) = 0.e0 107 zu0(:,jpj+1) = 0.e0 ; zv0(:,jpj+1) = 0.e0 108 zztms(:,1:jpj) = tms(:,:) ; zzfrld(:,1:jpj) = frld(:,:) 109 zu0(:,1:jpj) = u_ice(:,:) ; zv0(:,1:jpj) = v_ice(:,:) 110 111 zu_a(:,:) = zu0(:,:) ; zv_a(:,:) = zv0(:,:) 112 zu_n(:,:) = zu0(:,:) ; zv_n(:,:) = zv0(:,:) 95 zztms(:,0 ) = 0._wp ; zzfrld(:,0 ) = 0._wp 96 zztms(:,jpj+1) = 0._wp ; zzfrld(:,jpj+1) = 0._wp 97 zu0 (:,0 ) = 0._wp ; zv0 (:,0 ) = 0._wp 98 zu0 (:,jpj+1) = 0._wp ; zv0 (:,jpj+1) = 0._wp 99 zztms(:,1:jpj) = tms (:,:) ; zzfrld(:,1:jpj) = frld (:,:) 100 zu0 (:,1:jpj) = u_ice(:,:) ; zv0 (:,1:jpj) = v_ice(:,:) 101 zu_a (:, : ) = zu0 (:,:) ; zv_a (:, : ) = zv0 (:,:) 102 zu_n (:, : ) = zu0 (:,:) ; zv_n (:, : ) = zv0 (:,:) 113 103 114 104 !i 115 zi1 (:,:) = 0. e0116 zi2 (:,:) = 0. e0117 zpresh(:,:) = 0. e0118 zmasst(:,:) = 0. e0105 zi1 (:,:) = 0._wp 106 zi2 (:,:) = 0._wp 107 zpresh(:,:) = 0._wp 108 zmasst(:,:) = 0._wp 119 109 !i 120 110 !!gm violant 121 zfrld(:,:) =0. e0122 zcorl(:,:) =0. e0123 zmass(:,:) =0. e0124 za1ct(:,:) =0. e0125 za2ct(:,:) =0. e0111 zfrld(:,:) =0._wp 112 zcorl(:,:) =0._wp 113 zmass(:,:) =0._wp 114 za1ct(:,:) =0._wp 115 za2ct(:,:) =0._wp 126 116 !!gm end 127 117 128 zviszeta(:,:) = 0. e0129 zviseta (:,:) = 0. e0130 131 !i zviszeta(:,0 ) = 0. e0 ; zviseta(:,0 ) = 0.e0132 !i zviszeta(:,jpj ) = 0. e0 ; zviseta(:,jpj ) = 0.e0133 !i zviszeta(:,jpj+1) = 0. e0 ; zviseta(:,jpj+1) = 0.e0118 zviszeta(:,:) = 0._wp 119 zviseta (:,:) = 0._wp 120 121 !i zviszeta(:,0 ) = 0._wp ; zviseta(:,0 ) = 0._wp 122 !i zviszeta(:,jpj ) = 0._wp ; zviseta(:,jpj ) = 0._wp 123 !i zviszeta(:,jpj+1) = 0._wp ; zviseta(:,jpj+1) = 0._wp 134 124 135 125 … … 143 133 DO ji = 1 , jpi 144 134 ! only the sinus changes its sign with the hemisphere 145 zsang(ji,jj) = SIGN( 1. e0, fcor(ji,jj) ) * sangvg ! only the sinus changes its sign with the hemisphere135 zsang(ji,jj) = SIGN( 1._wp, fcor(ji,jj) ) * sangvg ! only the sinus changes its sign with the hemisphere 146 136 ! 147 137 zmasst(ji,jj) = tms(ji,jj) * ( rhosn * hsnm(ji,jj) + rhoic * hicm(ji,jj) ) 148 138 zpresh(ji,jj) = tms(ji,jj) * pstarh * hicm(ji,jj) * EXP( -c_rhg * frld(ji,jj) ) 149 139 !!gm :: stress given at I-point (F-point for the ocean) only compute the ponderation with the ice fraction (1-frld) 150 zi1(ji,jj) = tms(ji,jj) * ( 1. 0- frld(ji,jj) )151 zi2(ji,jj) = tms(ji,jj) * ( 1. 0- frld(ji,jj) )140 zi1(ji,jj) = tms(ji,jj) * ( 1._wp - frld(ji,jj) ) 141 zi2(ji,jj) = tms(ji,jj) * ( 1._wp - frld(ji,jj) ) 152 142 END DO 153 143 END DO … … 163 153 zstms = zztms(ji,jj ) * wght(ji,jj,2,2) + zztms(ji-1,jj ) * wght(ji,jj,1,2) & 164 154 & + zztms(ji,jj-1) * wght(ji,jj,2,1) + zztms(ji-1,jj-1) * wght(ji,jj,1,1) 165 zusw = 1. 0/ MAX( zstms, epsd )155 zusw = 1._wp / MAX( zstms, epsd ) 166 156 167 157 zt11 = zztms(ji ,jj ) * zzfrld(ji ,jj ) … … 201 191 ! Gradient of the sea surface height 202 192 zgsshx = ( (ssh_m(ji ,jj ) - ssh_m(ji-1,jj ))/e1u(ji-1,jj ) & 203 & + (ssh_m(ji ,jj-1) - ssh_m(ji-1,jj-1))/e1u(ji-1,jj-1) ) * 0.5 193 & + (ssh_m(ji ,jj-1) - ssh_m(ji-1,jj-1))/e1u(ji-1,jj-1) ) * 0.5_wp 204 194 zgsshy = ( (ssh_m(ji ,jj ) - ssh_m(ji ,jj-1))/e2v(ji ,jj-1) & 205 & + (ssh_m(ji-1,jj ) - ssh_m(ji-1,jj-1))/e2v(ji-1,jj-1) ) * 0.5 195 & + (ssh_m(ji-1,jj ) - ssh_m(ji-1,jj-1))/e2v(ji-1,jj-1) ) * 0.5_wp 206 196 207 197 ! Computation of the velocity field taking into account the ice-ice interaction. … … 219 209 ! ! ==================== ! 220 210 zindu = MOD( iter , 2 ) 221 zusdtp = ( zindu * 2. 0 + ( 1.0 - zindu ) * 1.0) * REAL( nbiter ) / rdt_ice211 zusdtp = ( zindu * 2._wp + ( 1._wp - zindu ) * 1._wp ) * REAL( nbiter ) / rdt_ice 222 212 223 213 ! Computation of free drift field for free slip boundary conditions. … … 241 231 zdgi = zt12 + zt21 242 232 ztrace2 = zdgp * zdgp 243 zdeter = zt11 * zt22 - 0.25 * zdgi * zdgi233 zdeter = zt11 * zt22 - 0.25_wp * zdgi * zdgi 244 234 245 235 ! Creep limit depends on the size of the grid. 246 zdelta = MAX( SQRT( ztrace2 + ( ztrace2 - 4. 0* zdeter ) * usecc2 ), creepl)236 zdelta = MAX( SQRT( ztrace2 + ( ztrace2 - 4._wp * zdeter ) * usecc2 ), creepl) 247 237 248 238 !- Computation of viscosities. … … 256 246 DO ji = 2, fs_jpim1 ! NO vector opt. 257 247 !* zc1u , zc2v 258 zvis11 = 2. 0* zviseta (ji-1,jj-1) + dm259 zvis12 = zviseta (ji-1,jj-1) + dm260 zvis21 = zviseta (ji-1,jj-1)261 zvis22 = zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1)248 zvis11 = 2._wp * zviseta (ji-1,jj-1) + dm 249 zvis12 = zviseta (ji-1,jj-1) + dm 250 zvis21 = zviseta (ji-1,jj-1) 251 zvis22 = zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 262 252 zdiag = zvis22 * ( akappa(ji-1,jj-1,1,1) + akappa(ji-1,jj-1,2,1) ) 263 253 zs11_11 = zvis11 * akappa(ji-1,jj-1,1,1) + zdiag … … 266 256 zs22_11 = zvis11 * akappa(ji-1,jj-1,2,1) + zdiag 267 257 268 zvis11 = 2. 0* zviseta (ji,jj-1) + dm269 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1)270 zvis12 = zviseta (ji,jj-1) + dm271 zvis21 = zviseta (ji,jj-1)258 zvis11 = 2._wp * zviseta (ji,jj-1) + dm 259 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1) 260 zvis12 = zviseta (ji,jj-1) + dm 261 zvis21 = zviseta (ji,jj-1) 272 262 zdiag = zvis22 * ( -akappa(ji,jj-1,1,1) + akappa(ji,jj-1,2,1) ) 273 263 zs11_21 = -zvis11 * akappa(ji,jj-1,1,1) + zdiag … … 276 266 zs21_21 = -zvis12 * akappa(ji,jj-1,1,2) + zvis21 * akappa(ji,jj-1,2,2) 277 267 278 zvis11 = 2. 0* zviseta (ji-1,jj) + dm279 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj)280 zvis12 = zviseta (ji-1,jj) + dm281 zvis21 = zviseta (ji-1,jj)282 zdiag = zvis22 * ( akappa(ji-1,jj,1,1) + akappa(ji-1,jj,2,1) )268 zvis11 = 2._wp * zviseta (ji-1,jj) + dm 269 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj) 270 zvis12 = zviseta (ji-1,jj) + dm 271 zvis21 = zviseta (ji-1,jj) 272 zdiag = zvis22 * ( akappa(ji-1,jj,1,1) + akappa(ji-1,jj,2,1) ) 283 273 zs11_12 = zvis11 * akappa(ji-1,jj,1,1) + zdiag 284 274 zs12_12 = -zvis12 * akappa(ji-1,jj,2,2) - zvis21 * akappa(ji-1,jj,1,2) … … 286 276 zs21_12 = -zvis12 * akappa(ji-1,jj,1,2) - zvis21 * akappa(ji-1,jj,2,2) 287 277 288 zvis11 = 2. 0* zviseta (ji,jj) + dm289 zvis22 = zviszeta(ji,jj) - zviseta(ji,jj)290 zvis12 = zviseta (ji,jj) + dm291 zvis21 = zviseta (ji,jj)278 zvis11 = 2._wp * zviseta (ji,jj) + dm 279 zvis22 = zviszeta(ji,jj) - zviseta(ji,jj) 280 zvis12 = zviseta (ji,jj) + dm 281 zvis21 = zviseta (ji,jj) 292 282 zdiag = zvis22 * ( -akappa(ji,jj,1,1) + akappa(ji,jj,2,1) ) 293 283 zs11_22 = -zvis11 * akappa(ji,jj,1,1) + zdiag … … 315 305 316 306 !* zc1v , zc2v. 317 zvis11 = 2. 0* zviseta (ji-1,jj-1) + dm318 zvis22 = zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1)319 zvis12 = zviseta (ji-1,jj-1) + dm320 zvis21 = zviseta (ji-1,jj-1)307 zvis11 = 2._wp * zviseta (ji-1,jj-1) + dm 308 zvis22 = zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 309 zvis12 = zviseta (ji-1,jj-1) + dm 310 zvis21 = zviseta (ji-1,jj-1) 321 311 zdiag = zvis22 * ( akappa(ji-1,jj-1,1,2) + akappa(ji-1,jj-1,2,2) ) 322 312 zs11_11 = zvis11 * akappa(ji-1,jj-1,1,2) + zdiag … … 325 315 zs21_11 = zvis12 * akappa(ji-1,jj-1,1,1) - zvis21 * akappa(ji-1,jj-1,2,1) 326 316 327 zvis11 = 2. 0* zviseta (ji,jj-1) + dm328 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1)329 zvis12 = zviseta (ji,jj-1) + dm330 zvis21 = zviseta (ji,jj-1)317 zvis11 = 2._wp * zviseta (ji,jj-1) + dm 318 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1) 319 zvis12 = zviseta (ji,jj-1) + dm 320 zvis21 = zviseta (ji,jj-1) 331 321 zdiag = zvis22 * ( akappa(ji,jj-1,1,2) + akappa(ji,jj-1,2,2) ) 332 322 zs11_21 = zvis11 * akappa(ji,jj-1,1,2) + zdiag … … 335 325 zs21_21 = -zvis12 * akappa(ji,jj-1,1,1) - zvis21 * akappa(ji,jj-1,2,1) 336 326 337 zvis11 = 2. 0* zviseta (ji-1,jj) + dm338 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj)339 zvis12 = zviseta (ji-1,jj) + dm340 zvis21 = zviseta (ji-1,jj)327 zvis11 = 2._wp * zviseta (ji-1,jj) + dm 328 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj) 329 zvis12 = zviseta (ji-1,jj) + dm 330 zvis21 = zviseta (ji-1,jj) 341 331 zdiag = zvis22 * ( akappa(ji-1,jj,1,2) - akappa(ji-1,jj,2,2) ) 342 332 zs11_12 = zvis11 * akappa(ji-1,jj,1,2) + zdiag … … 345 335 zs21_12 = zvis12 * akappa(ji-1,jj,1,1) - zvis21 * akappa(ji-1,jj,2,1) 346 336 347 zvis11 = 2. 0* zviseta (ji,jj) + dm348 zvis22 = zviszeta(ji,jj) - zviseta(ji,jj)349 zvis12 = zviseta (ji,jj) + dm350 zvis21 = zviseta (ji,jj)337 zvis11 = 2._wp * zviseta (ji,jj) + dm 338 zvis22 = zviszeta(ji,jj) - zviseta(ji,jj) 339 zvis12 = zviseta (ji,jj) + dm 340 zvis21 = zviseta (ji,jj) 351 341 zdiag = zvis22 * ( akappa(ji,jj,1,2) - akappa(ji,jj,2,2) ) 352 342 zs11_22 = zvis11 * akappa(ji,jj,1,2) + zdiag … … 388 378 ze22 = + akappa(ji,jj-1,2,2) * zv_a(ji+1,jj) + akappa(ji,jj-1,2,1) * zu_a(ji+1,jj) 389 379 ze21 = akappa(ji,jj-1,1,1) * zv_a(ji+1,jj) - akappa(ji,jj-1,1,2) * zu_a(ji+1,jj) 390 zvis11 = 2. 0* zviseta (ji,jj-1) + dm391 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1)392 zvis12 = zviseta (ji,jj-1) + dm393 zvis21 = zviseta (ji,jj-1)380 zvis11 = 2._wp * zviseta (ji,jj-1) + dm 381 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1) 382 zvis12 = zviseta (ji,jj-1) + dm 383 zvis21 = zviseta (ji,jj-1) 394 384 zdiag = zvis22 * ( ze11 + ze22 ) 395 385 zs11_21 = zvis11 * ze11 + zdiag … … 406 396 ze21 = akappa(ji-1,jj,1,1) * ( zv_a(ji ,jj+1) - zv_a(ji-1,jj+1) ) & 407 397 & - akappa(ji-1,jj,1,2) * ( zu_a(ji ,jj+1) + zu_a(ji-1,jj+1) ) 408 zvis11 = 2. 0* zviseta (ji-1,jj) + dm409 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj)410 zvis12 = zviseta (ji-1,jj) + dm411 zvis21 = zviseta (ji-1,jj)398 zvis11 = 2._wp * zviseta (ji-1,jj) + dm 399 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj) 400 zvis12 = zviseta (ji-1,jj) + dm 401 zvis21 = zviseta (ji-1,jj) 412 402 zdiag = zvis22 * ( ze11 + ze22 ) 413 403 zs11_12 = zvis11 * ze11 + zdiag … … 424 414 ze21 = akappa(ji,jj,1,1) * ( zv_a(ji+1,jj) + zv_a(ji+1,jj+1) - zv_a(ji ,jj+1) ) & 425 415 & - akappa(ji,jj,1,2) * ( zu_a(ji+1,jj) + zu_a(ji+1,jj+1) + zu_a(ji ,jj+1) ) 426 zvis11 = 2. 0* zviseta (ji,jj) + dm427 zvis22 = zviszeta(ji,jj) - zviseta(ji,jj)428 zvis12 = zviseta (ji,jj) + dm429 zvis21 = zviseta (ji,jj)416 zvis11 = 2._wp * zviseta (ji,jj) + dm 417 zvis22 = zviszeta(ji,jj) - zviseta(ji,jj) 418 zvis12 = zviseta (ji,jj) + dm 419 zvis21 = zviseta (ji,jj) 430 420 zdiag = zvis22 * ( ze11 + ze22 ) 431 421 zs11_22 = zvis11 * ze11 + zdiag … … 443 433 ze21 = akappa(ji-1,jj-1,1,1) * ( zv_a(ji ,jj-1) - zv_a(ji-1,jj-1) - zv_a(ji-1,jj) ) & 444 434 & - akappa(ji-1,jj-1,1,2) * ( zu_a(ji ,jj-1) + zu_a(ji-1,jj-1) + zu_a(ji-1,jj) ) 445 zvis11 = 2. 0* zviseta (ji-1,jj-1) + dm446 zvis22 = zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1)447 zvis12 = zviseta (ji-1,jj-1) + dm448 zvis21 = zviseta (ji-1,jj-1)435 zvis11 = 2._wp * zviseta (ji-1,jj-1) + dm 436 zvis22 = zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 437 zvis12 = zviseta (ji-1,jj-1) + dm 438 zvis21 = zviseta (ji-1,jj-1) 449 439 zdiag = zvis22 * ( ze11 + ze22 ) 450 440 zs11_11 = zvis11 * ze11 + zdiag … … 461 451 ze21 = akappa(ji,jj-1,1,1) * ( zv_a(ji+1,jj-1) - zv_a(ji ,jj-1) ) & 462 452 & - akappa(ji,jj-1,1,2) * ( zu_a(ji+1,jj-1) + zu_a(ji ,jj-1) ) 463 zvis11 = 2. 0* zviseta (ji,jj-1) + dm464 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1)465 zvis12 = zviseta (ji,jj-1) + dm466 zvis21 = zviseta (ji,jj-1)453 zvis11 = 2._wp * zviseta (ji,jj-1) + dm 454 zvis22 = zviszeta(ji,jj-1) - zviseta(ji,jj-1) 455 zvis12 = zviseta (ji,jj-1) + dm 456 zvis21 = zviseta (ji,jj-1) 467 457 zdiag = zvis22 * ( ze11 + ze22 ) 468 458 zs11_21 = zs11_21 + zvis11 * ze11 + zdiag … … 475 465 ze22 = - akappa(ji-1,jj,2,2) * zv_a(ji-1,jj) + akappa(ji-1,jj,2,1) * zu_a(ji-1,jj) 476 466 ze21 = - akappa(ji-1,jj,1,1) * zv_a(ji-1,jj) - akappa(ji-1,jj,1,2) * zu_a(ji-1,jj) 477 zvis11 = 2. 0* zviseta (ji-1,jj) + dm478 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj)479 zvis12 = zviseta (ji-1,jj) + dm480 zvis21 = zviseta (ji-1,jj)467 zvis11 = 2._wp * zviseta (ji-1,jj) + dm 468 zvis22 = zviszeta(ji-1,jj) - zviseta(ji-1,jj) 469 zvis12 = zviseta (ji-1,jj) + dm 470 zvis21 = zviseta (ji-1,jj) 481 471 zdiag = zvis22 * ( ze11 + ze22 ) 482 472 zs11_12 = zs11_12 + zvis11 * ze11 + zdiag … … 506 496 zvr = zv_a(ji,jj) - v_oce(ji,jj) 507 497 !!!! 508 zmod = SQRT( zur*zur + zvr*zvr ) * ( 1. 0- zfrld(ji,jj) )498 zmod = SQRT( zur*zur + zvr*zvr ) * ( 1._wp - zfrld(ji,jj) ) 509 499 za = rhoco * zmod 510 500 !!!! 511 !!gm chg resul za = rhoco * SQRT( zur*zur + zvr*zvr ) * ( 1. 0- zfrld(ji,jj) )501 !!gm chg resul za = rhoco * SQRT( zur*zur + zvr*zvr ) * ( 1._wp - zfrld(ji,jj) ) 512 502 zac = za * cangvg 513 503 zmpzas = alpha * zcorl(ji,jj) + za * zsang(ji,jj) 514 504 zmassdt = zusdtp * zmass(ji,jj) 515 zcorlal = ( 1. 0- alpha ) * zcorl(ji,jj)505 zcorlal = ( 1._wp - alpha ) * zcorl(ji,jj) 516 506 517 507 za1 = zmassdt * zu0(ji,jj) + zcorlal * zv0(ji,jj) + za1ct(ji,jj) & … … 527 517 zunw = ( ( za1 + zd1 ) * zc2 + ( za2 + zd2 ) * zc1 ) * zden 528 518 zvnw = ( ( za2 + zd2 ) * zb1 - ( za1 + zd1 ) * zb2 ) * zden 529 zmask = ( 1. 0 - MAX( rzero, SIGN( rone , 1.0- zmass(ji,jj) ) ) ) * tmu(ji,jj)519 zmask = ( 1._wp - MAX( rzero, SIGN( rone , 1._wp - zmass(ji,jj) ) ) ) * tmu(ji,jj) 530 520 531 521 zu_n(ji,jj) = ( zu_a(ji,jj) + om * ( zunw - zu_a(ji,jj) ) * tmu(ji,jj) ) * zmask … … 582 572 #else 583 573 !!---------------------------------------------------------------------- 584 !! Default option Dummy module NO 2.0 LIMsea-ice model574 !! Default option Dummy module NO VP & LIM-2 sea-ice model 585 575 !!---------------------------------------------------------------------- 586 576 CONTAINS -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limrst_2.F90
r2319 r2370 121 121 CALL iom_rstput( iter, nitrst, numriw, 'fsbbq' , fsbbq (:,:) ) 122 122 #if ! defined key_lim2_vp 123 CALL iom_rstput( iter, nitrst, numriw, 'stress1_i' , stress1_i (:,:) )124 CALL iom_rstput( iter, nitrst, numriw, 'stress2_i' , stress2_i (:,:) 125 CALL iom_rstput( iter, nitrst, numriw, 'stress12_i' , stress12_i(:,:) 123 CALL iom_rstput( iter, nitrst, numriw, 'stress1_i' , stress1_i (:,:) ) ! EVP rheology 124 CALL iom_rstput( iter, nitrst, numriw, 'stress2_i' , stress2_i (:,:) ) 125 CALL iom_rstput( iter, nitrst, numriw, 'stress12_i' , stress12_i(:,:) ) 126 126 #endif 127 127 CALL iom_rstput( iter, nitrst, numriw, 'sxice' , sxice (:,:) ) -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90
r2319 r2370 2 2 !!====================================================================== 3 3 !! *** MODULE limsbc_2 *** 4 !! computation of the flux at the sea ice/ocean interface4 !! LIM-2 : updates the fluxes at the ocean surface with ice-ocean fluxes 5 5 !!====================================================================== 6 6 !! History : LIM ! 2000-01 (H. Goosse) Original code 7 7 !! 1.0 ! 2002-07 (C. Ethe, G. Madec) re-writing F90 8 8 !! 3.0 ! 2006-07 (G. Madec) surface module 9 !! 3.3 ! 2009-05 (G.Garric, C. Bricaud) addition of the lim2_evp case 9 !! 3.3 ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 10 !! - ! 2010-11 (G. Madec) ice-ocean stress computed at each ocean time-step 10 11 !!---------------------------------------------------------------------- 11 12 #if defined key_lim2 … … 13 14 !! 'key_lim2' LIM 2.0 sea-ice model 14 15 !!---------------------------------------------------------------------- 15 !! ----------------------------------------------------------------------16 !! lim_sbc_ 2 : flux at the ice / ocean interface16 !! lim_sbc_flx_2 : update mass, heat and salt fluxes at the ocean surface 17 !! lim_sbc_tau_2 : update i- and j-stresses, and its modulus at the ocean surface 17 18 !!---------------------------------------------------------------------- 18 19 USE par_oce ! ocean parameters 20 USE phycst ! physical constants 19 21 USE dom_oce ! ocean domain 22 USE dom_ice_2 ! LIM-2: ice domain 23 USE ice_2 ! LIM-2: ice variables 20 24 USE sbc_ice ! surface boundary condition: ice 21 25 USE sbc_oce ! surface boundary condition: ocean 22 USE phycst ! physical constants23 USE ice_2 ! LIM2: ice variables24 26 25 27 USE lbclnk ! ocean lateral boundary condition - MPP exchanges 26 28 USE in_out_manager ! I/O manager 27 29 USE diaar5, ONLY : lk_diaar5 28 USE iom ! 30 USE iom ! I/O library 29 31 USE albedo ! albedo parameters 30 32 USE prtctl ! Print control … … 34 36 PRIVATE 35 37 36 PUBLIC lim_sbc_2 ! called by sbc_ice_lim_2 37 38 REAL(wp) :: r1_rdtice ! constant values 39 REAL(wp) :: epsi16 = 1.e-16 ! - - 40 REAL(wp) :: rzero = 0.e0 ! - - 41 REAL(wp) :: rone = 1.e0 ! - - 38 PUBLIC lim_sbc_flx_2 ! called by sbc_ice_lim_2 39 PUBLIC lim_sbc_tau_2 ! called by sbc_ice_lim_2 40 41 REAL(wp) :: r1_rdtice ! = 1. / rdt_ice 42 REAL(wp) :: epsi16 = 1.e-16_wp ! constant values 43 REAL(wp) :: rzero = 0._wp ! - - 44 REAL(wp) :: rone = 1._wp ! - - 42 45 ! 43 REAL(wp), DIMENSION(jpi,jpj) :: soce_r, sice_r ! constant SSS and ice salinity used in levitating sea-ice case 46 REAL(wp), DIMENSION(jpi,jpj) :: soce_0, sice_0 ! constant SSS and ice salinity used in levitating sea-ice case 47 48 REAL(wp), DIMENSION(jpi,jpj) :: utau_oce, vtau_oce ! air-ocean surface i- & j-stress [N/m2] 49 REAL(wp), DIMENSION(jpi,jpj) :: tmod_io ! modulus of the ice-ocean relative velocity [m/s] 44 50 45 51 !! * Substitutions … … 48 54 !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 49 55 !! $Id$ 50 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 51 !!---------------------------------------------------------------------- 52 56 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 57 !!---------------------------------------------------------------------- 53 58 CONTAINS 54 59 55 SUBROUTINE lim_sbc_ 2( kt )60 SUBROUTINE lim_sbc_flx_2( kt ) 56 61 !!------------------------------------------------------------------- 57 62 !! *** ROUTINE lim_sbc_2 *** … … 77 82 !! Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 78 83 !!--------------------------------------------------------------------- 79 INTEGER :: kt ! number of iteration84 INTEGER, INTENT(in) :: kt ! number of iteration 80 85 !! 81 INTEGER :: ji, jj 86 INTEGER :: ji, jj ! dummy loop indices 82 87 INTEGER :: ii0, ii1, ij0, ij1 ! local integers 83 88 INTEGER :: ifvt, i1mfr, idfr, iflt ! - - 84 89 INTEGER :: ial, iadv, ifral, ifrdv ! - - 85 REAL(wp) :: zqsr, zqns, zsang, zmod, zfm ! local scalars 86 REAL(wp) :: zinda, zfons, zemp, zztmp ! - - 87 REAL(wp) :: zfrldu, zutau, zu_io ! - - 88 REAL(wp) :: zfrldv, zvtau, zv_io ! - - 89 REAL(wp), DIMENSION(jpi,jpj) :: ztio_u, ztio_v ! 2D workspace 90 REAL(wp), DIMENSION(jpi,jpj) :: ztiomi, zqnsoce ! - - 90 REAL(wp) :: zqsr, zqns, zfm ! local scalars 91 REAL(wp) :: zinda, zfons, zemp ! - - 92 REAL(wp), DIMENSION(jpi,jpj) :: zqnsoce ! 2D workspace 91 93 REAL(wp), DIMENSION(jpi,jpj,1) :: zalb, zalbp ! 2D/3D workspace 92 93 94 !!--------------------------------------------------------------------- 94 95 95 96 96 IF( kt == nit000 ) THEN 97 97 IF(lwp) WRITE(numout,*) 98 #if defined key_lim2_vp 99 IF(lwp) WRITE(numout,*) 'lim_sbc_2 : LIM 2.0 sea-ice (VP rheology) - surface boundary condition' 100 #else 101 IF(lwp) WRITE(numout,*) 'lim_sbc_2 : LIM 2.0 sea-ice (EVP rheology) - surface boundary condition' 102 #endif 103 IF(lwp) WRITE(numout,*) '~~~~~~~~~ ' 98 IF(lwp) WRITE(numout,*) 'lim_sbc_flx_2 : LIM-2 sea-ice - surface boundary condition - Mass, heat & salt fluxes' 99 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ ' 104 100 ! 105 101 r1_rdtice = 1. / rdt_ice 106 102 ! 107 soce_r(:,:) = soce 108 sice_r(:,:) = sice 109 ! 110 IF( cp_cfg == "orca" ) THEN 111 ! ocean/ice salinity in the Baltic sea 112 DO jj = 1, jpj 113 DO ji = 1, jpi 114 IF( glamt(ji,jj) >= 14. .AND. glamt(ji,jj) <= 32. .AND. gphit(ji,jj) >= 54. .AND. gphit(ji,jj) <= 66. ) THEN 115 soce_r(ji,jj) = 4.e0 116 sice_r(ji,jj) = 2.e0 117 END IF 118 END DO 119 END DO 120 ! 121 END IF 122 END IF 103 soce_0(:,:) = soce ! constant SSS and ice salinity used in levitating sea-ice case 104 sice_0(:,:) = sice 105 ! 106 IF( cp_cfg == "orca" ) THEN ! decrease ocean & ice reference salinities in the Baltic sea 107 WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND. & 108 & 54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp ) 109 soce_0(:,:) = 4._wp 110 sice_0(:,:) = 2._wp 111 END WHERE 112 ENDIF 113 ! 114 ENDIF 123 115 124 116 !------------------------------------------! … … 126 118 !------------------------------------------! 127 119 128 !!gm129 !!gm CAUTION130 !!gm re-verifies the non solar expression, especially over open ocen131 !!gm132 120 zqnsoce(:,:) = qns(:,:) 133 121 DO jj = 1, jpj … … 190 178 191 179 fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj) ! ??? 192 180 ! 193 181 qsr (ji,jj) = zqsr ! solar heat flux 194 182 qns (ji,jj) = zqns - fdtcn(ji,jj) ! non solar heat flux … … 203 191 ! mass flux at the ocean surface ! 204 192 !------------------------------------------! 205 206 !!gm207 !!gm CAUTION208 !!gm re-verifies the emp & emps expression, especially the absence of 1-frld on zfm209 !!gm210 193 DO jj = 1, jpj 211 194 DO ji = 1, jpi 212 195 ! 213 196 #if defined key_coupled 214 197 ! freshwater exchanges at the ice-atmosphere / ocean interface (coupled mode) … … 223 206 ! ! ice-covered fraction: 224 207 #endif 225 208 ! 226 209 ! computing salt exchanges at the ice/ocean interface 227 zfons = ( soce_ r(ji,jj) - sice_r(ji,jj) ) * ( rdmicif(ji,jj) * r1_rdtice )228 210 zfons = ( soce_0(ji,jj) - sice_0(ji,jj) ) * ( rdmicif(ji,jj) * r1_rdtice ) 211 ! 229 212 ! converting the salt flux from ice to a freshwater flux from ocean 230 213 zfm = zfons / ( sss_m(ji,jj) + epsi16 ) 231 214 ! 232 215 emps(ji,jj) = zemp + zfm ! surface ocean concentration/dilution effect (use on SSS evolution) 233 216 emp (ji,jj) = zemp ! surface ocean volume flux (use on sea-surface height evolution) 234 217 ! 235 218 END DO 236 219 END DO 237 220 238 IF( lk_diaar5 ) THEN 221 IF( lk_diaar5 ) THEN ! AR5 diagnostics 239 222 CALL iom_put( 'isnwmlt_cea' , rdmsnif(:,:) * r1_rdtice ) 240 CALL iom_put( 'fsal_virt_cea', soce_r(:,:) * rdmicif(:,:) * r1_rdtice ) 241 CALL iom_put( 'fsal_real_cea', - sice_r(:,:) * rdmicif(:,:) * r1_rdtice ) 242 ENDIF 243 244 !------------------------------------------! 245 ! momentum flux at the ocean surface ! 246 !------------------------------------------! 247 248 IF ( ln_limdyn ) THEN ! Update the stress over ice-over area (only in ice-dynamic case) 249 ! ! otherwise the atmosphere-ocean stress is used everywhere 250 251 ! ... ice stress over ocean with a ice-ocean rotation angle (at I-point) 252 !CDIR NOVERRCHK 253 DO jj = 1, jpj 254 !CDIR NOVERRCHK 255 DO ji = 1, jpi 256 ! ... change the cosinus angle sign in the south hemisphere 257 zsang = SIGN(1.e0, gphif(ji,jj) ) * sangvg 258 ! ... ice velocity relative to the ocean at I-point 259 zu_io = u_ice(ji,jj) - u_oce(ji,jj) 260 zv_io = v_ice(ji,jj) - v_oce(ji,jj) 261 zmod = SQRT( zu_io * zu_io + zv_io * zv_io ) 262 zztmp = rhoco * zmod 263 ! ... components of ice stress over ocean with a ice-ocean rotation angle (at I-point) 264 ztio_u(ji,jj) = zztmp * ( cangvg * zu_io - zsang * zv_io ) 265 ztio_v(ji,jj) = zztmp * ( cangvg * zv_io + zsang * zu_io ) 266 ! ... module of ice stress over ocean (at I-point) 267 ztiomi(ji,jj) = zztmp * zmod 268 ! 269 END DO 270 END DO 271 272 DO jj = 2, jpjm1 273 DO ji = 2, jpim1 ! NO vector opt. 274 ! ... components of ice-ocean stress at U and V-points (from I-point values) 275 #if defined key_lim2_vp 276 zutau = 0.5 * ( ztio_u(ji+1,jj) + ztio_u(ji+1,jj+1) ) ! VP rheology 277 zvtau = 0.5 * ( ztio_v(ji,jj+1) + ztio_v(ji+1,jj+1) ) 278 #else 279 zutau = ztio_u(ji,jj) ! EVP rheology 280 zvtau = ztio_v(ji,jj) 281 #endif 282 ! ... open-ocean (lead) fraction at U- & V-points (from T-point values) 283 zfrldu = 0.5 * ( frld(ji,jj) + frld(ji+1,jj ) ) 284 zfrldv = 0.5 * ( frld(ji,jj) + frld(ji ,jj+1) ) 285 ! ... update components of surface ocean stress (ice-cover wheighted) 286 utau(ji,jj) = zfrldu * utau(ji,jj) + ( 1. - zfrldu ) * zutau 287 vtau(ji,jj) = zfrldv * vtau(ji,jj) + ( 1. - zfrldv ) * zvtau 288 ! ... module of ice-ocean stress at T-points (from I-point values) 289 zztmp = 0.25 * ( ztiomi(ji,jj) + ztiomi(ji+1,jj) + ztiomi(ji,jj+1) + ztiomi(ji+1,jj+1) ) 290 ! ... update module of surface ocean stress (ice-cover wheighted) 291 taum(ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1. - frld(ji,jj) ) * zztmp 292 ! 293 END DO 294 END DO 295 296 ! boundary condition on the stress (utau,vtau,taum) 297 CALL lbc_lnk( utau, 'U', -1. ) 298 CALL lbc_lnk( vtau, 'V', -1. ) 299 CALL lbc_lnk( taum, 'T', 1. ) 300 223 CALL iom_put( 'fsal_virt_cea', soce_0(:,:) * rdmicif(:,:) * r1_rdtice ) 224 CALL iom_put( 'fsal_real_cea', - sice_0(:,:) * rdmicif(:,:) * r1_rdtice ) 301 225 ENDIF 302 226 … … 305 229 !-----------------------------------------------! 306 230 307 IF ( lk_cpl ) THEN231 IF( lk_cpl ) THEN ! coupled case 308 232 ! Ice surface temperature 309 233 tn_ice(:,:,1) = sist(:,:) ! sea-ice surface temperature … … 314 238 ENDIF 315 239 316 IF(ln_ctl) THEN 240 IF(ln_ctl) THEN ! control print 317 241 CALL prt_ctl(tab2d_1=qsr , clinfo1=' lim_sbc: qsr : ', tab2d_2=qns , clinfo2=' qns : ') 318 242 CALL prt_ctl(tab2d_1=emp , clinfo1=' lim_sbc: emp : ', tab2d_2=emps , clinfo2=' emps : ') … … 321 245 CALL prt_ctl(tab2d_1=fr_i , clinfo1=' lim_sbc: fr_i : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice : ') 322 246 ENDIF 323 324 END SUBROUTINE lim_sbc_2 247 ! 248 END SUBROUTINE lim_sbc_flx_2 249 250 251 SUBROUTINE lim_sbc_tau_2( kt , pu_oce, pv_oce ) 252 !!------------------------------------------------------------------- 253 !! *** ROUTINE lim_sbc_tau *** 254 !! 255 !! ** Purpose : Update the ocean surface stresses due to the ice 256 !! 257 !! ** Action : * at each ice time step (every nn_fsbc time step): 258 !! - compute the modulus of ice-ocean relative velocity 259 !! at T-point (C-grid) or I-point (B-grid) 260 !! tmod_io = rhoco * | U_ice-U_oce | 261 !! - update the modulus of stress at ocean surface 262 !! taum = frld * taum + (1-frld) * tmod_io * | U_ice-U_oce | 263 !! * at each ocean time step (each kt): 264 !! compute linearized ice-ocean stresses as 265 !! Utau = tmod_io * | U_ice - pU_oce | 266 !! using instantaneous current ocean velocity (usually before) 267 !! 268 !! NB: - the averaging operator used depends on the ice dynamics grid (cp_ice_msh='I' or 'C') 269 !! - ice-ocean rotation angle only allowed in cp_ice_msh='I' case 270 !! - here we make an approximation: taum is only computed every ice time step 271 !! This avoids mutiple average to pass from T -> U,V grids and next from U,V grids 272 !! to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... 273 !! 274 !! ** Outputs : - utau, vtau : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes 275 !! - taum : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes 276 !!--------------------------------------------------------------------- 277 INTEGER , INTENT(in) :: kt ! ocean time-step index 278 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents 279 !! 280 INTEGER :: ji, jj ! dummy loop indices 281 REAL(wp) :: zfrldu, zat_u, zu_i, zutau_ice, zu_t, zmodt ! local scalar 282 REAL(wp) :: zfrldv, zat_v, zv_i, zvtau_ice, zv_t, zmodi ! - - 283 REAL(wp) :: zsang, zumt ! - - 284 REAL(wp), DIMENSION(jpi,jpj) :: ztio_u, ztio_v ! ocean stress below sea-ice 285 !!--------------------------------------------------------------------- 286 ! 287 IF( kt == nit000 .AND. lwp ) THEN ! control print 288 WRITE(numout,*) 289 WRITE(numout,*) 'lim_sbc_tau_2 : LIM 2.0 sea-ice - surface ocean momentum fluxes' 290 WRITE(numout,*) '~~~~~~~~~~~~~ ' 291 IF( lk_lim2_vp ) THEN ; WRITE(numout,*) ' VP rheology - B-grid case' 292 ELSE ; WRITE(numout,*) ' EVP rheology - C-grid case' 293 ENDIF 294 ENDIF 295 ! 296 SELECT CASE( cp_ice_msh ) 297 ! !-----------------------! 298 CASE( 'I' ) ! B-grid ice dynamics ! I-point (i.e. F-point with sea-ice indexation) 299 ! !--=--------------------! 300 ! 301 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step) 302 !CDIR NOVERRCHK 303 DO jj = 1, jpj !* modulus of ice-ocean relative velocity at I-point 304 !CDIR NOVERRCHK 305 DO ji = 1, jpi 306 zu_i = u_ice(ji,jj) - u_oce(ji,jj) ! ice-ocean relative velocity at I-point 307 zv_i = v_ice(ji,jj) - v_oce(ji,jj) 308 tmod_io(ji,jj) = SQRT( zu_i * zu_i + zv_i * zv_i ) ! modulus of this velocity (at I-point) 309 END DO 310 END DO 311 !CDIR NOVERRCHK 312 DO jj = 1, jpjm1 !* update the modulus of stress at ocean surface (T-point) 313 !CDIR NOVERRCHK 314 DO ji = 1, jpim1 ! NO vector opt. 315 ! ! modulus of U_ice-U_oce at T-point 316 zumt = 0.25_wp * ( tmod_io(ji+1,jj) + tmod_io(ji+1,jj+1) & 317 & + tmod_io(ji,jj+1) + tmod_io(ji+1,jj+1) ) 318 ! ! update the modulus of stress at ocean surface 319 taum(ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1._wp - frld(ji,jj) ) * rhoco * zumt * zumt 320 END DO 321 END DO 322 CALL lbc_lnk( taum, 'T', 1. ) 323 ! 324 utau_oce(:,:) = utau(:,:) !* save the air-ocean stresses at ice time-step 325 vtau_oce(:,:) = vtau(:,:) 326 ! 327 ENDIF 328 ! 329 ! !== at each ocean time-step ==! 330 ! 331 ! !* ice/ocean stress WITH a ice-ocean rotation angle at I-point 332 DO jj = 2, jpj 333 zsang = SIGN( 1._wp, gphif(1,jj) ) * sangvg ! change the cosine angle sign in the SH 334 DO ji = 2, jpi ! NO vect. opt. possible 335 ! ... ice-ocean relative velocity at I-point using instantaneous surface ocean current at u- & v-pts 336 zu_i = u_ice(ji,jj) - 0.5_wp * ( pu_oce(ji-1,jj ) + pu_oce(ji-1,jj-1) ) * tmu(ji,jj) 337 zv_i = v_ice(ji,jj) - 0.5_wp * ( pv_oce(ji ,jj-1) + pv_oce(ji-1,jj-1) ) * tmu(ji,jj) 338 ! ... components of stress with a ice-ocean rotation angle 339 zmodi = rhoco * tmod_io(ji,jj) 340 ztio_u(ji,jj) = zmodi * ( cangvg * zu_i - zsang * zv_i ) 341 ztio_v(ji,jj) = zmodi * ( cangvg * zv_i + zsang * zu_i ) 342 END DO 343 END DO 344 ! !* surface ocean stresses at u- and v-points 345 DO jj = 2, jpjm1 346 DO ji = 2, jpim1 ! NO vector opt. 347 ! ! ice-ocean stress at U and V-points (from I-point values) 348 zutau_ice = 0.5_wp * ( ztio_u(ji+1,jj) + ztio_u(ji+1,jj+1) ) 349 zvtau_ice = 0.5_wp * ( ztio_v(ji,jj+1) + ztio_v(ji+1,jj+1) ) 350 ! ! open-ocean (lead) fraction at U- & V-points (from T-point values) 351 zfrldu = 0.5_wp * ( frld(ji,jj) + frld(ji+1,jj ) ) 352 zfrldv = 0.5_wp * ( frld(ji,jj) + frld(ji ,jj+1) ) 353 ! ! update the surface ocean stress (ice-cover wheighted) 354 utau(ji,jj) = zfrldu * utau_oce(ji,jj) + ( 1._wp - zfrldu ) * zutau_ice 355 vtau(ji,jj) = zfrldv * vtau_oce(ji,jj) + ( 1._wp - zfrldv ) * zvtau_ice 356 END DO 357 END DO 358 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition 359 ! 360 ! 361 ! !-----------------------! 362 CASE( 'C' ) ! C-grid ice dynamics ! U & V-points (same as in the ocean) 363 ! !--=--------------------! 364 ! 365 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step) 366 !CDIR NOVERRCHK 367 DO jj = 2, jpjm1 !* modulus of the ice-ocean velocity at T-point 368 !CDIR NOVERRCHK 369 DO ji = fs_2, fs_jpim1 370 zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj) ! 2*(U_ice-U_oce) at T-point 371 zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1) 372 zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t ) ! |U_ice-U_oce|^2 373 ! ! update the modulus of stress at ocean surface 374 taum (ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1._wp - frld(ji,jj) ) * rhoco * zmodt 375 tmod_io(ji,jj) = SQRT( zmodt ) * rhoco ! rhoco*|Uice-Uoce| 376 END DO 377 END DO 378 CALL lbc_lnk( taum, 'T', 1. ) ; CALL lbc_lnk( tmod_io, 'T', 1. ) 379 ! 380 utau_oce(:,:) = utau(:,:) !* save the air-ocean stresses at ice time-step 381 vtau_oce(:,:) = vtau(:,:) 382 ! 383 ENDIF 384 ! 385 ! !== at each ocean time-step ==! 386 ! 387 DO jj = 2, jpjm1 !* ice stress over ocean WITHOUT a ice-ocean rotation angle 388 DO ji = fs_2, fs_jpim1 389 ! ! ocean area at u- & v-points 390 zfrldu = 0.5_wp * ( frld(ji,jj) + frld(ji+1,jj ) ) 391 zfrldv = 0.5_wp * ( frld(ji,jj) + frld(ji ,jj+1) ) 392 ! ! quadratic drag formulation without rotation 393 ! ! using instantaneous surface ocean current 394 zutau_ice = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) ) 395 zvtau_ice = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) ) 396 ! ! update the surface ocean stress (ice-cover wheighted) 397 utau(ji,jj) = zfrldu * utau_oce(ji,jj) + ( 1._wp - zfrldu ) * zutau_ice 398 vtau(ji,jj) = zfrldv * vtau_oce(ji,jj) + ( 1._wp - zfrldv ) * zvtau_ice 399 END DO 400 END DO 401 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition 402 ! 403 END SELECT 404 405 IF(ln_ctl) CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau : ', mask1=umask, & 406 & tab2d_2=vtau, clinfo2=' vtau : ' , mask2=vmask ) 407 ! 408 END SUBROUTINE lim_sbc_tau_2 325 409 326 410 #else 327 411 !!---------------------------------------------------------------------- 328 !! Default option : Dummy module NO LIM 2.0 sea-ice model 329 !!---------------------------------------------------------------------- 330 CONTAINS 331 SUBROUTINE lim_sbc_2 ! Dummy routine 332 END SUBROUTINE lim_sbc_2 412 !! Default option Empty module NO LIM 2.0 sea-ice model 413 !!---------------------------------------------------------------------- 333 414 #endif 334 415 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limtrp_2.F90
r2319 r2370 7 7 !! 2.0 ! 2001-05 (G. Madec, R. Hordoir) opa norm 8 8 !! - ! 2004-01 (G. Madec, C. Ethe) F90, mpp 9 !! 3.3 ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case9 !! 3.3 ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 10 10 !!---------------------------------------------------------------------- 11 11 #if defined key_lim2 … … 41 41 REAL(wp) :: rone = 1.e0 42 42 43 44 43 !! * Substitution 45 44 # include "vectopt_loop_substitute.h90" … … 88 87 ! ice velocities at ocean U- and V-points (zui_u,zvi_v) 89 88 ! --------------------------------------- 90 zvbord = 1.0 + ( 1.0 - bound ) ! zvbord=2 no-slip, =0 free slip boundary conditions 91 #if defined key_lim2_vp 92 DO jj = 1, jpjm1 93 DO ji = 1, jpim1 ! NO vector opt. 94 zui_u(ji,jj) = ( u_ice(ji+1,jj ) + u_ice(ji+1,jj+1) ) / ( MAX( tmu(ji+1,jj ) + tmu(ji+1,jj+1), zvbord ) ) 95 zvi_v(ji,jj) = ( v_ice(ji ,jj+1) + v_ice(ji+1,jj+1) ) / ( MAX( tmu(ji ,jj+1) + tmu(ji+1,jj+1), zvbord ) ) 96 END DO 97 END DO 98 CALL lbc_lnk( zui_u, 'U', -1. ) ; CALL lbc_lnk( zvi_v, 'V', -1. ) ! Lateral boundary conditions 99 #else 100 zui_u(:,:) = u_ice(:,:) ! EVP rheology: ice (u,v) at u- and v-points 101 zvi_v(:,:) = v_ice(:,:) 102 #endif 89 IF( lk_lim2_vp ) THEN ! VP rheology : B-grid sea-ice dynamics (I-point ice velocity) 90 zvbord = 1._wp + ( 1._wp - bound ) ! zvbord=2 no-slip, =0 free slip boundary conditions 91 DO jj = 1, jpjm1 92 DO ji = 1, jpim1 ! NO vector opt. 93 zui_u(ji,jj) = ( u_ice(ji+1,jj) + u_ice(ji+1,jj+1) ) / ( MAX( tmu(ji+1,jj)+tmu(ji+1,jj+1), zvbord ) ) 94 zvi_v(ji,jj) = ( v_ice(ji,jj+1) + v_ice(ji+1,jj+1) ) / ( MAX( tmu(ji,jj+1)+tmu(ji+1,jj+1), zvbord ) ) 95 END DO 96 END DO 97 CALL lbc_lnk( zui_u, 'U', -1. ) ; CALL lbc_lnk( zvi_v, 'V', -1. ) ! Lateral boundary conditions 98 ! 99 ELSE ! EVP rheology : C-grid sea-ice dynamics (u- & v-points ice velocity) 100 zui_u(:,:) = u_ice(:,:) ! EVP rheology: ice (u,v) at u- and v-points 101 zvi_v(:,:) = v_ice(:,:) 102 ENDIF 103 103 104 104 ! CFL test for stability 105 105 ! ---------------------- 106 zcfl = 0. e0106 zcfl = 0._wp 107 107 zcfl = MAX( zcfl, MAXVAL( ABS( zui_u(1:jpim1, : ) ) * rdt_ice / e1u(1:jpim1, : ) ) ) 108 108 zcfl = MAX( zcfl, MAXVAL( ABS( zvi_v( : ,1:jpjm1) ) * rdt_ice / e2v( : ,1:jpjm1) ) ) … … 193 193 ! DO jj = 1, jpjm1 ! NB: has not to be defined on jpj line and jpi row 194 194 ! DO ji = 1 , fs_jpim1 ! vector opt. 195 ! IF( MIN( zs0a(ji,jj) , zs0a(ji+1,jj) ) == 0.e0 ) pahu(ji,jj) = 0. e0196 ! IF( MIN( zs0a(ji,jj) , zs0a(ji,jj+1) ) == 0.e0 ) pahv(ji,jj) = 0. e0195 ! IF( MIN( zs0a(ji,jj) , zs0a(ji+1,jj) ) == 0.e0 ) pahu(ji,jj) = 0._wp 196 ! IF( MIN( zs0a(ji,jj) , zs0a(ji,jj+1) ) == 0.e0 ) pahv(ji,jj) = 0._wp 197 197 ! END DO 198 198 ! END DO -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90
r2287 r2370 4 4 !! Sea-Ice dynamics : 5 5 !!====================================================================== 6 !! history : 1.0 ! 2002-08 (C. Ethe, G. Madec) original VP code 7 !! 3.0 ! 2007-03 (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle) LIM3: EVP-Cgrid 8 !!---------------------------------------------------------------------- 6 9 #if defined key_lim3 7 10 !!---------------------------------------------------------------------- … … 11 14 !! lim_dyn_init : initialization and namelist read 12 15 !!---------------------------------------------------------------------- 13 !! * Modules used14 16 USE phycst 15 17 USE in_out_manager ! I/O manager … … 30 32 PRIVATE 31 33 32 !! * Accessibility 33 PUBLIC lim_dyn ! routine called by ice_step 34 PUBLIC lim_dyn ! routine called by ice_step 34 35 35 36 !! * Substitutions 36 37 # include "vectopt_loop_substitute.h90" 37 38 !! * Module variables39 REAL(wp) :: rone = 1.e0 ! constant value40 41 38 !!---------------------------------------------------------------------- 42 39 !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 43 40 !! $Id$ 44 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 45 !!---------------------------------------------------------------------- 46 41 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 42 !!---------------------------------------------------------------------- 47 43 CONTAINS 48 44 … … 59 55 !! - computation of the stress at the ocean surface 60 56 !! - treatment of the case if no ice dynamic 61 !! History :62 !! 1.0 ! 01-04 (LIM) Original code63 !! 2.0 ! 02-08 (C. Ethe, G. Madec) F90, mpp64 !! 3.0 ! 2007-03 (M.A. Morales Maqueda, S. Bouillon, M. Vancoppenolle)65 !! LIM3, EVP, C-grid66 57 !!------------------------------------------------------------------------------------ 67 58 INTEGER, INTENT(in) :: kt ! number of iteration 68 !! * Local variables59 !! 69 60 INTEGER :: ji, jj, jl, ja ! dummy loop indices 70 61 INTEGER :: i_j1, i_jpj ! Starting/ending j-indices for rheology 71 REAL(wp) :: zcoef ! temporaryscalar62 REAL(wp) :: zcoef ! local scalar 72 63 REAL(wp), DIMENSION(jpj) :: zind ! i-averaged indicator of sea-ice 73 64 REAL(wp), DIMENSION(jpj) :: zmsk ! i-averaged of tmask … … 158 149 zv_io(:,:) = v_ice(:,:) - ssv_m(:,:) 159 150 ! frictional velocity at T-point 151 zcoef = 0.5 * cw 160 152 DO jj = 2, jpjm1 161 153 DO ji = fs_2, fs_jpim1 ! vector opt. 162 ust2s(ji,jj) = 0.5 * cw & 163 & * ( zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj) & 164 & + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1) ) * tms(ji,jj) 154 ust2s(ji,jj) = zcoef * ( zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj) & 155 & + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1) ) * tms(ji,jj) 165 156 END DO 166 157 END DO … … 171 162 DO jj = 2, jpjm1 172 163 DO ji = fs_2, fs_jpim1 ! vector opt. 173 ust2s(ji,jj) = zcoef * tms(ji,jj) *SQRT( utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj) &174 & + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1))164 ust2s(ji,jj) = zcoef * SQRT( utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj) & 165 & + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) * tms(ji,jj) 175 166 END DO 176 167 END DO 177 168 ! 178 169 ENDIF 179 180 170 CALL lbc_lnk( ust2s, 'T', 1. ) ! T-point 181 171 … … 218 208 END DO 219 209 ENDIF 220 210 ! 221 211 END SUBROUTINE lim_dyn 212 222 213 223 214 SUBROUTINE lim_dyn_init … … 232 223 !! 233 224 !! ** input : Namelist namicedyn 234 !!235 !! history :236 !! 8.5 ! 03-08 (C. Ethe) original code237 !! 9.0 ! 07-03 (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle)238 !! EVP-Cgrid-LIM3239 225 !!------------------------------------------------------------------- 240 226 NAMELIST/namicedyn/ epsd, alpha, & … … 244 230 !!------------------------------------------------------------------- 245 231 246 ! Define the initial parameters 247 ! ------------------------- 248 249 ! Read Namelist namicedyn 250 REWIND ( numnam_ice ) 251 READ ( numnam_ice , namicedyn ) 252 IF(lwp) THEN 232 REWIND( numnam_ice ) ! Read Namelist namicedyn 233 READ ( numnam_ice , namicedyn ) 234 235 IF(lwp) THEN ! control print 253 236 WRITE(numout,*) 254 237 WRITE(numout,*) 'lim_dyn_init : ice parameters for ice dynamics ' … … 272 255 WRITE(numout,*) ' timescale for elastic waves telast = ', telast 273 256 WRITE(numout,*) ' coefficient for the solution of int. stresses alphaevp = ', alphaevp 274 275 ENDIF 276 277 usecc2 = 1.0 / ( ecc * ecc ) 278 rhoco = rau0 * cw 257 ENDIF 258 ! 259 IF( angvg /= 0._wp ) THEN 260 CALL ctl_warn( 'lim_dyn_init: turning angle for oceanic stress not properly coded for EVP ', & 261 & '(see limsbc module). We force angvg = 0._wp' ) 262 angvg = 0._wp 263 ENDIF 264 265 usecc2 = 1._wp / ( ecc * ecc ) 266 rhoco = rau0 * cw 279 267 angvg = angvg * rad 280 268 sangvg = SIN( angvg ) 281 269 cangvg = COS( angvg ) 282 pstarh = pstar / 2.0270 pstarh = pstar * 0.5_wp 283 271 284 272 ! Diffusion coefficients. -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90
r2357 r2370 11 11 #if defined key_lim3 || ( defined key_lim2 && ! defined key_lim2_vp ) 12 12 !!---------------------------------------------------------------------- 13 !! 'key_lim3' LIM3 sea-ice model 13 !! 'key_lim3' OR LIM-3 sea-ice model 14 !! 'key_lim2' AND NOT 'key_lim2_vp' VP LIM-2 sea-ice model 14 15 !!---------------------------------------------------------------------- 15 16 !! lim_rhg : computes ice velocities 16 17 !!---------------------------------------------------------------------- 17 !! * Modules used 18 USE phycst 19 USE par_oce 20 USE dom_oce 21 USE sbc_oce ! Surface boundary condition: ocean fields 22 USE sbc_ice ! Surface boundary condition: ice fields 23 USE lbclnk 24 USE lib_mpp 25 USE in_out_manager ! I/O manager 26 USE prtctl ! Print control 18 USE phycst ! Physical constant 19 USE par_oce ! Ocean parameters 20 USE dom_oce ! Ocean domain 21 USE sbc_oce ! Surface boundary condition: ocean fields 22 USE sbc_ice ! Surface boundary condition: ice fields 23 USE lbclnk ! Lateral Boundary Condition / MPP link 24 USE lib_mpp ! MPP library 25 USE in_out_manager ! I/O manager 26 USE prtctl ! Print control 27 27 #if defined key_lim3 28 USE ice 29 USE dom_ice 30 USE iceini 31 USE limitd_me 32 #endif 33 #if defined key_lim2 && ! defined key_lim2_vp 28 USE ice ! LIM-3: ice variables 29 USE dom_ice ! LIM-3: ice domain 30 USE iceini ! LIM-3: ice initialisation 31 USE limitd_me ! LIM-3: 32 #else 34 33 USE ice_2 ! LIM2: ice variables 35 34 USE dom_ice_2 ! LIM2: ice domain … … 37 36 #endif 38 37 39 40 38 IMPLICIT NONE 41 39 PRIVATE 42 40 43 !! * Routine accessibility 44 PUBLIC lim_rhg ! routine called by lim_dyn 45 41 PUBLIC lim_rhg ! routine called by lim_dyn (or lim_dyn_2) 42 43 REAL(wp) :: rzero = 0._wp ! constant values 44 REAL(wp) :: rone = 1._wp ! constant values 45 46 46 !! * Substitutions 47 47 # include "vectopt_loop_substitute.h90" 48 49 !! * Module variables50 REAL(wp) :: & ! constant values51 rzero = 0.e0 , &52 rone = 1.e053 48 !!---------------------------------------------------------------------- 54 49 !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 55 50 !! $Id$ 56 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)51 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 57 52 !!---------------------------------------------------------------------- 58 59 53 CONTAINS 60 54 61 55 SUBROUTINE lim_rhg( k_j1, k_jpj ) 62 63 56 !!------------------------------------------------------------------- 64 57 !! *** SUBROUTINE lim_rhg *** … … 107 100 !! e.g. in the Canadian Archipelago 108 101 !! 109 !! ** References : Hunke and Dukowicz, JPO97 110 !! Bouillon et al., 08, in prep (update this when 111 !! published) 112 !! Vancoppenolle et al., OM08 102 !! References : Hunke and Dukowicz, JPO97 103 !! Bouillon et al., Ocean Modelling 2009 104 !! Vancoppenolle et al., Ocean Modelling 2008 105 !!------------------------------------------------------------------- 106 INTEGER, INTENT(in) :: k_j1 ! southern j-index for ice computation 107 INTEGER, INTENT(in) :: k_jpj ! northern j-index for ice computation 113 108 !! 114 !!------------------------------------------------------------------- 115 ! * Arguments 116 ! 117 INTEGER, INTENT(in) :: & 118 k_j1 , & !: southern j-index for ice computation 119 k_jpj !: northern j-index for ice computation 120 121 ! * Local variables 122 INTEGER :: ji, jj !: dummy loop indices 123 124 INTEGER :: & 125 jter !: temporary integers 126 109 INTEGER :: ji, jj ! dummy loop indices 110 INTEGER :: jter ! local integers 127 111 CHARACTER (len=50) :: charout 128 129 REAL(wp) :: & 130 zt11, zt12, zt21, zt22, & !: temporary scalars 131 ztagnx, ztagny, & !: wind stress on U/V points 132 delta ! 133 134 REAL(wp) :: & 135 za, & !: 136 zstms, & !: temporary scalar for ice strength 137 zsang, & !: temporary scalar for coriolis term 138 zmask !: mask for the computation of ice mass 112 REAL(wp) :: zt11, zt12, zt21, zt22, ztagnx, ztagny, delta ! 113 REAL(wp) :: za, zstms, zsang, zmask ! local scalars 139 114 140 115 REAL(wp),DIMENSION(jpi,jpj) :: & … … 154 129 155 130 REAL(wp) :: & 156 dtevp, & !: time step for subcycling 157 dtotel, & !: 158 ecc2, & !: square of yield ellipse eccenticity 159 z0, & !: temporary scalar 160 zr, & !: temporary scalar 161 zcca, zccb, & !: temporary scalars 162 zu_ice2, & !: 163 zv_ice1, & !: 164 zddc, zdtc, & !: temporary array for delta on corners 165 zdst, & !: temporary array for delta on centre 166 zdsshx, zdsshy, & !: term for the gradient of ocean surface 167 sigma1, sigma2 !: internal ice stress 131 dtevp, & ! time step for subcycling 132 dtotel, & ! 133 ecc2, & ! square of yield ellipse eccenticity 134 z0, & ! temporary scalar 135 zr, & ! temporary scalar 136 zcca, zccb, & ! temporary scalars 137 zu_ice2, & ! 138 zv_ice1, & ! 139 zddc, zdtc, & ! temporary array for delta on corners 140 zdst, & ! temporary array for delta on centre 141 zdsshx, zdsshy, & ! term for the gradient of ocean surface 142 sigma1, sigma2 ! internal ice stress 143 144 REAL(wp),DIMENSION(jpi,jpj) :: zf1, zf2 ! arrays for internal stresses 168 145 169 146 REAL(wp),DIMENSION(jpi,jpj) :: & 170 zf1, zf2 !: arrays for internal stresses 171 172 REAL(wp),DIMENSION(jpi,jpj) :: & 173 zdd, zdt, & !: Divergence and tension at centre of grid cells 174 zds, & !: Shear on northeast corner of grid cells 175 deltat, & !: Delta at centre of grid cells 176 deltac, & !: Delta on corners 177 zs1, zs2, & !: Diagonal stress tensor components zs1 and zs2 178 zs12 !: Non-diagonal stress tensor component zs12 147 zdd, zdt, & ! Divergence and tension at centre of grid cells 148 zds, & ! Shear on northeast corner of grid cells 149 deltat, & ! Delta at centre of grid cells 150 deltac, & ! Delta on corners 151 zs1, zs2, & ! Diagonal stress tensor components zs1 and zs2 152 zs12 ! Non-diagonal stress tensor component zs12 179 153 180 154 REAL(wp) :: & 181 zresm , & !: Maximal error on ice velocity 182 zindb , & !: ice (1) or not (0) 183 zdummy !: dummy argument 184 185 REAL(wp),DIMENSION(jpi,jpj) :: & 186 zu_ice , & !: Ice velocity on previous time step 187 zv_ice , & 188 zresr !: Local error on velocity 155 zresm , & ! Maximal error on ice velocity 156 zindb , & ! ice (1) or not (0) 157 zdummy ! dummy argument 158 159 REAL(wp),DIMENSION(jpi,jpj) :: zu_ice, zv_ice, zresr ! Local error on velocity 189 160 !!------------------------------------------------------------------- 190 191 161 #if defined key_lim2 && ! defined key_lim2_vp 192 162 # if defined key_agrif … … 205 175 ! 206 176 ! Put every vector to 0 207 zpresh (:,:) = 0.0 ; zc1(:,:) = 0.0208 zpreshc(:,:) = 0. 0209 u_ice2 (:,:) = 0.0 ; v_ice1(:,:) = 0.0210 zdd (:,:) = 0.0 ; zdt(:,:) = 0.0 ; zds(:,:) = 0.0177 zpresh (:,:) = 0._wp ; zc1 (:,:) = 0._wp 178 zpreshc(:,:) = 0._wp 179 u_ice2 (:,:) = 0._wp ; v_ice1(:,:) = 0._wp 180 zdd (:,:) = 0._wp ; zdt (:,:) = 0._wp ; zds(:,:) = 0._wp 211 181 212 182 #if defined key_lim3 213 ! Ice strength on T-points 214 CALL lim_itd_me_icestrength(ridge_scheme_swi) 183 CALL lim_itd_me_icestrength( ridge_scheme_swi ) ! LIM-3: Ice strength on T-points 215 184 #endif 216 185 217 ! Ice mass and temp variables 218 !CDIR NOVERRCHK 219 DO jj = k_j1 , k_jpj 186 !CDIR NOVERRCHK 187 DO jj = k_j1 , k_jpj ! Ice mass and temp variables 220 188 !CDIR NOVERRCHK 221 189 DO ji = 1 , jpi 222 190 zc1(ji,jj) = tms(ji,jj) * ( rhosn * vt_s(ji,jj) + rhoic * vt_i(ji,jj) ) 223 191 #if defined key_lim3 224 zpresh(ji,jj) = tms(ji,jj) * strength(ji,jj) / 2.192 zpresh(ji,jj) = tms(ji,jj) * strength(ji,jj) * 0.5_wp 225 193 #endif 226 194 ! tmi = 1 where there is ice or on land 227 tmi(ji,jj) = 1.0 - ( 1.0 - MAX( 0.0 , SIGN ( 1.0 , vt_i(ji,jj) - & 228 epsd ) ) ) * tms(ji,jj) 195 tmi(ji,jj) = 1._wp - ( 1._wp - MAX( 0._wp , SIGN ( 1._wp , vt_i(ji,jj) - epsd ) ) ) * tms(ji,jj) 229 196 END DO 230 197 END DO … … 271 238 DO ji = fs_2, fs_jpim1 272 239 273 zt11 = tms(ji ,jj)*e1t(ji,jj)274 zt12 = tms(ji+1,jj) *e1t(ji+1,jj)275 zt21 = tms(ji,jj )*e2t(ji,jj)276 zt22 = tms(ji,jj+1) *e2t(ji,jj+1)240 zt11 = tms(ji ,jj) * e1t(ji ,jj) 241 zt12 = tms(ji+1,jj) * e1t(ji+1,jj) 242 zt21 = tms(ji,jj ) * e2t(ji,jj ) 243 zt22 = tms(ji,jj+1) * e2t(ji,jj+1) 277 244 278 245 ! Leads area. 279 zfrld1(ji,jj) = ( zt12 * ( 1.0 - at_i(ji,jj) ) + & 280 & zt11 * ( 1.0 - at_i(ji+1,jj) ) ) / ( zt11 + zt12 + epsd ) 281 zfrld2(ji,jj) = ( zt22 * ( 1.0 - at_i(ji,jj) ) + & 282 & zt21 * ( 1.0 - at_i(ji,jj+1) ) ) / ( zt21 + zt22 + epsd ) 246 zfrld1(ji,jj) = ( zt12 * ( 1.0 - at_i(ji,jj) ) + zt11 * ( 1.0 - at_i(ji+1,jj) ) ) / ( zt11 + zt12 + epsd ) 247 zfrld2(ji,jj) = ( zt22 * ( 1.0 - at_i(ji,jj) ) + zt21 * ( 1.0 - at_i(ji,jj+1) ) ) / ( zt21 + zt22 + epsd ) 283 248 284 249 ! Mass, coriolis coeff. and currents 285 250 zmass1(ji,jj) = ( zt12*zc1(ji,jj) + zt11*zc1(ji+1,jj) ) / (zt11+zt12+epsd) 286 251 zmass2(ji,jj) = ( zt22*zc1(ji,jj) + zt21*zc1(ji,jj+1) ) / (zt21+zt22+epsd) 287 zcorl1(ji,jj) = zmass1(ji,jj) * ( e1t(ji+1,jj)*fcor(ji,jj) + & 288 e1t(ji,jj)*fcor(ji+1,jj) ) & 289 / (e1t(ji,jj) + e1t(ji+1,jj) + epsd ) 290 zcorl2(ji,jj) = zmass2(ji,jj) * ( e2t(ji,jj+1)*fcor(ji,jj) + & 291 e2t(ji,jj)*fcor(ji,jj+1) ) & 292 / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 252 zcorl1(ji,jj) = zmass1(ji,jj) * ( e1t(ji+1,jj)*fcor(ji,jj) + e1t(ji,jj)*fcor(ji+1,jj) ) & 253 & / ( e1t(ji,jj) + e1t(ji+1,jj) + epsd ) 254 zcorl2(ji,jj) = zmass2(ji,jj) * ( e2t(ji,jj+1)*fcor(ji,jj) + e2t(ji,jj)*fcor(ji,jj+1) ) & 255 & / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 293 256 ! 294 257 u_oce1(ji,jj) = u_oce(ji,jj) … … 314 277 ! include it later 315 278 316 zdsshx = ( ssh_m(ji+1,jj) - ssh_m(ji,jj))/e1u(ji,jj)317 zdsshy = ( ssh_m(ji,jj+1) - ssh_m(ji,jj))/e2v(ji,jj)279 zdsshx = ( ssh_m(ji+1,jj) - ssh_m(ji,jj) ) / e1u(ji,jj) 280 zdsshy = ( ssh_m(ji,jj+1) - ssh_m(ji,jj) ) / e2v(ji,jj) 318 281 319 282 za1ct(ji,jj) = ztagnx - zmass1(ji,jj) * grav * zdsshx … … 330 293 ! Time step for subcycling 331 294 dtevp = rdt_ice / nevp 332 dtotel = dtevp / ( 2. 0* telast )295 dtotel = dtevp / ( 2._wp * telast ) 333 296 334 297 !-ecc2: square of yield ellipse eccenticrity (reminder: must become a namelist parameter) 335 ecc2 = ecc *ecc298 ecc2 = ecc * ecc 336 299 337 300 !-Initialise stress tensor 338 zs1 (:,:) = stress1_i(:,:)339 zs2 (:,:) = stress2_i(:,:)301 zs1 (:,:) = stress1_i (:,:) 302 zs2 (:,:) = stress2_i (:,:) 340 303 zs12(:,:) = stress12_i(:,:) 341 304 342 ! ----------------------!305 ! !----------------------! 343 306 DO jter = 1 , nevp ! loop over jter ! 344 ! ----------------------!307 ! !----------------------! 345 308 DO jj = k_j1, k_jpj-1 346 zu_ice(:,jj) = u_ice(:,jj) ! velocity at previous time step309 zu_ice(:,jj) = u_ice(:,jj) ! velocity at previous time step 347 310 zv_ice(:,jj) = v_ice(:,jj) 348 311 END DO … … 409 372 END DO 410 373 END DO 411 CALL lbc_lnk( v_ice1(:,:), 'U', -1. ) 412 CALL lbc_lnk( u_ice2(:,:), 'V', -1. ) 374 CALL lbc_lnk( v_ice1, 'U', -1. ) ; CALL lbc_lnk( u_ice2, 'V', -1. ) ! lateral boundary cond. 413 375 414 376 !CDIR NOVERRCHK … … 418 380 419 381 !- Calculate Delta at centre of grid cells 420 zdst = ( e2u( ji , jj ) * v_ice1(ji,jj)&421 & - e2u( ji-1, jj) * v_ice1(ji-1,jj) &422 & + e1v( ji , jj ) * u_ice2(ji,jj)&423 & - e1v( ji , jj-1) * u_ice2(ji,jj-1) &424 & ) 382 zdst = ( e2u(ji , jj) * v_ice1(ji ,jj) & 383 & - e2u(ji-1, jj) * v_ice1(ji-1,jj) & 384 & + e1v(ji, jj ) * u_ice2(ji,jj ) & 385 & - e1v(ji, jj-1) * u_ice2(ji,jj-1) & 386 & ) & 425 387 & / area(ji,jj) 426 388 427 delta = SQRT( zdd(ji,jj)*zdd(ji,jj) + & 428 & ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 ) 429 deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 + & 430 (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 389 delta = SQRT( zdd(ji,jj)*zdd(ji,jj) + ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 ) 390 deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 + (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 431 391 432 392 !-Calculate stress tensor components zs1 and zs2 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90
r2287 r2370 4 4 !! computation of the flux at the sea ice/ocean interface 5 5 !!====================================================================== 6 !! History : - ! 2006-07 (M. Vancoppelle) LIM3 original code 7 !! 3.0 ! 2008-03 (C. Tallandier) surface module 8 !! - ! 2008-04 (C. Tallandier) split in 2 + new ice-ocean coupling 6 !! History : - ! 2006-07 (M. Vancoppelle) LIM3 original code 7 !! 3.0 ! 2008-03 (C. Tallandier) surface module 8 !! - ! 2008-04 (C. Tallandier) split in 2 + new ice-ocean coupling 9 !! 3.3 ! 2010-05 (G. Madec) decrease ocean & ice reference salinities in the Baltic sea 10 !! ! + simplification of the ice-ocean stress calculation 9 11 !!---------------------------------------------------------------------- 10 12 #if defined key_lim3 … … 12 14 !! 'key_lim3' LIM 3.0 sea-ice model 13 15 !!---------------------------------------------------------------------- 14 !! lim_sbc : flux at the ice / ocean interface15 !! ----------------------------------------------------------------------16 USE oce16 !! lim_sbc_flx : updates mass, heat and salt fluxes at the ocean surface 17 !! lim_sbc_tau : update i- and j-stresses, and its modulus at the ocean surface 18 !!---------------------------------------------------------------------- 17 19 USE par_oce ! ocean parameters 18 20 USE par_ice ! ice parameters … … 35 37 PUBLIC lim_sbc_tau ! called by sbc_ice_lim 36 38 37 REAL(wp) :: epsi16 = 1.e-16 ! constant values 38 REAL(wp) :: rzero = 0.e0 39 REAL(wp) :: rone = 1.e0 40 41 REAL(wp), DIMENSION(jpi,jpj) :: utau_oce, vtau_oce !: air-ocean surface i- & j-stress [N/m2] 42 REAL(wp), DIMENSION(jpi,jpj) :: tmod_io !: modulus of the ice-ocean relative velocity [m/s] 43 REAL(wp), DIMENSION(jpi,jpj) :: ssu_mb , ssv_mb !: before mean ocean surface currents [m/s] 39 REAL(wp) :: r1_rdtice ! = 1. / rdt_ice 40 REAL(wp) :: epsi16 = 1.e-16_wp ! constant values 41 REAL(wp) :: rzero = 0._wp 42 REAL(wp) :: rone = 1._wp 43 44 REAL(wp), DIMENSION(jpi,jpj) :: utau_oce, vtau_oce ! air-ocean surface i- & j-stress [N/m2] 45 REAL(wp), DIMENSION(jpi,jpj) :: tmod_io ! modulus of the ice-ocean relative velocity [m/s] 46 47 REAL(wp), DIMENSION(jpi,jpj) :: soce_0, sice_0 ! constant SSS and ice salinity used in levitating sea-ice case 44 48 45 49 !! * Substitutions … … 48 52 !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 49 53 !! $Id$ 50 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 51 !!---------------------------------------------------------------------- 52 54 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 55 !!---------------------------------------------------------------------- 53 56 CONTAINS 54 55 SUBROUTINE lim_sbc_tau( kt, kcpl )56 !!-------------------------------------------------------------------57 !! *** ROUTINE lim_sbc_tau ***58 !!59 !! ** Purpose : Update the ocean surface stresses due to the ice60 !!61 !! ** Action : - compute the ice-ocean stress depending on kcpl:62 !! fluxes at the ice-ocean interface.63 !! Case 0 : old LIM-3 way, computed at ice time-step only64 !! Case 1 : at each ocean time step the stresses are computed65 !! using the current ocean velocity (now)66 !! Case 2 : combination of half case 0 + half case 167 !!68 !! ** Outputs : - utau : sea surface i-stress (ocean referential)69 !! - vtau : sea surface j-stress (ocean referential)70 !!71 !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90.72 !! Tartinville et al. 2001 Ocean Modelling, 3, 95-108.73 !!---------------------------------------------------------------------74 INTEGER :: kt ! number of ocean iteration75 INTEGER :: kcpl ! ice-ocean coupling indicator: =0 LIM-3 old case76 ! ! =1 stresses computed using now ocean velocity77 ! ! =2 combination of 0 and 1 cases78 !!79 INTEGER :: ji, jj ! dummy loop indices80 REAL(wp) :: zfrldu, zat_u, zu_ico, zutaui, zu_u, zu_ij, zu_im1j ! temporary scalar81 REAL(wp) :: zfrldv, zat_v, zv_ico, zvtaui, zv_v, zv_ij, zv_ijm1 ! - -82 REAL(wp) :: zsang, zztmp ! - -83 REAL(wp), DIMENSION(jpi,jpj) :: ztio_u, ztio_v ! ocean stress below sea-ice84 #if defined key_coupled85 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zalb ! albedo of ice under overcast sky86 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zalbp ! albedo of ice under clear sky87 #endif88 !!---------------------------------------------------------------------89 90 IF( kt == nit000 ) THEN91 IF(lwp) WRITE(numout,*)92 IF(lwp) WRITE(numout,*) 'lim_sbc_tau : LIM 3.0 sea-ice - surface ocean momentum fluxes'93 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ '94 ENDIF95 96 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN97 !CDIR NOVERRCHK98 DO jj = 2, jpjm1 !* modulus of the ice-ocean velocity99 !CDIR NOVERRCHK100 DO ji = fs_2, fs_jpim1101 zu_ij = u_ice(ji ,jj) - ssu_m(ji ,jj) ! (i ,j)102 zu_im1j = u_ice(ji-1,jj) - ssu_m(ji-1,jj) ! (i-1,j)103 zv_ij = v_ice(ji,jj ) - ssv_m(ji,jj ) ! (i,j )104 zv_ijm1 = v_ice(ji,jj-1) - ssv_m(ji,jj-1) ! (i,j-1)105 zztmp = 0.5 * ( zu_ij * zu_ij + zu_im1j * zu_im1j & ! mean of the square values instead106 & + zv_ij * zv_ij + zv_ijm1 * zv_ijm1 ) ! of the square of the mean values...107 ! WARNING, here we make an approximation for case 1 and 2: taum is not computed at each time108 ! step but only every nn_fsbc time step. This avoid mutiple avarage to pass from T -> U,V grids109 ! and next from U,V grids to T grid. Taum is used in tke, which should not be too sensitive to110 ! this approximaton...111 taum(ji,jj) = ( 1. - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zztmp112 tmod_io(ji,jj) = SQRT( zztmp )113 END DO114 END DO115 CALL lbc_lnk( taum, 'T', 1. ) ; CALL lbc_lnk( tmod_io, 'T', 1. )116 ENDIF117 118 SELECT CASE( kcpl )119 ! !--------------------------------!120 CASE( 0 ) ! LIM 3 old stress computation ! (at ice timestep only)121 ! !--------------------------------!122 ! !* ice stress over ocean with a ice-ocean rotation angle123 DO jj = 1, jpjm1124 zsang = SIGN( 1.e0, gphif(1,jj) ) * sangvg ! change the sinus angle sign in the south hemisphere125 DO ji = 1, fs_jpim1126 zu_u = u_ice(ji,jj) - u_oce(ji,jj) ! ice velocity relative to the ocean127 zv_v = v_ice(ji,jj) - v_oce(ji,jj)128 ! ! quadratic drag formulation with rotation129 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1)130 zutaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_u - zsang * zv_v )131 zvtaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_v + zsang * zu_u )132 ! ! bound for too large stress values133 ! IMPORTANT: the exponential below prevents numerical oscillations in the ice-ocean stress134 ! This is not physically based. A cleaner solution is offer in CASE kcpl=2135 ztio_u(ji,jj) = zutaui * EXP( - ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) )136 ztio_v(ji,jj) = zvtaui * EXP( - ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) )137 END DO138 END DO139 DO jj = 2, jpjm1 ! stress at the surface of the ocean140 DO ji = fs_2, fs_jpim1 ! vertor opt.141 zfrldu = 0.5 * ( ato_i(ji,jj) + ato_i(ji+1,jj ) ) ! open-ocean fraction at U- & V-points (from T-point values)142 zfrldv = 0.5 * ( ato_i(ji,jj) + ato_i(ji ,jj+1) )143 ! ! update surface ocean stress144 utau(ji,jj) = zfrldu * utau(ji,jj) + ( 1. - zfrldu ) * ztio_u(ji,jj)145 vtau(ji,jj) = zfrldv * vtau(ji,jj) + ( 1. - zfrldv ) * ztio_v(ji,jj)146 END DO147 END DO148 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition149 150 !151 ! !--------------------------------!152 CASE( 1 ) ! Use the now velocity ! (at each ocean timestep)153 ! !--------------------------------!154 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN155 utau_oce(:,:) = utau(:,:) !* save the air-ocean stresses at ice time-step156 vtau_oce(:,:) = vtau(:,:)157 ENDIF158 !159 DO jj = 2, jpjm1 !* ice stress over ocean with a ice-ocean rotation angle160 zsang = SIGN(1.e0, gphif(1,jj-1) ) * sangvg ! change the sinus angle sign in the south hemisphere161 DO ji = fs_2, fs_jpim1162 zat_u = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5 ! ice area at u and V-points163 zat_v = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5164 ! ! (u,v) ice-ocean velocity at (U,V)-point, resp.165 zu_u = u_ice(ji,jj) - ub(ji,jj,1)166 zv_v = v_ice(ji,jj) - vb(ji,jj,1)167 ! ! quadratic drag formulation with rotation168 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1)169 zutaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_u - zsang * zv_v )170 zvtaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_v + zsang * zu_u )171 ! ! stress at the ocean surface172 utau(ji,jj) = ( 1.- zat_u ) * utau_oce(ji,jj) + zat_u * zutaui173 vtau(ji,jj) = ( 1.- zat_v ) * vtau_oce(ji,jj) + zat_v * zvtaui174 END DO175 END DO176 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition177 178 !179 ! !--------------------------------!180 CASE( 2 ) ! mixed 0 and 2 cases ! (at each ocean timestep)181 ! !--------------------------------!182 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN183 utau_oce(:,:) = utau (:,:) !* save the air-ocean stresses at ice time-step184 vtau_oce(:,:) = vtau (:,:)185 ssu_mb (:,:) = ssu_m(:,:) !* save the ice-ocean velocity at ice time-step186 ssv_mb (:,:) = ssv_m(:,:)187 ENDIF188 DO jj = 2, jpjm1 !* ice stress over ocean with a ice-ocean rotation angle189 zsang = SIGN(1.e0, gphif(1,jj-1) ) * sangvg190 DO ji = fs_2, fs_jpim1191 zat_u = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5 ! ice area at u and V-points192 zat_v = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5193 !194 zu_ico = u_ice(ji,jj) - 0.5 * ( ub(ji,jj,1) + ssu_mb(ji,jj) ) ! ice-oce velocity using ub and ssu_mb195 zv_ico = v_ice(ji,jj) - 0.5 * ( vb(ji,jj,1) + ssv_mb(ji,jj) )196 ! ! quadratic drag formulation with rotation197 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1)198 zutaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_ico - zsang * zv_ico )199 zvtaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_ico + zsang * zu_ico )200 !201 utau(ji,jj) = ( 1.-zat_u ) * utau_oce(ji,jj) + zat_u * zutaui ! stress at the ocean surface202 vtau(ji,jj) = ( 1.-zat_v ) * vtau_oce(ji,jj) + zat_v * zvtaui203 END DO204 END DO205 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition206 !207 END SELECT208 209 IF(ln_ctl) CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau : ', mask1=umask, &210 & tab2d_2=vtau, clinfo2=' vtau : ' , mask2=vmask )211 !212 END SUBROUTINE lim_sbc_tau213 214 57 215 58 SUBROUTINE lim_sbc_flx( kt ) … … 235 78 !! Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 236 79 !!--------------------------------------------------------------------- 237 INTEGER :: kt ! number of iteration80 INTEGER, INTENT(in) :: kt ! number of iteration 238 81 !! 239 82 INTEGER :: ji, jj ! dummy loop indices … … 254 97 IF(lwp) WRITE(numout,*) 'lim_sbc_flx : LIM 3.0 sea-ice - heat salt and mass ocean surface fluxes' 255 98 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 99 ! 100 r1_rdtice = 1. / rdt_ice 101 ! 102 soce_0(:,:) = soce 103 sice_0(:,:) = sice 104 ! 105 IF( cp_cfg == "orca" ) THEN ! decrease ocean & ice reference salinities in the Baltic sea 106 WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND. & 107 & 54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp ) 108 soce_0(:,:) = 4._wp 109 sice_0(:,:) = 2._wp 110 END WHERE 111 ENDIF 112 ! 256 113 ENDIF 257 114 … … 298 155 ! fscmbq and ffltbif are obsolete 299 156 ! & + iflt * ffltbif(ji,jj) !!! only if one category is used 300 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) / rdt_ice &301 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) / rdt_ice&157 & + ifral * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice & 158 & + ifrdv * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) * r1_rdtice & 302 159 & + fhmec(ji,jj) & ! new contribution due to snow melt in ridging!! 303 160 & + fheat_rpo(ji,jj) & ! contribution from ridge formation … … 340 197 WRITE(numout,*) ' qcmif : ', qcmif(jiindx,jjindx) 341 198 WRITE(numout,*) ' qldif : ', qldif(jiindx,jjindx) 342 WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindx,jjindx) / rdt_ice343 WRITE(numout,*) ' qldif / dt: ', qldif(jiindx,jjindx) / rdt_ice199 WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindx,jjindx) * r1_rdtice 200 WRITE(numout,*) ' qldif / dt: ', qldif(jiindx,jjindx) * r1_rdtice 344 201 WRITE(numout,*) ' ifrdv : ', ifrdv 345 202 WRITE(numout,*) ' qfvbq : ', qfvbq(jiindx,jjindx) 346 203 WRITE(numout,*) ' qdtcn : ', qdtcn(jiindx,jjindx) 347 WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) / rdt_ice348 WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) / rdt_ice204 WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) * r1_rdtice 205 WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) * r1_rdtice 349 206 WRITE(numout,*) ' ' 350 207 WRITE(numout,*) ' fdtcn : ', fdtcn(jiindx,jjindx) … … 380 237 ! & - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) ) & ! remov. snow precip over ice 381 238 & - sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) ) & ! remov. snow precip over ice 382 & - rdmsnif(ji,jj) / rdt_ice& ! freshwaterflux due to snow melting239 & - rdmsnif(ji,jj) * r1_rdtice & ! freshwaterflux due to snow melting 383 240 ! new contribution from snow falling when ridging 384 241 & + fmmec(ji,jj) … … 386 243 ! computing salt exchanges at the ice/ocean interface 387 244 ! sice should be the same as computed with the ice model 388 zfons = ( soce - sice ) * ( rdmicif(ji,jj) / rdt_ice )245 zfons = ( soce_0(ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice 389 246 ! SOCE 390 zfons = ( sss_m (ji,jj) - sice ) * ( rdmicif(ji,jj) / rdt_ice )247 zfons = ( sss_m (ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice 391 248 392 249 !CT useless ! salt flux for constant salinity … … 446 303 END SUBROUTINE lim_sbc_flx 447 304 305 306 SUBROUTINE lim_sbc_tau( kt , pu_oce, pv_oce ) 307 !!------------------------------------------------------------------- 308 !! *** ROUTINE lim_sbc_tau *** 309 !! 310 !! ** Purpose : Update the ocean surface stresses due to the ice 311 !! 312 !! ** Action : * at each ice time step (every nn_fsbc time step): 313 !! - compute the modulus of ice-ocean relative velocity 314 !! (*rho*Cd) at T-point (C-grid) or I-point (B-grid) 315 !! tmod_io = rhoco * | U_ice-U_oce | 316 !! - update the modulus of stress at ocean surface 317 !! taum = frld * taum + (1-frld) * tmod_io * | U_ice-U_oce | 318 !! * at each ocean time step (every kt): 319 !! compute linearized ice-ocean stresses as 320 !! Utau = tmod_io * | U_ice - pU_oce | 321 !! using instantaneous current ocean velocity (usually before) 322 !! 323 !! NB: - ice-ocean rotation angle no more allowed 324 !! - here we make an approximation: taum is only computed every ice time step 325 !! This avoids mutiple average to pass from T -> U,V grids and next from U,V grids 326 !! to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... 327 !! 328 !! ** Outputs : - utau, vtau : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes 329 !! - taum : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes 330 !!--------------------------------------------------------------------- 331 INTEGER , INTENT(in) :: kt ! ocean time-step index 332 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pu_oce, pv_oce ! surface ocean currents 333 !! 334 INTEGER :: ji, jj ! dummy loop indices 335 REAL(wp) :: zat_u, zutau_ice, zu_t, zmodt ! local scalar 336 REAL(wp) :: zat_v, zvtau_ice, zv_t ! - - 337 !!--------------------------------------------------------------------- 338 339 IF( kt == nit000 ) THEN 340 IF(lwp) WRITE(numout,*) 341 IF(lwp) WRITE(numout,*) 'lim_sbc_tau : LIM-3 sea-ice - surface ocean momentum fluxes' 342 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 343 ENDIF 344 345 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN !== Ice time-step only ==! (i.e. surface module time-step) 346 !CDIR NOVERRCHK 347 DO jj = 2, jpjm1 !* update the modulus of stress at ocean surface (T-point) 348 !CDIR NOVERRCHK 349 DO ji = fs_2, fs_jpim1 350 ! ! 2*(U_ice-U_oce) at T-point 351 zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj) 352 zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1) 353 ! ! |U_ice-U_oce|^2 354 zmodt = 0.25_wp * ( zu_t * zu_t + zv_t * zv_t ) 355 ! ! update the ocean stress modulus 356 taum(ji,jj) = ( 1._wp - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zmodt 357 tmod_io(ji,jj) = rhoco * SQRT( zmodt ) ! rhoco * |U_ice-U_oce| at T-point 358 END DO 359 END DO 360 CALL lbc_lnk( taum, 'T', 1. ) ; CALL lbc_lnk( tmod_io, 'T', 1. ) 361 ! 362 utau_oce(:,:) = utau(:,:) !* save the air-ocean stresses at ice time-step 363 vtau_oce(:,:) = vtau(:,:) 364 ! 365 ENDIF 366 ! 367 ! !== every ocean time-step ==! 368 ! 369 DO jj = 2, jpjm1 !* update the stress WITHOUT a ice-ocean rotation angle 370 DO ji = fs_2, fs_jpim1 ! Vect. Opt. 371 zat_u = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5_wp ! ice area at u and V-points 372 zat_v = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5_wp 373 ! ! linearized quadratic drag formulation 374 zutau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) ) 375 zvtau_ice = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) ) 376 ! ! stresses at the ocean surface 377 utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_ice 378 vtau(ji,jj) = ( 1._wp - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice 379 END DO 380 END DO 381 CALL lbc_lnk( utau, 'U', -1. ) ; CALL lbc_lnk( vtau, 'V', -1. ) ! lateral boundary condition 382 ! 383 IF(ln_ctl) CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau : ', mask1=umask, & 384 & tab2d_2=vtau, clinfo2=' vtau : ' , mask2=vmask ) 385 ! 386 END SUBROUTINE lim_sbc_tau 387 448 388 #else 449 389 !!---------------------------------------------------------------------- -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90
r2287 r2370 9 9 #if defined key_lim3 || defined key_lim2 10 10 !!---------------------------------------------------------------------- 11 !! 'key_lim2' or 'key_lim3' : LIM 2.0 or 3.0sea-ice model11 !! 'key_lim2' or 'key_lim3' : LIM-2 or LIM-3 sea-ice model 12 12 !!---------------------------------------------------------------------- 13 13 USE par_oce ! ocean parameters … … 23 23 24 24 # if defined key_lim2 25 LOGICAL , PUBLIC, PARAMETER :: lk_lim2 = .TRUE. !: LIM-2 ice model 26 LOGICAL , PUBLIC, PARAMETER :: lk_lim3 = .FALSE. !: no LIM-3 27 CHARACTER(len=1), PUBLIC :: cigr_type = 'I' !: 'I'-grid ice-velocity (B-grid lower left corner) 25 LOGICAL , PUBLIC, PARAMETER :: lk_lim2 = .TRUE. !: LIM-2 ice model 26 LOGICAL , PUBLIC, PARAMETER :: lk_lim3 = .FALSE. !: no LIM-3 27 # if defined key_lim2_vp 28 CHARACTER(len=1), PUBLIC, PARAMETER :: cp_ice_msh = 'I' !: VP : 'I'-grid ice-velocity (B-grid lower left corner) 29 # else 30 CHARACTER(len=1), PUBLIC, PARAMETER :: cp_ice_msh = 'C' !: EVP: 'C'-grid ice-velocity 31 # endif 28 32 # endif 29 33 # if defined key_lim3 30 LOGICAL , PUBLIC, PARAMETER :: lk_lim2 = .FALSE.!: no LIM-231 LOGICAL , PUBLIC, PARAMETER :: lk_lim3 = .TRUE.!: LIM-3 ice model32 CHARACTER(len=1), PUBLIC :: cigr_type = 'C'!: 'C'-grid ice-velocity34 LOGICAL , PUBLIC, PARAMETER :: lk_lim2 = .FALSE. !: no LIM-2 35 LOGICAL , PUBLIC, PARAMETER :: lk_lim3 = .TRUE. !: LIM-3 ice model 36 CHARACTER(len=1), PUBLIC, PARAMETER :: cp_ice_msh = 'C' !: 'C'-grid ice-velocity 33 37 # endif 34 38 … … 41 45 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpl) :: alb_ice !: albedo of ice 42 46 43 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: utau_ice !: u-stress over ice (I-p oint for LIM2 or U,V-point for LIM3) [N/m2]44 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: vtau_ice !: v-stress over ice (I-p oint for LIM2 or U,V-point for LIM3) [N/m2]45 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fr1_i0 !: 1st fraction of sol. rad. which penetrateinside the ice cover46 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fr2_i0 !: 2nd fraction of sol. rad. which penetrateinside the ice cover47 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: utau_ice !: u-stress over ice (I-pt for VP or U,V-pts for EVP) [N/m2] 48 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: vtau_ice !: v-stress over ice (I-pt for VP or U,V-pts for EVP) [N/m2] 49 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fr1_i0 !: 1st fraction of Qsr which penetrates inside the ice cover 50 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: fr2_i0 !: 2nd fraction of Qsr which penetrates inside the ice cover 47 51 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: emp_ice !: solid freshwater budget over ice: sublivation - snow 48 52 … … 55 59 !! Default option NO LIM 2.0 or 3.0 sea-ice model 56 60 !!---------------------------------------------------------------------- 57 LOGICAL , PUBLIC, PARAMETER :: lk_lim2 58 LOGICAL , PUBLIC, PARAMETER :: lk_lim3 59 CHARACTER(len=1), PUBLIC :: cigr_type= '-' !: no grid ice-velocity61 LOGICAL , PUBLIC, PARAMETER :: lk_lim2 = .FALSE. !: no LIM-2 ice model 62 LOGICAL , PUBLIC, PARAMETER :: lk_lim3 = .FALSE. !: no LIM-3 ice model 63 CHARACTER(len=1), PUBLIC, PARAMETER :: cp_ice_msh = '-' !: no grid ice-velocity 60 64 #endif 61 65 … … 63 67 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 64 68 !! $Id$ 65 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 66 !!---------------------------------------------------------------------- 67 69 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 68 70 !!====================================================================== 69 71 END MODULE sbc_ice -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90
r2287 r2370 4 4 !! Surface module : variables defined in core memory 5 5 !!====================================================================== 6 !! History : 3.0 ! 2006-06 (G. Madec) Original code 7 !! - ! 2008-08 (G. Madec) namsbc moved from sbcmod 8 !! 3.3 ! 2010-04 (M. Leclair, G. Madec) Forcing averaged over 2 time steps 6 !! History : 3.0 ! 2006-06 (G. Madec) Original code 7 !! - ! 2008-08 (G. Madec) namsbc moved from sbcmod 8 !! 3.3 ! 2010-04 (M. Leclair, G. Madec) Forcing averaged over 2 time steps 9 !! - ! 2010-11 (G. Madec) ice-ocean stress always computed at each ocean time-step 9 10 !!---------------------------------------------------------------------- 10 11 USE par_oce ! ocean parameters … … 30 31 ! !: = 1 global mean of e-p-r set to zero at each nn_fsbc time step 31 32 ! !: = 2 annual global mean of e-p-r set to zero 32 INTEGER , PUBLIC :: nn_ico_cpl = 0 !: ice-ocean coupling indicator33 ! !: = 0 LIM-3 old case34 ! !: = 1 stresses computed using now ocean velocity35 ! !: = 2 combination of 0 and 1 cases36 33 37 34 !!---------------------------------------------------------------------- 38 35 !! Ocean Surface Boundary Condition fields 39 36 !!---------------------------------------------------------------------- 40 LOGICAL , PUBLIC :: lhftau = .FALSE. 37 LOGICAL , PUBLIC :: lhftau = .FALSE. !: HF tau used in TKE: mean(stress module) - module(mean stress) 41 38 !! !! now ! before !! 42 39 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: utau , utau_b !: sea surface i-stress (ocean referential) [N/m2] … … 53 50 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: emp_tot !: total E-P over ocean and ice [Kg/m2/s] 54 51 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: rnf , rnf_b !: river runoff [Kg/m2/s] 55 ! - ML - begin52 !! 56 53 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpts) :: sbc_tsc, sbc_tsc_b !: sbc content trend [K.m/s] 57 54 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: qsr_hc , qsr_hc_b !: heat content trend due to qsr flux [K.m/s] 58 ! - ML - end55 !! 59 56 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: tprecip !: total precipitation [Kg/m2/s] 60 57 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: sprecip !: solid precipitation [Kg/m2/s] … … 63 60 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: atm_co2 !: atmospheric pCO2 [ppm] 64 61 #endif 65 !!$ REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: rrunoff !: runoff66 !!$ REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: calving !: calving67 62 68 63 !!---------------------------------------------------------------------- … … 79 74 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 80 75 !! $Id$ 81 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)76 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 82 77 !!====================================================================== 83 78 END MODULE sbc_oce -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90
r2287 r2370 83 83 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 84 84 !! $Id$ 85 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)85 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 86 86 !!---------------------------------------------------------------------- 87 88 87 CONTAINS 89 88 … … 137 136 138 137 ! (NB: frequency positive => hours, negative => months) 139 ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation!140 ! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs!141 sn_utau = FLD_N( 'utau' , 24 , 'utau' , .true. , .false. , 'yearly' , '' , '')142 sn_vtau = FLD_N( 'vtau' , 24 , 'vtau' , .true. , .false. , 'yearly' , '' , '')143 sn_wndm = FLD_N( 'mwnd10m' , 24 , 'm_10' , .true. , .false. , 'yearly' , '' , '')144 sn_tair = FLD_N( 'tair10m' , 24 , 't_10' , .false. , .false. , 'yearly' , '' , '')145 sn_humi = FLD_N( 'humi10m' , 24 , 'q_10' , .false. , .false. , 'yearly' , '' , '')146 sn_ccov = FLD_N( 'ccover' , -1 , 'cloud' , .true. , .false. , 'yearly' , '' , '')147 sn_prec = FLD_N( 'precip' , -1 , 'precip' , .true. , .false. , 'yearly' , '' , '')138 ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation ! 139 ! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs ! 140 sn_utau = FLD_N( 'utau' , 24 , 'utau' , .true. , .false. , 'yearly' , '' , '' ) 141 sn_vtau = FLD_N( 'vtau' , 24 , 'vtau' , .true. , .false. , 'yearly' , '' , '' ) 142 sn_wndm = FLD_N( 'mwnd10m', 24 , 'm_10' , .true. , .false. , 'yearly' , '' , '' ) 143 sn_tair = FLD_N( 'tair10m', 24 , 't_10' , .false. , .false. , 'yearly' , '' , '' ) 144 sn_humi = FLD_N( 'humi10m', 24 , 'q_10' , .false. , .false. , 'yearly' , '' , '' ) 145 sn_ccov = FLD_N( 'ccover' , -1 , 'cloud' , .true. , .false. , 'yearly' , '' , '' ) 146 sn_prec = FLD_N( 'precip' , -1 , 'precip' , .true. , .false. , 'yearly' , '' , '' ) 148 147 149 148 REWIND( numnam ) ! ... read in namlist namsbc_clio … … 473 472 INTEGER :: ijpl ! number of ice categories (size of 3rd dim of input arrays) 474 473 !! 475 REAL(wp) :: zcoef, zmt1, zmt2, zmt3, ztatm3 474 REAL(wp) :: zcoef, zmt1, zmt2, zmt3, ztatm3 ! temporary scalars 476 475 REAL(wp) :: ztaevbk, zind1, zind2, zind3, ztamr ! - - 477 476 REAL(wp) :: zesi, zqsati, zdesidt ! - - … … 496 495 SELECT CASE( cd_grid ) 497 496 CASE( 'C' ) ! C-grid ice dynamics 498 ! Change from wind speed to wind stress over OCEAN (cao is used) 499 zcoef = cai / cao 497 zcoef = cai / cao ! Change from air-sea stress to air-ice stress 500 498 !CDIR COLLAPSE 501 499 DO jj = 1 , jpj … … 505 503 END DO 506 504 END DO 507 CASE( 'B' ) ! B-grid ice dynamics 508 ! stress from ocean U- and V-points to ice U,V point 509 !CDIR COLLAPSE 510 DO jj = 2, jpj 511 DO ji = 2, jpi ! B grid : no vector opt. 512 p_taui(ji,jj) = 0.5 * ( utau(ji-1,jj ) + utau(ji-1,jj-1) ) 513 p_tauj(ji,jj) = 0.5 * ( vtau(ji ,jj-1) + vtau(ji-1,jj-1) ) 505 CASE( 'I' ) ! I-grid ice dynamics: I-point (i.e. F-point lower-left corner) 506 zcoef = 0.5_wp * cai / cao ! Change from air-sea stress to air-ice stress!CDIR COLLAPSE 507 DO jj = 2, jpj ! stress from ocean U- and V-points to ice U,V point 508 DO ji = 2, jpi ! I-grid : no vector opt. 509 p_taui(ji,jj) = zcoef * ( utau(ji-1,jj ) + utau(ji-1,jj-1) ) 510 p_tauj(ji,jj) = zcoef * ( vtau(ji ,jj-1) + vtau(ji-1,jj-1) ) 514 511 END DO 515 512 END DO 516 CALL lbc_lnk( p_taui(:,:), 'I', -1. ) ! I-point (i.e. ice U-V point) 517 CALL lbc_lnk( p_tauj(:,:), 'I', -1. ) ! I-point (i.e. ice U-V point) 513 CALL lbc_lnk( p_taui(:,:), 'I', -1. ) ; CALL lbc_lnk( p_tauj(:,:), 'I', -1. ) ! I-point 518 514 END SELECT 519 515 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r2291 r2370 75 75 !! NEMO/OPA 3.3 , NEMO-consortium (2010) 76 76 !! $Id$ 77 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)77 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 78 78 !!---------------------------------------------------------------------- 79 79 CONTAINS … … 133 133 ! 134 134 ! (NB: frequency positive => hours, negative => months) 135 ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation!136 ! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs!137 sn_wndi = FLD_N( 'uwnd10m' , 24 , 'u_10' , .false. , .false. , 'yearly' , '' , '')138 sn_wndj = FLD_N( 'vwnd10m' , 24 , 'v_10' , .false. , .false. , 'yearly' , '' , '')139 sn_qsr = FLD_N( 'qsw' , 24 , 'qsw' , .false. , .false. , 'yearly' , '' , '')140 sn_qlw = FLD_N( 'qlw' , 24 , 'qlw' , .false. , .false. , 'yearly' , '' , '')141 sn_tair = FLD_N( 'tair10m' , 24 , 't_10' , .false. , .false. , 'yearly' , '' , '')142 sn_humi = FLD_N( 'humi10m' , 24 , 'q_10' , .false. , .false. , 'yearly' , '' , '')143 sn_prec = FLD_N( 'precip' , -1 , 'precip' , .true. , .false. , 'yearly' , '' , '')144 sn_snow = FLD_N( 'snow' , -1 , 'snow' , .true. , .false. , 'yearly' , '' , '')145 sn_tdif = FLD_N( 'taudif' , 24 , 'taudif' , .true. , .false. , 'yearly' , '' , '')135 ! ! file ! frequency ! variable ! time intep ! clim ! 'yearly' or ! weights ! rotation ! 136 ! ! name ! (hours) ! name ! (T/F) ! (T/F) ! 'monthly' ! filename ! pairs ! 137 sn_wndi = FLD_N( 'uwnd10m', 24 , 'u_10' , .false. , .false. , 'yearly' , '' , '' ) 138 sn_wndj = FLD_N( 'vwnd10m', 24 , 'v_10' , .false. , .false. , 'yearly' , '' , '' ) 139 sn_qsr = FLD_N( 'qsw' , 24 , 'qsw' , .false. , .false. , 'yearly' , '' , '' ) 140 sn_qlw = FLD_N( 'qlw' , 24 , 'qlw' , .false. , .false. , 'yearly' , '' , '' ) 141 sn_tair = FLD_N( 'tair10m', 24 , 't_10' , .false. , .false. , 'yearly' , '' , '' ) 142 sn_humi = FLD_N( 'humi10m', 24 , 'q_10' , .false. , .false. , 'yearly' , '' , '' ) 143 sn_prec = FLD_N( 'precip' , -1 , 'precip' , .true. , .false. , 'yearly' , '' , '' ) 144 sn_snow = FLD_N( 'snow' , -1 , 'snow' , .true. , .false. , 'yearly' , '' , '' ) 145 sn_tdif = FLD_N( 'taudif' , 24 , 'taudif' , .true. , .false. , 'yearly' , '' , '' ) 146 146 ! 147 147 REWIND( numnam ) ! read in namlist namsbc_core … … 396 396 !! caution : the net upward water flux has with mm/day unit 397 397 !!--------------------------------------------------------------------- 398 REAL(wp), INTENT(in ), DIMENSION(:,:,:) 399 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) 400 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) 401 REAL(wp), INTENT(in ), DIMENSION(:,:,:) 402 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) 403 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) 404 REAL(wp), INTENT( out), DIMENSION(:,:,:) 405 REAL(wp), INTENT( out), DIMENSION(:,:,:) 406 REAL(wp), INTENT( out), DIMENSION(:,:,:) 407 REAL(wp), INTENT( out), DIMENSION(:,:,:) 408 REAL(wp), INTENT( out), DIMENSION(:,:,:) 409 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) 410 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) 411 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) 412 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) 413 CHARACTER(len=1), INTENT(in ) 414 INTEGER, INTENT(in ) 398 REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: pst ! ice surface temperature (>0, =rt0 over land) [Kelvin] 399 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) :: pui ! ice surface velocity (i- and i- components [m/s] 400 REAL(wp), INTENT(in ), DIMENSION(jpi,jpj) :: pvi ! at I-point (B-grid) or U & V-point (C-grid) 401 REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: palb ! ice albedo (clear sky) (alb_ice_cs) [%] 402 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: p_taui ! i- & j-components of surface ice stress [N/m2] 403 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: p_tauj ! at I-point (B-grid) or U & V-point (C-grid) 404 REAL(wp), INTENT( out), DIMENSION(:,:,:) :: p_qns ! non solar heat flux over ice (T-point) [W/m2] 405 REAL(wp), INTENT( out), DIMENSION(:,:,:) :: p_qsr ! solar heat flux over ice (T-point) [W/m2] 406 REAL(wp), INTENT( out), DIMENSION(:,:,:) :: p_qla ! latent heat flux over ice (T-point) [W/m2] 407 REAL(wp), INTENT( out), DIMENSION(:,:,:) :: p_dqns ! non solar heat sensistivity (T-point) [W/m2] 408 REAL(wp), INTENT( out), DIMENSION(:,:,:) :: p_dqla ! latent heat sensistivity (T-point) [W/m2] 409 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: p_tpr ! total precipitation (T-point) [Kg/m2/s] 410 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: p_spr ! solid precipitation (T-point) [Kg/m2/s] 411 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: p_fr1 ! 1sr fraction of qsr penetration in ice (T-point) [%] 412 REAL(wp), INTENT( out), DIMENSION(jpi,jpj) :: p_fr2 ! 2nd fraction of qsr penetration in ice (T-point) [%] 413 CHARACTER(len=1), INTENT(in ) :: cd_grid ! ice grid ( C or B-grid) 414 INTEGER, INTENT(in ) :: pdim ! number of ice categories 415 415 !! 416 416 INTEGER :: ji, jj, jl ! dummy loop indices … … 449 449 ! ----------------------------------------------------------------------------- ! 450 450 SELECT CASE( cd_grid ) 451 CASE( ' B' ) ! B-grid ice dynamics : I-point (i.e. F-point with sea-ice indexation)451 CASE( 'I' ) ! B-grid ice dynamics : I-point (i.e. F-point with sea-ice indexation) 452 452 ! and scalar wind at T-point ( = | U10m - U_ice | ) (masked) 453 453 #if defined key_vectopt_loop … … 605 605 !! 9.0 ! 05-08 (L. Brodeau) Rewriting and optimization 606 606 !!---------------------------------------------------------------------- 607 !! * Arguments608 609 607 REAL(wp), INTENT(in) :: zu ! altitude of wind measurement [m] 610 608 REAL(wp), INTENT(in), DIMENSION(jpi,jpj) :: & … … 646 644 grav = 9.8, & ! gravity 647 645 kappa = 0.4 ! von Karman s constant 648 646 !!---------------------------------------------------------------------- 649 647 !! * Start 650 648 !! Air/sea differences … … 770 768 grav = 9.8, & ! gravity 771 769 kappa = 0.4 ! von Karman's constant 772 770 !!---------------------------------------------------------------------- 773 771 !! * Start 774 772 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
r2287 r2370 4 4 !! Surface Boundary Condition : momentum, heat and freshwater fluxes in coupled mode 5 5 !!====================================================================== 6 !! History : 2.0 ! 06-2007(R. Redler, N. Keenlyside, W. Park) Original code split into flxmod & taumod7 !! 3.0 ! 02-2008(G. Madec, C Talandier) surface module8 !! 3.1 ! 02-2009(G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface6 !! History : 2.0 ! 2007-06 (R. Redler, N. Keenlyside, W. Park) Original code split into flxmod & taumod 7 !! 3.0 ! 2008-02 (G. Madec, C Talandier) surface module 8 !! 3.1 ! 2009_02 (G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface 9 9 !!---------------------------------------------------------------------- 10 10 #if defined key_oasis3 || defined key_oasis4 … … 788 788 !! ** Method : transform the received stress from the atmosphere into 789 789 !! an atmosphere-ice stress in the (i,j) ocean referencial 790 !! and at the velocity point of the sea-ice model (c igr_type):790 !! and at the velocity point of the sea-ice model (cp_ice_msh): 791 791 !! 'C'-grid : i- (j-) components given at U- (V-) point 792 !! ' B'-grid: both components given at I-point792 !! 'I'-grid : B-grid lower-left corner: both components given at I-point 793 793 !! 794 794 !! The received stress are : … … 803 803 !! first as 2 components on the sphere 804 804 !! second as 2 components oriented along the local grid 805 !! third as 2 components on the c igr_typepoint805 !! third as 2 components on the cp_ice_msh point 806 806 !! 807 807 !! In 'oce and ice' case, only one vector stress field … … 809 809 !! so that it is now defined as (i,j) components given at U- 810 810 !! and V-points, respectively. Therefore, here only the third 811 !! transformation is done and only if the ice-grid is a ' B'-grid.812 !! 813 !! ** Action : return ptau_i, ptau_j, the stress over the ice at c igr_typepoint811 !! transformation is done and only if the ice-grid is a 'I'-grid. 812 !! 813 !! ** Action : return ptau_i, ptau_j, the stress over the ice at cp_ice_msh point 814 814 !!---------------------------------------------------------------------- 815 815 REAL(wp), INTENT(out), DIMENSION(jpi,jpj) :: p_taui ! i- & j-components of atmos-ice stress [N/m2] … … 872 872 ! 873 873 ! j+1 j -----V---F 874 ! ice stress on ice velocity point (c igr_type)! |874 ! ice stress on ice velocity point (cp_ice_msh) ! | 875 875 ! (C-grid ==>(U,V) or B-grid ==> I or F) j | T U 876 876 ! | | … … 879 879 ! i-1 i i 880 880 ! i i+1 (for I) 881 SELECT CASE ( c igr_type)881 SELECT CASE ( cp_ice_msh ) 882 882 ! 883 883 CASE( 'I' ) ! B-grid ==> I … … 1258 1258 END DO 1259 1259 CASE( 'weighted oce and ice' ) 1260 SELECT CASE ( c igr_type)1260 SELECT CASE ( cp_ice_msh ) 1261 1261 CASE( 'C' ) ! Ocean and Ice on C-grid ==> T 1262 1262 DO jj = 2, jpjm1 … … 1293 1293 CALL lbc_lnk( zitx1, 'T', -1. ) ; CALL lbc_lnk( zity1, 'T', -1. ) 1294 1294 CASE( 'mixed oce-ice' ) 1295 SELECT CASE ( c igr_type)1295 SELECT CASE ( cp_ice_msh ) 1296 1296 CASE( 'C' ) ! Ocean and Ice on C-grid ==> T 1297 1297 DO jj = 2, jpjm1 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r2287 r2370 5 5 !! & covered area using LIM sea-ice model 6 6 !! Sea-Ice model : LIM 3.0 Sea ice model time-stepping 7 !!====================================================================== 8 !! History : 2.0 ! 2006-12 (M. Vancoppenolle) Original code 9 !! 3.0 ! 2008-02 (C. Talandier) Surface module from icestp.F90 10 !! 9.0 ! 2008-04 (G. Madec) sltyle and lim_ctl routine 7 !!===================================================================== 8 !! History : 2.0 ! 2006-12 (M. Vancoppenolle) Original code 9 !! 3.0 ! 2008-02 (C. Talandier) Surface module from icestp.F90 10 !! - ! 2008-04 (G. Madec) sltyle and lim_ctl routine 11 !! 3.3 ! 2010-11 (G. Madec) ice-ocean stress always computed at each ocean time-step 11 12 !!---------------------------------------------------------------------- 12 13 #if defined key_lim3 … … 20 21 USE oce ! ocean dynamics and tracers 21 22 USE dom_oce ! ocean space and time domain 22 USE lib_mpp 23 USE lib_mpp ! MPP library 23 24 USE par_ice ! sea-ice parameters 24 USE ice 25 USE iceini 26 USE dom_ice 25 USE ice ! LIM-3: ice variables 26 USE iceini ! LIM-3: ice initialisation 27 USE dom_ice ! LIM-3: ice domain 27 28 28 29 USE sbc_oce ! Surface boundary condition: ocean fields … … 30 31 USE sbcblk_core ! Surface boundary condition: CORE bulk 31 32 USE sbcblk_clio ! Surface boundary condition: CLIO bulk 32 USE albedo 33 USE albedo ! ocean & ice albedo 33 34 34 35 USE phycst ! Define parameters for the routines … … 46 47 USE limvar ! Ice variables switch 47 48 48 USE lbclnk 49 USE lbclnk ! lateral boundary condition - MPP link 49 50 USE iom ! I/O manager library 50 51 USE in_out_manager ! I/O manager … … 56 57 PUBLIC sbc_ice_lim ! routine called by sbcmod.F90 57 58 58 CHARACTER(len=1) :: cl_grid = 'C' ! type of grid used in ice dynamics59 60 59 !! * Substitutions 61 60 # include "domzgr_substitute.h90" … … 64 63 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 65 64 !! $Id$ 66 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 67 !!---------------------------------------------------------------------- 68 65 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 66 !!---------------------------------------------------------------------- 69 67 CONTAINS 70 68 71 SUBROUTINE sbc_ice_lim( kt, kblk , kico)69 SUBROUTINE sbc_ice_lim( kt, kblk ) 72 70 !!--------------------------------------------------------------------- 73 71 !! *** ROUTINE sbc_ice_lim *** … … 91 89 INTEGER, INTENT(in) :: kt ! ocean time step 92 90 INTEGER, INTENT(in) :: kblk ! type of bulk (=3 CLIO, =4 CORE) 93 INTEGER, INTENT(in) :: kico ! ice-ocean stress treatment94 91 !! 95 92 INTEGER :: jl ! loop index … … 142 139 & qla_ice , dqns_ice , dqla_ice , & 143 140 & tprecip , sprecip , & 144 & fr1_i0 , fr2_i0 , c l_grid, jpl )141 & fr1_i0 , fr2_i0 , cp_ice_msh, jpl ) 145 142 ! 146 143 CASE( 4 ) ! CORE bulk formulation … … 149 146 & qla_ice , dqns_ice , dqla_ice , & 150 147 & tprecip , sprecip , & 151 & fr1_i0 , fr2_i0 , c l_grid, jpl )148 & fr1_i0 , fr2_i0 , cp_ice_msh, jpl ) 152 149 END SELECT 153 150 … … 160 157 ! ! Store previous ice values 161 158 !!gm : remark old_... should becomes ...b as tn versus tb 162 old_a_i (:,:,:) = a_i(:,:,:) ! ice area163 old_e_i (:,:,:,:) = e_i(:,:,:,:) ! ice thermal energy164 old_v_i (:,:,:) = v_i(:,:,:) ! ice volume165 old_v_s (:,:,:) = v_s(:,:,:) ! snow volume166 old_e_s (:,:,:,:) = e_s(:,:,:,:) ! snow thermal energy167 old_smv_i(:,:,:) = smv_i(:,:,:)! salt content168 old_oa_i (:,:,:) = oa_i(:,:,:)! areal age content159 old_a_i (:,:,:) = a_i (:,:,:) ! ice area 160 old_e_i (:,:,:,:) = e_i (:,:,:,:) ! ice thermal energy 161 old_v_i (:,:,:) = v_i (:,:,:) ! ice volume 162 old_v_s (:,:,:) = v_s (:,:,:) ! snow volume 163 old_e_s (:,:,:,:) = e_s (:,:,:,:) ! snow thermal energy 164 old_smv_i(:,:,:) = smv_i(:,:,:) ! salt content 165 old_oa_i (:,:,:) = oa_i (:,:,:) ! areal age content 169 166 170 167 ! ! intialisation to zero !!gm is it truly necessary ??? 171 d_a_i_thd (:,:,:) = 0.e0 ; d_a_i_trp(:,:,:) = 0.e0172 d_v_i_thd (:,:,:) = 0.e0 ; d_v_i_trp(:,:,:) = 0.e0173 d_e_i_thd (:,:,:,:) = 0.e0 ; d_e_i_trp(:,:,:,:) = 0.e0174 d_v_s_thd (:,:,:) = 0.e0 ; d_v_s_trp(:,:,:) = 0.e0175 d_e_s_thd (:,:,:,:) = 0.e0 ; d_e_s_trp(:,:,:,:) = 0.e0176 d_smv_i_thd(:,:,:) = 0.e0 ; d_smv_i_trp(:,:,:)= 0.e0177 d_oa_i_thd (:,:,:) = 0.e0 ; d_oa_i_trp(:,:,:)= 0.e0178 ! 179 fseqv (:,:)= 0.e0180 fsbri (:,:) = 0.e0 ;fsalt_res(:,:) = 0.e0168 d_a_i_thd (:,:,:) = 0.e0 ; d_a_i_trp (:,:,:) = 0.e0 169 d_v_i_thd (:,:,:) = 0.e0 ; d_v_i_trp (:,:,:) = 0.e0 170 d_e_i_thd (:,:,:,:) = 0.e0 ; d_e_i_trp (:,:,:,:) = 0.e0 171 d_v_s_thd (:,:,:) = 0.e0 ; d_v_s_trp (:,:,:) = 0.e0 172 d_e_s_thd (:,:,:,:) = 0.e0 ; d_e_s_trp (:,:,:,:) = 0.e0 173 d_smv_i_thd(:,:,:) = 0.e0 ; d_smv_i_trp(:,:,:) = 0.e0 174 d_oa_i_thd (:,:,:) = 0.e0 ; d_oa_i_trp (:,:,:) = 0.e0 175 ! 176 fseqv (:,:) = 0.e0 177 fsbri (:,:) = 0.e0 ; fsalt_res(:,:) = 0.e0 181 178 fsalt_rpo(:,:) = 0.e0 182 fhmec (:,:) = 0.e0 ; fhbri(:,:)= 0.e0183 fmmec (:,:) = 0.e0 ;fheat_res(:,:) = 0.e0184 fheat_rpo(:,:) = 0.e0 ; focea2D(:,:)= 0.e0185 fsup2D (:,:)= 0.e0179 fhmec (:,:) = 0.e0 ; fhbri (:,:) = 0.e0 180 fmmec (:,:) = 0.e0 ; fheat_res(:,:) = 0.e0 181 fheat_rpo(:,:) = 0.e0 ; focea2D (:,:) = 0.e0 182 fsup2D (:,:) = 0.e0 186 183 ! 187 diag_sni_gr(:,:) = 0.e0 ; diag_lat_gr(:,:) = 0.e0188 diag_bot_gr(:,:) = 0.e0 ; diag_dyn_gr(:,:) = 0.e0189 diag_bot_me(:,:) = 0.e0 ; diag_sur_me(:,:) = 0.e0184 diag_sni_gr(:,:) = 0.e0 ; diag_lat_gr(:,:) = 0.e0 185 diag_bot_gr(:,:) = 0.e0 ; diag_dyn_gr(:,:) = 0.e0 186 diag_bot_me(:,:) = 0.e0 ; diag_sur_me(:,:) = 0.e0 190 187 ! dynamical invariants 191 delta_i(:,:) = 0.e0 ; divu_i (:,:) = 0.e0 ;shear_i(:,:) = 0.e0188 delta_i(:,:) = 0.e0 ; divu_i(:,:) = 0.e0 ; shear_i(:,:) = 0.e0 192 189 193 190 CALL lim_rst_opn( kt ) ! Open Ice restart file … … 196 193 ! 197 194 #if ! defined key_c1d 198 ! Ice dynamics & transport (not in 1D case)195 ! Ice dynamics & transport (not in 1D case) 199 196 CALL lim_dyn( kt ) ! Ice dynamics ( rheology/dynamics ) 200 197 CALL lim_trp( kt ) ! Ice transport ( Advection/diffusion ) … … 204 201 CALL lim_itd_me ! Mechanical redistribution ! (ridging/rafting) 205 202 #endif 206 !207 203 ! ! Ice thermodynamics 208 204 CALL lim_var_glo2eqv ! equivalent variables … … 216 212 CALL lim_itd_th( kt ) ! Remap ice categories, lateral accretion ! 217 213 ! 218 ! ! Global variables update |214 ! ! Global variables update 219 215 CALL lim_var_agg( 1 ) ! requested by limupdate 220 216 CALL lim_update ! Global variables update … … 223 219 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx, 2, ' - Final state - ' ) ! control print 224 220 ! 225 ! ! Fluxes of mass and heat to the ocean | 226 CALL lim_sbc_flx( kt ) ! Ice/Ocean heat freshwater/salt fluxes 227 IF( ln_limdyn .AND. kico == 0 ) & ! Ice/Ocean stresses (only in ice-dynamic case) 228 & CALL lim_sbc_tau( kt, kico ) ! otherwise the atm.-ocean stresses are used everywhere 221 CALL lim_sbc_flx( kt ) ! Update surface ocean mass, heat and salt fluxes 229 222 ! 230 223 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx, 3, ' - Final state lim_sbc - ' ) ! control print … … 239 232 IF( ln_nicep ) CALL lim_ctl ! alerts in case of model crash 240 233 ! 241 ENDIF ! End sea-ice time step only 242 243 ! !--------------------------! 244 ! Ice/Ocean stresses (nn_ico_cpl=1 or 2 cases) ! at all ocean time step ! 245 ! !--------------------------! 246 IF( ln_limdyn .AND. kico /= 0 ) & 247 & CALL lim_sbc_tau( kt, kico ) 248 !!gm remark, in this case the ocean-ice stress is not saved in diag call above ..... find a solution!!! 234 ENDIF ! End sea-ice time step only 235 236 ! !--------------------------! 237 ! ! at all ocean time step ! 238 ! !--------------------------! 239 ! 240 ! ! Update surface ocean stresses (only in ice-dynamic case) 241 ! ! otherwise the atm.-ocean stresses are used everywhere 242 IF( ln_limdyn ) CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) ) ! using before instantaneous surf. currents 243 244 !!gm remark, the ocean-ice stress is not saved in ice diag call above ..... find a solution!!! 249 245 ! 250 246 END SUBROUTINE sbc_ice_lim … … 664 660 !!---------------------------------------------------------------------- 665 661 CONTAINS 666 SUBROUTINE sbc_ice_lim ( kt, kblk , kico) ! Dummy routine667 WRITE(*,*) 'sbc_ice_lim: You should not have seen this print! error?', kt, kblk , kico662 SUBROUTINE sbc_ice_lim ( kt, kblk ) ! Dummy routine 663 WRITE(*,*) 'sbc_ice_lim: You should not have seen this print! error?', kt, kblk 668 664 END SUBROUTINE sbc_ice_lim 669 665 #endif -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90
r2319 r2370 2 2 !!====================================================================== 3 3 !! *** MODULE sbcice_lim_2 *** 4 !! Surface module : update surface ocean boundary condition over ice 5 !! covered area using LIM sea-ice model 6 !! Sea-Ice model : LIM 2.0 Sea ice model time-stepping 4 !! Surface module : update surface ocean boundary condition over ice covered area using LIM sea-ice model 5 !! Sea-Ice model : LIM-2 Sea ice model time-stepping 7 6 !!====================================================================== 8 7 !! History : 1.0 ! 06-2006 (G. Madec) from icestp_2.F90 … … 12 11 #if defined key_lim2 13 12 !!---------------------------------------------------------------------- 14 !! 'key_lim2' : LIM 2.0 sea-ice model 15 !!---------------------------------------------------------------------- 16 !! sbc_ice_lim_2 : sea-ice model time-stepping and 17 !! update ocean sbc over ice-covered area 18 !!---------------------------------------------------------------------- 19 USE oce ! ocean dynamics and tracers 20 USE dom_oce ! ocean space and time domain 21 USE lib_mpp 13 !! 'key_lim2' : LIM-2 sea-ice model 14 !!---------------------------------------------------------------------- 15 !! sbc_ice_lim_2 : sea-ice model time-stepping and update ocean sbc over ice-covered area 16 !!---------------------------------------------------------------------- 17 USE oce ! ocean dynamics and tracers 18 USE dom_oce ! ocean space and time domain 22 19 USE ice_2 23 20 USE par_ice_2 … … 25 22 USE dom_ice_2 26 23 27 USE sbc_oce ! Surface boundary condition: ocean fields28 USE sbc_ice ! Surface boundary condition: ice fields29 USE sbcblk_core ! Surface boundary condition: CORE bulk30 USE sbcblk_clio ! Surface boundary condition: CLIO bulk31 USE sbccpl ! Surface boundary condition: coupled interface24 USE sbc_oce ! Surface boundary condition: ocean fields 25 USE sbc_ice ! Surface boundary condition: ice fields 26 USE sbcblk_core ! Surface boundary condition: CORE bulk 27 USE sbcblk_clio ! Surface boundary condition: CLIO bulk 28 USE sbccpl ! Surface boundary condition: coupled interface 32 29 USE albedo 33 30 34 USE phycst ! Define parameters for the routines35 USE eosbn2 ! equation of state31 USE phycst ! Define parameters for the routines 32 USE eosbn2 ! equation of state 36 33 USE limdyn_2 37 34 USE limtrp_2 38 35 USE limdmp_2 39 36 USE limthd_2 40 USE limsbc_2 ! sea surface boundary condition37 USE limsbc_2 ! sea surface boundary condition 41 38 USE limdia_2 42 39 USE limwri_2 43 40 USE limrst_2 44 41 45 USE lbclnk 46 USE iom ! I/O manager library 47 USE in_out_manager ! I/O manager 48 USE prtctl ! Print control 42 USE lbclnk ! lateral boundary condition - MPP link 43 USE lib_mpp ! MPP library 44 USE iom ! I/O manager library 45 USE in_out_manager ! I/O manager 46 USE prtctl ! Print control 49 47 50 48 IMPLICIT NONE … … 59 57 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 60 58 !! $Id$ 61 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 62 !!---------------------------------------------------------------------- 63 59 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 60 !!---------------------------------------------------------------------- 64 61 CONTAINS 65 62 … … 97 94 IF(lwp) WRITE(numout,*) 'sbc_ice_lim_2 : update ocean surface boudary condition' 98 95 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ via Louvain la Neuve Ice Model (LIM) time stepping' 99 96 ! 100 97 CALL ice_init_2 101 102 98 ENDIF 103 99 104 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 105 ! 100 ! !----------------------! 101 IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN ! Ice time-step only ! 102 ! !----------------------! 103 ! Bulk Formulea ! 104 !----------------! 106 105 ! ... mean surface ocean current at ice dynamics point 107 ! B-grid dynamics : I-point 108 DO jj = 2, jpj 109 DO ji = 2, jpi ! B grid : no vector opt. 110 u_oce(ji,jj) = 0.5 * ( ssu_m(ji-1,jj ) + ssu_m(ji-1,jj-1) ) * tmu(ji,jj) 111 v_oce(ji,jj) = 0.5 * ( ssv_m(ji ,jj-1) + ssv_m(ji-1,jj-1) ) * tmu(ji,jj) 106 SELECT CASE( cp_ice_msh ) 107 CASE( 'I' ) !== B-grid ice dynamics : I-point (i.e. F-point with sea-ice indexation) 108 DO jj = 2, jpj 109 DO ji = 2, jpi ! NO vector opt. possible 110 u_oce(ji,jj) = 0.5_wp * ( ssu_m(ji-1,jj ) + ssu_m(ji-1,jj-1) ) * tmu(ji,jj) 111 v_oce(ji,jj) = 0.5_wp * ( ssv_m(ji ,jj-1) + ssv_m(ji-1,jj-1) ) * tmu(ji,jj) 112 END DO 112 113 END DO 113 END DO 114 CALL lbc_lnk( u_oce, 'I', -1. ) ! I-point (i.e. F-point with ice indices) 115 CALL lbc_lnk( v_oce, 'I', -1. ) ! I-point (i.e. F-point with ice indices) 114 CALL lbc_lnk( u_oce, 'I', -1. ) ! I-point (i.e. F-point with ice indices) 115 CALL lbc_lnk( v_oce, 'I', -1. ) ! I-point (i.e. F-point with ice indices) 116 ! 117 CASE( 'C' ) !== C-grid ice dynamics : U & V-points (same as ocean) 118 u_oce(:,:) = ssu_m(:,:) ! mean surface ocean current at ice velocity point 119 v_oce(:,:) = ssv_m(:,:) 120 ! 121 END SELECT 116 122 117 123 ! ... masked sea surface freezing temperature [Kelvin] (set to rt0 over land) … … 142 148 & qla_ice , dqns_ice , dqla_ice , & 143 149 & tprecip , sprecip , & 144 & fr1_i0 , fr2_i0 , c l_grid, jpl )150 & fr1_i0 , fr2_i0 , cp_ice_msh , jpl ) 145 151 146 152 CASE( 4 ) ! CORE bulk formulation … … 149 155 & qla_ice , dqns_ice , dqla_ice , & 150 156 & tprecip , sprecip , & 151 & fr1_i0 , fr2_i0 , c l_grid, jpl )157 & fr1_i0 , fr2_i0 , cp_ice_msh , jpl ) 152 158 CASE( 5 ) ! Coupled formulation : atmosphere-ice stress only (fluxes provided after ice dynamics) 153 159 CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) … … 170 176 ! Ice model step ! 171 177 ! ---------------- ! 172 numit = numit + nn_fsbc !Ice model time step173 174 CALL lim_rst_opn_2 ( kt )! Open Ice restart file178 numit = numit + nn_fsbc ! Ice model time step 179 180 CALL lim_rst_opn_2 ( kt ) ! Open Ice restart file 175 181 #if ! defined key_c1d 176 ! Ice dynamics & transport (not in 1D case)177 CALL lim_dyn_2 ( kt )! Ice dynamics ( rheology/dynamics )178 CALL lim_trp_2 ( kt )! Ice transport ( Advection/diffusion )179 IF( ln_limdmp ) CALL lim_dmp_2 ( kt )! Ice damping182 ! Ice dynamics & transport (except in 1D case) 183 CALL lim_dyn_2 ( kt ) ! Ice dynamics ( rheology/dynamics ) 184 CALL lim_trp_2 ( kt ) ! Ice transport ( Advection/diffusion ) 185 IF( ln_limdmp ) CALL lim_dmp_2 ( kt ) ! Ice damping 180 186 #endif 181 187 #if defined key_coupled 182 IF( ksbc == 5 ) CALL sbc_cpl_ice_flx( reshape( frld, (/jpi,jpj,1/) ), & 183 & qns_tot, qns_ice, qsr_tot , qsr_ice, & 184 & emp_tot, emp_ice, dqns_ice, sprecip, & 188 ! ! Ice surface fluxes in coupled mode 189 IF( ksbc == 5 ) CALL sbc_cpl_ice_flx( reshape( frld, (/jpi,jpj,1/) ), & 190 & qns_tot, qns_ice, qsr_tot , qsr_ice, & 191 & emp_tot, emp_ice, dqns_ice, sprecip, & 185 192 ! optional arguments, used only in 'mixed oce-ice' case 186 & 193 & palbi = zalb_ice_cs, psst = sst_m, pist = sist ) 187 194 #endif 188 189 CALL lim_sbc_2 ( kt ) ! Ice/Ocean Mass & Heat fluxes195 CALL lim_thd_2 ( kt ) ! Ice thermodynamics 196 CALL lim_sbc_flx_2 ( kt ) ! update surface ocean mass, heat & salt fluxes 190 197 191 198 IF( ( MOD( kt+nn_fsbc-1, ninfo ) == 0 .OR. ntmoy == 1 ) .AND. .NOT. lk_mpp ) & 192 & 199 & CALL lim_dia_2 ( kt ) ! Ice Diagnostics 193 200 # if ! defined key_iomput 194 201 CALL lim_wri_2 ( kt ) ! Ice outputs 195 202 # endif 196 IF( lrst_ice )CALL lim_rst_write_2( kt ) ! Ice restart file203 IF( lrst_ice ) CALL lim_rst_write_2( kt ) ! Ice restart file 197 204 ! 198 ENDIF 205 ENDIF ! End sea-ice time step only 206 ! 207 ! !--------------------------! 208 ! ! at all ocean time step ! 209 ! !--------------------------! 210 ! 211 ! ! Update surface ocean stresses (only in ice-dynamic case) 212 ! ! otherwise the atm.-ocean stresses are used everywhere 213 IF( ln_limdyn ) CALL lim_sbc_tau_2( kt, ub(:,:,1), vb(:,:,1) ) ! using before instantaneous surf. currents 199 214 ! 200 215 END SUBROUTINE sbc_ice_lim_2 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
r2305 r2370 4 4 !! Surface module : provide to the ocean its surface boundary condition 5 5 !!====================================================================== 6 !! History : 3.0 ! 2006-07 (G. Madec) Original code 7 !! 3.1 ! 2008-08 (S. Masson, A. Caubel, E. Maisonnave, G. Madec) coupled interface 8 !! 3.3 ! 2010-04 (M. Leclair, G. Madec) Forcing averaged over 2 time steps 9 !! 3.3 ! 2010-10 (S. Masson) add diurnal cycle 10 !! 3.3 ! 09-2010 (D. Storkey) add ice boundary conditions (BDY) 6 !! History : 3.0 ! 2006-07 (G. Madec) Original code 7 !! 3.1 ! 2008-08 (S. Masson, A. Caubel, E. Maisonnave, G. Madec) coupled interface 8 !! 3.3 ! 2010-04 (M. Leclair, G. Madec) Forcing averaged over 2 time steps 9 !! 3.3 ! 2010-10 (S. Masson) add diurnal cycle 10 !! 3.3 ! 2010-09 (D. Storkey) add ice boundary conditions (BDY) 11 !! - ! 2010-11 (G. Madec) ice-ocean stress always computed at each ocean time-step 11 12 !!---------------------------------------------------------------------- 12 13 … … 74 75 !! 75 76 NAMELIST/namsbc/ nn_fsbc, ln_ana , ln_flx, ln_blk_clio, ln_blk_core, ln_cpl , & 76 & nn_ice , ln_dm2dc, ln_rnf, ln_ssr , nn_fwb , nn_ico_cpl77 & nn_ice , ln_dm2dc, ln_rnf, ln_ssr , nn_fwb 77 78 !!---------------------------------------------------------------------- 78 79 … … 107 108 WRITE(numout,*) ' Misc. options of sbc : ' 108 109 WRITE(numout,*) ' ice management in the sbc (=0/1/2/3) nn_ice = ', nn_ice 109 WRITE(numout,*) ' ice-ocean stress computation (=0/1/2) nn_ico_cpl = ', nn_ico_cpl110 110 WRITE(numout,*) ' daily mean to diurnal cycle qsr ln_dm2dc = ', ln_dm2dc 111 111 WRITE(numout,*) ' runoff / runoff mouths ln_rnf = ', ln_rnf … … 238 238 239 239 SELECT CASE( nn_ice ) ! Update heat and freshwater fluxes over sea-ice areas 240 CASE( 1 ) ; CALL sbc_ice_if ( kt ) 240 CASE( 1 ) ; CALL sbc_ice_if ( kt ) ! Ice-cover climatology ("Ice-if" model) 241 241 ! 242 CASE( 2 ) ; CALL sbc_ice_lim_2( kt, nsbc ) ! LIM 2.0ice model243 IF( lk_bdy ) CALL bdy_ice_frs( kt )242 CASE( 2 ) ; CALL sbc_ice_lim_2( kt, nsbc ) ! LIM-2 ice model 243 IF( lk_bdy ) CALL bdy_ice_frs ( kt ) ! BDY boundary condition 244 244 ! 245 CASE( 3 ) ; CALL sbc_ice_lim ( kt, nsbc , nn_ico_cpl) ! LIM 3.0ice model245 CASE( 3 ) ; CALL sbc_ice_lim ( kt, nsbc ) ! LIM-3 ice model 246 246 END SELECT 247 247 -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90
r2350 r2370 4 4 !! Ocean forcing: river runoff 5 5 !!===================================================================== 6 !! History : OPA ! 2000-11 (R. Hordoir, E. Durand) NetCDF FORMAT 7 !! NEMO 1.0 ! 2002-09 (G. Madec) F90: Free form and module 8 !! 3.0 ! 2006-07 (G. Madec) Surface module 9 !! 3.2 ! 2009-04 (B. Lemaire) Introduce iom_put 6 !! History : OPA ! 2000-11 (R. Hordoir, E. Durand) NetCDF FORMAT 7 !! NEMO 1.0 ! 2002-09 (G. Madec) F90: Free form and module 8 !! 3.0 ! 2006-07 (G. Madec) Surface module 9 !! 3.2 ! 2009-04 (B. Lemaire) Introduce iom_put 10 !! 3.3 ! 2010-10 (R. Furner, G. Madec) runoff distributed over ocean levels 10 11 !!---------------------------------------------------------------------- 11 12 … … 49 50 REAL(wp), PUBLIC, DIMENSION(jpk) :: rnfmsk_z !: river mouth mask (vert.) 50 51 51 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_rnf !: structure of inputriver runoff (file information, fields read)52 53 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_s_rnf !: structure of inputriver runoff salinity (file information, fields read)54 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_t_rnf !: structure of inputriver runoff temperature (file information, fields read)52 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_rnf !: structure: river runoff (file information, fields read) 53 54 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_s_rnf !: structure: river runoff salinity (file information, fields read) 55 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_t_rnf !: structure: river runoff temperature (file information, fields read) 55 56 56 57 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: h_rnf !: depth of runoff in m 57 58 INTEGER, PUBLIC, DIMENSION(jpi,jpj) :: nk_rnf !: depth of runoff in model levels 58 59 59 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpts) :: rnf_tsc_b, rnf_tsc !: before and now 60 ! !: temp. & sal. content of river runoffs [K.m/s & PSU.m/s] 60 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpts) :: rnf_tsc_b, rnf_tsc !: before and now T & S contents of runoffs [K.m/s & PSU.m/s] 61 61 62 62 !! * Substitutions … … 65 65 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 66 66 !! $Id$ 67 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)67 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 68 68 !!---------------------------------------------------------------------- 69 70 69 CONTAINS 71 70 … … 171 170 CALL iom_rstput( kt, nitrst, numrow, 'rnf_sc_b', rnf_tsc(:,:,jp_sal) ) 172 171 ENDIF 173 174 172 ! 175 173 END SUBROUTINE sbc_rnf 174 176 175 177 176 SUBROUTINE sbc_rnf_div( phdivn ) -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/step.F90
r2348 r2370 26 26 27 27 !!---------------------------------------------------------------------- 28 !! stp : OPA system time-stepping28 !! stp : OPA system time-stepping 29 29 !!---------------------------------------------------------------------- 30 30 USE step_oce ! time stepping definition modules 31 31 #if defined key_top 32 USE trcstp ! passive tracer time-stepping (trc_stp routine)33 #endif 34 USE asminc ! assimilation increments (tra_asm_inc, dyn_asm_inc routines)35 USE stpctl ! time stepping control (stp_ctl routine)36 USE restart ! ocean restart (rst_wri routine)37 USE prtctl ! Print control (prt_ctl routine)32 USE trcstp ! passive tracer time-stepping (trc_stp routine) 33 #endif 34 USE asminc ! assimilation increments (tra_asm_inc, dyn_asm_inc routines) 35 USE stpctl ! time stepping control (stp_ctl routine) 36 USE restart ! ocean restart (rst_wri routine) 37 USE prtctl ! Print control (prt_ctl routine) 38 38 39 39 #if defined key_agrif … … 44 44 PRIVATE 45 45 46 PUBLIC stp 46 PUBLIC stp ! called by opa.F90 47 47 48 48 !! * Substitutions … … 52 52 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 53 53 !! $Id$ 54 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 55 !!---------------------------------------------------------------------- 56 54 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 55 !!---------------------------------------------------------------------- 57 56 CONTAINS 58 57 … … 149 148 IF( ln_zps ) CALL zps_hde( kstp, jpts, tsb, gtsu, gtsv, & ! Partial steps: before horizontal gradient 150 149 & rhd, gru , grv ) ! of t, s, rd at the last ocean level 151 IF ( .NOT. ln_traldf_grif ) THEN152 CALL ldf_slp ( kstp, rhd, rn2b ) ! before slope of the lateral mixing150 IF( ln_traldf_grif ) THEN ! before slope for Griffies operator 151 CALL ldf_slp_grif( kstp ) 153 152 ELSE 154 CALL ldf_slp_grif( kstp ) 155 IF ( ln_dynldf_bilap .OR. ln_dynldf_iso ) CALL ldf_slp( kstp, rhd, rn2b ) 156 ENDIF 153 CALL ldf_slp( kstp, rhd, rn2b ) ! before slope for Madec operator 154 ENDIF 157 155 ENDIF 158 156 #if defined key_traldf_c2d -
branches/nemo_v3_3_beta/NEMOGCM/TOOLS/COMPILE/cfg.txt
r2281 r2370 4 4 ORCA2_OFF_PISCES OPA_SRC OFF_SRC TOP_SRC 5 5 POMME OPA_SRC NST_SRC 6 ORCA2_LIM3 OPA_SRC LIM_SRC_3 6 7 ORCA2_LIM OPA_SRC LIM_SRC_2 NST_SRC C1D_SRC
Note: See TracChangeset
for help on using the changeset viewer.