Constrained Large-Displacement Thermal Analysis

Temperature has significant effect on a wide range of physics and engineering applications that include manufacturing, automotive, machine, and aerospace systems. While expansion due to thermal loads in absence of boundary and motion constraints is considered stress-free process, such an expansion can lead to significant change in the reference-configuration geometry as well as material properties. Temperature effect cannot be ignored as evidenced by large number of investigations in the areas of solid and fluid mechanics [1–11]. In conventional computational approaches for solution of thermo-elasticity problems, thermal expansion has been assumed, for the most part, to produce volumetric changes without distortion or shear [12]. Nonetheless, volumetric change and thermal expansion can lack uniformity because of constraints or restrictions imposed on displacements of boundary points.

Therefore, in case of thermal-load applications, two different cases of displacements can be encountered. In the first case, motion of the boundary points of the continuum is not restricted. In this case, under assumption of homogeneous and isotropic material and constant thermal-expansion coefficients, all material points experience the same displacement-gradient change as the result of thermal-load application. This case of homogeneous continuum motion is referred to as an unconstrained thermal expansion. In a second and different scenario, the motion of material points is restricted due to boundary and/or motion constraints, and therefore, such points do not move freely when the thermal load is applied, as in case of a rod with its endpoints are not allowed to move during thermal-expansion process. This case is referred to in this paper as constrained expansion because displacements of material points are subjected to boundary and/or motion constraints when thermal loads are applied. Therefore, “constrained thermal expansion” refers to thermal expansion that is not stress-free because of kinematic restrictions.

In both cases of unconstrained and constrained thermal expansions, however, reference-configuration geometry needs to be accurately described. Consequently, geometry plays an important role in formulation of continuum governing equations. In computational geometry methods, control points are used with knot multiplicity and knot vectors to define shapes and degree of continuity of curves, surfaces, and volumes [13–18]. However, as discussed in the literature, use of computational geometry methods as mechanics analysis tools has limitations, which include rigid structure, use of nonmaterial points, and use of the concepts of knot vector and knot multiplicity that restrict the generality of formulating the multibody system (MBS) constraint equations [19].

Study of articulated mechanical systems (AMS) has been considered in mechanical-design literature separately from the vast field of thermal analysis. Articulated mechanical systems are subject to nonlinear motion constraints that must be enforced at position, velocity, and acceleration levels in the process of numerical integration of motion equations. Articulated mechanical systems, whose dynamics is governed by Lagrange–D'Alembert principle, often operate in high-temperature environment, as in the case of slider-crank mechanisms, commonly used in engines. Despite such common examples, the literature lacks an approach that allows solving thermo-elasticity problems of articulated mechanical systems. To address this limitation, a new approach is proposed in this investigation for solving thermo-elasticity problems, in which systems subjected to boundary and motion constraints may undergo arbitrarily large displacements. This approach is based on integration of Lagrange–D'Alembert principle and large-displacement thermal analysis based on multiplicative decomposition of position-gradient matrix.

In the computational multiplicative-decomposition approach proposed in this investigation, four configurations are used to define continuum displacements and solve thermo-elasticity problems [20–22]. The four configurations are conveniently described using FE absolute nodal coordinate formulation (ANCF), which allows for large displacements including rigid-body displacements, and employs position-gradient vectors with unique geometric meaning as tangents to coordinate lines [23]. Therefore, the computational approach proposed in this study allows applying thermal loads during constrained large displacements, does not impose restrictions on choice of thermal coefficients, accurately captures reference-configuration geometry, and allows for systematic integration with multibody system (MBS) algorithms. While ANCF finite elements have been widely used in many applications [24–52], such elements have not previously been used to provide objective evaluation for and comparison between unconstrained and constrained thermal expansions. The specific contributions of this study can be summarized as follows:

