# Force: Contact

Model ElementForce_Contact defines a 2-D or 3-D contact force between geometries on two rigid bodies.

## Description

Each body is characterized by a set of one or more geometries which can be a 3D mesh, an analytically defined solid, or a 2D curve. Whenever any geometry on the first body penetrates any geometry on the second body, contact normal and frictional forces are generated. The normal force tends to repulse motion along the common normal at the contact point. The frictional force tends to oppose relative sliding velocity at the contact point. The contact force vanishes when the geometries no longer intersect.

## Format

<Force_Contact
id                             = "integer"
[  label                          = "integer"   ]
[  full_label                     = "string"    ]
num_i_graphics                 = "integer"
i_graphics_id                  = "integer_list"
num_j_graphics                 = "integer"
j_graphics_id                  = "integer_list"
ignore_penetration_larger_than = "real"
master_surface                 = { "I" | "J" | "I_and_J" | "AUTO" }
enable_node_mode               = {"TRUE" | "FALSE" }
enable_analytical              = {"TRUE" | "FALSE" }
contact_stability_2d           = "real"
{
cnf_type                       = "POISSON"
penalty                        = "real"
restitution_coef               = "real"
normal_trans_vel               = "real"
|
cnf_type                       = "IMPACT"
stiffness                      = "real"
exponent                       = "real"
damping                        = "real"
dmax                           = "real"
|
cnf_type                       = "VOLUME"
i_elastic_modulus              = "real"
j_elastic_modulus              = "real"
i_layer_depth                  = "real"
j_layer_depth                  = "real"
exponent                       = "real"
damping                        = "real"
|
cnf_type                       = "USERCNF"
cnf_param_string               = "USER( [[par_1 [,...][,par_n]] )"
cnf_fnc_name                   = "string"
{ cnf_dll_name                   = "string"
|
script_name                    = "string"
interpreter                    = { "PYTHON" | "MATLAB" }
}
}
{
cff_type                       = "COULOMB_ON"
mu_static                      = "real"
mu_dynamic                     = "real"
stiction_trans_vel             = "real"
friction_trans_vel             = "real"
|
cff_type                       = "COULOMB_DYNAMICONLY"
mu_dynamic                     = "real"
friction_trans_vel             = "real"
|
cff_type                       = "COULOMB_OFF"
|
cff_type                       = "USERCFF"
cff_param_string               = "USER( [[par_1 [,...][,par_n]] )"
cff_fnc_name                   = "string"
{ cff_dll_name                   = "string"
|
script_name                    = "string"
interpreter                    = { "PYTHON" | "MATLAB" }
}
}
{
cpost_param_string             = "USER( [[par_1 [,...][,par_n]] )"
cpost_fnc_name                 = "string"
{ cpost_dll_name                 = "string"
|
cpost_script_name              = "string"
cpost_interpreter              = "PYTHON"
}
}

## Attributes

id
Element identification number (integer>0).

This number is unique among all Force_Contact elements and uniquely identifies the element.

Note: id is used only to identify the Force_Contact element to be modified. By itself, id cannot be modified.
label
The name of the Force_Contact element. This argument is optional.
full_label
The full label of the Force_Contact element. This attribute is typically populated by MotionView as per the model hierarchy. This argument is optional.
num_i_graphics
Specifies the number of Post_Graphic elements on the first body that are to be considered in evaluating the contact force.

num_i_graphics > 0.

i_graphics_id
This is the list of the Post_Graphic element IDs on the first body to be considered for contact. The number of IDs in this list is specified by num_i_graphics.
Note: All the Post_Graphic entities must belong to the same body.
num_j_graphics
Specifies the number of Post_Graphic elements on the second body that are to be considered in evaluating the contact force.
num_j_graphics > 0
j_graphics_id
This is the list of the Post_Graphic element IDs on the second body to be considered for contact. The number of IDs in this list is specified by num_j_graphics.
Note: All the Post_Graphic entities must belong to the same body.
ignore_penetration_larger_than
Use this attribute to change the maximum penetration depth that MotionSolve considers for evaluating contact forces. Any value above this threshold is ignored.
master_surface
Specifies which surface should be used as the master surface while calculating the penetration depth, point of contact, contact normal, and so on. Options include:

