Frictional contact of flexible and rigid bodies

 This paper is about planar frictional contact problems of both flexible and rigid bodies. For the flexible case a nonlinear finite element formulation is presented, which is based on a modified Coulomb friction law. Stick-slip motion is incorporated into the formulation through a radial return mapping scheme. Linearly interpolating four node elements and three node contact elements are utilized for the finite element discretization. The corresponding tangent stiffness matrices and residual vectors of the equations of motion are presented. In the rigid body case the contact problem is divided into impact and continual contact, which are mathematically described by linear complementarity problems. The impact in normal direction is modeled by a modified Poisson hypothesis, which is adapted to allow multiple impacts. The formulation of the tangential impact is grounded on Coulombs law of friction. The normal contact forces of the continual contact are such that colliding bodies are prevented from penetration and the corresponding tangential forces are expressed by Coulombs law of friction. Examples and comparisions between the different methods are presented. Received: 10 January 2001     

Direct computation of instability points for contact problems

 The extended system is known as a reliable algorithm for the direct computation of instability points on the equilibrium path of mechanical structures. This article describes the application of the extended system as critical point computation method to mechanical contact problems. In this type of problems inequality constraints have to be considered. Moreover a prediction method based on the extended system algorithm is presented which allows the detection of favorable starting values for a critical point computation on the equilibrium path. Dedicated to the memory of Prof. Mike Crisfield, for his cheerfulness and cooperation as a colleague and friend over many years.     

Adaptive multi-analysis strategy for contact problems with friction

The objective of the work presented here is to develop an efficient strategy for the parametric analysis of bolted joints designed for aerospace applications. These joints are used in elastic structural assemblies with local nonlinearities (such as unilateral contact with friction) under quasi-static loading. Our approach is based on a decomposition of an assembly into substructures (representing the parts) and interfaces (representing the connections). The problem within each substructure is solved by the finite element method, while an iterative scheme based on the LATIN method (Ladevèze in Nonlinear computational structural mechanics—new approaches and non-incremental methods of calculation, 1999) is used for the global resolution. The proposed strategy consists in calculating response surfaces (Rajashekhar and Ellingwood in Struct Saf 12:205–220, 1993) such that each point of a surface is associated with a design configuration. Each design configuration corresponds to a set of values of all the variable parameters (friction coefficients, prestresses) which are introduced into the mechanical analysis. Here, instead of carrying out a full calculation for each point of the surface, we propose to use the capabilities of the LATIN method and reutilize the solution of one problem (for one set of parameters) in order to solve similar problems (for the other sets of parameters) (Boucard and Champaney in Int J Numer Methods Eng 57:1259–1281, 2003). The strategy is adaptive in the sense that it takes into account the results of the previous calculations. The method presented can be used for several types of nonlinear problems requiring multiple analyses: for example, it has already been used for structural identification (Allix and Vidal in Comput Methods Appl Mech Eng 191:2727–2758, 2001).     

BEM and FEM analysis of Signorini contact problems with friction

In this paper, we present a comparative study of the boundary element method (BEM) and the finite element method (FEM) for analysis of Signorini contact problems in elastostatics with Coulomb's friction law. Particularities of each method and comparison with the penalty method are discussed. Numerical examples are included to demonstrate the present formulations and to highlight its performance.     

A mixed least squares method for solving problems in linear elasticity: theoretical study

 In a previous paper we proposed a mixed least squares method for solving problems in linear elasticity. The solution to the equations of linear elasticity was obtained via minimization of a least squares functional depending on displacements and stresses. The performance of the method was tested numerically for low order elements for classical examples with well known analytical solutions. In this paper we derive a condition for the existence and uniqueness of the solution of the discrete problem for both compressible and incompressible cases, and verify the uniqueness of the solution analytically for two low order piece-wise polynomial FEM spaces. Received: 20 January 2001 / Accepted: 14 June 2002 The authors gratefully acknowledge the financial support provided by NASA George C. Marshall Space Flight Centre under contract number NAS8-38779.     

An integral equation approach to the inclusion-crack interactions in three-dimensional infinite elastic domain

 In this paper, an integral equation method to the inclusion-crack interaction problem in three-dimensional elastic medium is presented. The method is implemented following the idea that displacement integral equation is used at the source points situated in the inclusions, whereas stress integral equation is applied to source points along crack surfaces. The displacement and stress integral equations only contain unknowns in displacement (in inclusions) and displacement discontinuity (along cracks). The hypersingular integrals appearing in stress integral equation are analytically transferred to line integrals (for plane cracks) which are at most weakly singular. Finite elements are adopted to discretize the inclusions into isoparametric quadratic 10-node tetrahedral or 20-node hexahedral elements and the crack surfaces are decomposed into discontinuous quadratic quadrilateral elements. Special crack tip elements are used to simulate the variation of displacements near the crack front. The stress intensity factors along the crack front are calculated. Numerical results are compared with other available methods. Received: 28 January 2002 / Accepted: 4 June 2002 The work described in this paper was partially supported by a grant from the Research Grant Council of the Hong Kong Special Administration Region, China (Project No.: HKU 7101/99 E).     