1. A new approach is proposed for the solution of thermo-elasticity problems, in which system components subjected to boundary and motion constraints may undergo arbitrarily large displacements. This approach, applicable to articulated mechanical systems, requires integration of Lagrange–D'Alembert principle and large-displacement thermal analysis based on multiplicative decomposition of position-gradient matrix.

2. Differences in kinematic description used in two different cases of unconstrained and constrained thermal expansions are explained. It is shown that in case of constrained thermal expansion, displacements of different material points as the result of applying thermal loads are not the same because of boundary constraints.

3. Configuration of constrained thermal expansion is introduced and used in the formulation of continuum kinematic equations. Use of strain additive decomposition is avoided to generalize the method to large-displacement and thermal analysis of MBS applications. To this end, use of Lagrange–D'Alembert principle is necessary for proper treatment of algebraic constraint equations enforced during system motion.

4. In the study of constrained thermal expansion considered in this investigation, four configurations are employed to develop kinematic equations: straight configuration, reference configuration, thermal-expansion configuration, and current configuration. The four-configuration approach allows changing temperature during system large displacement while capturing accurate reference-configuration geometry.

5. Fundamental differences between the approach used in this paper and conventional FE approaches used in thermal analyses are identified. As discussed in the paper, use of ANCF position-gradient vectors as nodal coordinates allows for describing reference and thermal-expansion configurations using nodal coordinates. Furthermore, fully parameterized ANCF beam elements capture Poisson effect, required for accurate modeling of cross section deformations.

6. Numerical results are presented to demonstrate using the approach presented in this investigation, compare between results obtained using assumptions of unconstrained and constrained thermal expansions, and compare results with results obtained using conventional FE approaches. MBS mechanism subjected to large displacements including finite rotations is considered to demonstrate integration of the proposed thermal-analysis approach with Lagrange–D'Alembert principle.

The transient heat conduction equation is a first-order partial differential equation, while the motion equations of the continuum are second-order ordinary differential equations obtained from the continuum-mechanics partial differential equations using approximation methods. In this study, the temperature profile is assumed known and the two equations are assumed decoupled, therefore, the solution of the first-order partial differential equation of the heat conduction is assumed known. Nonetheless, the heat equation can be solved using an interpolation similar to the ANCF interpolation that ensures continuity of the spatial derivatives. The use of such an interpolation for the temperature to solve the heat equation ensures the continuity of the temperature gradients at the nodal points when the heat equation is numerically solved.

In case of unconstrained thermal expansion, changes in gradients throughout the continuum are assumed the same for all material points; while in case of constrained expansion, material points may experience different displacements. These two different cases, unconstrained and constrained expansions, are further explained using the fully parameterized ANCF elements briefly reviewed in this section.

The basic difference between unconstrained and constrained expansions can be explained using displacement field of the fully parameterized ANCF finite elements. In general, ANCF displacement field is written as $r(x,t)=S(x)e(t)$⁠, where $r$ is the global position vector of an arbitrary FE point, $S$ is the FE shape-function matrix, $e$ is the FE vector of coordinates, $x=[x1x2x3]T$ is the vector of element spatial coordinates, and $t$ is time [40]. In case of ANCF elements, polynomial coefficients are replaced by position and position-gradient coordinates. Use of position-gradient coordinates eliminates need for using, as nodal coordinates, noncommutative finite rotations which may lack clear physical interpretations and are not directly related to curvature and torsion in case of large-displacement analysis [53].

For fully parameterized ANCF elements, the vector of coordinates at given node $k$ can be written in the spatial analysis as $ek=[rkTrx1kTrx2kTrx3kT]T, k=1,2,…,nn$⁠, where $nn$ is the number of element nodes, and $rxl=∂r/∂xl, l=1,2,3$⁠. The fact that independent polynomial coefficients can be replaced by ANCF position and position-gradient coordinates without encountering singularities implies that all ANCF coordinates are independent. Because nine position gradients define nine independent modes of displacements; three rotation and six deformation modes, ANCF kinematic description can be used to explain independence of shear and bending in case of beam problems [54].

