- added implementation of turbulence method: Liu et al., "Turbulent Details Simulation for SPH Fluids via Vorticity Refinement", Computer Graphics Forum, 2021

- Scene files now have a namespace for each pressure and non-pressure solver. This means existing scene files have to be adapted (see doc/file_format.md).
- fixed some bugs in the elasticity solvers

Author: janbender

CMake/NeighborhoodSearch.cmake

@@ -33,11 +33,14 @@ else(USE_GPU_NEIGHBORHOOD_SEARCH)
 	   Ext_NeighborhoodSearch
 	   PREFIX "${CMAKE_BINARY_DIR}/extern/CompactNSearch"
 	   GIT_REPOSITORY https://github.com/InteractiveComputerGraphics/CompactNSearch.git
-	   GIT_TAG "f6c0f32f3946554c91561cbaec4cfe8be0c23430"
+	   GIT_TAG "a9ab7c71ce264487660ecbaf81b5060bda462722"
 	   INSTALL_DIR ${ExternalInstallDir}/NeighborhoodSearch
 	   CMAKE_ARGS -DCMAKE_BUILD_TYPE=${EXT_CMAKE_BUILD_TYPE} -DCMAKE_CXX_FLAGS=${CMAKE_CXX_FLAGS} -DCMAKE_CXX_FLAGS_RELEASE=${CMAKE_CXX_FLAGS_RELEASE} -DCMAKE_INSTALL_PREFIX:PATH=${ExternalInstallDir}/NeighborhoodSearch -DUSE_DOUBLE_PRECISION:BOOL=${USE_DOUBLE_PRECISION} -DBUILD_DEMO:BOOL=OFF -DCMAKE_POLICY_VERSION_MINIMUM=3.10 -DCMAKE_PREFIX_PATH=${CMAKE_PREFIX_PATH}
 	)
 	set(NEIGHBORHOOD_ASSEMBLY_NAME CompactNSearch)
+	if(WIN32)
+		add_definitions( -D_SILENCE_STDEXT_ARR_ITERS_DEPRECATION_WARNING)
+	endif() 
 
 endif(USE_GPU_NEIGHBORHOOD_SEARCH)
 

Changelog.txt

@@ -1,3 +1,8 @@
+2.16.0
+	- added implementation of turbulence method: Liu et al., "Turbulent Details Simulation for SPH Fluids via Vorticity Refinement", Computer Graphics Forum, 2021
+	- Scene files now have a namespace for each pressure and non-pressure solver. This means existing scene files have to be adapted (see doc/file_format.md).
+	- fixed some bugs in the elasticity solvers
+
 2.15.1
 	- deactivated avx in Linux wheels due to problems 
 

SPlisHSPlasH/CMakeLists.txt

