Multigrid and domain decomposition techniques in ...



Finite Element Stress Analysis of FRP Composites Using Multigrid and Domain Decomposition Techniques

B.B.Mahanta1, M.S.Shah1, T. Kant2

1Scientific and Engineering Computing Group

Centre for Development of Advanced Computing

Pune 411 007

INDIA

2Department of Civil Engineering

Indian Institute of Technology, Bombay

Mumbai 400 076

INDIA

Abstract: - This paper deals with the stress analysis of FRP composites in structural mechanics using numerical solutions. Thereby, a finite element method (FEM) is used. The cost of solving the set of linear equations arising from the finite element discretization dominates the computational cost in this kind of problems, when using implicit solution schemes. For the solution of the obtained matrix equation system we have implemented a multigrid method (MG). Further, we have also implemented a domain decomposition technique (DDM) using the Schur complement method. We have experimented these methods in COMPOSIT, our own developed code and some case studies are presented.

Key-Words: - Stress Analysis, FRP Composites Materials, Finite Element Method, Multigrid Method, Domain Decomposition Method

1 Introduction

Structural components made up of fibre reinforced composite materials are being extensively used in high-to-low technology areas in recent years. Their industrial applications are multiplying rapidly because of their superior mechanical properties. The laminated fibre reinforced polymer composites (FRPC) are the combination of two or more materials on a macroscopic scale to form a useful material. These are essentially the heterogeneous materials in which the lamination is used to combine the best aspects of the constituent layers. In composite laminates, the mechanical properties can vary from layer-to-layer. There is also a very thin layer-to-layer interface consisting of matrix dominating material. Usually the properties of this interface affect the interlaminar strength of the composites. A potential problem in the construction of laminates is the presence and importance of interlaminar stresses between layers. Various approaches have been used for the analysis of composite and sandwich laminates. The analytical and experimental methods are not generally applicable, in real life applications, due to intricate modeling and sophisticated requirements and prohibitive costs.

In most of the practical situations composite laminates involve complicated geometry, loading and boundary conditions. Thus some approximate numerical methods, especially finite element methods (FEM), are to be used to tackle these problems [1]. The finite element method (FEM) has proven to be a versatile method for the simulation of continuous physical systems in many areas of science and engineering [2] FEM requires a discretization of a domain into a finite element mesh. This discretization, along with FE formulations, is used to assemble a set of simultaneous equations, which when solved provide an approximate solution of the system under consideration. The cost of solving the linear equations Au=f (A is a sparse stiffness matrix and f is a residual vector) dominates the computational cost of FEM, when using implicit solution schemes. Direct methods, based on Gauss elimination, are commonly used for the solution of these systems, as they are easy to use and efficient for moderately sized matrices. Accurate models of complex systems, however, often require that large matrices be solved; the asymptotic costs of direct methods are high in comparison to iterative methods. Iterative methods in turn benefit greatly from preconditioning of matrix A.

Multigrid (MG) is known to be one of the most effective solution methods for matrices that arise from discretized PDEs.[3,4,5] These methods are based on the multilevel or multi-scale paradigm. Our work involves using multigrid method to solve the equations arising from the finite element stress analysis of FRP composites. Further we have also made a study of non-overlapping domain decomposition method based on Schur complement method. Domain decomposition techniques are efficient tools for solving large problems in engineering analysis. [6,7] Despite their usefulness in dividing large computational tasks in the framework of traditional computer architecture, these methods are particularly attractive regarding parallel computations.

2 Mathematical Formulation

2.1 Finite Element Formulation

The total potential energy P of the laminate of middle surface area A enclosing a space of volume V can be written as [8]

P = U - W (1)

where U is the strain energy stored in the laminate, W is the work done by externally applied load.

The standard strain-displacement and the constitutive stress-strain relations are used to relate displacements and laminate stresses. The basic equation to be solved in static analysis can be written as,

Au=f (2)

where A is global, symmetric, positive-definite stiffness matrix, u is displacement vector and f is load vector.

2.2 Multigrid Methods

