diff --git a/.github/pull_request_template.md b/.github/pull_request_template.md new file mode 100644 index 0000000..86619d6 --- /dev/null +++ b/.github/pull_request_template.md @@ -0,0 +1,7 @@ +## This PR updates + +(Add description) + +## PR Checklist +- [ ] Did you include comments/documentation? +- [ ] Did you include and pass unit tests? \ No newline at end of file diff --git a/CMakeLists.txt b/CMakeLists.txt index a3e8e72..7e24ed4 100755 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -154,3 +154,16 @@ if(TEST) target_link_libraries(simulatorTest PUBLIC gtest testLib) gtest_discover_tests(simulatorTest) endif() + +# create benchmarks +option(BENCHMARK "Configure for building benchmarks") +if(BENCHMARK) + find_package(benchmark REQUIRED) + set(TARGET_NAME simulatorBenchmark) + add_executable(simulatorBenchmark) + target_sources(simulatorBenchmark PUBLIC benchmarks/benchmark.cpp) + target_link_libraries(simulatorBenchmark PUBLIC gtest gtest_main) + target_link_libraries(simulatorBenchmark PUBLIC gtest lbmLib) + target_link_libraries(simulatorBenchmark PUBLIC gtest simLib) + target_link_libraries(simulatorBenchmark PUBLIC benchmark::benchmark) +endif() diff --git a/benchmarks/CMakeLists.txt b/benchmarks/CMakeLists.txt new file mode 100644 index 0000000..712438e --- /dev/null +++ b/benchmarks/CMakeLists.txt @@ -0,0 +1,5 @@ +target_sources(${TARGET_NAME} PRIVATE ${SOURCE_LIST}) +target_link_libraries(${TARGET_NAME} PUBLIC gtest gtest_main) +target_link_libraries(${TARGET_NAME} PUBLIC nlohmann_json::nlohmann_json) +target_link_libraries(${TARGET_NAME} PUBLIC lbmLib) +target_link_libraries(${TARGET_NAME} PUBLIC simLib) diff --git a/benchmarks/benchmark.cpp b/benchmarks/benchmark.cpp new file mode 100755 index 0000000..dd105ad --- /dev/null +++ b/benchmarks/benchmark.cpp @@ -0,0 +1,24 @@ +#include "gtest/gtest.h" +#include "benchmark/benchmark.h" + +#include "../src/baseSimulator.h" +#include "../src/baseSimulator.hh" + +using T = double; + +void BM_simRun(benchmark::State& state) { + + std::string file = "../examples/Hybrid/Network4a.JSON"; + + // Load and set the network from a JSON file + arch::Network network = porting::networkFromJSON(file); + + // Load and set the simulation from a JSON file + sim::Simulation testSimulation = porting::simulationFromJSON(file, &network); + for (auto _ : state) { + testSimulation.simulate(); + } +} +BENCHMARK(BM_simRun); + +BENCHMARK_MAIN(); \ No newline at end of file diff --git a/examples/Hybrid/Network1a.py b/examples/Hybrid/Network1a.py index 57cdd83..fe751dc 100644 --- a/examples/Hybrid/Network1a.py +++ b/examples/Hybrid/Network1a.py @@ -51,9 +51,11 @@ def hybridContinuous(): simulation.setContinuousPhase(f0) simulation.setPoiseuilleResistanceModel() - simulation.addLbmSimulator("1a", "../STL/cross.stl", m0, [n5, n7, n8, n9], \ - [[1.0, 0.0], [0.0, -1.0], [0.0, 1.0], [-1.0, 0.0]], [1e-4, 1e-4, 1e-4, 1e-4], \ - 1e-4, 1e-1, 0.1, 20, 1e-1, 0.55) + s1 = simulation.addLbmSimulator("1a", "../STL/cross.stl", m0, [n5, n7, n8, n9], \ + [[1.0, 0.0], [0.0, -1.0], [0.0, 1.0], [-1.0, 0.0]], [1e-4, 1e-4, 1e-4, 1e-4], \ + 1e-4, 1e-1, 20, 1e-1, 0.55) + + s1.setNaiveScheme(0.1, 0.5, 10) network.sort() diff --git a/python/mmft/simulator/bindings.cpp b/python/mmft/simulator/bindings.cpp index 4f01cec..dc43371 100644 --- a/python/mmft/simulator/bindings.cpp +++ b/python/mmft/simulator/bindings.cpp @@ -89,7 +89,6 @@ PYBIND11_MODULE(pysimulator, m) { std::vector widths, T charPhysLength, T charPhysVelocity, - T alpha, T resolution, T epsilon, T tau) { @@ -111,7 +110,7 @@ PYBIND11_MODULE(pysimulator, m) { } return simulation.addLbmSimulator( name, stlFile, simulation.getNetwork()->getModule(moduleId), openings, charPhysLength, - charPhysVelocity, alpha, resolution, epsilon, tau)->getId(); + charPhysVelocity, resolution, epsilon, tau)->getId(); }, "Add a LBM simulator to the simulation.") .def("injectDroplet", [](sim::Simulation &simulation, int dropletId, T injectionTime, int channelId, T injectionPosition) { @@ -126,6 +125,9 @@ PYBIND11_MODULE(pysimulator, m) { sim::ResistanceModelPoiseuille* resistanceModel = new sim::ResistanceModelPoiseuille(simulation.getContinuousPhase()->getViscosity()); simulation.setResistanceModel(resistanceModel); }) + .def("setNaiveScheme", [](sim::Simulation & simulation, T alpha, T beta, int theta) { + simulation.setNaiveHybridScheme(alpha, beta, theta); + }) .def("simulate", &sim::Simulation::simulate) .def("print", &sim::Simulation::printResults) .def("loadSimulation", [](sim::Simulation &simulation, arch::Network &network, std::string file) { diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 59a7056..31f5adc 100755 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -9,6 +9,7 @@ set(HEADER_LIST # add subdirectories add_subdirectory(nodalAnalysis) add_subdirectory(architecture) +add_subdirectory(hybridDynamics) add_subdirectory(simulation) add_subdirectory(olbProcessors) add_subdirectory(porting) diff --git a/src/baseSimulator.h b/src/baseSimulator.h index 91df00c..1446e65 100644 --- a/src/baseSimulator.h +++ b/src/baseSimulator.h @@ -26,6 +26,9 @@ #include "nodalAnalysis/NodalAnalysis.h" +#include "hybridDynamics/Scheme.h" +#include "hybridDynamics/Naive.h" + #include "architecture/Channel.h" #include "architecture/ChannelPosition.h" #include "architecture/Edge.h" diff --git a/src/baseSimulator.hh b/src/baseSimulator.hh index b2e5192..0420699 100644 --- a/src/baseSimulator.hh +++ b/src/baseSimulator.hh @@ -20,6 +20,9 @@ #include "nodalAnalysis/NodalAnalysis.hh" +#include "hybridDynamics/Scheme.hh" +#include "hybridDynamics/Naive.hh" + #include "architecture/Channel.hh" #include "architecture/ChannelPosition.hh" #include "architecture/Edge.hh" diff --git a/src/hybridDynamics/CMakeLists.txt b/src/hybridDynamics/CMakeLists.txt new file mode 100644 index 0000000..69e655e --- /dev/null +++ b/src/hybridDynamics/CMakeLists.txt @@ -0,0 +1,13 @@ +set(SOURCE_LIST + Scheme.hh + Naive.hh +) + +set(HEADER_LIST + Scheme.h + Naive.h +) + +target_sources(${TARGET_NAME} PUBLIC ${SOURCE_LIST} ${HEADER_LIST}) +target_include_directories(${TARGET_NAME} PUBLIC ${CMAKE_CURRENT_LIST_DIR}) +target_link_libraries(${TARGET_NAME} PUBLIC eigen) \ No newline at end of file diff --git a/src/hybridDynamics/Naive.h b/src/hybridDynamics/Naive.h new file mode 100644 index 0000000..be18f41 --- /dev/null +++ b/src/hybridDynamics/Naive.h @@ -0,0 +1,77 @@ +/** + * @file Naive.h + */ + +#pragma once + +#include +#include +#include + +namespace arch { + +// Forward declared dependencies +template +class Module; + +} + +namespace mmft{ + +/** + * @brief The Naive Scheme is an update scheme that sets the relaxation factor of an iterative update scheme to a + * given constant for all nodes. + */ +template +class NaiveScheme : public Scheme { + +public: + /** + * @brief Constructor of the Naive Scheme with provided constants. + * @param[in] modules The map of modules with boundary nodes upon which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + NaiveScheme(const std::unordered_map>>& modules, T alpha, T beta, int theta); + + /** + * @brief Constructor of the Naive Scheme with provided constants. + * @param[in] modules The map of modules with boundary nodes upon which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + NaiveScheme(const std::shared_ptr> module, T alpha, T beta, int theta); + + /** + * @brief Constructor of the Naive Scheme with provided constants. + * @param[in] modules The map of modules with boundary nodes upon which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + NaiveScheme(const std::shared_ptr> module, std::unordered_map alpha, std::unordered_map beta, int theta); + + /** + * @brief Constructor of the Naive Scheme with provided constants. + * @param[in] nodeIds The Ids of the set of nodes on which this scheme acts + * @param[in] moduleIds The Ids of the set of modules on which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + NaiveScheme(std::vector nodeIds, std::vector moduleIds, T alpha, T beta, int theta); + + /** + * @brief Constructor of the Naive Scheme with provided constants. The ids of the nodes and modules + * are passed in the maps. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + NaiveScheme(std::unordered_map alpha, std::unordered_map beta, std::unordered_map theta); + +}; + +} // namespace mmft \ No newline at end of file diff --git a/src/hybridDynamics/Naive.hh b/src/hybridDynamics/Naive.hh new file mode 100644 index 0000000..8a5a897 --- /dev/null +++ b/src/hybridDynamics/Naive.hh @@ -0,0 +1,36 @@ +#include "Naive.h" + +namespace mmft { + +template +NaiveScheme::NaiveScheme(const std::unordered_map>>& modules, T alpha, T beta, int theta) : + Scheme(modules, alpha, beta, theta) { } + +template +NaiveScheme::NaiveScheme(const std::shared_ptr> module, T alpha, T beta, int theta) : + Scheme(module, alpha, beta, theta) { } + +template +NaiveScheme::NaiveScheme(const std::shared_ptr> module, std::unordered_map alpha, std::unordered_map beta, int theta) : + Scheme(module, alpha, beta, theta) { } + +template +NaiveScheme::NaiveScheme(std::vector nodeIds, std::vector moduleIds, T alpha, T beta, int theta) : + Scheme() +{ + for (auto nodeId : nodeIds) { + this->alpha.try_emplace(nodeId, alpha); + } + for (auto nodeId : nodeIds) { + this->beta.try_emplace(nodeId, beta); + } + for (auto moduleId : moduleIds) { + this->alpha.try_emplace(moduleId, theta); + } +} + +template +NaiveScheme::NaiveScheme(std::unordered_map alpha_, std::unordered_map beta_, std::unordered_map theta_) : + Scheme(alpha_, beta_, theta_) { } + +} // namespace mmft \ No newline at end of file diff --git a/src/hybridDynamics/Scheme.h b/src/hybridDynamics/Scheme.h new file mode 100644 index 0000000..77f4811 --- /dev/null +++ b/src/hybridDynamics/Scheme.h @@ -0,0 +1,175 @@ +/** + * @file Scheme.h + */ + +#pragma once + +#include +#include +#include + +namespace arch { + +// Forward declared dependencies +template +class Module; + +} + +namespace mmft{ + +/** + * @brief "Virtual" definition of a general update scheme that functions as the base definition for other + * update schemes. An update scheme defines the update rules between Abstract and CFD for Hybrid simulation. + * The "virtuality" stems from the protected specifier for all constructors. + */ +template +class Scheme { + +protected: + + std::unordered_map alpha; // relaxation value for pressure value updates on the Abstract-CFD interface nodes + std::unordered_map beta; // relaxation value for pressure value flow rate on the Abstract-CFD interface nodes + std::unordered_map theta; // Amount of LBM collide and stream iterations between update steps + + /** + * @brief Default constructor of a Scheme. + */ + Scheme(); + + /** + * @brief Constructor of a Scheme with the required maps. The ids of the nodes and modules are passed in the maps. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + Scheme(std::unordered_map alpha, std::unordered_map beta, std::unordered_map theta); + + /** + * @brief Constructor of a Scheme with provided constants. The constants are set for all provided modules (and, hence, nodes). + * @param[in] modules The map of modules with boundary nodes upon which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + Scheme(const std::unordered_map>>& modules, T alpha, T beta, int theta); + + /** + * @brief Constructor of a Scheme with provided constants. The constants are set for all provided modules (and, hence, nodes). + * @param[in] modules The map of modules with boundary nodes upon which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + Scheme(const std::shared_ptr> module, T alpha, T beta, int theta); + + /** + * @brief Constructor of a Scheme with provided constants. The constants are set for all provided modules (and, hence, nodes). + * @param[in] modules The map of modules with boundary nodes upon which this scheme acts. + * @param[in] alpha The relaxation value for the pressure value update. + * @param[in] beta The relaxation value of the flow rate value update. + * @param[in] theta The amount of LBM stream and collide cycles between updates for a module. + */ + Scheme(const std::shared_ptr> module, std::unordered_map alpha, std::unordered_map beta, int theta); + +public: + + /** + * @brief Sets the relaxation value for pressure updates for all nodes to the provided value. + * @param[in] alpha The relaxation value for pressure updates. + */ + void setAlpha(T alpha); + + /** + * @brief Sets the relaxation value for pressure updates for node with nodeId to the provided value. + * @param[in] nodeId The id of the node that is to be updated. + * @param[in] alpha The relaxation value for pressure updates. + */ + void setAlpha(int nodeId, T alpha); + + /** + * @brief Sets the relaxation value for pressure updates for all nodes to the provided value. + * @param[in] alpha The map of relaxation values and node Ids. + */ + void setAlpha(std::unordered_map alpha); + + /** + * @brief Sets the relaxation value for flow rate updates for all nodes to the provided value. + * @param[in] beta The relaxation value for flow rate updates. + */ + void setBeta(T beta); + + /** + * @brief Sets the relaxation value for flow rate updates for node with nodeId to the provided value. + * @param[in] nodeId The id of the node that is to be updated. + * @param[in] beta The relaxation value for flow rate updates. + */ + void setBeta(int nodeId, T beta); + + /** + * @brief Sets the relaxation value for flow rate updates for all nodes to the provided value. + * @param[in] beta The map of relaxation values and node Ids. + */ + void setBeta(std::unordered_map beta); + + /** + * @brief Sets the number of LBM iterations between update steps for all modules to the provided value. + * @param[in] theta The number of LBM iterations between update steps. + */ + void setTheta(int theta); + + /** + * @brief Sets the number of LBM iterations between update steps for module with moduleId to the provided value. + * @param[in] module The id of the module that is to be updated. + * @param[in] theta The number of LBM iterations between update steps. + */ + void setTheta(int moduleId, int theta); + + /** + * @brief Sets the number of LBM iterations between update steps for all modules to the provided value. + * @param[in] beta The map of number of LBM iterations and moduleIds. + */ + void setTheta(std::unordered_map theta); + + /** + * @brief Returns the relaxation value for pressure updates for node with nodeId. + * @param[in] nodeId The node for which alpha is returned. + * @returns alpha. + */ + T getAlpha(int nodeId) const; + + /** + * @brief Returns the relaxation value for pressure updates for all nodes. + * @returns Map of alpha values. + */ + const std::unordered_map& getAlpha() const; + + /** + * @brief Returns the relaxation value for flow rate updates for node with nodeId. + * @param[in] nodeId The node for which beta is returned. + * @returns beta. + */ + T getBeta(int nodeId) const; + + /** + * @brief Returns the relaxation value for flow rate updates for all nodes. + * @returns Map of beta values. + */ + const std::unordered_map& getBeta() const; + + /** + * @brief Returns the number of LBM iterations between update steps for module with moduleId. + * @param[in] moduleId The module for which theta is returned. + * @returns theta. + */ + int getTheta(int moduleId) const; + + /** + * @brief Returns the number of LBM iterations between update steps for all modules. + * @returns Map of theta values. + */ + const std::unordered_map& getTheta() const; + +}; + +} // namespace mmft \ No newline at end of file diff --git a/src/hybridDynamics/Scheme.hh b/src/hybridDynamics/Scheme.hh new file mode 100644 index 0000000..9f34417 --- /dev/null +++ b/src/hybridDynamics/Scheme.hh @@ -0,0 +1,122 @@ +#include "Scheme.h" + +namespace mmft { + +template +Scheme::Scheme() { } + +template +Scheme::Scheme(std::unordered_map alpha_, std::unordered_map beta_, std::unordered_map theta_) : + alpha(alpha_), beta(beta_), theta(theta_) { } + +template +Scheme::Scheme(const std::unordered_map>>& modules_, T alpha_, T beta_, int theta_) +{ + for (auto& [key, module] : modules_) { + for (auto& [nodeId, node] : module->getNodes()) { + alpha.try_emplace(nodeId, alpha_); + beta.try_emplace(nodeId, beta_); + } + theta.try_emplace(module->getId(), theta_); + } +} + +template +Scheme::Scheme(const std::shared_ptr> module_, T alpha_, T beta_, int theta_) +{ + for (auto& [nodeId, node] : module_->getNodes()) { + alpha.try_emplace(nodeId, alpha_); + beta.try_emplace(nodeId, beta_); + } + theta.try_emplace(module_->getId(), theta_); +} + +template +Scheme::Scheme(const std::shared_ptr> module_, std::unordered_map alpha_, std::unordered_map beta_, int theta_) : + alpha(alpha_), beta(beta_) +{ + theta.try_emplace(module_->getId(), theta_); +} + +template +void Scheme::setAlpha(T alpha_) { + for (auto a : alpha) { + a = alpha_; + } +} + +template +void Scheme::setAlpha(int nodeId_, T alpha_) { + alpha.at(nodeId_) = alpha_; +} + +template +void Scheme::setAlpha(std::unordered_map alpha_) { + alpha = alpha_; +} + +template +void Scheme::setBeta(T beta_) { + for (auto b : beta) { + b = beta_; + } +} + +template +void Scheme::setBeta(int nodeId_, T beta_) { + beta.at(nodeId_) = beta_; +} + +template +void Scheme::setBeta(std::unordered_map beta_) { + beta = beta_; +} + +template +void Scheme::setTheta(int theta_) { + for (auto t : theta) { + t = theta_; + } +} + +template +void Scheme::setTheta(int moduleId_, int theta_) { + theta.at(moduleId_) = theta_; +} + +template +void Scheme::setTheta(std::unordered_map theta_) { + theta = theta_; +} + +template +T Scheme::getAlpha(int nodeId_) const { + return alpha.at(nodeId_); +} + +template +const std::unordered_map& Scheme::getAlpha() const { + return alpha; +} + +template +T Scheme::getBeta(int nodeId_) const { + return beta.at(nodeId_); +} + +template +const std::unordered_map& Scheme::getBeta() const { + return beta; +} + +template +int Scheme::getTheta(int moduleId_) const { + return theta.at(moduleId_); +} + +template +const std::unordered_map& Scheme::getTheta() const { + return theta; +} + +} // namespace mmft \ No newline at end of file diff --git a/src/nodalAnalysis/NodalAnalysis.hh b/src/nodalAnalysis/NodalAnalysis.hh index 5527c92..3eff6c9 100644 --- a/src/nodalAnalysis/NodalAnalysis.hh +++ b/src/nodalAnalysis/NodalAnalysis.hh @@ -371,7 +371,7 @@ void NodalAnalysis::writeCfdSimulators(std::unordered_mapgetPressure(); T set_pressure = 0.0; if (old_pressure > 0 ) { - set_pressure = old_pressure + cfdSimulator.second->getAlpha() * ( new_pressure - old_pressure ); + set_pressure = old_pressure + cfdSimulator.second->getAlpha(key) * ( new_pressure - old_pressure ); } else { set_pressure = new_pressure; } @@ -387,7 +387,7 @@ void NodalAnalysis::writeCfdSimulators(std::unordered_mapgetOpenings().at(key).width; T set_flowRate = 0.0; if (old_flowRate > 0 ) { - set_flowRate = old_flowRate + 5 * cfdSimulator.second->getAlpha() * ( new_flowRate - old_flowRate ); + set_flowRate = old_flowRate + 5 * cfdSimulator.second->getAlpha(key) * ( new_flowRate - old_flowRate ); } else { set_flowRate = new_flowRate; } @@ -432,4 +432,4 @@ void NodalAnalysis::printSystem() { std::cout << "Vector x:\n" << x << "\n\n" << std::endl; } -} // namespace nodal +} // namespace nodal \ No newline at end of file diff --git a/src/porting/jsonPorter.hh b/src/porting/jsonPorter.hh index 1f7dc25..96d3cc9 100644 --- a/src/porting/jsonPorter.hh +++ b/src/porting/jsonPorter.hh @@ -112,6 +112,7 @@ void simulationFromJSON(json jsonString, arch::Network* network_, sim::Simula if (platform == sim::Platform::Continuous) { readSimulators(jsonString, simulation, network_); + readUpdateScheme(jsonString, simulation); network_->sortGroups(); } else if (platform == sim::Platform::Mixing) { readMixingModel(jsonString, simulation); @@ -119,6 +120,7 @@ void simulationFromJSON(json jsonString, arch::Network* network_, sim::Simula readMixtures(jsonString, simulation); readMixtureInjections(jsonString, simulation, activeFixture); readSimulators(jsonString, simulation, network_); + readUpdateScheme(jsonString, simulation); network_->sortGroups(); } else if (platform == sim::Platform::Ooc) { readMixingModel(jsonString, simulation); @@ -127,6 +129,7 @@ void simulationFromJSON(json jsonString, arch::Network* network_, sim::Simula readMixtureInjections(jsonString, simulation, activeFixture); readTissues(jsonString, simulation); readSimulators(jsonString, simulation, network_); + readUpdateScheme(jsonString, simulation); network_->sortGroups(); } else if (platform == sim::Platform::BigDroplet) { throw std::invalid_argument("Droplet simulations are currently only supported for Abstract simulations."); diff --git a/src/porting/jsonReaders.h b/src/porting/jsonReaders.h index 4d00b14..7eebfec 100644 --- a/src/porting/jsonReaders.h +++ b/src/porting/jsonReaders.h @@ -141,12 +141,20 @@ void readContinuousPhase (json jsonString, sim::Simulation& simulation, int a /** * @brief Construct and stores the CFD modules simulators that are included in the network as defined by the json string * @param[in] jsonString json string + * @param[in] simulation simulation object * @param[in] network pointer to the network - * @param[in] activeFixture active fixture */ template void readSimulators (json jsonString, sim::Simulation& simulation, arch::Network* network); +/** + * @brief Construct and stores the update scheme that is used for the Abstract-CFD coupling. + * @param[in] jsonString json string + * @param[in] simulation simulation object +*/ +template +void readUpdateScheme (json jsonString, sim::Simulation& simulation); + /** * @brief Sets channels in the network to pressure or flow rate pump, as defined by the json string * @param[in] jsonString json string diff --git a/src/porting/jsonReaders.hh b/src/porting/jsonReaders.hh index d24e32d..6141b67 100644 --- a/src/porting/jsonReaders.hh +++ b/src/porting/jsonReaders.hh @@ -254,7 +254,6 @@ void readSimulators(json jsonString, sim::Simulation& simulation, arch::Netwo std::string stlFile = simulator["stlFile"]; T charPhysLength = simulator["charPhysLength"]; T charPhysVelocity = simulator["charPhysVelocity"]; - T alpha = simulator["alpha"]; T resolution = simulator["resolution"]; T epsilon = simulator["epsilon"]; T tau = simulator["tau"]; @@ -270,7 +269,7 @@ void readSimulators(json jsonString, sim::Simulation& simulation, arch::Netwo if(simulator["Type"] == "LBM") { auto simulator = simulation.addLbmSimulator(name, stlFile, network->getModule(moduleId), Openings, charPhysLength, - charPhysVelocity, alpha, resolution, epsilon, tau); + charPhysVelocity, resolution, epsilon, tau); simulator->setVtkFolder(vtkFolder); } else if (simulator["Type"] == "Mixing") @@ -280,7 +279,7 @@ void readSimulators(json jsonString, sim::Simulation& simulation, arch::Netwo species.try_emplace(specieId, speciePtr.get()); } auto simulator = simulation.addLbmMixingSimulator(name, stlFile, network->getModule(moduleId), species, - Openings, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + Openings, charPhysLength, charPhysVelocity, resolution, epsilon, tau); simulator->setVtkFolder(vtkFolder); } else if (simulator["Type"] == "Organ") @@ -292,14 +291,14 @@ void readSimulators(json jsonString, sim::Simulation& simulation, arch::Netwo std::string organStlFile = simulator["organStlFile"]; int tissueId = simulator["tissue"]; auto simulator = simulation.addLbmOocSimulator(name, stlFile, tissueId, organStlFile, network->getModule(moduleId), species, - Openings, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + Openings, charPhysLength, charPhysVelocity, resolution, epsilon, tau); simulator->setVtkFolder(vtkFolder); } else if(simulator["Type"] == "ESS_LBM") { #ifdef USE_ESSLBM auto simulator = simulation.addEssLbmSimulator(name, stlFile, network->getModule(moduleId), Openings, charPhysLength, - charPhysVelocity, alpha, resolution, epsilon, tau); + charPhysVelocity, resolution, epsilon, tau); simulator->setVtkFolder(vtkFolder); #else throw std::invalid_argument("The simulator was not build using the ESS library."); @@ -308,6 +307,58 @@ void readSimulators(json jsonString, sim::Simulation& simulation, arch::Netwo } } +template +void readUpdateScheme(json jsonString, sim::Simulation& simulation) { + + /* Legacy definition of update scheme for hybrid simulation + Will be deprecated in next release. */ + if (!jsonString["simulation"].contains("updateScheme")) { + T alpha = jsonString["simulation"]["settings"]["simulators"][0]["alpha"]; + simulation.setNaiveHybridScheme(alpha, 5*alpha, 10); + } + else { + if (jsonString["simulation"]["updateScheme"]["scheme"] == "Naive") { + if (jsonString["simulation"]["updateScheme"]["scheme"].contains("alpha") && + jsonString["simulation"]["updateScheme"]["scheme"].contains("beta") && + jsonString["simulation"]["updateScheme"]["scheme"].contains("theta")) + { + T alpha = jsonString["simulation"]["updateScheme"]["scheme"]["alpha"]; + T beta = jsonString["simulation"]["updateScheme"]["scheme"]["beta"]; + int theta = jsonString["simulation"]["updateScheme"]["scheme"]["theta"]; + simulation.setNaiveHybridScheme(alpha, beta, theta); + return; + } + else { + for (auto& simulator : jsonString["simulation"]["settings"]["simulators"]) { + if (simulator.contains("alpha") && simulator.contains("beta") && simulator.contains("theta")) { + int moduleCounter = 0; + if (simulator["alpha"].is_number() && simulator["beta"].is_number()) { + T alpha = simulator["alpha"]; + T beta = simulator["beta"]; + int theta = simulator["theta"]; + simulation.setNaiveHybridScheme(moduleCounter, alpha, beta, theta); + } else if (simulator["alpha"].is_array() && simulator["beta"].is_array()) { + int nodeCounter = 0; + std::unordered_map alpha; + std::unordered_map beta; + int theta = simulator["theta"]; + for (auto& opening : simulator["Openings"]) { + alpha.try_emplace(opening["node"], simulator["alpha"][nodeCounter]); + alpha.try_emplace(opening["node"], simulator["beta"][nodeCounter]); + nodeCounter++; + } + simulation.setNaiveHybridScheme(moduleCounter, alpha, beta, theta); + } + moduleCounter++; + } else { + throw std::invalid_argument("alpha, beta or theta values are either not or ill-defined for Naive update scheme."); + } + } + } + } + } +} + template void readBoundaryConditions(json jsonString, sim::Simulation& simulation, int activeFixture) { if (jsonString["simulation"]["fixtures"][activeFixture].contains("boundaryConditions")) { diff --git a/src/simulation/CFDSim.h b/src/simulation/CFDSim.h index 54922f1..722c0db 100644 --- a/src/simulation/CFDSim.h +++ b/src/simulation/CFDSim.h @@ -13,17 +13,26 @@ namespace sim { template class CFDSimulator; - /** - * @brief Conduct theta iterations of the CFD simulation on the network. - * @param[in] network The network on which the CFD simulations are conducted. - */ - template - bool conductCFDSimulation(const std::unordered_map>>& cfdSimulators, int iteration); - - template - void coupleNsAdLattices(const std::unordered_map>>& cfdSimulators); - - template - bool conductADSimulation(const std::unordered_map>>& cfdSimulators); +/** + * @brief Conduct the NS simulation step (of theta iterations) for all CFD simulators on the simulation. + * The amount of collide and stream iterations, theta, is obtained from the update scheme of the simulator. + * @returns A boolean for whether all simulators have converged, or not. + */ +template +bool conductCFDSimulation(const std::unordered_map>>& cfdSimulators); + +/** + * @brief Conduct the coupling step between the AD lattice and the NS lattice for all simulators. + */ +template +void coupleNsAdLattices(const std::unordered_map>>& cfdSimulators); + +/** + * @brief Conduct the AD simulation step (of theta iterations) for all CFD simulators on the simulation. + * The amount of collide and stream iterations, theta, is obtained from the update scheme of the simulator. + * @returns A boolean for whether all simulators have converged, or not. + */ +template +bool conductADSimulation(const std::unordered_map>>& cfdSimulators); } // namespace sim \ No newline at end of file diff --git a/src/simulation/CFDSim.hh b/src/simulation/CFDSim.hh index d855f32..79b546c 100644 --- a/src/simulation/CFDSim.hh +++ b/src/simulation/CFDSim.hh @@ -3,7 +3,7 @@ namespace sim { template - bool conductCFDSimulation(const std::unordered_map>>& cfdSimulators, int iteration) { + bool conductCFDSimulation(const std::unordered_map>>& cfdSimulators) { bool allConverge = true; diff --git a/src/simulation/Simulation.h b/src/simulation/Simulation.h index c4a2b62..2319a84 100644 --- a/src/simulation/Simulation.h +++ b/src/simulation/Simulation.h @@ -30,6 +30,17 @@ class Opening; } +namespace mmft { + +// Forward declared dependencies +template +class Scheme; + +template +class NaiveScheme; + +} + namespace nodal { // Forward declared dependencies @@ -120,7 +131,7 @@ class Simulation { Type simType = Type::Abstract; ///< The type of simulation that is being done. Platform platform = Platform::Continuous; ///< The microfluidic platform that is simulated in this simulation. arch::Network* network; ///< Network for which the simulation should be conducted. - std::shared_ptr> nodalAnalysis; + std::shared_ptr> nodalAnalysis; ///< The nodal analysis object, used to conduct abstract simulation. std::unordered_map>> fluids; ///< Fluids specified for the simulation. std::unordered_map>> droplets; ///< Droplets which are simulated in droplet simulation. std::unordered_map>> species; ///< Species specified for the simulation. @@ -128,10 +139,11 @@ class Simulation { std::unordered_map>> dropletInjections; ///< Injections of droplets that should take place during a droplet simulation. std::unordered_map>> mixtures; ///< Mixtures present in the simulation. std::unordered_map>> mixtureInjections; ///< Injections of fluids that should take place during the simulation. - std::unordered_map>> cfdSimulators; + std::unordered_map>> cfdSimulators; ///< The set of CFD simulators, that conduct CFD simulations on . ResistanceModel* resistanceModel; ///< The resistance model used for the simulation. MembraneModel* membraneModel; ///< The membrane model used for an OoC simulation. MixingModel* mixingModel; ///< The mixing model used for a mixing simulation. + std::unordered_map>> updateSchemes; ///< The update scheme for Abstract-CFD coupling int continuousPhase = 0; ///< Fluid of the continuous phase. int iteration = 0; int maxIterations = 1e5; @@ -277,6 +289,33 @@ class Simulation { */ Mixture* addDiffusiveMixture(std::unordered_map*> species, std::unordered_map, std::vector, T>> specieDistributions); + /** + * @brief Define and set the naive update scheme for a hybrid simulation. + * @param[in] alpha The relaxation value for the pressure value update for all nodes. + * @param[in] beta The relaxation value for the flow rate value update for all nodes. + * @param[in] theta The amount of LBM stream and collide cycles between updates for all modules. + * @returns A shared_ptr to the created naive update scheme. + */ + std::shared_ptr> setNaiveHybridScheme(T alpha, T beta, int theta); + + /** + * @brief Define and set the naive update scheme for a hybrid simulation. + * @param[in] alpha The relaxation value for the pressure value update for all nodes of the module. + * @param[in] beta The relaxation value for the flow rate value update for all nodes of the module. + * @param[in] theta The amount of LBM stream and collide cycles between updates for the modules. + * @returns A shared_ptr to the created naive update scheme. + */ + std::shared_ptr> setNaiveHybridScheme(int moduleId, T alpha, T beta, int theta); + + /** + * @brief Define and set the naive update scheme for a hybrid simulation. + * @param[in] alpha The relaxation value for the pressure value update for all nodes of the module. + * @param[in] beta The relaxation value for the flow rate value update for all nodes of the module. + * @param[in] theta The amount of LBM stream and collide cycles between updates for the modules. + * @returns A shared_ptr to the created naive update scheme. + */ + std::shared_ptr> setNaiveHybridScheme(int moduleId, std::unordered_map alpha, std::unordered_map beta, int theta); + /** * @brief Create injection. * @param[in] dropletId Id of the droplet that should be injected. @@ -297,44 +336,42 @@ class Simulation { MixtureInjection* addMixtureInjection(int mixtureId, int channelId, T injectionTime); /** - * @brief Adds a new module to the network. - * @param[in] name Name of the module. + * @brief Adds a new simulator to the network. + * @param[in] name Name of the simulator. * @param[in] stlFile Location of the stl file that gives the geometry of the domain. * @param[in] module Shared pointer to the module on which this solver acts. * @param[in] openings Map of openings corresponding to the nodes. * @param[in] charPhysLength Characteristic physical length of this simulator. * @param[in] charPhysVelocity Characteristic physical velocity of this simulator. - * @param[in] alpha Relaxation parameter for this simulator. * @param[in] resolution Resolution of this simulator. * @param[in] epsilon Error tolerance for convergence criterion of this simulator. * @param[in] tau Relaxation time of this simulator (0.5 < tau < 2.0). - * @return Pointer to the newly created module. + * @return Pointer to the newly created simulator. */ lbmSimulator* addLbmSimulator(std::string name, std::string stlFile, std::shared_ptr> module, std::unordered_map> openings, - T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau); + T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau); /** - * @brief Adds a new module to the network. - * @param[in] name Name of the module. + * @brief Adds a new simulator to the network. + * @param[in] name Name of the simulator. * @param[in] stlFile Location of the stl file that gives the geometry of the domain. * @param[in] module Shared pointer to the module on which this solver acts. * @param[in] species Map of specieIds and speciePtrs of the species simulated in the AD fields of this simulator. * @param[in] openings Map of openings corresponding to the nodes. * @param[in] charPhysLength Characteristic physical length of this simulator. * @param[in] charPhysVelocity Characteristic physical velocity of this simulator. - * @param[in] alpha Relaxation parameter for this simulator. * @param[in] resolution Resolution of this simulator. * @param[in] epsilon Error tolerance for convergence criterion of this simulator. * @param[in] tau Relaxation time of this simulator (0.5 < tau < 2.0). - * @return Pointer to the newly created module. + * @return Pointer to the newly created simulator. */ lbmMixingSimulator* addLbmMixingSimulator(std::string name, std::string stlFile, std::shared_ptr> module, std::unordered_map*> species, - std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau); + std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau); /** - * @brief Adds a new module to the network. - * @param[in] name Name of the module. + * @brief Adds a new simulator to the network. + * @param[in] name Name of the simulator. * @param[in] stlFile Location of the stl file that gives the geometry of the domain. * @param[in] tissueId The Id of the tissue that the organ in this nodule consists of. * @param[in] organStlFile The location of the stl file describing the geometry of the organ. @@ -343,23 +380,22 @@ class Simulation { * @param[in] openings Map of openings corresponding to the nodes. * @param[in] charPhysLength Characteristic physical length of this simulator. * @param[in] charPhysVelocity Characteristic physical velocity of this simulator. - * @param[in] alpha Relaxation parameter for this simulator. * @param[in] resolution Resolution of this simulator. * @param[in] epsilon Error tolerance for convergence criterion of this simulator. * @param[in] tau Relaxation time of this simulator (0.5 < tau < 2.0). - * @return Pointer to the newly created module. + * @return Pointer to the newly created simulator. */ lbmOocSimulator* addLbmOocSimulator(std::string name, std::string stlFile, int tissueId, std::string organStlFile, std::shared_ptr> module, std::unordered_map*> species, - std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau); + std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau); /** - * @brief Adds a new module to the network. - * @param[in] name Name of the module. + * @brief Adds a new simulator to the network. + * @param[in] name Name of the simulator. * @param[in] module Shared pointer to the module on which this solver acts. * @param[in] openings Map of openings corresponding to the nodes. */ essLbmSimulator* addEssLbmSimulator(std::string name, std::string stlFile, std::shared_ptr> module, std::unordered_map> openings, - T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau); + T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau); /** * @brief Set the platform of the simulation. @@ -423,8 +459,7 @@ class Simulation { /** * @brief Calculate and set new state of the continuous fluid simulation. Move mixture positions and create new mixtures if necessary. - * - * @param timeStep Time step in s for which the new mixtures state should be calculated. + * @param[in] timeStep Time step in s for which the new mixtures state should be calculated. */ void calculateNewMixtures(double timeStep); @@ -575,12 +610,13 @@ class Simulation { Fluid* mixFluids(int fluid0Id, T volume0, int fluid1Id, T volume1); /** - * @brief Creates a new droplet from two existing droplets. - Please note that this only creates a new droplet inside the simulation, but the actual boundaries have to be set separately, which is usually done inside the corresponding merge events. - * @param droplet0Id Id of the first droplet. - * @param droplet1Id Id of the second droplet. - * @return Pointer to new droplet. - */ + * @brief Creates a new droplet from two existing droplets. Please note that this only creates a new + * droplet inside the simulation, but the actual boundaries have to be set separately, which is usually + * done inside the corresponding merge events. + * @param droplet0Id Id of the first droplet. + * @param droplet1Id Id of the second droplet. + * @return Pointer to new droplet. + */ Droplet* mergeDroplets(int droplet0Id, int droplet1Id); /** diff --git a/src/simulation/Simulation.hh b/src/simulation/Simulation.hh index 8fb8553..b22e436 100644 --- a/src/simulation/Simulation.hh +++ b/src/simulation/Simulation.hh @@ -198,14 +198,40 @@ namespace sim { return nullptr; } + template + std::shared_ptr> Simulation::setNaiveHybridScheme(T alpha, T beta, int theta) { + auto naiveScheme = std::make_shared>(network->getModules(), alpha, beta, theta); + for (auto& [key, simulator] : cfdSimulators) { + updateSchemes.try_emplace(simulator->getId(), naiveScheme); + simulator->setUpdateScheme(updateSchemes.at(simulator->getId())); + } + return naiveScheme; + } + + template + std::shared_ptr> Simulation::setNaiveHybridScheme(int moduleId, T alpha, T beta, int theta) { + auto naiveScheme = std::make_shared>(network->getModule(moduleId), alpha, beta, theta); + updateSchemes.try_emplace(moduleId, naiveScheme); + cfdSimulators.at(moduleId)->setUpdateScheme(updateSchemes.at(moduleId)); + return naiveScheme; + } + + template + std::shared_ptr> Simulation::setNaiveHybridScheme(int moduleId, std::unordered_map alpha, std::unordered_map beta, int theta) { + auto naiveScheme = std::make_shared>(network->getModule(moduleId), alpha, beta, theta); + updateSchemes.try_emplace(moduleId, naiveScheme); + cfdSimulators.at(moduleId)->setUpdateScheme(updateSchemes.at(moduleId)); + return naiveScheme; + } + template lbmSimulator* Simulation::addLbmSimulator(std::string name, std::string stlFile, std::shared_ptr> module, std::unordered_map> openings, - T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau) + T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau) { if (resistanceModel != nullptr) { // create Simulator auto id = cfdSimulators.size(); - auto addCfdSimulator = new lbmSimulator(id, name, stlFile, module, openings, resistanceModel, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + auto addCfdSimulator = new lbmSimulator(id, name, stlFile, module, openings, resistanceModel, charPhysLength, charPhysVelocity, resolution, epsilon, tau); // add Simulator cfdSimulators.try_emplace(id, addCfdSimulator); @@ -218,13 +244,13 @@ namespace sim { template lbmMixingSimulator* Simulation::addLbmMixingSimulator(std::string name, std::string stlFile, std::shared_ptr> module, std::unordered_map*> species, - std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau) + std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau) { std::cout << "Trying to add a mixing simulator" << std::endl; if (resistanceModel != nullptr) { // create Simulator auto id = cfdSimulators.size(); - auto addCfdSimulator = new lbmMixingSimulator(id, name, stlFile, module, species, openings, resistanceModel, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + auto addCfdSimulator = new lbmMixingSimulator(id, name, stlFile, module, species, openings, resistanceModel, charPhysLength, charPhysVelocity, resolution, epsilon, tau); // add Simulator cfdSimulators.try_emplace(id, addCfdSimulator); @@ -237,12 +263,12 @@ namespace sim { template lbmOocSimulator* Simulation::addLbmOocSimulator(std::string name, std::string stlFile, int tissueId, std::string organStlFile, std::shared_ptr> module, std::unordered_map*> species, - std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau) + std::unordered_map> openings, T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau) { if (resistanceModel != nullptr) { // create Simulator auto id = cfdSimulators.size(); - auto addCfdSimulator = new lbmOocSimulator(id, name, stlFile, tissues.at(tissueId), organStlFile, module, species, openings, resistanceModel, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + auto addCfdSimulator = new lbmOocSimulator(id, name, stlFile, tissues.at(tissueId), organStlFile, module, species, openings, resistanceModel, charPhysLength, charPhysVelocity, resolution, epsilon, tau); // add Simulator cfdSimulators.try_emplace(id, addCfdSimulator); @@ -255,13 +281,13 @@ namespace sim { template essLbmSimulator* Simulation::addEssLbmSimulator(std::string name, std::string stlFile, std::shared_ptr> module, std::unordered_map> openings, - T charPhysLength, T charPhysVelocity, T alpha, T resolution, T epsilon, T tau) + T charPhysLength, T charPhysVelocity, T resolution, T epsilon, T tau) { #ifdef USE_ESSLBM if (resistanceModel != nullptr) { // create Simulator auto id = cfdSimulators.size(); - auto addCfdSimulator = new essLbmSimulator(id, name, stlFile, module, openings, resistanceModel, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + auto addCfdSimulator = new essLbmSimulator(id, name, stlFile, module, openings, resistanceModel, charPhysLength, charPhysVelocity, resolution, epsilon, tau); // add Simulator cfdSimulators.try_emplace(id, addCfdSimulator); @@ -556,12 +582,12 @@ namespace sim { // Initialization of CFD domains while (! allConverged) { - allConverged = conductCFDSimulation(cfdSimulators, 1); + allConverged = conductCFDSimulation(cfdSimulators); } while (! allConverged || !pressureConverged) { // conduct CFD simulations - allConverged = conductCFDSimulation(cfdSimulators, 10); + allConverged = conductCFDSimulation(cfdSimulators); // compute nodal analysis again pressureConverged = nodalAnalysis->conductNodalAnalysis(cfdSimulators); @@ -590,13 +616,13 @@ namespace sim { // Initialization of NS CFD domains while (! allConverged) { - allConverged = conductCFDSimulation(cfdSimulators, 1); + allConverged = conductCFDSimulation(cfdSimulators); } // Obtain overal steady-state flow result while (! allConverged || !pressureConverged) { // conduct CFD simulations - allConverged = conductCFDSimulation(cfdSimulators, 10); + allConverged = conductCFDSimulation(cfdSimulators); // compute nodal analysis again pressureConverged = nodalAnalysis->conductNodalAnalysis(cfdSimulators); } @@ -631,12 +657,12 @@ namespace sim { // Initialization of CFD domains while (! allConverged) { - allConverged = conductCFDSimulation(cfdSimulators, 1); + allConverged = conductCFDSimulation(cfdSimulators); } while (! allConverged || !pressureConverged) { // conduct CFD simulations - allConverged = conductCFDSimulation(cfdSimulators, 10); + allConverged = conductCFDSimulation(cfdSimulators); // compute nodal analysis again pressureConverged = nodalAnalysis->conductNodalAnalysis(cfdSimulators); } diff --git a/src/simulation/simulators/cfdSimulator.h b/src/simulation/simulators/cfdSimulator.h index 9c5b885..fadb0ce 100755 --- a/src/simulation/simulators/cfdSimulator.h +++ b/src/simulation/simulators/cfdSimulator.h @@ -22,6 +22,13 @@ struct Opening; } +namespace mmft { + +template +class Scheme; + +} + namespace sim { /** @@ -30,30 +37,54 @@ namespace sim { template class CFDSimulator { protected: - int const id; ///< Id of the simulator. - - std::string name; ///< Name of the module. + int const id; ///< Id of the simulator. + std::string name; ///< Name of the simulator. std::string vtkFolder = "./tmp/"; ///< Folder in which vtk files will be saved. std::string vtkFile = "."; ///< File in which last file was saved. - bool initialized = false; ///< Is the module initialized? + bool initialized = false; ///< Is the simulator initialized. std::string stlFile; ///< The STL file of the CFD domain. - std::shared_ptr> cfdModule; - std::shared_ptr> moduleNetwork; ///< Fully connected graph as network for the initial approximation. - std::unordered_map> moduleOpenings; ///< Map of openings. - std::unordered_map groundNodes; ///< Map of nodes that communicate the pressure to the 1D solver. + std::shared_ptr> cfdModule; ///< A pointer to the module, upon which this simulator acts. + std::shared_ptr> moduleNetwork; ///< Fully connected graph as network for the initial approximation. + std::unordered_map> moduleOpenings; ///< Map of openings. + std::unordered_map groundNodes; ///< Map of nodes that communicate the pressure to the 1D solver. + std::shared_ptr> updateScheme; ///< The update scheme for Abstract-CFD coupling + + /** + * @brief Constructor of a CFDSimulator, which acts as a base definition for other simulators. + * @param[in] id Id of the simulator. + * @param[in] name The name of the simulator. + * @param[in] stlFile Location of the stl file that gives the geometry of the domain. + * @param[in] cfdModule Shared pointer to the module on which this solver acts. + * @param[in] openings Map of openings corresponding to the nodes. + * @param[in] ResistanceModel The resistance model used for the simulation, necessary to set the initial condition. + */ + CFDSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, + std::unordered_map> openings, ResistanceModel* ResistanceModel); + + /** + * @brief Constructor of a CFDSimulator, which acts as a base definition for other simulators. + * @param[in] id Id of the simulator. + * @param[in] name The name of the simulator. + * @param[in] stlFile Location of the stl file that gives the geometry of the domain. + * @param[in] cfdModule Shared pointer to the module on which this solver acts. + * @param[in] openings Map of openings corresponding to the nodes. + * @param[in] updateScheme Shared pointer to the update scheme for Abstract-CFD coupling. + * @param[in] ResistanceModel The resistance model used for the simulation, necessary to set the initial condition. + */ + CFDSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, + std::unordered_map> openings, std::shared_ptr> updateScheme, ResistanceModel* ResistanceModel); - T alpha; ///< Relaxation factor for convergence between 1D and CFD simulation. /** * @brief Define and prepare the coupling of the NS lattice with the AD lattices. */ - virtual void executeCoupling() { }; + virtual void executeCoupling() + { + throw std::runtime_error("The function executeCoupling is undefined for this CFD simulator."); + }; public: - - CFDSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, std::unordered_map> openings, T alpha, ResistanceModel* ResistanceModel); - /** * @brief Get id of the simulator. * @returns id. @@ -61,7 +92,7 @@ class CFDSimulator { int getId() const; /** - * @brief Get the fully connected graph of this module, that is used for the initial approximation. + * @brief Get a shared ptr to the module, upon which this simulator acts. * @return Network of the fully connected graph. */ std::shared_ptr> getModule() const; @@ -102,20 +133,45 @@ class CFDSimulator { */ void setInitialized(bool initialization); - void setVtkFolder(std::string vtkFolder_); + /** + * @brief Set the update scheme for Abstract-CFD coupling for this simulator. + */ + void setUpdateScheme(const std::shared_ptr>& updateScheme); + + /** + * @brief Set the path, where vtk output from the simulator should be stored. + * @param[in] vtkFolder A string containing the path to the vtk folder. + */ + void setVtkFolder(std::string vtkFolder); + /** + * @brief Get the location of the last created vtk file. + * @returns vtkFile. + */ std::string getVtkFile(); /** - * @brief Get the relaxation factor alpha. + * @brief Get the relaxation factor for pressure update values, alpha. * @returns alpha. */ - T getAlpha(); + T getAlpha(int nodeId); - // fully virtual functions + /** + * @brief Get the relaxation factor for flow rate update values, beta. + * @returns beta. + */ + T getBeta(int nodeId); + /** + * @brief Initialize an instance of the LBM solver for this simulator. + * @param[in] dynViscosity Dynamic viscosity of the simulated fluid in _kg / m s_. + * @param[in] density Density of the simulated fluid in _kg / m^3_. + */ virtual void lbmInit(T dynViscosity, T density) = 0; + /** + * @brief Conducts the collide and stream operations of the lattice. + */ virtual void solve() = 0; /** @@ -124,40 +180,108 @@ class CFDSimulator { */ virtual bool hasConverged() const = 0; + /** TODO: + * + */ virtual void storePressures(std::unordered_map pressure) = 0; + /** + * @brief Store the abstract pressures at the nodes on the module boundary in the simulator. + * @param[in] pressure Map of pressures and node ids. + */ virtual std::unordered_map getPressures() const = 0; + /** + * @brief Store the abstract flow rates at the nodes on the module boundary in the simulator. + * @param[in] flowRate Map of flow rates and node ids. + */ virtual void storeFlowRates(std::unordered_map flowRate) = 0; + /** + * @brief Get the flow rates at the boundary nodes. + * @returns Flow rates in m^3/s. + */ virtual std::unordered_map getFlowRates() const = 0; - virtual void storeConcentrations(std::unordered_map> concentrations) { } + /** + * @brief Store the abstract concentrations at the nodes on the module boundary in the simulator. + * @param[in] concentrations Map of concentrations and node ids. + */ + virtual void storeConcentrations(std::unordered_map> concentrations) + { + throw std::runtime_error("The function storeConcentrations is undefined for this CFD simulator."); + } - virtual std::unordered_map> getConcentrations() const { return std::unordered_map>(); } + /** + * @brief Get the concentrations at the boundary nodes. + * @returns Concentrations + */ + virtual std::unordered_map> getConcentrations() const + { + throw std::runtime_error("The function storeConcentrations is undefined for this CFD simulator."); + return std::unordered_map>(); + } + /** + * @brief Set the boundary values on the lattice at the module nodes. + * @param[in] iT Iteration step. + */ virtual void setBoundaryValues(int iT) = 0; - // virtual functions - - virtual void prepareGeometry() {} + /** + * @brief Prepare the LBM geometry of this simulator. + */ + virtual void prepareGeometry() + { + throw std::runtime_error("The function prepareGeomtry is undefined for this CFD simulator."); + } - virtual void prepareLattice() {} + /** + * @brief Prepare the LBM lattice on the LBM geometry. + */ + virtual void prepareLattice() + { + throw std::runtime_error("The function prepareLattice is undefined for this CFD simulator."); + } /** * @brief Conducts the collide and stream operations of the NS lattice. */ - virtual void nsSolve() {} + virtual void nsSolve() + { + throw std::runtime_error("The function prepareLattice is undefined for this CFD simulator."); + } /** * @brief Conducts the collide and stream operations of the AD lattice(s). */ - virtual void adSolve() {} + virtual void adSolve() + { + throw std::runtime_error("The function prepareLattice is undefined for this CFD simulator."); + } - virtual void writeVTK (int iT) {} + /** + * @brief Write the vtk file with results of the CFD simulation to file system. + * @param[in] iT Iteration step. + */ + virtual void writeVTK (int iT) + { + throw std::runtime_error("The function prepareLattice is undefined for this CFD simulator."); + } - virtual void storeCfdResults (int iT) {} + /** + * @brief Update the values at the module nodes based on the simulation result after stepIter iterations. + * @param[in] iT Iteration step. + */ + virtual void storeCfdResults (int iT) + { + throw std::runtime_error("The function prepareLattice is undefined for this CFD simulator."); + } + /** + * @brief Returns whether the AD lattices have converged or not. + * @returns Boolean for AD convergence. + */ virtual bool hasAdConverged() const { return false; } friend void coupleNsAdLattices(const std::unordered_map>>& cfdSimulators); diff --git a/src/simulation/simulators/cfdSimulator.hh b/src/simulation/simulators/cfdSimulator.hh index dc58015..552436d 100755 --- a/src/simulation/simulators/cfdSimulator.hh +++ b/src/simulation/simulators/cfdSimulator.hh @@ -3,8 +3,8 @@ namespace sim{ template -CFDSimulator::CFDSimulator (int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map> openings_, T alpha_, ResistanceModel* resistanceModel_) : - id(id_), name(name_), stlFile(stlFile_), cfdModule(cfdModule_), moduleOpenings(openings_), alpha(alpha_) +CFDSimulator::CFDSimulator (int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map> openings_, ResistanceModel* resistanceModel_) : + id(id_), name(name_), stlFile(stlFile_), cfdModule(cfdModule_), moduleOpenings(openings_) { // Create this module's network, required for initial condition moduleNetwork = std::make_shared> (cfdModule_->getNodes()); @@ -23,6 +23,13 @@ CFDSimulator::CFDSimulator (int id_, std::string name_, std::string stlFile_, } } +template +CFDSimulator::CFDSimulator (int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map> openings_, std::shared_ptr> updateScheme_, ResistanceModel* resistanceModel_) : + CFDSimulator (id_, name_, stlFile_, cfdModule_, openings_, resistanceModel_) +{ + updateScheme = updateScheme_; +} + template int CFDSimulator::getId() const { return id; @@ -49,6 +56,11 @@ void CFDSimulator::setInitialized(bool initialization_) { this->initialized = initialization_; } +template +void CFDSimulator::setUpdateScheme(const std::shared_ptr>& updateScheme_) { + this->updateScheme = updateScheme_; +} + template void CFDSimulator::setVtkFolder(std::string vtkFolder_) { this->vtkFolder = vtkFolder_; @@ -60,8 +72,13 @@ std::string CFDSimulator::getVtkFile() { } template -T CFDSimulator::getAlpha() { - return alpha; +T CFDSimulator::getAlpha(int nodeId_) { + return updateScheme->getAlpha(nodeId_); +} + +template +T CFDSimulator::getBeta(int nodeId_) { + return updateScheme->getBeta(nodeId_); } template diff --git a/src/simulation/simulators/olbContinuous.h b/src/simulation/simulators/olbContinuous.h index 8ba28c4..03d355e 100644 --- a/src/simulation/simulators/olbContinuous.h +++ b/src/simulation/simulators/olbContinuous.h @@ -49,7 +49,6 @@ using BounceBack = olb::BounceBack; int step = 0; ///< Iteration step of this module. int stepIter = 1000; ///< Number of iterations for the value tracer. int maxIter = 1e7; ///< Maximum total iterations. - int theta = 10; ///< Number of OLB iterations per communication iteration. std::unordered_map pressures; ///< Vector of pressure values at module nodes. std::unordered_map flowRates; ///< Vector of flowRate values at module nodes. @@ -132,17 +131,36 @@ using BounceBack = olb::BounceBack; * @param[in] relaxationTime Relaxation time tau for the LBM solver. */ lbmSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, std::unordered_map> openings, - ResistanceModel* resistanceModel, T charPhysLenth, T charPhysVelocity, T alpha, T resolution, T epsilon, T relaxationTime=0.932); + ResistanceModel* resistanceModel, T charPhysLenth, T charPhysVelocity, T resolution, T epsilon, T relaxationTime=0.932); /** - * @brief Initialize an instance of the LBM solver for this module. + * @brief Constructor of an lbm module. + * @param[in] id Id of the module. + * @param[in] name Name of the module. + * @param[in] pos Absolute position of the module in _m_, from the bottom left corner of the microfluidic device. + * @param[in] size Size of the module in _m_. + * @param[in] nodes Map of nodes that are on the boundary of the module. + * @param[in] openings Map of the in-/outlets of the module. + * @param[in] stlFile STL file that describes the geometry of the CFD domain. + * @param[in] charPhysLength Characteristic physical length of the geometry of the module in _m_. + * @param[in] charPhysVelocity Characteristic physical velocity of the flow in the module in _m/s_. + * @param[in] alpha Relaxation factor for the iterative updates between the 1D and CFD solvers. + * @param[in] resolution Resolution of the CFD mesh in gridpoints per charPhysLength. + * @param[in] epsilon Convergence criterion for the pressure values at nodes on the boundary of the module. + * @param[in] relaxationTime Relaxation time tau for the LBM solver. + */ + lbmSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, std::unordered_map> openings, + std::shared_ptr> updateScheme, ResistanceModel* resistanceModel, T charPhysLenth, T charPhysVelocity, T resolution, T epsilon, T relaxationTime=0.932); + + /** + * @brief Initialize an instance of the LBM solver for this simulator. * @param[in] dynViscosity Dynamic viscosity of the simulated fluid in _kg / m s_. * @param[in] density Density of the simulated fluid in _kg / m^3_. */ void lbmInit(T dynViscosity, T density) override; /** - * @brief Prepare the LBM geometry of this instance. + * @brief Prepare the LBM geometry of this simulator. */ void prepareGeometry() override; diff --git a/src/simulation/simulators/olbContinuous.hh b/src/simulation/simulators/olbContinuous.hh index f958f82..46b0370 100644 --- a/src/simulation/simulators/olbContinuous.hh +++ b/src/simulation/simulators/olbContinuous.hh @@ -6,14 +6,23 @@ namespace sim{ template lbmSimulator::lbmSimulator ( int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map> openings_, - ResistanceModel* resistanceModel_, T charPhysLength_, T charPhysVelocity_, T alpha_, T resolution_, T epsilon_, T relaxationTime_) : - CFDSimulator(id_, name_, stlFile_, cfdModule_, openings_, alpha_, resistanceModel_), + ResistanceModel* resistanceModel_, T charPhysLength_, T charPhysVelocity_, T resolution_, T epsilon_, T relaxationTime_) : + CFDSimulator(id_, name_, stlFile_, cfdModule_, openings_, resistanceModel_), charPhysLength(charPhysLength_), charPhysVelocity(charPhysVelocity_), resolution(resolution_), epsilon(epsilon_), relaxationTime(relaxationTime_) { this->cfdModule->setModuleTypeLbm(); } +template +lbmSimulator::lbmSimulator ( + int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map> openings_, + std::shared_ptr> updateScheme_, ResistanceModel* resistanceModel_, T charPhysLength_, T charPhysVelocity_, T resolution_, T epsilon_, T relaxationTime_) : + lbmSimulator(id_, name_, stlFile_, cfdModule_, openings_, resistanceModel_, charPhysLength_, charPhysVelocity_, resolution_, epsilon_, relaxationTime_) +{ + this->updateScheme = updateScheme_; +} + template void lbmSimulator::lbmInit (T dynViscosity, T density) { @@ -129,9 +138,9 @@ void lbmSimulator::writeVTK (int iT) { template void lbmSimulator::solve() { - // theta = 10 + int theta = this->updateScheme->getTheta(this->cfdModule->getId()); this->setBoundaryValues(step); - for (int iT = 0; iT < 10; ++iT){ + for (int iT = 0; iT < theta; ++iT){ writeVTK(step); lattice->collideAndStream(); step += 1; diff --git a/src/simulation/simulators/olbMixing.h b/src/simulation/simulators/olbMixing.h index 77e2d10..368c4e1 100644 --- a/src/simulation/simulators/olbMixing.h +++ b/src/simulation/simulators/olbMixing.h @@ -119,8 +119,28 @@ using NoADDynamics = olb::NoDynamics; * @param[in] relaxationTime Relaxation time tau for the LBM solver. */ lbmMixingSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, std::unordered_map*> species, - std::unordered_map> openings, ResistanceModel* resistanceModel, T charPhysLenth, T charPhysVelocity, - T alpha, T resolution, T epsilon, T relaxationTime=0.932, T adRelaxationTime=0.932); + std::unordered_map> openings, ResistanceModel* resistanceModel, T charPhysLenth, + T charPhysVelocity, T resolution, T epsilon, T relaxationTime=0.932, T adRelaxationTime=0.932); + + /** + * @brief Constructor of an lbm module. + * @param[in] id Id of the module. + * @param[in] name Name of the module. + * @param[in] pos Absolute position of the module in _m_, from the bottom left corner of the microfluidic device. + * @param[in] size Size of the module in _m_. + * @param[in] nodes Map of nodes that are on the boundary of the module. + * @param[in] openings Map of the in-/outlets of the module. + * @param[in] stlFile STL file that describes the geometry of the CFD domain. + * @param[in] charPhysLength Characteristic physical length of the geometry of the module in _m_. + * @param[in] charPhysVelocity Characteristic physical velocity of the flow in the module in _m/s_. + * @param[in] alpha Relaxation factor for the iterative updates between the 1D and CFD solvers. + * @param[in] resolution Resolution of the CFD mesh in gridpoints per charPhysLength. + * @param[in] epsilon Convergence criterion for the pressure values at nodes on the boundary of the module. + * @param[in] relaxationTime Relaxation time tau for the LBM solver. + */ + lbmMixingSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> cfdModule, std::unordered_map*> species, + std::unordered_map> openings, std::shared_ptr> updateScheme, ResistanceModel* resistanceModel, T charPhysLenth, + T charPhysVelocity, T resolution, T epsilon, T relaxationTime=0.932, T adRelaxationTime=0.932); /** * @brief Initialize an instance of the LBM solver for this module. diff --git a/src/simulation/simulators/olbMixing.hh b/src/simulation/simulators/olbMixing.hh index 86f33f5..188ddcb 100644 --- a/src/simulation/simulators/olbMixing.hh +++ b/src/simulation/simulators/olbMixing.hh @@ -6,15 +6,26 @@ namespace sim{ template lbmMixingSimulator::lbmMixingSimulator ( int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map*> species_, - std::unordered_map> openings_, ResistanceModel* resistanceModel_, T charPhysLength_, T charPhysVelocity_, - T alpha_, T resolution_, T epsilon_, T relaxationTime_, T adRelaxationTime_) : - lbmSimulator(id_, name_, stlFile_, cfdModule_, openings_, resistanceModel_, charPhysLength_, charPhysVelocity_, alpha_, resolution_, epsilon_, relaxationTime_), + std::unordered_map> openings_, ResistanceModel* resistanceModel_, T charPhysLength_, + T charPhysVelocity_, T resolution_, T epsilon_, T relaxationTime_, T adRelaxationTime_) : + lbmSimulator(id_, name_, stlFile_, cfdModule_, openings_, resistanceModel_, charPhysLength_, charPhysVelocity_, resolution_, epsilon_, relaxationTime_), species(species_), adRelaxationTime(adRelaxationTime_) - { - std::cout << "Creating module and setting its type to lbm" << std::endl; - this->cfdModule->setModuleTypeLbm(); - fluxWall.try_emplace(int(0), &zeroFlux); - } +{ + std::cout << "Creating module and setting its type to lbm" << std::endl; + this->cfdModule->setModuleTypeLbm(); + fluxWall.try_emplace(int(0), &zeroFlux); +} + +template +lbmMixingSimulator::lbmMixingSimulator ( + int id_, std::string name_, std::string stlFile_, std::shared_ptr> cfdModule_, std::unordered_map*> species_, + std::unordered_map> openings_, std::shared_ptr> updateScheme_, ResistanceModel* resistanceModel_, T charPhysLength_, + T charPhysVelocity_, T resolution_, T epsilon_, T relaxationTime_, T adRelaxationTime_) : + lbmMixingSimulator(id_, name_, stlFile_, cfdModule_, openings_, updateScheme_, resistanceModel_, charPhysLength_, charPhysVelocity_, resolution_, epsilon_, + relaxationTime_, species_, adRelaxationTime_) +{ + this->updateScheme = updateScheme_; +} template void lbmMixingSimulator::lbmInit (T dynViscosity, T density) { diff --git a/src/simulation/simulators/olbOoc.h b/src/simulation/simulators/olbOoc.h index 7266e4b..ec77395 100644 --- a/src/simulation/simulators/olbOoc.h +++ b/src/simulation/simulators/olbOoc.h @@ -81,8 +81,28 @@ using NoADDynamics = olb::NoDynamics; * @param[in] relaxationTime Relaxation time tau for the LBM solver. */ lbmOocSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> tissue, std::string organStlFile, std::shared_ptr> cfdModule, - std::unordered_map*> species, std::unordered_map> openings, ResistanceModel* resistanceModel, T charPhysLenth, T charPhysVelocity, - T alpha, T resolution, T epsilon, T relaxationTime=0.932, T adRelaxationTime=0.932); + std::unordered_map*> species, std::unordered_map> openings, ResistanceModel* resistanceModel, + T charPhysLenth, T charPhysVelocity, T resolution, T epsilon, T relaxationTime=0.932, T adRelaxationTime=0.932); + + /** + * @brief Constructor of an lbm module. + * @param[in] id Id of the module. + * @param[in] name Name of the module. + * @param[in] pos Absolute position of the module in _m_, from the bottom left corner of the microfluidic device. + * @param[in] size Size of the module in _m_. + * @param[in] nodes Map of nodes that are on the boundary of the module. + * @param[in] openings Map of the in-/outlets of the module. + * @param[in] stlFile STL file that describes the geometry of the CFD domain. + * @param[in] charPhysLength Characteristic physical length of the geometry of the module in _m_. + * @param[in] charPhysVelocity Characteristic physical velocity of the flow in the module in _m/s_. + * @param[in] alpha Relaxation factor for the iterative updates between the 1D and CFD solvers. + * @param[in] resolution Resolution of the CFD mesh in gridpoints per charPhysLength. + * @param[in] epsilon Convergence criterion for the pressure values at nodes on the boundary of the module. + * @param[in] relaxationTime Relaxation time tau for the LBM solver. + */ + lbmOocSimulator(int id, std::string name, std::string stlFile, std::shared_ptr> tissue, std::string organStlFile, std::shared_ptr> cfdModule, + std::unordered_map*> species, std::unordered_map> openings, std::shared_ptr> updateScheme, ResistanceModel* resistanceModel, + T charPhysLenth, T charPhysVelocity, T resolution, T epsilon, T relaxationTime=0.932, T adRelaxationTime=0.932); /** * @brief Initialize an instance of the LBM solver for this module. diff --git a/src/simulation/simulators/olbOoc.hh b/src/simulation/simulators/olbOoc.hh index 566db19..5c58573 100644 --- a/src/simulation/simulators/olbOoc.hh +++ b/src/simulation/simulators/olbOoc.hh @@ -6,15 +6,26 @@ namespace sim{ template lbmOocSimulator::lbmOocSimulator ( int id_, std::string name_, std::string stlFile_, std::shared_ptr> tissue_, std::string organStlFile_, std::shared_ptr> cfdModule_, - std::unordered_map*> species_, std::unordered_map> openings_, ResistanceModel* resistanceModel_, T charPhysLength_, - T charPhysVelocity_, T alpha_, T resolution_, T epsilon_, T relaxationTime_, T adRelaxationTime_) : + std::unordered_map*> species_, std::unordered_map> openings_, + ResistanceModel* resistanceModel_, T charPhysLength_, T charPhysVelocity_, T resolution_, T epsilon_, T relaxationTime_, T adRelaxationTime_) : lbmMixingSimulator(id_, name_, stlFile_, cfdModule_, species_, openings_, resistanceModel_, charPhysLength_, charPhysVelocity_, - alpha_, resolution_, epsilon_, relaxationTime_, adRelaxationTime_), + resolution_, epsilon_, relaxationTime_, adRelaxationTime_), tissue(tissue_), organStlFile(organStlFile_) - { - this->cfdModule->setModuleTypeLbm(); - this->fluxWall.try_emplace(int(0), &this->zeroFlux); - } +{ + this->cfdModule->setModuleTypeLbm(); + this->fluxWall.try_emplace(int(0), &this->zeroFlux); +} + +template +lbmOocSimulator::lbmOocSimulator ( + int id_, std::string name_, std::string stlFile_, std::shared_ptr> tissue_, std::string organStlFile_, std::shared_ptr> cfdModule_, + std::unordered_map*> species_, std::unordered_map> openings_, std::shared_ptr> updateScheme_, + ResistanceModel* resistanceModel_, T charPhysLength_, T charPhysVelocity_, T resolution_, T epsilon_, T relaxationTime_, T adRelaxationTime_) : + lbmOocSimulator(id_, name_, stlFile_, tissue_, organStlFile_, cfdModule_, species_, openings_, updateScheme_, resistanceModel_, charPhysLength_, charPhysVelocity_, + resolution_, epsilon_, relaxationTime_, adRelaxationTime_) +{ + this->updateScheme = updateScheme_; +} template void lbmOocSimulator::lbmInit (T dynViscosity, T density) { diff --git a/tests/hybrid/Hybrid.test.cpp b/tests/hybrid/Hybrid.test.cpp index e8fa1ef..652ada1 100644 --- a/tests/hybrid/Hybrid.test.cpp +++ b/tests/hybrid/Hybrid.test.cpp @@ -14,6 +14,7 @@ TEST(Hybrid, Case1a) { // define network arch::Network network; + testSimulation.setNetwork(&network); // nodes auto node0 = network.addNode(0.0, 0.0, true); @@ -68,7 +69,6 @@ TEST(Hybrid, Case1a) { std::string stlFile = "../examples/STL/cross.stl"; T charPhysLength = 1e-4; T charPhysVelocity = 1e-1; - T alpha = 0.1; T resolution = 20; T epsilon = 1e-1; T tau = 0.55; @@ -78,7 +78,8 @@ TEST(Hybrid, Case1a) { Openings.try_emplace(8, arch::Opening(network.getNode(8), std::vector({0.0, 1.0}), 1e-4)); Openings.try_emplace(9, arch::Opening(network.getNode(9), std::vector({-1.0, 0.0}), 1e-4)); - testSimulation.addLbmSimulator(name, stlFile, network.getModule(m0->getId()), Openings, charPhysLength, charPhysVelocity, alpha, resolution, epsilon, tau); + testSimulation.addLbmSimulator(name, stlFile, network.getModule(m0->getId()), Openings, charPhysLength, charPhysVelocity, resolution, epsilon, tau); + testSimulation.setNaiveHybridScheme(0.1, 0.5, 10); network.sortGroups(); // pressure pump @@ -90,7 +91,6 @@ TEST(Hybrid, Case1a) { network.isNetworkValid(); // Simulate - testSimulation.setNetwork(&network); testSimulation.simulate(); EXPECT_NEAR(network.getNodes().at(0)->getPressure(), 0, 1e-2);