A general algorithm for the numerical evaluation of nearly singular boundary integrals of various orders for two- and three-dimensional elasticity

 A general algorithm of the distance transformation type is presented in this paper for the accurate numerical evaluation of nearly singular boundary integrals encountered in elasticity, which, next to the singular ones, has long been an issue of major concern in computational mechanics with boundary element methods. The distance transformation is realized by making use of the distance functions, defined in the local intrinsic coordinate systems, which plays the role of damping-out the near singularity of integrands resulting from the very small distance between the source and the integration points. By taking advantage of the divergence-free property of the integrals with the nearly hypersingular kernels in the 3D case, a technique of geometric conversion over the auxiliary cone surfaces of the boundary element is designed, which is suitable also for the numerical evaluation of the hypersingular boundary integrals. The effects of the distance transformations are studied and compared numerically for different orders in the 2D case and in the different local systems in the 3D case using quadratic boundary elements. It is shown that the proposed algorithm works very well, by using standard Gaussian quadrature formulae, for both the 2D and 3D elastic problems. Received: 20 November 2001 / Accepted: 4 June 2002 The work was supported by the Science Foundation of Shanghai Municipal Commission of Education.     

A three-dimensional FE analysis of large deformations for impact loadings using tetrahedral elements

 A three-dimensional dynamic program for the anaysis of large deformations in contact-penetration problems is developed using the finite element Lagrangian method with explicit time integration. By incorporating a tetrahedral element, which allows a single-point integration without a special hourglass control scheme, this program can be more effective to the present problem. The position code algorithm is used to search contact surface. Eroding surfaces are also considered. The defense node algorithm was slightly modified for the calculation of contact forces. A study of obliquity effects on metallic plate perforation and ricochet processes in thin plates impacted by a sphere was conducted. It is well simulated that on separation of two parts of the sphere, the portion still within the crater tends to perforate, while the portion in contact with the plate surface ricochets. This deformation pattern is observed in experiments, especially at high obliquities. A long rod that impacts an oblique steel plate at high impact velocity was also simulated in order to study the dynamics of the rod caused by the three dimensional asymmetric contact. The agreement between simulated and experimental results is quite good. Fracture phenomena occuring at high obliquity deserves further investigations. Received: 20 February 2002 / Accepted: 20 September 2002     

A new boundary element approach for contact problems with friction

A new boundary element solution algorithm for two-dimensional and axisymmetric contact problems with friction, based on an independent discretization of the contacting surfaces and under static and proportional loading conditions, is presented. The solution procedure uses the element shape functions to distribute the geometry, tractions and displacements on each contact element. The contact constraints are then applied between each contacting node and the opposite contact segment. The overall boundary element matrix equations for the contacting bodies are coupled using the contact conditions at the interface without introducing any additional variables into the solution matrix. The algorithm is applied to several two-dimensional and axisymmetric frictional contact examples and the results obtained are in very good agreement with finite element and analytical solutions.     

A symmetric Galerkin BEM implementation for 3D elastostatic problems with an extension to curved elements

 Like the finite element method (FEM), the symmetric Galerkin boundary element method (SGBEM) can produce symmetric system matrices. While widely developed for two dimensional problems, the 3D-applications of the SGBEM are very rare. This paper deals with the regularization of the singular integrals in the case of 3D elastostatic problems. It is shown that the integration formulas can be extended to curved elements. In contrast to other techniques, the Kelvin fundamental solutions are used with no need to introduce the new kernel functions. The accuracy of the developed integration formulas is verified on a problem with known analytical solution. Received 6 November 2000     

Mesh matching and contact patch test

 Mesh matching discussed in this paper refers to a situation when two adjacent lines or surfaces with different meshes on them joint together. A mesh matching formulation is given by enforcing sticking conditions on the jointed lines and surfaces. The preservation of uniform displacements (stress distributions) across the matched meshes is investigated. This, in the contact context, is usually called patch test. The conditions to ensure the formulation pass the patch test are derived, that include partition of unity and compatible conditions of shape functions, as well as the inf-sup condition of the augmented equilibrium equations. The proposed formulation is a general one from which some special cases, such as one-pass and two-pass node-to-surface contact procedures, can be recovered. Various examples are given to show the patch test performances of lower and higher-order mesh matching (contact) formulations for both straight and curved inter faces in 1-D. For a general 2-D problem, an example demonstrates that, unlike its behaviours in 1-D, two-pass node-to-surface contact procedure using linear elements could not pass the patch test. Dedicated to the memory of Prof. Mike Crisfield, for his cheerfulness and cooperation as a colleague and friend over many years.     

