1. Nonlinear Behavior

1.1. Introduction

We begin our study of nonlinear computational solid mechanics in this chapter by surveying some frequently encountered sources of nonlinearity in engineering mechanics. This will be done in a rather elementary way, by discussing perhaps the simplest structural idealization, the truss member, which is assumed to transmit loads in the axial direction only. By introducing various nonlinearities into this system one at a time, we will motivate the more general discussion of nonlinear continuum mechanics, constitutive modeling, and numerical treatments to follow. This model system will serve as a template throughout the text as new continuum mechanical and computational ideas are introduced.

Following this motivation will be an introduction to the prescription of initial/boundary value problems in solid mechanics. This introduction will be provided by discussing a completely linear system, namely linear elastic behavior in a continuum subject to infinitesimal displacements. This treatment will include presentation of the relevant field equations, boundary conditions, and initial conditions, encompassing both dynamic and quasistatic problems in the discussion. Also featured is a brief discussion of the weak or integral form of the governing equations, providing a starting point for application of the finite element method. Examination of these aspects of problem formulation in the comparatively simple setting of linear elasticity allows one to concentrate on the ideas and concepts involved in problem description without the need for an overly burdensome notational structure.

In anticipation of nonlinear solid mechanics applications, however, we will find it necessary to expand this notational framework so that large deformation of solids can be accommodated. Fortunately, provided certain interpretations are kept in mind, the form of the governing equations is largely unchanged by the generalization of the linear elastic system. This chapter therefore provides an introduction to how this generalization can be made. However, it will be seen that the continuum description and constitutive modeling of solids undergoing large deformations are complex topics that should be understood in detail before formulating and implementing numerical strategies. The closely related topics of nonlinear continuum mechanics and constitutive modeling will therefore be the subjects of subsequent chapters, followed by significant discussions of numerical strategies.

We conclude with a short list of references the reader may find useful as background material. Throughout the text, we assume little or no familiarity with either the finite element method or nonlinear solid mechanics, but we do assume a basic level of familiarity with the mechanics of materials, linear continuum mechanics, and linear elasticity. The last section of this chapter provides some basic references in these areas for those wishing to fill gaps in knowledge.

1.2. Linear Structural Component

We consider the simple axial (or in structural terms, truss) member shown schematically in Fig. 1.1. We can think of this member as a straight bar of material, whose transverse dimensions are small compared to its overall length, and which can only transmit loads in the axial direction. Real-world examples include taut cables in tension, truss members, and similar rod-like objects.

../_images/nonlinear_behavior-fig01.png

Fig. 1.1 Axial model problem: schematic and local coordinate system.

We index the material with coordinates \(\mathbf{x}\) with values between \(0\) and \(L_{0}\). Assuming that all displacement of the rod occurs in the axial direction, we write this displacement as \(u(\mathbf{x},t)\), with \(t\) signifying time. The infinitesimal strain or engineering strain at any point \(\mathbf{x} \in (0,L_{0})\) is given by

(1.1)\[\epsilon _{E}(\mathbf{x},t)=\frac{\partial }{\partial\mathbf{x}} u(\mathbf{x},t) .\]

The true stress \(\sigma_{T}\) at any point in the bar and at any instant is described via

(1.2)\[\sigma_{T}(\mathbf{x},t)=\frac{P(\mathbf{x},t)}{A(\mathbf{x},t)},\]

where \(P\) is the total axial force acting at location \(\mathbf{x}\) and \(A\) is the current cross-sectional area at that location. If the cross-sectional area does not change very much as a result of the deformation, it is appropriate to define the nominal stress or engineering stress as

(1.3)\[\sigma_{E}=\frac{P(\mathbf{x},t)}{A_{0}(\mathbf{x})},\]

where \(A_{0}(\mathbf{x})\) is the initial cross-sectional area at point \(\mathbf{x}\). If the material behaves in a linear elastic manner then \(\sigma_{E}\) and \(\epsilon_{E}\) are related via

(1.4)\[\sigma_{E}=E\epsilon_{E} ,\]

where \(E\) is the elastic modulus, or Young’s modulus, for the material.

