CALCS PROC072
Pressure between two Bodies in Contact 180725

CONTENTS

GUIDANCE
NOTATION
INPUT DATA FORMAT
PROJECT EXAMPLE - INPUT FOR GENERAL CASE
PROJECT EXAMPLE - OUTPUT FOR GENERAL CASE
PROJECT EXAMPLE - CYLINDER ON FLAT SURFACE


GUIDANCE

APPLICATION
General case of contact between two bodies. This allows the specification of two mutually orthogonal principle axes of curvature 1/radius for each body.
That is, the specification of 2x2=4 radii, continuous in the area of contact.
There are then limits to the scope of contacting surfaces shapes. Even so, it should be possible to obtain reasonable results for a wide variety of contacting shapes.
Important Special Cases Include:
Wheel-Rail with wheel crown & rail head radii;
Contact between cylinder & flat surface;
Contact between convex & concave cylindrical surfaces (rotating bearing);
Contact between sphere & flat surface;
Contact between convex & concave spherical surfaces (spherical bearing);
Contact between two spheres.
See also: PROC033
INPUT TO PROCEDURE
Generally, all four radii are dissimilar: R11<>R12<>R21<>R22.
Curvatures 1/R are used, then input iR=0 for a flat surface (infinite radius)
Consider the Wheel-Rail with wheel crown & rail head radii.
The radii may then be specified as follows: R11 wheel radius, R12 wheel crown radius (R12<>0 or R12=0 for flat crown), R21 rail head radius (R21<>0 or R12=0 for flat head).
Example of psi angle: wheel cylindrical rim radius R11 & rail head radius R21 then psi=90
Contact between two spheres. The radii may then be specified as follows: R11=R12 & R21=R22.

Analysis for a cylinder is invoked by the following radii: R11=0 for flat OR R11>0 convex OR R11<0 concave contacting surface; R22=0; R21 cylinder and R22=0.

OUTPUT FROM PROCEDURE
Plot showing:
Elliptic shape of contact zone;
Lengths of ellipse principle axes;
Orientation angle of ellipse;
Rotated edges of continuous surfaces.
Calculations solving the contact problem with results:
Orientation angle of contact ellipse;
Semi-axis lengths of contact ellipse;
Maximum contact stress;
Strength assessment.
THEORETICAL BASIS
S.P.Timoshenko. Theory of Elasticity Ed.3. McGraw-Hill (1970). 'Pressure between two bodies in contact - More General Case'.
R.J.Roark & W.C.Young. Formulas for Stress & Strain Ed.5. McGraw-Hill (1975). Table33 Item2 Cylinder & Item4 'General case of two bodies in contact'.
BS EN 13001-3-3:2014 Cl. 5.3.

NOTATION

limiting width of surface (bodies 1 & 2)
A Bf(R11,R12,R21,R22) geometric functions
bwidth of contact line (cylinder analysis)
b1 b2
c dlengths of contact semi-axes
dh1 dh2depth of hardening from contact surface (bodies 1 & 2)
iRcurvature (ex. iR11=1/R11)
k1 k2f(E,v) material constants
m nf(A,B) geometric functions
OceOrientation angle of contact ellipse
Pcontact force
psiAngle (deg) between pr. axes (parallel to curvatures 1/R11 & 1/R21) of bodies 1 &
psi1 psi2Orientation angles of both bodies (pr. axes 1/R11 & 1/R21)
q0max contact pressure
R11Body 1 Principal Radius 1 Minumum (+ve for convex curvature)
R12Body 1 Principal Radius 2 Maximum
R21Body 2 Principal Radius 1 Minumum
R22Body 2 Principal Radius 2 Maximum
thetageometric constant

INPUT DATA FORMAT

LET E1=:nu1=:Fy1=:dh1=0:E2=:nu2=:Fy2=:dh2=0
LET P=:b1=:b2=
LET R11=:R12=0:R21=:R22=0
LET psi1=RAD(0):psi2=RAD(0)
PROC072

PROJECT INPUT


Twelve related test cases are presented.
Each group of four cases considers Body2 being orientated 0, 30, 60 & 90deg from Body1 which is not rotated.
The first group has equal min & max radii for the two bodies. That is: R11=R21 & R12=R22.
The second group (4 cases) has R11 less than R21 & R12=R22.
The third group (4 cases) has R11 greater than R21 & R12=R22.

