adding in verlet

index.html

@@ -105,6 +105,14 @@
             gl_FragColor = vec4(position-u_center, 0.0);
         }
     </script>
+    <script id="copyTexture" type="x-shader/x-fragment">
+        precision mediump float;
+        uniform sampler2D u_orig;
+        uniform vec2 u_textureDim;
+        void main(){
+            gl_FragColor = texture2D(u_orig, gl_FragCoord.xy/u_textureDim);
+        }
+    </script>
 
     <script id="positionCalcShader" type="x-shader/x-fragment">
         precision mediump float;
@@ -132,6 +140,30 @@
         }
     </script>
 
+    <script id="velocityCalcVerletShader" type="x-shader/x-fragment">
+        precision mediump float;
+        uniform vec2 u_textureDim;
+        uniform float u_dt;
+        uniform sampler2D u_position;
+        uniform sampler2D u_lastPosition;
+        uniform sampler2D u_mass;
+
+        void main(){
+            vec2 fragCoord = gl_FragCoord.xy;
+            vec2 scaledFragCoord = fragCoord/u_textureDim;
+
+            float isFixed = texture2D(u_mass, scaledFragCoord).y;
+            if (isFixed == 1.0){
+                gl_FragColor = vec4(0.0);
+                return;
+            }
+
+            vec3 position = texture2D(u_position, scaledFragCoord).xyz;
+            vec3 lastPosition = texture2D(u_position, scaledFragCoord).xyz;
+            gl_FragColor = vec4((position-lastPosition)/u_dt,0.0);
+        }
+    </script>
+
     <script id="velocityCalcShader" type="x-shader/x-fragment">
         precision mediump float;
         uniform vec2 u_textureDim;
@@ -336,6 +368,212 @@
         }
     </script>
 
