IBFaceSet Class Reference

#include <ERF_IBFaceSet.H>

Collaboration diagram for IBFaceSet:

Classes
struct	SunState
	Sun and irradiances of the last compute_shortwave() call, for reports. More...

Public Member Functions
	IBFaceSet (const IBSEBParams &params, int lev)
void	build (const amrex::MultiFab &blanking, const amrex::Geometry &geom)
int	n_faces () const
int	n_buildings () const
int	n_layers () const
int	state_ncomp () const
int	n_slots () const
const amrex::BoxArray &	state_boxarray () const
bool	has_state () const
amrex::MultiFab	make_state () const
void	compute_view_fractions ()
void	compute_longwave (const amrex::MultiFab &cons)
void	assign_materials ()
void	compute_ground (amrex::Real dt)
void	solve_balance (amrex::Real dt)
void	compute_sensible (const amrex::MultiFab &cons, const amrex::MultiFab &xvel, const amrex::MultiFab &yvel, const amrex::MultiFab &zvel, amrex::Real c_p, const amrex::MultiFab *olen_ground=nullptr, const amrex::MultiFab *pblh_ground=nullptr, amrex::Real z_i_bulk=0.0)
void	add_heat_flux_to_source (amrex::MultiFab &source, const amrex::MultiFab &cons, const amrex::Geometry &geom, amrex::Real c_p, amrex::Real rdOcp) const
void	compute_shortwave (amrex::Real time)
void	scatter_field (const amrex::Gpu::DeviceVector< amrex::Real > &v, amrex::MultiFab &out) const
void	dump_faces (const std::string &prefix) const
void	scatter_diagnostics (amrex::MultiFab &nfaces, amrex::MultiFab &tskin) const
void	save_state (amrex::MultiFab &state) const
void	load_state (const amrex::MultiFab &state)
void	report (amrex::Real time, int step, bool write_csv) const
const SunState &	sun () const
void	print_debug_summary () const
void	add_cost (double seconds)
void	set_init_cost (double seconds)
const std::vector< int > &	fab_start () const

Public Attributes
amrex::Gpu::DeviceVector< int >	d_i
amrex::Gpu::DeviceVector< int >	d_j
amrex::Gpu::DeviceVector< int >	d_k
	Fluid cell of the face. More...
amrex::Gpu::DeviceVector< int >	d_dir
	Face direction: 0 x, 1 y, 2 z. More...
amrex::Gpu::DeviceVector< int >	d_side
	Solid neighbour at -1 (low) or +1 (high) side. More...
amrex::Gpu::DeviceVector< int >	d_slot
	Slot of the face within its cell's checkpoint record. More...
amrex::Gpu::DeviceVector< int >	d_bid
	Building id (1-based) More...
amrex::Gpu::DeviceVector< int >	d_mat
	Material id (0 = uniform inputs) More...
amrex::Gpu::DeviceVector< amrex::Real >	d_albedo
amrex::Gpu::DeviceVector< amrex::Real >	d_emis
	Optical properties of the face. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_kth
amrex::Gpu::DeviceVector< amrex::Real >	d_rhocp
amrex::Gpu::DeviceVector< amrex::Real >	d_thick
	Slab conductivity [W/m/K], heat capacity [J/m3/K], thickness [m]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_area
	Face area [m2]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_xf
amrex::Gpu::DeviceVector< amrex::Real >	d_yf
amrex::Gpu::DeviceVector< amrex::Real >	d_zf
	Face centre [m]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_T_skin
	Skin temperature [K]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_T_slab
	Slab layers [K], face-major: [face * n_layers + layer]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_f_sky
amrex::Gpu::DeviceVector< amrex::Real >	d_f_ground
amrex::Gpu::DeviceVector< amrex::Real >	d_f_bldg
	View fractions, sum to 1. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_SW_abs
amrex::Gpu::DeviceVector< amrex::Real >	d_LW_net
amrex::Gpu::DeviceVector< amrex::Real >	d_H
amrex::Gpu::DeviceVector< amrex::Real >	d_LE
amrex::Gpu::DeviceVector< amrex::Real >	d_G
amrex::Gpu::DeviceVector< amrex::Real >	d_Q_ext
	External incident flux [W/m2], e.g. fire radiation. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_shadow
	1 when the ray to the sun hits a building, else 0 More...
amrex::Gpu::DeviceVector< amrex::Real >	d_SW_direct_in
	Direct beam incident on the face. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_SW_diffuse_in
	Sky and ground-reflected diffuse incident. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_LW_down_in
	Incoming longwave on the face (sky, ground and walls) More...
amrex::Gpu::DeviceVector< amrex::Real >	d_T_air
	Air temperature of the fluid cell. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_theta_air
	Potential temperature of the fluid cell [K]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_rho
	Density of the fluid cell [kg/m3]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_U_tan
	Tangential wind speed at the fluid cell [m/s]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_ustar
	Friction velocity of the wall function [m/s]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_H_coeff
	H = H_coeff (theta_skin - theta_air) [W/m2/K]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_hbld
	Height of the face's building column [m]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_w_star
	Convective velocity scale of the step [m/s]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_olen
	Obukhov length used [m] (1e30 when neutral) More...
amrex::Gpu::DeviceVector< amrex::Real >	d_z_i
	Depth in w* [m]: mixed layer above a roof, building height for a wall. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_LW_ext
	Incoming longwave from sky and ground only [W/m2]. More...
amrex::Gpu::DeviceVector< amrex::Real >	d_resid
	|balance| after the solve [W/m2] More...
amrex::Gpu::DeviceVector< int >	d_niter
	Newton iterations of the last solve. More...
amrex::Gpu::DeviceVector< int >	d_col_top
int	m_col_i0 = 0
int	m_col_j0 = 0
int	m_col_nx = 0
int	m_col_ny = 0
amrex::BoxArray	m_state_ba
amrex::BoxArray	m_xfer_ba
amrex::DistributionMapping	m_xfer_dm
std::vector< int >	m_xfer_src
amrex::Long	m_domain_cells = 0
	Cells of the level, for the debug summary of the field. More...

Private Attributes
IBSEBParams	m_params
	Copy of the inputs. More...
int	m_lev = 0
	AMR level, for reports. More...
int	m_nface = 0
	Faces owned by this rank. More...
int	m_nslots = 0
	Slots per cell of the checkpoint field (largest face count on a cell) More...
int	m_nbld = 0
	Buildings of the level (global) More...
std::vector< int >	m_fab_start
	See fab_start() More...
std::vector< amrex::Long >	m_bld_nface
	Faces per building. More...
std::vector< amrex::Real >	m_bld_area
	Face area per building [m2]. More...
amrex::Long	m_nface_dir [3] = {0, 0, 0}
	Faces per direction. More...
amrex::Real	m_area_total = 0.0
	Total face area [m2]. More...
std::vector< amrex::Long >	m_bld_ncol
std::vector< int >	m_bld_ilo
std::vector< int >	m_bld_ihi
std::vector< int >	m_bld_jlo
std::vector< int >	m_bld_jhi
int	m_nx = 0
int	m_ny = 0
	Columns of the domain. More...
amrex::Real	m_x_lo = 0.0
amrex::Real	m_y_lo = 0.0
	Domain origin [m]. More...
amrex::Real	m_dx [3] = {0.0, 0.0, 0.0}
	Cell sizes [m]. More...
bool	m_per_x = false
bool	m_per_y = false
	Periodicity in x and y. More...
amrex::Real	m_col_top_max = 0.0
	Tallest column [m]. More...
amrex::Real	m_z_ground = 0.0
	Ground height [m] (flat, the domain bottom) More...
amrex::Real	m_max_path = 0.0
	Ray path cap [m]. More...
SunState	m_sun
double	m_cost_s = 0.0
	Wall-clock seconds spent in the per-step update on this rank. More...
amrex::Long	m_cost_n = 0
	Number of per-step updates timed. More...
double	m_init_cost_s = 0.0
	Wall-clock seconds of build() and compute_view_fractions() More...

Constructor & Destructor Documentation

IBFaceSet()

IBFaceSet::IBFaceSet ( const IBSEBParams &  params, int  lev  )

inline

Construct an empty face set for one level; build() fills it.

Parameters

params	Inputs of the balance (copied; the set never reads ParmParse).
lev	AMR level the set belongs to, for reports only.

: m_params(params), m_lev(lev) {}

Member Function Documentation

add_cost()

void IBFaceSet::add_cost ( double  seconds )

inline

Add the wall-clock time of one per-step update of this rank, for the cost line of report().

{ m_cost_s += seconds; ++m_cost_n; }

add_heat_flux_to_source()

void IBFaceSet::add_heat_flux_to_source	(	amrex::MultiFab &	source,
		const amrex::MultiFab &	cons,
		const amrex::Geometry &	geom,
		amrex::Real	c_p,
		amrex::Real	rdOcp
	)		const

Add the face heat flux to the temperature equation: every face deposits H A / (c_p V Pi) into the rho-theta source of its fluid cell (Pi the Exner function, so the flux heats the cell by H A / (c_p V) in temperature). Called after the sources are rebuilt at every slow stage; adds, never overwrites, and uses atomic adds because a corner cell receives several faces.

Face flux into the temperature equation. The rho-theta tendency of the fluid cell is H A / (c_p V Pi); with several faces on one cell the adds are atomic.

{
    const auto dx = geom.CellSizeArray();
    const Real vol = dx[0] * dx[1] * dx[2];
    const int*  pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const Real* pH = d_H.data();  const Real* pA = d_area.data();
    for (MFIter mfi(source); mfi.isValid(); ++mfi) {
        const int f0 = m_fab_start[mfi.LocalIndex()];
        const int f1 = m_fab_start[mfi.LocalIndex() + 1];
        auto const& src = source.array(mfi);
        auto const& c   = cons.const_array(mfi);
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f = f0 + m;
            const int i = pi[f], j = pj[f], k = pk[f];
            const Real exner = getExnergivenRTh(c(i, j, k, RhoTheta_comp), rdOcp);
            const Real tend  = pH[f] * pA[f] / (c_p * vol * exner);
            Gpu::Atomic::AddNoRet(&src(i, j, k, RhoTheta_comp), tend);
        });
    }
    Gpu::streamSynchronize();
}

Here is the call graph for this function:

assign_materials()

void IBFaceSet::assign_materials

(

)

Assign every face its material properties: from the material library by building id, or the uniform inputs. Called once after build(); fills d_mat, d_albedo, d_emis, d_k, d_rhocp, d_thick and initialises the slab layers to the skin temperature.

Materials per face from the library by building id, or the uniform inputs. Host work on the static geometry, once.

{
    if (m_params.material_file.empty()) { return; }
    const Vector<IBSEBMaterial> table = read_ibseb_materials(m_params.material_file);
    auto find = [&](int id) -> const IBSEBMaterial& {
        for (const auto& m : table) { if (m.mat_id == id) return m; }
        Abort("erf.ibseb: material id " + std::to_string(id) + " is not in " + m_params.material_file);
        return table[0];
    };
    std::vector<int> h_b(m_nface);
    Gpu::copy(Gpu::deviceToHost, d_bid.begin(), d_bid.end(), h_b.begin());
    Gpu::streamSynchronize();
    std::vector<int>  h_mat(m_nface);
    std::vector<Real> h_alb(m_nface), h_emi(m_nface), h_k(m_nface), h_rc(m_nface), h_th(m_nface);
    for (int n = 0; n < m_nface; ++n) {
        const int b = h_b[n];
        int id = m_params.material_default;
        if (b >= 1 && b <= static_cast<int>(m_params.material_by_building.size())) {
            id = m_params.material_by_building[b - 1];
        }
        const IBSEBMaterial& m = find(id);
        h_mat[n] = m.mat_id; h_alb[n] = m.albedo; h_emi[n] = m.emissivity;
        h_k[n] = m.k_therm; h_rc[n] = m.rho_cp; h_th[n] = m.thickness;
    }
    upload(d_mat, h_mat); upload(d_albedo, h_alb); upload(d_emis, h_emi);
    upload(d_kth, h_k); upload(d_rhocp, h_rc); upload(d_thick, h_th);
    Gpu::streamSynchronize();
    if (m_params.debug) {
        Print() << "[IBSEB DEBUG] lev=" << m_lev << " materials: " << table.size() << " in "
                << m_params.material_file << ", default id " << m_params.material_default;
        for (size_t b = 0; b < m_params.material_by_building.size(); ++b) {
            Print() << ", building " << b + 1 << " -> " << m_params.material_by_building[b];
        }
        Print() << "\n";
        for (const auto& m : table) {
            Print() << "[IBSEB DEBUG]   material " << m.mat_id << " (" << m.name << "): albedo=" << m.albedo
                    << " emissivity=" << m.emissivity << " k=" << m.k_therm << " W/m/K rho_cp=" << m.rho_cp
                    << " J/m3/K thickness=" << m.thickness << " m\n";
        }
    }
}

Here is the call graph for this function:

build()

