Stresses at the integration points in ParaView

[FEAnalyst]

One of the limitations of CalculiX is that it can only output averaged nodal stresses in .frd files. If you need them at the integration points (where they are originally calculated and most accurate) you can use .dat file output but that’s not very handy with larger models. Fortunately, ParaView can help with this. It features a whole group of filters called Quadrature Points. Here’s what you can do with them, in order:

1. Generate Quadrature Scheme Dictionary - uses a built-in library of Gauss points (unfortunately, only full integration linear and quadratic elements of the following types are supported at the moment: triangle, quadrilateral, tetrahedron*), necessary for subsequent filters
2. Interpolate to Quadrature Points - as the name suggests, this filter interpolates the nodal data to integration points
3. Generate Quadrature Points - generates integration points so that you can visualize them and make measurements using point selection tools like Hover Points On

For better visibility, increase the Point Size a bit (you can also render them as spheres) and show the results of the Interpolate to Quadrature Points Filter with Representation set to Wireframe or some Opacity.

*Linear tetrahedron elements in CalculiX have a single integration point while ParaView assumes 4 integration points for this type of element. Of course, ParaView doesn’t account for CalculiX-specific expansion of 2D elements to 3D ones. So it might be best to use this functionality for quadratic tetrahedrons for now. Fortunately, they are the most commonly used elements in PrePoMax.

[mkraska]

Showing the integration point values is an interesting feature.
Yet, I guess that the integration point values displayed as described above aren’t the true IP values, because they are back-interpolated from averaged nodal values if you don’t suppress the averaging.

You will most likely have different values in ParaView and in the dat file. In particular, in ParaView you still can have IPs with the equivalent stress exceeding the yield stress, whereas in the dat file this should not happen.

To avoid this, you need to output the nodal values for each element separately. How to do this has been discussed in other threads here. I made a python script for this purpose.

CalculiX-Examples/Scripts at master · calculix/CalculiX-Examples - separate.py

[FEAnalyst]

Of course, this is just an interpolation of averaged nodal stresses from the .frd file to integration points. It’s a simple workaround if someone wants to generate those points (they are actually created and displayed) and make measurements at them. Their coordinates can be displayed as well. More accurate solutions require scripts (also in ParaView). It would be best if a user could choose the way stresses are written to the .frd file by CalculiX using keywords. For example, Abaqus provides the following options:

*ELEMENT OUTPUT, POSITION=AVERAGED AT NODES / CENTROIDAL / INTEGRATION POINTS / NODES

Calculix could have nodal averaged, nodal and integration point output.

[mkraska]

I understand that the FRD file format generally isn’t capable of handling integration point variables. Yet, with the separation.py approach you get the precise IP values in Paraview assuming that the interpolations on ccx and ParaView side are equivalent. No need for additional scripting inside ParaView.