+    <script id="positionCalcVerletShader" type="x-shader/x-fragment">
+        precision mediump float;
+        uniform vec2 u_textureDim;
+        uniform vec2 u_textureDimEdges;
+        uniform vec2 u_textureDimFaces;
+        uniform vec2 u_textureDimCreases;
+        uniform vec2 u_textureDimNodeCreases;
+        uniform vec2 u_textureDimNodeFaces;
+        uniform float u_creasePercent;
+        uniform float u_dt;
+        uniform float u_axialStiffness;
+        uniform sampler2D u_lastPosition;
+        uniform sampler2D u_lastLastPosition;
+        uniform sampler2D u_lastVelocity;
+        uniform sampler2D u_originalPosition;
+        uniform sampler2D u_externalForces;
+        uniform sampler2D u_mass;
+        uniform sampler2D u_meta;//[beamsIndex, numBeam, nodeCreaseMetaIndex, numCreases]
+        uniform sampler2D u_beamMeta;//[k, d, length, otherNodeIndex]
+        uniform sampler2D u_creaseMeta;//[k, d, targetTheta]
+        uniform sampler2D u_nodeCreaseMeta;//[creaseIndex, nodeIndex, -, -]
+        uniform sampler2D u_normals;
+        uniform sampler2D u_theta;//[theta, z, normal1Index, normal2Index]
+        uniform sampler2D u_creaseGeo;//[h1, h2, coef1, coef2]
+        uniform sampler2D u_meta2;//[nodesFaceIndex, numFaces]
+        uniform sampler2D u_nodeFaceMeta;//[faceIndex, a, b, c]
+        uniform sampler2D u_nominalTriangles;//[angleA, angleB, angleC]
+
+        vec4 getFromArray(float index1D, vec2 dimensions, sampler2D tex){
+            vec2 index = vec2(mod(index1D, dimensions.x)+0.5, floor(index1D/dimensions.x)+0.5);
+            vec2 scaledIndex = index/dimensions;
+            return texture2D(tex, scaledIndex);
+        }
+
+        vec3 getPosition(float index1D){
+            vec2 index = vec2(mod(index1D, u_textureDim.x)+0.5, floor(index1D/u_textureDim.x)+0.5);
+            vec2 scaledIndex = index/u_textureDim;
+            return texture2D(u_lastPosition, scaledIndex).xyz + texture2D(u_originalPosition, scaledIndex).xyz;
+        }
+
+        void main(){
+            vec2 fragCoord = gl_FragCoord.xy;
+            vec2 scaledFragCoord = fragCoord/u_textureDim;
+
+            vec2 mass = texture2D(u_mass, scaledFragCoord).xy;
+            if (mass[1] == 1.0){//fixed
+                gl_FragColor = vec4(0.0);
+                return;
+            }
+            vec3 force = texture2D(u_externalForces, scaledFragCoord).xyz;
+            vec3 lastPosition = texture2D(u_lastPosition, scaledFragCoord).xyz;
+            vec3 lastLastPosition = texture2D(u_lastLastPosition, scaledFragCoord).xyz;
+            vec3 lastVelocity = texture2D(u_lastVelocity, scaledFragCoord).xyz;
+            vec3 originalPosition = texture2D(u_originalPosition, scaledFragCoord).xyz;
+
+            vec4 neighborIndices = texture2D(u_meta, scaledFragCoord);
+            vec4 meta = texture2D(u_meta, scaledFragCoord);
+            vec2 meta2 = texture2D(u_meta2, scaledFragCoord).xy;
+
+            float nodeError = 0.0;
+
+            for (int j=0;j<100;j++){//for all beams (up to 100, had to put a const int in here)
+                if (j >= int(meta[1])) break;
+
+                vec4 beamMeta = getFromArray(meta[0]+float(j), u_textureDimEdges, u_beamMeta);
+
+                float neighborIndex1D = beamMeta[3];
+                vec2 neighborIndex = vec2(mod(neighborIndex1D, u_textureDim.x)+0.5, floor(neighborIndex1D/u_textureDim.x)+0.5);
+                vec2 scaledNeighborIndex = neighborIndex/u_textureDim;
+                vec3 neighborLastPosition = texture2D(u_lastPosition, scaledNeighborIndex).xyz;
+                vec3 neighborLastVelocity = texture2D(u_lastVelocity, scaledNeighborIndex).xyz;
+                vec3 neighborOriginalPosition = texture2D(u_originalPosition, scaledNeighborIndex).xyz;
+
+                vec3 nominalDist = neighborOriginalPosition-originalPosition;
+                vec3 deltaP = neighborLastPosition-lastPosition+nominalDist;
+                float deltaPLength = length(deltaP);
+                deltaP -= deltaP*(beamMeta[2]/deltaPLength);
+                nodeError += abs(deltaPLength/length(nominalDist) - 1.0);
+                vec3 deltaV = neighborLastVelocity-lastVelocity;
+
+                vec3 _force = deltaP*beamMeta[0] + deltaV*beamMeta[1];
+                force += _force;
+            }
+            nodeError /= meta[1];
+
+            for (int j=0;j<100;j++){//for all creases (up to 100, had to put a const int in here)
+                if (j >= int(meta[3])) break;
+
+                vec4 nodeCreaseMeta = getFromArray(meta[2]+float(j), u_textureDimNodeCreases, u_nodeCreaseMeta);
+
+                float creaseIndex1D = nodeCreaseMeta[0];
+                vec2 creaseIndex = vec2(mod(creaseIndex1D, u_textureDimCreases.x)+0.5, floor(creaseIndex1D/u_textureDimCreases.x)+0.5);
+                vec2 scaledCreaseIndex = creaseIndex/u_textureDimCreases;
+
+                vec4 thetas = texture2D(u_theta, scaledCreaseIndex);
+                vec3 creaseMeta = texture2D(u_creaseMeta, scaledCreaseIndex).xyz;//[k, d, targetTheta]
+                vec4 creaseGeo = texture2D(u_creaseGeo, scaledCreaseIndex);//[h1, h2, coef1, coef2]
+                if (creaseGeo[0]< 0.0) continue;//crease disabled bc it has collapsed too much
+
+                float targetTheta = creaseMeta[2] * u_creasePercent;
+                float angForce = creaseMeta[0]*(targetTheta-thetas[0]);// + creaseMeta[1]*thetas[1];
+
+                float nodeNum = nodeCreaseMeta[1];//1, 2, 3, 4
+
+                if (nodeNum > 2.0){//crease reaction, node is on a crease
+
+                    //node #1
+                    vec3 normal1 = getFromArray(thetas[2], u_textureDimFaces, u_normals).xyz;
+
+                    //node #2
+                    vec3 normal2 = getFromArray(thetas[3], u_textureDimFaces, u_normals).xyz;
+
+                    float coef1 = creaseGeo[2];
+                    float coef2 = creaseGeo[3];
+
+                    if (nodeNum == 3.0){
+                        coef1 = 1.0-coef1;
+                        coef2 = 1.0-coef2;
+                    }
+
+                    vec3 _force = -angForce*(coef1/creaseGeo[0]*normal1 + coef2/creaseGeo[1]*normal2);
+                    force += _force;
+
+                } else {
+
+                    float normalIndex1D = thetas[2];//node #1
+                    float momentArm = creaseGeo[0];//node #1
+                    if (nodeNum == 2.0) {
+                        normalIndex1D = thetas[3];//node #2
+                        momentArm = creaseGeo[1];//node #2
+                    }
+
+                    vec3 normal = getFromArray(normalIndex1D, u_textureDimFaces, u_normals).xyz;
+
+                    vec3 _force = angForce/momentArm*normal;
+                    force += _force;
+                }
+            }
+
+            for (int j=0;j<100;j++){//for all faces (up to 100, had to put a const int in here)
+                if (j >= int(meta2[1])) break;
+
+                vec4 faceMeta = getFromArray(meta2[0]+float(j), u_textureDimNodeFaces, u_nodeFaceMeta);//[face index, a, b, c]
+                vec3 nominalAngles = getFromArray(faceMeta[0], u_textureDimFaces, u_nominalTriangles).xyz;//[angA, angB, angC]
+
+                int faceIndex = 0;
+                if (faceMeta[2] < 0.0) faceIndex = 1;
+                if (faceMeta[3] < 0.0) faceIndex = 2;
+
+                //get node positions
+                vec3 a = faceIndex == 0 ? lastPosition+originalPosition : getPosition(faceMeta[1]);
+                vec3 b = faceIndex == 1 ? lastPosition+originalPosition : getPosition(faceMeta[2]);
+                vec3 c = faceIndex == 2 ? lastPosition+originalPosition : getPosition(faceMeta[3]);
+
+                //calc angles
+                vec3 ab = b-a;
+                vec3 ac = c-a;
+                vec3 bc = c-b;
+
+                float lengthAB = length(ab);
+                float lengthAC = length(ac);
+                float lengthBC = length(bc);
+
+                float tol = 0.0000001;
+                if (abs(lengthAB) < tol || abs(lengthBC) < tol || abs(lengthAC) < tol) continue;
+
+                ab /= lengthAB;
+                ac /= lengthAC;
+                bc /= lengthBC;
+
+                vec3 angles = vec3(acos(dot(ab, ac)),
+                    acos(-1.0*dot(ab, bc)),
+                    acos(dot(ac, bc)));
+                vec3 anglesDiff = nominalAngles-angles;
+
+                vec3 normal = getFromArray(faceMeta[0], u_textureDimFaces, u_normals).xyz;
+
+                //calc forces
+                anglesDiff *= u_axialStiffness/100.0;
+                if (faceIndex == 0){//a
+                    vec3 normalCrossAC = cross(normal, ac)/lengthAC;
+                    vec3 normalCrossAB = cross(normal, ab)/lengthAB;
+                    force -= anglesDiff[0]*(normalCrossAC - normalCrossAB);
+                    force -= anglesDiff[1]*normalCrossAB;
+                    force += anglesDiff[2]*normalCrossAC;
+                } else if (faceIndex == 1){
+                    vec3 normalCrossAB = cross(normal, ab)/lengthAB;
+                    vec3 normalCrossBC = cross(normal, bc)/lengthBC;
+                    force -= anglesDiff[0]*normalCrossAB;
+                    force += anglesDiff[1]*(normalCrossAB + normalCrossBC);
+                    force -= anglesDiff[2]*normalCrossBC;
+                } else if (faceIndex == 2){
+                    vec3 normalCrossAC = cross(normal, ac)/lengthAC;
+                    vec3 normalCrossBC = cross(normal, bc)/lengthBC;
+                    force += anglesDiff[0]*normalCrossAC;
+                    force -= anglesDiff[1]*normalCrossBC;
+                    force += anglesDiff[2]*(normalCrossBC - normalCrossAC);
+                }
+
+            }
+
+            vec3 nextPosition = force*u_dt*u_dt/mass[0] + 2.0*lastPosition - lastLastPosition;
+            gl_FragColor = vec4(nextPosition,nodeError);//position.a has error info
+        }
+    </script>
+
     <script id="thetaCalcShader" type="x-shader/x-fragment">
         #define TWO_PI 6.283185307179586476925286766559
         precision mediump float;
