CN105313336A - Thin-shell 3D printing optimizing method - Google Patents

Abstract

The invention discloses a thin-shell 3D printing optimizing method. The thin-shell 3D printing optimizing method comprises the following three steps that firstly, the maximum thickness of each vertex position is determined; secondly, a model is segmented, intersections of zone boundaries are extended to obtain transition areas, a thickness parameter is adopted for each independent segment, and a thickness parameter is adopted for each vertex; and thirdly, according to needs, the thickness parameters are optimized, and an optimized three-dimensional model is obtained. Through the framework provided in the method, the 3D printable model which has enough structural strength and consumes few materials under specific usage scenarios can be worked out rapidly and automatically. The method is clear, the speed is high, and the result is robust. The method can be used for the fields of computer aided design, 3D printing model generation and the like.

Description

A kind of thin walled shell 3D prints optimization method

Technical field

The present invention relates to computer graphical and Computer Aided Design field, 3D printing technique field, particularly relate to thin walled shell 3D and print the method that optimization thickness parameter obtains printable entity.

Background technology

It is increase material manufacturing technology (AM that 3D prints (3 D-printing), being commonly called as additiveManufacturing), by Design of digital model, use dusty material (plastics or metal), adopt the shaping mode of superposition to generate 3D solid by hot-melt adhesion techniques.Intelligent digitalized technology plays the part of important function in 3D print procedure.Computer realize, can the intention realizing user and generate printable entity.Along with technology development, 3D printer can adopt multiple material, comprises powder, expected, metal etc., and the chocolate under special applications scene, human stem cell etc.Processed by layering, superpose shaping mode and realize manufacturing complex model.Based on above advantage, 3D printing technique is in industry manufacture, and Aero-Space, robot, the fields such as biomedicine are fast-developing.Also due to the reduction of 3D printer manufacturing cost, this technology is allowed to be promoted rapidly in vast domestic consumer.

In the process of modeling, user sets up grid model by 3D modeling software (as Maya, 3DSMax etc.).By this kind of method establishment model normally bounding surface model out, and normally very thin can not being directly used for of the entity generated like this prints.In order to address this problem, most of 3D print software adopts following two kinds of methods: 1, filled in curved surface hole, and overall curved surface is closed and generated entity; 2, a thickness is set to block mold, generate after equidistant surface is closed and form entity.But for method 1, significantly can change model surface structure, inner space Print All can be become entity by printing like this, can increase the duration of printing, consume more material in print procedure simultaneously.And for method 2, the too small meeting that thickness is determined makes the object printed be difficult to practical requirement, bear certain external force; Thickness crosses the overall outward appearance of conference impact.

Had some researchers to be absorbed in recent years structural strength that how research strengthens printing objects.STAVA and VANEK proposes, by improving body form and increasing the structural strength supporting and increase object, to reduce the material printing and consume, see STAVA simultaneously, O., VANEK, J., BENES, B., CARR, N., ANDM ˇ ECH, R.:Stressrelief:Improvingstructuralstrengthof3dprintable, 566,2012, ACMTrans.Graph.31,4,48:1 – 48:11.Wang etc. support by construction framework the balance that shape structure realizes robust construction and saving material.See WANG, W., WANG, T.Y., YANG, Z., LIU, L., TONG, X., TONG, W., DENG, J., CHEN, F., ANDLIU, X.2013.Cost-effectiveprintingof3dobjectswithskin-framest ructures.ACMTrans.Graph.32,6,177:1 – 177:10.But above method is inapplicable for the object of thin shell-like structure.Because shell shape object surfaces externally and internally is all visible, increasing the way supported can affect face shaping, is difficult to meet user demand when user designs simultaneously.

Summary of the invention

The invention provides a kind of thin walled shell 3D and print optimization method, the structural strength reached according to model and the needs of input can be realized, automatically calculate every dot thickness parameter, and generate printable airtight grid model.

In order to achieve the above object, the present invention is by the following technical solutions: a kind of thin walled shell 3D prints optimization method, comprises the following steps:

(1) according to three-dimensional grid model geometry, the maximum gauge that each vertex position can reach and minimum thickness is determined;

(2) according to the geometry of model, Fragmentation is carried out to model, obtain multiple points of panel region, expansion is carried out to intersection, burst zone boundary and obtains transitional region.The maximum gauge obtained according to step 1 and minimum thickness, in point panel region each independently sheet adopt a thickness parameter, a thickness parameter is adopted for each summit in transitional region;

(3) the stress intensity requirement in the use that reaches is needed according to model, bear the situation of external force and print the material used, use the emulation of particular finite element modeling and sensitivity analysis technology, by alternating iteration, the thickness parameter in step 2 is optimized, try to achieve the optimum thickness on each summit on grid model, and generate the threedimensional model after optimizing.

Further, determine each summit maximum gauge method in step (1), be specially: according to input three-dimensional grid model, represent with M, try to achieve its bounding box at x, y, z-axis is of a size of (a, b, c), then model each some maximum gauge is set as:

