Simple beam bending FEM vs. theory comparison error 4

[rcsingh723]

I am having some trouble getting accurate values for max. stress and max. deflection for an FEM model of a simply supported I-beam in bending. I am using I-DEAS 12 NX to solve the model. I'm pretty new to finite element analysis as a whole, but I know that getting any more than 10% error (relative to Euler beam theory) in max deflection or stress for a 60x10x8 in. should not happen. Strangely enough, using the same boundary conditions and meshing techniques for a 60x10x8 in. rectangular solid, I find the results I was expecting. Here are some more details:

-For each geometry, there are 3 models. One uses beam elements, one uses solid (hexa) elements, and one uses shell (quad) elements.

-The boundary conditions (pin on one side, roller on the other) are applied to the mesh, using constraint elements on each face, with the node in the middle of the face being the independent node (where the b.c.'s are applied). Similarly, the applied force of 1000 lbs was modeled as a series of point forces applied along the neutral axis of the beam at its center.

-As for mesh density, I modeled the rectangular prism using an element size of about 2.5 in x 2.5 in x 2 in for the hex elements, just to give an idea. The I-beam model required a finer mesh for model "convergence" (there seem to be singularities at the b.c.'s and loads for the I-beam, regardless), so we're talking element sizes of about 1 in x .5 in x .5 in.

Any thoughts on what I (or possibly I-DEAS) is doing wrong? I'm thinking Saint-Venant's principle, but the magnitude of the max stress error and max deflection error in the I-beam model is much bigger than the errors in the rectangular prism.

Hope this wasn't too ridiculous to read.

Thanks,

Raymond C. Singh

 

[GBor]

If a problem this simple is off this far, check your units to make sure your inputs are consistent.

 

[rb1957]

not clear on your "beam" model ... the solid (hexa) model is easy enough to understand, your shell model presumably has elements at the mid-thickness planes ... be aware that you need to support the out-of-plane freedoms or use plate elements (with bending stiffness). but for your "beam" model are you modelling the flanges with beam elements and the web with a shell, or are you modelling the upper cap as a beam and connecting the two caps with a shear-only shell ?

either way your results should match theory very well, and any model with singularities has be judged very suspiciously (ie find out what the singularities are and mix them !)

 

[rcsingh723]

Thanks for replying. I looked at the .dat file of the model from NASTRAN, and everything checks out...geometry, forces, DOF, moments of inertia, etc. I don't think that's it.

 

[rcsingh723]

rb1957,

For the beam element model of the I-beam, I am using the beam section property tool in I-DEAS to create an I cross-section with the proper dimensions. The beam elements are oriented along the length of the beam. Also, what do you mean by 'mixing' singularities? Thanks for the help.

 

[fetraining]

Raymond,

I worked in tech support for several nastran providers and this came up every once in a while.

The finite element models will include bending stiffness terms and shear stiffness terms, whether they are created using beams, shells or solids. (more later on relative accuracies).

Euler beam theory does not take account of shear stiffness. So if you used the equations out of Roark or similar you will get a stiffer result. For example for a cantilever beam Roark has the term:

d = PL^3/(3EI) for bending stiffness, but the additional term is PL/(GA) to account for shear stiffness.

Full solution is d = PL^3/(3EI) + PL/(GA)

I don't have the shear term for simply supported beam to hand, but I'm sure a good textbook, web seerch or other posts will reveal it.

with regards to variations in element types:

Beams should give excellent results as long as there are sufficient to follow the required curvature (high order buckling or vibration modes need more).

Shells will give good results, but I recommend two or three shells minimum in the web to get a good through depth shear distribution. The end constraints are more complicated however. A method here is to link all grids at each end with a rigid element. The single independent grid is at the centroid of the cross section. This is simply supported, or fixed dependent on the test. This is an attempt to model the 'theory' that all sections remain plane.

Solid elements require many more elements through the web depth to get a good answer. the performance varies amongst types, but assuming you have a standard Nastran type brick ( 8 noded Hexa) then at least 8 elements through web of the I beam. If you have a higher order brick then it would be a good idea. End constraint issues are similar to the shell. You want to strap all the grids at each end face to a single supported grid with rigid elements. Again it keeps sections planar and replicates what the theory is assuming.

Don't use Tetrahedral elements unless you are really forced to, and then never used 4 noded Tets. These will give bad results. If you have to use 10 noded Tets then make the mesh in the beam very fine and as regular as possible. Avoid 'jaggies' which destroy accuracy (high aspect ratio icicle type shapes).

I hope this helps. I think the fundamental reason for the discrepancy is in the hand calculation. But the end conditions are alo tricky in shells/solids.

There always comes a point where trying to benchmark against 'theory' is a diminishing return because of the limitations inherent in the theoretical assuptions versus real world applications.

regards,

Tony

 

[rcsingh723]

Thanks for the reply. This is where things get strange. With respect to stiffness, the max stress and max deflection values are all greater in magnitude than what is given by theory. During mesh convergence, I saw that though the models were converging, rather than diverging such that max. stress goes to infinity, they were approaching their asymptotes from lower stresses, or rather, they were converging, but converged to a more pliant model solution.

The beam element model gave spot-on results for stress across the board, but strangely enough, the deflection results were just as bad as the shell and solid mesh grids (i.e. max deflection error was roughly 36% relative to theory in all three for an I-beam 60 inches long, while for the same beam but with a rectangular cross-section, the error was about 5%).

I figured theory would introduce some error (which is actually the problem that I was originally trying to answer...at what length does Euler beam theory no longer apply), but my manager as well as coworkers have said that they never had this problem. I'll use Roark and see what I come up with, but I think there's something else going on. Even using basic Euler beam theory, I shouldn't be getting 40% error in max. deflection or stress values.

 

[fetraining]

Raymond,

if the beams are showing greater deflection than Euler theory, then I still suspect the missing shear term in the theory. One way to test this is to manually increase the shear stiffness factor on the beams. The idea is that if you set this to a very high value it replicates the missing (i.e. infinite) shear term in the Euler theory. i have done this in the past to match euler theory and demonstrate thispoint.

You mention a Nastran input file - is this what I-DEAS is writing out? If so look at the K1, K2 terms on the PBEAM cards (forgive my age). They default to 1.0, set them to about 20 and see if you get closer to the Euler theory with its limitations.

regards,

Tony

 

[rb1957]

i don't think the problem is shear deflections (since he's putting a pretty small shear load onto a pretty big section).

i wonder about the support constraints ... i understand that you applied x,y,z restraint at one end and y,z and the other (x being axial). at the end with the x- restraint, did you apply it to every node on the surface ? (if so, this would also create a moment support).

for the beam model ... ok, you modelled the section as a single beam element; did you restrain the torsion (Rx) freedom ?

 

[GBor]

I have a more fundamental question. Your beam is 60 inches long and 8 deep? I believe beam equations generally assume that your length is much greater than your depth. By "much greater", I think that would require at least a 10x factor.

Application of beam theory to this situation may be stretching the limits of the equations a bit...

 

[rcsingh723]

fetraining:

I got I-DEAS to export the model file as a .dat file to use with MSC/Nastran. I will fiddle around with that later today and see what I get.

rb1957:

What I did was use a constraint element to tie the node at the end face's midpoint to all the other nodes on the face. Then, I applied the x,y,z displacement restraints (and x,y rotational restraints) to the independent node (midpoint). This is different than directly applying the restraints to every node on the face. I think the latter is what you're referring to, which isn't what I'm doing. It seems like the deformation shape is sound.

For the beam element model, I'm not sure I follow. Perhaps I said it wrong, but I used about 10-15 elements, with each beam element being that I-beam cross section w/ its respective properties.

GBor:

That's what I thought as well, when I first conducted this study. An aspect ratio of 10 or greater is what I usually use as well. The error shouldn't be this large though, regardless of when beam theory starts to fall apart. I mean, using these values, beam theory isn't "accurate" until an aspect ratio of at least 20. I will ask my manager and see if he can set up a quick model and compare the results.

 

[btrueblood]

I've seen the same thing, but given as 15x L/d, when the section has thin webs.

But, going thru some numbers for an I-section with 1" uniform webs, 10" deep, 8" wide, 60" long, and a rect. section 10" deep, 8" wide, and 60" long, the Roark shear formulas give ys/yb (ratio of shear deflection to bending deflection) of about 0.1% for the rectangular section, and 0.04% for the I-section.

Agree with Rb, the problem is likely in the boundary conditions. Check and see if there are any out-of-plane deflections or rotations. Plot the y-deflection vs. span and see if you have a smooth curve, or if perhaps there is a big jump near the ends due to stress concentrations at your assumed fixed conditions at the midplane of the beam.

 

[rcsingh723]

Yeah, checked those. While there are out-of-plane deflections, they are a couple orders of magnitude smaller than the bending deflections, as predicted by btrueblood's data. There is in fact a large jump in the % error for the stress data around the loads. I guess what I was asking was why the magnitude of the errors was so large for the I-beam, and acceptable for the rectangular section. Another thing I noticed just now, was that if I increase the thicknesses of the flanges/web (from 1 in. to 2 in.), the models produce less error. Perhaps the accuracy of the model at the b.c.'s has something to do with the moment of inertia of the beam at hand?

 

[GBor]

Boundary conditions are usually a good place to go...

 

[fetraining]

Raymond,

I am still convinced the shear stiffness term missing from Roark is the culprit. Particularly as you have a short stubby beam where this is an important term.

I set up a nastran beam model with length L = 60 ins, center load W = 10,000 lbs. I beam 10 ins deep , 8 ins wide. E = 1.07e7 psi.

central deflection in model = .01573 ins.

Checked with Roark, d max = Wl^3/(48EI) = 0.011428 ins
(I is 386 ins^4)

So FE model appears to be 37% too flexible.

But Roark ignores shear stiffness term. We can simulate this in Nastran Beam by putting K shear correction term to a large number (1e6). This effectively knocks out the shear stiffness term and simulates Roark.

Rerun and the nastran central deflection = .01143 ins. this exactly matches Roark calulation.

So the Roark formula misses the shear stiffness term which accounts for 37% of the delection.

I did the same thing with the rectangular 10 by 10 beam. I is now 833 in^4.

Roark gives d = 0.00505 ins.
Nastran with default K term = 0.00549 ins
Nastran with high K term = 0.00505 ins.

So in this case the effect of the missing shear stiffness in Roark is less (10%)because the flexural stiffness term dominates.

The percentages tie up with what you originally noted as the discrepancies.

I can send you the Nastran decks if you want.

The shell and solid FE models will exhibit the same trend, they will have a shear stifness term inherent. If the mesh density is poor then they may struggle to get the term accuratly, but it will not be infinite as in the case of Roark calculation.

regards,

Tony

 

[rb1957]

going back aways, raymond, you said that you tied all the nodes one the end surface to a single node and contrained this in x, y, and z. i think this is constraining each node on the face; check the displacements of extreme nodes, check the stress on the extreme elements. if i'm right you'll see a significant moment coming out, the displacement of extreme nodes will be zero, and the stresses of extreme elements will have different signs.

 

[rcsingh723]

I went over my model set thoroughly with my manager, and after talking to him about some of the suggestions made by you guys, he deduced that the varying errors in deflection/stress are in fact a by-product of the inaccuracy of beam theory, and not any eccentricities or mistakes in the model itself. I believe you are on to something, Tony....I'm in the process of tabulating data, but I'll go back and change K1 and K2 and see what happens. I don't necessarily need the NASTRAN model files, but I guess they'd be a good means for comparison, so if you could send those over, that would be fine. Thanks very much for your help, guys!

 

[fetraining]

Raymond,

I'm new to the forum so this is my attempt to upload the Nastran file. it is called 'beam_shear_test.nas'.

Both beam property types (I and rectangular) are in the deck, plus annotated comment/uncomment lines for k.

regards,

Tony

 

[rb1957]

not that it matters much, but your PBEAM 3 is for a 10x10 section

 

[fetraining]

rb1957,

thanks for spotting that - yes I mislabeled it as

$ rectangular solid 10 by 8 select PBEAM 3

It should read ' 10 by 10 '

just to confirm the model file and hand calcs are for a 10 by 10.

regards,

Tony