@@ -824,6 +1062,17 @@
     <div class="fullWidth">
         <a href="#" class="seeMore closed" data-id="simulationSettings"><span class="fui-triangle-down"></span><b>Simulation Settings:</b></a><a class="about" href="#" id="aboutStiffness"><span class="fui-question-circle"></span></a><br/>
         <div id="simulationSettings" class="hide">
+            <div class="indent">Integration:
+            <div class="indent">
+                <label class="radio">
+                    <input name="integrationType" value="verlet" data-toggle="radio" class="custom-radio" type="radio"><span class="icons"><span class="icon-unchecked"></span><span class="icon-checked"></span></span>
+                    Verlet (second order)
+                </label>
+                <label class="radio">
+                    <input name="integrationType" value="euler" data-toggle="radio" class="custom-radio" type="radio"><span class="icons"><span class="icon-unchecked"></span><span class="icon-checked"></span></span>
+                    Euler (first order)
+                </label>
+            </div></div><br/>
             <div class="sliderInput paddingBottom" id="axialStiffness">
                 <span class="label-slider">Axial Stiffness : </span><div class="flat-slider ui-slider ui-corner-all ui-slider-horizontal ui-widget ui-widget-content"></div>
                 <input value="" placeholder="" class="float form-control" type="text">

js/controls.js

