Advances in Engineering Software Vol. 29, No. 3–6, pp. 249–263, 1998 q 1998 Published by Elsevier Science Ltd Printed in Great Britain. All rights reserved PII: S 0 9 6 5 - 9 9 7 8 ( 9 8 ) 0 0 0 0 9 - X 0965-9978/98/$19.00 + 0.00
Advanced manufacturing of large-scale composite structures: process modeling, manufacturing simulations and massively parallel computing platforms R. V. Mohan a, K. K. Tamma a,*, D. R. Shires b & A. Mark b a
Department of Mechanical Engineering, Institute of Technology, 111 Church Street S.E., University of Minnesota, Minneapolis, MN 55455, USA b US Army Research Laboratory, Aberdeen Proving Ground, MD 21005, USA Newly emerging composite manufacturing processes, where there exist only limited industrial experience, demonstrate a definite need for process simulations to reduce the time and cost associated with the product and process developments. The physical and geometric complexity of the net-shape composite structural parts produced by resin transfer molding (RTM) lead to larger computational problem sizes and computational complexity demanding more memory and computational time which can be provided by high-performance and massively parallel computing platforms. The associated computational costs and times for large-scale simulations is a major factor for the effective use of massively parallel computing platforms. However, the limited availability of these massively parallel computing platforms and the demand for their computing resources make it imperative to have improved computational methodologies and algorithms for physical problems for effective use of the massively parallel computing platforms. An effective pure finite element methodology for the process simulation analysis of mold filling in RTM is employed for large-scale process simulation of composite structures and the computational effectiveness of the methodology is demonstrated for large-scale process simulations. The Connection Machine CM-5 is the massively parallel computing platform employed in this study. Issues regarding the implementation and software development of these manufacturing process simulations, behavior of the linear system solution methodologies during largescale composite manufacturing process simulations, computation and communication loads, performance and scalability are specifically addressed. Demonstrative process simulation applications involving practical large-scale composite structures such as the Comanche helicopter keel beam and an aircraft composite wing box section are also presented. q 1998 Published by Elsevier Science Ltd. All rights reserved Key words: large-scale composite structures, process modeling, massively parallel computing platforms.
platforms is clearly seen in traditional computationally intensive physical engineering application areas such as computational fluid dynamics and computational structural mechanics. These new computational platforms also have initiated research in many areas of computational sciences including applied areas of linear and non-linear algebra, iterative solvers, graph partioning and domain decompositions, as well as in fundamental computer science areas such as compiler research, high-speed networking, computer communications, etc., involving developments both in theoretical and software aspects. However,
1 INTRODUCTION The advent of high-performance and massively parallel computing platforms has a significant impact in various aspects of engineering, mathematical and physical sciences. The enormous computational power, memory and speed of these computing platforms now permit realistic modeling and simulation of large complex physical problems which were once impossible. The impact of the computing *Author to whom all correspondence should be addressed. 249
250
R. V. Mohan et al.
Fig. 1. RAH-66 comanche helicopter: new composite weapon system.
high-performance computing is yet to make a definitive impact in application areas relevant to manufacturing. The use of high-performance and massively parallel computing in manufacturing process simulations is very limited and has been attempted in very special applications. From an industrial perspective, new manufacturing processes have a limited industrial experience and data base, and can benefit significantly from computer modeling and simulation of the manufacturing processes. One of the emerging manufacturing processes employed in the production of high-performance composite structural configurations is resin transfer molding (RTM). The process involves placing a fiber preform in a matched mold cavity and injecting a reactive liquid resin into the mold cavity through the porous fiber medium via one or more injection ports. The liquid resin flows through the fiber preform impregnating the fiber layers. As the resin fills and impregnates the fiber preform in the mold cavity, the mold is heated and the part is cured with the liquid resin undergoing polymerization reactions to become rigid. After
curing, the net-shape high-performance composite part with the desired fiber orientations is removed from the mold. RTM has a significant potential for cost-effective manufacturing of net-shape composite structural parts. Large and complex parts can be manufactured by RTM. One such illustrative aerospace application is the RAH-66 Comanche (Fig. 1) with its composite keel beam which is targeted to be manufactured as a single piece. Another illustrative application is the wing box structure (Fig. 2). Other large complex structural applications include the use of the process to manufacture thin and thick composite structural parts both in the defense and in the commercial industry. However, manufacturing of these large complex parts poses significant challenges to obtain a good process maturation and a working production methodology. Manufacturing process simulations involving the understanding and analysis of the physical process can be employed to obtain improved process conditions based on parametric study optimizations. Such process simulation studies reduce the time and the cost associated with repeated trial and errors
Fig. 2. Composite wing section.36
Advanced manufacturing of large-scale composite structures involving mold and tool design changes which are cost- and time-prohibitive, especially in large complex parts. The acceptability of the manufactured composite part and the part strength with its desired fiber orientation is influenced by the complete impregnation of the fiber preform during the filling stage without any void formation, resin-rich areas due to race tracking, dry spots, etc., affecting both the part strength and quality. A clear understanding of the mold filling process and process simulation studies based on physical models for the RTM mold filling process are crucial for the successful manufacture of RTM structural components with high quality. Numerical methods and, in particular, the finite element method with its ability to model complex geometries and material variations form the basis of the process simulations. Complex large-scale geometric structures involve finite element discretized meshes with large numbers of elements and nodes to accurately represent the geometric and material variations thereby, increasing the computational problem sizes during the simulations. The geometrical and computational complexity is significantly increased when large threedimensional complex geometries are involved. The larger computational problem sizes and computational complexity demand more memory and computational time. For larger computational domains involving large problem sizes (larger number of nodes and elements in the case of the finite element method), massively parallel computing platforms have shown great promise and have been employed in finite element calculations involving structural, thermal and fluid dynamics.1–8 Large-scale composite manufacturing process simulations involving large-scale complex composite structures are focused upon in this paper. Attention is restricted to isothermal situations merely for illustration of the present developments. The overall efforts emanate from a general purpose simulation code entitled ‘On Composites Technology Of Polymeric Useful Structures’ (OCTOPUS) for applications to isothermal/non-isothermal situations in 2-D/3-D thin and thick composite parts, and with capabilities on different computing platforms.
2 RESIN IMPREGNATION AND PROCESS MODELING IN RTM The resin impregnation and mold filling process in RTM for composite structures involves the flow of a polymeric resin through a fiber network until the mold is completely filled. The process involves interaction between the polymeric resin in the liquid phase and the fiber porous medium at a number of different length-scales. The manufacturing process simulations are based on macroscopic flow models which are described by Darcy’s law9 relating the fluid flowrate to the pressure gradient, fluid viscosity and the permeability of the porous medium. Physically, resin impregnation and flow permeating through a porous media is a free surface moving boundary problem whereby the
251
field equations and the free surface have to be solved and tracked. Physical-based numerical simulations can be carried out either with fixed mesh or moving mesh methods. In mold filling and resin flow impregnation situations such as those seen in RTM, the primary interest is in the temporal progression of the resin inside a complex fixed mold cavity filled with fiber preform. Considering the geometric, material complexities of composite structural components leading to diverging and merging flowfronts and race tracking effects, Eulerian fixed mesh approaches are effective in RTM flow impregnation process simulations.
3 LARGE-SCALE MANUFACTURING SIMULATIONS AND MODELING ALGORITHMS Large complex structures need models which actually capture accurately both the geometric and physical phenomena. This may involve the use of three-dimensional flow simulations as warranted by the geometry, thickness and fiber preforms employed. There also exists a need for accuracy improvements by refining the discretization of the computational domain. All these have a serious impact on the computational time and power requirements. Massively parallel computational platforms have a potential to keep the computation times within reasonable limits. Physical modeling and computational algorithmic methodologies play an important role in computational times. For large-scale computations, it becomes critical to have algorithms that are physically accurate and permit a faster solution of the computational physical problem. It is not only essential to have optimal data structures and communication strategies for massively parallel computing architectures, but it is also very important to have improved computational algorithms and methodologies to further improve the computational performance of large-scale simulations. RTM flow modeling and mold filling process simulations involving a pressure-driven flow through a complex anisotropic fiber preform have traditionally employed a commonly known approach termed the finite element– control volume10–14 method to solve for the field equations and track the progressing resin front. The finite element–control volume method treats the transient mold filling problem as a series of quasi-steady state problems and is an explicit filling approach. The time step increment between each of the quasi-steady states is restricted by a stability condition that only the smallest control volume regions are completely filled during any time increment based on the local flowrates and the unfilled mesh control volumes. This stability condition restricts the time step increments strictly when large structures and refined meshes are involved as the time step increment is limited by the geometric size of the control volume regions. The time step increment restriction increases the number of quasi-steady state analysis steps significantly, thereby, increasing the computing cost and time involved in the
252
R. V. Mohan et al.
large-scale process simulations. Further, front progression updates involve computation of velocity fields, the associated communications between element level computations and nodal level computations further adding to the computational cost. New effective methodologies with improved physical accuracies and reduced computational cost are necessary to tackle large-scale manufacturing process simulations. These significantly improve the affordability of large-scale process simulations. In this context, a pure finite element methodology which is an implicit filling technique and is based on a transient resin mass balance equation was originally developed by Mohan et al.15–17. This methodology has proven to provide a physically accurate, computationally faster and algorithmically better solution strategy for process flow modeling simulations involving a pressure-driven flow through a complex anisotropic preform which can be either two-dimensional in the case of thin composite sections and three-dimensional in the case of thick composite sections. The pure finite element methodology does not have the time step restrictions of the finite element–control volume methodology. Furthermore, with fill domain solution updates based on efficient matrix– vector products, we have improved computational times by a multitude of factors for unstructured meshes involving large-scale complex geometries. This paper discusses the effectiveness of the pure finite element methodology for large-scale composite manufacturing process simulations employing a massively parallel computing platform. Of the several massively parallel computing platforms available, the Connection Machine system CM-5 was available for our use and was used as a target architecture. Discussions involve the computational efficiency of the pure finite element methodology for large-scale composite process simulations, linear system solution methodologies and their behavior during largescale composite manufacturing process simulations, computation and communication loads, performance and scalability. Demonstrative process simulation applications will involve illustrative large-scale structural composite components. Though the discussions are based on the massively parallel computing platform CM-5, the computational strategies can be extended to other massively parallel computing systems with appropriate programming models and implementation constructs. The physical and computational effectiveness of the pure finite element methodology over the finite element–control volume approaches will be seen in any computational architecture and has been demonstrated in computing architectures which includes workstations, traditional supercomputer systems and other massively parallel computing platforms. This is due to the fact that the physical and computational effectiveness is due to the improved physical and computational methodology employed and not to the implementation. This clearly emphasizes the need for new physically accurate and computationally efficient algorithms and methodologies for the
physical modeling of engineering problems to obtain optimal performance of the numerical simulations.
4 FLUID FLOW BEHAVIOR AND MODELING IN RESIN TRANSFER MOLDING The macroscopic flow behavior in LCM processes through the fiber preform is modeled as a flow through a porous media. The flow through a porous media has traditionally been described by Darcy’s law.9 Darcy’s law relates the fluid flow rate to the pressure gradient, fluid viscosity and the permeability of the porous medium: u¯ ¼ ¹
K =P m
(1)
where u¯ is the volume averaged velocity of the resin at a macroscopic scale, m is the viscosity of the resin and K is the permeability of the porous fiber preform. The permeability is a measure of the ease of flow through the porous medium and is a function of the fiber architecture and porosity of the medium. This macroscopic flow velocity model is employed in understanding the mold fill behavior such as the front progression, flow front location, mold filling time, pressure distribution, effect of location of injectors, vents, injection flow rates and pressures. From a computational perspective, assuming isothermal conditions, the main issue in RTM fluid flow modeling is the tracking of the resin front surface and the pressure field as they progress through the porous media. The geometrical complexity of the parts and the need to analyze the flow progression behavior inside a closed mold cavity makes Eulerian fixed mesh approaches effective for such simulations. Most of the RTM flow modeling simulations employ a finite element–control volume approach and philosophy to track the flow field and the filled regions inside the mold cavity.10–14 The method involves associating a parameter called the fill factor with each node that is associated with a control volume region. The fill factor indicates whether a given control volume region is full or empty. The transient mold filling problem is treated as a quasi-steady state problem. The flow front progression is based on time increments between each of the quasi-steady states. The time increment employed allows only the smallest unfilled/partially filled control volume mesh region to be completely filled in a single time increment step. This technique is normally termed an explicit filing technique. The restriction of the time step increments ensures the stability of the quasisteady state approximations and the satisfaction of Courant conditions. However, the time step increment restriction increases the number of quasi-steady analysis steps significantly for practical process simulation applications involving complex parts, thereby increasing the computing time and cost associated with the process simulations. For very large-scale problems of complex geometrical shapes, such approaches may potentially hinder the completion of an analysis.
Advanced manufacturing of large-scale composite structures Other approaches that are cited in the literature for tracking of the flow fronts in Eulerian fixed meshes such as those in RTM simulations include: (1) the marker cell (MAC) approach, which is a common technique employed in metal casting simulations;18 and (2) the volume of fluid (VOF) method, which is another method that is commonly employed in mold flow simulations in metal castings.19–22 In view of the difficulties and deficiencies cited earlier, new effective methodologies with improved physical accuracies and reduced computational costs are necessary for affordable process simulations. This significantly improves the affordability of process simulations and permits quicker simulation runs for various real-time parametric studies to obtain optimal process conditions. In this context, a pure finite element methodology originally developed by Mohan et al.15–17 and proposed as part of the present efforts is an implicit method of computation and is based on the transient resin mass balance equation and is described next. Some significant advantages of this method are that the computed location of the flow front at any instant of time does not depend on the time step size employed to reach that stage and there exists no need to form control volumes from typical generated finite element meshes. 4.1 Pure finite element methodology In what we originally termed ‘the pure finite element methodology’ for applications to large-scale practical problems, the objective of the resin transfer mold filling analysis is the conservation of the resin mass at any instant of time. The fill factor implicitly defines the amount of resin present at any time inside the mold cavity and its distribution. Considering a general mold cavity forming the Eulerian domain Q shown in Fig. 3, at any instant of time part of the mold cavity is filled with the resin and the fill factor W
253
denotes the status of the mold cavity. In the filled regions, the fill factor W is greater than zero with a maximum value of the fill factor being equal to one. In the unfilled regions, the fill factor is equal to zero. For a constant density resin based on incompressible flow conditions, and Darcy’s flow field, the modified mass conservation equation for the resin mass can be represented as Z Z ] K =P dQ (2) W dQ ¼ W=· Q ]t Q m The modified mass balance equation involves the fill factor and the pressure field. In the case of non-isothermal situations involving polymerization curing conditions, the resin viscosity m is a function of the temperature and the degree of cure and will continually vary during the mold filling. The resin flow inside a mold cavity based on Darcy’s flow approximation is a pressure-driven flow and the pressure gradients are negligible (=P < 0) in unfilled and partially filled regions, W ¼ 1 in the completely filled regions, and the variation in the fill factor is similar to those in the VOF method.19–21 These variations are similar to those in the so-called finite element–control volume method. The modified form of the mass balance equation involving the fill factor and the pressure field represents the governing model equation for the present finite element formulations and is given by Z Z ] K =P dQ (3) W dQ ¼ =· Q ]t Q m A differential form of eqn (3) is given by ]W K ¼ =· =P ]t m
(4)
The boundary and initial conditions for the mass balance
Fig. 3. General partially filled mold cavity.
254
R. V. Mohan et al.
equation can be stated as follows: At mold surface
]P ¼0 ]n
At resin front P ¼ 0 At mold inlet P ¼ P0 (prescribed pressure), or
un ¼ ¹
K (=P) ¼ u0 (prescribed flow rate) m
At injection ports W(t $ 0) ¼ 1:0
ð5Þ
These forms of the mass balance equation involving the fill factor and the pressure field and the associated boundary conditions are discretized next using the finite element method. Finite element discretizations for the pressure field and the fill factor field are introduced and application of Galerkin finite element method leads to a discretized system of equations for the pressure and the fill factor which are solved in an iterative manner. The pure finite element method for mold filling analysis in RTM has significant physical, algorithmic and computational advantages which are described elsewhere.15–17 The pure finite element method briefly described does not involve the restrictions of the time step increments seen in the explicit finite element–control volume method and computes the position of the flow front at each of the discrete time steps which are selected by the analyst. This leads to significant advantages in process simulation times, improving the affordability of process simulations and thus further advancing efforts towards making real-time interactive simulations and virtual process manufacturing a reality.
5 HARDWARE ARCHITECTURE: THE CONNECTION MACHINE CM-5 SYSTEM The massively parallel computing platform employed in the present process simulation developments is the connection machine system CM-5,23–26 which is a distributed memory computing architecture supporting both data parallel and message passing architectures in a single integrated universal architecture. The simulation implementations are based on a data parallel paradigm involving single instruction multiple data (SIMD) programming model and the language constructs based on CM Fortran (CMF) involving Fortran 90 constructs and additional data layouts to provide the necessary programming implementations. The CMF programming language syntax allows for the distribution of Fortran array elements within the distributed processing nodes. 5.1 Data structures and implementation issues Data structures which are representative of the massively parallel computing architectures are necessary to achieve good performance for the physical problem on massively
parallel processing systems. Johnson and Mathur2 discuss a description of possible data structures for the finite element method along with the storage and arithmetic requirements. In general, for the data parallel programming models employed in the simulations, two levels of parallelism are possible in the finite element computations: (1) parallelism at the element level with each finite element forming the data parallel object; (2) nodal level parallelism with each finite element node forming the parallel data object. RTM process simulations require both these primary parallel level data objects. Traditional explicit time integration implementations of transient thermal and structural problems27–30 employ both element level and nodal level data parallel objects. For the process simulation implementations, data structures are suitably selected based on both the physical and computational aspects in various kernels and the overall objective is to reduce the amount of any associated communications between the element level and nodal level data structures. 5.2 Communication issues In the data parallel programming implementations with the element level and the nodal level data structures each residing in their own processing nodes, and the exchange of information between the two, which in the finite element implementations depend upon the element connectivities, inter-processor communications across the processing nodes are necessary. Such communications can be a substantial part of the total run-time, thus affecting the overall efficiency of the implementations. The current CM-5 compiler has no means of recognizing the communications that are needed across the processing nodes. Element mapping techniques which preserve the locality of the element connections and the interconnected nodes to the adjacent processing nodes to reduce the communication distances have been considered.31,32 Such efficient mappings involve a substantial processing time. The computational mesh currently employed in process modeling simulations is Eulerian and the optimal routing algorithms for the routing of the data between processing nodes can be used. Due to the Eulerian nature of the mesh, once the communication pattern is set up, additional set-ups are not necessary and the created communication patterns can be repeatedly used. Connection Machine System Scientific Library (CMSSL) routines to generate the dual mesh connectivity to layout the elements in the nearby processors based on the hardware architecture and the processor layout are available. These routines based on Lanczos type of algorithms improve the communication pattern. However, the current versions of these library routines have limitations. For example, 2-D thin shell elements have no more than two elements to an edge. But, for 2-D thin shell elements used in 3-D space made of built-up structures, there can be more than two elements to an edge. Mold filling finite element implementations of the Connection Machine CM-5 necessitate two different types
Advanced manufacturing of large-scale composite structures of communication operations: (1) gather—the information available at the global nodal level is gathered to the element level which is local to the nodes of a given element; (2) scatter—the information available at the element level for each local node of the element is scattered back to the nodes. Such an operation involves sending the element level data to the global nodes which are shared by different elements and adding the entries at the global nodal level. CMSSL communication primitives have been employed to perform these required communications. 5.3 System solver and distributed memory data parallel RTM process modeling simulations Traditional finite element computations have normally employed an assembly operation of the element level stiffness matrices into a global system of linear equations which are normally sparse. Direct linear system solvers have been normally employed for RTM process modeling simulations involving an incompressible flow of the resin through a mold cavity. In a parallel computing environment, both computational and inter-processor communication aspects have to be taken into account. In the programming model and hardware architecture environment employed in the Connection Machine CM-5, the direct solvers involve fill-ins across the matrix entries requiring a lot of interprocessor communication, making it expensive. Direct sparse matrix solution techniques involving highly unstructured meshes have not yet been implemented effectively on these massively parallel computing platforms, except in some special cases involving narrow bandwidths such as those in tri-diagonal systems. Instead, alternate techniques which do not require the communication-bound assembly operations as well as the fill-ins and the associated communications of the direct solvers are employed. One such solution technique that is commonly advocated and has been employed, the iterative solution technique for large systems of unstructured equations arising from the finite element computations, is described next. Resin transfer process modeling simulations require the solution of an equation system of the form Ax ¼ b. A general iterative conjugate gradient residual approach to solve the above system of equations is given by:33 Start: Set p0 ¼ r0 ¼ b ¹ Ax0 For i ¼ 0, 1, 2, …, Until convergence do Compute ai ¼ (ri , Api )=(Api , Api ) x i þ 1 ¼ x i þ a i pi ri þ 1 ¼ ri ¹ ai Api X pi þ 1 ¼ ri þ 1 þ ij ¼ 0 Bij pj Bij are chosen such that (Api þ 1 , Apj ) ¼ 0 for 0 # j # i: In the above, the operation (·,·) denotes an inner product, and an operation of the form Aui indicates a matrix–vector product involving the matrix A and the vector u. Thus, the conjugate gradient-based iterative solvers involve evaluation of matrix–vector products for the evaluation of the
255
residuals. In the data parallel programming and hardware architecture model of CM-5, these residuals can be computed at the element level in parallel. The required matrix–vector products are performed for each of the elements. These element level products can be scattered back to the nodes in the computation of residuals, which are employed to check on convergence. These techniques are similar to the element by element techniques.34,35 The matrix–vector product Au can be considered as: ! nelem X A e ue (6) Au ¼ e¼1
where Ae is the dense element level stiffness matrix and the summation is the assembly over all the elements. Thus, the matrix–vector products can be evaluated at the element level thereby avoiding an explicit assembly operation of the stiffness matrices. The programming model and the data structures employ elements as data level parallel objects and these matrix–vector products can be done in parallel. The element level products and the residuals are then scattered into global nodal vectors for use in the inner product evaluations of the conjugate gradient method. The scatter communication set-up and the scatter communication routines are employed for this and the inner products are computed as nodal level parallel operations. Preconditioning techniques are employed to improve the convergence of the linear system of equations and the computations are performed at the element level. The above iterative solution strategy does not involve an explicit assembly of the stiffness matrices and the interprocessor communication associated with such operations.
6 RTM PROCESS MODELING AND CM-5 EQUATION SYSTEM SOLVER PARALLEL IMPLEMENTATIONS The major significant cost in the numerical computations involving the flow modeling and process modeling simulations is the solution of the linear system of equations. The linear system of equations which are diagonally dominant continually changes during different stages of mold filling based on the fill factors associated with the nodes. For the empty and partially filled regions of the mold, the linear system of equations arising from the finite element discretizations of the physical model equations are modified so that the diagonal entry corresponding to the node which is in the unfilled or partially filled region is one and the non-diagonal entries are zero. The right-side vector is also modified so that a zero gauge pressure is obtained in the regions which are partially filled or empty. Thus the linear system of equations continually changes during the transient progression of the front. In a typical RTM mold filling simulation, the ratio of the total number of injection port nodes to the total number of nodes in the mesh is typically small. In addition, the linear system tends
256
R. V. Mohan et al.
to be nearly an identity matrix with the condition number being nearly one during the beginning and early stages of the mold filling. As the mold filling progresses, the nature and the condition number of the matrix system change continually. The effect of these changes are clearly evident in the iteration counts for the system solution convergence as the filling progresses. 6.1 Preconditioning To improve upon the convergence of the iterative solution methods, preconditioning and scaling methods are employed in which a modified system with an improved condition number is employed during the solution procedure. Since the resulting finite element matrix equations in mold filling simulations are positive definite and diagonally dominant, a diagonal preconditioner which has shown good convergence in many applications is employed. Pre-conditioned iterative methods involve a solution of the modified system of equations. For a linear system of equations given by, Ax ¼ b
(7)
the modified system of equations is given by A˜ x˜ ¼ b˜
(8)
where A˜ ¼ W L¹ 1 AW R¹ 1
(9)
x˜ ¼ W R x
(10)
b˜ ¼ W L¹ 1 b
(11)
The scaling in these preconditioning methods is performed by using a symmetric, positive definite matrix and setting the left and right preconditioners to be the same so as to set W L ¹ W R ¼ W 1=2
(12)
The matrix system involved in the mold filling simulations is positive definite and involves positive entries on the diagonal. In this case, the choice of the preconditioner,
which is also a choice in incompressible flow problems, is given by W ¼ diagA
(13)
Studies indicate that the choice of this preconditioner improved the convergence of the system. Studies during the mold filling simulations of a geometrically complex component indicated that the number of iterations for a convergence tolerance of about 1e ¹ 06 3 initial residual increases with the increasing number of filled nodes and filled regions in the mold. An illustrative application of a large composite structure flow simulation is the modeling of the 10-ft composite keel beam section. The mesh geometry employed in the simulations is shown in Fig. 4 and the time progression of the fill contours for representative injection locations is shown in Fig. 5. For this 10-ft keel beam section mold filling simulation involving 29,171 nodal points, the change in number of iterations for the solution convergence plotted against the percentage of filled nodes of the mold is shown in Fig. 6. Effective data structures and solution strategies have been implemented to develop an RTM flow process simulation software code on the massively parallel computing platform of CM-5. The implementations have been verified, validated and benchmarked first with simple geometries where analytical solutions are available for code debugging. Similar solutions were obtained for the same physical problem and the physical conditions providing the confidence with the computer implementations and the element level residual iterative solver. The accuracy and behavior of the solutions obtained from the CM-5 implementations are similar and are not described here.15–17 For complex geometries involving moderately small problem sizes for the CM-5 environment, the solution behavior and the computational advantages, when comparing the explicit finite element–control volume methodologies and the pure finite element methodologies, are similar to those obtained earlier.15–17 Illustrative large-scale applications and the performance issues are discussed next.
Fig. 4. 10-ft keel beam section: computational mesh.
Advanced manufacturing of large-scale composite structures
257
Fig. 5. 10-ft keel beam section: fill pattern.
7 ILLUSTRATIVE APPLICATIONS Illustrative examples are presented demonstrating the applications of the process simulations to large-scale composite structures. First, a simple illustrative example demonstrating the physical accuracy of the pure finite element methodology is presented. The numerical front solutions are compared with the analytical solution in the case of a simple circular plate geometry. Subsequently, large-scale practical applications are demonstrated.
forming the inlet gate of the mold is considered. The inner radius of the plate at radius R0 is subjected to a constant pressure of P, and the flow is driven by this injection pressure. The cavity is of thickness H, and the permeability of the fiber preform is taken to be K, with a porosity of fiber compaction denoted by F. Assuming a perfect complete radial flow, the time taken by the flow front to reach a radius R f is given by
t¼ 7.1 Circular plate with constant injection pressure at the center—comparison with analytical solution A circular plate with a hole of radius R0 at the center
HF KP
ZRf r ln r R0 dr R0 ZH 1 dz h
(14)
0
For constant viscosity h0 , the time to reach a radius of R f is
Fig. 6. Solver iteration variations: 10-ft keel beam section.
258
R. V. Mohan et al.
Fig. 7. Circular plate geometry; constant pressure injection. (a) Mold plate mesh geometry; (b) fill contours in the mold.
Advanced manufacturing of large-scale composite structures
259
compared with both the present pure finite element formulations and the traditional explicit finite–element– control volume formulations as shown in Fig. 8(a). From the figure, it is clear that the pure finite element methodology has a better resolution of the flow fronts compared to the explicit finite element–control volume formulations. This clearly demonstrates the physical advantage of the pure finite element methodology in preserving the mass conservation. The comparison of flow front locations computed using different time step sizes employing the pure finite formulation is shown in Fig. 8(b). The numerical simulations are based on employing a two-dimensional quadrilateral mesh of 600 elements, with a twodimensional isotropic permeability tensor, and an outer radius taken to be 0.1 m. The results compare well with the analytical solution which is based on a perfect radial flow. The pure finite element method is an implicit approach and computes the flow front locations at a given time with good accuracy. The flow front location does not depend on the time step size employed, and the computed resin front advancement and impregnated region at any given time remain the same irrespective of the time step size employed to reach that stage. The numerical accuracy of the pressure solution, however, depends on the time step size. Improved flow front progression resolution can be obtained with smaller time step sizes. 7.2 Computational costs and computational methodologies
Fig. 8. Circular plate mold-flow comparisons: constant pressure injection. (a) Comparison with explicit FE/CV; (b) pure FE comparisons.
given by t¼
h0 F R 2 Rf R2 R2 ln ¹ fþ 0 R0 4 4 KP 2
(15)
A circular plate model with the inner radius R0 ¼ 0:0015 m, with permeability K ¼ 44:0E ¹ 12 m 2, porosity F ¼ 0.805, viscosity m ¼ 0.02 poise, and subjected to an injection pressure of 69,000 Pa is considered. The geometric mesh employed and the flow front contours are shown in Fig. 7. The flow front as computed by the analytical solution is
Computational cost and time of large-scale simulations are major factors for effective use of massively parallel computing platforms. The computing platforms provide the required computational power and the memory to perform large-scale simulations. However, the limited availability of these computing systems and the demand for their computational resources make it imperative to have improved computational methodologies and algorithms for physical problems for the effective use of these massively parallel computing platforms. The pure finite element methodology has been shown to be computationally efficient compared to the traditional finite element–control volume methodology. This is more evident in large-scale RTM flow simulations. The finer resolved meshes employed to model the different physical regions and the geometrical complexities restrict the time increments that can be employed thereby increasing the number of quasi-steady steps for complete filling during the analysis. The pure finite element method does not have such restrictions and the computational advantage is very clearly seen even while employing time step increments providing good accurate resolution of the resin time progression and distribution. For the 10-ft keel beam section involving 29,171 nodes and 58,187 elements the computational time ratios on a 64 processor node partition are shown in Table 1. The computational advantage of the pure finite element–control volume methodology is clearly seen from the ratio of the computational times presented.
260
R. V. Mohan et al.
Table 1. Comparison of computational times: 10-ft keel beam section
Table 2. Mesh configurations for 24-ft beam example Mesh geometry
Method Explicit FE–CV Pure implicit FE (Dt ¼ 0.5 s)
Computational time* 31.28 1.0
*Relative to the actual computational time corresponding to Dt ¼ 0.4 s.
The number of iterations taken by the explicit finite element–control volume methodology is 20,567 compared to the 363 iterations taken by the pure finite element methodology for the time step size of 0.5 s. These numbers clearly demonstrate the effectiveness of the pure finite element methodology in large-scale RTM simulations and the need for effective computational methodologies while employing massively parallel computing platforms. The mesh geometry employed in the simulations is shown in Fig. 4 and the time progression of the fill contours or representative injection locations are shown in Fig. 5. In the case of other very large-scale flow simulations involving large-scale complex geometric models and finite element meshes, whose performance results are presented next, our experiences indicated that the time step increments were very small when employing the finite element–control volume method and the complete filling will definitely require the use of large computational resources. Since the computational performance and the advantage of the pure finite element methodology have been very well demonstrated, explicit finite element–control volume methodology calculations were not attempted for these large-scale models.
8 COMPUTATIONAL SCALABILITY AND LARGESCALE RTM SIMULATIONS A geometrical complex large-scale composite structure of a
Mesh Mesh Mesh Mesh Mesh
1 2 3 4 5
Number of nodes
Number of elements
45 547 135 492 297 576 405 327 892 332
89 945 269 835 594 756 809 505 1 784 268
24-ft keel beam is also considered for RTM flow process simulation studies. The geometric model is meshed with unstructured finite elements. For computational scalability and performance studies, five different unstructured mesh sizes involving increasing number of elements and nodes are considered for the same physical geometry. The flow progression contours for a representative injection condition for this 24-ft keel beam geometry are shown in Fig. 9. The details of the five different mesh configurations employed are shown in Table 2. The mesh geometries for the thin section process simulations involved two-dimensional triangular elements built in a three-dimensional space. The process simulations involve the computation of the element local permeability values which are computationally intensive. The kernels involve arithmetic operations as well as mathematical trigonometric intrinsic functions making the kernel to be highly computational intensive. The computational speed and processing power provided by the massively parallel computing platforms such as the CM-5 can be understood by considering the computing time and their behavior. The computational time was based on a 512 CM-5 partition (this partition allows dedicated processing thereby reducing the influence of other work loads on these processors) for this kernel involving the computation of element level local permeabilities. The computational operations are based on element level data parallel objects and the number of elements employed is used as the measure of the problem size. The timing comparisons are shown in Table 3.
Fig. 9. Temporal flow progression contours in 24-ft keel beam section.
Advanced manufacturing of large-scale composite structures Table 3. Computational performance and problem size: based on local permeability kernel Problem size
Time (s)
Size ratio
Time ratio
269 835 809 505 1 784 268
3.328 10.0605 24.040
1.0 3.0 6.61
1.0 3.022 7.22
The variation of the size and time ratios when the same number of processing nodes are employed for various problem sizes clearly demonstrates the computational scalability of the data structures, parallel programming model and the hardware architecture employed. Experience indicates that a significant amount of time during the process simulations is spent on the iterative solution of the linear system arising from the unstructured finite element mesh computations. The conjugate gradient algorithm theoretically converges in N iterations, where N is the order of the system. The number of iterations for convergence were very much less than the order of the system involved. An analysis of the time break out for one pass of the conjugate gradient iterative solver for the mesh 5 configuration (a highly unstructured mesh) indicated that nearly 99 per cent of the time is spent in the two communication operations. Thus the communication time dominates the computational time of the iterative solver. The actual communication time, however, is within the bounds of the physical machine configuration. An improvement in communication time can be obtained through a reorder mesh connectivity based on the machine architecture where the dual mesh connectivity lays out the elements in nearby processors. However, the current version of the TMC–CMSSL library routine generate_dual has the restriction that when two-dimensional elements such as the triangular and quadrilateral elements are involved, no more than two elements can share a single edge. In our case involving twodimensional thin shell elements built in a three-dimensional space, more than two finite elements can share a common
261
edge resulting in an error code returned by the TMC– CMSSL routine. A newer version of the routine may overcome the problem and may be developed in close collaboration with the original developers. One of the main time-consuming parts in large-scale RTM flow simulations is the generation of large-scale meshes. Currently, the mesh geometry for RTM process simulations are created based on CAD geometry model of the RTM part based on commercial programs such as ProEngineer. These parametric-based mesh generation programs normally work on graphical workstation platforms and are limited by their computing power and memory. In our computations involving large-scale meshes, we obtained the maximum possible refined mesh configurations from these mesh generation software tools and obtained fine large-scale meshes by the refinement of these original meshes. The process is time-consuming and tedious. This clearly shows a need for good parametric-based CAD geometry automatic mesh generators which would significantly enhance the time involved in performing large-scale simulations. Likewise, an improvement in the development and availability of large-scale data visualization tools improves the time involved in performing large-scale simulations. All these clearly demonstrate the need for collaborative research and development efforts involved for large-scale simulations employing massively parallel computing platforms. These computing platforms make the large-scale complex simulations a reality while requiring careful development to reach that stage.
9 COMPOSITE WING STUB BOX An illustrative NASA application involving a new generation aircraft wing schematic is shown in Fig. 2. The wing structure has a composite wing stub box which is to be potentially manufactured by resin transfer molding. For purposes of illustration, process simulations were performed
Fig. 10. Temporal flow progression contours in composite wing stub box.
262
R. V. Mohan et al.
on a representative wing box geometry based on a thin shell geometry. The net-shape fabrication of the wing box involve laying and pressing of the fiber to the required preform shape configurations. The bending of the fibers around the edges cause deformation and permeability variations leading to race tracking effects along the edges and are represented and taken into account by the appropriate increase in the permeability in these corner regions during the process simulations. The fill progression contours for illustrative injections at the narrower end are shown in Fig. 10. The effect of race tracking along the edges and the variations in flow progression due to these race tracking edges are clearly seen.
Carolina, under grant number DAAH04-96-1-0172. Special thanks are due to Mr Walter Roy, Dr Shawn Walsh and Dr Dennis Viechnicki of the Materials Division of the US Army Research Laboratory, MD. Thanks are also due to Mr Bill Mermagen of the US Army Research Laboratory and the IMT activities at ARL for the encouragement and support. Special thanks are also due to Mr C. Nietubicz and ARL/MSRC computing facilities. Additional support in the form of computer grants from the Minnesota Supercomputer Institute (MSI), University of Minnesota, is also acknowledged.
REFERENCES 10 CONCLUDING REMARKS High-performance computing is yet to make a definite impact in the area of advanced manufacturing which is vital to US industry needs. The use of high-performance and massively parallel computing in manufacturing process simulations is very limited. Manufacturing processes such as the resin transfer molding with its capacity to manufacture large, geometrically complex composite parts justify the use of the high-performance computing resources to obtain optimal process solutions. Large-scale process simulations need computational algorithms which are physically accurate and permit faster solution of the computational physical problem. The success of the large-scale composite process simulations not only depends on the optimal data structures and communication strategies that are employed for the massively parallel computing architectures, but also on the development and enhancement of the computational algorithms and methodologies employed for the solution of the physical problem. In this framework, an effective computational methodology for the resin impregnation in the composite manufacturing process by resin transfer molding is presented and the effectiveness is demonstrated for large-scale process simulations. Implementation, solution strategies, communication issues have been discussed with reference to the massively parallel computing architecture of Connection Machine CM-5 that was available for our use. Our experiences and observations of the computational behavior of the linear system solvers during the large-scale composite manufacturing process simulations are discussed. Practical illustrative applications have been presented. The computational scalability of the data structures, parallel programming models and hardware architecture employed have been illustrated with large-scale computational problem sizes.
ACKNOWLEDGEMENTS The authors from the University of Minnesota are very pleased to acknowledge the support from the US Army Research Office (ARO), Research Triangle Park, North
1. Plaskacz, E. J., Kennedy, J. M., Belyschko, T. and Greenwell, D. L., Finite element analysis on the Connection Machine. Comp. Meth. Appl. Mech. Eng., 1990, 81, 229–254. 2. Johnson, S. L. and Mathur, K. K., Data structures and algorithms for the finite element method on a data parallel supercomputer. Int. J. Num. Meth. Eng., 1990, 29, 881–908. 3. Mathur, K. K. and Johnson, S. L., The finite element method on a data parallel computing system. Int. J. High Speed Computing, 1989, 1, 29–44. 4. Farhat, C., Fezoui, L. and Lanteri, S., Computational fluid dynamics with irregular grids on the Connection Machine. Comp. Meth. Appl. Mech. Eng., 1992. 5. Johan, Z., Hughes, T. J. R., Mathur, K. K. and Johnsson, S. L., A data parallel finite element method for computational fluid dynamics on the Connection Machine system. Computer Methods in Applied Mechanics and Engineering, 1992, 99, 113–134. 6. Mohan, R. V., Avila, A. F., Tamma, K. K. and Namburu, R. R., Three-dimensional transient thermal analysis with explicit finite element representations—parallel implementations and performance studies: the Connection Machine (CM-5). Numerical Heat Transfer, Part B, 1995, 28, 277–291. 7. Namburu, R. R., Mohan, R. V., Rostam Abadi, F., Beck, R. and Garcia, R., Partitioned and adaptive self-starting explicit strategies for CSD: data parallel implementation. In AIAA 34th Structural Dynamics, Materials Conference, La Jolla, April 1993. 8. Namburu, R. R., Turner, D. A. and Tamma, K. K., Effective data parallel self-starting explicit methodology for computational structural dynamics on the connection machine CM-5. International Journal for Numerical Methods in Engineering, 1995, 38(19), 3211–3226. 9. Darcy, H., Les Fontaines Publiques de la Ville de Dijon. Delmont, Paris, 1856. 10. Fracchia, C. A., Numerical simulation of resin transfer mold filling. Master’s thesis, University of Illinois at UrbanaChampaign, IL, 1990. 11. Tucker, C. L., Fracchia, C. A. and Castro, J., A finite element/control volume simulation of resin transfer mold filling. In Proc. of the American Society for Composites, 4th Technical Conference, Lancaster, PA, 1957. 12. Bruschke, M. V. and Advani, S. G., A numerical simulation of the resin transfer mold filling process. In ANTEC, 1990. 13. Bruschke, M. V. and Advani, S. G., A finite element/control volume approach to mold filling in anisotropic porous media. Poly. Comp., 1990, 11(6). 14. Trouchu, F., Gauvin, R. and Gao, D. M., Numerical analysis of the resin transfer molding process by the finite element method. Advances in Polymer Technology, 1993, 12(4), 329– 342.
Advanced manufacturing of large-scale composite structures 15. Mohan, R. V., Ngo, N. D., Tamma, K. K. and Fickie, K. D., On a pure finite element based methodology for resin transfer mold filling simulations. In Numerical Methods for Thermal Problems, vol. IX, ed. R. W. Lewis and P. Durbetaki, Atlanta, GA, Pineridge Press, July 1995, pp. 1287–1310. 16. Mohan, R. V., Ngo, N. D., Tamma, K. K. and Fickie, K. D., A pure finite element based methodology for resin transfer mold filling simulations. Technical Report ARL-TR-975, US Army Research Laboratory, March 1996. 17. Mohan, R. V., Ngo, N. D., Tamma, K. K. and Fickie, K. D., Process modeling and implicit tracking of moving fronts for three-dimensional thick composites manufacturing. In AIAA96-0725, 34th Aerospace Sciences Meeting, Reno, NV, January 1996. 18. Hwan, W.-S. and Stoehr, R. A., Molten metal flow pattern prediction for complete solidification analysis of near net shape castings. Materials Science and Technology, 1988, 4, 240–250. 19. Dhatt, G., Gao, D. M. and Ben Cheikh, A., A finite element simulation of mold flow in moulds. Int. J. Num. Meth. Eng., 1990, 30, 821–831. 20. Hirt, C. W. and Nichols, B. D., Volume of Fluid method for the dynamics of free boundaries. J. Computational Physics, 1981, 39, 201–225. 21. Stoehr, R. A. and Wang, C., Coupled heat transfer and fluid flow in the filling of castings. AFS Transactions, March 1988, 733–740. 22. Swaminathan, C. R. and Voller, V. R., A time-implicit filling algorithm. Applied Math. Modelling, 1994, 18, 101–108. 23. TMC. CMSSL For CM Fortran, 1992. 24. TMC. CMSSL Release Notes for the CM-5, version 3.0 beta edition, 1992. 25. Thinking Machines Corporation, 245, First Street, Cambridge, MA 02142. The Connection Machine CM-5 Technical Summary, 1991. 26. TMC. CM Fortran Manual, versions 1.0 and 1.1 edition, 1991.
263
27. Mohan, R. V., Avila, A. F. and Tamma, K. K., Three-dimensional transient thermal analysis with explicit finite element representations—parallel implementations and performance studies: the Connection Machine (CM-5). In AIAA 29th Thermophysics Conference, Colarado Springs, June 1994. 28. Namburu, R. R., Turner, D. A. and Tamma, K. K., A data parallel self-starting explicit method for computational structural dynamics on the connection machine CM-5. In 1992 ASME Winter Annual Meeting, Anaheim, CA, 1992. 29. Namburu, R. R. and Rostam-Abadi, F., A data parallel finite element explicit method for computational heat transfer. In 31st Aerospace Sciences Meeting and Exhibit, Reno, NV, 1993. 30. Mohan, R. V., Namburu Raju R. and Tamma, K. K., A data parallel explicit subcycling time integration strategy for computational heat transfer. In ASME Winter Annual Meeting, December 1993. 31. Farhat, C. and Wilson, E., A new finite element concurrent computer program architecture. Int. J. Num. Meth. Eng., 1987, 24, 1771–1792. 32. Dahl, E. D., Mapping and compiled communication on the Connection Machine system. In Proc. 5th Distributed Memory Computing Conference, Los Alamitos, CA, IEEE Computer Society Press, 1990. 33. Saad, Y. and Schultz, M. H., GMRES: A generalized minimal residual algorithm for solving nonsymmetric linear systems. SIAM J. Sci. Statis. Comput., 1986, 7, 856–869. 34. Carey, G. F., Barragy, E., Mclay, R. and Sharma, M., Element-by-element vector and parallel computations. Comm. Appl. Num. Methods, 1988, 4, 299–308. 35. Carey, G. F., Parallelism in finite element modeling. Comm. Appl. Num. Methods, 1986, 2, 281–287. 36. Wang, J. T. and Ransom, J. B., Application of interface technology in nonlinear analysis of a stitched/RFI composite wing stub box. AIAA Paper 97–1190, 38th AIAA/ASME/ ASCI/AHS/ASC Structures, Structural Dynamics and Materials Conference, 1997, Part 3, 2295–2310.