are derivatives of the shape functions with respect to the FE spatial coordinates.

Unconstrained expansion refers in this study to displacement mode in which all element material points have same change in displacement gradients during thermal-expansion process. This mode of displacement is observed in case of stress-free thermal expansion without any constraints imposed on the motion of the continuum boundary points. In process of thermal expansion, the orientation of position-gradient vectors is assumed to remain unchanged, and therefore, there is no shearing. That is, gradient vectors experience only stretch and such a stretch leads to only volumetric change.

where $ro=[x1x2]T$ is the arbitrary-point position before expansion. This equation shows that all material points have the same change in displacement gradients. It is clear that in this case of unconstrained expansion, beam length and height change linearly as functions of $α1$ and $α2$⁠, respectively. Figure 1 shows an example of such unconstrained stretch. The element nodal positions after expansion are determined by integrating differential equation of Eq. (6).

In case of constrained expansions, motion of continuum boundary points during thermal-load application is constrained. That is, position of these boundary points cannot be determined from the differential relationship presented in the preceding section alone since independent coordinates must be identified first and used to determine dependent gradients. Material-point translations can be independent of gradients because an infinitesimal volume in planar analysis has six modes of displacements, stretch of gradients at a point does not always imply displacement of this particular point. That is strains at a point that has zero displacements can be different from zero; a fact that can be easily demonstrated using ANCF displacement field.

Clearly, this equation leads to definitions of positions of material points different from the definitions obtained in case of unconstrained expansion. Using this equation, one can also show the effect of the FE interpolation on position gradients within the element. Figure 2 shows an example of such constrained thermal expansion in case of the planar beam element. The results of this figure show that in case of constrained thermal expansion, nodal positions do not change, causing displacements of continuum interior points depicted by dashed lines. Effect of restrictions caused by boundary constraints becomes clear by comparing with results of interior displacements in case of unconstrained thermal expansion shown by dashed lines shown in Fig. 1. Figure 3 shows plots of positions and displacements of points along beam centerline $(x2=0)$ whose overall length does not change as result of the expansion. Figure 4 shows norm of $rx1$ as function of spatial coordinate $x1$ when $x2=0$⁠. Using Eq. (9), it can be shown that the norm of $rx2$ remains constant along the beam and is equal to $1+α2$⁠. Position-gradient vectors can be also defined using Eq. (9) as $rx1=[1+α1(1−6(ξ−ξ2))0]T$ and $rx2=[01+α2]T$⁠.

Conditions $r1=[00]T$ and $r2=[l0]T$ for nodal positions at thermal-analysis preprocessing stage do not imply that motion of nodes are constrained in any subsequent dynamic analysis. That is, the conditions $r1=[00]T$ and $r2=[l0]T$ are used to determine the displacement field after constant initial thermal expansion, but nodes are free to move in subsequent dynamic analysis as demonstrated in the section of numerical examples.

