Skip to content

Commit 5aaa034

Browse files
janbenderjanbender
committed
- added citation file
- fixed rigid body colors in visualization - fixed export of static rigid bodies
1 parent a081124 commit 5aaa034

File tree

12 files changed

+64
-189
lines changed

12 files changed

+64
-189
lines changed

‎CITATION.cff‎

Lines changed: 22 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,22 @@
1+
cff-version: 1.2.0
2+
title: SPlisHSPlasH Library
3+
message: >-
4+
If you use this software, please cite it using the
5+
metadata from this file.
6+
type: software
7+
authors:
8+
- given-names: Jan
9+
family-names: Bender
10+
identifiers:
11+
- type: url
12+
value: >-
13+
https://github.com/InteractiveComputerGraphics/SPlisHSPlasH
14+
repository-code: >-
15+
https://github.com/InteractiveComputerGraphics/SPlisHSPlasH
16+
keywords:
17+
- Smoothed Particle Hydrodynamics
18+
- SPH
19+
- fluids
20+
- simulation
21+
license: MIT
22+
version: 2.12.0

‎Changelog.txt‎

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,8 @@
1+
2.12.4
2+
- added citation file
3+
- fixed rigid body colors in visualization
4+
- fixed export of static rigid bodies
5+
16
2.12.3
27
- added missing include
38
- fixed addKeyFunc in Python interface

‎README.md‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -262,3 +262,12 @@ https://link.springer.com/chapter/10.1007/978-981-16-1781-2_38
262262
https://openreview.net/forum?id=FEBFJ98FKx
263263
* Muzaffer Akbay, Nicholas Nobles, Victor Zordan, Tamar Shinar. An extended partitioned method for conservative solid-fluid coupling. ACM Transactions on Graphics 37, 4, 2018
264264
https://dl.acm.org/doi/10.1145/3197517.3201345
265+
* Samuel Carensac, Nicolas Pronost & Saïda Bouakaz. Optimizations for predictive–corrective particle-based fluid simulation on GPU. The Visual Computer volume 39, 2023
266+
* Iliya Starodubtsev, Pavel Vasev, Yuliya Starodubtseva, Igor Tsepelev. Numerical Simulation and Visualization of Lava Flows. Scientific Visualization, 2022, volume 14, number 5
267+
https://doi.org/10.26583/sv.14.5.05
268+
* Arjun Mani, Ishaan Preetam Chandratreya, Elliot Creager, Carl Vondrick, Richard Zemel. SurfsUp: Learning Fluid Simulation for Novel Surfaces. arXiv 2023
269+
https://doi.org/10.48550/arXiv.2304.06197
270+
* Jian Chen. Incompressible Projective Smooth Particle Hydrodynamic. Proc. CNCI 2019
271+
https://doi.org/10.2991/cnci-19.2019.80
272+
* Min Hyung Kee, Kiwon Um, HyunMo Kang, and JungHyun Han. An Optimization-based SPH Solver for Simulation of Hyperelastic Solids. Computer Graphics Forum, 42, 2, 2023
273+
https://doi.org/10.1111/cgf.14756

‎SPlisHSPlasH/DFSPH/TimeStepDFSPH.cpp‎

