-J/o 9s 

ANALYSIS OF SHELL-TYPE STRUCTURES SUBJECTED TO 

TIME-DEPENDENT MECHANICAL AND THERMAL LOADING 

A SEMI-ANNUAL STATUS REPORT 

SUBMIlTED TO 

NASA-LEWIS RESEARCH CENTER 

CLEVELAND, OHIO 

6. J. Simitses 

School of Aerospace Engineering 

Georgia Institute of Technology 

Atlanta, Georgia 

May 1989 

(NASA GRANT: NAG 3-534) 

INTRODUCTION 

The objective of the present research is to develop a general mathematical model and sol uti on methodol ogi es for analyzing structural response of 

thin, metal 1 i c she1 1 -type structures under 1 arge transient , cycl i c or static 

thermomechanical loads. Among the system responses, which are associated with 

these load conditions, are thermal buckling, creep buckling and racheting. 

Thus, geometric as well as material-type nonlinearities (of high order) can be 

anticipated and must be considered in the development of the mathematical 

model. Furthermore, this must also be accommodated in the solution procedures. 

SUmARY OF PROGRESS 

The progress to date has been elaborated upon in an interim scientific 

report submitted to the sponsor during the summer of 1986, and in a series of 

semiannual progress reports. The most recent of these is dated April, 1988. 

A complete true ab-initio rate theory of kinematics and kinetics for . 

continuum and curved thin structures, without any restriction on the magnitude 

of the strains or the deformation, was formulated. The time dependence and 

large strain behavior are incorporated through the introduction of the time 

rates of the metric and curvature in two coordinate systems; a fixed (spatial) 

one and a convected (material) coordinate system. The relations between the 

time derivative and the covariant derivatives (gradients) have been developed 

for curved space and motion, so that the velocity components supply the 

{tiASA-CR- 184989) ANA LYSI S OF SHELL- TYPE N89-24669 

STRUCTURES SUBJECTED TO TIHE-DEPENDEWT 

MECHANICAL AND THERFiAL LOADING Semiannual 

Status Report (Georqia Inst. of Tech.) 9 p Unclas 

CSCL 20K G3/39 0210906 
connection between the equations of motion and the time rate of change of the 

metric and curvature tensors. 

The metric tensor (time rate of change) in the convected material coordinate system is linearly decomposed into elastic and plastic parts. In this 

formulation, a yield function is assumed, which is dependent on the rate of 

change of stress, metric, temperature, and a set of internal variables. 

Moreover, a hypoelastic law was chosen to describe the thermoelastic part of 

the deformation. 

A time and temperature dependent viscoplastic model was formulated in 

this convected material system to account for finite strains and rotations. 

The history and temperature dependence were incorporated through the introduction of internal variables. The choice of these variables, as well as their 

evolution, was motivated by phenomenological thermodynamic considerations. 

The noni sothermal el asti c-viscopl astic deformation process was described 

completely by "thermodynamic state'' equations. Most investigators (in the 

area of viscoplasticity) employ plastic strains as state variables. Our study 

shows that, in general, use of plastic strains as state variables may lead to 

inconsistencies with regard to thermodynamic considerations. Furthermore, the 

approach and formulation employed by all previous investigators lead to the 

condition that all plastic work is completely dissipated. This, however, is 

in contradiction with experimental evidence, from which it emerges that part 

of the plastic work is used for producing residual stresses in the lattice, 

which, when phenomenologically considered, causes hardening. Both limitations 

are not present in our formulation, because of the inclusion of the "thermodynami c state" equations. 

The obtained complete rate field equations consist of the principles of 

the rate of the virtual power and the rate of conservation of energy, of the 

constitutive relations, and of boundary and initial conditions. These formulations provide a sound basis for the formulation of the adopted finite 

element solution procedures. 

The derived shell theory, in the least restricted form, before any 

The two-dimensional, impulse-integral form of the equations of 

motion and the Second Law of Thermodynamics (Clausius-Duhem inequality) for a shell follow naturally and exactly from their three-disimplifying assumptions are imposed, has the following desirable features: 

mensional counterparts. 

Unique and concrete definitions of shell variables such as stress 

