- Timestamp:
- 2016-12-19T13:15:59+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_r6519_HPC_4/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r7037 r7508 141 141 INTEGER :: ifpr ! dummy loop indice 142 142 INTEGER :: jfld ! dummy loop arguments 143 INTEGER :: ji, jj ! dummy loop indices 143 144 INTEGER :: ios ! Local integer output status for namelist read 144 145 ! … … 194 195 CALL fld_fill( sf, slf_i, cn_dir, 'sbc_blk_core', 'flux formulation for ocean surface boundary condition', 'namsbc_core' ) 195 196 ! 196 sfx(:,:) = 0._wp ! salt flux; zero unless ice is present (computed in limsbc(_2).F90) 197 ! 197 !$OMP PARALLEL DO schedule(static) private(jj, ji) 198 DO jj = 1, jpj 199 DO ji = 1, jpi 200 sfx(ji,jj) = 0._wp ! salt flux; zero unless ice is present (computed in limsbc(_2).F90) 201 END DO 202 END DO 198 203 ENDIF 199 204 … … 205 210 #if defined key_cice 206 211 IF( MOD( kt - 1, nn_fsbc ) == 0 ) THEN 207 qlw_ice(:,:,1) = sf(jp_qlw)%fnow(:,:,1) 208 qsr_ice(:,:,1) = sf(jp_qsr)%fnow(:,:,1) 209 tatm_ice(:,:) = sf(jp_tair)%fnow(:,:,1) 210 qatm_ice(:,:) = sf(jp_humi)%fnow(:,:,1) 211 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac 212 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac 213 wndi_ice(:,:) = sf(jp_wndi)%fnow(:,:,1) 214 wndj_ice(:,:) = sf(jp_wndj)%fnow(:,:,1) 212 !$OMP PARALLEL DO schedule(static) private(jj, ji) 213 DO jj = 1, jpj 214 DO ji = 1, jpi 215 qlw_ice(ji,jj,1) = sf(jp_qlw)%fnow(ji,jj,1) 216 qsr_ice(ji,jj,1) = sf(jp_qsr)%fnow(ji,jj,1) 217 tatm_ice(ji,jj) = sf(jp_tair)%fnow(ji,jj,1) 218 qatm_ice(ji,jj) = sf(jp_humi)%fnow(ji,jj,1) 219 tprecip(ji,jj) = sf(jp_prec)%fnow(ji,jj,1) * rn_pfac 220 sprecip(ji,jj) = sf(jp_snow)%fnow(ji,jj,1) * rn_pfac 221 wndi_ice(ji,jj) = sf(jp_wndi)%fnow(ji,jj,1) 222 wndj_ice(ji,jj) = sf(jp_wndj)%fnow(ji,jj,1) 223 END DO 224 END DO 215 225 ENDIF 216 226 #endif … … 267 277 zcoef_qsatw = 0.98 * 640380. / rhoa 268 278 269 !$OMP PARALLEL WORKSHARE 270 zst(:,:) = pst(:,:) + rt0 ! convert SST from Celcius to Kelvin (and set minimum value far above 0 K) 279 !$OMP PARALLEL DO schedule(static) private(jj, ji) 280 DO jj = 1, jpj 281 DO ji = 1, jpi 282 zst(ji,jj) = pst(ji,jj) + rt0 ! convert SST from Celcius to Kelvin (and set minimum value far above 0 K) 271 283 272 284 ! ----------------------------------------------------------------------------- ! … … 275 287 276 288 ! ... components ( U10m - U_oce ) at T-point (unmasked) 277 zwnd_i(:,:) = 0.e0 278 zwnd_j(:,:) = 0.e0 279 !$OMP END PARALLEL WORKSHARE 289 zwnd_i(ji,jj) = 0.e0 290 zwnd_j(ji,jj) = 0.e0 291 END DO 292 END DO 280 293 #if defined key_cyclone 281 294 CALL wnd_cyc( kt, zwnd_i, zwnd_j ) ! add analytical tropical cyclone (Vincent et al. JGR 2012) … … 312 325 ! ocean albedo assumed to be constant + modify now Qsr to include the diurnal cycle ! Short Wave 313 326 zztmp = 1. - albo 314 IF( ln_dm2dc ) THEN ; qsr(:,:) = zztmp * sbc_dcy( sf(jp_qsr)%fnow(:,:,1) ) * tmask(:,:,1) 315 ELSE ; qsr(:,:) = zztmp * sf(jp_qsr)%fnow(:,:,1) * tmask(:,:,1) 327 IF( ln_dm2dc ) THEN 328 qsr(:,:) = zztmp * sbc_dcy( sf(jp_qsr)%fnow(:,:,1) ) * tmask(:,:,1) 329 ELSE 330 qsr(:,:) = zztmp * sf(jp_qsr)%fnow(:,:,1) * tmask(:,:,1) 316 331 ENDIF 317 332 … … 319 334 DO jj = 1, jpj 320 335 DO ji = 1, jpi 321 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 Wave322 ! ----------------------------------------------------------------------------- !323 ! II Turbulent FLUXES !324 ! ----------------------------------------------------------------------------- !325 326 ! ... specific humidity at SST and IST327 zqsatw(ji,jj) = zcoef_qsatw * EXP( -5107.4 / zst(ji,jj) )336 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 337 ! ----------------------------------------------------------------------------- ! 338 ! II Turbulent FLUXES ! 339 ! ----------------------------------------------------------------------------- ! 340 341 ! ... specific humidity at SST and IST 342 zqsatw(ji,jj) = zcoef_qsatw * EXP( -5107.4 / zst(ji,jj) ) 328 343 329 344 END DO … … 346 361 ! ... add the HF tau contribution to the wind stress module? 347 362 IF( lhftau ) THEN 348 taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:,1) 363 !$OMP PARALLEL DO schedule(static) private(jj, ji) 364 DO jj = 1, jpj 365 DO ji = 1, jpi 366 taum(ji,jj) = taum(ji,jj) + sf(jp_tdif)%fnow(ji,jj,1) 367 END DO 368 END DO 349 369 ENDIF 350 370 CALL iom_put( "taum_oce", taum ) ! output wind stress module … … 371 391 DO jj = 1, jpj 372 392 DO ji = 1, jpi 373 IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN 374 !! q_air and t_air are (or "are almost") given at 10m (wind reference height) 375 zevap(ji,jj) = rn_efac*MAX( 0._wp, rhoa*Ce(ji,jj)*( zqsatw(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) )*wndm(ji,jj) ) ! Evaporation 376 zqsb (ji,jj) = cpa*rhoa*Ch(ji,jj)*( zst (ji,jj) - sf(jp_tair)%fnow(ji,jj,1) )*wndm(ji,jj) ! Sensible Heat 377 ELSE 378 !! q_air and t_air are not given at 10m (wind reference height) 379 ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!! 380 zevap(ji,jj) = rn_efac*MAX( 0._wp, rhoa*Ce(ji,jj)*( zqsatw(ji,jj) - zq_zu(ji,jj) )*wndm(ji,jj) ) ! Evaporation 381 zqsb (ji,jj) = cpa*rhoa*Ch(ji,jj)*( zst (ji,jj) - zt_zu(ji,jj) )*wndm(ji,jj) ! Sensible Heat 382 ENDIF 383 zqla (ji,jj) = Lv * zevap(ji,jj) ! Latent Heat 384 393 IF( ABS( rn_zu - rn_zqt) < 0.01_wp ) THEN 394 !! q_air and t_air are (or "are almost") given at 10m (wind reference height) 395 zevap(ji,jj) = rn_efac*MAX( 0._wp, rhoa*Ce(ji,jj)*( zqsatw(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) )*wndm(ji,jj) ) ! Evaporation 396 zqsb (ji,jj) = cpa*rhoa*Ch(ji,jj)*( zst (ji,jj) - sf(jp_tair)%fnow(ji,jj,1) )*wndm(ji,jj) ! Sensible Heat 397 ELSE 398 !! q_air and t_air are not given at 10m (wind reference height) 399 ! Values of temp. and hum. adjusted to height of wind during bulk algorithm iteration must be used!!! 400 zevap(ji,jj) = rn_efac*MAX( 0._wp, rhoa*Ce(ji,jj)*( zqsatw(ji,jj) - zq_zu(ji,jj) )*wndm(ji,jj) ) ! Evaporation 401 zqsb (ji,jj) = cpa*rhoa*Ch(ji,jj)*( zst (ji,jj) - zt_zu(ji,jj) )*wndm(ji,jj) ! Sensible Heat 402 ENDIF 403 zqla (ji,jj) = Lv * zevap(ji,jj) ! Latent Heat 385 404 END DO 386 405 END DO … … 417 436 ! 418 437 #if defined key_lim3 419 !$OMP PARALLEL WORKSHARE 420 qns_oce(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) ! non solar without emp (only needed by LIM3) 421 qsr_oce(:,:) = qsr(:,:) 422 !$OMP END PARALLEL WORKSHARE 438 !$OMP PARALLEL DO schedule(static) private(jj, ji) 439 DO jj = 1, jpj 440 DO ji = 1, jpi 441 qns_oce(ji,jj) = zqlw(ji,jj) - zqsb(ji,jj) - zqla(ji,jj) ! non solar without emp (only needed by LIM3) 442 qsr_oce(ji,jj) = qsr(ji,jj) 443 END DO 444 END DO 423 445 #endif 424 446 ! … … 431 453 CALL iom_put( "qsr_oce" , qsr ) ! output downward solar heat over the ocean 432 454 CALL iom_put( "qt_oce" , qns+qsr ) ! output total downward heat over the ocean 433 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac ! output total precipitation [kg/m2/s] 434 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! output solid precipitation [kg/m2/s] 455 !$OMP PARALLEL DO schedule(static) private(jj, ji) 456 DO jj = 1, jpj 457 DO ji = 1, jpi 458 tprecip(ji,jj) = sf(jp_prec)%fnow(ji,jj,1) * rn_pfac ! output total precipitation [kg/m2/s] 459 sprecip(ji,jj) = sf(jp_snow)%fnow(ji,jj,1) * rn_pfac ! output solid precipitation [kg/m2/s] 460 END DO 461 END DO 435 462 CALL iom_put( 'snowpre', sprecip * 86400. ) ! Snow 436 463 CALL iom_put( 'precip' , tprecip * 86400. ) ! Total precipitation … … 477 504 478 505 !!gm brutal.... 479 !$OMP PARALLEL WORKSHARE 480 utau_ice (:,:) = 0._wp 481 vtau_ice (:,:) = 0._wp 482 wndm_ice (:,:) = 0._wp 483 !$OMP END PARALLEL WORKSHARE 506 !$OMP PARALLEL DO schedule(static) private(jj, ji) 507 DO jj = 1, jpj 508 DO ji = 1, jpi 509 utau_ice (ji,jj) = 0._wp 510 vtau_ice (ji,jj) = 0._wp 511 wndm_ice (ji,jj) = 0._wp 512 END DO 513 END DO 484 514 !!gm end 485 515 … … 515 545 ! 516 546 CASE( 'C' ) ! C-grid ice dynamics : U & V-points (same as ocean) 517 !$OMP PARALLEL DO schedule(static) private(jj,ji,zwndi_t,zwndj_t) 547 !$OMP PARALLEL 548 !$OMP DO schedule(static) private(jj,ji,zwndi_t,zwndj_t) 518 549 DO jj = 2, jpj 519 550 DO ji = fs_2, jpi ! vect. opt. … … 523 554 END DO 524 555 END DO 525 !$OMP PARALLELDO schedule(static) private(jj,ji)556 !$OMP DO schedule(static) private(jj,ji) 526 557 DO jj = 2, jpjm1 527 558 DO ji = fs_2, fs_jpim1 ! vect. opt. … … 532 563 END DO 533 564 END DO 565 !$OMP END PARALLEL 534 566 CALL lbc_lnk( utau_ice, 'U', -1. ) 535 567 CALL lbc_lnk( vtau_ice, 'V', -1. ) … … 637 669 END DO 638 670 ! 639 !$OMP WORKSHARE 640 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac ! total precipitation [kg/m2/s] 641 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! solid precipitation [kg/m2/s] 642 !$OMP END WORKSHARE 671 !$OMP DO schedule(static) private(jj, ji) 672 DO jj = 1, jpj 673 DO ji = 1, jpi 674 tprecip(ji,jj) = sf(jp_prec)%fnow(ji,jj,1) * rn_pfac ! total precipitation [kg/m2/s] 675 sprecip(ji,jj) = sf(jp_snow)%fnow(ji,jj,1) * rn_pfac ! solid precipitation [kg/m2/s] 676 END DO 677 END DO 643 678 !$OMP END PARALLEL 644 679 CALL iom_put( 'snowpre', sprecip * 86400. ) ! Snow precipitation … … 650 685 ! --- evaporation --- ! 651 686 z1_lsub = 1._wp / Lsub 652 !$OMP PARALLEL WORKSHARE 653 evap_ice (:,:,:) = rn_efac * qla_ice (:,:,:) * z1_lsub ! sublimation 654 devap_ice(:,:,:) = rn_efac * dqla_ice(:,:,:) * z1_lsub ! d(sublimation)/dT 655 zevap (:,:) = rn_efac * ( emp(:,:) + tprecip(:,:) ) ! evaporation over ocean 656 657 ! --- evaporation minus precipitation --- ! 658 zsnw(:,:) = 0._wp 659 !$OMP END PARALLEL WORKSHARE 687 688 !$OMP PARALLEL 689 !$OMP DO schedule(static) private(jl, jj, ji) 690 DO jl = 1, jpl 691 DO jj = 1, jpj 692 DO ji = 1, jpi 693 evap_ice (ji,jj,jl) = rn_efac * qla_ice (ji,jj,jl) * z1_lsub ! sublimation 694 devap_ice(ji,jj,jl) = rn_efac * dqla_ice(ji,jj,jl) * z1_lsub ! d(sublimation)/dT 695 END DO 696 END DO 697 END DO 698 !$OMP DO schedule(static) private(jj, ji) 699 DO jj = 1, jpj 700 DO ji = 1, jpi 701 zevap (ji,jj) = rn_efac * ( emp(ji,jj) + tprecip(ji,jj) ) ! evaporation over ocean 702 ! --- evaporation minus precipitation --- ! 703 zsnw(ji,jj) = 0._wp 704 END DO 705 END DO 706 !$OMP END PARALLEL 660 707 CALL lim_thd_snwblow( pfrld, zsnw ) ! snow distribution over ice after wind blowing 708 661 709 emp_oce(:,:) = pfrld(:,:) * zevap(:,:) - ( tprecip(:,:) - sprecip(:,:) ) - sprecip(:,:) * (1._wp - zsnw ) 662 710 emp_ice(:,:) = SUM( a_i_b(:,:,:) * evap_ice(:,:,:), dim=3 ) - sprecip(:,:) * zsnw … … 664 712 665 713 ! --- heat flux associated with emp --- ! 666 qemp_oce(:,:) = - pfrld(:,:) * zevap(:,:) * sst_m(:,:) * rcp 667 668 & + sprecip(:,:) * ( 1._wp - zsnw ) *& ! solid precip at min(Tair,Tsnow)669 670 qemp_ice(:,:) = sprecip(:,:) * zsnw * 671 714 qemp_oce(:,:) = - pfrld(:,:) * zevap(:,:) * sst_m(:,:) * rcp & ! evap at sst 715 & + ( tprecip(:,:) - sprecip(:,:) ) * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & ! liquid precip at Tair 716 & + sprecip(:,:) * ( 1._wp - zsnw ) * & ! solid precip at min(Tair,Tsnow) 717 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 718 qemp_ice(:,:) = sprecip(:,:) * zsnw * & ! solid precip (only) 719 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 672 720 673 721 ! --- total solar and non solar fluxes --- ! … … 675 723 qsr_tot(:,:) = pfrld(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 ) 676 724 677 ! --- heat content of precip over ice in J/m3 (to be used in 1D-thermo) --- ! 725 ! --- heat content of precip over ice in J/m3 (to be used in 1D-thermo) 726 ! --- ! 678 727 qprec_ice(:,:) = rhosn * ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 679 728 … … 684 733 ! But we do not have Tice => consider it at 0°C => evap=0 685 734 END DO 686 687 735 CALL wrk_dealloc( jpi,jpj, zevap, zsnw ) 688 736 #endif … … 693 741 ! ( Maykut and Untersteiner, 1971 ; Ebert and Curry, 1993 ) 694 742 ! 695 !$OMP PARALLEL WORKSHARE 696 fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 697 fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 698 !$OMP END PARALLEL WORKSHARE 743 !$OMP PARALLEL DO schedule(static) private(jj, ji) 744 DO jj = 1, jpj 745 DO ji = 1, jpi 746 fr1_i0(ji,jj) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 747 fr2_i0(ji,jj) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 748 END DO 749 END DO 699 750 ! 700 751 ! … … 753 804 ! 754 805 INTEGER :: j_itt 806 INTEGER :: ji, jj ! dummy loop indices 755 807 INTEGER , PARAMETER :: nb_itt = 5 ! number of itterations 756 808 LOGICAL :: l_zt_equal_zu = .FALSE. ! if q and t are given at different height than U … … 787 839 !! Neutral coefficients at 10m: 788 840 IF( ln_cdgw ) THEN ! wave drag case 789 !$OMP PARALLEL WORKSHARE 790 cdn_wave(:,:) = cdn_wave(:,:) + rsmall * ( 1._wp - tmask(:,:,1) ) 791 ztmp0 (:,:) = cdn_wave(:,:) 792 !$OMP END PARALLEL WORKSHARE 841 !$OMP PARALLEL DO schedule(static) private(jj, ji) 842 DO jj = 1, jpj 843 DO ji = 1, jpi 844 cdn_wave(ji,jj) = cdn_wave(ji,jj) + rsmall * ( 1._wp - tmask(ji,jj,1) ) 845 ztmp0 (ji,jj) = cdn_wave(ji,jj) 846 END DO 847 END DO 793 848 ELSE 794 849 ztmp0 = cd_neutral_10m( U_zu )
Note: See TracChangeset
for help on using the changeset viewer.