Lines changed: 4 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -574,10 +574,9 @@ void TimeStepDFSPH::pressureSolveIteration(const unsigned int fluidModelIndex, R
574574
//////////////////////////////////////////////////////////////////////////
575575
Real& p_rho2_i = m_simulationData.getPressureRho2(fluidModelIndex, i);
576576
const Real residuum = min(s_i - aij_pj, static_cast<Real>(0.0)); // r = b - A*p
577-
p_rho2_i -= residuum * m_simulationData.getFactor(fluidModelIndex, i);
578-
579-
//p_rho2_i = max(p_rho2_i - 0.5*(s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), 0.0);
577+
//p_rho2_i -= residuum * m_simulationData.getFactor(fluidModelIndex, i);
580578

579+
p_rho2_i = max(p_rho2_i - 0.5 * (s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), 0.0);
581580

582581
//////////////////////////////////////////////////////////////////////////
583582
// Compute the sum of the density errors
@@ -662,9 +661,8 @@ void TimeStepDFSPH::divergenceSolveIteration(const unsigned int fluidModelIndex,
662661
if (numNeighbors < 7)
663662
residuum = 0.0;
664663
}
665-
pv_rho2_i -= residuum * m_simulationData.getFactor(fluidModelIndex, i);
666-
667-
//pv_rho2_i = max(pv_rho2_i - 0.5*(s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), 0.0);
664+
//pv_rho2_i -= residuum * m_simulationData.getFactor(fluidModelIndex, i);
665+
pv_rho2_i = max(pv_rho2_i - 0.5*(s_i - aij_pj) * m_simulationData.getFactor(fluidModelIndex, i), 0.0);
668666

669667

670668
//////////////////////////////////////////////////////////////////////////

‎SPlisHSPlasH/Utilities/SceneParameterObjects.h‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -327,7 +327,7 @@ namespace Utilities
327327
Vector3r scale;
328328
bool dynamic;
329329
bool isWall;
330-
Eigen::Matrix<float, 4, 1, Eigen::DontAlign> color;
330+
Vector4r color;
331331
std::string mapFile;
332332
bool mapInvert;
333333
Real mapThickness;
@@ -346,7 +346,7 @@ namespace Utilities
346346
scale = Vector3r::Ones();
347347
dynamic = false;
348348
isWall = false;
349-
color = Eigen::Vector4f(1.0f, 0.0f, 0.0f, 0.0f);
349+
color = Vector4r(1.0, 0.0, 0.0, 0.0);
350350
samplingMode = 0;
351351
isAnimated = false;
352352
// Maps
@@ -357,7 +357,7 @@ namespace Utilities
357357
}
358358

359359
BoundaryParameterObject(std::string samplesFile_, std::string meshFile_, Vector3r translation_, Vector3r axis_, Real angle_, Vector3r scale_,
360-
bool dynamic_, bool isWall_, Eigen::Matrix<float, 4, 1, Eigen::DontAlign> color_, std::string mapFile_, bool mapInvert_,
360+
bool dynamic_, bool isWall_, Vector4r color_, std::string mapFile_, bool mapInvert_,
361361
Real mapThickness_, Eigen::Matrix<unsigned int, 3, 1, Eigen::DontAlign> mapResolution_, unsigned int samplingMode_, bool isAnimated_)
362362
{
363363
samplesFile = samplesFile_;

‎SPlisHSPlasH/Vorticity/MicropolarModel_Bender2017.cpp‎

Lines changed: 0 additions & 132 deletions
Original file line numberDiff line numberDiff line change
@@ -132,71 +132,6 @@ void MicropolarModel_Bender2017::step()
132132
delta_angAcceli_avx += ((vi_avx - vj_avx) % V_gradW) * (nut_density_i_intertiaInverse * density0_avx);
133133
);
134134