In case of constrained thermal expansion, boundary and motion constraints can in general be written in a vector form as $C(x,e,t)=0$⁠, where $C$ is vector function of $nc$ algebraic constraint equations that restrict motion of boundary or interior nodes, and $e$ is vector of FE mesh nodal coordinates which has dimension $n$⁠. While simple examples are used in this study to focus on constrained-expansion concepts, vector $C$ of constraint functions can be nonlinear and can depend on large-displacement coordinates. Motion constraints on boundary of the FE mesh can be used to eliminate degrees-of-freedom using Lagrange–D'Alembert procedure, which is the foundation of MBS algorithms. A virtual change in constraint functions $C$ leads to $Ceδe=0$⁠, where $Ce=∂C/∂e$ is constraint-Jacobian matrix. When ANCF elements are used, the vector of mesh coordinates $e$ can be written as $e=eo+ed$⁠, where $eo$ is constant vector of nodal coordinates in the reference configuration, and $ed$ is the displacement vector. Therefore, one can write $Ceδe=Ceδed=0$⁠. Because of the constraints, the vector of nodal displacements can be partitioned as independent coordinates $edi$ and dependent coordinates $edd$⁠, that is, $ed=[ediTeddT]T$⁠. Using this partitioning and assuming that constraint equations are linearly independent, one can write $Ceδed=Ceddδedd+Cediδedi=0$⁠, where $Cedd=∂C/∂edd$ is $nc×nc$ square nonsingular constraint Jacobian matrix associated with the dependent coordinates $edd$⁠, and $Cedi=∂C/∂edi$ is $nc×(n−nc)$ constraint Jacobian matrix associated with the independent coordinates $edi$⁠. One can, therefore, write $δedd=−Cedd−1Cediδedi$⁠. This equation, leads to $δe=[δediTδeddT]T=[I−(Cedd−1Cedi)T]Tδedi$⁠, where $I$ is $(n−nc)×(n−nc)$ identity matrix. It follows that $δe=Bdiδedi$⁠, where $Bdi=[I−(Cedd−1Cedi)T]T$ is velocity-transformation matrix. In case of linear constraints, equation $δe=Bdiδedi$ leads to $e=Bdiedi$⁠, as demonstrated in the appendix of the paper. The partitioning of the coordinates as independent and dependent is accomplished using a Gaussian forward elimination procedure applied to the nonsquare Jacobian matrix of the constraint equations. The linearly independent columns after the Gaussian elimination define the dependent coordinates. This procedure of identifying the independent coordinates has been used in the MBS literature and there are no numerical problems associated with such a procedure as long as the kinematic constraint equations are linearly independent.

A procedure based on Lagrange–D'Alembert principle for treatment of boundary and motion constraints requires satisfying the algebraic constraint equations at the position, velocity, and acceleration levels during computer simulations. For systems subjected to large displacements, spinning motion, and nonlinear kinematic constraints, thermal loads can be applied during system motion with thermal and material coefficients that are temperature-dependent.

In case of large strains, use of assumptions of strain additive decomposition often adopted for small-deformation thermal analysis is not justified [55,56]. In the strain additive decomposition, based on super-position principle, displacement gradients (not position gradients) are used to define thermal expansion. This simplified approach, which has been used with conventional FE displacement fields, does not allow for properly capturing complex reference-configuration geometries. In this conventional thermal-analysis approach, normal strains are approximated using gradients of displacement vector $u=[u1u2u3]T$ as $εii=∂ui/∂xi, i=1,2,3$⁠. Because such linearization of Green–Lagrange strains is not applicable to large-displacement problems and neglects geometric nonlinearities, in this section, assumptions underlying strain additive decomposition are highlighted before discussing a gradient multiplicative decomposition. ANCF finite elements have been also used in thermo-elasticity analysis in different applications [57–61]. However, these investigations are not based on the multiplicative decomposition adopted in this study [62].

This equation shows the strain additive decomposition, based on principle of super-position, and consequently, is not applicable to nonlinear large-strain problems.

where in this case, the matrix of position-gradient vectors $J$ used to define Green–Lagrange strain tensor is $J=∂r/∂XΘ$⁠, $Jo+Θ=JΘXJo$⁠, and $dv=|JJo+Θ|dV=|J||Jo+Θ|dV$⁠. Effect of constant temperature can be accounted for at a preprocessing stage to reduce number of mathematical operations during numerical solution of system equations.

The four-configuration approach discussed in this section allows for systematically accounting for time- and space-dependent temperatures. The position-gradient matrix $J=∂r/∂X=JrΘJΘX$ can be iteratively updated during computer simulations depending on changes in the temperature value. The method also allows fluctuations in temperatures and for treatment of nonlinear constraints on ANCF gradients using Lagrange–D'Alembert principle.

