[[Imported from SVN: r8972]]

[[Imported from SVN: r8862]]

[[Imported from SVN: r8861]]

[[Imported from SVN: r8687]]

[[Imported from SVN: r8630]]

[[Imported from SVN: r8629]]

[[Imported from SVN: r8627]]

There were two issues fixed by this patch:

a) The bending energy only involves the third director.
   One must look at it as a function in S^2, not as
   a subset of coefficients of a function to SO(3).
   This means that we have to use a different projection
   in one place.
b) A certain subset of the coefficients (see the code)
   must always be zero, because they are parts of
   3-vectors which are known to be perpendicular to (0,0,1).
   For some strange reason they aren't.  I now hardwire
   them to be 0 anyway.  This makes my test pass.
   Maybe one day I'll find the true reason :-)

[[Imported from SVN: r8626]]

of Neumann boundary forces.  But I am not sure and I want to
postpone the issue.  Therefore this patch adds a cpp warning.

[[Imported from SVN: r8625]]

Don't use a loop with 4 iteration to clear an array,
if the array has only length 3.

Boy, this did take a while to find...

[[Imported from SVN: r8593]]

Another easy patch which makes the code a lot faster.

Also, this patch removes the quasi-Newton code.  So far that code
hasn't helped at all.  Quasi-Newton methods will need more work
in the future.

[[Imported from SVN: r8588]]

the ones with i<=j.  The rest follows from symmetry of
the Hessian.

Needless to say, this makes the Hessian assembly quite
a bit faster.

[[Imported from SVN: r8587]]

[[Imported from SVN: r8572]]

[[Imported from SVN: r8571]]

[[Imported from SVN: r8570]]

[[Imported from SVN: r8569]]

Because trace(sym(A)) == trace(A)

[[Imported from SVN: r8568]]

[[Imported from SVN: r8567]]

[[Imported from SVN: r8566]]

[[Imported from SVN: r8565]]

[[Imported from SVN: r8562]]

[[Imported from SVN: r8559]]

[[Imported from SVN: r8554]]

[[Imported from SVN: r8553]]

[[Imported from SVN: r8543]]

[[Imported from SVN: r8508]]

[[Imported from SVN: r8507]]

This change makes the localgeodesicfefunction test fail again.
So there must be a bug in the corresponding non-finite-difference
code.  Looking at simplex cases confirms me in this opinion.
But I don't find the bug...

[[Imported from SVN: r8463]]

[[Imported from SVN: r8462]]

[[Imported from SVN: r8409]]

implement the remaining curvature and bending terms for the energy gradient.  These are not working yet, though.

[[Imported from SVN: r8408]]

with respect to the quaternion coordinates.  Or, in other words,
the second derivative of the map from quaternions to 
orthogonal matrices.

[[Imported from SVN: r8407]]

So far, the long form of the quadratic bending energy has been
implemented and tested, and it appears to do the right thing.

This is a snapshot commit.  The code remains out-commented until
it is finished.

[[Imported from SVN: r8403]]

[[Imported from SVN: r8385]]

[[Imported from SVN: r8384]]

[[Imported from SVN: r8376]]

New implementation properly handles the two sheets.

Now the TargetSpace test passes again.  Yay!

[[Imported from SVN: r8375]]

The new implementation should take proper care of the two-sheet structure.

I am writing 'should' here: it still fails the test sometimes,
but then the result appears to be close enough to think rounding errors
are responsible.

[[Imported from SVN: r8373]]

The new implementation properly handles the double-sheet structure

[[Imported from SVN: r8372]]

don't keep a private implementation of derivativeOfArcCosSquared -- use the one from UnitVector instead

[[Imported from SVN: r8371]]