void IBFaceSet::build	(	const amrex::MultiFab &	blanking,
		const amrex::Geometry &	geom
	)

Detect this rank's faces from the cell-centred blanking and allocate every per-face array.

Runs on the host once at initialisation: the blanking is copied fab by fab (a pinned copy on GPU builds), the valid cells are scanned, and the arrays are uploaded to the device at the end. Must be called after the blanking has been built and its ghost cells filled, since the solid neighbour of a fluid cell on a box edge lives in a ghost cell.

Parameters

blanking	Cell-centred blanking of the level (needs at least one ghost cell, filled).
geom	Geometry of the level: cell sizes for the face areas, domain and periodicity for the neighbour tests.

Detect the faces owned by this rank and allocate the per-face arrays.

The scan visits every valid cell of every local fab in MFIter order. A fluid cell (blanking < 0.5) contributes one face for each of its six neighbours that is solid (blanking >= 0.5); the neighbour may sit in a ghost cell, which is why the blanking must have its ghost cells filled. Neighbours outside the domain in a non-periodic direction are skipped, so a building against a non-periodic boundary has no face there. Solid cells contribute nothing but the column mask used for the building ids.

Faces are appended in scan order, which makes them contiguous per fab; m_fab_start records where each fab's faces begin. The per-direction, per-building and total counts are reduced over all ranks here, once, so report() does not have to reduce static numbers every time.