Kinematic equations developed in this section based on the four-configuration multiplicative decomposition demonstrate unique geometric interpretation of position gradients as tangent to coordinate lines (parameters) [40,63–65]. Displacement gradients do not have the same geometric meaning [23]. Matrix of displacement-gradient vectors can be singular or even null matrix, while matrix of position-gradient vectors is always nonsingular.

Three configurations are used in existing multiplicative-decomposition thermal-stress-analysis approaches: (1) Initial stress-free configuration at uniform reference temperature; (2) Deformed current configuration characterized by nonuniform stress and temperature fields; and (3) Intermediated material configuration determined by release of isothermal stresses of the current configuration to allow thermal stress-free deformation to occur. The significance of intermediate material configuration was emphasized in the literature [66–68]. Using this approach, the total deformation gradient is decomposed into product of purely elastic and thermal parts. To characterize material behavior in thermo-elasticity context, the constitutive equations are developed using this decomposition. Multiplicative decomposition of deformation-gradient matrix was introduced in thermo-elasticity literatures several decades ago [69–72]. Several researchers used multiplicative decomposition of deformation gradients for modeling thermo-elastic continua [20,22,73–79]. Recently, the multiplicative decomposition of deformation gradients was used with finite–discrete element method for thermal fracturing analysis of rocks [80]. Nonetheless, when FE method is used in thermo-elastic large-displacement analysis, an incremental-rotation approach is often used when structural beam, plate, and shell FE meshes are used. This approach is fundamentally different from the nonincremental approach used in this investigation. Furthermore, in a previous investigation, the stress-free thermally loaded configuration is defined by $JΘ=(1+∫ΘoΘα(Θ)dΘ)I$⁠, where $α(Θ)$ is temperature-dependent thermal-expansion coefficient, and $Θo$ is the initial temperature [21]. However, in the ANCF four-configuration multiplicative-decomposition approach, the stress-free thermally loaded configuration is defined by the matrix $JΘX$ when the temperature remains constant or varies with time, demonstrating generality of the ANCF multiplicative-decomposition approach in accounting for the reference-configuration geometry.

The thermal-analysis approach described in this section is applicable to both unconstrained and constrained thermal expansions. The difference between the two cases is boundary conditions and motion constraints used in determining the displacements as the result of thermal-load application. Using the multiplicative decomposition with ANCF finite elements, which impose no restrictions on the amount of FE deformation or rotation, general nonincremental thermal-analysis formulations and computational algorithms can be developed to accurately account for the reference-configuration geometry.

This strain equation ensures that the reference-configuration geometry and stress-free thermal expansion do not enter into the elastic-stress formulation. Nonsingular thermal-gradient matrix $JΘX$ can be defined as $JΘX=αΘΔΘ+I$⁠, where $αΘ$ is diagonal matrix with thermal-expansion coefficients $αΘk, k=1,2,3$⁠, as its diagonal elements. The approach used in this study allows using temperature- and space-dependent thermal-expansion and material coefficients, that is, one can write $αΘk=αΘk(Θ,x),k=1,2,3$⁠. Figure 6 compares magnitude of two gradient vectors $rx1$ and $rx2$ as function of spatial coordinate $x1$ at $x2=0$ in cases of unconstrained and constrained thermal expansions. The results show that, in case of constrained expansion, norm of $rx1$ is quadratic due to boundary restrictions, and minimum gradient stretch occurs at $ξ=0.5$⁠.