T _Max＝0.05*(a+b+c-min{a,b,c}-max{a,b,c})

The most minimal thickness that minimum thickness can print according to 3D printer is determined:

T _min＝0.5mm

According to model meshes geometric properties, thickness correction is carried out to every bit on M, adopt the distance field function d (X of averageddistancefunctions function computing grid model M, M), the distance of representation space mid point X to grid model M, point x is the point be positioned in model M, F (x, t) represents the integral curve by an x represent that some x points to the vector of F (x, t), N (x) represents x point normal orientation in model M.

In curved surface and normal direction equidirectional side, maximum gauge for:

In like manner with normal direction opposite side, maximum gauge for:

For some x, revised maximum gauge T on M _gx () is expressed as

$T_G\left(x\right)=2\min \{T_G^{-}\left(x\right),T_G^{+}\left(x\right),0.5T_{Max}\}$

Further, according to model geometric feature in step (2), adopt fuzzycut method to carry out burst to model, be specially:

Border between expansion burst: the number of vertex of input model represents with n, the iterations N of setting extended boundary _expand(at least carrying out once), in each iterative process, adds the summit in the monocycle neighborhood on summits all in transitional region to transitional region.Through N _expandafter iteration, the transitional region obtained is expressed as B _t;

For in point panel region each independently in the thickness parameter of sheet and transitional region the thickness parameter on each summit adopt following way to obtain:

In transitional region, number of vertices is n _t, a thickness parameter is given on each summit, represents with vectorial β.In β, the thickness parameter of the corresponding respective vertices of each component, is expressed as

After carrying out border extension, remaining part is divided into part not adjacent to each other, is called segmented areas, is expressed as B _s, non-conterminous block count is n _s, each independent segmented areas adopts same burst thickness, represents with vectorial α.The thickness in the corresponding corresponding sub-block region of each component of α, is expressed as

represent the upper thickness limit of α, represent the lower thickness limit of α.For α _iall vertex sets of being comprised by its point of panel region of thickness parameter excursion determine, namely

$T_{\min }<{\mathrm{\&alpha;}}_i<\min \left\{T_G\left(x\right)\right|x\&Element;{\mathrm{sx}}_{{\mathrm{\&alpha;}}_i}\}$

represent the upper thickness limit of β, represent the lower thickness limit of β.For β _ithickness parameter excursion by the model vertices of its correspondence determine, namely

$T_{\min }<{\mathrm{\&beta;}}_i<T_G\left(x_{{\mathrm{\&beta;}}_i}\right)$

Thus, can the thickness t on each summit on Confirming model M according to each summit thickness parameter β of burst thickness alpha and precomputation.Namely can in the hope of t for given α, β.

Further, step (3) comprises step:

(4.1) model M is split as secondary housing unit and line slab unit, accelerates computational speed Bulk stiffness matrix, thickness parameter α can be expressed as, the multinomial of β.

(4.2) thickness parameter α is set up by finite element equation, the impact that β change changes each vertex v onMises power σ of model

According to finite element equation, under external force F effect, according to Bulk stiffness matrix K _syscan the displacement U on each summit of solving model:

K _sysU＝F

For the thickness t at summit i place in model M _ithe impact changed global displacement produces can be expressed as:

$\frac{\&part;F}{\&part;t_i}=\frac{\&part;K_{sys}}{\&part;t_i}U+K_{sys}\frac{\&part;U}{\&part;t_i}$

Because external force remains unchanged, not along with thickness t _ichange, above formula can be reduced to:

$\frac{\&part;U}{\&part;t_i}=-K_{sys}^{-1}\frac{\&part;K_{sys}}{\&part;t_i}U$

In conjunction with vonMises computing formula, can try to achieve according to the change of displacement to thickness parameter

(4.3) set up optimization object function, by alternating iteration optimizing process, optimize α respectively, β.

Optimization aim is under the condition of external force is specified in model carrying, ensures that material does not rupture, meets structural strength, and the printed material required for maintenance is minimum and maintenance generation curved surface is smooth.Namely available following form is expressed as:

$\underset{{\mathrm{\&Delta;t}}_i^k}{\min }\underset{i=1}{\overset{n}{\&Sigma;}}(t_i^k+{\mathrm{\&Delta;t}}_i^k)s_i$