{
    const Box& domain = geom.Domain();
    const int nx  = domain.length(0);
    const int ny  = domain.length(1);
    const int ilo = domain.smallEnd(0);
    const int jlo = domain.smallEnd(1);
    const auto dx = geom.CellSizeArray();
    const auto plo = geom.ProbLoArray();
    const Real face_area[3] = { dx[1] * dx[2], dx[0] * dx[2], dx[0] * dx[1] };
    m_nx = nx; m_ny = ny;
    m_x_lo = plo[0]; m_y_lo = plo[1];
    for (int d = 0; d < 3; ++d) { m_dx[d] = dx[d]; }
    m_per_x = geom.isPeriodic(0); m_per_y = geom.isPeriodic(1);
    m_max_path = 4.0 * (geom.ProbLength(0) + geom.ProbLength(1));
    m_z_ground = plo[2];
 
    // The face list of a cell, and with it the slot numbering of the
    // checkpoint record (see state_ncomp()), must not depend on the
    // decomposition: every valid cell has to see its six neighbours, which
    // needs a filled ghost layer on the blanking.
    AMREX_ALWAYS_ASSERT_WITH_MESSAGE(blanking.nGrowVect().min() >= 1,
                                     "erf.ibseb: the blanking needs at least one ghost cell");
 
    std::vector<int>  h_i, h_j, h_k, h_dir, h_side, h_nbi, h_nbj, h_slot;
    std::vector<Real> h_area, h_xf, h_yf, h_zf;
    m_domain_cells = domain.numPts();
 
    // Bounding box of the solid columns (global), so the column arrays below
    // cover the built area only: 8 bytes per built column per rank.
    int bi0 = std::numeric_limits<int>::max(), bi1 = std::numeric_limits<int>::min();
    int bj0 = bi0, bj1 = bi1;
    for (MFIter mfi(blanking); mfi.isValid(); ++mfi) {
        const Box& bx = mfi.validbox();
        const HostFab hf(blanking[mfi]);
        auto const& b = hf.array();
        const auto lo = lbound(bx);
        const auto hi = ubound(bx);
        for (int k = lo.z; k <= hi.z; ++k) {
        for (int j = lo.y; j <= hi.y; ++j) {
        for (int i = lo.x; i <= hi.x; ++i) {
            if (b(i, j, k) >= 0.5) {
                bi0 = std::min(bi0, i); bi1 = std::max(bi1, i);
                bj0 = std::min(bj0, j); bj1 = std::max(bj1, j);
            }
        }}}
    }
    ParallelDescriptor::ReduceIntMin(bi0); ParallelDescriptor::ReduceIntMin(bj0);
    ParallelDescriptor::ReduceIntMax(bi1); ParallelDescriptor::ReduceIntMax(bj1);
    const bool any_solid = (bi1 >= bi0) && (bj1 >= bj0);
    m_col_i0 = any_solid ? bi0 - ilo : 0;  m_col_j0 = any_solid ? bj0 - jlo : 0;
    m_col_nx = any_solid ? bi1 - bi0 + 1 : 0;  m_col_ny = any_solid ? bj1 - bj0 + 1 : 0;
    const int bw = m_col_nx, bh = m_col_ny;
    const size_t ncol = static_cast<size_t>(bw) * bh;
    // Highest solid cell of each built column of this rank, as a domain-
    // relative k index, -1 where the column carries no solid cell; reduced
    // below so every rank sees the same columns. This one array is both the
    // solid-column mask of the labelling (a column is solid exactly when its
    // top is not -1) and the column top of the ray cast, and it is an index
    // rather than a height so that it costs 4 bytes per column instead of the
    // 24 a separate mask, index and height in Real took.
    std::vector<int> coltop_k(ncol, -1);
 
    m_fab_start.clear();
    for (MFIter mfi(blanking); mfi.isValid(); ++mfi) {
        m_fab_start.push_back(static_cast<int>(h_i.size()));
        const Box& bx = mfi.validbox();
        const HostFab hf(blanking[mfi]);
        auto const& b = hf.array();
        const Box& fbox = Box(b);
        const auto lo = lbound(bx);
        const auto hi = ubound(bx);
        for (int k = lo.z; k <= hi.z; ++k) {
        for (int j = lo.y; j <= hi.y; ++j) {
        for (int i = lo.x; i <= hi.x; ++i) {
            if (b(i, j, k) >= 0.5) {
                const size_t c = static_cast<size_t>(i - bi0) * bh + (j - bj0);
                coltop_k[c] = std::max(coltop_k[c], k - domain.smallEnd(2));
                continue;
            }
            // Slot of each face within its cell's checkpoint record: the
            // faces of a cell are numbered in the (dir, side) order of the
            // loops below, which the rebuild from the same blanking
            // reproduces on restart whatever the rank count.
            int cell_slot = 0;
            for (int d = 0; d < 3; ++d) {
                for (int s = -1; s <= 1; s += 2) {
                    IntVect nb(i, j, k);
                    nb[d] += s;
                    if (!fbox.contains(nb)) { continue; }
                    if (!geom.isPeriodic(d) && !domain.contains(nb)) { continue; }
                    if (b(nb) < 0.5) { continue; }
                    h_i.push_back(i); h_j.push_back(j); h_k.push_back(k);
                    h_dir.push_back(d); h_side.push_back(s);
                    h_slot.push_back(cell_slot++);
                    // Face centre: the fluid cell centre moved half a cell toward the solid.
                    Real c[3] = { Real(plo[0] + (i + 0.5) * dx[0]), Real(plo[1] + (j + 0.5) * dx[1]), Real(plo[2] + (k + 0.5) * dx[2]) };
                    c[d] += 0.5 * s * dx[d];
                    h_xf.push_back(c[0]); h_yf.push_back(c[1]); h_zf.push_back(c[2]);
                    // Column of the solid neighbour, wrapped for periodic directions.
                    int ci = nb[0] - ilo, cj = nb[1] - jlo;
                    ci = (ci % nx + nx) % nx;
                    cj = (cj % ny + ny) % ny;
                    h_nbi.push_back(ci); h_nbj.push_back(cj);
                    h_area.push_back(face_area[d]);
                }
            }
        }}}
    }
    m_fab_start.push_back(static_cast<int>(h_i.size()));
    m_nface = static_cast<int>(h_i.size());
 
    // Width of the checkpoint record: the largest number of faces on one cell
    // anywhere on the level (see state_ncomp()), reduced so every rank lays
    // the field out the same way.
    m_nslots = 0;
    for (int n = 0; n < m_nface; ++n) { m_nslots = std::max(m_nslots, h_slot[n] + 1); }
    ParallelDescriptor::ReduceIntMax(m_nslots);
 
    // Checkpoint layout: 4 x 4 column blocks clipped to the k-range that owns
    // faces (see state_boxarray()), the same on every rank from a reduced
    // block k-range map; and the transfer layer, this rank's grids cut by
    // those blocks, gathered into one global BoxArray owned by the rank whose
    // grid each piece came from.
    {
        constexpr int BS = 4;
        const int klo = domain.smallEnd(2);
        const int nbx = (nx + BS - 1) / BS, nby = (ny + BS - 1) / BS;
        // Lowest and highest face-owning cell of each block. Face-owning cells
        // are the fluid cells against solid, a one-cell shell, so a block over
        // the interior of a wide footprint owns faces only in the layer above
        // the roof and its box collapses to that layer instead of carrying the
        // whole column from the ground.
        std::vector<int> blk_lo(static_cast<size_t>(nbx) * nby, std::numeric_limits<int>::max());
        std::vector<int> blk_hi(static_cast<size_t>(nbx) * nby, -1);
        for (int n = 0; n < m_nface; ++n) {
            const size_t p = static_cast<size_t>((h_i[n] - ilo) / BS) * nby + (h_j[n] - jlo) / BS;
            blk_lo[p] = std::min(blk_lo[p], h_k[n] - klo);
            blk_hi[p] = std::max(blk_hi[p], h_k[n] - klo);
        }
        ParallelDescriptor::ReduceIntMin(blk_lo.data(), static_cast<int>(blk_lo.size()));
        ParallelDescriptor::ReduceIntMax(blk_hi.data(), static_cast<int>(blk_hi.size()));
        BoxList bl;
        for (int bi = 0; bi < nbx; ++bi) {
            for (int bj = 0; bj < nby; ++bj) {
                const size_t p = static_cast<size_t>(bi) * nby + bj;
                const int kmax = blk_hi[p];
                if (kmax < 0) { continue; }
                const int kmin = blk_lo[p];
                const IntVect blo(ilo + bi * BS, jlo + bj * BS, klo + kmin);
                const IntVect bhi(std::min(ilo + (bi + 1) * BS - 1, domain.bigEnd(0)),
                                  std::min(jlo + (bj + 1) * BS - 1, domain.bigEnd(1)), klo + kmax);
                bl.push_back(Box(blo, bhi));
            }
        }
        m_state_ba = BoxArray(bl);
        Vector<Box> mybx;
        m_xfer_src.clear();
        for (MFIter mfi(blanking); mfi.isValid(); ++mfi) {
            for (const auto& pr : m_state_ba.intersections(mfi.validbox())) {
                mybx.push_back(pr.second);
                m_xfer_src.push_back(mfi.LocalIndex());
            }
        }
        const int nranks = ParallelDescriptor::NProcs();
        const int nmine = static_cast<int>(mybx.size());
        Vector<int> counts(nranks, 0);
        if (nranks == 1) {
            // The non-MPI implementation of AllGather leaves its output
            // untouched, so populate the single-rank count explicitly.
            counts[0] = nmine;
        } else {
            ParallelAllGather::AllGather(nmine, counts.data(), ParallelDescriptor::Communicator());
        }
        Vector<Box> allbx(mybx);
        AllGatherBoxes(allbx);
        Vector<int> owner;
        owner.reserve(allbx.size());
        for (int r = 0; r < nranks; ++r) { for (int n = 0; n < counts[r]; ++n) { owner.push_back(r); } }
        AMREX_ALWAYS_ASSERT(owner.size() == allbx.size());
        BoxList xl;
        for (const Box& b : allbx) { xl.push_back(b); }
        m_xfer_ba = BoxArray(xl);
        m_xfer_dm = DistributionMapping(owner);
    }
 
    // Buildings: 4-connected solid columns, numbered in scan order. The column
    // mask is reduced so every rank labels the same columns with the same ids;
    // the labelling itself is a plain depth-first flood fill on the host.
    if (ncol > 0) {
        ParallelDescriptor::ReduceIntMax(coltop_k.data(), static_cast<int>(ncol));
    }
    // Column index in the built box of a face's solid neighbour (a 0-based,
    // wrapped domain column; the neighbour is solid, so it lies in the box).
    auto col_index = [&] (int n) {
        const size_t c = static_cast<size_t>(h_nbi[n] - m_col_i0) * bh + (h_nbj[n] - m_col_j0);
        AMREX_ALWAYS_ASSERT(h_nbi[n] >= m_col_i0 && h_nbi[n] < m_col_i0 + bw && h_nbj[n] >= m_col_j0 && h_nbj[n] < m_col_j0 + bh);
        return c;
    };
    // Column tops for the ray cast: the index map itself, which
    // ibseb::column_top() turns into a height one cell above the highest
    // solid cell, the ground where the column is fluid.
    upload(d_col_top, coltop_k);
    m_col_top_max = plo[2];
    for (size_t c = 0; c < ncol; ++c) {
        if (coltop_k[c] >= 0) {
            m_col_top_max = std::max(m_col_top_max, Real(plo[2] + (coltop_k[c] + 1.0) * dx[2]));
        }
    }
    // Height of each face's building column, for the wall-function depth.
    {
        std::vector<Real> h_hb(h_i.size());
        for (size_t n = 0; n < h_i.size(); ++n) {
            const int kt = coltop_k[col_index(n)];
            h_hb[n] = (kt >= 0) ? (kt + 1.0) * dx[2] : 0.0;
        }
        upload(d_hbld, h_hb);
    }
    std::vector<int> h_bid(m_nface);
    {
    std::vector<int> label(ncol, 0);
    std::vector<int> stack;
    m_nbld = 0;
    for (int ci = 0; ci < bw; ++ci) {
        for (int cj = 0; cj < bh; ++cj) {
            const size_t p = static_cast<size_t>(ci) * bh + cj;
            if (coltop_k[p] < 0 || label[p] > 0) { continue; }
            ++m_nbld;
            label[p] = m_nbld;
            stack.clear();
            stack.push_back(static_cast<int>(p));
            while (!stack.empty()) {
                const int q = stack.back(); stack.pop_back();
                const int qi = q / bh, qj = q % bh;
                const int ni[4] = {qi - 1, qi + 1, qi, qi};
                const int nj[4] = {qj, qj, qj - 1, qj + 1};
                for (int n = 0; n < 4; ++n) {
                    if (ni[n] < 0 || ni[n] >= bw || nj[n] < 0 || nj[n] >= bh) { continue; }
                    const size_t r = static_cast<size_t>(ni[n]) * bh + nj[n];
                    if (coltop_k[r] >= 0 && label[r] == 0) {
                        label[r] = m_nbld;
                        stack.push_back(static_cast<int>(r));
                    }
                }
            }
        }
    }
    for (int n = 0; n < m_nface; ++n) {
        h_bid[n] = label[col_index(n)];
    }
    // Footprints, for the debug summary: columns and bounding box per building.
    m_bld_ncol.assign(m_nbld + 1, 0);
    m_bld_ilo.assign(m_nbld + 1, nx); m_bld_ihi.assign(m_nbld + 1, -1);
    m_bld_jlo.assign(m_nbld + 1, ny); m_bld_jhi.assign(m_nbld + 1, -1);
    for (int ci = 0; ci < bw; ++ci) {
        for (int cj = 0; cj < bh; ++cj) {
            const int b = label[static_cast<size_t>(ci) * bh + cj];
            if (b == 0) { continue; }
            const int gi = ci + m_col_i0 + ilo, gj = cj + m_col_j0 + jlo;
            m_bld_ncol[b] += 1;
            m_bld_ilo[b] = std::min(m_bld_ilo[b], gi); m_bld_ihi[b] = std::max(m_bld_ihi[b], gi);
            m_bld_jlo[b] = std::min(m_bld_jlo[b], gj); m_bld_jhi[b] = std::max(m_bld_jhi[b], gj);
        }
    }
    }   // label and stack, one int per column each, are freed here
    // The column map has done its work too; free it before the per-face
    // device uploads below.
    std::vector<int>().swap(coltop_k);
 
    // Static totals: per direction, per building, and the summed area.
    std::vector<Long> nd(3, 0);
    std::vector<Long> bn(m_nbld + 1, 0);
    std::vector<Real> ba(m_nbld + 1, 0.0);
    Real area = 0.0;
    for (int n = 0; n < m_nface; ++n) {
        nd[h_dir[n]] += 1;
        bn[h_bid[n]] += 1;
        ba[h_bid[n]] += h_area[n];
        area += h_area[n];
    }
    ParallelDescriptor::ReduceLongSum(nd.data(), 3);
    ParallelDescriptor::ReduceLongSum(bn.data(), m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(ba.data(), m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(area);
    for (int d = 0; d < 3; ++d) { m_nface_dir[d] = nd[d]; }
    m_bld_nface = bn;
    m_bld_area  = ba;
    m_area_total = area;
 
    // Device arrays. The state starts uniform at T_skin_init (the slab too;
    // the interior boundary is applied when the slab is advanced); the view
    // fractions and fluxes start at zero and are filled by
    // compute_view_fractions() and the per-step routines.
    const size_t nf = static_cast<size_t>(m_nface);
    upload(d_i, h_i); upload(d_j, h_j); upload(d_k, h_k);
    upload(d_dir, h_dir); upload(d_side, h_side); upload(d_bid, h_bid);
    upload(d_slot, h_slot);
    upload(d_area, h_area);
    upload(d_xf, h_xf); upload(d_yf, h_yf); upload(d_zf, h_zf);
    fill(d_mat,      nf, 0);
    fill(d_albedo,   nf, m_params.albedo);
    fill(d_emis,     nf, m_params.emissivity);
    fill(d_kth,      nf, m_params.k_therm);
    fill(d_rhocp,    nf, m_params.rho_cp);
    fill(d_thick,    nf, m_params.thickness);
    fill(d_T_skin,   nf, m_params.T_skin_init);
    fill(d_T_slab,   nf * static_cast<size_t>(n_layers()), m_params.T_skin_init);
    // View fractions: placeholders (a roof sees the whole sky, a wall half sky
    // and half ground) that compute_view_fractions() replaces at initialisation.
    {
        std::vector<Real> fs(nf), fg(nf), fb(nf, 0.0);
        for (size_t n = 0; n < nf; ++n) {
            const bool roof = (h_dir[n] == 2);
            fs[n] = roof ? 1.0 : 0.5;
            fg[n] = roof ? 0.0 : 0.5;
        }
        upload(d_f_sky, fs); upload(d_f_ground, fg); upload(d_f_bldg, fb);
    }
    fill(d_shadow,        nf, Real(0.0));
    fill(d_SW_direct_in,  nf, Real(0.0));
    fill(d_SW_diffuse_in, nf, Real(0.0));
    fill(d_LW_down_in,    nf, Real(0.0));
    fill(d_T_air,         nf, m_params.T_skin_init);
    fill(d_theta_air,     nf, m_params.T_skin_init);
    fill(d_rho,           nf, Real(1.0));
    fill(d_U_tan,         nf, Real(0.0));
    fill(d_ustar,         nf, Real(0.0));
    fill(d_SW_abs,   nf, Real(0.0));
    fill(d_LW_net,   nf, Real(0.0));
    fill(d_H,        nf, Real(0.0));
    fill(d_LE,       nf, Real(0.0));
    fill(d_G,        nf, Real(0.0));
    fill(d_Q_ext,    nf, m_params.Q_ext_uniform);
    fill(d_H_coeff,  nf, Real(0.0));
    fill(d_w_star,   nf, Real(0.0));
    fill(d_olen,     nf, Real(1.0e30));
    fill(d_z_i,      nf, Real(0.0));
    fill(d_LW_ext,   nf, Real(0.0));
    fill(d_resid,    nf, Real(0.0));
    fill(d_niter,    nf, 0);
    Gpu::streamSynchronize();
 
    if (m_params.debug) { print_debug_summary(); }
}

Here is the call graph for this function:

compute_ground()

void IBFaceSet::compute_ground

(

amrex::Real

dt

)

Ground heat flux of the current step: advance every face's slab by dt with the skin temperature at the top and the interior temperature at the bottom (implicit, ERF_IBSEBSlab.H) and store the conduction into the slab, G, positive into the wall.

Slab conduction on every face: one kernel over the faces, each solving its own tridiagonal system with the skin at the top and the interior at the bottom.

{
    const int  nl   = n_layers();
    const Real Tint = m_params.T_interior;
    const Real* pT  = d_T_skin.data();
    const Real* pkt = d_kth.data(); const Real* prc = d_rhocp.data(); const Real* pth = d_thick.data();
    Real* pS = d_T_slab.data();   Real* pG = d_G.data();
    ParallelFor(m_nface, [=] AMREX_GPU_DEVICE (int f) noexcept {
        Real T[ibseb::SLAB_MAX_LAYERS];
        for (int l = 0; l < nl; ++l) { T[l] = pS[f * nl + l]; }
        const Real dz = pth[f] / nl;
        pG[f] = ibseb::advance_slab_dirichlet(T, pT[f], Tint, pkt[f], prc[f], dz, dt, nl);
        for (int l = 0; l < nl; ++l) { pS[f * nl + l] = T[l]; }
    });
    Gpu::streamSynchronize();
}

Here is the call graph for this function:

compute_longwave()

void IBFaceSet::compute_longwave

(

const amrex::MultiFab &

cons

)

Longwave of the current step. Incoming on a face: f_sky LW_sky + f_ground eps_g sigma T_g^4 + f_bldg sigma T_skin^4, the sky term fixed or gray (sky_emissivity sigma T_air^4 with the air temperature of the face's fluid cell), the building term the isothermal-surroundings approximation. Net: eps (LW_in - sigma T_skin^4), positive into the face.

Parameters

cons	Conserved state of the level, for the air temperature.

Longwave of the current step. One kernel per fab over that fab's faces, reading the air temperature of the fluid cell from the conserved state (potential temperature and density through the equation of state).

{
    const Real  sigma  = ibseb::SIGMA_SB;
    const bool  gray   = (m_params.lw_mode == "gray");
    const Real  lw_fix = m_params.lw_down, eps_sky = m_params.sky_emissivity;
    const Real  eps_g  = m_params.emissivity_ground;
    const Real* peps   = d_emis.data();
    const Real  Tg4    = std::pow(m_params.T_ground, 4);
    const int*  pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const Real* pfs = d_f_sky.data(); const Real* pfg = d_f_ground.data(); const Real* pfb = d_f_bldg.data();
    const Real* pT = d_T_skin.data();
    Real* pTa = d_T_air.data(); Real* pin = d_LW_down_in.data(); Real* pnet = d_LW_net.data();
    Real* pext = d_LW_ext.data();
    for (MFIter mfi(cons); mfi.isValid(); ++mfi) {
        const int f0 = m_fab_start[mfi.LocalIndex()];
        const int f1 = m_fab_start[mfi.LocalIndex() + 1];
        auto const& c = cons.const_array(mfi);
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f = f0 + m;
            const Real rho = c(pi[f], pj[f], pk[f], Rho_comp);
            const Real Ta  = getTgivenRandRTh(rho, c(pi[f], pj[f], pk[f], RhoTheta_comp));
            const Real lw_sky = gray ? eps_sky * sigma * Ta * Ta * Ta * Ta : lw_fix;
            const Real Ts4 = pT[f] * pT[f] * pT[f] * pT[f];
            // Sky and ground do not depend on the face; the wall term does
            // (isothermal surroundings), so the balance keeps them apart.
            const Real lw_ext = pfs[f] * lw_sky + pfg[f] * eps_g * sigma * Tg4;
            const Real lw_in  = lw_ext + pfb[f] * sigma * Ts4;
            pTa[f]  = Ta;
            pext[f] = lw_ext;
            pin[f]  = lw_in;
            pnet[f] = peps[f] * (lw_in - sigma * Ts4);
        });
    }
    Gpu::streamSynchronize();
}

Here is the call graph for this function:

compute_sensible()

void IBFaceSet::compute_sensible	(	const amrex::MultiFab &	cons,
		const amrex::MultiFab &	xvel,
		const amrex::MultiFab &	yvel,
		const amrex::MultiFab &	zvel,
		amrex::Real	c_p,
		const amrex::MultiFab *	olen_ground = nullptr,
		const amrex::MultiFab *	pblh_ground = nullptr,
		amrex::Real	z_i_bulk = 0.0
	)

Sensible heat of the current step through a wall function on every face: the tangential wind of the fluid cell at half a cell from the wall gives u* with the roughness z0_wall, the skin-to-air potential-temperature difference gives theta* with z0h_wall, and H = rho c_p u* theta*, positive out of the face. Neutral, or with the surface layer's stability functions on roofs when asked. The latent flux is not modelled and stays zero.

Beyond neutral: with convective_velocity = deardorff the wind is sqrt(U_tan^2 + (beta w*)^2) with w* from the previous step's H (zero when H is into the face) and a depth that is the mixed layer above a roof or the building height for a wall; with stability_correction the roofs iterate the similarity functions on the face's own Obukhov length, seeded from the ground surface layer's 2D field at the face's column when given.

Fills d_H, d_ustar, d_U_tan, d_rho, d_theta_air, d_H_coeff (the coefficient of the skin-to-air potential-temperature difference that the balance keeps frozen), d_w_star, d_olen and d_z_i; T_air was read by compute_longwave(), which must run first in a step.

Parameters

cons	Conserved state of the level (density, rho theta).
xvel,yvel,zvel	Face-centred velocities of the level.
c_p	Specific heat of air [J/kg/K].
olen_ground	Ground surface layer's Obukhov length (2D field on the level's box array collapsed to k = 0), or null.
pblh_ground	Its boundary-layer height field, or null.
z_i_bulk	Mixed-layer depth of the level from the bulk Richardson diagnostic or the fixed input [m].