resultants and couples, rate of deformation, spin and entropy 

resultants can be obtained in terms of weighted integrals of the . 

three-dimensional quantities through the thickness. 

There are no series expansions in the thickness direction. 

There is no need for making use of the Kirchhoff Hypotheses in the 

kinematics. 

- All approximations can be postponed until the descretization process 

of the integral forms of the First Law of Thermodynamics 

A by-product of the descent from three-dimensional theory is that 

the two-dimensional temperature field (that emerges) is not a 

through-the-thi ckness average, but a true point by point 

distribution. This is contrary to what one finds in the literature 

concerning thermal stresses in the shell. 
To develop geometrical ly nonl i near , doubly curved finite shel 1 elements 

the basic equations of nonlinear shell theories have to be transferred into 

the finite element model. As these equations in general are written in tensor 

notation, their implementation into the finite element matrix formulation 

requires considerable effort. 

The nonlinear element matrices are directly derived from the incrementally formulated nonlinear shell equations, by using a tensor-oriented procedure. 

The cl assi cal thin shel 1 theory based on the Ki rchoff-Love hypotheses (Formulation D in Appendix A) was employed for this purpose. For this formulation, 

we are using the "natural" degrees of freedom per mid-surface shell node: 

three incremental velocities and the rates of rotations about the material 

coordinates in a mixed form. 

A description of the developed element and related finite element code 

are given in Appendix B. This exposition provides 'information concerning the 

formulation, the finite element and how it is employed in the solution of 

shel 1-1 i ke configurations. Moreover a compl ete description incl udi ng program 

flow chart, listing, input instructions to the user and explanation of output 

are also included in Appendix B. 

The quasi-linear nature of the principle of the rate of virtual power 

suggests the adoption of an incremental approach to numerical integration with 

respect to time. The availability of the field formulation provides assurance 

of the completeness of the incremental equations and allows the use of any 

convenient procedure for spatial integration over the domain V. In the 

present instance, the choice has been made in favor of a simple first order 

expansion in time for the construction of incremental solutions from the 

results of finite element spatial integration of the governing equations. 

The procedure employed permits the rates of the field formulation to be 

interpreted as increments in the numerical solution. This is particularly 

convenient for the construction of incremental boundary condition histories. 

Even under the condition of static external loads and slowly growing 

creep effects, the presence of snap-through buck1 ing makes the inertial 

effects significant. In dynamic analyses, the applied body forces include 

inertial forces. Assuming that the mass of the body considered is preserved, 

the mass matrix can be evaluated prior to the time integration using the 

ini ti a1 configuration. 

Finite element solution of any boundary-value problem involves the 

solution of the equilibrium equations (global) together with the constitutive 

equations (local). Both sets of equations are solved simultaneously in a step 

by step manner. The incremental form of the global and local equations can be 

achieved by taking the integration over the incremental time step t=t. -t.. 

The rectangular rule has been applied to execute the resulting time idkigrdtion. 
Clearly, the numerical solution involves iteration. A simp1 ified version 

of the Ri ks-Wempner constant-arc-length method has been uti 1 ized. This iteration procedure which is a generalization of the displacement control method 

also allows to trace the nonlinear response beyond bifurcation points. In 

contrast to the conventional Newton-Raphson techniques, the iteration of the 

method takes place in the velocity and load rate space. The load step of the 

first solution in each increment is limited by controlling the length ds of 

the tangent. Either the length is kept constant in each step or it is adapted 

to the characteristics of the solution. In each step the triangular-size 

stiffness matrix has to be checked for negative diagonal terms, indicating 

that a critical point is reached. 

One of the most challenging aspects of finite strain formulations is to 

locate analytical solution with which to compare the proposed formulation. 

Typically, as a first problem, a large strain uniaxial test case was analyzed. 

The case considered examines the rate-dependent plastic response of a bar to a 

deformation history that includes segments of loading, unloading, and reloading, each occurring at varying strain and temperature rates. Moreover, it was 

shown that the proposed formulation generates no strain energy under a pure 

rigid body rotation. These are surely important demonstrations but they only 

represent a partial test because the principal stretch directions remain 