@@ -71,23 +71,20 @@ set(ICSPH_SOURCE_FILES
 	)	
 
 set(SURFACETENSION_HEADER_FILES
-	SurfaceTension/SurfaceTensionBase.h
 	SurfaceTension/SurfaceTension_Becker2007.h
 	SurfaceTension/SurfaceTension_Akinci2013.h
 	SurfaceTension/SurfaceTension_He2014.h
 	SurfaceTension/SurfaceTension_Jeske2023.h
 	)
 
 set(SURFACETENSION_SOURCE_FILES
-	SurfaceTension/SurfaceTensionBase.cpp
 	SurfaceTension/SurfaceTension_Becker2007.cpp
 	SurfaceTension/SurfaceTension_Akinci2013.cpp
 	SurfaceTension/SurfaceTension_He2014.cpp
 	SurfaceTension/SurfaceTension_Jeske2023.cpp
 	)	
 
 set(VISCOSITY_HEADER_FILES
-	Viscosity/ViscosityBase.h
 	Viscosity/Viscosity_Standard.h
 	Viscosity/Viscosity_Bender2017.h
 	Viscosity/Viscosity_Peer2015.h
@@ -97,7 +94,6 @@ set(VISCOSITY_HEADER_FILES
 	)
 
 set(VISCOSITY_SOURCE_FILES
-	Viscosity/ViscosityBase.cpp
 	Viscosity/Viscosity_Standard.cpp
 	Viscosity/Viscosity_Bender2017.cpp
 	Viscosity/Viscosity_Peer2015.cpp
@@ -107,38 +103,34 @@ set(VISCOSITY_SOURCE_FILES
 	)	
 
 set(VORTICITY_HEADER_FILES
-	Vorticity/VorticityBase.h
 	Vorticity/MicropolarModel_Bender2017.h
 	Vorticity/VorticityConfinement.h
+	Vorticity/VorticityRefinement_Liu2021.h
 	)
 
 set(VORTICITY_SOURCE_FILES
-	Vorticity/VorticityBase.cpp
 	Vorticity/MicropolarModel_Bender2017.cpp
 	Vorticity/VorticityConfinement.cpp
+	Vorticity/VorticityRefinement_Liu2021.cpp
 	)	
 
 set(DRAG_HEADER_FILES
-	Drag/DragBase.h
 	Drag/DragForce_Gissler2017.h
 	Drag/DragForce_Macklin2014.h
 	)
 
 set(DRAG_SOURCE_FILES
-	Drag/DragBase.cpp
 	Drag/DragForce_Gissler2017.cpp
 	Drag/DragForce_Macklin2014.cpp
 	)	
 
 set(ELASTICITY_HEADER_FILES
-	Elasticity/ElasticityBase.h
 	Elasticity/Elasticity_Becker2009.h	
 	Elasticity/Elasticity_Peer2018.h
 	Elasticity/Elasticity_Kugelstadt2021.h
 	)
 
 set(ELASTICITY_SOURCE_FILES
-	Elasticity/ElasticityBase.cpp
 	Elasticity/Elasticity_Becker2009.cpp
 	Elasticity/Elasticity_Peer2018.cpp
 	Elasticity/Elasticity_Kugelstadt2021.cpp

SPlisHSPlasH/DFSPH/TimeStepDFSPH.cpp

@@ -15,7 +15,11 @@ using namespace SPH;
 using namespace std;
 using namespace GenParam;
 
-
+std::string TimeStepDFSPH::METHOD_NAME = "DFSPH";
+int TimeStepDFSPH::SOLVER_ITERATIONS = -1;
+int TimeStepDFSPH::MIN_ITERATIONS = -1;
+int TimeStepDFSPH::MAX_ITERATIONS = -1;
+int TimeStepDFSPH::MAX_ERROR = -1;
 int TimeStepDFSPH::SOLVER_ITERATIONS_V = -1;
 int TimeStepDFSPH::MAX_ITERATIONS_V = -1;
 int TimeStepDFSPH::MAX_ERROR_V = -1;
@@ -27,6 +31,10 @@ TimeStepDFSPH::TimeStepDFSPH() :
 	m_simulationData()
 {
 	m_simulationData.init();
+	m_iterations = 0;
+	m_minIterations = 2;
+	m_maxIterations = 100;
+	m_maxError = static_cast<Real>(0.01);
 	m_iterationsV = 0;
 	m_enableDivergenceSolver = true;
 	m_maxIterationsV = 100;
@@ -38,11 +46,11 @@ TimeStepDFSPH::TimeStepDFSPH() :
 	for (unsigned int fluidModelIndex = 0; fluidModelIndex < nModels; fluidModelIndex++)
 	{
 		FluidModel *model = sim->getFluidModel(fluidModelIndex);
-		model->addField({ "factor", FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getFactor(fluidModelIndex, i); } });
-		model->addField({ "advected density", FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getDensityAdv(fluidModelIndex, i); } });
-		model->addField({ "p / rho^2", FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getPressureRho2(fluidModelIndex, i); }, true });
-		model->addField({ "p_v / rho^2", FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getPressureRho2_V(fluidModelIndex, i); }, true });
-		model->addField({ "pressure acceleration", FieldType::Vector3, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getPressureAccel(fluidModelIndex, i)[0]; } });
+		model->addField({ "factor", METHOD_NAME, FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getFactor(fluidModelIndex, i); } });
+		model->addField({ "advected density", METHOD_NAME, FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getDensityAdv(fluidModelIndex, i); } });
+		model->addField({ "p / rho^2", METHOD_NAME, FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getPressureRho2(fluidModelIndex, i); }, true });
+		model->addField({ "p_v / rho^2", METHOD_NAME, FieldType::Scalar, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getPressureRho2_V(fluidModelIndex, i); }, true });
+		model->addField({ "pressure acceleration", METHOD_NAME, FieldType::Vector3, [this, fluidModelIndex](const unsigned int i) -> Real* { return &m_simulationData.getPressureAccel(fluidModelIndex, i)[0]; } });
 	}
 }
 
