Changeset 4333 for branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC
- Timestamp:
- 2013-12-11T18:34:00+01:00 (10 years ago)
- Location:
- branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 10 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdy_oce.F90
r4292 r4333 104 104 REAL(wp), DIMENSION(jp_bdy) :: rn_time_dmp_out !: Damping time scale in days at radiation outflow points 105 105 106 #if defined key_lim2107 CHARACTER(len=20), DIMENSION(jp_bdy) :: nn_ice_lim 2! Choice of boundary condition for sea ice variables108 INTEGER, DIMENSION(jp_bdy) :: nn_ice_lim 2_dta!: = 0 use the initial state as bdy dta ;106 #if ( defined key_lim2 || defined key_lim3 ) 107 CHARACTER(len=20), DIMENSION(jp_bdy) :: nn_ice_lim ! Choice of boundary condition for sea ice variables 108 INTEGER, DIMENSION(jp_bdy) :: nn_ice_lim_dta !: = 0 use the initial state as bdy dta ; 109 109 !: = 1 read it in a NetCDF file 110 110 #endif -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90
r4292 r4333 197 197 198 198 #if defined key_lim2 199 IF( nn_ice_lim 2_dta(ib_bdy) .eq. 0 ) THEN199 IF( nn_ice_lim_dta(ib_bdy) .eq. 0 ) THEN 200 200 ilen1(:) = nblen(:) 201 201 IF( dta%ll_frld ) THEN … … 649 649 #if defined key_lim2 650 650 ! sea ice 651 IF( nn_ice_lim 2_dta(ib_bdy) .eq. 1 ) THEN651 IF( nn_ice_lim_dta(ib_bdy) .eq. 1 ) THEN 652 652 653 653 IF( dta%ll_frld ) THEN … … 724 724 ENDIF 725 725 726 ENDIF 726 727 #endif 727 728 ! Recalculate field counts … … 850 851 #if defined key_lim2 851 852 IF (nn_ice_lim(ib_bdy) .gt. 0) THEN 852 IF( nn_ice_lim 2_dta(ib_bdy) .eq. 0 ) THEN853 IF( nn_ice_lim_dta(ib_bdy) .eq. 0 ) THEN 853 854 ALLOCATE( dta_bdy(ib_bdy)%frld(nblen(1)) ) 854 855 ALLOCATE( dta_bdy(ib_bdy)%hicif(nblen(1)) ) -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90
r4292 r4333 45 45 46 46 REAL(wp) :: epsi20 = 1.e-20_wp ! module constants 47 REAL(wp) :: epsi10 = 1.e-10_wp ! min area allowed by ice model 47 48 !!---------------------------------------------------------------------- 48 49 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 91 92 INTEGER, INTENT(in) :: ib_bdy ! BDY set index !! 92 93 94 INTEGER :: jpbound ! 0 = incoming ice 95 ! 1 = outgoing ice 93 96 INTEGER :: jb, jk, jgrd, jl ! dummy loop indices 94 INTEGER :: ji, jj 97 INTEGER :: ji, jj, ii, ij ! local scalar 95 98 REAL(wp) :: zwgt, zwgt1 ! local scalar 96 REAL(wp) :: zinda, ztmelts 99 REAL(wp) :: zinda, ztmelts, zdh 100 101 REAL(wp), PARAMETER :: zsal = 6.3 ! arbitrary salinity for incoming ice 102 REAL(wp), PARAMETER :: ztem = 270.0 ! arbitrary temperature for incoming ice 103 REAL(wp), PARAMETER :: zage = 30.0 ! arbitrary age for incoming ice 97 104 !!------------------------------------------------------------------------------ 98 105 ! 99 106 IF( nn_timing == 1 ) CALL timing_start('bdy_ice_frs') 100 107 ! 101 !IF( kt /= nit000 ) THEN102 108 jgrd = 1 ! Everything is at T-points here 103 109 ! … … 170 176 ht_i(ji,jj,jl) = ( ht_i(ji,jj,jl) * zwgt1 + dta%ht_i(jb,jl) * zwgt ) * tmask(ji,jj,1) ! Ice depth 171 177 ht_s(ji,jj,jl) = ( ht_s(ji,jj,jl) * zwgt1 + dta%ht_s(jb,jl) * zwgt ) * tmask(ji,jj,1) ! Snow depth 178 179 ! ----------------- 180 ! Pathological case 181 ! ----------------- 182 ! In case a) snow load would be in excess or b) ice is coming into a warmer environment that would lead to 183 ! very large transformation from snow to ice (see limthd_dh.F90) 184 185 ! Then, a) transfer the snow excess into the ice (different from limthd_dh) 186 zdh = MAX( 0._wp, ( rhosn * ht_s(ji,jj,jl) + ( rhoic - rau0 ) * ht_i(ji,jj,jl) ) * r1_rau0 ) 187 ! Or, b) transfer all the snow into ice (if incoming ice is likely to melt as it comes into a warmer environment) 188 !zdh = MAX( 0._wp, ht_s(ji,jj,jl) * rhosn / rhoic ) 189 190 ! recompute ht_i, ht_s 191 ht_i(ji,jj,jl) = MIN( hi_max(jl), ht_i(ji,jj,jl) + zdh ) 192 ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic / rhosn ) 193 172 194 ENDDO 173 195 CALL lbc_bdy_lnk( a_i(:,:,jl), 'T', 1., ib_bdy ) 174 196 CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy ) 175 197 CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy ) 176 198 ENDDO 199 ! retrieve at_i 200 at_i(:,:) = 0._wp 201 DO jl = 1, jpl 202 at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 203 END DO 204 205 DO jl = 1, jpl 177 206 DO jb = 1, idx%nblen(jgrd) 178 207 ji = idx%nbi(jb,jgrd) 179 208 jj = idx%nbj(jb,jgrd) 180 ! clem : condition on ice thickness depends on the ice velocity 181 ! if velocity is outward (strictly), then ice thickness, volume... must be equal to adjacent values 182 IF ( u_ice(ji+1,jj ) < 0. .AND. umask(ji-1,jj ,1) == 0. ) THEN 183 ht_i(ji,jj,jl) = ht_i(ji+1,jj ,jl) 184 a_i (ji,jj,jl) = a_i (ji+1,jj ,jl) 185 ht_s(ji,jj,jl) = ht_s(ji+1,jj ,jl) 186 ENDIF 187 IF ( u_ice(ji-1,jj ) > 0. .AND. umask(ji+1,jj ,1) == 0. ) THEN 188 ht_i(ji,jj,jl) = ht_i(ji-1,jj ,jl) 189 a_i (ji,jj,jl) = a_i (ji-1,jj ,jl) 190 ht_s(ji,jj,jl) = ht_s(ji-1,jj ,jl) 191 ENDIF 192 IF ( v_ice(ji ,jj+1) < 0. .AND. vmask(ji ,jj-1,1) == 0. ) THEN 193 ht_i(ji,jj,jl) = ht_i(ji ,jj+1,jl) 194 a_i (ji,jj,jl) = a_i (ji ,jj+1,jl) 195 ht_s(ji,jj,jl) = ht_s(ji ,jj+1,jl) 196 ENDIF 197 IF ( v_ice(ji ,jj-1) > 0. .AND. vmask(ji ,jj+1,1) == 0. ) THEN 198 ht_i(ji,jj,jl) = ht_i(ji ,jj-1,jl) 199 a_i (ji,jj,jl) = a_i (ji ,jj-1,jl) 200 ht_s(ji,jj,jl) = ht_s(ji ,jj-1,jl) 201 ENDIF 202 203 zinda = 1.0 - MAX( 0.0_wp , SIGN ( 1.0_wp , - a_i(ji,jj,jl) ) ) ! 0 if no ice 204 ! -------------------- 209 210 ! condition on ice thickness depends on the ice velocity 211 ! if velocity is outward (strictly), then ice thickness, volume... must be equal to adjacent values 212 jpbound = 0; ii = ji; ij = jj; 213 214 IF ( u_ice(ji+1,jj ) < 0. .AND. umask(ji-1,jj ,1) == 0. ) jpbound = 1; ii = ji+1; ij = jj 215 IF ( u_ice(ji-1,jj ) > 0. .AND. umask(ji+1,jj ,1) == 0. ) jpbound = 1; ii = ji-1; ij = jj 216 IF ( v_ice(ji ,jj+1) < 0. .AND. vmask(ji ,jj-1,1) == 0. ) jpbound = 1; ii = ji ; ij = jj+1 217 IF ( v_ice(ji ,jj-1) > 0. .AND. vmask(ji ,jj+1,1) == 0. ) jpbound = 1; ii = ji ; ij = jj-1 218 219 zinda = 1.0 - MAX( 0.0_wp , SIGN ( 1.0_wp , - at_i(ii,ij) + 0.01 ) ) ! 0 if no ice 220 221 ! concentration and thickness 222 a_i (ji,jj,jl) = a_i (ii,ij,jl) * zinda 223 ht_i(ji,jj,jl) = ht_i(ii,ij,jl) * zinda 224 ht_s(ji,jj,jl) = ht_s(ii,ij,jl) * zinda 225 205 226 ! Ice and snow volumes 206 ! --------------------207 227 v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl) 208 228 v_s(ji,jj,jl) = ht_s(ji,jj,jl) * a_i(ji,jj,jl) 209 229 210 ! ------------- 211 ! Ice salinity 212 !--------------- 213 sm_i(ji,jj,jl) = zinda * bulk_sal + ( 1.0 - zinda ) * s_i_min 230 SELECT CASE( jpbound ) 231 232 CASE( 0 ) ! velocity is inward 233 234 ! Ice salinity, age, temperature 235 sm_i(ji,jj,jl) = zinda * zsal + ( 1.0 - zinda ) * s_i_min 236 o_i(ji,jj,jl) = zinda * zage + ( 1.0 - zinda ) 237 t_su(ji,jj,jl) = zinda * ztem + ( 1.0 - zinda ) * ztem 238 DO jk = 1, nlay_s 239 t_s(ji,jj,jk,jl) = zinda * ztem + ( 1.0 - zinda ) * rtt 240 END DO 241 DO jk = 1, nlay_i 242 t_i(ji,jj,jk,jl) = zinda * ztem + ( 1.0 - zinda ) * rtt 243 s_i(ji,jj,jk,jl) = zinda * zsal + ( 1.0 - zinda ) * s_i_min 244 END DO 245 246 CASE( 1 ) ! velocity is outward 247 248 ! Ice salinity, age, temperature 249 sm_i(ji,jj,jl) = zinda * sm_i(ii,ij,jl) + ( 1.0 - zinda ) * s_i_min 250 o_i(ji,jj,jl) = zinda * o_i(ii,ij,jl) + ( 1.0 - zinda ) 251 t_su(ji,jj,jl) = zinda * t_su(ii,ij,jl) + ( 1.0 - zinda ) * rtt 252 DO jk = 1, nlay_s 253 t_s(ji,jj,jk,jl) = zinda * t_s(ii,ij,jk,jl) + ( 1.0 - zinda ) * rtt 254 END DO 255 DO jk = 1, nlay_i 256 t_i(ji,jj,jk,jl) = zinda * t_i(ii,ij,jk,jl) + ( 1.0 - zinda ) * rtt 257 s_i(ji,jj,jk,jl) = zinda * s_i(ii,ij,jk,jl) + ( 1.0 - zinda ) * s_i_min 258 END DO 259 260 END SELECT 261 262 ! contents 214 263 smv_i(ji,jj,jl) = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) 215 216 !----------217 ! Ice age218 !----------219 o_i(ji,jj,jl) = zinda * 1.0 + ( 1.0 - zinda )220 264 oa_i(ji,jj,jl) = o_i(ji,jj,jl) * a_i(ji,jj,jl) 221 222 !------------------------------223 ! Sea ice surface temperature224 !------------------------------225 t_su(ji,jj,jl) = zinda * 270.0 + ( 1.0 - zinda ) * t_bo(ji,jj)226 227 !------------------------------------228 ! Snow temperature and heat content229 !------------------------------------230 265 DO jk = 1, nlay_s 231 t_s(ji,jj,jk,jl) = zinda * 270.00 + ( 1.0 - zinda ) * rtt232 266 ! Snow energy of melting 233 267 e_s(ji,jj,jk,jl) = zinda * rhosn * ( cpic * ( rtt - t_s(ji,jj,jk,jl) ) + lfus ) … … 235 269 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac 236 270 ! Multiply by volume, so that heat content in 10^9 Joules 237 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * & 238 v_s(ji,jj,jl) / nlay_s 239 END DO !jk 240 241 !----------------------------------------------- 242 ! Ice salinities, temperature and heat content 243 !----------------------------------------------- 271 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * v_s(ji,jj,jl) / nlay_s 272 END DO 244 273 DO jk = 1, nlay_i 245 t_i(ji,jj,jk,jl) = zinda * 270.00 + ( 1.0 - zinda ) * rtt 246 s_i(ji,jj,jk,jl) = zinda * bulk_sal + ( 1.0 - zinda ) * s_i_min 247 ztmelts = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K 248 274 ztmelts = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K 249 275 ! heat content per unit volume 250 276 e_i(ji,jj,jk,jl) = zinda * rhoic * & 251 ( cpic * ( ztmelts - t_i(ji,jj,jk,jl) ) & 252 + lfus * ( 1.0 - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-epsi20) ) & 253 - rcp * ( ztmelts - rtt ) ) 254 277 ( cpic * ( ztmelts - t_i(ji,jj,jk,jl) ) & 278 + lfus * ( 1.0 - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-epsi20) ) & 279 - rcp * ( ztmelts - rtt ) ) 255 280 ! Correct dimensions to avoid big values 256 281 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac 257 258 282 ! Mutliply by ice volume, and divide by number of layers to get heat content in 10^9 J 259 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * &260 area(ji,jj) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) / nlay_i261 END DO ! jk 283 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * area(ji,jj) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) / nlay_i 284 END DO 285 262 286 263 287 END DO !jb 264 288 265 266 289 CALL lbc_bdy_lnk( a_i(:,:,jl), 'T', 1., ib_bdy ) ! lateral boundary conditions 267 290 CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy ) … … 288 311 #endif 289 312 ! 290 !ENDIF !nit000/=0291 313 IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_frs') 292 314 ! … … 294 316 295 317 296 SUBROUTINE bdy_ice_lim_dyn( kn, pvar1, pvar2, pvar12)318 SUBROUTINE bdy_ice_lim_dyn( cd_type ) 297 319 !!------------------------------------------------------------------------------ 298 320 !! *** SUBROUTINE bdy_ice_lim_dyn *** 299 321 !! 300 322 !! ** Purpose : Apply dynamics boundary conditions for sea-ice in the cas of unstructured open boundaries. 301 !! kn = 1:u_ice and v_ice are equal to the value of the adjacent grid point if this latter is not ice free323 !! u_ice and v_ice are equal to the value of the adjacent grid point if this latter is not ice free 302 324 !! if adjacent grid point is ice free, then u_ice and v_ice are equal to ocean velocities 303 !! kn = 2: the stress tensor is set to 0 (i.e. pvar1, pvar2, pvar12)304 325 !! 305 326 !! 2013-06 : C. Rousset 306 327 !!------------------------------------------------------------------------------ 307 328 !! 308 INTEGER, INTENT( in ) :: kn ! set up of ice vel (kn=1) or stress tensor (kn=2) 309 REAL(wp), INTENT( inout ), DIMENSION(:,:), OPTIONAL :: pvar1, pvar2, pvar12 ! stress tensor components (zs1,zs2,zs12) 329 CHARACTER(len=1), INTENT(in) :: cd_type ! nature of velocity grid-points 310 330 INTEGER :: jb, jgrd ! dummy loop indices 311 331 INTEGER :: ji, jj ! local scalar … … 317 337 ! 318 338 DO ib_bdy=1, nb_bdy 319 ! 320 SELECT CASE( nn_ice_lim(ib_bdy) ) 321 322 CASE(jp_none) 323 CYCLE 324 325 CASE(jp_frs) 326 327 IF ( kn == 1 ) THEN ! set up of u_ice and v_ice 328 329 jgrd = 2 ! u velocity 330 DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 331 ji = idx_bdy(ib_bdy)%nbi(jb,jgrd) 332 jj = idx_bdy(ib_bdy)%nbj(jb,jgrd) 333 zflag = idx_bdy(ib_bdy)%flagu(jb) 334 335 IF ( ABS( zflag ) == 1. ) THEN ! eastern and western boundaries 336 ! one of the two zmsk is always 0 (because of zflag) 337 zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji+1,jj) ) ) ! 0 if no ice 338 zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji-1,jj) ) ) ! 0 if no ice 339 ! 340 SELECT CASE( nn_ice_lim(ib_bdy) ) 341 342 CASE('none') 343 344 CYCLE 345 346 CASE('frs') 347 348 349 SELECT CASE ( cd_type ) 350 351 CASE ( 'U' ) 352 353 jgrd = 2 ! u velocity 354 DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 355 ji = idx_bdy(ib_bdy)%nbi(jb,jgrd) 356 jj = idx_bdy(ib_bdy)%nbj(jb,jgrd) 357 zflag = idx_bdy(ib_bdy)%flagu(jb) 339 358 340 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce) 341 u_ice (ji,jj) = u_ice(ji+1,jj) * 0.5 * ABS( zflag + 1._wp ) * zmsk1 + & 342 & u_ice(ji-1,jj) * 0.5 * ABS( zflag - 1._wp ) * zmsk2 + & 343 & u_oce(ji ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 344 ELSE ! everywhere else 345 u_ice(ji,jj) = u_oce(ji,jj) 346 ENDIF 347 ! mask ice velocities 348 zinda = 1.0 - MAX( 0.0_wp, SIGN( 1.0_wp , - at_i(ji,jj) ) ) ! 0 if no ice 349 u_ice(ji,jj) = zinda * u_ice(ji,jj) 350 351 ENDDO 359 IF ( ABS( zflag ) == 1. ) THEN ! eastern and western boundaries 360 ! one of the two zmsk is always 0 (because of zflag) 361 zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji+1,jj) ) ) ! 0 if no ice 362 zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji-1,jj) ) ) ! 0 if no ice 363 364 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce) 365 u_ice (ji,jj) = u_ice(ji+1,jj) * 0.5 * ABS( zflag + 1._wp ) * zmsk1 + & 366 & u_ice(ji-1,jj) * 0.5 * ABS( zflag - 1._wp ) * zmsk2 + & 367 & u_oce(ji ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 368 ELSE ! everywhere else 369 !u_ice(ji,jj) = u_oce(ji,jj) 370 u_ice(ji,jj) = 0._wp 371 ENDIF 372 ! mask ice velocities 373 zinda = 1.0 - MAX( 0.0_wp , SIGN ( 1.0_wp , - at_i(ji,jj) + 0.01 ) ) ! 0 if no ice 374 u_ice(ji,jj) = zinda * u_ice(ji,jj) 375 376 ENDDO 377 378 CALL lbc_bdy_lnk( u_ice(:,:), 'U', -1., ib_bdy ) 379 380 CASE ( 'V' ) 381 382 jgrd = 3 ! v velocity 383 DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 384 ji = idx_bdy(ib_bdy)%nbi(jb,jgrd) 385 jj = idx_bdy(ib_bdy)%nbj(jb,jgrd) 386 zflag = idx_bdy(ib_bdy)%flagv(jb) 387 388 IF ( ABS( zflag ) == 1. ) THEN ! northern and southern boundaries 389 ! one of the two zmsk is always 0 (because of zflag) 390 zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj+1) ) ) ! 0 if no ice 391 zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj-1) ) ) ! 0 if no ice 392 393 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce) 394 v_ice (ji,jj) = v_ice(ji,jj+1) * 0.5 * ABS( zflag + 1._wp ) * zmsk1 + & 395 & v_ice(ji,jj-1) * 0.5 * ABS( zflag - 1._wp ) * zmsk2 + & 396 & v_oce(ji,jj ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 397 ELSE ! everywhere else 398 !v_ice(ji,jj) = v_oce(ji,jj) 399 v_ice(ji,jj) = 0._wp 400 ENDIF 401 ! mask ice velocities 402 zinda = 1.0 - MAX( 0.0_wp , SIGN ( 1.0_wp , - at_i(ji,jj) + 0.01 ) ) ! 0 if no ice 403 v_ice(ji,jj) = zinda * v_ice(ji,jj) 404 405 ENDDO 406 407 CALL lbc_bdy_lnk( v_ice(:,:), 'V', -1., ib_bdy ) 408 409 END SELECT 352 410 353 jgrd = 3 ! v velocity 354 DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd) 355 ji = idx_bdy(ib_bdy)%nbi(jb,jgrd) 356 jj = idx_bdy(ib_bdy)%nbj(jb,jgrd) 357 zflag = idx_bdy(ib_bdy)%flagv(jb) 358 359 IF ( ABS( zflag ) == 1. ) THEN ! northern and southern boundaries 360 ! one of the two zmsk is always 0 (because of zflag) 361 zmsk1 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj+1) ) ) ! 0 if no ice 362 zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj-1) ) ) ! 0 if no ice 363 364 ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce) 365 v_ice (ji,jj) = v_ice(ji,jj+1) * 0.5 * ABS( zflag + 1._wp ) * zmsk1 + & 366 & v_ice(ji,jj-1) * 0.5 * ABS( zflag - 1._wp ) * zmsk2 + & 367 & v_oce(ji,jj ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 368 ELSE ! everywhere else 369 v_ice(ji,jj) = v_oce(ji,jj) 370 ENDIF 371 ! mask ice velocities 372 zinda = 1.0 - MAX( 0.0_wp, SIGN( 1.0_wp , - at_i(ji,jj) ) ) ! 0 if no ice 373 v_ice(ji,jj) = zinda * v_ice(ji,jj) 374 375 ENDDO 376 377 CALL lbc_bdy_lnk( u_ice(:,:), 'U', -1., ib_bdy ) 378 CALL lbc_bdy_lnk( v_ice(:,:), 'V', -1., ib_bdy ) 379 380 ELSEIF ( kn == 2 ) THEN ! set up of stress tensor (not sure it works) 411 CASE DEFAULT 412 CALL ctl_stop( 'bdy_ice_lim_dyn : unrecognised option for open boundaries for ice fields' ) 413 END SELECT 381 414 382 jgrd = 1 ! T grid383 DO jb = 1, idx_bdy(ib_bdy)%nblen(jgrd)384 ji = idx_bdy(ib_bdy)%nbi(jb,jgrd)385 jj = idx_bdy(ib_bdy)%nbj(jb,jgrd)386 387 pvar1 (ji,jj) = 0._wp388 pvar2 (ji,jj) = 0._wp389 pvar12(ji,jj) = 0._wp390 391 ENDDO392 393 ENDIF394 395 CASE DEFAULT396 CALL ctl_stop( 'bdy_ice_lim_dyn : unrecognised option for open boundaries for ice fields' )397 END SELECT398 399 415 ENDDO 400 416 401 417 IF( nn_timing == 1 ) CALL timing_stop('bdy_ice_lim_dyn') 402 418 403 419 END SUBROUTINE bdy_ice_lim_dyn 404 420 -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90
r4294 r4333 488 488 DO igrd = 1, jpbgrd 489 489 id_dummy = iom_varid( inum, 'nbi'//cgrid(igrd), kdimsz=kdimsz ) 490 nblendta(igrd,ib_bdy) = kdimsz(1) 491 jpbdtau = MAX(jpbdtau, kdimsz(1)) 490 !clem nblendta(igrd,ib_bdy) = kdimsz(1) 491 !clem jpbdtau = MAX(jpbdtau, kdimsz(1)) 492 nblendta(igrd,ib_bdy) = MAXVAL(kdimsz) 493 jpbdtau = MAX(jpbdtau, MAXVAL(kdimsz)) 492 494 ENDDO 493 495 CALL iom_close( inum ) -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90
r4328 r4333 36 36 LOGICAL, PUBLIC :: ln_diahsb !: check the heat and salt budgets 37 37 38 REAL( dp), SAVE :: frc_t , frc_s , frc_v ! global forcing trends39 REAL( dp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ssh_ini !40 REAL( dp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hc_loc_ini, sc_loc_ini, e3t_ini !41 REAL( dp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hcssh_loc_ini, scssh_loc_ini !38 REAL(wp), SAVE :: frc_t , frc_s , frc_v ! global forcing trends 39 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ssh_ini ! 40 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hc_loc_ini, sc_loc_ini, e3t_ini ! 41 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hcssh_loc_ini, scssh_loc_ini ! 42 42 43 43 !! * Substitutions … … 67 67 !! 68 68 INTEGER :: jk ! dummy loop indice 69 REAL( dp) :: zdiff_hc , zdiff_sc ! heat and salt content variations70 REAL( dp) :: zdiff_v1 , zdiff_v2 ! volume variation71 REAL( dp) :: z_hc , z_sc ! heat and salt content72 REAL( dp) :: z_v1 , z_v2 ! volume73 REAL( dp) :: zdeltat ! - -74 REAL( dp) :: z_frc_trd_t , z_frc_trd_s ! - -75 REAL( dp) :: z_frc_trd_v ! - -76 REAL( dp), POINTER, DIMENSION(:,:) :: zsurf !69 REAL(wp) :: zdiff_hc , zdiff_sc ! heat and salt content variations 70 REAL(wp) :: zdiff_v1 , zdiff_v2 ! volume variation 71 REAL(wp) :: z_hc , z_sc ! heat and salt content 72 REAL(wp) :: z_v1 , z_v2 ! volume 73 REAL(wp) :: zdeltat ! - - 74 REAL(wp) :: z_frc_trd_t , z_frc_trd_s ! - - 75 REAL(wp) :: z_frc_trd_v ! - - 76 REAL(wp), POINTER, DIMENSION(:,:) :: zsurf ! 77 77 !!--------------------------------------------------------------------------- 78 78 IF( nn_timing == 1 ) CALL timing_start('dia_hsb') … … 104 104 ! 2a - Content variations ! 105 105 ! ------------------------ ! 106 zdiff_v2 = 0._ dp107 zdiff_hc = 0._ dp108 zdiff_sc = 0._ dp106 zdiff_v2 = 0._wp 107 zdiff_hc = 0._wp 108 zdiff_sc = 0._wp 109 109 ! volume variation (calculated with ssh) 110 110 zdiff_v1 = glob_sum( zsurf(:,:) * ( sshn(:,:) - ssh_ini(:,:) ) ) 111 111 DO jk = 1, jpkm1 112 112 ! volume variation (calculated with scale factors) 113 zdiff_v2 = zdiff_v2 & 114 & + glob_sum( zsurf(:,:) * tmask(:,:,jk) * ( fse3t_n(:,:,jk) - e3t_ini(:,:,jk) ) ) 113 zdiff_v2 = zdiff_v2 + glob_sum( zsurf(:,:) * tmask(:,:,jk) * ( fse3t_n(:,:,jk) - e3t_ini(:,:,jk) ) ) 115 114 ! heat content variation 116 zdiff_hc = zdiff_hc & 117 & + glob_sum( zsurf(:,:) * tmask(:,:,jk) * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) & 115 zdiff_hc = zdiff_hc + glob_sum( zsurf(:,:) * tmask(:,:,jk) * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) & 118 116 & - hc_loc_ini(:,:,jk) ) ) 119 117 ! salt content variation 120 zdiff_sc = zdiff_sc & 121 & + glob_sum( zsurf(:,:) * tmask(:,:,jk) * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_sal) & 118 zdiff_sc = zdiff_sc + glob_sum( zsurf(:,:) * tmask(:,:,jk) * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_sal) & 122 119 & - sc_loc_ini(:,:,jk) ) ) 123 120 ENDDO … … 140 137 ! 2b - Content ! 141 138 ! ----------------------- ! 142 z_v2 = 0._ dp143 z_hc = 0._ dp144 z_sc = 0._ dp139 z_v2 = 0._wp 140 z_hc = 0._wp 141 z_sc = 0._wp 145 142 ! volume (calculated with ssh) 146 143 z_v1 = glob_sum( zsurf(:,:) * sshn(:,:) ) 147 144 DO jk = 1, jpkm1 148 145 ! volume (calculated with scale factors) 149 z_v2 = z_v2 & 150 & + glob_sum( zsurf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) ) 146 z_v2 = z_v2 + glob_sum( zsurf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) ) 151 147 ! heat content 152 z_hc = z_hc & 153 & + glob_sum( zsurf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) ) 148 z_hc = z_hc + glob_sum( zsurf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) ) 154 149 ! salt content 155 z_sc = z_sc & 156 & + glob_sum( zsurf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) * tsn(:,:,jk,jp_sal) ) 150 z_sc = z_sc + glob_sum( zsurf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) * tsn(:,:,jk,jp_sal) ) 157 151 ENDDO 158 152 ! add ssh if not vvl … … 170 164 CALL iom_put( 'bgtemper' , z_hc / z_v2 ) ! Temperature (C) 171 165 CALL iom_put( 'bgsaline' , z_sc / z_v2 ) ! Salinity (psu) 172 CALL iom_put( 'bgheatco' , zdiff_hc * rau0 * rcp * 1.e-9_ dp ) ! Heat content variation (10^9 J)166 CALL iom_put( 'bgheatco' , zdiff_hc * rau0 * rcp * 1.e-9_wp ) ! Heat content variation (10^9 J) 173 167 CALL iom_put( 'bgsaltco' , zdiff_sc * 1.e-9 ) ! Salt content variation (psu*km3) 174 168 CALL iom_put( 'bgvolssh' , zdiff_v1 * 1.e-9 ) ! volume ssh (km3) … … 176 170 CALL iom_put( 'bgvoltot' , zdiff_v2 * 1.e-9 ) ! volume total (km3) 177 171 CALL iom_put( 'bgfrcvol' , frc_v * 1.e-9 ) ! vol - surface forcing (volume) 178 CALL iom_put( 'bgfrctem' , frc_t * rau0 * rcp * 1.e-9_ dp ) ! hc - surface forcing (heat content)172 CALL iom_put( 'bgfrctem' , frc_t * rau0 * rcp * 1.e-9_wp ) ! hc - surface forcing (heat content) 179 173 CALL iom_put( 'bgfrcsal' , frc_s * 1.e-9 ) ! sc - surface forcing (salt content) 180 174 ! … … 224 218 IF(lwp) THEN ! Control print 225 219 WRITE(numout,*) 220 WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 221 WRITE(numout,*) '~~~~~~~~~~~~' 226 222 WRITE(numout,*) ' Namelist namhsb : set hsb parameters' 227 223 WRITE(numout,*) ' Switch for hsb diagnostic (T) or not (F) ln_diahsb = ', ln_diahsb … … 308 304 hcssh_loc_ini(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) ! initial heat content in ssh 309 305 scssh_loc_ini(:,:) = tsn(:,:,1,jp_sal) * sshn(:,:) ! initial salt content in ssh 310 frc_v = 0._ dp311 frc_t = 0._ dp312 frc_s = 0._ dp306 frc_v = 0._wp 307 frc_t = 0._wp 308 frc_s = 0._wp 313 309 ENDIF 314 310 -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/IOM/restart.F90
r4313 r4333 23 23 USE eosbn2 ! equation of state (eos bn2 routine) 24 24 USE trdmld_oce ! ocean active mixed layer tracers trends variables 25 USE domvvl ! variable volume 25 26 USE divcur ! hor. divergence and curl (div & cur routines) 26 27 USE sbc_ice, ONLY : lk_lim3 … … 210 211 CALL iom_get( numror, jpdom_autoglo, 'hdivb' , hdivb ) 211 212 CALL iom_get( numror, jpdom_autoglo, 'sshb' , sshb ) 213 IF( lk_vvl ) CALL iom_get( numror, jpdom_autoglo, 'fse3t_b', fse3t_b(:,:,:) ) 212 214 ELSE 213 215 neuler = 0 -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r4306 r4333 793 793 Cd_n10(:,:) = cdn_wave 794 794 ELSE 795 Cd_n10 = 1 E-3 * ( 2.7/dU10 + 0.142 + dU10/13.09 ) ! L & Y eq. (6a)795 Cd_n10 = 1.e-3 * ( 2.7/dU10 + 0.142 + dU10/13.09 ) ! L & Y eq. (6a) 796 796 ENDIF 797 797 sqrt_Cd_n10 = sqrt(Cd_n10) 798 Ce_n10 = 1 E-3 * ( 34.6 * sqrt_Cd_n10 ) ! L & Y eq. (6b)799 Ch_n10 = 1 E-3*sqrt_Cd_n10*(18*stab + 32.7*(1-stab)) ! L & Y eq. (6c), (6d)798 Ce_n10 = 1.e-3 * ( 34.6 * sqrt_Cd_n10 ) ! L & Y eq. (6b) 799 Ch_n10 = 1.e-3*sqrt_Cd_n10*(18.*stab + 32.7*(1.-stab)) ! L & Y eq. (6c), (6d) 800 800 !! 801 801 !! Initializing transfert coefficients with their first guess neutral equivalents : … … 824 824 825 825 !! Updating the neutral 10m transfer coefficients : 826 Cd_n10 = 1 E-3 * (2.7/U_n10 + 0.142 + U_n10/13.09) ! L & Y eq. (6a)826 Cd_n10 = 1.e-3 * (2.7/U_n10 + 0.142 + U_n10/13.09) ! L & Y eq. (6a) 827 827 sqrt_Cd_n10 = sqrt(Cd_n10) 828 Ce_n10 = 1 E-3 * (34.6 * sqrt_Cd_n10) ! L & Y eq. (6b)828 Ce_n10 = 1.e-3 * (34.6 * sqrt_Cd_n10) ! L & Y eq. (6b) 829 829 stab = 0.5 + sign(0.5,zeta) 830 Ch_n10 = 1 E-3*sqrt_Cd_n10*(18.*stab + 32.7*(1-stab)) ! L & Y eq. (6c), (6d)830 Ch_n10 = 1.e-3*sqrt_Cd_n10*(18.*stab + 32.7*(1.-stab)) ! L & Y eq. (6c), (6d) 831 831 832 832 !! Shifting the neutral 10m transfer coefficients to ( zu , zeta ) : … … 927 927 Cd_n10(:,:) = cdn_wave 928 928 ELSE 929 Cd_n10 = 1 E-3*( 2.7/dU10 + 0.142 + dU10/13.09 )929 Cd_n10 = 1.e-3*( 2.7/dU10 + 0.142 + dU10/13.09 ) 930 930 ENDIF 931 931 sqrt_Cd_n10 = sqrt(Cd_n10) 932 Ce_n10 = 1 E-3*( 34.6 * sqrt_Cd_n10 )933 Ch_n10 = 1 E-3*sqrt_Cd_n10*(18*stab + 32.7*(1- stab))932 Ce_n10 = 1.e-3*( 34.6 * sqrt_Cd_n10 ) 933 Ch_n10 = 1.e-3*sqrt_Cd_n10*(18.*stab + 32.7*(1. - stab)) 934 934 935 935 !! Initializing transf. coeff. with their first guess neutral equivalents : … … 973 973 ELSE 974 974 !! Updating the neutral 10m transfer coefficients : 975 Cd_n10 = 1 E-3 * (2.7/U_n10 + 0.142 + U_n10/13.09) ! L & Y eq. (6a)975 Cd_n10 = 1.e-3 * (2.7/U_n10 + 0.142 + U_n10/13.09) ! L & Y eq. (6a) 976 976 sqrt_Cd_n10 = sqrt(Cd_n10) 977 Ce_n10 = 1 E-3 * (34.6 * sqrt_Cd_n10) ! L & Y eq. (6b)977 Ce_n10 = 1.e-3 * (34.6 * sqrt_Cd_n10) ! L & Y eq. (6b) 978 978 stab = 0.5 + sign(0.5,zeta_u) 979 Ch_n10 = 1 E-3*sqrt_Cd_n10*(18.*stab + 32.7*(1-stab)) ! L & Y eq. (6c-6d)979 Ch_n10 = 1.e-3*sqrt_Cd_n10*(18.*stab + 32.7*(1.-stab)) ! L & Y eq. (6c-6d) 980 980 !! 981 981 !! -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r4269 r4333 58 58 USE in_out_manager ! I/O manager 59 59 USE prtctl ! Print control 60 USE lib_fortran ! 60 61 61 62 #if defined key_bdy … … 167 168 ! 168 169 IF( ln_nicep ) THEN ! control print at a given point 169 jiindx = 6 ; jjindx = 47170 WRITE(numout,*) ' The debugging point is : jiindx : ',jiindx, ' jjindx : ',jjindx170 jiindx = 177 ; jjindx = 112 171 IF(lwp) WRITE(numout,*) ' The debugging point is : jiindx : ',jiindx, ' jjindx : ',jjindx 171 172 ENDIF 172 173 ENDIF … … 303 304 fhbri (:,:) = 0._wp ; fheat_mec(:,:) = 0._wp ; fheat_res(:,:) = 0._wp 304 305 fhmec (:,:) = 0._wp ; 305 fmmec (:,:) = 0._wp 306 fmmec (:,:) = 0._wp 307 fmmflx (:,:) = 0._wp 306 308 focea2D(:,:) = 0._wp 307 309 fsup2D (:,:) = 0._wp … … 328 330 CALL lim_rst_opn( kt ) ! Open Ice restart file 329 331 ! 330 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx, 1, ' - Beginning the time step - ' ) ! control print 331 ! 332 IF( .NOT. lk_c1d ) THEN ! Ice dynamics & transport (except in 1D case) 332 IF( ln_nicep ) CALL lim_prt_state( kt, jiindx, jjindx, 1, ' - Beginning the time step - ' ) ! control print 333 ! ---------------------------------------------- 334 ! ice dynamics and transport (except in 1D case) 335 ! ---------------------------------------------- 336 IF( .NOT. lk_c1d ) THEN 333 337 CALL lim_dyn( kt ) ! Ice dynamics ( rheology/dynamics ) 334 338 CALL lim_trp( kt ) ! Ice transport ( Advection/diffusion ) 335 339 CALL lim_var_glo2eqv ! equivalent variables, requested for rafting 336 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx,-1, ' - ice dyn & trp - ' ) ! control print340 IF( ln_nicep ) CALL lim_prt_state( kt, jiindx, jjindx,-1, ' - ice dyn & trp - ' ) ! control print 337 341 CALL lim_itd_me ! Mechanical redistribution ! (ridging/rafting) 338 342 CALL lim_var_agg( 1 ) 343 #if defined key_bdy 344 ! bdy ice thermo 345 CALL lim_var_glo2eqv ! equivalent variables 346 CALL bdy_ice_lim( kt ) 347 CALL lim_itd_me_zapsmall 348 CALL lim_var_agg(1) 349 IF( ln_nicep ) CALL lim_prt_state( kt, jiindx, jjindx, 1, ' - ice thermo bdy - ' ) ! control print 350 #endif 339 351 CALL lim_update1 340 352 ENDIF … … 349 361 old_oa_i(:,:,:) = oa_i(:,:,:) 350 362 old_smv_i(:,:,:) = smv_i (:,:,:) 351 ! ! Ice thermodynamics 363 364 ! ---------------------------------------------- 365 ! ice thermodynamic 366 ! ---------------------------------------------- 352 367 CALL lim_var_glo2eqv ! equivalent variables 353 368 CALL lim_var_agg(1) ! aggregate ice categories 369 ! previous lead fraction and ice volume for flux calculations 370 pfrld(:,:) = 1._wp - at_i(:,:) 371 phicif(:,:) = vt_i(:,:) 372 ! 354 373 CALL lim_var_bv ! bulk brine volume (diag) 355 374 CALL lim_thd( kt ) ! Ice thermodynamics … … 357 376 oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * zcoef 358 377 CALL lim_var_glo2eqv ! this CALL is maybe not necessary (Martin) 359 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx, 1, ' - ice thermodyn. - ' ) ! control print378 IF( ln_nicep ) CALL lim_prt_state( kt, jiindx, jjindx, 1, ' - ice thermodyn. - ' ) ! control print 360 379 CALL lim_itd_th( kt ) ! Remap ice categories, lateral accretion ! 361 !362 ! ! Global variables update363 380 CALL lim_var_agg( 1 ) ! requested by limupdate 364 CALL lim_update2 ! Global variables update 365 #if defined key_bdy 366 CALL lim_var_glo2eqv ! 367 CALL bdy_ice_lim( kt ) ! clem modif: bdy ice 368 CALL bdy_ice_lim_dyn( 1 ) !?? repeated from limrhg.F90 369 #endif 381 CALL lim_update2 ! Global variables update 382 370 383 CALL lim_var_glo2eqv ! equivalent variables (outputs) 371 384 CALL lim_var_agg(2) ! aggregate ice thickness categories 372 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx, 2, ' - Final state - ' ) ! control print385 IF( ln_nicep ) CALL lim_prt_state( kt, jiindx, jjindx, 2, ' - Final state - ' ) ! control print 373 386 ! 374 387 CALL lim_sbc_flx( kt ) ! Update surface ocean mass, heat and salt fluxes 375 388 ! 376 IF( ln_nicep ) CALL lim_prt_state( jiindx, jjindx, 3, ' - Final state lim_sbc - ' ) ! control print389 IF( ln_nicep ) CALL lim_prt_state( kt, jiindx, jjindx, 3, ' - Final state lim_sbc - ' ) ! control print 377 390 ! 378 391 ! ! Diagnostics and outputs … … 387 400 CALL lim_var_glo2eqv ! ??? 388 401 ! 389 IF( ln_nicep ) CALL lim_ctl 402 IF( ln_nicep ) CALL lim_ctl( kt ) ! alerts in case of model crash 390 403 ! 391 404 ENDIF ! End sea-ice time step only … … 414 427 415 428 416 SUBROUTINE lim_ctl 429 SUBROUTINE lim_ctl( kt ) 417 430 !!----------------------------------------------------------------------- 418 431 !! *** ROUTINE lim_ctl *** … … 420 433 !! ** Purpose : Alerts in case of model crash 421 434 !!------------------------------------------------------------------- 435 INTEGER, INTENT(in) :: kt ! ocean time step 422 436 INTEGER :: ji, jj, jk, jl ! dummy loop indices 423 437 INTEGER :: inb_altests ! number of alert tests (max 20) … … 439 453 DO ji = 1, jpi 440 454 IF( v_i(ji,jj,jl) /= 0._wp .AND. a_i(ji,jj,jl) == 0._wp ) THEN 441 WRITE(numout,*) ' ALERTE 2 : Incompatible volume and concentration '442 WRITE(numout,*) ' at_i ', at_i(ji,jj)443 WRITE(numout,*) ' Point - category', ji, jj, jl444 WRITE(numout,*) ' a_i *** a_i_old ', a_i (ji,jj,jl), old_a_i (ji,jj,jl)445 WRITE(numout,*) ' v_i *** v_i_old ', v_i (ji,jj,jl), old_v_i (ji,jj,jl)446 WRITE(numout,*) ' d_a_i_thd/trp ', d_a_i_thd(ji,jj,jl), d_a_i_trp(ji,jj,jl)447 WRITE(numout,*) ' d_v_i_thd/trp ', d_v_i_thd(ji,jj,jl), d_v_i_trp(ji,jj,jl)455 !WRITE(numout,*) ' ALERTE 2 : Incompatible volume and concentration ' 456 !WRITE(numout,*) ' at_i ', at_i(ji,jj) 457 !WRITE(numout,*) ' Point - category', ji, jj, jl 458 !WRITE(numout,*) ' a_i *** a_i_old ', a_i (ji,jj,jl), old_a_i (ji,jj,jl) 459 !WRITE(numout,*) ' v_i *** v_i_old ', v_i (ji,jj,jl), old_v_i (ji,jj,jl) 460 !WRITE(numout,*) ' d_a_i_thd/trp ', d_a_i_thd(ji,jj,jl), d_a_i_trp(ji,jj,jl) 461 !WRITE(numout,*) ' d_v_i_thd/trp ', d_v_i_thd(ji,jj,jl), d_v_i_trp(ji,jj,jl) 448 462 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 449 463 ENDIF … … 459 473 DO ji = 1, jpi 460 474 IF( ht_i(ji,jj,jl) > 50._wp ) THEN 461 CALL lim_prt_state(ji, jj, 2, ' ALERTE 3 : Very thick ice ' )475 !CALL lim_prt_state( kt, ji, jj, 2, ' ALERTE 3 : Very thick ice ' ) 462 476 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 463 477 ENDIF … … 470 484 DO jj = 1, jpj 471 485 DO ji = 1, jpi 472 IF( MAX( ABS( u_ice(ji,jj) ), ABS( v_ice(ji,jj) ) ) > 0.5 .AND. &486 IF( MAX( ABS( u_ice(ji,jj) ), ABS( v_ice(ji,jj) ) ) > 1.5 .AND. & 473 487 & at_i(ji,jj) > 0._wp ) THEN 474 CALL lim_prt_state(ji, jj, 1, ' ALERTE 4 : Very fast ice ' )475 WRITE(numout,*) ' ice strength : ', strength(ji,jj)476 WRITE(numout,*) ' oceanic stress utau : ', utau(ji,jj)477 WRITE(numout,*) ' oceanic stress vtau : ', vtau(ji,jj)478 WRITE(numout,*) ' sea-ice stress utau_ice : ', utau_ice(ji,jj)479 WRITE(numout,*) ' sea-ice stress vtau_ice : ', vtau_ice(ji,jj)480 WRITE(numout,*) ' oceanic speed u : ', u_oce(ji,jj)481 WRITE(numout,*) ' oceanic speed v : ', v_oce(ji,jj)482 WRITE(numout,*) ' sst : ', sst_m(ji,jj)483 WRITE(numout,*) ' sss : ', sss_m(ji,jj)484 WRITE(numout,*)488 !CALL lim_prt_state( kt, ji, jj, 1, ' ALERTE 4 : Very fast ice ' ) 489 !WRITE(numout,*) ' ice strength : ', strength(ji,jj) 490 !WRITE(numout,*) ' oceanic stress utau : ', utau(ji,jj) 491 !WRITE(numout,*) ' oceanic stress vtau : ', vtau(ji,jj) 492 !WRITE(numout,*) ' sea-ice stress utau_ice : ', utau_ice(ji,jj) 493 !WRITE(numout,*) ' sea-ice stress vtau_ice : ', vtau_ice(ji,jj) 494 !WRITE(numout,*) ' oceanic speed u : ', u_oce(ji,jj) 495 !WRITE(numout,*) ' oceanic speed v : ', v_oce(ji,jj) 496 !WRITE(numout,*) ' sst : ', sst_m(ji,jj) 497 !WRITE(numout,*) ' sss : ', sss_m(ji,jj) 498 !WRITE(numout,*) 485 499 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 486 500 ENDIF … … 494 508 DO ji = 1, jpi 495 509 IF( tms(ji,jj) <= 0._wp .AND. at_i(ji,jj) > 0._wp ) THEN 496 CALL lim_prt_state(ji, jj, 1, ' ALERTE 6 : Ice on continents ' )497 WRITE(numout,*) ' masks s, u, v : ', tms(ji,jj), tmu(ji,jj), tmv(ji,jj)498 WRITE(numout,*) ' sst : ', sst_m(ji,jj)499 WRITE(numout,*) ' sss : ', sss_m(ji,jj)500 WRITE(numout,*) ' at_i(ji,jj) : ', at_i(ji,jj)501 WRITE(numout,*) ' v_ice(ji,jj) : ', v_ice(ji,jj)502 WRITE(numout,*) ' v_ice(ji,jj-1) : ', v_ice(ji,jj-1)503 WRITE(numout,*) ' u_ice(ji-1,jj) : ', u_ice(ji-1,jj)504 WRITE(numout,*) ' u_ice(ji,jj) : ', v_ice(ji,jj)510 !CALL lim_prt_state( kt, ji, jj, 1, ' ALERTE 6 : Ice on continents ' ) 511 !WRITE(numout,*) ' masks s, u, v : ', tms(ji,jj), tmu(ji,jj), tmv(ji,jj) 512 !WRITE(numout,*) ' sst : ', sst_m(ji,jj) 513 !WRITE(numout,*) ' sss : ', sss_m(ji,jj) 514 !WRITE(numout,*) ' at_i(ji,jj) : ', at_i(ji,jj) 515 !WRITE(numout,*) ' v_ice(ji,jj) : ', v_ice(ji,jj) 516 !WRITE(numout,*) ' v_ice(ji,jj-1) : ', v_ice(ji,jj-1) 517 !WRITE(numout,*) ' u_ice(ji-1,jj) : ', u_ice(ji-1,jj) 518 !WRITE(numout,*) ' u_ice(ji,jj) : ', v_ice(ji,jj) 505 519 ! 506 520 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 … … 516 530 DO jj = 1, jpj 517 531 DO ji = 1, jpi 518 !!gm test twice sm_i ... ???? bug? 519 IF( ( ( ABS( sm_i(ji,jj,jl) ) < 0.5 ) .OR. & 520 ( ABS( sm_i(ji,jj,jl) ) < 0.5 ) ) .AND. & 521 ( a_i(ji,jj,jl) > 0._wp ) ) THEN 522 ! CALL lim_prt_state(ji,jj,1, ' ALERTE 7 : Very fresh ice ' ) 532 IF( sm_i(ji,jj,jl) < 0.1 .AND. a_i(ji,jj,jl) > 0._wp ) THEN 533 ! CALL lim_prt_state(kt,ji,jj,1, ' ALERTE 7 : Very fresh ice ' ) 523 534 ! WRITE(numout,*) ' sst : ', sst_m(ji,jj) 524 535 ! WRITE(numout,*) ' sss : ', sss_m(ji,jj) … … 541 552 ( ABS( o_i(ji,jj,jl) ) < 0._wp) ) .AND. & 542 553 ( a_i(ji,jj,jl) > 0._wp ) ) THEN 543 CALL lim_prt_state(ji, jj, 1, ' ALERTE 9 : Wrong ice age ')554 !CALL lim_prt_state( kt, ji, jj, 1, ' ALERTE 9 : Wrong ice age ') 544 555 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 545 556 ENDIF … … 553 564 DO jj = 1, jpj 554 565 DO ji = 1, jpi 555 IF( ABS( sfx (ji,jj) ) .GT. 1.0e-2 ) THEN 556 CALL lim_prt_state(ji, jj, 3, ' ALERTE 5 : High salt flux ' )557 DO jl = 1, jpl558 WRITE(numout,*) ' Category no: ', jl559 WRITE(numout,*) ' a_i : ', a_i (ji,jj,jl) , ' old_a_i : ', old_a_i (ji,jj,jl)560 WRITE(numout,*) ' d_a_i_trp : ', d_a_i_trp(ji,jj,jl) , ' d_a_i_thd : ', d_a_i_thd(ji,jj,jl)561 WRITE(numout,*) ' v_i : ', v_i (ji,jj,jl) , ' old_v_i : ', old_v_i (ji,jj,jl)562 WRITE(numout,*) ' d_v_i_trp : ', d_v_i_trp(ji,jj,jl) , ' d_v_i_thd : ', d_v_i_thd(ji,jj,jl)563 WRITE(numout,*) ' '564 END DO566 IF( ABS( sfx (ji,jj) ) .GT. 1.0e-2 ) THEN ! = 1 psu/day for 1m ocean depth 567 !CALL lim_prt_state( kt, ji, jj, 3, ' ALERTE 5 : High salt flux ' ) 568 !DO jl = 1, jpl 569 !WRITE(numout,*) ' Category no: ', jl 570 !WRITE(numout,*) ' a_i : ', a_i (ji,jj,jl) , ' old_a_i : ', old_a_i (ji,jj,jl) 571 !WRITE(numout,*) ' d_a_i_trp : ', d_a_i_trp(ji,jj,jl) , ' d_a_i_thd : ', d_a_i_thd(ji,jj,jl) 572 !WRITE(numout,*) ' v_i : ', v_i (ji,jj,jl) , ' old_v_i : ', old_v_i (ji,jj,jl) 573 !WRITE(numout,*) ' d_v_i_trp : ', d_v_i_trp(ji,jj,jl) , ' d_v_i_thd : ', d_v_i_thd(ji,jj,jl) 574 !WRITE(numout,*) ' ' 575 !END DO 565 576 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 566 577 ENDIF … … 573 584 DO jj = 1, jpj 574 585 DO ji = 1, jpi 575 IF( ABS( qns(ji,jj) ) > 1500._wp .AND. at_i(ji,jj) > 0._wp )THEN586 IF( ABS( qns(ji,jj) ) > 1500._wp .AND. at_i(ji,jj) > 0._wp ) THEN 576 587 ! 577 WRITE(numout,*) ' ALERTE 8 : Very high non-solar heat flux'578 WRITE(numout,*) ' ji, jj : ', ji, jj579 WRITE(numout,*) ' qns : ', qns(ji,jj)580 WRITE(numout,*) ' sst : ', sst_m(ji,jj)581 WRITE(numout,*) ' sss : ', sss_m(ji,jj)582 WRITE(numout,*) ' qcmif : ', qcmif(ji,jj)583 WRITE(numout,*) ' qldif : ', qldif(ji,jj)584 WRITE(numout,*) ' qcmif : ', qcmif(ji,jj) / rdt_ice585 WRITE(numout,*) ' qldif : ', qldif(ji,jj) / rdt_ice586 WRITE(numout,*) ' qfvbq : ', qfvbq(ji,jj)587 WRITE(numout,*) ' qdtcn : ', qdtcn(ji,jj)588 WRITE(numout,*) ' qfvbq / dt: ', qfvbq(ji,jj) / rdt_ice589 WRITE(numout,*) ' qdtcn / dt: ', qdtcn(ji,jj) / rdt_ice590 WRITE(numout,*) ' fdtcn : ', fdtcn(ji,jj)591 WRITE(numout,*) ' fhmec : ', fhmec(ji,jj)592 WRITE(numout,*) ' fheat_mec : ', fheat_mec(ji,jj)593 WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj)594 WRITE(numout,*) ' fhbri : ', fhbri(ji,jj)588 !WRITE(numout,*) ' ALERTE 8 : Very high non-solar heat flux' 589 !WRITE(numout,*) ' ji, jj : ', ji, jj 590 !WRITE(numout,*) ' qns : ', qns(ji,jj) 591 !WRITE(numout,*) ' sst : ', sst_m(ji,jj) 592 !WRITE(numout,*) ' sss : ', sss_m(ji,jj) 593 !WRITE(numout,*) ' qcmif : ', qcmif(ji,jj) 594 !WRITE(numout,*) ' qldif : ', qldif(ji,jj) 595 !WRITE(numout,*) ' qcmif : ', qcmif(ji,jj) / rdt_ice 596 !WRITE(numout,*) ' qldif : ', qldif(ji,jj) / rdt_ice 597 !WRITE(numout,*) ' qfvbq : ', qfvbq(ji,jj) 598 !WRITE(numout,*) ' qdtcn : ', qdtcn(ji,jj) 599 !WRITE(numout,*) ' qfvbq / dt: ', qfvbq(ji,jj) / rdt_ice 600 !WRITE(numout,*) ' qdtcn / dt: ', qdtcn(ji,jj) / rdt_ice 601 !WRITE(numout,*) ' fdtcn : ', fdtcn(ji,jj) 602 !WRITE(numout,*) ' fhmec : ', fhmec(ji,jj) 603 !WRITE(numout,*) ' fheat_mec : ', fheat_mec(ji,jj) 604 !WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj) 605 !WRITE(numout,*) ' fhbri : ', fhbri(ji,jj) 595 606 ! 596 CALL lim_prt_state(ji, jj, 2, ' ')607 !CALL lim_prt_state( kt, ji, jj, 2, ' ') 597 608 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 598 609 ! … … 611 622 DO ji = 1, jpi 612 623 ztmelts = -tmut * s_i(ji,jj,jk,jl) + rtt 613 IF( t_i(ji,jj,jk,jl) >= ztmelts .AND. v_i(ji,jj,jl) > 1.e- 6&624 IF( t_i(ji,jj,jk,jl) >= ztmelts .AND. v_i(ji,jj,jl) > 1.e-10 & 614 625 & .AND. a_i(ji,jj,jl) > 0._wp ) THEN 615 WRITE(numout,*) ' ALERTE 10 : Very warm ice'616 WRITE(numout,*) ' ji, jj, jk, jl : ', ji, jj, jk, jl617 WRITE(numout,*) ' t_i : ', t_i(ji,jj,jk,jl)618 WRITE(numout,*) ' e_i : ', e_i(ji,jj,jk,jl)619 WRITE(numout,*) ' s_i : ', s_i(ji,jj,jk,jl)620 WRITE(numout,*) ' ztmelts : ', ztmelts626 !WRITE(numout,*) ' ALERTE 10 : Very warm ice' 627 !WRITE(numout,*) ' ji, jj, jk, jl : ', ji, jj, jk, jl 628 !WRITE(numout,*) ' t_i : ', t_i(ji,jj,jk,jl) 629 !WRITE(numout,*) ' e_i : ', e_i(ji,jj,jk,jl) 630 !WRITE(numout,*) ' s_i : ', s_i(ji,jj,jk,jl) 631 !WRITE(numout,*) ' ztmelts : ', ztmelts 621 632 inb_alp(ialert_id) = inb_alp(ialert_id) + 1 622 633 ENDIF … … 626 637 END DO 627 638 628 ialert_id = 1 ! reference number of this alert 629 cl_alname(ialert_id) = ' NO alerte 1 ' ! name of the alert 630 WRITE(numout,*) 631 WRITE(numout,*) ' All alerts at the end of ice model ' 632 DO ialert_id = 1, inb_altests 633 WRITE(numout,*) ialert_id, cl_alname(ialert_id)//' : ', inb_alp(ialert_id), ' times ! ' 634 END DO 635 ! 639 ! sum of the alerts on all processors 640 IF( lk_mpp ) THEN 641 DO ialert_id = 1, inb_altests 642 CALL mpp_sum(inb_alp(ialert_id)) 643 END DO 644 ENDIF 645 646 ! print alerts 647 IF( lwp ) THEN 648 ialert_id = 1 ! reference number of this alert 649 cl_alname(ialert_id) = ' NO alerte 1 ' ! name of the alert 650 WRITE(numout,*) ' time step ',kt 651 WRITE(numout,*) ' All alerts at the end of ice model ' 652 DO ialert_id = 1, inb_altests 653 WRITE(numout,*) ialert_id, cl_alname(ialert_id)//' : ', inb_alp(ialert_id), ' times ! ' 654 END DO 655 ENDIF 656 ! 636 657 END SUBROUTINE lim_ctl 637 658 638 659 639 SUBROUTINE lim_prt_state( k i, kj, kn, cd1 )660 SUBROUTINE lim_prt_state( kt, ki, kj, kn, cd1 ) 640 661 !!----------------------------------------------------------------------- 641 662 !! *** ROUTINE lim_prt_state *** … … 651 672 !! n : number of the option 652 673 !!------------------------------------------------------------------- 674 INTEGER , INTENT(in) :: kt ! ocean time step 653 675 INTEGER , INTENT(in) :: ki, kj, kn ! ocean gridpoint indices 654 676 CHARACTER(len=*), INTENT(in) :: cd1 ! 655 677 !! 656 INTEGER :: jl 678 INTEGER :: jl, ji, jj 657 679 !!------------------------------------------------------------------- 658 680 659 WRITE(numout,*) cd1 ! print title 660 661 !---------------- 662 ! Simple state 663 !---------------- 664 665 IF ( kn == 1 .OR. kn == -1 ) THEN 666 WRITE(numout,*) ' lim_prt_state - Point : ',ki,kj 667 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 668 WRITE(numout,*) ' Simple state ' 669 WRITE(numout,*) ' masks s,u,v : ', tms(ki,kj), tmu(ki,kj), tmv(ki,kj) 670 WRITE(numout,*) ' lat - long : ', gphit(ki,kj), glamt(ki,kj) 671 WRITE(numout,*) ' Time step : ', numit 672 WRITE(numout,*) ' - Ice drift ' 673 WRITE(numout,*) ' ~~~~~~~~~~~ ' 674 WRITE(numout,*) ' u_ice(i-1,j) : ', u_ice(ki-1,kj) 675 WRITE(numout,*) ' u_ice(i ,j) : ', u_ice(ki,kj) 676 WRITE(numout,*) ' v_ice(i ,j-1): ', v_ice(ki,kj-1) 677 WRITE(numout,*) ' v_ice(i ,j) : ', v_ice(ki,kj) 678 WRITE(numout,*) ' strength : ', strength(ki,kj) 679 WRITE(numout,*) 680 WRITE(numout,*) ' - Cell values ' 681 WRITE(numout,*) ' ~~~~~~~~~~~ ' 682 WRITE(numout,*) ' cell area : ', area(ki,kj) 683 WRITE(numout,*) ' at_i : ', at_i(ki,kj) 684 WRITE(numout,*) ' vt_i : ', vt_i(ki,kj) 685 WRITE(numout,*) ' vt_s : ', vt_s(ki,kj) 686 DO jl = 1, jpl 687 WRITE(numout,*) ' - Category (', jl,')' 688 WRITE(numout,*) ' a_i : ', a_i(ki,kj,jl) 689 WRITE(numout,*) ' ht_i : ', ht_i(ki,kj,jl) 690 WRITE(numout,*) ' ht_s : ', ht_s(ki,kj,jl) 691 WRITE(numout,*) ' v_i : ', v_i(ki,kj,jl) 692 WRITE(numout,*) ' v_s : ', v_s(ki,kj,jl) 693 WRITE(numout,*) ' e_s : ', e_s(ki,kj,1,jl)/1.0e9 694 WRITE(numout,*) ' e_i : ', e_i(ki,kj,1:nlay_i,jl)/1.0e9 695 WRITE(numout,*) ' t_su : ', t_su(ki,kj,jl) 696 WRITE(numout,*) ' t_snow : ', t_s(ki,kj,1,jl) 697 WRITE(numout,*) ' t_i : ', t_i(ki,kj,1:nlay_i,jl) 698 WRITE(numout,*) ' sm_i : ', sm_i(ki,kj,jl) 699 WRITE(numout,*) ' smv_i : ', smv_i(ki,kj,jl) 700 WRITE(numout,*) 701 WRITE(numout,*) ' Pathological case : ', patho_case(ki,kj,jl) 702 END DO 703 ENDIF 704 IF( kn == -1 ) THEN 705 WRITE(numout,*) ' Mechanical Check ************** ' 706 WRITE(numout,*) ' Check what means ice divergence ' 707 WRITE(numout,*) ' Total ice concentration ', at_i (ki,kj) 708 WRITE(numout,*) ' Total lead fraction ', ato_i(ki,kj) 709 WRITE(numout,*) ' Sum of both ', ato_i(ki,kj) + at_i(ki,kj) 710 WRITE(numout,*) ' Sum of both minus 1 ', ato_i(ki,kj) + at_i(ki,kj) - 1.00 711 ENDIF 712 713 714 !-------------------- 715 ! Exhaustive state 716 !-------------------- 717 718 IF ( kn .EQ. 2 ) THEN 719 WRITE(numout,*) ' lim_prt_state - Point : ',ki,kj 720 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 721 WRITE(numout,*) ' Exhaustive state ' 722 WRITE(numout,*) ' lat - long ', gphit(ki,kj), glamt(ki,kj) 723 WRITE(numout,*) ' Time step ', numit 724 WRITE(numout,*) 725 WRITE(numout,*) ' - Cell values ' 726 WRITE(numout,*) ' ~~~~~~~~~~~ ' 727 WRITE(numout,*) ' cell area : ', area(ki,kj) 728 WRITE(numout,*) ' at_i : ', at_i(ki,kj) 729 WRITE(numout,*) ' vt_i : ', vt_i(ki,kj) 730 WRITE(numout,*) ' vt_s : ', vt_s(ki,kj) 731 WRITE(numout,*) ' u_ice(i-1,j) : ', u_ice(ki-1,kj) 732 WRITE(numout,*) ' u_ice(i ,j) : ', u_ice(ki,kj) 733 WRITE(numout,*) ' v_ice(i ,j-1): ', v_ice(ki,kj-1) 734 WRITE(numout,*) ' v_ice(i ,j) : ', v_ice(ki,kj) 735 WRITE(numout,*) ' strength : ', strength(ki,kj) 736 WRITE(numout,*) ' d_u_ice_dyn : ', d_u_ice_dyn(ki,kj), ' d_v_ice_dyn : ', d_v_ice_dyn(ki,kj) 737 WRITE(numout,*) ' old_u_ice : ', old_u_ice(ki,kj) , ' old_v_ice : ', old_v_ice(ki,kj) 738 WRITE(numout,*) 739 740 DO jl = 1, jpl 741 WRITE(numout,*) ' - Category (',jl,')' 742 WRITE(numout,*) ' ~~~~~~~~ ' 743 WRITE(numout,*) ' ht_i : ', ht_i(ki,kj,jl) , ' ht_s : ', ht_s(ki,kj,jl) 744 WRITE(numout,*) ' t_i : ', t_i(ki,kj,1:nlay_i,jl) 745 WRITE(numout,*) ' t_su : ', t_su(ki,kj,jl) , ' t_s : ', t_s(ki,kj,1,jl) 746 WRITE(numout,*) ' sm_i : ', sm_i(ki,kj,jl) , ' o_i : ', o_i(ki,kj,jl) 747 WRITE(numout,*) ' a_i : ', a_i(ki,kj,jl) , ' old_a_i : ', old_a_i(ki,kj,jl) 748 WRITE(numout,*) ' d_a_i_trp : ', d_a_i_trp(ki,kj,jl) , ' d_a_i_thd : ', d_a_i_thd(ki,kj,jl) 749 WRITE(numout,*) ' v_i : ', v_i(ki,kj,jl) , ' old_v_i : ', old_v_i(ki,kj,jl) 750 WRITE(numout,*) ' d_v_i_trp : ', d_v_i_trp(ki,kj,jl) , ' d_v_i_thd : ', d_v_i_thd(ki,kj,jl) 751 WRITE(numout,*) ' v_s : ', v_s(ki,kj,jl) , ' old_v_s : ', old_v_s(ki,kj,jl) 752 WRITE(numout,*) ' d_v_s_trp : ', d_v_s_trp(ki,kj,jl) , ' d_v_s_thd : ', d_v_s_thd(ki,kj,jl) 753 WRITE(numout,*) ' e_i1 : ', e_i(ki,kj,1,jl)/1.0e9 , ' old_ei1 : ', old_e_i(ki,kj,1,jl)/1.0e9 754 WRITE(numout,*) ' de_i1_trp : ', d_e_i_trp(ki,kj,1,jl)/1.0e9, ' de_i1_thd : ', d_e_i_thd(ki,kj,1,jl)/1.0e9 755 WRITE(numout,*) ' e_i2 : ', e_i(ki,kj,2,jl)/1.0e9 , ' old_ei2 : ', old_e_i(ki,kj,2,jl)/1.0e9 756 WRITE(numout,*) ' de_i2_trp : ', d_e_i_trp(ki,kj,2,jl)/1.0e9, ' de_i2_thd : ', d_e_i_thd(ki,kj,2,jl)/1.0e9 757 WRITE(numout,*) ' e_snow : ', e_s(ki,kj,1,jl) , ' old_e_snow : ', old_e_s(ki,kj,1,jl) 758 WRITE(numout,*) ' d_e_s_trp : ', d_e_s_trp(ki,kj,1,jl) , ' d_e_s_thd : ', d_e_s_thd(ki,kj,1,jl) 759 WRITE(numout,*) ' smv_i : ', smv_i(ki,kj,jl) , ' old_smv_i : ', old_smv_i(ki,kj,jl) 760 WRITE(numout,*) ' d_smv_i_trp: ', d_smv_i_trp(ki,kj,jl) , ' d_smv_i_thd: ', d_smv_i_thd(ki,kj,jl) 761 WRITE(numout,*) ' oa_i : ', oa_i(ki,kj,jl) , ' old_oa_i : ', old_oa_i(ki,kj,jl) 762 WRITE(numout,*) ' d_oa_i_trp : ', d_oa_i_trp(ki,kj,jl) , ' d_oa_i_thd : ', d_oa_i_thd(ki,kj,jl) 763 WRITE(numout,*) ' Path. case : ', patho_case(ki,kj,jl) 764 END DO !jl 765 766 WRITE(numout,*) 767 WRITE(numout,*) ' - Heat / FW fluxes ' 768 WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' 769 ! WRITE(numout,*) ' sfx_bri : ', sfx_bri (ki,kj) 770 ! WRITE(numout,*) ' sfx : ', sfx (ki,kj) 771 ! WRITE(numout,*) ' sfx_res : ', sfx_res(ki,kj) 772 WRITE(numout,*) ' fmmec : ', fmmec (ki,kj) 773 WRITE(numout,*) ' fhmec : ', fhmec (ki,kj) 774 WRITE(numout,*) ' fhbri : ', fhbri (ki,kj) 775 WRITE(numout,*) ' fheat_mec : ', fheat_mec(ki,kj) 776 WRITE(numout,*) 777 WRITE(numout,*) ' sst : ', sst_m(ki,kj) 778 WRITE(numout,*) ' sss : ', sss_m(ki,kj) 779 WRITE(numout,*) 780 WRITE(numout,*) ' - Stresses ' 781 WRITE(numout,*) ' ~~~~~~~~ ' 782 WRITE(numout,*) ' utau_ice : ', utau_ice(ki,kj) 783 WRITE(numout,*) ' vtau_ice : ', vtau_ice(ki,kj) 784 WRITE(numout,*) ' utau : ', utau (ki,kj) 785 WRITE(numout,*) ' vtau : ', vtau (ki,kj) 786 WRITE(numout,*) ' oc. vel. u : ', u_oce (ki,kj) 787 WRITE(numout,*) ' oc. vel. v : ', v_oce (ki,kj) 788 ENDIF 789 790 !--------------------- 791 ! Salt / heat fluxes 792 !--------------------- 793 794 IF ( kn .EQ. 3 ) THEN 795 WRITE(numout,*) ' lim_prt_state - Point : ',ki,kj 796 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 797 WRITE(numout,*) ' - Salt / Heat Fluxes ' 798 WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' 799 WRITE(numout,*) ' lat - long ', gphit(ki,kj), glamt(ki,kj) 800 WRITE(numout,*) ' Time step ', numit 801 WRITE(numout,*) 802 WRITE(numout,*) ' - Heat fluxes at bottom interface ***' 803 WRITE(numout,*) ' qsr : ', qsr(ki,kj) 804 WRITE(numout,*) ' qns : ', qns(ki,kj) 805 WRITE(numout,*) 806 WRITE(numout,*) ' - Salt fluxes at bottom interface ***' 807 WRITE(numout,*) ' emp : ', emp (ki,kj) 808 WRITE(numout,*) ' sfx_bri : ', sfx_bri(ki,kj) 809 WRITE(numout,*) ' sfx : ', sfx (ki,kj) 810 WRITE(numout,*) ' sfx_res : ', sfx_res(ki,kj) 811 WRITE(numout,*) ' sfx_mec : ', sfx_mec(ki,kj) 812 WRITE(numout,*) ' - Heat fluxes at bottom interface ***' 813 WRITE(numout,*) ' fheat_res : ', fheat_res(ki,kj) 814 WRITE(numout,*) 815 WRITE(numout,*) ' - Momentum fluxes ' 816 WRITE(numout,*) ' utau : ', utau(ki,kj) 817 WRITE(numout,*) ' vtau : ', vtau(ki,kj) 818 ENDIF 819 WRITE(numout,*) ' ' 820 ! 681 DO ji = mi0(ki), mi1(ki) 682 DO jj = mj0(kj), mj1(kj) 683 684 WRITE(numout,*) ' time step ',kt,' ',cd1 ! print title 685 686 !---------------- 687 ! Simple state 688 !---------------- 689 690 IF ( kn == 1 .OR. kn == -1 ) THEN 691 WRITE(numout,*) ' lim_prt_state - Point : ',ji,jj 692 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 693 WRITE(numout,*) ' Simple state ' 694 WRITE(numout,*) ' masks s,u,v : ', tms(ji,jj), tmu(ji,jj), tmv(ji,jj) 695 WRITE(numout,*) ' lat - long : ', gphit(ji,jj), glamt(ji,jj) 696 WRITE(numout,*) ' Time step : ', numit 697 WRITE(numout,*) ' - Ice drift ' 698 WRITE(numout,*) ' ~~~~~~~~~~~ ' 699 WRITE(numout,*) ' u_ice(i-1,j) : ', u_ice(ji-1,jj) 700 WRITE(numout,*) ' u_ice(i ,j) : ', u_ice(ji,jj) 701 WRITE(numout,*) ' v_ice(i ,j-1): ', v_ice(ji,jj-1) 702 WRITE(numout,*) ' v_ice(i ,j) : ', v_ice(ji,jj) 703 WRITE(numout,*) ' strength : ', strength(ji,jj) 704 WRITE(numout,*) 705 WRITE(numout,*) ' - Cell values ' 706 WRITE(numout,*) ' ~~~~~~~~~~~ ' 707 WRITE(numout,*) ' cell area : ', area(ji,jj) 708 WRITE(numout,*) ' at_i : ', at_i(ji,jj) 709 WRITE(numout,*) ' vt_i : ', vt_i(ji,jj) 710 WRITE(numout,*) ' vt_s : ', vt_s(ji,jj) 711 DO jl = 1, jpl 712 WRITE(numout,*) ' - Category (', jl,')' 713 WRITE(numout,*) ' a_i : ', a_i(ji,jj,jl) 714 WRITE(numout,*) ' ht_i : ', ht_i(ji,jj,jl) 715 WRITE(numout,*) ' ht_s : ', ht_s(ji,jj,jl) 716 WRITE(numout,*) ' v_i : ', v_i(ji,jj,jl) 717 WRITE(numout,*) ' v_s : ', v_s(ji,jj,jl) 718 WRITE(numout,*) ' e_s : ', e_s(ji,jj,1,jl)/1.0e9 719 WRITE(numout,*) ' e_i : ', e_i(ji,jj,1:nlay_i,jl)/1.0e9 720 WRITE(numout,*) ' t_su : ', t_su(ji,jj,jl) 721 WRITE(numout,*) ' t_snow : ', t_s(ji,jj,1,jl) 722 WRITE(numout,*) ' t_i : ', t_i(ji,jj,1:nlay_i,jl) 723 WRITE(numout,*) ' sm_i : ', sm_i(ji,jj,jl) 724 WRITE(numout,*) ' smv_i : ', smv_i(ji,jj,jl) 725 WRITE(numout,*) 726 END DO 727 ENDIF 728 IF( kn == -1 ) THEN 729 WRITE(numout,*) ' Mechanical Check ************** ' 730 WRITE(numout,*) ' Check what means ice divergence ' 731 WRITE(numout,*) ' Total ice concentration ', at_i (ji,jj) 732 WRITE(numout,*) ' Total lead fraction ', ato_i(ji,jj) 733 WRITE(numout,*) ' Sum of both ', ato_i(ji,jj) + at_i(ji,jj) 734 WRITE(numout,*) ' Sum of both minus 1 ', ato_i(ji,jj) + at_i(ji,jj) - 1.00 735 ENDIF 736 737 738 !-------------------- 739 ! Exhaustive state 740 !-------------------- 741 742 IF ( kn .EQ. 2 ) THEN 743 WRITE(numout,*) ' lim_prt_state - Point : ',ji,jj 744 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 745 WRITE(numout,*) ' Exhaustive state ' 746 WRITE(numout,*) ' lat - long ', gphit(ji,jj), glamt(ji,jj) 747 WRITE(numout,*) ' Time step ', numit 748 WRITE(numout,*) 749 WRITE(numout,*) ' - Cell values ' 750 WRITE(numout,*) ' ~~~~~~~~~~~ ' 751 WRITE(numout,*) ' cell area : ', area(ji,jj) 752 WRITE(numout,*) ' at_i : ', at_i(ji,jj) 753 WRITE(numout,*) ' vt_i : ', vt_i(ji,jj) 754 WRITE(numout,*) ' vt_s : ', vt_s(ji,jj) 755 WRITE(numout,*) ' u_ice(i-1,j) : ', u_ice(ji-1,jj) 756 WRITE(numout,*) ' u_ice(i ,j) : ', u_ice(ji,jj) 757 WRITE(numout,*) ' v_ice(i ,j-1): ', v_ice(ji,jj-1) 758 WRITE(numout,*) ' v_ice(i ,j) : ', v_ice(ji,jj) 759 WRITE(numout,*) ' strength : ', strength(ji,jj) 760 WRITE(numout,*) ' d_u_ice_dyn : ', d_u_ice_dyn(ji,jj), ' d_v_ice_dyn : ', d_v_ice_dyn(ji,jj) 761 WRITE(numout,*) ' old_u_ice : ', old_u_ice(ji,jj) , ' old_v_ice : ', old_v_ice(ji,jj) 762 WRITE(numout,*) 763 764 DO jl = 1, jpl 765 WRITE(numout,*) ' - Category (',jl,')' 766 WRITE(numout,*) ' ~~~~~~~~ ' 767 WRITE(numout,*) ' ht_i : ', ht_i(ji,jj,jl) , ' ht_s : ', ht_s(ji,jj,jl) 768 WRITE(numout,*) ' t_i : ', t_i(ji,jj,1:nlay_i,jl) 769 WRITE(numout,*) ' t_su : ', t_su(ji,jj,jl) , ' t_s : ', t_s(ji,jj,1,jl) 770 WRITE(numout,*) ' sm_i : ', sm_i(ji,jj,jl) , ' o_i : ', o_i(ji,jj,jl) 771 WRITE(numout,*) ' a_i : ', a_i(ji,jj,jl) , ' old_a_i : ', old_a_i(ji,jj,jl) 772 WRITE(numout,*) ' d_a_i_trp : ', d_a_i_trp(ji,jj,jl) , ' d_a_i_thd : ', d_a_i_thd(ji,jj,jl) 773 WRITE(numout,*) ' v_i : ', v_i(ji,jj,jl) , ' old_v_i : ', old_v_i(ji,jj,jl) 774 WRITE(numout,*) ' d_v_i_trp : ', d_v_i_trp(ji,jj,jl) , ' d_v_i_thd : ', d_v_i_thd(ji,jj,jl) 775 WRITE(numout,*) ' v_s : ', v_s(ji,jj,jl) , ' old_v_s : ', old_v_s(ji,jj,jl) 776 WRITE(numout,*) ' d_v_s_trp : ', d_v_s_trp(ji,jj,jl) , ' d_v_s_thd : ', d_v_s_thd(ji,jj,jl) 777 WRITE(numout,*) ' e_i1 : ', e_i(ji,jj,1,jl)/1.0e9 , ' old_ei1 : ', old_e_i(ji,jj,1,jl)/1.0e9 778 WRITE(numout,*) ' de_i1_trp : ', d_e_i_trp(ji,jj,1,jl)/1.0e9, ' de_i1_thd : ', d_e_i_thd(ji,jj,1,jl)/1.0e9 779 WRITE(numout,*) ' e_i2 : ', e_i(ji,jj,2,jl)/1.0e9 , ' old_ei2 : ', old_e_i(ji,jj,2,jl)/1.0e9 780 WRITE(numout,*) ' de_i2_trp : ', d_e_i_trp(ji,jj,2,jl)/1.0e9, ' de_i2_thd : ', d_e_i_thd(ji,jj,2,jl)/1.0e9 781 WRITE(numout,*) ' e_snow : ', e_s(ji,jj,1,jl) , ' old_e_snow : ', old_e_s(ji,jj,1,jl) 782 WRITE(numout,*) ' d_e_s_trp : ', d_e_s_trp(ji,jj,1,jl) , ' d_e_s_thd : ', d_e_s_thd(ji,jj,1,jl) 783 WRITE(numout,*) ' smv_i : ', smv_i(ji,jj,jl) , ' old_smv_i : ', old_smv_i(ji,jj,jl) 784 WRITE(numout,*) ' d_smv_i_trp: ', d_smv_i_trp(ji,jj,jl) , ' d_smv_i_thd: ', d_smv_i_thd(ji,jj,jl) 785 WRITE(numout,*) ' oa_i : ', oa_i(ji,jj,jl) , ' old_oa_i : ', old_oa_i(ji,jj,jl) 786 WRITE(numout,*) ' d_oa_i_trp : ', d_oa_i_trp(ji,jj,jl) , ' d_oa_i_thd : ', d_oa_i_thd(ji,jj,jl) 787 END DO !jl 788 789 WRITE(numout,*) 790 WRITE(numout,*) ' - Heat / FW fluxes ' 791 WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' 792 WRITE(numout,*) ' emp : ', emp (ji,jj) 793 WRITE(numout,*) ' sfx : ', sfx (ji,jj) 794 WRITE(numout,*) ' sfx_thd : ', sfx_thd(ji,jj) 795 WRITE(numout,*) ' sfx_bri : ', sfx_bri (ji,jj) 796 WRITE(numout,*) ' sfx_mec : ', sfx_mec (ji,jj) 797 WRITE(numout,*) ' sfx_res : ', sfx_res(ji,jj) 798 WRITE(numout,*) ' fmmec : ', fmmec (ji,jj) 799 WRITE(numout,*) ' fhmec : ', fhmec (ji,jj) 800 WRITE(numout,*) ' fhbri : ', fhbri (ji,jj) 801 WRITE(numout,*) ' fheat_mec : ', fheat_mec(ji,jj) 802 WRITE(numout,*) 803 WRITE(numout,*) ' sst : ', sst_m(ji,jj) 804 WRITE(numout,*) ' sss : ', sss_m(ji,jj) 805 WRITE(numout,*) 806 WRITE(numout,*) ' - Stresses ' 807 WRITE(numout,*) ' ~~~~~~~~ ' 808 WRITE(numout,*) ' utau_ice : ', utau_ice(ji,jj) 809 WRITE(numout,*) ' vtau_ice : ', vtau_ice(ji,jj) 810 WRITE(numout,*) ' utau : ', utau (ji,jj) 811 WRITE(numout,*) ' vtau : ', vtau (ji,jj) 812 WRITE(numout,*) ' oc. vel. u : ', u_oce (ji,jj) 813 WRITE(numout,*) ' oc. vel. v : ', v_oce (ji,jj) 814 ENDIF 815 816 !--------------------- 817 ! Salt / heat fluxes 818 !--------------------- 819 820 IF ( kn .EQ. 3 ) THEN 821 WRITE(numout,*) ' lim_prt_state - Point : ',ji,jj 822 WRITE(numout,*) ' ~~~~~~~~~~~~~~ ' 823 WRITE(numout,*) ' - Salt / Heat Fluxes ' 824 WRITE(numout,*) ' ~~~~~~~~~~~~~~~~ ' 825 WRITE(numout,*) ' lat - long ', gphit(ji,jj), glamt(ji,jj) 826 WRITE(numout,*) ' Time step ', numit 827 WRITE(numout,*) 828 WRITE(numout,*) ' - Heat fluxes at bottom interface ***' 829 WRITE(numout,*) ' qsr : ', qsr(ji,jj) 830 WRITE(numout,*) ' qns : ', qns(ji,jj) 831 WRITE(numout,*) ' fdtcn : ', fdtcn(ji,jj) 832 WRITE(numout,*) ' qcmif : ', qcmif(ji,jj) * r1_rdtice 833 WRITE(numout,*) ' qldif : ', qldif(ji,jj) * r1_rdtice 834 WRITE(numout,*) 835 WRITE(numout,*) ' - Salt fluxes at bottom interface ***' 836 WRITE(numout,*) ' emp : ', emp (ji,jj) 837 WRITE(numout,*) ' sfx_bri : ', sfx_bri(ji,jj) 838 WRITE(numout,*) ' sfx : ', sfx (ji,jj) 839 WRITE(numout,*) ' sfx_res : ', sfx_res(ji,jj) 840 WRITE(numout,*) ' sfx_mec : ', sfx_mec(ji,jj) 841 WRITE(numout,*) ' - Heat fluxes at bottom interface ***' 842 WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj) 843 WRITE(numout,*) 844 WRITE(numout,*) ' - Momentum fluxes ' 845 WRITE(numout,*) ' utau : ', utau(ji,jj) 846 WRITE(numout,*) ' vtau : ', vtau(ji,jj) 847 ENDIF 848 WRITE(numout,*) ' ' 849 ! 850 END DO 851 END DO 852 821 853 END SUBROUTINE lim_prt_state 822 854 -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
r4230 r4333 343 343 ! ! (update freshwater fluxes) 344 344 !RBbug do not understand why see ticket 667 345 CALL lbc_lnk( emp, 'T', 1. )345 !clem-bugsal CALL lbc_lnk( emp, 'T', 1. ) 346 346 ! 347 347 IF( kt == nit000 ) THEN ! set the forcing field at nit000 - 1 ! -
branches/2013/dev_MERGE_2013/NEMOGCM/NEMO/OPA_SRC/TRA/traqsr.F90
r4292 r4333 33 33 USE timing ! Timing 34 34 USE sbc_ice, ONLY : lk_lim3 35 36 35 37 36 IMPLICIT NONE
Note: See TracChangeset
for help on using the changeset viewer.