Matrix $Jo+Θ=JΘXJo$⁠, which accounts for combined effects of reference-configuration geometry and thermal volumetric change, can be formulated in straightforward manner using ANCF finite elements. To this end, the equation $ro+Θ(x)=S(x)eo+Θ$ is used, where the vector of element coordinates $eo+Θ$ has position-gradient vectors that properly account for the reference-configuration geometry and thermal-load effect by using stretch coefficients $χgk=|(rXk)o+Θ|=dok(1+αΘkΔΘ), k=1,2,3$⁠. These stretch coefficients ensure that gradient-vector orientations are not affected by the thermal expansion in absence of boundary and motion constraints. In case of constant temperature, the matrix $Jo+Θ=JΘXJo$ can be evaluated at a preprocessing stage since such a matrix remains constant. In case of time-varying temperature, $ΔΘ=ΔΘ(x,t)$⁠; one can write $J=JrΘ=∂r/∂XΘ$ and $dr=JJo+Θdx$⁠, where $JΘX=∂XΘ/∂X$ and $JΘX=αΘΔΘ+I$⁠. By providing time-dependent temperature profile, properly updating matrix $JΘX$ at integration points, and following a procedure similar to the one previously discussed; one can show that thermal expansion does not contribute to elastic stresses.

The Green–Lagrange strain tensor defined in this section can be used with the constitutive material model to define virtual work of elastic or viscous forces. Inertia forces can be defined using virtual work or kinetic energy. The details of formulating ANCF elastic and inertia forces and of MBS equations of motion are presented in the literature [40].

In this section, the procedure proposed in this paper for the constrained thermal expansion is first applied to case of unconstrained thermal expansion for verification purpose and to obtain a solution that serves as a reference to measure the effect of the thermal load in case of constrained thermal expansion. To this end, two examples are considered; the first is structural problem of cantilever beam, while the second is simple MBS slider crank mechanism with flexible connecting rod. To focus on demonstrating the main concepts introduced in this study, case of constant and uniform temperature and constant thermal-expansion coefficients is considered. In all simulations performed, the Poisson ratio is assumed to be 0.3.

The first example used to obtain unconstrained thermal expansion results is a beam with fixed-free boundary conditions, as shown in Fig. 7. The beam is assumed made of soft material with modulus of elasticity $E=1.2×108 Pa$⁠, mass density $ρ=1500 kg/m3$⁠, and coefficient of thermal expansion $α=1×10−4 (1/°C)$⁠. The beam length $l=1 m$⁠, and cross section area $A=0.3×0.3 m2$⁠. The beam is modeled using 10 two-dimensional ANCF fully parameterized planar beam elements. Elastic forces are formulated using the general continuum-mechanics approach using plane stress assumptions and neglecting damping effect. The beam is subjected to two types of loads; namely, constant axial load of 1000 N applied at the beam-free end, and uniform thermal load. Three simulations are performed with thermal loads $ΔΘ=0 °C, 100 °C,and 200 °C$⁠. Results due to application of these three different thermal loads are recorded. Figure 8 shows axial deformation at the beam-free end. The results shown demonstrate that, in case of unconstrained thermal expansion, thermal-load effect on the deformation is not significant.

The second example shown in Fig. 9 and used to study effect of unconstrained thermal expansion is a slider crank mechanism with a flexible connecting rod subjected to a uniform thermal load. The point, connecting the slider block with the connecting rod, is permitted to extend when the thermal load is applied. This is a typical unconstrained thermal expansion problem. The crankshaft, assumed rigid, rotates with a constant angular velocity $ω=π rad/s$⁠. At the initial configuration, the connecting rod and the crankshaft are assumed to be in horizontal position. A linear elastic model is used for the connecting rod which is assumed to be made of soft material with modulus of elasticity $E=2×108 Pa$⁠, mass density $ρ=7200 kg/m3$⁠, and coefficient of thermal expansion $α=80×10−6 (1/°C)$⁠. The rod undeformed length is $l=3.048×10−1m$ and cross-sectional area is $A=(5.5×10−3)×(5.76×10−3) m2$⁠. Temperature changes are assumed $ΔΘ=50 °C, 100 °C, and 200 °C$⁠. General continuum-mechanics approach is used to formulate the elastic forces using plane-stress assumptions. Strain split method (SSM) is used to alleviate locking problems [81]. Gravity effect is not considered in this example to focus on temperature effect, and total simulation time is assumed 2 s. The connecting rod is modeled using six ANCF fully parameterized planar beam elements, and the slider block is assumed massless. To verify the model results, the connecting-rod midpoint transverse deformation is compared with results obtained using the floating frame of reference (FFR) formulation. FFR connecting-rod model consists of four planar elements and FFR simply supported reference conditions are used. The FFR deformation is obtained using first four mode shapes. Figure 10 shows that thermal-load effect is negligible. The results presented also demonstrate that temperature effect is not significant in case of unconstrained thermal expansion when MBS joint constraints do not restrict the thermal expansion.