constant. Finally, a problem which was discussed by Nagtegaal and de Jong, 

and others too, as a problem which demonstrates limitations of the constitutive models in many strain formulation,, is the Couette flow problem. This 

problem is solved as a third example. The results of these test problems show 

that: 

- The formulation can accommodate very large strains and rotations. 

- The formulation does not display the oscillatory behavior in the 

stresses of the Couette flow problem. 

- The model incorporates the simplification associated with 

rate-insensitive elastic response without losing the ability to 

model rate temperature dependent yield strength and plasticity. 

The problem of buckling of shallow arches under transient 

thermomechanical load was investigated next. The analysis was performed with 

the aid of 24 paralinear isoparametric elements. The paralinear isoparametric 

element is such that the thickness is small compared to other dimensions. The 

characteristics of the el ement are defined by the geometry and i nterpol ati on 

functions, which are linear in the thickness direction and parabolic in the 

longitudinal direction. Consequently, the element allows for shear strain 

energy since normal s to a mi d-surface before deformation remain straight , but 

not necessarily normal to the midsurface after deformation. 

The developed solution scheme is capable of predicting response which 

includes pre- and post-buckling with thermal creep and plastic effects. The 

solution procedure was demonstrated through several examples which include 

both creep and snap-through behavior. 

The last set of problems which are under investigation consists of creep 

or thermal buckling, with plastic effects, of shells of revolution. 
In addition, following a more traditional approach, a method was developed for bounding the response (solution) of bars and beams of (linear) 

visoelastic material behavior, based on nonlinear kinematic relations. 

In connection with the progress to date, two papers were published by the 

AIAA Journal in 1986 and 1987. Moreover, a paper entitled, 'Non-Isothermal 

El astovi scoplasti c Analysis of P1 anar Curved Beams" was presented at the 3rd 

Symposium on Nonlinear Constitutive Relations for High-Temperature Appl ications, held at the University of Akron, on June 11-13, 1986. A descriptive 

abstract of this paper was published in the meeting proceedings and the full 

paper appeared in a special publication (NASA CP 10010). Copies of the above 

have been sent to the sponsor. 

In addition, the two papers presented at the 28th AIAA/ASME/ASCE/AHS SDM 

Conference and published in the proceedings of this conference, have been 

accepted for publication; both by the AIAA Journal. These two papers deal 

with applications to snap-through and creep buckling of bars and arches. Most 

of this work was also presented at the NASA-Lewis Conference on Structural 

Integrity and Durability of Reusable Space Propulsion Systems on May 1987 in 

Cleveland, Copies of these papers will be forwarded to the Sponsor, as soon as 

they appear in print. 

In connection with the more traditional approach a paper accepted for 

presentation and for pub1 i cati on in the .Proceedings of the speci a1 Symposi um 

on Constitutive Equations at the ASME Winter Annual Meeting, Chicago, IL., 

November 28 - December 2, 1988. The title of the paper is "Creep Analysis of 

Beams and Arches Based on a Hereditary Visco- Elastic- Plastic Constitutive 

Law". Copies of this paper have been forwarded to the sponsor in December, 

1988 1. 

Moreover, one paper was published in the the AIAA/ASME/ .../ AHS 30th SDM 

Conference. A copy of this paper is attached, herewith. 

Finally, a finite element has been developed and it is currently being 

tested for a less restricted shell formulation (Formulation C; see Appendix 

A). A description of this newly developed element and the related finite 

element code will be submitted to the sponsor as soon as the testing is 

completed and reliable results have been obtained. 
FUTURE TASKS 

The main thrust of the additional tasks is to develop a finite element 

and select a code, which will be made available to all users and which will be 

based on the most general (but practical) nonlinear shell formulation possible 

and nonlinear constitutive relations to predict the response of shell-like 

structures, when subjected to time-dependent thermomechanical 1 oads with large 

excursions. 

This should be completed in the next six months and the complete package 

will be delivered to the sponsor. 
APPENDIX A 

The various shell theory approximations (formulations) are based on the 

use of certai n simp1 i fyi ng assumpti ons regardi ng the geometry and ki nemati cs 