Wall function on every face. One kernel per fab over that fab's faces: the cell-centred velocity of the fluid cell from its face values, its tangential part with respect to the wall, the neutral log law for u* and theta*, and on request two additions: a convective velocity scale in the wind (Beljaars' gustiness form with Deardorff's w* from the previous step's flux) and, on roofs, the surface layer's similarity functions iterated on the face's own Obukhov length. The skin temperature is converted to a potential temperature with the fluid cell's Exner function before the difference is taken.

The Obukhov length of a face is L = u*^2 theta / (kappa g theta*) (negative over a hot face). Walls stay on the log law even with the stability correction: the similarity functions assume a horizontal surface, and the convective scale is what carries free convection on a wall. The ground's 2D Obukhov field only seeds the roof iteration.

{
    const Real z0   = m_params.z0_wall, z0h = m_params.z0h_wall;
    const bool stab = m_params.stability_correction;
    const bool conv = (m_params.convective_velocity == "deardorff");
    const bool seed = stab && (olen_ground != nullptr) && (m_params.obukhov_seed == "ground");
    const bool use_pblh = conv && (pblh_ground != nullptr) && (m_params.z_i_mode == "pblh");
    const Real beta = m_params.beta_conv;
    const Real relax = m_params.obukhov_relax;
    const Real zi_fallback = (z_i_bulk > 0.0) ? z_i_bulk : m_params.z_i;
    const Real dxa[3] = {m_dx[0], m_dx[1], m_dx[2]};
    const Real z_lo = m_z_ground;
    const int*  pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const int*  pd = d_dir.data();
    const Real* pT = d_T_skin.data(); const Real* phb = d_hbld.data(); const Real* pzf = d_zf.data();
    Real* pth = d_theta_air.data(); Real* prho = d_rho.data(); Real* pU = d_U_tan.data();
    Real* pus = d_ustar.data();     Real* pH = d_H.data();     Real* pTa = d_T_air.data();
    Real* pHc = d_H_coeff.data();   Real* pws = d_w_star.data(); Real* pol = d_olen.data(); Real* pzi = d_z_i.data();
    for (MFIter mfi(cons); mfi.isValid(); ++mfi) {
        const int f0 = m_fab_start[mfi.LocalIndex()];
        const int f1 = m_fab_start[mfi.LocalIndex() + 1];
        auto const& c = cons.const_array(mfi);
        auto const& u = xvel.const_array(mfi);
        auto const& v = yvel.const_array(mfi);
        auto const& w = zvel.const_array(mfi);
        // The surface layer's 2D fields share the level's distribution with
        // the boxes collapsed to k = 0, so the same iterator indexes them.
        Array4<const Real> ol = seed     ? olen_ground->const_array(mfi) : Array4<const Real>{};
        Array4<const Real> pb = use_pblh ? pblh_ground->const_array(mfi) : Array4<const Real>{};
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f = f0 + m;
            const int i = pi[f], j = pj[f], k = pk[f], d = pd[f];
            const Real rho = c(i, j, k, Rho_comp);
            const Real rth = c(i, j, k, RhoTheta_comp);
            const Real th  = rth / rho;
            const Real Ta  = getTgivenRandRTh(rho, rth);
            const Real th_skin = pT[f] * th / Ta;            // skin temperature as a potential temperature
            // Cell-centred velocity and its tangential part with respect to the wall.
            Real U[3] = { Real(0.5) * (u(i, j, k) + u(i + 1, j, k)),
                          Real(0.5) * (v(i, j, k) + v(i, j + 1, k)),
                          Real(0.5) * (w(i, j, k) + w(i, j, k + 1)) };
            U[d] = 0.0;
            const Real Ut = amrex::max(std::sqrt(U[0] * U[0] + U[1] * U[1] + U[2] * U[2]), Real(1.0e-3));
            // Convective velocity scale from the previous step's flux out of the face.
            Real wstar = 0.0, depth = 0.0;
            if (conv && pH[f] > 0.0) {
                if (d == 2) {
                    Real zi = use_pblh ? pb(i, j, 0) : zi_fallback;
                    if (!(zi > 0.0 && zi < 1.0e5)) { zi = zi_fallback; }
                    depth = amrex::max(zi - (pzf[f] - z_lo), phb[f]);
                } else {
                    depth = phb[f];
                }
                wstar = std::cbrt(CONST_GRAV / th * pH[f] / (rho * c_p) * depth);
            }
            const Real Ut_eff = std::sqrt(Ut * Ut + beta * beta * wstar * wstar);
            const Real delta = 0.5 * dxa[d];
            const Real lnm = std::log(delta / z0), lnh = std::log(delta / z0h);
            Real ustar = KAPPA * Ut_eff / lnm;
            Real thstar = KAPPA * (th - th_skin) / lnh;
            Real lnh_eff = lnh;
            Real olen = 1.0e30;
            if (stab && d == 2) {
                // Roofs: the surface layer's similarity functions on the
                // face's own Obukhov length, a few fixed-point passes, seeded
                // from the ground's field at this column when available.
                similarity_funs sf;
                Real L = 1.0e30;
                if (seed) {
                    const Real Lg = ol(i, j, 0);
                    if (std::abs(Lg) > 0.0 && std::abs(Lg) < 1.0e5) { L = Lg; }
                }
                if (!(std::abs(L) < 1.0e5)) {
                    L = ustar * ustar * th / (KAPPA * CONST_GRAV * thstar + 1.0e-20);
                }
                // Fixed point on L until u*, theta* and L agree to 1e-6, so
                // the stored values are mutually consistent (a fixed number
                // of passes leaves L one pass behind). u* is under-relaxed
                // between passes, as the surface layer's iteration is in
                // erf-model #3486, which keeps a strongly unstable roof from
                // oscillating at low wind; the converged values do not
                // depend on the factor.
                Real us_it = ustar;
                for (int it = 0; it < 50; ++it) {
                    const Real zeta = delta / L;
                    const Real psi_m = sf.calc_psi_m(zeta), psi_h = sf.calc_psi_h(zeta);
                    const Real us_new = KAPPA * Ut_eff / amrex::max(lnm - psi_m, Real(0.1));
                    us_it   = (it == 0) ? us_new : (1.0 - relax) * us_it + relax * us_new;
                    lnh_eff = amrex::max(lnh - psi_h, Real(0.1));
                    thstar  = KAPPA * (th - th_skin) / lnh_eff;
                    const Real L_new = us_it * us_it * th / (KAPPA * CONST_GRAV * thstar + 1.0e-20);
                    const bool done = std::abs(L_new - L) <= 1.0e-6 * std::abs(L_new) &&
                                      std::abs(us_new - us_it) <= 1.0e-6 * us_new;
                    L = L_new;
                    if (done) { us_it = us_new; break; }
                }
                ustar = us_it;
                olen  = L;
            }
            pth[f]  = th;
            pTa[f]  = Ta;
            prho[f] = rho;
            pU[f]   = Ut;
            pus[f]  = ustar;
            pws[f]  = wstar;
            pol[f]  = olen;
            pzi[f]  = depth;
            pHc[f]  = rho * c_p * KAPPA * ustar / lnh_eff;   // H = H_coeff (theta_skin - theta_air)
            pH[f]   = -rho * c_p * ustar * thstar;           // positive out of the face
        });
    }
    Gpu::streamSynchronize();
}

compute_shortwave()

void IBFaceSet::compute_shortwave

(

amrex::Real

time

)

Shortwave of the current step: sun position and irradiances from the prescribed provider, then per face the direct beam on the face (zero when the ray toward the sun hits a building), the diffuse light through the view fractions, and the absorbed sum.

Fills d_shadow, d_SW_direct_in, d_SW_diffuse_in and d_SW_abs; stores the sun vector and irradiances for the reports.

Parameters

time	Simulation time [s]; with sun_mode = solar it is added to time_zero_utc_s.

Shortwave of the current step. The sun and the irradiances are scalars computed on the host from the provider inputs; the per-face work is one kernel: incidence cosine, ray cast against the column tops, direct and diffuse incident, absorbed.

With sun_mode = fixed the direct-normal irradiance and the horizontal diffuse are the inputs as given. With sun_mode = solar the sun follows the site and time and the clear-sky formulas give both irradiances; the sun below the horizon gives zero everywhere.

Diffuse on a face: f_sky * diffuse_h + f_ground * albedo_ground * (direct_h + diffuse_h), the second term being the ground-reflected part, with the direct on a horizontal surface direct_h = dni * cos z. The view fractions are the placeholders of build() until compute_view_fractions() has run.

{
    // ---- Sun and irradiances (host scalars) ----
    SunState s;
    if (m_params.sun_mode == "fixed") {
        s.zenith  = m_params.sun_zenith_deg  * PI / 180.0;
        s.azimuth = m_params.sun_azimuth_deg * PI / 180.0;
        const Real cz = std::cos(s.zenith);
        s.dni       = (cz > 0.0) ? m_params.sw_direct_normal : 0.0;
        s.diffuse_h = (cz > 0.0) ? m_params.sw_diffuse : 0.0;
    } else {
        const Real t_utc = m_params.time_zero_utc_s + time;
        const Real decl  = ibseb::solar_declination(m_params.day_of_year);
        // Hour angle from UTC and the longitude (a time-zone offset would cancel).
        const Real ha    = ibseb::solar_hour_angle(t_utc, m_params.longitude_deg, m_params.day_of_year);
        s.zenith  = ibseb::solar_zenith(m_params.latitude_deg, decl, ha);
        s.azimuth = ibseb::solar_azimuth(m_params.latitude_deg, decl, ha, s.zenith);
        const Real cz = std::cos(s.zenith);
        const Real e0 = ibseb::earth_sun_distance_factor(m_params.day_of_year);
        s.dni       = ibseb::clear_sky_dni(cz, m_params.solar_constant, m_params.sw_transmission, e0);
        s.diffuse_h = ibseb::clear_sky_diffuse_h(cz, m_params.solar_constant, m_params.sw_transmission,
                                                 e0, m_params.sw_diffuse_coeff);
    }
    ibseb::sun_vector(s.zenith, s.azimuth, s.sx, s.sy, s.sz);
    m_sun = s;
 
    // ---- Per face ----
    const Real sx = s.sx, sy = s.sy, sz = s.sz;
    const Real dni = s.dni, dif_h = s.diffuse_h;
    const Real dir_h = dni * amrex::max(Real(0.0), sz);
    const Real alb_g = m_params.albedo_ground;
    const Real* palb = d_albedo.data();
    const int*  col_top = d_col_top.data();
    const int   nx = m_nx, ny = m_ny;
    const int   ci0 = m_col_i0, cj0 = m_col_j0, cbw = m_col_nx, cbh = m_col_ny;
    const Real  x_lo = m_x_lo, y_lo = m_y_lo, dx = m_dx[0], dy = m_dx[1], dz = m_dx[2];
    const bool  per_x = m_per_x, per_y = m_per_y;
    const Real  z_ground = m_z_ground, z_max = m_col_top_max, max_path = m_max_path;
    const int*  pd  = d_dir.data();  const int* ps = d_side.data();
    const Real* pxf = d_xf.data();   const Real* pyf = d_yf.data();  const Real* pzf = d_zf.data();
    const Real* pfs = d_f_sky.data(); const Real* pfg = d_f_ground.data();
    Real* psh = d_shadow.data();  Real* pdir = d_SW_direct_in.data();
    Real* pdif = d_SW_diffuse_in.data();  Real* pabs = d_SW_abs.data();
    ParallelFor(m_nface, [=] AMREX_GPU_DEVICE (int f) noexcept {
        // Outward normal: opposite to the side the solid is on.
        Real n[3] = {0.0, 0.0, 0.0};
        n[pd[f]] = -static_cast<Real>(ps[f]);
        const Real cosi = n[0] * sx + n[1] * sy + n[2] * sz;
        Real shadow = 0.0, direct = 0.0;
        if (sz > 0.0 && cosi > 0.0 && dni > 0.0) {
            shadow = ibseb::ray_blocked(pxf[f], pyf[f], pzf[f], sx, sy, sz, col_top, nx, ny, ci0, cj0, cbw, cbh,
                                        x_lo, y_lo, dx, dy, dz, per_x, per_y, z_ground, z_max, max_path) ? 1.0 : 0.0;
            direct = dni * cosi * (1.0 - shadow);
        }
        const Real diffuse = pfs[f] * dif_h + pfg[f] * alb_g * (dir_h + dif_h);
        psh[f]  = shadow;
        pdir[f] = direct;
        pdif[f] = diffuse;
        pabs[f] = (1.0 - palb[f]) * (direct + diffuse);
    });
    Gpu::streamSynchronize();
 
    if (m_params.debug) {
        Print() << "[IBSEB DEBUG] lev=" << m_lev << " sun: zenith=" << s.zenith * 180.0 / PI
                << " deg azimuth=" << s.azimuth * 180.0 / PI << " deg s=(" << s.sx << "," << s.sy << "," << s.sz
                << ") DNI=" << s.dni << " W/m2 diffuse_h=" << s.diffuse_h << " W/m2\n";
    }
}

Here is the call graph for this function:

compute_view_fractions()

void IBFaceSet::compute_view_fractions

