cn345-7949/88 $3.00 + 0.00 0 1988 Pergamon Press plc
Compurers & Sfrucrures Vol. 29. No. 2, PP. X39-315, 1988 Printed in Great Britain.
A GLOBAL-LOCAL FINITE ELEMENT METHOD SUITABLE FOR PARALLEL COMPUTATIONS C. T. SUN and K. M. MAO School of Aeronautics & Astronautics, Purdue University, West Lafayette, IN 47907, U.S.A. (Received 20 July 1987)
Abstract-A
global-local finite element analysis procedure is developed based on the fast convergent nature of FEM in displacement. A special scheme is used to utilize the global displacement solution as boundary conditions for local regions of interest. In the local region, a refined mesh is used for further stress analysis. This global-local procedure can be easily programmed in parallel on MIMD multiprocessor computers for significant time savings. A Sequent Balance 21000 system is used for demonstrating the parallel programming.
INTRODUCTION Finite element analysis of a structure with steep stress gradients requires very fine meshes. Structures with many cutouts may need a large number of degrees of freedom and, consequently, a great deal of computing time. For global responses of composite laminates, effective stiffnesses can yield excellent results. However, fine finite element meshes are needed in order to obtain interlaminar stresses. It is well known that displacement-formulated finite elements converge in displacement much faster than in stress. In other words, coarser meshes can achieve a good convergence in displacement but not in stress. In the present study, a global-local method is developed based on the fast convergent nature of FEM in displacement. A special scheme is used to utilize the global displacement solution as boundary conditions for a local region of interest. In the local region, a refined mesh is used for further analysis to improve the accuracy of the finite element solution in stress. The proposed procedures are particularly suitable for computation on multiprocessor computers. In order to obtain the stress fields in numerous subregions of a structure, many local regions of refined meshes must be considered. These computations can be performed concurrently on different processors in the computer, thereby saving much computing time. Recently, interest in parallel computation for structural mechanics increased dramatically because of availability of a number of multiprocessor computers with parallel architectures [l]. So far, most of the effort has been placed in developing parallel algorithms for solving systems of equations and differential equations. A few authors have begun to bring mechanics into the parallel programs. One notable example is the Finite Element Machine developed at NASA Langley Research Center [2,3]. This multiprocessing system has interconnected microprocessors, each of which is assigned to a node or a
set of nodes.
Formulation of the finite elements and solutions for the resulting global equilibrium equations on the Finite Element Machine was discussed in [3]. In a recent paper, Law [4] presented an iterative procedure for parallel computation with finite elements. Each finite element is mapped onto a processor of a multiprocessor system where the stiffness matrix for the element is generated. The global equilibrium equations are satisfied by minimizing the unbalanced nodal forces using the conjugate gradient method. The main difference between the present method and those of [3,4] is that in the present method the finite element model of a structure is broken up into local regions which are solved concurrently on separate processors in a multiprocessor machine. Each individual process is no different from the process on a sequential machine. Of course, parallel computation can be invoked to solve the systems of equations for each local region, thus further increasing the time saving. Two illustrative examples are considered in this paper: a cantilever beam of composite laminate, and a center-cracked isotropic panel under uniform tension. Application of parallel computation is demonstrated using the Balance 21000 multiprocessor computer by Sequent.
THE GLOBALLOCALPROCEDURE The present global-local method is based on the fact that the displacement field calculated using displacement-formulated finite elements converges much more rapidly than the stress field. In other words, a relatively coarse mesh may yield a reasonably accurate solution in displacement. This displacement solution can then be used as a displacement boundary condition for a local region in which the stress field is of interest. The local region is analyzed with a refined mesh. 309
C. T. SUNand K. M.
310
MAO
Fig. 1. Global and local finite element meshes. Fig. 2. Geometry and loading of the composite beam. this procedure consists of two steps of computation. The first step is a global analysis of the displacement field in the whole structure using a coarse mesh. The second step is the stress analysis of local regions with refined meshes using the global displacements as boundary conditions for the local regions. Avoiding errors in the global solution, local analysis requires regions larger than the region of interest. As illustrated in Fig. 1, if the stress field in the shaded area is to be calculated, the local region should be taken at least one element away (shown by the region ABCD). After the mesh in the local region is refined, additional nodal degrees of freedom are introduced. The global solution only yields displacements at the nodes of the original coarse mesh; the displacement boundary conditions for the refined mesh along the edges ABCD (see Fig. 1) must be derived from the global solution. One method is to prescribe the displacements at the refined nodes employing the shape function of the original finite element used in the global analysis. Thus, between two adjacent nodes of the original global mesh, displacement along the boundary is completely determined following the global analysis. If traction is prescribed on part of the boundary of the local region, then the actual force condition should be incorporated in the local analysis. One important feature in the procedure described above is that the local analysis for each region can be performed independently. This leads naturally to the use of multipr~ssor computers. In theory, if there are enough processors, computing time should not increase as the number of local regions expands. In this paper, two numerical examples are given. One involves stress analysis of a thick laminate, and the other deals with calculating the stress field near a crack tip in a center-cracked panel. The global-local method is first demonstrated with a conventional sequential computer and then using the Sequent multiprocessor computer. All numerical examples employ eight-node isoparametric plane finite elements. Thus,
ANALYSIS
OF A THICK
LAMINATE
Advanced fiber composites have gained wide application in various kinds of structures. They are usually
used in the form of laminates. The material properties are not uniform through the thickness. Consequently, very fine finite element meshes are needed to analyze interlaminar stresses. The global method provides an efficient way to obtain these localized stress fields. Consider a cantilever beam of a ~aphite/epoxy composite laminate consisting of 40 plies. The laminate is formed by repeating the O/90 pair 20 times and is denoted by [O/90],. The beam is subjected to a uniform transverse load of 1.2 lb/in., as shown in Fig. 2. A state of plane strain parallel to the x-y plane is assumed in the analysis. The basic uni~rectional fiber composite is assumed to have the following material constants: E, = 20 x lo6 psi E2 = E, = 1.43 x lo6 psi Gt2 = G,, = G,, = 0.7 x lo6 psi vu = vu = 0.3 v23= 0.1. The definitions of the elastic moduli posites can be found in [5]. Three solutions are obtained using coarse mesh, and the global-local spectively. Eight-noded isoparametric are used with consistent loads.
for fiber coma fine mesh, a procedure, refinite elements
1 Fine mesh solution-the whole beam is divided uniformly into 40 x 80 rectangular eight-noded isoparamet~c elements (40 elements in the longitudinal direction). Thus, each lamina is modeled by two elements over its thickness. The solution obtained with this mesh is considered converged and is used for comparison with the other two solutions. 2. Coarse mesh solution-a 10 x 10 uniform rectangular mesh is used to model the beam. Thus, each element covers four laminae. The stiffness matrix is generated by integrating the energies from lamina to lamina in each element. 3. Global-local solution-by taking the coarse mesh solution as the global solution, local regions are selected and meshes are refined. A local region contains a 3 x 3 coarse mesh which is further refined into a 12 x 24 mesh for the local analysis. Note that the grid size in the local
Global-local
311
finite element method for parallel computations
A
‘+OsSO MESH
0
10~10
MESH
q GLOBAL-LOCAL
m
Fig. 3. Transverse displacement at the bottom surface of the beam.
analysis is identical to that in the fine mesh solution. Depending of the purpose of analysis, many local regions may be needed to cover the area of interest.
Comparisons of the three solutions are made in Figs 3-5. Figure 3 shows the transverse displacements of the bottom surface of the beam obtained with the fine mesh and the coarse mesh. Both solutions
-Am0
-.mln
A 0 [II
‘+OrSO
MESH MESH GLOBAL-LOCAL
10~10
-.mD-
c:
E *.Sm-
$ z -.w-
normal stress uv along y-axis x =0.1375in.
Method
-.mm-
normal stress o,, along x-axis y = 0.1475 in.
at
at
coincide everywhere, indicating rapid convergence in displacement. Figure 4 shows the transverse normal stress (a,) at y = 0.1475 in. along the x-axis. The results clearly indicate that the global-local solution has converged to the fine mesh solution. Appreciable errors are noted in the coarse mesh solution. In Fig. 5, transverse normal stress is plotted along the y-axis at x = 0.1375 in. from the fixed end. Again, the globalqocal solution agrees very well with the fine mesh solution. The coarse mesh solution exhibits a zigzag behavior, with the average values falling near the converged solution. All computation is performed on a VAX 11/780 computer. The CPU times corresponding to the three methods are listed in Table 1. The CPU time for the global-local method includes 131 set for the global analysis (using the 10 x 10 mesh) and 179 set for the local analysis of one local region. Using a sequential machine, each additional local analysis requires an additional 179 sec. It is evident that the global-local procedure may achieve great savings in computing time if a small number of local regions is considered. Table 1. Computing time on VAX 1l/780 computer for the composite beam
-Am-
Fig. 4. Transverse
Fig. 5. Transverse
Fine mesh (40 x 80) Coarse mesh (10 x 10) Global-local (one local region)
Time (set) 20551 131 310
312
C.
ANALYSIS
OF A CRACKED
T. SUN and K.
M. MAO
PANEL
q-6.0
(lb/in)
(a)
The second evaluative example is a square panel of an isotropic material with a center crack. The panel is subjected to a uniform tension as shown in Fig. 6. Due to symmetry, only a quadrant of the plate is taken for analysis. Similar to the previous example, three different methods are used to analyze the stress distribution near the crack tip. The coarse mesh contains 8 x 5 = 40 rectangular eight-noded isoparametric finite elements as shown in Fig. 7. The fine mesh model is obtained by further dividing each coarse element into a uniform mesh of 5 x 5 = 25 elements. In the global-local analysis of the crack tip stress field in the region JKGH, the local region is selected to be ACFI. For stresses in the region KLFG, the local region is BDEH (see Fig. 7(b)). The local regions have the same mesh size as the fine mesh model. In the analysis of region ACFI, nodal displacements are prescribed along the edges AC, CF, and IA, using the global solution in conjunction with the quadratic interpolation between two global nodes. At other nodes along edge IF (but not including nodes I and F) the exact boundary conditions are imposed. The elastic constants used in the analysis are E = 2.4 x lo6 psi v = 0.2. Stresses near the crack tip along the x-direction at y = 0.025 in. are calculated using the element displacement function. The results are presented in Figs 8-10. The comparison in Figs 8-10 indicates that while the global solution is inaccurate, the global-local solution agrees well with the fine mesh solution up to the crack tip.
(b)
Fig. 7. Global finite element mesh and local regions for the cracked panel.
However, we note slightly larger errors in the global-local solution in the region between x = 0.25-0.65 in., especially with stress component a,. This is due to the use of the crack tip node (node H) as a boundary node for the local region BDEH. Due to the stress singularity nature at the crack tip, the error in this node’s displacement obtained from global analysis may be quite significant. This error may affect the accuracy of stress analysis for the
-m mm-
A 0 [II
FINE MESH CORRSE MESH GLOBAL-LOCAL
pl.ca-
Fig. 6. Geometry and loading of the cracked panel.
Fig. 8. Stress (uJ distribution near the crack tip along x-axis at y = 0.025 in.
Global-local
finite element method for parallel computations
A FINE MESH 0 CORRSE MESH ITI GLOBRL-LOCAL
@.I-
313
A Q
FINE IIESH GLOBAL-LOCAL
l.m-
1962 rr x
rm-
8 8 l.76-
km-
IS-
I-
m
0
Fig. 9. Stress (a,) distribution near the crack tip along x-axis at y = 0.025 in.
region KLFG. Indeed, examining the nodal displacement at node H and the nearby nodes, we noted that the coarse mesh displacement at the crack tip (node H) deviates substantially from the fine mesh solution (see Table 2). At other nodes in the neighborhood of the crack tip, the two solutions agree reasonably well. To avoid the aforementioned problem, the local region ADEI is chosen for analysis of the stress field
.
!Stress (a,) distribution
near the crack tip using an enlarged local region.
in the region JLFH. The new results for Us are presented in Fig. 11. Great improvement in the global-local solution is evident. The CPU times on the VAX 11/780 computer for the three methods are given in Table 3. The time for the global-local analysis includes two local regions. The speedup due to the global-local procedure is about 17 times.
,.a
A FINE MESH 0 CIYRRSE MESH II1 GLOBAL-LOCAL
s.mo
pP.ca-
b.mo
\ ,e.m-
Y.om
Z
=
rr c
r.mo
H
rc.m-
$ E
e.mo
IS&Q-
1o.m -
,.a0
0
-I
B >
.mo
Fig. 10. Stress (a,) distribution near the crack tip along x-axis at y = 0.025 in.
,.ca-
!a
Fig. 12. Stress (a,) distribution near the crack tip using undersized local regions.
C.
314
T. SLJNand K. M.
Table 2. Comparison of fine mesh and coarse mesh horizontal displacement components at nodes near the crack tip Fine mesh
Node Solution (10e6in.) H J B I
Coarse mesh Solution (10e6in.) Error (%)
-0.2365 -0.9669 -3.521 -0.1528
-0.2121 -0.9712 - 3.582 -0.1512
10.3 0.45 1.74 1.06
It is important to bear in mind that the local region must be selected so that its boundary is at least one coarse element away from the region where the stress field is of interest. The purpose is to minimize the effect of error in the global displacement solution at the boundary nodes. For example, to calculate the stress field in region JKGH, the minimum size of the local region should be ACFI. If region JKGH is used as the local region, then errors may result. To illustrate this point, we use local regions JKGH and KLFG to calculate the near tip stress or along y = 0.025 in. The results are presented in Fig. 12, which also includes the fine mesh solution for comparison. USE OF PARALLEL
PROGRAMMING
In global-local analysis, local analyses can be performed concurrently. Thus, this global-local procedure is natural for the use of parallel programming on multiprocessor computers. The computer used in this study is the MIMD Balance 21000 system by Sequent. It has 12 identical processors (CPUs) and a single common memory. The Balance CPUs are general-purpose, 32-bit microprocessors. This computer runs on the DYNIX operating system. DYNIX balances the system load
Table 3. computing time on VAX 11/780 computer for the cracked panel Method
Time (set)
Fine mesh (1000 elements) Coarse mesh (40 elements) Global-local (2 local regions)
4276 27 249
Table 4. Computing time on Balance 21000 for the cracked panel Method
Time (set)
Fine mesh (1000 elements) Coarse mesh (40 elements) Global-local (2 local regions)
4951 29 176
MAO
among the available processors, keeping all processors busy as long as work is available. In global-local analysis, one processor is employed to perform global analysis with a coarse mesh with the resulting nodal displacements stored in the common memory. Subsequently, a number of processors (depending on how many local regions are selected) are called upon to concurrently perform the local analyses. In each process, conventional Gaussian elimination procedure is used to solve the finite element system of equations. Although intercommunication among these processors is allowed, it is not needed in this analysis. Thus, no time is spent on interprocessor communications. Consequently, if the finite element meshes of the selected local regions are identical, the number of local regions should not alter the computing time. As an example, we use the cracked panel considered in the previous section. The mesh designs remain the same as described previously. The CPU times for the three methods are listed in Table 4. Note that in the fine mesh (1000 elements) solution, only one processor of the Balance system is used. The results shown in Table 4 indicate a speedup of about 30 times if the global-local procedure is employed in conjunction with parallel programming. The speedup will be even greater for larger structures or for structures with multiple cutouts or cracks, where large numbers of finite elements are needed.
CONCLUSION
The present global-local finite element analysis procedure can be used to achieve great time savings in both sequential computers and MIMD multiprocessor computers. In parallel programming, this procedure does not require interprocessor communications and, consequently, the computing time does not depend on how many local regions are considered. Even greater time savings are possible by employing parallel programming in the global analysis and in each local analysis. The accuracy of the global-local procedure is affected by the choice of local regions. One rule is to avoid using global nodes of questionable accuracies (such as the crack tip) as boundary nodes for the local region. Acknowledgement-This research was supported by NASA Langley Research Center under Grant No. NAG-l-581. Helpful discussions with Mr W. J. Stroud are gratefully acknowledged.
REFERENCES
1. A. K. Noor (ed.), Impact of New Computing Systems on Comoutational Mechanics. ASME, Boston, MA (1983). 2. H. F Jordan and P. L. Sawyer, Amulti-m~cropr&essdr system for finite element structural analysis. Comput. Struct. 10, 21-29 (1979).
Global-local
finite element method for parallel computations
3. L. M. Adams, Iterative algorithms for large sparse linear systems on parallel computers. NASA CR166027, NASA Langley Research Center, Hampton, VA (1982).
315
4. K. H. Law, A parallel finite element solution. Compur. Srrucr. 23, 845-858 (1986). 5. R. M. Jones, Mechanics of Composite Materials. Scripta, Washington, DC (1975).