I: The surface of body I is used as the master surface.

J: The surface of body J is used as the master surface.

I_AND_J: Surfaces of both bodies, I and J, are used as master surface alternatively and the result is summed up.

AUTO: Let MotionSolve decide which geometry to use as the master surface. This is the default and recommended option.

enable_node_mode
For mesh based geometries, specifies whether the normal forces are calculated at the centers of the interconnected triangular elements (FALSE), or at their connection nodes (TRUE).
The default is FALSE.
enable_analytical
Specifies whether MotionSolve should use an analytical description, where possible. Select TRUE or FALSE.
The default is FALSE.
See Comment 1 for more details on usage.
contact_stability_2d
Specifies a stability value that is used when calculating the contact point and penetration depth in 2D curve contact.
The default value for this parameter is 0.1. 9
cnf_type
Specifies whether the POISSON, IMPACT, VOLUME or user defined (USERCNF) force model is to be used for calculating the contact normal force.
A contact force is computed as a function of the penetration between the intersecting geometries defined for the contact. In other words, as the geometries intersect, a repulsive force is generated which is based on the type of contact selected. 2
penalty (POISSON)
Specifies the stiffness parameter that is used to calculate the spring force. A large value for penalty permits only a small penetration between the two contacting geometries; a small value permits a larger penetration. Consider using a penalty value that is as small as possible, but still captures realistic deformation, to improve solver performance.
penalty ≥ 0
restitution_coef (POISSON)
Defines the coefficient of restitution (COR) between the geometries in contact. This is used in the computation of the contact force. A value of zero specifies perfectly plastic contact meaning that the two bodies coalesce after contact. A value of one specifies perfectly elastic contact. No energy is lost in the collision and the relative velocity of separation equals the relative velocity of approach.
COR = (relative speed after collision)/(relative speed before collision)
1 ≥ restitution_coef ≥ 0
normal_trans_vel (POISSON)
Defines the velocity limit between the two bodies at which full damping is applied in the contact force. If unset by the user, MotionSolve initializes it to 1. See Comment 4 for more details on this parameter’s use in the damping force.
normal_trans_vel > 0
stiffness (IMPACT)
The stiffness parameter of the contact; this is the same parameter as found for the IMPACT function. A large value for stiffness permits only a small penetration between the two contacting geometries; a small value permits a larger penetration. This parameter must be non-negative. Consider using a stiffness value that is as small as possible, but still captures accurate deformation, to improve solver performance.
stiffness ≥ 0
exponent (IMPACT)
The exponent of the force deformation characteristic of the contact interface. For a stiffening spring characteristic, it must be greater than 1.0 and for a softening spring characteristic, it must be less than 1.0.
exponent > 0
damping (IMPACT)
The maximum damping coefficient. It must be non-negative.
damping ≥ 0
dmax (IMPACT)
The penetration at which full damping is applied. It must be a positive number.
dmax > 0
i_elastic_modulus (VOLUME)
Specifies the elastic modulus for the geometries belonging to Body I for the Volume force model. The value of this elastic modulus can be derived from the Bulk and Shear Modulus of a material. 5
i_elastic_modulus > 0
j_elastic_modulus (VOLUME)
Specifies the elastic modulus for the geometries belonging to Body J for the Volume force model. The value of this elastic modulus can be derived from the Bulk and Shear Modulus of a material. 5
j_elastic_modulus > 0
i_layer_depth (VOLUME)
Specifies the size of the layer depth for the geometries belonging to body I. 5
i_layer_depth > 0
j_layer_depth (VOLUME)
Specifies the size of the layer depth for the geometries belonging to body J. 5
j_layer_depth > 0
exponent (VOLUME)
The exponent of the force deformation characteristic of the contact interface. For a stiffening spring characteristic, it must be greater than 1.0 and for a softening spring characteristic, it must be less than 1.0.
exponent > 0
damping (VOLUME)
The coefficient of damping used to calculate the damping force for the VOLUME force model
damping ≥ 0
cnf_param_string (USERCNF)
Specify a list of parameters that are passed from the data file to the user written subroutine, CNFSUB.
cnf_dll_name (USERCNF)
Specifies the path and name of the shared library containing the contact force calculation subroutine. MotionSolve uses this information to load the user subroutine in the library at run time.
cnf_fnc_name (USERCNF)
Specifies an alternate name for the contact normal force subroutine, that is other than CNFSUB.
script_name (USERCNF, USERCFF)
Specifies the path and name of the user written script that contains the routine specified by cnf_fnc_name or cff_fnc_name.
interpreter (USERCNF, USERCFF)
Specifies the interpreted language that the user script is written in. Valid choices are MATLAB or PYTHON.
cff_type