To begin we consider the case of static equilibrium where inertial effects are either negligible or nonexistent and the response is therefore independent of time. One can in this case suppress the time argument in (1.2) and (1.4). The balance of linear momentum for the static case is expressed at each point \(\mathbf{x}\) by

(1.5)\[\frac{\mathrm d}{\mathrm d x} \left( A_{0}(\mathbf{x})\sigma_{E}(\mathbf{x}) \right) = f(\mathbf{x}),\]

where \(f\) is the applied external body loading, assumed to be axial, with units of force per unit length. Substitution of (1.4) into (1.5) gives the following ordinary differential equation for \(u(\mathbf{x})\) on the domain \((0,L_{0})\):

(1.6)\[\frac{\mathrm d}{\mathrm d x} \left( EA_{0} \frac{\mathrm d}{\mathrm d x}(u) \right) = f(\mathbf{x}).\]

If we assume that the cross-section is uniform so that \(A_{0}\) does not vary with \(\mathbf{x}\), and that the material is homogeneous so that \(E\) does not vary throughout the rod, then

(1.7)\[EA_{0} \frac{\mathrm d^{2}}{\mathrm dx^{2}} u(\mathbf{x}) = f(\mathbf{x}),\]

We note that (1.7) is a linear, second order differential equation for the unknown displacement field \(u\). In order to pose a mathematical problem that can be uniquely solved it is necessary to pose two boundary conditions on the unknown \(u\). We will be interested primarily in two types, corresponding to prescribed displacement and prescribed force (or stress) boundary conditions. An example of a displacement boundary condition would be

(1.8)\[u(0) = \bar{u},\]

while an example of a force boundary condition is

(1.9)\[\sigma_{E}(L_{0}) = E \frac{\mathrm d u}{\mathrm d x} \left( L_{0} \right) = \bar{\sigma},\]

where \(\bar{u}\) and \(\bar{\sigma}\) are prescribed values for the displacement and axial stress at the left and right bar ends, respectively. In mathematics parlance, the boundary condition in (1.8) is called a Dirichlet boundary condition while the boundary condition represented by (1.9) is a Neumann boundary condition. Dirichlet boundary conditions involve the unknown independent variable itself, while Neumann boundary conditions are expressed in terms of its derivatives.

Virtually any combination of such boundary conditions can be applied to our problem, but only one boundary condition (either a Neumann or Dirichlet condition) can be applied at each endpoint. In the case where Neumann (stress) conditions are applied at both ends of the bar, the solution \(u(\mathbf{x})\) is only determinable up to an arbitrary constant (this fact can be verified by applying separation of variables to (1.7)).

We now consider a particular case of this linear problem that will be useful in considering some of the various nonlinearities to be discussed below. In particular, suppose \(f=0\) on the domain \((0,L_{0})\), and furthermore consider the boundary conditions

(1.10)\[u=0 \text{ at } x=0,\]

and

(1.11)\[\sigma _{E} = \frac{ F^{ext} }{ A_{0} } \text{ at } x=L,\]

where \(F^{ext}\) is an applied force on the right end of the rod.

In this case, examination of (1.5) yields

(1.12)\[A_{0}\frac{\mathrm d}{\mathrm d x} \left( \sigma _{E}(\mathbf{x}) \right) = 0,\]

meaning that \(\sigma _{E}\) does not vary along the length of the rod. Since \(\sigma _{E}\) is proportional to \(\epsilon _{E}\) (see (1.4)), the strain must also be a constant value along the rod length.

Finally, in view of (1.1) we conclude that \(u(\mathbf{x}\)) must vary linearly with \(\mathbf{x}\). In other words, we know that the solution \(u(\mathbf{x}\)) must take the form

(1.13)\[u(\mathbf{x}) = u(0) + \delta \left( \frac{x}{L} \right) = \delta \left( \frac{x}{L} \right) ,\]

where \(\delta\) is the elongation, or difference between the left and right end displacement. The problem therefore reduces to finding the elongation produced by the applied force \(F^{ext}\). This problem is trivially solved and leads to the familiar linear relationship between \(F^{ext}\) and \(\delta\):