of the shell configuration. 

These are: 

Assumption I: The material points which are on the normal to the reference surface before deformation will be on the same normal after deformation. 

Assumption 11: The shell is sufficiently thin so that we can assume 

linear dependence of the position of any material point (in the deformed 

state) to the normal (to the reference surface) coordinate (in the 

deformed state). The linear dependence can easily be changed to parabolic, cubic, or any desired degree of approximation. 

Assumption 111: The rate of change of the velocity gradients with 

respect to in-plan coordinates on the two boundary shell surfaces is 

negligibly small. -- Assumption IV: 

reference surface is negligibly small. 

The rate of change of the distance of a material from the 

On the basis of the above four simplifying assumptions, several formulations result, for the analysis of thin shells. 

Formulation A: This formulation makes use of Assumption I, only. 

Formulation B: This formulation employs Assumptions I and 11. 

Formulation C: This formulation employs Assumption I, I1 and 111. 

Formulation 0: This is the classical thin shell theory based on the 

Kirchhoff-Love hypotheses of Assumptions I, 11, 111, IV, 

as applied to large deformation theory. 

These formulations are arranged in such a manner that we move from the least 

restrictive (A) to the most restrictive (D). 

In addition to this a fifth formulation (E) can easily be devised and 

this formulation in terms of order of restriction is similar to Formulation A. 

Formulation E makes use of Assumption I1 only. 

APPENDIX B FINITE ELEMENT AND RELATED CODE FOR FORMULATION D 

In this Appendix, a description of the developed finite element and the 

related finite element code is presented. First, the essentials of the 

element are described and then the complete solution procedure is presented 

with sufficient detail. This includes a flow chart, a line by line listing of 

the computer program, input data information, and explanation of the output. 
Btl. THE SHELL ELEMENT 

A brief description highlighting the essential features of 

the shell element development and the related code used in this 

work is given here. 

In order to derive discrete algorithm based on the finite 

element displacement method we approximate the velocity field by 

index-oriented notation, which allows the separate representation 

of the shape functions (the specific expression depends on the 

decided upon degress - of - freedom, Lagrangian, Hermitian, etc.) 

for the tangential velocities uk and for the normal velocityJi . 

Upper indices imply the columns, lower indices the rows of a 

matrix expression, and the summation is carried out spanning the 

number of degrees-of-freedom. VM and VM represent, therefore, the 

vector of nodal velocities by the row and by the column respectively. We get the shape functions for the partial derivatives 

of the velocity shape functionsvz jQT : 

The main idea of this formulation is the development of 

shape functions for further mechanical and thermal variables by 

the application of well-known tensor procedure on the basic shape 

functions (3) and (4). Taking, for example, the covariant 

derivative 

’ 

and inserting. (3) and (l), we can define the shape function of 

the covariant- derivative: 

vd; P 

In the same way we receive the shape function of the 

internal variables, for example the rate of deformation and the 

spin tensors: 
The same procedure is now applied to the shift tensor 

which is responsible for the llexacttt distributions over the 

thickness. 

Finally the shape functions of the internal forces and 

temperature variables can be derived from the shape functions of . the rate of deformation and spin tensor via the constitutive 

relations; for example: 

All these expression are now substituted into the rate of 

the first law of thermodynamics to obtain the element ltstiffnesstl 

equations. The developed element matrices are implemented into 

the global relation of the complete shell structure by standard 

assemblage process considering incidence and boundary conditions. 

The present curvilinear formulation of the element enables 

the precise description of geometry, external loads and temperatures and the fulfillment of the convergence criteria, while the 

rigid body motion condition can only be satisfied in an approximate manner. The tensor -oriented formulation renders the . 

optional use of various shape functions for the tangential 

velocities v and the normal velocity v3. 

'The shell element which have been used up to today is based 

on the bicubic Hermite polynomial with 4x12 generalized velocities and 4 temperatures. Numerical integration spanning the 

element domain was applied (16 points of integration), whereby 

area and boundary integrals were replaced by double integration 

with respect to the curvilinear B4 - coordinates. 

dA= d0' d.b2 (11) 

dS= rip dBd d-0') (12)