PROCtitle("1000mm Wheel Test","","Pmax=500t")
LET E1=2070:nu1=.3:Fy1=8:dh1=0:E2=2070:nu2=.3:Fy2=8:dh2=0:REM 180727
LET P=500:b1=50:b2=60
LET R11=100:R12=20000:R21=100:R22=20000:REM 180725
LET psi1=RAD(0):psi2=RAD(0):PROC072
LET psi1=RAD(0):psi2=RAD(30):PROC072
LET psi1=RAD(0):psi2=RAD(60):PROC072
LET psi1=RAD(0):psi2=RAD(90):PROC072
LET R11=100:R12=20000:R21=1000:R22=20000
LET psi1=RAD(0):psi2=RAD(0):PROC072
LET psi1=RAD(0):psi2=RAD(30):PROC072
LET psi1=RAD(0):psi2=RAD(60):PROC072
LET psi1=RAD(0):psi2=RAD(90):PROC072
LET R11=1000:R12=20000:R21=100:R22=20000
LET psi1=RAD(0):psi2=RAD(0):PROC072
LET psi1=RAD(0):psi2=RAD(30):PROC072
LET psi1=RAD(0):psi2=RAD(60):PROC072
LET psi1=RAD(0):psi2=RAD(90):PROC072

PROJECT OUTPUT

The calculated results & plots allow interesting comparisions to be made between these cases.
For example, for the three radii groups and four orentation set, note the variation in aspect ratio, stresses & orientation of the contact ellipse.
These are intuitively correct and confirm the importance of radii & orientation specification on shape & area of contact zone and stress magnitude.


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 0
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.02
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 1E-4
Oce= atn_YiRdivXiR= 5E-3
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 74.63
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.99
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 24.89
d= beta*(P*KD*CE)^(1/3)= 0.9188
fc= 1.5*P/(c*d)= 32.8
RF_001= (3*Fy1)/(fc)= 0.7317*
RF_002= (3*Fy2)/(fc)= 0.7317*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 0.5236
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.01864
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 5.093E-3
Oce= atn_YiRdivXiR= 0.2668
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 74.63
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.8574
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 8.557
d= beta*(P*KD*CE)^(1/3)= 1.601
fc= 1.5*P/(c*d)= 54.74
RF_003= (3*Fy1)/(fc)= 0.4385*
RF_004= (3*Fy2)/(fc)= 0.4385*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 1.047
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.01496
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 8.735E-3
Oce= atn_YiRdivXiR= 0.5286
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 74.63
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.495
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 4.736
d= beta*(P*KD*CE)^(1/3)= 2.304
fc= 1.5*P/(c*d)= 68.75
RF_005= (3*Fy1)/(fc)= 0.3491*
RF_006= (3*Fy2)/(fc)= 0.3491*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 1.571
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 9.95E-3
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 0.01005
Oce= atn_YiRdivXiR= 0.7904
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 74.63
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 3.201
d= beta*(P*KD*CE)^(1/3)= 3.201
fc= 1.5*P/(c*d)= 73.18
RF_007= (3*Fy1)/(fc)= 0.3279*
RF_008= (3*Fy2)/(fc)= 0.3279*
END072
LET R11=100:R12=20000:R21=1000:R22=20000


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 0
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.011
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 1E-4
Oce= atn_YiRdivXiR= 9.091E-3
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.982
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 24.59
d= beta*(P*KD*CE)^(1/3)= 1.248
fc= 1.5*P/(c*d)= 24.44
RF_009= (3*Fy1)/(fc)= 0.982*
RF_010= (3*Fy2)/(fc)= 0.982*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 0.5236
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.01084
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 5.933E-4
Oce= atn_YiRdivXiR= 0.05467
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.9421
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 15.2
d= beta*(P*KD*CE)^(1/3)= 1.594
fc= 1.5*P/(c*d)= 30.96
RF_011= (3*Fy1)/(fc)= 0.7753*
RF_012= (3*Fy2)/(fc)= 0.7753*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 1.047
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.01046
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 9.41E-4
Oce= atn_YiRdivXiR= 0.08975
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.8568
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 10.41
d= beta*(P*KD*CE)^(1/3)= 1.954
fc= 1.5*P/(c*d)= 36.88
RF_013= (3*Fy1)/(fc)= 0.6507*
RF_014= (3*Fy2)/(fc)= 0.6507*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 1.571
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 9.95E-3
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 1.05E-3
Oce= atn_YiRdivXiR= 0.1051
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.8108
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 9.203
d= beta*(P*KD*CE)^(1/3)= 2.091
fc= 1.5*P/(c*d)= 38.97
RF_015= (3*Fy1)/(fc)= 0.6159*
RF_016= (3*Fy2)/(fc)= 0.6159*
END072
LET R11=1000:R12=20000:R21=100:R22=20000


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 0
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 0.011
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 1E-4
Oce= atn_YiRdivXiR= 9.091E-3
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.982
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 24.59
d= beta*(P*KD*CE)^(1/3)= 1.248
fc= 1.5*P/(c*d)= 24.44
RF_017= (3*Fy1)/(fc)= 0.982*
RF_018= (3*Fy2)/(fc)= 0.982*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 0.5236
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 9.635E-3
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 5.093E-3
Oce= atn_YiRdivXiR= 0.4863
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.9421
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 15.2
d= beta*(P*KD*CE)^(1/3)= 1.594
fc= 1.5*P/(c*d)= 30.96
RF_019= (3*Fy1)/(fc)= 0.7753*
RF_020= (3*Fy2)/(fc)= 0.7753*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 1.047
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 5.957E-3
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 8.735E-3
Oce= atn_YiRdivXiR= 0.9723
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.8568
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 10.41
d= beta*(P*KD*CE)^(1/3)= 1.954
fc= 1.5*P/(c*d)= 36.88
RF_021= (3*Fy1)/(fc)= 0.6507*
RF_022= (3*Fy2)/(fc)= 0.6507*
END072