(1.14)\[\frac{ E A_{0} }{ L_{0} } \delta = F^{ext} ;\]

in other words, we have a simple linear spring with stiffness \(E A_{0} / L_{0}\). After solving for \(\delta\) one may merely substitute (1.13) to obtain the desired expression for \(u(\mathbf{x})\).

1.3. Material Nonlinearity

We examine the case of a material nonlinearity by replacing (1.4) with generic relationship between \(\sigma_{E}\) and \(\epsilon_{E}\),

(1.15)\[\sigma _{E} = \hat{\sigma} (\epsilon_{E}) ,\]

where \(\hat{\sigma}\) is a smooth and generally nonlinear function, see Fig. 1.2.

../_images/nonlinear_behavior-fig02.png

Fig. 1.2 Schematic of a nonlinear one-dimensional stress-strain relation.

We make few restrictions on the specific form of \(\hat{\sigma}\), other than to assume that \(\frac{ \mathrm d }{\mathrm d \epsilon_{E} } \hat{\sigma} > 0\) for all values of \(\epsilon_{E}\). If we retain the assumption that \(f=0\) and impose boundary conditions (1.10) and (1.11) then (1.12) is still valid, i.e.,

(1.16)\[\sigma _{E} = \frac{F^{ext}}{A_{0}}\]

throughout the rod. Furthermore, since we assume that a one-to-one relation exists between \(\sigma _{E}\) and \(\epsilon_{E}\), we conclude that, just as in the linear material case, the strain is a constant value in the rod given by

(1.17)\[\epsilon_{E} = \frac{\delta}{L_{0}}.\]

We can solve the problem by finding \(\delta\) as before, but now we must solve the nonlinear equation

(1.18)\[A_{0} \hat{\sigma} \left( \frac{\delta}{L_{0}} \right) = F^{ext}.\]

We can express (1.18) as an equation for the displacement at the right end which we denote as \(d_{L} = u(L)\). We can write

(1.19)\[N \left( d_{0} \right) = F^{ext},\]

where \(N \left( d_{0} \right)\) is a nonlinear function of the unknown \(d_{L}\) defined in this case as

(1.20)\[N \left( d_{0} \right) = A_{0} \hat{\sigma} \left( \frac{d_{L}}{L_{0}} \right).\]

In general, (1.20) will not have a closed-form solution and some sort of iterative procedure is necessary. Nonlinear equation solving is discussed at length in Chapter Section 13. Here we resort to one of the more recognized and widely-used procedures, Newton-Raphson iteration. In this method one introduces a set of indices k corresponding to the iterations, and given a current iterate \(d^{k}_{L}\), a first-order Taylor series expansion of (1.20) is utilized to generate the next iterate \(d^{k+1}_{L}\) via

(1.21)\[0 = F^{ext} - N \left( d^{k+1}_{L} \right) \approx F^{ext} - \left( N \left( d^{k}_{L} \right) + \frac{\mathrm d}{\mathrm d d_{L}} N \left( d^{k}_{L} \right) \Delta d_{L} \right),\]

where

(1.22)\[d^{k+1}_{L} = d^{k}_{L} + \Delta d_{L} .\]

(1.21) can be expressed more compactly as

(1.23)\[K \left( d^{k}_{L} \right) \Delta d_{L} = R \left( d^{k}_{L} \right),\]

where \(R \left( d^{k}_{L} \right)\), the residual or out-of-balance force, is given by

(1.24)\[R \left( d^{k}_{L} \right) := F^{ext} - N \left( d^{k}_{L} \right),\]

and \(K \left( d^{k}_{L} \right)\), the incremental or tangent stiffness, is written as

(1.25)\[K \left( d^{k}_{L} \right) := \frac{\mathrm d}{\mathrm d d_{L}} N \left( d^{k}_{L} \right).\]

The Newton-Raphson procedure is then carried out by recursively solving (1.23) and (1.22).

1.4. Geometric Nonlinearity

Geometric nonlinearities are induced by nonlinearities in the kinematic description of the system at hand. We will identify and work with several nonlinearities of this general type in great detail in Section 4, Section 5, and Section 6, but to begin we consider two particular cases in the context of our simple model problem.

