erf_sbm Namespace Reference

Namespaces
	remap_detail

Classes
struct	AdvectionFacePolicyInput
struct	AdvectionFacePolicy
struct	SBMAMRStateView
struct	BulkProjection
class	SBMBulkProjection
struct	ConstraintTerm
struct	LinearConstraint
struct	ConstraintGroup
struct	ConstraintDescriptor
struct	FlattenedConstraintSet
struct	AttachedPropertySupportDescriptor
struct	ConstraintClosureChunk
struct	FixtureInitialStateCanonicality
struct	SpectralPopulationSpec
struct	LiquidProjectionSpec
struct	AttachedPropertyDescriptor
struct	PopulationLayout
struct	SBMLayoutSpec
struct	LayoutValidation
class	SBMLayout
struct	PopulationRemapView
struct	RouteDestination
struct	PacketRoutingPlan
struct	PacketApplicationCoreResult
struct	ReconstructionDeltaView
struct	IntegratedMomentsCoreResult
struct	PacketApplicationResult
struct	ReconstructionDelta
struct	IntegratedMoments
struct	RepresentationIdentity
class	SBMStateManager
class	SBMTransport
struct	SpectralGridSpec
struct	GridValidation
struct	EndpointTransform
class	SpectralGrid

Enumerations
enum class	AdvectionFaceKind : int { Periodic , ImpermeableWall , AdvectiveOutflow , PrescribedInflow }
enum class	AdvectionFaceAction : int {   PeriodicSharedFace , ZeroNormalFlux , InteriorDonor , ExplicitSpectralInflow ,   RejectInwardOutflow , RejectMissingInflowState , RejectIncompleteInflowState }
enum class	PropertyKind { ExtensiveMass , NumberCarried , MassBoundedSubset }
enum class	PropertyRemapPolicy { CarrierBinConservative , Independent }
enum class	SupportRequirement { None , PositiveMass , PositiveNumber }
enum class	PopulationPhase { Liquid , Ice , Aerosol , Generic }
enum class	HostWritePath {   Advection , Diffusion , Source , Positivity ,   Microphysics , Wall , Boundary }
enum class	RemapStatus {   Ok , ZeroWaterResidual , BelowSupportedGrid , Overflow ,   NumericalUnderflow , Invalid }
enum class	ReconstructionStatus { Populated , Empty , Invalid }
enum class	CoordinateKind { Mass , Radius }
enum class	MomentMode { OneMoment , TwoMoment }

Functions
constexpr AdvectionFacePolicy	make_advection_face_policy (const AdvectionFacePolicyInput &input) noexcept
bool	RestrictMappedSpectrum (const SBMLayout &layout, const SBMAMRStateView &fine, const SBMAMRStateView &coarse, const amrex::IntVect &ratio, const int coarse_level, MultiFab &coarse_candidate, std::string &diagnostic)
bool	ProlongCarrierRelativeSpectrum (const SBMLayout &layout, const SBMAMRStateView &coarse, const SBMAMRStateView &fine_target, const amrex::Geometry &coarse_geometry, const amrex::Geometry &fine_geometry, const amrex::IntVect &ratio, const int fine_level, MultiFab &fine_candidate, std::string &diagnostic)
bool	RestrictMappedSpectrum (const SBMLayout &layout, const SBMAMRStateView &fine, const SBMAMRStateView &coarse, const amrex::IntVect &ratio, int coarse_level, amrex::MultiFab &coarse_candidate, std::string &diagnostic)
bool	ProlongCarrierRelativeSpectrum (const SBMLayout &layout, const SBMAMRStateView &coarse, const SBMAMRStateView &fine_target, const amrex::Geometry &coarse_geometry, const amrex::Geometry &fine_geometry, const amrex::IntVect &ratio, int fine_level, amrex::MultiFab &fine_candidate, std::string &diagnostic)
std::string	canonical_real (const amrex::Real value, const bool allow_unbounded_nan=false)
std::vector< ConstraintGroup >	make_constraint_groups (const SBMLayout &layout)
FlattenedConstraintSet	flatten_constraint_groups (const std::vector< ConstraintGroup > &groups, const int state_components)
FlattenedConstraintSet	make_constraint_descriptors (const SBMLayout &layout)
std::vector< AttachedPropertySupportDescriptor >	make_attached_property_support_descriptors (const SBMLayout &layout)
std::vector< ConstraintClosureChunk >	make_constraint_closure_chunks (const SBMLayout &layout, const int max_groups)
EndpointTransform	transform_two_moment (const amrex::Real C, const amrex::Real M, const amrex::Real lower, const amrex::Real upper)
std::pair< amrex::Real, amrex::Real >	inverse_two_moment (const amrex::Real L, const amrex::Real H, const amrex::Real lower, const amrex::Real upper)
bool	property_support_is_admissible (const amrex::Real property, const amrex::Real carrier_number, const amrex::Real lower, const amrex::Real upper, const amrex::Real tolerance) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real	evaluate_constraint_descriptor (const ConstraintDescriptor &constraint, const ConstraintTerm *terms, const amrex::Real *state) noexcept
FixtureInitialStateCanonicality	validate_fixture_initial_state_canonicality (const SBMLayout &layout, const std::vector< amrex::Real > &state) noexcept
const char *	host_write_path_name (const HostWritePath path) noexcept
bool	host_write_allowed (const bool sbm_active, const int component) noexcept
bool	host_write_range_allowed (const bool sbm_active, const int start_component, const int num_components) noexcept
void	require_host_write_allowed (const bool sbm_active, const int component, const HostWritePath path)
void	require_host_write_range_allowed (const bool sbm_active, const int start_component, const int num_components, const HostWritePath path)
PopulationRemapView	population_remap_view (const SBMLayout &layout, const int population_id) noexcept
PacketApplicationResult	apply_packet_routing (const SBMLayout &layout, const PacketRoutingPlan &plan, const std::vector< amrex::Real > &property_packet_amounts, std::vector< amrex::Real > &candidate_state)
ReconstructionStatus	reconstruct_bin (const SBMLayout &layout, const int population_id, const int bin, const std::vector< amrex::Real > &persisted_state, ReconstructionDelta &reconstruction)
bool	integrate_interval (const ReconstructionDelta &reconstruction, const amrex::Real lower, const amrex::Real upper, const bool include_upper, IntegratedMoments &integral)
bool	project_reconstruction_bin (const SBMLayout &layout, const ReconstructionDelta &reconstruction, std::vector< amrex::Real > &candidate_state)
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	valid_population_remap_context (const PopulationRemapView &population) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	valid_population_remap_view (const PopulationRemapView &population) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	same_population_context (const PopulationRemapView &left, const PopulationRemapView &right) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE PacketRoutingPlan	plan_packet_routing (const PopulationRemapView population, const amrex::Real packet_number, const amrex::Real particle_mass) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE PacketApplicationCoreResult	apply_packet_routing_core (const PopulationRemapView population, const PacketRoutingPlan plan, const amrex::Real *property_packet_amounts, const int property_count, amrex::Real *residual_properties, const int residual_capacity, amrex::Real *candidate_state, const int state_count) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ReconstructionStatus	reconstruct_bin_core (const PopulationRemapView population, const int bin, const amrex::Real *persisted_state, const int state_count, amrex::Real *property_per_particle, const int property_capacity, ReconstructionDeltaView &reconstruction) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	integrate_interval_core (const ReconstructionDeltaView reconstruction, const amrex::Real lower, const amrex::Real upper, const bool include_upper, amrex::Real *attached_properties, const int property_capacity, IntegratedMomentsCoreResult &integral) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	project_reconstruction_bin_core (const PopulationRemapView population, const ReconstructionDeltaView reconstruction, amrex::Real *candidate_state, const int state_count) noexcept
constexpr RepresentationIdentity	representation_identity (const MomentMode mode) noexcept
std::string	restart_schema (const SBMLayout &layout)
bool	restart_schema_matches (const SBMLayout &layout, const std::string &persisted)
bool	authoritative_state_admissible (const amrex::MultiFab &spectrum, const SBMLayout &layout, const int level, std::string *diagnostic)
bool	restart_projection_matches (const amrex::MultiFab &spectrum, const amrex::MultiFab &persisted_core, const SBMBulkProjection &projection, const int qc_component, const int qr_component, const amrex::Real tolerance_scale)
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void	copy_host_slow_component (const amrex::Array4< amrex::Real > &current, const amrex::Array4< const amrex::Real > &predictor, const int i, const int j, const int k, const int component, const bool sbm_active, const int qc_component, const int qr_component) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	try_two_moment_to_endpoints (const amrex::Real C, const amrex::Real M, const amrex::Real lower, const amrex::Real upper, EndpointTransform &result) noexcept
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool	try_endpoints_to_two_moment (const amrex::Real L, const amrex::Real H, const amrex::Real lower, const amrex::Real upper, amrex::Real &C, amrex::Real &M) noexcept

Enumeration Type Documentation

AdvectionFaceAction

enum erf_sbm::AdvectionFaceAction : int

strong

Result of applying the M3 spectral state policy at one oriented face.

Enumerator
PeriodicSharedFace
ZeroNormalFlux
InteriorDonor
ExplicitSpectralInflow
RejectInwardOutflow
RejectMissingInflowState
RejectIncompleteInflowState

                               : int {
    PeriodicSharedFace,
    ZeroNormalFlux,
    InteriorDonor,
    ExplicitSpectralInflow,
    RejectInwardOutflow,
    RejectMissingInflowState,
    RejectIncompleteInflowState
};

AdvectionFaceKind

enum erf_sbm::AdvectionFaceKind : int

strong

Physical face categories understood by the bounded M3 advection contract.

Enumerator
Periodic
ImpermeableWall
AdvectiveOutflow
PrescribedInflow

                             : int {
    Periodic,
    ImpermeableWall,
    AdvectiveOutflow,
    PrescribedInflow
};

CoordinateKind

enum erf_sbm::CoordinateKind

strong

Enumerator
Mass
Radius

{ Mass, Radius };

HostWritePath

enum erf_sbm::HostWritePath

strong

Enumerator
Advection
Diffusion
Source
Positivity
Microphysics
Wall
Boundary

                         {
    Advection,
    Diffusion,
    Source,
    Positivity,
    Microphysics,
    Wall,
    Boundary
};

MomentMode

enum erf_sbm::MomentMode

strong

Enumerator
OneMoment
TwoMoment

{ OneMoment, TwoMoment };

PopulationPhase

enum erf_sbm::PopulationPhase

strong

Enumerator
Liquid
Ice
Aerosol
Generic

{ Liquid, Ice, Aerosol, Generic };

PropertyKind

enum erf_sbm::PropertyKind

strong

Semantic class of an attached carrier-bin extensive quantity.

ExtensiveMass is nonnegative additional mass and is not part of the water projection. NumberCarried is a nonnegative count-like quantity. A MassBoundedSubset is nonnegative and may be constrained by its carrier mass when that future physics is implemented.

Enumerator
ExtensiveMass
NumberCarried
MassBoundedSubset

{ ExtensiveMass, NumberCarried, MassBoundedSubset };

PropertyRemapPolicy

enum erf_sbm::PropertyRemapPolicy

strong

P1/P2 remap contract; independent remapping is reserved for future work.

Enumerator
CarrierBinConservative
Independent

{ CarrierBinConservative, Independent };

ReconstructionStatus

enum erf_sbm::ReconstructionStatus

strong

Enumerator
Populated
Empty
Invalid

{ Populated, Empty, Invalid };

RemapStatus

enum erf_sbm::RemapStatus

strong

Enumerator
Ok
ZeroWaterResidual
BelowSupportedGrid
Overflow
NumericalUnderflow
Invalid

                       {
    Ok, ZeroWaterResidual, BelowSupportedGrid, Overflow, NumericalUnderflow, Invalid
};

SupportRequirement

enum erf_sbm::SupportRequirement

strong

Enumerator
None
PositiveMass
PositiveNumber

{ None, PositiveMass, PositiveNumber };

Function Documentation

apply_packet_routing()

PacketApplicationResult erf_sbm::apply_packet_routing	(	const SBMLayout &	layout,
		const PacketRoutingPlan &	plan,
		const std::vector< amrex::Real > &	property_packet_amounts,
		std::vector< amrex::Real > &	candidate_state
	)