(

)

View fractions: sample a cosine-weighted hemisphere around every face's outward normal with view_n_az x view_n_el rays through the column walk and count where they end: sky, ground or a building. The three fractions sum to one and are view factors. Called once at initialisation; the geometry is static.

Hemisphere sampling for the view fractions. One kernel over the faces, each looping over its rays; the count of rays ending on the sky, the ground and a building over the total gives the three fractions. Roofs point up and never see the ground; a wall on flat open ground sees half sky and half ground, which the regression test checks.

{
    const int n_az = m_params.view_n_az, n_el = m_params.view_n_el;
    const int*  col_top = d_col_top.data();
    const int   nx = m_nx, ny = m_ny;
    const int   ci0 = m_col_i0, cj0 = m_col_j0, cbw = m_col_nx, cbh = m_col_ny;
    const Real  x_lo = m_x_lo, y_lo = m_y_lo, dx = m_dx[0], dy = m_dx[1], dz = m_dx[2];
    const bool  per_x = m_per_x, per_y = m_per_y;
    const Real  z_ground = m_z_ground, z_max = m_col_top_max, max_path = m_max_path;
    const int*  pd = d_dir.data();  const int* ps = d_side.data();
    const Real* pxf = d_xf.data();  const Real* pyf = d_yf.data();  const Real* pzf = d_zf.data();
    Real* pfs = d_f_sky.data(); Real* pfg = d_f_ground.data(); Real* pfb = d_f_bldg.data();
    ParallelFor(m_nface, [=] AMREX_GPU_DEVICE (int f) noexcept {
        int n_sky = 0, n_gnd = 0, n_bld = 0;
        for (int ie = 0; ie < n_el; ++ie) {
            for (int ia = 0; ia < n_az; ++ia) {
                Real sx, sy, sz;
                ibseb::hemisphere_direction(pd[f], -ps[f], ia, ie, n_az, n_el, sx, sy, sz);
                const int hit = ibseb::ray_hit(pxf[f], pyf[f], pzf[f], sx, sy, sz, col_top, nx, ny, ci0, cj0, cbw, cbh,
                                               x_lo, y_lo, dx, dy, dz, per_x, per_y, z_ground, z_max, max_path);
                if (hit == ibseb::RAY_SKY) ++n_sky; else if (hit == ibseb::RAY_GROUND) ++n_gnd; else ++n_bld;
            }
        }
        const Real inv = 1.0 / static_cast<Real>(n_az * n_el);
        pfs[f] = n_sky * inv; pfg[f] = n_gnd * inv; pfb[f] = n_bld * inv;
    });
    Gpu::streamSynchronize();
 
    if (m_params.debug) {
        std::vector<Real> fs(m_nface), fg(m_nface), fb(m_nface);
        Gpu::copy(Gpu::deviceToHost, d_f_sky.begin(), d_f_sky.end(), fs.begin());
        Gpu::copy(Gpu::deviceToHost, d_f_ground.begin(), d_f_ground.end(), fg.begin());
        Gpu::copy(Gpu::deviceToHost, d_f_bldg.begin(), d_f_bldg.end(), fb.begin());
        Gpu::streamSynchronize();
        Real s[3] = {0.0, 0.0, 0.0};
        for (int n = 0; n < m_nface; ++n) { s[0] += fs[n]; s[1] += fg[n]; s[2] += fb[n]; }
        ParallelDescriptor::ReduceRealSum(s, 3);
        const Real nf_all = static_cast<Real>(m_nface_dir[0] + m_nface_dir[1] + m_nface_dir[2]);
        Print() << "[IBSEB DEBUG] lev=" << m_lev << " view fractions: " << n_az * n_el
                << " rays per face, mean f_sky=" << s[0] / nf_all << " f_ground=" << s[1] / nf_all
                << " f_bldg=" << s[2] / nf_all << "\n";
    }
}

Here is the call graph for this function:

dump_faces()

void IBFaceSet::dump_faces

(

const std::string &

prefix

)

const

Write every face of this rank to <prefix>.rank<N>.csv with its geometry, view fractions, shadow flag, shortwave terms and skin temperature, for the regression tests and for users who want the raw list.

Per-rank face dump. Every rank writes its own file so no gather is needed; readers concatenate <prefix>.rank*.csv.

{
    const int nf = m_nface;
    auto host = [&](const Gpu::DeviceVector<Real>& d) {
        std::vector<Real> h(nf);
        Gpu::copy(Gpu::deviceToHost, d.begin(), d.end(), h.begin());
        return h;
    };
    auto hosti = [&](const Gpu::DeviceVector<int>& d) {
        std::vector<int> h(nf);
        Gpu::copy(Gpu::deviceToHost, d.begin(), d.end(), h.begin());
        return h;
    };
    const auto i = hosti(d_i), j = hosti(d_j), k = hosti(d_k), dir = hosti(d_dir), side = hosti(d_side), bid = hosti(d_bid);
    const auto xf = host(d_xf), yf = host(d_yf), zf = host(d_zf), area = host(d_area);
    const auto fs = host(d_f_sky), fg = host(d_f_ground), fb = host(d_f_bldg);
    const auto sh = host(d_shadow), sdir = host(d_SW_direct_in), sdif = host(d_SW_diffuse_in), sabs = host(d_SW_abs);
    const auto T = host(d_T_skin), Ta = host(d_T_air), lwin = host(d_LW_down_in), lwn = host(d_LW_net);
    const auto th = host(d_theta_air), rho = host(d_rho), Ut = host(d_U_tan), us = host(d_ustar), H = host(d_H);
    const auto mat = hosti(d_mat);
    const auto alb = host(d_albedo), emi = host(d_emis), kk = host(d_kth), rc = host(d_rhocp), thk = host(d_thick), G = host(d_G);
    const auto Hc = host(d_H_coeff), Qx = host(d_Q_ext), LE = host(d_LE), lwe = host(d_LW_ext), res = host(d_resid);
    const auto nit = hosti(d_niter);
    const auto ws = host(d_w_star), ol = host(d_olen), zi = host(d_z_i), hb = host(d_hbld);
    const int nl = n_layers();
    std::vector<Real> slab(static_cast<size_t>(nf) * nl);
    Gpu::copy(Gpu::deviceToHost, d_T_slab.begin(), d_T_slab.end(), slab.begin());
    Gpu::streamSynchronize();
    std::ofstream f(prefix + ".rank" + std::to_string(ParallelDescriptor::MyProc()) + ".csv");
    f << "i,j,k,dir,side,bid,x_m,y_m,z_m,area_m2,f_sky,f_ground,f_bldg,shadow,SW_direct_in,SW_diffuse_in,SW_abs,T_skin,T_air,LW_in,LW_net,"
         "theta_air,rho,U_tan,ustar,H,mat,albedo,emissivity,k_therm,rho_cp,thickness,G,T_slab_top,T_slab_bottom,"
         "H_coeff,Q_ext,LE,LW_ext,resid,n_iter,w_star,Olen,z_i,h_bld";
    for (int l = 0; l < nl; ++l) { f << ",T_slab_" << l; }
    f << "\n";
    f << std::setprecision(12);
    for (int n = 0; n < nf; ++n) {
        f << i[n] << "," << j[n] << "," << k[n] << "," << dir[n] << "," << side[n] << "," << bid[n] << ","
          << xf[n] << "," << yf[n] << "," << zf[n] << "," << area[n] << ","
          << fs[n] << "," << fg[n] << "," << fb[n] << "," << sh[n] << ","
          << sdir[n] << "," << sdif[n] << "," << sabs[n] << "," << T[n] << ","
          << Ta[n] << "," << lwin[n] << "," << lwn[n] << ","
          << th[n] << "," << rho[n] << "," << Ut[n] << "," << us[n] << "," << H[n] << ","
          << mat[n] << "," << alb[n] << "," << emi[n] << "," << kk[n] << "," << rc[n] << "," << thk[n] << ","
          << G[n] << "," << slab[static_cast<size_t>(n) * nl] << "," << slab[static_cast<size_t>(n) * nl + nl - 1] << ","
          << Hc[n] << "," << Qx[n] << "," << LE[n] << "," << lwe[n] << "," << res[n] << "," << nit[n]
          << "," << ws[n] << "," << ol[n] << "," << zi[n] << "," << hb[n];
        for (int l = 0; l < nl; ++l) { f << "," << slab[static_cast<size_t>(n) * nl + l]; }
        f << "\n";
    }
}

Referenced by report().

Here is the call graph for this function:

Here is the caller graph for this function:

fab_start()

const std::vector<int>& IBFaceSet::fab_start ( ) const

inline

First face of each local fab in MFIter order; size number of fabs + 1. Faces fab_start()[n] to fab_start()[n+1]-1 belong to the fab with local index n, which is what the scatter kernels rely on.

{ return m_fab_start; }

has_state()

bool IBFaceSet::has_state ( ) const

inline

Whether the level has a checkpoint field to write (any faces at all).

{ return !m_state_ba.empty(); }

load_state()

void IBFaceSet::load_state

(

const amrex::MultiFab &

state

)

Read the face state (skin, flux, slab) back from a field written by save_state(), after build() has recreated the list from the blanking: the reverse ParallelCopy onto the transfer field, then each face reads its slot. Only the slots that correspond to a face of this list are read.

Read the skin temperature, the sensible flux and the slab temperatures back from the state field, the inverse of save_state(), into a list that build() has just recreated with the same slot numbering.

{
    AMREX_ALWAYS_ASSERT(state.nComp() == state_ncomp() && state.boxArray() == m_state_ba);
    const int nl = n_layers();
    const int* pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const int* pslot = d_slot.data();
    Real* pT = d_T_skin.data(); Real* pS = d_T_slab.data(); Real* pH = d_H.data();
    MultiFab mine(m_xfer_ba, m_xfer_dm, state_ncomp(), 0);
    mine.ParallelCopy(state);
    for (MFIter mfi(mine); mfi.isValid(); ++mfi) {
        const int li = m_xfer_src[mfi.LocalIndex()];
        const int f0 = m_fab_start[li];
        const int f1 = m_fab_start[li + 1];
        const Box bx = mfi.validbox();
        auto const& st = mine.const_array(mfi);
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f  = f0 + m;
            if (!bx.contains(IntVect(pi[f], pj[f], pk[f]))) { return; }
            const int c0 = pslot[f] * (2 + nl);
            pT[f] = st(pi[f], pj[f], pk[f], c0);
            pH[f] = st(pi[f], pj[f], pk[f], c0 + 1);
            for (int l = 0; l < nl; ++l) { pS[f * nl + l] = st(pi[f], pj[f], pk[f], c0 + 2 + l); }
        });
    }
    Gpu::streamSynchronize();
    if (m_params.debug) {
        std::vector<Real> h_T(m_nface);
        Gpu::copy(Gpu::deviceToHost, d_T_skin.begin(), d_T_skin.end(), h_T.begin());
        Gpu::streamSynchronize();
        Real tmin = 1.0e30, tmax = -1.0e30;
        for (int n = 0; n < m_nface; ++n) { tmin = std::min(tmin, h_T[n]); tmax = std::max(tmax, h_T[n]); }
        ParallelDescriptor::ReduceRealMin(tmin);
        ParallelDescriptor::ReduceRealMax(tmax);
        Print() << "[IBSEB DEBUG] lev=" << m_lev << " face state loaded: T_skin in ["
                << tmin << ", " << tmax << "] K\n";
    }
}

Here is the call graph for this function:

make_state()

MultiFab IBFaceSet::make_state

(

)

const

A field on state_boxarray() with state_ncomp() components, for the checkpoint.

Write the skin temperature, the sensible flux and the slab temperatures into the state field. The record of a face starts at component slot * (2 + n_layers), the slot being the face's number within its own cell (see state_ncomp()); no two faces of a cell share a slot, so the writes are plain stores.

{
    if (m_state_ba.empty()) { return MultiFab(); }
    return MultiFab(m_state_ba, DistributionMapping(m_state_ba), state_ncomp(), 0);
}

Here is the call graph for this function:

n_buildings()

int IBFaceSet::n_buildings ( ) const

inline

Number of buildings of the level, the same on every rank.

{ return m_nbld; }

n_faces()

int IBFaceSet::n_faces ( ) const

inline

Number of faces owned by this rank (not the global count; see report()).

{ return m_nface; }

n_layers()

int IBFaceSet::n_layers ( ) const

inline

Conduction layers per face, from erf.ibseb.n_slab_layers.

{ return m_params.n_slab_layers; }

Referenced by dump_faces(), load_state(), print_debug_summary(), save_state(), and state_ncomp().

Here is the caller graph for this function:

n_slots()

int IBFaceSet::n_slots ( ) const

inline

Slots per cell of the checkpoint field: the largest face count on a cell of the level.

{ return m_nslots; }

print_debug_summary()

void IBFaceSet::print_debug_summary

(

)

const

