From c67d33c25e0e191f69a3aa5b20aab73b8e2fd8fe Mon Sep 17 00:00:00 2001 From: Benjamin Kraft Date: Sun, 29 Sep 2024 22:02:46 +0200 Subject: [PATCH] soft tetrahedrons --- include/vulkan/vertex.hpp | 2 +- shaders/pbd.comp | 68 +++++++++++++++++++++++++++++++++------ src/simulation.cpp | 4 +-- src/soft_body.cpp | 27 ++++++++++++++-- 4 files changed, 85 insertions(+), 16 deletions(-) diff --git a/include/vulkan/vertex.hpp b/include/vulkan/vertex.hpp index 2df0d1d..21a591d 100644 --- a/include/vulkan/vertex.hpp +++ b/include/vulkan/vertex.hpp @@ -9,7 +9,7 @@ struct Vertex { alignas(16) glm::vec3 color; alignas(16) glm::vec3 normal; alignas(16) glm::vec3 velocity; - glm::vec3 prevPosition; + alignas(16) glm::vec3 prevPosition; float inverseMass; static VkVertexInputBindingDescription getBindingDescription(); diff --git a/shaders/pbd.comp b/shaders/pbd.comp index 1dec462..14791b8 100644 --- a/shaders/pbd.comp +++ b/shaders/pbd.comp @@ -17,7 +17,13 @@ struct Edge { float restLength; }; - +struct Tetrahedron { + uint a; + uint b; + uint c; + uint d; + float restVolume; +}; layout (std430, set = 0, binding = 0) buffer VertexBuffer { @@ -28,6 +34,10 @@ layout (std430, set = 0, binding = 2) buffer EdgeBuffer { Edge edges[]; }; +layout (std430, set = 0, binding = 4) buffer TetrahedronBuffer { + Tetrahedron tetrahedra[]; +}; + layout (std140, set = 0, binding = 5) uniform Sizes { uint vertexCount; uint faceCount; @@ -57,30 +67,68 @@ void preSolve(uint vID){ } void solveEdge(uint eID){ - Edge e = edges[eID]; + Edge edge = edges[eID]; - Vertex v1 = vertices[e.a]; - Vertex v2 = vertices[e.b]; + Vertex v1 = vertices[edge.a]; + Vertex v2 = vertices[edge.b]; vec3 diff = v1.position - v2.position; float currentLength = length(diff); - float alpha = 0.0001 / (dt * dt); + float alpha = 0.3 / dt / dt; - float s = (currentLength - e.restLength) / (1 + 1 + alpha); + float s = -(currentLength - edge.restLength) / (v1.w + v2.w + alpha); vec3 g1 = normalize(diff); vec3 g2 = -g1; - vec3 delta1 = -s * 1 * g1; - vec3 delta2 = -s * 1 * g2; + vec3 delta1 = s * v1.w * g1; + vec3 delta2 = s * v2.w * g2; - vertices[e.a].position += delta1; - vertices[e.b].position += delta2; + vertices[edge.a].position += delta1; + vertices[edge.b].position += delta2; } void solveTetrahedron(uint tetID){ + Tetrahedron tetrahedron = tetrahedra[tetID]; + + Vertex va = vertices[tetrahedron.a]; + Vertex vb = vertices[tetrahedron.b]; + Vertex vc = vertices[tetrahedron.c]; + Vertex vd = vertices[tetrahedron.d]; + + vec3 a = va.position; + vec3 b = vb.position; + vec3 c = vc.position; + vec3 d = vd.position; + + float volumeError = dot(d - a, cross(b - a, c - a)) / 6 - tetrahedron.restVolume; + + vec3 ga = cross(d - b, c - b) / 6; + vec3 gb = cross(c - a, d - a) / 6; + vec3 gc = cross(d - a, b - a) / 6; + vec3 gd = cross(b - a, c - a) / 6; + + float w = + va.w * dot(ga, ga) + + vb.w * dot(gb, gb) + + vc.w * dot(gc, gc) + + vd.w * dot(gd, gd); + + if (w == 0) return; + + float alpha = 0 / dt / dt; + float s = -volumeError / (w + alpha); + + vec3 addA = s * ga * va.w; + vec3 addB = s * gb * vb.w; + vec3 addC = s * gc * vc.w; + vec3 addD = s * gd * vd.w; + vertices[tetrahedron.a].position += addA; + vertices[tetrahedron.b].position += addB; + vertices[tetrahedron.c].position += addC; + vertices[tetrahedron.d].position += addD; } void postSolve(uint vID){ diff --git a/src/simulation.cpp b/src/simulation.cpp index 8df125e..259e82b 100644 --- a/src/simulation.cpp +++ b/src/simulation.cpp @@ -23,7 +23,7 @@ Simulation::Simulation() { properties.gravity = {0, -9.81, 0}; properties.k = 10; - properties.dt = 1.f / 144.f / static_cast(properties.k); + properties.dt = 1.f / 60.f / static_cast(properties.k); propertiesBuffer = make_unique( sizeof(Properties), &properties, sizeof(properties), @@ -148,8 +148,6 @@ void Simulation::createMeshBuffers() { : Buffer(size, data, size, VK_BUFFER_USAGE_STORAGE_BUFFER_BIT | additionalUsageFlags, VMA_MEMORY_USAGE_AUTO_PREFER_DEVICE, 0) {} }; - vertices[0].position += glm::vec3(0, 0.2, 0); - vertexBuffer = make_unique(vertices.data(), vertices.size() * sizeof(Vertex), VK_BUFFER_USAGE_VERTEX_BUFFER_BIT); faceBuffer = make_unique(faces.data(), faces.size() * sizeof(Face), VK_BUFFER_USAGE_INDEX_BUFFER_BIT); edgeBuffer = make_unique(constraintData.edges.data(), constraintData.edges.size() * sizeof(Edge)); diff --git a/src/soft_body.cpp b/src/soft_body.cpp index 86150a6..2b75d38 100644 --- a/src/soft_body.cpp +++ b/src/soft_body.cpp @@ -9,7 +9,7 @@ SoftBody::SoftBody(Mesh* mesh, float compliance) : compliance(compliance) { tetgenbehavior behavior; - behavior.parse_commandline(std::string("pYa0.001Qfez").data()); + behavior.parse_commandline(std::string("pYQfeza0.1").data()); tetgenio in; in.numberofpoints = static_cast(mesh->vertices.size()); @@ -79,13 +79,36 @@ SoftBody::SoftBody(Mesh* mesh, float compliance) : compliance(compliance) { constraintData.edges.emplace_back(Edge(a, b, length)); } + float density = 1; + vector masses(vertices.size()); + constraintData.tetrahedra.reserve(out.numberoftetrahedra); for (size_t i = 0; i < out.numberoftetrahedra; i++){ uint32_t a = out.tetrahedronlist[i * 4 + 0]; uint32_t b = out.tetrahedronlist[i * 4 + 1]; uint32_t c = out.tetrahedronlist[i * 4 + 2]; uint32_t d = out.tetrahedronlist[i * 4 + 3]; - constraintData.tetrahedra.emplace_back(Tetrahedron(a, b, c, d, 0)); + + glm::vec3 p1 = vertices[a].position; + glm::vec3 p2 = vertices[b].position; + glm::vec3 p3 = vertices[c].position; + glm::vec3 p4 = vertices[d].position; + + float volume = glm::dot(p4 - p1, glm::cross(p2 - p1, p3 - p1)) / 6; + assert(volume > 0); + + float pointMass = volume * density / 4; + + masses[a] += pointMass; + masses[b] += pointMass; + masses[c] += pointMass; + masses[d] += pointMass; + + constraintData.tetrahedra.emplace_back(Tetrahedron(a, b, c, d, volume)); + } + + for (size_t i = 0; i < masses.size(); i++){ + vertices[i].inverseMass = 1.0f / masses[i]; } splitConstraints();