Example for 1D analysis

[bernd]

Resently theres was a Fem analysis example of the standard CalculiX Cantilever added to the Start workbench of FreeCAD. The 1D example was not added because I encoutered problems with the CalculiX 2nd order beam element. Atached the CalculiX cantilever with B32 2nd order beam elements. It gives 110 mm instead off around 87 mm deformation. I tried a few things but it seams B32 give wrong results. I made a post on CalculiX mailing list. Since this mailing list is readable only for registered users and I do not would like to send files to all mailing list members I uses this thread to post the files. With this the files are stored visible for everyone interested.

EDIT: updated example on the 25th of September 2018 (problem resists with ccx 2.11)

Link to example tutorial:
http://www.freecadweb.org/wiki/index.ph ... tilever_3D

Link to example in ccx manual (the manual is outdated, but the example is still up to date):
http://web.mit.edu/calculix_v2.7/Calcul ... node7.html



EDIT: the file provided here has wrong results in the result object. Rerun with ccx version 2.17 and higher to obtain correct results. Find a FreeCAD file with correct results in this post: https://forum.freecadweb.org/viewtopic. ... 66#p518866

[bernd]

attached my post on the mailing list and the files

Hi,

1) would it be possible to do a parameter study? Increase the number of B31,
B32 and B32R elements and check what value it is converging to?

2) can you decrease the cross section of the beam and repeat the calculation
(the beam is not really slender)?

Best Greetings,

Guido

Am Freitag, 3. Juni 2016 19:36:42 schrieb Bernd Hahnebach calculix@b75.ch
[CALCULIX]:
> Hi Calculix user and developer,
>
> Calculation where made with ccx 2.10
>
> I encountered a problem with beam b32 elements which I'm not able to
> explain. It is about the CalculiX cantilever example from the manual.
> The very first example. I modeled the beam with B31 and B32 elements.
> By the use of 10 B31 elements (11 nodes) I get a very good result off
> 88 mm deflection. I get the same with shell and solid elements.
>
> But with 5 B32 elements (11 nodes) I get nearly 20 % to more
> deformation than I should get (102mm). It gets even worse, if I
> increase the B32 element count the deflection gets highter means I do
> not converge to the real result but to another one?
>
> Could someone bring light into this?
>
> Kind regards bernd

beam_b31.inp

Code: Select all

*Node, NSET=Nall
1,     0, 500, 500
2,   800, 500, 500
3,  1600, 500, 500
4,  2400, 500, 500
5,  3200, 500, 500
6,  4000, 500, 500
7,  4800, 500, 500
8,  5600, 500, 500
9,  6400, 500, 500
10, 7200, 500, 500
11, 8000, 500, 500

*Element, TYPE=B31, ELSET=Eall
21, 1, 2
22, 2, 3
23, 3, 4
24, 4, 5
25, 5, 6
26, 6, 7
27, 7, 8
28, 8, 9
29, 9, 10
30, 10, 11


*ELSET,ELSET=Mat0Beam0
Eall

*NSET,NSET=FemConstraintFixed
1,

*MATERIAL, NAME=MechanicalMaterial
*ELASTIC
210000 ,  0.300

*BEAM SECTION, ELSET=Mat0Beam0, MATERIAL=MechanicalMaterial, SECTION=RECT
1000 , 1000 


******************************************
*STEP
*STATIC

*BOUNDARY
FemConstraintFixed,1
FemConstraintFixed,2
FemConstraintFixed,3
FemConstraintFixed,4
FemConstraintFixed,5
FemConstraintFixed,6

*CLOAD
11,3,-9.0E+06


*NODE FILE, OUTPUT=2d
U
*EL FILE
S, E
*NODE PRINT , NSET=Nall
U 
*EL PRINT , ELSET=Eall
S 

*END STEP 

beam_b32.inp

Code: Select all

*Node, NSET=Nall
1,     0, 500, 500
2,   800, 500, 500
3,  1600, 500, 500
4,  2400, 500, 500
5,  3200, 500, 500
6,  4000, 500, 500
7,  4800, 500, 500
8,  5600, 500, 500
9,  6400, 500, 500
10, 7200, 500, 500
11, 8000, 500, 500