Specifies the friction force model that will be used to compute the contact friction force. Choose from "COULOMB_ON", "COULOMB_OFF", "DYNAMIC_ONLY" or "USERCFF".

cff_type = "COULOMB_ON" specifies that friction forces are to be calculated and applied at the contact location. MotionSolve uses the Coulomb model for friction force. See Comment 7 below.

cff_type = "COULOMB_OFF" specifies that frictional forces are not calculated or applied at the contact location. Friction is turned off.

cff_type = "COULOMB_DYNAMICONLY" specifies that only the dynamic friction is active. See Comment 7.

cff_type = "USERCFF" specifies that a user defined friction force model is to be used.

mu_static

Defines the coefficient of static friction when the friction is in the static regime.

MotionSolve uses a step function to transition between the static and dynamic friction regimes. 7

mu_static > 0
mu_dynamic
Defines the coefficient of dynamic friction when the friction is in the dynamic regime.

MotionSolve uses a step function to transition between the static and dynamic friction regimes. The coefficient of dynamic friction is smaller than or equal to the coefficient of static friction. 7

mu_staticmu_dynamic > 0
stiction_trans_vel
The velocity at which the friction regime transitions to stiction. 7
friction_trans_vel
The velocity at which the friction regime transitions to friction. 7
|friction_trans_vel| ≥ |stiction_trans_vel|
cff_dll_name (USERCFF)
Specifies the path and name of the shared library containing the contact friction sub subroutine. MotionSolve uses this information to load the user subroutine in the library at run time.
cff_fnc_name (USERCFF)
Specifies an alternate name for the contact normal force subroutine. In other words, other than CFFSUB.
cff_param_string (USERCFF)
Specify a list of parameters that are passed from the data file to the user written subroutine, CFFSUB.
cpost_dll_name
Specifies the path and name of the shared library containing the contact post sub subroutine. MotionSolve uses this information to load the user subroutine in the library at run time.
cpost_fnc_name
Specifies an alternate name for the contact post subroutine. In other words, other than CONTACTPOST.
cpost_param_string
Specify a list of parameters that are passed from the data file to the user written subroutine, CONTACTPOST.
cpost_interpreter
Specifies the interpreted language that the user script is written in. Currently, only Python is supported as an interpreter for the CONTACTPOST subroutine.
cpost_script_name
Specifies the path and name of the user written script that contains the routine specified by cpost_fnc_name.

## Modeling and Debugging Tips