@@ -66,6 +74,26 @@ void TimeStepDFSPH::initParameters()
 {
 	TimeStep::initParameters();
 
+	SOLVER_ITERATIONS = createNumericParameter("iterations", "Iterations", &m_iterations);
+	setGroup(SOLVER_ITERATIONS, "Simulation|DFSPH");
+	setDescription(SOLVER_ITERATIONS, "Iterations required by the pressure solver.");
+	getParameter(SOLVER_ITERATIONS)->setReadOnly(true);
+
+	MIN_ITERATIONS = createNumericParameter("minIterations", "Min. iterations", &m_minIterations);
+	setGroup(MIN_ITERATIONS, "Simulation|DFSPH");
+	setDescription(MIN_ITERATIONS, "Minimal number of iterations of the pressure solver.");
+	static_cast<NumericParameter<unsigned int>*>(getParameter(MIN_ITERATIONS))->setMinValue(0);
+
+	MAX_ITERATIONS = createNumericParameter("maxIterations", "Max. iterations", &m_maxIterations);
+	setGroup(MAX_ITERATIONS, "Simulation|DFSPH");
+	setDescription(MAX_ITERATIONS, "Maximal number of iterations of the pressure solver.");
+	static_cast<NumericParameter<unsigned int>*>(getParameter(MAX_ITERATIONS))->setMinValue(1);
+
+	MAX_ERROR = createNumericParameter("maxError", "Max. density error(%)", &m_maxError);
+	setGroup(MAX_ERROR, "Simulation|DFSPH");
+	setDescription(MAX_ERROR, "Maximal density error (%).");
+	static_cast<RealParameter*>(getParameter(MAX_ERROR))->setMinValue(static_cast<Real>(1e-6));
+
 	SOLVER_ITERATIONS_V = createNumericParameter("iterationsV", "Iterations (divergence)", &m_iterationsV);
 	setGroup(SOLVER_ITERATIONS_V, "Simulation|DFSPH");
 	setDescription(SOLVER_ITERATIONS_V, "Iterations required by the divergence solver.");
@@ -575,7 +603,7 @@ void TimeStepDFSPH::pressureSolveIteration(const unsigned int fluidModelIndex, R
 			const Real residuum = min(s_i - aij_pj, static_cast<Real>(0.0));     // r = b - A*p
 			//p_rho2_i -= residuum * m_simulationData.getFactor(fluidModelIndex, i);
 
-			p_rho2_i = max(p_rho2_i - 0.5 * (s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), 0.0);
+			p_rho2_i = max(p_rho2_i - static_cast<Real>(0.5) * (s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), static_cast<Real>(0.0));
 
 			//////////////////////////////////////////////////////////////////////////
 			// Compute the sum of the density errors
@@ -661,7 +689,7 @@ void TimeStepDFSPH::divergenceSolveIteration(const unsigned int fluidModelIndex,
 					residuum = 0.0;
 			}
 			//pv_rho2_i -= residuum * m_simulationData.getFactor(fluidModelIndex, i);
-			pv_rho2_i = max(pv_rho2_i - 0.5*(s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), 0.0);
+			pv_rho2_i = max(pv_rho2_i - static_cast<Real>(0.5)*(s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), static_cast<Real>(0.0));
 
 
 			//////////////////////////////////////////////////////////////////////////

SPlisHSPlasH/DFSPH/TimeStepDFSPH.h

@@ -26,6 +26,10 @@ namespace SPH
 	protected:
 		SimulationDataDFSPH m_simulationData;
 		const Real m_eps = static_cast<Real>(1.0e-5);
+		unsigned int m_iterations;
+		Real m_maxError;
+		unsigned int m_minIterations;
+		unsigned int m_maxIterations;
 		bool m_enableDivergenceSolver;
 		unsigned int m_iterationsV;
 		Real m_maxErrorV;
@@ -49,6 +53,11 @@ namespace SPH
 		virtual void initParameters();
 
 	public:
+		static std::string METHOD_NAME;
+		static int SOLVER_ITERATIONS;
+		static int MIN_ITERATIONS;
+		static int MAX_ITERATIONS;
+		static int MAX_ERROR;
 		static int SOLVER_ITERATIONS_V;
 		static int MAX_ITERATIONS_V;
 		static int MAX_ERROR_V;
@@ -62,6 +71,8 @@ namespace SPH
 		virtual void reset();
 
 		virtual void resize();
+		virtual std::string getMethodName() { return METHOD_NAME; }
+		virtual int getNumIterations() { return m_iterations; }
 	};
 }
 

SPlisHSPlasH/Drag/DragBase.cpp

@@ -7,20 +7,36 @@
 #include "../BoundaryModel_Bender2019.h"
 
 using namespace SPH;
+using namespace GenParam;
+
+std::string DragForce_Gissler2017::METHOD_NAME = "Gissler et al. 2017";
+int DragForce_Gissler2017::DRAG_COEFFICIENT = -1;
 
 // C.Gissler, S.Band, A.Peer, M.Ihmsen, M.Teschner,
 // "Approximate Air-Fluid Interactions for SPH,"
 // VRIPHYS 2017
 
 DragForce_Gissler2017::DragForce_Gissler2017(FluidModel *model) :
