Nonlinear Axial Compression (Buckling) of a Cylinder not Converging

[fethio]

Dear All,

My and my student have been trying to run a nonlinear buckling analysis on a reactor skirt, basically a cylindrical shell with holes, under axial compression in CalculiX 2.23, but I am having convergence problems (running the stock solver).

The workflow based on a previous discussion Nonlinear buckling of a cylindrical shell is as follows:

First, I worked on linear buckling. Created an INP file containing 1 buckling mode, and applied the load as compression with *CLOAD. After solving, opened the FRD file and wrote the following commands to obtain the deformed mesh:

1. datasets-displacement-entity:4 All
2. viewing-toggle add displacement
3. scal d 10000
   1. prnt se
4. send all abq
5. send all (nodeset_groups) nam

In the proceeding nonlinear buckling analysis the newly created files are used.

I have tried both load control and displacement control, with and without plasticity defined, but in all cases the analysis failed after a downwards displacement of about 1% of the cylinder height.

The geometry and group sets are in the following link: https://limewire.com/d/GrE9z#FTTW4TZXr9

The nonlinear buckling analysis is run with the below INP file. But it fails to converge at an axial compressive displacement of about 1% of the cylinders height. The culprits in my opinion are

1. the mesh being too coarse
2. the stock solver not being able to handle buckling

Any advice, comments by experienced CalculiX Users will be greatly appreciated.

Fethi

*Heading
nonlinear buckling analysis, skirt with openings and holes under load, shell mesh
*INCLUDE, INPUT=all.msh
*INCLUDE, INPUT=base_plate_center_YFIX.nam
*INCLUDE, INPUT=skir_sless_top_LOAD.nam

*BOUNDARY
Nbase_plate_center_YFIX,1,3

*MATERIAL,NAME=STEEL
*ELASTIC
210000,0.3
*DENSITY
7.8E-9

*SHELL SECTION, ELSET=Eall, MATERIAL=STEEL
25
**
****step 1, nonlinear displacement control buckling analysis
*STEP,NLGEOM,INC=1000
*STATIC
0.02,1,1e-4,0.1
*BOUNDARY
Nskir_sless_top_LOAD,3,3,100
*NODE FILE, OUTPUT=2D
U
*EL FILE
S
*NODE PRINT, NSET=Nskir_sless_top_LOAD
U
*NODE PRINT, NSET=Nbase_plate_center_YFIX, TOTALS=ONLY
RF
*END STEP

[fgr]

@fethio I’m not sure the Limewire link is working, when I try I get the following error message

[rsmith]

Looking at your displacement image, the contours follow the element boundaries. That is suspicious. How can you be sure that this displacement is not an artefact of your mesh?

I would suggest creating a regular mesh using the geometry commands in cgx.

[fethio]

@fgr Thanks for pointing out the broken link, I fixed it, should work now.

@rsmith

Looking at your displacement image, the contours follow the element boundaries. That is suspicious. How can you be sure that this displacement is not an artefact of your mesh?

The displacement contours are from the first buckling mode. I did not understand what you mean by a mesh artefect, it was produced by gmsh from a STEP file (it is not a straight cylinder with holes on it).

[jbr]

@fethio

The culprits in my opinion are

1. the mesh being too coarse
2. the stock solver not being able to handle buckling

Any advice, comments by experienced CalculiX Users will be greatly appreciated.

I reviewed the model and I agree with the mesh being too coarse. Not sure if I agree with the second point.

I would start with remeshing the model with 1st order quads. Make sure to split the geometry in a way that makes it easier for further refinement if needed. I noticed that the geometry is not fully symmetrical, and it deviates from a cylinder.

I would start first with just a simple cylinder, then add the cutouts. Follow the approach from Dr. Kraska - Non-linear Buckling - CalculiX-Examples - Created by Martin Kraska

[rsmith]

Looking at the purple color in the image, it seems to follow the edges of the elements to a degree.

In other words, if you were to do the same calculation with a different mesh, would it result in the same displacement field? If not, the displacement field could be said to be an artefact of the mesh.

[fgr]

@fethio I just deleted this content because it mostly could be misunderstood as nonsense

[johanngil]

I added some plasticity data to your model and then it converged. My file ( Unique Download Link | WeTransfer ) includes 2 alternative data sets for the *PLASTIC card to ensure convergence stability vs. the details of the *PLASTIC card. The resulting displacement U3 in axial direction at the end of the analysis shows a significant crunch of the geometry where you apply the boundary conditions. I would recommend to refine the mesh in axial direction in that area to further analyze the zone of plasticity:

The plot of the total axial force (vertical axis) vs. the U3 displacement (horizontal axis) shows that the analysis ends behind the buckling point:

[fethio]

@jbr

I would start first with just a simple cylinder, then add the cutouts. Follow the approach from Dr. Kraska - Non-linear Buckling - CalculiX-Examples - Created by Martin Kraska

@jbr thanks for the lead, we are using Prof Kraska’s examples heavily in our learning.

@rsmith thanks for pointing out the meaning of what you called an artefact of the displacement field.

@fgr I am afraid I couldn’t follow the element pass/fail test you suggested by taking a cutout from the mesh

@johanngil

I added some plasticity data to your model and then it converged.

@johanngil thanks a lot for trying out plasticity. We will follow the recommendation to refine the mesh around the cylinder mouth. As a matter of fact, the image you shared confirms very well qualitatively with the first buckling deformation mode we observed (the sort of triangulation of the opening)

