How to visualise the difference between two fields defined over the same mesh?

First save the data into two ASCII VTK files. For example,

nmagpp --vtk=m.vtk --vtkascii --fields=m simulation_name

Note the option "--vtkascii" to force the creation of a ASCII file.
Let's say this command created the two files m-000000.vtk and m-000001.vtk.
You can now use the library pyvtk to load the two files, compute the difference and save it back to a third file.

import numpy, pyvtk
a = pyvtk.VtkData("m-000000.vtk")
b = pyvtk.VtkData("m-000001.vtk")
va = a.point_data.data[0].vectors
vb = b.point_data.data[0].vectors
for i in range(len(va)):
    va[i] = list(numpy.array(va[i]) - numpy.array(vb[i]))
a.tofile("difference.vtk")

Save this text to a file named diff.py and run it as:

python diff.py

You'll get a third file with name difference.vtk containing the difference of the two fields.

A step forward

If you are repeating this operation many times, it may become annoying to open again and again the diff.py file to change the names of the input files. You can then modify the script as follows:

import sys, numpy, pyvtk
a = pyvtk.VtkData(sys.argv[1])
b = pyvtk.VtkData(sys.argv[2])
va = a.point_data.data[0].vectors
vb = b.point_data.data[0].vectors
for i in range(len(va)):
    assert a.structure.points[i] == b.structure.points[i]
    va[i] = list(numpy.array(va[i]) - numpy.array(vb[i]))
a.tofile(sys.argv[3])

The name of the files are taken from the command line. You can then compute the difference using:

python diff.py a.vtk b.vtk a_minus_b.vtk 

Notice that in the last version of the script we also added the line,

assert a.structure.points[i] == b.structure.points[i]

which does just check that the two files are using the same set of points (i.e. the same mesh).

Also available in: HTML TXT