[Libmesh-devel] [Libmesh-users] A possible bug in 3D?

[John P. <pet...@cf...>]

dx...@my... writes:
 > Hi, I'd like to report a possible bug in the 3D element matrix calculation.
 > 
 > The Hamiltonian in my test problem is:
 > 
 > H =  -1/2*laplacian
 >        + 1/2*x^2+2*x^4+1/2*x^6
 >        + 1/2*y^2+2*y^4+1/2*y^6
 >        + 1/2*z^2+2*z^4+1/2*z^6
 >        + x*y + x*z + y*z
 > 
 > So, the code for the element stiffness matrix is
 > 
 > const Real x = q_point[qp](0);
 > const Real y = q_point[qp](1);
 > const Real z = q_point[qp](2);
 > 
 > Ke(i,j) += JxW[qp]*(.5*dphi[i][qp]*dphi[j][qp]
 >                               + (      .5*pow(x, 2.)+2.*pow(x, 4.)+.5*pow(x, 6.)
 >                                       +.5*pow(y, 2.)+2.*pow(y, 4.)+.5*pow(y, 6.)
 >                                       +.5*pow(z, 2.)+2.*pow(z, 4.)+.5*pow(z, 6.)
 >                                       +x*y+x*z+y*z)
 >                               *phi[i][qp]*phi[j][qp]);
 > 
 > The problem is that the matrix assembly stalled when I had 6th power
 > functions in my code, I pin-pointed this problem by simply replacing
 > pow(x, 6.) with x*x*x*x*x*x and the code ran just fine. The odd thing
 > is that pow(x, 6) works for a 2D test case.

"Stalled"?  What does that mean, exactly?  And why is a problem with "pow"
a problem with libmesh?

-John

 > I'd appreciate it if anyone will look into this and let me know.
 > 
 > Thanks a lot in advance,
 > 
 > -DX