Merge lp:~asc/fluidity/regularisespatialaspect into lp:fluidity

[Adam Candy]

Added an option '/io/regularise_spatial_aspect', which when enabled will regularise the spatial aspect of the output mesh. This rescales
spatial extent of output position mesh so all spatial extents are of the same length (the smallest of the spatial extents), i.e. scale so
the domain has an aspect ratio of 1:1, or 1:1:1. Note this option should not be used in conjunction with checkpointing (at present).

I certainly find this useful for testing purposes when working with large-aspect ratio domains. It modifies the spatial properties of the
mesh, so is strictly for visualisation and test simulations.

Is this the best place to put the option in the schema?

(now with VTK_interfaces.F90 properly commited)

[Stephan Kramer (s-kramer)]

Hm, I feel this option shouldn't be implemented inside vtk_write_fields(). vtk_write_fields is a low level routine that should just do what it's told to in the code, dump the fields that it's passed and not be dependent on options tree behaviour. Why don't you put this in vtk_write_state_new_option() and pass down the rescaled coordinate to vtk_write_fields() from there? That would also fix the checkpointing issue (which is a definite no-go for me).

I'm pondering how to make this option a bit more generally applicable:
- in general if I rescale my output for viewing (this is actually fairly trivial in paraview), I actually like to control the exact aspect ratio myself. Blowing it up to an exact 1:1 aspect ratio typically heavily distorts the features/mesh elements.
- it should also be made to work on the sphere, in which case it's hard to work out the rescaling parameter automatically
- if you have a free surface with mesh movement, automatically rescaling will cause the rescaling to be different every dump, which is probably not desired

Maybe have this option with a user specified rescaling? Or alternatively, just have a script that rescales the vtus afterwards. I'd like to hear other people's opinion on this...

[Adam Candy (asc)]

Thanks for the feedback Stephan - some useful comments, particularly with regard to free surface simulations. This can be considered work-in-progress, and currently I propose the following:

 1. Move the implementation to a higher level, as suggested,
 2. Add options to provide scaling manually,

which deals with points 1 (fully) and 3 (adequately).