• Confirm that models run without contact before adding them.
• Start modeling with contact friction off, and add it later after getting the normal force running.
• Some experimentation is usually required for tuning the contact parameters to get realistic and stable results.
• Small values for stiction_trans_vel and friction_trans_vel cause the integrator to take too small of steps when the coefficient of friction is large. Increase the value for stiction_trans_vel if you are experiencing integration difficulties. In our experience, this does not change the results significantly.
• Large values for the penalty or stiffness decreases the penetration allowed between geometries; however they also tend to decrease the integrator step size. Consider allowing as much penetration of the geometries as will still capture the behavior of the system that you're interested in measuring. This will be model-dependent. The contact force computation will be much smoother in general if a moderate amount of penetration (due to a lower stiffness/penalty) is allowed. Ask yourself, "Will 0.1 mm of penetration change the behavior of my system? Will 1 mm?"
• In general, a finer mesh at the contact surface will help obtain more accurate results. However, a fine mesh also increases the computational effort of the collision-detection algorithm, which may slow down the simulation.
• Check the normals of the geometries to make sure that they have been defined properly for the contact. In the Tria mesh geometry, the normal is defined by a right-hand rule of the mesh nodes. For example, you can view the normals using the check normals tool in the Contact panel in MotionView.
• If the contacts are not detected or pass through each other during the simulation, try one of the following: decreasing the integrator step size (for example, h_max), decreasing the error tolerance, increasing the penalty/stiffness (if too soft), or checking the element normals.

Alternatively, you can turn on the "zero crossing" action for the contact sensor in the Contact panel in MotionView. By doing this, MotionSolve will attempt to determine the moment of first contact more accurately which may alleviate this problem. For more information on this sensor, please refer to the documentation for Sensor_Event.

• If the integrator fails, use DebugOutput to determine if the contact is causing problems. If it is, then try one of the following: decreasing the contact penalty/stiffness and/or damping (if too large), decreasing h_max, increasing the integrator error tolerance, adding damping to other stiff entities in the model, or tuning the contact friction parameter as described previously.

## Limitations

The current implementation of Force_Contact has some limitations.
• It is only supported for rigid body and point mass objects. Flexible bodies are not supported.
• Stiction or stick-slip phenomena are not supported at this time.
• For analytical Cylinder to Box contact, a contact between the Edge or Vertex of a Box to the flat side of the Cylinder is not supported at this time.

## Example-1: Mesh-to-Mesh Contact

The image below demonstrates an epicyclic gear train that is used in a child's toy. When the child cranks the sun gear, the other gears rotate.

The gear train consists of two contacts: one between the Sun gear and the Planet gear, and a second between the Planet gear and the Ring gear. The gears are made out of copper, and some oil is used to lubricate the system.

The Force_Contact statements for this example are shown below. The system is modeled in SI units.

The first contact is between the Sun gear and the Planet gear.
<Force_Contact
id                    = "1"
num_i_graphics        = "1"
i_graphics_id         = "82"
num_j_graphics        = "1"
j_graphics_id         = "92"
cnf_type              = "POISSON
penalty               = "1E6"
restitution_coef      = "0.52"
cff_type              = "COULOMB_ON"
mu_static             = "0.08"
mu_dynamic            = "0.05"
stiction_trans_vel    = "0.005"
friction_trans_vel    = "0.01"
/>
The second contact, shown below, is between the Planet and the Ring gears.
<Force_Contact
id                    = "2"
num_i_graphics        = "1"
i_graphics_id         = "93"
num_j_graphics        = "1"
j_graphics_id         = "73"
cnf_type              = "POISSON"
penalty               = "1E6"
restitution_coef      = "0.52"
cff_type              = "COULOMB_ON"
mu_static             = "0.08"
mu_dynamic            = "0.06"
stiction_trans_vel    = "0.005"
friction_trans_vel    = "0.0l"
/>

For a more detailed example, see the tutorial MV-1010: Contact Simulation using MotionSolve.

## Example-2: Analytical 3D Contact

This example illustrates the use of analytical 3D contact. There is no change to the user actions in MotionView. Behind the scenes, the solver recognizes the geometry pair in contact as being “simple” and it switches to using a specialized analytical contact formulation.

