Merge pull request #2241 from peastman/cross

Fixed bug involving nested cross products

platforms/cuda/src/CudaExpressionUtilities.cpp

@@ -97,6 +97,7 @@ void CudaExpressionUtilities::processExpression(stringstream& out, const Express
             // at once, so check to see if both are needed.
 
             vector<const ExpressionTreeNode*> nodes;
+            nodes.push_back(&node);
             for (int j = 0; j < (int) allExpressions.size(); j++)
                 findRelatedCustomFunctions(node, allExpressions[j].getRootNode(), nodes);
             vector<string> nodeNames;
@@ -155,7 +156,7 @@ void CudaExpressionUtilities::processExpression(stringstream& out, const Express
                     if (derivOrder[0] == 0 && derivOrder[1] == 0)
                         out << nodeNames[j] << " = make_" << tempType << "(dot(" << child1 << ", " << child2 << "));\n";
                     else
-                        throw OpenMMException("Unsupported derivative order for cross()");
+                        throw OpenMMException("Unsupported derivative order for dot()");
                 }
             }
             else if (node.getOperation().getName() == "cross") {

platforms/opencl/src/OpenCLExpressionUtilities.cpp

@@ -88,6 +88,7 @@ void OpenCLExpressionUtilities::processExpression(stringstream& out, const Expre
             // at once, so check to see if both are needed.
 
             vector<const ExpressionTreeNode*> nodes;
+            nodes.push_back(&node);
             for (int j = 0; j < (int) allExpressions.size(); j++)
                 findRelatedCustomFunctions(node, allExpressions[j].getRootNode(), nodes);
             vector<string> nodeNames;
@@ -146,7 +147,7 @@ void OpenCLExpressionUtilities::processExpression(stringstream& out, const Expre
                     if (derivOrder[0] == 0 && derivOrder[1] == 0)
                         out << nodeNames[j] << " = dot(" << child1 << ", " << child2 << ");\n";
                     else
-                        throw OpenMMException("Unsupported derivative order for cross()");
+                        throw OpenMMException("Unsupported derivative order for dot()");
                 }
             }
             else if (node.getOperation().getName() == "cross") {

tests/TestCustomIntegrator.h

@@ -1041,9 +1041,11 @@ void testVectorFunctions() {
     integrator.addGlobalVariable("sumy", 0.0);
     integrator.addPerDofVariable("angular", 0.0);
     integrator.addPerDofVariable("shuffle", 0.0);
+    integrator.addPerDofVariable("multicross", 0.0);
     integrator.addComputeSum("sumy", "x*vector(0, 1, 0)");
     integrator.addComputePerDof("angular", "cross(v, x)");
     integrator.addComputePerDof("shuffle", "dot(vector(_z(x), _x(x), _y(x)), v)");
+    integrator.addComputePerDof("multicross", "cross(vector(1, 0, 0), cross(vector(0, 0, 1), vector(1, 0, 0)))");
     OpenMM_SFMT::SFMT sfmt;
     init_gen_rand(0, sfmt);
     vector<Vec3> positions(numParticles);
@@ -1061,12 +1063,14 @@ void testVectorFunctions() {
     // See if the expressions were computed correctly.
     
     double sumy = 0;
-    vector<Vec3> angular, shuffle;
+    vector<Vec3> angular, shuffle, multicross;
     integrator.getPerDofVariable(0, angular);
     integrator.getPerDofVariable(1, shuffle);
+    integrator.getPerDofVariable(2, multicross);
     for (int i = 0; i < numParticles; i++) {
         ASSERT_EQUAL_VEC(velocities[i].cross(positions[i]), angular[i], 1e-5);
         ASSERT_EQUAL_VEC(Vec3(1, 1, 1)*velocities[i].dot(Vec3(positions[i][2], positions[i][0], positions[i][1])), shuffle[i], 1e-5);
+        ASSERT_EQUAL_VEC(Vec3(0, 0, 1), multicross[i], 1e-5);
         sumy += positions[i][1];
     }
     ASSERT_EQUAL_TOL(sumy, integrator.getGlobalVariable(0), 1e-5);