Changeset 921 for trunk/NEMO/LIM_SRC_3/limthd_lac.F90
- Timestamp:
- 2008-05-13T10:28:52+02:00 (16 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/LIM_SRC_3/limthd_lac.F90 (modified) (10 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/LIM_SRC_3/limthd_lac.F90
r888 r921 25 25 USE limtab 26 26 USE limcons 27 27 28 28 IMPLICIT NONE 29 29 PRIVATE … … 50 50 51 51 CONTAINS 52 52 53 53 SUBROUTINE lim_thd_lac 54 54 !!------------------------------------------------------------------- … … 146 146 zalphai , & !: factor describing how old and new layers overlap each other [m] 147 147 zindb 148 148 149 149 REAL(wp), DIMENSION(jpij,jkmax+1,jpl) :: & 150 150 zqm0 , & !: old layer-system heat content … … 188 188 189 189 !!-----------------------------------------------------------------------! 190 190 191 191 et_i_init(:,:) = 0.0 192 192 et_s_init(:,:) = 0.0 … … 195 195 zeps6 = 1.0e-6 196 196 197 !------------------------------------------------------------------------------!198 ! 1) Conservation check and changes in each ice category199 !------------------------------------------------------------------------------!197 !------------------------------------------------------------------------------! 198 ! 1) Conservation check and changes in each ice category 199 !------------------------------------------------------------------------------! 200 200 IF ( con_i ) THEN 201 201 CALL lim_column_sum (jpl, v_i, vt_i_init) … … 205 205 ENDIF 206 206 207 !------------------------------------------------------------------------------|208 ! 2) Convert units for ice internal energy209 !------------------------------------------------------------------------------|207 !------------------------------------------------------------------------------| 208 ! 2) Convert units for ice internal energy 209 !------------------------------------------------------------------------------| 210 210 DO jl = 1, jpl 211 DO jk = 1, nlay_i 212 DO jj = 1, jpj 213 DO ji = 1, jpi 214 !Energy of melting q(S,T) [J.m-3] 215 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / & 216 MAX( area(ji,jj) * v_i(ji,jj,jl) , zeps ) * & 217 nlay_i 218 zindb = 1.0-MAX(0.0,SIGN(1.0,-v_i(ji,jj,jl))) !0 if no ice and 1 if yes 219 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl)*unit_fac*zindb 211 DO jk = 1, nlay_i 212 DO jj = 1, jpj 213 DO ji = 1, jpi 214 !Energy of melting q(S,T) [J.m-3] 215 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / & 216 MAX( area(ji,jj) * v_i(ji,jj,jl) , zeps ) * & 217 nlay_i 218 zindb = 1.0-MAX(0.0,SIGN(1.0,-v_i(ji,jj,jl))) !0 if no ice and 1 if yes 219 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl)*unit_fac*zindb 220 END DO 220 221 END DO 221 END DO 222 END DO 222 END DO 223 223 END DO 224 224 225 !------------------------------------------------------------------------------!226 ! 3) Collection thickness of ice formed in leads and polynyas227 !------------------------------------------------------------------------------!225 !------------------------------------------------------------------------------! 226 ! 3) Collection thickness of ice formed in leads and polynyas 227 !------------------------------------------------------------------------------! 228 228 ! hicol is the thickness of new ice formed in open water 229 229 ! hicol can be either prescribed (frazswi = 0) … … 248 248 IF (fraz_swi.eq.1.0) THEN 249 249 250 251 252 253 254 255 256 257 258 259 260 261 DO ji = 1, jpi262 263 IF ( tms(ji,jj) * ( qcmif(ji,jj) - qldif(ji,jj) ) > 0.e0 ) THEN264 !-------------265 ! Wind stress266 !-------------267 ! C-grid wind stress components268 ztaux = ( utaui_ice(ji-1,jj ) * tmu(ji-1,jj ) &269 270 ztauy = ( vtaui_ice(ji ,jj-1) * tmv(ji ,jj-1) &271 272 ! Square root of wind stress273 ztenagm = SQRT( SQRT( ztaux * ztaux + ztauy * ztauy ) )274 275 !---------------------276 ! Frazil ice velocity277 !---------------------278 zvfrx = zgamafr * zsqcd * ztaux / MAX(ztenagm,zeps)279 zvfry = zgamafr * zsqcd * ztauy / MAX(ztenagm,zeps)280 281 !-------------------282 ! Pack ice velocity283 !-------------------284 ! C-grid ice velocity285 zindb = MAX(0.0, SIGN(1.0, at_i(ji,jj) ))286 zvgx = zindb * ( u_ice(ji-1,jj ) * tmu(ji-1,jj ) &287 288 zvgy = zindb * ( v_ice(ji ,jj-1) * tmv(ji ,jj-1) &289 290 291 !-----------------------------------292 ! Relative frazil/pack ice velocity293 !-----------------------------------294 ! absolute relative velocity295 zvrel2 = MAX( ( zvfrx - zvgx ) * ( zvfrx - zvgx ) + &296 297 298 zvrel(ji,jj) = SQRT(zvrel2)299 300 !---------------------301 ! Iterative procedure302 !---------------------303 hicol(ji,jj) = zhicrit + 0.1304 hicol(ji,jj) = zhicrit + hicol(ji,jj) / &305 306 307 308 iter = 1309 iterate_frazil = .true.310 311 DO WHILE ( iter .LT. 100 .AND. iterate_frazil )312 zf = ( hicol(ji,jj) - zhicrit ) * ( hicol(ji,jj)**2 - zhicrit**2 ) &313 314 zfp = ( hicol(ji,jj) - zhicrit ) * ( 3.0*hicol(ji,jj) + zhicrit ) &315 316 zhicol_new = hicol(ji,jj) - zf/zfp317 hicol(ji,jj) = zhicol_new318 319 iter = iter + 1320 321 END DO ! do while322 323 ENDIF ! end of selection of pixels where ice forms324 325 END DO ! loop on ji ends326 END DO ! loop on jj ends250 !-------------------- 251 ! Physical constants 252 !-------------------- 253 hicol(:,:) = 0.0 254 255 zhicrit = 0.04 ! frazil ice thickness 256 ztwogp = 2. * rau0 / ( grav * 0.3 * ( rau0 - rhoic ) ) ! reduced grav 257 zsqcd = 1.0 / SQRT( 1.3 * cai ) ! 1/SQRT(airdensity*drag) 258 zgamafr = 0.03 259 260 DO jj = 1, jpj 261 DO ji = 1, jpi 262 263 IF ( tms(ji,jj) * ( qcmif(ji,jj) - qldif(ji,jj) ) > 0.e0 ) THEN 264 !------------- 265 ! Wind stress 266 !------------- 267 ! C-grid wind stress components 268 ztaux = ( utaui_ice(ji-1,jj ) * tmu(ji-1,jj ) & 269 + utaui_ice(ji ,jj ) * tmu(ji ,jj ) ) / 2.0 270 ztauy = ( vtaui_ice(ji ,jj-1) * tmv(ji ,jj-1) & 271 + vtaui_ice(ji ,jj ) * tmv(ji ,jj ) ) / 2.0 272 ! Square root of wind stress 273 ztenagm = SQRT( SQRT( ztaux * ztaux + ztauy * ztauy ) ) 274 275 !--------------------- 276 ! Frazil ice velocity 277 !--------------------- 278 zvfrx = zgamafr * zsqcd * ztaux / MAX(ztenagm,zeps) 279 zvfry = zgamafr * zsqcd * ztauy / MAX(ztenagm,zeps) 280 281 !------------------- 282 ! Pack ice velocity 283 !------------------- 284 ! C-grid ice velocity 285 zindb = MAX(0.0, SIGN(1.0, at_i(ji,jj) )) 286 zvgx = zindb * ( u_ice(ji-1,jj ) * tmu(ji-1,jj ) & 287 + u_ice(ji,jj ) * tmu(ji ,jj ) ) / 2.0 288 zvgy = zindb * ( v_ice(ji ,jj-1) * tmv(ji ,jj-1) & 289 + v_ice(ji,jj ) * tmv(ji ,jj ) ) / 2.0 290 291 !----------------------------------- 292 ! Relative frazil/pack ice velocity 293 !----------------------------------- 294 ! absolute relative velocity 295 zvrel2 = MAX( ( zvfrx - zvgx ) * ( zvfrx - zvgx ) + & 296 ( zvfry - zvgy ) * ( zvfry - zvgy ) & 297 , 0.15 * 0.15 ) 298 zvrel(ji,jj) = SQRT(zvrel2) 299 300 !--------------------- 301 ! Iterative procedure 302 !--------------------- 303 hicol(ji,jj) = zhicrit + 0.1 304 hicol(ji,jj) = zhicrit + hicol(ji,jj) / & 305 ( hicol(ji,jj) * hicol(ji,jj) - & 306 zhicrit * zhicrit ) * ztwogp * zvrel2 307 308 iter = 1 309 iterate_frazil = .true. 310 311 DO WHILE ( iter .LT. 100 .AND. iterate_frazil ) 312 zf = ( hicol(ji,jj) - zhicrit ) * ( hicol(ji,jj)**2 - zhicrit**2 ) & 313 - hicol(ji,jj) * zhicrit * ztwogp * zvrel2 314 zfp = ( hicol(ji,jj) - zhicrit ) * ( 3.0*hicol(ji,jj) + zhicrit ) & 315 - zhicrit * ztwogp * zvrel2 316 zhicol_new = hicol(ji,jj) - zf/zfp 317 hicol(ji,jj) = zhicol_new 318 319 iter = iter + 1 320 321 END DO ! do while 322 323 ENDIF ! end of selection of pixels where ice forms 324 325 END DO ! loop on ji ends 326 END DO ! loop on jj ends 327 327 328 328 ENDIF ! End of computation of frazil ice collection thickness 329 329 330 !------------------------------------------------------------------------------!331 ! 4) Identify grid points where new ice forms332 !------------------------------------------------------------------------------!330 !------------------------------------------------------------------------------! 331 ! 4) Identify grid points where new ice forms 332 !------------------------------------------------------------------------------! 333 333 334 334 !------------------------------------- … … 349 349 END DO 350 350 351 IF( lwp) THEN351 IF( ln_nicep ) THEN 352 352 WRITE(numout,*) 'lim_thd_lac : nbpac = ', nbpac 353 353 ENDIF … … 360 360 361 361 IF ( nbpac > 0 ) THEN 362 363 CALL tab_2d_1d( nbpac, zat_i_ac (1:nbpac) , at_i , &364 365 DO jl = 1, jpl366 CALL tab_2d_1d( nbpac, za_i_ac(1:nbpac,jl) , a_i(:,:,jl) , &367 368 CALL tab_2d_1d( nbpac, zv_i_ac(1:nbpac,jl) , v_i(:,:,jl) , &369 370 CALL tab_2d_1d( nbpac, zoa_i_ac(1:nbpac,jl) , oa_i(:,:,jl) , &371 372 CALL tab_2d_1d( nbpac, zsmv_i_ac(1:nbpac,jl), smv_i(:,:,jl), &373 374 DO jk = 1, nlay_i375 CALL tab_2d_1d( nbpac, ze_i_ac(1:nbpac,jk,jl), e_i(:,:,jk,jl) , &376 377 END DO ! jk378 END DO ! jl379 380 CALL tab_2d_1d( nbpac, qldif_1d (1:nbpac) , qldif , &381 382 CALL tab_2d_1d( nbpac, qcmif_1d (1:nbpac) , qcmif , &383 384 CALL tab_2d_1d( nbpac, t_bo_b (1:nbpac) , t_bo , &385 386 CALL tab_2d_1d( nbpac, fseqv_1d (1:nbpac) , fseqv , &387 388 CALL tab_2d_1d( nbpac, hicol_b (1:nbpac) , hicol , &389 390 CALL tab_2d_1d( nbpac, zvrel_ac (1:nbpac) , zvrel , &391 392 393 !------------------------------------------------------------------------------!394 ! 5) Compute thickness, salinity, enthalpy, age, area and volume of new ice395 !------------------------------------------------------------------------------!396 397 !----------------------398 ! Thickness of new ice399 !----------------------400 DO ji = 1, nbpac401 zh_newice(ji) = hiccrit(1)402 END DO403 IF ( fraz_swi .EQ. 1.0 ) zh_newice(:) = hicol_b(:)404 405 !----------------------406 ! Salinity of new ice407 !----------------------408 409 IF ( num_sal .EQ. 1 ) THEN410 zs_newice(:) = bulk_sal411 ENDIF ! num_sal412 413 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN414 415 DO ji = 1, nbpac416 zs_newice(ji) = MIN( 4.606 + 0.91 / zh_newice(ji) , s_i_max )417 zji = MOD( npac(ji) - 1, jpi ) + 1418 zjj = ( npac(ji) - 1 ) / jpi + 1419 zs_newice(ji) = MIN( 0.5*sss_m(zji,zjj) , zs_newice(ji) )420 END DO ! jl421 422 ENDIF ! num_sal423 424 IF ( num_sal .EQ. 3 ) THEN425 zs_newice(:) = 2.3426 ENDIF ! num_sal427 428 !-------------------------429 ! Heat content of new ice430 !-------------------------431 ! We assume that new ice is formed at the seawater freezing point432 DO ji = 1, nbpac433 ztmelts = - tmut * zs_newice(ji) + rtt ! Melting point (K)434 ze_newice(ji) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) &435 436 437 438 ze_newice(ji) = MAX( ze_newice(ji) , 0.0 ) + &439 440 441 END DO ! ji442 !----------------443 ! Age of new ice444 !----------------445 DO ji = 1, nbpac446 zo_newice(ji) = 0.0447 END DO ! ji448 449 !--------------------------450 ! Open water energy budget451 !--------------------------452 DO ji = 1, nbpac453 zqbgow(ji) = qldif_1d(ji) - qcmif_1d(ji) !<0454 END DO ! ji455 456 !-------------------457 ! Volume of new ice458 !-------------------459 DO ji = 1, nbpac460 zv_newice(ji) = - zqbgow(ji) / ze_newice(ji)461 462 ! A fraction zfrazb of frazil ice is accreted at the ice bottom463 zfrazb = ( TANH ( Cfrazb * ( zvrel_ac(ji) - vfrazb ) ) &464 465 zdh_frazb(ji) = zfrazb*zv_newice(ji)466 zv_newice(ji) = ( 1.0 - zfrazb ) * zv_newice(ji)467 END DO468 469 !---------------------------------470 ! Salt flux due to new ice growth471 !---------------------------------472 IF ( ( num_sal .EQ. 4 ) ) THEN473 DO ji = 1, nbpac474 zji = MOD( npac(ji) - 1, jpi ) + 1475 zjj = ( npac(ji) - 1 ) / jpi + 1476 fseqv_1d(ji) = fseqv_1d(ji) + &477 478 479 END DO480 ELSE481 DO ji = 1, nbpac482 zji = MOD( npac(ji) - 1, jpi ) + 1483 zjj = ( npac(ji) - 1 ) / jpi + 1484 fseqv_1d(ji) = fseqv_1d(ji) + &485 486 487 END DO ! ji488 ENDIF489 490 !------------------------------------491 ! Diags for energy conservation test492 !------------------------------------493 DO ji = 1, nbpac494 ! Volume495 zji = MOD( npac(ji) - 1, jpi ) + 1496 zjj = ( npac(ji) - 1 ) / jpi + 1497 vt_i_init(zji,zjj) = vt_i_init(zji,zjj) + zv_newice(ji)498 ! Energy499 zde = ze_newice(ji) / unit_fac500 zde = zde * area(zji,zjj) * zv_newice(ji)501 et_i_init(zji,zjj) = et_i_init(zji,zjj) + zde502 END DO503 504 ! keep new ice volume in memory505 CALL tab_1d_2d( nbpac, v_newice , npac(1:nbpac), zv_newice(1:nbpac) , &506 507 508 !-----------------509 ! Area of new ice510 !-----------------511 DO ji = 1, nbpac512 za_newice(ji) = zv_newice(ji) / zh_newice(ji)513 ! diagnostic514 zji = MOD( npac(ji) - 1, jpi ) + 1515 zjj = ( npac(ji) - 1 ) / jpi + 1516 diag_lat_gr(zji,zjj) = zv_newice(ji) / rdt_ice517 END DO !ji518 519 !------------------------------------------------------------------------------!520 ! 6) Redistribute new ice area and volume into ice categories !521 !------------------------------------------------------------------------------!522 523 !-----------------------------------------524 ! Keep old ice areas and volume in memory525 !-----------------------------------------526 zv_old(:,:) = zv_i_ac(:,:)527 za_old(:,:) = za_i_ac(:,:)528 529 !-------------------------------------------530 ! Compute excessive new ice area and volume531 !-------------------------------------------532 ! If lateral ice growth gives an ice concentration gt 1, then533 ! we keep the excessive volume in memory and attribute it later534 ! to bottom accretion535 DO ji = 1, nbpac536 ! vectorize537 IF ( za_newice(ji) .GT. ( 1.0 - zat_i_ac(ji) ) ) THEN538 zda_res(ji) = za_newice(ji) - (1.0 - zat_i_ac(ji) )539 zdv_res(ji) = zda_res(ji) * zh_newice(ji)540 za_newice(ji) = za_newice(ji) - zda_res(ji)541 zv_newice(ji) = zv_newice(ji) - zdv_res(ji)542 ELSE543 zda_res(ji) = 0.0544 zdv_res(ji) = 0.0545 ENDIF546 END DO ! ji547 548 !------------------------------------------------549 ! Laterally redistribute new ice volume and area550 !------------------------------------------------551 zat_i_ac(:) = 0.0552 553 DO jl = 1, jpl554 DO ji = 1, nbpac555 ! vectorize556 IF ( ( hi_max(jl-1) .LT. zh_newice(ji) ) &557 558 za_i_ac(ji,jl) = za_i_ac(ji,jl) + za_newice(ji)559 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + zv_newice(ji)560 zat_i_ac(ji) = zat_i_ac(ji) + za_i_ac(ji,jl)561 zcatac(ji) = jl562 ENDIF563 END DO ! ji564 END DO ! jl565 566 !----------------------------------567 ! Heat content - lateral accretion568 !----------------------------------569 DO ji = 1, nbpac570 jl = zcatac(ji) ! categroy in which new ice is put571 ! zindb = 0 if no ice and 1 if yes572 zindb = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , -za_old(ji,jl) ) )573 ! old ice thickness574 zhice_old(ji,jl) = zv_old(ji,jl) &575 576 ! difference in thickness577 zdhex(ji) = MAX( 0.0, zh_newice(ji) - zhice_old(ji,jl) )578 ! is ice totally new in category jl ?579 zswinew(ji) = MAX( 0.0, SIGN( 1.0 , - za_old(ji,jl) + epsi11 ) )580 END DO581 582 DO jk = 1, nlay_i583 DO ji = 1, nbpac584 jl = zcatac(ji)585 zqold = ze_i_ac(ji,jk,jl) ! [ J.m-3 ]586 zalphai = MIN( zhice_old(ji,jl) * jk / nlay_i , &587 588 589 590 ze_i_ac(ji,jk,jl) = &591 zswinew(ji) * ze_newice(ji) &592 + ( 1.0 - zswinew(ji) ) * &593 ( za_old(ji,jl) * zqold * zhice_old(ji,jl) / nlay_i &594 + za_newice(ji) * ze_newice(ji) * zalphai &595 + za_newice(ji) * ze_newice(ji) * zdhex(ji) / nlay_i ) / &596 ( ( zv_i_ac(ji,jl) ) / nlay_i )597 598 END DO !ji599 END DO !jl600 601 !-----------------------------------------------602 ! Add excessive volume of new ice at the bottom603 !-----------------------------------------------604 ! If the ice concentration exceeds 1, the remaining volume of new ice605 ! is equally redistributed among all ice categories in which there is606 ! ice607 608 ! Fraction of level ice609 jm = 1610 zat_i_lev(:) = 0.0611 612 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)613 DO ji = 1, nbpac614 zat_i_lev(ji) = zat_i_lev(ji) + za_i_ac(ji,jl)615 END DO616 END DO617 618 WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindx, 1:jpl)619 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)620 DO ji = 1, nbpac621 zindb = MAX( 0.0, SIGN( 1.0, zdv_res(ji) ) )622 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + &623 624 625 END DO ! ji626 END DO ! jl627 WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindx, 1:jpl)628 629 !---------------------------------630 ! Heat content - bottom accretion631 !---------------------------------632 jm = 1633 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)634 DO ji = 1, nbpac635 ! zindb = 0 if no ice and 1 if yes636 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 &637 638 zhice_old(ji,jl) = zv_i_ac(ji,jl) / &639 640 zdhicbot(ji,jl) = zdv_res(ji) / MAX( za_i_ac(ji,jl) , zeps ) &641 642 643 644 ! thickness of residual ice645 zdummy(ji,jl) = zv_i_ac(ji,jl)/MAX(za_i_ac(ji,jl),zeps)*zindb646 END DO !ji647 END DO !jl648 649 ! old layers thicknesses and enthalpies650 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)651 DO jk = 1, nlay_i652 DO ji = 1, nbpac653 zthick0(ji,jk,jl)= zhice_old(ji,jl) / nlay_i654 zqm0 (ji,jk,jl)= ze_i_ac(ji,jk,jl) * zthick0(ji,jk,jl)655 END DO !ji656 END DO !jk657 END DO !jl658 659 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)660 DO ji = 1, nbpac661 zthick0(ji,nlay_i+1,jl) = zdhicbot(ji,jl)662 zqm0 (ji,nlay_i+1,jl) = ze_newice(ji)*zdhicbot(ji,jl)663 END DO ! ji664 END DO ! jl665 666 ! Redistributing energy on the new grid667 ze_i_ac(:,:,:) = 0.0668 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)669 DO jk = 1, nlay_i670 DO layer = 1, nlay_i + 1671 DO ji = 1, nbpac672 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , &673 674 ! Redistributing energy on the new grid675 zweight = MAX ( &676 MIN( zhice_old(ji,jl) * layer , zdummy(ji,jl) * jk ) - &677 MAX( zhice_old(ji,jl) * ( layer - 1 ) , zdummy(ji,jl) * &362 363 CALL tab_2d_1d( nbpac, zat_i_ac (1:nbpac) , at_i , & 364 jpi, jpj, npac(1:nbpac) ) 365 DO jl = 1, jpl 366 CALL tab_2d_1d( nbpac, za_i_ac(1:nbpac,jl) , a_i(:,:,jl) , & 367 jpi, jpj, npac(1:nbpac) ) 368 CALL tab_2d_1d( nbpac, zv_i_ac(1:nbpac,jl) , v_i(:,:,jl) , & 369 jpi, jpj, npac(1:nbpac) ) 370 CALL tab_2d_1d( nbpac, zoa_i_ac(1:nbpac,jl) , oa_i(:,:,jl) , & 371 jpi, jpj, npac(1:nbpac) ) 372 CALL tab_2d_1d( nbpac, zsmv_i_ac(1:nbpac,jl), smv_i(:,:,jl), & 373 jpi, jpj, npac(1:nbpac) ) 374 DO jk = 1, nlay_i 375 CALL tab_2d_1d( nbpac, ze_i_ac(1:nbpac,jk,jl), e_i(:,:,jk,jl) , & 376 jpi, jpj, npac(1:nbpac) ) 377 END DO ! jk 378 END DO ! jl 379 380 CALL tab_2d_1d( nbpac, qldif_1d (1:nbpac) , qldif , & 381 jpi, jpj, npac(1:nbpac) ) 382 CALL tab_2d_1d( nbpac, qcmif_1d (1:nbpac) , qcmif , & 383 jpi, jpj, npac(1:nbpac) ) 384 CALL tab_2d_1d( nbpac, t_bo_b (1:nbpac) , t_bo , & 385 jpi, jpj, npac(1:nbpac) ) 386 CALL tab_2d_1d( nbpac, fseqv_1d (1:nbpac) , fseqv , & 387 jpi, jpj, npac(1:nbpac) ) 388 CALL tab_2d_1d( nbpac, hicol_b (1:nbpac) , hicol , & 389 jpi, jpj, npac(1:nbpac) ) 390 CALL tab_2d_1d( nbpac, zvrel_ac (1:nbpac) , zvrel , & 391 jpi, jpj, npac(1:nbpac) ) 392 393 !------------------------------------------------------------------------------! 394 ! 5) Compute thickness, salinity, enthalpy, age, area and volume of new ice 395 !------------------------------------------------------------------------------! 396 397 !---------------------- 398 ! Thickness of new ice 399 !---------------------- 400 DO ji = 1, nbpac 401 zh_newice(ji) = hiccrit(1) 402 END DO 403 IF ( fraz_swi .EQ. 1.0 ) zh_newice(:) = hicol_b(:) 404 405 !---------------------- 406 ! Salinity of new ice 407 !---------------------- 408 409 IF ( num_sal .EQ. 1 ) THEN 410 zs_newice(:) = bulk_sal 411 ENDIF ! num_sal 412 413 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN 414 415 DO ji = 1, nbpac 416 zs_newice(ji) = MIN( 4.606 + 0.91 / zh_newice(ji) , s_i_max ) 417 zji = MOD( npac(ji) - 1, jpi ) + 1 418 zjj = ( npac(ji) - 1 ) / jpi + 1 419 zs_newice(ji) = MIN( 0.5*sss_m(zji,zjj) , zs_newice(ji) ) 420 END DO ! jl 421 422 ENDIF ! num_sal 423 424 IF ( num_sal .EQ. 3 ) THEN 425 zs_newice(:) = 2.3 426 ENDIF ! num_sal 427 428 !------------------------- 429 ! Heat content of new ice 430 !------------------------- 431 ! We assume that new ice is formed at the seawater freezing point 432 DO ji = 1, nbpac 433 ztmelts = - tmut * zs_newice(ji) + rtt ! Melting point (K) 434 ze_newice(ji) = rhoic * ( cpic * ( ztmelts - t_bo_b(ji) ) & 435 + lfus * ( 1.0 - ( ztmelts - rtt ) & 436 / ( t_bo_b(ji) - rtt ) ) & 437 - rcp * ( ztmelts-rtt ) ) 438 ze_newice(ji) = MAX( ze_newice(ji) , 0.0 ) + & 439 MAX( 0.0 , SIGN( 1.0 , - ze_newice(ji) ) ) & 440 * rhoic * lfus 441 END DO ! ji 442 !---------------- 443 ! Age of new ice 444 !---------------- 445 DO ji = 1, nbpac 446 zo_newice(ji) = 0.0 447 END DO ! ji 448 449 !-------------------------- 450 ! Open water energy budget 451 !-------------------------- 452 DO ji = 1, nbpac 453 zqbgow(ji) = qldif_1d(ji) - qcmif_1d(ji) !<0 454 END DO ! ji 455 456 !------------------- 457 ! Volume of new ice 458 !------------------- 459 DO ji = 1, nbpac 460 zv_newice(ji) = - zqbgow(ji) / ze_newice(ji) 461 462 ! A fraction zfrazb of frazil ice is accreted at the ice bottom 463 zfrazb = ( TANH ( Cfrazb * ( zvrel_ac(ji) - vfrazb ) ) & 464 + 1.0 ) / 2.0 * maxfrazb 465 zdh_frazb(ji) = zfrazb*zv_newice(ji) 466 zv_newice(ji) = ( 1.0 - zfrazb ) * zv_newice(ji) 467 END DO 468 469 !--------------------------------- 470 ! Salt flux due to new ice growth 471 !--------------------------------- 472 IF ( ( num_sal .EQ. 4 ) ) THEN 473 DO ji = 1, nbpac 474 zji = MOD( npac(ji) - 1, jpi ) + 1 475 zjj = ( npac(ji) - 1 ) / jpi + 1 476 fseqv_1d(ji) = fseqv_1d(ji) + & 477 ( sss_m(zji,zjj) - bulk_sal ) * rhoic * & 478 zv_newice(ji) / rdt_ice 479 END DO 480 ELSE 481 DO ji = 1, nbpac 482 zji = MOD( npac(ji) - 1, jpi ) + 1 483 zjj = ( npac(ji) - 1 ) / jpi + 1 484 fseqv_1d(ji) = fseqv_1d(ji) + & 485 ( sss_m(zji,zjj) - zs_newice(ji) ) * rhoic * & 486 zv_newice(ji) / rdt_ice 487 END DO ! ji 488 ENDIF 489 490 !------------------------------------ 491 ! Diags for energy conservation test 492 !------------------------------------ 493 DO ji = 1, nbpac 494 ! Volume 495 zji = MOD( npac(ji) - 1, jpi ) + 1 496 zjj = ( npac(ji) - 1 ) / jpi + 1 497 vt_i_init(zji,zjj) = vt_i_init(zji,zjj) + zv_newice(ji) 498 ! Energy 499 zde = ze_newice(ji) / unit_fac 500 zde = zde * area(zji,zjj) * zv_newice(ji) 501 et_i_init(zji,zjj) = et_i_init(zji,zjj) + zde 502 END DO 503 504 ! keep new ice volume in memory 505 CALL tab_1d_2d( nbpac, v_newice , npac(1:nbpac), zv_newice(1:nbpac) , & 506 jpi, jpj ) 507 508 !----------------- 509 ! Area of new ice 510 !----------------- 511 DO ji = 1, nbpac 512 za_newice(ji) = zv_newice(ji) / zh_newice(ji) 513 ! diagnostic 514 zji = MOD( npac(ji) - 1, jpi ) + 1 515 zjj = ( npac(ji) - 1 ) / jpi + 1 516 diag_lat_gr(zji,zjj) = zv_newice(ji) / rdt_ice 517 END DO !ji 518 519 !------------------------------------------------------------------------------! 520 ! 6) Redistribute new ice area and volume into ice categories ! 521 !------------------------------------------------------------------------------! 522 523 !----------------------------------------- 524 ! Keep old ice areas and volume in memory 525 !----------------------------------------- 526 zv_old(:,:) = zv_i_ac(:,:) 527 za_old(:,:) = za_i_ac(:,:) 528 529 !------------------------------------------- 530 ! Compute excessive new ice area and volume 531 !------------------------------------------- 532 ! If lateral ice growth gives an ice concentration gt 1, then 533 ! we keep the excessive volume in memory and attribute it later 534 ! to bottom accretion 535 DO ji = 1, nbpac 536 ! vectorize 537 IF ( za_newice(ji) .GT. ( 1.0 - zat_i_ac(ji) ) ) THEN 538 zda_res(ji) = za_newice(ji) - (1.0 - zat_i_ac(ji) ) 539 zdv_res(ji) = zda_res(ji) * zh_newice(ji) 540 za_newice(ji) = za_newice(ji) - zda_res(ji) 541 zv_newice(ji) = zv_newice(ji) - zdv_res(ji) 542 ELSE 543 zda_res(ji) = 0.0 544 zdv_res(ji) = 0.0 545 ENDIF 546 END DO ! ji 547 548 !------------------------------------------------ 549 ! Laterally redistribute new ice volume and area 550 !------------------------------------------------ 551 zat_i_ac(:) = 0.0 552 553 DO jl = 1, jpl 554 DO ji = 1, nbpac 555 ! vectorize 556 IF ( ( hi_max(jl-1) .LT. zh_newice(ji) ) & 557 .AND. ( zh_newice(ji) .LE. hi_max(jl) ) ) THEN 558 za_i_ac(ji,jl) = za_i_ac(ji,jl) + za_newice(ji) 559 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + zv_newice(ji) 560 zat_i_ac(ji) = zat_i_ac(ji) + za_i_ac(ji,jl) 561 zcatac(ji) = jl 562 ENDIF 563 END DO ! ji 564 END DO ! jl 565 566 !---------------------------------- 567 ! Heat content - lateral accretion 568 !---------------------------------- 569 DO ji = 1, nbpac 570 jl = zcatac(ji) ! categroy in which new ice is put 571 ! zindb = 0 if no ice and 1 if yes 572 zindb = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , -za_old(ji,jl) ) ) 573 ! old ice thickness 574 zhice_old(ji,jl) = zv_old(ji,jl) & 575 / MAX ( za_old(ji,jl) , zeps ) * zindb 576 ! difference in thickness 577 zdhex(ji) = MAX( 0.0, zh_newice(ji) - zhice_old(ji,jl) ) 578 ! is ice totally new in category jl ? 579 zswinew(ji) = MAX( 0.0, SIGN( 1.0 , - za_old(ji,jl) + epsi11 ) ) 580 END DO 581 582 DO jk = 1, nlay_i 583 DO ji = 1, nbpac 584 jl = zcatac(ji) 585 zqold = ze_i_ac(ji,jk,jl) ! [ J.m-3 ] 586 zalphai = MIN( zhice_old(ji,jl) * jk / nlay_i , & 587 zh_newice(ji) ) & 588 - MIN( zhice_old(ji,jl) * ( jk - 1 ) & 589 / nlay_i , zh_newice(ji) ) 590 ze_i_ac(ji,jk,jl) = & 591 zswinew(ji) * ze_newice(ji) & 592 + ( 1.0 - zswinew(ji) ) * & 593 ( za_old(ji,jl) * zqold * zhice_old(ji,jl) / nlay_i & 594 + za_newice(ji) * ze_newice(ji) * zalphai & 595 + za_newice(ji) * ze_newice(ji) * zdhex(ji) / nlay_i ) / & 596 ( ( zv_i_ac(ji,jl) ) / nlay_i ) 597 598 END DO !ji 599 END DO !jl 600 601 !----------------------------------------------- 602 ! Add excessive volume of new ice at the bottom 603 !----------------------------------------------- 604 ! If the ice concentration exceeds 1, the remaining volume of new ice 605 ! is equally redistributed among all ice categories in which there is 606 ! ice 607 608 ! Fraction of level ice 609 jm = 1 610 zat_i_lev(:) = 0.0 611 612 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 613 DO ji = 1, nbpac 614 zat_i_lev(ji) = zat_i_lev(ji) + za_i_ac(ji,jl) 615 END DO 616 END DO 617 618 IF( ln_nicep ) WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindx, 1:jpl) 619 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 620 DO ji = 1, nbpac 621 zindb = MAX( 0.0, SIGN( 1.0, zdv_res(ji) ) ) 622 zv_i_ac(ji,jl) = zv_i_ac(ji,jl) + & 623 zindb * zdv_res(ji) * za_i_ac(ji,jl) / & 624 MAX( zat_i_lev(ji) , epsi06 ) 625 END DO ! ji 626 END DO ! jl 627 IF( ln_nicep ) WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindx, 1:jpl) 628 629 !--------------------------------- 630 ! Heat content - bottom accretion 631 !--------------------------------- 632 jm = 1 633 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 634 DO ji = 1, nbpac 635 ! zindb = 0 if no ice and 1 if yes 636 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 & 637 , - za_i_ac(ji,jl ) ) ) 638 zhice_old(ji,jl) = zv_i_ac(ji,jl) / & 639 MAX( za_i_ac(ji,jl) , zeps ) * zindb 640 zdhicbot(ji,jl) = zdv_res(ji) / MAX( za_i_ac(ji,jl) , zeps ) & 641 * zindb & 642 + zindb * zdh_frazb(ji) ! frazil ice 643 ! may coalesce 644 ! thickness of residual ice 645 zdummy(ji,jl) = zv_i_ac(ji,jl)/MAX(za_i_ac(ji,jl),zeps)*zindb 646 END DO !ji 647 END DO !jl 648 649 ! old layers thicknesses and enthalpies 650 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 651 DO jk = 1, nlay_i 652 DO ji = 1, nbpac 653 zthick0(ji,jk,jl)= zhice_old(ji,jl) / nlay_i 654 zqm0 (ji,jk,jl)= ze_i_ac(ji,jk,jl) * zthick0(ji,jk,jl) 655 END DO !ji 656 END DO !jk 657 END DO !jl 658 659 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 660 DO ji = 1, nbpac 661 zthick0(ji,nlay_i+1,jl) = zdhicbot(ji,jl) 662 zqm0 (ji,nlay_i+1,jl) = ze_newice(ji)*zdhicbot(ji,jl) 663 END DO ! ji 664 END DO ! jl 665 666 ! Redistributing energy on the new grid 667 ze_i_ac(:,:,:) = 0.0 668 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 669 DO jk = 1, nlay_i 670 DO layer = 1, nlay_i + 1 671 DO ji = 1, nbpac 672 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , & 673 - za_i_ac(ji,jl ) ) ) 674 ! Redistributing energy on the new grid 675 zweight = MAX ( & 676 MIN( zhice_old(ji,jl) * layer , zdummy(ji,jl) * jk ) - & 677 MAX( zhice_old(ji,jl) * ( layer - 1 ) , zdummy(ji,jl) * & 678 678 ( jk - 1 ) ) , 0.0 ) & 679 / ( MAX(nlay_i * zthick0(ji,layer,jl),zeps) ) * zindb680 ze_i_ac(ji,jk,jl) = ze_i_ac(ji,jk,jl) + &681 682 END DO ! ji683 END DO ! layer684 END DO ! jk685 END DO ! jl686 687 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2)688 DO jk = 1, nlay_i689 DO ji = 1, nbpac690 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 &691 692 ze_i_ac(ji,jk,jl) = ze_i_ac(ji,jk,jl) / &693 694 695 END DO696 END DO697 END DO698 699 !------------700 ! Update age701 !------------702 DO jl = 1, jpl703 DO ji = 1, nbpac704 !--ice age705 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - &706 707 zoa_i_ac(ji,jl) = za_old(ji,jl) * zoa_i_ac(ji,jl) / &708 709 END DO ! ji710 END DO ! jl711 712 !-----------------713 ! Update salinity714 !-----------------715 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN716 717 DO jl = 1, jpl718 DO ji = 1, nbpac719 !zindb = 0 if no ice and 1 if yes720 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - &721 722 zdv = zv_i_ac(ji,jl) - zv_old(ji,jl)723 zsmv_i_ac(ji,jl) = ( zsmv_i_ac(ji,jl) + zdv * zs_newice(ji) ) * &724 725 END DO ! ji726 END DO ! jl727 728 ENDIF ! num_sal729 730 731 !------------------------------------------------------------------------------!732 ! 8) Change 2D vectors to 1D vectors733 !------------------------------------------------------------------------------!734 735 DO jl = 1, jpl736 CALL tab_1d_2d( nbpac, a_i(:,:,jl) , npac(1:nbpac) , &737 738 CALL tab_1d_2d( nbpac, v_i(:,:,jl) , npac(1:nbpac) , &739 740 CALL tab_1d_2d( nbpac, oa_i(:,:,jl), npac(1:nbpac) , &741 742 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) &743 CALL tab_1d_2d( nbpac, smv_i(:,:,jl) , npac(1:nbpac) , &744 745 DO jk = 1, nlay_i746 CALL tab_1d_2d( nbpac, e_i(:,:,jk,jl) , npac(1:nbpac), &747 748 END DO ! jk749 END DO !jl750 CALL tab_1d_2d( nbpac, fseqv , npac(1:nbpac), fseqv_1d (1:nbpac) , &751 752 753 ENDIF ! nbpac > 0754 755 !------------------------------------------------------------------------------!756 ! 9) Change units for e_i757 !------------------------------------------------------------------------------!679 / ( MAX(nlay_i * zthick0(ji,layer,jl),zeps) ) * zindb 680 ze_i_ac(ji,jk,jl) = ze_i_ac(ji,jk,jl) + & 681 zweight * zqm0(ji,layer,jl) 682 END DO ! ji 683 END DO ! layer 684 END DO ! jk 685 END DO ! jl 686 687 DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 688 DO jk = 1, nlay_i 689 DO ji = 1, nbpac 690 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 & 691 , - zv_i_ac(ji,jl) ) ) !0 if no ice 692 ze_i_ac(ji,jk,jl) = ze_i_ac(ji,jk,jl) / & 693 MAX( zv_i_ac(ji,jl) , zeps) & 694 * za_i_ac(ji,jl) * nlay_i * zindb 695 END DO 696 END DO 697 END DO 698 699 !------------ 700 ! Update age 701 !------------ 702 DO jl = 1, jpl 703 DO ji = 1, nbpac 704 !--ice age 705 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - & 706 za_i_ac(ji,jl) ) ) ! 0 if no ice and 1 if yes 707 zoa_i_ac(ji,jl) = za_old(ji,jl) * zoa_i_ac(ji,jl) / & 708 MAX( za_i_ac(ji,jl) , zeps ) * zindb 709 END DO ! ji 710 END DO ! jl 711 712 !----------------- 713 ! Update salinity 714 !----------------- 715 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) THEN 716 717 DO jl = 1, jpl 718 DO ji = 1, nbpac 719 !zindb = 0 if no ice and 1 if yes 720 zindb = 1.0 - MAX( 0.0 , SIGN( 1.0 , - & 721 zv_i_ac(ji,jl) ) ) ! 0 if no ice and 1 if yes 722 zdv = zv_i_ac(ji,jl) - zv_old(ji,jl) 723 zsmv_i_ac(ji,jl) = ( zsmv_i_ac(ji,jl) + zdv * zs_newice(ji) ) * & 724 zindb 725 END DO ! ji 726 END DO ! jl 727 728 ENDIF ! num_sal 729 730 731 !------------------------------------------------------------------------------! 732 ! 8) Change 2D vectors to 1D vectors 733 !------------------------------------------------------------------------------! 734 735 DO jl = 1, jpl 736 CALL tab_1d_2d( nbpac, a_i(:,:,jl) , npac(1:nbpac) , & 737 za_i_ac(1:nbpac,jl) , jpi, jpj ) 738 CALL tab_1d_2d( nbpac, v_i(:,:,jl) , npac(1:nbpac) , & 739 zv_i_ac(1:nbpac,jl) , jpi, jpj ) 740 CALL tab_1d_2d( nbpac, oa_i(:,:,jl), npac(1:nbpac) , & 741 zoa_i_ac(1:nbpac,jl), jpi, jpj ) 742 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) ) & 743 CALL tab_1d_2d( nbpac, smv_i(:,:,jl) , npac(1:nbpac) , & 744 zsmv_i_ac(1:nbpac,jl) , jpi, jpj ) 745 DO jk = 1, nlay_i 746 CALL tab_1d_2d( nbpac, e_i(:,:,jk,jl) , npac(1:nbpac), & 747 ze_i_ac(1:nbpac,jk,jl), jpi, jpj ) 748 END DO ! jk 749 END DO !jl 750 CALL tab_1d_2d( nbpac, fseqv , npac(1:nbpac), fseqv_1d (1:nbpac) , & 751 jpi, jpj ) 752 753 ENDIF ! nbpac > 0 754 755 !------------------------------------------------------------------------------! 756 ! 9) Change units for e_i 757 !------------------------------------------------------------------------------! 758 758 759 759 DO jl = 1, jpl … … 767 767 ! of layers to get heat content in 10^9 Joules 768 768 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * & 769 770 769 area(ji,jj) * v_i(ji,jj,jl) / & 770 nlay_i 771 771 END DO 772 772 END DO 773 773 END DO 774 774 END DO 775 776 !------------------------------------------------------------------------------|777 ! 10) Conservation check and changes in each ice category778 !------------------------------------------------------------------------------|775 776 !------------------------------------------------------------------------------| 777 ! 10) Conservation check and changes in each ice category 778 !------------------------------------------------------------------------------| 779 779 780 780 IF ( con_i ) THEN 781 CALL lim_column_sum (jpl, v_i, vt_i_final) 782 fieldid = 'v_i, limthd_lac' 783 CALL lim_cons_check (vt_i_init, vt_i_final, 1.0e-6, fieldid) 784 785 CALL lim_column_sum_energy(jpl, nlay_i, e_i, et_i_final) 786 fieldid = 'e_i, limthd_lac' 787 CALL lim_cons_check (et_i_final, et_i_final, 1.0e-3, fieldid) 788 789 CALL lim_column_sum (jpl, v_s, vt_s_final) 790 fieldid = 'v_s, limthd_lac' 791 CALL lim_cons_check (vt_s_init, vt_s_final, 1.0e-6, fieldid) 792 793 ! CALL lim_column_sum (jpl, e_s(:,:,1,:) , et_s_init) 794 ! fieldid = 'e_s, limthd_lac' 795 ! CALL lim_cons_check (et_s_init, et_s_final, 1.0e-3, fieldid) 796 797 WRITE(numout,*) ' vt_i_init : ', vt_i_init(jiindx,jjindx) 798 WRITE(numout,*) ' vt_i_final: ', vt_i_final(jiindx,jjindx) 799 WRITE(numout,*) ' et_i_init : ', et_i_init(jiindx,jjindx) 800 WRITE(numout,*) ' et_i_final: ', et_i_final(jiindx,jjindx) 781 CALL lim_column_sum (jpl, v_i, vt_i_final) 782 fieldid = 'v_i, limthd_lac' 783 CALL lim_cons_check (vt_i_init, vt_i_final, 1.0e-6, fieldid) 784 785 CALL lim_column_sum_energy(jpl, nlay_i, e_i, et_i_final) 786 fieldid = 'e_i, limthd_lac' 787 CALL lim_cons_check (et_i_final, et_i_final, 1.0e-3, fieldid) 788 789 CALL lim_column_sum (jpl, v_s, vt_s_final) 790 fieldid = 'v_s, limthd_lac' 791 CALL lim_cons_check (vt_s_init, vt_s_final, 1.0e-6, fieldid) 792 793 ! CALL lim_column_sum (jpl, e_s(:,:,1,:) , et_s_init) 794 ! fieldid = 'e_s, limthd_lac' 795 ! CALL lim_cons_check (et_s_init, et_s_final, 1.0e-3, fieldid) 796 797 IF( ln_nicep ) THEN 798 WRITE(numout,*) ' vt_i_init : ', vt_i_init(jiindx,jjindx) 799 WRITE(numout,*) ' vt_i_final: ', vt_i_final(jiindx,jjindx) 800 WRITE(numout,*) ' et_i_init : ', et_i_init(jiindx,jjindx) 801 WRITE(numout,*) ' et_i_final: ', et_i_final(jiindx,jjindx) 802 ENDIF 801 803 802 804 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
