The ultimate way to postprocess OpenFoam data in Python (updated to Pyvista)

Abstract

In this post, I use foamToVTK in OpenFoam to convert OpenFoam data into legacy VTK (The Visualization ToolKit)) format, then use vtkInterface for data manipulation in Python under Ubuntu.

Introduction

OpenFoam is a popular open source code for computational fluid dynamics (CFD). Although it contains various helpful postpocessing modules in command line such as postProcess, it is still designed for convenience but not flexibility. For example, it only provides operations that are very common in the context of fluid mechanics or vector mathematics and it hides the details of operation such as the numerical scheme to approximate the derivatives. Most of the time, OpenFoam saves the data in folder named by the current time and in each folder contains a special OpenFoam format-txt like data, which is also designed for convenience such that one can directly read the result in the field data. However, if one wants to manipulate the data in a more flexible sense in the modern data-driven era, a Python-script-driven manipulation of data is extremely favorable. Also, to avoid dealing with the mesh in the script, if would be great if one can simply add the modified field on the original mesh.

Fortunately, with the help of an awesome Python package on Github: currently Pyvista previously vtkInterface , originally by Alex Kaszynski, one can easily leverage the powerful libraries in Python environment to postprocessing traditional, mature, standard and specialized scientific computing data and immediately put them back in to again, leverage the existing powerful visualization software in scientific computing community.

Using pip to install vtkInterface

Note that it supports well in Python 3.5+.

sudo pip install pyvista

Tutorial: 2D Flow past cylinder

The material can be obtained from Wolf Dynamics at this link.

Prepare data

1. untar the .tar file

   tar -zxvf vortex_shedding.tar.gz ./
   

2. go to c1 directory for running a standard case

   cd c1
   blockMesh
   checkMesh
   icoFoam > log &
   

Convert OpenFoam default format to VTK

foamToVTK

Using Python to manipulate VTK data

import pyvista as vtki
import numpy as np

## grid is the central object in VTK where every field is added on to grid
grid = vtki.UnstructuredGrid('./VTK/c1_1000.vtk')

## point-wise information of geometry is contained
print grid.points

## get a dictionary contains all cell/point information
print grid.cell_arrays # note that cell-based and point-based are in different size
print grid.point_arrays # 

## get a field in numpy array
p_cell = grid.cell_arrays['p']

## create a new cell field of pressure^2
p2_cell = p_cell**2
grid._add_cell_scalar(p2_cell, 'p2')

## remember to save the modified vtk
grid.save('./VTK/c1_1000_shaowu.vtk')

Visualize the new field in ParaView

paraview