To examine case of constrained thermal expansion, a beam with fixed-fixed boundary conditions, shown in Fig. 11, is used as an example. The beam is assumed made of soft material with modulus of elasticity $E=2×108 Pa$⁠, mass density $ρ=7200 kg/m3$⁠, and coefficient of thermal expansion $α=80×10−6 (1/ °C)$⁠. The beam length is $l=1 m$⁠, and cross section area is $A=0.05×0.05 m2$⁠. The beam is modeled using 10 two-dimensional ANCF fully parameterized planar beam elements. Elastic forces are formulated using general continuum mechanics approach using plane-stress assumptions and neglecting effect of damping and gravity. SSM approach is used to alleviate locking problem. The beam is subjected to two types of loads; namely, constant vertical load of 20 N applied at beam midpoint, and uniform thermal load. Three simulation scenarios that correspond to thermal loads $ΔΘ=0 °C, 50 °C, and 100 °C$ are considered. To verify ANCF results, an ansys workbench model was developed to solve the same problem using transient structural solver where thermal and mechanical loads are applied simultaneously. To capture accurately thermal-load effect on beam cross section, the beam is modeled in ANSYS using three-dimensional solid elements (SOLID 186) as recommended in the literature [82]. The results, presented in Figs. 12–14 for midpoint vertical deformation for the two models, show that as the thermal load increases, deformation amplitude increases. Deformation increase is due to increase in the restoring forces required to prevent beam axial expansion. Consequently, these forces increase as thermal load increases [83–88]. Furthermore, the results show that ANCF fully parameterized two-dimensional elements can model large deformations and changes in cross section dimensions, which cannot be accurately modeled using conventional beam elements. Figures 15–17 show normal strain components for the two models at different thermal loads. It is shown that the initial axial compressive stresses lead to initial jump in axial normal stresses. Furthermore, transverse strain initially increases as result of the stretching due to Poisson effect.

In the analysis of the thermal-load effect, distinction is made between two scenarios. In the first scenario, stress-free expansion is allowed by assuming that space is adjusted to accommodate for the stress-free thermal expansion [89–95]. In this case of unconstrained thermal expansion, there is no initial stress due to thermal loads because the expansion in some directions is not prevented by joints or boundary conditions. The second problem is the constrained expansion considered in this study, where boundary constraints are enforced during thermal expansion. This case leads to initial stresses as previously mentioned.

In order to examine the effect of constrained thermal expansion on the solution of MBS problems, the same slider crank mechanism example used in Sec. 7.1 is considered. In case of constrained thermal expansion, the connecting-rod endpoints are not allowed to move when the thermal load is applied before the motion starts. Figure 18 shows the deviation in the slider position at different thermal loads from its position in case of no thermal load. Figures 19–21 show the connecting rod midpoint transverse deformation for different temperature values, while Figs. 22–25 shows the effect of thermal loads on normal strain components. The results presented in these figures demonstrate that the thermal-load effect on the dynamic response increases as the temperature increases. These results also show significant effect of the initial deformation in case of constrained thermal expansion.

