diff --git a/src/Makefile b/src/Makefile index 3cfa00b..c1ae188 100644 --- a/src/Makefile +++ b/src/Makefile @@ -1,5 +1,5 @@ CXX=llvm-g++ -COMPILE=-c -std=c++17 \ +COMPILE=-g -c -std=c++17 \ -I/usr/local/Cellar/glew/2.1.0/include \ -I/usr/local/Cellar/glfw/3.3/include \ -I/usr/local/Cellar/eigen/3.3.7/include/eigen3 \ diff --git a/src/constants.h b/src/constants.h index 735f1b0..67fd5ae 100755 --- a/src/constants.h +++ b/src/constants.h @@ -40,6 +40,7 @@ const static glm::vec2 worldOrigin(0.f, 0.f); const static float sdfCellSize = H_INV; const static int sdfWidth = X_GRID; const static int sdfHeight = Y_GRID; +const static float NARROW_BAND = 1.f; /* ----- TRANSFER ----- */ #if INTERPOLATION == 1 diff --git a/src/main b/src/main index e6d2930..903983e 100755 Binary files a/src/main and b/src/main differ diff --git a/src/main.cpp b/src/main.cpp index 9935fe6..0340279 100755 --- a/src/main.cpp +++ b/src/main.cpp @@ -101,7 +101,7 @@ int main(int argc, char **argv) { // AddParticles(); // P2G, compute grid forces, etc. - // Update(); + Update(); // Display frame at desired rate if (t_count % (int)(DT_render / DT) == 0) { @@ -112,43 +112,42 @@ int main(int argc, char **argv) { glfwSwapBuffers(window); - // #if RECORD_VIDEO - // std::string dir = "../result/sim"; - // - // // zero padding - // std::string num = std::to_string(frameNumber); - // num = std::string(4 - num.length(), '0') + num; - // std::string output = dir + num + ".png"; - // - // // use opencv to save a frame - // // the frame size becomes twice from frame 2, - // // thus multiplying 2 is needed. - // // don't know why - // cv::Mat outImg(Y_WINDOW * 2, X_WINDOW * 2, CV_8UC3); - // glPixelStorei(GL_PACK_ALIGNMENT, 4); - // glReadPixels(0, 0, X_WINDOW * 2, Y_WINDOW * 2, GL_BGR, - // GL_UNSIGNED_BYTE, - // (GLvoid *)outImg.data); - // cv::flip(outImg, outImg, 0); - // cv::imwrite(output, outImg); - // - // std::cout << "Frame " << frameNumber << " saved." << '\n'; - // - // frameNumber++; - // #endif +#if RECORD_VIDEO + std::string dir = "../result/sim"; + + // zero padding + std::string num = std::to_string(frameNumber); + num = std::string(4 - num.length(), '0') + num; + std::string output = dir + num + ".png"; + + // use opencv to save a frame + // the frame size becomes twice from frame 2, + // thus multiplying 2 is needed. + // don't know why + cv::Mat outImg(Y_WINDOW * 2, X_WINDOW * 2, CV_8UC3); + glPixelStorei(GL_PACK_ALIGNMENT, 4); + glReadPixels(0, 0, X_WINDOW * 2, Y_WINDOW * 2, GL_BGR, GL_UNSIGNED_BYTE, + (GLvoid *)outImg.data); + cv::flip(outImg, outImg, 0); + cv::imwrite(output, outImg); + + std::cout << "Frame " << frameNumber << " saved." << '\n'; + + frameNumber++; +#endif glfwPollEvents(); } // end outer if // Don't forget to reset grid every frame !! - // Simulation->ResetGrid(); + Simulation->ResetGrid(); t_count++; } // end while - // #if RECORD_VIDEO - // images2video(); - // #endif +#if RECORD_VIDEO + images2video(); +#endif glfwTerminate(); #endif diff --git a/src/node.h b/src/node.h index d6422bf..f9d68f1 100755 --- a/src/node.h +++ b/src/node.h @@ -11,7 +11,7 @@ class Node { /* Data */ float Mi; // Node mass - Vector2f Xi; // Node position + Vector2f Xi; // Node position (index, also world space position) Vector2f Vi; // Node velocity, before update and after force Vector2f Vi_col; // Node velocity, after collision Vector2f Vi_fri; // Node velocity, after friction diff --git a/src/particle.h b/src/particle.h index 401ff82..ca21321 100755 --- a/src/particle.h +++ b/src/particle.h @@ -191,7 +191,7 @@ class Snow : public Particle { float VOL = 2 * PI * R_BALL * R_BALL / static_cast(NP); float MASS = VOL * RHO_snow / 100.0f; - Vector2f v = Vector2f(0.f, -9.8f); // Initial velocity + Vector2f v = Vector2f(0.f, -9.8f * 10.f); // Initial velocity Matrix2f a = Matrix2f(0); // get position from image diff --git a/src/result.mp4 b/src/result.mp4 index 16f5692..3bcf68d 100644 Binary files a/src/result.mp4 and b/src/result.mp4 differ diff --git a/src/sdf.h b/src/sdf.h index c086889..c958548 100644 --- a/src/sdf.h +++ b/src/sdf.h @@ -59,7 +59,7 @@ class Polygon { polys.push_back(square); polys.push_back(object2); - polys.push_back(object3); + // polys.push_back(object3); return polys; } diff --git a/src/solver.cpp b/src/solver.cpp index 5293763..bdf63ec 100755 --- a/src/solver.cpp +++ b/src/solver.cpp @@ -125,6 +125,7 @@ float Solver::getDistance(glm::vec2 p) { return 9999.f; } else { int idx = idx_x + idx_y * (sdfWidth - 1); + // std::cout << idx << '\n'; return nodes[idx].sdfDist; } } @@ -159,6 +160,84 @@ glm::vec2 Solver::getGradient(glm::vec2 p) { return gd; } +void Solver::applySdfCollision(Node &node) { + // node world space position + glm::vec2 pos(node.Xi[0], node.Xi[1]); + glm::vec2 v(node.Vi[0], node.Vi[1]); + + // std::cout << glm::to_string(pos) << '\n'; + + // distance + float dist = getDistance(pos); + glm::vec2 n; + // Polygon *co; // collision object + bool isCollisionOn = false; + + glm::vec2 vco(0.f, 0.f); + + // 基于 sdf 的碰撞检测 + // 和边界处的碰撞检测 + // 唯一不同的是法向量 n 的计算方式 + // 而碰撞后的速度,采用同样��处理 + // for sdf collision detection + // narrow band threshold + if (glm::abs(dist) < NARROW_BAND) { + n = -getGradient(pos); + // co = getPolygon(pos); // get collision object + isCollisionOn = true; + + // object velocity + // linear + // vec2 vlin = co->v; + glm::vec2 vlin(0.f, 0.f); + + // rotational + // vec2 center = (co->lb + co->rt) * 0.5f; + // vec2 r = pos - center; + // + // vec2 vrot; + // vrot.x = -r.y / length(r); + // vrot.y = r.x / length(r); + glm::vec2 vrot(0.f, 0.f); + + // for visualization convenience + float scale = 0.1f; + + // vrot *= co->omega * length(r) * scale; + + vco = vlin + vrot; + + // relative velocity + glm::vec2 vrel = v - vco; + + float vnLength = glm::dot(n, vrel); + // std::cout << vnLength << ", " << isCollisionOn << '\n'; + + // if not separating + if (vnLength < 0.f && isCollisionOn) { + glm::vec2 vt = vrel - vnLength * n; + float mu = 0.55f; + + if (length(vt) <= -mu * vnLength) { + vrel = glm::vec2(0.f, 0.f); + } else { + vrel = vt + mu * vnLength * normalize(vt); + } + + // back to particle velocity + v = vrel + vco; + + // std::cout << v.x << ", " << v.y << '\n'; + + node.Vi[0] = v.x; + node.Vi[1] = v.y; + + // std::cout << "node velocity: " << node.Vi[0] << ", " << node.Vi[1] + // << '\n'; + } + } +} + /* ----------------------------------------------------------------------- | MATERIAL POINT METHOD ALGORITHM | @@ -233,7 +312,9 @@ void Solver::UpdateNodes() { nodes[i].Vi += nodes[i].Fi; // Apply collisions and frictions - nodes[i].NodeCollisions(); + applySdfCollision(nodes[i]); // sdf-based + + nodes[i].NodeCollisions(); // borders #if FRICTION nodes[i].NodeFrictions(); #else diff --git a/src/solver.h b/src/solver.h index 9a2b86d..de68c4a 100755 --- a/src/solver.h +++ b/src/solver.h @@ -43,6 +43,7 @@ class Solver { float nearest_distance(glm::vec2, Polygon &); float getDistance(glm::vec2); glm::vec2 getGradient(glm::vec2); + void applySdfCollision(Node &); void Draw(); // Draw particles, border and nodes (if selected) void WriteToFile(