Convergence issue in 2D contact problem

1

Dear CalculiX community

During Covid lockdown in May 2020, I found the opportunity to explore CalculiX. I performed several basic analyses and shared them from my GitHub account here. Currently, I am working on a two-dimensional problem with contact where I have two bodies, one modelled with axisymmetric elements and the other using plane stress elements. Only half of the geometry is modelled due to symmetry. Both bodies are under centrifugal loading. I tried several things but I could not manage getting the model to converge. I am using bconverged distribution for Windows upgraded to version 2.17.

I copied my files below. Do you have any comments or suggestions? I appreciate any help.

Regards
Süleyman

.fbd

#
## Geometry
#

pnt p1 0.000000 60.000000 0
pnt p2 2.000000 60.000000 0
pnt p3 2.000000 59.000000 0
pnt p4 3.000000 59.000000 0
pnt p5 2.233956 58.357212 0
pnt p6 3.519531 56.825124 0
pnt p7 3.136509 56.503730 0
pnt p8 3.569521 56.253730 0
pnt p9 2.609641 54.591168 0
pnt p10 0.000000 54.591168 0
pnt p11 0.000000 56.253730 0
pnt p12 0.000000 56.825124 0
pnt p13 0.000000 58.357212 0
pnt p14 0.000000 59.000000 0
pnt p15 2.350826 58.217932 0
pnt p16 2.870357 56.043055 0
pnt p17 3.636401 56.685843 0
pnt p18 3.736382 55.543055 0
pnt p19 2.436255 55.291168 0
pnt p20 3.302280 54.791168 0
pnt p21 2.436255 54.291168 0
pnt p22 2.000000 54.291168 0
pnt p23 0.000000 54.291168 0
pnt p24 2.200000 58.541379 0
pnt p25 2.548909 59.500000 0
pnt p26 4.548909 59.500000 0
pnt p27 8.013011 56.685843 0
pnt p28 7.811507 55.543055 0
pnt p29 7.238111 52.291168 0
pnt p30 3.736382 52.291168 0
pnt p31 3.302280 52.291168 0
pnt p32 2.436255 52.291168 0
pnt p33 2.000000 51.291168 0
pnt p34 2.000000 52.291168 0
pnt p35 1.292893 51.998275 0
pnt p36 1.000000 51.291168 0
pnt p37 3.000000 27.828427 0
pnt p38 1.000000 27.828427 0
pnt p39 1.585786 26.414214 0
pnt p40 5.000000 23.000000 0
pnt p41 5.000000 20.000000 0
pnt p42 0.000000 20.000000 0
pnt p43 0.000000 23.000000 0
pnt p44 0.000000 26.414214 0
pnt p45 0.000000 27.828427 0
pnt p46 0.000000 51.291168 0


# Element type parameterization.
# valu axisEtyp qu8cr # Eight-node axisymmetric element with reduced integration.
# valu psEtyp qu8sr # Eight-node plane stress element with reduced integration.
valu axisEtyp qu4c # Four-node axisymmetric element.
valu psEtyp qu4s # Four-node plane stress element.


# Line Divisions
valu div1 20
valu div2 6
valu div3 10
valu div4 60
valu div5 20
valu div6 15
valu div7 30
valu div8 15
valu div9 60
valu div10 30
valu div11 6
valu div12 12
valu div13 6
valu div14 20
valu div15 80
valu div16 10
valu div17 20
valu div18 12


# Blade lines
line l1 p1 p2 div1
line l2 p2 p3 div2
line l3 p3 p5 p4
line l4 p5 p6 div4
line l5 p6 p8 p7
line l6 p8 p9 div6
line l7 p9 p10 div1
line l8 p10 p11 div6
line l9 p11 p12 div5
line l10 p12 p13 div4
line l11 p13 p14 div3
line l12 p14 p1 div2
line l13 p14 p3 div1
line l14 p13 p5 div1
line l15 p12 p6 div1
line l16 p11 p8 div1

# Disk lines
line l17 p15 p17 div9	
line l18 p17 p18 p16
line l19 p18 p20 div11	
line l20 p20 p21 p19
line l21 p21 p22 div13
line l22 p22 p23 div14
line l23 p15 p24 div8	
line l24 p24 p25 div7	
line l25 p25 p26 div8
line l26 p26 p27 div9	
line l27 p27 p28 div10	
line l28 p28 p29 div7	
line l29 p29 p30 div7
line l30 p30 p31 div11	
line l31 p31 p32 div12	
line l32 p32 p34 div13 
line l33 p34 p35 p33
line l34 p35 p36 p33
line l35 p36 p38 div15	
line l36 p38 p39 p37
line l37 p39 p40 div17	
line l38 p40 p41 div18	
line l39 p41 p42 div7	
line l40 p42 p43 div18	
line l41 p43 p44 div17	
line l42 p44 p45 div16	
line l43 p45 p46 div15	
line l44 p46 p23 div14	
line l45 p15 p26 div7
line l46 p17 p27 div7
line l47 p18 p28 div7
line l48 p18 p30 div7 
line l49 p20 p31 div7
line l50 p21 p32 div7
line l51 p22 p34 div7
line l52 p23 p35 div7
line l53 p46 p36 div7
line l54 p45 p38 div7
line l55 p44 p39 div7
line l56 p43 p40 div7


# Line divisions.
div l3 div3
div l5 div5
div l18 div10
div l20 div12
div l33 div14
div l34 div14
div l36 div16

