CN102837435A - Flow front prediction method of non-isothermal resin transfer molding based on mid-plane model - Google Patents

Abstract

The invention provides a flow front prediction method of non-isothermal resin transfer molding based on a mid-plane model, comprising the following steps: (1) constructing a node three-dimensional control volume by a die cavity mid-plane grid; (2) setting molding technological parameters; (3) setting material parameters; (4) calculating resin flowing pressure and speed; (5) calculating resin solidifying reaction; (6) calculating the temperatures of resin and fiber reinforcements; (7) updating flow front and material physical properties so as to ensure time step; (8) ensuring whether the die cavity is filled, and entering the step (9) if the die cavity has already been filled, or entering the step (5); and (9) predicting the distribution of dry spots. Based on the mid-plane grid, the three-dimensional control volume is constructed and the finite element method of the control volume is used, the changes of viscidity, speed, temperature, curing degree and the like in the thickness direction and the influence to the calculation of pressure field are considered in the calculation, and the flow front of resin transfer molding and the distribution of various physical quantity fields can be predicted more accurately.

Description

Flow front Forecasting Methodology based on the non-isothermal resin transfer moulding of middle surface model

Technical field

The invention belongs to numerical simulation method; Be particularly related to a kind of flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model, be used for the non-isothermal resin transfer molding of in-plane size much larger than the composite product of thickness direction size.

Background technology

Resin transfer moulding (Resin Transfer Molding; RTM) being a kind of resin to be injected the technology that close die is soaked into reinforcing material and curing molding, also is one of important method of producing in high quality fibers reinforced composite goods, compares with hot press molding technology with prepreg technology; The equipment investment of RTM technological forming is few, the production environmental protection; Stick with paste the technology ratio with traditional hand, the RTM quality of item is stable, it is few to discharge unwanted volatile matter, and RTM goods mechanical property tool designability, stowing pressure is low, cost is low; Be fit to the production and the mechanization production of large complicated part, be widely used in fields such as Aero-Space, boats and ships, automobile, power electronics.

Resin transfer moulding mainly comprises operations such as shop layer, blanking, dipping, moulding, curing.Wherein the resin dipping process is a vital ring, and dipping process directly is reflected on the stowing operation flow front.Injection pressure, injection speed, inject time, temperature, curing reaction, gate location, gate size, cast gate quantity and many factors such as shape of product and thickness all have material impact to the resin flows forward position; Flow front develops inhomogeneous or flow front converges the quality problems that all can cause composite product; As: the excessive distortion that causes washing away fibre reinforcement and cause fibre reinforcement of injection pressure; The too fast dipping that causes that flows is bad; The flow front meet can form dry spot, long meeting inject time make resin in flow process, a large amount of curing reactions just take place and production efficiency low, or the like.If the solution of these problems relies on traditional experience method, trial-and-error method; Promptly constantly evaluate product design and technological parameter again after the die trial repeatedly; Can make that cost height, the construction cycle of resin transfer molding (RTM) process are long; And the percent defective of product is high, production is unstable, therefore, can waste lot of manpower and material resources.And, a kind of means of economic, practical and science are provided for RTM product design and technological design then based on the flow front of the method for numerical simulation prediction resin transfer moulding.

Have the scholar to adopt the flow front of The Numerical Simulation Methods resin transfer moulding, like the RTM resin simulation softward of BJ University of Aeronautics & Astronautics's exploitation, in the filling process not account temperature and curing reaction to the influence of flow front.In addition; Also have and consider that in the resin fill die cavity process temperature and curing degree influence flow front; See " execute and fly; Dong Xianghuai. the numerical simulation of non-isothermal RTM technology. composite journal, 26 (4): 146-150,2009 " and " Chen Ren-liang; Li Ming-cheng.Study of flow and temperature in resin transfer molding process by numerical model.Transactions of Nanjing University of Aeronautics&Astronautics.20 (2): 165-171; 2003 ", and be Green's (Green) formula or divergence theorem but it adopts on middle veil lattice in-plane, on thickness direction, adopt finite difference method; And other the variation of physical quantity (like speed, viscosity, curing degree, flow rate etc.) on thickness direction do not considered in the only variation of account temperature on thickness direction.In the pressure solution procedure; Said method does not consider that temperature and the curing degree on the thickness direction changes the viscosity variation that causes; And then cause the variation of the flowing velocity of diverse location on the thickness direction, in the temperature solution procedure, said method adopts 1 relatively poor rank upstreame scheme of precision in the processing of convective term; Flowing velocity does not change on the thickness direction; Be coupled to flood partly and conduct with the heat of not flooding the heat transfer of part and ignore on the middle veil lattice in-plane, in curing reaction was found the solution, said method was not considered the variation of the resin flows speed of diverse location on the thickness direction, different curing reaction speed and curing degree.

Along with global market competition is growing more intense, require as far as possible science, design RTM product and technological parameter is set exactly, improve product quality and production efficiency, therefore, need a kind of method of prediction flow front of comprehensive science.

Summary of the invention

The object of the invention is to provide a kind of flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model; Method for numerical simulation to solve non-isothermal RTM technology of the prior art adopts finite difference method on thickness direction; And the only variation of account temperature on thickness direction; Do not consider other the variation of physical quantity (like speed, viscosity, curing degree etc.) on thickness direction, do not consider that impregnation of fibers strengthens body and the impregnation of fibers not high technical matters of precision that strengthen the coupled heat transfer between the body and on middle veil lattice in-plane, only adopt the numerical computations that 1 rank upstreame scheme etc. causes not.

The object of the invention realizes through following technical scheme:

A kind of flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model may further comprise the steps:

(1) makes the three-dimensional control of node volume by veil lattice in the die cavity;

(2) molding technique parameter is set;

(3) material parameter is set;

(4) calculate resin flowing pressure and speed;

(5) calculate the resin curing reaction;

(6) temperature of calculating resin and fibre reinforcement;

(7) upgrade flow front and material property, confirm time step;

(8) confirm whether die cavity is full of, full if mold cavity has been filled, then get into step (9), otherwise get into step (5);

(9) the prediction dry spot distributes.

Preferably; The making the three-dimensional control of node volume step by veil lattice in the die cavity and further comprise of said step (1): import veil lattice in the die cavity; According to the grid number of plies that thickness direction divided; Connect the corresponding node of adjacent two layers unit on the thickness direction, obtain the multilayer prism elements, the body-centered, the center of area, the limit mid point that connect prism elements can the three-dimensional control of structure node volumes.

Preferably; In the calculating resin flowing pressure and speed step of said step (4); Ignore products thickness side and upward pressure and change and along the velocity component of thickness direction, the resin plane flowing velocity on the thickness direction changes, the calculation procedure of pressure and speed further comprises:

(1) according to the viscosity of diverse location on the Viscosity Model calculating products thickness direction, resin viscosity is the function of temperature and curing degree:

$\mathrm{\&eta;}(T,\mathrm{\&alpha;})=\mathrm{Bexp}\left(\frac{T_b}{T}\right)\&CenterDot;{\left(\frac{{\mathrm{\&alpha;}}_{\mathrm{gel}}}{{\mathrm{\&alpha;}}_{\mathrm{gel}}-\mathrm{\&alpha;}}\right)}^{{\mathrm{<figure-callout\; id="1"\; label="C"\; filenames="HDA00002016184000031.png,HDA00002016184000032.png"\; state="\{\{state\}\}">C</figure-callout>}}_1+{\mathrm{<figure-callout\; id="2"\; label="C"\; filenames="HDA00002016184000032.png"\; state="\{\{state\}\}">C</figure-callout>}}_2\mathrm{\&alpha;}}$

In the formula, η is a dynamic viscosity, and T is a temperature, and α is a curing degree, B, T _b, C ₁, C ₂, α _GelBe material parameter;

(2) governing equation that adopts the conduct of Darcy formula and mass-conservation equation to find the solution pressure and speed:

$u_i=-\frac{{\mathrm{\&lambda;}}_{\mathrm{ij}}}{\mathrm{\&eta;}}{p}_{,j}$

u _i，i=0

In above-mentioned two formulas, following formula is a Darcy formula, and following formula is a mass-conservation equation, subscript ", " expression local derviation, u _iBe the superficial velocity of resin fluid, i, j=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, λ _Ij, η is respectively the permeability tensor component and the dynamic viscosity of resin, p is a resin flows pressure;

Consider viscosity change on the thickness direction, be based on the pressure equation of middle surface model:

${\&Integral;}_{\mathrm{\&Gamma;}}({\&Integral;}_0^h-\frac{{\mathrm{\&lambda;}}_{\mathrm{ij}}}{\mathrm{\&eta;}\left(z\right)}{p}_{,\mathrm{ij}}\mathrm{dz})n_i\mathrm{dS}=0$

In the formula, h is a products thickness, and η (z) expression viscosities il is the function of thickness direction coordinate z, and Γ is the node control volume surface, n _i(i=1,2) be node control volume surface along in method vector outside the unit of veil lattice in-plane;

(3) utilize control volume finite element model for solving pressure equation;

(4) by node pressure, according to the resin flows speed of diverse location on the Darcy formula calculated thickness direction.

Preferably, in the calculating resin curing reaction step of said step (5), utilize control volume finite element model for solving curing reaction equation, be divided into multilayer on the thickness, every layer of resin flows speed is different, and convective term adopts quality weighting upstreame scheme; Because the temperature and the resin flows speed of diverse location are inconsistent on the thickness direction, so the curing degree of diverse location is also inconsistent on the thickness direction;

$\mathrm{\&phi;}{\mathrm{\&alpha;}}_{,t}+u_i{\mathrm{\&alpha;}}_{,i}=\mathrm{\&phi;}\stackrel{\&CenterDot;}{G}(T,\mathrm{\&alpha;})$

