feat(geometry/surface): System for prescripting surface condition and building computational node for surface condition
The surface needs some defined condition (e.g. Isobaric) This module defines how the surface is defined. Again, much like the deformation system this does not actually implement much of the code to deform the surface to that condition; rather, it implements code to allow the prescription to be described robustly at compile time
This commit is contained in:
83
src/include/serif/geometry/surface/automatic.hpp
Normal file
83
src/include/serif/geometry/surface/automatic.hpp
Normal file
@@ -0,0 +1,83 @@
|
||||
#pragma once
|
||||
|
||||
#include "serif/geometry/surface/detail/automatic_relations.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
enum class PressureSurfaceSelection {
|
||||
supported,
|
||||
unavailable,
|
||||
ambiguous
|
||||
};
|
||||
|
||||
namespace detail {
|
||||
/** @brief Strong boundary enforcement requires an actual solved, scalar H1 equation row. */
|
||||
template <typename Field, typename Rows, typename SolvedTerms, typename PreferredField>
|
||||
concept EligibleSurfaceRow = discretization::blocks::SpatialTerm<Field> && discretization::blocks::HasResidualBlock<Field> &&
|
||||
utils::dependencies::ContainsKey<Field, Rows> && utils::dependencies::ContainsKey<Field, SolvedTerms> &&
|
||||
(std::same_as<PreferredField, void> || std::same_as<PreferredField, Field>);
|
||||
|
||||
/** @brief Search one precedence level, without asking the EOS for a relation when pressure is already an unknown. */
|
||||
template <int Priority, typename EOS, typename Field, typename Fields, typename Rows, typename SolvedTerms, typename PreferredField>
|
||||
[[nodiscard]] consteval auto field_candidates() {
|
||||
if constexpr (!EligibleSurfaceRow<Field, Rows, SolvedTerms, PreferredField>) {
|
||||
return std::tuple<>{};
|
||||
}
|
||||
|
||||
else if constexpr (Priority == 0) {
|
||||
if constexpr (std::same_as<typename Field::Quantity, dimensions::Pressure>) {
|
||||
using Bindings = utils::dependencies::Bindings<utils::dependencies::Binding<dimensions::Pressure, Field>>;
|
||||
return std::tuple<std::type_identity<PressureFieldSurfaceFormulation<EOS, Bindings>>>{};
|
||||
}
|
||||
else {
|
||||
return std::tuple<>{};
|
||||
}
|
||||
}
|
||||
else if constexpr (eos::models::EOSModel<EOS>) {
|
||||
return relation_candidates<Priority == 1, EOS, Field, Fields>(typename EOS::Relations{});
|
||||
}
|
||||
else {
|
||||
return std::tuple<>{};
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Prefer a pressure unknown, then an explicit inverse, then a direct pressure residual. */
|
||||
template <int Priority, typename EOS, typename Rows, typename SolvedTerms, typename PreferredField, typename... Fields>
|
||||
[[nodiscard]] consteval auto select_candidates(utils::dependencies::KeyList<Fields...> fields) {
|
||||
using FieldsList = utils::dependencies::KeyList<Fields...>;
|
||||
using Candidates = decltype(std::tuple_cat(field_candidates<Priority, EOS, Fields, FieldsList, Rows, SolvedTerms, PreferredField>()...));
|
||||
|
||||
if constexpr (std::tuple_size_v<Candidates> != 0 || Priority == 2) {
|
||||
return Candidates{};
|
||||
}
|
||||
else {
|
||||
return select_candidates<Priority + 1, EOS, Rows, SolvedTerms, PreferredField>(fields);
|
||||
}
|
||||
}
|
||||
|
||||
template <typename Candidates>
|
||||
struct SelectedPressureSurface {
|
||||
static constexpr PressureSurfaceSelection status = std::tuple_size_v<Candidates> == 0 ? PressureSurfaceSelection::unavailable : PressureSurfaceSelection::ambiguous;
|
||||
};
|
||||
|
||||
template <typename Candidate>
|
||||
struct SelectedPressureSurface<std::tuple<Candidate>> {
|
||||
static constexpr PressureSurfaceSelection status = PressureSurfaceSelection::supported;
|
||||
using Formulation = typename Candidate::type;
|
||||
};
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief Compile an isobaric condition from available physical fields, actual equation rows, and the EOS catalog.
|
||||
* @note This searches declared algebraic relations, never performs inversion, and rejects equal-priority choices.
|
||||
* Global thermodynamic inputs may opt in with thermodynamic_state = true; units alone do not establish that role.
|
||||
* An optional preferred field narrows the row search without choosing an EOS relation for the model user.
|
||||
*/
|
||||
template <typename EOS, utils::dependencies::IsKeyList AvailableTerms, utils::dependencies::IsKeyList ActualRows,
|
||||
utils::dependencies::IsKeyList SolvedTerms = AvailableTerms, typename PreferredField = void>
|
||||
struct AutomaticPressureSurface : detail::SelectedPressureSurface<decltype(detail::select_candidates<0, std::remove_cvref_t<EOS>, ActualRows, SolvedTerms, PreferredField>(detail::thermodynamic_fields(AvailableTerms{})))> {
|
||||
};
|
||||
|
||||
/** @brief Whether the available equations and numerical backend select one implementable isobaric residual. */
|
||||
template <typename EOS, typename AvailableTerms, typename ActualRows, typename SolvedTerms = AvailableTerms, typename PreferredField = void>
|
||||
concept AutomaticPressureSurfaceCompilable = requires { typename AutomaticPressureSurface<EOS, AvailableTerms, ActualRows, SolvedTerms, PreferredField>::Formulation; };
|
||||
}
|
||||
108
src/include/serif/geometry/surface/compilation.hpp
Normal file
108
src/include/serif/geometry/surface/compilation.hpp
Normal file
@@ -0,0 +1,108 @@
|
||||
#pragma once
|
||||
|
||||
#include "serif/geometry/surface/detail/relations.hpp"
|
||||
#include "serif/geometry/surface/isobaric.hpp"
|
||||
|
||||
namespace serif::geometry::surface
|
||||
{
|
||||
/**
|
||||
* @brief Whether the available state selects exactly one pressure-to-carrier EOS relation.
|
||||
* @tparam Carrier Optional physical quantity whose equation enforces the boundary condition.
|
||||
* @note An ambiguous catalog never selects its first entry. Specifying a carrier narrows the same search.
|
||||
*/
|
||||
template <typename EOS, typename Bindings, typename Carrier = void>
|
||||
concept PressureSurfaceCompilable = eos::models::EOSModel<EOS> && SurfaceStateBindings<Bindings> &&
|
||||
(detail::MatchingPressureRelations<std::remove_cvref_t<EOS>, std::remove_cvref_t<Bindings>, Carrier, typename std::remove_cvref_t<EOS>::Relations>::count == 1);
|
||||
|
||||
template <typename EOS, typename Bindings, typename Carrier = void>
|
||||
requires PressureSurfaceCompilable<EOS, Bindings, Carrier>
|
||||
using pressure_relation_t = std::tuple_element_t<0, typename detail::MatchingPressureRelations<std::remove_cvref_t<EOS>, std::remove_cvref_t<Bindings>, Carrier, typename std::remove_cvref_t<EOS>::Relations>::Matches>;
|
||||
|
||||
namespace detail
|
||||
{
|
||||
/** @brief Common physical metadata; numerical values and preparation belong to the dependency context and cache. */
|
||||
template <typename EOS, typename Bindings, typename RowQuantity, typename RelationT, bool ResolveCarrier>
|
||||
struct PressureSurfaceDescription
|
||||
{
|
||||
using EquationOfState = std::remove_cvref_t<EOS>;
|
||||
using StateBindings = std::remove_cvref_t<Bindings>;
|
||||
using Relation = RelationT;
|
||||
using PhysicalQuantity = dimensions::Pressure;
|
||||
using ResidualQuantity = std::conditional_t<ResolveCarrier, RowQuantity, dimensions::Pressure>;
|
||||
using CarrierQuantity = RowQuantity;
|
||||
using CarrierTerm = utils::dependencies::bound_type_t<StateBindings, RowQuantity>;
|
||||
using Test = CarrierTerm;
|
||||
using Support = discretization::domain::StellarSurfaceBoundary;
|
||||
using EOSStateQuantities = std::conditional_t<ResolveCarrier, typename RelationQuantities<Relation>::NonPressureInputs, typename RelationQuantities<Relation>::AllInputs>;
|
||||
using StateQuantities = utils::dependencies::CombinedKeyList<utils::dependencies::KeyList<RowQuantity>, EOSStateQuantities>;
|
||||
using StateTerms = BoundTerms<StateBindings, StateQuantities>;
|
||||
using EOSStateTerms = BoundTerms<StateBindings, EOSStateQuantities>;
|
||||
using Operands = typename SurfaceOperands<StateTerms>::Type;
|
||||
|
||||
static constexpr bool resolves_carrier = ResolveCarrier;
|
||||
static constexpr bool uses_eos = true;
|
||||
};
|
||||
}
|
||||
|
||||
/** @brief The boundary residual is carrier minus the EOS value at prescribed pressure; no carrier quantity is privileged. */
|
||||
template <typename EOS, typename Bindings, typename Carrier = void>
|
||||
requires PressureSurfaceCompilable<EOS, Bindings, Carrier>
|
||||
struct PressureSurfaceFormulation final : detail::PressureSurfaceDescription<EOS, Bindings, typename pressure_relation_t<EOS, Bindings, Carrier>::OutputQuantity, pressure_relation_t<EOS, Bindings, Carrier>, true>
|
||||
{
|
||||
};
|
||||
|
||||
/** @brief A selected explicit inverse relation; choosing it requires no numerical EOS inversion. */
|
||||
template <eos::models::EOSModel EOS, SurfaceStateBindings Bindings, typename Relation>
|
||||
requires eos::relations::CatalogContainingRelation<typename EOS::Relations, Relation> && (detail::matches_pressure_relation<EOS, Bindings, void, Relation>())
|
||||
struct InversePressureSurfaceFormulation final : detail::PressureSurfaceDescription<EOS, Bindings, typename Relation::OutputQuantity, Relation, true>
|
||||
{
|
||||
};
|
||||
|
||||
/** @brief A direct pressure residual requires an explicitly chosen pressure relation and an input quantity whose equation row it replaces. */
|
||||
template <typename EOS, typename Bindings, typename RowQuantity, typename Relation>
|
||||
concept DirectPressureSurfaceCompilable = eos::models::EOSModel<EOS> && SurfaceStateBindings<Bindings> && eos::relations::ValidThermodynamicRelation<Relation> &&
|
||||
eos::relations::CatalogContainingRelation<typename std::remove_cvref_t<EOS>::Relations, Relation> &&
|
||||
std::same_as<typename Relation::OutputQuantity, dimensions::Pressure> && eos::relations::RelationContainingInput<Relation, RowQuantity> &&
|
||||
SurfaceResidualRow<Bindings, RowQuantity> &&
|
||||
utils::dependencies::SatisfiesKeyList<Bindings, typename detail::RelationQuantities<Relation>::AllInputs> &&
|
||||
detail::supports_partials<EOS, Relation>(typename detail::RelationQuantities<Relation>::AllInputs{});
|
||||
|
||||
/**
|
||||
* @brief Enforce P_EOS(state) - P_surface = 0 without inventing an inverse EOS relation.
|
||||
* @note The explicit row and relation make the choice of residual, units, and Newton scaling visible.
|
||||
*/
|
||||
template <typename EOS, typename Bindings, typename RowQuantity, typename Relation>
|
||||
requires DirectPressureSurfaceCompilable<EOS, Bindings, RowQuantity, Relation>
|
||||
struct DirectPressureSurfaceFormulation final : detail::PressureSurfaceDescription<EOS, Bindings, RowQuantity, Relation, false>
|
||||
{
|
||||
};
|
||||
|
||||
/** @brief An existing pressure unknown enforces P - P_surface directly, without evaluating the EOS. */
|
||||
template <typename EOS, SurfaceStateBindings Bindings>
|
||||
requires SurfaceResidualRow<Bindings, dimensions::Pressure>
|
||||
struct PressureFieldSurfaceFormulation final : detail::PressureSurfaceDescription<EOS, Bindings, dimensions::Pressure, eos::relations::Relation<dimensions::Pressure, dimensions::Pressure>, true>
|
||||
{
|
||||
static constexpr bool uses_eos = false;
|
||||
};
|
||||
|
||||
namespace detail
|
||||
{
|
||||
template <typename Candidate>
|
||||
inline constexpr bool is_pressure_surface_formulation = false;
|
||||
|
||||
template <typename EOS, typename Bindings, typename Carrier>
|
||||
inline constexpr bool is_pressure_surface_formulation<PressureSurfaceFormulation<EOS, Bindings, Carrier>> = true;
|
||||
|
||||
template <typename EOS, typename Bindings, typename Relation>
|
||||
inline constexpr bool is_pressure_surface_formulation<InversePressureSurfaceFormulation<EOS, Bindings, Relation>> = true;
|
||||
|
||||
template <typename EOS, typename Bindings, typename RowQuantity, typename Relation>
|
||||
inline constexpr bool is_pressure_surface_formulation<DirectPressureSurfaceFormulation<EOS, Bindings, RowQuantity, Relation>> = true;
|
||||
|
||||
template <typename EOS, typename Bindings>
|
||||
inline constexpr bool is_pressure_surface_formulation<PressureFieldSurfaceFormulation<EOS, Bindings>> = true;
|
||||
}
|
||||
|
||||
template <typename Candidate>
|
||||
concept PressureSurfaceFormulationType = detail::is_pressure_surface_formulation<Candidate>;
|
||||
}
|
||||
108
src/include/serif/geometry/surface/computations.hpp
Normal file
108
src/include/serif/geometry/surface/computations.hpp
Normal file
@@ -0,0 +1,108 @@
|
||||
#pragma once
|
||||
|
||||
#include <tuple>
|
||||
#include <type_traits>
|
||||
#include <utility>
|
||||
|
||||
#include "serif/geometry/surface/detail/batches.hpp"
|
||||
#include "serif/utils/dependencies/computation_graph.hpp"
|
||||
|
||||
namespace serif::geometry::surface
|
||||
{
|
||||
namespace detail
|
||||
{
|
||||
/** @brief Declare the exact EOS inputs from the formulation, so derivative metadata and cache dependencies cannot diverge. */
|
||||
template <PressureSurfaceFormulationType Formulation, template <typename> typename SampleKey, typename... Terms>
|
||||
[[nodiscard]] auto make_surface_eos_computation(utils::dependencies::KeyList<Terms...>)
|
||||
{
|
||||
using namespace utils::dependencies;
|
||||
using EOS = typename Formulation::EquationOfState;
|
||||
|
||||
if constexpr (!Formulation::uses_eos)
|
||||
{
|
||||
// Reuse the carrier residual stage with one constant target; pressure unknowns need no EOS dependency.
|
||||
return make_computation<keys::EOSResponse, InputList<keys::PressureCondition>>(
|
||||
[](const IsobaricSurface& condition)
|
||||
{
|
||||
return std::vector<SurfaceEOSResponse<Formulation>>{{.value = condition.target_pressure()}};
|
||||
});
|
||||
}
|
||||
else if constexpr (Formulation::resolves_carrier)
|
||||
{
|
||||
return make_computation<keys::EOSResponse, InputList<keys::EquationOfState, keys::PressureCondition, SampleKey<Terms>...>>(
|
||||
[](const EOS& equation_of_state, const IsobaricSurface& condition, const BoundaryValues<Terms>&... values)
|
||||
{
|
||||
return evaluate_response_samples<Formulation>(equation_of_state, condition.target_pressure(), values...);
|
||||
});
|
||||
}
|
||||
else
|
||||
{
|
||||
return make_computation<keys::EOSResponse, InputList<keys::EquationOfState, SampleKey<Terms>...>>(
|
||||
[](const EOS& equation_of_state, const BoundaryValues<Terms>&... values)
|
||||
{
|
||||
return evaluate_response_samples<Formulation>(equation_of_state, dimensions::PressureValue{0.0}, values...);
|
||||
});
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Keep the inexpensive boundary residual separate from the reusable EOS response and its partials. */
|
||||
template <PressureSurfaceFormulationType Formulation, template <typename> typename SampleKey>
|
||||
[[nodiscard]] auto make_surface_linearization_computation()
|
||||
{
|
||||
using namespace utils::dependencies;
|
||||
using Responses = std::vector<SurfaceEOSResponse<Formulation>>;
|
||||
|
||||
if constexpr (Formulation::resolves_carrier)
|
||||
{
|
||||
using Carrier = typename Formulation::CarrierTerm;
|
||||
return make_computation<keys::Linearization, InputList<keys::EOSResponse, SampleKey<Carrier>>>(
|
||||
[](const Responses& responses, const BoundaryValues<Carrier>& carrier)
|
||||
{
|
||||
return linearize_carrier_samples<Formulation>(responses, carrier);
|
||||
});
|
||||
}
|
||||
else
|
||||
{
|
||||
return make_computation<keys::Linearization, InputList<keys::EOSResponse, keys::PressureCondition>>(
|
||||
[](const Responses& responses, const IsobaricSurface& condition)
|
||||
{
|
||||
return linearize_pressure_samples<Formulation>(condition, responses);
|
||||
});
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Expose the existing EOS and residual stages for composition with backend sampling or other physics. */
|
||||
template <PressureSurfaceFormulationType Formulation, template <typename> typename SampleKey = std::type_identity_t>
|
||||
[[nodiscard]] auto make_pressure_surface_computations()
|
||||
{
|
||||
return std::tuple{
|
||||
detail::make_surface_eos_computation<Formulation, SampleKey>(typename Formulation::EOSStateTerms{}),
|
||||
detail::make_surface_linearization_computation<Formulation, SampleKey>()};
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief Describe the isobaric calculation as an ordinary dependency graph, without evaluating the EOS.
|
||||
* @tparam Formulation A carrier-resolved or explicitly selected direct-pressure surface formulation.
|
||||
* @return Two computations: the EOS response and the boundary residual with analytical derivatives.
|
||||
* @note Supply pressure_condition and boundary_values, and bind keys::EquationOfState with a tracking policy when the formulation uses the EOS.
|
||||
* A pressure-field formulation needs no EOS source. Use make_computation_cache<SurfaceFailure> and its normal try_prepare/get interface;
|
||||
* the cache owns the complete preparation lifecycle.
|
||||
*/
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] auto make_pressure_surface_graph()
|
||||
{
|
||||
return std::apply([](auto... computations)
|
||||
{
|
||||
return utils::dependencies::make_computation_graph(std::move(computations)...);
|
||||
}, make_pressure_surface_computations<Formulation>());
|
||||
}
|
||||
|
||||
/** @brief Infer a unique carrier-resolved formulation directly from an EOS and the available physical state. */
|
||||
template <typename EOS, typename Bindings, typename Carrier = void>
|
||||
requires PressureSurfaceCompilable<EOS, Bindings, Carrier>
|
||||
[[nodiscard]] auto make_pressure_surface_graph()
|
||||
{
|
||||
return make_pressure_surface_graph<PressureSurfaceFormulation<EOS, Bindings, Carrier>>();
|
||||
}
|
||||
}
|
||||
@@ -0,0 +1,87 @@
|
||||
#pragma once
|
||||
|
||||
#include "serif/geometry/surface/compilation.hpp"
|
||||
#include "serif/discretization/blocks/traits.hpp"
|
||||
|
||||
namespace serif::geometry::surface::detail {
|
||||
/** @brief Globals require an explicit physical role: matching units alone would mistake C or a multiplier for enthalpy. */
|
||||
template <typename Term>
|
||||
concept ThermodynamicBoundaryField = discretization::blocks::ScalarTraceTerm<Term> &&
|
||||
discretization::blocks::HasQuantity<Term> && dimensions::ThermodynamicQuantityType<typename Term::Quantity> &&
|
||||
(discretization::blocks::SpatialTerm<Term> || requires { requires Term::thermodynamic_state; });
|
||||
|
||||
/** @brief Collect physical fields, retaining distinct terms with equal quantities so ambiguity remains visible. */
|
||||
template <typename... Terms>
|
||||
[[nodiscard]] consteval auto thermodynamic_fields(utils::dependencies::KeyList<Terms...>) {
|
||||
return utils::dependencies::CombinedKeyList<std::conditional_t<ThermodynamicBoundaryField<Terms>, utils::dependencies::KeyList<Terms>, utils::dependencies::KeyList<>>...>{};
|
||||
}
|
||||
|
||||
/** @brief Return a unique field for a quantity; the chosen residual carrier supplies its own quantity unambiguously. */
|
||||
template <typename Quantity, typename Carrier, typename... Fields>
|
||||
[[nodiscard]] consteval auto quantity_field(utils::dependencies::KeyList<Fields...>) {
|
||||
using Matches = decltype(std::tuple_cat(std::declval<std::conditional_t<std::same_as<Quantity, typename Fields::Quantity>, std::tuple<Fields>, std::tuple<>>>()...));
|
||||
if constexpr (std::same_as<Quantity, typename Carrier::Quantity>) {
|
||||
return std::type_identity<Carrier>{};
|
||||
}
|
||||
else if constexpr (std::tuple_size_v<Matches> == 1) {
|
||||
return std::type_identity<std::tuple_element_t<0, Matches>>{};
|
||||
}
|
||||
else {
|
||||
return std::type_identity<void>{};
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Bind only the quantities read by this candidate, ignoring unrelated fields and quantities. */
|
||||
template <typename Carrier, typename Fields, typename... Quantities>
|
||||
[[nodiscard]] consteval auto relation_bindings(utils::dependencies::KeyList<Quantities...>) {
|
||||
using namespace utils::dependencies;
|
||||
|
||||
if constexpr ((!std::same_as<void, typename decltype(quantity_field<Quantities, Carrier>(Fields{}))::type> && ...)) {
|
||||
return std::type_identity<Bindings<Binding<Quantities, typename decltype(quantity_field<Quantities, Carrier>(Fields{}))::type>...>>{};
|
||||
}
|
||||
else {
|
||||
return std::type_identity<void>{};
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Test a declared inverse and every required analytic partial before publishing its formulation. */
|
||||
template <typename EOS, typename Carrier, typename Fields, typename Relation>
|
||||
[[nodiscard]] consteval auto inverse_candidate() {
|
||||
using namespace utils::dependencies;
|
||||
using Quantities = CombinedKeyList<KeyList<typename Carrier::Quantity>, typename surface::detail::RelationQuantities<Relation>::NonPressureInputs>;
|
||||
using Bindings = typename decltype(relation_bindings<Carrier, Fields>(Quantities{}))::type;
|
||||
|
||||
if constexpr (SurfaceStateBindings<Bindings> && surface::detail::matches_pressure_relation<EOS, Bindings, typename Carrier::Quantity, Relation>()) {
|
||||
return std::tuple<std::type_identity<InversePressureSurfaceFormulation<EOS, Bindings, Relation>>>{};
|
||||
}
|
||||
else {
|
||||
return std::tuple<>{};
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief A direct pressure relation may replace only a row belonging to one of its solved inputs. */
|
||||
template <typename EOS, typename Carrier, typename Fields, typename Relation>
|
||||
[[nodiscard]] consteval auto direct_candidate() {
|
||||
using namespace utils::dependencies;
|
||||
using Quantities = CombinedKeyList<KeyList<typename Carrier::Quantity>, typename surface::detail::RelationQuantities<Relation>::AllInputs>;
|
||||
using Bindings = typename decltype(relation_bindings<Carrier, Fields>(Quantities{}))::type;
|
||||
|
||||
if constexpr (DirectPressureSurfaceCompilable<EOS, Bindings, typename Carrier::Quantity, Relation>) {
|
||||
return std::tuple<std::type_identity<DirectPressureSurfaceFormulation<EOS, Bindings, typename Carrier::Quantity, Relation>>>{};
|
||||
}
|
||||
else {
|
||||
return std::tuple<>{};
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Keep every admissible relation; catalog order never selects a physical implementation. */
|
||||
template <bool Inverse, typename EOS, typename Carrier, typename Fields, typename... Relations>
|
||||
[[nodiscard]] consteval auto relation_candidates(eos::relations::RelationCatalog<Relations...>) {
|
||||
if constexpr (Inverse) {
|
||||
return std::tuple_cat(inverse_candidate<EOS, Carrier, Fields, Relations>()...);
|
||||
}
|
||||
else {
|
||||
return std::tuple_cat(direct_candidate<EOS, Carrier, Fields, Relations>()...);
|
||||
}
|
||||
}
|
||||
}
|
||||
114
src/include/serif/geometry/surface/detail/batches.hpp
Normal file
114
src/include/serif/geometry/surface/detail/batches.hpp
Normal file
@@ -0,0 +1,114 @@
|
||||
#pragma once
|
||||
|
||||
#include <array>
|
||||
|
||||
#include "serif/geometry/surface/detail/evaluation.hpp"
|
||||
#include "serif/geometry/surface/samples.hpp"
|
||||
|
||||
namespace serif::geometry::surface::detail
|
||||
{
|
||||
/** @brief All EOS columns must describe the same ordered, nonempty set of boundary points. */
|
||||
template <typename... Terms>
|
||||
[[nodiscard]] std::expected<std::size_t, SurfaceFailure> surface_sample_count(const BoundaryValues<Terms>&... columns)
|
||||
{
|
||||
if constexpr (sizeof...(Terms) == 0)
|
||||
{
|
||||
// A pressure-only relation has one constant response, broadcast over the carrier samples later.
|
||||
return 1;
|
||||
}
|
||||
else
|
||||
{
|
||||
const std::array<std::size_t, sizeof...(Terms)> counts{columns.values.size()...};
|
||||
if (counts.front() == 0)
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::empty_surface, .message = "A surface calculation requires at least one boundary sample."});
|
||||
}
|
||||
for (const std::size_t count : counts)
|
||||
{
|
||||
if (count != counts.front())
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::inconsistent_sample_count, .message = "Surface terms must be sampled at the same ordered boundary points; their sample counts differ."});
|
||||
}
|
||||
}
|
||||
return counts.front();
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Assemble the conventional equation point from typed columns; no mesh assumptions enter the physical calculation. */
|
||||
template <PressureSurfaceFormulationType Formulation, typename... Terms>
|
||||
[[nodiscard]] SurfacePoint<Formulation> surface_point(const std::size_t sample, const BoundaryValues<Terms>&... columns)
|
||||
{
|
||||
SurfacePoint<Formulation> point{};
|
||||
((point.fields.template get<Terms>() = columns.values[sample]), ...);
|
||||
return point;
|
||||
}
|
||||
|
||||
/** @brief Collect one ordered batch, attaching the first failure's sample index before any result is published. */
|
||||
template <typename Value, typename Evaluate>
|
||||
[[nodiscard]] std::expected<std::vector<Value>, SurfaceFailure> collect_surface_samples(std::size_t count, Evaluate evaluate)
|
||||
{
|
||||
std::vector<Value> result;
|
||||
result.reserve(count);
|
||||
for (std::size_t sample = 0; sample < count; ++sample)
|
||||
{
|
||||
auto value = evaluate(sample);
|
||||
if (!value)
|
||||
{
|
||||
SurfaceFailure failure = value.error();
|
||||
failure.sample = sample;
|
||||
return std::unexpected(failure);
|
||||
}
|
||||
result.push_back(std::move(*value));
|
||||
}
|
||||
return result;
|
||||
}
|
||||
|
||||
/** @brief Evaluate the EOS once per state sample, or once for a genuinely pressure-only relation. */
|
||||
template <PressureSurfaceFormulationType Formulation, typename... Terms>
|
||||
[[nodiscard]] std::expected<std::vector<SurfaceEOSResponse<Formulation>>, SurfaceFailure> evaluate_response_samples(const typename Formulation::EquationOfState& equation_of_state, const dimensions::PressureValue pressure, const BoundaryValues<Terms>&... columns)
|
||||
{
|
||||
const auto count = surface_sample_count(columns...);
|
||||
if (!count)
|
||||
{
|
||||
return std::unexpected(count.error());
|
||||
}
|
||||
return collect_surface_samples<SurfaceEOSResponse<Formulation>>(*count, [&](std::size_t sample)
|
||||
{
|
||||
return evaluate_surface_response<Formulation>(equation_of_state, pressure, surface_point<Formulation>(sample, columns...));
|
||||
});
|
||||
}
|
||||
|
||||
/** @brief Broadcast only a pressure-only relation; every state-dependent carrier sample needs its own response. */
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] std::expected<std::vector<SurfaceLinearization<Formulation>>, SurfaceFailure> linearize_carrier_samples(const std::vector<SurfaceEOSResponse<Formulation>>& responses, const BoundaryValues<typename Formulation::CarrierTerm>& carrier)
|
||||
{
|
||||
constexpr bool constant_response = std::same_as<typename Formulation::EOSStateTerms, utils::dependencies::KeyList<>>;
|
||||
const auto count = surface_sample_count(carrier);
|
||||
if (!count)
|
||||
{
|
||||
return std::unexpected(count.error());
|
||||
}
|
||||
if (responses.size() != (constant_response ? 1 : *count))
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::inconsistent_sample_count, .message = "The carrier and EOS inputs must refer to the same boundary points."});
|
||||
}
|
||||
return collect_surface_samples<SurfaceLinearization<Formulation>>(*count, [&](std::size_t sample) -> std::expected<SurfaceLinearization<Formulation>, SurfaceFailure>
|
||||
{
|
||||
if (!std::isfinite(carrier.values[sample]))
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::non_finite_state, .message = "A boundary carrier value is non-finite."});
|
||||
}
|
||||
return checked_linearization<Formulation>(dimensions::PressureValue{0.0}, surface_point<Formulation>(sample, carrier), responses[constant_response ? 0 : sample]);
|
||||
});
|
||||
}
|
||||
|
||||
/** @brief Changing prescribed pressure reuses all values and partials of a direct-pressure EOS relation. */
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] std::expected<std::vector<SurfaceLinearization<Formulation>>, SurfaceFailure> linearize_pressure_samples(const IsobaricSurface& condition, const std::vector<SurfaceEOSResponse<Formulation>>& responses)
|
||||
{
|
||||
return collect_surface_samples<SurfaceLinearization<Formulation>>(responses.size(), [&](std::size_t sample)
|
||||
{
|
||||
return checked_linearization<Formulation>(condition.target_pressure(), SurfacePoint<Formulation>{}, responses[sample]);
|
||||
});
|
||||
}
|
||||
} // namespace serif::geometry::surface::detail
|
||||
135
src/include/serif/geometry/surface/detail/evaluation.hpp
Normal file
135
src/include/serif/geometry/surface/detail/evaluation.hpp
Normal file
@@ -0,0 +1,135 @@
|
||||
#pragma once
|
||||
|
||||
#include <cmath>
|
||||
#include <expected>
|
||||
|
||||
#include "serif/eos/evaluation.hpp"
|
||||
#include "serif/geometry/surface/linearization.hpp"
|
||||
|
||||
namespace serif::geometry::surface::detail
|
||||
{
|
||||
/** @brief Read one thermodynamic role from its bound term, preserving the EOS argument's physical type. */
|
||||
template <typename Quantity, PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] dimensions::QuantityValue<Quantity> state_value(const SurfacePoint<Formulation>& point)
|
||||
{
|
||||
using Term = utils::dependencies::bound_type_t<typename Formulation::StateBindings, Quantity>;
|
||||
return dimensions::QuantityValue<Quantity>{point.fields.template get<Term>()};
|
||||
}
|
||||
|
||||
/** @brief Only the carrier-resolved formulation substitutes prescribed pressure for an EOS input. */
|
||||
template <typename Quantity, PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] dimensions::QuantityValue<Quantity> relation_input(const dimensions::PressureValue pressure, const SurfacePoint<Formulation>& point)
|
||||
{
|
||||
if constexpr (Formulation::resolves_carrier && std::same_as<Quantity, dimensions::Pressure>)
|
||||
{
|
||||
return pressure;
|
||||
}
|
||||
else
|
||||
{
|
||||
return state_value<Quantity, Formulation>(point);
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Validate precisely the state read by this calculation, leaving unrelated state outside its dependencies. */
|
||||
template <PressureSurfaceFormulationType Formulation, typename... Quantities>
|
||||
[[nodiscard]] bool finite_state(const SurfacePoint<Formulation>& point, utils::dependencies::KeyList<Quantities...>)
|
||||
{
|
||||
return (std::isfinite(state_value<Quantities, Formulation>(point).value()) && ...);
|
||||
}
|
||||
|
||||
/** @brief Add each quantity's chain-rule contribution, including distinct roles bound to the same unknown. */
|
||||
template <typename Quantity, PressureSurfaceFormulationType Formulation, typename... Inputs>
|
||||
void accumulate_partial(const typename Formulation::EquationOfState& equation_of_state, const dimensions::PressureValue pressure, const SurfacePoint<Formulation>& point, SurfaceEOSResponse<Formulation>& response, dimensions::QuantityList<Inputs...>)
|
||||
{
|
||||
using Output = typename Formulation::Relation::OutputQuantity;
|
||||
using Term = utils::dependencies::bound_type_t<typename Formulation::StateBindings, Quantity>;
|
||||
|
||||
const dimensions::PartialDerivative<Output, Quantity> derivative = eos::partial_derivative<Output, Quantity>(equation_of_state, relation_input<Inputs, Formulation>(pressure, point)...);
|
||||
response.field_derivatives.template get<Term>() += derivative.value();
|
||||
}
|
||||
|
||||
/** @brief Evaluate one selected relation and all required partials without any hidden EOS inversion. */
|
||||
template <PressureSurfaceFormulationType Formulation, typename... Inputs, typename... StateQuantities>
|
||||
[[nodiscard]] SurfaceEOSResponse<Formulation> evaluate_relation(const typename Formulation::EquationOfState& equation_of_state, const dimensions::PressureValue pressure, const SurfacePoint<Formulation>& point, dimensions::QuantityList<Inputs...> inputs, utils::dependencies::KeyList<StateQuantities...>)
|
||||
{
|
||||
using Output = typename Formulation::Relation::OutputQuantity;
|
||||
|
||||
SurfaceEOSResponse<Formulation> response{};
|
||||
response.value = eos::evaluate<Output>(equation_of_state, relation_input<Inputs, Formulation>(pressure, point)...);
|
||||
(accumulate_partial<StateQuantities, Formulation>(equation_of_state, pressure, point, response, inputs), ...);
|
||||
return response;
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief Reject non-finite thermodynamic inputs, EOS outputs, and partials as values.
|
||||
* @note Existing EOS implementations can still throw their documented domain errors; these propagate to the caller.
|
||||
*/
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] std::expected<SurfaceEOSResponse<Formulation>, SurfaceFailure> evaluate_surface_response(const typename Formulation::EquationOfState& equation_of_state, const dimensions::PressureValue pressure, const SurfacePoint<Formulation>& point)
|
||||
{
|
||||
if constexpr (!Formulation::uses_eos)
|
||||
{
|
||||
// A pressure unknown uses the prescribed pressure itself as its target, without an EOS relation.
|
||||
return SurfaceEOSResponse<Formulation>{.value = pressure};
|
||||
}
|
||||
else
|
||||
{
|
||||
if (!finite_state<Formulation>(point, typename Formulation::EOSStateQuantities{}))
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::non_finite_state, .message = "The surface EOS inputs must be finite; supply an admissible boundary state."});
|
||||
}
|
||||
|
||||
const SurfaceEOSResponse<Formulation> response = evaluate_relation<Formulation>(equation_of_state, pressure, point, typename Formulation::Relation::InputQuantities{}, typename Formulation::EOSStateQuantities{});
|
||||
bool finite = std::isfinite(response.value.value());
|
||||
for (const double derivative : response.field_derivatives.values)
|
||||
{
|
||||
finite = finite && std::isfinite(derivative);
|
||||
}
|
||||
if (!finite)
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::non_finite_eos_response, .message = "The surface EOS returned a non-finite value or partial derivative; check its thermodynamic domain."});
|
||||
}
|
||||
return response;
|
||||
}
|
||||
}
|
||||
|
||||
/** @brief Form the declared residual from an already evaluated EOS response, keeping the carrier identity explicit. */
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] SurfaceLinearization<Formulation> linearize_surface(const dimensions::PressureValue pressure, const SurfacePoint<Formulation>& point, const SurfaceEOSResponse<Formulation>& response)
|
||||
{
|
||||
SurfaceLinearization<Formulation> result{};
|
||||
result.field_derivatives = response.field_derivatives;
|
||||
|
||||
if constexpr (Formulation::resolves_carrier)
|
||||
{
|
||||
result.residual = state_value<typename Formulation::CarrierQuantity, Formulation>(point).value() - response.value.value();
|
||||
for (double& derivative : result.field_derivatives.values)
|
||||
{
|
||||
derivative = -derivative;
|
||||
}
|
||||
result.field_derivatives.template get<typename Formulation::CarrierTerm>() += 1.0;
|
||||
}
|
||||
else
|
||||
{
|
||||
result.residual = response.value.value() - pressure.value();
|
||||
}
|
||||
return result;
|
||||
}
|
||||
|
||||
/** @brief Catch overflow in residual subtraction or derivative accumulation before publishing a linearization. */
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] std::expected<SurfaceLinearization<Formulation>, SurfaceFailure> checked_linearization(const dimensions::PressureValue pressure, const SurfacePoint<Formulation>& point, const SurfaceEOSResponse<Formulation>& response)
|
||||
{
|
||||
const SurfaceLinearization<Formulation> result = linearize_surface<Formulation>(pressure, point, response);
|
||||
bool finite = std::isfinite(result.residual);
|
||||
for (const double derivative : result.field_derivatives.values)
|
||||
{
|
||||
finite = finite && std::isfinite(derivative);
|
||||
}
|
||||
if (!finite)
|
||||
{
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::non_finite_linearization, .message = "The surface residual or derivative is non-finite; reduce the trial state before using this linearization."});
|
||||
}
|
||||
return result;
|
||||
}
|
||||
}
|
||||
97
src/include/serif/geometry/surface/detail/relations.hpp
Normal file
97
src/include/serif/geometry/surface/detail/relations.hpp
Normal file
@@ -0,0 +1,97 @@
|
||||
#pragma once
|
||||
|
||||
#include <tuple>
|
||||
#include <type_traits>
|
||||
|
||||
#include "serif/discretization/forms/value.hpp"
|
||||
#include "serif/eos/models/concepts.hpp"
|
||||
#include "serif/utils/dependencies/key_list_operations.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
/** @brief Scalar thermodynamic data may come from any declared scalar space; boundary enforcement is a separate backend capability. */
|
||||
template <typename Term>
|
||||
concept SurfaceScalarTerm = discretization::blocks::Term<Term> && requires {
|
||||
typename discretization::forms::OperandValue<discretization::forms::Operand<Term>>;
|
||||
requires std::same_as<discretization::forms::OperandValue<discretization::forms::Operand<Term>>, double>;
|
||||
};
|
||||
|
||||
namespace detail {
|
||||
template <typename Bindings>
|
||||
struct AreSurfaceBindings : std::false_type {
|
||||
};
|
||||
|
||||
template <typename... Entries>
|
||||
struct AreSurfaceBindings<utils::dependencies::Bindings<Entries...>> : std::bool_constant<
|
||||
sizeof...(Entries) != 0 && (dimensions::ThermodynamicQuantityType<typename Entries::Key> && ...) &&
|
||||
(SurfaceScalarTerm<typename Entries::Target> && ...)> {
|
||||
};
|
||||
}
|
||||
|
||||
/** @brief Reuse dependency bindings, retaining distinct quantity roles even when they refer to one term. */
|
||||
template <typename Bindings>
|
||||
concept SurfaceStateBindings = utils::dependencies::IsBindings<Bindings> && detail::AreSurfaceBindings<std::remove_cvref_t<Bindings>>::value;
|
||||
|
||||
/** @brief The quantity enforcing the boundary must have an actual equation row, not just sampled values. */
|
||||
template <typename Bindings, typename Quantity>
|
||||
concept SurfaceResidualRow = utils::dependencies::BindingsWithBinding<Bindings, Quantity> &&
|
||||
discretization::blocks::HasResidualBlock<utils::dependencies::bound_type_t<Bindings, Quantity>>;
|
||||
|
||||
namespace detail {
|
||||
template <typename Relation>
|
||||
struct RelationQuantities;
|
||||
|
||||
template <typename Output, typename... Inputs>
|
||||
struct RelationQuantities<eos::relations::Relation<Output, Inputs...>> {
|
||||
using AllInputs = utils::dependencies::KeyList<Inputs...>;
|
||||
using NonPressureInputs = utils::dependencies::CombinedKeyList<
|
||||
std::conditional_t<std::same_as<Inputs, dimensions::Pressure>, utils::dependencies::KeyList<>, utils::dependencies::KeyList<Inputs>>...>;
|
||||
};
|
||||
|
||||
/** @brief Check only the derivatives of unknown state; prescribed pressure has no Newton column. */
|
||||
template <typename EOS, typename Relation, typename... Quantities>
|
||||
[[nodiscard]] consteval bool supports_partials(utils::dependencies::KeyList<Quantities...>) {
|
||||
return (eos::models::EOSSupportingPartialDerivative<EOS, Relation, Quantities> && ...);
|
||||
}
|
||||
|
||||
/** @brief A resolved carrier must be supplied by one pressure-input relation, with all remaining inputs and partials available. */
|
||||
template <typename EOS, typename Bindings, typename Carrier, typename Relation>
|
||||
[[nodiscard]] consteval bool matches_pressure_relation() {
|
||||
using Output = typename Relation::OutputQuantity;
|
||||
using Inputs = typename RelationQuantities<Relation>::NonPressureInputs;
|
||||
|
||||
return !std::same_as<Output, dimensions::Pressure> &&
|
||||
(std::same_as<Carrier, void> || std::same_as<Carrier, Output>) &&
|
||||
eos::relations::RelationContainingInput<Relation, dimensions::Pressure> &&
|
||||
SurfaceResidualRow<Bindings, Output> &&
|
||||
utils::dependencies::SatisfiesKeyList<Bindings, Inputs> && supports_partials<EOS, Relation>(Inputs{});
|
||||
}
|
||||
|
||||
template <typename EOS, typename Bindings, typename Carrier, typename Catalog>
|
||||
struct MatchingPressureRelations;
|
||||
|
||||
template <typename EOS, typename Bindings, typename Carrier, typename... Relations>
|
||||
struct MatchingPressureRelations<EOS, Bindings, Carrier, eos::relations::RelationCatalog<Relations...>> {
|
||||
using Matches = decltype(std::tuple_cat(std::declval<std::conditional_t<matches_pressure_relation<EOS, Bindings, Carrier, Relations>(), std::tuple<Relations>, std::tuple<>>>()...));
|
||||
static constexpr std::size_t count = std::tuple_size_v<Matches>;
|
||||
};
|
||||
|
||||
template <typename Terms>
|
||||
struct SurfaceOperands;
|
||||
|
||||
template <typename... Terms>
|
||||
struct SurfaceOperands<utils::dependencies::KeyList<Terms...>> {
|
||||
using Type = utils::types::TypeList<discretization::forms::Operand<Terms>...>;
|
||||
};
|
||||
|
||||
template <typename Terms>
|
||||
struct UniqueTerms;
|
||||
|
||||
template <typename... Terms>
|
||||
struct UniqueTerms<utils::types::TypeList<Terms...>> {
|
||||
using Type = utils::dependencies::CombinedKeyList<utils::dependencies::KeyList<Terms>...>;
|
||||
};
|
||||
|
||||
template <typename Bindings, typename Quantities>
|
||||
using BoundTerms = typename UniqueTerms<utils::dependencies::ResolvedTargets<Bindings, Quantities>>::Type;
|
||||
}
|
||||
}
|
||||
26
src/include/serif/geometry/surface/evaluation.hpp
Normal file
26
src/include/serif/geometry/surface/evaluation.hpp
Normal file
@@ -0,0 +1,26 @@
|
||||
#pragma once
|
||||
|
||||
#include "serif/geometry/surface/detail/evaluation.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
/**
|
||||
* @brief Evaluate the chosen isobaric residual and its analytical derivatives at one boundary point.
|
||||
* @param equation_of_state The EOS declared by Formulation.
|
||||
* @param condition The prescribed surface pressure, independent of the choice of carrier.
|
||||
* @param point Values of the formulation's bound scalar terms, already sampled at a common boundary point.
|
||||
* @return A conventional equation linearization, or a numerical rejection. No preparation or borrowed state is retained.
|
||||
* @note Domain exceptions from an existing throwing EOS propagate; no exception is used for surface validation.
|
||||
*/
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
[[nodiscard]] std::expected<SurfaceLinearization<Formulation>, SurfaceFailure> evaluate_pressure_surface(const typename Formulation::EquationOfState& equation_of_state, const IsobaricSurface& condition, const SurfacePoint<Formulation>& point) {
|
||||
if (!detail::finite_state<Formulation>(point, typename Formulation::StateQuantities{})) {
|
||||
return std::unexpected(SurfaceFailure{.reason = SurfaceFailureReason::non_finite_state, .message = "The surface state must be finite; supply finite values for every required thermodynamic quantity."});
|
||||
}
|
||||
|
||||
const auto response = detail::evaluate_surface_response<Formulation>(equation_of_state, condition.target_pressure(), point);
|
||||
if (!response) {
|
||||
return std::unexpected(response.error());
|
||||
}
|
||||
return detail::checked_linearization<Formulation>(condition.target_pressure(), point, *response);
|
||||
}
|
||||
}
|
||||
@@ -0,0 +1,40 @@
|
||||
#pragma once
|
||||
|
||||
#include <array>
|
||||
#include <type_traits>
|
||||
|
||||
#include "serif/dimensions/dimensions.hpp"
|
||||
#include "serif/dimensions/type_alias.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
struct PressureSurfaceDescriptor final {
|
||||
double target_pressure;
|
||||
};
|
||||
|
||||
class IsobaricSurface final {
|
||||
public:
|
||||
using PhysicalQuantity = dimensions::Pressure;
|
||||
using TargetValue = dimensions::PressureValue;
|
||||
|
||||
public:
|
||||
struct Parameters final {
|
||||
dimensions::PressureValue PSurf;
|
||||
};
|
||||
|
||||
public:
|
||||
explicit IsobaricSurface(const Parameters parameters);
|
||||
explicit IsobaricSurface (const TargetValue targetPressure);
|
||||
|
||||
[[nodiscard]] TargetValue target_pressure() const noexcept;
|
||||
[[nodiscard]] std::array<double, 1> continuation_identity() const noexcept;
|
||||
[[nodiscard]] PressureSurfaceDescriptor descriptor() const noexcept;
|
||||
|
||||
private:
|
||||
TargetValue m_targetPressure;
|
||||
};
|
||||
|
||||
template <typename Candidate>
|
||||
concept IsobaricSurfaceType = std::same_as<std::remove_cvref_t<Candidate>, IsobaricSurface>;
|
||||
|
||||
using Isobaric = IsobaricSurface;
|
||||
}
|
||||
41
src/include/serif/geometry/surface/linearization.hpp
Normal file
41
src/include/serif/geometry/surface/linearization.hpp
Normal file
@@ -0,0 +1,41 @@
|
||||
#pragma once
|
||||
|
||||
#include <cstddef>
|
||||
#include <limits>
|
||||
#include <string_view>
|
||||
|
||||
#include "serif/equations/point.hpp"
|
||||
#include "serif/geometry/surface/compilation.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
enum class SurfaceFailureReason {
|
||||
non_finite_state,
|
||||
non_finite_eos_response,
|
||||
non_finite_linearization,
|
||||
inconsistent_sample_count,
|
||||
empty_surface
|
||||
};
|
||||
|
||||
/** @brief A rejected thermodynamic candidate, with a sample index when evaluating a boundary batch. */
|
||||
struct SurfaceFailure {
|
||||
SurfaceFailureReason reason;
|
||||
std::string_view message;
|
||||
std::size_t sample = std::numeric_limits<std::size_t>::max();
|
||||
};
|
||||
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
using SurfacePoint = equations::EquationPoint<typename Formulation::Operands>;
|
||||
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
using SurfaceLinearization = equations::PointLinearization<typename Formulation::Operands>;
|
||||
|
||||
/**
|
||||
* @brief An EOS value and all of its local derivatives, evaluated together at one state.
|
||||
* @note Derivatives sharing a term are accumulated; the later residual calculation supplies the sign and carrier identity.
|
||||
*/
|
||||
template <PressureSurfaceFormulationType Formulation>
|
||||
struct SurfaceEOSResponse {
|
||||
dimensions::QuantityValue<typename Formulation::Relation::OutputQuantity> value{0.0};
|
||||
equations::OperandValues<typename Formulation::Operands, double> field_derivatives{};
|
||||
};
|
||||
}
|
||||
70
src/include/serif/geometry/surface/samples.hpp
Normal file
70
src/include/serif/geometry/surface/samples.hpp
Normal file
@@ -0,0 +1,70 @@
|
||||
#pragma once
|
||||
|
||||
#include <utility>
|
||||
#include <vector>
|
||||
|
||||
#include "serif/geometry/surface/linearization.hpp"
|
||||
#include "serif/utils/dependencies/resource.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
namespace keys {
|
||||
struct EquationOfState {
|
||||
};
|
||||
struct PressureCondition {
|
||||
};
|
||||
struct EOSResponse {
|
||||
};
|
||||
struct Linearization {
|
||||
};
|
||||
}
|
||||
|
||||
/**
|
||||
* @brief Values of one scalar term at an ordered collection of boundary points.
|
||||
* @note Every term in a preparation must use the same point ordering. Sampling and finite-element trace choices belong to the backend.
|
||||
*/
|
||||
template <SurfaceScalarTerm Term>
|
||||
struct BoundaryValues {
|
||||
std::vector<double> values;
|
||||
|
||||
/// Compare the actual sampled values, including their count, without a lossy checksum.
|
||||
[[nodiscard]] bool operator==(const BoundaryValues&) const = default;
|
||||
};
|
||||
|
||||
namespace detail {
|
||||
struct BoundaryValueTracking : utils::dependencies::ValueSnapshot {
|
||||
/// Numerical samples have value semantics; rebuilding an equal context need not reevaluate the physics.
|
||||
template <SurfaceScalarTerm Term>
|
||||
[[nodiscard]] const void* identity(const BoundaryValues<Term>&) const {
|
||||
return nullptr;
|
||||
}
|
||||
};
|
||||
|
||||
struct PressureConditionTracking {
|
||||
/// Only prescribed pressure can affect the isobaric calculation.
|
||||
[[nodiscard]] dimensions::PressureValue capture(const IsobaricSurface& condition) const {
|
||||
return condition.target_pressure();
|
||||
}
|
||||
|
||||
/// Compare pressure exactly, retaining the dependency engine's normal change semantics.
|
||||
[[nodiscard]] bool matches(const dimensions::PressureValue captured, const IsobaricSurface& condition) const {
|
||||
return captured == condition.target_pressure();
|
||||
}
|
||||
|
||||
/// An equivalent pressure condition may be supplied by a newly constructed context.
|
||||
[[nodiscard]] const void* identity(const IsobaricSurface&) const {
|
||||
return nullptr;
|
||||
}
|
||||
};
|
||||
}
|
||||
|
||||
/** @brief Supply owned boundary samples to an ordinary dependency context, keyed by their physical term. */
|
||||
template <SurfaceScalarTerm Term>
|
||||
[[nodiscard]] auto boundary_values(std::vector<double> values) {
|
||||
return utils::dependencies::own<Term>(BoundaryValues<Term>{.values = std::move(values)}, detail::BoundaryValueTracking{});
|
||||
}
|
||||
|
||||
/** @brief Supply the isobaric condition as an owned dependency source, with pressure tracked by value. */
|
||||
[[nodiscard]] inline auto pressure_condition(IsobaricSurface condition) {
|
||||
return utils::dependencies::own<keys::PressureCondition>(std::move(condition), detail::PressureConditionTracking{});
|
||||
}
|
||||
}
|
||||
17
src/include/serif/geometry/surface/surface.hpp
Normal file
17
src/include/serif/geometry/surface/surface.hpp
Normal file
@@ -0,0 +1,17 @@
|
||||
#pragma once
|
||||
|
||||
/**
|
||||
* @file surface.hpp
|
||||
* @brief Isobaric conditions, EOS-selected formulations, and reusable boundary calculations.
|
||||
*
|
||||
* IsobaricSurface describes the prescribed pressure. PressureSurfaceFormulation selects one compatible EOS relation and its carrier;
|
||||
* DirectPressureSurfaceFormulation explicitly chooses a forward pressure residual and its equation row. Neither selects a finite-element trace.
|
||||
*
|
||||
* Use evaluate_pressure_surface for one sampled point, or make_pressure_surface_graph for an ordered batch. The latter returns an ordinary
|
||||
* dependency graph: the caller supplies the EOS with its tracking policy, the pressure condition, and typed boundary values. Results use the
|
||||
* same PointLinearization representation as other equations. Finite-element row enforcement remains a separate numerical capability.
|
||||
*/
|
||||
|
||||
#include "serif/geometry/surface/compilation.hpp"
|
||||
#include "serif/geometry/surface/computations.hpp"
|
||||
#include "serif/geometry/surface/evaluation.hpp"
|
||||
@@ -0,0 +1,31 @@
|
||||
#include <cmath>
|
||||
#include <stdexcept>
|
||||
#include <array>
|
||||
|
||||
#include "serif/geometry/surface/isobaric.hpp"
|
||||
#include "serif/utils/error/errors.hpp"
|
||||
|
||||
namespace serif::geometry::surface {
|
||||
using namespace utils::errors;
|
||||
|
||||
IsobaricSurface::IsobaricSurface(const Parameters parameters) : IsobaricSurface(parameters.PSurf) {}
|
||||
|
||||
IsobaricSurface::IsobaricSurface(const TargetValue targetPressure) : m_targetPressure(targetPressure) {
|
||||
serif_ensure<SERiFErrorCode::NON_FINITE_VALUE_ERROR, std::invalid_argument>(std::isfinite(targetPressure.value()), "The target surface pressure must be finite, got {}", targetPressure.value());
|
||||
serif_ensure<SERiFErrorCode::NEGATIVE_VALUE_ERROR, std::invalid_argument>(targetPressure.value() >= 0, "The target surface pressure must be non-negative, got {}", targetPressure.value());
|
||||
}
|
||||
|
||||
IsobaricSurface::TargetValue IsobaricSurface::target_pressure() const noexcept {
|
||||
return m_targetPressure;
|
||||
}
|
||||
|
||||
std::array<double, 1> IsobaricSurface::continuation_identity() const noexcept {
|
||||
return {m_targetPressure.value()};
|
||||
}
|
||||
|
||||
PressureSurfaceDescriptor IsobaricSurface::descriptor() const noexcept {
|
||||
return PressureSurfaceDescriptor{
|
||||
.target_pressure = m_targetPressure.value()
|
||||
};
|
||||
}
|
||||
}
|
||||
|
||||
Reference in New Issue
Block a user