- Timestamp:
- 2015-02-11T18:27:52+01:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90
r5078 r5080 123 123 !! and Elizabeth C. Hunke, LANL are gratefully acknowledged 124 124 !!--------------------------------------------------------------------! 125 INTEGER :: ji, jj, jk, jl ! dummy loop index 126 INTEGER :: niter, nitermax = 20 ! local integer 127 LOGICAL :: asum_error ! flag for asum .ne. 1 125 INTEGER :: ji, jj, jk, jl ! dummy loop index 126 INTEGER :: niter ! local integer 128 127 INTEGER :: iterate_ridging ! if true, repeat the ridging 129 REAL(wp) :: w1, tmpfac! local scalar128 REAL(wp) :: za, zfac ! local scalar 130 129 CHARACTER (len = 15) :: fieldid 131 130 REAL(wp), POINTER, DIMENSION(:,:) :: closing_net ! net rate at which area is removed (1/s) … … 137 136 REAL(wp), POINTER, DIMENSION(:,:) :: esnow_mlt ! energy needed to melt snow in ocean (J m-2) 138 137 REAL(wp), POINTER, DIMENSION(:,:) :: vt_i_init, vt_i_final ! ice volume summed over categories 138 ! 139 INTEGER, PARAMETER :: nitermax = 20 139 140 ! 140 141 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b … … 236 237 ! would be removed. Reduce the opening rate proportionately. 237 238 IF ( ato_i(ji,jj) > epsi10 .AND. athorn(ji,jj,0) > 0.0 ) THEN 238 w1= athorn(ji,jj,0) * closing_gross(ji,jj) * rdt_ice239 IF ( w1> ato_i(ji,jj)) THEN240 tmpfac = ato_i(ji,jj) / w1241 closing_gross(ji,jj) = closing_gross(ji,jj) * tmpfac242 opning(ji,jj) = opning(ji,jj) * tmpfac243 ENDIF !w1244 ENDIF !at0i and athorn245 246 END DO ! ji247 END DO ! jj239 za = athorn(ji,jj,0) * closing_gross(ji,jj) * rdt_ice 240 IF ( za > ato_i(ji,jj)) THEN 241 zfac = ato_i(ji,jj) / za 242 closing_gross(ji,jj) = closing_gross(ji,jj) * zfac 243 opning(ji,jj) = opning(ji,jj) * zfac 244 ENDIF 245 ENDIF 246 247 END DO 248 END DO 248 249 249 250 ! correction to closing rate / opening if excessive ice removal … … 256 257 DO ji = 1, jpi 257 258 IF ( a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0._wp )THEN 258 w1= athorn(ji,jj,jl) * closing_gross(ji,jj) * rdt_ice259 IF ( w1> a_i(ji,jj,jl) ) THEN260 tmpfac = a_i(ji,jj,jl) / w1261 closing_gross(ji,jj) = closing_gross(ji,jj) * tmpfac262 opning (ji,jj) = opning (ji,jj) * tmpfac259 za = athorn(ji,jj,jl) * closing_gross(ji,jj) * rdt_ice 260 IF ( za > a_i(ji,jj,jl) ) THEN 261 zfac = a_i(ji,jj,jl) / za 262 closing_gross(ji,jj) = closing_gross(ji,jj) * zfac 263 opning (ji,jj) = opning (ji,jj) * zfac 263 264 ENDIF 264 265 ENDIF 265 END DO !ji266 END DO ! jj267 END DO !jl266 END DO 267 END DO 268 END DO 268 269 269 270 ! 3.3 Redistribute area, volume, and energy. … … 322 323 ! Convert ridging rate diagnostics to correct units. 323 324 ! Update fresh water and heat fluxes due to snow melt. 324 325 asum_error = .false.326 327 325 DO jj = 1, jpj 328 326 DO ji = 1, jpi 329 330 IF(ABS(asum(ji,jj) - kamax) > epsi10 ) asum_error = .true.331 327 332 328 dardg1dt(ji,jj) = dardg1dt(ji,jj) * r1_rdtice … … 434 430 !!---------------------------------------------------------------------- 435 431 INTEGER, INTENT(in) :: kstrngth ! = 1 for Rothrock formulation, 0 for Hibler (1979) 436 437 INTEGER :: ji,jj, jl ! dummy loop indices 438 INTEGER :: ksmooth ! smoothing the resistance to deformation 439 INTEGER :: numts_rm ! number of time steps for the P smoothing 440 REAL(wp) :: hi, zw1, zp, zdummy, zzc, z1_3 ! local scalars 432 INTEGER :: ji,jj, jl ! dummy loop indices 433 INTEGER :: ksmooth ! smoothing the resistance to deformation 434 INTEGER :: numts_rm ! number of time steps for the P smoothing 435 REAL(wp) :: zhi, zp, z1_3 ! local scalars 441 436 REAL(wp), POINTER, DIMENSION(:,:) :: zworka ! temporary array used here 442 437 !!---------------------------------------------------------------------- … … 464 459 ! 465 460 IF( a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0._wp ) THEN 466 hi = v_i(ji,jj,jl) / a_i(ji,jj,jl)461 zhi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 467 462 !---------------------------- 468 463 ! PE loss from deforming ice 469 464 !---------------------------- 470 strength(ji,jj) = strength(ji,jj) - athorn(ji,jj,jl) * hi *hi465 strength(ji,jj) = strength(ji,jj) - athorn(ji,jj,jl) * zhi * zhi 471 466 472 467 !-------------------------- 473 468 ! PE gain from rafting ice 474 469 !-------------------------- 475 strength(ji,jj) = strength(ji,jj) + 2._wp * araft(ji,jj,jl) * hi *hi470 strength(ji,jj) = strength(ji,jj) + 2._wp * araft(ji,jj,jl) * zhi * zhi 476 471 477 472 !---------------------------- … … 479 474 !---------------------------- 480 475 strength(ji,jj) = strength(ji,jj) + aridge(ji,jj,jl)/krdg(ji,jj,jl) & 481 * z1_3 * ( hrmax(ji,jj,jl)**3 - hrmin(ji,jj,jl)**3) / ( hrmax(ji,jj,jl)-hrmin(ji,jj,jl) )482 !!gm Optimization: (a**3-b**3)/(a-b) = a*a+ab+b*b ==> less costly operations even if a**3 is replaced by a*a*a...476 * z1_3 * ( hrmax(ji,jj,jl)**2 + hrmin(ji,jj,jl)**2 + hrmax(ji,jj,jl) * hrmin(ji,jj,jl) ) 477 !!(a**3-b**3)/(a-b) = a*a+ab+b*b 483 478 ENDIF 484 479 ! … … 486 481 END DO 487 482 END DO 488 489 zzc = rn_pe_rdg * Cp ! where Cp = (g/2)*(rhow-rhoi)*(rhoi/rhow) and rn_pe_rdg accounts for frictional dissipation 490 strength(:,:) = zzc * strength(:,:) / aksum(:,:) 491 483 484 strength(:,:) = rn_pe_rdg * Cp * strength(:,:) / aksum(:,:) 485 ! where Cp = (g/2)*(rhow-rhoi)*(rhoi/rhow) and rn_pe_rdg accounts for frictional dissipation 492 486 ksmooth = 1 493 487 … … 513 507 DO jj = 1, jpj 514 508 DO ji = 1, jpi 515 IF ( bv_i(ji,jj) > 0.0 ) THEN516 zdummy = MIN ( bv_i(ji,jj), 0.10 ) * MIN( bv_i(ji,jj), 0.10 )517 ELSE518 zdummy = 0.0519 ENDIF520 509 strength(ji,jj) = strength(ji,jj) * exp(-5.88*SQRT(MAX(bv_i(ji,jj),0.0))) 521 510 END DO … … 536 525 CALL lbc_lnk( strength, 'T', 1. ) 537 526 538 DO jj = 2, jpj -1539 DO ji = 2, jpi -1527 DO jj = 2, jpjm1 528 DO ji = 2, jpim1 540 529 IF ( ( asum(ji,jj) - ato_i(ji,jj) ) > epsi10) THEN ! ice is present 541 530 zworka(ji,jj) = 4.0 * strength(ji,jj) & … … 545 534 & + strength(ji,jj+1) * tmask(ji,jj+1,1) 546 535 547 zw 1 = 4.0 + tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1)548 zworka(ji,jj) = zworka(ji,jj) / zw1536 zworka(ji,jj) = zworka(ji,jj) / & 537 & ( 4.0 + tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1) ) 549 538 ELSE 550 539 zworka(ji,jj) = 0._wp … … 553 542 END DO 554 543 555 DO jj = 2, jpj -1556 DO ji = 2, jpi -1544 DO jj = 2, jpjm1 545 DO ji = 2, jpim1 557 546 strength(ji,jj) = zworka(ji,jj) 558 547 END DO … … 609 598 !!---------------------------------------------------------------------! 610 599 INTEGER :: ji,jj, jl ! dummy loop indices 611 REAL(wp) :: Gstari, astari, hi, hrmean, zdummy ! local scalar600 REAL(wp) :: Gstari, astari, zhi, hrmean, zdummy ! local scalar 612 601 REAL(wp), POINTER, DIMENSION(:,:) :: zworka ! temporary array used here 613 602 REAL(wp), POINTER, DIMENSION(:,:,:) :: Gsum ! Gsum(n) = sum of areas in categories 0 to n … … 697 686 athorn(ji,jj,jl) = 0.0 698 687 ENDIF 699 END DO ! ji700 END DO ! jj701 END DO ! jl688 END DO 689 END DO 690 END DO 702 691 703 692 ELSE !--- Exponential, more stable formulation (Lipscomb et al, 2007) … … 791 780 792 781 IF (a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0.0 ) THEN 793 hi = v_i(ji,jj,jl) / a_i(ji,jj,jl)794 hrmean = MAX(SQRT(rn_hstar* hi),hi*krdgmin)795 hrmin(ji,jj,jl) = MIN(2.0* hi, 0.5*(hrmean +hi))782 zhi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 783 hrmean = MAX(SQRT(rn_hstar*zhi), zhi*krdgmin) 784 hrmin(ji,jj,jl) = MIN(2.0*zhi, 0.5*(hrmean + zhi)) 796 785 hrmax(ji,jj,jl) = 2.0*hrmean - hrmin(ji,jj,jl) 797 hraft(ji,jj,jl) = kraft* hi798 krdg(ji,jj,jl) = hrmean / hi786 hraft(ji,jj,jl) = kraft*zhi 787 krdg(ji,jj,jl) = hrmean / zhi 799 788 ELSE 800 789 hraft(ji,jj,jl) = 0.0 … … 844 833 INTEGER :: ij ! horizontal index, combines i and j loops 845 834 INTEGER :: icells ! number of cells with aicen > puny 846 REAL(wp) :: hL, hR, farea, z dummy, zdummy0, ztmelts ! left and right limits of integration835 REAL(wp) :: hL, hR, farea, ztmelts ! left and right limits of integration 847 836 REAL(wp) :: zsstK ! SST in Kelvin 848 837
Note: See TracChangeset
for help on using the changeset viewer.