Print the [IBSEB DEBUG] description of the set: the inputs, every rank's face count and per-fab ranges, the buildings with their footprint columns and bounding boxes, and the device memory of the arrays. Called once from build() when erf.ibseb.debug is set.

Debug description of the set, in the style of the fire module's [FIRE DEBUG] lines: what was read, what every rank holds, what the buildings look like, and what the arrays cost. Per-rank lines use AllPrint() so every rank reports; the rest is printed by the I/O rank.

{
    const int nl = n_layers();
    Print() << "[IBSEB DEBUG] lev=" << m_lev << " inputs: n_slab_layers=" << nl
            << " T_skin_init=" << m_params.T_skin_init << " K"
            << " T_interior=" << m_params.T_interior << " K"
            << " csv_file=" << m_params.csv_file << " csv_int=" << m_params.csv_int << "\n";
    Print() << "[IBSEB DEBUG] lev=" << m_lev << " global faces=" << (m_nface_dir[0] + m_nface_dir[1] + m_nface_dir[2])
            << " (x=" << m_nface_dir[0] << " y=" << m_nface_dir[1] << " z=" << m_nface_dir[2] << ")"
            << " area=" << m_area_total << " m2 buildings=" << m_nbld << "\n";
    for (int b = 1; b <= m_nbld; ++b) {
        Print() << "[IBSEB DEBUG]   building " << b << ": columns=" << m_bld_ncol[b]
                << " i=[" << m_bld_ilo[b] << "," << m_bld_ihi[b] << "]"
                << " j=[" << m_bld_jlo[b] << "," << m_bld_jhi[b] << "]"
                << " faces=" << m_bld_nface[b] << " area=" << m_bld_area[b] << " m2\n";
    }
    // Per-rank ownership: one line per rank with its fab ranges.
    const int nfab = static_cast<int>(m_fab_start.size()) - 1;
    const size_t bytes = static_cast<size_t>(m_nface) *
        (7 * sizeof(int) + (11 + static_cast<size_t>(nl)) * sizeof(Real));
    std::ostringstream os;
    os << "[IBSEB DEBUG] lev=" << m_lev << " rank " << ParallelDescriptor::MyProc()
       << ": faces=" << m_nface << " over " << nfab << " fab(s), device memory "
       << bytes / 1024.0 << " kB; fab ranges:";
    for (int f = 0; f < nfab; ++f) { os << " [" << m_fab_start[f] << "," << m_fab_start[f + 1] << ")"; }
    AllPrint() << os.str() << "\n";
}

Here is the call graph for this function:

report()

void IBFaceSet::report	(	amrex::Real	time,
		int	step,
		bool	write_csv
	)		const

Print the [IBSEB] summary line and, when asked, append one row per building to erf.ibseb.csv_file.

The summary carries the global face counts per direction, the number of buildings, the total face area and the range of the skin temperature; the CSV rows carry per building the face count, the area, the area-weighted mean fluxes and skin temperature with its range, the largest balance residual, and the sun of the step (zenith, azimuth, direct-normal and diffuse irradiance). All numbers are reduced over the ranks; the file is written by the I/O rank only, with a header only when the file does not exist yet, so a restarted run appends.

Parameters

time	Simulation time [s] written to the rows.
step	Step number written to the rows.
write_csv	Whether to append rows, or only print the summary.

Summary line and per-building CSV rows. The dynamic quantities (skin temperature range and area-weighted mean per building) are gathered on the host from the device arrays and reduced; the static counts and areas were reduced in build().

{
    // Skin temperature: min, max and per-building mean over all ranks.
    std::vector<Real> h_T(m_nface), h_A(m_nface);
    std::vector<int>  h_b(m_nface);
    Gpu::copy(Gpu::deviceToHost, d_T_skin.begin(), d_T_skin.end(), h_T.begin());
    Gpu::copy(Gpu::deviceToHost, d_area.begin(),   d_area.end(),   h_A.begin());
    Gpu::copy(Gpu::deviceToHost, d_bid.begin(),    d_bid.end(),    h_b.begin());
    Gpu::streamSynchronize();
 
    std::vector<Real> h_S(m_nface), h_sh(m_nface), h_L(m_nface), h_H(m_nface), h_G(m_nface), h_r(m_nface), h_Q(m_nface), h_w(m_nface);
    Gpu::copy(Gpu::deviceToHost, d_w_star.begin(), d_w_star.end(), h_w.begin());
    Gpu::copy(Gpu::deviceToHost, d_resid.begin(), d_resid.end(), h_r.begin());
    Gpu::copy(Gpu::deviceToHost, d_Q_ext.begin(), d_Q_ext.end(), h_Q.begin());
    Gpu::copy(Gpu::deviceToHost, d_SW_abs.begin(), d_SW_abs.end(), h_S.begin());
    Gpu::copy(Gpu::deviceToHost, d_shadow.begin(), d_shadow.end(), h_sh.begin());
    Gpu::copy(Gpu::deviceToHost, d_LW_net.begin(), d_LW_net.end(), h_L.begin());
    Gpu::copy(Gpu::deviceToHost, d_H.begin(), d_H.end(), h_H.begin());
    Gpu::copy(Gpu::deviceToHost, d_G.begin(), d_G.end(), h_G.begin());
    Gpu::streamSynchronize();
 
    std::vector<Real> bsum(m_nbld + 1, 0.0), bsw(m_nbld + 1, 0.0), bsh(m_nbld + 1, 0.0), blw(m_nbld + 1, 0.0), bH(m_nbld + 1, 0.0);
    std::vector<Real> bG(m_nbld + 1, 0.0), bQ(m_nbld + 1, 0.0), btmin(m_nbld + 1, 1.0e30), btmax(m_nbld + 1, -1.0e30), bres(m_nbld + 1, 0.0);
    Real tmin = 1.0e30, tmax = -1.0e30, sw_sum = 0.0, sw_max = 0.0, sh_area = 0.0, lw_sum = 0.0, H_sum = 0.0, G_sum = 0.0;
    Real res_max = 0.0, Q_sum = 0.0, w_max = 0.0;
    for (int n = 0; n < m_nface; ++n) {
        w_max = std::max(w_max, h_w[n]);
        G_sum += h_G[n] * h_A[n];
        Q_sum += h_Q[n] * h_A[n];
        res_max = std::max(res_max, h_r[n]);
        bG[h_b[n]]   += h_G[n] * h_A[n];
        bQ[h_b[n]]   += h_Q[n] * h_A[n];
        btmin[h_b[n]] = std::min(btmin[h_b[n]], h_T[n]);
        btmax[h_b[n]] = std::max(btmax[h_b[n]], h_T[n]);
        bres[h_b[n]]  = std::max(bres[h_b[n]], h_r[n]);
        bsum[h_b[n]] += h_T[n] * h_A[n];
        bsw[h_b[n]]  += h_S[n] * h_A[n];
        bsh[h_b[n]]  += h_sh[n] * h_A[n];
        blw[h_b[n]]  += h_L[n] * h_A[n];
        bH[h_b[n]]   += h_H[n] * h_A[n];
        tmin = std::min(tmin, h_T[n]);
        tmax = std::max(tmax, h_T[n]);
        sw_sum  += h_S[n] * h_A[n];
        sw_max   = std::max(sw_max, h_S[n]);
        sh_area += h_sh[n] * h_A[n];
        lw_sum  += h_L[n] * h_A[n];
        H_sum   += h_H[n] * h_A[n];
    }
    ParallelDescriptor::ReduceRealSum(bsum.data(), m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(bsw.data(),  m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(bsh.data(),  m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(blw.data(),  m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(bH.data(),   m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(bG.data(),   m_nbld + 1);
    ParallelDescriptor::ReduceRealSum(bQ.data(),   m_nbld + 1);
    ParallelDescriptor::ReduceRealMin(btmin.data(), m_nbld + 1);
    ParallelDescriptor::ReduceRealMax(btmax.data(), m_nbld + 1);
    ParallelDescriptor::ReduceRealMax(bres.data(),  m_nbld + 1);
    ParallelDescriptor::ReduceRealMax(res_max);
    ParallelDescriptor::ReduceRealMax(w_max);
    ParallelDescriptor::ReduceRealSum(Q_sum);
    ParallelDescriptor::ReduceRealMin(tmin);
    ParallelDescriptor::ReduceRealMax(tmax);
    ParallelDescriptor::ReduceRealSum(sw_sum);
    ParallelDescriptor::ReduceRealMax(sw_max);
    ParallelDescriptor::ReduceRealSum(sh_area);
    ParallelDescriptor::ReduceRealSum(lw_sum);
    ParallelDescriptor::ReduceRealSum(H_sum);
    ParallelDescriptor::ReduceRealSum(G_sum);
    const Real G_mean   = (m_area_total > 0.0) ? G_sum / m_area_total : 0.0;
    const Real sw_mean  = (m_area_total > 0.0) ? sw_sum / m_area_total : 0.0;
    const Real sh_frac  = (m_area_total > 0.0) ? sh_area / m_area_total : 0.0;
    const Real lw_mean  = (m_area_total > 0.0) ? lw_sum / m_area_total : 0.0;
    const Real H_mean   = (m_area_total > 0.0) ? H_sum / m_area_total : 0.0;
 
    const Long nf_all = m_nface_dir[0] + m_nface_dir[1] + m_nface_dir[2];
    Print() << "[IBSEB] lev=" << m_lev << " step=" << step << " t=" << time
            << " faces=" << nf_all
            << " (x=" << m_nface_dir[0] << " y=" << m_nface_dir[1] << " z=" << m_nface_dir[2] << ")"
            << " buildings=" << m_nbld
            << " area=" << std::setprecision(10) << m_area_total << " m2"
            << " T_skin_min=" << ((nf_all > 0) ? tmin : Real(0.0))
            << " T_skin_max=" << ((nf_all > 0) ? tmax : Real(0.0))
            << " SW_abs_mean=" << sw_mean << " SW_abs_max=" << sw_max
            << " shadow_frac=" << sh_frac << " LW_net_mean=" << lw_mean
            << " H_mean=" << H_mean << " G_mean=" << G_mean
            << " Q_ext_mean=" << ((m_area_total > 0.0) ? Q_sum / m_area_total : Real(0.0))
            << " resid_max=" << res_max
            << " w_star_max=" << w_max
            << " H_total_W=" << std::setprecision(12) << H_sum << "\n";
    // Cost line: the slowest rank's per-step time of the balance and its face
    // count, so the cost of the face list on a larger case can be estimated.
    {
        Real ms = (m_cost_n > 0) ? 1000.0 * m_cost_s / m_cost_n : 0.0;
        Real ms_max = ms; int nf_max = m_nface;
        ParallelDescriptor::ReduceRealMax(ms_max);
        ParallelDescriptor::ReduceIntMax(nf_max);
        Print() << "[IBSEB] cost: lev=" << m_lev << " advance_ms_per_step_max=" << std::setprecision(4) << ms_max
                << " faces_per_rank_max=" << nf_max << " ranks=" << ParallelDescriptor::NProcs()
                << " init_s=" << m_init_cost_s << "\n";
    }
 
    if (m_params.debug) {
        for (int b = 1; b <= m_nbld; ++b) {
            const Real tmean = (m_bld_area[b] > 0.0) ? bsum[b] / m_bld_area[b] : Real(0.0);
            const Real swm = (m_bld_area[b] > 0.0) ? bsw[b] / m_bld_area[b] : Real(0.0);
            const Real shf = (m_bld_area[b] > 0.0) ? bsh[b] / m_bld_area[b] : Real(0.0);
            const Real lwm = (m_bld_area[b] > 0.0) ? blw[b] / m_bld_area[b] : Real(0.0);
            const Real Hm  = (m_bld_area[b] > 0.0) ? bH[b]  / m_bld_area[b] : Real(0.0);
            const Real Gm  = (m_bld_area[b] > 0.0) ? bG[b]  / m_bld_area[b] : Real(0.0);
            Print() << "[IBSEB DEBUG]   building " << b << ": faces=" << m_bld_nface[b]
                    << " area=" << m_bld_area[b] << " m2 T_skin_mean=" << tmean << " K"
                    << " [" << btmin[b] << ", " << btmax[b] << "]"
                    << " SW_abs_mean=" << swm << " W/m2 shadow_frac=" << shf
                    << " LW_net_mean=" << lwm << " W/m2 H_mean=" << Hm << " W/m2 G_mean=" << Gm
                    << " W/m2 resid_max=" << bres[b] << " W/m2\n";
        }
    }
 
    if (write_csv && !m_params.dump_faces_file.empty()) {
        std::string prefix = m_params.dump_faces_file;
        if (m_params.dump_faces_tag_step) {
            std::ostringstream os;
            os << prefix << ".step" << std::setw(6) << std::setfill('0') << step;
            prefix = os.str();
        }
        dump_faces(prefix);
    }
 
    if (!write_csv || !ParallelDescriptor::IOProcessor()) { return; }
    bool need_header = true;
    {
        std::ifstream probe(m_params.csv_file, std::ios::ate);
        if (probe.good() && probe.tellg() > 0) { need_header = false; }
    }
    std::ofstream csv(m_params.csv_file, std::ios::app);
    if (need_header) {
        csv << "time_s,step,level,building,n_faces,area_m2,T_skin_mean_K,SW_abs_mean_Wm2,shadow_frac,LW_net_mean_Wm2,H_mean_Wm2,"
               "G_mean_Wm2,Q_ext_mean_Wm2,T_skin_min_K,T_skin_max_K,resid_max_Wm2,"
               "sun_zenith_deg,sun_azimuth_deg,dni_Wm2,diffuse_h_Wm2\n";
    }
    for (int b = 1; b <= m_nbld; ++b) {
        const Real tmean = (m_bld_area[b] > 0.0) ? bsum[b] / m_bld_area[b] : Real(0.0);
        const Real swm   = (m_bld_area[b] > 0.0) ? bsw[b]  / m_bld_area[b] : Real(0.0);
        const Real shf   = (m_bld_area[b] > 0.0) ? bsh[b]  / m_bld_area[b] : Real(0.0);
        const Real lwm   = (m_bld_area[b] > 0.0) ? blw[b]  / m_bld_area[b] : Real(0.0);
        const Real Hm    = (m_bld_area[b] > 0.0) ? bH[b]   / m_bld_area[b] : Real(0.0);
        const Real Gm    = (m_bld_area[b] > 0.0) ? bG[b]   / m_bld_area[b] : Real(0.0);
        const Real Qm    = (m_bld_area[b] > 0.0) ? bQ[b]   / m_bld_area[b] : Real(0.0);
        csv << std::setprecision(10) << time << "," << step << "," << m_lev << "," << b << ","
            << m_bld_nface[b] << "," << m_bld_area[b] << "," << tmean << "," << swm << "," << shf << "," << lwm << "," << Hm << ","
            << Gm << "," << Qm << "," << btmin[b] << "," << btmax[b] << "," << bres[b] << ","
            << m_sun.zenith * 180.0 / PI << "," << m_sun.azimuth * 180.0 / PI << "," << m_sun.dni << "," << m_sun.diffuse_h << "\n";
    }
}

Here is the call graph for this function:

save_state()

void IBFaceSet::save_state

(

amrex::MultiFab &

state

)

const

Write the face state (skin temperature, sensible flux and slab temperatures) into a field from make_state(), for the checkpoint.

The faces of this rank are written into a transfer field on this rank's grids cut by the state boxes, which a ParallelCopy moves onto the state's own distribution. Each face writes its own slot (see state_ncomp()), so no two faces touch the same component of the same cell and no atomics are needed; cells without a face in a slot keep zero there.

{
    AMREX_ALWAYS_ASSERT(state.nComp() == state_ncomp() && state.boxArray() == m_state_ba);
    const int nl = n_layers();
    const int* pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const int* pslot = d_slot.data();
    const Real* pT = d_T_skin.data(); const Real* pS = d_T_slab.data(); const Real* pH = d_H.data();
    MultiFab mine(m_xfer_ba, m_xfer_dm, state_ncomp(), 0);
    mine.setVal(0.0);
    for (MFIter mfi(mine); mfi.isValid(); ++mfi) {
        const int li = m_xfer_src[mfi.LocalIndex()];
        const int f0 = m_fab_start[li];
        const int f1 = m_fab_start[li + 1];
        const Box bx = mfi.validbox();
        auto const& st = mine.array(mfi);
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f  = f0 + m;
            if (!bx.contains(IntVect(pi[f], pj[f], pk[f]))) { return; }
            const int c0 = pslot[f] * (2 + nl);
            st(pi[f], pj[f], pk[f], c0)     = pT[f];
            st(pi[f], pj[f], pk[f], c0 + 1) = pH[f];
            for (int l = 0; l < nl; ++l) { st(pi[f], pj[f], pk[f], c0 + 2 + l) = pS[f * nl + l]; }
        });
    }
    Gpu::streamSynchronize();
    state.setVal(0.0);
    state.ParallelCopy(mine);
    if (m_params.debug) {
        Print() << "[IBSEB DEBUG] lev=" << m_lev << " face state saved: " << state_ncomp()
                << " components (" << m_nslots << " slots x " << (2 + n_layers()) << ") on " << m_state_ba.size()
                << " boxes, " << m_state_ba.numPts() << " of the level's " << m_domain_cells << " cells\n";
    }
}

Here is the call graph for this function:

scatter_diagnostics()

void IBFaceSet::scatter_diagnostics	(	amrex::MultiFab &	nfaces,
		amrex::MultiFab &	tskin
	)		const

Scatter the faces into two cell-centred diagnostics for the plotfile.

Parameters

[out]	nfaces	Number of wall faces touching each fluid cell (0 away from buildings; up to 3 at an outside corner, 6 in a one-cell slot).
[out]	tskin	Mean skin temperature of those faces, 0 where there are none.

Both fields must be on the level's box array with no ghost cells. Faces of one fab are accumulated with atomic adds because a corner cell receives several faces.

Scatter the face count and the mean skin temperature into cell-centred fields. One kernel per fab over that fab's faces (contiguous, see fab_start()), with atomic adds since a corner cell collects several faces, then a kernel over the box divides the temperature sum by the count.

{
    nfaces.setVal(0.0);
    tskin.setVal(0.0);
    const int* pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const Real* pT = d_T_skin.data();
    for (MFIter mfi(nfaces); mfi.isValid(); ++mfi) {
        const int f0 = m_fab_start[mfi.LocalIndex()];
        const int f1 = m_fab_start[mfi.LocalIndex() + 1];
        auto const& n = nfaces.array(mfi);
        auto const& t = tskin.array(mfi);
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f = f0 + m;
            Gpu::Atomic::AddNoRet(&n(pi[f], pj[f], pk[f]), Real(1.0));
            Gpu::Atomic::AddNoRet(&t(pi[f], pj[f], pk[f]), pT[f]);
        });
        const Box& bx = mfi.validbox();
        ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
            if (n(i, j, k) > 0.0) { t(i, j, k) /= n(i, j, k); }
        });
    }
    Gpu::streamSynchronize();
}