$\begin{array}{cc}s.t.& {\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{\max },\&ForAll;v\&Element;M,\end{array}$

$t^{\min }<t_i^k+{\mathrm{\&Delta;t}}_i^k<t_i^{\max },i=\mathrm{1...}n,$

Wherein s _irepresent summit i adjoin triangle surface area and control curve smooth degree, t ^minprint most minimal thickness by 3D printer to determine, w _irepresent the weight on its neighborhood summit, if the number of vertex adjacent with an i is N _ring, then σ _maxmaterial selected by printing determines.

be equivalent to the volume V of closed mould _m.

Solved by iterative optimization procedure, in kth time iterative process, stage 1 and solving of stage 2 are respectively:

Stage 1:

$\underset{{\mathrm{\&Delta;t}}_i^k}{\min }\underset{i=1}{\overset{n}{\&Sigma;}}({\mathrm{\&alpha;}}_i^k+{\mathrm{\&Delta;\&alpha;}}_i^k)s_i$

$\begin{array}{cc}s.t.& {\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{\max },\&ForAll;v\&Element;M,\end{array}$

${\mathrm{\&alpha;}}_i^{\min }<{\mathrm{\&alpha;}}_i^k+{\mathrm{\&Delta;\&alpha;}}_i^k<{\mathrm{\&alpha;}}_i^{\max },i=\mathrm{1...}{n}_s,$

In stage 1 solution procedure, fixing β value is constant, and the solution solving linear optimization problem is α ^k+ Δ α ^k, by solving Laplace's equation, obtain β ^k.

Stage 2:

$\underset{{\mathrm{\&Delta;t}}_i^k}{\min }\underset{i=1}{\overset{n}{\&Sigma;}}({\mathrm{\&beta;}}_i^k+{\mathrm{\&Delta;\&beta;}}_i^k)s_i$

$\begin{array}{cc}s.t.& {\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{\max },\&ForAll;v\&Element;M,\end{array}$

$t^{\min }<{\mathrm{\&beta;}}_i^k+{\mathrm{\&Delta;\&beta;}}_i^k<{\mathrm{\&beta;}}_i^{\max },i=\mathrm{1...}{n}_t,$

After the solving of stage 2, obtain β ^k+ Δ β ^k, same by solving Laplace's equation herein, try to achieve α ^k+1, as the initial value of k+1 iteration.

Stage 1 and stage 2 carry out successively, and circulation solves, until constraints meets (ensure that material does not rupture, meet structural strength, the printed material required for maintenance is minimum and maintenance generation curved surface is smooth), and optimization aim V _mnot in decline, then stop optimizing process, obtain final result t ^f.

According to the thickness parameter solving each summit on the grid M that obtains according to the position on integral curve, each to both sides place obtains generation surface location, the threedimensional model after being namely optimized.

Beneficial effect of the present invention is:

1, block mold is divided into point panel region and transitional region, decreases optimized variable, make optimization more easy;

2, adopt special secondary housing unit and the combination of linear plane stress element to carry out analog simulation, accelerate computational speed;

3, use sensitivity analysis technology that optimization problem is reduced to linear programming, adopt iterative optimization techniques further, obtain the feasible solution of this optimization problem; Under the prerequisite not affecting face shaping, material consumption is minimized, meet user demand when user designs simultaneously, comprise stressing conditions etc.

Accompanying drawing explanation

Fig. 1 is technical scheme flow chart of the present invention.

Fig. 2 is initial model schematic diagram of the present invention;

Fig. 3 is zoning schematic diagram;

Fig. 4 is for optimizing rear model result figure.

Detailed description of the invention

As shown in Figure 1, one, for thin walled shell 3D print thickness parameter optimization method, comprises following three steps, and the maximum gauge that each vertex position can reach determined by the three-dimensional grid model geometry that 1, as required 3D print; 2, according to the geometry of model, Fragmentation is carried out to model, obtains a point panel region, expansion is carried out to intersection, zone boundary and obtains transitional region, for in point panel region each independently sheet adopt a thickness parameter, a thickness parameter is adopted for each summit in transitional region; 3, the user demand reached is needed according to model, bear the situation of external force and print the material used, use the emulation of particular finite element modeling and sensitivity analysis technology, by alternating iteration technology, the thickness that optimizing each summit on grid model needs to reach also generates the threedimensional model closed.

Below in conjunction with embodiment and accompanying drawing 2-4, the present invention is described in detail.

Now specifically introduce three steps of this method:

1) geometry of the three-dimensional grid model of 3D printing as required, determine that the method for the maximum gauge that each vertex position can reach is as follows:

According to user input three-dimensional grid, represent with M, try to achieve its bounding box at x, y, z-axis is of a size of (a, b, c).Due to the distortion that excessive thickness can cause FEM model to emulate, according to the maximum gauge that the size estimation shell shape object of block mold can reach, then model each some maximum gauge is set as:

T _Max＝0.05*(a+b+c-min{a,b,c}-max{a,b,c})

The most minimal thickness that minimum thickness can print according to 3D printer is determined:

T _min＝0.5mm

According to model meshes geometric properties, thickness correction is carried out to every bit on M, adopt the distance field function d (X of averageddistancefunctions function computing grid model M, M), the distance of representation space mid point X to grid model M, point x is the point be positioned in model M, F (x, t) represents the integral curve by an x represent that some x points to F (x, t) vector, N (x) represents x point normal orientation in model M, see PengJianbo, Kristjansson, DanielZorin, Denis (2004) .Interactivemodelingoftopologicallycomplexgeometricdetai l.ACMTransactionsonGraphics (TOG).

In curved surface and normal direction equidirectional side, maximum gauge for:

In like manner with normal direction opposite side, maximum gauge for:

For some x on M according to the maximum gauge T of geometry computations _gx () is expressed as

$T_G\left(x\right)=2\min \{T_G^{-}\left(x\right),T_G^{+}\left(x\right),0.5T_{Max}\}$

Adopt averageddistancefunctions can change integration direction according to the geometric properties of grid due to the integral curve of this function herein, the integral curve between difference can not selfing, effectively avoids the impact that the concavo-convex change of curved surface causes.Method according to normal orientation structure equidistant surface then cannot avoid the selfing situation at concave recess.T _gx the value of () is to ensure that the correct of limit element artificial module is suitable for.

2) according to model geometric feature, fuzzycut method is adopted to carry out burst to model, with reference to KatzSagi, TalAyellet, 2003, Hierarchicalmeshdecompositionusingfuzzyclusteringandcuts, ACMTransactionsonGraphics (TOG).Burst result as shown in Figure 3.Composition graphs describes

The number of the cut zone obtained can be specified by user or program judges automatically.On the basis of zoning, user from main regulation, by change default parameters, can expand or reduce respective regions.

Border between expansion burst, the number of vertex of input model represents with n, the iterations N of setting extended boundary _expand(at least carrying out once), in each iterative process, adds the monocycle neighborhood inner vertex on summits all in transitional region to transitional region.Through N _expandafter iteration, the transitional region obtained is expressed as B _t, in transitional region, number of vertices is n _t, a thickness parameter is given on each summit, represents with vectorial β.The thickness parameter of the corresponding respective vertices of each component of β, is expressed as after carrying out border extension, remaining part is divided into part not adjacent to each other, is called that segmented areas is expressed as B _s, non-conterminous block count is n _s, each independent segmented areas adopts same thickness, represents with vectorial α.The thickness in the corresponding corresponding sub-block region of each component of α, is expressed as

represent α upper thickness limit, represent the lower limit of α thickness parameter.For α _iall vertex sets of being comprised by its point of panel region of thickness parameter excursion determine, namely

$T_{\min }<{\mathrm{\&alpha;}}_i<\min \left\{T_G\left(x\right)\right|x\&Element;{\mathrm{sx}}_{{\mathrm{\&alpha;}}_i}\}$

represent β upper thickness limit, represent the lower limit of β thickness parameter.For β _ithickness parameter excursion by the model vertices of its correspondence determine, namely

$T_{\min }<{\mathrm{\&beta;}}_i<T_G\left(x_{{\mathrm{\&beta;}}_i}\right)$

This is had to set up α, the correspondence of the thickness t on each summit on β to model M.Namely can in the hope of t for given α, β, vice versa.

Adopt fuzzycut to divide surface methodology and effectively can utilize geological information.The effect of segmented areas is to keep the complete of model global feature, reduce the free degree optimizing part, accelerate algorithmic statement, the effect of transitional region is the segmented areas that smooth link is different, can not observe the notable difference because varied in thickness causes on mode shape.

3) Fig. 2 is initial model schematic diagram, and as can be seen from the figure, the external force that model need bear is 5N, and in the present embodiment, the elastic modelling quantity printing the material used represents with E, and Poisson's ratio ν represents.

