If you set aside for a moment the constraint on the vector x, then the problem of maximizing
(a sum over repeated indices is implied) translates, upon taking the derivative with respect to xi,
into the equation
Dijxj = 0.
This is a special form of the
eigenvalue equation for the matrix D, with the
eigenvalue being 0. In order for non-trivial solutions, the matrix D cannot have an inverse, so it's
determinant must vanish. Finding the eigenvectors can
be done with PDL; see
Mastering Algorithms with Perl for a discussion. The
Math::Cephes::Matrix module can also do this for real symmetric matrices. Generally, the components of the
eigenvectors found in this way are not completely determined; by convention, the normalization xixi = 1 is imposed.
Update: The preceding needs to be augmented by a
test on the
Hessian of the matrix to determine what type of
critical point is present.