Merge branch 'update_testKernel' into 'master'

Update test kernel

Removes `doself1_density_vec_2` function. Updates `testKernel.c` to include tests for `kernel_deval_1_vec` and `kernel_deval_2_vec`.

See merge request !338

tests/difffloat.py

@@ -35,13 +35,18 @@ file1 = sys.argv[1]
 file2 = sys.argv[2]
 number_to_check = -1
 
-if len(sys.argv) == 5:
-    number_to_check = int(sys.argv[4])
-
 fileTol = ""
 if len(sys.argv) >= 4:
     fileTol = sys.argv[3]
 
+if len(sys.argv) >= 5:
+    number_to_check = int(sys.argv[4])
+
+if len(sys.argv) == 6:
+    ignoreSmallRhoDh = int(sys.argv[5])
+else:
+    ignoreSmallRhoDh = 0
+    
 data1 = loadtxt(file1)
 data2 = loadtxt(file2)
 if fileTol != "":
@@ -102,8 +107,11 @@ for i in range(n_lines_to_check):
             print ""
             error = True
 
-        if abs(data1[i,j]) < 1e-6 and + abs(data2[i,j]) < 1e-6 : continue
-            
+        if abs(data1[i,j]) + abs(data2[i,j]) < 1e-6 : continue
+
+        # Ignore pathological cases with rho_dh
+        if ignoreSmallRhoDh and j == 8 and abs(data1[i,j]) < 2e-4: continue
+        
         if( rel_diff > 1.1*relTol[j]):
             print "Relative difference larger than tolerance (%e) for particle %d, column %d:"%(relTol[j], i,j)
             print "%10s:           a = %e"%("File 1", data1[i,j])

tests/test27cellsPerturbed.sh.in

@@ -10,7 +10,7 @@ do
 
     if [ -e brute_force_27_perturbed.dat ]
     then
-	python @srcdir@/difffloat.py brute_force_27_perturbed.dat swift_dopair_27_perturbed.dat @srcdir@/tolerance_27_perturbed.dat 6
+	python @srcdir@/difffloat.py brute_force_27_perturbed.dat swift_dopair_27_perturbed.dat @srcdir@/tolerance_27_perturbed.dat 6 1
     else
 	exit 1
     fi

tests/testKernel.c

@@ -39,7 +39,7 @@ int main() {
   const float numPoints_inv = 1. / numPoints;
 
   for (int i = 0; i < numPoints; ++i) {
-    u[i] = i * 2.5f * numPoints_inv / h;
+    u[i] = i * 2.25f * numPoints_inv / h;
   }
 
   for (int i = 0; i < numPoints; ++i) {
@@ -55,18 +55,22 @@ int main() {
 
 #ifdef WITH_VECTORIZATION
 
+  printf("\nVector Output for kernel_deval_1_vec\n");
+  printf("-------------\n");
+
+  /* Test vectorised kernel that uses one vector. */
   for (int i = 0; i < numPoints; i += VEC_SIZE) {
 
     vector vx, vx_h;
     vector W_vec, dW_vec;
 
     for (int j = 0; j < VEC_SIZE; j++) {
-      vx.f[j] = (i + j) * 2.5f / numPoints;
+      vx.f[j] = (i + j) * 2.25f / numPoints;
     }
 
     vx_h.v = vx.v / vec_set1(h);
 
-    kernel_deval_vec(&vx_h, &W_vec, &dW_vec);
+    kernel_deval_1_vec(&vx_h, &W_vec, &dW_vec);
 
     for (int j = 0; j < VEC_SIZE; j++) {
       printf("%2d: h= %f H= %f x=%f W(x,h)=%f dW(x,h)=%f\n", i + j, h,
@@ -85,6 +89,63 @@ int main() {
     }
   }
 
+  printf("\nVector Output for kernel_deval_2_vec\n");
+  printf("-------------\n");
+  
+  /* Test vectorised kernel that uses two vectors. */
+  for (int i = 0; i < numPoints; i += VEC_SIZE) {
+
+    vector vx, vx_h;
+    vector W_vec, dW_vec;
+
+    vector vx_2, vx_h_2;
+    vector W_vec_2, dW_vec_2;
+    
+    for (int j = 0; j < VEC_SIZE; j++) {
+      vx.f[j] = (i + j) * 2.25f / numPoints;
+      vx_2.f[j] = (i + j) * 2.25f / numPoints;
+    }
+
+    vx_h.v = vx.v / vec_set1(h);
+    vx_h_2.v = vx_2.v / vec_set1(h);
+
+    kernel_deval_2_vec(&vx_h, &W_vec, &dW_vec, &vx_h_2, &W_vec_2, &dW_vec_2);
+
+    /* Check first vector results. */
+    for (int j = 0; j < VEC_SIZE; j++) {
+      printf("%2d: h= %f H= %f x=%f W(x,h)=%f dW(x,h)=%f\n", i + j, h,
+             h * kernel_gamma, vx.f[j], W_vec.f[j], dW_vec.f[j]);
+
+      if (fabsf(W_vec.f[j] - W[i + j]) > 2e-7) {
+        printf("Invalid value ! scalar= %e, vector= %e\n", W[i + j],
+               W_vec.f[j]);
+        return 1;
+      }
+      if (fabsf(dW_vec.f[j] - dW[i + j]) > 2e-7) {
+        printf("Invalid value ! scalar= %e, vector= %e\n", dW[i + j],
+               dW_vec.f[j]);
+        return 1;
+      }
+    }
+    
+    /* Check second vector results. */
+    for (int j = 0; j < VEC_SIZE; j++) {
+      printf("%2d: h= %f H= %f x=%f W(x,h)=%f dW(x,h)=%f\n", i + j, h,
+             h * kernel_gamma, vx_2.f[j], W_vec_2.f[j], dW_vec_2.f[j]);
+
+      if (fabsf(W_vec_2.f[j] - W[i + j]) > 2e-7) {
+        printf("Invalid value ! scalar= %e, vector= %e\n", W[i + j],
+               W_vec_2.f[j]);
+        return 1;
+      }
+      if (fabsf(dW_vec_2.f[j] - dW[i + j]) > 2e-7) {
+        printf("Invalid value ! scalar= %e, vector= %e\n", dW[i + j],
+               dW_vec_2.f[j]);
+        return 1;
+      }
+    }
+  }
+
   printf("\nAll values are consistent\n");
 
 #endif

tests/tolerance_27_perturbed.dat

 #   ID      pos_x      pos_y      pos_z        v_x        v_y        v_z           rho        rho_dh        wcount     wcount_dh         div_v       curl_vx       curl_vy       curl_vz
     0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1.2e-6     1e-4	    5e-5       2e-3		 3.1e-6	     3e-6	   3e-6		 3e-6
-    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1.2e-2	    1e-5       1e-4		 2e-5	     2e-3	   2e-3	 	 2e-3
+    0	    1e-6       1e-6	  1e-6 	       1e-6 	  1e-6	     1e-6	   1e-6	      1e-3	    1e-5       1e-4		 2e-5	     2e-3	   2e-3	 	 2e-3