High frequencies in the displacement solutions presented in this section are attributed to the thermal expansion considered at preprocessing stage in this study. When assembling the mechanism based on the original dimensions before thermal expansion, the thermally expanded connecting rod has initial deformation that produces initial stresses and high frequencies in the solution. Future investigations will consider gradual increase in the temperature.

The approach presented in this paper and used to solve the examples and obtain results presented in this section is based on using the multiplicative decomposition of the matrix of position-gradient vectors and four different configurations previously described to account for the thermal-load effect. Using the position-gradient vectors as nodal coordinates for ANCF finite elements allows for systematic treatment of constraints imposed on the boundary of the deformable bodies. It also allows for systematic implementation of the proposed multiplicative decomposition approach in general large displacement MBS algorithms. When position gradients are used as nodal coordinates, different joints can be defined. For example, in the planar analysis, the use of the conventional finite elements that employ two translation coordinates and one rotation coordinate for a node imposes restrictions on the joint type and deformation modes that can be used. Fixing these three coordinates defines conventional clamped joint, which eliminates translation and rotation at nodal point. However, such a conventional clamped joint does not impose any constraints on the position gradients and does not describe properly the geometry in case of thermal expansion. ANCF finite elements, on the other hand, allow imposing different clamped joints, including partially clamped joint, which eliminates translations and rotations at a point without imposing any constraints on the position gradients. For such partially clamped joint, thermal expansion is the same in all direction. Another joint type is fully clamped joint, which imposes constraints on translation, rotation, and position gradients. Figure 26 shows the effect of thermal load on the geometry when both types of joints are used. In case of fully clamped joint at first beam end, which is fixed, cross section does not stretch while the free end stretches, producing tapered geometry. In case of partially clamped joint, on the hand, cross section stretches with the same amount at both ends of the beam since no constraints are imposed on position gradients.

In the thermal analysis, two different situations can be encountered. In the first case, continuum boundary is not constrained and all material points experience same displacement-gradient change as result of thermal-load application. This case, referred to in this paper as unconstrained thermal expansion, produces uniform motion within the continuum. In the second case, referred to as constrained thermal expansion, motion of boundary points is restricted, and therefore, such points are not allowed to move freely when the thermal load is applied. In this second case, displacements of material points are not the same, regardless of whether temperature is constant throughout the homogeneous continuum. In this study, a comparison between two scenarios is made using ANCF finite elements. ANCF position-gradient vectors, which have a unique geometric meaning as tangent to coordinate lines, allow for systematically describing two different unconstrained and constrained thermal expansions using multiplicative decomposition of matrix of position-gradient vectors. Numerical results are presented to demonstrate using the approach proposed in this study and compare between the two different expansion patterns. The approach proposed in this study can be used in future investigations for different reference geometries and practical applications [96,97].

• National Science Foundation (Funder ID: 10.13039/100000001).

All the data used in this investigation are presented in the paper.

This equation shows that, in case of one element, the position-gradient vectors can be written in terms of independent coordinates, which can freely vary for different loading conditions.

The procedure described in this appendix can be generalized to the case of FE mesh that consists of arbitrary number of elements. In this case, the vector $e$ represents vector of mesh nodal coordinates which can have dimension $n$⁠. Matrix $Bdi$ in this case is an $n×nc$ matrix, where $nc$ is number of constraint conditions imposed on the coordinates $e$⁠. In this case, the vector of independent coordinates $edi$ has dimension $ndi=n−nc$⁠.

If constraints are imposed on position coordinates and no constraints are imposed on position gradients, one can identify independent gradients, which define the expansion. For example, in the case of a single element mesh, if two FE endpoints are not allowed to move, vector $edi$ can be written as $edi=[(rx11)diT(rx21)diT(rx12)diT(rx22)diT]T$⁠. These independent gradient vectors can be written in terms of product of thermal-expansion coefficient and temperature change, $αΘΔΘ$⁠. The original gradients can be obtained using the velocity transformation as $ed=Bdiedi$⁠.