The basic iteration or relaxation methods are good in damping the oscillatory components of the error.[9] They are not effective on the smooth modes. So the idea behind multigrid is to use relaxation on the fine grid until oscillatory modes are sufficiently reduced and convergence slows down. Then we move to the coarse grid where the (fine grid) smooth modes are more oscillatory and hence more effectively reduced by relaxation.

Simple (and inexpensive) iterative methods like Gauss-Seidel, Jacobi, and block Jacobi are effective at reducing the high frequency error but are ineffectual in reducing the low frequency error. These simple solvers are called smoothers as they “smooth” the error (actually they reduce high energy components, leaving the low energy components, which are “smooth” in, for example, Poisson’s equation with constant coefficients). The ineffectiveness of simple iterative methods can be ameliorated by projecting the solution onto a smaller (coarse) space, that can resolve the low frequency content of the solution, in exactly the same manner as the finite element method projects the continuous solution

onto a finite dimensional subspace to compute an approximate solution. This coarse grid correction is then added to the current solution. Thus, the goal of a multigrid method is to construct, and compose, a series of function spaces in which iterative solvers and small direct solves can work together to economically reduce the entire spectrum of the error [10]

Figure 1 shows the multigrid V-cycle, using a smoother, restriction operator R that maps residuals from the fine grid i to the next coarse grid i + 1, and prolongation operator P to map corrections from coarse grid i + 1 to fine grid i. The basic algorithm is as follows

procedure MG(level,Ah,uh,fh)

if level=coarsest then

solve coarse grid equation

else

Pre-smooth Ahuh = fh (m1 times)

Compute residual rh = fh - Ahuh

Restrict fH = Rrh

Call MG(level - 1, AH,vH=0,fH) (mc times)

Interpolate vh = PvH

Correct uhnew = uh + vh

Post-smooth Ahuh = fh (m2 times)

[pic]

Fig 1. Multigrid ‘V Cycle’

2.3 Domain Decomposition Methods

