Penalty factor abaqus vs calculix

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.

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.

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.

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.

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 ).
Check â€śkscalemaxâ€ť keyword in Manual.

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

``````** 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
Slave1,Master1
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.

Please let us know what you get from ABAQUS.

1 Like