Remove ReturnData::nx, ReturnData::nxtrue (#2966)

dweindl · web-flow · commit 3d46b6a6e472 · 2025-10-06T19:00:31.000+02:00
Replacements `ReturnData::nx_rdata` and `ReturnData::nxtrue_rdata` have been installed long ago. Let's get rid of the aliases.
diff --git a/include/amici/rdata.h b/include/amici/rdata.h
@@ -468,19 +468,6 @@ class ReturnData : public ModelDimensions {
      */
     int status = AMICI_NOT_RUN;
 
-    /**
-     * @brief Number of state variables.
-     *
-     * (alias `nx_rdata`, kept for backward compatibility)
-     */
-    int nx{0};
-
-    /**
-     * Number of state variables in the unaugmented system
-     * (alias nxtrue_rdata, kept for backward compatibility)
-     */
-    int nxtrue{0};
-
     /** Number of parameters w.r.t. which sensitivities were requested */
     int nplist{0};
 
diff --git a/include/amici/serialization.h b/include/amici/serialization.h
@@ -208,8 +208,8 @@ void serialize(
 ) {
     ar& dynamic_cast<amici::ModelDimensions&>(r);
     ar & r.id;
-    ar & r.nx;
-    ar & r.nxtrue;
+    ar & r.nx_rdata;
+    ar & r.nxtrue_rdata;
     ar & r.nplist;
     ar & r.nmaxevent;
     ar & r.nt;
diff --git a/models/model_calvetti_py/main.cpp b/models/model_calvetti_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_dirac_py/main.cpp b/models/model_dirac_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_events_py/main.cpp b/models/model_events_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_jakstat_adjoint_py/main.cpp b/models/model_jakstat_adjoint_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_nested_events_py/main.cpp b/models/model_nested_events_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_neuron_py/main.cpp b/models/model_neuron_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_robertson_py/main.cpp b/models/model_robertson_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/models/model_steadystate_py/main.cpp b/models/model_steadystate_py/main.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/python/sdist/amici/numpy.py b/python/sdist/amici/numpy.py
@@ -388,12 +388,12 @@ def __init__(self, rdata: ReturnDataPtr | ReturnData):
             )
         self._field_dimensions = {
             "ts": [rdata.nt],
-            "x": [rdata.nt, rdata.nx],
-            "x0": [rdata.nx],
-            "x_ss": [rdata.nx],
-            "sx": [rdata.nt, rdata.nplist, rdata.nx],
-            "sx0": [rdata.nplist, rdata.nx],
-            "sx_ss": [rdata.nplist, rdata.nx],
+            "x": [rdata.nt, rdata.nx_rdata],
+            "x0": [rdata.nx_rdata],
+            "x_ss": [rdata.nx_rdata],
+            "sx": [rdata.nt, rdata.nplist, rdata.nx_rdata],
+            "sx0": [rdata.nplist, rdata.nx_rdata],
+            "sx_ss": [rdata.nplist, rdata.nx_rdata],
             # observables
             "y": [rdata.nt, rdata.ny],
             "sigmay": [rdata.nt, rdata.ny],
diff --git a/src/hdf5.cpp b/src/hdf5.cpp
@@ -421,7 +421,7 @@ void writeReturnData(
 
     if (!rdata.x.empty())
         createAndWriteDouble2DDataset(
-            file, hdf5Location + "/x", rdata.x, rdata.nt, rdata.nx
+            file, hdf5Location + "/x", rdata.x, rdata.nt, rdata.nx_rdata
         );
 
     if (!rdata.y.empty())
@@ -463,13 +463,13 @@ void writeReturnData(
 
     if (!rdata.sx0.empty())
         createAndWriteDouble2DDataset(
-            file, hdf5Location + "/sx0", rdata.sx0, rdata.nplist, rdata.nx
+            file, hdf5Location + "/sx0", rdata.sx0, rdata.nplist, rdata.nx_rdata
         );
 
     if (!rdata.sx.empty())
         createAndWriteDouble3DDataset(
             file, hdf5Location + "/sx", rdata.sx, rdata.nt, rdata.nplist,
-            rdata.nx
+            rdata.nx_rdata
         );
 
     if (!rdata.sy.empty())
@@ -699,7 +699,7 @@ void writeReturnDataDiagnosis(
 
     if (!rdata.J.empty())
         createAndWriteDouble2DDataset(
-            file, hdf5Location + "/J", rdata.J, rdata.nx, rdata.nx
+            file, hdf5Location + "/J", rdata.J, rdata.nx_rdata, rdata.nx_rdata
         );
 
     if (!rdata.x_ss.empty())
diff --git a/src/main.template.cpp b/src/main.template.cpp
@@ -78,14 +78,14 @@ int main() {
 
     std::cout << "State sensitivities for timepoint " << rdata->ts[i_time]
               << std::endl; // nt x nplist x nx
-    for (int i_state = 0; i_state < rdata->nx; ++i_state) {
+    for (int i_state = 0; i_state < rdata->nx_rdata; ++i_state) {
         std::cout << "\td(" << state_ids[i_state] << ")/d("
                   << parameter_ids[model->plist(i_nplist)] << ") = ";
 
         // rdata->sx is a flat 3D array in row-major ordering
         std::cout << rdata->sx
-                         [i_time * rdata->nplist * rdata->nx
-                          + i_nplist * rdata->nx + i_state];
+                         [i_time * rdata->nplist * rdata->nx_rdata
+                          + i_nplist * rdata->nx_rdata + i_state];
         std::cout << std::endl;
     }
 
diff --git a/src/rdata.cpp b/src/rdata.cpp
@@ -37,8 +37,6 @@ ReturnData::ReturnData(
 )
     : ModelDimensions(model_dimensions_)
     , ts(std::move(ts_))
-    , nx(nx_rdata)
-    , nxtrue(nxtrue_rdata)
     , nplist(plist_.size())
     , nmaxevent(nmaxevent_)
     , nt(ts.size())
@@ -135,7 +133,7 @@ void ReturnData::initializeFullReporting(bool quadratic_llh) {
     sigmaz.resize(nmaxevent * nz, 0.0);
 
     rz.resize(nmaxevent * nz, 0.0);
-    x.resize(nt * nx, 0.0);
+    x.resize(nt * nx_rdata, 0.0);
     w.resize(nt * nw, 0.0);
 
     preeq_numsteps.resize(3, 0);
@@ -159,17 +157,17 @@ void ReturnData::initializeFullReporting(bool quadratic_llh) {
         }
     }
 
-    x0.resize(nx, getNaN());
-    x_ss.resize(nx, getNaN());
+    x0.resize(nx_rdata, getNaN());
+    x_ss.resize(nx_rdata, getNaN());
 
     if (sensi >= SensitivityOrder::first) {
-        sx0.resize(nx * nplist, getNaN());
-        sx_ss.resize(nx * nplist, getNaN());
+        sx0.resize(nx_rdata * nplist, getNaN());
+        sx_ss.resize(nx_rdata * nplist, getNaN());
 
         if (sensi_meth == SensitivityMethod::forward
             || sensi >= SensitivityOrder::second) {
             // for second order we can fill in from the augmented states
-            sx.resize(nt * nx * nplist, 0.0);
+            sx.resize(nt * nx_rdata * nplist, 0.0);
             sz.resize(nmaxevent * nz * nplist, 0.0);
             srz.resize(nmaxevent * nz * nplist, 0.0);
         }
@@ -339,7 +337,7 @@ void ReturnData::getDataOutput(
     int it, Model& model, SolutionState const& sol, ExpData const* edata
 ) {
     if (!x.empty()) {
-        model.fx_rdata(slice(x, it, nx), sol.x);
+        model.fx_rdata(slice(x, it, nx_rdata), sol.x);
     }
     if (!w.empty())
         model.getExpression(slice(w, it, nw), ts[it], sol.x);
@@ -378,7 +376,7 @@ void ReturnData::getDataSensisFSA(
 ) {
     if (!sx.empty()) {
         model.fsx_rdata(
-            gsl::span<realtype>(&sx.at(it * nplist * nx), nplist * nx), sol.sx,
+            gsl::span<realtype>(&sx.at(it * nplist * nx_rdata), nplist * nx_rdata), sol.sx,
             sol.x
         );
     }
@@ -652,14 +650,14 @@ void ReturnData::invalidate(int const it_start) {
     invalidateSLLH();
 
     if (!x.empty())
-        std::fill(x.begin() + nx * it_start, x.end(), getNaN());
+        std::fill(x.begin() + nx_rdata * it_start, x.end(), getNaN());
     if (!y.empty())
         std::fill(y.begin() + ny * it_start, y.end(), getNaN());
     if (!w.empty())
         std::fill(w.begin() + nw * it_start, w.end(), getNaN());
 
     if (!sx.empty())
-        std::fill(sx.begin() + nx * nplist * it_start, sx.end(), getNaN());
+        std::fill(sx.begin() + nx_rdata * nplist * it_start, sx.end(), getNaN());
     if (!sy.empty())
         std::fill(sy.begin() + ny * nplist * it_start, sy.end(), getNaN());
 }
@@ -729,10 +727,10 @@ void ReturnData::applyChainRuleFactorToSimulationResults(Model const& model) {
             && sensi_meth == SensitivityMethod::adjoint) {
             if (!sx.empty() && !x.empty())
                 for (int ip = 0; ip < nplist; ++ip)
-                    for (int ix = 0; ix < nxtrue; ++ix)
+                    for (int ix = 0; ix < nxtrue_rdata; ++ix)
                         for (int it = 0; it < nt; ++it)
-                            sx.at(ix + nxtrue * (ip + it * nplist))
-                                = x.at(it * nx + nxtrue + ip * nxtrue + ix);
+                            sx.at(ix + nxtrue_rdata * (ip + it * nplist))
+                                = x.at(it * nx_rdata + nxtrue_rdata + ip * nxtrue_rdata + ix);
 
             if (!sy.empty() && !y.empty())
                 for (int ip = 0; ip < nplist; ++ip)
@@ -789,14 +787,14 @@ void ReturnData::applyChainRuleFactorToSimulationResults(Model const& model) {
                             val *= pcoefficient[ip];
         };
 
-        chain_rule(sx, nxtrue, nx, nt);
+        chain_rule(sx, nxtrue_rdata, nx_rdata, nt);
         chain_rule(sy, nytrue, ny, nt);
         chain_rule(ssigmay, nytrue, ny, nt);
         chain_rule(sz, nztrue, nz, nmaxevent);
         chain_rule(ssigmaz, nztrue, nz, nmaxevent);
         chain_rule(srz, nztrue, nz, nmaxevent);
-        chain_rule(sx0, nxtrue, nx, 1);
-        chain_rule(sx_ss, nxtrue, nx, 1);
+        chain_rule(sx0, nxtrue_rdata, nx_rdata, 1);
+        chain_rule(sx_ss, nxtrue_rdata, nx_rdata, 1);
     }
 
     if (o2mode == SecondOrderMode::full) {
@@ -843,7 +841,7 @@ void ReturnData::applyChainRuleFactorToSimulationResults(Model const& model) {
                         }
         };
 
-        chain_rule(sx, nxtrue, nx, nt);
+        chain_rule(sx, nxtrue_rdata, nx_rdata, nt);
         chain_rule(sy, nytrue, ny, nt);
         chain_rule(ssigmay, nytrue, ny, nt);
         chain_rule(sz, nztrue, nz, nmaxevent);
@@ -883,7 +881,7 @@ void ReturnData::applyChainRuleFactorToSimulationResults(Model const& model) {
                     }
         };
 
-        chain_rule(sx, nxtrue, nx, nt);
+        chain_rule(sx, nxtrue_rdata, nx_rdata, nt);
         chain_rule(sy, nytrue, ny, nt);
         chain_rule(ssigmay, nytrue, ny, nt);
         chain_rule(sz, nztrue, nz, nmaxevent);
diff --git a/tests/cpp/testfunctions.cpp b/tests/cpp/testfunctions.cpp
@@ -307,7 +307,7 @@ void verifyReturnDataSensitivities(H5::H5File const& file, std::string const& re
                 expected = hdf5::getDoubleDataset3D(file, resultPath + "/sx", m, n, o);
                 for(int ip = 0; ip < model->nplist(); ++ip)
                     checkEqualArray(&expected[ip * model->nt() * model->nxtrue_rdata],
-                            &rdata->sx[ip * model->nt() * rdata->nx],
+                            &rdata->sx[ip * model->nt() * rdata->nx_rdata],
                             model->nt() * model->nxtrue_rdata, atol, rtol, "sx");
             } else {
                 ASSERT_TRUE(rdata->sx.empty());
diff --git a/tests/cpp/unittests/testSerialization.cpp b/tests/cpp/unittests/testSerialization.cpp
@@ -14,8 +14,8 @@ void checkReturnDataEqual(
     ASSERT_EQ(r.id, s.id);
     ASSERT_EQ(r.np, s.np);
     ASSERT_EQ(r.nk, s.nk);
-    ASSERT_EQ(r.nx, s.nx);
-    ASSERT_EQ(r.nxtrue, s.nxtrue);
+    ASSERT_EQ(r.nx_rdata, s.nx_rdata);
+    ASSERT_EQ(r.nxtrue_rdata, s.nxtrue_rdata);
     ASSERT_EQ(r.nx_solver, s.nx_solver);
     ASSERT_EQ(r.nx_solver_reinit, s.nx_solver_reinit);
     ASSERT_EQ(r.ny, s.ny);