For Prepomax, it might be an option to use the dat file values (with assumed locations of the IPs. Based on the dat file, all the options you listed from ABAQUS would be available in the postprocessor.

[Matej]

I wonder if the Python script can affect the computational time, compatibility with other SPCs or MPCs like rigid bodies, and whether it works for shell elements. I think this approach would be nice when two materials share the same mesh. Using the procedure, the material domains could be split and correct nodal values obtained.

So a universal solution would be to use .dat file to get IP values and then use the post-processing set by the user. I am only afraid that this would work only for smaller meshes since the .dat file will be huge.

[mkraska]

In principle, the computation time should not be affected too much, because the number of dofs does not change. However, the amount of data changes massively, because internal nodes of the mesh are split into 8 nodes in a structured hex mesh.

As to the materials: It is on my todo list to create the split not by element but by element set. That would prevent averaging over set boundariesm much like ABAQUS doesn’t average over different materials. Should not be too complex to do. Perhaps one could call separate.py with the mesh file and some additional set definition files (e.g. .nam files from CGX). Inside Prepomax you could offer this for materials by default.

Recently I tested the approach for plane strain elements and it seems to work, except for the node numbering depencence, where separate.py is not to blame. So it might work also for shells but I didn’t test it.

[jaumesagues]

In my case, I have results (stresses) at gauss points and I need to write vtk files with these results and later open this vtk file in paraview.

[ANYS]

Hello,

I have a .dat file for a unique element with coordinates and stresses according to the standard ccx output format.
¿Is there any chance/trick to import those int.Points as new nodes with their stress values assigned to them just for visualization purposes?.
Thanks

S T E P 1


INCREMENT 1


stresses (elem, integ.pnt.,sxx,syy,szz,sxy,sxz,syz) for set INT_E and time 0.1000000E+01

84 1 -4.424639E+05 -5.033502E+07 5.008302E+07 -8.359406E+05 6.832757E+05 1.130077E+04
84 2 -5.184794E+05 -4.977215E+07 5.028015E+07 -8.324547E+05 6.869696E+05 1.125365E+04
84 3 -5.947606E+05 -4.921010E+07 5.047727E+07 -8.293110E+05 6.906657E+05 1.121116E+04
84 4 -4.916443E+05 -5.056893E+07 4.994240E+07 -5.128480E+05 6.821303E+05 6.934875E+03
84 5 -5.673642E+05 -4.998341E+07 5.014623E+07 -4.988423E+05 6.859105E+05 6.745522E+03
84 6 -6.433571E+05 -4.939875E+07 5.035005E+07 -4.851960E+05 6.896928E+05 6.561028E+03
84 7 1.124758E+04 -5.046714E+07 5.009931E+07 -1.106555E+05 6.773608E+05 1.497355E+03
84 8 -6.498273E+04 -4.985944E+07 5.030927E+07 -8.629843E+04 6.812311E+05 1.168041E+03
84 9 -1.414924E+05 -4.925265E+07 5.051919E+07 -6.231801E+04 6.851035E+05 8.438203E+02
84 10 -4.517438E+05 -5.033504E+07 5.009223E+07 -8.361759E+05 -5.254211E-08 0.000000E+00
84 11 -5.277718E+05 -4.977216E+07 5.028943E+07 -8.326889E+05 7.005614E-08 -4.104852E-09
84 12 -6.040655E+05 -4.921009E+07 5.048662E+07 -8.295440E+05 8.757018E-09 2.531325E-09
84 13 -5.008798E+05 -5.056893E+07 4.995161E+07 -5.128949E+05 -5.254211E-08 9.577988E-10
84 14 -5.766124E+05 -4.998340E+07 5.015552E+07 -4.988880E+05 -1.401123E-07 2.189254E-09
84 15 -6.526179E+05 -4.939873E+07 5.035940E+07 -4.852406E+05 -1.751404E-08 -4.104852E-10
84 16 2.062261E+03 -5.046715E+07 5.010846E+07 -1.105067E+05 -8.757018E-08 -4.378509E-09
84 17 -7.418082E+04 -4.985944E+07 5.031848E+07 -8.614846E+04 -7.005614E-08 4.104852E-09
84 18 -1.507033E+05 -4.925262E+07 5.052848E+07 -6.216694E+04 -1.313553E-08 2.736568E-10
84 19 -4.424639E+05 -5.033502E+07 5.008302E+07 -8.359406E+05 -6.832757E+05 -1.130077E+04
84 20 -5.184794E+05 -4.977215E+07 5.028015E+07 -8.324547E+05 -6.869696E+05 -1.125365E+04
84 21 -5.947606E+05 -4.921010E+07 5.047727E+07 -8.293110E+05 -6.906657E+05 -1.121116E+04
84 22 -4.916443E+05 -5.056893E+07 4.994240E+07 -5.128480E+05 -6.821303E+05 -6.934875E+03
84 23 -5.673642E+05 -4.998341E+07 5.014623E+07 -4.988423E+05 -6.859105E+05 -6.745522E+03
84 24 -6.433571E+05 -4.939875E+07 5.035005E+07 -4.851960E+05 -6.896928E+05 -6.561028E+03
84 25 1.124758E+04 -5.046714E+07 5.009931E+07 -1.106555E+05 -6.773608E+05 -1.497355E+03
84 26 -6.498273E+04 -4.985944E+07 5.030927E+07 -8.629843E+04 -6.812311E+05 -1.168041E+03
84 27 -1.414924E+05 -4.925265E+07 5.051919E+07 -6.231801E+04 -6.851035E+05 -8.438203E+02

global coordinates (elem, integ.pnt.,x,y,z) for set INT_E and time 0.1000000E+01

84 1 1.000855E+00 1.773224E-02 -1.353137E-02
84 2 1.000081E+00 1.774597E-02 -1.352090E-02
84 3 9.993063E-01 1.775969E-02 -1.351043E-02
84 4 1.000742E+00 9.992254E-03 -1.352984E-02
84 5 9.999676E-01 1.000000E-02 -1.351937E-02
84 6 9.991931E-01 1.000775E-02 -1.350890E-02
84 7 1.000688E+00 2.252285E-03 -1.352911E-02
84 8 9.999134E-01 2.254033E-03 -1.351864E-02
84 9 9.991388E-01 2.255781E-03 -1.350817E-02
84 10 1.000947E+00 1.773224E-02 0.000000E+00
84 11 1.000172E+00 1.774597E-02 0.000000E+00
84 12 9.993976E-01 1.775969E-02 -1.029992E-18
84 13 1.000834E+00 9.992254E-03 1.734723E-18
84 14 1.000059E+00 1.000000E-02 2.602085E-18
84 15 9.992845E-01 1.000775E-02 4.336809E-19
84 16 1.000779E+00 2.252285E-03 0.000000E+00
84 17 1.000005E+00 2.254033E-03 1.517883E-18
84 18 9.992302E-01 2.255781E-03 -3.794708E-19
84 19 1.000855E+00 1.773224E-02 1.353137E-02
84 20 1.000081E+00 1.774597E-02 1.352090E-02
84 21 9.993063E-01 1.775969E-02 1.351043E-02
84 22 1.000742E+00 9.992254E-03 1.352984E-02
84 23 9.999676E-01 1.000000E-02 1.351937E-02
84 24 9.991931E-01 1.000775E-02 1.350890E-02
84 25 1.000688E+00 2.252285E-03 1.352911E-02
84 26 9.999134E-01 2.254033E-03 1.351864E-02
84 27 9.991388E-01 2.255781E-03 1.350817E-02

[NorbertH]

Can you share your .inp and the results? When the model got no offsets in the node and element ids this should be doable.

[ANYS]

Analisys is described as LE4: Axisymmetric hyperbolic shell under uniform internal pressure.

LE4: Axisymmetric hyperbolic shell under uniform internal pressure

Target solution: On the midsurface at point K the meridional stress, σϕϕ, is −50.0 MPa and the hoop stress, σθθ, is 50.0 MPa.

My nodal values are σϕϕ, −49.76 MPa and the hoop stress, σθθ, is 50.42 MPa which is perfectly fine but I was tring to look at the closest integration point just to confirm.

I end up ploting the integration points as new nodes into the *inp. Later I selected them graphically and identify the coordinate for later search on the list (not very elegant). Here is not so difficult but it can be cumbersome in the midle of a solid and with arbitrary element orientation. It would be nice if prepomax could read the dat with that information and plot it.

LE4_Axisymmetric_shell_internal_pressure_fine.inp (42.2 KB)

[NorbertH]

The extra nodes prevents me from linking the results with the mesh. I commented them out.
I can only do it for the 2d output otherwise the linking fails.

https://transfer.coreform.com/Xg6PLXQo89/LE4_Axisymmetric_shell_internal_pressure_fine.zip

First find your desired element. And then take a look at the integrations points. As calculix expands 2d elements into 3d internally the integration points are overlapping in the thickness. Just select them all and view them in the table.

[ANYS]

Thanks NorbertH. I will look at that procedure.

[NorbertH]

I would have added a detailed description with screenshots but i recently joined the forum and can’t add images -.-

[FEAnalyst]

Please try again now.

[NorbertH]

test image

[NorbertH]

okay this works.
should i just show the paraview part or the whole workflow?
I would need to reference a commercial mesher and a open source plugin for the complete workflow.

[ANYS]

Wow. Looks impressive. That seems exactly the kind of tool I would like to see. Maybe just limiting to 1 or two elements at the same time to avoid overloading the postprocess. Looks like the integration points also contain the stresses values as they are changing color ¿isn’t it? Is that change the stress or the displacement?
Do you think this could be added to Prepomax?

[NorbertH]

those were displacements. But you could also visualize moving stress values.
Doesn’t really matter as its only limited by the capabilities of paraview.

As long as it can be coded it can be added but i think that this will be a lot of work.
I don’t know how many hours i put just into the converter. And i didn’t had to care to create the gui for visualization.

As there is no rejection i will create a post for the whole workflow. just need some freetime for it.

[NorbertH]

Here a short visualization example

[NorbertH]

as for the workflow.
you will need cubit and the opensource calculix plugin.

Before importing the .inp i first need to comment some nodes out, switch to 2d output and change the element set for the .dat output.

Cubit got a gui but can also be used headless as it basically works with the command line.

So you just need to copy paste the lines here into cubits command line. Don’t forget to fix your filepath. This will import your mesh. Then i can choose the correct calculix element type and create a job to run the file.
After the run is completed the results gets loaded and converted.

reset
import abaqus mesh geometry "/home/user/coreform/LE4_Axisymmetric_shell_internal_pressure_fine.inp" feature_angle 135
ccx block 1 element_type cax8
ccx create job name "convert" filepath "/home/user/coreform/LE4_Axisymmetric_shell_internal_pressure_fine.inp"
ccx run job 1

Now we can already take a look at the results through the job monitor. Just double click on the job in the model tree and choose the frd or dat tab to query, plot or export the results to a csv.

To view the results in paraview just click on the “result paraview” button.
Activate the multiblock inspector in paraview. This way you can choose which element, node or sidesets you want to view.

Now i will select the interesting set for the integration point and its values. To also have some contrast i will also show the default element set.

I select the points with the point trough method.
With the FindData Panel you can view the values of the integration points.

When i view the Attribute “point data” the column “ids” is the integration point id.

While the Attribute “cell data” the column “ids” should show the element id. But here its not…i need to take a look at that.

Normally this will look like this.

So now you got a way to visualize the integration points results and i got a bug more to fix

Edit:
No bug.
Just the wrong way of selection. You will need to select the integration point as elements to get the cells data for the element id.