Changeset 14280
- Timestamp:
- 2021-01-08T13:21:47+01:00 (3 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r14122_ENHANCE-14_smueller_OSMOSIS_streamlining/src/OCE/ZDF/zdfosm.F90
r14264 r14280 39 39 !! 'ln_zdfosm' OSMOSIS scheme 40 40 !!---------------------------------------------------------------------- 41 !! zdf_osm : update momentum and tracer Kz from osm scheme 42 !! zdf_osm_init : initialization, namelist read, and parameters control 43 !! osm_rst : read (or initialize) and write osmosis restart fields 44 !! tra_osm : compute and add to the T & S trend the non-local flux 45 !! trc_osm : compute and add to the passive tracer trend the non-local flux (TBD) 46 !! dyn_osm : compute and add to u & v trensd the non-local flux 41 !! zdf_osm : update momentum and tracer Kz from osm scheme 42 !! zdf_osm_vertical_average : compute vertical averages over boundary layers 43 !! zdf_osm_init : initialization, namelist read, and parameters control 44 !! osm_rst : read (or initialize) and write osmosis restart fields 45 !! tra_osm : compute and add to the T & S trend the non-local flux 46 !! trc_osm : compute and add to the passive tracer trend the non-local flux (TBD) 47 !! dyn_osm : compute and add to u & v trensd the non-local flux 47 48 !! 48 49 !! Subroutines in revised code. … … 271 272 INTEGER, DIMENSION(jpi,jpj) :: ibld ! level of boundary layer base 272 273 INTEGER, DIMENSION(jpi,jpj) :: imld ! level of mixed-layer depth (pycnocline top) 273 INTEGER, DIMENSION(jpi,jpj) :: jp_ext , jp_ext_mle! offset for external level274 INTEGER, DIMENSION(jpi,jpj) :: jp_ext ! offset for external level 274 275 INTEGER, DIMENSION(jpi, jpj) :: j_ddh ! Type of shear layer 275 276 … … 295 296 REAL(wp), DIMENSION(jpi,jpj) :: zdt_bl,zds_bl,zdu_bl,zdv_bl,zdb_bl ! difference between blayer average and parameter at base of blayer 296 297 REAL(wp), DIMENSION(jpi,jpj) :: zdt_ml,zds_ml,zdu_ml,zdv_ml,zdb_ml ! difference between mixed layer average and parameter at base of blayer 297 REAL(wp), DIMENSION(jpi,jpj) :: zdt_mle,zds_mle,zdu_mle,zdv_mle,zdb_mle ! difference between MLE layer average and parameter at base of blayer298 298 ! REAL(wp), DIMENSION(jpi,jpj) :: zwth_ent,zws_ent ! heat and salinity fluxes at the top of the pycnocline 299 299 REAL(wp) :: zwth_ent,zws_ent ! heat and salinity fluxes at the top of the pycnocline … … 368 368 zdt_ml(:,:) = 0._wp ; zds_ml(:,:) = 0._wp ; zdu_ml(:,:) = 0._wp ; zdv_ml(:,:) = 0._wp 369 369 zdb_ml(:,:) = 0._wp 370 zdt_mle(:,:) = 0._wp ; zds_mle(:,:) = 0._wp ; zdu_mle(:,:) = 0._wp371 zdv_mle(:,:) = 0._wp ; zdb_mle(:,:) = 0._wp372 370 zwth_ent = 0._wp ; zws_ent = 0._wp 373 371 ! … … 580 578 ! Averages over well-mixed and boundary layer 581 579 jp_ext(:,:) = 2 582 CALL zdf_osm_vertical_average(ibld, jp_ext, zt_bl, zs_bl, zb_bl, zu_bl, zv_bl, zdt_bl, zds_bl, zdb_bl, zdu_bl, zdv_bl) 580 CALL zdf_osm_vertical_average( Kbb, Kmm, & 581 & ibld, zt_bl, zs_bl, zb_bl, zu_bl, zv_bl, & 582 & jp_ext, zdt_bl, zds_bl, zdb_bl, zdu_bl, zdv_bl ) 583 583 ! jp_ext(:,:) = ibld(:,:) - imld(:,:) + 1 584 CALL zdf_osm_vertical_average(ibld, ibld-imld+1, zt_ml, zs_ml, zb_ml, zu_ml, zv_ml, zdt_ml, zds_ml, zdb_ml, zdu_ml, zdv_ml) 584 CALL zdf_osm_vertical_average( Kbb, Kmm, & 585 & ibld, zt_ml, zs_ml, zb_ml, zu_ml, zv_ml, ibld-imld+1, & 586 & zdt_ml, zds_ml, zdb_ml, zdu_ml, zdv_ml ) 585 587 ! Velocity components in frame aligned with surface stress. 586 588 CALL zdf_osm_velocity_rotation( zcos_wind, zsin_wind, zu_ml, zv_ml, zdu_ml, zdv_ml ) … … 595 597 IF ( hmle(ji,jj) >= gdepw(ji,jj,jk,Kmm) ) mld_prof(ji,jj) = MIN(mbkt(ji,jj), jk) 596 598 END_3D 597 jp_ext_mle(:,:) = 2598 CALL zdf_osm_vertical_average(mld_prof, jp_ext_mle, zt_mle, zs_mle, zb_mle, zu_mle, zv_mle, zdt_mle, zds_mle, zdb_mle, zdu_mle, zdv_mle)599 CALL zdf_osm_vertical_average( Kbb, Kmm, & 600 & mld_prof, zt_mle, zs_mle, zb_mle, zu_mle, zv_mle ) 599 601 600 602 DO_2D( 0, 0, 0, 0 ) … … 635 637 END_2D 636 638 637 CALL zdf_osm_vertical_average(ibld, jp_ext, zt_bl, zs_bl, zb_bl, zu_bl, zv_bl, zdt_bl, zds_bl, zdb_bl, zdu_bl, zdv_bl ) 639 CALL zdf_osm_vertical_average( Kbb, Kmm, & 640 & ibld, zt_bl, zs_bl, zb_bl, zu_bl, zv_bl, & 641 & jp_ext, zdt_bl, zds_bl, zdb_bl, zdu_bl, zdv_bl ) 638 642 ! jp_ext = ibld-imld+1 639 CALL zdf_osm_vertical_average(imld-1, ibld-imld+1, zt_ml, zs_ml, zb_ml, zu_ml, zv_ml, zdt_ml, zds_ml, zdb_ml, zdu_ml, zdv_ml) 643 CALL zdf_osm_vertical_average( Kbb, Kmm, & 644 & imld-1, zt_ml, zs_ml, zb_ml, zu_ml, zv_ml, & 645 & ibld-imld+1, zdt_ml, zds_ml, zdb_ml, zdu_ml, zdv_ml ) 640 646 ! Rate of change of hbl 641 647 CALL zdf_osm_calculate_dhdt( zdhdt, zddhdt ) … … 664 670 ! is external level in bounds? 665 671 666 CALL zdf_osm_vertical_average( ibld, jp_ext, zt_bl, zs_bl, zb_bl, zu_bl, zv_bl, zdt_bl, zds_bl, zdb_bl, zdu_bl, zdv_bl ) 672 CALL zdf_osm_vertical_average( Kbb, Kmm, & 673 & ibld, zt_bl, zs_bl, zb_bl, zu_bl, zv_bl, & 674 & jp_ext, zdt_bl, zds_bl, zdb_bl, zdu_bl, zdv_bl ) 667 675 ! 668 676 ! … … 679 687 ! Average over the depth of the mixed layer in the convective boundary layer 680 688 ! jp_ext = ibld - imld +1 681 CALL zdf_osm_vertical_average( imld-1, ibld-imld+1, zt_ml, zs_ml, zb_ml, zu_ml, zv_ml, zdt_ml, zds_ml, zdb_ml, zdu_ml, zdv_ml ) 689 CALL zdf_osm_vertical_average( Kbb, Kmm, & 690 & imld-1, zt_ml, zs_ml, zb_ml, zu_ml, zv_ml, & 691 & ibld-imld+1, zdt_ml, zds_ml, zdb_ml, zdu_ml, zdv_ml ) 682 692 ! rotate mean currents and changes onto wind align co-ordinates 683 693 ! … … 1540 1550 END SUBROUTINE zdf_osm_osbl_state 1541 1551 1542 1543 SUBROUTINE zdf_osm_vertical_average( jnlev_av, jp_ext, zt, zs, zb, zu, zv, zdt, zds, zdb, zdu, zdv )1544 !!---------------------------------------------------------------------1545 !! *** ROUTINE zdf_vertical_average ***1546 !!1547 !! ** Purpose : Determines vertical averages from surface to jnlev.1548 !!1549 !! ** Method : Averages are calculated from the surface to jnlev.1550 !! The external level used to calculate differences is ibld+ibld_ext1551 !!1552 !!----------------------------------------------------------------------1553 1554 INTEGER, DIMENSION(jpi,jpj) :: jnlev_av ! Number of levels to average over.1555 INTEGER, DIMENSION(jpi,jpj) :: jp_ext1556 1557 ! Alan: do we need zb?1558 REAL(wp), DIMENSION(jpi,jpj) :: zt, zs, zb ! Average temperature and salinity1559 REAL(wp), DIMENSION(jpi,jpj) :: zu,zv ! Average current components1560 REAL(wp), DIMENSION(jpi,jpj) :: zdt, zds, zdb ! Difference between average and value at base of OSBL1561 REAL(wp), DIMENSION(jpi,jpj) :: zdu, zdv ! Difference for velocity components.1562 1563 INTEGER :: jk, ji, jj, ibld_ext1564 REAL(wp) :: zthick, zthermal, zbeta1565 1566 1567 IF( ln_timing ) CALL timing_start('zdf_osm_va')1568 zt = 0._wp1569 zs = 0._wp1570 zu = 0._wp1571 zv = 0._wp1572 DO_2D( 0, 0, 0, 0 )1573 zthermal = rab_n(ji,jj,1,jp_tem) !ideally use ibld not 1??1574 zbeta = rab_n(ji,jj,1,jp_sal)1575 ! average over depth of boundary layer1576 zthick = epsln1577 DO jk = 2, jnlev_av(ji,jj)1578 zthick = zthick + e3t(ji,jj,jk,Kmm)1579 zt(ji,jj) = zt(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm)1580 zs(ji,jj) = zs(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm)1581 zu(ji,jj) = zu(ji,jj) + e3t(ji,jj,jk,Kmm) &1582 & * ( uu(ji,jj,jk,Kbb) + uu(ji - 1,jj,jk,Kbb) ) &1583 & / MAX( 1. , umask(ji,jj,jk) + umask(ji - 1,jj,jk) )1584 zv(ji,jj) = zv(ji,jj) + e3t(ji,jj,jk,Kmm) &1585 & * ( vv(ji,jj,jk,Kbb) + vv(ji,jj - 1,jk,Kbb) ) &1586 & / MAX( 1. , vmask(ji,jj,jk) + vmask(ji,jj - 1,jk) )1587 END DO1588 zt(ji,jj) = zt(ji,jj) / zthick1589 zs(ji,jj) = zs(ji,jj) / zthick1590 zu(ji,jj) = zu(ji,jj) / zthick1591 zv(ji,jj) = zv(ji,jj) / zthick1592 zb(ji,jj) = grav * zthermal * zt(ji,jj) - grav * zbeta * zs(ji,jj)1593 ibld_ext = jnlev_av(ji,jj) + jp_ext(ji,jj)1594 IF ( ibld_ext < mbkt(ji,jj) ) THEN1595 zdt(ji,jj) = zt(ji,jj) - ts(ji,jj,ibld_ext,jp_tem,Kmm)1596 zds(ji,jj) = zs(ji,jj) - ts(ji,jj,ibld_ext,jp_sal,Kmm)1597 zdu(ji,jj) = zu(ji,jj) - ( uu(ji,jj,ibld_ext,Kbb) + uu(ji-1,jj,ibld_ext,Kbb ) ) &1598 & / MAX(1. , umask(ji,jj,ibld_ext ) + umask(ji-1,jj,ibld_ext ) )1599 zdv(ji,jj) = zv(ji,jj) - ( vv(ji,jj,ibld_ext,Kbb) + vv(ji,jj-1,ibld_ext,Kbb ) ) &1600 & / MAX(1. , vmask(ji,jj,ibld_ext ) + vmask(ji,jj-1,ibld_ext ) )1601 zdb(ji,jj) = grav * zthermal * zdt(ji,jj) - grav * zbeta * zds(ji,jj)1602 ELSE1603 zdt(ji,jj) = 0._wp1604 zds(ji,jj) = 0._wp1605 zdu(ji,jj) = 0._wp1606 zdv(ji,jj) = 0._wp1607 zdb(ji,jj) = 0._wp1608 ENDIF1609 END_2D1610 IF( ln_timing ) CALL timing_stop('zdf_osm_va')1611 END SUBROUTINE zdf_osm_vertical_average1612 1552 1613 1553 SUBROUTINE zdf_osm_velocity_rotation( zcos_w, zsin_w, zu, zv, zdu, zdv ) … … 2478 2418 END SUBROUTINE zdf_osm 2479 2419 2420 SUBROUTINE zdf_osm_vertical_average( Kbb, Kmm, & 2421 & knlev, pt, ps, pb, pu, pv, & 2422 & kp_ext, pdt, pds, pdb, pdu, pdv ) 2423 !!--------------------------------------------------------------------- 2424 !! *** ROUTINE zdf_vertical_average *** 2425 !! 2426 !! ** Purpose : Determines vertical averages from surface to knlev, 2427 !! and optionally the differences between these vertical 2428 !! averages and values at an external level 2429 !! 2430 !! ** Method : Averages are calculated from the surface to knlev. 2431 !! The external level used to calculate differences is 2432 !! knlev+kp_ext 2433 !!---------------------------------------------------------------------- 2434 INTEGER, INTENT(in ) :: Kbb, Kmm ! Ocean time-level indices 2435 INTEGER, DIMENSION(jpi,jpj), INTENT(in ) :: knlev ! Number of levels to average over. 2436 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pt, ps ! Average temperature and salinity 2437 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pb ! Average buoyancy 2438 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pu, pv ! Average current components 2439 INTEGER, DIMENSION(jpi,jpj), INTENT(in ), OPTIONAL :: kp_ext ! External-level offsets 2440 REAL(wp), DIMENSION(jpi,jpj), INTENT( out), OPTIONAL :: pdt, pds ! Difference between average temperature, salinity, 2441 REAL(wp), DIMENSION(jpi,jpj), INTENT( out), OPTIONAL :: pdb ! buoyancy, 2442 REAL(wp), DIMENSION(jpi,jpj), INTENT( out), OPTIONAL :: pdu, pdv ! velocity components and the OSBL 2443 ! 2444 INTEGER :: jk, jkflt, jkmax, ji, jj ! Loop indices 2445 INTEGER :: ibld_ext ! External-layer index 2446 REAL(wp), DIMENSION(jpi,jpj) :: zthick ! Layer thickness 2447 REAL(wp) :: zthermal, zbeta ! Thermal/haline expansion/contraction coefficients 2448 !!---------------------------------------------------------------------- 2449 ! 2450 IF( ln_timing ) CALL timing_start('zdf_osm_va') 2451 ! 2452 ! Averages over depth of boundary layer 2453 pt(:,:) = 0.0_wp 2454 ps(:,:) = 0.0_wp 2455 pu(:,:) = 0.0_wp 2456 pv(:,:) = 0.0_wp 2457 zthick(:,:) = epsln 2458 jkflt = jpk 2459 jkmax = 0 2460 DO_2D( 0, 0, 0, 0 ) 2461 IF ( knlev(ji,jj) < jkflt ) jkflt = knlev(ji,jj) 2462 IF ( knlev(ji,jj) > jkmax ) jkmax = knlev(ji,jj) 2463 END_2D 2464 DO_3D( 0, 0, 0, 0, 2, jkflt ) ! Upper, flat part of layer 2465 zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm) 2466 pt(ji,jj) = pt(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm) 2467 ps(ji,jj) = ps(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) 2468 pu(ji,jj) = pu(ji,jj) + e3t(ji,jj,jk,Kmm) * & 2469 & ( uu(ji,jj,jk,Kbb) + uu(ji - 1,jj,jk,Kbb) ) / & 2470 & MAX( 1.0_wp , umask(ji,jj,jk) + umask(ji - 1,jj,jk) ) 2471 pv(ji,jj) = pv(ji,jj) + e3t(ji,jj,jk,Kmm) * & 2472 & ( vv(ji,jj,jk,Kbb) + vv(ji,jj - 1,jk,Kbb) ) / & 2473 & MAX( 1.0_wp , vmask(ji,jj,jk) + vmask(ji,jj - 1,jk) ) 2474 END_3D 2475 DO_3D( 0, 0, 0, 0, jkflt+1, jkmax ) ! Lower, non-flat part of layer 2476 IF ( knlev(ji,jj) >= jk ) THEN 2477 zthick(ji,jj) = zthick(ji,jj) + e3t(ji,jj,jk,Kmm) 2478 pt(ji,jj) = pt(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_tem,Kmm) 2479 ps(ji,jj) = ps(ji,jj) + e3t(ji,jj,jk,Kmm) * ts(ji,jj,jk,jp_sal,Kmm) 2480 pu(ji,jj) = pu(ji,jj) + e3t(ji,jj,jk,Kmm) * & 2481 & ( uu(ji,jj,jk,Kbb) + uu(ji - 1,jj,jk,Kbb) ) / & 2482 & MAX( 1.0_wp , umask(ji,jj,jk) + umask(ji - 1,jj,jk) ) 2483 pv(ji,jj) = pv(ji,jj) + e3t(ji,jj,jk,Kmm) * & 2484 & ( vv(ji,jj,jk,Kbb) + vv(ji,jj - 1,jk,Kbb) ) / & 2485 & MAX( 1.0_wp , vmask(ji,jj,jk) + vmask(ji,jj - 1,jk) ) 2486 END IF 2487 END_3D 2488 DO_2D( 0, 0, 0, 0 ) 2489 pt(ji,jj) = pt(ji,jj) / zthick(ji,jj) 2490 ps(ji,jj) = ps(ji,jj) / zthick(ji,jj) 2491 pu(ji,jj) = pu(ji,jj) / zthick(ji,jj) 2492 pv(ji,jj) = pv(ji,jj) / zthick(ji,jj) 2493 zthermal = rab_n(ji,jj,1,jp_tem) ! ideally use ibld not 1?? 2494 zbeta = rab_n(ji,jj,1,jp_sal) 2495 pb(ji,jj) = grav * zthermal * pt(ji,jj) - grav * zbeta * ps(ji,jj) 2496 END_2D 2497 ! 2498 ! Differences between vertical averages and values at an external layer 2499 IF ( PRESENT( kp_ext ) ) THEN 2500 DO_2D( 0, 0, 0, 0 ) 2501 ibld_ext = knlev(ji,jj) + kp_ext(ji,jj) 2502 IF ( ibld_ext < mbkt(ji,jj) ) THEN 2503 pdt(ji,jj) = pt(ji,jj) - ts(ji,jj,ibld_ext,jp_tem,Kmm) 2504 pds(ji,jj) = ps(ji,jj) - ts(ji,jj,ibld_ext,jp_sal,Kmm) 2505 pdu(ji,jj) = pu(ji,jj) - ( uu(ji,jj,ibld_ext,Kbb) + uu(ji-1,jj,ibld_ext,Kbb ) ) / & 2506 & MAX(1.0_wp , umask(ji,jj,ibld_ext ) + umask(ji-1,jj,ibld_ext ) ) 2507 pdv(ji,jj) = pv(ji,jj) - ( vv(ji,jj,ibld_ext,Kbb) + vv(ji,jj-1,ibld_ext,Kbb ) ) / & 2508 & MAX(1.0_wp , vmask(ji,jj,ibld_ext ) + vmask(ji,jj-1,ibld_ext ) ) 2509 zthermal = rab_n(ji,jj,1,jp_tem) ! ideally use ibld not 1?? 2510 zbeta = rab_n(ji,jj,1,jp_sal) 2511 pdb(ji,jj) = grav * zthermal * pdt(ji,jj) - grav * zbeta * pds(ji,jj) 2512 ELSE 2513 pdt(ji,jj) = 0.0_wp 2514 pds(ji,jj) = 0.0_wp 2515 pdu(ji,jj) = 0.0_wp 2516 pdv(ji,jj) = 0.0_wp 2517 pdb(ji,jj) = 0.0_wp 2518 ENDIF 2519 END_2D 2520 END IF 2521 ! 2522 IF( ln_timing ) CALL timing_stop('zdf_osm_va') 2523 ! 2524 END SUBROUTINE zdf_osm_vertical_average 2480 2525 2481 2526 SUBROUTINE zdf_osm_init( Kmm )
Note: See TracChangeset
for help on using the changeset viewer.