135-
//////////////////////////////////////////////////////////////////////////
136-
// Boundary
137-
//////////////////////////////////////////////////////////////////////////
138-
if (sim->getBoundaryHandlingMethod() == BoundaryHandlingMethods::Akinci2012)
139-
{
140-
forall_boundary_neighbors_avx(
141-
const Vector3f8 vj_avx = convertVec_zero(&sim->getNeighborList(fluidModelIndex, pid, i)[j], &bm_neighbor->getVelocity(0), count);
142-
143-
// Viscosity
144-
const Vector3f8 xixj = xi_avx - xj_avx;
145-
const Vector3f8 &omegaij = omegai_avx; // ToDo: omega_j
146-
const Scalarf8 Vj_avx = convert_zero(&sim->getNeighborList(fluidModelIndex, pid, i)[j], &bm_neighbor->getVolume(0), count);
147-
const Vector3f8 gradW = CubicKernel_AVX::gradW(xixj);
148-
const Scalarf8 mj_avx = density0_avx * Vj_avx;
149-
150-
// XSPH for angular velocity field
151-
//delta_angAcceli_avx -= omegaij * factor_avx * (mj_avx / density_i_avx) * CubicKernel_AVX::W(xixj);
152-
153-
// difference curl
154-
const Vector3f8 a = (omegaij % gradW) * (nut_density_i * mj_avx);
155-
delta_ai_avx += a;
156-
delta_angAcceli_avx += ((vi_avx - vj_avx) % gradW) * (nut_density_i_intertiaInverse * mj_avx);
157-
bm_neighbor->addForce(xj_avx, -a * mi_avx, count);
158-
);
159-
}
160-
else if (sim->getBoundaryHandlingMethod() == BoundaryHandlingMethods::Koschier2017)
161-
{
162-
forall_density_maps(
163-
Vector3r vj;
164-
bm_neighbor->getPointVelocity(xi, vj);
165-
const Vector3r &omegaij = omegai; // ToDo: omega_j
166-
167-
// XSPH for angular velocity field
168-
//angAcceli -= invDt * m_inertiaInverse * zeta * (density0 / density_i) * omegaij * rho;
169-
170-
// difference curl
171-
const Vector3r a = nu_t * density0 / density_i * (omegaij.cross(gradRho));
172-
ai += a;
173-
angAcceli += nu_t * density0 / density_i * m_inertiaInverse * ((vi - vj).cross(gradRho));
174-
175-
bm_neighbor->addForce(xj, -m_model->getMass(i) * a);
176-
);
177-
}
178-
else if (sim->getBoundaryHandlingMethod() == BoundaryHandlingMethods::Bender2019)
179-
{
180-
forall_volume_maps(
181-
Vector3r vj;
182-
bm_neighbor->getPointVelocity(xj, vj);
183-
const Vector3r omegaij = omegai; // ToDo: omega_j
184-
185-
const Vector3r xij = xi - xj;
186-
const Vector3r gradW = sim->gradW(xij);
187-
188-
// XSPH for angular velocity field
189-
//angAcceli -= (invDt * m_inertiaInverse * zeta * (density0 * Vj / density_i) * sim->W(xij)) * omegaij;
190-
191-
// difference curl
192-
const Vector3r a = (nu_t * 1.0 / density_i * density0 * Vj) * (omegaij.cross(gradW));
193-
ai += a;
194-
angAcceli += (nu_t * 1.0 / density_i * m_inertiaInverse * density0 * Vj) * (vi - vj).cross(gradW);
195-
196-
bm_neighbor->addForce(xj, -m_model->getMass(i) * a);
197-
);
198-
}
199-
200135
ai[0] += delta_ai_avx.x().reduce();
201136
ai[1] += delta_ai_avx.y().reduce();
202137
ai[2] += delta_ai_avx.z().reduce();
@@ -283,73 +218,6 @@ void MicropolarModel_Bender2017::step()
283218
angAcceli += nu_t * 1.0/density_i * m_inertiaInverse * (m_model->getMass(neighborIndex) * (vi - vj).cross(gradW));
284219
);
285220