The model consists of a cylinder that is placed on a flat box. The cylinder rotates about its axis, which is parallel to the box top surface. A dynamic simulation is performed to evaluate the contact behavior for this simple case. The rotation of the cylinders simply causes triangles in the meshes to move in and out of contact. There is no change in the contact kinematics. Thus, one would expect the contact algorithm to report a constant force.

The problem is solved using 3D analytical contact as well as mesh-to-mesh contact. The contact forces are shown in the plot below.

The analytical contact reports a constant force as expected. The mesh-to-mesh contact representation, in contrast, reports a lower force. Moreover, as triangles move in and out of contact, in mesh-to-mesh contact, the force changes to reflect this. The variation in the force can be decreased by using finer tessellations, however this would increase the cost of the simulation.

The image below compares the performance of the analytical contact against the mesh-to-mesh contact. The analytical contact (left) is about 21x faster than the mesh-to-mesh contact (right) for this simple model.

1. MotionSolve use a geometric representation of the contacting bodies to compute the contact force. Geometric representations can be 3D solid meshes, 3D analytical shapes, and 2D curves.

Many geometric combinations are possible for 3D contact. When enable_analytical is set to “TRUE”, MotionSolve uses analytical representations to compute accurate, fast, and smooth contact forces. When this is not possible, MotionSolve uses a general 3D mesh-to-mesh-contact.

Support for the different combinations is shown in the table below.

Table 1: Contact kinematics calculation methods for supported geometries
In the above table, the entries in each cell have the following meaning:
Analytical
Each of the contacting geometries is represented by an analytical expression. Contact between the geometries is determined analytically. Contact point locations and normal are computed without using any tessellation representations. Curved geometries are represented exactly.
Semi
One of the geometries is represented analytically, i.e., it is not tessellated. However, the other contacting geometry is tessellated and represented as a triangular mesh. The accuracy is dependent on the accuracy of the tessellation.
Mesh
Both geometries are tessellated. A general mesh-to-mesh contact detection method is used. Contact accuracy is dependent on the accuracy of the mesh. The contact detection speed is inversely dependent on the mesh density. In other words, it will take more time to compute the contact kinematics between fine meshes.
Using a 3D mesh representation: This geometrical representation of the body is termed a graphic for that body. MotionSolve relies on a collision engine to detect intersection or penetration between two graphics. To do this, the collision engine requires the graphics to be meshed using triangular elements. A meshed graphic is a set of interconnected triangular elements that describes a closed solid. This is illustrated in the figures below:
MotionSolve prescribes certain conditions for the meshed representation of a body's graphic(s). These conditions must be met to generate realistic contact forces. These are:
• Only triangular shell elements are allowed in the meshed representation. Quadrilateral or other shell elements are not permitted.
• The mesh must be a surface mesh i.e. a mesh comprised of 2D triangular shell elements with no thickness. Solid meshes that make use of tetrahedral, hexahedral or other element types are not required.
• The mesh must enclose a volume; it must be closed. In other words, there must not be any free edges or cracks in the mesh. For example, a flat surface with zero volume is not allowed.
• T-connections are not allowed in the mesh - this means that an edge of a triangular element must be shared by no more than two elements.
• The surface normal for a meshed representation of a body’s graphic(s) must point outward from the surface that can come into contact. The normal direction must be defined in the corresponding graphic component (see Post: Graphic) .

For more information on these rules, please refer to Best Practices for Running 3D Contact Models in MotionSolve.

Graphics for rigid bodies may be defined in many ways:
• For spherical geometries, MotionSolve can use an analytical description of the geometry (based on a radius and a location for the center of the sphere) while computing the contact point, contact normal, and other quantities. This approach is followed whenever MotionSolve detects a primitive sphere geometry (of type = "SPHERE" in Post_Graphic) for Body i, Body j, or both. There are two advantages of using an analytical description for a spherical body:
• Since there is no surface mesh, there are no errors due to the approximation of the surface by a triangulated mesh. This means the contact forces are more accurate;
• The solution time is typically much smaller for contact involving a sphere and a generalized 3D mesh.