A finite element formulation of three-dimensional contact problems with slip and friction

The paper describes a special finite element for three-dimensional, large displacement analysis of contact problems with slip and friction. This element may be used to model contact between several finite element bodies or contact between a finite element body and a flexible or rigid geometrical surface fixed in space or moving with time. The contact formulation is based on the concept of a spring-supported, moving disk that transfers normal contact forces and Coulomb friction forces. The contact surface has a finite, prescribed boundary.The contact element has been incorporated into the general-purpose, nonlinear, finite element program FENRIS. Three examples of its application are described in the paper.J. W. Simons was previously NTNF Fellow, Division of Technology, Trondheim, Norway     

A mode II crack problem with sliding contact in an orthotropic strip

The problem of a centrally cracked, linear elastic orthotropic strip loaded in bending by three point forces is analyzed and discussed. Coulomb friction is assumed between the crack faces to study the influence of the friction coefficient on the strain energy release rate. Under certain simplifying assumptions the problem is reduced to the solution of a singular integral equation which is evaluated numerically. The results are compared with the solution of the same problem obtained using the beam theory; limits of the application of beam theory for the reduction of experimental data are discussed.     

Solution of large deformation contact problems with friction between Blatz-Ko hyperelastic bodies

The present paper is devoted to the analysis of the contact problems with Coulomb friction and large deformation between two hyperelastic bodies. One approach to separate the material nonlinearity and contact nonlinearity is presented. The total Lagrangian formulation is adopted to describe the geometrically nonlinear behavior. Nondifferentiable contact potentials are regularized by means of the augmented Lagrangian method. Numerical examples are carried out in two cases: rigid-deformable contact and deformable-deformable contact with large slips. The numerical results prove that the proposed approach is robust and efficient concerning numerical stability.     

Efficient and accurate approach for powder compaction problems

 In this paper, a new approach for powder cold compaction simulations is presented. A density-dependent plastic model within the framework of finite strain multiplicative hyperelastoplasticity is used to describe the highly nonlinear material behaviour; the Coulomb dry friction model is used to capture friction effects at die-powder contact; and an Arbitrary Lagrangian–Eulerian (ALE) formulation is used to avoid the (usual) excessive distortion of Lagrangian meshes caused by large mass fluxes. Several representative examples, involving structured and unstructured meshes are simulated. The results obtained agree with the experimental data and other numerical results reported in the literature. It is shown that, contrary to other Lagrangian and adaptive h-remeshing approaches recently reported for this type of problems, the present approach verifies the mass conservation principle with very low relative errors (less than 1% in all ALE examples and exactly in the pure Lagrangian examples). Moreover, thanks to the use of an ALE formulation and in contrast with other simulations, the presented density distributions do not present spurious oscillations. Received: 20 March 2002 / Accepted: 15 October 2002 The partial financial support of the Ministerio de Ciencia y Tecnología (grant number DPI 2001-2204) is gratefully acknowledged.     

Solving thermal and phase change problems with the eXtended finite element method

 The application of the eXtended finite element method (X-FEM) to thermal problems with moving heat sources and phase boundaries is presented. Of particular interest is the ability of the method to capture the highly localized, transient solution in the vicinity of a heat source or material interface. This is effected through the use of a time-dependent basis formed from the union of traditional shape functions with a set of evolving enrichment functions. The enrichment is constructed through the partition of unity framework, so that the system of equations remains sparse and the resulting approximation is conforming. In this manner, local solutions and arbitrary discontinuities that cannot be represented by the standard shape functions are captured with the enrichment functions. A standard time-projection algorithm is employed to account for the time-dependence of the enrichment, and an iterative strategy is adopted to satisfy local interface conditions. The separation of the approximation into classical shape functions that remain fixed in time and the evolving enrichment leads to a very efficient solution strategy. The robustness and utility of the method is demonstrated with several benchmark problems involving moving heat sources and phase transformations. Received 20 May 2001 / Accepted 19 December 2001     

A point interpolation mesh free method for static and frequency analysis of two-dimensional piezoelectric structures

 A mesh free method called point interpolation method (PIM) is presented for static and mode-frequency analysis of two-dimensional piezoelectric structures. In the present method, the problem domain and its boundaries are represented by a set of properly scattered nodes. The displacements and the electric potential of a point are interpolated by the values of nodes in its local support domain using shape functions derived based on a point interpolation scheme. Techniques are discussed to surmount the singularity of the moment matrix. Variational principle together with linear constitutive piezoelectric equations is used to establish a set of system equations for arbitrary-shaped piezoelectric structures. These equations are assembled for all quadrature points and solved for displacements and electric potentials. A polynomial PIM program has been developed in MATLAB with matrix triangularization algorithm (MTA), which automatically performs a proper node enclosure and a proper basis selection. Examples are also presented to demonstrate the accuracy and stability of the present method and their results are compared with the conventional FEM results from ABAQUS as well as the analytical or experimental ones. Received: 6 February 2002 / Accepted: 5 August 2002     

