Skip to content

Commit 727d4df

Browse files
committed
Fix unit test registration and fitter order-1 bug
Register test_radiation_ss_model and test_radiation_ss_fitter in CMakeLists.txt (source files existed but were never wired into CTest). Fix fitter rejecting order-1 fits - the loop started at order=2 and FitFromSVD explicitly rejected O<2, so rank-1 kernels (single exponential) always returned invalid. Allow order=1 in both paths. Switch TwoModeExponential test from pointwise relative error to peak-normalized absolute error (relative error is a poor metric for decaying kernels - can blow up in the tail where |K| is small).
1 parent 5aeba2a commit 727d4df

4 files changed

Lines changed: 15 additions & 14 deletions

File tree

include/hydroc/radiation/radiation_types.h

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -92,7 +92,8 @@ struct StateSpaceOptions {
9292
/**
9393
* @brief Minimum R² (coefficient of determination) threshold for fit acceptance.
9494
*
95-
* The fitter adds modes until either R² >= r2_threshold or max_order is reached.
95+
* The fitter tries all orders up to max_order and keeps the best R².
96+
* Higher orders can add fast-decaying poles that improve transient accuracy.
9697
* Values closer to 1.0 require better fits but may need more modes.
9798
* Typical values: 0.90-0.99. Default: 0.95.
9899
*/

src/hydro/radiation/radiation_ss_fitter.cpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -121,7 +121,7 @@ StateSpaceFitResult RadiationStateSpaceFitter::FitKernel(
121121
// enough" R² prevents these fast poles from getting through, causing a
122122
// characteristic slight lag after the first response lobe.
123123
// =========================================================================
124-
for (int order = 2; order <= max_possible_order; ++order) {
124+
for (int order = 1; order <= max_possible_order; ++order) {
125125
StateSpaceFitResult candidate = FitFromSVD(K, dt, order, y, U, V, svh, hankel_size);
126126

127127
if (candidate.IsValid() && candidate.r2 > result.r2) {
@@ -150,7 +150,7 @@ StateSpaceFitResult RadiationStateSpaceFitter::FitFromSVD(
150150
const int N = static_cast<int>(K.size());
151151
const int O = order;
152152

153-
if (svh.size() < O || O < 2) {
153+
if (svh.size() < O || O < 1) {
154154
return result;
155155
}
156156

tests/unit/CMakeLists.txt

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -32,3 +32,5 @@ endfunction()
3232
# --- Register unit tests ---
3333
message(STATUS "Add HydroChrono unit tests")
3434
build_unit_test(test_added_mass_determinism)
35+
build_unit_test(test_radiation_ss_model)
36+
build_unit_test(test_radiation_ss_fitter)

tests/unit/test_radiation_ss_fitter.cpp

Lines changed: 9 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -165,22 +165,20 @@ bool TestTwoModeExponential() {
165165
TEST_ASSERT(result.order >= 2, "Should find at least 2 modes");
166166
TEST_ASSERT(result.r2 >= 0.98, "R² should be >= 0.98");
167167

168-
// Verify reconstruction
168+
// Verify reconstruction using peak-normalized absolute error.
169+
// Relative error is a poor metric for decaying kernels (blows up in the tail).
169170
Eigen::VectorXd K_fit = RadiationStateSpaceFitter::ReconstructKernel(result, dt, N);
170171

171-
double max_rel_error = 0.0;
172-
for (int k = 0; k < N; ++k) {
173-
if (std::abs(K(k)) > 0.01) {
174-
double rel_error = std::abs(K(k) - K_fit(k)) / std::abs(K(k));
175-
max_rel_error = std::max(max_rel_error, rel_error);
176-
}
177-
}
178-
179-
TEST_ASSERT(max_rel_error < 0.1, "Max relative error should be < 10%");
172+
const double K_peak = K.lpNorm<Eigen::Infinity>();
173+
double max_abs_error = (K - K_fit).lpNorm<Eigen::Infinity>();
174+
double normalized_error = max_abs_error / K_peak;
180175

181176
std::cout << " R² = " << std::fixed << std::setprecision(6) << result.r2 << std::endl;
182177
std::cout << " Order = " << result.order << std::endl;
183-
std::cout << " Max relative error = " << max_rel_error << std::endl;
178+
std::cout << " Max abs error = " << max_abs_error
179+
<< " (" << std::setprecision(2) << 100.0 * normalized_error << "% of peak)" << std::endl;
180+
181+
TEST_ASSERT(normalized_error < 0.05, "Peak-normalized error should be < 5%");
184182
std::cout << " PASSED" << std::endl;
185183
return true;
186184
}

0 commit comments

Comments
 (0)