Behavioral Non-Portability in Scientific Numeric Computing

The precise semantics of floating-point arithmetic programs depends on the execution platform, including the compiler and the target hardware, even if these are fully compliant with the widely adopted IEEE 754 floating-point standard. Platform dependencies infringe on the highly desirable goal of software portability (which is in fact promised by heterogeneous computing frameworks like OpenCL): the same program run on the same inputs on different platforms can produce different results. Serious doubts on the portability of numeric applications arise when these differences are behavioral, i.e. when they lead to changes in the control flow of a program.

In this research we design an algorithm that takes a numeric procedure and determines an input that may lead to different branching decisions depending on how the arithmetic in the procedure is compiled. Our algorithm proceeds in two phases. The first phase determines a candidate input, such that minor variations of the value of exactly one of the input variables causes different branching decisions. In the second phase, the algorithm efficiently searches though all possible platform parameters such as decisions made by the compiler about evaluation order, and use or non-use of FMA hardware features. The goal is to find two sets of platform parameters that lead to different branching decisions.

Our implementation of the algorithm requires minimal intervention by the user. We illustrate its operation on a diverse set of examples, characteristic of scientific numeric computing, where control flow divergence actually occurs across different execution platforms.

Experiments

We test our algorithm on three procedures: (1) ray tracing, in the raySphere function there is a branching decision over the D variable (D> 0). The inputs are s, r, which represent the coordinates of the sphere origin and the ray, and radiusSq, a scalar value giving the square of the radius of sphere; (2) molecular dynamics, this open source code calculates the Lennard Jones potential of a molecular system. The energy is calculated based on pairwise distances of atoms and their velocity and acceleration. In the ComputeAccel function a decision branch occurs over (rr , rrCut). rr is the pairwise distance between atoms j1 and j2 and rrCut is the threshold for pairwise distances. The calculation of rr is directly effected by a density variable of the atom using a face centered cubic model; (3) long summation of 32 numbers, we compare serial C code which accumulates a value to a register, and a reduction kernel, written in OpenCL which is the common way to implement long summation on a parallel architecture. The result of the sum is compared to zero.

1. Generating Unstable Inputs

To test our algorithm, we first apply transformations to each procedure: we attach a main function that calls the tested procedure; we annotate the symbolic variables and volatile expressions by calling klee_make_symbolic_with_sort and klee_tag_reorderable. Following shows the ray tracing code after transformation.

//Dot Product 3-Vectors
float dot3(float *a, float *b){     
  float r = a[0] * b[0] + a[1] * b[1] 
           + a[2] * b[2];
  klee_tag_reorderable(&r, 0, 2);
  return r;
}

int raySphere(float *r, float *s, float radiusSq) {
  float A = dot3(r,r);                       
  float B = -2.0 * dot3(s,r);               
  float C = dot3(s,s) - radiusSq;          
  float D = B*B - 4*A*C;                       
  if (D > 0)
    return 0;
} 

int main() {
  float r[3], s[3];
  float radiusSq;
  klee_make_symbolic_with_sort(&r, sizeof(r), "r", 8, 32);
  klee_make_symbolic_with_sort(&s, sizeof(s), "s", 8, 32);
  klee_make_symbolic_with_sort(&radiusSq, sizeof(radiusSq), "radiusSq", 0, 32);
  return raySphere(r, s, radius);
} 

For a conditional q in the execution path, there may exist multiple inputs that exhibit control flow instability. In our experiments we split the domian of one input variable into subintervals, and run our algorithm on each of these subintervals. In the table below it shows the experiment results.

2. Control Flow Divergence

To show control flow divergence, we took the generated sets of inputs and tested them on six different architectures. We convert the serial c code into OpenCL code so that we can run the same code and compare the results on these heterogeneous architectures. Below is the list of these architectures.

2.1 Ray Tracing Results

For the generated 408 sets of inputs. 45 of these produce results on different sides of zero for D on different architectures. In the table below some of these inputs are shown. D1, D2 and D3 show the value of D for D1: machine 1,2, D2: machine 3,4,5 and D3: machine 6.

2.2 Molecular Dynamics Results

Following table shows some of the results. R is the difference in pairwise distances of the atoms j1 and j2 to rrCut (rr-rrCut). R1,R2 shows result for machine 1,2,6 and machine 3,4,5 respectively.

2.3 Long summation and Reduction Results

The 32 input numbers have a range between [-1.05,1]. Sum (1,2,3,5,6) shows the reduction results for machine 1,2,3,5,6 and sum 4 shows the reduction results for machine 4. Our code compared the result of the summation to 0. The reduction kernel and serial summation results on the inputs generated produce results on either side of zero.

People:

• Faculty: Miriam Leeser, Thomas Wahl
• Research assistant: Mahsa Bayati, Jaideep Ramachandran, Yijia Gu

Publications:

• Miriam Leeser, Jaideep Ramachandran, Thomas Wahl, and Devon Yablonski. OpenCL floating point software on heterogeneous architectures- portable or not? In Workshop on Numerical Software Verification (NSV), 2012. [pdf ]
• Miriam Leeser, Saoni Mukherjee, Jaideep Ramachandran, and Thomas Wahl. Make it real: Effective floating-point reasoning via exact arithmetic. In Design Automation and Testing in Europe (DATE), 2014.
• Yijia Gu, Thomas Wahl, Mahsa Bayati, and Miriam Leeser. Behavioral non-portability in scientific numeric computing. In Parallel and Distributed Computing (EURO-PAR), pages 558--569, 2015. [pdf]
  

Sponsorship: