CN102184298A - Method for storing and generating stiffness matrix for finite-element analysis in metal bulk plastic forming - Google Patents
Method for storing and generating stiffness matrix for finite-element analysis in metal bulk plastic forming Download PDFInfo
- Publication number
- CN102184298A CN102184298A CN2011101300739A CN201110130073A CN102184298A CN 102184298 A CN102184298 A CN 102184298A CN 2011101300739 A CN2011101300739 A CN 2011101300739A CN 201110130073 A CN201110130073 A CN 201110130073A CN 102184298 A CN102184298 A CN 102184298A
- Authority
- CN
- China
- Prior art keywords
- grid model
- stiffness matrix
- matrix
- node
- zero submatrices
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Landscapes
- Complex Calculations (AREA)
Abstract
The invention relates to a method for storing and generating a stiffness matrix for finite-element analysis in metal bulk plastic forming. The method capable of remarkably saving a storage space of the stiffness matrix and improving the solving efficiency of a stiffness equation comprises the following steps of: (1.1) determining a generalized adjacent node pair of a grid node and a correspondence relationship between the generalized adjacent node pair and non-zero sub-matrixes in a total stiffness matrix of a grid model according to a finite-element grid model of a bulk plastic forming workpiece; (1.2) determining the positions and quantity of the non-zero sub-matrixes in the total stiffness matrix of the grid model based on the generalized adjacent node pair, and storing the non-zero sub-matrixes; (1.3) determining the non-zero sub-matrix of the total stiffness matrix in which each unit stiffness matrix element of the grid model is located and the position of each unit stiffness matrix element in a total stiffness matrix storage array, and assembling each unit stiffness matrix element into the total stiffness matrix storage array of the grid model; and (1.4) solving a total stiffness equation of the grid model by using a conjugate gradient method for relaxation pre-treatment of symmetrical positive definite coefficient matrixes.
Description
Technical field
The present invention relates in the metal volume Plastic Forming storage of finite element analysis stiffness matrix and generation method, be suitable for forging, the finite element numerical simulation of extruding, rolling, swaging equal-volume forming process.
Background technology
In commercial production, forging, extruding, metal volume Plastic Forming manufacturing technology such as rolling support one of major technique of the national economic development and national defense construction, have extensive use in industries such as automobile, space flight and aviation, equipment manufacturing, weapons, the energy, shipbuildings.Plastic yield takes place with the metal material (blank) of (cold, warm, thermoforming) under room temperature or the high temperature in volume plastic forming process exactly under forging equipment and the mould action that is installed on it thereof, and be configured as the process of the part of required form and performance, volume Plastic Forming material generally is assumed to be and just (glues) plastic material, therefore, big and violent plastic yield takes place in blank in mold cavity, have material nonlinearity, geometrical non-linearity and complicated boundary condition, finding the solution of forming process is very complicated.Can forming technology and mould design for the quality of drip molding and be shaped smoothly has material impact, also is the key that influences energy-saving and cost-reducing and Accurate Shaping.And grasp plastic yield behavior and the flowing law of metal material in mold cavity is the key theory foundation of design forming technology and mould.Finite Element Method is as the important method of engineering problem numerical analysis, brought into play vital role in the engineering field, metal volume plastic forming process numerical simulation method is the main method of research metal material plastic yield behavior and flowing law, can obtain metal material flowing law and various field variable distribution results, the feasibility of checking technological design, parameter selection and die design schemes, and then the optimization optimization of realization forming technology and mould, improve speed, quality and the reduction cost of development of product development.
Finite element analysis for metal volume Plastic Forming problem, at first by the plastic yield governing equation, set up the mathematical model of metal volume plastic yield mechanics problem, geometric model (mould and sex change material) is carried out finite element method to disperse, set up the finite element stiffness equations, realize the numerical simulation of forming process, numerical analysis methods such as finite element have particularly been brought into play vital role in fields such as metal forming manufacturing and structure analyses in the engineering field.
The step that finite element numerical is found the solution metal volume Plastic Forming problem is: the CAD geometric model (geometric configuration of deformable material) of setting up problem, the CAD geometric model is carried out finite elements to disperse, obtain the finite element grid model of geometric model, this grid model is made up of the grid cell of limited quantity; At first set up the stiffness equations of each unit in the grid model, the stiffness equations of all unit is assembled, form non-linear global stiffness system of equations about all node speeds of grid model; By to linearization process, obtain linear rigidity system of equations about the node speed increment about the non-linear rigidity equation of node speed; Adopt the Newton-Raphson process of iteration that velocity field is carried out iterative.Because the node speed increment is found the solution, therefore, for a metal volume plastic forming process, need be decomposed into curring time up to a hundred even a hundreds of time step (length), in each time step, all need computing unit stiffness matrix and nodal force vector, be assembled into the global stiffness equation of node speed increment, apply velocity boundary conditions, the global stiffness equation of solution node speed increment (large-scale system of linear equations), obtain the speed correction, with the speed correction velocity field is revised again, repeat above-mentioned steps then, velocity field is carried out iterative, restrain up to velocity field, generally speaking, need 7 to 10 iteration can obtain the convergent velocity field at least; According to the speed of convergence field that obtains, refresh distortion workpiece shape and die location, just can finish the numerical evaluation in the time step.And then carry out the interior numerical evaluation of next step-length, until reaching desired deformation extent.
This shows,, need the frequent global stiffness equation (large-scale system of linear equations) of finding the solution for the numerical simulation of metal volume plastic forming process.For example for the 3-dimensional metal volume Plastic Forming problem of medium complexity, generally need the simulation step number more than 200 times, each roughly needs 9 velocity fields of iteration in the simulation step, and each velocity field iteration just need be found the solution once large-scale system of linear equations, if the finite element interstitial content of a metal volume Plastic Forming problem is 50,000 (three-dimensional problems of medium complexity), just need finds the solution the large-scale linearity side that constitutes by 150,000 equations and organize 1800 times.If analyze more challenge, its amount of calculation is huger.
In recent years, needs along with heavy constructions such as great equipment, boats and ships, automobile and wind-powered electricity generations, increasing to demand large-scale, heavy duty metal volume Plastic Forming part, need the shaping problem of analysis complicated day by day, physical dimension is big, the node and the element number of finite element model are many, data volume is huge, the FEM (finite element) calculation inefficiency, required calculator memory increases severely, even when analyzing large-scale and ultra-large type problem, the problem of memory space inadequate often appears, this has limited the application of finite element numerical method in METHOD IN METAL FORMING PROCESSES is analyzed greatly.
Metal volume Plastic Forming problem not only, and much the finite element analysis of other engineering problem finally all can be summed up as and finds the solution large-scale system of linear equations, finite element find the solution efficient depend primarily on the system of linear equations stiffness matrix storage, generation and stiffness equations thereof find the solution efficient.
Metal volume Plastic Forming problem finite element matrix is large-scale, banded sparse matrix, nonzeros generally concentrates on the diagonal line zone, at this situation with in order to improve counting yield, people have proposed the storage means of many conserve storage, for example, the half-band width storage means, the variable bandwidth storage means, sparse system of equations pretreatment iteration method for solving etc., although these methods have reduced the storage space of stiffness matrix, but stored the part nonzero element also unavoidablely, wasted the part storage area, these methods all depend on node serial number simultaneously, need be optimized the node serial number order, thereby also waste computing time.Existing one-dimension variable bandwidth compression and storage method adopts one-dimension array to preserve the nonzero element of stiffness matrix by the order of row, adopt two supplementary number group records row number and the initial element of every row position in the stiffness matrix array, saved storage space, avoided the node numbering optimization problem, but auxiliary array that is adopted and stiffness matrix storage array have same length, still have the problem that consumes storage space and waste computing time.
Therefore, during engineering problems such as employing finite element analysis metal volume Plastic Forming, the efficient compression storage and the generation method of finite element matrix element become the problem that presses for solution.
Summary of the invention
The objective of the invention is the big and not high problem of stiffness equations counting yield of stiffness matrix data space for existing in the metal volume Plastic Forming problem finite element numerical simulation process, storage of finite element analysis stiffness matrix and generation method in a kind of new metal volume Plastic Forming are provided, it not only has remarkable saving stiffness matrix storage space and improves the characteristics of counting yield, and having a increase with metal volume Plastic Forming problem finite element grid model node and element number, its conserve space is remarkable more and find the solution the high more characteristics of efficient.
The present invention realizes by following technical scheme:
Storage of finite element analysis stiffness matrix and generation method in a kind of metal volume Plastic Forming may further comprise the steps:
(1.1) at first adopt the counter 3-D geometric models of asking equipment to obtain metal volume Plastic Forming problem such as existing CAD design software or three-coordinates measuring machine, adopt existing grid generation method to generate the finite element grid model of Forming Workpiece, next is according to the finite element grid model topology structure of Plastic Forming workpiece, determine " the broad sense adjacent node to " of its grid node, set up the corresponding relation between " broad sense adjacent node to " and the workpiece finite element grid model global stiffness matrix non-zero submatrices;
(1.2) based on " broad sense adjacent node to ", store non-zero submatrices, opening up non-zero submatrices row number storage array J is the integer array of indexes, the row of the global stiffness matrix non-zero submatrices that record is corresponding with " broad sense adjacent node to " in the workpiece finite element grid model number, open up another integer array of indexes K, the position of the first non-zero submatrices of every row in array J in the record workpiece finite element grid rigidity of model matrix, determine the position and the quantity of non-zero submatrices in the workpiece finite element grid rigidity of model matrix, according to the corresponding relation of non-zero submatrices and nonzero element, determine the position and the quantity of workpiece finite element grid rigidity of model matrix nonzero element;
(1.3) non-zero submatrices at each element stiffness matrix element place in workpiece finite element grid model global stiffness matrix in the calculating workpiece finite element grid model, utilize the data among array of indexes J and the K, the position (concrete determining step see 4.1～4.9) of determining unit stiffness matrix element in grid model global stiffness matrix stores array A; According to the create-rule of grid model global stiffness matrix, each element stiffness matrix element of grid model is assembled among the global stiffness matrix stores array A of grid model with sitting in the right seat; Behind all unit in the traversal grid model, finally generate the global stiffness matrix stores array of metal volume shaping problem grid model;
(1.4) adopt the lax pre-service of existing symmetric positive definite matrix of coefficients volume gradient method altogether, global stiffness equation (large-scale system of linear equations) to grid model is found the solution, and obtains the numerical solution of the speed increment field of all nodes in the metal volume forming process grid model.
In the described step (1.1), determine that the method for " the broad sense adjacent node to " of grid node is:
(2.1) in the finite element grid model of bulk forming problem, the node that any two nodes in each unit constitute is to the node of node self and self formation (comprise to), is called " broad sense adjacent node to "; With the sequential node of adjacent node to being called " in order broad sense adjacent node to "; Otherwise, be called " unordered broad sense adjacent node to ", (I J), then constitutes two " broad sense adjacent node to " in order if the global node of certain two node in unit is numbered in the grid model
With
If do not consider the order of node, then constitute " unordered broad sense adjacent node a to " B
IJ, promptly think
With
Be that same node is right, be designated as B
IJOr B
JIFor the convenience of using and representing, can take the ascending form of subscript to represent " unordered broad sense adjacent node to ", the promptly unified B that is designated as
IJ(I≤J);
(2.2) in grid model, in the global stiffness matrix of each " in order broad sense adjacent node to " corresponding grid model one with the identical non-zero submatrices of institute's problem analysis dimension, the position of non-zero submatrices is by the global node sequence number decision of " in order broad sense adjacent node to ", for example, " in order broad sense adjacent node to " in grid model
The non-zero submatrices that I is capable in the global stiffness matrix of corresponding grid model, J is listed as;
(2.3) in grid model, in the upper triangular matrix of the global stiffness matrix of each " unordered broad sense adjacent node to " corresponding grid model one with the identical non-zero submatrices of institute's problem analysis dimension, the position of non-zero submatrices is by the global node sequence number decision of " unordered broad sense adjacent node to "; " unordered broad sense adjacent node to " B in grid model for example
IJ(I≤J), the non-zero submatrices that I is capable in the upper triangular matrix of the global stiffness matrix of corresponding grid model, J is listed as;
(2.4) because the global stiffness matrix of metal volume plastic forming process finite element analysis is a symmetric matrix, for conserve space can only be stored its upper triangular matrix, therefore the global stiffness matrix stores of the grid model of mentioning among the present invention is meant the storage of its upper triangular matrix, non-zero submatrices because in the grid model in the upper triangular matrix of " unordered broad sense adjacent node to " and the global stiffness matrix of grid model is corresponding one by one, and the row of the right global node of node numbering and non-zero submatrices number, the row correspondence, therefore by this corresponding relation, number can determine the quantity of non-zero submatrices in the upper triangular matrix of grid model stiffness matrix and position (row number with row number) by the quantity of all in the grid model " unordered broad sense adjacent nodes to " and global node, and can further determine the quantity and the position of the nonzero element of composition non-zero submatrices.
In the described step (1.2), further comprising the steps of:
(3.1) adopt a floating type array A by each nonzero element in the row sequential storage grid model global stiffness matrix upper triangular matrix;
(3.2) adopt the row number of the non-zero submatrices in shaping array J record and " unordered broad sense adjacent node to " corresponding global stiffness matrix upper triangular matrix;
(3.3) adopt the numbering of first non-zero submatrices of every row in array J in the shaping array K record global stiffness matrix upper triangular matrix, and determine that with this non-zero submatrices is at the row of global stiffness matrix upper triangular matrix number, because K[i] be the numbering of capable first non-zero submatrices of i in array J in the global stiffness matrix upper triangular matrix, and K[i+1] be the numbering of capable first non-zero submatrices of i+1 in array J in the global stiffness matrix upper triangular matrix, numbering that hence one can see that in array J is from K[i] to K[i+1] the sub-square of all non-zeros of-1 correspondence is all capable at i.
Finite element grid model after the geometric model of institute's problem analysis is divided be the two a mixed cell model of the hexahedral element, tetrahedron element of three-dimensional or its; Finite element grid model after perhaps the geometric model of institute's problem analysis is divided be the two a mixed cell model of two-dimentional quadrilateral units, triangular element or its, can adopt method of the present invention to carry out the storage of grid model global stiffness matrix.
In the described step (1.3),, be m, remember that the integral body of n the node in this unit is numbered I[i] if the node number of unit is the dimension of n, physical quantity for a unit in the grid model, i=1,2 ..., n, this element stiffness matrix is that (m * n, m * n), then the arbitrary element of this element stiffness matrix is U
U
ijkl＝U[(i-1)×m+k][(j-1)×m+l]
Wherein, i, j=1,2 ..., n is the part numbering of this cell node, k, and l=1,2 ..., m is physical quantity dimension numbering;
Determine U
IjklThe process of sequence number is as follows in grid model global stiffness matrix stores array A:
(4.1) determine U
IjklIntegral node numbering: I=I[i], J=I[j], wherein I, J are U
IjklThe non-zero submatrices A at place
IJIn row in grid model global stiffness matrix upper triangular matrix number and row number; The position of non-zero submatrices as shown in Figure 1.
(4.2) determine non-zero submatrices A
IIThe quantity of non-zero submatrices: N before
I=K[I]-1;
(4.3) determine non-zero submatrices A
IIThe quantity of nonzero element before:
NN
I＝N
I×m
2-(I-1)×m×(m-1)/2
Wherein, m
2Be the quantity of each non-zero submatrices nonzero element, the minimizing item of the nonzero element that m (m-1)/2 causes for the diagonal angle non-zero submatrices of every row; As the non-zero submatrices A among Fig. 2
IIBe the three-dimensional submatrix of m=3, only stored 6 nonzero elements, compare with other submatrix and deposit 3 less.
(4.4) determine the capable quantity of nonzero element before of k in the capable non-zero submatrices of I:
NN
IK＝(K[I+1]-K[I])×m×(k-1)
(4.5) determine that I is capable, the quantity of the non-zero submatrices before the J row.Search array J is from label K[I] to K[I+1] J[II]=the array label II of J, during I was capable, it was N that J is listed as non-zero submatrices quantity before
J=II-K[I];
(4.6) determine U
IjklNumbering in the global stiffness matrix, being expert at:
NN
Jl＝S
J×m-TT+l
Wherein, TT is for because the minimizing item of the nonzero element that causes of diagonal angle non-zero submatrices, TT=k (k-1)/2;
(4.7) determine U
IjklNumbering in grid model global stiffness matrix storage arrays A:
NN
IJkl＝NN
I+NN
Ik+NN
Jl
(4.8) according to the create-rule and the U of grid model global stiffness matrix
IjklPosition Number NN in grid model global stiffness matrix storage arrays A
IJkl, with U
IjklAssemble or be added among the global stiffness matrix storage arrays A of grid model;
(4.9) repeating step (4.1)～(4.8), each element in all element stiffness matrixs in the traversal grid model, the global stiffness matrix one dimension that can obtain grid model is stored array A.
Find the solution numerical analysis method that grid model adopts and finally be summed up as storage, generation and stiffness equations thereof (being system of linear equations) sparse, the symmetric population stiffness matrix and find the solution, the parameter of the global stiffness equation solution of the grid model of being analyzed is full dose or increment.
The invention has the beneficial effects as follows: propose according to the present invention " the broad sense adjacent node to " and and grid model global stiffness matrix non-zero submatrices between corresponding relation, non-zero submatrices and position thereof in the save mesh model global stiffness matrix have changed the method for nonzero element in the save mesh model global stiffness matrix in the past and position thereof; The required storage space of compression and storage method of the present invention is linear function with the increase of grid model number of nodes to be increased, being different from the required storage space of existing compression and storage method is quadratic function with the increase of grid model number of nodes and increases, therefore, method of the present invention has not only significantly been saved the storage space of grid model global stiffness matrix and has been improved operation efficiency, and with the increase of grid model node and element number, the ability of its conserve storage is remarkable more; Lax pre-service based on the grid model global stiffness matrix equation of storing non-zero submatrices is total to the volume gradient method, need not the division matrix of progressively overrelaxation iteration of symmetry is inverted and opened up array in addition and deposit the division matrix, can directly carry out computing, find the solution and the operation efficiency height the division nonzeros; Therefore, be particularly suitable for the metal volume this increment that need carry out repeatedly that is shaped and find the solution numerical analysis with the large complicated engineering problem of iterative.
Description of drawings
Fig. 1 is the distribution schematic diagram of non-zero submatrices in the global stiffness matrix;
Fig. 2 is the position view of nonzero element in the global stiffness matrix;
Fig. 3 is metal storage of (gluing) Plastic Forming grid model global stiffness matrix equation compression just and generating algorithm process flow diagram;
The grid model of Fig. 4 for constituting by one eight node hexahedral element;
The grid model of Fig. 5 for constituting by two eight node hexahedral elements;
Fig. 6 a is the forging process finite element numerical simulation;
Fig. 6 b is the forging process finite element numerical simulation;
Fig. 6 c is the forging process finite element numerical simulation;
Fig. 6 d is the forging process finite element numerical simulation.
Embodiment
Below in conjunction with accompanying drawing and embodiment, the present invention will be further described.
Fig. 3 is metal storage of (gluing) Plastic Forming grid model global stiffness matrix equation compression just and generating algorithm process flow diagram.According to shown in Figure 3, metal (glue) Plastic Forming grid model global stiffness matrix equation compression just store and the generating algorithm flow process as follows: according to the geometric model and the finite element grid model of the metal forming problem that has generated, determine " the broad sense adjacent node to " of whole grid model and calculate its quantity; Open up grid model global stiffness matrix stores array A, open up non-zero submatrices row number storage array J, open up the first non-zero submatrices location storage of every row array K; All nodes in the grid model are circulated,, comprise all unit of node i in the search grid model for arbitrary node i; Search comprises all nodes that all unit comprised of node i, and determines " the broad sense adjacent node to " of node i; According to the numbering that comprises node i " broad sense adjacent node to ", generate array of indexes J; According to the quantity that comprises node i " broad sense adjacent node to ", generate array of indexes K; Whether the circulation of judging all nodes in the grid model finishes, if also do not finish, then next node is repeated above-mentioned steps, until all node processing in the grid model are finished.All unit in the grid model are circulated, generate the stiffness matrix U of each unit, the position of each element of judging unit stiffness matrix non-zero submatrices at place in grid model global stiffness matrix upper triangular matrix, the non-zero submatrices before calculating place non-zero submatrices is expert at and the quantity of nonzero element; The quantity of all nonzero elements before the element of the grid model global stiffness matrix upper triangular matrix of computing unit stiffness matrix element correspondence is expert at; The serial number that the grid model global stiffness matrix upper triangular matrix element of computing unit stiffness matrix element correspondence is expert at; The Position Number of grid model global stiffness matrix upper triangular matrix element in grid model global stiffness matrix stores array A of computing unit stiffness matrix element correspondence, according to the create-rule of grid model global stiffness matrix, with this element stiffness matrix element grid model global stiffness matrix upper triangular matrix storage array of packing into.Judging whether all elements calculates in this cell matrix finishes, if do not finish as yet, then the next element in this unit is carried out above-mentioned processing, and all elements disposes in to this unit.Judge whether to still have the unit that does not generate stiffness matrix then, if also have, then carry out the cycle calculations of next unit, the stiffness matrix until all unit finishes according to the above-mentioned steps generation, thereby obtains the upper triangular matrix and the stiffness equations thereof of whole grid model global stiffness matrix.At last, adopt lax pre-service volume gradient method altogether, find the solution the global stiffness equation of whole grid model, the increment that obtains velocity field is separated, and utilizes this increment to separate then and refreshes the velocity field of grid model in this moment.
The quantity of " broad sense adjacent node to " is by the topological structure decision of grid model in the grid model, any two nodes that belong to same unit can be formed one " broad sense adjacent node to ", therefore " broad sense adjacent node to " quantity determines, calculate exactly in the whole grid model, the unduplicated node of same unit that belongs to is to quantity, for the grid model of different units type, determine that the principle of " broad sense adjacent node to " is identical, but concrete grammar is different.
Because non-zero submatrices is corresponding one by one in the global stiffness matrix upper triangular matrix of " broad sense adjacent node to " and grid model, can determine the quantity of non-zero submatrices in the global stiffness matrix upper triangular matrix by the quantity of " broad sense adjacent node to ", and then the quantity of definite nonzero element.Yet because the global stiffness matrix of grid model adopts the upper triangular matrix storage, so each non-zero submatrices on the diagonal line also only needs storage to go up the triangle submatrix.If establish the dimension of m for the problem of finding the solution, the nonzero element of diagonal line non-zero submatrices is m (m+1)/2, and off-diagonal non-zero submatrices data are m
2, can determine the quantity of the nonzero element of grid model global stiffness matrix upper triangular matrix thus.
The grid model that constitutes with eight node hexahedral elements by limited quantity is that example is illustrated the aforementioned calculation process.
For the grid model of being made up of eight node hexahedral elements, its " broad sense adjacent node to " is divided into three types, and the one, the limit adjacent node is right; The 2nd, the face adjacent node is right; The 3rd, the body adjacent node is right.The limit adjacent node is to by the decision of the quantity of adjacent node line in the unit, corresponding one of every line " broad sense adjacent node to "; The right quantity of face adjacent node is by the quantity decision of elemental area, and each face has two diagonal line, corresponding one of every face diagonal " broad sense adjacent node to "; Right for the body adjacent node, there are four body diagonals each unit, corresponding one of every body diagonal " broad sense adjacent node to ".If N
1Be number of nodes (node that is node and self generation is to quantity), N
2Be the quantity (being the right quantity of limit adjacent node) of adjacent node line in the grid model, N
3(each face has two face adjacent nodes right, 2N for the face number of unit in the grid model
3The quantity that the presentation surface adjacent node is right), N
4(each unit has two 4 individual adjacent nodes right, 4N for the quantity of unit in the grid model
4Expression body adjacent node right quantity), then the quantity of grid model " broad sense adjacent node to " is:
S
SM＝N
1+N
2+2×N
3+4×N
4
N wherein
1Be the right quantity of node that node self constitutes, the type node to the non-zero submatrices row of correspondence number be listed as number consistently, be positioned on the diagonal line, so N
1Quantity for the diagonal line non-zero submatrices.
For three-dimensional problem (m=3), then the data that need store of each diagonal line non-zero submatrices are m (m+1)/2=6, and the data that the off-diagonal non-zero submatrices need be stored are m
2=9.Therefore the quantity of the global stiffness matrix upper triangular matrix nonzero element of grid model is:
S
NZ＝6×N
1+9×(N
2+2×N
3+4×N
4)
Be that example is described further the aforementioned calculation process with two concrete grid models below.
For the grid model that constitutes by a hexahedral element shown in the accompanying drawing 4, N
1=8, N
2=12, N
3=6, N
4=1, when adopting the upper triangular matrix storage, its non-zero submatrices quantity is:
S
SM＝N
1+(N
2+2×N
3+4×N
4)＝8+(12+2×6+1×4)＝36
The global stiffness matrix upper triangular matrix nonzero element quantity that grid model shown in Figure 4 need be preserved is:
S
NZ＝6×N
1+9×(N
2+2×N
3+4×N
4)＝6×8+9×28＝300
In like manner, for the grid model that constitutes by two hexahedral elements shown in the accompanying drawing 5, N
1=12, N
2=20, N
3=11, N
4=2.When adopting the upper triangular matrix storage, its non-zero submatrices quantity is:
S
SM＝N
1+(N
2+2×N
3+4×N
4)＝12+(20+2×11+2×4)＝62
The global stiffness matrix upper triangular matrix nonzero element quantity that grid model shown in Figure 5 need be preserved is:
S
NZ＝6×N
1+9×(N
2+2×N
3+4×N
4)＝6×12+9×50＝522
Determined the quantity of non-zero submatrices and nonzero element, the storage of the global stiffness matrix upper triangular matrix of the grid model of the storage space of the global stiffness matrix upper triangular matrix of distribution network lattice model, and then realization exactly.
Global stiffness matrix with the grid model that is made of two hexahedral elements shown in Figure 5 is an example below, and the compression storage advantage of global stiffness matrix of the present invention is described.
The data length of the J array in the global stiffness matrix storage method involved in the present invention is the quantity of global stiffness matrix non-zero submatrices, i.e. S
SM=62, and the data length of K array is an element number 12, then the global stiffness matrix stores data volume of this grid model is 62+12=74.If adopt CSR (Compressed Sparse Row format) form storage means commonly used at present, then the data length of its J array is the quantity of all nonzero elements, i.e. S
NZ=522, the data length of K array is the dimension of global stiffness equation, is 12 * 3=36, and therefore, total amount of data is 558.If only compare two auxiliary array J and K, global stiffness matrix storage method then involved in the present invention is compared with the CSR method, and the storage space of its saving is: (558-74)/and 558=87%.If floating type array A adopts single precision floating datum to preserve, global stiffness matrix storage method then involved in the present invention is compared with the CSR method, and the space that the global stiffness matrix stores is saved is: 1-(74+552 * 2)/(558+552 * 2)=29%.As seen, for example shown in Figure 5, the storage space that the present invention saves auxiliary array J and K reaches 87%, saves global stiffness matrix stores space 29%, and the space-saving effect is obvious.
Below by a concrete forging example, beneficial effect of the present invention is described further.
Accompanying drawing 6a-Fig. 6 d is the forging finite element result that a certain moment adopts different dividing elements quantity to obtain in forging process, wherein the unit of all grid models all is eight node hexahedral elements, forging grid model among Fig. 6 a contains 1565 unit and 2045 nodes, forging grid model among Fig. 6 b contains 3284 unit and 4119 nodes, forging grid model among Fig. 6 c contains 5339 unit and 6557 nodes, and the forging grid model among Fig. 6 d contains 8530 unit and 10217 nodes.From EMS memory occupation and two aspects of counting yield, the algorithm that compression storage iterative algorithm of the present invention and half band storage are directly found the solution compares below.Wherein the computing machine that is adopted is: CPU:Intel (R) Core (TM) 2, and Duo, CPU E8400,3.00GHz, internal memory: 3.25, operating system: Windows XP.
Aspect EMS memory occupation, the element of the global stiffness equation of forging grid model adopts single precision floating datum to preserve, i.e. 4 bytes of each data occupancy; Auxiliary array adopts the shaping data to preserve, i.e. 2 bytes of each data occupancy.Under different nodes and element number, the shared storage space of storing in different ways of forging grid model global stiffness matrix is compared, the result is as shown in table 1.
The comparison in table 1 forging grid model global stiffness matrix stores space under different storage meanss
By table 1 as seen, the shared storage space of compression and storage method of the present invention is far smaller than the storage space of half band storage means, when the element number of forging grid model is 1565, and conserve storage 93.2%, when the element number of forging grid model is 8530, conserve storage 97.4%.As seen, with the increase of grid model node or element number, it is increasing that storage space is saved ratio.This is main because compression and storage method of the present invention has only been preserved the global stiffness nonzeros, the quantity of each row nonzero element is only relevant with the quantity of the broad sense adjacent node of this row corresponding node in the global stiffness matrix, and irrelevant with the number of nodes of integral grid model.In other words, for compression and storage method of the present invention, storage space is only relevant with the line number increase of global stiffness matrix with the proportional relation of number of nodes, is linear trends of change.And for partly being with storage means, the storage space of each row is increased with number of nodes by the half-band width decision of global stiffness matrix, and half-band width generally also increases thereupon, and the line number of global stiffness matrix and every capable all increases to some extent of data that need storage, thereby present the para-curve variation tendency.Therefore, along with the increase of number of nodes, partly being quadratic function with the required storage space of storage means increases, and the required storage space of compression and storage method of the present invention is the linear function increase, and its space-saving advantage is more remarkable.
Aspect counting yield, be example only with the global stiffness equation single calculation efficient of finding the solution a certain moment in the forging process herein, half-and-half the band storage is directly found the solution with two kinds of methods of compression storage iterative of the present invention and is compared, and the result is shown in subordinate list 2.
Table 2 is found the solution efficiency ratio for the grid model of different units and number of nodes when different storage means
As can be seen from the table, when the element number of forging grid model was 1565, counting yield improved 69.4%, when the element number of forging grid model is 8530, improved counting yield 88.3%.As seen, the efficient that the efficient of compression storage iterative of the present invention is directly found the solution far above half band storage, and also number of nodes is many more, finds the solution efficient and improves obvious more.Therefore, the method that proposes of the present invention is particularly suitable for the grid model unit and interstitial content is big, the finite element numerical analysis of the forging of finding the solution of iterating, extruding, metal volume shaping problem such as rolling.
Claims (6)
1. storage of finite element analysis stiffness matrix and generation method in the metal volume Plastic Forming is characterized in that may further comprise the steps:
(1.1) at first adopt existing CAD software or three-coordinates measuring machine to obtain the 3-D geometric model of metal volume Plastic Forming problem; Adopt existing grid generation method to generate the finite element grid model of Forming Workpiece; Next is according to the finite element grid model topology structure of Plastic Forming workpiece, determine " the broad sense adjacent node to " of its grid node, set up the corresponding relation between " broad sense adjacent node to " and the workpiece finite element grid model global stiffness matrix non-zero submatrices;
(1.2) based on " broad sense adjacent node to ", store non-zero submatrices, opening up non-zero submatrices row number storage array J is the integer array of indexes, the row of the global stiffness matrix non-zero submatrices that record is corresponding with " broad sense adjacent node to " in the workpiece finite element grid model number, open up another integer array of indexes K, the position of the first non-zero submatrices of every row in array J in the record workpiece finite element grid rigidity of model matrix upper triangular matrix, determine the position and the quantity of non-zero submatrices in the workpiece finite element grid rigidity of model matrix, according to the corresponding relation of non-zero submatrices and nonzero element, determine the position and the quantity of workpiece finite element grid rigidity of model matrix nonzero element;
(1.3) non-zero submatrices at each element stiffness matrix element place in workpiece finite element grid model global stiffness matrix in the calculating workpiece finite element grid model, utilize the data among array of indexes J and the K, the position of determining unit stiffness matrix element in grid model global stiffness matrix stores array A; According to the create-rule of grid model global stiffness matrix, each element stiffness matrix element of grid model is assembled among the global stiffness matrix stores array A of grid model with sitting in the right seat; Behind all unit in the traversal grid model, finally generate the global stiffness matrix stores array of metal volume shaping problem grid model;
(1.4) adopt the lax pre-service of existing symmetric positive definite matrix of coefficients volume gradient method altogether, global stiffness equation to grid model is that large-scale system of linear equations is found the solution, and obtains the numerical solution of the speed increment field of all nodes in the metal volume forming process grid model.
2. in the metal volume Plastic Forming according to claim 1 finite element analysis stiffness matrix storage and generation method is characterized in that, in the described step (1.1), determine grid node " the broad sense adjacent node to " method be:
(2.1) in the finite element grid model of bulk forming problem, the node that any two nodes in each unit constitute is right, comprises that promptly the node of node self and self formation is right, is called " broad sense adjacent node to "; With the sequential node of adjacent node to being called " in order broad sense adjacent node to "; Otherwise, be called " unordered broad sense adjacent node to ", (I J), then constitutes two " broad sense adjacent node to " in order if the global node of certain two node in unit is numbered in the grid model
With
If do not consider the order of node, then constitute " unordered broad sense adjacent node a to " B
IJ, promptly think
With
Be that same node is right, be designated as B
IJOr B
JIFor the convenience of using and representing, can take the ascending form of subscript to represent " unordered broad sense adjacent node to ", the promptly unified B that is designated as
IJ(I≤J);
(2.2) in grid model, in the global stiffness matrix of each " in order broad sense adjacent node to " corresponding grid model one with the identical non-zero submatrices of institute's problem analysis dimension, the position of non-zero submatrices is determined by the global node sequence number of " broad sense adjacent node to " in order;
(2.3) in grid model, in the upper triangular matrix of the global stiffness matrix of each " unordered broad sense adjacent node to " corresponding grid model one with the identical non-zero submatrices of institute's problem analysis dimension, the position of non-zero submatrices is by the global node sequence number decision of " unordered broad sense adjacent node to ";
(2.4) the global stiffness matrix of metal volume plastic forming process finite element analysis is a symmetric matrix, only stores the triangular matrix on it; Non-zero submatrices in the grid model in the upper triangular matrix of " unordered broad sense adjacent node to " and the global stiffness matrix of grid model is corresponding one by one, and the row of the right global node of node numbering and non-zero submatrices number, be listed as number corresponding, therefore by this corresponding relation, quantity and position by non-zero submatrices in the upper triangular matrix of the quantity of all in the grid model " unordered broad sense adjacent nodes to " and global node number definite grid model stiffness matrix, at once number with row number, thereby further determine to form the quantity and the position of the nonzero element of non-zero submatrices.
3. storage of finite element analysis stiffness matrix and generation method is characterized in that in the metal volume Plastic Forming according to claim 1, and be in the described step (1.2), further comprising the steps of:
(3.1) adopt a floating type array A by each nonzero element in the row sequential storage grid model global stiffness matrix upper triangular matrix;
(3.2) adopt the row number of the sub-square of non-zero in shaping array J record and " unordered broad sense adjacent node to " corresponding global stiffness matrix upper triangular matrix;
(3.3) adopt the numbering of first non-zero submatrices of every row in array J in the shaping array K record global stiffness matrix upper triangular matrix, and determine that with this non-zero submatrices is at the row of global stiffness matrix upper triangular matrix number, because K[i] be the numbering of capable first non-zero submatrices of i in array J in the global stiffness matrix upper triangular matrix, and K[i+1] be the numbering of capable first non-zero submatrices of i+1 in array J in the global stiffness matrix upper triangular matrix, thus in array J numbering from K[i] to K[i+1] the sub-square of all non-zeros of-1 correspondence is all capable at i.
4. storage of finite element analysis stiffness matrix and generation method in the metal volume Plastic Forming according to claim 1, it is characterized in that the finite element grid model after the geometric model of institute's problem analysis is divided be the two a mixed cell model of the hexahedral element, tetrahedron element of three-dimensional or its; Finite element grid model after perhaps the geometric model of institute's problem analysis is divided be the two a mixed cell model of two-dimentional quadrilateral units, triangular element or its.
5. according to storage of finite element analysis stiffness matrix and generation method in claim 1 or the 4 described metal volume Plastic Forming, it is characterized in that, in the described step (1.3), for a unit in the grid model, if the node number of this unit is that n, dimension are m, remember that the integral body of n the node in this unit is numbered I[i], i=1,2,, n, element stiffness matrix are U (m * n, m * n), then the arbitrary element of this element stiffness matrix is
U
ijkl＝U[(i-1)×m+k][(j-1)×m+l]
Wherein, i, j=1,2 ..., n is the part numbering of this cell node, k, and l=1,2 ..., m is physical quantity dimension numbering;
Determine U
IjklThe process of sequence number is as follows in grid model global stiffness matrix stores array A:
(4.1) determine U
IjklIntegral node numbering: I=I[i], J=I[j], wherein I, J are U
IjklThe non-zero submatrices A at place
IJIn row in grid model global stiffness matrix upper triangular matrix number and row number;
(4.2) determine non-zero submatrices A
IIThe quantity of non-zero submatrices: N before
I=K[I]-1;
(4.3) determine non-zero submatrices A
IIThe quantity of nonzero element before:
NN
I＝N
I×m
2-(I-1)×m×(m-1)/2
Wherein, m
2Be the quantity of each non-zero submatrices nonzero element, the minimizing item of the nonzero element that m (m-1)/2 causes for the diagonal angle non-zero submatrices of every row;
(4.4) determine the capable quantity of nonzero element before of k in the capable non-zero submatrices of I:
NN
IK＝(K[I+1]-K[I])×m×(k-1)
(4.5) determine that I is capable, the quantity of the non-zero submatrices before the J row; Search array J is from label K[I] to K[I+1] J[II]=the array label II of J, during I was capable, it was N that J is listed as non-zero submatrices quantity before
J=II-K[I];
(4.6) determine U
IjklNumbering in grid model global stiffness matrix, being expert at:
NN
Jl＝S
J×m-TT+l
Wherein, TT is for because the minimizing item of the nonzero element that causes of diagonal angle non-zero submatrices, TT=k (k-1)/2;
(4.7) determine U
IjklNumbering in grid model global stiffness matrix storage arrays A:
NN
IJkl＝NN
I+NN
Ik+NN
Jl
(4.8) according to the create-rule and the U of grid model global stiffness matrix
IjklPosition Number NN in grid model global stiffness matrix storage arrays A
IJkl, with U
IjklAssemble or be added among the global stiffness matrix storage arrays A of grid model;
(4.9) repeating step (4.1)～(4.8), each element in all element stiffness matrixs in the traversal grid model, the global stiffness matrix one dimension that can obtain grid model is stored array A.
6. storage of finite element analysis stiffness matrix and generation method in the metal volume Plastic Forming according to claim 1, it is characterized in that: it is Solving Linear that the analysis numerical analysis method that grid model adopted finally is summed up as storage, generation and stiffness equations thereof sparse, the symmetric population stiffness matrix, and the global stiffness equation solution parameter of the grid model of being analyzed is full dose or increment.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN 201110130073 CN102184298B (en) | 2011-05-18 | 2011-05-18 | Method for storing and generating stiffness matrix for finite-element analysis in metal bulk plastic forming |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN 201110130073 CN102184298B (en) | 2011-05-18 | 2011-05-18 | Method for storing and generating stiffness matrix for finite-element analysis in metal bulk plastic forming |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN102184298A true CN102184298A (en) | 2011-09-14 |
| CN102184298B CN102184298B (en) | 2013-06-19 |
Family
ID=44570475
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN 201110130073 Active CN102184298B (en) | 2011-05-18 | 2011-05-18 | Method for storing and generating stiffness matrix for finite-element analysis in metal bulk plastic forming |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN102184298B (en) |
Cited By (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN107908844A (en) * | 2017-11-07 | 2018-04-13 | 武汉科技大学 | Metamorphic mechanisms kinematics model auto-creating method |
| CN108724817A (en) * | 2018-05-04 | 2018-11-02 | 河南工业大学 | A kind of intelligence manufacture method of thermoplastic composite metal sandwich slab products |
| CN109388840A (en) * | 2017-08-10 | 2019-02-26 | 利弗莫尔软件技术公司 | The dedicated programmed computer of numerical simulation for the Metal Forming Process with predefined load paths and corresponding grid Adjusted Option |
| CN109948279A (en) * | 2019-03-29 | 2019-06-28 | 江苏精研科技股份有限公司 | A kind of emulation design method of metalwork shaping |
| CN111008497A (en) * | 2019-12-06 | 2020-04-14 | 集美大学 | Generation method of finite element total stiffness matrix and terminal |
| US10713400B2 (en) | 2017-04-23 | 2020-07-14 | Cmlabs Simulations Inc. | System and method for executing a simulation of a constrained multi-body system |
| CN111965694A (en) * | 2020-08-03 | 2020-11-20 | 中国石油天然气集团有限公司 | Method and device for determining position of seismic physical point and seismic observation system |
| CN112434451A (en) * | 2020-10-28 | 2021-03-02 | 西安交通大学 | Finite element analysis method based on block parallel computation |
| CN115081297A (en) * | 2022-08-23 | 2022-09-20 | 天津大学 | Overall stiffness matrix solving method in composite material plane elastic finite element analysis |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1139242A (en) * | 1996-02-02 | 1997-01-01 | 曾攀 | Numerical composite unit method for obtaining static and dynamic characteristics of engineering structure |
| CN1386201A (en) * | 2000-06-29 | 2002-12-18 | 目标储油层公司 | Method and system for solving finite element models using multi-phase physics |
| CN1609351A (en) * | 2004-11-24 | 2005-04-27 | 东南大学 | Set-up method for generalized Bingham soft soil rheological deformation analogue body |
-
2011
- 2011-05-18 CN CN 201110130073 patent/CN102184298B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1139242A (en) * | 1996-02-02 | 1997-01-01 | 曾攀 | Numerical composite unit method for obtaining static and dynamic characteristics of engineering structure |
| CN1386201A (en) * | 2000-06-29 | 2002-12-18 | 目标储油层公司 | Method and system for solving finite element models using multi-phase physics |
| CN1609351A (en) * | 2004-11-24 | 2005-04-27 | 东南大学 | Set-up method for generalized Bingham soft soil rheological deformation analogue body |
Cited By (15)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US10713400B2 (en) | 2017-04-23 | 2020-07-14 | Cmlabs Simulations Inc. | System and method for executing a simulation of a constrained multi-body system |
| CN109388840A (en) * | 2017-08-10 | 2019-02-26 | 利弗莫尔软件技术公司 | The dedicated programmed computer of numerical simulation for the Metal Forming Process with predefined load paths and corresponding grid Adjusted Option |
| CN109388840B (en) * | 2017-08-10 | 2023-02-10 | 利弗莫尔软件技术公司 | Special programming computer for numerical simulation of metal forming processes with predefined load paths and corresponding grid adjustment schemes |
| CN107908844A (en) * | 2017-11-07 | 2018-04-13 | 武汉科技大学 | Metamorphic mechanisms kinematics model auto-creating method |
| CN108724817A (en) * | 2018-05-04 | 2018-11-02 | 河南工业大学 | A kind of intelligence manufacture method of thermoplastic composite metal sandwich slab products |
| CN109948279A (en) * | 2019-03-29 | 2019-06-28 | 江苏精研科技股份有限公司 | A kind of emulation design method of metalwork shaping |
| CN109948279B (en) * | 2019-03-29 | 2022-12-20 | 江苏精研科技股份有限公司 | Simulation design method for shaping metal piece |
| CN111008497B (en) * | 2019-12-06 | 2022-11-11 | 集美大学 | Method for generating finite element total stiffness matrix and terminal |
| CN111008497A (en) * | 2019-12-06 | 2020-04-14 | 集美大学 | Generation method of finite element total stiffness matrix and terminal |
| CN111965694A (en) * | 2020-08-03 | 2020-11-20 | 中国石油天然气集团有限公司 | Method and device for determining position of seismic physical point and seismic observation system |
| CN111965694B (en) * | 2020-08-03 | 2024-01-30 | 中国石油天然气集团有限公司 | Method and device for determining position of physical point of earthquake and earthquake observation system |
| CN112434451A (en) * | 2020-10-28 | 2021-03-02 | 西安交通大学 | Finite element analysis method based on block parallel computation |
| CN112434451B (en) * | 2020-10-28 | 2024-04-02 | 西安交通大学 | Finite element analysis method based on block parallel computing |
| CN115081297B (en) * | 2022-08-23 | 2022-11-04 | 天津大学 | Overall stiffness matrix solving method in composite material plane elastic finite element analysis |
| CN115081297A (en) * | 2022-08-23 | 2022-09-20 | 天津大学 | Overall stiffness matrix solving method in composite material plane elastic finite element analysis |
Also Published As
| Publication number | Publication date |
|---|---|
| CN102184298B (en) | 2013-06-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN102184298B (en) | Method for storing and generating stiffness matrix for finite-element analysis in metal bulk plastic forming | |
| CN110110413B (en) | Structural topology optimization method based on material field reduction progression expansion | |
| CN109145427B (en) | Porous structure design and optimization method based on three-cycle minimum curved surface | |
| CN111709171B (en) | Isogeometric solving and heat dissipation topology generation method for heat flow strong coupling problem | |
| CN107609141B (en) | Method for performing rapid probabilistic modeling on large-scale renewable energy data | |
| CN105550744A (en) | Nerve network clustering method based on iteration | |
| CN108846212B (en) | Rigid frame pile internal force and displacement design calculation method | |
| CN103034766B (en) | A kind of laying angular direction of definite Test of Laminate Composites and the method for thickness | |
| CN106372347A (en) | Dynamic response topological optimization method implemented by application of improved bi-directional evolutionary structural optimization (BESO) to equivalent static load method | |
| CN104573281A (en) | Complex space curved surface thin plate forming die face designing method taking springback compensation | |
| CN113204906B (en) | Multiphase material topology optimization design method and system considering structural stability | |
| Gao et al. | An exact block-based reanalysis method for local modifications | |
| CN111507026A (en) | Dual-grid multi-scale finite element method for simulating node Darcy permeation flow rate | |
| CN103902504A (en) | Method for calculating inherent frequency of Euler-Bernoulli beam through improved differential transformation method | |
| Shi et al. | Analysis of price Stackelberg duopoly game with bounded rationality | |
| CN114347029A (en) | Model order reduction method for rapid simulation of pneumatic soft robot | |
| CN102063525A (en) | Method for generating practical multidisciplinary design optimization model automatically | |
| CN111539138B (en) | Method for solving time domain response sensitivity of structural dynamics peak based on step function | |
| CN111680441B (en) | Gradient lattice sandwich board structure suitable for thermal working condition | |
| CN112163385A (en) | Parallel non-grid method and system for solving strong heat fluid-solid coupling problem | |
| CN102799750B (en) | Method for quickly generating common side and non-common sides of geometry surface triangle | |
| CN111274624A (en) | Multi-working-condition special-shaped node topology optimization design method based on RBF proxy model | |
| CN114417681B (en) | Two-dimensional structure deformation monitoring method and device based on dynamic decision and neural network | |
| CN113486512B (en) | Flutter analysis method for functional gradient variable-thickness blade model | |
| CN113343512B (en) | Mobile-U-Net-based multi-scale topology optimization design method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant |