WIP: Add PowerLawParameterParameterDependence

.github/workflows/ci_bug.yml

@@ -0,0 +1,26 @@
+name: CI
+
+on:
+  push:
+    branches:
+      - master
+  pull_request:
+    branches:
+      - master
+
+jobs:
+  test:
+    runs-on: ubuntu-latest
+    strategy:
+      fail-fast: true
+    defaults:
+      run:
+        shell: bash -l {0}
+    env:
+      BUILD_DIR: ${{ github.workspace }}/build
+    steps:
+      - name: Test
+        run: |
+          echo "${{ github.workspace }}"
+          echo "${BUILD_DIR}"
+

src/libcadet/CMakeLists.txt

@@ -148,6 +148,7 @@ set(LIBCADET_SOURCES
 	${CMAKE_SOURCE_DIR}/src/libcadet/model/paramdep/ParameterDependenceBase.cpp
 	${CMAKE_SOURCE_DIR}/src/libcadet/model/paramdep/LiquidSaltSolidParameterDependence.cpp
 	${CMAKE_SOURCE_DIR}/src/libcadet/model/paramdep/DummyParameterDependence.cpp
+	${CMAKE_SOURCE_DIR}/src/libcadet/model/paramdep/PowerLawParameterDependence.cpp
 )
 
 # LIBCADET_NONLINALG_SOURCES holds all source files for LIBCADET_NONLINALG target

src/libcadet/ParameterDependenceFactory.cpp

@@ -22,6 +22,7 @@ namespace cadet
 			void registerLiquidSaltSolidParamDependence(std::unordered_map<std::string, std::function<model::IParameterStateDependence*()>>& paramDeps);
 			void registerDummyParamDependence(std::unordered_map<std::string, std::function<model::IParameterStateDependence*()>>& paramDeps);
 			void registerDummyParamDependence(std::unordered_map<std::string, std::function<model::IParameterParameterDependence*()>>& paramDeps);
+			void registerPowerLawParamDependence(std::unordered_map<std::string, std::function<model::IParameterParameterDependence*()>>& paramDeps);
 		}
 	}
 
@@ -33,6 +34,7 @@ namespace cadet
 
 		// Register all ParamParam dependencies
 		model::paramdep::registerDummyParamDependence(_paramParamDeps);
+		model::paramdep::registerPowerLawParamDependence(_paramParamDeps);
 	}
 
 	ParameterDependenceFactory::~ParameterDependenceFactory() { }

src/libcadet/model/paramdep/DummyParameterDependence.cpp

@@ -107,13 +107,13 @@ class DummyParameterParameterDependence : public ParameterParameterDependenceBas
 	}
 
 	template <typename ParamType>
-	ParamType getValue(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd) const
+	ParamType getValueImpl(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd) const
 	{
 		return 0.0;
 	}
 
 	template <typename ParamType>
-	ParamType getValue(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd, const ParamType& val) const
+	ParamType getValueImpl(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd, const ParamType& val) const
 	{
 		return 0.0;
 	}

src/libcadet/model/paramdep/LiquidSaltSolidParameterDependence.cpp

@@ -307,7 +307,7 @@ namespace paramdep
 		paramDeps[PowerLiquidSaltSolidParameterStateDependence::identifier()] = []() { return new PowerLiquidSaltSolidParameterStateDependence(); };
 		paramDeps[ColloidalAffinityLiquidSaltSolidParameterStateDependence::identifier()] = []() { return new ColloidalAffinityLiquidSaltSolidParameterStateDependence(); };
 	}
-}  // namespace reaction
+}  // namespace paramdep
 
 }  // namespace model
 

src/libcadet/model/paramdep/ParameterDependenceBase.hpp

@@ -241,25 +241,25 @@ class ParameterParameterDependenceBase : public IParameterParameterDependence
  *          
  *          The implementation is inserted inline in the class declaration.
  */
-#define CADET_PARAMETERPARAMETERDEPENDENCE_BOILERPLATE                                                                                                                                              \
+#define CADET_PARAMETERPARAMETERDEPENDENCE_BOILERPLATE                                                                                                                                          \
 	virtual double getValue(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd, double val) const        \
 	{                                                                                                                                                                                           \
-		return getValue<double>(unitOpIdx, params, colPos, comp, parType, bnd, val);                                                                                                            \
+		return getValueImpl<double>(unitOpIdx, params, colPos, comp, parType, bnd, val);                                                                                                        \
 	}                                                                                                                                                                                           \
                                                                                                                                                                                                 \
 	virtual active getValue(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd, const active& val) const \
 	{                                                                                                                                                                                           \
-		return getValue<active>(unitOpIdx, params, colPos, comp, parType, bnd, val);                                                                                                            \
+		return getValueImpl<active>(unitOpIdx, params, colPos, comp, parType, bnd, val);                                                                                                        \
 	}                                                                                                                                                                                           \
                                                                                                                                                                                                 \
 	virtual double getValue(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd) const                    \
 	{                                                                                                                                                                                           \
-		return getValue<double>(unitOpIdx, params, colPos, comp, parType, bnd);                                                                                                                 \
+		return getValueImpl<double>(unitOpIdx, params, colPos, comp, parType, bnd);                                                                                                             \
 	}                                                                                                                                                                                           \
                                                                                                                                                                                                 \
 	virtual active getValueActive(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd) const              \
 	{                                                                                                                                                                                           \
-		return getValue<active>(unitOpIdx, params, colPos, comp, parType, bnd);                                                                                                                 \
+		return getValueImpl<active>(unitOpIdx, params, colPos, comp, parType, bnd);                                                                                                             \
 	}
 
 

