Clang gives error, GCC and Visual Studio ok: dimension errors in a custom matrix class

Question

Test case is here https://gist.github.com/8875125 (the test itself is at the bottom of the page, the rest is class definition and operator overloads used by the test - suggestions for how to reduce this any further are welcome)

Expected behavior, given by g++ and Visual Studio, is that the for loop at the end of dmtest.cxx should run successfully, L.Print("") should output the entries of a 20-by-20 matrix.

For some reason building with clang (from package clang-3.3 on Ubuntu 13.10, or with Xcode CLT in OS X 10.9, or with the clang 3.1 package in x86 Cygwin), the DMatrix colon() operator results in different dimensions at an intermediate iteration of the for loop at line 67 of dmtest.cxx. All dimensions are expected to remain constant, but according to lldb/gdb they are not. Is there a bug in one of the DMatrix class operators that's resulting in incorrect or undefined behavior that just gets lucky with g++, or what? Appreciate any help.

Edit: Okay, the absolutely amazing C-Reduce tool (http://embed.cs.utah.edu/creduce/) reduced my code down to the following almost-meaningless testcase:

class A
{
public:
    A (const A &);
    ~A ();
};

A::A (const A &)
{
}

A & fn1 ()
{
}

A::~A ()
{
}

A a = fn1 ();
int
main ()
{
}

g++ compiles this to an executable which runs and successfully does nothing, but clang++ gives a compile-time warning ("control reaches end of non-void function") and "Illegal instruction (core dumped)" when run. I'll just send this to the clang mailing list, since I have no idea what's going on now.

Answer 1

Thanks to some help on the Clang mailing list, I now have a solution to the original problem. The tiny test case above was reduced too far by creduce (I should have told it not to accept code that gives compiler warnings) and is not the cause of the original issue.

If I change line 67 of dmtest.cxx from

L(i,colon()) = ( L(i,colon())&(  x-pointx(j)*ones(1,length(x)) )  )/(pointx(i)-pointx(j));

to

L(i,colon()) = elemProduct( L(i,colon()),(  x-pointx(j)*ones(1,length(x)) )  )/(pointx(i)-pointx(j));

then Clang and GCC both compile and run successfully, giving the same output.

There may be some undefined behavior triggered by the allocation of temporary matrices, since

DMatrix& DMatrix::operator &(const DMatrix B) const

is by-value whereas

DMatrix& elemProduct( const DMatrix& A, const DMatrix& B )

is by-reference. Otherwise operator& just calls elemProduct, so they should result in the same numerical output.