To force MotionSolve to use a meshed representation of the sphere geometry, you can specify enable_sphere_to_mesh to be FALSE. This way, MotionSolve will mesh the primitive sphere geometry and use the meshed representation of the sphere to compute the contact related quantities.

• For simple shapes, you may use built-in primitive graphics provided by MotionSolve such as Sphere, Frustum, Box, and so on. For these primitive shapes, MotionSolve meshes the surface internally (using the attribute refinement_level in Post: Graphic) to determine the density of the mesh.

More complicated shapes can be defined in a CAD package and imported to MotionSolve using the Import CAD Tool in MotionView. Please see Import CAD or FE for more information on how CAD geometry can be meshed using this tool. Alternatively, you may also use HyperMesh to mesh your CAD graphics.

Using a 3D analytical representation: When the two contacting geometries are “simple”, the general mesh-to-mesh contact method is not used. Simple geometries include sphere, box, cylinder and frustrum. In this case, simpler analytical computations are employed to detect the contact point locations and the directions of the normal. Tessellation of these geometries is not required, and it is not used.

Using a 2D representation: MotionSolve also supports contact between two curves that represent two rigid bodies. 2D curve graphics can be specified by referencing a Post_Graphic of type = "CURVE". Using 2D curve-based contact can be advantageous over 3D mesh-based contact when:
1. It is known a prior that the contact occurs only along the curves that are defined. Further, contact is expected to occur within the plane in which the two curves are defined. In other words, there are no out-of-plane contact forces that are expected.
2. The curves between which the contact is to be calculated are smooth and represent the curvature of the 3D geometry well.
There may be several advantages of simulating rigid body contact using 2D curves over 3D tessellated geometry:
1. More accurate results: Results are more accurate since there is less error due to discretization. Discretization of 3D geometry into triangles occurs over three dimensions as opposed to two in curves.
2. Easily improve accuracy: Result accuracy can be improved rapidly by introducing more points in the curve within MotionView. Improving the accuracy of results when using 3D contact requires you to re-import the graphic and re-mesh, or go outside the MotionView environment into a CAD tool to re-mesh. Re-meshing a 3D geometry takes more time than refining a 2D curve.
3. Faster simulation times: For geometries where contact occurs over a curve, using a curve-curve contact is much faster than using a 3D tessellated geometry since the solver has to work harder to determine contact for a 3D tessellated geometry which increases the overall simulation time.
As with 3D mesh-based contact, certain restrictions apply:
1. Curves have to be flat and defined in the x-y plane relative to the reference marker specified in Post: Graphic.
2. CURSUB routines are currently not supported.
3. All i_graphic_id's and all j_graphic_id's have to point to curve graphics, which are defined in the same plane.
4. The z-direction of each graphic reference marker has to point in the same direction.
5. You must ensure that the graphics will remain in one plane during simulation. The easiest way to achieve this is by adding an inplane joint to the model.
6. There can be no contact between two bodies that are represented only by straight lines. This is because the curve contact algorithm requires there to be at least two points of intersection.
Note: 2D curve contact uses a smooth representation of the specified curves which is independent of the nseg specified in Post: Graphic.
2. Three pre-defined contact force models are available for computing the normal force: POISSON, IMPACT and VOLUME. If none of these is suitable, you may write your own force model as a subroutine or script and use that to compute the normal contact force.
The normal force computed in the POISSON, IMPACT and VOLUME methods consists of two components: an elastic spring force and a dissipative damping force. Each method has its own way of computing these forces, described in the subsequent comments.
${F}_{contact}={F}_{spring}+{F}_{damping}$
3. In the IMPACT model, the contact normal force is modeled with the IMPACT function. The penetration depth and rate are computed by the collision detection algorithm and fed into the IMPACT function.
${F}_{spring}=K{z}^{\mathrm{exp}}$
${F}_{damping}=-STEP\left(z,-{d}_{\mathrm{max}},c,0,0\right)\frac{dz}{dt}$
$K$
Stiffness
$z$
Contact penetration depth.
$exp$
Exponent for the force deformation characteristic.
${d}_{\mathrm{max}}$
Penetration depth at which full damping is applied.
$c$
Damping coefficient.
4. In the POISSON model, the elastic spring force component of the normal force is modeled using the Hertz contact theory, while the damping force uses the coefficient of restitution and a user defined velocity limit as part of the calculation. The equations that define this force are listed below.
${F}_{spring}=K\sqrt{{z}^{3}}$
${F}_{damping}=-{F}_{spring}\left(\left(\frac{1-{C}_{R}{}^{2}}{1+{C}_{R}{}^{2}}\right)\ast STEP\left(\frac{dz}{dt},-{v}_{norm,trans},-1.0,{v}_{norm_trans},1.0\right)\right)$
$K$
Penalty parameter (stiffness).
$z$
Contact penetration depth.
${C}_{R}$
The coefficient of restitution. This is equal to 0 for a perfectly plastic collision and 1 for a perfectly elastic collision.
${v}_{norm,trans}$
Normal velocity at which full damping is applied.
Note: For a coefficient of restitution = 0, when the two bodies are colliding head to head, and dz/dt > vnorm_trans, then the damping force is equal to the spring force, and the total contact normal force is zero.