src/libcadet/model/paramdep/PowerLawParameterDependence.cpp

@@ -0,0 +1,87 @@
+// =============================================================================
+//  CADET
+//  
+//  Copyright © 2008-2022: The CADET Authors
+//            Please see the AUTHORS and CONTRIBUTORS file.
+//  
+//  All rights reserved. This program and the accompanying materials
+//  are made available under the terms of the GNU Public License v3.0 (or, at
+//  your option, any later version) which accompanies this distribution, and
+//  is available at http://www.gnu.org/licenses/gpl.html
+// =============================================================================
+
+#include "model/paramdep/ParameterDependenceBase.hpp"
+#include "SimulationTypes.hpp"
+#include "cadet/ParameterId.hpp"
+
+#include <cmath>
+
+namespace cadet
+{
+namespace model
+{
+
+/**
+ * @brief Defines a power law parameter parameter dependence
+ * @details The parameter dependence is defined as
+ * @f[\begin{align}
+	p_{new}(p) = \alpha p^k,
+   \end{align} @f]
+   where @f$ \alpha @f$ and @f$ k @f$ are (meta-)parameters.
+ */
+class PowerLawParameterParameterDependence : public ParameterParameterDependenceBase
+{
+public:
+
+	PowerLawParameterParameterDependence() { }
+	virtual ~PowerLawParameterParameterDependence() CADET_NOEXCEPT { }
+
+	static const char* identifier() { return "POWER_LAW"; }
+	virtual const char* name() const CADET_NOEXCEPT { return PowerLawParameterParameterDependence::identifier(); }
+
+	CADET_PARAMETERPARAMETERDEPENDENCE_BOILERPLATE
+
+protected:
+	active _base;
+	active _exponent;
+
+	virtual bool configureImpl(IParameterProvider& paramProvider, UnitOpIdx unitOpIdx, ParticleTypeIdx parTypeIdx, BoundStateIdx bndIdx, const std::string& name)
+	{
+		const std::string baseName = name + "_BASE";
+		const std::string expName = name + "_EXPONENT";
+		_base = paramProvider.getDouble(baseName);
+		_exponent = paramProvider.getDouble(expName);
+
+		_parameters[makeParamId(hashStringRuntime(baseName), unitOpIdx, CompIndep, parTypeIdx, bndIdx, ReactionIndep, SectionIndep)] = &_base;
+		_parameters[makeParamId(hashStringRuntime(expName), unitOpIdx, CompIndep, parTypeIdx, bndIdx, ReactionIndep, SectionIndep)] = &_exponent;
+
+		return true;
+	}
+
+	template <typename ParamType>
+	ParamType getValueImpl(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd) const
+	{
+		return 0.0;
+	}
+
+	template <typename ParamType>
+	ParamType getValueImpl(UnitOpIdx unitOpIdx, const std::unordered_map<ParameterId, active*>& params, const ColumnPosition& colPos, int comp, int parType, int bnd, ParamType val) const
+	{
+		using std::pow;
+		return static_cast<ParamType>(_base) * pow(val, static_cast<ParamType>(_exponent));
+	}
+
+};
+
+
+namespace paramdep
+{
+	void registerPowerLawParamDependence(std::unordered_map<std::string, std::function<model::IParameterParameterDependence*()>>& paramDeps)
+	{
+		paramDeps[PowerLawParameterParameterDependence::identifier()] = []() { return new PowerLawParameterParameterDependence(); };
+	}
+}  // namespace paramdep
+
+}  // namespace model
+
+}  // namespace cadet

src/libcadet/model/parts/ConvectionDispersionOperator.cpp