# Line set definitions.
seta lset1 l l1 l2 l13 l12
seta lset2 l l13 l3 l14 l11
seta lset3 l l14 l4 l15 l10
seta lset4 l l15 l5 l16 l9
seta lset5 l l16 l6 l7 l8

seta lset6 l l23 l24 l25 l45
seta lset7 l l17 l45 l26 l46
seta lset8 l l18 l46 l27 l47
seta lset9 l l48 l47 l28 l29
seta lset10 l l19 l48 l30 l49
seta lset11 l l20 l49 l31 l50
seta lset12 l l21 l51 l32 l50
seta lset13 l l22 l51 l33 l52
seta lset14 l l44 l52 l34 l53
seta lset15 l l43 l53 l35 l54
seta lset16 l l42 l54 l36 l55
seta lset17 l l41 l55 l37 l56
seta lset18 l l40 l56 l38 l39

seta bcon l l3 l4 l5
seta dcon l l23 l17 l18

# Create surfaces based on line sets.
surf s1 lset1
surf s2 lset2
surf s3 lset3
surf s4 lset4
surf s5 lset5
surf s6 lset6
surf s7 lset7
surf s8 lset8
surf s9 lset9
surf s10 lset10
surf s11 lset11
surf s12 lset12
surf s13 lset13
surf s14 lset14
surf s15 lset15
surf s16 lset16
surf s17 lset17
surf s18 lset18

# Make sure that all surfaces share a common direction.
flip s1 auto
flip s18 auto
flip all

# Surface set definitions for meshhing.
seta psurf s s1 s2 s3 s4 s5
seta axisurf s s6 s7 s8 s9 s10 s11 s12 s13 s14 s15 s16 s17 s18


#
## Meshing.
#

# Element type is assigned to surfaces.
elty psurf psEtyp
elty axisurf axisEtyp

# Mesh the geometry.
mesh all
# Mesh is plotted.
plot m all

# Contact Surface-to-surface
seta sbcon l bcon
comp sbcon do
comp sbcon do
send sbcon abq sur
seta sdcon l dcon
comp sdcon do
comp sdcon do
send sdcon abq sur

# Generate midside nodes.
mids all gen

# All nodes are stored in a set.
seta nodes n all
# Boundary condition nodes are selected from the all nodes and stored in a set.
enq nodes symm_BC rec 0 _ 0 0.001 a


# Node set containing boundary condition nodes are written to a file.
send symm_BC abq nam

# Exported nodal information to a file.
send nodes abq

# Exported elemental information to a file.
seta axi e axisurf
send axi abq
seta ps e psurf
send ps abq

# Set view orientation.
rot -z



# Save image of the mesh model with BCs.
ulin Geometry.
view elem
plot l all
frame
hcpy png geom

# Save image of the mesh model with BCs.
ulin Finite element model. BCs (red).
view elem
plot f all g
plus n symm_BC r 6
frame
hcpy png pre



# sets plot
ulin Disk contact (blue), blade contact (red)
frame
zoom 5
tra d 16
plot f sbcon r 6
plus f sdcon b 6
plus f sdcon b 6

view elem
hcpy png conta

.inp

*INCLUDE, INPUT=nodes.msh
*INCLUDE, INPUT=axi.msh
*INCLUDE, INPUT=ps.msh
*INCLUDE, INPUT=symm_BC.nam 
*INCLUDE, INPUT=sbcon.sur
*INCLUDE, INPUT=sdcon.sur


*BOUNDARY
Nsymm_BC,1,1,0

*ELSET,ELSET=Eall
Eaxi,Eps


*MATERIAL,NAME=MatAxi
*ELASTIC
1.7E5,.3
*DENSITY
7.8E-9

*MATERIAL,NAME=MatPs
*ELASTIC
1.85E5,.3
*DENSITY
8.22E-9

*SOLID SECTION,ELSET=Eaxi,MATERIAL=MatAxi
*SOLID SECTION,ELSET=Eps,MATERIAL=MatPs
9500

*SURFACE INTERACTION, NAME=contact
*SURFACE BEHAVIOR, PRESSURE-OVERCLOSURE=LINEAR
10000000,0.1,0.01
*FRICTION
0.3,1000000

*CONTACT PAIR, INTERACTION=contact, TYPE=SURFACE TO SURFACE
Ssbcon,Ssdcon

*STEP, NLGEOM, INC=200
*STATIC
*DLOAD
Eall,CENTRIF,2020,0.,0.,0.,1.,0.,0.
*NODE PRINT,NSET=Nall
U
*EL PRINT,ELSET=Eall
S
*EL FILE
U,S
*CONTACT FILE
CDIS,cstr
*END STEP
2

I can’t see how the geometry looks because I’m not using CGX but I’m curious what you intend the physical meaning of contact between 2D and axisymmetric elements to be. It seems kind of hacky.

Maybe cyclic symmetry would be more appropriate than axial symmetry?

Since axisymmetric elements are expanded to a thin wedge, are those expanded element in the appropriate location to contact with the plane stress elements? Maybe they’re missing each other?

3

Hi @suleymanmuti, I imported what you shared for both files, and ccx errors out about missing node set S1 and NALL. Check if those were properly reproduced from the fbd file.

4

No fixing at least one node in the y-axis, why?

hcpy_1
hcpy_11920×1053 17 KB

5

Hi. Thank you for the replies.

y-axis was the symmetry axis and x-axis was the axis of rotation. I figured out that y-axis needs to be the axis of rotation. I swapped the axes and it works fine.



cdt
cdt779×615 9.68 KB

\

1 Like