The first type of nonlinearity we consider is introduced by the use of nonlinear strain and stress measures in definition of the stress-strain relation. As an example, let us consider alternatives to (1.1) and (1.3), which defined the engineering stain \(\epsilon _{E}\) and engineering stress \(\sigma _{E}\) that we have utilized to this point. When used in our model problem with \(f=0\) and boundary conditions (1.10) and (1.11), we have seen that the engineering strain does not vary over the rod’s length, having a constant value \(\delta/L_{0}\). For this strain measure to be appropriate, the deformation \(\delta\) should be infinitesimal. In the presence of larger deformations, the true strain or logarithmic strain is often used,

(1.26)\[\epsilon _{T} = \int_{L_0}^L \frac{\mathrm d \gamma}{\gamma} = \log \left( \frac{L}{L_{0}} \right) = \log \left( 1 + \epsilon _{E} \right).\]

Similarly, if the cross-sectional area \(A\) changes appreciably during the process, it is likely that the engineering stress \(\sigma _{E}\) should be replaced by the true stress \(\sigma _{T}\) defined in (1.2). In the case of our model problem, this would imply

(1.27)\[\sigma _{T} = \frac{F^{ext}}{A} ,\]

where \(A\) is to be interpreted as the cross-sectional area in the final (deformed) configuration.

Relating this area to the elongation \(\delta\) requires a constitutive assumption to be made. For example, if we assume the rod consists of an elastic material, we could approximate this variation by considering the area to vary according to Poisson’s effect. This would require that for each differential increment \(\mathrm d \epsilon _{T}\) in the axial true strain, each lateral dimension should change by a factor of \((1-\nu \mathrm d \epsilon _{T})\), where \(\nu\) is Poisson’s ratio for the material.

At a given instant of the loading process, therefore, an incremental change in the area \(A\) can be approximated via

(1.28)\[\begin{split}A + \mathrm d A & = \left(1-\nu \mathrm d \epsilon _{T} \right)^{2} A \notag \\ & \approx \left(1-2\nu \mathrm d \epsilon _{T} \right)A.\end{split}\]

(1.28) implies that

(1.29)\[\begin{split}\frac{1}{A} \mathrm d A & = -2\nu \mathrm d \epsilon _{T} \notag \\ & = -2 \nu \frac{\mathrm d \epsilon _{T}}{\mathrm d L} \mathrm d L \notag \\ & = -2 \nu \left( \frac{1}{L} \right) \mathrm d L.\end{split}\]

Integrating (1.29) between the initial and the final configurations gives

(1.30)\[A = A_{0} \left( \frac{L_{0}}{L} \right)^{2\nu} = A_{0} \left( \frac{L_{0}}{L_{0}+\delta} \right)^{2\nu} .\]

If we assume Hooke’s Law,

(1.31)\[\sigma _{T} = E \epsilon _{T},\]

we can use (1.26), (1.27), and (1.30) to conclude that

(1.32)\[EA_{0} \log \left( \frac{L_{0} + \delta}{L_{0}} \right) \left( \frac{L_{0}}{L_{0}+ \delta} \right)^{2\nu} = F^{ext},\]

which is a nonlinear equation governing the elongation \(\delta\). Note that this nonlinearity is not caused by any sort of nonlinear stress-strain relation, but instead results from the observation that the amount of deformation may not be small, necessitating more general representations of stress and strain.

The second sort of nonlinearity we wish to consider is that caused by large superimposed rigid body rotations and translations that introduce nonlinearities into many problems even when the strains in the material are well-approximated by infinitesimal measures. Toward this end we refer to Fig. 1.3, in which we embed our one-dimensional truss element in a two-dimensional frame. We locate one end of the rod at the origin and consider this end to be pinned so that it is free to rotate but not translate. The other end of the rod, initially located at coordinates \(x^{0}_{1},\ x^{0}_{2}\), is subjected to a (vector valued) force \(\mathbf{F}^{ext}\), which need not be directed along the axis of the rod.