1000mm Wheel - PROC072 Contact General Case - Pmax=500t

LET iR11=1/R11:iR12=1/R12:iR21=1/R21:iR22=1/R22
psi= ABS(psi2-psi1)= 1.571
XiR= iR11*COS(psi1)-iR12*SIN(psi1)+iR21*COS(psi2)-iR22*SIN(psi2)= 9.5E-4
YiR= iR11*SIN(psi1)+iR12*COS(psi1)+iR21*SIN(psi2)+iR22*COS(psi2)= 0.01005
Oce= atn_YiRdivXiR= 1.477
CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
KD= 1.5/(iR11+iR12+iR21+iR22)= 135.1
Ctheta= (KD/1.5)*SQR((iR11-iR12)^2+(iR21-iR22)^2+2*(iR11-iR12)*(iR21-iR22)*COS(2*psi))= 0.8108
LET theta=ACS(Ctheta):thetaD=DEG(theta)
alpha=FN_LIP(Ctheta,alpha())
beta=FN_LIP(Ctheta,beta())
c= alpha*(P*KD*CE)^(1/3)= 9.203
d= beta*(P*KD*CE)^(1/3)= 2.091
fc= 1.5*P/(c*d)= 38.97
RF_023= (3*Fy1)/(fc)= 0.6159*
RF_024= (3*Fy2)/(fc)= 0.6159*
END072

120mm Cylinder - PROC072 Contact General Case - P=20t

LET E1=2070:nu1=.3:Fy1=9:dh1=.5:E2=2070:nu2=.3:Fy2=9:dh2=0
LET P=20:b1=10:b2=8
LET R11=0:R12=0:R21=6:R22=0
LET psi1=RAD(0):psi2=RAD(0):P=20


120mm Cylinder - PROC072 Contact by Cylinder 180725 - P=20t

CE= (1-nu1^2)/E1+(1-nu2^2)/E2= 8.792E-4
LET D1=2*R11:D2=2*R21
KD= D2= 12
p= P/b2= 2.5
b= 1.6*SQR(p*KD*CE)= 0.2599
fc= .798*SQR(p/(KD*CE))= 12.28
RF_025= (dh1)/(.4*b)= 4.81
RF_026= (3*Fy1)/(fc)= 2.198
RF_027= (3*Fy2)/(fc)= 2.198
END072

PROJECT END

END072