*Element, TYPE=B32, ELSET=Eall
21, 1, 2, 3
22, 3, 4, 5
23, 5, 6, 7
24, 7, 8, 9
25, 9, 10, 11


*ELSET,ELSET=Mat0Beam0
Eall

*NSET,NSET=FemConstraintFixed
1,

*MATERIAL, NAME=MechanicalMaterial
*ELASTIC
210000 ,  0.300

*BEAM SECTION, ELSET=Mat0Beam0, MATERIAL=MechanicalMaterial, SECTION=RECT
1000 , 1000 


******************************************
*STEP
*STATIC

*BOUNDARY
FemConstraintFixed,1
FemConstraintFixed,2
FemConstraintFixed,3
FemConstraintFixed,4
FemConstraintFixed,5
FemConstraintFixed,6

*CLOAD
11,3,-9.0E+06


*NODE FILE, OUTPUT=2d
U
*EL FILE
S, E
*NODE PRINT , NSET=Nall
U 
*EL PRINT , ELSET=Eall
S 

*END STEP 

[bernd]

attached my results I'm gone show the CalculiX devloper ...

summary screen

screen.jpg (113.48 KiB) Viewed 6203 times

summary text

summary.txt

(1.43 KiB) Downloaded 124 times

calculix input files zipped

beam-test-calculiations.zip

(29.84 KiB) Downloaded 142 times

The FreeCAD file the CalculiX files are created with (for each mesh an analysis)
The 2D and 3D example could be loaded from Start workbench

beam-various-line_meshes.fcstd

(58.95 KiB) Downloaded 136 times

[fandaL]

2) can you decrease the cross section of the beam and repeat the calculation
(the beam is not really slender)?

When changing only the cross section in input files (and B32 to B32R), one get same results for 10x10 cross section for all element types as an analytical solution gives. When the slenderness increases quadratic types B32 and B32R seem to be more sensitive and become more flexible than B31.

Deformations D3

deformations_depending_on_cross_section.png (3.55 KiB) Viewed 6153 times

[bernd]

What version of ccx have you used?

I still think it is an error in the B32 element. Even if it only apply on not slender beams. It should not be that far away. B32 should be not more stiff than B31 for the same node count?! What was found on mailing list is.

- ccx 2.8.2, 2.9 and 2.10 seams to give wrong results
- ccx 2.6 and 2.7 seam to give correct results

I'm still waiting for a response from the developer.

bernd

[fandaL]

What version of ccx have you used?

I forgot.
ccx 2.8.2 on Windows 7

[bernd]

I will copy the post from CalculiX mailing list here since it is not readable for public ... One was in the second topic post ...

Hi Guido,

thanks for your reply. Find my answers inside your post.

Kind regards bernd

1) would it be possible to do a parameter study? Increase the number of B31,
B32 and B32R elements and check what value it is converging to?

sure. I do not want to spam the mailing list with lots of files. Find
them in a forum thread here
viewtopic.php?f=18&t=16044&p=127089#p127089

2) can you decrease the cross section of the beam and repeat the calculation
(the beam is not really slender)?

What would be a good cross section in your opinion. I could provide
results for a slender beam too.

Hi,
it seams that there is an issue with the newer versions of calculix.
If I compute the beam_b32 example with ccx 2.10 I get a deflection of
102mm. However, if I use ccx 2.7 (or 2.61) I get a deflection of 85mm.
I cannot test with 2.8 or 2.9 but it seams, that the internals of the b32
have been changed somewhere between these versions.
br
Jakob

Just tested,

ccx 2.10, ccx 2.9, ccx 2.8.2 gives wrong results for me.

Damn I deleted other old binaries of ccx just a week ago.

Bernd

Found a ccx 2.6.1 in my trash

With this binary I get results I would expect. Means I can confirm what Jakob wrote. Between 2.7 and 2.8.2 seams to happen some changes which breaks the B32 elements.

bernd

[bernd]

and Guidos answer ...

Hi Bernd,

I found the reason for the different behavior:

up to version 2.7 rotational boundary conditions were simulated with knot
MPC's, in version 2.8 and higher they are similated with mean rotation MPC's.