-	DragBase(model)
+	NonPressureForceBase(model)
 {
+	m_dragCoefficient = static_cast<Real>(0.01);
 }
 
 DragForce_Gissler2017::~DragForce_Gissler2017(void)
 {
 }
 
+void DragForce_Gissler2017::initParameters()
+{
+	NonPressureForceBase::initParameters();
+
+	DRAG_COEFFICIENT = createNumericParameter("drag", "Drag coefficient", &m_dragCoefficient);
+	setGroup(DRAG_COEFFICIENT, "Fluid Model|Drag force");
+	setDescription(DRAG_COEFFICIENT, "Coefficient for the drag force computation");
+	RealParameter* rparam = static_cast<RealParameter*>(getParameter(DRAG_COEFFICIENT));
+	rparam->setMinValue(0.0);
+}
+
 #ifdef USE_AVX
 
 void DragForce_Gissler2017::step()
@@ -69,7 +85,7 @@ void DragForce_Gissler2017::step()
 	const Real y_coeff = (C_F * We_i_wo_v * c_def) / (C_k * C_b);
 
 	const Real n_full = 38;
-	const Real n_full_23 = n_full * 2.0/3.0;
+	const Real n_full_23 = n_full * static_cast<Real>(2.0/3.0);
 
 	#pragma omp parallel default(shared)
 	{

SPlisHSPlasH/Drag/DragBase.h

@@ -3,7 +3,7 @@
 
 #include "SPlisHSPlasH/Common.h"
 #include "SPlisHSPlasH/FluidModel.h"
-#include "DragBase.h"
+#include "SPlisHSPlasH/NonPressureForceBase.h"
 
 namespace SPH
 {
@@ -13,9 +13,10 @@ namespace SPH
 	* References:
 	* - [GPB+17] Christoph Gissler, Stefan Band, Andreas Peer, Markus Ihmsen, and Matthias Teschner. Approximate air-fluid interactions for SPH. In Virtual Reality Interactions and Physical Simulations, 1-10. April 2017. URL: http://dx.doi.org/10.2312/vriphys.20171081
 	*/
-	class DragForce_Gissler2017 : public DragBase
+	class DragForce_Gissler2017 : public NonPressureForceBase
 	{
 	protected:
+		Real m_dragCoefficient;
 		const Real rho_a = static_cast<Real>(1.2041);
 		const Real sigma = static_cast<Real>(0.0724);
 		const Real mu_l = static_cast<Real>(0.00102);
@@ -25,11 +26,17 @@ namespace SPH
 		const Real C_b = static_cast<Real>(0.5);
 		const Real mu_a = static_cast<Real>(0.00001845);
 
+		virtual void initParameters();
+
 	public:
+		static std::string METHOD_NAME;
+		static int DRAG_COEFFICIENT;
+
 		DragForce_Gissler2017(FluidModel *model);
 		virtual ~DragForce_Gissler2017(void);
 
 		static NonPressureForceBase* creator(FluidModel* model) { return new DragForce_Gissler2017(model); }
+		virtual std::string getMethodName() { return METHOD_NAME; }
 
 		virtual void step();
 		virtual void reset();

SPlisHSPlasH/Drag/DragForce_Gissler2017.cpp

@@ -2,16 +2,32 @@
 #include "SPlisHSPlasH/TimeManager.h"
 
 using namespace SPH;
+using namespace GenParam;
+
+std::string DragForce_Macklin2014::METHOD_NAME = "Macklin et al. 2014";
+int DragForce_Macklin2014::DRAG_COEFFICIENT = -1;
 
 DragForce_Macklin2014::DragForce_Macklin2014(FluidModel *model) :
-	DragBase(model)
+	NonPressureForceBase(model)
 {
+	m_dragCoefficient = static_cast<Real>(0.01);
 }
 
 DragForce_Macklin2014::~DragForce_Macklin2014(void)
 {
 }
 
+void DragForce_Macklin2014::initParameters()
+{
+	NonPressureForceBase::initParameters();
+
+	DRAG_COEFFICIENT = createNumericParameter("drag", "Drag coefficient", &m_dragCoefficient);
+	setGroup(DRAG_COEFFICIENT, "Fluid Model|Drag force");
+	setDescription(DRAG_COEFFICIENT, "Coefficient for the drag force computation");
+	RealParameter* rparam = static_cast<RealParameter*>(getParameter(DRAG_COEFFICIENT));
+	rparam->setMinValue(0.0);
+}
+
 void DragForce_Macklin2014::step()
 {
 	const Real density0 = m_model->getDensity0();