This article is a summary of a white paper written by Ben Names (Sr. Stress Engineer) on The State of Composite Buckling in NX Nastran, as described in the attached white paper.
For a stress analyst, there are typically two important factors when designing and sizing a structure. The first is the strength capability captured by stress, the second is the static stability, which is typically handled with a buckling analysis. At Structures.Aero, where we specialize in the analysis of composite structures, we have come across a few oddities that generally only come up when running buckling on composite structures that we felt needed documenting. One particularly troubling behavior we have found is the appearance of short wave-length buckles on cored panels. The summary of my findings is that as of version 12.2 of NX Nastran (release date September 2018), the best way to run a buckling analysis when using plate elements is to use CQUAD4/CTRIA3 elements with no offsets applied to the elements.
Buckling with offsets
Although it may not be very well known, it is documented in the NX NASTRAN QRG that the offset term (ZOFFS) should not be used when differential stiffness is calculated in solution 103 (linear buckling) as “the differential stiffness calculation does not include the offset vectors” for CQUAD4 and CQUAD8 elements. Interestingly, this note is omitted from the QRG for CQUADR elements which were originally introduced to add a drilling DOF to the CQUAD4 element. Also interesting is within the realm of linear buckling, it is a common exercise to show that eccentricity does not impact the buckling eigenvalue. With that exercise informing me that the offset should matter, and the QRG informing me that the adding an offset with result in incorrect buckling eigenvalues for CQUAD4 and CQUAD8 elements, I had to investigate. To test this idea with a FEM, I created six different rectangular plate models in which offsets and element types were varied:
Metallic plate test fem buckling models with plate boundary conditions.
Each plate has the dimensions of 1”x5” with a plate thickness of t = 0.05”, and all plates are made of a generic Aluminum with E = 10.3 Msi. One of the cases I investigated was a simple Euler-column buckling solution. To that end, a unit load was applied to the plate at one end, and each plate was simply supported. Using equation 5, we can calculate the Pcr=42.4 lbs:
Analytical Euler column buckling results
The results from Nastran paint a clear picture, and quickly enumerate why using CQUADR’s are the preferred element. All three element types produce almost exact results when no offsets are applied to the model, but only CQUADR elements correctly produce an unchanged eigenvalue when using offsets.
Plate buckling results with Euler column boundary conditions
I went on to test the same plate in several different loading and boundary condition configurations and always came to the same conclusion: it appears the only way to accurately analyze buckling with plate elements is to model a structure with CQUAD4 elements at the midplane or model them with CQUADR elements. Unfortunately, that isn’t the end of the story.
Buckling with core
Another common feature within composite structures is the inclusion of core within a layup to increase a structure’s bending stiffness EI and help it to resist buckling. Unfortunately, because core is very soft this can lead a panel to have a dramatically lower transverse shear stiffness. To investigate how the effects of a soft transverse shear stiffness impact plate/beam buckling, I returned to the simply supported plate which resulted in an Euler column bucking load of Pcr = 42.4 lbs. If we model the beam as a shear deformable (Timoshenko) beam, we get a very similar buckling eigenvalue of Pcr = 42.5. What is more interesting however is how the critical buckling eigenvalue changes as the transverse shears stiffness is forced to be softer. Below is a plot of the critical buckling load vs the ratio of the transverse shear stiffness to bending stiffness:
Timoshenko beam critical buckling eigenvalue vs transverse shear bending stiffness ratio
From this simplified investigation, two things became clear:
- As the GAK decreases, Pcr also decreases
From this, it would seem that when analyzing a structure with a soft transverse shear stiffness, such as a cored panel, it is very important to accurately quantify the transverse shear stiffness as it will likely dominate the magnitude of the critical buckling eigenvalue for the structure. So the open question becomes: does NASTRAN calculate the transverse shear stiffness for CQUAD4 and CQUADR elements correctly? Unfortunately, the answer is no. Using several different methods, I found that NASTRAN under-predicted the transverse shear stiffness for CQUADR elements, often biasing their stiffness towards the softest material in the laminate stack:
FEM derived cross-sectional stiffness values
As one would expect then, the critical buckling eigenvalues for the panels are also under-predicted for cored panels using CQUADR elements:
Plate critical buckling eigenvalue using Euler column boundary conditions
The findings of this paper are a little disheartening for any stress engineer that analyzes composite structures in NX Nastran. From the first section, it is clear that only CQUADR elements handle offsets correctly. However, the second section shows that CQUADR elements have an incorrect transverse shear stiffness, so the only accurate way to run a buckling analysis on composite structures is to do so using CQUAD4 elements with no offsets. Though the results also shed a lot of light on some subtleties of bucking, it doesn’t answer the question as to whether the short wavelength buckles produced by NX Nastran on cored panels are real, and if they are what their significance is. Often when these short wavelength buckles appear, they are also very mesh dependent. For example, if one were to see these buckles on a short wavelength panel, and then refine the mesh by doubling the nodes, one would see a similar mode shape with double the number of peaks in that mode at approximately the same eigenvalue. This mesh dependent oscillating solution that does not converge is a very typical sign of numerical convergence issues. It is my recommendation that these eigenvalues be ignored, though it would be better to rule these out by addressing the validity by running a nonlinear analysis.
Should there be an opportunity to investigate the short wavelength buckles in the future, the goal of the future study would be to determine the validity of the short wavelength buckles using nonlinear analysis, and possibly models using 3D core, as well as outlining workarounds to get accurate buckling solutions without having to sacrifice the presence of offsets within the model.