Ok, this is the short version of the explanation. Now the long version:

- all information is actually contained in sections 6.2.13 and 6.7.3 of the
User's Manual.

- beam elements are expanded into 3D-elements. 3D elements only have
translational DOFs. Consequently, the rotational DOFs of the beam elements
have to be simulated using MPCs.

- up to version 2.7 rotational boundary conditions (and moments) were
simulated by applying expandable rigid body MPCs (= knots). All nodes
belonging to a knot MPC behave as if they were part of a monolytic structure
characterized by 3 translations, 3 rotations and 3 stretch parameters.

- in version 2.8 and higher the rotational boundary conditions (and moments)
are simulated by mean rotation MPCs. These conditions are not as strong as a
knot MPC. For instance, if the rotation about the x-axis is to be zero, the
mean rotation MPC specifies that the mean rotation about the x-axis of all
nodes belonging to the MPC has to be zero. The rotation of any individual
node does not have to be zero.

- using B31 elements the constrained section is expanded into 4 nodes. The
section is constrained for all 6 DOFs. Up to version 2.7 a rigid body is
created, i.e. all 12 DOFs (4 nodes x 3 DOFs) are reduced to 9 DOFs (the three
translations, three rotations and three stretch parameters of the body), of
which six are fixed to zero (translations and rotations): three DOFs (the
stretches) are left free. In version 2.8 and higher mean rotation MPCs are
introduced for the rotations (and mean translation for the translations),
i.e. 6 MPCs are generated: consequently 6 DOFs are left undetermined. So 6
DOFs left free (Version 2.8) vs. 3 DOFs left free( Version 2.7). This,
however, does not seem to have much effect on the results.

- using B32 elements the constrained section is expanded into 8 nodes, i.e. 24
DOFs. Up to version 2.7 all these DOFs are reduced to 9 of which 6 are fixed:
three DOFs are left free. In version 2.8 and higher 3 mean translation and 3
mean rotation MPCs are generated, leaving 18 DOFs undetermined. This leads to
a much more flexible behavior of the beam.

- if you model the beam in 3D from the beginning, the mean rotation approach
corresponds to removing only the minimum number of DOFs at the fixed section
in order to avoid any rigid body motion.

- the mean rotation approach has increased the agreement of the results
between CalculiX and ABAQUS, especially for shells (I do not recall whether a
lot of calculations were performed for beams). Rigid bodies introduced too
many constraints, leading to too small deformations.

- The theory of the mean rotation MPC is in my book (Wiley 2004).

- Knots are still introduced in the case of sharp angles in shells and/or
beams. Consequently, you may be able to reproduce the results of version 2.7
and older by adding a beam under 90 degrees at the fixation.

Hope this helps,

Guido

Hi all,

a slight addition: it may be that in version 2.7 the knot was even more
restricted, it had probably only 7 DOFs (3 translations, 3 rotations and one
isotropic expansion) instead of 9 now.

Best Greetings,

Guido

[bernd]

and Martin Kraska did say something in this regard too.

Hi Guido,

this topic is a good occasion to revisit the question, how the drastic differences in the convergence behaviour of the cantilever example

https://github.com/mkraska/CalculiX-Exa ... r/Streifen

might be explained.

The objective is to deform a cantilever beam by a 90° prescribed rotation at the tip. In the example, I compare beam and shell models with *boundary to apply the rotation. The beam model doesn't even reach 12% while the shell model succeeds.

Another two models use solid elements with the rotation applied via mean rotation (reaches 28%) and rigid body constraint (reaches 74%).

It seems that the împlicit mean rotation mpc's as used in the shell model are more efficient than the ones exposed by the *rigid body keyword.

Mit freundlichen Grüßen/Best regards

Prof. Dr.-Ing. Martin Kraska

Hi Martin,

I guess the differences are due to the different boundary conditions.
Nonlinear problems can react quite strongly to the way boundary conditions
are applied.

By the way, you can print the MPC's generated in CalculiX (all of them: the
user's and the internal ones) by de-commenting the call to writempc at the
beginning of cascade.c and recompile.

Best Greetings,

Guido

[bernd]

issue #3056