ERF_TwoStreamRadiation.cpp File Reference

TwoStreamRadiation: the two-stream radiation model and its surface state. More...

#include <ERF_TwoStreamRadiation.H>
#include <ERF_RadStruct.H>
#include <AMReX_VisMF.H>
#include <AMReX_PlotFileUtil.H>
#include <AMReX_Utility.H>
#include <AMReX_Math.H>
#include <ERF_RadiationDiagnostics.H>
#include <ERF_TwoStreamColumn.H>
#include <ERF_TwoStreamCanopyForcing.H>
#include <ERF_Constants.H>
#include <ERF_PrognosticCloudFraction.H>
#include <ERF_AerosolOpticalDepth.H>
#include <ERF_SimplifiedSEB.H>
#include <ERF_SurfaceMoisture.H>
#include <ERF_SEBTurbulentFlux.H>
#include <ERF_OrbCosZenith.H>
#include <AMReX_Print.H>
#include <AMReX_ParallelDescriptor.H>
#include <AMReX_Gpu.H>
#include <ERF_IndexDefines.H>
#include <ERF_EOS.H>
#include <algorithm>
#include <cmath>
#include <ctime>
#include <limits>

Include dependency graph for ERF_TwoStreamRadiation.cpp:

Functions
TwoStreamSunDate	two_stream_sun_date (const RadChoice &rc, TwoStreamRadiation::OrbitalCache &orbit, double epoch_time)
	The calendar sun of the two-stream model at a UTC epoch time (start_time + t). More...
void	copy_surface_plane (const MultiFab &src2d, MultiFab &dst)

Detailed Description

TwoStreamRadiation: the two-stream radiation model and its surface state.

Computes SW/LW fluxes and per-level heating rates using real per-column vertical sweeps over the atmospheric grid. Reads temperature and density from the state and properly accumulates optical depth through all vertical levels.

Capabilities include:

• Height-varying optical depth: optional "cloud_layer" tau_profile_type adds cloud_tau_per_layer on top of the clear-sky background within [cloud_base_height_m, cloud_top_height_m].
• Cloud fraction masking: blends clear-sky and cloudy-column fluxes via F = (1 - cloud_fraction) * F_clear + cloud_fraction) * F_cloudy.
• Diffuse (scattered) SW flux: two-stream reflectance/transmittance per layer combined with the surface albedo by the adding method, giving upward and downward diffuse streams (see ERF_TwoStreamSW.H).
• Per-level heating rate output: writes the SW/LW radiative tendencies of potential temperature, (dT/dt)/pi, to a 2-component MultiFab (component 0 = SW, component 1 = LW), mirroring the RRTMGP convention expected by the RhoTheta source term.

Includes support for:

• Height-varying surface properties (albedo, emissivity, temperature)
• Dynamic cloud fraction from relative humidity and cloud water
• Prescribed bulk aerosol optical depth
• Dynamic solar geometry (time-varying solar position)
• Surface energy balance diagnostics and prognostic updates

Vertical orientation follows ERF: k = kmin is the surface layer and k = kmax the top layer. SW sweeps downward from kmax to kmin; LW sweeps downward (TOA -> surface) and then upward (surface -> TOA). Layer temperature is obtained from rho*theta through the Exner function.

CSV diagnostics (domain means unless noted): SW_surface is the SW absorbed by the surface, SW_TOA the incident SW at the top of the atmosphere, SW_up_TOA the reflected SW leaving the top, LW_net_surface the net (up - down) LW at the surface, LW_up_TOA the outgoing LW at the top, and heating_rate_max the max(|Q_sw|+|Q_lw|) over the column evaluations.

Function Documentation

copy_surface_plane()

void copy_surface_plane	(	const MultiFab &	src2d,
		MultiFab &	dst
	)

{
    BoxList flat = dst.boxArray().boxList();
    for (Box& b : flat) { b.setRange(2, 0); }
    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(
        BoxArray(std::move(flat)) == src2d.boxArray() &&
        dst.DistributionMap() == src2d.DistributionMap(),
        "copy_surface_plane: the destination's horizontal grids differ from the source's");
    for (MFIter mfi(dst); mfi.isValid(); ++mfi) {
        const Box& fab_box = dst[mfi].box();
        AMREX_ALWAYS_ASSERT_WITH_MESSAGE(fab_box.smallEnd(2) <= 0 && fab_box.bigEnd(2) >= 0,
            "copy_surface_plane: the destination must hold the k = 0 plane in its valid or ghost region");
        const Box plane = makeSlab(mfi.validbox(), 2, 0);
        dst[mfi].template copy<RunOn::Device>(src2d[mfi], plane, 0, plane, 0, 1);
    }
}

Referenced by TwoStreamRadiation::write_land_forcing().

Here is the call graph for this function:

Here is the caller graph for this function:

two_stream_sun_date()

TwoStreamSunDate two_stream_sun_date	(	const RadChoice &	rc,
		TwoStreamRadiation::OrbitalCache &	orbit,
		double	epoch_time
	)

The calendar sun of the two-stream model at a UTC epoch time (start_time + t).

The calendar sun of the two-stream model at a UTC epoch time: the calendar day with the time of day as its fraction, the declination and the Earth-Sun distance factor, from the date as the RRTMGP interface forms it and the orbital code of ERF_OrbCosZenith.H (orbital_params once per year into orbit, with the erf.rad_orbital_* overrides passed to it, then orbital_decl per call). set_solar_state() takes the sun of every sweep from here, and the immersed-boundary balance its erf.ibseb.sun_mode = two_stream sun, so the two cannot drift apart.

{
    // Calendar date of this call (UTC), as the RRTMGP interface forms it.
    time_t timestamp = time_t(epoch_time);
    struct tm timeinfo{};
#if defined(_WIN32)
    gmtime_s(&timeinfo, &timestamp);
#else
    gmtime_r(&timestamp, &timeinfo);
#endif
    const int year = (rc.rad_orbital_year >= 0) ? rc.rad_orbital_year : timeinfo.tm_year + 1900;
    const int mon = timeinfo.tm_mon + 1;
    const int day = timeinfo.tm_mday;
    const int sec = timeinfo.tm_hour*3600 + timeinfo.tm_min*60 + timeinfo.tm_sec;
 
    // Orbital parameters of the year (Berger 1978) unless overridden: a few
    // hundred series terms and six vectors, so once per year, not per step.
    if (orbit.year != year) {
        int  iyear = year;
        double eccen = rc.rad_orbital_eccentricity;
        double obliq = rc.rad_orbital_obliquity;
        double mvelp = rc.rad_orbital_mvelp;
        double obliqr = 0.0, lambm0 = 0.0, mvelpp = 0.0;
        orbital_params(iyear, eccen, obliq, mvelp, obliqr, lambm0, mvelpp);
        orbit.year = year;
        orbit.eccen = eccen; orbit.obliqr = obliqr; orbit.lambm0 = lambm0; orbit.mvelpp = mvelpp;
    }
 
    // Declination and Earth-Sun distance factor of the day.
    double calday = orbital_calday(year, mon, day, sec);
    double eccen = orbit.eccen, mvelpp = orbit.mvelpp, lambm0 = orbit.lambm0, obliqr = orbit.obliqr;
    double delta = 0.0, eccf = 1.0;
    orbital_decl(calday, eccen, mvelpp, lambm0, obliqr, delta, eccf);
 
    TwoStreamSunDate out;
    out.calday = calday;
    out.declin = delta;
    out.eccf   = eccf;
    return out;
}

Referenced by ERF::ibseb_set_two_stream_sun().

Here is the call graph for this function:

Here is the caller graph for this function: