[ET Trac] #2944: CCE_Export Issues Ticket
Zach Etienne
trac-noreply at einsteintoolkit.org
Thu May 14 13:36:01 CDT 2026
#2944: CCE_Export Issues Ticket
Reporter: Zach Etienne
Status: new
Milestone:
Version:
Type: bug
Priority: major
Component:
Changes (by Zach Etienne):
@{5a73cd832d61371e861f3270}
The below review of CCE_Export was the product of a number of LLM worker agents, including Qwen3.6-27b, Claude Sonnet-4.7-XHigh, as well as a supervisor GPT-5.5-XHigh agent.
```
* Issue 1
- Severity high
- Description of the error
The interpolation helpers ignore too many error/status paths. The ADMBase
variable names used here are fixed by the thorn and inherited from ADMBase,
so "bad user-supplied variable name" is not the main risk. However,
CCTK_VarIndex, CCTK_InterpHandle, Util_TableCreate, the table setters,
CCTK_CoordSystemHandle, CCTK_InterpGridArrays, and Util_TableDestroy can all
report failure, and the output buffers are consumed immediately afterward.
This is especially important for out-of-domain extraction radii and missing
or misconfigured interpolation support.
- Filename + line number(s)
repos/CCE_Export/src/interpolate.cc:11, 33-60, 68, 85-112
repos/einsteinanalysis/Multipole/src/interpolate.cc:13-19, 69-79
repos/einsteinanalysis/ET_BHaHAHA/src/BHaHAHA_interpolate_metric_data.c:142-151
- Original code snippet
CCTK_INT variable_index = CCTK_VarIndex(name.c_str());
const int operator_handle =
CCTK_InterpHandle("Hermite polynomial interpolation");
int param_table_handle = Util_TableCreate(UTIL_TABLE_FLAGS_DEFAULT);
Util_TableSetFromString(param_table_handle, "order=3 ...");
Util_TableSetIntArray(param_table_handle, 4, operand_indices,
"operand_indices");
Util_TableSetIntArray(param_table_handle, 4, operation_codes,
"operation_codes");
const int coord_system_handle = CCTK_CoordSystemHandle("cart3d");
CCTK_InterpGridArrays(cctkGH, num_dims, operator_handle,
param_table_handle, coord_system_handle,
CCTK_MyProc(cctkGH) == 0 ? array_size : 0,
CCTK_VARIABLE_REAL, interp_coords, num_input_arrays,
input_array_indices, num_output_arrays,
output_array_types, output_arrays);
- Suggested patch
Add abort-level checks before any interpolated data are consumed:
CCTK_INT variable_index = CCTK_VarIndex(name.c_str());
if (variable_index < 0)
CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING,
"Could not find variable '%s'", name.c_str());
const int operator_handle =
CCTK_InterpHandle("Hermite polynomial interpolation");
if (operator_handle < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not find Hermite interpolator");
int param_table_handle = Util_TableCreate(UTIL_TABLE_FLAGS_DEFAULT);
if (param_table_handle < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not create interpolation table");
int ierr = Util_TableSetFromString(param_table_handle, "order=3 ...");
if (ierr < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not populate interpolation table");
ierr = Util_TableSetIntArray(param_table_handle, 4, operand_indices,
"operand_indices");
if (ierr < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not set operand_indices");
ierr = Util_TableSetIntArray(param_table_handle, 4, operation_codes,
"operation_codes");
if (ierr < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not set operation_codes");
const int coord_system_handle = CCTK_CoordSystemHandle("cart3d");
if (coord_system_handle < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not find cart3d coordinate system");
ierr = CCTK_InterpGridArrays(...);
if (ierr < 0)
CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING,
"Interpolation failed for '%s' with status %d",
name.c_str(), ierr);
ierr = Util_TableDestroy(param_table_handle);
if (ierr < 0)
CCTK_WARN(CCTK_WARN_ABORT, "Could not destroy interpolation table");
* Issue 2
- Severity medium
- Description of the error
The interpolation call must remain collective, but the spherical-harmonic
work after interpolation is only useful on rank 0. Cactus schedule option
GLOBAL means once per grid hierarchy/component scheduling level, not only on
MPI rank 0. The existing rank-0 guard around output confirms that all ranks
enter the routine. Nonzero ranks request zero interpolation points and then
perform decomposition work whose results are discarded.
- Filename + line number(s)
repos/CCE_Export/schedule.ccl:3-8
repos/CCE_Export/src/cce_export.cc:135-183
repos/flesh/doc/UsersGuide/Appendices.tex:1034-1085
- Original code snippet
Extract_Metric_Shift_Lapse_On_Sphere(...);
Compute_Ylms(th, ph, re_ylms, im_ylms, lmax, array_size);
for (int i = 0; i < 3; i++) {
for (int j = i; j < 3; j++) {
Decompose_Spherical_Harmonics(...);
}
}
...
if (CCTK_MyProc(cctkGH) == 0) {
Output_Decomposed_Metric_Data(...);
}
- Suggested patch
Keep the collective interpolation outside the rank guard, then perform mode
decomposition and output only on rank 0:
Extract_Metric_Shift_Lapse_On_Sphere(...);
if (CCTK_MyProc(cctkGH) == 0) {
Compute_Ylms(th, ph, re_ylms, im_ylms, lmax, array_size);
for (int i = 0; i < 3; i++) {
for (int j = i; j < 3; j++) {
Decompose_Spherical_Harmonics(...);
}
}
for (int i = 0; i < 3; i++) {
Decompose_Spherical_Harmonics(... beta ...);
Decompose_Spherical_Harmonics(... dr_beta ...);
Decompose_Spherical_Harmonics(... dt_beta ...);
}
Decompose_Spherical_Harmonics(... alpha ...);
Decompose_Spherical_Harmonics(... dr_alpha ...);
Decompose_Spherical_Harmonics(... dt_alpha ...);
Output_Decomposed_Metric_Data(...);
}
* Issue 3
- Severity high
- Description of the error
CCTK_REAL precision is configurable in Cactus/ET. The writer stores data in
vector<CCTK_REAL> but tells HDF5 that the memory buffer is H5T_NATIVE_DOUBLE.
This is only safe for CCTK_REAL_PRECISION_8. CarpetIOHDF5 maps CCTK_REAL to
H5T_NATIVE_FLOAT, H5T_NATIVE_DOUBLE, or H5T_NATIVE_LDOUBLE according to the
configured precision, which is the relevant ET convention. The utility reader
also needs a deliberate memory datatype if the on-disk datatype changes.
- Filename + line number(s)
repos/CCE_Export/src/h5_export.cc:65, 132
repos/CCE_Export/src/util/cce_export_hdf5_to_ascii.c:50, 76-77
repos/Carpet/CarpetIOHDF5/src/CarpetIOHDF5.hh:22-27
repos/Carpet/CarpetIOHDF5/src/CarpetIOHDF5.cc:249-302
- Original code snippet
dataset = H5Dcreate(file, datasetname.c_str(), H5T_NATIVE_DOUBLE, dataspace,
cparms);
...
HDF5_ERROR(H5Dwrite(dataset, H5T_NATIVE_DOUBLE, memdataspace, filespace,
H5P_DEFAULT, data));
- Suggested patch
Use a CCTK_REAL-aware memory datatype for writes. If SpECTRE compatibility
requires double-precision datasets on disk, keep the dataset datatype double
and let HDF5 convert from the correct memory datatype:
#if defined(CCTK_REAL_PRECISION_16)
#define CCE_HDF5_REAL_MEMORY H5T_NATIVE_LDOUBLE
#elif defined(CCTK_REAL_PRECISION_8)
#define CCE_HDF5_REAL_MEMORY H5T_NATIVE_DOUBLE
#elif defined(CCTK_REAL_PRECISION_4)
#define CCE_HDF5_REAL_MEMORY H5T_NATIVE_FLOAT
#else
#error "Unsupported CCTK_REAL precision"
#endif
const hid_t file_real_type = H5T_NATIVE_DOUBLE;
HDF5_ERROR(dataset = H5Dcreate(file, datasetname.c_str(), file_real_type,
dataspace, cparms));
HDF5_ERROR(H5Dwrite(dataset, CCE_HDF5_REAL_MEMORY, memdataspace, filespace,
H5P_DEFAULT, data));
In the utility reader, allocate the requested output type deliberately and
pass the matching HDF5 memory datatype, for example:
double *data = NULL;
status = H5Dread(dataset_id, H5T_NATIVE_DOUBLE, H5S_ALL, H5S_ALL,
H5P_DEFAULT, data);
* Issue 4
- Severity medium
- Description of the error
This thorn has C++14/C++17 dependencies without an explicit local contract.
Full ET context does not make this disappear: Cactus documentation still says
C++11 is required, and several simfactory optionlists still select C++11 or
C++14, even though newer generic and machine optionlists commonly use C++17.
The code uses std::make_unique, comma-separated using-declarations, and
std::filesystem. The filesystem feature-test check is also unreliable because
it tests __cpp_lib_filesystem before including a header that would define it.
- Filename + line number(s)
repos/CCE_Export/src/h5_export.cc:3-8, 13
repos/CCE_Export/src/cce_export.cc:12
repos/CCE_Export/src/cce_export.hh:15
repos/CCE_Export/src/spherical_harmonic_decomposition.cc:5, 81-82
repos/flesh/doc/UsersGuide/Notes.tex:47-49
repos/simfactory2/mdb/optionlists/generic.cfg:30-32
repos/simfactory2/mdb/optionlists/debian-cuda.cfg:44
repos/simfactory2/mdb/optionlists/wheeler.cfg:19
- Original code snippet
#if defined __cpp_lib_filesystem && __cpp_lib_filesystem < 201703L
#include <experimental/filesystem>
namespace fs = std::experimental::filesystem;
#else
#include <filesystem>
namespace fs = std::filesystem;
#endif
using std::string, std::ostringstream, std::map, std::ios, std::setprecision;
static auto fr = std::make_unique<CCTK_REAL []>(array_size);
- Suggested patch
Either make C++17 an explicit requirement, or remove the C++17-only code. A
minimal explicit-contract patch is:
#if __cplusplus < 201703L
#error "CCE_Export requires C++17 or later"
#endif
#include <filesystem>
namespace fs = std::filesystem;
using std::string;
using std::ostringstream;
using std::map;
using std::ios;
using std::setprecision;
Also document the C++17 requirement in README/ThornGuide and ensure the
build configuration used for this thorn selects a C++17-capable compiler
mode. If C++11 compatibility is desired instead, replace std::filesystem
with Cactus/local path handling and replace comma-separated using-declarations
with individual declarations.
* Issue 5
- Severity high
- Description of the error
For gravitational-wave extraction/CCE, r = 0 is not a valid extraction
radius or worldtube. The default configuration activates one radius at
radius[0] = 0.0, so a user who merely activates the thorn with defaults is
requesting mathematically invalid data. The implementation then compounds
this by computing radial derivatives through x^i/r partial_i, which divides
by radius[rad_index]. Unlike Multipole, CCE_Export explicitly exports radial
derivatives of worldtube data, so zero radius must be rejected, not treated
as an allowed lower endpoint.
- Filename + line number(s)
repos/CCE_Export/param.ccl:9-17
repos/CCE_Export/src/interpolate.cc:173-179, 197-202, 215-218
- Original code snippet
CCTK_INT nradii "How many radii to export the data on" STEERABLE=recover
{
0:100 :: "number of radii"
} 1
CCTK_REAL radius[101] "Radii to export the metric data on" STEERABLE=recover
{
0.0:* :: "Extraction radius"
} 0.0
dr_g.at(i).at(j).at(array_index) =
(xs.at(array_index) / radius[rad_index]) *
dx_g.at(i).at(j).at(array_index) +
(ys.at(array_index) / radius[rad_index]) *
dy_g.at(i).at(j).at(array_index) +
(zs.at(array_index) / radius[rad_index]) *
dz_g.at(i).at(j).at(array_index);
- Suggested patch
Make the default inactive or positive, and reject every active non-positive
extraction radius before interpolation:
CCTK_INT nradii "How many radii to export the data on" STEERABLE=recover
{
0:100 :: "number of radii"
} 0
CCTK_REAL radius[101] "Radii to export the metric data on" STEERABLE=recover
{
(0.0:* :: "Positive extraction radius"
} 1.0
for (int ir = 0; ir < nradii; ++ir) {
if (radius[ir] <= 0.0) {
CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING,
"CCE_Export::radius[%d] = %.17g is invalid: CCE/worldtube "
"extraction radii must be strictly positive",
ir, static_cast<double>(radius[ir]));
}
}
* Issue 6
- Severity medium
- Description of the error
Output is documented and implemented as one file per extraction radius, but
the filename uses a float argument and fixed two-decimal formatting. Distinct
CCTK_REAL radii can therefore collide in the same output filename. The
parameter documentation mentions two decimal places, but that does not
protect the data from accidental collisions and is inconsistent with "a
separate h5 file for each desired extraction radius."
- Filename + line number(s)
repos/CCE_Export/src/h5_export.hh:25-39
repos/CCE_Export/src/h5_export.cc:154, 168-171
repos/CCE_Export/src/cce_export.cc:178-182
repos/CCE_Export/param.ccl:25-29
repos/CCE_Export/doc/documentation.tex:49-50
- Original code snippet
void Output_Decomposed_Metric_Data(..., float rad, int lmax) {
...
basename << base_file_name << "R" << setiosflags(ios::fixed)
<< setprecision(2) << rad << "." << extension;
string output_name = (fs::path(my_out_dir) / basename.str()).string();
}
- Suggested patch
Preserve CCTK_REAL precision through the output path, increase filename
precision, and validate that active radii produce unique names:
void Output_Decomposed_Metric_Data(..., CCTK_REAL rad, int lmax) {
...
basename << base_file_name << "R" << setiosflags(ios::fixed)
<< setprecision(12) << static_cast<double>(rad)
<< "." << extension;
string output_name = (fs::path(my_out_dir) / basename.str()).string();
}
std::set<std::string> output_names;
for (int ir = 0; ir < nradii; ++ir) {
std::ostringstream name;
name << base_file_name << "R" << std::fixed << std::setprecision(12)
<< static_cast<double>(radius[ir]) << "." << extension;
if (!output_names.insert(name.str()).second) {
CCTK_VWarn(CCTK_WARN_ABORT, __LINE__, __FILE__, CCTK_THORNSTRING,
"Multiple CCE_Export radii map to output file '%s'",
name.str().c_str());
}
}
* Issue 7
- Severity medium
- Description of the error
The documentation is wrong in two places. The metric components are
interpolated directly from ADMBase; the code uses the ADM extrinsic-curvature
relation to compute partial_t gamma_ij, not to compute gamma_ij itself. The
spherical-harmonic coefficient integral implemented here is over dOmega using
sin(theta) dtheta dphi and does not multiply by r^2.
- Filename + line number(s)
repos/CCE_Export/doc/documentation.tex:76-84, 93-100
repos/CCE_Export/src/interpolate.cc:166-170, 221-311
repos/CCE_Export/src/spherical_harmonic_decomposition.cc:84-96
- Original code snippet
The metric is computed from the extrinsic curvature, lapse, and shift as:
\begin{equation}
K_{ij} = \frac{1}{2 \alpha} \left ( -\partial_0 \gamma_{ij}
+ \beta^k \partial_k \gamma_{ij}
+ \gamma_{ki} \partial_j \beta^k
+ \gamma_{kj} \partial_i \beta^k \right ) \, .
\end{equation}
C^{lm}(t, r) = \int {}_0 Y_{lm}^* u(t, r, \theta, \varphi) r^2 d \Omega
- Suggested patch
Update the ThornGuide to describe the implemented mathematics:
The metric components are interpolated directly from ADMBase. The metric
time derivative is computed by inverting the ADM definition of extrinsic
curvature:
\begin{equation}
\partial_t \gamma_{ij} =
-2\alpha K_{ij}
+ \beta^k \partial_k \gamma_{ij}
+ \gamma_{ki} \partial_j \beta^k
+ \gamma_{kj} \partial_i \beta^k .
\end{equation}
The spherical-harmonic coefficients are
\begin{eqnarray}
C^{lm}(t, r) =
\int {}_0 Y_{lm}^* u(t, r, \theta, \varphi) d\Omega .
\end{eqnarray}
* Issue 8
- Severity medium
- Description o
--
Ticket URL: https://bitbucket.org/einsteintoolkit/tickets/issues/2944/cce_export-issues-ticket
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.einsteintoolkit.org/pipermail/trac/attachments/20260514/5ddb803b/attachment-0001.htm>
More information about the Trac
mailing list