Implement missing CAPI methods (#115)

* Add various missing CAPI methods

* Add write last sensitivity state for units

* Do not write out particle dimension for LRMP

* Add method to return CAPI version number

* Update documentation

Co-authored-by: Johannes Schmoelder <[email protected]>
Co-authored-by: Jan Breuer <[email protected]>

doc/interface/return_data.rst

@@ -16,15 +16,15 @@ Group /input/return
 
 ``WRITE_SOLUTION_LAST``
 
-   Write full solution state vector at last time point (optional, defaults to 0)
+   Write full solution state vector at last time point, including the state derivative vector (optional, defaults to 0)
 
    =============  ==========================
    **Type:** int  **Range:** :math:`\{0,1\}`
    =============  ==========================
 
 ``WRITE_SENS_LAST``
 
-   Write full sensitivity state vectors at last time point (optional, defaults to 0)
+   Write full sensitivity state vectors at last time point, including the state derivative vector (optional, defaults to 0)
 
    =============  ==========================
    **Type:** int  **Range:** :math:`\{0,1\}`
@@ -292,9 +292,17 @@ Group /input/return/unit_XXX
 
 ``WRITE_SOLUTION_LAST_UNIT``
 
-   Write solution state vector of this unit at last time point
+   Write solution state vector of this unit at last time point, including the state derivative vector (optional, defaults to 0)
 
    =============  ==========================
    **Type:** int  **Range:** :math:`\{0,1\}`
    =============  ==========================
 
+``WRITE_SENS_LAST_UNIT``
+
+   Write sensitivity state vector of this unit at last time point, including the state derivative vector (optional, defaults to 0)
+
+   =============  ==========================
+   **Type:** int  **Range:** :math:`\{0,1\}`
+   =============  ==========================
+

doc/interface/system.rst

@@ -32,23 +32,15 @@ Group /input/model
 
 ``INIT_STATE_SENSY_XXX``
 
-   Number of unit operations in the system
+   Initial full time derivative state vector of the :math:`\texttt{XXX}` th sensitivity system (optional, can currently not be specified on unit operation level)
 
    ================  ==================================
    **Type:** double  **Length:** :math:`\texttt{NDOF}`
    ================  ==================================
 
 ``INIT_STATE_SENSYDOT_XXX``
 
-   Initial full state vector of the :math:`\texttt{XXX}` th sensitivity system (optional, unit operation specific initial data is ignored)
-
-   ================  ==================================
-   **Type:** double  **Length:** :math:`\texttt{NDOF}`
-   ================  ==================================
-
-``NUNITS``
-
-   Initial full time derivative state vector of the :math:`\texttt{XXX}` th sensitivity system (optional, unit operation specific initial data is ignored)
+   Initial full state vector of the :math:`\texttt{XXX}` th sensitivity system (optional, can currently not be specified on unit operation level)
 
    ================  ==================================
    **Type:** double  **Length:** :math:`\texttt{NDOF}`

include/cadet/LibVersionInfo.hpp

@@ -81,6 +81,13 @@ namespace cadet
 	 * @return Git refspec
 	 */
 	CADET_API const char* getLibraryBuildHost() CADET_NOEXCEPT;
+
+	/**
+	 * @brief Returns the latest C-API version implemented by CADET
+	 * @sa cadetGetLatestCAPIVersion()
+	 * @return C-API Version number
+	 */
+	CADET_API const char* getLatestCAPIVersion() CADET_NOEXCEPT;
 } // namespace cadet
 
 #endif  // LIBCADET_LIBVERSIONINFO_HPP_

include/common/Driver.hpp

@@ -455,7 +455,6 @@ class Driver
 		else
 			_writeLastState = false;
 