A mixed finite element method for beam and frame problems

 In this work we consider solutions for the Euler-Bernoulli and Timoshenko theories of beams in which material behavior may be elastic or inelastic. The formulation relies on the integration of the local constitutive equation over the beam cross section to develop the relations for beam resultants. For this case we include axial, bending and shear effects. This permits consideration in a direct manner of elastic and inelastic behavior with or without shear deformation. A finite element solution method is presented from a three-field variational form based on an extension of the Hu–Washizu principle to permit inelastic material behavior. The approximation for beams uses equilibrium satisfying axial force and bending moments in each element combined with discontinuous strain approximations. Shear forces are computed as derivative of bending moment and, thus, also satisfy equilibrium. For quasi-static applications no interpolation is needed for the displacement fields, these are merely expressed in terms of nodal values. The development results in a straight forward, variationally consistent formulation which shares all the properties of so-called flexibility methods. Moreover, the approach leads to a shear deformable formulation which is free of locking effects – identical to the behavior of flexibility based elements. The advantages of the approach are illustrated with a few numerical examples. Dedicated to the memory of Prof. Mike Crisfield, for his cheerfulness and cooperation as a colleague and friend over many years.     

Stability analysis of an explicit finite element scheme for plane wave motions in elastic solids

 Expressions for critical timesteps are provided for an explicit finite element method for plane elastodynamic problems in isotropic, linear elastic solids. Both 4-node and 8-node quadrilateral elements are considered. The method involves solving for the eigenvalues directly from the eigenvalue problem at the element level. The characteristic polynomial is of order 8 for 4-node elements and 16 for 8-node elements. Due to the complexity of these equations, direct solution of these polynomials had not been attempted previously. The commonly used critical time-step estimates in the literature were obtained by reducing the characteristic equation for 4-node elements to a second-order equation involving only the normal strain modes of deformation. Furthermore, the results appear to be valid only for lumped-mass 4-node elements. In this paper, the characteristic equations are solved directly for the eigenvalues using <ty>Mathematica<ty> and critical time-step estimates are provided for both lumped and consistent mass matrix formulations. For lumped-mass method, both full and reduced integration are considered. In each case, the natural modes of deformation are obtained and it is shown that when Poisson's ratio is below a certain transition value, either shear-mode or hourglass mode of deformation dominates depending on the formulation. And when Poisson's ratio is above the transition value, in all the cases, the uniform normal strain mode dominates. Consequently, depending on Poisson's ratio the critical time-step also assumes two different expressions. The approach used in this work also has a definite pedagogical merit as the same approach is used in obtaining time-step estimates for simpler problems such as rod and beam elements. Received: 8 January 2002 / Accepted: 12 July 2002 The support of NSF under grant number DMI-9820880 is gratefully acknowledged.     

Parallelized FVM algorithm for three-dimensional viscoelastic flows

 A parallel implementation for the finite volume method (FVM) for three-dimensional (3D) viscoelastic flows is developed on a distributed computing environment through Parallel Virtual Machine (PVM). The numerical procedure is based on the SIMPLEST algorithm using a staggered FVM discretization in Cartesian coordinates. The final discretized algebraic equations are solved with the TDMA method. The parallelisation of the program is implemented by a domain decomposition strategy, with a master/slave style programming paradigm, and a message passing through PVM. A load balancing strategy is proposed to reduce the communications between processors. The three-dimensional viscoelastic flow in a rectangular duct is computed with this program. The modified Phan-Thien–Tanner (MPTT) constitutive model is employed for the equation system closure. Computing results are validated on the secondary flow problem due to non-zero second normal stress difference N ₂. Three sets of meshes are used, and the effect of domain decomposition strategies on the performance is discussed. It is found that parallel efficiency is strongly dependent on the grid size and the number of processors for a given block number. The convergence rate as well as the total efficiency of domain decomposition depends upon the flow problem and the boundary conditions. The parallel efficiency increases with increasing problem size for given block number. Comparing to two-dimensional flow problems, 3D parallelized algorithm has a lower efficiency owing to largely overlapped block interfaces, but the parallel algorithm is indeed a powerful means for large scale flow simulations. Received: 2 July 2002 / Accepted: 15 November 2002 This research is supported by an A^⋆STAR Grant EMT/00/011.