# _swx_

Member

166

1138 Excellent

• Rank
Member

• Interests
Programming
1. ## rigid body fixed constraint

I have implemented something similar to http://matthias-mueller-fischer.ch/publications/fractureSG2013.pdf In that paper they write "Also, in nature, objects do not necessarily fracture at the impact location only. Both problems can be solved by performing a simplified stress analysis on the connectivity graph of a compound". It also seems that Nvidia Blast does something very similar (in Blast\sdk\extensions\stress\source\NvBlastExtStressSolver.cpp). My first approach was to voxelize the compounds and then compute deformation using FEM, but it takes forever for an iterative solver to converge when dealing with very stiff materials.
2. ## rigid body fixed constraint

I don't know whether the orientation constraint derivation is correct or not, but it looks like the same Jacobian is used by other people as well, e.g http://www.randygaul.net/2013/11/10/orientation-constraint-derivation/ Using the updated velocities when calculating the rotation impulse greatly improves convergence, so I think I can use it for my intended purpose. I don't think that I will need any stabilization since I will be using it for a simplified structural integrity solver where nothing actually moves.
3. ## rigid body fixed constraint

I'm trying to implement a constraint that fixes the relative orientation and position of two rigid bodies. My implementation is based on this document http://danielchappuis.ch/download/ConstraintsDerivationRigidBody3D.pdf The resulting velocity adjustments do not satisfy the constraint and I don't know that I'm doing wrong. #include <Eigen/Eigen> using namespace Eigen; typedef double Real; typedef Eigen::Matrix<Real, 3, 3> Matrix3r; typedef Eigen::Matrix<Real, 3, 1> Vector3r; typedef Eigen::Quaternion<Real> Quaternionr; typedef Eigen::AlignedBox<Real, 3> AlignedBox3r; template <typename Scalar> inline Eigen::Matrix<Scalar, 3, 3> crossProductMatrix(const Eigen::Matrix<Scalar, 3, 1> &v) { Eigen::Matrix<Scalar, 3, 3> v_hat; v_hat << 0, -v(2), v(1), v(2), 0, -v(0), -v(1), v(0), 0; return v_hat; } struct RigidBody { Matrix3r R; // world-space rotation of the rigid body Vector3r x; // world-space position of the rigid body Vector3r v; // linear velocity Vector3r w; // angular velocity Matrix3r Iinv; // local-space inertia tensor Real invmass; }; Matrix3r compute_inertia_tensor(const AlignedBox3r& box, Real m) { const Real a = box.sizes().x(); const Real b = box.sizes().y(); const Real c = box.sizes().z(); return (Real(1.0 / 3.0) * m * Vector3r { b * b + c * c, a * a + c * c, a * a + b * b }).asDiagonal(); } void computeMatrixK( const Vector3r &connector, const Real invMass, const Vector3r &x, const Matrix3r &inertiaInverseW, Matrix3r &K) { const Vector3r r = connector - x; K.setZero(); K += invMass * Matrix3r::Identity(); K += crossProductMatrix(r) * inertiaInverseW * crossProductMatrix(r).transpose(); } bool init_RigidBodyFixedConstraint( const Real invMass0, // inverse mass is zero if body is static const Vector3r &x0, // center of mass of body 0 const Vector3r &v0, // velocity of body 0 const Matrix3r &inertiaInverseW0, // inverse inertia tensor (world space) of body 0 const Vector3r &omega0, // angular velocity of body 0 const Real invMass1, // inverse mass is zero if body is static const Vector3r &x1, // center of mass of body 1 const Vector3r &v1, // velocity of body 1 const Matrix3r &inertiaInverseW1, // inverse inertia tensor (world space) of body 1 const Vector3r &omega1, // angular velocity of body 1 const Vector3r &connector, // world-space coordinates of anchor Eigen::Matrix<Real, 3, 7> &constraintInfo) { Matrix3r K1, K2; computeMatrixK(connector, invMass0, x0, inertiaInverseW0, K1); computeMatrixK(connector, invMass1, x1, inertiaInverseW1, K2); const Matrix3r inv_K_trans = (K1 + K2).inverse(); const Matrix3r inv_K_rot = (inertiaInverseW0 + inertiaInverseW1).inverse(); constraintInfo.block<3, 1>(0, 0) = connector; constraintInfo.block<3, 3>(0, 1) = inv_K_trans; constraintInfo.block<3, 3>(0, 4) = inv_K_rot; return true; } bool velocitySolve_RigidBodyFixedConstraint( const Real invMass0, // inverse mass is zero if body is static const Vector3r &x0, // center of mass of body 0 const Vector3r &v0, // velocity of body 0 const Matrix3r &inertiaInverseW0, // inverse inertia tensor (world space) of body 0 const Vector3r &omega0, // angular velocity of body 0 const Real invMass1, // inverse mass is zero if body is static const Vector3r &x1, // center of mass of body 1 const Vector3r &v1, // velocity of body 1 const Matrix3r &inertiaInverseW1, // inverse inertia tensor (world space) of body 1 const Vector3r &omega1, // angular velocity of body 1 Eigen::Matrix<Real, 3, 7> &constraintInfo, // precomputed constraint info Vector3r &corr_v0, Vector3r &corr_omega0, Vector3r &corr_v1, Vector3r &corr_omega1) { const Vector3r connector = constraintInfo.block<3, 1>(0, 0); const Matrix3r invKt = constraintInfo.block<3, 3>(0, 1); const Matrix3r invKr = constraintInfo.block<3, 3>(0, 4); const Vector3r r0 = connector - x0; const Vector3r r1 = connector - x1; const Vector3r J_trans_v = (v1 + omega1.cross(r1)) - (v0 + omega0.cross(r0)); const Vector3r J_rot_v = omega1 - omega0; const Vector3r lambda_trans = -invKt * J_trans_v; const Vector3r lambda_rot = -invKr * J_rot_v; corr_v0.setZero(); corr_omega0.setZero(); corr_v1.setZero(); corr_omega1.setZero(); if (invMass0 != 0.0f) { corr_v0 -= invMass0 * lambda_trans; corr_omega0 -= inertiaInverseW0 * r0.cross(lambda_trans); corr_omega0 -= inertiaInverseW0 * lambda_rot; } if (invMass1 != 0.0f) { corr_v1 += invMass1 * lambda_trans; corr_omega1 += inertiaInverseW1 * r1.cross(lambda_trans); corr_omega1 += inertiaInverseW1 * lambda_rot; } return true; } Real evalError( const Vector3r &x0, // center of mass of body 0 const Vector3r &v0, // velocity of body 0 const Vector3r &omega0, // angular velocity of body 0 const Vector3r &x1, // center of mass of body 1 const Vector3r &v1, // velocity of body 1 const Vector3r &omega1, // angular velocity of body 1 const Vector3r &connector) { const Vector3r r0 = connector - x0; const Vector3r r1 = connector - x1; const Vector3r J_trans_v = v1 + omega1.cross(r1) - v0 - omega0.cross(r0); const Vector3r J_rot_v = -omega0 + omega1; return J_trans_v.norm() + J_rot_v.norm(); } int main(int argc, char ** argv) { const Real extents = 2.0f; // each cube is 2x2x2 RigidBody A, B; A.x = Vector3r{ 0.0f, 1.0f, 0.0f }; B.x = A.x + Vector3r::UnitY() * (1.0f * extents); const AlignedBox3r box{ Vector3r{ -1.0f, -1.0f, -1.0f }, Vector3r{ +1.0f, +1.0f, +1.0f } }; const Real massA = 20.0f; const Real massB = 40.0f; A.Iinv = compute_inertia_tensor(box, massA).inverse(); B.Iinv = compute_inertia_tensor(box, massB).inverse(); A.invmass = massA == 0.0 ? 0.0 : 1.0 / massA; B.invmass = massB == 0.0 ? 0.0 : 1.0f / massB; B.v = Vector3r::UnitX() * +1.0f; A.v = Vector3r::UnitY() * -1.0f; //A.w.setZero(); //B.w.setZero(); B.w = { 0.1f, 0.2f, 0.3f }; A.w = { 0.4f, 0.1f, 0.7f }; A.R.setIdentity(); B.R.setIdentity(); // world-space inverse inertia tensors const Matrix3r inertiaInverseW1 = A.R * A.Iinv * A.R.transpose(); const Matrix3r inertiaInverseW2 = B.R * B.Iinv * B.R.transpose(); const Vector3r connector{ 0.0f, 2.0f, 0.0f }; Eigen::Matrix<Real, 3, 7> constraintInfo; init_RigidBodyFixedConstraint( A.invmass, A.x, A.v, inertiaInverseW1, A.w, B.invmass, B.x, B.v, inertiaInverseW2, B.w, connector, constraintInfo); for (int i = 0; i < 1; ++i) { Vector3r corr_v0, corr_omega0; Vector3r corr_v1, corr_omega1; velocitySolve_RigidBodyFixedConstraint( A.invmass, A.x, A.v, inertiaInverseW1, A.w, B.invmass, B.x, B.v, inertiaInverseW2, B.w, constraintInfo, corr_v0, corr_omega0, corr_v1, corr_omega1); A.v += corr_v0; A.w += corr_omega0; B.v += corr_v1; B.w += corr_omega1; } const Real error = evalError(A.x, A.v, A.w, B.x, B.v, B.w, connector); return 0; }
4. ## Getting bounding box of convex polyhedron without the points?

