For visualization check emails with Christoph Geuzaine.
Also, Coherence; and Mesh.Algorithm = 5; are two very powerful options.

Solve Phi_u for choosen initial conditions (Full frap, half frap) and for choosen dimension (3D,2D or like 1D). Comparison with matlab results for full FRAP, Root mean square, computing time. Video possible

Create many non-overlapping spheres in a box, write to file for easy copying into
.geo file. After creation of .msh file this can be converted to .xdmf which is suitable
for parallel use.

Code can be run for nu=10^-15 or 10^-6, significant deviations for precision=1.

Works reasonably well, but deviations around 0.5%. Might need sphere
inside sphere geometry, instead of box. Potential influence of corners.

Also introduce check for mass conservation in prob_laplace. Currently
does not sum to 1, a factor of 2 is missing somewhere.

This is done with a loop that iteratively creates array 'x' by starting
at zero and continuing until gradient is hit, at which point the mesh
starts becoming finer and then coarser again.

Whether or not two meshes agree well numerically seems to depend on
how well connected two differently spaced mesh domains are.