View source

<div align="center">
<font size="+2" color="darkblue">
A Customized Python Module for CFD Flow Analysis within VisTrails
</font><br />
'''''by Joel E. Tohline, Jinghua Ge, Wesley Even, &amp; Erik Anderson'''''
</div>

'''<font color="#007700">A relatively simple, customized Python module that plugs smoothly into an otherwise standard workflow within VisTrails facilitates a quantitative analysis of complex fluid flows in simulations of merging binary stars.</font>'''

==<font color="#0000DD">Introduction</font>==

Researchers in the open source community are steadily improving scientific visualization tools.  These tools are providing a wider array of sophisticated probes for data analysis and a wider assortment of effective user-friendly interfaces.  They're also making it easier for researchers in the computational science community &#150; across many disciplines &#150; to effectively analyze huge datasets by drawing on the human brain's acute ability to sort through complex and time-varying visual patterns.  The astrophysics group at Louisiana State University (LSU), for example, routinely uses volume-rendering and ray-tracing algorithms in conjunction with animation techniques to examine the time-varying behavior of isodensity surfaces that arise in computational fluid dynamic (CFD) simulations of mass-transferring and merging binary star systems.  Although such analyses generally provide only a qualitative identification and assessment of structure within a given dataset, the insight gained from visual inspection can nevertheless be extremely valuable.  For example, it was through visual inspection that researchers at LSU initially spotted the nonlinear development of triangular-, square-, and pentagonal-shaped tidal resonances in recent simulations.

LSU's astrophysics group has begun to incorporate VisTrails into its arsenal of scientific visualization and data analysis tools.  VisTrails primarily interested the group a few years ago because it provides a user-friendly workflow interface to the extensive VTK software library.  It also automatically tracks the provenance of data analysis efforts.  However, what most impresses us now is the ease with which VisTrails facilitates the insertion of home-grown analysis modules into an otherwise VTK-based workflow.  Taking advantage of this additional programming versatility, we have gained a greater appreciation of the role that visualization tools can play in the quantitative assessment of results from large-scale simulations.  In this article, we first describe the VTK-based workflow that we initially constructed in VisTrails to view streamlines within each binary mass-transfer simulation.  We then describe the Python module, whose insertion into this workflow has permitted us to identify values of key rotational frequencies associated with such flows.

==Aside==

Before diving into a discussion of the VisTrails workflows that have been developed for this study, it will be useful to download the relevant .vt module from the VisTrails database.
<div align="center">

<!-- vistrail host="vistrails.sci.utah.edu" db="vistrails" vtid="11" version="1826" tag="" showspreadsheetonly="True"/ -->

<!-- vistrail host="vistrails.sci.utah.edu" db="vistrails" vtid="18" version="1833" tag="" showspreadsheetonly="True"/ -->

<vistrail host="vistrails.sci.utah.edu" db="vistrails" vtid="19" version="1834" tag="" showspreadsheetonly="True"/>
</div>

==<font color="#0000DD">Base Workflow</font>==