@@ -676,9 +676,11 @@ bool RadialConvectionDispersionOperatorBase::notifyDiscontinuousSectionTransitio
  */
 void RadialConvectionDispersionOperatorBase::setFlowRates(const active& in, const active& out, const active& colPorosity) CADET_NOEXCEPT
 {
+	const double pi = 3.1415926535897932384626434;
+
 	// If we have cross section area, interstitial velocity is given by network flow rates
 	if (_colLength > 0.0)
-		_curVelocity = _dir * in / (2.0 * M_PI * _colLength * colPorosity);
+		_curVelocity = _dir * in / (2.0 * pi * _colLength * colPorosity);
 }
 
 /**

test/CMakeLists.txt

@@ -31,22 +31,24 @@ endif()
 add_executable(testLogging testLogging.cpp)
 target_link_libraries(testLogging PRIVATE CADET::CompileOptions)
 
-add_executable(testRadialKernel testRadialKernel.cpp TimeIntegrator.cpp ${CMAKE_SOURCE_DIR}/src/libcadet/Logging.cpp ${CMAKE_SOURCE_DIR}/src/libcadet/AutoDiff.cpp)
-target_link_libraries(testRadialKernel PRIVATE CADET::CompileOptions CADET::AD SUNDIALS::sundials_idas ${SUNDIALS_NVEC_TARGET} )
-target_include_directories(testRadialKernel PRIVATE ${CMAKE_SOURCE_DIR}/test)
-
-if (ENABLE_GRM_2D)
-	target_include_directories(testRadialKernel PRIVATE ${CMAKE_BINARY_DIR}/src/libcadet)
-	if (SUPERLU_FOUND)
-		target_link_libraries(testRadialKernel PRIVATE SuperLU::SuperLU)
-	endif()
-	if (UMFPACK_FOUND)
-		target_link_libraries(testRadialKernel PRIVATE UMFPACK::UMFPACK)
+if (NOT WIN32)
+	add_executable(testRadialKernel testRadialKernel.cpp TimeIntegrator.cpp ${CMAKE_SOURCE_DIR}/src/libcadet/Logging.cpp ${CMAKE_SOURCE_DIR}/src/libcadet/AutoDiff.cpp)
+	target_link_libraries(testRadialKernel PRIVATE CADET::CompileOptions CADET::AD SUNDIALS::sundials_idas ${SUNDIALS_NVEC_TARGET} )
+	target_include_directories(testRadialKernel PRIVATE ${CMAKE_SOURCE_DIR}/test ${CMAKE_BINARY_DIR}/src/libcadet)
+
+	if (ENABLE_GRM_2D)
+		target_include_directories(testRadialKernel PRIVATE ${CMAKE_BINARY_DIR}/src/libcadet)
+		if (SUPERLU_FOUND)
+			target_link_libraries(testRadialKernel PRIVATE SuperLU::SuperLU)
+		endif()
+		if (UMFPACK_FOUND)
+			target_link_libraries(testRadialKernel PRIVATE UMFPACK::UMFPACK)
+		endif()
 	endif()
-endif()
 
-list(APPEND TEST_NONLINALG_TARGETS testRadialKernel)
-list(APPEND TEST_HDF5_TARGETS testRadialKernel)
+	list(APPEND TEST_NONLINALG_TARGETS testRadialKernel)
+	list(APPEND TEST_HDF5_TARGETS testRadialKernel)
+endif()
 
 
 # CATCH unit tests
@@ -74,8 +76,8 @@ add_executable(testRunner testRunner.cpp JsonTestModels.cpp ColumnTests.cpp Unit
 	$<TARGET_OBJECTS:libcadet_object>)
 
 target_link_libraries(testRunner PRIVATE CADET::CompileOptions CADET::AD SUNDIALS::sundials_idas ${SUNDIALS_NVEC_TARGET} ${TBB_TARGET})
+target_include_directories(testRunner PRIVATE ${CMAKE_BINARY_DIR}/src/libcadet)
 if (ENABLE_GRM_2D)
-	target_include_directories(testRunner PRIVATE ${CMAKE_BINARY_DIR}/src/libcadet)
 	if (SUPERLU_FOUND)
 		target_link_libraries(testRunner PRIVATE SuperLU::SuperLU)
 	endif()
@@ -112,9 +114,10 @@ endforeach()
 
 # ---------------------------------------------------
 
-
-add_executable(testCAPIv1 "${CMAKE_CURRENT_BINARY_DIR}/Paths_$<CONFIG>.cpp" testCAPIv1.cpp)
-target_include_directories(testCAPIv1 PRIVATE ${CMAKE_SOURCE_DIR}/include ${CMAKE_SOURCE_DIR}/ThirdParty/json)
+if (WIN32)
+	add_executable(testCAPIv1 "${CMAKE_CURRENT_BINARY_DIR}/Paths_$<CONFIG>.cpp" testCAPIv1.cpp)
+	target_include_directories(testCAPIv1 PRIVATE ${CMAKE_SOURCE_DIR}/include ${CMAKE_SOURCE_DIR}/ThirdParty/json)
+endif()
 
 # ---------------------------------------------------
 

test/testRadialKernel.cpp

@@ -30,7 +30,7 @@
 #define TEST_MANUFACTURED_TEXPT 1
 
 // Uncomment the next line to enable logging output of CADET in unit tests
-#define CADETTEST_ENABLE_LOG
+//#define CADETTEST_ENABLE_LOG
 
 #ifdef CADETTEST_ENABLE_LOG
 	#include "cadet/Logging.hpp"