../_images/nonlinear_behavior-fig03.png

Fig. 1.3 Model problem with infinitesimal motions superposed on large rigid body motions.

We note that under the restriction of small motions this problem is ill-posed because the rod is incapable of transmitting anything but axial force (\(\mathbf{F}^{ext}\) would need to act in the axial direction). However, in the current context we allow unlimited rotation with the result that the rod will rotate until it aligns with \(\mathbf{F}^{ext}\) in its equilibrium condition. In fact this observation allows us to guess the solution to the problem. Since we assume that the axial response of the rod is completely linear, we may deduce that the final elongation is given by

(1.33)\[\delta = \frac{L_{0} \| \mathbf{F}^{ext} \| }{EA_{0}},\]

where \(\| \mathbf{F}^{ext} \|\) denotes the Euclidean length of the vector \(\mathbf{F}^{ext}\). The final orientation of the rod must coincide with the direction \(\mathbf{F}^{ext}\), so we can write the final position of the end of the rod, using the coordinates \(\left(x^{f}_{1}, x^{f}_{2}\right)\), as

(1.34)\[\begin{split}\begin{bmatrix} x^{f}_{1} \\ x^{f}_{2} \end{bmatrix} = \frac{L_{0}}{\| \mathbf{F}^{ext} \|} \left( 1+ \frac{\| \mathbf{F}^{ext} \|}{EA_{0}} \right) \begin{bmatrix} F^{ext}_{1} \\ F^{ext}_{2} \end{bmatrix} ,\end{split}\]

or, writing the solution in terms of the rod end displacements \(d_{1}\) and \(d_{2}\),

(1.35)\[\begin{split}\begin{bmatrix} {d_{1}} \\ {d_{2}} \end{bmatrix} = \frac{L_{0}}{\| \mathbf{F}^{ext} \|} \left( 1+ \frac{\| \mathbf{F}^{ext} \|}{EA_{0}} \right) \begin{bmatrix} {F^{ext}_{1}} \\ {F^{ext}_{2}} \end{bmatrix} - \begin{bmatrix} {x^{0}_{1}} \\ {x^{0}_{2}} \end{bmatrix} .\end{split}\]

It is instructive to proceed as though we do not know the solution summarized in (1.35) and formulate the equilibrium equations governing \(d_{1}\) and \(d_{2}\).

If we observe that the elongation \(\delta\) of the rod can be written as

(1.36)\[\delta = \sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^{0}_{2} \right)^{2} } \, -\, L_{0},\]

then (1.33) gives the relationship between \(\| \mathbf{F}^{ext} \|\) and the unknown displacements. Furthermore, as noted above, the direction of \(\mathbf{F}^{ext}\) is given by

(1.37)\[\begin{split}\frac{\mathbf{F}^{ext}}{ \| \mathbf{F}^{ext} \|} = \frac{1}{\sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^{0}_{2} \right)^{2} }} \begin{bmatrix} d_1 + x^{0}_{1} \\ d_2 + x^{0}_{2} \end{bmatrix} .\end{split}\]

Combining these facts gives the equation that governs \(d_1\) and \(d_2\),

(1.38)\[\begin{split}\begin{bmatrix} F^{ext}_1 \\ F^{ext}_2 \end{bmatrix} = EA_0 \frac{\sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^{0}_{2} \right)^{2} } -L_0}{L_0 \sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^ {0}_{2} \right)^{2} }} \begin{bmatrix} d_1 + x^{0}_{1} \\ d_2 + x^{0}_{2} \end{bmatrix} .\end{split}\]

The reader may wish to verify this equation by substituting the solution (1.35) into (1.38).

(1.38) is a nonlinear, vector-valued equation for the unknowns \(d_1\) and \(d_2\). Recalling the generic form for nonlinear equations we introduced in the one dimensional case in (1.19), we could write this generically as

(1.39)\[\mathbf{N} ( \mathbf{d} ) = \mathbf{F}^{ext} ,\]

where

(1.40)\[\begin{split}\mathbf{d} = \begin{bmatrix} d_1 \\ d_2 \end{bmatrix}\end{split}\]

and