$\stackrel{\&CenterDot;}{G}(T,\mathrm{\&alpha;})=[{\mathrm{<figure-callout\; id="1"\; label="A"\; filenames="HDA00002016184000031.png,HDA00002016184000032.png"\; state="\{\{state\}\}">A</figure-callout>}}_1\exp (-\frac{{\mathrm{<figure-callout\; id="1"\; label="E"\; filenames="HDA00002016184000031.png,HDA00002016184000032.png"\; state="\{\{state\}\}">E</figure-callout>}}_1}{\mathrm{RT}})+{\mathrm{<figure-callout\; id="2"\; label="A"\; filenames="HDA00002016184000032.png"\; state="\{\{state\}\}">A</figure-callout>}}_2\exp (-\frac{{\mathrm{<figure-callout\; id="2"\; label="E"\; filenames="HDA00002016184000032.png"\; state="\{\{state\}\}">E</figure-callout>}}_2}{\mathrm{RT}}){\mathrm{\&alpha;}}^m]{({\mathrm{\&alpha;}}_f-\mathrm{\&alpha;})}^n$

α _f=a ₀T+b ₀

In the formula, Being curing reaction speed, is the function of temperature and curing degree, subscript ", " expression local derviation, and φ is the fibre reinforced materials porosity, and t is the time, and T is an absolute temperature, u _iThe velocity component of veil lattice in-plane in the expression, i=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, α _fThe final conversion degree of isothermal under the expression assigned temperature, A ₁, A ₂, E ₁, E ₂, m, n, a ₀, b ₀Be material constant, R=8.314J/ (molK) is a universal gas constant.

Preferably, in the temperature step of the calculating resin of said step (6) and fibre reinforcement, utilize control volume finite element model for solving energy conservation equation; It is the temperature equation; Be divided into multilayer on the thickness, every layer of resin flows speed is different, and convective term adopts quality weighting upstreame scheme; Consider not impregnation of fibers enhancing body and the coupled heat transfer of impregnation of fibers precast body in the calculating, heat conduction in considering simultaneously on the veil lattice in-plane and the heat conduction on the thickness direction;

$[{\mathrm{\&phi;\&rho;}}_r{C}_{\mathrm{pr}}+(1-\mathrm{\&phi;}){\mathrm{\&rho;}}_f{C}_{\mathrm{pf}}]T_{,t}+{\mathrm{\&rho;}}_r{C}_{\mathrm{pr}}\left(u_j{T}_{,j}\right)={\mathrm{kT}}_{,\mathrm{ii}}+\mathrm{\&phi;\&Delta;H}\stackrel{\&CenterDot;}{G}$

(1-φ)ρ _fC _pfT _，t=kT _，ii

In above-mentioned two formulas, following formula is that impregnation of fibers strengthens the energy conservation equation of body, and following formula be the energy conservation equation of impregnation of fibers enhancing body not; This step is with above-mentioned two formula Solving Coupled, and subscript ", " is represented local derviation; T is the temperature of fibre reinforcement and resin, and t is the time, u _jThe velocity component of veil lattice in-plane in the expression, i, j=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, and i=3 is the coordinate components on the thickness direction, and φ is the porosity of fibre reinforcement, ρ _r, C _PrBe respectively resin density, specific heat at constant pressure, ρ _f, C _PfBe respectively fiber reinforcement volume density, specific heat at constant pressure, Be the thermal source item of mold filling process, wherein Δ H, Be respectively curing reaction heat, curing reaction speed, k is the equivalent coefficient of heat conduction;

$k=\frac{k_r{k}_f}{k_r{w}_f+k_f{w}_r},$ $w_r=\frac{\mathrm{\&phi;}/{\mathrm{\&rho;}}_f}{\mathrm{\&phi;}/{\mathrm{\&rho;}}_f+(1-\mathrm{\&phi;})/{\mathrm{\&rho;}}_r},$ w _f=1-w _r

In the formula, w _r, w _fBe respectively the mass ratio of resin and fibre reinforcement, k _r, k _fBe respectively the coefficient of heat conduction of resin and fiber.

Preferably, the said molding technique parameter in the said step (2) comprises mold temperature, resin temperature and fiber reinforcement temperature, and said molding technique parameter also comprises a kind of in injection pressure, inject time or the injection speed.

Preferably, the said material parameter in the said step (3) comprises density, specific heat capacity, the coefficient of heat conduction, porosity and the permeability of resin and fibre reinforcement, and said material parameter also comprises the material constant of Viscosity Model and curing reaction model.

Preferably, the material constant of Viscosity Model comprises A ₁, A ₂, E ₁, E ₂, m, n, a ₀, b ₀, the material constant of said curing reaction model comprises B, T _b, C ₁, C ₂, α _Gel

Preferably; Said step (7) further comprises: according to the velocity field that calculates; Adopt the flow network analytic approach; Select the shortest required filling time of the control volume of each node on whole thickness direction in the current flow front as time step, upgrade the filling rate of other node control volume in the flow front and the material property of diverse location, said material property comprises viscosity, specific heat, density and the coefficient of heat conduction.

Preferably, said step (9) further comprises: according to the loading time of each node on the grid, the morning and evening of contrast filling, the prediction dry spot distributes, i.e. the dry spot position that possibly occur.

Compared with prior art, the present invention has following advantage:

1, compares with traditional middle surface model; Construct three-dimensional control volume among the present invention and used the control volume finite-element method; The variation of viscosity, speed, temperature, curing degree etc. on the thickness direction and the influence that pressure field is calculated thereof have been considered in the calculating; In energy conservation equation, adopt thermal balance model; Considered that impregnation of fibers does not strengthen body and the coupled heat transfer and the conduction of the heat on the middle veil lattice in-plane of impregnation of fibers precast body, the convective term of temperature and curing degree is handled and has been adopted quality weighting upstreame scheme, has improved the stability and the computational accuracy of numerical computations; More, can predict the distribution of resin transfer moulding flow front and various Physical Quantity Field more accurately near actual formulations of resin transfer molding processes;

2, the present invention is the basis with the fundamental equation of porous media hydrodynamics; Introduce reasonably hypothesis and simplification; By the three-dimensional control of middle veil lattice model structure node volume, adopt control volume finite-element method simulation resin packed type chamber process, calculate pressure, speed, curing degree and temperature in the resin flow process; The flow front of prediction resin fill; Can be and rational gate location, gate size, cast gate quantity and flash mouth position are set scientific basis is provided, avoid the dry spot generation of defects, thereby improve product quality and production efficiency.

Description of drawings

Fig. 1 is the middle veil lattice model of certain goods and the sketch map of cast gate;

Fig. 2 is a FB(flow block) of the present invention;

Fig. 3 is the three-dimensional control of structure node of the present invention volume sketch map;

Fig. 4 is a quality weighting upstreame scheme sketch map of the present invention;

The curing degree distribution schematic diagram of products thickness central core when Fig. 5 is the die cavity sand off;

The Temperature Distribution sketch map of products thickness central core when Fig. 6 is the die cavity sand off;

Fig. 7 is a flow front result of calculation sketch map of the present invention;

Fig. 8 is a dry spot position prediction sketch map of the present invention.

The specific embodiment

Below in conjunction with accompanying drawing, specify the present invention.

See also Fig. 1, be the sketch map of veil lattice and gate location in the geometrical model of certain goods, this products thickness is 5mm.The step of prediction flow front is as shown in Figure 2, and concrete steps are following.

(1) makes the three-dimensional control of node volume by veil lattice in the die cavity.See also Fig. 3; Subscript i among the figure, i+1 (l >=i>0; L is the unit number of plies of middle surface model thickness direction, l=13 in this example) i, i+1 node layer on the expression thickness direction, connect i, i+1 layer unit corresponding node is configured to tri-prism element; The son control volume that connects tri-prism element body-centered, the center of area, limit mid point structure node again, the node control volume in all tri-prism elements is merged into the three-dimensional control of node volume.

(2) molding technique parameter is set.Molding technique parameter comprises injection pressure or inject time or injection speed (3 optional 1), mold temperature, resin temperature, fiber reinforcement temperature etc.; It is 30 ℃ that the resin implantation temperature is set, and the fiber reinforcement temperature is 60 ℃, and the die cavity wall temperature is 60 ℃.

When the technological parameter that is provided with is injection pressure, with the boundary condition of injection pressure as pressure equation calculating;

When the technological parameter that is provided with is injection rate,, the flow of injection is equaled or the injection rate that approaches to be provided with satisfies the injection rate requirement through continuous adjustment injection pressure;

When the technological parameter that is provided with is inject time, be translated into injection rate according to the volume of goods, find the solution pressure equation according to the mode of injection rate again; In this example, injection pressure being set is 100KPa.

(3) material parameter is set.In this example; Adopt parameter, Viscosity Model and the curing reaction model of certain material; When selecting other material for use; Only need to revise the flow front prediction that relevant material parameter can realize the non-isothermal resin transfer moulding of other material, the fibre reinforcement porosity is 60%, and material parameter and resin viscosity parameter, resin solidification reaction rate parameter are shown in table 1-table 3.

Table 1, material parameter

Table 2, resin viscosity parameter

Table 3, resin solidification reaction rate parameter

In the resin viscosity of material parameter and the curing reaction speed, T is a temperature, and α is a curing degree;

(4) calculate resin flowing pressure and speed.According to the viscosity of diverse location on the Viscosity Model calculating products thickness direction, resin viscosity is the function of temperature and curing degree, adopts the resin viscosity model in the step (3).

The governing equation that adopts the conduct of Darcy formula and mass-conservation equation to find the solution pressure and speed:

$u_i=-\frac{{\mathrm{\&lambda;}}_{\mathrm{ij}}}{\mathrm{\&eta;}}{\mathrm{\&rho;}}_{,j}$

u _i,i=0

In above-mentioned two formulas, following formula is a Darcy formula, and following formula is a mass-conservation equation, subscript ", " expression local derviation, u _iBe the superficial velocity of resin fluid, i, j=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, λ _Ij, η is respectively the permeability tensor component and the viscosity of resin, p is a resin flows pressure;

