Penalty factor abaqus vs calculix

1

Does anyone know how calculate penalty stiffness in abaqus HARD contact?
I tried to do CNORM/CDISP*CAREA but I didnt get good results.
I would like to know to perform a correlation between calculix and abaqus penalty factor.

2

https://docs.software.vt.edu/abaqusv2022/English/SIMACAEITNRefMap/simaitn-c-contactdiagnostics.htm

3

Unfortunately, this shows everything except for penalty stiffness unless it’s non-default. The details of some calculations done internally by Abaqus (like the representative underlying element stiffness) are not documented and only briefly mentioned in the manuals. Even support engineers don’t know them.

4

Yes I agreee. The manual does not talk deeply about how calculate It.
I could calculate It like CNORMF/(COPENĂ— CNAREA).
But I am not sure.

5

I don’t think this is the right way. You should try a variation study in Abaqus for a well-known problem, say the 2D Hertz contact problem (CalculiX-Examples/Contact/Hertz_2D at master · calculix/CalculiX-Examples · GitHub). Then you can establish the values by checking certain conditions and then compare that to Calculix if you want.

6

1D: (Hooks law) : F=K*e

2D: F/A=K/A*e

F/A =K’*e : [N/m2] = [Pa] = [ Pa/m] * [m].

K’= Contact Stiffnes= F/A/e=P/e

That simplification (your formula) is only valid for a uniform displacement of the whole surface in the normal direction and all the undelaying springs in the contact having the same stiffness.

There is a way to satisfy that condition. Set up a contact between a Rigid Master Surface and a Slave which should be part of a unique HEX8 element.

Impose a uniaxial compression test to assure all the springs get the same load.

I have imposed the compression by means of a displacement driven on the master surface. In my previous test with load driven, stiffness was artificially relaxing, I guess to help the convergence (Thanks for the hint Victor :ok_hand:).
Check “kscalemax” keyword in Manual.

Then your formula is ready to be valid.

If we do that, I think you could extract the Stiffness used by Abaqus for a Hard Contact. Do not assume the Stiffness is constant (linear overclosure response). Plot the whole path (CSTR,CDIS) and Stiffness should be the slope.

As an example, I have set up two cubes as descrived. One contact is HARD and the second is Linear (with the expected 50E for the Stiffnes).
Cube of 1mx1mx1m

imagen
imagen1153Ă—640 11.1 KB

** Generated by Mecway 24
*NODE
1,-1.2,0.6,0
2,1.3,-0.6,0
3,1.3,0.6,0
4,-1.2,0.6,-0.2
5,1.3,0.6,-0.2
6,1.3,-0.6,-0.2
7,0,0,-0.15
8,-1.2,-0.6,0
9,0.1,0.6,-0.2
10,-1.1,-0.5,0
11,1.2,-0.5,0
12,1.2,0.5,0
13,-1.1,0.5,0
14,-0.1,0.5,0
15,0.2,0.5,0
16,-0.1,-0.5,0
17,0.2,-0.5,0
18,-1.1,0.5,1
19,-1.1,-0.5,1
20,-0.1,-0.5,1
21,-0.1,0.5,1
22,0.2,0.5,1
23,0.2,-0.5,1
24,1.2,-0.5,1
25,1.2,0.5,1
26,-1.2,-0.6,-0.2
27,0,0.6,0
28,0,0.6,-0.2
29,0,-0.6,0
30,0,-0.6,-0.2
31,0.1,-0.6,0
32,0.1,0.6,0
33,0.1,-0.6,-0.2
*ELEMENT,TYPE=C3D8R
1,17,23,24,11,15,22,25,12
2,27,29,8,1,28,30,26,4
3,13,18,19,10,14,21,20,16
4,3,2,31,32,5,6,33,9
*NSET,NSET=NSET_CONTACT
11
12
15
17
*ELSET,ELSET=C3D8
1
3
*ELSET,ELSET=Bases
2
4
*ELSET,ELSET=Hard
3
*ELSET,ELSET=Linear
1
*SURFACE,NAME=Master2
2,S1
*SURFACE,NAME=Slave1
1,S6
*SURFACE,NAME=Slave2
3,S6
*SURFACE,NAME=Master1
4,S1
*MATERIAL,NAME=Steel
*ELASTIC,TYPE=ISOTROPIC
100000000000,0.3
*DENSITY
7850
*SOLID SECTION,ELSET=Bases,MATERIAL=Steel
*SOLID SECTION,ELSET=Hard,MATERIAL=Steel
*SOLID SECTION,ELSET=Linear,MATERIAL=Steel
*BOUNDARY
2,2,,0
8,1,,0
8,2,,0
18,3,,0
19,1,,0
19,2,,0
19,3,,0
20,2,,0
20,3,,0
21,3,,0
22,3,,0
23,1,,0
23,2,,0
23,3,,0
24,2,,0
24,3,,0
25,3,,0
29,2,,0
31,1,,0
31,2,,0
*AMPLITUDE,NAME=A_1
0,0
1,5E-05
*TIME POINTS,NAME=timepointsname,GENERATE
0,1,0.01
*SURFACE INTERACTION,NAME=SI_6
*SURFACE BEHAVIOR,PRESSURE-OVERCLOSURE=LINEAR
5E12
*SURFACE INTERACTION,NAME=SI_7
*SURFACE BEHAVIOR,PRESSURE-OVERCLOSURE=HARD
*CONTACT PAIR,INTERACTION=SI_6,TYPE=SURFACE TO SURFACE,ADJUST=1
Slave1,Master1
*CONTACT PAIR,INTERACTION=SI_7,TYPE=SURFACE TO SURFACE,ADJUST=1
Slave2,Master2
*STEP,NLGEOM=YES,INC=110,AMPLITUDE=STEP
*STATIC,SOLVER=PARDISO
0.01,1,0,0.01
*BOUNDARY,AMPLITUDE=A_1
1,3,,1
*BOUNDARY,AMPLITUDE=A_1
2,3,,1
*BOUNDARY,AMPLITUDE=A_1
3,3,,1
*BOUNDARY,AMPLITUDE=A_1
8,3,,1
*BOUNDARY,AMPLITUDE=A_1
27,3,,1
*BOUNDARY,AMPLITUDE=A_1
29,3,,1
*BOUNDARY,AMPLITUDE=A_1
31,3,,1
*BOUNDARY,AMPLITUDE=A_1
32,3,,1
*NODE FILE,GLOBAL=YES,TIME POINTS=timepointsname
U,RF
*EL FILE,TIME POINTS=timepointsname
S,NOE,ENER
*CONTACT FILE,TIME POINTS=timepointsname
CSTR,CDIS
*NODE FILE, NSET=NSET_CONTACT
CELS
*END STEP

This is my Pressure Overclosure response. Linear, with constant slope (Stiffness), and its value is the exact input value 5E12 for both contacts.

imagen

Please let us know what you get from ABAQUS.

1 Like