feat(beam-reference-qualification): step 0 — comparison-metric-and-entity-matching
This commit is contained in:
KOKO\Mimi
@@ -0,0 +1,442 @@
|
||||
#include <fesa/validation/comparison.hpp>
|
||||
|
||||
#include <algorithm>
|
||||
#include <array>
|
||||
#include <cmath>
|
||||
#include <cstddef>
|
||||
#include <iomanip>
|
||||
#include <limits>
|
||||
#include <optional>
|
||||
#include <set>
|
||||
#include <sstream>
|
||||
#include <string>
|
||||
#include <string_view>
|
||||
#include <tuple>
|
||||
#include <utility>
|
||||
#include <vector>
|
||||
|
||||
namespace fesa {
|
||||
namespace {
|
||||
|
||||
using PositionKey = std::tuple<
|
||||
ReferenceQuantity,
|
||||
std::string,
|
||||
std::int64_t,
|
||||
std::optional<std::int64_t>>;
|
||||
|
||||
void add_failure(
|
||||
std::vector<Diagnostic>& failures,
|
||||
std::string code,
|
||||
std::string message) {
|
||||
failures.push_back({
|
||||
DiagnosticStage::validation,
|
||||
Severity::error,
|
||||
std::move(code),
|
||||
std::move(message),
|
||||
std::nullopt,
|
||||
});
|
||||
}
|
||||
|
||||
std::string_view quantity_name(const ReferenceQuantity quantity) {
|
||||
switch (quantity) {
|
||||
case ReferenceQuantity::displacement:
|
||||
return "displacement";
|
||||
case ReferenceQuantity::reaction:
|
||||
return "reaction";
|
||||
case ReferenceQuantity::internal_force:
|
||||
return "internal_force";
|
||||
case ReferenceQuantity::centroid_stress:
|
||||
return "centroid_stress";
|
||||
}
|
||||
return "unknown";
|
||||
}
|
||||
|
||||
std::size_t expected_component_count(
|
||||
const ReferenceQuantity quantity) {
|
||||
switch (quantity) {
|
||||
case ReferenceQuantity::displacement:
|
||||
case ReferenceQuantity::reaction:
|
||||
case ReferenceQuantity::internal_force:
|
||||
return 6;
|
||||
case ReferenceQuantity::centroid_stress:
|
||||
return 1;
|
||||
}
|
||||
return 0;
|
||||
}
|
||||
|
||||
std::string component_name(
|
||||
const ReferenceQuantity quantity,
|
||||
const std::size_t index) {
|
||||
switch (quantity) {
|
||||
case ReferenceQuantity::displacement:
|
||||
return std::string{NodalFrame::displacement_components[index]};
|
||||
case ReferenceQuantity::reaction:
|
||||
return std::string{NodalFrame::reaction_components[index]};
|
||||
case ReferenceQuantity::internal_force:
|
||||
return std::string{
|
||||
BeamElementFrame::section_force_components[index]};
|
||||
case ReferenceQuantity::centroid_stress:
|
||||
return std::string{BeamElementFrame::axial_stress_component};
|
||||
}
|
||||
return "unknown";
|
||||
}
|
||||
|
||||
std::string number_text(const double value) {
|
||||
std::ostringstream stream;
|
||||
stream << std::setprecision(std::numeric_limits<double>::max_digits10)
|
||||
<< value;
|
||||
return stream.str();
|
||||
}
|
||||
|
||||
std::string tolerance_text(const Tolerance tolerance) {
|
||||
return " relative_tolerance=" + number_text(tolerance.relative) +
|
||||
" absolute_scale=" + number_text(tolerance.absolute_scale);
|
||||
}
|
||||
|
||||
std::string position_text(
|
||||
const ReferenceQuantity quantity,
|
||||
const ResultPosition& position) {
|
||||
std::string text = "quantity=" + std::string{quantity_name(quantity)} +
|
||||
" instance=" + position.instance_name +
|
||||
" entity=" + std::to_string(position.entity_label);
|
||||
if (position.end_node_label.has_value()) {
|
||||
text += " end_node=" +
|
||||
std::to_string(*position.end_node_label);
|
||||
} else {
|
||||
text += " end_node=n/a";
|
||||
}
|
||||
return text;
|
||||
}
|
||||
|
||||
std::string scalar_failure_text(
|
||||
const ComparisonSample& sample,
|
||||
const std::size_t component,
|
||||
const double normalized_error) {
|
||||
return position_text(sample.quantity, sample.position) +
|
||||
" component=" + component_name(sample.quantity, component) +
|
||||
" reference=" + number_text(sample.reference[component]) +
|
||||
" actual=" + number_text(sample.actual[component]) +
|
||||
" normalized_error=" + number_text(normalized_error) +
|
||||
tolerance_text(sample.tolerance);
|
||||
}
|
||||
|
||||
std::string unevaluable_failure_text(
|
||||
const ReferenceQuantity quantity,
|
||||
const ResultPosition& position,
|
||||
const Tolerance tolerance,
|
||||
const std::string_view reason) {
|
||||
return position_text(quantity, position) +
|
||||
" component=n/a reference=n/a actual=n/a "
|
||||
"normalized_error=inf" + tolerance_text(tolerance) +
|
||||
" reason=" + std::string{reason};
|
||||
}
|
||||
|
||||
bool origin_matches(
|
||||
const EntityOrigin& origin,
|
||||
const std::string& instance_name,
|
||||
const std::int64_t local_label) {
|
||||
return origin.instance_name == instance_name &&
|
||||
origin.local_label == local_label;
|
||||
}
|
||||
|
||||
ComparisonSampleMatch matching_failure(
|
||||
const ReferenceQuantity quantity,
|
||||
const ResultPosition& position,
|
||||
const Tolerance tolerance,
|
||||
std::string code,
|
||||
const std::string_view reason) {
|
||||
std::vector<Diagnostic> failures;
|
||||
add_failure(
|
||||
failures,
|
||||
std::move(code),
|
||||
unevaluable_failure_text(
|
||||
quantity, position, tolerance, reason));
|
||||
return {std::nullopt, std::move(failures)};
|
||||
}
|
||||
|
||||
std::vector<double> as_vector(const std::array<double, 6>& values) {
|
||||
return {values.begin(), values.end()};
|
||||
}
|
||||
|
||||
} // namespace
|
||||
|
||||
ComparisonSampleMatch make_comparison_sample(
|
||||
const ResultFrame& frame,
|
||||
const ReferenceQuantity quantity,
|
||||
const ResultPosition& position,
|
||||
const std::span<const double> reference,
|
||||
const Tolerance tolerance) {
|
||||
std::vector<double> actual;
|
||||
|
||||
if (
|
||||
quantity == ReferenceQuantity::displacement ||
|
||||
quantity == ReferenceQuantity::reaction) {
|
||||
if (position.end_node_label.has_value()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.invalid_result_position",
|
||||
"nodal_position_has_end_node");
|
||||
}
|
||||
|
||||
std::optional<std::size_t> matched_index;
|
||||
for (
|
||||
std::size_t index = 0;
|
||||
index < frame.nodal.origins.size();
|
||||
++index) {
|
||||
if (origin_matches(
|
||||
frame.nodal.origins[index],
|
||||
position.instance_name,
|
||||
position.entity_label)) {
|
||||
if (matched_index.has_value()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.unknown_result_origin",
|
||||
"ambiguous_result_origin");
|
||||
}
|
||||
matched_index = index;
|
||||
}
|
||||
}
|
||||
if (!matched_index.has_value()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.unknown_result_origin",
|
||||
"unknown_result_origin");
|
||||
}
|
||||
|
||||
const auto& field =
|
||||
quantity == ReferenceQuantity::displacement
|
||||
? frame.nodal.displacement
|
||||
: frame.nodal.reaction;
|
||||
if (*matched_index >= field.size()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.component_count_mismatch",
|
||||
"missing_actual_components");
|
||||
}
|
||||
actual = as_vector(field[*matched_index]);
|
||||
} else {
|
||||
if (!position.end_node_label.has_value()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.invalid_element_node_pair",
|
||||
"missing_end_node");
|
||||
}
|
||||
|
||||
const BeamElementFrame* matched_beam = nullptr;
|
||||
for (const BeamElementFrame& beam : frame.element.beams) {
|
||||
if (origin_matches(
|
||||
beam.origin,
|
||||
position.instance_name,
|
||||
position.entity_label)) {
|
||||
if (matched_beam != nullptr) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.unknown_result_origin",
|
||||
"ambiguous_result_origin");
|
||||
}
|
||||
matched_beam = &beam;
|
||||
}
|
||||
}
|
||||
if (matched_beam == nullptr) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.unknown_result_origin",
|
||||
"unknown_result_origin");
|
||||
}
|
||||
|
||||
std::optional<NodeId> end_node;
|
||||
for (
|
||||
std::size_t index = 0;
|
||||
index < frame.nodal.origins.size() &&
|
||||
index < frame.nodal.node_ids.size();
|
||||
++index) {
|
||||
if (origin_matches(
|
||||
frame.nodal.origins[index],
|
||||
position.instance_name,
|
||||
*position.end_node_label)) {
|
||||
if (end_node.has_value()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.unknown_result_origin",
|
||||
"ambiguous_end_node_origin");
|
||||
}
|
||||
end_node = frame.nodal.node_ids[index];
|
||||
}
|
||||
}
|
||||
if (!end_node.has_value()) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.unknown_result_origin",
|
||||
"unknown_end_node_origin");
|
||||
}
|
||||
|
||||
const BeamSectionResult* matched_end = nullptr;
|
||||
for (const BeamSectionResult& end : matched_beam->end_results) {
|
||||
if (end.end_node == *end_node) {
|
||||
matched_end = &end;
|
||||
break;
|
||||
}
|
||||
}
|
||||
if (matched_end == nullptr) {
|
||||
return matching_failure(
|
||||
quantity,
|
||||
position,
|
||||
tolerance,
|
||||
"validation.invalid_element_node_pair",
|
||||
"node_is_not_element_end");
|
||||
}
|
||||
|
||||
if (quantity == ReferenceQuantity::internal_force) {
|
||||
actual = as_vector(matched_end->section_force);
|
||||
} else {
|
||||
actual = {matched_end->centroid_sigma_xx};
|
||||
}
|
||||
}
|
||||
|
||||
ComparisonSample matched{
|
||||
quantity,
|
||||
position,
|
||||
{reference.begin(), reference.end()},
|
||||
std::move(actual),
|
||||
tolerance,
|
||||
};
|
||||
return {std::move(matched), {}};
|
||||
}
|
||||
|
||||
ComparisonReport compare_samples(
|
||||
const std::span<const ComparisonSample> samples) {
|
||||
ComparisonReport report{true, 0.0, {}};
|
||||
std::set<PositionKey> positions;
|
||||
|
||||
for (const ComparisonSample& sample : samples) {
|
||||
const PositionKey key{
|
||||
sample.quantity,
|
||||
sample.position.instance_name,
|
||||
sample.position.entity_label,
|
||||
sample.position.end_node_label,
|
||||
};
|
||||
if (!positions.insert(key).second) {
|
||||
add_failure(
|
||||
report.failures,
|
||||
"validation.duplicate_result_position",
|
||||
unevaluable_failure_text(
|
||||
sample.quantity,
|
||||
sample.position,
|
||||
sample.tolerance,
|
||||
"duplicate_result_position"));
|
||||
report.maximum_normalized_error =
|
||||
std::numeric_limits<double>::infinity();
|
||||
continue;
|
||||
}
|
||||
|
||||
const std::size_t expected =
|
||||
expected_component_count(sample.quantity);
|
||||
if (
|
||||
sample.reference.size() != expected ||
|
||||
sample.actual.size() != expected) {
|
||||
add_failure(
|
||||
report.failures,
|
||||
"validation.component_count_mismatch",
|
||||
unevaluable_failure_text(
|
||||
sample.quantity,
|
||||
sample.position,
|
||||
sample.tolerance,
|
||||
"component_count_mismatch") +
|
||||
" expected=" + std::to_string(expected) +
|
||||
" reference_count=" +
|
||||
std::to_string(sample.reference.size()) +
|
||||
" actual_count=" +
|
||||
std::to_string(sample.actual.size()));
|
||||
report.maximum_normalized_error =
|
||||
std::numeric_limits<double>::infinity();
|
||||
continue;
|
||||
}
|
||||
|
||||
if (
|
||||
!std::isfinite(sample.tolerance.relative) ||
|
||||
!std::isfinite(sample.tolerance.absolute_scale) ||
|
||||
sample.tolerance.relative < 0.0 ||
|
||||
sample.tolerance.absolute_scale < 0.0) {
|
||||
add_failure(
|
||||
report.failures,
|
||||
"validation.invalid_tolerance",
|
||||
unevaluable_failure_text(
|
||||
sample.quantity,
|
||||
sample.position,
|
||||
sample.tolerance,
|
||||
"invalid_tolerance"));
|
||||
report.maximum_normalized_error =
|
||||
std::numeric_limits<double>::infinity();
|
||||
continue;
|
||||
}
|
||||
|
||||
for (std::size_t component = 0; component < expected; ++component) {
|
||||
const double reference = sample.reference[component];
|
||||
const double actual = sample.actual[component];
|
||||
if (!std::isfinite(reference) || !std::isfinite(actual)) {
|
||||
add_failure(
|
||||
report.failures,
|
||||
"validation.nonfinite_comparison_value",
|
||||
scalar_failure_text(
|
||||
sample,
|
||||
component,
|
||||
std::numeric_limits<double>::infinity()));
|
||||
report.maximum_normalized_error =
|
||||
std::numeric_limits<double>::infinity();
|
||||
continue;
|
||||
}
|
||||
|
||||
const double denominator =
|
||||
sample.tolerance.absolute_scale +
|
||||
sample.tolerance.relative * std::abs(reference);
|
||||
if (!(denominator > 0.0) || !std::isfinite(denominator)) {
|
||||
add_failure(
|
||||
report.failures,
|
||||
"validation.invalid_tolerance",
|
||||
scalar_failure_text(
|
||||
sample,
|
||||
component,
|
||||
std::numeric_limits<double>::infinity()));
|
||||
report.maximum_normalized_error =
|
||||
std::numeric_limits<double>::infinity();
|
||||
continue;
|
||||
}
|
||||
|
||||
const double normalized_error =
|
||||
std::abs(actual - reference) / denominator;
|
||||
report.maximum_normalized_error = std::max(
|
||||
report.maximum_normalized_error,
|
||||
normalized_error);
|
||||
if (!(normalized_error <= 1.0)) {
|
||||
add_failure(
|
||||
report.failures,
|
||||
"validation.tolerance_exceeded",
|
||||
scalar_failure_text(
|
||||
sample, component, normalized_error));
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
report.passed = report.failures.empty();
|
||||
return report;
|
||||
}
|
||||
|
||||
} // namespace fesa
|
||||
Reference in New Issue
Block a user