{
    PacketApplicationResult result;
    result.status = plan.status;
    if (plan.status != RemapStatus::Ok && plan.status != RemapStatus::ZeroWaterResidual) {
        return result;
    }
    if (&property_packet_amounts == &candidate_state) {
        result.status = RemapStatus::Invalid;
        return result;
    }
    const auto population = population_remap_view(layout, plan.population.population_id);
    if (population.property_count >= 0) {
        result.residual_properties.resize(static_cast<std::size_t>(population.property_count),
                                          amrex::Real(0.0));
    }
    const auto core = apply_packet_routing_core(
        population, plan,
        property_packet_amounts.data(), static_cast<int>(property_packet_amounts.size()),
        result.residual_properties.data(), static_cast<int>(result.residual_properties.size()),
        candidate_state.data(), static_cast<int>(candidate_state.size()));
    result.status = core.status;
    result.residual_number = core.residual_number;
    result.residual_water_mass = core.residual_water_mass;
    result.roundoff_water_mass_correction = core.roundoff_water_mass_correction;
    result.normalized_roundoff = core.normalized_roundoff;
    return result;
}

Here is the call graph for this function:

apply_packet_routing_core()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE PacketApplicationCoreResult erf_sbm::apply_packet_routing_core ( const PopulationRemapView  population, const PacketRoutingPlan  plan, const amrex::Real *  property_packet_amounts, const int  property_count, amrex::Real *  residual_properties, const int  residual_capacity, amrex::Real *  candidate_state, const int  state_count  )

noexcept

Allocation-free packet application. The first pass validates the complete plan, residual outputs, and every candidate increment. The second pass commits only after all failure paths have been ruled out. Input/output and property/state storage supplied by the caller must not overlap.