In this model, the energy loss is modeled with a coefficient of restitution. Table 2 lists the coefficient of restitution for some materials. Contact is assumed to be between like materials, for example, Brass-Brass, Steel-Steel, and so on.

Table 2: Coefficient of Restitution for Common Materials
Material Coefficient of Restitution
Brass 0.30
Bronze 0.52
Copper 0.22
Cork 0.60
Elm 0.48
Glass 0.96
Iron 0.67
Ivory 0.78
Maple 0.65
Rubber 0.75
Steel 0.90
5. In the VOLUME model, the elastic spring force component of the normal force is modeled similar to the IMPACT model, while the damping force uses a damping coefficient to compute the damping force. The equations that define this force are listed below.
${F}_{spring}={A}_{con}K{z}^{\mathrm{exp}}$
${F}_{damping}=-c\frac{dz}{dt}$
where
$K$
Contact stiffness.
$z$
Contact penetration depth.
$exp$
Exponent for the force deformation characteristic.
$c$
Damping coefficient.
${A}_{con}$
Area of contact.

The VOLUME model assumes that both the colliding bodies are surrounded by a layer of springs whose stiffness is determined by the material's elastic modulus properties and the depth of this layer. The individual stiffness for each body (I and J) is calculated as:

$\begin{array}{l}{c}_{i}=\frac{{M}_{i}}{{d}_{i}}\\ {c}_{j}=\frac{{M}_{j}}{{d}_{j}}\end{array}$
where,
${c}_{i},{c}_{j}$
Stiffness parameter for I and J bodies, respectively.
${M}_{i},{M}_{j}$
Elastic Modulus for I and J bodies, respectively.
${d}_{i},{d}_{j}$
Layer depth for I and J bodies, respectively.

The equivalent contact stiffness is calculated assuming that the springs for body I and J are in series:

$K=\frac{{c}_{i}{c}_{j}}{{c}_{i}+{c}_{j}}$
Note: The equivalent stiffness calculated above has units of Force per cubic length. The VOLUME force method utilizes the contact area and penetration depth to retrieve the contact force in Newton (or the model unit for force).

The elastic modulus M (sometimes also called P-wave modulus) for both the bodies in contact is calculated using the Bulk and Shear Modulus of their respective materials. These can be specified in the Contact panel in MotionView. The elastic modulus is calculated as

$M=K+\frac{4G}{3}$
where,
$M$
Elastic or P-wave modulus.
$K$
Bulk Modulus.
$G$
Shear Modulus.

