Changeset 7753 for trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90
- Timestamp:
- 2017-03-03T12:46:59+01:00 (7 years ago)
- File:
-
- 1 edited
-
trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90 (modified) (25 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90
r7698 r7753 316 316 #if defined key_cice 317 317 IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN 318 !$OMP PARALLEL DO schedule(static) private(jj, ji) 319 DO jj = 1, jpj 320 DO ji = 1, jpi 321 qlw_ice(ji,jj,1) = sf(jp_qlw)%fnow(ji,jj,1) 322 END DO 323 END DO 324 IF( ln_dm2dc ) THEN ; qsr_ice(:,:,1) = sbc_dcy( sf(jp_qsr)%fnow(:,:,1) ) 325 ELSE 326 !$OMP PARALLEL DO schedule(static) private(jj, ji) 327 DO jj = 1, jpj 328 DO ji = 1, jpi 329 qsr_ice(ji,jj,1) = sf(jp_qsr)%fnow(ji,jj,1) 330 END DO 331 END DO 318 qlw_ice(:,:,1) = sf(jp_qlw )%fnow(:,:,1) 319 IF( ln_dm2dc ) THEN ; qsr_ice(:,:,1) = sbc_dcy( sf(jp_qsr)%fnow(:,:,1) ) 320 ELSE ; qsr_ice(:,:,1) = sf(jp_qsr)%fnow(:,:,1) 332 321 ENDIF 333 !$OMP PARALLEL DO schedule(static) private(jj, ji) 334 DO jj = 1, jpj 335 DO ji = 1, jpi 336 tatm_ice(ji,jj) = sf(jp_tair)%fnow(ji,jj,1) 337 qatm_ice(ji,jj) = sf(jp_humi)%fnow(ji,jj,1) 338 tprecip(ji,jj) = sf(jp_prec)%fnow(ji,jj,1) * rn_pfac 339 sprecip(ji,jj) = sf(jp_snow)%fnow(ji,jj,1) * rn_pfac 340 wndi_ice(ji,jj) = sf(jp_wndi)%fnow(ji,jj,1) 341 wndj_ice(ji,jj) = sf(jp_wndj)%fnow(ji,jj,1) 342 END DO 343 END DO 322 tatm_ice(:,:) = sf(jp_tair)%fnow(:,:,1) 323 qatm_ice(:,:) = sf(jp_humi)%fnow(:,:,1) 324 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac 325 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac 326 wndi_ice(:,:) = sf(jp_wndi)%fnow(:,:,1) 327 wndj_ice(:,:) = sf(jp_wndj)%fnow(:,:,1) 344 328 ENDIF 345 329 #endif … … 398 382 ! 399 383 400 !$OMP PARALLEL DO schedule(static) private(jj, ji) 401 DO jj = 1, jpj 402 DO ji = 1, jpi 403 ! local scalars ( place there for vector optimisation purposes) 404 zst(ji,jj) = pst(ji,jj) + rt0 ! convert SST from Celcius to Kelvin (and set minimum value far above 0 K) 405 406 ! ... components ( U10m - U_oce ) at T-point (unmasked) 407 !!gm move zwnd_i (_j) set to zero inside the key_cyclone ??? 408 zwnd_i(ji,jj) = 0._wp 409 zwnd_j(ji,jj) = 0._wp 410 END DO 411 END DO 384 ! local scalars ( place there for vector optimisation purposes) 385 zst(:,:) = pst(:,:) + rt0 ! convert SST from Celcius to Kelvin (and set minimum value far above 0 K) 386 412 387 ! ----------------------------------------------------------------------------- ! 413 388 ! 0 Wind components and module at T-point relative to the moving ocean ! 414 389 ! ----------------------------------------------------------------------------- ! 415 390 391 ! ... components ( U10m - U_oce ) at T-point (unmasked) 392 !!gm move zwnd_i (_j) set to zero inside the key_cyclone ??? 393 zwnd_i(:,:) = 0._wp 394 zwnd_j(:,:) = 0._wp 416 395 #if defined key_cyclone 417 396 CALL wnd_cyc( kt, zwnd_i, zwnd_j ) ! add analytical tropical cyclone (Vincent et al. JGR 2012) 418 !$OMP PARALLEL DO schedule(static) private(jj, ji)419 397 DO jj = 2, jpjm1 420 398 DO ji = fs_2, fs_jpim1 ! vect. opt. … … 424 402 END DO 425 403 #endif 426 !$OMP PARALLEL DO schedule(static) private(jj, ji)427 404 DO jj = 2, jpjm1 428 405 DO ji = fs_2, fs_jpim1 ! vect. opt. … … 434 411 CALL lbc_lnk( zwnd_j(:,:) , 'T', -1. ) 435 412 ! ... scalar wind ( = | U10m - U_oce | ) at T-point (masked) 436 !$OMP PARALLEL DO schedule(static) private(jj, ji) 437 DO jj = 1, jpj 438 DO ji = 1, jpi 439 wndm(ji,jj) = SQRT( zwnd_i(ji,jj) * zwnd_i(ji,jj) & 440 & + zwnd_j(ji,jj) * zwnd_j(ji,jj) ) * tmask(ji,jj,1) 441 442 END DO 443 END DO 413 wndm(:,:) = SQRT( zwnd_i(:,:) * zwnd_i(:,:) & 414 & + zwnd_j(:,:) * zwnd_j(:,:) ) * tmask(:,:,1) 415 444 416 ! ----------------------------------------------------------------------------- ! 445 417 ! I Radiative FLUXES ! … … 449 421 zztmp = 1. - albo 450 422 IF( ln_dm2dc ) THEN ; qsr(:,:) = zztmp * sbc_dcy( sf(jp_qsr)%fnow(:,:,1) ) * tmask(:,:,1) 451 ELSE 452 !$OMP PARALLEL DO schedule(static) private(jj, ji) 453 DO jj = 1, jpj 454 DO ji = 1, jpi 455 qsr(ji,jj) = zztmp * sf(jp_qsr)%fnow(ji,jj,1) * tmask(ji,jj,1) 456 END DO 457 END DO 458 ENDIF 459 460 !$OMP PARALLEL DO schedule(static) private(jj, ji) 461 DO jj = 1, jpj 462 DO ji = 1, jpi 463 zqlw(ji,jj) = ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * zst(ji,jj)*zst(ji,jj)*zst(ji,jj)*zst(ji,jj) ) * tmask(ji,jj,1) ! Long Wave 464 END DO 465 END DO 423 ELSE ; qsr(:,:) = zztmp * sf(jp_qsr)%fnow(:,:,1) * tmask(:,:,1) 424 ENDIF 425 426 zqlw(:,:) = ( sf(jp_qlw)%fnow(:,:,1) - Stef * zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:) ) * tmask(:,:,1) ! Long Wave 466 427 467 428 … … 500 461 END IF 501 462 502 !$OMP PARALLEL 503 !$OMP DO schedule(static) private(jj, ji) 504 DO jj = 1, jpj 505 DO ji = 1, jpi 506 Cd_oce(ji,jj) = Cd(ji,jj) ! record value of pure ocean-atm. drag (clem) 507 END DO 508 END DO 509 510 !$OMP DO schedule(static) private(jj, ji) 463 Cd_oce(:,:) = Cd(:,:) ! record value of pure ocean-atm. drag (clem) 464 511 465 DO jj = 1, jpj ! tau module, i and j component 512 466 DO ji = 1, jpi … … 517 471 END DO 518 472 END DO 519 !$OMP END PARALLEL520 473 521 474 ! ! add the HF tau contribution to the wind stress module 522 IF( lhftau ) THEN 523 !$OMP PARALLEL DO schedule(static) private(jj, ji) 524 DO jj = 1, jpj 525 DO ji = 1, jpi 526 taum(ji,jj) = taum(ji,jj) + sf(jp_tdif)%fnow(ji,jj,1) 527 END DO 528 END DO 529 END IF 475 IF( lhftau ) taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:,1) 530 476 531 477 CALL iom_put( "taum_oce", taum ) ! output wind stress module … … 534 480 ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines 535 481 ! Note the use of MAX(tmask(i,j),tmask(i+1,j) is to mask tau over ice shelves 536 !$OMP PARALLEL DO schedule(static) private(jj, ji)537 482 DO jj = 1, jpjm1 538 483 DO ji = 1, fs_jpim1 … … 551 496 552 497 ! zqla used as temporary array, for rho*U (common term of bulk formulae): 553 !$OMP PARALLEL DO schedule(static) private(jj, ji) 554 DO jj = 1, jpj 555 DO ji = 1, jpi 556 zqla(ji,jj) = zrhoa(ji,jj) * zU_zu(ji,jj) 557 END DO 558 END DO 498 zqla(:,:) = zrhoa(:,:) * zU_zu(:,:) 559 499 560 500 IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN 561 501 !! q_air and t_air are given at 10m (wind reference height) 562 !$OMP PARALLEL DO schedule(static) private(jj, ji) 563 DO jj = 1, jpj 564 DO ji = 1, jpi 565 zevap(ji,jj) = rn_efac*MAX( 0._wp, zqla(ji,jj)*Ce(ji,jj)*(zsq(ji,jj) - sf(jp_humi)%fnow(ji,jj,1)) ) ! Evaporation, using bulk wind speed 566 END DO 567 END DO 568 zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch(:,:)*(zst(:,:) - ztpot(:,:) ) ! Sensible Heat, using bulk wind speed 502 zevap(:,:) = rn_efac*MAX( 0._wp, zqla(:,:)*Ce(:,:)*(zsq(:,:) - sf(jp_humi)%fnow(:,:,1)) ) ! Evaporation, using bulk wind speed 503 zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch(:,:)*(zst(:,:) - ztpot(:,:) ) ! Sensible Heat, using bulk wind speed 569 504 ELSE 570 505 !! q_air and t_air are not given at 10m (wind reference height) 571 506 ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!! 572 !$OMP PARALLEL DO schedule(static) private(jj, ji) 573 DO jj = 1, jpj 574 DO ji = 1, jpi 575 zevap(ji,jj) = rn_efac*MAX( 0._wp, zqla(ji,jj)*Ce(ji,jj)*(zsq(ji,jj) - zq_zu(ji,jj) ) ) ! Evaporation ! using bulk wind speed 576 END DO 577 END DO 507 zevap(:,:) = rn_efac*MAX( 0._wp, zqla(:,:)*Ce(:,:)*(zsq(:,:) - zq_zu(:,:) ) ) ! Evaporation ! using bulk wind speed 578 508 zqsb (:,:) = cp_air(sf(jp_humi)%fnow(:,:,1))*zqla(:,:)*Ch(:,:)*(zst(:,:) - zt_zu(:,:) ) ! Sensible Heat ! using bulk wind speed 579 509 ENDIF … … 597 527 ! ----------------------------------------------------------------------------- ! 598 528 ! 599 !$OMP PARALLEL DO schedule(static) private(jj, ji) 600 DO jj = 1, jpj 601 DO ji = 1, jpi 602 emp (ji,jj) = ( zevap(ji,jj) & ! mass flux (evap. - precip.) 603 & - sf(jp_prec)%fnow(ji,jj,1) * rn_pfac ) * tmask(ji,jj,1) 604 ! 605 qns(ji,jj) = zqlw(ji,jj) - zqsb(ji,jj) - zqla(ji,jj) & ! Downward Non Solar 606 & - sf(jp_snow)%fnow(ji,jj,1) * rn_pfac * lfus & ! remove latent melting heat for solid precip 607 & - zevap(ji,jj) * pst(ji,jj) * rcp & ! remove evap heat content at SST 608 & + ( sf(jp_prec)%fnow(ji,jj,1) - sf(jp_snow)%fnow(ji,jj,1) ) * rn_pfac & ! add liquid precip heat content at Tair 609 & * ( sf(jp_tair)%fnow(ji,jj,1) - rt0 ) * rcp & 610 & + sf(jp_snow)%fnow(ji,jj,1) * rn_pfac & ! add solid precip heat content at min(Tair,Tsnow) 611 & * ( MIN( sf(jp_tair)%fnow(ji,jj,1), rt0_snow ) - rt0 ) * cpic * tmask(ji,jj,1) 612 ! 529 emp (:,:) = ( zevap(:,:) & ! mass flux (evap. - precip.) 530 & - sf(jp_prec)%fnow(:,:,1) * rn_pfac ) * tmask(:,:,1) 531 ! 532 qns(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) & ! Downward Non Solar 533 & - sf(jp_snow)%fnow(:,:,1) * rn_pfac * lfus & ! remove latent melting heat for solid precip 534 & - zevap(:,:) * pst(:,:) * rcp & ! remove evap heat content at SST 535 & + ( sf(jp_prec)%fnow(:,:,1) - sf(jp_snow)%fnow(:,:,1) ) * rn_pfac & ! add liquid precip heat content at Tair 536 & * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & 537 & + sf(jp_snow)%fnow(:,:,1) * rn_pfac & ! add solid precip heat content at min(Tair,Tsnow) 538 & * ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) 539 ! 613 540 #if defined key_lim3 614 qns_oce(ji,jj) = zqlw(ji,jj) - zqsb(ji,jj) - zqla(ji,jj) ! non solar without emp (only needed by LIM3)615 qsr_oce(ji,jj) = qsr(ji,jj)541 qns_oce(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) ! non solar without emp (only needed by LIM3) 542 qsr_oce(:,:) = qsr(:,:) 616 543 #endif 617 END DO618 END DO619 544 ! 620 545 IF ( nn_ice == 0 ) THEN … … 626 551 CALL iom_put( "qsr_oce" , qsr ) ! output downward solar heat over the ocean 627 552 CALL iom_put( "qt_oce" , qns+qsr ) ! output total downward heat over the ocean 628 !$OMP PARALLEL DO schedule(static) private(jj, ji) 629 DO jj = 1, jpj 630 DO ji = 1, jpi 631 tprecip(ji,jj) = sf(jp_prec)%fnow(ji,jj,1) * rn_pfac ! output total precipitation [kg/m2/s] 632 sprecip(ji,jj) = sf(jp_snow)%fnow(ji,jj,1) * rn_pfac ! output solid precipitation [kg/m2/s] 633 END DO 634 END DO 553 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac ! output total precipitation [kg/m2/s] 554 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! output solid precipitation [kg/m2/s] 635 555 CALL iom_put( 'snowpre', sprecip * 86400. ) ! Snow 636 556 CALL iom_put( 'precip' , tprecip * 86400. ) ! Total precipitation … … 679 599 CALL wrk_alloc( jpi,jpj, Cd ) 680 600 681 !$OMP PARALLEL DO schedule(static) private(jj, ji) 682 DO jj = 1, jpj 683 DO ji = 1, jpi 684 Cd(ji,jj) = Cd_ice 685 END DO 686 END DO 601 Cd(:,:) = Cd_ice 687 602 688 603 ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) … … 698 613 zrhoa (:,:) = rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 699 614 700 !$OMP PARALLEL DO schedule(static) private(jj, ji) 701 DO jj = 1, jpj 702 DO ji = 1, jpi 703 !!gm brutal.... 704 utau_ice (ji,jj) = 0._wp 705 vtau_ice (ji,jj) = 0._wp 706 wndm_ice (ji,jj) = 0._wp 707 !!gm end 708 END DO 709 END DO 615 !!gm brutal.... 616 utau_ice (:,:) = 0._wp 617 vtau_ice (:,:) = 0._wp 618 wndm_ice (:,:) = 0._wp 619 !!gm end 710 620 711 621 ! ----------------------------------------------------------------------------- ! … … 715 625 CASE( 'I' ) ! B-grid ice dynamics : I-point (i.e. F-point with sea-ice indexation) 716 626 ! and scalar wind at T-point ( = | U10m - U_ice | ) (masked) 717 !$OMP PARALLEL DO schedule(static) private(jj,ji,zwndi_f,zwndj_f,zwnorm_f,zwndi_t,zwndj_t)718 627 DO jj = 2, jpjm1 719 628 DO ji = 2, jpim1 ! B grid : NO vector opt … … 740 649 ! 741 650 CASE( 'C' ) ! C-grid ice dynamics : U & V-points (same as ocean) 742 !$OMP PARALLEL DO schedule(static) private(jj,ji,zwndi_t,zwndj_t)743 651 DO jj = 2, jpj 744 652 DO ji = fs_2, jpi ! vect. opt. … … 748 656 END DO 749 657 END DO 750 !$OMP PARALLEL DO schedule(static) private(jj,ji)751 658 DO jj = 2, jpjm1 752 659 DO ji = fs_2, fs_jpim1 ! vect. opt. … … 793 700 REAL(wp) :: zztmp, z1_lsub ! - - 794 701 REAL(wp), DIMENSION(:,:,:), POINTER :: z_qlw ! long wave heat flux over ice 795 REAL(wp), DIMENSION(:,:,:), POINTER :: zevap_ice3d, zqns_ice3d, zqsr_ice3d796 702 REAL(wp), DIMENSION(:,:,:), POINTER :: z_qsb ! sensible heat flux over ice 797 703 REAL(wp), DIMENSION(:,:,:), POINTER :: z_dqlw ! long wave heat sensitivity over ice 798 704 REAL(wp), DIMENSION(:,:,:), POINTER :: z_dqsb ! sensible heat sensitivity over ice 799 705 REAL(wp), DIMENSION(:,:) , POINTER :: zevap, zsnw ! evaporation and snw distribution after wind blowing (LIM3) 800 REAL(wp), DIMENSION(:,:) , POINTER :: zevap_ice2d, zqns_ice2d, zqsr_ice2d801 706 REAL(wp), DIMENSION(:,:) , POINTER :: zrhoa 802 707 REAL(wp), DIMENSION(:,:) , POINTER :: Cd ! transfer coefficient for momentum (tau) … … 805 710 IF( nn_timing == 1 ) CALL timing_start('blk_ice_flx') 806 711 ! 807 CALL wrk_alloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb , zevap_ice3d, zqns_ice3d, zqsr_ice3d)808 CALL wrk_alloc( jpi,jpj, zrhoa , zevap_ice2d, zqns_ice2d, zqsr_ice2d)712 CALL wrk_alloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 713 CALL wrk_alloc( jpi,jpj, zrhoa) 809 714 CALL wrk_alloc( jpi,jpj, Cd ) 810 715 811 !$OMP PARALLEL DO schedule(static) private(jj, ji) 812 DO jj = 1, jpj 813 DO ji = 1, jpi 814 Cd(ji,jj) = Cd_ice 815 END DO 816 END DO 716 Cd(:,:) = Cd_ice 817 717 818 718 ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) … … 831 731 ! 832 732 zztmp = 1. / ( 1. - albo ) 833 !$OMP PARALLEL 834 !$OMP DO schedule(static) private(jl,jj,ji,zst2,zst3) ! ========================== ! 835 DO jl = 1, jpl ! Loop over ice categories ! 836 ! ! ========================== ! 733 ! ! ========================== ! 734 DO jl = 1, jpl ! Loop over ice categories ! 735 ! ! ========================== ! 837 736 DO jj = 1 , jpj 838 737 DO ji = 1, jpi … … 882 781 END DO 883 782 ! 884 !$OMP DO schedule(static) private(jj, ji) 885 DO jj = 1, jpj 886 DO ji = 1, jpi 887 tprecip(ji,jj) = sf(jp_prec)%fnow(ji,jj,1) * rn_pfac ! total precipitation [kg/m2/s] 888 sprecip(ji,jj) = sf(jp_snow)%fnow(ji,jj,1) * rn_pfac ! solid precipitation [kg/m2/s] 889 END DO 890 END DO 891 !$OMP END PARALLEL 783 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac ! total precipitation [kg/m2/s] 784 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! solid precipitation [kg/m2/s] 892 785 CALL iom_put( 'snowpre', sprecip * 86400. ) ! Snow precipitation 893 786 CALL iom_put( 'precip' , tprecip * 86400. ) ! Total precipitation … … 898 791 ! --- evaporation --- ! 899 792 z1_lsub = 1._wp / Lsub 900 !$OMP PARALLEL 901 !$OMP DO schedule(static) private(jl,jj,ji) 793 evap_ice (:,:,:) = rn_efac * qla_ice (:,:,:) * z1_lsub ! sublimation 794 devap_ice(:,:,:) = rn_efac * dqla_ice(:,:,:) * z1_lsub ! d(sublimation)/dT 795 zevap (:,:) = rn_efac * ( emp(:,:) + tprecip(:,:) ) ! evaporation over ocean 796 797 ! --- evaporation minus precipitation --- ! 798 zsnw(:,:) = 0._wp 799 CALL lim_thd_snwblow( pfrld, zsnw ) ! snow distribution over ice after wind blowing 800 emp_oce(:,:) = pfrld(:,:) * zevap(:,:) - ( tprecip(:,:) - sprecip(:,:) ) - sprecip(:,:) * (1._wp - zsnw ) 801 emp_ice(:,:) = SUM( a_i_b(:,:,:) * evap_ice(:,:,:), dim=3 ) - sprecip(:,:) * zsnw 802 emp_tot(:,:) = emp_oce(:,:) + emp_ice(:,:) 803 804 ! --- heat flux associated with emp --- ! 805 qemp_oce(:,:) = - pfrld(:,:) * zevap(:,:) * sst_m(:,:) * rcp & ! evap at sst 806 & + ( tprecip(:,:) - sprecip(:,:) ) * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & ! liquid precip at Tair 807 & + sprecip(:,:) * ( 1._wp - zsnw ) * & ! solid precip at min(Tair,Tsnow) 808 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 809 qemp_ice(:,:) = sprecip(:,:) * zsnw * & ! solid precip (only) 810 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 811 812 ! --- total solar and non solar fluxes --- ! 813 qns_tot(:,:) = pfrld(:,:) * qns_oce(:,:) + SUM( a_i_b(:,:,:) * qns_ice(:,:,:), dim=3 ) + qemp_ice(:,:) + qemp_oce(:,:) 814 qsr_tot(:,:) = pfrld(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 ) 815 816 ! --- heat content of precip over ice in J/m3 (to be used in 1D-thermo) --- ! 817 qprec_ice(:,:) = rhosn * ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 818 819 ! --- heat content of evap over ice in W/m2 (to be used in 1D-thermo) --- 902 820 DO jl = 1, jpl 903 DO jj = 1 , jpj 904 DO ji = 1, jpi 905 evap_ice (ji,jj,jl) = rn_efac * qla_ice (ji,jj,jl) * z1_lsub ! sublimation 906 devap_ice(ji,jj,jl) = rn_efac * dqla_ice(ji,jj,jl) * z1_lsub ! d(sublimation)/dT 907 END DO 908 END DO 821 qevap_ice(:,:,jl) = 0._wp ! should be -evap_ice(:,:,jl)*( ( Tice - rt0 ) * cpic * tmask(:,:,1) ) 822 ! But we do not have Tice => consider it at 0degC => evap=0 909 823 END DO 910 !911 !$OMP DO schedule(static) private(jj, ji)912 DO jj = 1, jpj913 DO ji = 1, jpi914 zevap (ji,jj) = rn_efac * ( emp(ji,jj) + tprecip(ji,jj) ) ! evaporation over ocean915 916 ! --- evaporation minus precipitation --- !917 zsnw(ji,jj) = 0._wp918 END DO919 END DO920 !$OMP END PARALLEL921 CALL lim_thd_snwblow( pfrld, zsnw ) ! snow distribution over ice after wind blowing922 !$OMP PARALLEL923 !$OMP DO schedule(static) private(jj,ji)924 DO jj = 1, jpj925 DO ji = 1, jpi926 emp_oce(ji,jj) = pfrld(ji,jj) * zevap(ji,jj) - ( tprecip(ji,jj) - sprecip(ji,jj) ) - sprecip(ji,jj) * (1._wp - zsnw(ji,jj))927 END DO928 END DO929 !$OMP END DO NOWAIT930 !$OMP DO schedule(static) private(jl,jj,ji)931 DO jl = 1, jpl932 DO jj = 1 , jpj933 DO ji = 1, jpi934 zevap_ice3d(ji,jj,jl) = a_i_b(ji,jj,jl) * evap_ice(ji,jj,jl)935 zqns_ice3d(ji,jj,jl) = a_i_b(ji,jj,jl) * qns_ice(ji,jj,jl)936 zqsr_ice3d(ji,jj,jl) = a_i_b(ji,jj,jl) * qsr_ice(ji,jj,jl)937 END DO938 END DO939 END DO940 !$OMP END DO NOWAIT941 !$OMP DO schedule(static) private(jj,ji)942 DO jj = 1, jpj943 DO ji = 1, jpi944 zevap_ice2d(ji,jj) = 0._wp945 zqns_ice2d(ji,jj) = 0._wp946 zqsr_ice2d(ji,jj) = 0._wp947 END DO948 END DO949 DO jl = 1, jpl950 !$OMP DO schedule(static) private(jj,ji)951 DO jj = 1 , jpj952 DO ji = 1, jpi953 zevap_ice2d(ji,jj) = zevap_ice2d(ji,jj) + zevap_ice3d(ji,jj,jl)954 zqns_ice2d(ji,jj) = zqns_ice2d(ji,jj) + zqns_ice3d(ji,jj,jl)955 zqsr_ice2d(ji,jj) = zqsr_ice2d(ji,jj) + zqsr_ice3d(ji,jj,jl)956 END DO957 END DO958 END DO959 !$OMP DO schedule(static) private(jj,ji)960 DO jj = 1 , jpj961 DO ji = 1, jpi962 emp_ice(ji,jj) = zevap_ice2d(ji,jj) - sprecip(ji,jj) * zsnw(ji,jj)963 emp_tot(ji,jj) = emp_oce(ji,jj) + emp_ice(ji,jj)964 965 ! --- heat flux associated with emp --- !966 qemp_oce(ji,jj) = - pfrld(ji,jj) * zevap(ji,jj) * sst_m(ji,jj) * rcp & ! evap at sst967 & + ( tprecip(ji,jj) - sprecip(ji,jj) ) * ( sf(jp_tair)%fnow(ji,jj,1) - rt0 ) * rcp & ! liquid precip at Tair968 & + sprecip(ji,jj) * ( 1._wp - zsnw(ji,jj) ) * & ! solid precip at min(Tair,Tsnow)969 & ( ( MIN( sf(jp_tair)%fnow(ji,jj,1), rt0_snow ) - rt0 ) * cpic * tmask(ji,jj,1) - lfus )970 qemp_ice(ji,jj) = sprecip(ji,jj) * zsnw(ji,jj) * & ! solid precip (only)971 & ( ( MIN( sf(jp_tair)%fnow(ji,jj,1), rt0_snow ) - rt0 ) * cpic * tmask(ji,jj,1) - lfus )972 973 ! --- total solar and non solar fluxes --- !974 qns_tot(ji,jj) = pfrld(ji,jj) * qns_oce(ji,jj) + zqns_ice2d(ji,jj) + qemp_ice(ji,jj) + qemp_oce(ji,jj)975 qsr_tot(ji,jj) = pfrld(ji,jj) * qsr_oce(ji,jj) + zqsr_ice2d(ji,jj)976 977 ! --- heat content of precip over ice in J/m3 (to be used in 1D-thermo) --- !978 qprec_ice(ji,jj) = rhosn * ( ( MIN( sf(jp_tair)%fnow(ji,jj,1), rt0_snow ) - rt0 ) * cpic * tmask(ji,jj,1) - lfus )979 END DO980 END DO981 !$OMP END DO NOWAIT982 983 ! --- heat content of evap over ice in W/m2 (to be used in 1D-thermo) ---984 !$OMP DO schedule(static) private(jl,jj,ji)985 DO jl = 1, jpl986 DO jj = 1, jpj987 DO ji = 1, jpi988 qevap_ice(ji,jj,jl) = 0._wp ! should be -evap_ice(:,:,jl)*( ( Tice - rt0 ) * cpic * tmask(:,:,1) )989 ! But we do not have Tice => consider it at 0degC => evap=0990 END DO991 END DO992 END DO993 !$OMP END PARALLEL994 824 995 825 CALL wrk_dealloc( jpi,jpj, zevap, zsnw ) … … 1001 831 ! ( Maykut and Untersteiner, 1971 ; Ebert and Curry, 1993 ) 1002 832 ! 1003 !$OMP PARALLEL DO schedule(static) private(jj,ji) 1004 DO jj = 1, jpj 1005 DO ji = 1, jpi 1006 fr1_i0(ji,jj) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 1007 fr2_i0(ji,jj) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 1008 END DO 1009 END DO 833 fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 834 fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 1010 835 ! 1011 836 ! … … 1019 844 ENDIF 1020 845 1021 CALL wrk_dealloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb , zevap_ice3d, zqns_ice3d, zqsr_ice3d)846 CALL wrk_dealloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 1022 847 CALL wrk_dealloc( jpi,jpj, zrhoa ) 1023 CALL wrk_dealloc( jpi,jpj, Cd , zevap_ice2d, zqns_ice2d, zqsr_ice2d)848 CALL wrk_dealloc( jpi,jpj, Cd ) 1024 849 ! 1025 850 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_flx') … … 1083 908 !!---------------------------------------------------------------------------------- 1084 909 ! 1085 !$OMP PARALLEL DO schedule(static) private(jj,ji,ztmp,ze_sat)1086 910 DO jj = 1, jpj 1087 911 DO ji = 1, jpi … … 1120 944 !!---------------------------------------------------------------------------------- 1121 945 ! 1122 !$OMP PARALLEL DO schedule(static) private(jj,ji,zrv,ziRT)1123 946 DO jj = 1, jpj 1124 947 DO ji = 1, jpi
Note: See TracChangeset
for help on using the changeset viewer.