{
    PacketApplicationCoreResult result;
    result.status = plan.status;
    if (plan.status != RemapStatus::Ok && plan.status != RemapStatus::ZeroWaterResidual) {
        return result;
    }
    if (!valid_population_remap_context(population) ||
        !same_population_context(population, plan.population) ||
        !amrex::Math::isfinite(plan.packet_number) ||
        !amrex::Math::isfinite(plan.packet_particle_mass) ||
        plan.packet_number < amrex::Real(0.0) || plan.packet_particle_mass < amrex::Real(0.0) ||
        property_count != population.property_count ||
        (property_count > 0 && (property_packet_amounts == nullptr || residual_properties == nullptr)) ||
        residual_capacity < property_count || candidate_state == nullptr ||
        state_count != population.state_components || plan.destination_count < 0 ||
        plan.destination_count > 2 || !amrex::Math::isfinite(plan.residual_number_weight) ||
        plan.residual_number_weight < amrex::Real(0.0) ||
        plan.residual_number_weight > amrex::Real(1.0)) {
        result.status = RemapStatus::Invalid;
        return result;
    }
    if (!remap_detail::route_plan_is_consistent(population, plan)) {
        result.status = RemapStatus::Invalid;
        return result;
    }
    for (int property = 0; property < property_count; ++property) {
        if (!remap_detail::finite_nonnegative(property_packet_amounts[property])) {
            result.status = RemapStatus::Invalid;
            return result;
        }
    }
    if (plan.packet_number == amrex::Real(0.0)) {
        for (int property = 0; property < property_count; ++property) {
            if (property_packet_amounts[property] != amrex::Real(0.0)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
        }
    }
    amrex::Real weight_sum = plan.residual_number_weight;
    for (int destination = 0; destination < plan.destination_count; ++destination) {
        const auto& route = plan.destinations[destination];
        if (route.bin < 0 || route.bin >= population.nbins ||
            !amrex::Math::isfinite(route.number_weight) || route.number_weight < amrex::Real(0.0) ||
            route.number_weight > amrex::Real(1.0) ||
            !remap_detail::finite_nonnegative(route.particle_mass)) {
            result.status = RemapStatus::Invalid;
            return result;
        }
        for (int previous = 0; previous < destination; ++previous) {
            if (plan.destinations[previous].bin == route.bin) {
                result.status = RemapStatus::Invalid;
                return result;
            }
        }
        weight_sum += route.number_weight;
    }
    if (plan.packet_number > amrex::Real(0.0)) {
        const amrex::Real weight_tolerance = amrex::Real(16.0) *
            std::numeric_limits<amrex::Real>::epsilon();
        if (!amrex::Math::isfinite(weight_sum) ||
            std::abs(weight_sum - amrex::Real(1.0)) > weight_tolerance ||
            ((plan.status == RemapStatus::ZeroWaterResidual) !=
             (plan.residual_number_weight > amrex::Real(0.0)))) {
            result.status = RemapStatus::Invalid;
            return result;
        }
    } else if (plan.destination_count != 0 ||
               plan.residual_number_weight != amrex::Real(0.0)) {
        result.status = RemapStatus::Invalid;
        return result;
    }
 
    amrex::Real residual_number = amrex::Real(0.0);
    auto product_status = remap_detail::checked_product(
        plan.packet_number, plan.residual_number_weight, residual_number);
    if (product_status != remap_detail::ProductStatus::Ok) {
        result.status = remap_detail::remap_product_status(product_status);
        return result;
    }
    amrex::Real roundoff_water_mass_correction = amrex::Real(0.0);
    // Validate all outputs and candidate increments before changing either.
    for (int property = 0; property < property_count; ++property) {
        amrex::Real residual = amrex::Real(0.0);
        product_status = remap_detail::checked_product(
            property_packet_amounts[property], plan.residual_number_weight, residual);
        if (product_status != remap_detail::ProductStatus::Ok) {
            result.status = remap_detail::remap_product_status(product_status);
            return result;
        }
    }
    for (int destination = 0; destination < plan.destination_count; ++destination) {
        const auto& route = plan.destinations[destination];
        amrex::Real number_increment = amrex::Real(0.0);
        amrex::Real mass_increment = amrex::Real(0.0);
        const amrex::Real deposited_mass = route.particle_mass;
        product_status = remap_detail::checked_product(
            plan.packet_number, route.number_weight, number_increment);
        if (product_status == remap_detail::ProductStatus::Ok) {
            product_status = remap_detail::checked_product(
                number_increment, deposited_mass, mass_increment);
        }
        if (product_status != remap_detail::ProductStatus::Ok) {
            result.status = remap_detail::remap_product_status(product_status);
            return result;
        }
        if (population.moment_mode == MomentMode::TwoMoment && plan.normalized_roundoff) {
            const amrex::Real mass_delta =
                route.particle_mass - plan.packet_particle_mass;
            if (!amrex::Math::isfinite(mass_delta) ||
                mass_delta == amrex::Real(0.0)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
            const amrex::Real mass_delta_magnitude = mass_delta < amrex::Real(0.0) ?
                -mass_delta : mass_delta;
            amrex::Real correction_magnitude = amrex::Real(0.0);
            product_status = remap_detail::checked_product(
                plan.packet_number, mass_delta_magnitude, correction_magnitude);
            if (product_status != remap_detail::ProductStatus::Ok) {
                result.status = remap_detail::remap_product_status(product_status);
                return result;
            }
            const amrex::Real correction = mass_delta < amrex::Real(0.0) ?
                -correction_magnitude : correction_magnitude;
            if (!amrex::Math::isfinite(correction) ||
                !amrex::Math::isfinite(roundoff_water_mass_correction + correction)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
            roundoff_water_mass_correction += correction;
        }
        const int mass_component = population.mass_offset + route.bin;
        const amrex::Real old_mass = candidate_state[mass_component];
        if (!remap_detail::finite_nonnegative(old_mass) ||
            !remap_detail::finite_nonnegative(old_mass + mass_increment)) {
            result.status = RemapStatus::Invalid;
            return result;
        }
        const amrex::Real next_mass = old_mass + mass_increment;
        amrex::Real next_number = amrex::Real(0.0);
        if (population.number_offset >= 0) {
            const int number_component = population.number_offset + route.bin;
            const amrex::Real old_number = candidate_state[number_component];
            if (!remap_detail::finite_nonnegative(old_number)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
            next_number = old_number + number_increment;
            if (!remap_detail::finite_nonnegative(next_number)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
            if (next_number == amrex::Real(0.0)) {
                if (next_mass != amrex::Real(0.0)) {
                    result.status = RemapStatus::Invalid;
                    return result;
                }
            } else {
                amrex::Real mean = amrex::Real(0.0);
                const auto quotient_status = remap_detail::checked_quotient(
                    next_mass, next_number, mean);
                if (quotient_status != remap_detail::QuotientStatus::Ok) {
                    result.status = remap_detail::remap_quotient_status(quotient_status);
                    return result;
                }
                if (!remap_detail::canonical_two_moment_moments(
                        population, route.bin, next_number, next_mass)) {
                    result.status = RemapStatus::Invalid;
                    return result;
                }
            }
        } else {
            const auto quotient_status = remap_detail::checked_quotient(
                next_mass, population.pivots[route.bin], next_number);
            if (quotient_status != remap_detail::QuotientStatus::Ok) {
                result.status = remap_detail::remap_quotient_status(quotient_status);
                return result;
            }
        }
        for (int property = 0; property < property_count; ++property) {
            amrex::Real increment = amrex::Real(0.0);
            const int component = population.property_component_offsets[property] + route.bin;
            const amrex::Real old_value = candidate_state[component];
            product_status = remap_detail::checked_product(
                property_packet_amounts[property], route.number_weight, increment);
            if (product_status != remap_detail::ProductStatus::Ok) {
                result.status = remap_detail::remap_product_status(product_status);
                return result;
            }
            if (!remap_detail::finite_nonnegative(old_value)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
            const amrex::Real next_value = old_value + increment;
            if (!remap_detail::finite_nonnegative(next_value)) {
                result.status = RemapStatus::Invalid;
                return result;
            }
            if (next_number == amrex::Real(0.0)) {
                if (next_value != amrex::Real(0.0)) {
                    result.status = RemapStatus::Invalid;
                    return result;
                }
            } else {
                amrex::Real per_particle = amrex::Real(0.0);
                const auto quotient_status = remap_detail::checked_quotient(
                    next_value, next_number, per_particle);
                if (quotient_status != remap_detail::QuotientStatus::Ok) {
                    result.status = remap_detail::remap_quotient_status(quotient_status);
                    return result;
                }
            }
        }
    }
 
    for (int destination = 0; destination < plan.destination_count; ++destination) {
        const auto& route = plan.destinations[destination];
        const amrex::Real number_increment = plan.packet_number * route.number_weight;
        const amrex::Real deposited_mass = route.particle_mass;
        candidate_state[population.mass_offset + route.bin] += number_increment * deposited_mass;
        if (population.number_offset >= 0) {
            candidate_state[population.number_offset + route.bin] += number_increment;
        }
        for (int property = 0; property < property_count; ++property) {
            candidate_state[population.property_component_offsets[property] + route.bin] +=
                property_packet_amounts[property] * route.number_weight;
        }
    }
    for (int property = 0; property < property_count; ++property) {
        residual_properties[property] =
            property_packet_amounts[property] * plan.residual_number_weight;
    }
    result.residual_number = residual_number;
    result.residual_water_mass = amrex::Real(0.0);
    result.roundoff_water_mass_correction = roundoff_water_mass_correction;
    result.normalized_roundoff = plan.normalized_roundoff;
    return result;
}

Referenced by apply_packet_routing().

Here is the call graph for this function:

Here is the caller graph for this function:

authoritative_state_admissible()

bool erf_sbm::authoritative_state_admissible	(	const amrex::MultiFab &	spectrum,
		const SBMLayout &	layout,
		int	level,
		std::string *	diagnostic = nullptr
	)

Validate every valid-cell value against the layout's runtime constraints.

{
    auto reject = [diagnostic](std::string message) {
        if (diagnostic) *diagnostic = std::move(message);
        return false;
    };
    if (diagnostic) diagnostic->clear();
    if (level < 0) return reject("SBM authoritative restart state has an invalid level");
    if (spectrum.nComp() != layout.ncomp()) {
        return reject("SBM authoritative restart state component count does not match its layout");
    }
 
    const auto groups = make_constraint_groups(layout);
    std::vector<PopulationRemapView> remap_views;
    remap_views.reserve(layout.populations().size());
    for (const auto& population : layout.populations()) {
        remap_views.push_back(population_remap_view(layout, population.population_id));
    }
    std::vector<amrex::Real> state(static_cast<std::size_t>(layout.ncomp()));
    const int precision = std::numeric_limits<amrex::Real>::max_digits10;
 
    for (amrex::MFIter mfi(spectrum); mfi.isValid(); ++mfi) {
        const amrex::FArrayBox& source_fab = spectrum[mfi];
        const amrex::FArrayBox* host_source = &source_fab;
        std::unique_ptr<amrex::FArrayBox> host_fab;
#ifdef AMREX_USE_GPU
        if (source_fab.arena()->isManaged() || source_fab.arena()->isDevice()) {
            host_fab = std::make_unique<amrex::FArrayBox>(
                source_fab.box(), source_fab.nComp(), amrex::The_Pinned_Arena());
            amrex::Gpu::dtoh_memcpy_async(
                host_fab->dataPtr(), source_fab.dataPtr(),
                static_cast<std::size_t>(source_fab.size()) * sizeof(amrex::Real));
            amrex::Gpu::streamSynchronize();
            host_source = host_fab.get();
        }
#endif
        const auto values = host_source->const_array();
        for (amrex::BoxIterator iterator(mfi.validbox()); iterator.ok(); ++iterator) {
            const amrex::IntVect cell = iterator();
            for (int component = 0; component < layout.ncomp(); ++component) {
                const amrex::Real value = values(cell, component);
                if (!std::isfinite(value)) {
                    std::ostringstream message;
                    message << "SBM authoritative restart state is inadmissible"
                            << ": level=" << level << ", cell=(";
                    for (int direction = 0; direction < AMREX_SPACEDIM; ++direction) {
                        if (direction != 0) message << ',';
                        message << cell[direction];
                    }
                    message << "), component=" << component
                            << ", constraint=finite, value="
                            << std::setprecision(precision) << value;
                    const auto group = std::find_if(groups.begin(), groups.end(),
                        [component](const ConstraintGroup& candidate) {
                            return candidate.contains(component);
                        });
                    if (group != groups.end()) {
                        message << ", population=" << group->population_id
                                << " (" << group->semantic_id << ")"
                                << ", bin=" << group->bin;
                    }
                    return reject(message.str());
                }
                state[static_cast<std::size_t>(component)] = value;
            }
 
            for (std::size_t population_index = 0;
                 population_index < remap_views.size(); ++population_index) {
                const auto& view = remap_views[population_index];
                if (view.moment_mode != MomentMode::TwoMoment) continue;
                for (int bin = 0; bin < view.nbins; ++bin) {
                    if (remap_detail::canonical_two_moment_bin_state(
                            view, bin, state.data(), static_cast<int>(state.size()))) continue;
                    std::ostringstream message;
                    message << "SBM authoritative restart state is inadmissible"
                            << ": level=" << level << ", cell=(";
                    for (int direction = 0; direction < AMREX_SPACEDIM; ++direction) {
                        if (direction != 0) message << ',';
                        message << cell[direction];
                    }
                    message << "), population=" << view.population_id
                            << ", bin=" << bin
                            << ", constraint=canonical-two-moment-bin-state";
                    return reject(message.str());
                }
            }
 
            for (const auto& group : groups) {
                amrex::Real margin = amrex::Real(0.0);
                std::string failed_constraint;
                if (group.admissible(state, &margin, &failed_constraint)) continue;
 
                std::ostringstream message;
                message << "SBM authoritative restart state is inadmissible"
                        << ": level=" << level << ", cell=(";
                for (int direction = 0; direction < AMREX_SPACEDIM; ++direction) {
                    if (direction != 0) message << ',';
                    message << cell[direction];
                }
                message << "), population=" << group.population_id
                        << " (" << group.semantic_id << ")"
                        << ", bin=" << group.bin
                        << ", constraint=" << failed_constraint
                        << ", margin=" << std::setprecision(precision) << margin
                        << ", components={";
                for (std::size_t member = 0; member < group.members.size(); ++member) {
                    if (member != 0) message << ',';
                    const int component = group.members[member];
                    message << component << ':'
                            << std::setprecision(precision)
                            << state[static_cast<std::size_t>(component)];
                }
                message << '}';
                return reject(message.str());
            }
 
            for (std::size_t population_index = 0;
                 population_index < remap_views.size(); ++population_index) {
                const auto& view = remap_views[population_index];
                if (view.moment_mode != MomentMode::OneMoment) continue;
                for (int bin = 0; bin < view.nbins; ++bin) {
                    if (remap_detail::canonical_persisted_bin_state(
                            view, bin, state.data(), static_cast<int>(state.size()))) continue;
                    std::ostringstream message;
                    message << "SBM authoritative restart state is inadmissible"
                            << ": level=" << level << ", cell=(";
                    for (int direction = 0; direction < AMREX_SPACEDIM; ++direction) {
                        if (direction != 0) message << ',';
                        message << cell[direction];
                    }
                    message << "), population=" << view.population_id
                            << ", bin=" << bin
                            << ", constraint=canonical-one-moment-bin-state";
                    return reject(message.str());
                }
            }
        }
    }
    return true;
}

Referenced by erf_sbm::SBMTransport::advance_stage(), ProlongCarrierRelativeSpectrum(), ERF::ReadCheckpointFile(), and RestrictMappedSpectrum().

Here is the call graph for this function:

Here is the caller graph for this function:

canonical_real()

std::string erf_sbm::canonical_real ( const amrex::Real  value, const bool  allow_unbounded_nan = false  )

inline

{
    if (allow_unbounded_nan && std::isnan(value)) return "unbounded";
    std::ostringstream out;
    out.imbue(std::locale::classic());
    out.setf(std::ios_base::scientific, std::ios_base::floatfield);
    out << std::setprecision(std::numeric_limits<amrex::Real>::max_digits10) << value;
    return out.str();
}

Referenced by erf_sbm::SpectralGrid::identity(), and erf_sbm::SBMLayout::SBMLayout().

Here is the caller graph for this function:

copy_host_slow_component()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void erf_sbm::copy_host_slow_component ( const amrex::Array4< amrex::Real > &  current, const amrex::Array4< const amrex::Real > &  predictor, const int  i, const int  j, const int  k, const int  component, const bool  sbm_active, const int  qc_component, const int  qr_component  )

noexcept

Perform ERF's predictor-to-current copy for one slow state component.

{
    const bool sbm_projected_liquid =
        sbm_active &&
        (component == qc_component || component == qr_component);
    if (!sbm_projected_liquid) {
        current(i, j, k, component) = predictor(i, j, k, component);
    }
}

Referenced by erf_slow_rhs_post().

Here is the caller graph for this function:

evaluate_constraint_descriptor()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real erf_sbm::evaluate_constraint_descriptor ( const ConstraintDescriptor &  constraint, const ConstraintTerm *  terms, const amrex::Real *  state  )

noexcept

Evaluate a validated flattened constraint in host or device code.

{
    amrex::Real result = amrex::Real(0.0);
    for (int i = 0; i < constraint.term_count; ++i) {
        const auto& term = terms[constraint.term_offset + i];
        result = std::fma(term.coefficient, state[term.component], result);
    }
    return result;
}

flatten_constraint_groups()

FlattenedConstraintSet erf_sbm::flatten_constraint_groups	(	const std::vector< ConstraintGroup > &	groups,
		int	state_components
	)

Flatten arbitrary semantic groups through the production sparse format.

{
    if (state_components < 0) {
        throw std::invalid_argument("SBM production constraints require a nonnegative state size");
    }
    FlattenedConstraintSet flattened;
    for (std::size_t gi = 0; gi < groups.size(); ++gi) {
        const auto& group = groups[gi];
        for (std::size_t ci = 0; ci < group.constraints.size(); ++ci) {
            const auto& constraint = group.constraints[ci];
            if (constraint.terms.empty() ||
                constraint.terms.size() > static_cast<std::size_t>(std::numeric_limits<int>::max()) ||
                flattened.terms.size() > static_cast<std::size_t>(std::numeric_limits<int>::max()) -
                    constraint.terms.size()) {
                throw std::invalid_argument(
                    "SBM production constraints require a nonempty representable sparse term range: " +
                    constraint.semantic_id);
            }
            ConstraintDescriptor descriptor;
            descriptor.group_index = static_cast<int>(gi);
            descriptor.population_id = group.population_id;
            descriptor.bin = group.bin;
            descriptor.constraint_index = static_cast<int>(ci);
            descriptor.term_offset = static_cast<int>(flattened.terms.size());
            descriptor.term_count = static_cast<int>(constraint.terms.size());
            for (const auto& term : constraint.terms) {
                if (term.component < 0 || term.component >= state_components ||
                    !finite(term.coefficient)) {
                    throw std::invalid_argument(
                        "SBM production constraint has an invalid sparse term: " +
                        constraint.semantic_id);
                }
                flattened.terms.push_back(term);
            }
            flattened.constraints.push_back(descriptor);
        }
    }
    return flattened;
}

Referenced by make_constraint_descriptors().

Here is the call graph for this function:

Here is the caller graph for this function:

host_write_allowed()

bool erf_sbm::host_write_allowed ( bool  sbm_active, int  component  )

noexcept

The native host retains qv; SBM exclusively owns the compact liquid lanes.

{
    return host_write_range_allowed(sbm_active, component, 1);
}

Here is the call graph for this function:

host_write_path_name()

const char * erf_sbm::host_write_path_name ( const HostWritePath  path )

noexcept

{
    switch (path) {
    case HostWritePath::Advection: return "advection";
    case HostWritePath::Diffusion: return "diffusion";
    case HostWritePath::Source: return "source";
    case HostWritePath::Positivity: return "positivity";
    case HostWritePath::Microphysics: return "microphysics";
    case HostWritePath::Wall: return "wall";
    case HostWritePath::Boundary: return "boundary";
    }
    return "unknown";
}

Referenced by require_host_write_range_allowed().

Here is the caller graph for this function:

host_write_range_allowed()

bool erf_sbm::host_write_range_allowed ( bool  sbm_active, int  start_component, int  num_components  )

noexcept

Check every component in a native writer's complete, non-empty write range.

{
    if (start_component < 0 || num_components <= 0 ||
        num_components > std::numeric_limits<int>::max() - start_component) {
        return false;
    }
    if (!sbm_active) return true;
    const int end_component = start_component + num_components;
    return !(start_component <= RhoQ2_comp && RhoQ2_comp < end_component) &&
           !(start_component <= RhoQ3_comp && RhoQ3_comp < end_component);
}

Referenced by host_write_allowed(), and require_host_write_range_allowed().

Here is the caller graph for this function:

integrate_interval()

bool erf_sbm::integrate_interval	(	const ReconstructionDelta &	reconstruction,
		const amrex::Real	lower,
		const amrex::Real	upper,
		const bool	include_upper,
		IntegratedMoments &	integral
	)

{
    std::vector<amrex::Real> properties(reconstruction.property_per_particle.size(),
                                        amrex::Real(0.0));
    ReconstructionDeltaView view;
    view.bin = reconstruction.bin;
    view.empty = reconstruction.empty;
    view.normalized_roundoff = reconstruction.normalized_roundoff;
    view.number = reconstruction.number;
    view.particle_mass = reconstruction.particle_mass;
    view.property_per_particle = reconstruction.property_per_particle.data();
    view.property_count = static_cast<int>(reconstruction.property_per_particle.size());
    IntegratedMomentsCoreResult core;
    if (!integrate_interval_core(view, lower, upper, include_upper,
                                 properties.data(), static_cast<int>(properties.size()), core)) {
        return false;
    }
    IntegratedMoments next;
    next.number = core.number;
    next.water_mass = core.water_mass;
    next.attached_properties = std::move(properties);
    integral = std::move(next);
    return true;
}

Here is the call graph for this function:

integrate_interval_core()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::integrate_interval_core ( const ReconstructionDeltaView  reconstruction, const amrex::Real  lower, const amrex::Real  upper, const bool  include_upper, amrex::Real *  attached_properties, const int  property_capacity, IntegratedMomentsCoreResult &  integral  )

noexcept

Integrate one reconstructed delta over [lower, upper), optionally including upper.

{
    if (!amrex::Math::isfinite(lower) || !amrex::Math::isfinite(upper) || !(lower < upper) ||
        !remap_detail::finite_nonnegative(reconstruction.number) ||
        !remap_detail::finite_nonnegative(reconstruction.particle_mass) ||
        reconstruction.property_count < 0 || property_capacity < reconstruction.property_count ||
        (reconstruction.property_count > 0 &&
         (reconstruction.property_per_particle == nullptr || attached_properties == nullptr)) ||
        (!reconstruction.empty && reconstruction.number <= amrex::Real(0.0)) ||
        (reconstruction.empty && (reconstruction.number != amrex::Real(0.0) ||
                                  reconstruction.particle_mass != amrex::Real(0.0)))) {
        return false;
    }
    for (int property = 0; property < reconstruction.property_count; ++property) {
        const amrex::Real value = reconstruction.property_per_particle[property];
        if (!remap_detail::finite_nonnegative(value) ||
            (reconstruction.empty && value != amrex::Real(0.0))) {
            return false;
        }
    }
    const bool contains = !reconstruction.empty && reconstruction.particle_mass >= lower &&
        (reconstruction.particle_mass < upper ||
         (include_upper && reconstruction.particle_mass == upper));
    IntegratedMomentsCoreResult next;
    if (contains) {
        next.number = reconstruction.number;
        if (!remap_detail::finite_product(reconstruction.number,
                                           reconstruction.particle_mass, next.water_mass)) {
            return false;
        }
        for (int property = 0; property < reconstruction.property_count; ++property) {
            amrex::Real amount = amrex::Real(0.0);
            if (!remap_detail::finite_product(reconstruction.number,
                    reconstruction.property_per_particle[property], amount)) {
                return false;
            }
        }
    }
    for (int property = 0; property < reconstruction.property_count; ++property) {
        amrex::Real amount = amrex::Real(0.0);
        if (contains && !remap_detail::finite_product(reconstruction.number,
                reconstruction.property_per_particle[property], amount)) {
            return false;
        }
        attached_properties[property] = amount;
    }
    integral = next;
    return true;
}

Referenced by integrate_interval().

Here is the call graph for this function:

Here is the caller graph for this function:

inverse_two_moment()

std::pair< amrex::Real, amrex::Real > erf_sbm::inverse_two_moment	(	amrex::Real	L,
		amrex::Real	H,
		amrex::Real	lower,
		amrex::Real	upper
	)

Inverse transform with finite-state and nonnegative-endpoint checks.

{
    return SpectralGrid::endpoints_to_two_moment(L, H, lower, upper);
}

make_advection_face_policy()

constexpr AdvectionFacePolicy erf_sbm::make_advection_face_policy ( const AdvectionFacePolicyInput &  input )

constexprnoexcept

Decide the face-level spectral advection semantics without performing ERF's deferred physical-boundary ghost-fill lifecycle. At a physical face whose high-order stencil is not qualified, the admissible donor proposal remains available and high-order transport falls back to it.

{
    using Action = AdvectionFaceAction;
    if (input.kind == AdvectionFaceKind::Periodic) {
        return {Action::PeriodicSharedFace, false, false, true};
    }
    if (input.kind == AdvectionFaceKind::ImpermeableWall ||
        input.outward_carrier == amrex::Real(0.0)) {
        return {Action::ZeroNormalFlux, false, false, false};
    }
 
    if (input.outward_carrier > amrex::Real(0.0)) {
        // Outward flow uses the interior spectrum. The donor path is always
        // defined; high order is used only with a qualified physical stencil.
        return {Action::InteriorDonor, true, false,
                input.high_order_stencil_available};
    }
 
    // Reverse flow contradicts a face declared as advective outflow. External
    // spectral data must not silently change that declared boundary policy.
    if (input.kind == AdvectionFaceKind::AdvectiveOutflow) {
        return {Action::RejectInwardOutflow, false, false, false};
    }
 
    // The remaining negative-carrier boundary category is prescribed inflow.
    // A complete explicit atomic-group spectrum is required; qc/qr cannot
    // supply missing bin or attached-property composition.
    if (!input.explicit_spectral_state_present) {
        return {Action::RejectMissingInflowState, false, false, false};
    }
    if (!input.explicit_spectral_group_complete) {
        return {Action::RejectIncompleteInflowState, false, false, false};
    }
    return {Action::ExplicitSpectralInflow, false, true,
            input.high_order_stencil_available};
}

make_attached_property_support_descriptors()

std::vector< AttachedPropertySupportDescriptor > erf_sbm::make_attached_property_support_descriptors

(

const SBMLayout &

layout

)

{
    std::vector<AttachedPropertySupportDescriptor> result;
    const auto groups = make_constraint_groups(layout);
    for (std::size_t group_index = 0; group_index < groups.size(); ++group_index) {
        const auto& group = groups[group_index];
        const auto population = std::find_if(layout.populations().begin(), layout.populations().end(),
            [&](const PopulationLayout& candidate) {
                return candidate.population_id == group.population_id;
            });
        if (population == layout.populations().end()) {
            throw std::logic_error("SBM attached-property group references an unknown population");
        }
        for (const int property_index : group.attached_property_indices) {
            const auto& property = layout.attached_properties()[static_cast<std::size_t>(property_index)];
            AttachedPropertySupportDescriptor descriptor;
            descriptor.group_index = static_cast<int>(group_index);
            descriptor.property_index = property_index;
            descriptor.property_component = layout.property_offset(property_index) + group.bin;
            descriptor.mass_component = population->mass_offset + group.bin;
            descriptor.number_component = population->number_offset >= 0 ?
                population->number_offset + group.bin : -1;
            descriptor.two_moment = descriptor.number_component >= 0 ? 1 : 0;
            descriptor.carrier_component = descriptor.two_moment ?
                descriptor.number_component : descriptor.mass_component;
            descriptor.pivot = population->grid.pivot(group.bin);
            // The transport envelope is expressed as S/carrier.  For a
            // two-moment carrier this is S/C, while a one-moment property
            // uses C=M/m_pivot and therefore has normalized bounds
            // (support_min/m_pivot, support_max/m_pivot).
            const amrex::Real support_scale = descriptor.two_moment ?
                amrex::Real(1.0) : descriptor.pivot;
            descriptor.hard_min = property.support_min / support_scale;
            descriptor.hard_max = property.support_max / support_scale;
            descriptor.has_hard_min = finite(property.support_min) ? 1 : 0;
            descriptor.has_hard_max = finite(property.support_max) ? 1 : 0;
            result.push_back(descriptor);
        }
    }
    return result;
}

Here is the call graph for this function:

make_constraint_closure_chunks()

std::vector< ConstraintClosureChunk > erf_sbm::make_constraint_closure_chunks	(	const SBMLayout &	layout,
		const int	max_groups
	)

{
    if (max_groups <= 0) {
        throw std::invalid_argument("SBM scratch chunk size must be a positive number of complete groups");
    }
    const auto groups = make_constraint_groups(layout);
    std::vector<ConstraintClosureChunk> chunks;
    for (std::size_t group_index = 0; group_index < groups.size(); ++group_index) {
        if (chunks.empty() ||
            static_cast<int>(chunks.back().group_indices.size()) >= max_groups) {
            chunks.push_back({});
        }
        auto& chunk = chunks.back();
        chunk.group_indices.push_back(static_cast<int>(group_index));
        for (const int component : groups[group_index].members) {
            if (std::find(chunk.components.begin(), chunk.components.end(), component) == chunk.components.end()) {
                chunk.components.push_back(component);
            }
        }
    }
    return chunks;
}

Here is the call graph for this function:

make_constraint_descriptors()

FlattenedConstraintSet erf_sbm::make_constraint_descriptors

(

const SBMLayout &

layout

)

Flatten complete semantic groups into bounded device descriptors.

{
    return flatten_constraint_groups(make_constraint_groups(layout), layout.ncomp());
}

Here is the call graph for this function:

make_constraint_groups()

std::vector< ConstraintGroup > erf_sbm::make_constraint_groups

(

const SBMLayout &

layout

)

Build deterministic groups from the immutable runtime layout.

{
    std::vector<ConstraintGroup> groups;
    for (const auto& population : layout.populations()) {
        const int nbins = population.grid.nbins();
        for (int bin = 0; bin < nbins; ++bin) {
            ConstraintGroup group;
            group.population_id = population.population_id;
            group.bin = bin;
            group.semantic_id = population.semantic_id + ".bin." + std::to_string(bin);
            group.moment_mode = population.moment_mode;
            const int mass = population.mass_offset + bin;
            group.members.push_back(mass);
            group.transport_members.push_back(mass);
            const int number = population.number_offset >= 0 ? population.number_offset + bin : -1;
            if (number >= 0) {
                group.members.push_back(number);
                group.transport_members.push_back(number);
                const amrex::Real lower = population.grid.edges()[static_cast<std::size_t>(bin)];
                const amrex::Real upper = population.grid.edges()[static_cast<std::size_t>(bin + 1)];
                const amrex::Real denominator = upper - lower;
                add_constraint(group, "number_nonnegative", {{number, amrex::Real(1.0)}});
                add_constraint(group, "endpoint_low",
                               {{number, upper / denominator},
                                {mass, -amrex::Real(1.0) / denominator}});
                add_constraint(group, "endpoint_high",
                               {{mass, amrex::Real(1.0) / denominator},
                                {number, -lower / denominator}});
            } else {
                add_constraint(group, "mass_nonnegative", {{mass, amrex::Real(1.0)}});
            }
 
            for (std::size_t property_index = 0;
                 property_index < layout.attached_properties().size(); ++property_index) {
                const auto& property = layout.attached_properties()[property_index];
                if (property.carrier_population != population.population_id) continue;
                const int property_component = layout.property_offset(static_cast<int>(property_index)) + bin;
                group.members.push_back(property_component);
                group.transport_members.push_back(property_component);
                group.attached_property_indices.push_back(static_cast<int>(property_index));
                add_constraint(group, property.semantic_id + ".nonnegative", {{property_component, amrex::Real(1.0)}});
 
                const bool has_lower = finite(property.support_min) && property.support_min > amrex::Real(0.0);
                const bool has_upper = finite(property.support_max) &&
                    property.support_max >= property.support_min;
                if (has_upper) {
                    if (number >= 0) {
                        // Two-moment storage is physical (M,C), so support
                        // bounds are bounds on the carrier number C.
                        add_constraint(group, property.semantic_id + ".support_upper",
                                       {{number, property.support_max},
                                        {property_component, -amrex::Real(1.0)}});
                    } else {
                        const amrex::Real pivot = population.grid.pivot(bin);
                        add_constraint(group, property.semantic_id + ".support_upper",
                                       {{mass, property.support_max / pivot},
                                        {property_component, -1.0}});
                    }
                }
                if (has_lower) {
                    if (number >= 0) {
                        add_constraint(group, property.semantic_id + ".support_lower",
                                       {{number, -property.support_min},
                                        {property_component, amrex::Real(1.0)}});
                    } else {
                        const amrex::Real pivot = population.grid.pivot(bin);
                        add_constraint(group, property.semantic_id + ".support_lower",
                                       {{mass, -property.support_min / pivot},
                                        {property_component, amrex::Real(1.0)}});
                    }
                }
 
                if (property.kind == PropertyKind::MassBoundedSubset) {
                    if (number >= 0) {
                        // The carrier mass is already the authoritative M
                        // component in two-moment storage.
                        add_constraint(group, property.semantic_id + ".carrier_mass_minus_subset",
                                       {{mass, amrex::Real(1.0)}, {property_component, -amrex::Real(1.0)}});
                    } else {
                        add_constraint(group, property.semantic_id + ".carrier_mass_minus_subset",
                                       {{mass, amrex::Real(1.0)}, {property_component, -amrex::Real(1.0)}});
                    }
                }
            }
            groups.push_back(std::move(group));
        }
    }
    return groups;
}

Referenced by authoritative_state_admissible(), ERF::ERF_shared(), make_attached_property_support_descriptors(), make_constraint_closure_chunks(), and make_constraint_descriptors().

Here is the call graph for this function:

Here is the caller graph for this function:

plan_packet_routing()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE PacketRoutingPlan erf_sbm::plan_packet_routing ( const PopulationRemapView  population, const amrex::Real  packet_number, const amrex::Real  particle_mass  )

noexcept

Classify one packet without allocation or state mutation. Two-moment routing uses half-open intervals, includes the global top edge, and admits only four representable steps of outer-boundary roundoff. The plan keeps the original packet mass while its destination records the exact boundary mass. The immutable metadata must have been fully qualified once before this call; the per-packet path performs only constant-time context checks and a logarithmic interval search.

{
    PacketRoutingPlan plan;
    plan.population = population;
    plan.packet_number = packet_number;
    plan.packet_particle_mass = particle_mass;
    if (!valid_population_remap_context(population) ||
        !amrex::Math::isfinite(packet_number) || !amrex::Math::isfinite(particle_mass) ||
        packet_number < amrex::Real(0.0) || particle_mass < amrex::Real(0.0)) {
        return plan;
    }
    if (packet_number == amrex::Real(0.0)) {
        plan.status = RemapStatus::Ok;
        return plan;
    }
 
    if (population.moment_mode == MomentMode::TwoMoment) {
        // A dry positive-number packet is residual material, not a populated
        // zero-water liquid interval, even when the first edge is exactly zero.
        if (particle_mass == amrex::Real(0.0)) {
            plan.status = RemapStatus::ZeroWaterResidual;
            plan.residual_number_weight = amrex::Real(1.0);
            return plan;
        }
        if (particle_mass < population.edges[0]) {
            if (!remap_detail::within_boundary_roundoff(
                    particle_mass, population.edges[0], true)) {
                plan.status = RemapStatus::BelowSupportedGrid;
                return plan;
            }
            plan.status = RemapStatus::Ok;
            plan.destination_count = 1;
            plan.destinations[0] = {0, amrex::Real(1.0), population.edges[0]};
            plan.normalized_roundoff = true;
            return plan;
        }
        if (particle_mass > population.edges[population.nbins]) {
            if (!remap_detail::within_boundary_roundoff(
                    particle_mass, population.edges[population.nbins], false)) {
                plan.status = RemapStatus::Overflow;
                return plan;
            }
            plan.status = RemapStatus::Ok;
            plan.destination_count = 1;
            plan.destinations[0] = {population.nbins - 1, amrex::Real(1.0),
                                    population.edges[population.nbins]};
            plan.normalized_roundoff = true;
            return plan;
        }
        const int bin = particle_mass == population.edges[population.nbins] ?
            population.nbins - 1 :
            amrex::bisect(population.edges, 0, population.nbins, particle_mass);
        plan.status = RemapStatus::Ok;
        plan.destination_count = 1;
        plan.destinations[0] = {bin, amrex::Real(1.0), particle_mass};
        return plan;
    }
 
    if (particle_mass == amrex::Real(0.0)) {
        plan.status = RemapStatus::ZeroWaterResidual;
        plan.residual_number_weight = amrex::Real(1.0);
        return plan;
    }
    if (particle_mass < population.pivots[0]) {
        amrex::Real wet_weight = amrex::Real(0.0);
        const auto quotient_status = remap_detail::checked_quotient(
            particle_mass, population.pivots[0], wet_weight);
        if (quotient_status != remap_detail::QuotientStatus::Ok) {
            plan.status = remap_detail::remap_quotient_status(quotient_status);
            return plan;
        }
        plan.status = RemapStatus::ZeroWaterResidual;
        plan.destination_count = 1;
        plan.destinations[0] = {0, wet_weight, population.pivots[0]};
        plan.residual_number_weight = amrex::Real(1.0) - wet_weight;
        return plan;
    }
    if (particle_mass > population.pivots[population.nbins - 1]) {
        plan.status = RemapStatus::Overflow;
        return plan;
    }
    const int lower_bin = particle_mass == population.pivots[population.nbins - 1] ?
        population.nbins - 1 :
        amrex::bisect(population.pivots, 0, population.nbins - 1, particle_mass);
    if (particle_mass == population.pivots[lower_bin]) {
        plan.status = RemapStatus::Ok;
        plan.destination_count = 1;
        plan.destinations[0] = {lower_bin, amrex::Real(1.0), population.pivots[lower_bin]};
        return plan;
    }
    if (lower_bin + 1 < population.nbins) {
        const amrex::Real lower = population.pivots[lower_bin];
        const amrex::Real upper = population.pivots[lower_bin + 1];
        const amrex::Real width = upper - lower;
        amrex::Real lower_weight = amrex::Real(0.0);
        amrex::Real upper_weight = amrex::Real(0.0);
        const auto lower_status = remap_detail::checked_quotient(
            upper - particle_mass, width, lower_weight);
        const auto upper_status = remap_detail::checked_quotient(
            particle_mass - lower, width, upper_weight);
        if (lower_status != remap_detail::QuotientStatus::Ok ||
            upper_status != remap_detail::QuotientStatus::Ok) {
            const auto quotient_status =
                lower_status != remap_detail::QuotientStatus::Ok ? lower_status : upper_status;
            plan.status = remap_detail::remap_quotient_status(quotient_status);
            return plan;
        }
        plan.status = RemapStatus::Ok;
        plan.destination_count = 2;
        plan.destinations[0] = {lower_bin, lower_weight, lower};
        plan.destinations[1] = {lower_bin + 1, upper_weight, upper};
        return plan;
    }
    return plan;
}

Here is the call graph for this function:

population_remap_view()

PopulationRemapView erf_sbm::population_remap_view ( const SBMLayout &  layout, int  population_id  )

noexcept

Host adapter building a stable view into the immutable SBMLayout.

{
    const auto& populations = layout.populations();
    const auto found = std::find_if(populations.begin(), populations.end(),
        [population_id](const PopulationLayout& population) {
            return population.population_id == population_id;
        });
    if (found == populations.end()) return {};
 
    PopulationRemapView view;
    view.population_id = found->population_id;
    view.moment_mode = found->moment_mode;
    view.coordinate_kind = found->grid.coordinate_kind();
    view.edges = found->grid.edges().data();
    view.pivots = found->grid.pivots().data();
    view.nbins = found->grid.nbins();
    view.mass_offset = found->mass_offset;
    view.number_offset = found->number_offset;
    view.state_components = layout.ncomp();
    view.property_component_offsets = found->property_component_offsets.data();
    view.property_count = static_cast<int>(found->property_component_offsets.size());
    // SBMLayout validates these immutable arrays and constructs the component
    // offsets. This adapter is used on per-packet paths, so do not repeat the
    // O(nbins + properties^2) full qualification here. Keep the full validator
    // for arbitrary raw views; this layout-derived view needs only a constant-
    // time context check before entering the device/core routines.
    if (!valid_population_remap_context(view)) return {};
    return view;
}

Referenced by apply_packet_routing(), authoritative_state_admissible(), project_reconstruction_bin(), reconstruct_bin(), erf_sbm::SBMTransport::SBMTransport(), and validate_fixture_initial_state_canonicality().

Here is the call graph for this function:

Here is the caller graph for this function:

project_reconstruction_bin()

bool erf_sbm::project_reconstruction_bin	(	const SBMLayout &	layout,
		const ReconstructionDelta &	reconstruction,
		std::vector< amrex::Real > &	candidate_state
	)

{
    const auto population = population_remap_view(layout, reconstruction.population_id);
    ReconstructionDeltaView view;
    view.population = population;
    view.bin = reconstruction.bin;
    view.empty = reconstruction.empty;
    view.normalized_roundoff = reconstruction.normalized_roundoff;
    view.number = reconstruction.number;
    view.particle_mass = reconstruction.particle_mass;
    view.property_per_particle = reconstruction.property_per_particle.data();
    view.property_count = static_cast<int>(reconstruction.property_per_particle.size());
    return project_reconstruction_bin_core(population, view, candidate_state.data(),
                                            static_cast<int>(candidate_state.size()));
}

Here is the call graph for this function:

project_reconstruction_bin_core()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::project_reconstruction_bin_core ( const PopulationRemapView  population, const ReconstructionDeltaView  reconstruction, amrex::Real *  candidate_state, const int  state_count  )

noexcept

Project one reconstructed bin into caller-owned state after full validation.

{
    if (!valid_population_remap_context(population) ||
        !same_population_context(population, reconstruction.population) ||
        reconstruction.bin < 0 || reconstruction.bin >= population.nbins ||
        candidate_state == nullptr || state_count != population.state_components ||
        reconstruction.property_count != population.property_count ||
        (reconstruction.property_count > 0 && reconstruction.property_per_particle == nullptr) ||
        !remap_detail::finite_nonnegative(reconstruction.number) ||
        !remap_detail::finite_nonnegative(reconstruction.particle_mass) ||
        (!reconstruction.empty && reconstruction.number <= amrex::Real(0.0)) ||
        (reconstruction.empty && (reconstruction.number != amrex::Real(0.0) ||
                                  reconstruction.particle_mass != amrex::Real(0.0)))) {
        return false;
    }
    if (!reconstruction.empty) {
        if (population.moment_mode == MomentMode::OneMoment) {
            if (reconstruction.particle_mass != population.pivots[reconstruction.bin]) return false;
        }
    }
    amrex::Real mass = amrex::Real(0.0);
    if (!remap_detail::finite_product(reconstruction.number,
                                      reconstruction.particle_mass, mass)) {
        return false;
    }
    amrex::Real projected_number = reconstruction.number;
    if (population.moment_mode == MomentMode::OneMoment) {
        if (remap_detail::checked_quotient(
                mass, population.pivots[reconstruction.bin], projected_number) !=
            remap_detail::QuotientStatus::Ok) return false;
    } else if (!remap_detail::canonical_two_moment_moments(
                   population, reconstruction.bin, projected_number, mass)) {
        return false;
    }
    for (int property = 0; property < population.property_count; ++property) {
        const amrex::Real value = reconstruction.property_per_particle[property];
        amrex::Real amount = amrex::Real(0.0);
        if (!remap_detail::finite_nonnegative(value) ||
            (reconstruction.empty && value != amrex::Real(0.0)) ||
            !remap_detail::finite_product(reconstruction.number, value, amount)) {
            return false;
        }
        if (projected_number == amrex::Real(0.0)) {
            if (amount != amrex::Real(0.0)) return false;
        } else {
            amrex::Real per_particle = amrex::Real(0.0);
            if (remap_detail::checked_quotient(
                    amount, projected_number, per_particle) !=
                remap_detail::QuotientStatus::Ok) return false;
        }
    }
    const int bin = reconstruction.bin;
    candidate_state[population.mass_offset + bin] = mass;
    if (population.number_offset >= 0) {
        candidate_state[population.number_offset + bin] = reconstruction.number;
    }
    for (int property = 0; property < population.property_count; ++property) {
        candidate_state[population.property_component_offsets[property] + bin] =
            reconstruction.number * reconstruction.property_per_particle[property];
    }
    return true;
}

Referenced by project_reconstruction_bin().

Here is the call graph for this function:

Here is the caller graph for this function:

ProlongCarrierRelativeSpectrum() [1/2]

bool erf_sbm::ProlongCarrierRelativeSpectrum	(	const SBMLayout &	layout,
		const SBMAMRStateView &	coarse,
		const SBMAMRStateView &	fine_target,
		const amrex::Geometry &	coarse_geometry,
		const amrex::Geometry &	fine_geometry,
		const amrex::IntVect &	ratio,
		const int	fine_level,
		MultiFab &	fine_candidate,
		std::string &	diagnostic
	)

{
    diagnostic.clear();
    std::string ratio_diagnostic;
    const bool ratio_ok = valid_ratio(ratio, ratio_diagnostic);
    const bool level_ok = fine_level > 0;
 
    std::string coarse_structure_diagnostic;
    std::string fine_structure_diagnostic;
    const bool coarse_structure_ok = validate_view_structure(
        coarse, layout, "coarse", coarse_structure_diagnostic);
    const bool fine_structure_ok = validate_view_structure(
        fine_target, layout, "fine target", fine_structure_diagnostic);
    const bool same_time = coarse.spectrum_time == fine_target.spectrum_time;
    const bool local_structure_ok = ratio_ok && level_ok && coarse_structure_ok &&
                                    fine_structure_ok && same_time;
    const std::string structure_diagnostic =
        !ratio_ok ? ratio_diagnostic
        : !level_ok ? "prolongation fine level must be greater than zero"
        : !coarse_structure_ok ? coarse_structure_diagnostic
        : !fine_structure_ok ? fine_structure_diagnostic
        : !same_time ? "prolongation requires same-time source and target tuples"
                     : std::string{};
    if (!collective_all_true(
            local_structure_ok, structure_diagnostic,
            "SBM prolongation input tuple is invalid on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    const bool candidate_layout_ok =
        fine_candidate.nComp() == layout.ncomp() &&
        erf_auxiliary::SameCellLayout(fine_candidate, *fine_target.spectrum);
    const bool separate_from_coarse = candidate_is_separate(fine_candidate, coarse);
    const bool separate_from_fine_target =
        candidate_is_separate(fine_candidate, fine_target);
    bool local_candidate_ok = candidate_layout_ok && separate_from_coarse &&
                              separate_from_fine_target;
    std::string local_candidate_diagnostic =
        "prolongation candidate must be separate from both input tuples and match the fine SBM layout";
    if (local_candidate_ok) {
        for (int direction = 0; direction < AMREX_SPACEDIM; ++direction) {
            if (!coarse_geometry.isPeriodic(direction) ||
                !fine_geometry.isPeriodic(direction)) {
                local_candidate_ok = false;
                local_candidate_diagnostic =
                    "M4a carrier-relative prolongation requires periodic "
                    "geometry";
                break;
            }
        }
    }
    if (local_candidate_ok) {
        amrex::Box refined_coarse_domain = coarse_geometry.Domain();
        refined_coarse_domain.refine(ratio);
        if (refined_coarse_domain != fine_geometry.Domain()) {
            local_candidate_ok = false;
            local_candidate_diagnostic =
                "prolongation refinement ratio does not map the coarse domain "
                "to the fine domain";
        }
    }
    if (!collective_all_true(
            local_candidate_ok, local_candidate_diagnostic,
            "SBM prolongation candidate or geometry is invalid on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    std::string coarse_carrier_diagnostic;
    std::string fine_carrier_diagnostic;
    const bool coarse_carriers_ok = validate_view_carriers(
        coarse, "coarse", coarse_carrier_diagnostic);
    const bool fine_carriers_ok = validate_view_carriers(
        fine_target, "fine target", fine_carrier_diagnostic);
    std::string coarse_spectrum_diagnostic;
    const bool coarse_spectrum_ok = validate_view_spectrum(
        coarse, layout, fine_level - 1, "coarse", coarse_spectrum_diagnostic);
    const bool local_sources_ok = coarse_carriers_ok && fine_carriers_ok &&
                                  coarse_spectrum_ok;
    const std::string source_diagnostic =
        !coarse_carriers_ok ? coarse_carrier_diagnostic
        : !fine_carriers_ok ? fine_carrier_diagnostic
        : !coarse_spectrum_ok ? coarse_spectrum_diagnostic
                              : std::string{};
    if (!collective_all_true(
            local_sources_ok, source_diagnostic,
            "SBM prolongation source is inadmissible on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    MultiFab coarse_z(coarse.spectrum->boxArray(),
                      coarse.spectrum->DistributionMap(), layout.ncomp(), 0);
    const bool coarse_z_ok = form_carrier_relative_state(
        *coarse.spectrum, *coarse.dry_air_density,
        coarse.density_component, coarse_z, diagnostic);
    if (!collective_all_true(
            coarse_z_ok, diagnostic,
            "forming coarse carrier-relative SBM state failed on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    MultiFab fine_z(fine_target.spectrum->boxArray(),
                    fine_target.spectrum->DistributionMap(), layout.ncomp(), 0);
    amrex::PhysBCFunctNoOp no_physical_boundary;
    amrex::Vector<amrex::BCRec> bcs(static_cast<std::size_t>(layout.ncomp()));
    // PCInterp currently ignores these records under the periodic-only M4a
    // gate. Keep every face explicit; this is not a nonperiodic BC policy.
    for (auto& bc : bcs) {
        for (int direction = 0; direction < AMREX_SPACEDIM; ++direction) {
            bc.setLo(direction, amrex::BCType::int_dir);
            bc.setHi(direction, amrex::BCType::int_dir);
        }
    }
    amrex::InterpFromCoarseLevel(
        fine_z, amrex::IntVect(0), static_cast<Real>(coarse.spectrum_time),
        coarse_z, 0, 0, layout.ncomp(), coarse_geometry, fine_geometry,
        no_physical_boundary, 0, no_physical_boundary, 0, ratio,
        &amrex::pc_interp, bcs, 0);
 
    MultiFab fine_candidate_scratch(fine_target.spectrum->boxArray(),
                                    fine_target.spectrum->DistributionMap(),
                                    layout.ncomp(), 0);
    const bool reconstructed_ok = reconstruct_carrier_relative_state(
        fine_z, *fine_target.dry_air_density, fine_target.density_component,
        *fine_target.mapped_measure, fine_target.measure_component,
        fine_candidate_scratch, diagnostic);
    std::string candidate_admission_diagnostic;
    const bool candidate_admissible = authoritative_state_admissible(
        fine_candidate_scratch, layout, fine_level,
        &candidate_admission_diagnostic);
    const bool local_candidate_admissible = reconstructed_ok &&
                                            candidate_admissible;
    const std::string candidate_diagnostic =
        !reconstructed_ok ? diagnostic
                          : !candidate_admissible
                                ? "prolonged SBM candidate is inadmissible: " +
                                      candidate_admission_diagnostic
                                : std::string{};
    if (!collective_all_true(
            local_candidate_admissible, candidate_diagnostic,
            "prolonged SBM candidate is inadmissible on another MPI rank",
            diagnostic)) {
        return false;
    }
    MultiFab::Copy(fine_candidate, fine_candidate_scratch, 0, 0, layout.ncomp(),
                   0);
    return true;
}

Here is the call graph for this function:

ProlongCarrierRelativeSpectrum() [2/2]

bool erf_sbm::ProlongCarrierRelativeSpectrum	(	const SBMLayout &	layout,
		const SBMAMRStateView &	coarse,
		const SBMAMRStateView &	fine_target,
		const amrex::Geometry &	coarse_geometry,
		const amrex::Geometry &	fine_geometry,
		const amrex::IntVect &	ratio,
		int	fine_level,
		amrex::MultiFab &	fine_candidate,
		std::string &	diagnostic
	)

Form a same-time piecewise-constant carrier-relative fine candidate. On success, overwrite the destination's valid cells with the complete, admitted result. On failure (including a rejection on another MPI rank), leave every destination valid cell and all source views unchanged. This operator does not itself commit authoritative SBM state.

property_support_is_admissible()

bool erf_sbm::property_support_is_admissible ( amrex::Real  property, amrex::Real  carrier_number, amrex::Real  lower, amrex::Real  upper, amrex::Real  tolerance = 0.0  )

noexcept

Host reference for the finite property-support constraints.

{
    if (!finite(property) || !finite(carrier_number) || !finite(lower) || !finite(upper) ||
        carrier_number < amrex::Real(0.0) || lower < amrex::Real(0.0) || upper < lower) return false;
    if (carrier_number == amrex::Real(0.0)) {
        return std::abs(property) <= tolerance;
    }
    return property >= lower * carrier_number - tolerance &&
           property <= upper * carrier_number + tolerance;
}

Here is the call graph for this function:

reconstruct_bin()

ReconstructionStatus erf_sbm::reconstruct_bin	(	const SBMLayout &	layout,
		const int	population_id,
		const int	bin,
		const std::vector< amrex::Real > &	persisted_state,
		ReconstructionDelta &	reconstruction
	)

{
    const auto population = population_remap_view(layout, population_id);
    std::vector<amrex::Real> properties(
        population.property_count >= 0 ? static_cast<std::size_t>(population.property_count) : 0u,
        amrex::Real(0.0));
    ReconstructionDeltaView view;
    const auto status = reconstruct_bin_core(
        population, bin, persisted_state.data(), static_cast<int>(persisted_state.size()),
        properties.data(), static_cast<int>(properties.size()), view);
    if (status == ReconstructionStatus::Invalid) return status;
 
    ReconstructionDelta next;
    next.population_id = population_id;
    next.bin = bin;
    next.empty = view.empty;
    next.normalized_roundoff = view.normalized_roundoff;
    next.number = view.number;
    next.particle_mass = view.particle_mass;
    next.property_per_particle = std::move(properties);
    reconstruction = std::move(next);
    return status;
}

Here is the call graph for this function:

reconstruct_bin_core()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE ReconstructionStatus erf_sbm::reconstruct_bin_core ( const PopulationRemapView  population, const int  bin, const amrex::Real *  persisted_state, const int  state_count, amrex::Real *  property_per_particle, const int  property_capacity, ReconstructionDeltaView &  reconstruction  )

noexcept

{
    if (!valid_population_remap_context(population) || bin < 0 || bin >= population.nbins ||
        persisted_state == nullptr || state_count != population.state_components ||
        property_capacity < population.property_count ||
        (population.property_count > 0 && property_per_particle == nullptr)) {
        return ReconstructionStatus::Invalid;
    }
    if (!remap_detail::canonical_persisted_bin_state(
            population, bin, persisted_state, state_count)) {
        return ReconstructionStatus::Invalid;
    }
    const amrex::Real mass = persisted_state[population.mass_offset + bin];
    if (!remap_detail::finite_nonnegative(mass)) return ReconstructionStatus::Invalid;
    ReconstructionDeltaView next;
    next.population = population;
    next.bin = bin;
    next.property_per_particle = property_per_particle;
    next.property_count = population.property_count;
    if (population.moment_mode == MomentMode::OneMoment) {
        next.particle_mass = population.pivots[bin];
        if (remap_detail::checked_quotient(
                mass, next.particle_mass, next.number) != remap_detail::QuotientStatus::Ok) {
            return ReconstructionStatus::Invalid;
        }
    } else {
        const amrex::Real number = persisted_state[population.number_offset + bin];
        if (!remap_detail::finite_nonnegative(number)) return ReconstructionStatus::Invalid;
        if (!remap_detail::canonical_two_moment_bin_state(
                population, bin, persisted_state, state_count)) {
            return ReconstructionStatus::Invalid;
        }
        const amrex::Real lower = population.edges[bin];
        const amrex::Real upper = population.edges[bin + 1];
        EndpointTransform endpoints;
        if (!try_two_moment_to_endpoints(number, mass, lower, upper, endpoints)) {
            return ReconstructionStatus::Invalid;
        }
        next.number = number;
        // Endpoint tolerance is only a realizability aid for the transform.
        // Persisted state has already passed strict canonical ownership and is
        // never normalized by no-process reconstruction.
        next.normalized_roundoff = false;
        if (number > amrex::Real(0.0)) {
            if (remap_detail::checked_quotient(
                    mass, number, next.particle_mass) != remap_detail::QuotientStatus::Ok) {
                return ReconstructionStatus::Invalid;
            }
        }
    }
    next.empty = next.number == amrex::Real(0.0);
    if (next.empty) next.particle_mass = amrex::Real(0.0);
    for (int property = 0; property < population.property_count; ++property) {
        const amrex::Real amount = persisted_state[
            population.property_component_offsets[property] + bin];
        if (!remap_detail::finite_nonnegative(amount)) return ReconstructionStatus::Invalid;
        if (next.empty) {
            if (amount != amrex::Real(0.0)) return ReconstructionStatus::Invalid;
        } else {
            amrex::Real per_particle = amrex::Real(0.0);
            if (remap_detail::checked_quotient(
                    amount, next.number, per_particle) != remap_detail::QuotientStatus::Ok) {
                return ReconstructionStatus::Invalid;
            }
        }
    }
    for (int property = 0; property < population.property_count; ++property) {
        const amrex::Real amount = persisted_state[
            population.property_component_offsets[property] + bin];
        if (next.empty) {
            property_per_particle[property] = amrex::Real(0.0);
        } else if (remap_detail::checked_quotient(
                       amount, next.number, property_per_particle[property]) !=
                   remap_detail::QuotientStatus::Ok) {
            return ReconstructionStatus::Invalid;
        }
    }
    reconstruction = next;
    return next.empty ? ReconstructionStatus::Empty : ReconstructionStatus::Populated;
}

Referenced by reconstruct_bin().

Here is the call graph for this function:

Here is the caller graph for this function:

representation_identity()

constexpr RepresentationIdentity erf_sbm::representation_identity ( const MomentMode  mode )

constexprnoexcept

{
    switch (mode) {
    case MomentMode::OneMoment:
        return {"fixed-pivot-1m-v1", "fixed-pivot-delta-v1", "fixed-pivot-two-center-v1"};
    case MomentMode::TwoMoment:
        return {"interval-2m-v1", "mean-delta-2m-v1", "actual-mass-interval-v1"};
    }
    return {"invalid-representation", "invalid-reconstruction", "invalid-packet-remap"};
}

Referenced by erf_sbm::SBMLayout::inspection(), and erf_sbm::SBMLayout::SBMLayout().

Here is the caller graph for this function:

require_host_write_allowed()

void erf_sbm::require_host_write_allowed	(	bool	sbm_active,
		int	component,
		HostWritePath	path
	)

Fail closed when a native host writer reaches a component owned by SBM.

{
    require_host_write_range_allowed(sbm_active, component, 1, path);
}

Here is the call graph for this function:

require_host_write_range_allowed()

void erf_sbm::require_host_write_range_allowed	(	bool	sbm_active,
		int	start_component,
		int	num_components,
		HostWritePath	path
	)

Fail closed when any component in a native writer's range is owned by SBM.

{
    if (!host_write_range_allowed(sbm_active, start_component, num_components)) {
        const std::string range = num_components > 0 && start_component >= 0 &&
            num_components <= std::numeric_limits<int>::max() - start_component
            ? " [" + std::to_string(start_component) + "," +
              std::to_string(start_component + num_components) + ")"
            : " (invalid range)";
        amrex::Abort(std::string("SBM owns compact qc/qr; native ") +
                     host_write_path_name(path) + " write denied for component range" + range);
    }
}

Referenced by erf_slow_rhs_post(), if(), and require_host_write_allowed().

Here is the call graph for this function:

Here is the caller graph for this function:

restart_projection_matches()

bool erf_sbm::restart_projection_matches	(	const amrex::MultiFab &	spectrum,
		const amrex::MultiFab &	persisted_core,
		const SBMBulkProjection &	projection,
		int	qc_component,
		int	qr_component,
		amrex::Real	tolerance_scale = amrex::Real(128.0)
	)

Compare the projected spectrum to checkpointed qc/qr without mutating core.

{
    if (spectrum.boxArray() != persisted_core.boxArray() ||
        spectrum.DistributionMap() != persisted_core.DistributionMap() ||
        qc_component < 0 || qr_component < 0 ||
        qc_component >= persisted_core.nComp() || qr_component >= persisted_core.nComp() ||
        !std::isfinite(tolerance_scale) || !(tolerance_scale >= amrex::Real(0.0)) ||
        !spectrum.is_finite(0, spectrum.nComp(), 0) ||
        !persisted_core.is_finite(qc_component, 1, 0) ||
        !persisted_core.is_finite(qr_component, 1, 0)) {
        return false;
    }
 
    amrex::MultiFab expected(persisted_core.boxArray(), persisted_core.DistributionMap(), 2, 0);
    amrex::MultiFab stored(persisted_core.boxArray(), persisted_core.DistributionMap(), 2, 0);
    expected.setVal(amrex::Real(0.0));
    for (amrex::MFIter mfi(expected); mfi.isValid(); ++mfi) {
        projection.apply_to_core(mfi.validbox(), spectrum.const_array(mfi),
                                 expected.array(mfi), 0, 1);
    }
    if (!expected.is_finite(0, 2, 0)) return false;
    amrex::MultiFab::Copy(stored, persisted_core, qc_component, 0, 1, 0);
    amrex::MultiFab::Copy(stored, persisted_core, qr_component, 1, 1, 0);
    if (!stored.is_finite(0, 2, 0)) return false;
 
    const amrex::Real scale = std::max({expected.norm0(0), expected.norm0(1),
                                        stored.norm0(0), stored.norm0(1)});
    amrex::MultiFab::Subtract(expected, stored, 0, 0, 2, 0);
    if (!expected.is_finite(0, 2, 0)) return false;
    const amrex::Real difference = std::max(expected.norm0(0), expected.norm0(1));
    const amrex::Real tolerance = tolerance_scale *
        std::numeric_limits<amrex::Real>::epsilon() *
        std::max(scale, std::numeric_limits<amrex::Real>::min());
    return std::isfinite(scale) && std::isfinite(difference) &&
           std::isfinite(tolerance) && difference <= tolerance;
}

Referenced by ERF::ReadCheckpointFile().

Here is the call graph for this function:

Here is the caller graph for this function:

restart_schema()

std::string erf_sbm::restart_schema

(

const SBMLayout &

layout

)

Exact semantic restart schema; no implicit checkpoint migration is performed.

{
    std::ostringstream schema;
    schema << "ERF-SBM-RESTART-v1\n"
           << "layout=" << layout.schema_identity() << '\n'
           << "constraint-policy=linear-groups-plus-canonical-persisted-v1\n"
           << "projection=liquid-mass-sum-to-qc-qr-v1\n"
           << "transport=sbm-m3-mapped-group-fct-v1\n"
           << "high-order=endpoint-number-wenoz3-v1\n";
    return schema.str();
}

Referenced by restart_schema_matches(), and ERF::WriteCheckpointFile().

Here is the call graph for this function:

Here is the caller graph for this function:

restart_schema_matches()

bool erf_sbm::restart_schema_matches	(	const SBMLayout &	layout,
		const std::string &	persisted
	)

Compare persisted and expected schemas byte-for-byte.

{
    return restart_schema(layout) == persisted;
}

Referenced by ERF::ReadCheckpointFile().

Here is the call graph for this function:

Here is the caller graph for this function:

RestrictMappedSpectrum() [1/2]

bool erf_sbm::RestrictMappedSpectrum	(	const SBMLayout &	layout,
		const SBMAMRStateView &	fine,
		const SBMAMRStateView &	coarse,
		const amrex::IntVect &	ratio,
		const int	coarse_level,
		MultiFab &	coarse_candidate,
		std::string &	diagnostic
	)

{
    diagnostic.clear();
    std::string ratio_diagnostic;
    const bool ratio_ok = valid_ratio(ratio, ratio_diagnostic);
    const bool level_ok = coarse_level >= 0;
    std::string fine_structure_diagnostic;
    std::string coarse_structure_diagnostic;
    const bool fine_structure_ok = validate_view_structure(
        fine, layout, "fine", fine_structure_diagnostic);
    const bool coarse_structure_ok = validate_view_structure(
        coarse, layout, "coarse", coarse_structure_diagnostic);
    const bool same_time = fine.spectrum_time == coarse.spectrum_time;
    // average_down requires complete fine children. Do not query coarsenability
    // until the ratio is positive and the fine spectrum pointer is valid.
    const bool fine_boxes_coarsenable =
        ratio_ok && fine_structure_ok &&
        fine.spectrum->boxArray().coarsenable(ratio);
    const bool local_structure_ok = ratio_ok && level_ok && fine_structure_ok &&
                                    coarse_structure_ok && same_time &&
                                    fine_boxes_coarsenable;
    const std::string structure_diagnostic =
        !ratio_ok ? ratio_diagnostic
        : !level_ok ? "restriction coarse level must be nonnegative"
        : !fine_structure_ok ? fine_structure_diagnostic
        : !coarse_structure_ok ? coarse_structure_diagnostic
        : !same_time ? "restriction requires fine and coarse state tuples at the same semantic time"
        : !fine_boxes_coarsenable
              ? "restriction fine BoxArray is not coarsenable by the refinement ratio; fine boxes must align to complete coarse cells"
                     : std::string{};
    if (!collective_all_true(
            local_structure_ok, structure_diagnostic,
            "SBM restriction input tuple is invalid on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    const bool candidate_layout_ok =
        coarse_candidate.nComp() == layout.ncomp() &&
        erf_auxiliary::SameCellLayout(coarse_candidate, *coarse.spectrum);
    const bool separate_from_fine = candidate_is_separate(coarse_candidate, fine);
    const bool separate_from_coarse = candidate_is_separate(coarse_candidate, coarse);
    const bool local_candidate_ok = candidate_layout_ok && separate_from_fine &&
                                    separate_from_coarse;
    if (!collective_all_true(
            local_candidate_ok,
            "restriction candidate must be separate from both input tuples and match the coarse SBM layout",
            "SBM restriction candidate aliases an input on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    std::string fine_carrier_diagnostic;
    std::string coarse_carrier_diagnostic;
    const bool fine_carriers_ok = validate_view_carriers(
        fine, "fine", fine_carrier_diagnostic);
    const bool coarse_carriers_ok = validate_view_carriers(
        coarse, "coarse", coarse_carrier_diagnostic);
    std::string fine_spectrum_diagnostic;
    std::string coarse_spectrum_diagnostic;
    const bool fine_spectrum_ok = validate_view_spectrum(
        fine, layout, coarse_level + 1, "fine", fine_spectrum_diagnostic);
    const bool coarse_spectrum_ok = validate_view_spectrum(
        coarse, layout, coarse_level, "coarse", coarse_spectrum_diagnostic);
    const bool local_sources_ok = fine_carriers_ok && coarse_carriers_ok &&
                                  fine_spectrum_ok && coarse_spectrum_ok;
    const std::string source_diagnostic =
        !fine_carriers_ok ? fine_carrier_diagnostic
        : !coarse_carriers_ok ? coarse_carrier_diagnostic
        : !fine_spectrum_ok ? fine_spectrum_diagnostic
        : !coarse_spectrum_ok ? coarse_spectrum_diagnostic
                              : std::string{};
    if (!collective_all_true(
            local_sources_ok, source_diagnostic,
            "SBM restriction source is inadmissible on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    MultiFab fine_h(fine.spectrum->boxArray(), fine.spectrum->DistributionMap(),
                    layout.ncomp(), 0);
    MultiFab coarse_h(coarse.spectrum->boxArray(), coarse.spectrum->DistributionMap(),
                      layout.ncomp(), 0);
    const int coverage_component = layout.ncomp();
    const int ncomp = layout.ncomp();
    MultiFab fine_support(fine.spectrum->boxArray(),
                          fine.spectrum->DistributionMap(), layout.ncomp() + 1,
                          0);
    MultiFab coarse_support(coarse.spectrum->boxArray(),
                            coarse.spectrum->DistributionMap(),
                            layout.ncomp() + 1, 0);
    std::string fine_mapped_diagnostic;
    const bool fine_mapped_ok = form_mapped_state(
        *fine.spectrum, *fine.mapped_measure, fine.measure_component, fine_h,
        fine_mapped_diagnostic);
    if (!collective_all_true(
            fine_mapped_ok, fine_mapped_diagnostic,
            "SBM restriction mapped-state formation failed on another MPI rank",
            diagnostic)) {
        return false;
    }
 
    fine_support.setVal(Real(0.0));
    for (amrex::MFIter mfi(fine_h, amrex::TilingIfNotGPU()); mfi.isValid();
         ++mfi) {
        const amrex::Box box = mfi.tilebox();
        const auto mapped = fine_h.const_array(mfi);
        const auto support = fine_support.array(mfi);
        amrex::ParallelFor(
            box, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
                for (int component = 0; component < ncomp; ++component) {
                    support(i, j, k, component) =
                        mapped(i, j, k, component) > Real(0.0) ? Real(1.0)
                                                               : Real(0.0);
                }
                support(i, j, k, coverage_component) = Real(1.0);
            });
    }
 
    // Both averages use the same fine coverage. The extra support components
    // remember whether any child had positive mapped inventory, even if the
    // inventory average itself rounds to zero.
    coarse_h.setVal(Real(0.0));
    coarse_support.setVal(Real(0.0));
    amrex::average_down(fine_h, coarse_h, 0, layout.ncomp(), ratio);
    amrex::average_down(fine_support, coarse_support, 0, layout.ncomp() + 1,
                        ratio);
 
    MultiFab lost_positive_support(coarse_h.boxArray(),
                                   coarse_h.DistributionMap(), 1, 0);
    lost_positive_support.setVal(Real(0.0));
    for (amrex::MFIter mfi(coarse_h, amrex::TilingIfNotGPU()); mfi.isValid();
         ++mfi) {
        const amrex::Box box = mfi.tilebox();
        const auto mapped = coarse_h.const_array(mfi);
        const auto support = coarse_support.const_array(mfi);
        const auto bad = lost_positive_support.array(mfi);
        amrex::ParallelFor(
            box, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
                int cell_bad = 0;
                if (support(i, j, k, coverage_component) > Real(0.0)) {
                    for (int component = 0; component < ncomp; ++component) {
                        if (support(i, j, k, component) > Real(0.0) &&
                            mapped(i, j, k, component) == Real(0.0)) {
                            cell_bad = 1;
                        }
                    }
                }
                bad(i, j, k, 0) = static_cast<Real>(cell_bad);
            });
    }
    // The nonnegative bad mask has max zero when clear; an empty rank's local
    // MultiFab maximum is negative, so both cases satisfy this predicate.
    const bool local_support_preserved =
        lost_positive_support.max(0, 0, true) <= Real(0.0);
    if (!collective_all_true(local_support_preserved,
                             "restriction averaging underflow erased positive "
                             "mapped support to zero",
                             "restriction averaging underflow erased positive "
                             "mapped support on another MPI rank",
                             diagnostic)) {
        return false;
    }
 
    // Build the complete candidate independently. Only cells represented by
    // fine coverage are reconstructed, so uncovered coarse values never take
    // an unnecessary multiply/divide round trip. The public destination is
    // untouched until every rank has admitted this scratch result.
    MultiFab coarse_candidate_scratch(coarse.spectrum->boxArray(),
                                      coarse.spectrum->DistributionMap(),
                                      layout.ncomp(), 0);
    MultiFab::Copy(coarse_candidate_scratch, *coarse.spectrum, 0, 0,
                   layout.ncomp(), 0);
    const bool divided_ok = divide_mapped_state(
        coarse_h, *coarse.mapped_measure, coarse.measure_component,
        coarse_support, coverage_component, coarse_candidate_scratch,
        diagnostic);
    std::string candidate_admission_diagnostic;
    const bool candidate_admissible = authoritative_state_admissible(
        coarse_candidate_scratch, layout, coarse_level,
        &candidate_admission_diagnostic);
    const bool local_candidate_admissible = divided_ok && candidate_admissible;
    const std::string candidate_diagnostic =
        !divided_ok ? diagnostic
                    : !candidate_admissible
                          ? "restricted SBM candidate is inadmissible: " +
                                candidate_admission_diagnostic
                          : std::string{};
    if (!collective_all_true(
            local_candidate_admissible, candidate_diagnostic,
            "restricted SBM candidate is inadmissible on another MPI rank",
            diagnostic)) {
        return false;
    }
    MultiFab::Copy(coarse_candidate, coarse_candidate_scratch, 0, 0,
                   layout.ncomp(), 0);
    return true;
}

Here is the call graph for this function:

RestrictMappedSpectrum() [2/2]

bool erf_sbm::RestrictMappedSpectrum	(	const SBMLayout &	layout,
		const SBMAMRStateView &	fine,
		const SBMAMRStateView &	coarse,
		const amrex::IntVect &	ratio,
		int	coarse_level,
		amrex::MultiFab &	coarse_candidate,
		std::string &	diagnostic
	)

Form an H-conservative coarse candidate from matching-time spectral views. On success, overwrite the destination's valid cells with the complete, admitted result. On failure (including a rejection on another MPI rank), leave every destination valid cell and all source views unchanged. This operator does not itself commit authoritative SBM state.

same_population_context()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::same_population_context ( const PopulationRemapView &  left, const PopulationRemapView &  right  )

noexcept

{
    return left.population_id == right.population_id &&
           left.moment_mode == right.moment_mode &&
           left.coordinate_kind == right.coordinate_kind &&
           left.edges == right.edges && left.pivots == right.pivots &&
           left.nbins == right.nbins && left.mass_offset == right.mass_offset &&
           left.number_offset == right.number_offset &&
           left.state_components == right.state_components &&
           left.property_component_offsets == right.property_component_offsets &&
           left.property_count == right.property_count;
}

Referenced by apply_packet_routing_core(), project_reconstruction_bin_core(), and erf_sbm::remap_detail::route_plan_is_consistent().

Here is the caller graph for this function:

transform_two_moment()

EndpointTransform erf_sbm::transform_two_moment	(	amrex::Real	C,
		amrex::Real	M,
		amrex::Real	lower,
		amrex::Real	upper
	)

Compatibility entry point delegating to the canonical grid-layer transform.

{
    return SpectralGrid::two_moment_to_endpoints(C, M, lower, upper);
}

try_endpoints_to_two_moment()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::try_endpoints_to_two_moment ( const amrex::Real  L, const amrex::Real  H, const amrex::Real  lower, const amrex::Real  upper, amrex::Real &  C, amrex::Real &  M  )

noexcept

Allocation-free canonical inverse transform for host and device callers.

{
    if (!amrex::Math::isfinite(L) || !amrex::Math::isfinite(H) || !amrex::Math::isfinite(lower) ||
        !amrex::Math::isfinite(upper) || L < amrex::Real(0.0) || H < amrex::Real(0.0) ||
        lower < amrex::Real(0.0) || !(lower < upper)) {
        return false;
    }
    const amrex::Real count = L + H;
    const amrex::Real mass = std::fma(lower, L, upper * H);
    if (!amrex::Math::isfinite(count) || !amrex::Math::isfinite(mass) ||
        count < amrex::Real(0.0) || mass < amrex::Real(0.0)) {
        return false;
    }
    C = count;
    M = mass;
    return true;
}

Referenced by erf_sbm::SpectralGrid::endpoints_to_two_moment().

Here is the caller graph for this function:

try_two_moment_to_endpoints()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::try_two_moment_to_endpoints ( const amrex::Real  C, const amrex::Real  M, const amrex::Real  lower, const amrex::Real  upper, EndpointTransform &  result  )

noexcept

Allocation-free canonical forward transform used by host and device code. Returns false instead of throwing when the state is invalid or materially unrealizable. The throwing host API below delegates to this function.

{
    if (!amrex::Math::isfinite(C) || !amrex::Math::isfinite(M) || !amrex::Math::isfinite(lower) ||
        !amrex::Math::isfinite(upper) || lower < amrex::Real(0.0) ||
        !(lower < upper) || C < amrex::Real(0.0)) {
        return false;
    }
    const amrex::Real scale = std::abs(M) + std::abs(lower * C) + std::abs(upper * C);
    const amrex::Real moment_tolerance = amrex::Real(128.0) *
        std::numeric_limits<amrex::Real>::epsilon() * scale;
    const amrex::Real low_numerator = std::fma(upper, C, -M);
    const amrex::Real high_numerator = std::fma(-lower, C, M);
    if (!amrex::Math::isfinite(scale) || !amrex::Math::isfinite(moment_tolerance) ||
        !amrex::Math::isfinite(low_numerator) || !amrex::Math::isfinite(high_numerator) ||
        low_numerator < -moment_tolerance || high_numerator < -moment_tolerance) {
        return false;
    }
    if (C == amrex::Real(0.0)) {
        if (std::abs(M) > moment_tolerance) return false;
        result = {amrex::Real(0.0), amrex::Real(0.0), moment_tolerance,
                  moment_tolerance / (upper - lower), M != amrex::Real(0.0)};
        return amrex::Math::isfinite(result.endpoint_tolerance);
    }
    const amrex::Real denominator = upper - lower;
    amrex::Real L = low_numerator / denominator;
    amrex::Real H = high_numerator / denominator;
    const amrex::Real endpoint_tolerance = moment_tolerance / denominator;
    if (!amrex::Math::isfinite(L) || !amrex::Math::isfinite(H) || !amrex::Math::isfinite(endpoint_tolerance)) {
        return false;
    }
    const bool normalized = L < amrex::Real(0.0) || H < amrex::Real(0.0);
    if (L < amrex::Real(0.0)) L = amrex::Real(0.0);
    if (H < amrex::Real(0.0)) H = amrex::Real(0.0);
    result = {L, H, moment_tolerance, endpoint_tolerance, normalized};
    return true;
}

Referenced by erf_sbm::SBMTransport::advance_stage(), reconstruct_bin_core(), erf_sbm::SpectralGrid::two_moment_realizable(), and erf_sbm::SpectralGrid::two_moment_to_endpoints().

Here is the caller graph for this function:

valid_population_remap_context()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::valid_population_remap_context ( const PopulationRemapView &  population )

noexcept

Constant-time structural/context checks for already qualified immutable metadata. Callers must derive views from a validated SBMLayout or call the full validator below once before using a raw view in the remapping cores.

{
    if (population.population_id < 0 ||
        population.coordinate_kind != CoordinateKind::Mass ||
        (population.moment_mode != MomentMode::OneMoment &&
         population.moment_mode != MomentMode::TwoMoment) ||
        population.nbins <= 0 || population.edges == nullptr ||
        population.state_components < population.nbins || population.mass_offset < 0 ||
        population.property_count < 0 ||
        population.property_count > population.state_components / population.nbins ||
        (population.property_count > 0 && population.property_component_offsets == nullptr)) {
        return false;
    }
    if (population.mass_offset > population.state_components - population.nbins) return false;
    if (population.moment_mode == MomentMode::OneMoment) {
        if (population.number_offset != -1 || population.pivots == nullptr) return false;
    } else if (population.number_offset < 0 ||
               population.number_offset > population.state_components - population.nbins) {
        return false;
    }
    if (population.moment_mode == MomentMode::TwoMoment &&
        population.number_offset < population.mass_offset + population.nbins &&
        population.mass_offset < population.number_offset + population.nbins) {
        return false;
    }
    return true;
}

Referenced by apply_packet_routing_core(), erf_sbm::remap_detail::canonical_one_moment_bin_state_view(), erf_sbm::remap_detail::canonical_two_moment_bin_state_view(), erf_sbm::remap_detail::canonical_two_moment_moments(), plan_packet_routing(), population_remap_view(), project_reconstruction_bin_core(), reconstruct_bin_core(), erf_sbm::remap_detail::route_plan_is_consistent(), and valid_population_remap_view().

Here is the caller graph for this function:

valid_population_remap_view()

AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool erf_sbm::valid_population_remap_view ( const PopulationRemapView &  population )

noexcept

Full O(nbins + properties^2) validator for untrusted/arbitrary views.

{
    if (!valid_population_remap_context(population)) return false;
    for (int edge = 0; edge <= population.nbins; ++edge) {
        const amrex::Real value = population.edges[edge];
        if (!amrex::Math::isfinite(value) || value < amrex::Real(0.0) ||
            (edge > 0 && !(population.edges[edge - 1] < value))) {
            return false;
        }
    }
    if (population.moment_mode == MomentMode::OneMoment) {
        for (int bin = 0; bin < population.nbins; ++bin) {
            const amrex::Real pivot = population.pivots[bin];
            if (!amrex::Math::isfinite(pivot) || pivot <= amrex::Real(0.0) ||
                pivot < population.edges[bin] || pivot > population.edges[bin + 1] ||
                (bin > 0 && !(population.pivots[bin - 1] < pivot))) {
                return false;
            }
        }
    }
    for (int property = 0; property < population.property_count; ++property) {
        const int offset = population.property_component_offsets[property];
        if (offset < 0 || offset > population.state_components - population.nbins ||
            (offset < population.mass_offset + population.nbins &&
             population.mass_offset < offset + population.nbins)) {
            return false;
        }
        if (population.moment_mode == MomentMode::TwoMoment &&
            offset < population.number_offset + population.nbins &&
            population.number_offset < offset + population.nbins) {
            return false;
        }
        for (int previous = 0; previous < property; ++previous) {
            const int previous_offset = population.property_component_offsets[previous];
            if (offset < previous_offset + population.nbins &&
                previous_offset < offset + population.nbins) {
                return false;
            }
        }
    }
    return true;
}

Here is the call graph for this function:

validate_fixture_initial_state_canonicality()

FixtureInitialStateCanonicality erf_sbm::validate_fixture_initial_state_canonicality ( const SBMLayout &  layout, const std::vector< amrex::Real > &  state  )

inlinenoexcept

Apply the persisted-state admission rule to each bin before ERF accepts a user/fixture initial spectrum. Kept as the constructor's shared validation seam so unit tests can exercise the exact candidate-layout check.

{
    FixtureInitialStateCanonicality result;
    if (state.size() != static_cast<std::size_t>(layout.ncomp())) return result;
 
    result.canonical = true;
    const int state_count = static_cast<int>(state.size());
    for (const auto& population : layout.populations()) {
        const auto view = population_remap_view(layout, population.population_id);
        for (int bin = 0; bin < population.grid.nbins(); ++bin) {
            if (!remap_detail::canonical_persisted_bin_state(
                    view, bin, state.data(), state_count)) {
                result.canonical = false;
                result.population_id = population.population_id;
                result.bin = bin;
                result.moment_mode = population.moment_mode;
                return result;
            }
        }
    }
    return result;
}

Referenced by ERF::ERF_shared().

Here is the call graph for this function:

Here is the caller graph for this function: