I’m looking for an answer regarding how ParaView (currently using v5.6.0) calculates the gradient on an unstructured dataset. I have a volume mesh filled with velocity data from my ANSYS Fluent CFD simulation; I use my scalar array/vector for CFD velocities and get the gradient from using the “Gradient Of Unstructured DataSet” filter within ParaView.
I’ve tried looking for the source code and all over the forums, but I’ve come up with nothing. If you are able to point me towards this specific filter’s source code and/or explain how ParaView numerically solves for the unstructured gradient, I would be extremely appreciative. It appears that ParaView will do what I want (eventually get Wall Shear Stress), but I need to be able to argue what order of accuracy/understand exactly what is going on with this calculation.
The implementation of ParaView’s Graident Of Unstructured DataSet filter is the VTK class vtkGradientFilter. You can get to the source code by downloading the ParaView source code, the VTK source code, or going to the VTK git repository. You can browse the code from the web here:
The operation of the filter is simple. It first iterates over all the cells and computes local gradients for each one. It does this by using the Derivatives method for the appropriate cell object. That is, each type of cell for an unstructured grid has an associated VTK class to manage that type of cell, and each of these has a Derivatives method to compute a local gradient based on parametric coordinates inside the cell. These are all located in the Common/DataModel directory of the VTK repository. Below are links to the class objects for tetrahedra and hexahedra, which I find are the most common cells in an unstructured grid.
Of course, when computing these local gradients, there are (typically) discontinuities at the boundaries of neighboring cells. To manage this and create a continuous gradient, the Gradient Of Unstructured DataSet filter averages out all gradients on cells incident to each point. There are actually 2 modes to do this. The first (default) mode computes for each cell a separate gradient for each of its points. (Cell types like a hexahedron can have different gradient values at each of its point.) When gradients are averaged over a point, it averages the gradient field for each incident cell at that point. However, if you turn on the Faster Approximation option, the filter will only compute one gradient per cell at the cell’s center and then use those to average at points (basically a standard cell-to-point operation).
Thank you for the response, Kenneth! I’ve been looking through the scripts trying to understand… I’m not a programmer, so I’m having a little difficulty identifying what math goes into finding the gradient. I’ve found the areas describing the Jacobian Inverse and Interpolation Derivatives, but I wanted to check-in again.
Note: I am using the cell type “Tetrahedron” with my velocity data saved as a vector within “Point Data.” Is the gradient calculation only considering points which touch the cell/the nodes of each tetrahedron-cell? I believe the Jacobian Inverse identifies an empty matrix and is then filled by the “functionDerivs()” line–could be way off with that. If you’re able to provide further information on how the local derivative works mathematically/that section of the code, I’d appreciate it once again!
If I’m unclear, I could potentially try to write-up by-hand.
Yes, the gradient computation within the cell only considers those points incident on that cell. A tetrahedron, which is a simplex, has a completely linear interpolation and so the gradient, which is the same everywhere in the tetrahedron, is essentially the slope of the scalar in that space.
Like I said, the gradients in adjacent tetrahedra will be different and thus will be averaged along the points.
Does this remain the same/very similar for the “Gradient” filter? I have the option to generate a structured grid. With either method, do you know of a way to dictate how many elements along the surface normal-vector should be used for calculating the gradient at the wall/boundary?
With the unstructured method, I’m now looking at using the “Point Dataset Interpolator” filter to bring the gradients you discussed onto the surface. If I could sample along the normal, this might work as well. It appears that I can only choose from a radius or total amount of nearest-neighbors.
No, the Gradient filter, which uses vtkImageGradient, does a significantly different computation than Gradient of Unstructured Grid. This filter uses finite differences to estimate the gradient. For each point, it takes the difference of the neighbors in the X direction for the X component of the gradient, and likewise for the Y and Z directions.
That’s really good to know… Is there a way to create a structured grid from my currently unstructured mesh with ParaView? I’d like to stay within the interface; however, if not, what form of imported geometry would you recommend?
The “Table To Structured Grid” doesn’t appear to be working very well. Although, I haven’t spent more than a day whittling away at it. The finite difference method should solve my issue.
The easiest way to create a structured grid from an unstructured mesh (with volumetric elements like tetrahedra) is to use the Resample to Image filter.
You mention using Table to Structured Grid. That sounds like maybe your data is actually points stored in something like a csv file. If that is the case, you might try the Table To Points filter and then the Gaussian Resampling filter.
I’ve been working on getting the structured grid and applying the Gradient filter once my Resample To Image is generated. However, I’m having issues getting Gradient to work. When I import my CFD results, I have:
Point Data: contains point 3D space orientation, velocity vector, velocity-component array(s) [x, y, z]
Cell Data: Block Number, Cell ID, Cell Type (Tetrahedron)
I can create the structured grid, and the vectors/scalars are transferred over well. I’m able to visualize nicely with the “Volume” rendering. I use Transform to scale the data back to the correct scale (m → mm) before creating the structured grid. I am able to apply Gradient to either version of the structured grid, but I’m having errors.
For the case where I’ve transformed > grid generation, I get:
No Scalar Field has been specified - assuming 1 component!
Which is strange, because I am feeding it a scalar. When applying to the NON-transformed grid, I don’t receive any errors. However, the bounds of the gradient-dataset are stretched/do not stay at the same location or scale. The top two images show the structured grid/shape of my volume render. The bottom two show how Gradient caused a massive expansion.
I think you can ignore error messages that you are getting. The Gradients still seem to be computed.
I don’t have a simple workaround for the rescaling of the data. You could use the transformation controls under the display parameters to rescale and move the data to where it is supposed to be.
Hi Kenneth and Jmuskat,
I’m stack on a similar issue and hope you can help me!
As J, I want to compute Wss from velocity data. My raw data is a little bit different though. I have a volume of velocities (x,y,z) and an .stl file which defines the region of interest I want to study. The geometry is an aorta model. The velocity data came from mri data and the stl from a segmentation of the aorta from images.
This is what I did: GradientOfUnstructuredData of the velocity field to get the velocity tensor–> resample it over the .stl model–> extracted surface normals from the stl geometry.
I would like to ask you:
if what I did so far is correct
How would you compute the WSS having the 9 components of the tensor and the 3 components of the normal to the surface?
This was my idea:
WSS_x= viscosity* (t1n_x + t2n_y + t3n_z);
WSS_y= viscosity (t4n_x + t5n_y + t6n_z);
WSS_y= viscosity (t7n_x + t8n_y + t9*n_z);
and then compute the magnitude
Thanks a lot!
So far what you describe (compute gradients on volume then resample on stl) is how I would recommend doing it.
For the second part, can you run a second computation after resampling on the stl data? If the stl data does not directly have the normals, you can run the Generate Surface Normals to create it. You can then use the Python Calculator to compute the WSS.
Hi Kenneth,
I was trying to validate the wss computation from paraview with results from cfd simulation.
Input in paraview=velocity field from cfd and then I computed gradients etc…
I realized that comparing du/dx (from cfd) and Gradient0 (from paraview) are different (one order less). And of course the wss it’s different too…Do you have any idea why this is happening?
Thanks a lot!
I’m not sure what the issue is. This question started with a scaling problem. Could the issue be that the CFD code and ParaView are computing the gradient with different scales of dx?
I’ve checked the scales of both velocity and size of the grid and they are both in meters (in paraview). Do you think that remeshing the velocity field in paraview could change something? My idea so far is that the computation of the derivatives in paraview is done by finite difference while fluent could work on finite volumes so probably is more accurate. Do you think is more correct computing the gradient on cell data or point data? Because also changing from point to cell the results are completely different…
Thanks a lot!
Sorry to get into the thread. But I want to compare measured velocity fields from medical imaging and velocity fields from CFX as a result of a CFD simulation. What formats did you use to export your velocity field from CFD? Was it in ANSYS or is it a custom code format.
I’m facing a similar problem. I developed a pipeline to calculate the wall shear stress from 4D Flow MRI data. To validate my pipeline I loaded a Computational Fluid Dynamic simulation I did on Ansys fluent (Tetahedral volumetric mesh) and I applied the “Gradient of Unstructured DataSet” Filter on the velocity vector obtained from the CFD simulation to compute the stress tensor. Nevertheless, the wss estimate I obtain is an order of magnitude lower than that obtained from Ansys. I computed the Vorticity using the filter and compared it with that from the CFD simulation and the value seems to be half that computed on Ansys. I’m pretty sure all the units are ok. Any help is more than welcomed!
Hi Xabier, I went trough the same streamwork as you (4d flow-cfd-paraview) and found the same results. I think the problem is in how the partial derivatives are computed. Let me know if you find a solution, I m curios to see what we can do with paraview
I struggled with this for 2 weeks tried to calculate WSS using Paraview. In the end it was not working due to the limitations in Paraview and it was not possible to calculate the slope near the wall with a predefined filter.
So, I decided to make a script in Python. Luckily the new PyVista library have some nice functionalities for 3D objects, mesh, and surfaces. You can find an example code on my github page.
Feel free to contribute and leave any comments. https://github.com/EdwardFerdian/wss_mri_calculator
