Dec 152012

Today, as I was debugging a molecular dynamics code I stumbled on a bug. After finding the bug and correcting it, I went to make a tea and while I was waiting for the microwave to heat up the water, I realized that although I could still remember the exasperating way in which the bug manifested itself, I already forgot the cause of the bug. I then concluded that it would be wise-and perhaps even useful-to other programmers to collect the bugs that I have stumbled upon so as to remember and avoid to make the same mistakes more than once. Here is the first one.

Forgetting to define the copy constructor. This is a pretty silly mistake, since once you create a class in C++ you should quickly write the copy constructor and the assignment operator (=). In my case when I ran my MD code in gdb I was getting an error in the line:

void FF_TullyTest1::compute_force(Snapshot& snapshot)
  // Build the electronic Hamiltonian
  double V11 = sign * A * (1.0 - std::exp(-B * std::fabs(x)));
  double V12 = C * std::exp(-D * x * x);
  H[0][0] = V11;  // GDB says: "H[0][0] is not a valid memory address"

After confirming that H (a matrix member of the class FF_TullyTest1) is allocated, I was a bit puzzled. It turns out that the bug originated from a different class:

/// Stores a trajectory
class Trajectory
    /// Initialize with no snapshot data
    Trajectory(int id, int natoms, int nstates, std::vector<double> mass, std::vector<int> Z);
    /// Initialize with snapshot data
    Trajectory(int id, int natoms, int nstates, std::vector<double> mass, std::vector<int> Z,
               std::vector<vec_double> r, std::vector<vec_double> v, int state);
    /// Position (au)
    vec_double* r_;
    /// Position in the image cell (au)
    vec_double* ri_;

The problem is that the class Trajectory contains a few pointers and does not define a copy constructor. This combination is pretty dangerous because the default copy constructor is not pointer friendly. So a when a Trajectory object is copied using the default copy constructor, like in the code below, you have a bug:

    // Copy the trajectory object and link the snapshot object
    trajectory_ = SharedTrajectory(new Trajectory(trajectory));

The problem is that instead of creating a copy of the data, the default copy constructor just performs a copy of all the class data (including pointers) and calls the default constructor for each of its member objects. Solution: write (and remember to update) a customized copy constructor. All this may be avoided if instead of storing the coordinates as arrays (vec_double* r_;) I were to use the standard template library vector (std::vector r_;). However, I doubt vector would be as fast as a raw array [In C++11 one can easily access the underlying array of a std::vector directly with a new function called data()].

 December 15, 2012  Programming

Sorry, the comment form is closed at this time.