Fenics Paraview Cannot Read Point Data Array
(updated 26. May)
15th May:
I was officially selected for the Google summertime of code 2017 under the project Develop XDMF format for visualisation and checkpointing. I am very happy to begin work on this project and to aid FEniCS community solve the problems related to the IO.
Nosotros have had the very get-go Google hangout with the mentors on Fri 12th. We've spoken virtually the general moving picture of this project, i.e. how electric current state of XDMF IO is not able to both save for visualization and save for checkpointing computed results. Electric current checkpointing format is binary HDF5 heavy data. FEniCS has methods that could read and write the (HDF5) data, but a priori data virtually the function space (how to interpret the data) must be provided. On the other hand, Paraview itself has no idea how to interpret raw FEniCS data - then XDMF XML structure must be provided in gild to visualize.
26th May:
My DEV enviroment
My Os is the LTS Xubuntu (Ubuntu + xfce) sixteen.04.2. I am using CLion IDE from JetBrains. HW: Intel(R) Core(TM) i5-6200U CPU @ 2.30GHz, 12 GB RAM.
Few words nearly my work
FEniCS is composed of several parts maintained by developers, namely
- dijitso (python)
- dolfin (C++)
- ffc (python)
- fiat (python)
- instant (python)
- mshr (python)
- ufl (python)
Python parts are installed as python modules - no need for C/C++ compilation. The only compiled part is the DOLFIN, umbrella which connects all the other parts and contains some important submodules (IO, visualisation, solvers handling, mesh handling, etc.).
My project is about IO (a dolfin function) so I won't edit any of python parts (hopefully :)).
The virtually of the work is concentrated to XDMFFile.cpp and XDMFFile.h class.
There are some college level user interface methods to write/read dolfin functions, meshes, mesh-value collections etc. Let me consider the bones example of writing a dolfin mesh to a XDMF file from python interface.
import dolfin as d mesh = d.UnitSquareMesh(16, xvi) file = d.XDMFFile("mesh.xdmf") file.write(mesh, d.XDMFFile.Encoding_ASCII) The constructor XDMFFile() takes std::string filename and optionally MPI communicator. If no communicator is provided then MPI_COMM_WORLD is used. XML handling in dolfin is via pugi library, so information technology instantiates new empty pugi xml document. XDMFFile likewise registers few parameters into global dolfin parameter variable.
Parameter rewrite_function_mesh takes boolean value (true by default) and applies when time series is written. If true then mesh is stored for each fourth dimension step, if false then only intial mesh is stored.
Boolean functions_share_mesh is past default false. If more functions are saved into the XDMFFile and so they all share the same mesh for the time step and this should optimize file size.
The concluding parameters is (bool, default false) flush_output, which allows for datasets to be saved after each time step. Could be used to inspect the result of computation while running. Performace cost is important here, plough to true only in debug cases.
Second line mesh = d.UnitSquareMesh ... constructs dolfin Mesh of triangles on \([0,1] \times [0,1]\). Somewhen, void XDMFFile::write(const Mesh &mesh, Encoding encoding) is chosen.
First coding period (30. May - 27. June)
(updated 22nd June)
13th June:
For the first coding period blog postal service, let me review some of XDMF current functionality in DOLFIN and relate it to what needs to be washed during my GSoC projection.
I have discussed in the previous post how, e.g., the XDMF interface to store a Mesh is used. The Mesh IO currently works fine in serial and parallel.
If you run the post-obit code in fenics
import dolfin equally d mesh = d.UnitSquareMesh(5, 5) file_ascii = d.XDMFFile("xdmf_ascii.xdmf") file_ascii.write(mesh, d.XDMFFile.Encoding_ASCII) information technology produces the following single XML file
<?xml version="two.0"?> <!DOCTYPE Xdmf Arrangement "Xdmf.dtd" []> <Xdmf Version= "3.0" xmlns:xi= "http://www.w3.org/2001/XInclude" > <Domain> <Grid Name= "mesh" GridType= "Uniform" > <Topology NumberOfElements= "8" TopologyType= "Triangle" NodesPerElement= "3" > <DataItem Dimensions= "8 3" NumberType= "UInt" Format= "XML" > 0 1 4 0 3 4 1 two v one four 5 iii 4 vii three 6 vii iv 5 viii four 7 8 </DataItem> </Topology> <Geometry GeometryType= "XY" > <DataItem Dimensions= "9 ii" Format= "XML" > 0 0 0.5 0 1 0 0 0.5 0.5 0.five 1 0.5 0 i 0.5 1 ane 1 </DataItem> </Geometry> </Grid> </Domain> </Xdmf> XDMF starts with the header definition and and then encapsulates <Grid> into <Domain> tag. <Topology> describes connectivity of each prison cell (triangles in the instance of UnitSquareMesh). Basically, the topology node contains n-tuple for each cell, which number the vertices of that cell. This numbering and then applies to the <Geometry> tag, where coordinates of points are present. For example, the start cell has vertices numbered 0,one,4 and these vertices take coordinates (0, 0), (0.five, 0), (0.five, 0.5).
If we read this XDMF3 file into ParaView, we go far at
mesh test
This was a simple idea, how the Mesh IO works.
My project is mostly about improving the IO for Functions. Current XDMF could save the Function of whatever chemical element type - but for higher order functions and discontinuous functions it makes some crude local interpolation which decimates the information - so it cannot be fully read back again to fenics. The reason for this is, that XDMF in fenics was meant to be used for visualization mainly and one can truly visualize only "linear" elements with given values at jail cell vertices.
Let me save the office
x y
on a unit of measurement foursquare triangular mesh to a current XDMF format. The code looks similar
import dolfin as d mesh = d.UnitSquareMesh(2, two) file_ascii = d.XDMFFile("xdmf_ascii.xdmf") V = d.FunctionSpace(mesh, "CG", 1) f = d.Role(Five) f.assign(d.project(d.Expression("10[0] * 10[1]", caste= two), Five)) file_ascii.write(f, d.XDMFFile.Encoding_ASCII) and resulting file contains autonomously from the geometry and topology data the following tag
<Attribute Name= "f_4" AttributeType= "Scalar" Center= "Node" > <DataItem Dimensions= "nine ane" Format= "XML" > 1.238011374155259e-fifteen 2.053478914687549e-xv 1.977584762613556e-15 ane.986258379993447e-xv 0.2500000000000034 0.5000000000000033 ane.951563910473904e-xv 0.5000000000000032 ane.000000000000005 </DataItem> </Attribute> These are values of function
10 y
at vertices and are easy to visualize. Doing and then in ParaView gives the following picture
function_test
Although this looks overnice, the higher order functions (piecewise quadratic) in current XDMF are stored only as piecewise linear functions (vertex values).
The master objective of the beginning coding period is to develop parallel write and read functions for XDMF, which store full part-related information = degrees of freedom mapping, values of degrees of freedom and some representation of FiniteElement used.
22nd June:
With the partial support of NumFOCUS I've attended FEniCS 2017 conference hosted in the capital of Grand duchy of luxembourg, 12th-15th June. I've had the pleasure to meet my mentors(C. Richardson and J. Unhurt) in person and discuss with them my GSoC projection.
fenics17
(left) Ivan Yaschuk, Chris Richardson and (right) me.
Back to piece of work now. I've added public methods into XDMFFile class, namely write_checkpoint and read_checkpoint. They are used to save and read Function. Underlying information structure is the same as in existing (and by user tested) HDF5File. Instead of storing vertex values (as described above) I store full information, i.e. cell ordering, number of degrees of freedom per cell, degrees of liberty mapping and bodily values of degrees of freedom.
When user wants to read back (checkpoint) a Role, he must specify an empty Part object where new values are be loaded.
I've too tried my mentor'south ParaView filters - have a look here. Since they are designed to a bit different format of XDMF output every bit I am preparing, I've fabricated slight change to it.
Afterall, functionality, which was not possible in FEniCS then far is e.g. to write in parallel (I did on 4 processes) DG2 part (piecewise polynomial of degree ii only discontinuous through elements), read information technology back, change it and visualize in ParaView.
A nice plot of a time-series of DG2 function due to new functionality:
dg2gif
Second coding period (30. June - 28. July)
(updated 25th July)
07th July:
The first coding period is over and I have my code merged to the master FEniCS branch! See here.
The current state of XDMF FEniCS side functionality is more than-or-less fine. We can relieve functions in parallel and read them back.
There is an ongoing discussion about calculation new attributes to XDMF format that would fit our electric current FEniCS output model. Y'all can track the discussion here. Our ultimate goal is to append the XDMF3 standard and modify the VTK's understanding of XDMF.
VTK (Visualization Toolkit) is a toolkit maintained by Kitware (world wide web.kitware.com) to read and write popular data file formats. It also provides some own file formats. Nice description of VTK's file formats is here http://www.vtk.org/wp-content/uploads/2015/04/file-formats.pdf.
They accept a gitlab repository at https://gitlab.kitware.com/vtk/vtk/ and ParaView, Kitware'southward visualization software is based on VTK data model. In addition, whatever change we brand to the XDMF standard it must map to the electric current VTK model.
25th July:
We take agreed on a format for (FEniCS'south) finite element role written to XDMF. Preliminarily the format is
<Filigree Name= "A filigree" > <Topology></Topology> <Geometry></Geometry> <Aspect Name= "Function_0" Eye= "Other" possibly_other_identifiers= "" ItemType= "FiniteElementFunction" ElementFamily= "DG" ElementDegree= "2" > <DataItem Name= "Alphabetize" > <!-- outset dataitem used for indices (dofmap)--> </DataItem> <DataItem Proper name= "Value" > <!-- second dataitem used for values--> </DataItem> </Attribute> <Attribute Name= "Function_1" Center= "Other" possible_other_identifiers= "" ItemType= "FiniteElementFunction" ElementFamily= "DG" ElementDegree= "2" > <DataItem Proper name= "Alphabetize" > <!-- first dataitem used for indices (dofmap)--> </DataItem> <DataItem Name= "Value" > <!-- second dataitem used for values--> </DataItem> </Attribute> </Grid> In <Topology> and <Geometry> tags the mesh connectivity and mesh point coordinates are saved. Novelty cames in ItemType aspect of Attribute, where you tin can specify FiniteElementFunction and consequently information technology's family and caste in ElementFamily, ElementDegree. Within the attribute there are at to the lowest degree ii DataItem. First - by the position - there is dataitem containing indices for values of freedom. These are just integer coordinates pointing to the second dataitem. The second dataitem contains values of degrees of liberty (bladder numbers). I volition call them dofmap and dofvalues dataitems.
Dofmap is partitioned into tuples for each cell, i.e. dofmap dataitem containing
ways there are two cells - each containing three degrees of freedom and together with dofvalues
nosotros can say, that 0-th jail cell has three dofs with values 0, 1, 0 and ane-st cell has three dofs with values 0, ane, 0.
My latest accomplishment is related to XDMF and VTK libraries. In order to brand ParaView read and visualize this new FiniteElementFunction format I had to make XDMF cpp library understand information technology. This work is in co-operative https://gitlab.kitware.com/michalhabera/xdmf/commits/michal/add together-finite-chemical element-function. Chief difficulty comes in reading Attribute with more than one DataItems. The problem is, that Attribute is in cpp library a child of DataItem - so information technology tin can naturally "contain" simply i dataitem. My workaround to this applies when the attribute is being populated with children. Instead of populating 2 dataitems (which is impossible), I parse dofmap and dofvalues dataitem into i. This is done simply constructing new dataitem with values symbolically written in array notation dofvalues[dofmap]. From the signal of view of VTK, which sees only the resulting populated structure, in that location is only one dataitem under finite element function aspect. This dataitem contains correctly arranged values of degrees of freedom.
In VTK cpp library there are at least 2 ways how to build unstructured grid. First way is to gear up connectivity and geometry equally a whole and aside and so ready is to vtkUnstructuredGrid. Second fashion is to append vtkUnstructuredGrid cell by jail cell using InsertNextCell and betoken by signal using InsertNextPoint. These two approaches differ. The start is more than optimized and stores merely one value per node, simply the second is capable of storing discontinuous functions (considering they accept more values for the same node).
In summary, I have made XDMF to parse multiple dataitems under i attribute and I have besides changed the way VTK builds unstructured filigree from XDMF grid according to python filters prepared by my mentor. Paraview tin can with my XDMF and VTK branches natively (no python filters) read our new XDMF checkpointing format for the simplest functions (CG1, DG1 on triangles).
Third coding menstruum (28. July - the end, 25th August)
(updated 18th Baronial)
18th August
The last coding period is in its meridian and this is my first blog postal service in August.
I accept finally (with the help of my mentor David DeMarle) addressed the "multiple dataitems" issue mentioned in previous postal service. In other words, the new logic is no more spread between both libraries - XDMF and VTK, but is lifted but to the VTK level.
I have as well implemented special type of finite element function - belonging to so called Raviart-Thomas infinite. This is, one tin can say, the first nontrivial role space visualized. Degrees of liberty for RT1 (of degree=1) functions on triangles are situated in midpoints of triangle edges and they specify normal component of the vector field. The vector field is continuous inside each triangle, normal components are constant on each edge and are continuous through edges. This defines uniquely a vector function and information technology is possible to compute its values at nodes. These values could be used to visualize a piecewise discontinuous vector function in VTK/ParaView. My results are shown the following flick. 
The y-component of vector field (x, 0) interpolated to RT1 space is shown. New plot on right is compared to the old manner of plotting (the one-time way consisted of interpolation into continuous function and its visualization).
Several other function spaces were added (continuous and discontinuous functions on quadrilaterals - this is possible with the work of another GSoC FEniCS'due south pupil Ivan Yashchuk).
Now I am testing the work and preparing terminal merge/pull requests to the XDMF and VTK libraries.
Michal Habera 2017, last modified: November 26 2017 15:47:01.
Source: https://www2.karlin.mff.cuni.cz/~habera/?p=gsoc17
0 Response to "Fenics Paraview Cannot Read Point Data Array"
Post a Comment