Here is the call graph for this function:

scatter_field()

void IBFaceSet::scatter_field	(	const amrex::Gpu::DeviceVector< amrex::Real > &	v,
		amrex::MultiFab &	out
	)		const

Scatter one per-face array into a cell-centred field as the mean over the faces touching each cell (zero where there are none).

Scatter one per-face array as a per-cell mean, the same accumulation as scatter_diagnostics() with the face count as the divisor.

{
    MultiFab cnt(out.boxArray(), out.DistributionMap(), 1, 0);
    cnt.setVal(0.0);
    out.setVal(0.0);
    const int* pi = d_i.data();  const int* pj = d_j.data();  const int* pk = d_k.data();
    const Real* pv = v.data();
    for (MFIter mfi(out); mfi.isValid(); ++mfi) {
        const int f0 = m_fab_start[mfi.LocalIndex()];
        const int f1 = m_fab_start[mfi.LocalIndex() + 1];
        auto const& n = cnt.array(mfi);
        auto const& o = out.array(mfi);
        ParallelFor(f1 - f0, [=] AMREX_GPU_DEVICE (int m) noexcept {
            const int f = f0 + m;
            Gpu::Atomic::AddNoRet(&n(pi[f], pj[f], pk[f]), Real(1.0));
            Gpu::Atomic::AddNoRet(&o(pi[f], pj[f], pk[f]), pv[f]);
        });
        ParallelFor(mfi.validbox(), [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
            if (n(i, j, k) > 0.0) { o(i, j, k) /= n(i, j, k); }
        });
    }
    Gpu::streamSynchronize();
}

Here is the call graph for this function:

set_init_cost()

void IBFaceSet::set_init_cost ( double  seconds )

inline

Record the wall-clock time of the initialisation (build and view fractions).

{ m_init_cost_s = seconds; }

solve_balance()

void IBFaceSet::solve_balance

(

amrex::Real

dt

)

Prognostic balance of the current step. Per face: the slab's linear response to the skin temperature for this step, the Newton solve of ERF_IBSEBBalance.H with the shortwave, the external flux, the sky and ground longwave, the frozen wall-function coefficient and the latent flux, the slab advanced with the solution, and the fluxes (LW_in, LW_net, H, G) rewritten at the solution. Fills d_resid (|balance| at the end, non-zero only at a bound or the iteration cap) and d_niter. Needs compute_shortwave(), compute_longwave() and compute_sensible() of the same step first.

Parameters

dt	Atmospheric time step [s], the slab's step.

Prognostic balance on every face: one kernel over the faces, each with its own slab response, Newton solve and slab advance. Everything the solve needs was stored by the three flux routines of the step, so no fab data is touched here. The fluxes are rewritten at the solution so that the dump, the reports and the heat deposited into the air all describe the same closed balance; d_resid holds what is left of it.

{
    const int  nl   = n_layers();
    const Real Tint = m_params.T_interior;
    const Real Tmin = m_params.T_skin_min, Tmax = m_params.T_skin_max;
    const Real tol  = m_params.newton_tol_K, max_step = m_params.newton_max_step_K;
    const int  max_iter = m_params.newton_max_iter;
    const Real* pkt = d_kth.data(); const Real* prc = d_rhocp.data(); const Real* pth = d_thick.data();
    const Real* peps = d_emis.data(); const Real* pfb = d_f_bldg.data();
    const Real* psw = d_SW_abs.data(); const Real* pq = d_Q_ext.data(); const Real* pext = d_LW_ext.data();
    const Real* pHc = d_H_coeff.data(); const Real* pTa = d_T_air.data(); const Real* ptha = d_theta_air.data();
    const Real* pLE = d_LE.data();
    Real* pT = d_T_skin.data(); Real* pS = d_T_slab.data(); Real* pG = d_G.data();
    Real* pH = d_H.data(); Real* pin = d_LW_down_in.data(); Real* pnet = d_LW_net.data();
    Real* pres = d_resid.data(); int* pnit = d_niter.data();
    ParallelFor(m_nface, [=] AMREX_GPU_DEVICE (int f) noexcept {
        Real T[ibseb::SLAB_MAX_LAYERS];
        for (int l = 0; l < nl; ++l) { T[l] = pS[f * nl + l]; }
        const Real dz = pth[f] / nl;
        // Conduction the implicit slab step will take for a skin temperature Ts: a Ts - b.
        Real a = 0.0, b = 0.0;
        ibseb::slab_skin_response(T, Tint, pkt[f], prc[f], dz, dt, nl, a, b);
        const Real exner = pTa[f] / ptha[f];
        const Real q_abs = peps[f] * pq[f];
        int n_it = 0; Real res = 0.0;
        const Real Ts = ibseb::solve_skin_balance(pT[f], psw[f], q_abs, pext[f], peps[f], pfb[f],
                                                  pHc[f], exner, ptha[f], pLE[f], a, b,
                                                  Tmin, Tmax, max_step, tol, max_iter, n_it, res);
        // Advance the slab with the solution and store every term at it.
        pG[f] = ibseb::advance_slab_dirichlet(T, Ts, Tint, pkt[f], prc[f], dz, dt, nl);
        for (int l = 0; l < nl; ++l) { pS[f * nl + l] = T[l]; }
        const Real Ts4 = Ts * Ts * Ts * Ts;
        pin[f]  = pext[f] + pfb[f] * ibseb::SIGMA_SB * Ts4;
        pnet[f] = peps[f] * (pin[f] - ibseb::SIGMA_SB * Ts4);
        pH[f]   = pHc[f] * (Ts / exner - ptha[f]);
        pT[f]   = Ts;
        // The residual with the stored terms (the slab's own G, not a Ts - b),
        // so it reports the closure of what is written out.
        pres[f] = std::abs(psw[f] + q_abs + pnet[f] - pH[f] - pLE[f] - pG[f]);
        pnit[f] = n_it;
    });
    Gpu::streamSynchronize();
 
    if (m_params.debug) {
        std::vector<Real> h_r(m_nface), h_T(m_nface);
        std::vector<int>  h_n(m_nface);
        Gpu::copy(Gpu::deviceToHost, d_resid.begin(), d_resid.end(), h_r.begin());
        Gpu::copy(Gpu::deviceToHost, d_T_skin.begin(), d_T_skin.end(), h_T.begin());
        Gpu::copy(Gpu::deviceToHost, d_niter.begin(), d_niter.end(), h_n.begin());
        Gpu::streamSynchronize();
        Real rmax = 0.0; int nmax = 0; Long nclamp = 0;
        for (int n = 0; n < m_nface; ++n) {
            rmax = std::max(rmax, h_r[n]); nmax = std::max(nmax, h_n[n]);
            if (h_T[n] <= Tmin || h_T[n] >= Tmax) { ++nclamp; }
        }
        ParallelDescriptor::ReduceRealMax(rmax);
        ParallelDescriptor::ReduceIntMax(nmax);
        ParallelDescriptor::ReduceLongSum(nclamp);
        Print() << "[IBSEB DEBUG] lev=" << m_lev << " balance: dt=" << dt << " s newton iterations max=" << nmax
                << " residual max=" << rmax << " W/m2 faces at a bound=" << nclamp << "\n";
    }
}