Consider viscosity change on the thickness direction, be based on the pressure equation of middle surface model:

${\&Integral;}_{\mathrm{\&Gamma;}}({\&Integral;}_0^h-\frac{{\mathrm{\&lambda;}}_{\mathrm{ij}}}{\mathrm{\&eta;}\left(z\right)}{p}_{,\mathrm{ij}}\mathrm{dz})n_i\mathrm{dS}=0$

In the formula, h is a products thickness, and η (z) expression viscosities il is the function of thickness direction coordinate z, and Γ is the node control volume surface, n _i(i=1,2) be node control volume surface along in method vector outside the unit of veil lattice in-plane; Establishing infiltration principal direction and fibre reinforcement consistent with resin flows speed principal direction in this example is isotropic material, and the permeability tensor component is:

${\mathrm{\&lambda;}}_{\mathrm{ij}}=\left\{\begin{array}{cc}1.0\&times;{10}^{-8}& i=j\\ 0& i\&NotEqual;j\end{array}\right.,$

Set up the pressure equation on the control volume of each node, apply boundary condition, all node pressure equation solutions of simultaneous can calculate node pressure.

Utilizing control volume finite element model for solving pressure equation, is 0 along thickness direction pressure partial derivative, and tri-prism element internal pressure partial derivative adopts the partial derivative of the shape function of plane triangle unit:

p _,i=N _β,ip _β

In the formula, N _β(β=1,2,3) are the shape function of β the node in plane triangle unit, p _βPressure for β node of triangular element;

By node pressure, according to the resin flows speed of diverse location on the Darcy formula calculated thickness direction.

(5) calculate the resin curing reaction.Adopt the control volume finite-element method, time difference adopts backward difference, and convective term adopts the upstreame scheme of quality weighting, each node control volume is set up the semi-implicit prediction-correction equation of node curing degree:

${\&Integral;}_V\mathrm{\&phi;}\frac{{\mathrm{\&alpha;}}^{t+\mathrm{\&theta;\&Delta;t}}-{\mathrm{\&alpha;}}^t}{\mathrm{\&theta;\&Delta;t}}\mathrm{dV}+{\&Integral;}_{\mathrm{\&Gamma;}}{m}_{\mathrm{\&gamma;}}{\mathrm{\&alpha;}}_{\mathrm{\&gamma;}}^{t+\mathrm{\&theta;\&Delta;t}}{n}_i\mathrm{dS}={\&Integral;}_V\mathrm{\&phi;}\stackrel{\&CenterDot;}{G}(T,{\mathrm{\&alpha;}}^t)\mathrm{dV}$

${\&Integral;}_V\mathrm{\&phi;}\frac{{\mathrm{\&alpha;}}^{t+\mathrm{\&Delta;t}}-{\mathrm{\&alpha;}}^t}{\mathrm{\&Delta;t}}\mathrm{dV}+{\&Integral;}_{\mathrm{\&Gamma;}}{m}_{\mathrm{\&gamma;}}{\mathrm{\&alpha;}}_{\mathrm{\&gamma;}}^{t+\mathrm{\&Delta;t}}{n}_i\mathrm{dS}={\&Integral;}_V\mathrm{\&phi;}\stackrel{\&CenterDot;}{G}(T,{\mathrm{\&alpha;}}^{t+\mathrm{\&theta;\&Delta;t}})\mathrm{dV}$

In the formula, following formula is a predictive equation, and following formula is a correction equation, and V is the volume of node control volume, α ^{T+ Δ t}Be t+ Δ t curing degree constantly, α ^tBe t curing degree constantly, Δ t is a time step, 0<θ≤1, employing prediction-alignment technique can improve the computational accuracy of curing degree, m _γ(γ=r, s t) are the flow on the control volume surface gamma owing to ignore the velocity component on the thickness direction, so convective term exist only in the in-plane of veil lattice, as shown in Figure 4,

Be the lip-deep curing degree of control volume:

${\stackrel{\&CenterDot;}{m}}_r={\&Integral;}_o^a{u}_i\&CenterDot;n_i\frac{h}{l}\mathrm{ds}$

${\stackrel{\&CenterDot;}{m}}_s={\&Integral;}_o^b{u}_i\&CenterDot;n_i\frac{h}{l}\mathrm{ds}$

${\stackrel{\&CenterDot;}{m}}_t={\&Integral;}_o^c{u}_i\&CenterDot;n_i\frac{h}{l}\mathrm{ds}$

${\mathrm{\&alpha;}}_r=\left\{\begin{array}{ccc}{f}_P{\mathrm{\&alpha;}}_t+(1-f_P){\mathrm{\&alpha;}}_1& f_P=\min [\max (-{\stackrel{\&CenterDot;}{m}}_t/{\stackrel{\&CenterDot;}{m}}_r,0),1]& {\stackrel{\&CenterDot;}{m}}_r>0\\ f_P{\mathrm{\&alpha;}}_s+(1-f_P){\mathrm{\&alpha;}}_2& f_P=\min [\max (-{\stackrel{\&CenterDot;}{m}}_s/{\stackrel{\&CenterDot;}{m}}_r,0),1]& {\stackrel{\&CenterDot;}{m}}_r<0\end{array}\right.$

${\mathrm{\&alpha;}}_s=\left\{\begin{array}{ccc}{f}_P{\mathrm{\&alpha;}}_t+(1-f_P){\mathrm{\&alpha;}}_3& f_P=\min [\max ({\stackrel{\&CenterDot;}{m}}_t/{\stackrel{\&CenterDot;}{m}}_s,0),1]& {\stackrel{\&CenterDot;}{m}}_s>0\\ f_P{\mathrm{\&alpha;}}_r+(1-f_P){\mathrm{\&alpha;}}_2& f_P=\min [\max (-{\stackrel{\&CenterDot;}{m}}_r/{\stackrel{\&CenterDot;}{m}}_s,0),1]& {\stackrel{\&CenterDot;}{m}}_s<0\end{array}\right.$

${\mathrm{\&alpha;}}_t=\left\{\begin{array}{ccc}{f}_P{\mathrm{\&alpha;}}_r+(1-f_P){\mathrm{\&alpha;}}_1& f_P=\min [\max (-{\stackrel{\&CenterDot;}{m}}_r/{\stackrel{\&CenterDot;}{m}}_t,0),1]& {\stackrel{\&CenterDot;}{m}}_s>0\\ f_P{\mathrm{\&alpha;}}_s+(1-f_P){\mathrm{\&alpha;}}_3& f_P=\min [\max ({\stackrel{\&CenterDot;}{m}}_s/{\stackrel{\&CenterDot;}{m}}_t,0),1]& {\stackrel{\&CenterDot;}{m}}_s<0\end{array}\right.$

In the formula, α _γ(γ=r, s, t) curing degree of an expression Atria node.The curing degree of the mid-depth layer during end-of-fill distributes as shown in Figure 5.

(6) temperature of calculating resin and fibre reinforcement.Adopt the control volume finite-element method, time difference adopts backward difference, and convective term adopts the upstreame scheme (similar curing degree) of quality weighting, and half implicit equation of each node control volume being set up node temperature is:

${\&Integral;}_V[{\mathrm{\&phi;\&rho;}}_r{C}_{\mathrm{pr}}+(1-\mathrm{\&phi;}){\mathrm{\&rho;}}_f{C}_{\mathrm{pf}}]\frac{T^{t+\mathrm{\&Delta;t}}-T^t}{\mathrm{\&Delta;t}}\mathrm{dV}+{\&Integral;}_{\mathrm{\&Gamma;}}{\mathrm{\&rho;}}_r{C}_{\mathrm{pr}}{m}_{\mathrm{\&gamma;}}{T}_{\mathrm{\&gamma;}}^{t+\mathrm{\&Delta;t}}{n}_i\mathrm{dS}$

$-{\&Integral;}_{\mathrm{\&Gamma;}}{\mathrm{kT}}_{,i}^{t+\mathrm{\&Delta;t}}{n}_i\mathrm{dS}={\&Integral;}_V\mathrm{\&phi;\&Delta;H}\stackrel{\&CenterDot;}{G}\mathrm{dV}$

${\&Integral;}_V\left[(1-\mathrm{\&phi;}){\mathrm{\&rho;}}_f{C}_{\mathrm{pf}}\right]\frac{T^{t+\mathrm{\&Delta;t}}-T^t}{\mathrm{\&Delta;t}}\mathrm{dV}-{\&Integral;}_{\mathrm{\&Gamma;}}{\mathrm{kT}}_{,i}^{t+\mathrm{\&Delta;t}}{n}_i\mathrm{dS}=0$

In the formula, T _{, i}=N _{β, i}T _β, on middle veil lattice in-plane, N _βFor the shape function of middle veil lattice unit, on thickness direction, N _βBe 2 dotted line unit shape functions.Above-mentioned two equations of simultaneous can be found the solution the temperature in the whole die cavity.The Temperature Distribution of the mid-depth layer during end-of-fill is as shown in Figure 6.

(7) upgrade flow front and material property, confirm the time step step.According to the velocity field that calculates; Adopt the flow network analytic approach; Select the shortest required filling time of each node control volume in the current flow front as time step, upgrade the filling rate of other node control volume in the flow front and the material property (viscosity, density, specific heat, the coefficient of heat conduction etc.) of diverse location.

(8) confirm whether die cavity is full of step.If mold cavity has been filled full, then get into (9), otherwise get into step (5).It is as shown in Figure 7 that die cavity is filled the flow front distribution that obtains when expiring, the moment that this place of the color showing at every place is full of on the goods among the figure.

(9) prediction dry spot distribution step.According to the loading time of each node on the grid, the morning and evening of contrast filling, the node place of all node loading times dry spot possibly occur, like red some position among Fig. 8 around loading time is later than.

Compare with traditional middle surface model; Construct three-dimensional control volume among the present invention and used the control volume finite-element method; The variation of viscosity, speed, temperature, curing degree etc. on the thickness direction and the influence that pressure field is calculated thereof have been considered in the calculating; In energy conservation equation, adopt thermal balance model; Considered that impregnation of fibers does not strengthen body and the coupled heat transfer and the conduction of the heat on the middle veil lattice in-plane of impregnation of fibers precast body, the convective term of temperature and curing degree is handled and has been adopted quality weighting upstreame scheme, has improved the stability and the computational accuracy of numerical computations; More, can predict the distribution of resin transfer moulding flow front and various Physical Quantity Field more accurately near actual formulations of resin transfer molding processes.The present invention is the basis with the fundamental equation of porous media hydrodynamics; Introduce reasonably hypothesis and simplification; By the three-dimensional control of middle veil lattice model structure node volume, adopt control volume finite-element method simulation resin packed type chamber process, calculate pressure, speed, curing degree and temperature in the resin flow process; The flow front of prediction resin fill; Can be and rational gate location, gate size, cast gate quantity and flash mouth position are set scientific basis is provided, avoid the dry spot generation of defects, thereby improve product quality and production efficiency.

More than the disclosed several specific embodiments that are merely the application, but the application is not limited thereto, any those skilled in the art can think variation, all should drop in the application's the protection domain.

Claims (10)

1. the flow front Forecasting Methodology based on the non-isothermal resin transfer moulding of middle surface model is characterized in that, may further comprise the steps:

(1) makes the three-dimensional control of node volume by veil lattice in the die cavity;

(2) molding technique parameter is set;

(3) material parameter is set;

(4) calculate resin flowing pressure and speed;

(5) calculate the resin curing reaction;

(6) temperature of calculating resin and fibre reinforcement;

(7) upgrade flow front and material property, confirm time step;

(8) confirm whether die cavity is full of, full if mold cavity has been filled, then get into step (9), otherwise get into step (5);

(9) the prediction dry spot distributes.

2. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; The making the three-dimensional control of node volume step by veil lattice in the die cavity and further comprise of said step (1): import veil lattice in the die cavity; According to the grid number of plies that thickness direction divided; Connect the corresponding node of adjacent two layers unit on the thickness direction, obtain the multilayer prism elements, the body-centered, the center of area, the limit mid point that connect prism elements can the three-dimensional control of structure node volumes.

3. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; In the calculating resin flowing pressure and speed step of said step (4); Ignore products thickness side and upward pressure and change and along the velocity component of thickness direction, the resin plane flowing velocity on the thickness direction changes, the calculation procedure of pressure and speed further comprises:

(1) according to the viscosity of diverse location on the Viscosity Model calculating products thickness direction, resin viscosity is the function of temperature and curing degree:

$\mathrm{\&eta;}(T,\mathrm{\&alpha;})=\mathrm{Bexp}\left(\frac{T_b}{T}\right)\&CenterDot;{\left(\frac{{\mathrm{\&alpha;}}_{\mathrm{gel}}}{{\mathrm{\&alpha;}}_{\mathrm{gel}}-\mathrm{\&alpha;}}\right)}^{C_1+C_2\mathrm{\&alpha;}}$

In the formula, η is a dynamic viscosity, and T is a temperature, and α is a curing degree, B, T _b, C ₁, C ₂, α _GelBe material parameter;

(2) governing equation that adopts the conduct of Darcy formula and mass-conservation equation to find the solution pressure and speed:

$u_i=-\frac{{\mathrm{\&lambda;}}_{\mathrm{ij}}}{\mathrm{\&eta;}}{p}_{,j}$

u _i,i=0

In above-mentioned two formulas, following formula is a Darcy formula, and following formula is a mass-conservation equation, subscript ", " expression local derviation, u _iBe the superficial velocity of resin fluid, i, j=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, λ _Ij, η is respectively the permeability tensor component and the dynamic viscosity of resin, p is a resin flows pressure;

Consider viscosity change on the thickness direction, be based on the pressure equation of middle surface model:

${\&Integral;}_{\mathrm{\&Gamma;}}({\&Integral;}_0^h-\frac{{\mathrm{\&lambda;}}_{\mathrm{ij}}}{\mathrm{\&eta;}\left(z\right)}{p}_{,\mathrm{ij}}\mathrm{dz})n_i\mathrm{dS}=0$

In the formula, h is a products thickness, and η (z) expression viscosities il is the function of thickness direction coordinate z, and Γ is the node control volume surface, n _i(i=1,2) be node control volume surface along in method vector outside the unit of veil lattice in-plane;

(3) utilize control volume finite element model for solving pressure equation;

(4) by node pressure, according to the resin flows speed of diverse location on the Darcy formula calculated thickness direction.

4. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; In the calculating resin curing reaction step of said step (5), utilize control volume finite element model for solving curing reaction equation, be divided into multilayer on the thickness; Every layer of resin flows speed is different, and convective term adopts quality weighting upstreame scheme; Because the temperature and the resin flows speed of diverse location are inconsistent on the thickness direction, so the curing degree of diverse location is also inconsistent on the thickness direction;

$\mathrm{\&phi;}{\mathrm{\&alpha;}}_{,t}+u_i{\mathrm{\&alpha;}}_{,i}=\mathrm{\&phi;}\stackrel{\&CenterDot;}{G}(T,\mathrm{\&alpha;})$

$\stackrel{\&CenterDot;}{G}(T,\mathrm{\&alpha;})=[A_1\exp (-\frac{E_1}{\mathrm{RT}})+A_2\exp (-\frac{E_2}{\mathrm{RT}}){\mathrm{\&alpha;}}^m]{({\mathrm{\&alpha;}}_f-\mathrm{\&alpha;})}^n$

α _f=a ₀T+b ₀

In the formula,

Being curing reaction speed, is the function of temperature and curing degree, subscript ", " expression local derviation, and φ is the fibre reinforced materials porosity, and t is the time, and T is an absolute temperature, u _iThe velocity component of veil lattice in-plane in the expression, i=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, α _fThe final conversion degree of isothermal under the expression assigned temperature, A ₁, A ₂, E ₁, E ₂, m, n, a ₀, b ₀Be material constant, R=8.314J/ (molK) is a universal gas constant.

5. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that, in the temperature step of the calculating resin of said step (6) and fibre reinforcement, utilize control volume finite element model for solving energy conservation equation; It is the temperature equation; Be divided into multilayer on the thickness, every layer of resin flows speed is different, and convective term adopts quality weighting upstreame scheme; Consider not impregnation of fibers enhancing body and the coupled heat transfer of impregnation of fibers precast body in the calculating, heat conduction in considering simultaneously on the veil lattice in-plane and the heat conduction on the thickness direction;

$[{\mathrm{\&phi;\&rho;}}_r{C}_{\mathrm{pr}}+(1-\mathrm{\&phi;}){\mathrm{\&rho;}}_f{C}_{\mathrm{pf}}]T_{,t}+{\mathrm{\&rho;}}_r{C}_{\mathrm{pr}}\left(u_j{T}_{,j}\right)={\mathrm{kT}}_{,\mathrm{ii}}+\mathrm{\&phi;\&Delta;H}\stackrel{\&CenterDot;}{G}$

(1-φ)ρ _fC _pfT _，t=kT _，ii

In above-mentioned two formulas, following formula is that impregnation of fibers strengthens the energy conservation equation of body, and following formula be the energy conservation equation of impregnation of fibers enhancing body not; This step is with above-mentioned two formula Solving Coupled, and subscript ", " is represented local derviation; T is the temperature of fibre reinforcement and resin, and t is the time, u _jThe velocity component of veil lattice in-plane in the expression, i, j=1,2 for the local Cartesian coordinates on the middle veil lattice in-plane is a coordinate components, and i=3 is the coordinate components on the thickness direction, and φ is the porosity of fibre reinforcement, ρ _r, C _PrBe respectively resin density, specific heat at constant pressure, ρ _f, C _PfBe respectively fiber reinforcement volume density, specific heat at constant pressure,

Be the thermal source item of mold filling process, wherein Δ H,

Be respectively curing reaction heat, curing reaction speed, k is the equivalent coefficient of heat conduction;

$k=\frac{k_r{k}_f}{k_r{w}_f+k_f{w}_r},$ $w_r=\frac{\mathrm{\&phi;}/{\mathrm{\&rho;}}_f}{\mathrm{\&phi;}/{\mathrm{\&rho;}}_f+(1-\mathrm{\&phi;})/{\mathrm{\&rho;}}_r},$ w _f=1-w _r

In the formula, w _r, w _fBe respectively the mass ratio of resin and fibre reinforcement, k _r, k _fBe respectively the coefficient of heat conduction of resin and fiber.

6. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; Said molding technique parameter in the said step (2) comprises mold temperature, resin temperature and fiber reinforcement temperature, and said molding technique parameter also comprises a kind of in injection pressure, inject time or the injection speed.

7. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; Said material parameter in the said step (3) comprises density, specific heat capacity, the coefficient of heat conduction, porosity and the permeability of resin and fibre reinforcement, and said material parameter also comprises the material constant of Viscosity Model and curing reaction model.

8. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1 is characterized in that the material constant of Viscosity Model comprises A ₁, A ₂, E ₁, E ₂, m, n, a ₀, b ₀, the material constant of said curing reaction model comprises B, T _b, C ₁, C ₂, α _Gel

9. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; Said step (7) further comprises: according to the velocity field that calculates; Adopt the flow network analytic approach; Select the shortest required filling time of the control volume of each node on whole thickness direction in the current flow front as time step, upgrade the filling rate of other node control volume in the flow front and the material property of diverse location, said material property comprises viscosity, specific heat, density and the coefficient of heat conduction.

10. the flow front Forecasting Methodology of the non-isothermal resin transfer moulding based on middle surface model as claimed in claim 1; It is characterized in that; Said step (9) further comprises: according to the loading time of each node on the grid; The morning and evening of contrast filling, the prediction dry spot distributes, i.e. the dry spot position that possibly occur.

Images