The Bulk and Shear Modulus for homogenous, isotropic materials can be calculated from the Young's Modulus and Poisson ratio as

$\begin{array}{l}G=\frac{E}{2\left(1+\upsilon \right)}\\ K=\frac{E}{3\left(1-2\upsilon \right)}\end{array}$
where,
$E$
Young's modulus
$\upsilon$
Poisson ratio

Both the elastic modulus and layer depth determine the value of the equivalent contact stiffness. Since for a given material, the elastic modulus is typically a constant value, you may change the layer depth to obtain an appropriate penetration in your model. For example, assume that two steel bodies are in contact. For steel (G = 160GPa, K = 79GPa), M = 292.33GPa. If you specify the layer depth to be 100mm, the contact stiffness is ~1.461E+03 N/mm3. If you reduce this layer depth by a factor of 10, i.e. 10mm, the contact stiffness goes up by a factor of 10 i.e. ~1.461E+04 N/mm3.

The VOLUME force model thus lets you derive the contact stiffness from physical material properties of the colliding bodies.

6. A good first estimate for the stiffness parameters (for both POISSON and IMPACT models) can be obtained by making use of the classical Hertzian contact theory for non-adhesive contact. The general case for contact between two spherical bodies of radius R1 and R2 is illustrated in the image below.

Hertz theory is developed using static load cases. The stiffness is a function of the material properties of the two bodies in contact as well as the geometry at the contact patch. The stiffness coefficient is shown in red.

Assume that two steel balls are in contact. The radii of the steel balls are 30mm; the Young's modulus for steel is 200 GPa; Poisson's ratio for steel is 0.3. Using the formula shown above, the stiffness coefficient for contact is: K = 1.7945E+10 N/m = 1.7945E+07 N/mm

This should be considered a starting point for the penalty or stiffness value. You may need to tune this parameter, as high values of stiffness can degrade the performance of the solver. However, you may find that you are able to allow a smaller penalty, and thus more penetration, and still capture the overall system behavior.

7. Two friction models are supported in MotionSolve
• Coulomb model of friction which calculates the friction force as the normal force multiplied by the coefficient of friction. This coefficient is calculated as a function of the slip velocity that is explained below
• A user defined friction force model that can be implemented as a subroutine or a script

For the first approach, the friction force at the contact patch is computed in the same manner for POISSON, IMPACT and VOLUME force models. The Coulomb friction law is used which models the friction force as a coefficient of friction multiplied by the normal force. The coefficient of friction is modeled as a function of the slip velocity. The frictional force is modeled as a viscous force according to the following equations:

(1)

where,
${V}_{t}$
is the current slip speed at the point of contact
${V}_{s}$

is the stiction transition slip speed at which the full value of ${\mu }_{s}$ is used for the coefficient of friction

${V}_{d}$
is the dynamic friction slip speed at which the full value of ${\mu }_{d}$ is used for the coefficient of friction
${\mu }_{s}$
is the coefficient of static friction
${\mu }_{d}$
is the coefficient of dynamic friction

Figure 8 illustrates the functional relation between coefficient of friction and slip velocity

The blue curve depicts the full friction model covering all three regimes: static, transition to sliding friction, and sliding friction. This is the friction model applied when you select cff_type = "COULOMB_ON".
• The static regime shown in light blue.
• The transition from static to dynamic friction show in white.
• The dynamic friction regime shown in light orange.

The red curve depicts the friction model covering only the sliding friction regime. This is the friction model applied when you select cff_type = "DYNAMICS_ONLY".

8. μ is the friction coefficient that is to be applied. As shown in table 3, it is a function of the slip velocity and dependent on the static and dynamic friction coefficients. The friction force opposes the direction of the slip velocity. The table below lists the static and dynamic friction coefficients for some common material pairs. The friction values are extremely sensitive to surface coatings, so both dry and greasy values are specified.