Here is the call graph for this function:

state_boxarray()

const amrex::BoxArray& IBFaceSet::state_boxarray ( ) const

inline

Boxes of the checkpoint field: 4 x 4 column blocks clipped to the k-range that actually owns faces, so the field follows the shell of face-owning cells rather than the built volume. A block over the interior of a wide footprint keeps only the layer above the roof, while a block on a wall still spans its height. Built identically on every rank by build() from a reduced block k-range map; the same buildings give the same boxes on restart, whatever the rank count. Empty when the level has no faces.

{ return m_state_ba; }

state_ncomp()

int IBFaceSet::state_ncomp ( ) const

inline

Number of components of the checkpoint state field.

The field carries n_slots() slots of 2 + n_layers() values per cell: component slot * (2 + n_layers()) is the skin temperature of that slot, the next its sensible flux (read as the previous step's by the convective velocity scale), and the following n_layers() components are its slab layers.

A cell carries one slot per face it owns, numbered in the (dir, side) order in which build() finds them, so the width is the largest face count on any one cell of the level rather than the six a cell could hold in principle. Around ordinary buildings that maximum is two or three (an outside edge or corner); six needs a one-cell slot. build() rebuilds the numbering from the same blanking on restart, so the slots line up whatever the rank count.

{ return m_nslots * (2 + n_layers()); }

Referenced by load_state(), make_state(), and save_state().

Here is the call graph for this function:

Here is the caller graph for this function:

sun()

const SunState& IBFaceSet::sun ( ) const

inline

{ return m_sun; }

Member Data Documentation

d_albedo

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_albedo

Referenced by compute_shortwave(), and dump_faces().

d_area

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_area

Face area [m2].

Referenced by add_heat_flux_to_source(), dump_faces(), and report().

d_bid

amrex::Gpu::DeviceVector<int> IBFaceSet::d_bid

Building id (1-based)

Referenced by dump_faces(), and report().

d_col_top

amrex::Gpu::DeviceVector<int> IBFaceSet::d_col_top

Referenced by compute_shortwave().

d_dir

amrex::Gpu::DeviceVector<int> IBFaceSet::d_dir

Face direction: 0 x, 1 y, 2 z.

Referenced by compute_shortwave(), and dump_faces().

d_emis

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_emis

Optical properties of the face.

Referenced by dump_faces().

d_f_bldg

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_f_bldg

View fractions, sum to 1.

Referenced by dump_faces().

d_f_ground

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_f_ground

Referenced by compute_shortwave(), and dump_faces().

d_f_sky

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_f_sky

Referenced by compute_shortwave(), and dump_faces().

d_G

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_G

Referenced by dump_faces(), and report().

d_H

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_H

Referenced by add_heat_flux_to_source(), dump_faces(), load_state(), report(), and save_state().

d_H_coeff

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_H_coeff

H = H_coeff (theta_skin - theta_air) [W/m2/K].

Referenced by dump_faces().

d_hbld

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_hbld

Height of the face's building column [m].

Referenced by dump_faces().

d_i

amrex::Gpu::DeviceVector<int> IBFaceSet::d_i

Referenced by add_heat_flux_to_source(), dump_faces(), load_state(), save_state(), scatter_diagnostics(), and scatter_field().

d_j

amrex::Gpu::DeviceVector<int> IBFaceSet::d_j

Referenced by add_heat_flux_to_source(), dump_faces(), load_state(), save_state(), scatter_diagnostics(), and scatter_field().

d_k

amrex::Gpu::DeviceVector<int> IBFaceSet::d_k

Fluid cell of the face.

Referenced by add_heat_flux_to_source(), dump_faces(), load_state(), save_state(), scatter_diagnostics(), and scatter_field().

d_kth

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_kth

Referenced by dump_faces().

d_LE

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_LE

Referenced by dump_faces().

d_LW_down_in

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_LW_down_in

Incoming longwave on the face (sky, ground and walls)

Referenced by dump_faces().

d_LW_ext

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_LW_ext

Incoming longwave from sky and ground only [W/m2].

Referenced by dump_faces().

d_LW_net

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_LW_net

Referenced by dump_faces(), and report().

d_mat

amrex::Gpu::DeviceVector<int> IBFaceSet::d_mat

Material id (0 = uniform inputs)

Referenced by dump_faces().

d_niter

amrex::Gpu::DeviceVector<int> IBFaceSet::d_niter

Newton iterations of the last solve.

Referenced by dump_faces().

d_olen

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_olen

Obukhov length used [m] (1e30 when neutral)

Referenced by dump_faces().

d_Q_ext

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_Q_ext

External incident flux [W/m2], e.g. fire radiation.

Referenced by dump_faces(), and report().

d_resid

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_resid

|balance| after the solve [W/m2]

Referenced by dump_faces(), and report().

d_rho

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_rho

Density of the fluid cell [kg/m3].

Referenced by dump_faces().

d_rhocp

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_rhocp

Referenced by dump_faces().

d_shadow

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_shadow

1 when the ray to the sun hits a building, else 0

Referenced by compute_shortwave(), dump_faces(), and report().

d_side

amrex::Gpu::DeviceVector<int> IBFaceSet::d_side

Solid neighbour at -1 (low) or +1 (high) side.

Referenced by compute_shortwave(), and dump_faces().

d_slot

amrex::Gpu::DeviceVector<int> IBFaceSet::d_slot

Slot of the face within its cell's checkpoint record.

Referenced by load_state(), and save_state().

d_SW_abs

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_SW_abs

Referenced by compute_shortwave(), dump_faces(), and report().

d_SW_diffuse_in

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_SW_diffuse_in

Sky and ground-reflected diffuse incident.

Referenced by compute_shortwave(), and dump_faces().

d_SW_direct_in

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_SW_direct_in

Direct beam incident on the face.

Referenced by compute_shortwave(), and dump_faces().

d_T_air

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_T_air

Air temperature of the fluid cell.

Referenced by dump_faces().

d_T_skin

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_T_skin

Skin temperature [K].

Referenced by dump_faces(), load_state(), report(), save_state(), and scatter_diagnostics().

d_T_slab

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_T_slab

Slab layers [K], face-major: [face * n_layers + layer].

Referenced by dump_faces(), load_state(), and save_state().

d_theta_air

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_theta_air

Potential temperature of the fluid cell [K].

Referenced by dump_faces().

d_thick

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_thick

Slab conductivity [W/m/K], heat capacity [J/m3/K], thickness [m].

Referenced by dump_faces().

d_U_tan

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_U_tan

Tangential wind speed at the fluid cell [m/s].

Referenced by dump_faces().

d_ustar

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_ustar

Friction velocity of the wall function [m/s].

Referenced by dump_faces().

d_w_star

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_w_star

Convective velocity scale of the step [m/s].

Referenced by dump_faces(), and report().

d_xf

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_xf

Referenced by compute_shortwave(), and dump_faces().

d_yf

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_yf

Referenced by compute_shortwave(), and dump_faces().

d_z_i

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_z_i

Depth in w* [m]: mixed layer above a roof, building height for a wall.

Referenced by dump_faces().

d_zf

amrex::Gpu::DeviceVector<amrex::Real> IBFaceSet::d_zf

Face centre [m].

Referenced by compute_shortwave(), and dump_faces().

m_area_total

amrex::Real IBFaceSet::m_area_total = 0.0

private

Total face area [m2].

Referenced by print_debug_summary(), and report().

m_bld_area

std::vector<amrex::Real> IBFaceSet::m_bld_area

private

Face area per building [m2].

Referenced by print_debug_summary(), and report().

m_bld_ihi

std::vector<int> IBFaceSet::m_bld_ihi

private

Referenced by print_debug_summary().

m_bld_ilo

std::vector<int> IBFaceSet::m_bld_ilo

private

Referenced by print_debug_summary().

m_bld_jhi

std::vector<int> IBFaceSet::m_bld_jhi

private

Referenced by print_debug_summary().

m_bld_jlo

std::vector<int> IBFaceSet::m_bld_jlo

private

Referenced by print_debug_summary().

m_bld_ncol

std::vector<amrex::Long> IBFaceSet::m_bld_ncol

private

Referenced by print_debug_summary().

m_bld_nface

std::vector<amrex::Long> IBFaceSet::m_bld_nface

private

Faces per building.

Referenced by print_debug_summary(), and report().

m_col_i0

int IBFaceSet::m_col_i0 = 0

Referenced by compute_shortwave().

m_col_j0

int IBFaceSet::m_col_j0 = 0

Referenced by compute_shortwave().

m_col_nx

int IBFaceSet::m_col_nx = 0

Referenced by compute_shortwave().

m_col_ny

int IBFaceSet::m_col_ny = 0

Referenced by compute_shortwave().

m_col_top_max

amrex::Real IBFaceSet::m_col_top_max = 0.0

private

Tallest column [m].

Referenced by compute_shortwave().

m_cost_n

amrex::Long IBFaceSet::m_cost_n = 0

private

Number of per-step updates timed.

Referenced by add_cost(), and report().

m_cost_s

double IBFaceSet::m_cost_s = 0.0

private

Wall-clock seconds spent in the per-step update on this rank.

Referenced by add_cost(), and report().

m_domain_cells

amrex::Long IBFaceSet::m_domain_cells = 0

Cells of the level, for the debug summary of the field.

Referenced by save_state().

m_dx

amrex::Real IBFaceSet::m_dx[3] = {0.0, 0.0, 0.0}

private

Cell sizes [m].

Referenced by compute_shortwave().

m_fab_start

std::vector<int> IBFaceSet::m_fab_start

private

See fab_start()

Referenced by add_heat_flux_to_source(), fab_start(), load_state(), print_debug_summary(), save_state(), scatter_diagnostics(), and scatter_field().

m_init_cost_s

double IBFaceSet::m_init_cost_s = 0.0

private

Wall-clock seconds of build() and compute_view_fractions()

Referenced by report(), and set_init_cost().

m_lev

int IBFaceSet::m_lev = 0

private

AMR level, for reports.

Referenced by compute_shortwave(), load_state(), print_debug_summary(), report(), and save_state().

m_max_path

amrex::Real IBFaceSet::m_max_path = 0.0

private

Ray path cap [m].

Referenced by compute_shortwave().

m_nbld

int IBFaceSet::m_nbld = 0

private

Buildings of the level (global)

Referenced by n_buildings(), print_debug_summary(), and report().

m_nface

int IBFaceSet::m_nface = 0

private

Faces owned by this rank.

Referenced by compute_shortwave(), dump_faces(), load_state(), n_faces(), print_debug_summary(), and report().

m_nface_dir

amrex::Long IBFaceSet::m_nface_dir[3] = {0, 0, 0}

private

Faces per direction.

Referenced by print_debug_summary(), and report().

m_nslots

int IBFaceSet::m_nslots = 0

private

Slots per cell of the checkpoint field (largest face count on a cell)

Referenced by n_slots(), save_state(), and state_ncomp().

m_nx

int IBFaceSet::m_nx = 0

private

Referenced by compute_shortwave().

m_ny

int IBFaceSet::m_ny = 0

private

Columns of the domain.

Referenced by compute_shortwave().

m_params

IBSEBParams IBFaceSet::m_params

private

Copy of the inputs.

Referenced by compute_shortwave(), load_state(), n_layers(), print_debug_summary(), report(), and save_state().

m_per_x

bool IBFaceSet::m_per_x = false

private

Referenced by compute_shortwave().

m_per_y

bool IBFaceSet::m_per_y = false

private

Periodicity in x and y.

Referenced by compute_shortwave().

m_state_ba

amrex::BoxArray IBFaceSet::m_state_ba

Referenced by has_state(), load_state(), make_state(), save_state(), and state_boxarray().

m_sun

SunState IBFaceSet::m_sun

private

Referenced by compute_shortwave(), report(), and sun().

m_x_lo

amrex::Real IBFaceSet::m_x_lo = 0.0

private

Referenced by compute_shortwave().

m_xfer_ba

amrex::BoxArray IBFaceSet::m_xfer_ba

Referenced by load_state(), and save_state().

m_xfer_dm

amrex::DistributionMapping IBFaceSet::m_xfer_dm

Referenced by load_state(), and save_state().

m_xfer_src

std::vector<int> IBFaceSet::m_xfer_src

Referenced by load_state(), and save_state().

m_y_lo

amrex::Real IBFaceSet::m_y_lo = 0.0

private

Domain origin [m].

Referenced by compute_shortwave().

m_z_ground

amrex::Real IBFaceSet::m_z_ground = 0.0

private

Ground height [m] (flat, the domain bottom)

Referenced by compute_shortwave().

---

The documentation for this class was generated from the following files:

• Source/ImmersedBoundarySEB/ERF_IBFaceSet.H
• Source/ImmersedBoundarySEB/ERF_IBFaceSet.cpp