(1.41)\[\begin{split}\mathbf{N} ( \mathbf{d} ) := EA_0 \frac{\sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^{0}_{2} \right)^{2} } -L_0}{L_0 \sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^ {0}_{2} \right)^{2} }} \begin{bmatrix} d_1 + x^{0}_{1} \\ d_2 + x^{0}_{2} \end{bmatrix} .\end{split}\]

Just as was done in the last section for the one degree of freedom case, we could introduce a Newton-Raphson strategy to solve (1.39) via

(1.42)\[\mathbf{K} ( \mathbf{d}^k ) \Delta d = **R** ( \mathbf{d}^k ) = \mathbf{F}^{ext} - \mathbf{N} ( \mathbf{d}^k ),\]

and

(1.43)\[\mathbf{d}^{k+1} = \mathbf{d}^k + \Delta \mathbf{d} ,\]

where

(1.44)\[\begin{split}\mathbf{K} ( \mathbf{d}^k ) := \frac{\partial \mathbf{N}} {\partial \mathbf{d}} ( \mathbf{d}^k ) = \left \lbrack \vphantom{\begin{array}{cc} \frac{\partial N_1}{\partial d_1} & \frac{\partial N_1}{\partial d_2} \\ \Big ( \Big ) \\ \frac{\partial N_2}{\partial d_1} & \frac{\partial N_2}{\partial d_2} %\\ \end{array}} \begin{array}{cc} \frac{\partial N_1}{\partial d_1} & \frac{\partial N_1}{\partial d_2} \\ \\ \frac{\partial N_2}{\partial d_1} & \frac{\partial N_2}{\partial d_2} %\\ \end{array} \right \rbrack _{\mathbf{d} = \mathbf{d}^k}\end{split}\]

Carrying out the calculation of \(\mathbf{K} ( \mathbf{d}^k )\) for the specific \(\mathbf{N} ( \mathbf{d} )\) at hand gives

(1.45)\[\mathbf{K} ( \mathbf{d}^k ) = \mathbf{K}_{direct} ( \mathbf{d}^k ) +\mathbf{K}_{geom} ( \mathbf{d}^k ).\]

\(\mathbf{K}_{direct} ( \mathbf{d}^k )\) is given by

(1.46)\[\begin{split}\mathbf{K}_{direct} ( \mathbf{d}^k ) := \frac{EA_0}{\sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^{0}_{2} \right)^{2} }} \begin{bmatrix} \left( d^k_{1}+x^{0}_{1} \right)^{2} & \left( d^k_{1}+x^{0}_{1} \right)\left( d^k_{2}+x^{0}_{2} \right) \\ \left( d^k_{1}+x^{0}_{1} \right)\left( d^k_{2}+x^{0}_{2} \right) & \left( d^k_{2}+x^{0}_{2} \right)^{2} \end{bmatrix}\end{split}\]

and \(\mathbf{K}_{geom} ( \mathbf{d}^k )\) by

(1.47)\[\begin{split}\mathbf{K}_{geom} ( \mathbf{d}^k ) := EA_0 \left( \frac{1}{L_0} - \frac{1}{\sqrt { \left( d_{1}+x^{0}_{1} \right)^{2} + \left( d_{2}+x^{0}_{2} \right)^{2} }} \right) \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} .\end{split}\]

As the notation suggests, \(\mathbf{K}_{direct}\) is sometimes referred to as the direct stiffness, or that part of the stiffness emanating directly from the material stiffness of the system at hand. \(\mathbf{K}_{geom}\), on the other hand, is sometimes called the geometric stiffness, and arises not from inherent stiffness of the material but by virtue of the large motions in the problem.

To gain insight into these issues in the current context, consider the case where \(\| \mathbf{d}^k \| \ll \| \mathbf{x}^0 \|\), the case where the motions are small in comparison to the rod’s length. In this case we find

(1.48)\[\mathbf{K}_{geom} ( \mathbf{d}^k ) \rightarrow 0,\]

and