-		// Adding a new output variable
 		std::ostringstream oss;
 		for (int i = 0; i <= _sim->model()->maxUnitOperationId(); ++i)
 		{
@@ -481,6 +480,25 @@ class Driver
 		else
 			_writeLastStateSens = false;
 
+		for (int i = 0; i <= _sim->model()->maxUnitOperationId(); ++i)
+		{
+			oss.str("");
+			oss << "unit_" << std::setfill('0') << std::setw(3) << std::setprecision(0) << i;
+			if (pp.exists(oss.str()))
+			{
+				pp.pushScope(oss.str());
+
+				if (pp.exists("WRITE_SENS_LAST_UNIT"))
+				{
+					const bool writeLastStateSensUnit = pp.getBool("WRITE_SENS_LAST_UNIT");
+					if (writeLastStateSensUnit)
+						_writeLastStateSensUnitId.push_back(i);
+				}
+
+				pp.popScope();
+			}
+		}
+
 		pp.popScope(); // scope return
 
 		if (applyInSimulator)
@@ -588,6 +606,40 @@ class Driver
 			}
 		}
 
+		for (int i = 0; i < _writeLastStateSensUnitId.size(); ++i)
+		{
+			unsigned int sliceStart;
+			unsigned int sliceEnd;
+			std::tie(sliceStart, sliceEnd) = _sim->model()->getModelStateOffsets(_writeLastStateSensUnitId[i]);
+
+			oss.str("");
+			oss << "unit_" << std::setfill('0') << std::setw(3) << std::setprecision(0) << _writeLastStateSensUnitId[i];
+			writer.pushGroup("sensitivity");
+			std::string unitID = oss.str();
+
+			unsigned int len = 0;
+
+			const std::vector<double const*> lastSens = _sim->getLastSensitivities(len);
+			const std::vector<double const*> lastSensdot = _sim->getLastSensitivityDerivatives(len);
+
+			const int nSens = lastSens.size();
+
+			for (int i = 0; i < nSens; i++)
+			{
+				oss.str("");
+				oss << std::setfill('0') << std::setw(3) << std::setprecision(0) << i;
+				writer.pushGroup("param_" + oss.str());
+				writer.pushGroup(unitID);
+				writer.template vector<double>("LAST_SENS_Y", sliceEnd - sliceStart, lastSens[i] + sliceStart);
+				writer.template vector<double>("LAST_SENS_YDOT", sliceEnd - sliceStart, lastSensdot[i] + sliceStart);
+				writer.popGroup();
+				writer.popGroup();
+			}
+
+			writer.popGroup();
+			writer.popGroup();
+		}
+
 		writer.popGroup();
 
 		if (writer.exists("meta"))
@@ -651,6 +703,27 @@ class Driver
 		}
 	}
 
+	inline void setWriteLastSensStateOfUnit(UnitOpIdx uid, bool writeLastSensStateUnit)
+	{
+		const std::vector<UnitOpIdx>::iterator it = std::find(_writeLastStateSensUnitId.begin(), _writeLastStateSensUnitId.end(), uid);
+		if (writeLastSensStateUnit)
+		{
+			if (it == _writeLastStateSensUnitId.end())
+			{
+				// Unit not in list, so add it
+				_writeLastStateSensUnitId.push_back(uid);
+			}
+		}
+		else
+		{
+			if (it != _writeLastStateSensUnitId.end())
+			{
+				// Unit is listed, so remove it
+				_writeLastStateSensUnitId.erase(it);
+			}
+		}
+	}
+
 	inline void setWriteLastState(bool writeLastState) CADET_NOEXCEPT { _writeLastState = writeLastState; }
 	inline void setWriteLastStateSens(bool writeLastState) CADET_NOEXCEPT { _writeLastStateSens = writeLastState; }
 	inline void setWriteSolutionTimes(bool solTimes) CADET_NOEXCEPT
@@ -670,6 +743,7 @@ class Driver
 	bool _writeLastState;
 	std::vector<UnitOpIdx> _writeLastStateUnitId;
 	bool _writeLastStateSens;