[[Image:Figure1_Tohline.jpg |border|center|600px|Figure 1.  Screenshots of the window within the VisTrails builder that displays user-designed visualization workflows. (a = left) The base workflow we constructed from standard VTK-based modules.  (b = inset) A segment of the workflow that's hidden inside the <font face="Courier">Draw_Streamlines</font> group module. (c = right) The customized workflow we used to create Figure 2, in which we inserted the <font face="Courier">SwitchCoord</font> module containing our customized Python script into the base workflow.]]

Within VisTrails, we initially selected various VTK-based modules to do the following, in sequence (see Figure 1a):
*''Read simulation data.'' We used <font face="Courier">vtkPLOT3DReader</font> to read in one file containing the <math>(x, y, z)</math> coordinate locations of every vertex on our 3D cylindrical coordinate mesh and a separate file containing the fluid's mass-density (scalar) and momentum-density (3D vector) at every grid vertex.
*''Outline cylindrical domain boundary.'' As shown, we enlisted <font face="Courier">vtkStructuredGridOutlineFilter</font>, <font face="Courier">vtkPolyDataMapper</font>, and <font face="Courier">vtkActor</font>.
*''Define isodensity surfaces.'' We rendered two nested isodensity surfaces to outline high- (red) and low-density (blue) flow regions.  The <font face="Courier">Red_contour</font> and <font face="Courier">Blue_contour</font> module groups each contain <font face="Courier">vtkContourFilter</font>, <font face="Courier">vtkDataSetMapper</font>, <font face="Courier">vtkProperty</font>, and <font face="Courier">vtkActor</font>.
*''Draw streamlines.'' As Figure 1b shows, each of the eight separate <font face="Courier">Draw_Streamlines</font> module groups uses <font face="Courier">vtkStreamLine</font>, <font face="Courier">vtkTubeFilter</font>, <font face="Courier">vtkDataSetMapper</font>, <font face="Courier">vtkProperty</font>, <font face="Courier">vtkActor</font>, <font face="Courier">vtkSphereSource</font>, <font face="Courier">vtkPolyDataMapper</font>, and <font face="Courier">vtkLODActor</font> to trace an individual streamline within the flow.  Streamline lengths are set by feeding a common <font face="Courier">Propagation_Time</font> into all eight module groups.

Vistrails renders the output from the various workflow actors in a composite scene using <font face="Courier">vtkRenderer</font> as viewed by an observer located at a position that <font face="Courier">vtkCamera</font> specifies.  Finally, the module <font face="Courier">vtkCell</font> directs this scene to the VisTrails interactive spreadsheet.

In this initially constructed ''base workflow,'' VisTrails pipes the 3D vector field representing the momentum density distribution from the <font face="Courier">vtkPLOT3DReader</font> module directly into each of the eight <font face="Courier">Draw_Streamline</font> module groups.  This base workflow &#150; which VisTrails assembles using generically available <font face="Courier">vtk</font> modules &#150; lets us examine the behavior of streamlines in our binary mass-transfer simulations, but only from the frame of reference, <math>\Omega_0</math>, in which we originally performed each simulation (see Figure 2b, labeled <math>\Delta\Omega = 0.00</math>).

==<font color="#0000DD">Customized Python Module</font>==
To make it possible for us to examine the properties of binary mass-transfer flows from reference frames that have a range of different angular frequencies of rotation, <math>\Omega_\mathrm{frame}  = (\Omega_0 + \Delta\Omega)</math>, we wrote a Python-based module &#150; <font face="Courier">SwitchCoord</font> &#150; for insertion into the base VisTrails workflow.  Figure 1c shows our resulting customized VisTrails workflow; it differs very little from the base workflow.  The [http://www.vistrails.org/index.php/User:Tohline/CustomizedPythonModule#Sidebar sidebar] shows the complete Python source code from our customized module.  The code segment that performs the required ''physics analysis'' is short and straightforward.  In particular, the <font face="Courier">SwitchCoord</font> module performs the following operations at each grid vertex:
* converts the <math>(x, y)</math> Cartesian to <math>(R, \phi)</math> cylindrical coordinates;
* divides the momentum components by density to obtain the velocity components if the density is greater than <font face="Courier">minp</font> (otherwise, it sets the velocity components to zero);
* shifts the azimuthal velocity component <math>v_\phi</math> to a new, rotating frame of reference by adding <math>R\times \Delta\Omega</math>;
* converts the cylindrical velocity components to Cartesian velocity components; and
* normalizes the velocities to the maximum velocity <font face="Courier">maxnorm</font> found across the domain where densities are greater than <font face="Courier">minp</font>.

We designed the output ports on <font face="Courier">SwitchCoord</font> to provide access to the same type of structured arrays that <font face="Courier">vtkPLOT3DReader</font> generates.  But, in our customized VisTrails workflow, which includes <font face="Courier">SwitchCoord</font> (see Figure 1c), the 3D vector field that VisTrails pipes into each of the eight <font face="Courier">Draw_Streamline</font> module groups represents the fluid's velocity distribution as viewed from the rotating frame of reference that the floating-point scalar, <font face="Courier">Omega_frame</font>, specifies.

==<font color="#0000DD">Interpretation of Results</font>==

[[Image:Figure2_Tohline.jpg |border|center|600px|Figure 2.  Screenshot of the VisTrails spreadsheet after we used five different values of <font face="Courier">Omega_frame</font> <math>(\Delta \Omega)</math> to execute our customized workflow.  Each 3D-rendered image displays eight equatorial-plane streamlines and a pair of isodensity surfaces (lower density surface colored blue; higher density surface colored red) that outline the structure of both stars as well as the connecting mass-transfer stream.  The propagation time is the same for all eight streamlines; along six streamlines, VisTrails carries out the integration in both directions from the location marked by a small colored sphere.]]
<center>'''Figure 2''' (scroll over image for caption)</center>

Figure 2 displays 3D renderings of the flow from one of our binary mass-transfer simulations as generated by our customized VisTrails workflow.  We have generated images assuming five different frame rotation frequencies, as specified by <math>\Delta\Omega</math>.  Aside from labeling <math>\Delta\Omega</math> values under each image, we produced Figure 2 by simply taking a screenshot of the VisTrails interactive spreadsheet.  The spreadsheet feature has proven to be extremely useful in this analysis because it facilitates the side-by-side comparison of scenes that VisTrails has rendered using different parameter values.  And, although we can't demonstrate it here in print, VisTrails lets users zoom, pan, and interactively rotate all 3D-rendered scenes simultaneously.

==<span id="Sidebar">Sidebar:</span> <font color="#0000DD">SwithCoord Python Module</font>==
Here we present a complete listing of the python code form our customized <font face="Courier">SwitchCoord</font> program module.  At the beginning of the "Physics Analysis" segment of the code, we assign names to the data arrays that have been acquired as input from <font face="Courier">vtkPLOT3DReader</font>: ''<font face="Courier">pcoords</font>'' is a tuple that identifies the Cartesian-based coordinate location of each grid vertex, ''<font face="Courier">density</font>'' is a scalar that specifies the mass density, and ''<font face="Courier">momentum</font>'' is a tuple that specifies the values of the cylindrical-coordinate-based vector momentum at each grid vertex.  We detail the remaining logic of the physics analysis in the main text.

<pre>
import core.modules.module_registry
from core.modules.vistrails_module import Module, ModuleError
import vtk, math
version="0.0.0"
name="SwitchCoord"
identifier="edu.lsu.switchcoord"

class SwitchCoord(Module):
    def compute(self):
        minp = self.getInputFromPort("min_density")
        Domega = self.getInputFromPort("Domega")
        dataset=self.getInputFromPort("dataset")
        output = self.create_instance_of_type('edu.utah.sci.vistrails.vtk', 'vtkStructuredGrid')
        output.vtkInstance = vtk.vtkStructuredGrid()
        mydata=output.vtkInstance
        mydata.DeepCopy(dataset.vtkInstance)
        self.op(mydata, minp, omega)
        self.setResult("changed_dataset", output)

    ###################################
    ##
    ## Begin:  Physics Analysis
    ##
    ###################################

    def op(self, mydata, minp, omega):
        extent=mydata.GetExtent()
        pcoords = mydata.GetPoints().GetData()
        density = mydata.GetPointData().GetScalars("Density")
        momentum = mydata.GetPointData().GetVectors("Momentum")

        maxnorm = 0.0
        for i in range(0, mydata.GetNumberOfPoints()):
            [x, y, z] = pcoords.GetTuple3(i)  
            [_v1, _v2, _v3] = momentum.GetTuple3(i)
            p = density.GetValue(i)
            r = math.sqrt(x*x + y*y)         
            phi = math.atan2(y, x)
            if p < minp: 
                vx=vy=vz=0
            else:
                vr = _v1 / p
                vphi = _v2 / (p) + r * omega
                vz = _v3 / p
                vx = vr * math.cos(phi) - vphi * math.sin(phi)
                vy = vr * math.sin(phi) + vphi * math.cos(phi)

                norm = math.sqrt(vx*vx + vy*vy + vz*vz)
                if norm > maxnorm:
                    maxnorm = norm
            momentum.SetTuple3(i, vx, vy, vz)

        for i in range(0, mydata.GetNumberOfPoints()):
            [vx, vy, vz] = momentum.GetTuple3(i)
            vx = vx/maxnorm
            vy = vy/maxnorm
            vz = vz/maxnorm
            momentum.SetTuple3(i, vx, vy, vz)

    ###################################
    ##
    ## End:  Physics Analysis
    ##
    ###################################

def initialize(*args, **keywords):
    reg=core.modules.module_registry.registry
    reg.add_module(SwitchCoord)
    reg.add_input_port(SwitchCoord, "scalar_range", 
                       [core.modules.basic_modules.Float, 
                        core.modules.basic_modules.Float]) 
    reg.add_input_port(SwitchCoord, "min_density", core.modules.basic_modules.Float)
    reg.add_input_port(SwitchCoord, "omega", core.modules.basic_modules.Float)
    reg.add_input_port(SwitchCoord, "dataset", 
                       (reg.get_descriptor_by_name(
                'edu.utah.sci.vistrails.vtk', 'vtkStructuredGrid').module) )
    reg.add_output_port(SwitchCoord, "changed_dataset", 
                        (reg.get_descriptor_by_name(
                'edu.utah.sci.vistrails.vtk', 'vtkStructuredGrid').module) )

def package_dependencies():
    import core.packagemanager
    manager = core.packagemanager.get_package_manager()
    if manager.has_package('edu.utah.sci.vistrails.vtk'):
        return ['edu.utah.sci.vistrails.vtk']
    else:
        return []
</pre>