Remove Neumann free BC type and directly set RHS for Neumann integrated
Reimplement Neumann point based conditions

The previous implementation did not integrate with other boundary
condition types and was difficult to understand to fix, so the
Neumann point conditions have been re-written from scratch.

This adds a new BoundaryConditions_NeumannIntegrate routine
that integrates point Neumann conditions over faces in 3D and
lines in 2D, and a BoundaryConditionsNeumannType for storing
information about the Neumann point conditions.

Tracker items 1809 and 2111.
Fix uninitialised value bug with boundary lines and faces

Lines and faces that weren't actually on a boundary could
have had boundary set to true, as it wasn't initialised.
Don't integrate faces/lines where an integrated value is set

Allows you to set up a problem like a cube or square where
only one face/side is integrated over.
Remove old Neumann condition routines

Also remove Pressure Poisson face integration routine, as this isn't
actually used anywhere and makes use of the old Neumann conditions
structure.
Remove old Neumann analytic functions

These don't seem to have ever worked properly, and don't
look easy to understand to fix and maintain them.
Deallocate Neumann information if error in initialise

Also make sure boundary conditions neumannBoundaryConditions pointer
is nullified in finalise routine.
Support 1D point Neumann conditions and remove unused array
Move integrated Neumann values to separate field parameter set

This will allow problems to use a mix of integrated and point
Neumann values, which are added together to give the force vector.
Add BC type for integrated Neumann with no terms from point values

Add a "NEUMANN_INTEGRATED_ONLY" boundary condition type. This sets
an integrated RHS value directly, and prevents the integration of
any point values on a line or face containing this DOF.
Fix distributed_matrix_by_vector_add to work with non-square matrices

In this case the product vector domain mapping matches the row
mapping, and the multiplied vector domain mapping matches the
column mapping.
Make sure integrated values parameter set is created on all nodes
Use distributed matrices and vectors for Neumann point integration

This is a step towards using a distributed sparse matrix to decrease
the memory usage, and also allows us to store the resulting product
vector directly in the field parameter set vector for the Neumann
integrated values.
Add bounds check for gfortran 4.6 with debug options

This would have helped me earlier...
Add new list routines

- List_AppendList for appending a list to another list
- List_ClearItems for removing all items without deallocating memory.
Compressed column matrix storage fixes

- Fix checking the row and column indices.
- Support multiplying a compressed column matrix by a vector in
  distribute_matrix_by_vector_add.
Set sparsity structure for Neumann integration matrix

This uses a compressed row format for storage, which takes
longer to calculate the sparsity structure in this case but
is faster when performing the matrix multiplication.
Fixes for face derivatives

Decomposition face calculations now allowed (previously commented out)
Fixes to DECOMPOSITION_TOPOLOGY_FACES_CALCULATE for calculation of domain faces
Update to DERIVATIVE_NUMBERS_IN_LOCAL_FACE in basis types and its calculation for Hermite bases
Merge Prasad's face_derivatives branch into neumann_conditions

This is required for getting derivative and version numbers on faces.
Correct derivative numbers and scale factors in Neumann conditions

For lines and faces, the element derivative number needs to be
obtained from the local face/line derivative numbers. And the scale
factors also need to be taken from the face/line.