+	std::vector<UnitOpIdx> _writeLastStateSensUnitId;
 
 	/**
 	 * @brief Sets section times and section continuity from the given parameter provider

include/common/SolutionRecorderImpl.hpp

@@ -461,6 +461,11 @@ class InternalStorageUnitOpRecorder : public ISolutionRecorder
 	inline unsigned int numComponents() const CADET_NOEXCEPT { return _nComp; }
 	inline unsigned int numInletPorts() const CADET_NOEXCEPT { return _nInletPorts; }
 	inline unsigned int numOutletPorts() const CADET_NOEXCEPT { return _nOutletPorts; }
+	inline unsigned int numAxialCells() const CADET_NOEXCEPT { return _nAxialCells; }
+	inline unsigned int numRadialCells() const CADET_NOEXCEPT { return _nRadialCells; }
+	inline unsigned int numParticleTypes() const CADET_NOEXCEPT { return _nParShells.size(); }
+	inline unsigned int numParticleShells(unsigned int parType = 0) const CADET_NOEXCEPT { return _nParShells[parType]; }
+	inline unsigned int numBoundStates(unsigned int parType = 0) const CADET_NOEXCEPT { return _nBoundStates[parType]; }
 
 	inline double const* time() const CADET_NOEXCEPT { return _time.data(); }
 	inline double const* inlet() const CADET_NOEXCEPT { return _data.inlet.data(); }
@@ -491,6 +496,9 @@ class InternalStorageUnitOpRecorder : public ISolutionRecorder
 	inline double const* sensSolidDot(unsigned int idx, unsigned int parType = 0) const CADET_NOEXCEPT { return _sensDot[idx].solid[parType].data(); }
 	inline double const* sensFluxDot(unsigned int idx) const CADET_NOEXCEPT { return _sensDot[idx].flux.data(); }
 	inline double const* sensVolumeDot(unsigned int idx) const CADET_NOEXCEPT { return _sensDot[idx].volume.data(); }
+	inline double const* primaryCoordinates() const CADET_NOEXCEPT { return _axialCoords.data(); }
+	inline double const* secondaryCoordinates() const CADET_NOEXCEPT { return _radialCoords.data(); }
+	inline double const* particleCoordinates() const CADET_NOEXCEPT { return _particleCoords.data(); }
 protected:
 
 	struct Storage
@@ -1181,6 +1189,11 @@ class InternalStorageSystemRecorder : public ISolutionRecorder
 	inline InternalStorageUnitOpRecorder* recorder(unsigned int idx) CADET_NOEXCEPT { return _recorders[idx]; }
 	inline InternalStorageUnitOpRecorder* recorder(unsigned int idx) const CADET_NOEXCEPT { return _recorders[idx]; }
 
+	inline int nUnitOperations() CADET_NOEXCEPT
+	{
+		return _recorders.size();
+	}
+
 	inline InternalStorageUnitOpRecorder* unitOperation(UnitOpIdx idx) CADET_NOEXCEPT
 	{
 		for (InternalStorageUnitOpRecorder* rec : _recorders)
@@ -1208,6 +1221,7 @@ class InternalStorageSystemRecorder : public ISolutionRecorder
 	}
 
 	inline double const* time() const CADET_NOEXCEPT { return _time.data(); }
+	inline unsigned int numSensitivites() const CADET_NOEXCEPT { return _numSens; }
 
 protected:
 

src/libcadet/VersionInfo.cpp.in

@@ -23,6 +23,7 @@ namespace cadet
 	const char COMPILER[] = "@CMAKE_CXX_COMPILER_ID@ @CMAKE_CXX_COMPILER_VERSION@";
 	const char BUILD_HOST[] = "@CMAKE_SYSTEM@";
 	const char COMPILER_FLAGS[] = "@CMAKE_CXX_FLAGS@";
+	const char LATEST_CAPI_VERSION[] = "1.0.0";
 
 	const char* getLibraryVersion() CADET_NOEXCEPT
 	{
@@ -63,6 +64,11 @@ namespace cadet
 	{
 		return cadet::BUILD_HOST;
 	}
+
+	const char* getLatestCAPIVersion() CADET_NOEXCEPT
+	{
+		return cadet::LATEST_CAPI_VERSION;
+	}
 }
 
 extern "C"
@@ -107,4 +113,9 @@ extern "C"
 	{
 		return cadet::getLibraryBuildHost();
 	}
+
+	CADET_API const char* cdtGetLatestCAPIVersion()
+	{
+		return cadet::getLatestCAPIVersion();
+	}
 }