Adopt triangle quadratic housing unit and linear plane stress element built-up pattern, Bulk stiffness matrix can be expressed as thickness parameter α, the multinomial of β.Adopt sensitivity analysis technology, set up thickness parameter α by finite element equation, the impact that β change changes each vertex v onMises power σ of model.

Vector t represents the thickness on each summit in model M, then for unit i, and the thickness t of j, m ^ecan be expressed as:

$t^e=\frac{t_i+t_j+t_m}{3}$

For typical triangular element, node i, j, m number counterclockwise, the expression of usable floor area coordinate, then its shape function N ^pcan be expressed as:

$N^P=\left[\begin{array}{cccccc}{N}_i^P& 0& N_j^P& 0& N_m^P& 0\\ 0& N_i^P& 0& N_j^P& 0& N_m^P\end{array}\right]$

Wherein, can be expressed as:

$N_i^P=\frac{a_i+b_{i}x+c_{i}y}{2\mathrm{\&Delta;}}$

Wherein, a _i, b _i, c _icomputing formula is as follows:

$\left\{\begin{array}{c}{a}_i=x_i{y}_m-x_m{y}_j\\ b_i=y_j-y_m\\ c_i=x_m-x_j\end{array}\right.$

Other coefficient can be tried to achieve by the sequential loop of subscript i, j, m, wherein

$2\mathrm{\&Delta;}=\det \left|\begin{array}{ccc}1& x_i& y_i\\ 1& x_j& y_j\\ 1& x_m& y_m\end{array}\right|$

Thus can compute matrix B ^p,

$B^P=\left[\begin{array}{ccc}{B}_i^P& B_j^P& B_m^P\end{array}\right]$

$B_i^P=\left[\begin{array}{cc}\frac{\&part;N_i^P}{\&part;x}& 0\\ 0& \frac{\&part;N_i^P}{\&part;y}\\ \frac{\&part;N_i^P}{\&part;y}& \frac{\&part;N_i^P}{\&part;x}\end{array}\right]=\frac{1}{2\mathrm{\&Delta;}}\left[\begin{array}{cc}{b}_i& 0\\ 0& c_i\\ c_i& b_i\end{array}\right]$

According to the character of material and plane stress element, matrix D ^pbe expressed as:

$D^P=\frac{E}{1-\mathrm{\&nu;}}\left[\begin{array}{ccc}1& \mathrm{\&nu;}& 0\\ \mathrm{\&nu;}& 1& 0\\ 0& 0& \frac{1-\mathrm{\&nu;}}{2}\end{array}\right]$

According to finite element theory, the thickness t of plane stress element ijm ^erepresent, then stiffness matrix be expressed as

$K_{ijm}^P=\&Integral;\&Integral;{\left(B^P\right)}^T{D}^P{B}^P{t}^{e}dxdy$

For secondary housing unit, its computing formula and plane stress element similar, its shape function N ^sbe expressed as:

$N^S=\left[\begin{array}{c}{P}_1-P_4+P_6+2(P_7-P_9)\\ -b_j(P_9-P_6)-b_m{P}_7\\ -c_j(P_9-P_6)-b_m{P}_7\\ P_2-P_5+P_4+2(P_8-P_7)\\ -b_m(P_7-P_4)-b_i{P}_8\\ -c_m(P_7-P_4)-c_i{P}_8\\ P_3-P_6+P_5+2(P_9-P_8)\\ -b_i(P_8-P_5)-b_j{P}_9\\ -c_i(P_8-P_5)-c_j{P}_9\end{array}\right]$

$P^T=\left[\begin{array}{c}{L}_i\\ L_j\\ L_m\\ L_i{L}_j\\ L_j{L}_m\\ L_m{L}_i\\ L_i^2{L}_j+\frac{1}{2}{L}_i{L}_j{L}_m\left(3\right(1-{\mathrm{\&mu;}}_m)L_i-(1+3{\mathrm{\&mu;}}_m)L_j+(1+{\mathrm{\&mu;}}_m\left)L_m\right)\\ L_j^2{L}_m+\frac{1}{2}{L}_i{L}_j{L}_m\left(3\right(1-{\mathrm{\&mu;}}_i)L_j-(1+3{\mathrm{\&mu;}}_i)L_m+(1+3{\mathrm{\&mu;}}_i\left)L_i\right)\\ L_m^2{L}_i+\frac{1}{2}{L}_i{L}_j{L}_m\left(3\right(1-{\mathrm{\&mu;}}_j)L_m-(1+3{\mathrm{\&mu;}}_j)L_j+(1+3{\mathrm{\&mu;}}_j\left)L_j\right)\end{array}\right]$

Wherein L _i, L _j, L _mbe expressed as,

$L_i=\frac{a_i+b_{i}x+c_{i}y}{2\mathrm{\&Delta;}},L_j=\frac{a_j+b_{j}x+c_{j}y}{2\mathrm{\&Delta;}},L_m=\frac{a_m+b_{m}x+c_{m}y}{2\mathrm{\&Delta;}}$

Wherein l _i, l _j, l _mrepresent limit i, the length of j, m, then μ in above formula _i, μ _j, μ _mbe expressed as:

${\mathrm{\&mu;}}_i=\frac{l_m^2-l_j^2}{l_i},{\mathrm{\&mu;}}_j=\frac{l_i^2-l_m^2}{l_j},{\mathrm{\&mu;}}_m=\frac{l_j^2-l_i^2}{l_m}$

Thus can compute matrix B ^s,

$B^S=\left[\begin{array}{cc}\frac{\&part;}{\&part;x}& 0\\ 0& \frac{\&part;}{\&part;y}\\ \frac{\&part;}{\&part;y}& \frac{\&part;}{\&part;x}\end{array}\right]\left[\begin{array}{c}\frac{\&part;}{\&part;x}\\ \frac{\&part;}{\&part;y}\end{array}\right]N^S$

According to the character of material and secondary housing unit, matrix D ^sbe expressed as:

$D^S=\frac{E{\left(t^e\right)}^3}{12(1-{\mathrm{\&nu;}}^2)}\left[\begin{array}{ccc}1& \mathrm{\&nu;}& 0\\ \mathrm{\&nu;}& 1& 0\\ 0& 0& \frac{1-\mathrm{\&nu;}}{2}\end{array}\right]$

According to finite element theory, the thickness t of plane stress element ijm ^erepresent, then stiffness matrix be expressed as

$K_{ijm}^S=\&Integral;\&Integral;\&Integral;{\left(B^S\right)}^T{D}^S{B}^{S}dxdydz$

Therefore at this, each unit is combined by secondary housing unit and linear plane stress element, the stiffness matrix K of its unit ijm _elebe expressed as:

$K_{ele}=\left(\begin{array}{ccc}{K}_{ijm}^P& 0& 0\\ 0& K_{ijm}^S& 0\\ 0& 0& K_z\end{array}\right)$

Wherein K _z=1e-8, makes whole equation to solve.K different herein _zvalue can affect the effect of finite element simulation significantly.Excessive K _zoverall unit can be made really up to the mark, too small K _zoverall equation intangibility can be made.

Above-mentioned stiffness matrix supposes that ijm is in xy plane, when 3D solid is applied, needs to solve spin matrix T _ijm, then under xyz coordinate system,

$K_{ele}^{global}=T_{ijm}{K}_{ele}{T}_{ijm}^T$

According to finite element theory, by all elements assemble and obtain model Bulk stiffness matrix K _sys.

In conjunction with above formula, K _syscan be explicit be expressed as t ^emultinomial.Due to t ^ethere are mapping relations with α, β, be namely equivalent to K _syscan be explicit be expressed as α, the multinomial of β.

According to finite element equation, under external force F effect, according to Bulk stiffness matrix K _syscan the displacement U on each summit of solving model:

K _sysU＝F

According to sensitivity analysis principle, for the thickness t at summit i place in model M _ithe impact changed global displacement produces can be expressed as:

$\frac{\&part;F}{\&part;t_i}=\frac{\&part;K_{sys}}{\&part;t_i}U+K_{sys}\frac{\&part;U}{\&part;t_i}$

Because external force remains unchanged, not along with thickness t _ichange, above formula can be reduced to:

$\frac{\&part;U}{\&part;t_i}=-K_{sys}^{-1}\frac{\&part;K_{sys}}{\&part;t_i}U$

In conjunction with vonMises computing formula, can try to achieve according to the change of displacement to thickness parameter

If do not adopt sensitivity analysis method, the problems referred to above are Solution of Nonlinear Optimal Problems, solve comparatively slow and are difficult to obtain convergence solution.If do not adopt alternating iteration optimizing process, then because variable weight difference in segmented areas and transitional region is comparatively large, and constraint excessively can cause Optimization Solution part unstable.

Optimization problem is converted into linear programming problem by employing sensitivity analysis method, is easy to solve accelerate solution procedure, and in finite element simulation process, reduces calculating, accelerate solving speed.

In implementation procedure, K _syscan be explicit be expressed as α, the polynomial character of β significantly improves computational speed.K _sysbe stored as α, the multinomial of β, need the α gone out each iterative at optimizing process, β recalculates K _sys, then only α, β need be multiplied by corresponding coefficient.In addition, in calculating time also can obtain exact solution fast.The key that this character is set up is to choose aforesaid secondary housing unit and the combination of linear plane stress element.If do not adopt this kind of unit, then cannot obtain this character in computational process, can computational speed be had a strong impact on.

Therefore set up optimization object function, by alternating iteration optimizing process, optimize α respectively, β.

Optimization aim is under the condition of external force is specified in model carrying, ensures that material does not rupture, and the printed material required for maintenance is minimum and maintenance generation curved surface is smooth.Namely available following form is expressed as:

$\underset{{\mathrm{\&Delta;t}}_i^k}{\min }\underset{i=1}{\overset{n}{\&Sigma;}}(t_i^k+{\mathrm{\&Delta;t}}_i^k)s_i$

$\begin{array}{cc}s.t.& {\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{\max },\&ForAll;v\&Element;M,\end{array}$

$t^{\min }<y_i^k+{\mathrm{\&Delta;t}}_i^k<t_i^{\max },i=\mathrm{1...}n,$

Wherein s _irepresent summit i adjoin triangle surface area and control curve smooth degree, t ^minprint most minimal thickness by 3D printer to determine, w _irepresent the weight on its neighborhood summit, if the number of vertex adjacent with an i is N _ring, σ _maxdetermined by printed material.Then be equivalent to the volume V of closed mould _m.

Try to achieve convergence solution in order to what realize fast and stable, solved by iterative optimization procedure, in kth time iterative process, stage 1 and solving of stage 2 are respectively:

Stage 1:

$\underset{{\mathrm{\&Delta;t}}_i^k}{\min }\underset{i=1}{\overset{n}{\&Sigma;}}({\mathrm{\&alpha;}}_i^k+{\mathrm{\&Delta;\&alpha;}}_i^k)s_i$

$\begin{array}{cc}s.t.& {\mathrm{\&sigma;}}_v^k\end{array}+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{\max },\&ForAll;v\&Element;M,$

${\mathrm{\&alpha;}}_i^{\min }<{\mathrm{\&alpha;}}_i^k+{\mathrm{\&Delta;\&alpha;}}_i^k<{\mathrm{\&alpha;}}_i^{\max },i=\mathrm{1...}{n}_s,$

In stage 1 solution procedure, fixing β value is constant, and the solution solving linear optimization problem is α ^k+ Δ α ^k, by solving Laplace's equation, obtain β ^k.

Stage 2:

$\underset{{\mathrm{\&Delta;t}}_i^k}{\min }\underset{i=1}{\overset{n}{\&Sigma;}}({\mathrm{\&beta;}}_i^k+{\mathrm{\&Delta;\&beta;}}_i^k)s_i$

$\begin{array}{cc}s.t.& {\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{\max },\&ForAll;v\&Element;M,\end{array}$

$t^{\min }<{\mathrm{\&beta;}}_i^k+{\mathrm{\&Delta;\&beta;}}_i^k<{\mathrm{\&beta;}}_i^{\max },i=\mathrm{1...}{n}_t,$

After the solving of stage 2, obtain β ^k+ Δ β ^k, same by solving Laplace's equation herein, try to achieve α ^k+1, as the initial value of k+1 iteration.

Stage 1 and stage 2 carry out successively, and circulation solves, until constraints meets, and optimization aim V _mnot in decline, then stop optimizing process, obtain final result t ^f.

If do not adopt iteration optimization to solve part, then solution procedure does not restrain, and cannot obtain feasible solution.

According to the thickness parameter solving each summit on the grid M that obtains according to the position on integral curve, each as both sides place obtains generation surface location, namely obtains the threedimensional model closed.

Claims (4)

1. thin walled shell 3D prints an optimization method, it is characterized in that: comprise the following steps:

(1) according to three-dimensional grid model geometry, the maximum gauge that each vertex position can reach and minimum thickness is determined;

(2) according to the geometry of model, Fragmentation is carried out to model, obtain multiple points of panel region, expansion is carried out to intersection, burst zone boundary and obtains transitional region.The maximum gauge obtained according to step 1 and minimum thickness, in point panel region each independently sheet adopt a thickness parameter, a thickness parameter is adopted for each summit in transitional region;

(3) the stress intensity requirement in the use that reaches is needed according to model, bear the situation of external force and print the material used, use the emulation of particular finite element modeling and sensitivity analysis technology, by alternating iteration, the thickness parameter in step 2 is optimized, try to achieve the optimum thickness on each summit on grid model, and generate the threedimensional model after optimizing.

2. method according to claim 1, it is characterized in that: in step (1), determine each summit maximum gauge method, be specially: according to input three-dimensional grid model, represent with M, try to achieve its bounding box at x, y, z-axis is of a size of (a, b, c), then model each some maximum gauge is set as:

T _Max＝0.05*(a+b+c-min{a,b,c}-max{a,b,c})

The most minimal thickness that minimum thickness can print according to 3D printer is determined:

T _min＝0.5mm

According to model meshes geometric properties, thickness correction is carried out to every bit on M, adopt the distance field function d (X of averageddistancefunctions function computing grid model M, M), the distance of representation space mid point X to grid model M, point x is the point be positioned in model M, F (x, t) represents the integral curve by an x represent that some x points to the vector of F (x, t), N (x) represents x point normal orientation in model M.

In curved surface and normal direction equidirectional side, maximum gauge for:

In like manner with normal direction opposite side, maximum gauge for:

For some x, revised maximum gauge T on M _gx () is expressed as

$T_G\left(x\right)=2min\{T_G^{-}\left(x\right),T_G^{+}\left(x\right),0.5T_{Max}\}$

3. method according to claim 1, is characterized in that: according to model geometric feature in step (2), adopts fuzzycut method to carry out burst to model, is specially:

Border between expansion burst: the number of vertex of input model represents with n, the iterations N of setting extended boundary _expand(at least carrying out once), in each iterative process, adds the summit in the monocycle neighborhood on summits all in transitional region to transitional region.Through N _expandafter iteration, the transitional region obtained is expressed as B _t;

For in point panel region each independently in the thickness parameter of sheet and transitional region the thickness parameter on each summit adopt following way to obtain:

In transitional region, number of vertices is n _t, a thickness parameter is given on each summit, represents with vectorial β.In β, the thickness parameter of the corresponding respective vertices of each component, is expressed as

After carrying out border extension, remaining part is divided into part not adjacent to each other, is called segmented areas, is expressed as B _s, non-conterminous block count is n _s, each independent segmented areas adopts same burst thickness, represents with vectorial α.The thickness in the corresponding corresponding sub-block region of each component of α, is expressed as

represent the upper thickness limit of α, represent the lower thickness limit of α.For α _iall vertex sets of being comprised by its point of panel region of thickness parameter excursion determine, namely

$T_{min}<{\mathrm{\&alpha;}}_i<min\left\{T_G\left(x\right)\right|x\&Element;{\mathrm{sx}}_{{\mathrm{\&alpha;}}_i}\}$

represent the upper thickness limit of β, represent the lower thickness limit of β.For β _ithickness parameter excursion by the model vertices of its correspondence determine, namely

$T_{min}<{\mathrm{\&beta;}}_i<T_G\left(x_{{\mathrm{\&beta;}}_i}\right)$

Thus, can the thickness t on each summit on Confirming model M according to each summit thickness parameter β of burst thickness alpha and precomputation.Namely can in the hope of t for given α, β.

4. method according to claim 1, is characterized in that: step (3) comprises step:

(4.1) model M is split as secondary housing unit and line slab unit, accelerates computational speed Bulk stiffness matrix, thickness parameter α can be expressed as, the multinomial of β.

(4.2) thickness parameter α is set up by finite element equation, the impact that β change changes each vertex v onMises power σ of model

According to finite element equation, under external force F effect, according to Bulk stiffness matrix K _syscan the displacement U on each summit of solving model:

K _sysU＝F

For the thickness t at summit i place in model M _ithe impact changed global displacement produces can be expressed as:

$\frac{\&part;F}{\&part;t_i}=\frac{\&part;K_{sys}}{\&part;t_i}U+K_{sys}\frac{\&part;U}{\&part;t_i}$

Because external force remains unchanged, not along with thickness t _ichange, above formula can be reduced to:

$\frac{\&part;U}{\&part;t_i}=-K_{sys}^{-1}\frac{\&part;K_{sys}}{\&part;t_i}U$

In conjunction with vonMises computing formula, can try to achieve according to the change of displacement to thickness parameter

(4.3) set up optimization object function, by alternating iteration optimizing process, optimize α respectively, β.

Optimization aim is under the condition of external force is specified in model carrying, ensures that material does not rupture, meets structural strength, and the printed material required for maintenance is minimum and maintenance generation curved surface is smooth.Namely available following form is expressed as:

$\underset{{\mathrm{\&Delta;t}}_i^k}{min}\underset{i=1}{\overset{n}{\&Sigma;}}(t_i^k+{\mathrm{\&Delta;t}}_i^k)s_i$

$s.t.{\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{max},\&ForAll;v\&Element;M,$

$t^{min}<t_i^k+{\mathrm{\&Delta;t}}_i^k<t_i^{\max },i=1...n,$

Wherein s _irepresent summit i adjoin triangle surface area and control curve smooth degree, t ^minprint most minimal thickness by 3D printer to determine, w _irepresent the weight on its neighborhood summit, if the number of vertex adjacent with an i is N _ring, then σ _maxmaterial selected by printing determines.

be equivalent to the volume V of closed mould _m.

Solved by iterative optimization procedure, in kth time iterative process, stage 1 and solving of stage 2 are respectively:

Stage 1:

$\underset{{\mathrm{\&Delta;t}}_i^k}{min}\underset{i=1}{\overset{n}{\&Sigma;}}({\mathrm{\&alpha;}}_i^k+{\mathrm{\&Delta;\&alpha;}}_i^k)s_i$

s.t. ${\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{max},\&ForAll;v\&Element;M,$

${\mathrm{\&alpha;}}_i^{\min }<{\mathrm{\&alpha;}}_i^k+{\mathrm{\&Delta;\&alpha;}}_i^k<{\mathrm{\&alpha;}}_i^{\max },i=1...n_s,$

In stage 1 solution procedure, fixing β value is constant, and the solution solving linear optimization problem is α ^k+ Δ α ^k, by solving Laplace's equation, obtain β ^k.

Stage 2:

$\underset{{\mathrm{\&Delta;t}}_i^k}{min}\underset{i=1}{\overset{n}{\&Sigma;}}({\mathrm{\&beta;}}_i^k+{\mathrm{\&Delta;\&beta;}}_i^k)s_i$

s.t. ${\mathrm{\&sigma;}}_v^k+\frac{\&part;{\mathrm{\&sigma;}}_v}{\&part;t}\left(t_i^k\right){\mathrm{\&Delta;t}}_i^k<{\mathrm{\&sigma;}}_{max},\&ForAll;v\&Element;M,$

$t^{min}<{\mathrm{\&beta;}}_i^k+{\mathrm{\&Delta;\&beta;}}_i^k<{\mathrm{\&beta;}}_i^{max},i=1...n_t,$

After the solving of stage 2, obtain β ^k+ Δ β ^k, same by solving Laplace's equation herein, try to achieve α ^k+1, as the initial value of k+1 iteration.

Stage 1 and stage 2 carry out successively, and circulation solves, until constraints meets (ensure that material does not rupture, meet structural strength, the printed material required for maintenance is minimum and maintenance generation curved surface is smooth), and optimization aim V _mnot in decline, then stop optimizing process, obtain final result t ^f.

According to the thickness parameter solving each summit on the grid M that obtains according to the position on integral curve, each to both sides place obtains generation surface location, the threedimensional model after being namely optimized.