286-
//////////////////////////////////////////////////////////////////////////
287-
// Boundary
288-
//////////////////////////////////////////////////////////////////////////
289-
if (sim->getBoundaryHandlingMethod() == BoundaryHandlingMethods::Akinci2012)
290-
{
291-
forall_boundary_neighbors(
292-
const Vector3r &vj = bm_neighbor->getVelocity(neighborIndex);
293-
const Vector3r &omegaj = Vector3r::Zero();//m_omega[neighborIndex];
294-
295-
// Viscosity
296-
const Vector3r xij = xi - xj;
297-
const Vector3r omegaij = omegai - omegaj;
298-
const Vector3r gradW = sim->gradW(xij);
299-
300-
// XSPH for angular velocity field
301-
//angAcceli -= invDt * m_inertiaInverse * zeta * (density0 * bm_neighbor->getVolume(neighborIndex) / density_i) * omegaij * sim->W(xij);
302-
303-
// Viscosity
304-
//angAcceli += d * m_inertiaInverse * zeta * (density0 * bm_neighbor->getVolume(neighborIndex) / density_i) * omegaij.dot(xij) / (xij.squaredNorm() + 0.01*h2) * gradW;
305-
306-
// difference curl
307-
const Vector3r a = nu_t * 1.0 / density_i * density0 * bm_neighbor->getVolume(neighborIndex) * (omegaij.cross(gradW));
308-
ai += a;
309-
angAcceli += nu_t * 1.0 / density_i * m_inertiaInverse * (density0 * bm_neighbor->getVolume(neighborIndex) * (vi - vj).cross(gradW));
310-
311-
bm_neighbor->addForce(xj, -model->getMass(i) * a);
312-
);
313-
}
314-
else if (sim->getBoundaryHandlingMethod() == BoundaryHandlingMethods::Koschier2017)
315-
{
316-
forall_density_maps(
317-
Vector3r vj;
318-
bm_neighbor->getPointVelocity(xi, vj);
319-
const Vector3r omegaij = omegai; // ToDo: omega_j
320-
321-
// // XSPH for angular velocity field
322-
// angAcceli -= invDt * m_inertiaInverse * zeta * (density0 / density_i) * omegaij * rho;
323-
324-
// difference curl
325-
const Vector3r a = nu_t * density0 / density_i * (omegaij.cross(gradRho));
326-
ai += a;
327-
angAcceli += nu_t * density0 / density_i * m_inertiaInverse * ((vi - vj).cross(gradRho));
328-
329-
bm_neighbor->addForce(xj, -model->getMass(i) * a);
330-
);
331-
}
332-
else if (sim->getBoundaryHandlingMethod() == BoundaryHandlingMethods::Bender2019)
333-
{
334-
forall_volume_maps(
335-
Vector3r vj;
336-
bm_neighbor->getPointVelocity(xj, vj);
337-
const Vector3r &omegaij = omegai; // ToDo: omega_j
338-
339-
const Vector3r xij = xi - xj;
340-
const Vector3r gradW = sim->gradW(xij);
341-
342-
// // XSPH for angular velocity field
343-
// angAcceli -= invDt * m_inertiaInverse * zeta * (density0 * Vj / density_i) * omegaij * sim->W(xij);
344-
345-
// difference curl
346-
const Vector3r a = nu_t * 1.0 / density_i * density0 * Vj * (omegaij.cross(gradW));
347-
ai += a;
348-
angAcceli += nu_t * 1.0 / density_i * m_inertiaInverse * (density0 * Vj * (vi - vj).cross(gradW));
349-
350-
bm_neighbor->addForce(xj, -model->getMass(i) * a);
351-
);
352-
}
353221
angAcceli -= 2.0 * m_inertiaInverse * nu_t * omegai;
354222
}
355223
}

‎Simulator/Exporter/RigidBodyExporter_OBJ.cpp‎

Lines changed: 4 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -51,23 +51,13 @@ void RigidBodyExporter_OBJ::writeRigidBodies(const unsigned int frame)
5151
Simulation* sim = Simulation::getCurrent();
5252
const unsigned int nBoundaryModels = sim->numberOfBoundaryModels();
5353

54-
// check if we have a static model
55-
bool isStatic = true;
5654
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
5755
{
5856
BoundaryModel* bm = sim->getBoundaryModel(i);
59-
if (bm->getRigidBodyObject()->isDynamic() || bm->getRigidBodyObject()->isAnimated())
60-
{
61-
isStatic = false;
62-
break;
63-
}
64-
}
6557

