close

Вход

Забыли?

вход по аккаунту

?

906

код для вставкиСкачать
INTERNATIONAL JOURNAL FOR NUMERICAL METHODS IN ENGINEERING, VOL.
39, 2837-2862 (1996)
SENSITIVITY ANALYSIS FORMULATION FOR
THREE-DIMENSIONAL CONDUCTION HEAT TRANSFER
WITH COMPLEX GEOMETRIES USING A
BOUNDARY ELEMENT METHOD
SEONG JIN PARK AND TAI HUN KWON
Department of Mechanical Engineering, Pohang University of Science and Technology, Pohang 790-784, South Korea
SUMMARY
A special boundary integral formulation had been proposed to analyse many engineering problems of conduction heat transfer in complex three-dimensional geometries (closely spaced surface and circular hole in
infinite domain or simple modification of it) by Rezayat and Burton. One example of such geometries is the
mold sets in the injection molding process. In this paper, an efficient and accurate approach for the design
sensitivity analysis (DSA) is presented for these kinds of problems in the similar complex geometries using
the direct differentiation approach (DDA) based on the above special boundary integral formulation. The
present approach utilizes the implicit differentiation of the boundary integral equations with respect to the
design variables (radii and locations of circular holes) to yield the sensitivity equations. A sample problem
(heat transfer of injection molding cooling system) is solved to demonstrate the accuracy of the present
sensitivity analysis formulation. Although the techniques introduced here are applied to a particular problem
in heat transfer of injection molding cooling system, their potential application is quite broad.
KEY WORDS: special boundary integral formulation; design sensitivity analysis (DSA); direct differentiation approach
(DDA); three-dimensional conduction heat transfer
1. INTRODUCTION
Many engineering problems involve heat transfer, stress analysis, fluid mechanics, electronics,
acoustics, etc., in regions with complex three-dimensional geometries (closely spaced surface and
circular hole in infinite domain or simple modification of it). Some examples of such geometries
are the mold sets for the process of injection molding of plastics and the region of fracture
mechanics problem for crack tip stress analysis, etc. The present study is concerned about a
particular heat transfer in the injection molding process. In this process, molten plastic is injected
into the cavity of the closed mold, solidified in place, and then mechanically ejected. Heat is
removed from molds mostly by water flowing in channels drilled in the mold.' Rezayat and
Burton2 proposed the special boundary integral formulation to analyse the heat transfer in injection
molds of these complex geometries. The proposed formulation consists of the boundary formula,
its gradient and the exterior formula to calculate the heat transfer between a thermoplastic part
and the corresponding mold.
Design sensitivity analysis (DSA) will ultimately play an important role in the design optimization process after the above analysis has been achieved. Design sensitivities represent the gradients
of variables of interest with respect to design variables (DVs). This information is required by
non-linear optimization algorithms with first-order methods to determine the search direction in
CCC 0029-598 11961162837-26
0 1996 by John Wiley & Sons, Ltd.
Received 10 April 1995
Revised 22 September 1995
2838
S. J. PARK AND T. H. KWON
which to proceed in altering the existing design of an object as it evolves towards an optimal
configuration. These algorithms change the configuration of the object incrementally and require a
number of such increments or iterations to obtain the optimal configuration. Inaccurate evaluations
of design sensitivities may lead to an increased number of iterations and may even render the effort
unsuccessful. The determination of accurate design sensitivities is thus of crucial significance and
has attracted the attention of numerous researchers. For example, there are some research works
on this subject of optimization problems in manufacturing processes of plastic materials: Barone
and Caulk3 and Matsumoto and c o - w o r k e r ~for
~ ~thermal
~
design of injection molds, Forcucci and
Kwon6 and Matsumoto and co-workers7~*for thermal design of compression molds, etc.
The Boundary Element Method (BEM) has special advantages when applied to design sensitivity
analysis as pointed out by, among others, Barone and Yang' and Saigal and co-workers.1° Two
methods, namely the direct differentiation approach (DDA) and the adjoint structure approach
(ASA), emerge as the most significant ones for the evaluation of design sensitivities. Both the
approaches have recently been developed in the context of boundary element method for elastic
systems: for examples, see the references by Barone and Yang' for two-dimensional systems, Saigal
and co-workers" for axisymmetric systems, etc. These approaches have more recently been applied
to several thermal system such as steady-state diffusion problems: for examples, Reference 11 for
two-dimensional and axisymmetric heat diffusion problems, Reference 12 for three-dimensional
heat conduction problems, etc.
In this paper, an efficient and accurate approach for the design sensitivity analysis is proposed for
the heat transfer of injection molding process using the direct differentiation approach based on the
special boundary integral formulation proposed by Rezayat and Burton.2 The major advantage of
this approach lies in the fact that the accurate sensitivity analysis results are obtained simultaneously
for all design variables. We can generally consider radii and locations of circular holes (cooling
channels in an injection mold) as design variables in many engineering problem^.'^.'^ Using an
example, the numerical results of the sensitivity analysis simulation developed in the present
study are amply discussed in terms of accuracy and computational efficiency. The accuracy of
the present formulation is demonstrated through favourable comparisons between the sensitivity
analysis results obtained from the present direct differentiation approach and those obtained by a
numerical derivative approach using the forward finite difference method (FDM) in this example.
2. BOUNDARY INTEGRAL FORMULATION FOR HEAT TRANSFER ANALYSIS
2.1. Governing equation
The configuration of a very simple injection mold is schematically shown in Figures 1 and 2.
0, Sp, I?, Sc, SE,and ii denote the mold region, the part surface, the mid-surface of part, the
cooling channel surfaces, the mold exterior surfaces, and the outward unit normal vector in part,
respectively. We can describe the part surface by only the mid-surface of the part I? instead of
the whole surface (positive and negative surface as indicated by Sp' and Si- in Figure 2.) because
plastic parts are usually thin (sheet-like). It is interesting to note that in most injection molding
applications, heat loss through exterior surface is very small (typically less than 5 per cent). In
such cases, even though the normal treatment of the mold exterior surface can be applied, the
mold exterior surface may be treated as an infinite adiabatic sphere without significant change of
re~u1ts.l~
It may be noted that this approximation does not require modelling of the mold exterior
surface and thus it leads to substantial savings in the memory requirement and CPU time. Thus,
the given geometry consists of closely spaced surfaces and circular holes in an infinite domain.
The governing differential equation for the heat transfer of injection mold cooling stage under the
2839
SENSITIVITY ANALYSIS FORMULATION
Figure 1. A simple example of mold bounded by
r, SC and SE
~moldaurfaas,
region n
Figure 2. Schematic diagram of a cross-sectional mold for notations
so-called cycle-averaged concept’3%
l 4 may be written as (steady-state heat conduction equation)
V 2 T = 0 in 0
(1)
+ +
where T is the temperature. The boundary conditions on the boundary S = Sp SC SE are given
as
-k-
aT
= h ( T - Tb) on SC
an
(2 )
aT
on SE
an
where h, k , h, and Tb denote the outward unit normal vector, the thermal conductivity of the mold
material, the heat transfer coefficient in the cooling channels and the coolant bulk temperature,
respectively, and (aT/an)o is assumed to be specified varying with the location in the part. (As
a matter of fact, (dT/an)o is to be determined as a part of solution of the system. However, the
-=0
2840
S. J. PARK AND T. H. KWON
following method for the sensitivity analysis is valid regardless of this assumption concerning
(aT/an)O.)
2.2. Boundary integral formulation for analysis
A standard boundary element formulation for the three-dimensional Laplace equation given by
equation (1 ) governing the conduction heat transfer in a mold leads to”
Here x and 5 are position vectors in space, r = 15 - X I , and a denotes a solid angle formed by
the boundary surface at the position x . Note that ci = 211 on the smooth boundary surface, and
a = 471 or 0 at internal points or external points of the domain, respectively.
For any two closely spaced surfaces such as the mold cavity surfaces, equation (3) leads to
redundancy in the final system of linear algebra equations; in such cases, a modified procedure
as described in Reference 2 needs to be used. According to this modification, a mid-surface is
considered rather than two closely spaced surfaces. For example, in the case of mold cavity surfaces
S p which consist of two closely spaced surfaces, Sp’ and S
,; it is represented by the cavity midsurface I? as shown in Figure 2. It might be noted that each mid-surface element has two degrees
of freedom (e.g. temperature on Sp’ and S;)
in order to derive the extra equation corresponding
to the additional degree of freedom, we also used an approach similar to that used for crack
analysis in fracture mechanics.2 For this purpose, a derivative of equation (3) with respect to the
normal direction, 3 ( i + in Figure 2), at the mid-surface element is taken. See Reference 2 for
further details of this modified approach for cavity surfaces.
For cooling channel (circular hole) surfaces (SC as shown in Figures 1 and 2), we use a special
formulation based on the line-sink approximation?. 3*16- I8 This formulation avoids discretization of
the circular channels along the circumferential direction and thus saves large amount of computer
memory and time.
When we approximate the exterior surface by an adiabatic sphere at infinite, the integral parts of
equation (3) and its gradient over SE do not need to be evaluated because the following integrals
hold:2,l9
[!g
-
$(!)TI
dS =41cT,
(4)
In equation ( 4 ) , T, denotes the temperature on the exterior surface (at infinite point).
The followings are the final BEM formulae for these modifications:
for a point x on the mid-surface of the part (I‘) we have the following pair of integral equations:
+ 4xT,
(6)
SENSITIVITY ANALYSIS FORMULATION
284 1
and for a point on the axis of the cylindrical segment of the cooling channels we have the following
integral equation:
In equations (6)-(8), N is the total number of cooling channels, uk is the radius of kth cooling
channel, 1 is the local co-ordinate in the axial direction of each cooling channel, and 0 is the local
circumferential co-ordinate of each cooling channel as shown in Figure 3. It is noted again that
( A + in Figure 2) in equation (7) is the outward unit normal vector in the plus plane at the source
point (x).
Figure 3. Co-ordinate transformation rule in one cooling channel element for sensitivity analysis with respect to
location of cooling channel
2842
S. J. PARK AND T. H. KWON
3. BOUNDARY INTEGRAL FORMULATION FOR SENSITIVITY ANALYSIS
3.1. Design variables
The boundary element method is also a strong tool in design sensitivity analyses when the
problem can be treated without the domain discretization and the design variables can be defined
on the boundary. In the mold cooling system design, the size and the arrangement of the channels
can be treated as the design variables defined only on the boundary. Therefore, the boundary
element method can be effectively applied to the sensitivity analyses of this kind of heat transfer
problems. The sensitivity values of the temperature and the heat flux with respect to these design
variables give us valuable information for the design optimization. Therefore, in this paper, we
consider the radii and the locations of the cooling channels as design variables. (The coolant
flow rate and the coolant bulk temperature can be also considered as design variables in injection
mold problems, but will not be discussed here.) In other engineering problems, the size and the
arrangement of the circular holes may be very interesting design variables for designers of those
systems.
3.2. Boundary integral formulation for sensitiviry analysis
The design sensitivity analysis consists of determination of the variation of an objective hnction
with respect to a variation of each design variable. In the present study, only the sizes and the
locations of the cooling channels are considered as design variables. To determine the sensitivity
of a heat transfer system with respect to each design variable, the implicit differentiation of equations (6)--(8) with respect to each design variable is performed. This leads to the following design
sensitivity equations for each design variable.
First, we consider the radius of jth cooling channel, u j , as a design variable:
for the part surface, the implicit derivatives of equations (6) and (7) with respect to the radius of
jth cooling channel, a,, give us the following pair of integral equations for the sensitivity analysis:
- aT+(x)
tl-
aaj
+ tl+-
aT-(x)
aa,
SENSITIVITY ANALYSIS FORMULATION
2843
(10)
for the ith cooling channel surface, the implicit derivative of equation (8) with respect to the radius
of jth cooling channel, u j , gives us the following integral equation for the sensitivity analysis:
+ 4K- a r m
dUj
(11)
In equations (9)-(11), the third integral terms appear because r is the function of the design
variable uj and the domain of integral (the surface of jth cooling channel) also changes.
Next, we consider the location of jth cooling channel as a design variable (a vector):
for the part surface, the implicit gradient of equations (6) and (7) with respect to the location of
the jth cooling channel gives us a pair of integral equations for the sensitivity analysis as follows
in a gradient form:
cr-VT+(x)
+ a+VT-(x)
+ 4aVTm
2844
S. J. PARK AND T. H. KWON
-
2
(a(:))
av an+
1
(VT' - V T - ) dS(Q
for the ith cooling channel surface, the implicit gradient of equation (8) with respect to the
location of the jth cooling channel gives us the following integral equation for the sensitivity
analysis, depending on the value of i:
when i # j ,
0=
A [:{ (g)
+ g)}
&
V
+ 4nV T ,
when i = j,
V(
-
1
( ! ) ( V T + - V T - ) dS(5)
SENSITIVITY ANALYSIS FORMULATION
2845
+ 4nVT,
(15)
In the above equations (12)-( 15), V t and 0, specifically denote gradient with respect to 5 and x,
respectively, depending on which is the design variable, where r = 15 - XI. Because the location is
a vector, the co-ordinate transformation rule should be introduced for the sensitivity analysis with
respect to the location of the cooling channel. The co-ordinate transformation rule between global
co-ordinate system and local co-ordinate system is illustrated in Figure 3 and briefly described in
Appendix I.
3.3. Solution procedure
Once the integrals on each element are calculated, the discretized boundary element formulae
for the analysis can be manipulated to obtain the following matrix equation:
where [Hi,] and [ Cij] are hnctions of geometry of boundary surfaces. Some of [Hi,] and [ Gi,]
are singular integral, which can be evaluated analytically as described in Appendix 11. Also, the
integrals over 8 in equations (6)-(8) can be evaluated in closed forms using the complete elliptic
integrals, See Appendix I1 for details. Other integrals can be evaluated by the Gaussian quadrature
integration method. Next, boundary conditions can be introduced into equation (16) to obtain a
system of linear algebraic equations:
[ ~ i j ]
={fiI
(17)
where [Aii] and { f i } reflect boundary conditions. For this particular problem, temperature T is
taken to be an unknown on each element with dT/dn being eliminated with the help of a mixed
boundary condition. This system of equations can be solved either by LU-decomposition through
Gaussian elimination or by an under-relaxation iterative method.
In similar manner, the discretized boundary element formulae for the sensitivity analysis can be
manipulated to obtain the following form of equation:
[ ~ i j ]{Tj,x) =
{fi.x) -
= {f;,x)
[ ~ i j , {~ T] j )
(18)
where X is a design variable and, 'X' stands for a partial derivative with respect to X. The
integrals in above equations can be evaluated in the similar manner to the above analysis case.
See Appendix I11 for the detailed integral formulae for various terms appearing in equations (9)(15). It should be noted that [Ai,] in equation (18) is exactly the same as that in equation (17)
for the analysis, and that only the forcing terms { f ; , x } are to be evaluated for each design
variable. Therefore, once { f i , x } are evaluated for all design variables, the sensitivity analysis
results, i.e. { T ~ J } can
, be determined simultaneously with multiple forcing matrices for all design variables, thus saving the computational time. This fact is the most fascinating point of the
present sensitivity analysis approach. In this regard, the major advantage of LU-decomposition
solver lies in the fact that the LU-factorization of the matrix [ A i j ] , formed during the analysis, can be multiply reused in this sensitivity analysis. On the other hand, an iterative method
2846
S. J. PARK AND T. H. KWON
which takes less computer memory than a direct solver can also efficiently evaluate {Tj,x}
for all design variables, thus would be very useful when one needs large number of boundary
elements.
4. RESULTS AND DISCUSSION
A simple injection mold cooling problem as shown in Figure 4 is presented here as an example
to demonstrate the accuracy and the computational efficiency of the sensitivity analysis method
formulated in the present study. The mold has two cooling channels, and the dimension of the
part is 200 x 200 x 2 (unit: mm). The global co-ordinate system was chosen such that the
x-co-ordinate is in the axial direction of the cooling channels, the y-co-ordinate is in the thickness
direction of the part, and the z-co-ordinate is in the spanwise direction of the part, as shown in
Figure 4. The sizes and the locations of cooling channel #1 and cooling channel #2 are made
intentionally unsymmetrically to show their effects on the cooling analysis. The discretization for
the boundary element analysis is shown in Figure 5. In Figure 5 , each asterisk represents a node.
We used triangular elements for the part surface and line elements for the cooling channel surface. Number of elements in the part and two cooling channels are 558 and 30, respectively. It
d i n g line H
(rsdlur: 0.6 cm)
3.5 cm
+ plane
t h i i w 0.2 an
oMin
COOliiline82
(radius: 0.4 c m
I
)
d
b
20 cm
Figure 4. Diagram for geometry of an example problem
SENSITIVITY ANALYSIS FORMULATION
2847
Figure 5. Boundary element mesh for the example
might be mentioned that each element has a constant temperature and a constant heat flux in
the boundary element analysis in the present study. The thermal conductivity of mold material
used in this example is 100 W/m .K.For this example, we used the following boundary conditions:
= 10 OC/m,
(a) for part surface: uniform temperature gradient (&'/an),
(b) for both cooling channels: Tb = 20 "C, h = lo4 W/m2 . K.
Figure 6 illustrates analysis results for the part surface temperature distribution in each plane.
The elapsed CPU time for this analysis is 0:28:19 (h:min:s) in the machine SUN SPARC 10 (22.9
MFLOPS/102 mips). (All computations were run on the same machine with a UNIX operating
system using the same Fortran 77 compiler and complier options.) We have used four significant
digits for each evaluated number (possible in single precision calculation). What we mean by
the number of significant digits, n, is that we used the following the maximum norm convergent
criteria in using an iterative method:
""I
TFew - T.Old
cold
1
'lo-"
where Tiis an unknown, i.e. temperature, in each element. The maximum surface temperature in
each plane is located in the position which is opposite to two cooling channels in the sense of
both y-co-ordinate and z-co-ordinate, respectively as indicated in Figure 6. Moreover, both the
minimum temperature and maximum temperature in the plus plane are higher than those of the
minus plane, thus cooling channel #1 has more dominant cooling effect than cooling channel #2,
as expected.
As far as the sensitivity analysis is concerned, we consider the radius, the y-co-ordinate and the
z-co-ordinate of each cooling channel as design variables because the effect of the x-co-ordinate is
negligible. Therefore, six design variables exist (three design variables for each cooling channel)
in this example. The elapsed CPU time for the sensitivity analysis using the direct differentiation
approach is 2:27:00 in the same machine, thus the total elapsed CPU time is 2:55:19 (one analysis
and one sensitivity analysis). In the sensitivity analysis, we have used three significant digits for
each evaluated number (possible in single precision calculation). In this kind of geometries, the
exact solution does not exist. Therefore, in order to check the accuracy of the sensitivity analysis
results, we should directly compare sensitivity analysis results obtained from the direct differentiation approach with those from another numerical method, namely a forward finite difference
method. Concerning the accuracy of the forward finite difference method, we have investigated the
effects of the number of significant digits used in the computation and the magnitude of difference
2848
S. J. PARK
AND T. H. KWON
max: 146.46
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
143.39
140.33
137.27
134.21
131.15
128.08
125.02
121.96
118.90
115.84
112.78
4: 109.71
3: 106.65
2: 103.59
1: 100.53
min: 97.47
max: 138.75
15: 136.06
14:
13:
12:
11:
10
9:
8:
7:
6:
5:
4:
3:
2:
1:
133.36
130.67
127.97
125.27
122.58
119.88
117.18
114.48
111.79
109.09
108.40
103.70
101.01
98.31
min: 95.61
Figure 6. Surface temperature distribution in part (unit: "C): (a) plus plane, (b) minus plane
in design variables (AX) on the sensitivity analysis results. Figure 7 illustrates those effects on
the sensitivity analysis with respect to a design variable, the radius of the cooling channel #1,
R1. In Figure 7, the ordinates of both plots (a) and (b) represent the root-mean-square value of
SENSITIVITY ANALYSIS FORMULATION
2849
Figure 7. Some criteria for finite difference method: (a) criterion on the number of the significant digits, (b) criterion
on the magnitude of design variable difference
d ( T + - T-)/dRI, i.e.
1
with the abscissa of plot (a) and (b) representing the number of significant digits and the magnitude
of difference, respectively. The magnitude of difference is set to be 0.01 percent (0.00006cm /
0.6cm) to obtain results in Figure 7(a), while the number of significant digits is set to be 10 for
Figure 7(b). According to Figure 7, it was found that at least six-digit significant figure and 0.01
per cent difference are required for accurate computations for the sensitivity analysis. Therefore,
2850
S. J. PARK AND T. H. KWON
max: -30.41
15: -31.54
14: -32.68
13: -33.81
12: -34.94
11: -36.08
10: -37.21
9: -38.35
8: -39.48
7: -40.62
6: -41.75
5: -42.89
4: -44.02
3: -45.15
2: -46.29
1: 47.42
min: -48.56
max: -38.70
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
-40.29
-41.89
-43.48
-45.08
46.67
-48.27
49.86
-51.46
-53.05
-54.66
-56.24
-57.04
-59.43
2: -61.03
1: -62.62
min: -64.22
Figure 8. Sensitivity results of part surface temperature by the direct differentiation approach with respect to radius of
cooling channel #1 (unit: OC/cm): (a) plus plane, (b) minus plane
we have used a double-precision calculation in order to satisfy the above requirements when the
forward finite difference method was used for the sensitivity analysis. The sensitivity analysis by
the forward finite difference method required the total elapsed CPU time of 10:14:23 (the elapsed
CPU time is 1:17:46 for each run of the analysis program and total 7 runs of analysis program
285 1
SENSITIVITY ANALYSIS FORMULATION
max: 30.63
15: -31.81
14: -32.98
13: -34.15
12: -35.33
11: -36.50
10: -37.68
9:
8:
7:
6:
38.85
-40.02
-41.20
-42.37
5: -43.55
4: 44.72
3: 45.90
2: 47.07
1: -48.24
min: -40.42
max: -39.20
15: -40.85
14: -42.50
13: -44.16
12: -45.81
11: -47.46
10: -40.11
9: -50.76
8: -52.41
7: -54.06
6: -55.71
5: 37.36
4: 4 . 0 1
3: -60.66
2: -62.31
1: -63.96
min: -65.61
Figure 9. Sensitivity results of part surface temperature by the forward finite difference method with respect to radius
of cooling channel #1 (unit: OCkm): (a) plus plane, (b) minus plane
are needed for 6 design variables), thus requiring about 3 times as much as that for sensitivity
analysis by the direct differentiation approach in this particular example. It should be noted that
the ratio of the total elapsed CPU time by the forward finite difference method to that by the direct
differentiation approach approximately linearly increases with the number of design variables (i.e.
cooling channels). As shown in Figure 7, sensitivity analysis results by the forward finite difference
2852
S. 1. PARK AND T. H. KWON
9
max: 3.49
L
1
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
2:
1:
3.32
3.14
2.97
2.80
2.62
2.45
2.27
2.10
1.92
1.75
1.57
1.40
1.23
1.05
0.88
min: 0.70
max:
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
2:
1:
min:
10.93
10.30
9.66
9.02
8.38
7.74
7.10
6.46
5.83
5.19
4.55
3.91
3.27
2.63
2.00
1.36
0.72
Figure 10. Sensitivity results of part surface temperature by the direct differentiation approach with respect to y-co-ordinate
of cooling channel #I (unit: OC/cm): (a) plus plane, (b) minus plane
method converge to those by the direct differentiation approach as the number of significant digits
increases and the magnitude of difference decreases, which in turn proves the accuracy of both
approaches. Thus, we can conclude that the direct differentiation approach is more accurate and
faster than the forward finite difference method.
2853
SENSITIVITY ANALYSIS FORMULATION
I
I - I
max: 3.68
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
2:
1:
mln:
max:
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
2:
3.50
3.31
3.13
2.94
2.76
2.58
2.30
2.21
2.03
1.84
1.66
1.48
1.29
1.11
0.92
0.74
11.51
10.84
10.17
9.50
8.82
8.15
7.48
6.81
6.14
5.47
4.79
4.12
3.45
2.78
2.11
1: 1.43
min: 0.76
Figure 1 1. Sensitivity results of part surface temperature by the forward finite difference method with respect to y-co-ordinate
of cooling channel #I (unit: OCkm): (a) plus plane, (b) minus plane
For each design variable of cooling channel #I, sensitivities of part surface temperature with
respect to the radius, the y-co-ordinate, the z-co-ordinate both by the direct differentiation approach
and by the forward finite difference method are shown in Figures 8 and 9, Figures 10 and 11, and
Figures 12 and 13, respectively. (Six-digit significant figure and 0.01 per cent difference in design
2854
S. J. PARK AND T. H. KWON
max:
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
v
I
4:
3:
2:
1:
1.26
1.08
0.91
0.73
0.56
0.38
0.21
0.03
-0.15
-0.32
-0.50
-0.67
-0.85
-1.02
-1.20
-1.37
min: -1.55
I
max:
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
2:
1:
5.50
4.92
4.33
3I75
3.17
2.59
2.00
1.42
0.84
0.26
-0.32
-0.91
-1.49
-2.07
-2.65
-3.24
min: -3.82
Figure 12. Sensitivity results of part surface temperature by the direct differentiation approach with respect to z-co-ordinate
of cooling channel # I (unit: OC/cm) (a) plus plane. (b) minus plane
variables are adopted in the finite difference calculations.) These figures reveal that the computed
sensitivity results are in an excellent agreement with each other. The maximum relative error found
in this example was of about 6-59 per cent. (The relative error on each element is defined as a ratio
of the absolute value of the difference between sensitivity result by the forward finite difference
2855
SENSITIVITY ANALYSIS FORMULATION
1
max:
15:
14:
13:
12:
11:
10:
9:
8:
7:
6:
5:
4:
3:
2:
1:
1.32
1.13
0.95
0.76
0.58
0.40
0.21
0.03
-0.15
-0.34
-0.52
-0.70
-0.89
-1.07
-1.26
-1.44
min: -1.62
max: 5.76
15: 5.15
14: 4.54
13: 3.93
12: 3.32
11: 2.71
10: 2.10
9: 1.49
8: 0.88
7: 0.27
6: -0.34
5: -0.95
4: -1.56
3: -2.17
2: -2.78
1: -3.39
min: -4.00
Figure 13. Sensitivity results of part surface temperature by the forward finite di..-rence method with respect to z-coordinate of cooling channel # 1 (unit: OC/cm) (a) plus plane, (b) minus plane
method and that by the direct differentiation approach to that by the direct differentiationapproach.)
The maximum error may be reduced by increasing the number of significant digits and reducing
the magnitude of difference in the forward finite difference method. From the above results, we
have found the following characteristics of each design variable:
2856
s. J.
PARK AND T.
n. KWON
(a) An increase of radius of cooling channel #1 makes the part surface temperature decrease
mainly due to the increase of surface area of the cooling channel.
(b) An increase of y-co-ordinate of cooling channel #1 make the part surface temperature
increase due to the increase of the distance between the part and cooling channel #l.
(c) As z-co-ordinate of cooling channel #1 increases (i.e. moving in the downward direction
in Figure 4), the part surface temperature increases in the portion of the part where the
distance between the part position and cooling channel #1 increases while it decreases in
the other portion. To no surprise, the line of no change of part surface temperature, that is,
zero sensitivity coincides exactly with the location of cooling channel #l.
It may just be mentioned that we have also compared sensitivity analysis results in the case of
each design variable of cooling channel #2 in the similar manner to the case of cooling channel
#1, and obtained the similar accuracy and results to the case of cooling channel #l.
5 . CONCLUDING REMARKS
An efficient and accurate approach for the design sensitivity analysis of steady-state heat conduction-type problem in complex three-dimensional geometries using boundary element method was
presented. The presented formulae are based on the implicit differentiation of the special boundary
integral equations with respect to the design variables (radii and locations of circular holes) to
obtain the sensitivity equations. Using an illustrative example problem, the numerical results of
the sensitivity analysis developed here is amply discussed in terms of accuracy and computational
efficiency. The accuracy of the present formulation is demonstrated through comparisons between
the present direct differentiation approach and the forward finite difference method. It was concluded
that the former is faster and more accurate than the latter.
These programs will be very useful for designers to obtain optimal configuration of radii and
locations of circular holes (cooling channels in the injection mold) to minimize an objective function in many engineering problems. This optimization problem may be automated and solved
numerically using any standard optimization software, for example, CONMINZowhich is based on
Davidon-Fletcher-Powell method2' treating the constraints via an Augmented Lagrangian Multiplier Method.
Finally, the procedures and techniques proposed in the present paper can be applied to similar
three-dimensional heat conduction problem. The extension of these techniques to engineering problems which could be represented by the Laplace operator should be straightforward, because the
kernels of the integral equations for heat transfer and those engineering problems are fundamentally
similar.
ACKNOWLEDGEMENTS
Portions of the research discussed herein have been supported by grants from R & D Center
in Samsung Electronics (Consumer Electronics Business) and Korea Science and Engineering
Foundation (KOSEF). The authors gratefully acknowledge these financial supports.
APPENDIX I
Co-ordinate transformation rule
The following orthogonal matix represents the co-ordinate transformation rule used for the sensitivity analysis with respect to the location of the cooling channel. The basis {i,j,i } represents
2857
SENSITIVITY ANALYSIS FORMULATION
the local co-ordinate system in each line element of the cooling channel and the basis { f , j , k }
represents the global co-ordinate system as shown in Figure 3. All symbols used in the following
orthogonal matrix are shown in Figure 3:
where
P I = P, - x ,
Pz = P, - y , and P3 = P, - z
In the special case that P is on the axis of the same element as Q, the above co-ordinate
transformation rule can be simplified to permit the following form (in representing v, unprimed
system { v,, vy,vz} for the global co-ordinate system and primed system { vi, v:, v:} for the local
co-ordinate system):
vii+ v:j = (vx - v:lx)i+ ( v y - v i l y ) j + (vz - V;I,)R
v; L = v:(lxi + lJ + 1,R)
where
L = i = 1,i + l y j + lzK
1
-
v, = 1 - 3 = lxvx
+ lYVY + l,v,
APPENDIX I1
Integral formulae for heat transfer analysis
Singular integrals in triangular element.2 The singular integrals over a flat area A in triangular
element can be evaluated analytically. The symbols used in the following closed forms are shown
in Figure 14:
1
+ -1nl
r3 I
tan ((02 -k a3)/2)
Wa3/2)
I
1
12
+ r
tan ((03 -k ai )/2}
tan ( aI /2)
In the above equations, A denotes the area of the triangular element.
I]
2858
S. 1. PARK AND T. H. KWON
Figure 14. Notations for geometric information of triangular element
Closed form of integrals in line element. The integrals over 6 can be evaluated in closed forms.
In the special case that P is on the axis of the same element as Q, the integrals over 0 and 1
can be also evaluated in closed forms. The symbols used in the following definitions and closed
forms are shown in Figure 3.
We first define some fkctions used in evaluating integrals over 8. To avoid conflicts with other
symbols in the paper, we shall use non-standard, calligraphic symbols for the following functions.
X and & denote the complete elliptic integrals of the first and second kind, respectively, defined
as
X ( m )=
de
JiYza
Lnj2
&(m) = J 0x ’ 2 d C Z & d ~
9,
9, 9 and 9’denote simple modifications of X and & and can be evaluated using the derivatives
of X and 8:
where
+
m = B ~ / ( A B~ ~ )
SENSITIVITY ANALYSIS FORMULATION
W , X,CY and 9’ are more general forms of
tl = L t3
=
~ p , t2 =
P,22, 9 and 9,
respectively,
d a 2 +z;
Jq,
+
A’ = S:
s;, B2 = 4aRp
m = B ~ / ( A ~B ~ )
+
In the analysis formulation, the following integral formulae are needed:
a d 0 = 4 a W ( m ; 1)
12 = ~ 2 n - a( - ) 1a d 0 = 4 a X ( m ; s ~ , 2 R p )
an
r
2859
2860
S. J. PARK AND T. H. KWON
I3
=
14
=
2na
I
(--)
ad0
1'" (%
a
a
= 4a%(m;s3, -2avi)
i
(--))ud0=
12aY(m;~1~3,-2~~,-B~v~)+4av~%(m;-1,2)
In the special case that P is on the axis of the same element as Q, the above formulae can
be simplified to obtain closed forms of integrations over 1. In such a case, the following integral
formulae are needed:
JI =
J2 =
/dL 1'"1
1"12" (-)
lLL2"av (--)
a d o d l = 27ta In
r
a 1 a d o d l = 2n
an r
a 1 a dB d l = 27tav;
J3
=
J4
= l L 1 2 ' - ( - (a- ) ) a d 0 id l
av an r
(2+ 2)
(i i)
-
=27ta2vi
APPENDIX 111
Integral formulae for sensitivity analysis
In the sensitivity analysis formulation with respect to the radius of the cooling channel, the
following integral formulae are needed:
all - 1 1
_
- - + 4a%(m; -SI, -2Rp)
aa
u
-=
'2
aa
a
313
13
aa
a
ar4
14
-
12aq(m;s~,4Rps1,4R~)
+ 4a%(m; 1,0)
- = - - 12aY(m;sls3, -2s4, -B2v:)
-=
aa
+ 4a v:%(m;
1, -2)
- - 6 0 a ~ ( m ; s ~ s 3 , 2 s ~ s ~ , 4 -2B'
R pRspv:)
g,
a
+12a9(m;2slv:
+ s3, -2(3sl
- Rp)v:,
-8Rpv:)
In the sensitivity analysis formulation with respect to the location of the cooling channel, the
following integral formulae are needed:
~ t ( 1 1 ) = -Vx(Ii) = 4 a ~ ( m ; - s 1 , 2 a ) 2 ^ - 4 a s ~ ~ (Im
, O;)
~ t ( 1 2 )=
-v,(I~)
=
{ 1 2 a ~ ( m-s:,~s:,B')
;
i
-I- 4 a ~ ( m1,-2)}
;
2^
- 12s2Y(m; s,, 2Rp, 0)
+
vt(I3) = -Vx(13) = - 1 2 a Y ( m ; s ~ s ~ , - 2 a ( s ~ v s3),4a2vi)zI
~
+ 4 8 a 3 v ~ ~ ( m ; 0 , - 1 , 1 ) j -12as2Y(m;s3,-2av~,O) i + 4 a O%(m; 1,O)
QU4) = -VXI4)
= { -6Oab(m;sfs3, -2s:(s3
+ av:),4as7,2aB2v:) + 12a?Y(m;slv: + s3, -2Fgr8av:)}
+ {-240a3v~Z(m;0,s~,3Rp- a , - 2 R p ) + 9 6 a 2 v ~ ~ ( m ; 01,, - l ) } j
2^
SENSITIVITY ANALYSIS FORMULATION
2861
Next, in the sensitivity analysis formulation with respect to the radius of the cooling channel for
the special case that P is on the axis of the same element as Q,the following integral formulae
are needed
aa
a
Finally, in the sensitivity analysis formulation with respect to the location of the cooling channel
for the same special case as above, the following integral formulae are also needed:
V c ( J l )= -V,(JI)
= -2za
Vt(J2)= -V,(J2) = -272a2
($ - $)
REFERENCES
I . K. J. Singh, ‘Design of mold cooling system’, in A. 1. Isayev (ed.), Injection and Compression Molding Fundamentals.
Marcel Dekker, New York, 1987, pp. 567-605.
2. M. Rezayat and T. Burton, ‘A boundary-integral formulation for complex three-dimensional geometries’, Inr. j . numer.
methods eng., 29, 263-273 (1990).
3. M. R. Barone and D. A. Caulk, ‘Optimal thermal design of injection molds for filled thermosets’, Polym. Eng. Sci.,
25, 608617 (1985).
4. T. Matsumoto and M.Tanaka, ‘Optimum design of cooling lines in injection moulds by using boundary element design
sensitivity analysis’, Finite Elements Anal. and Des., 14, 177-185 (1993).
5. T. Matsumoto, M. Tanaka and T. Takahashi, ‘Determination of cooling tempemture history in injection molding’, in
H. D. Bui, M. Tanaka, M. Bonnet, H. Maigre, and M. Reynier (eds), Inuerse Problems in Engineering Mechanics,
Bakema, Rotterdam, 1994, pp. 295-300.
6. S. J. Forcucci and T. H. Kwon. ‘A computer aided design system for three-dimensional compression mold heating’,
ASME J. Eng. Indusrry, 111, 361-368 (1989).
2862
S . J. PARK AND T. H. KWON
7. T. Matsumoto, M. Tanaka and M. Miyagawa, ‘Boundary element system for mold cooling/heating design’, in
C. A. Brebbia and J. J. Rencis (eds), Boundary Element XV-Vol. 2: Stress Analysis, Computational Mechanics
Publications, Boston 1993, pp. 461-475.
8. T. Matsumoto, M. Tanaka and N. Ishii, ‘Boundary element system for design sensitivity analysis of compression mold
heating’, in B. M. Kwak and M. Tanaka (eds), Computational Engineering, Elsevier Science Publishers, Amsterdam
1993, pp. 465-470.
9. M. Barone and R. J. Yang, ‘Boundary integral equations for recovery of design sensitivities in shape optimization’,
AIAA J., 26, 589-594 (1988).
10. S. Saigal, J. T. Borggaard and J. H. Kane, ‘Boundary element implicit differentiation equations for design sensitivities
of axisymmetric structures’, Inr. J. Solids Struct., 25, 527-538 (1989).
1 1. S. Saigal and A. Chandra, ‘Shape sensitivities and optimal configurations for heat diffusion problems: a BEM approach’,
J. Heat Transfer, 113, 287-295 (1991).
12. K. G. Prasad and J. H. Kane, ‘Three-dimensional boundary element thermal shape sensitivity analysis’, Inr. J. Heat
Mass Transfer, 35, 1427-1 439 (1 992).
13. K. Himasekhar, J. Lottey and K. K. Wang, ‘CAD of mold cooling in injection molding using a three-dimensional
numerical simulation’, ASME J. Eng. Industry, 144, 213-221 (1992).
14. T. H. Kwon, ‘Mold cooling system design using boundary element method’, A S M E J. Eng. Industry, 110, 384-394
(1988).
15. C. A. Brebbia, J. C. F. Telles and L. C. Wrobel, Boundary Element Techniques-Theory and Application in
Engineering, Springer, Berlin, 1984.
16. M. R. Barone and D. A. Caulk, ‘Special boundary integral equations for approximate solution of Laplace’s equation
in two-dimensional regions with circular holes’, Q. J. Mech. Appl. Marh., 34, 265-286 (1981).
17. J. C. F. Telles, L. C. Wrobel, W. J. Mansur and J. P. S. Azevedo, ‘Boundary elements for cathodic protection problems’,
in C. A. Brebbia and G. Maier (eds), Boundary Element VII-Proc. 7th Int. Con$, Springer, Berlin, 1985.
18. J. C. F. Telles, W. J. Mansur, L. C. M‘robel and L. C. Marinho, ‘Numerical simulation of a cathodically protected
semi-submersible platform using the PROCAT system’, Corrosion J., 46, 5 13-5 18 ( 1990).
19. J. C. F. Telles, W. J. Mansur and L. C. Wrobel, ‘On boundary elements for external potential problems’, Mech. Res.
Commun., 11, 373-377 (1984).
20. P. C. Haarhoff and J. D. Buys, ‘A new method for the optimization of a nonlinear function subject to nonlinear
constraint’, Comput. J . , 13, 178-184 (1970).
21. R. Fletcher and M. J. D. Powell, ‘A rapidly convergent descent method for minimization’, Comput. J., 6, 163-168
(1963).
Документ
Категория
Без категории
Просмотров
8
Размер файла
1 217 Кб
Теги
906
1/--страниц
Пожаловаться на содержимое документа