Future work (not critical to merging):

 3. Scale appropriately on the sphere (possibly just the option to extend in the radial direction - fixes Stephan's point 2),
 4. When the domain is scaled automatically, add the option to fix this scaling. This will allow the automatic scaling to be used in free surface simulations (assuming the free surface perturbations are small).

review: Needs Fixing

[David Ham (david-ham)]

As Stephan pointed out above, I'm not sure this belongs in Fluidity as opposed to in a post-processing stage.

In doing 4. Note that the calculated scale will need to be written back to the options tree, otherwise the information will get lost when restarting off a checkpoint.

[Adam Candy (asc)]

I agree that this should be a post-processing step, assuming that it does not incur a significant computational cost.

Unfortunately scaling the domain as described above is not well-supported by many of the unstructured mesh post-processing tools we use, and these transforms are made as an extra post-processing step - which are not so convenient for the user and take time to generate. In any case they do not natively support an automatic rescaling to a 1:1 or 1:1:1 aspect ratio - which I have found very useful running the ice shelf cavity simulations over a range of domain sizes.

I think building these features into Fluidity for diagnostic output is of benefit to users, particularly those regularly running in realistic and hence, non-idealised domains where more post-processing is required.

I can imagine that there could be a scenario where a user has mistakenly switched on this option and post-processes expecting the domain to be in a regular scaling - but I think we have to expect that the user knows what they are doing, and as long as these options are well documented, this should not be a problem.

With regard to David's note concerning point 4 - it's not the intention to apply transformations to checkpoint VTUs, only diagnostic output VTUs.

[David Ham (david-ham)]

On 5 September 2011 17:17, Adam Candy <email address hidden> wrote:
>
> With regard to David's note concerning point 4 - it's not the intention to apply transformations to checkpoint VTUs, only diagnostic output VTUs.

You misunderstand my point. The point is that if you are automatically
calculating a rescaling because the domain size is changing, you need
to record that rescaling because the rescaling you calculate on
resuming from a checkpoint won't be the same as the original one you
calculated.

--
Dr David Ham
Applied Modelling and Computation Group
Department of Earth Science and Engineering
Imperial College London

http://www.imperial.ac.uk/people/david.ham

[Adam Candy (asc)]

Ah yes, that is the case when the scaling is calculated at the start. In the other cases where the scaling is done on the fly, or prescribed explicitly in the flml, this won't be necessary. I'll add this as point 5 in the future work.

Unmerged revisions

3530. By Adam Candy on 2011-07-20

This (quite important) line got skipped in the previous commits to the domain regularisation option.

3529. By Adam Candy on 2011-07-20

Improved schema commet for new option /io/regularise_spatial_aspect.

3528. By Adam Candy on 2011-07-20

Added an option '/io/regularise_spatial_aspect', which when enabled will regularise the spatial aspect of the output mesh.
This rescales spatial extent of output position mesh so all spatial extents are of the same length (the smallest of the spatial extents), i.e. scale so the domain has an aspect ratio of 1:1, or 1:1:1.
Note this option should not be used in conjunction with checkpointing (at present).

=== modified file 'femtools/VTK_interfaces.F90'
--- femtools/VTK_interfaces.F90	2011-02-07 19:13:51 +0000
+++ femtools/VTK_interfaces.F90	2011-07-20 13:58:35 +0000
@@ -251,6 +251,7 @@
     type(scalar_field) :: l_model
     type(vector_field) :: v_model(3)
     type(tensor_field) :: t_model
+    real, dimension(position%dim) :: spatialextent
     logical :: dgify_fields ! should we DG-ify the fields -- make them discontinous?
     integer, allocatable, dimension(:)::ghost_levels
     real, allocatable, dimension(:,:) :: tempval
@@ -450,6 +451,24 @@
     do i=1, position%dim
       v_field_buffer(:,i)=v_model(position%dim)%val(i,:)
     end do
+
+    ! Regularise the domain so all spatial extents are of the same length
+    ! (the smallest of the spatial extents)
+    ! i.e. scale so the domain has an aspect ratio of 1:1, or 1:1:1
+    if (have_option("/io/regularise_spatial_aspect")) then
+      spatialextent = 1.0
+      do i=1, position%dim
+        spatialextent(i) = abs(  maxval( v_model(position%dim)%val(i,:) ) - minval( v_model(position%dim)%val(i,:) ) )
+        ewrite(3,*) "Spatial: ", minval( v_model(position%dim)%val(i,:) ),  maxval( v_model(position%dim)%val(i,:) )
+      end do
+      ewrite(3,*) "Regularised spatial aspect, spatial extents: ", spatialextent
+      spatialextent(:) = minval(spatialextent(:)) / spatialextent(:)
+      ewrite(3,*) "Regularised spatial aspect, spatial scalings: ", spatialextent
+      do i=1, position%dim
+        v_field_buffer(:,i)= spatialextent(i) * v_field_buffer(:,i)
+      end do
+    end if
+
     do i=position%dim+1, 3
       v_field_buffer(:,i)=0.0
     end do

=== modified file 'schemas/fluidity_options.rnc'
--- schemas/fluidity_options.rnc	2011-07-07 13:22:20 +0000
+++ schemas/fluidity_options.rnc	2011-07-20 13:58:35 +0000
@@ -132,6 +132,14 @@
                   attribute name { xsd:string }
                }
             ),
+            ## Regularise spatial aspect - rescales spatial extent of output position mesh.
+            ## Regularise the domain so all spatial extents are of the same length
+            ## (the smallest of the spatial extents),
+            ## i.e. scale so the domain has an aspect ratio of 1:1, or 1:1:1.
+            ## Note this option should not be used in conjunction with checkpointing (at present).
+            element regularise_spatial_aspect {
+               comment
+            }?,
             ## Options for convergence analysis.
             element convergence {
                ## Whether to enable the creation of a convergence

=== modified file 'schemas/fluidity_options.rng'
--- schemas/fluidity_options.rng	2011-07-14 15:54:58 +0000
+++ schemas/fluidity_options.rng	2011-07-20 13:58:35 +0000
@@ -164,6 +164,16 @@
           </element>
         </choice>
         <optional>
+          <element name="regularise_spatial_aspect">
+            <a:documentation>Regularise spatial aspect - rescales spatial extent of output position mesh.
+Regularise the domain so all spatial extents are of the same length
+(the smallest of the spatial extents),
+i.e. scale so the domain has an aspect ratio of 1:1, or 1:1:1.
+Note this option should not be used in conjunction with checkpointing (at present).</a:documentation>
+            <ref name="comment"/>
+          </element>
+        </optional>
+        <optional>
           <element name="convergence">
             <a:documentation>Options for convergence analysis.</a:documentation>
             <optional>