I think the standard approach is to have a really large "world cube" as the initial shape and clip it against the planes. This is the approach in most papers I've read and the one I'm using. I am however using a plane based representation for all the actual operations, and only compute points for rendering since this allows me to have 100% robust CSG operations.

6. ## OpenGL OpenGL supporting libs

One good reason not to use GLM is that Eigen is much more versatile and suitable for more than simple graphics calculations :) GLM doesn't have any of the advanced features that Eigen has and mixing different math libraries is messy as hell.
7. ## Baking high-poly to low-poly

www.cs.unc.edu/~olano/papers/SuccessiveMapping.IJCGA.pdf?   "In this paper we present a new simplification algorithm, which computes a piece-wise linear mapping between the original surface and the simplified surface"   Some images from an old implementation I did: Left: High poly Middle: Low poly Right: High poly geometry mapped onto low-poly mesh   Any property can be transferred from the high poly version.
8. ## sugar consumption

I'm on a ketogenic diet, so I eat almost no carbohydrates at all. No negative side effects as far as I can tell :)
9. ## When to Use Interpolation in Animations

Who exactly decides what's allowed, and how come Bitsquid is excempt from those rules?
10. ## When to Use Interpolation in Animations

I'm using a system similar to the one described here: http://bitsquid.blogspot.se/2009/11/bitsquid-low-level-animation-system.html http://bitsquid.blogspot.se/2011/10/low-level-animation-part-2.html
11. ## Question on conversion functions

struct Vector2D{ float x, y; }; struct Point2D{ float x, y; operator Vector2D() const { Vector2D r = { x, y }; return r; } }; Point2D operator-(const Point2D &lh, const Point2D &rh) { Point2D r = { lh.x - rh.x, lh.y - rh.y }; return r; } int main() { Point2D x = { 1, 1 }; Point2D y = { 2, 2 }; Vector2D z = x - y; }
12. ## Question on conversion functions

I don't recommend having different types for points and vectors though... it will make your code much messier with lots of conversions.
13. ## Question on conversion functions

Avoid conversion operators and find alternative solutions like: Vector2D operator-(const Position2D &lh, const Position2D &rh) { return Vector2D(lh.x - rh.x, lh.y - rh.y); }
14. ## Are there any standards to organizing game file and data chunks?

I use json for pretty much all my resource files except textures, models and other very large resource types. This allows my resource compiler to automatically find dependencies between the resources.
15. ## Looking for fast c++ math library

http://macton.smugmug.com/gallery/8936708_T6zQX#!i=593426709&k=BrHWXdJ is a good read if you want to write really efficient code, even if the author is a bit smug ;)