@@ -390,6 +390,10 @@ function initControls(globals){
         $("#aboutVRmodal").modal("show");
     });
 
+    setRadio("integrationType", globals.integrationType, function(val){
+        globals.integrationType = val;
+    });
+
 
     setRadio("simType", globals.simType, function(val){
         globals.simType = val;

js/dynamic/dynamicSolver.js

@@ -16,6 +16,7 @@ function initDynamicSolver(globals){
     var originalPosition;
     var position;
     var lastPosition;
+    var lastLastPosition;//for verlet integration
     var velocity;
     var lastVelocity;
     var externalForces;
@@ -62,6 +63,7 @@ function initDynamicSolver(globals){
     function reset(){
         globals.gpuMath.step("zeroTexture", [], "u_position");
         globals.gpuMath.step("zeroTexture", [], "u_lastPosition");
+        if (globals.integrationType == "verlet") globals.gpuMath.step("zeroTexture", [], "u_lastLastPosition");
         globals.gpuMath.step("zeroTexture", [], "u_velocity");
         globals.gpuMath.step("zeroTexture", [], "u_lastVelocity");
         globals.gpuMath.step("zeroThetaTexture", ["u_lastTheta"], "u_theta");
@@ -112,6 +114,7 @@ function initDynamicSolver(globals){
             globals.gpuMath.setProgram("centerTexture");
             globals.gpuMath.setUniformForProgram("centerTexture", "u_center", [avgPosition.x, avgPosition.y, avgPosition.z], "3f");
             globals.gpuMath.step("centerTexture", ["u_lastPosition"], "u_position");
+            if (globals.integrationType == "verlet") globals.gpuMath.step("copyTexture", ["u_position"], "u_lastLastPosition");
             globals.gpuMath.swapTextures("u_position", "u_lastPosition");
             globals.gpuMath.step("zeroTexture", [], "u_lastVelocity");
             globals.gpuMath.step("zeroTexture", [], "u_velocity");
@@ -142,12 +145,22 @@ function initDynamicSolver(globals){
         //already at textureDimCreasesxtextureDimCreases
         gpuMath.step("updateCreaseGeo", ["u_lastPosition", "u_originalPosition", "u_creaseMeta2"], "u_creaseGeo");
 
-        gpuMath.setProgram("velocityCalc");
-        gpuMath.setSize(textureDim, textureDim);
-        gpuMath.step("velocityCalc", ["u_lastPosition", "u_lastVelocity", "u_originalPosition", "u_externalForces",
-            "u_mass", "u_meta", "u_beamMeta", "u_creaseMeta", "u_nodeCreaseMeta", "u_normals", "u_theta", "u_creaseGeo",
-            "u_meta2", "u_nodeFaceMeta", "u_nominalTriangles"], "u_velocity");
-        gpuMath.step("positionCalc", ["u_velocity", "u_lastPosition", "u_mass"], "u_position");
+        if (globals.integrationType == "verlet"){
+            gpuMath.setProgram("positionCalcVerlet");
+            gpuMath.setSize(textureDim, textureDim);
+            gpuMath.step("positionCalcVerlet", ["u_lastPosition", "u_lastLastPosition", "u_lastVelocity", "u_originalPosition", "u_externalForces",
+                "u_mass", "u_meta", "u_beamMeta", "u_creaseMeta", "u_nodeCreaseMeta", "u_normals", "u_theta", "u_creaseGeo",
+                "u_meta2", "u_nodeFaceMeta", "u_nominalTriangles"], "u_position");
+            gpuMath.step("velocityCalcVerlet", ["u_position", "u_lastPosition", "u_mass"], "u_velocity");
+            gpuMath.swapTextures("u_lastPosition", "u_lastLastPosition");
+        } else {//euler
+            gpuMath.setProgram("velocityCalc");
+            gpuMath.setSize(textureDim, textureDim);
+            gpuMath.step("velocityCalc", ["u_lastPosition", "u_lastVelocity", "u_originalPosition", "u_externalForces",
+                "u_mass", "u_meta", "u_beamMeta", "u_creaseMeta", "u_nodeCreaseMeta", "u_normals", "u_theta", "u_creaseGeo",
+                "u_meta2", "u_nodeFaceMeta", "u_nominalTriangles"], "u_velocity");
+            gpuMath.step("positionCalc", ["u_velocity", "u_lastPosition", "u_mass"], "u_position");
+        }
 
         gpuMath.swapTextures("u_theta", "u_lastTheta");
         gpuMath.swapTextures("u_velocity", "u_lastVelocity");
@@ -261,6 +274,7 @@ function initDynamicSolver(globals){
 
         gpuMath.initTextureFromData("u_position", textureDim, textureDim, "FLOAT", position, true);
         gpuMath.initTextureFromData("u_lastPosition", textureDim, textureDim, "FLOAT", lastPosition, true);
+        gpuMath.initTextureFromData("u_lastLastPosition", textureDim, textureDim, "FLOAT", lastLastPosition, true);
         gpuMath.initTextureFromData("u_velocity", textureDim, textureDim, "FLOAT", velocity, true);
         gpuMath.initTextureFromData("u_lastVelocity", textureDim, textureDim, "FLOAT", lastVelocity, true);
         gpuMath.initTextureFromData("u_theta", textureDimCreases, textureDimCreases, "FLOAT", theta, true);
@@ -269,6 +283,7 @@ function initDynamicSolver(globals){
 
         gpuMath.initFrameBufferForTexture("u_position", true);
         gpuMath.initFrameBufferForTexture("u_lastPosition", true);
+        gpuMath.initFrameBufferForTexture("u_lastLastPosition", true);
         gpuMath.initFrameBufferForTexture("u_velocity", true);
         gpuMath.initFrameBufferForTexture("u_lastVelocity", true);
         gpuMath.initFrameBufferForTexture("u_theta", true);
@@ -292,6 +307,12 @@ function initDynamicSolver(globals){
         gpuMath.setUniformForProgram("positionCalc", "u_mass", 2, "1i");
         gpuMath.setUniformForProgram("positionCalc", "u_textureDim", [textureDim, textureDim], "2f");
 
+        gpuMath.createProgram("velocityCalcVerlet", vertexShader, document.getElementById("velocityCalcVerletShader").text);
+        gpuMath.setUniformForProgram("velocityCalcVerlet", "u_position", 0, "1i");
+        gpuMath.setUniformForProgram("velocityCalcVerlet", "u_lastPosition", 1, "1i");
+        gpuMath.setUniformForProgram("velocityCalcVerlet", "u_mass", 2, "1i");
+        gpuMath.setUniformForProgram("velocityCalcVerlet", "u_textureDim", [textureDim, textureDim], "2f");
+
         gpuMath.createProgram("velocityCalc", vertexShader, document.getElementById("velocityCalcShader").text);
         gpuMath.setUniformForProgram("velocityCalc", "u_lastPosition", 0, "1i");
         gpuMath.setUniformForProgram("velocityCalc", "u_lastVelocity", 1, "1i");
@@ -317,6 +338,32 @@ function initDynamicSolver(globals){
         gpuMath.setUniformForProgram("velocityCalc", "u_creasePercent", globals.creasePercent, "1f");
         gpuMath.setUniformForProgram("velocityCalc", "u_axialStiffness", globals.axialStiffness, "1f");
 
+        gpuMath.createProgram("positionCalcVerlet", vertexShader, document.getElementById("positionCalcVerletShader").text);
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_lastPosition", 0, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_lastLastPosition", 1, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_lastVelocity", 2, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_originalPosition", 3, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_externalForces", 4, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_mass", 5, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_meta", 6, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_beamMeta", 7, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_creaseMeta", 8, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_nodeCreaseMeta", 9, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_normals", 10, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_theta", 11, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_creaseGeo", 12, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_meta2", 13, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_nodeFaceMeta", 14, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_nominalTriangles", 15, "1i");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_textureDim", [textureDim, textureDim], "2f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_textureDimEdges", [textureDimEdges, textureDimEdges], "2f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_textureDimFaces", [textureDimFaces, textureDimFaces], "2f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_textureDimCreases", [textureDimCreases, textureDimCreases], "2f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_textureDimNodeCreases", [textureDimNodeCreases, textureDimNodeCreases], "2f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_textureDimNodeFaces", [textureDimNodeFaces, textureDimNodeFaces], "2f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_creasePercent", globals.creasePercent, "1f");
+        gpuMath.setUniformForProgram("positionCalcVerlet", "u_axialStiffness", globals.axialStiffness, "1f");
+
         gpuMath.createProgram("thetaCalc", vertexShader, document.getElementById("thetaCalcShader").text);
         gpuMath.setUniformForProgram("thetaCalc", "u_normals", 0, "1i");
         gpuMath.setUniformForProgram("thetaCalc", "u_lastTheta", 1, "1i");
@@ -349,6 +396,10 @@ function initDynamicSolver(globals){
         gpuMath.setUniformForProgram("centerTexture", "u_lastPosition", 0, "1i");
         gpuMath.setUniformForProgram("centerTexture", "u_textureDim", [textureDim, textureDim], "2f");
 
+        gpuMath.createProgram("copyTexture", vertexShader, document.getElementById("copyTexture").text);
+        gpuMath.setUniformForProgram("copyTexture", "u_orig", 0, "1i");
+        gpuMath.setUniformForProgram("copyTexture", "u_textureDim", [textureDim, textureDim], "2f");
+
         gpuMath.createProgram("updateCreaseGeo", vertexShader, document.getElementById("updateCreaseGeo").text);
         gpuMath.setUniformForProgram("updateCreaseGeo", "u_lastPosition", 0, "1i");
         gpuMath.setUniformForProgram("updateCreaseGeo", "u_originalPosition", 1, "1i");
@@ -457,6 +508,7 @@ function initDynamicSolver(globals){
         }
         globals.gpuMath.initTextureFromData("u_lastPosition", textureDim, textureDim, "FLOAT", lastPosition, true);
         globals.gpuMath.initFrameBufferForTexture("u_lastPosition", true);
+        if (globals.integrationType == "verlet") globals.gpuMath.step("copyTexture", ["u_lastPosition"], "u_lastLastPosition");
     }
 
     function setCreasePercent(percent){
@@ -504,6 +556,7 @@ function initDynamicSolver(globals){
         originalPosition = new Float32Array(textureDim*textureDim*4);
         position = new Float32Array(textureDim*textureDim*4);
         lastPosition = new Float32Array(textureDim*textureDim*4);
+        lastLastPosition = new Float32Array(textureDim*textureDim*4);
         velocity = new Float32Array(textureDim*textureDim*4);
         lastVelocity = new Float32Array(textureDim*textureDim*4);
         externalForces = new Float32Array(textureDim*textureDim*4);

js/globals.js

@@ -52,6 +52,7 @@ function initGlobals(){
         //dynamic sim settings
         percentDamping: 0.85,
         density: 1,
+        integrationType: "verlet",
 
         strainClip: 5.0,//for strain visualization, % strain that is drawn red