(1.49)\[\begin{split}\mathbf{K}_{direct} ( \mathbf{d}^k ) \rightarrow \frac{EA_0}{L_0} \begin{bmatrix} \cos \theta \cos \theta & \cos \theta \sin \theta \\ \cos \theta \sin \theta & \cos \theta \cos \theta \end{bmatrix} ,\end{split}\]

where \(\theta = \arctan \left( \frac{x^{0}_2} {x^{0}_1} \right)\) is the angle between the original axis of the rod and the positive \(x\)-axis. In other words, when the motions become small, the geometric stiffness vanishes and the direct stiffness reduces to the familiar stiffness matrix associated with a two-dimensional truss member.

1.5. Contact Nonlinearity

A final type of nonlinearity we wish to consider is that created due to contact with another deformable or rigid body. As a simple model problem for this case we refer to Fig. 1.4, where we consider a prescribed motion \(\bar{d}\) of the left end of our one-dimensional rod and solve for the static equilibrium of the unknown displacement \(d\) of the right end, subject to the constraint

(1.50)\[g(d) = d-g_0 \leq 0,\]

where \(g_0\) is the initial separation, or gap, between the right end of the rod and the rigid obstacle.

../_images/nonlinear_behavior-fig04.png

Fig. 1.4 Schematic of the rigid obstacle problem.

Even if we assume that the motions are small and the material response of the rod is elastic, the equations governing the response of our rod are nonlinear. To see this, let us choose \(d\) as our unknown and construct the following residual \(R(d)\) for our system:

(1.51)\[R(d) = \frac{EA_0}{L_0} (d- \bar{d} ) +F_c,\]

Here \(F_c\), the contact force between the obstacle and the rod (assumed positive in compression), is subject to the constraints

(1.52)\[F_c \geq 0; \: g(d) \leq 0 \: \text{ and } F_c g(d) = 0 \: .\]

Equations (1.52) are called Kuhn-Tucker complementary conditions in mathematical parlance and physically require that the contact force be compressive, that the rod end not penetrate the obstacle, and that the contact force only be nonzero when \(g=0\), i.e. when contact between the rod and obstacle occurs. In fact \(F_c\) is a Lagrange multiplier in this problem, enforcing the kinematic constraint (1.50). We see that the condition operating on the right end of the bar is neither a Dirichlet nor a Neumann boundary condition; in fact, both the stress and the displacement at this point are unknown but are related to each other through the constraints expressed in Equations (1.52).

Plots of the residual defined by Equations (1.51) and (1.52) are given in Fig. 1.5 for the two distinct cases of interest: where contact does not occur (when \(\bar d < g_0\)) and where contact does occur (when \(\bar d \geq g_0\)). The solutions (i.e. the zeros of \(R\)) are readily apparent. When no contact occurs \(d= \bar d\), while in the case of contact \(d=g_0\). The internal stresses generated in the bar are then readily deduced.

One may note from Fig. 1.5 some important practical features of this problem. First, in both cases the admissible region for \(d\) is restricted to be less than \(g_0\). Second, at the value \(d=g_0\), each diagram shows the residual be multiple valued, which is a direct consequence of the fact that in this condition (i.e., where \(g=0\)), \(F_c\) can be any positive number.

../_images/nonlinear_behavior-fig05.png

Fig. 1.5 Plots of residuals verses displacement for the rigid obstacle problem: (a) the case where \(\bar d < g_0\) (no contact); (b) the case where \(\bar d \geq g_0\).

Finally, although the solution to our simple model problem is readily guessed, we can see from both cases that the plot of \(R\) versus \(d\) is only piecewise linear; the kink in each diagram indicates the fact that a finite tangent stiffness operates when contact is not active, changing to an infinite effective stiffness imposed by Equations (1.52) when contact between the rod and obstacle is detected. This contact detection therefore becomes an important feature in general strategies for contact problems, and introduces both nonlinearities and non smoothness into the global equations as this rather simple example demonstrates.

The books [[1], [2], [3], [4], [5], [6], [7]] are suggested for those readers wishing to reinforce their knowledge of linear elasticity, elementary continuum mechanics, and/or fundamentals of solid mechanics. They are presented in alphabetical order, with no other significance to be attached to the order of presentation.