66-
// If we have a static model, write the data only for the first frame.
67-
// Otherwise each frame is exported.
68-
if (m_isFirstFrame || !isStatic)
69-
{
70-
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
58+
// If we have a static body, write the data only for the first frame.
59+
// Otherwise each frame is exported.
60+
if (m_isFirstFrame || bm->getRigidBodyObject()->isDynamic() || bm->getRigidBodyObject()->isAnimated())
7161
{
7262
std::string fileName = "rb_data_";
7363
fileName = fileName + std::to_string(i) + "_" + std::to_string(frame) + ".obj";
@@ -85,7 +75,6 @@ void RigidBodyExporter_OBJ::writeRigidBodies(const unsigned int frame)
8575
outfile << "# Created by SPlisHSPlasH version " << SPLISHSPLASH_VERSION << "\n";
8676
outfile << "g default\n";
8777

88-
BoundaryModel* bm = sim->getBoundaryModel(i);
8978
const std::vector<Vector3r>& vertices = bm->getRigidBodyObject()->getVertices();
9079
const std::vector<unsigned int>& faces = bm->getRigidBodyObject()->getFaces();
9180
int n_vertices = (int)vertices.size();
@@ -104,7 +93,7 @@ void RigidBodyExporter_OBJ::writeRigidBodies(const unsigned int frame)
10493
{
10594
for (int j = 0; j < n_triangles; j++)
10695
{
107-
outfile << "f " << faces[3 * j + 0] + 1 << " " << faces[3* j + 1] + 1 << " " << faces[3 * j + 2] + 1 << "\n";
96+
outfile << "f " << faces[3 * j + 0] + 1 << " " << faces[3 * j + 1] + 1 << " " << faces[3 * j + 2] + 1 << "\n";
10897
}
10998
}
11099
outfile.close();

‎Simulator/Exporter/RigidBodyExporter_PLY.cpp‎

Lines changed: 4 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -52,29 +52,18 @@ void RigidBodyExporter_PLY::writeRigidBodies(const unsigned int frame)
5252
Simulation* sim = Simulation::getCurrent();
5353
const unsigned int nBoundaryModels = sim->numberOfBoundaryModels();
5454

55-
// check if we have a static model
56-
bool isStatic = true;
5755
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
5856
{
5957
BoundaryModel* bm = sim->getBoundaryModel(i);
60-
if (bm->getRigidBodyObject()->isDynamic() || bm->getRigidBodyObject()->isAnimated())
61-
{
62-
isStatic = false;
63-
break;
64-
}
65-
}
6658

67-
// If we have a static model, write the data only for the first frame.
68-
// Otherwise each frame is exported.
69-
if (m_isFirstFrame || !isStatic)
70-
{
71-
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
59+
// If we have a static body, write the data only for the first frame.
60+
// Otherwise each frame is exported.
61+
if (m_isFirstFrame || bm->getRigidBodyObject()->isDynamic() || bm->getRigidBodyObject()->isAnimated())
7262
{
7363
std::string fileName = "rb_data_";
7464
fileName = fileName + std::to_string(i) + "_" + std::to_string(frame) + ".ply";
7565
std::string exportFileName = FileSystem::normalizePath(m_exportPath + "/" + fileName);
7666

77-
BoundaryModel* bm = sim->getBoundaryModel(i);
7867
const std::vector<Vector3r>& vertices = bm->getRigidBodyObject()->getVertices();
7968
const std::vector<unsigned int>& faces = bm->getRigidBodyObject()->getFaces();
8069
int n_vertices = (int)vertices.size();
@@ -83,7 +72,7 @@ void RigidBodyExporter_PLY::writeRigidBodies(const unsigned int frame)
8372
// Suppose these hold your data
8473
std::vector<std::array<double, 3>> meshVertexPositions;
8574
std::vector<std::vector<size_t>> meshFaceIndices;
86-
75+
8776
// vertices
8877
meshVertexPositions.resize(n_vertices);
8978
for (int j = 0u; j < n_vertices; j++)

‎Simulator/Exporter/RigidBodyExporter_VTK.cpp‎

Lines changed: 6 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -47,27 +47,19 @@ void RigidBodyExporter_VTK::writeRigidBodies(const unsigned int frame)
4747
Simulation* sim = Simulation::getCurrent();
4848
const unsigned int nBoundaryModels = sim->numberOfBoundaryModels();
4949

50-
// check if we have a static model
51-
bool isStatic = true;
52-
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
53-
{
54-
BoundaryModel* bm = sim->getBoundaryModel(i);
55-
if (bm->getRigidBodyObject()->isDynamic() || bm->getRigidBodyObject()->isAnimated())
56-
{
57-
isStatic = false;
58-
break;
59-
}
60-
}
61-
6250
#ifdef USE_DOUBLE
6351
const char* real_str = " double\n";
6452
#else
6553
const char* real_str = " float\n";
6654
#endif
6755

68-
if (m_isFirstFrame || !isStatic)
56+
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
6957
{
70-
for (unsigned int i = 0; i < sim->numberOfBoundaryModels(); i++)
58+
BoundaryModel* bm = sim->getBoundaryModel(i);
59+
60+
// If we have a static body, write the data only for the first frame.
61+
// Otherwise each frame is exported.
62+
if (m_isFirstFrame || bm->getRigidBodyObject()->isDynamic() || bm->getRigidBodyObject()->isAnimated())
7163
{
7264
std::string fileName = "rb_data_";
7365
fileName = fileName + std::to_string(i) + "_" + std::to_string(frame) + ".vtk";
@@ -87,7 +79,6 @@ void RigidBodyExporter_VTK::writeRigidBodies(const unsigned int frame)
8779
outfile << "BINARY\n";
8880
outfile << "DATASET UNSTRUCTURED_GRID\n";
8981

90-
BoundaryModel* bm = sim->getBoundaryModel(i);
9182
const std::vector<Vector3r>& vertices = bm->getRigidBodyObject()->getVertices();
9283
const std::vector<unsigned int>& faces = bm->getRigidBodyObject()->getFaces();
9384
int n_vertices = (int)vertices.size();

‎Simulator/GUI/imgui/Simulator_GUI_imgui.cpp‎

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -646,7 +646,8 @@ void Simulator_GUI_imgui::renderBoundary()
646646
if ((renderWalls == 1) || (!scene.boundaryModels[body]->isWall))
647647
{
648648
BoundaryModel_Akinci2012 *bm = static_cast<BoundaryModel_Akinci2012*>(sim->getBoundaryModel(body));
649-
Simulator_OpenGL::renderBoundaryParticles(bm, scene.boundaryModels[body]->color.data(), base->getRenderMinValue(0), base->getRenderMaxValue(0));
649+
Eigen::Matrix<float, 4, 1, Eigen::DontAlign> col = scene.boundaryModels[body]->color.cast<float>();
650+
Simulator_OpenGL::renderBoundaryParticles(bm, col.data(), base->getRenderMinValue(0), base->getRenderMaxValue(0));
650651
}
651652
}
652653
}
@@ -657,7 +658,8 @@ void Simulator_GUI_imgui::renderBoundary()
657658
if ((renderWalls == 3) || (!scene.boundaryModels[body]->isWall))
658659
{
659660
BoundaryModel *bm = sim->getBoundaryModel(body);
660-
Simulator_OpenGL::renderBoundary(bm, scene.boundaryModels[body]->color.data());
661+
Eigen::Matrix<float, 4, 1, Eigen::DontAlign> col = scene.boundaryModels[body]->color.cast<float>();
662+
Simulator_OpenGL::renderBoundary(bm, col.data());
661663
}
662664
}
663665
}

0 commit comments

Comments
 (0)