We assume that the domain ( with boundary (( is partitioned into N non-overlapping subdomains (1, (2, …, (N with boundaries ((1, ((2,…, ((N. After discretization, we obtain a linear system Au=f where the matrix A is sparse, unstructured, and symmetric positive definite. Let B be the set of all the indices of the discretized points, which belong to the interfaces between the subdomains. Grouping the points corresponding to B in the vector uB and the points corresponding to the interior I of the subdomains in uI, we get the reordered problem.[6,11]

[pic] (3)

The matrix AII can be reordered to a block diagonal matrix in which the i-th diagonal block corresponds to the internal variables of subdomain (I. Eliminating uI from the second block row of Equation (3) leads to the reduced problem

[pic] (4)

where

[pic] (5)

is the Schur complement matrix of the matrix AII in A The matrix S is symmetric positive definite.

Let ( be the interface between the subdomains. The local interface of a subdomain can be defined as ( i = ((i /((. Let Ri : (( ( i be the canonical point wise restriction which maps vectors defined on ( onto vectors defined on (I, and let RTI : (i( ( be its transpose. For a stiffness matrix A arising from a finite element discretization, the Schur complement matrix can be written as the sum of N (local) smaller Schur complement matrices

[pic] (6)where

[pic] (7)

where the local Schur complement matrix S(i) associated with subdomain (i is computed from the subdomain stiffness matrix A(i) , defined by

[pic] (8)

In a parallel (multi-processor) computing environment, each subdomain matrix A(i) is assigned to one processor. In this way, the local Schur complement matrices S(i) can be computed simultaneously. The (global) Schur complement matrix S is available through Equation (8) and is never assembled explicitly.

The implementation of this method usually consists of three steps. They are summarized in the algorithm as follows

Step 1: Factorize the matrix AII and compute the Schur complement matrix S.

Step 2: Solve the Schur complement system for uB

[pic] (9)

Step 3 : Solve the system for uI

[pic] (10)

Since the matrix AII is a block diagonal matrix, steps 1 and 3 each consist of N independent smaller steps (one for each subdomain). Furthermore, since each subdomain matrix A(i) is assigned to a different processor, these N steps can be performed in parallel. When the number of subdomains increases, the amount of parallelism naturally increases in steps 1 and 3 (but the total amount of work decreases). The solution of the Schur complement problem in step 2 becomes more complex when the number of subdomains gets larger; the Schur complement matrix S is not assembled and therefore efficient parallel direct or iterative schemes are required for step 2.

3 COMPOSIT Program

We have developed COMPOSIT program, which is a finite element code for the stress analysis of FRP composite structures. [12] A higher-order facet element has been used for this. Many commercial packages have been developed based on FEM, which can handle composites analysis. But these packages are based on either classical lamination theory (CLT) or the first-order shear deformation theory (FOSTs). These approaches either neglect the effects of transverse shear strains (as in CLT) or they assume a constant shear strain through a laminate thickness (as in FOST)[13].

As a general methodology C1 continuous FEMs can not be used, as these require first derivative of displacement also to be continuous, thus requiring complicated formulations. Hence, the approach of using simple C0 isoparametric FE formulation is generally attractive due to ease in software development. In most of the 2D theories, developed to predict the membrane-flexural analysis of composite laminates, the assumed cross-sectional deformation of a laminate appears to be reasonably unrealistic. The higher-order shear deformation theory is used here [14,15] which considers the realistic non-linear variation of displacement across the laminate thickness.

Essentially Taylor's series expansion method is used to deduce a 2D formulation of a 3D elasticity problem. The displacement model considered in the present software is given below. This model does not use complete material constitutive law.

[pic]

In the above relations, the terms 'u', 'v' and 'w' are the displacements of a general point (x,y,z) in the laminate domain in the 'x', 'y' and 'z' directions respectively. The parameters [pic] are the inplane displacements and [pic] is the transverse displacement of a point (x,y) on the laminate middle plane. The functions [pic] are the rotations of the normal to the laminate middle plane about x and y axes respectively. The parameters [pic] are the higher order terms in the Taylor's series expansion and they represent higher-order transverse cross sectional deformation modes (see Fig 2).

[pic]

Fig 2 : Laminate geometry and positive set of lamina/laminate reference axes, displacement components and fibre orientations.

The theories account for warping of transverse cross section and parabolic variation of transverse shear deformation through the plate/shell thickness. Nine noded Lagrangian and four and eight noded Serendipity quadrilateral elements with provision for selective numerical integration scheme based on Gauss-Legendre product rules are used. In post processing phase, the values of stress-resultants, in plane and transverse interlaminar stresses are computed in the local co-ordinate system within the element at the Gauss points. The inplane stresses are computed by using constitutive relations. Transverse interlaminar stresses are evaluated by using the 3-D/2-D elasticity equilibrium equations. Transverse interlaminar stresses are evaluated by using direct integration; finite difference and exact surface fitting methods and these are evaluated separately in each element. In general this stress analysis, for a real life problem, imposes great demands on computational requirements. The solver used in COMPOSIT for the solution of the linear set of equations is a direct Cholesky solver.

4 MG Implementation

The multigrid algorithm that has been implemented in the COMPOSIT code uses the following algorithmic components.

• Gauss Seidel iterative method has been used as the smoother.

• For restriction of the fine grid variables to coarse grid variables, full weighing restriction operator has been used.

• A direct cholesky solver has been used as the coarse grid solver at the lowest level.

• Prolongation of the coarse grid variables to the fine grid is done using bilinear interpolation operator.

• In case of slow convergence, acceleration is obtained using Successive Over Relaxation (SOR) method.

4.1 Case Studies and Numerical Validation

The numerical solution obtained by implementing the multigrid method in COMPOSIT is validated by solving the following stress analysis problem. A two layered angle-ply(-150/150) square plate is analyzed[16] The plate is of unit length and the thickness equal to 0.1 units. It is clamped on the edges with a uniformly distributed load of unity. The material properties of lamina are as follows:

[pic]

The problem has been solved using a 2, 3 and 4 level multigrid V- cycle. It was observed that there is no major difference in the performance of the multigrid algorithm with increasing levels. The performance depends on various parameters. More detailed studies have to be carried out to evaluate the performance of MG algorithms in this class of problems.

5 DDM Implementation

The sequence of steps constituting the domain decomposition implementation to the COMPOSIT code is as follows:

• The finite element domain is decomposed into N smaller subdomains.

• The Schur complement S(i) for each domain is computed.

• The Schur complement system is solved using a direct or iterative solver to obtain the boundary unknowns. However, explicit computation of the global Schur complement is an expensive operation thus iterative solvers such as the preconditioned conjugate gradient algorithm give better efficiency.

• Finally, the solution for the interior unknowns is obtained at each subdomain using the computed interface boundary values.

5.1 Case Studies and Numerical Validation

The numerical solution obtained by implementing the domain decomposition method in COMPOSIT is validated by solving the following stress analysis problem. A two layered angle-ply(-150/150) square plate is analyzed. [13] The plate is of unit length and the thickness equal to 0.1 units. It is clamped on the edges with a uniformly distributed load of unity. The material properties of lamina are as follows:

[pic]

This case has been solved by dividing the problem into 4 domains.

6 CONCLUSIONS

Multigrid and domain decomposition techniques have been developed and implemented in the COMPOSIT finite element code. The numerical results obtained were validated with existing results. The multigrid solver was converging slowly for the kinds of problems solved. Domain decomposition techniques are found more suitable for solving structural mechanics problems compared with MG methods. DD techniques are efficient for solving large problems. Despite their usefulness in dividing large computational tasks in the framework of traditional computer architecture, these methods are particularly attractive regarding parallel computations.

References: -

[1] Reddy J.N., On the generalization of

displacement based laminate theories, ASME

App. Mechanics Reviews, Vol 42(11), Part 2,

1989, pp. 213-222

[2] Zienkiewicz, O, C., The Finite Element

Method, McGraw-Hill, 1989

[3]Adams, F. M., A Distributed Memory

Unstructured Gauss-Seidel Algorithm for Multigrid Smoothers, Proceedings Supercomputing '0, 2001

[4] Adams, F. M, Multigrid Preconditioners

for Finite Element Method Matrices



reports_pre1998/algorithm_development

/../mg.htm

[5] Saha, A. K., Biswas G, Application of a Multigrid Method for Flows Behind Bluff

Bodies in a Channel, Journal of Institution of Engineers (Aerospace), Vol 77, 1996

[6] Al-Nasra, M and Nguyen, D. T., An

algorithm for domain decomposition in

finite element analysis, Computers and

Structures,Vol. 39, No. ¾, 1991, pp. 227-

289

[7] Nguyen, D. T., and Tungkahotara, S,

Parallel Finite Element Domain

Decomposition For Structural/Acoustic

Analysis., J. of Computational and Applied

Mechanics, Vol. 4, No 2, 2003, pp. 189-201

[8] Hinton, E. and Owen, D.J.R,, Finite Element Software for Plates and Shells, Pineridge Press Ltd., UK, 1977

[9] Hanson V. E., A Multigrid Tutorial, Virtual Proceedings, Copper Mountain Conference on Multigrid Methods, 2003

[10] Adams, F. M., and Demmel, J., W.,

Parallel Multigrid Solver for 3D

Unstructured Finite Element Problems,

Proceedings Supercomputing ’99,, 1999

[11] Martinez Aja, A., Sub-modelling techniques for static analysis, Proceedings MSC Software First South European Technology Conference, 2000

[12] Kant T., Shah M.S., Mahanta, B. B., Thermo-mechanical analysis of fiber reinforced composite plates and shells on supercomputers, Proceedings Structural Engineering Convention (SEC2003), IIT Kharagpur, India, 2003

[13] Kant, T., Owen D.R.J., Zienkiewic, O.C., A refined higher-order C plate bending element, Computers and Structures, Vol 15(2), 1982, pp. 177-183

[14] Kant T., Manjunatha, B.S., An unsymmetric FRC laminate C0 finite element model with twelve degrees of freedom per node, Int. J. Engg. Computing., Vol 5(4), 1988, pp. 300-308

[15] Kant T., Menon M.P., Higher-order theories for composite and sandwich cylindrical shells with C0 finite elements, Computers and. Structures, Vol 33, pp. 1191-1204, 1989

[16] Khare, R. K., Thermal Stresses in Fibre Reinforced Composite Plates and Shells- Some Studies, ,PhD Thesis, SGSITS, Indore, 1996

................
................

In order to avoid copyright disputes, this page is only a partial summary.

Google Online Preview   Download