[fgr]

@fethio I’m apologize for not explain myself in clear text, so I hope this will give more sense.

If I remove the fixation of your structure the data set is not even able to success the simple translation/movement load of the structure with the Nonlinear solver which tells me that something is wrong either with the dataset or the solver.

Then I got the idea that it could be a single element that caused the problem why I extracted a few elements from your structure. Trying to test the success of movement of just a single element.

In that case that the solver isn’t capable of success for a simple translation/movement load of a single element without other loads something is wrong either with the solver or the element.

Apparently, it seems to be both a single and a combination of elements, that will not succes a simple translation/movement load.

Personal for this reason I will never trust the simulation/solver/dataset before it is capable of success a simple translation/movement load.

[fethio]

@fgr could you clarify the notion of a simple translation/motion load?

Do you mean applying a load to an unconstrained structure so as to display the rigid body modes? But should not the solver exit with a failure to converge error if the structure is not properly constrained?

What is really the point in releasing all the boundary conditions and applying a load?

[fgr]

@fethio the solver can be stepwise force controlled or stepwise displacement controlled or both combined.
In your dataset you have chosen the solver to be stepwise displacement controlled.

When the solver is displacement controlled a small displacement is added to a specific part of the structure and when the forces are balanced the solver converged and moves further to next step and so on until end of the defined displacement.

So, without any external & internal forces the structure is just balanced from the beginning, the solver should just converge due to the balanced structure and move forward to the next step and so on until end of the defined displacement, so the normal result should be that the complete structure just moves from A to B as defined.
In case of this does not happen, something is wrong.

Personal I would also worry about the strange solver info message in your dataset about the flipped normal in the shown cross section

*INFO in gen3dnor: in node 15 opposite normals are defined

…

…
*INFO in gen3dnor: in node 333 opposite normals are defined

and last but not least the solver info message of MPC in common node 7705 which doesn’t exist in your dataset nor among the expanded nodes in the xxxxxx.12d file or admitted, I haven’t been able to find node 7705 anywhere.

*INFO in cascade: linear MPCs and
nonlinear MPCs depend on each other
common node: 7705 in direction 2

[fethio]

I preferred displacement control over load control for its better convergence property.

The skirt (cylinder) geometry consists of three adjacent axial parts. In the CAD modeling one of these parts somehow ended up with surface normals pointing in the direction opposite to the others. Thus the *INFO in gen3dnor. I am not sure if this should affect the accuracy of the solution, but maybe it adversely affects the model stability or solution convergence? Not sure

You’re right, the *INFO about MPCs common node: 7705 in direction 2 is strange indeed, as there is no node numbered 7705. Maybe the MPC is addressing the DOF numbering instead of the nodes.

[Disla]

1- Build a cylinder (without holes) to start.

2-Find the analytical solution of your problem and compute it.

15.4.1. Buckling of Thin Cylindrical Shells - Abbott Aerospace UK Ltd

3- Determine the right set up, minimum number of elements (In diameter and length) and type of element to achieve the best possible result compared to your analytical solution.

4- Add the holes and then face your problem with some confidence.

If you are not able to solve first three steps then I can help you.

EDITED: By the way, if you’re analyzing a skirt, you should consider the stiffness of the upper part of the structure it’s supporting. But that’s another topic. The first step is to model a basic cylinder so you have a minimally reliable starting point. In the process, you’ll see the problems inherent in this type of model—how boundary conditions, element type, and mesh density all affect the results. And that’s even before introducing imperfections or asymmetries into the geometry. I think your model is still a long way from capturing the behavior of a piece of equipment’s skirt.

[johanngil]

The “*INFO about MPCs common node: 7705 in direction 2” is a minor bug in ccx: If you harmonize the shell normals pointing for adjacent elements no longer in opposite directions then the INFO message for the non-existing node 7705 no longer appears.

[Calc_em]

@fethio

I preferred displacement control over load control for its better convergence property.

Yes, displacement control is much better for convergence. But, apart from refining the mesh, there are a few important things to consider here:

1. The shape of your initial imperfection:

It may have a large impact on how the nonlinear part goes, so you may consider changing it (e.g. replacing it with just a small dent in the middle of the shell - sort of loading imperfection) - linear buckling mode shapes are not always the best choice, especially with coarse meshes. Shells are particularly sensitive to imperfections.

Imperfections in buckling design

2. Boundary conditions - this is critical. As you can see in the thread linked in your first post, testing different BCs for the top and bottom of the shell (with or without rigid body constraints - just keep in mind that they typically cause non-convergence when used directly with shells in nonlinear analyses) can yield very different results.

3. CalculiX doesn’t have Riks solver (yet) or even stabilization for static steps, so to aid convergence, you would have to use dynamic implicit with damping or dynamic explicit, but these don’t work so well in CalculiX either. OpenRadioss would be a good choice for the latter analysis type.

4. I doubt plasticity will help here, it will only add complexity and the structure may yield too soon. It might be better to start with GNIA before proceeding to full-blown GMNIA: Process of nonlinear design on example of polygonal shell column - Enterfea

[fethio]

@Disla thanks for the guidance. Will follow the suggested path

@johanngil thanks for the tip. I was careless about this inconsistency of normals in our model, but will correct this before advancing.

[Calc_em]

@Calc_em

apart from refining the mesh

Also try different element types, but avoid S4R in buckling problems, even linear: 3DNLG-4: Lateral torsional buckling of an elastic cantilever subjected to a transverse end load