Identification of Thermal Conductivity of Modern Materials Using the Finite Element Method and Nelder-Mead's Optimization Algorithm

Due to the fact that modern materials are widely used in aerospace industry, automotive industry, medicine and many others it is very important to find out a way to meet all its needful parameters. There is a lot of composites, which arise from a combination of at least two different components on the macroscopic level and whose parameters are unknown. Therefore, it is difficult to get to know their possibilities and functionality. A kind of composite materials are Functionally Graded Materials (FGM). They are characterized by the fact that its composition and structure gradually change over the volume, which follows from changes in properties of material (Miyamoto et al., 1999). There are many works on this topic.


Heat transfer
The heat transfer can be defined as a movement of energy which is caused of temperature difference.It can be provided by the three mechanisms.First of them is a conduction, which can be described as diffusion, which is held in a stationary medium and occurs because of temperature gradient.The mentioned medium can be in form of solid or fluid.Next is convective, which appears as a result of fluid motion.The last one is radiation, which follows from electromagnetic waves between two surfaces, on which different temperatures are.Additionally those surfaces must comply with a condition that the first surface is visible to an infinitesimally small observer on the second surface.
The heat transfer by conduction can be defined by the heat equation where: T -is the temperature, ρ -is the density, C p -is the heat capacity at constant pressure, k -is the thermal conductivity and Q -is a heat source or heat sink.Taking into consideration a steady-state model the temperature does not change with time.
The thermal conductivity describes a relationship between the heat flux vector q and the gradient of temperature T (Bejan&Kraus, 2003), so it takes a form of The heat flux, mentioned above, is a kind of boundary condition, which can be described as where: q 0 -is the input heat flux, h• T −T -is used for convective heat transfer, where h is the heat transfer coefficient and T is the ambient bulk temperature, C • T −T -is used for radiation heat transfer, where T is the temperature of surrounding radiation environment and C is a product of surface emissivity and the Stefan-Boltzmann constant.

The Nelder-Mead algorithm
The Nelder-Mead algorithm is a method that does not require to determine a derivative of an objective function.This function is determined in few specific points, different in each iteration.The first simplex algorithm was defined by Spendley in 1962.In 1965 Nelder and Mead improved it and turned the simplex search into an optimization algorithm by adding options like: reflection, expansion, contraction and shrinking.Thanks those operations, which speed up the process of optimization, the simplex deforms in way that was suggested by Nelder and Mead to adapt better to the objective functions (Nelder & Mead, 1965).
An n-dimensional simplex with n+1 vertices p0,p1,p2,…,pn is the smallest convex set which contains these points, where set {pj -p0 : 1 ≤ j ≤ n} must consist of linearly independent vector.In the two-dimensional space the simplex can be created from any triangle and in three-dimensional space from any tetrahedra.
In this method selected initial simplex is modified by means of elementary geometric operations called: reflection, expansion, contraction and shrinking.As a result of each of them the vertex, where value of the objective function takes the highest value (the "worst" vertex), is replaced by another -"better".In this way the simplex is coming more and more to local minimum of examined function.
Finding the minimum of the objective function must be preceded by an analysis, where as a result the vertices, in which the objective function takes the smallest and the highest value, are marked in the following way (see Fig. 1a):  pmin -the vertex where the objective function takes the smallest value:  pmax -the vertex where the objective function takes the highest value:  ̅ -centroid of the points (the vertex pmax is excluded) see: After determining points pmin, p , pmax, a procedure of minimization of the objective function can begin.In each iteration the following stages can be specified: reflection, expansion, contraction and shrinking, which are described below (Weise, 2009).
Reflection -is based on determining a point which is symmetrical image of point pmax relative to point p .New point is marked as podb (see Fig. 1b) and its coordinates are designed by formula: Value of reflection coefficient is in the range of , , but usually it is assumed that = .
After reflection stage, depending on value of the objective function in reflection point f podb , we consider few excluded cases (8), ( 9) and ( 10), which determine further investigation in given iterations: If in calculated point podb the objective function takes value (8) then the reflection is accepted.The new simplex is designed by replacing the vertex pmax with podb.Next indexes min, max and location of point p are updated and if a stop condition, which is described later, is not fulfilled a new iteration begins with new reflection.
Expansion -Assume that reflection inequality (9) was fulfilled, which means that a vertex which was found in reflection stage is better point than pmin (it is closer to minimum of objective function f).
It suggests that next steps of finding the minimum should follow in this direction.Because of this the reflection is not accepted and the calculations are carried out by the expansion (see Fig. 2).A new point is calculated and marked as pe: where >1 is an expansion coefficient (usually = ).Next, a value of the objective function in new point is calculated f pe , and:  if f podb < f pmin then the expansion is successful, new simplex is designed by replacing pmax with pe (new simplex is designed by vertices pe, p2, pmin -Fig.2a); then indexes min and max and a location of point p are updated and after checking the stop condition next iteration begins;  else when f pe ≥ f podb , pmax is replaced by podb (new simplex is designed by vertices podb, p2, pmin -Fig.2b) and it follows as previous (indexes are updated, stop condition is checked and next iteration begins) Reflection cannot be accepted also in case when f odb ≥ f pmax , see (10).In this situation occurs contraction of a simplex, whose new vertex is counted according to the formula: where a coefficient of contraction β takes a value β , , usually β=0.5 (see Fig. 3a).If point pz leads to improvement, which means f pz < f pmax , then point pmax is replaced by point pz and a new simplex is created (designated by pz, p , pmin).Next indexes are updated, stop condition is checked and next iteration begins.
Shrinking.This stage takes place when after contraction an inequality ( 13) is fulfilled: In this situation point pmin remains unchanged, and the whole simplex is shrinking according formula ( 14): where , is a shrinking coefficient and usually = , (see Fig. 3b).A simplex which is build of a new obtained points p0, …, pn is used in next iteration (if the stop condition is not fulfilled).In this papers two stop conditions were used.The first when an absolute value of difference between f pmin and f pmax is smaller than accuracy solution abs f pmin -f pmax < (15) and the second when a number of iterations is bigger than maximum number of iterations

Reconstruction of thermal parameters in 1-D domain
As a first a reconstruction of thermal parameters was carried out.This simulation was made for heat transfer in 1-D space, in a domain which length was 1 m.The boundary condition was the heat flux at both ends of the domain.Basing on the temperature distribution in area thermal parameters of the issue were designated.Those parameters were: a thermal www.intechopen.com Finite Element Analysis -From Biomedical Applications to Industrial Developments 294 conductivity (k), a transversal convective heat transfer coefficient (h1 and h2) and an external temperature around both ends of the domain (Tinf1 and Tinf2).The reconstruction of mentioned factors was performed by the optimization.

Stage I
In the first stage desired temperature distribution was defined as and minimized integral have a form of The start simplex for particular parameters is presented in Table 1.

Stage II
In the second stage of calculation some restrictions were imposed to the optimal parameters, as it is presented below: Thereby, we did not have to minimize the roots of objective function.
The start simplex for particular parameters is presented in Table 4. Values which were determined are within established limits.

Stage III
Next stage of the research was optimization of the thermal parameters of material in which coefficient of thermal conductivity was dependent on spatial variable x, like in the Functionally Graded Materials.In those composites temperature distribution with given www.intechopen.com Finite Element Analysis -From Biomedical Applications to Industrial Developments 296 boundary conditions are usually nonlinear.In these paper it is assumed that parameter k(x) has polynomial form In this stage the heat transfer equation takes a form of where k(x) -thermal conductivity depends on spatial variable x.
The following boundary conditions (different temperature on ends) was assumed for calculations Now the vector of parameters have a form of: p=[p 0 , p 1 , p 2 ].
There were two tasks solved for one function of T x , one integral and for two different start simplexes.
The desired function T x has a form of: Minimized integral was defined as:

Task 1
The calculation began with the start simplex described in In Fig. 4 disparity between the expected and the obtained temperature distribution is presented.Distribution of the coefficient of thermal conductivity, for k x, p is presented in Fig. 5.In Fig. 6 disparity between the expected and the obtained temperature distribution is presented.Distribution of the coefficient of thermal conductivity, for k x, p is presented in Fig. 7.  Despite the fact that the thermal conductivity coefficient seems to look identical in task 1 and 2, there is some difference.In Fig. 8 the disparity in distribution of mentioned thermal conductivity is shown.Concluding, although in task 2 the start simplex was wider (bigger) than in task 1, the solution was found with the same accuracy.Finite Element Analysis -From Biomedical Applications to Industrial Developments 300 Summarizing, as it was shown above in this section, it is possible to provide a reconstruction of parameters using hybrid method (FEM with Neleder-Mead).Carrying out the simulation for 1-D domain with length 1 m and defined boundary condition the parameters (the thermal conductivity, the transversal convective heat transfer coefficient and the external temperature around both ends of the domain) can be designated.Values of those parameters can be calculated within some restrictions, which can be specified for material which is examined.It is also possible to designate the value of thermal conductivity parameter of FGM which has a polynomial form.

Reconstruction of thermal parameters in 2-D space
In this subsection calculations were made to designate the distribution of the thermal conductivity in 2-D domain.Cylinder with radius r=1m and height z=1m was analysed in 2D axial symmetry model.An axis of symmetry is designated as r=0.In this case it was also assumed that the distribution of the thermal conductivity has a form of polynomial as: A boundary condition, such that temperature at the top and at the bottom of the cylinder was equal in order that T 01 =400 K and T 02 =300 K.For axis r=0 axial symmetry was assumed and on the circumference of the cylinder was assumed a thermal insulation.In these calculations some restrictions have been imposed.An integral I, which contains sum of three integrals was minimized, as shown below where: I is minimized integral, I1 -is an absolute value of difference between expected and obtained temperature distribution I2 -determined part of domain where a relationship such that k(z)>kmin, where kmin is a minimum value, is satisfied, I3 -determined part of domain where a relationship such that k(z)<kmax,, where kmax is a maximum value, is satisfied.
There were two tasks computed, each for one function of expected temperature.Each task was calculated for two variants of values for kmin and kmax -each of them for three start simplexes, as it is presented in subsection below.

Task 1 -First function of temperature distribution
In this task the expected temperature distribution took a form of T 1 z =300.481+171.Three different start simplexes were assumed and collected in Table 10.As it was mentioned above there were two variants of calculations.Results and assumptions for them are presented in subsections 5.1.1 and 5.1.2.

Variant 1
For calculations below we defined restrictions as follows: kmin=20 and kmax=120.Which means that we were looking for k(z) distribution in this range: 20<k(z)<120 For the start simplex A, distribution of the thermal conductivity for minimized value k(pmin) is shown in Fig. 9. Disparity between the expected and the obtained temperature distribution was also examined (see Fig. 10) and it varies between -0.5 and 0.47.In all cases the results were achieved with solution accuracy of 1e-5.Temperature distribution was similar for all simplexes.

Variant 2
For calculations below we defined restrictions as follows: kmin=10 and kmax=320.Which means that we were looking for k(z) distribution in this range: 10<k(z)<320 For the start simplex A, distribution of the thermal conductivity for minimized value k(pmin) is shown in Fig. 15.Disparity between the expected and the obtained temperature distribution was also examined (see Fig. 16) and it varies between -3.8 and 2.8.For the start simplex B, distribution of the thermal conductivity for minimized value k(pmin) is shown in Fig. 17.Disparity between the expected and the obtained temperature distribution was also examined (see Fig. 18) and it varies between -0.56 and 0.46.In all cases the results were achieved with solution accuracy of 1e-5.Temperature distribution was similar for all simplexes.

Task 2 -Second function of temperature distribution
The form of the second expected temperature distribution looks as follows: T 1 z =298.794+149.For this function also two variants were calculated.Results and assumptions for them are presented in chapters 5.2.1 and 5.2.2.

Variant 1
The restrictions defined for these variants take a form of: kmin=20 and kmax=120.This means that we were looking for k(z) distribution in range like: 20<k(z)<120 For the start simplex A, distribution of the thermal conductivity for minimized value k(pmin) is plotted in Fig. 21.Between the expected and the obtained temperature distribution was some disparity which is from -1.48 and 1.2 (see Fig. 22).

Variant 2
The restrictions defined for these variants take a form of: kmin=10 and kmax=320.This means that we were looking for k(z) distribution in range like: 10<k(z)<320 In Fig. 27 the distribution of the thermal conductivity for minimized value k(pmin) is plotted for the start simplex A. It was some disparity, between the expected and the obtained temperature distribution, and value of it was from -1 to 1.2 (see Fig. 28).In all cases the results were achieved with solution accuracy of 1e-5.For all simplexes the temperature distribution was similar.
Summarizing, in this chapter some new possibilities of Identification of Thermal Conductivity of Modern Materials using the Finite Element Method and Nelder -Mead's Optimization Algorithm were proposed.Simulating heat transfer in 2-D axial symmetry model, made of the Functionally Graded Material, it is possible to designate its thermal conductivity distribution.This parameter can have the polynomial form.It is also possible to calculate its values within some restrictions.All solutions were found regardless of how far the start simplex was from the solution.

Conclusion
Because experimental evaluation of thermal parameters of composites is expensive and time consuming, computational methods have been found to be efficient alternatives for predicting the best parameters of composites.As it was presented in this chapter the Nelder-Mead algorithm connected with the Finite Element Method can be used to optimize many different issues.It has its applicable in problems where it is difficult or impossible to designate the gradient of the objective function.The developed hybrid method can be used for optimization of the heat transfer problems.
In the section 4 of this chapter reconstruction of parameters was provided.Some heat transfer parameters in one-dimensional domain with length 1m, for defined boundary conditions were designated using numerical calculations.The thermal conductivity, the transversal convective heat transfer coefficient and the external temperature around both ends of the domain were calculated within some defined restrictions.
Next, in section 5 possibility of designation of the thermal conductivity was shown.The 2-D axial symmetry model was considered where heat transfer was simulated.The thermal conductivity was in polynomial form.There was also possible to put some restrictions on the searched parameters.
The hybrid method, which was proposed here, can be very helpful in designating any parameters of modern materials like for example Functionally Graded Materials.Proposed method can be also used instead of destructive testing of materials.

Fig. 4 .
Fig. 4. Disparity between the expected and the obtained temperature distribution (task 1)

Fig. 6 .
Fig. 6.Disparity between the expected and the obtained temperature distribution (task 2)

Fig. 8 .
Fig. 8. Disparity between the thermal conductivity coefficient in task 2 and task 1

Fig. 26 .
Fig. 26.Variant 1, start simplex C -Disparity between the expected and the obtained temperature distribution

Fig. 29 .
Fig. 29.Variant 2, start simplex B-Distribution of the thermal conductivity Identification of Thermal Conductivity of Modern MaterialsUsing the Finite Element Method and Nelder -Mead's Optimization Algorithm www.intechopen.com

Table 1
Identification of Thermal Conductivity of Modern MaterialsUsing the Finite Element Method and Nelder -Mead's Optimization Algorithm 295 . Start simplex for the stage I Numbers, placed in Table1, are square roots of searched thermal parameters.This assumption guarantees that values will be positive.Results of the calculations are presented below in Table2 and Table 3. Required accuracy of solution was obtained after 88 steps with F=0.0105 for following set of parameters.www.intechopen.com

Table 4 .
Start simplex for the stage II Numbers placed in Table4are not square roots of searched thermal parameters, but there are values which minimize the objective function.

Table 5 .
Values of minimized parameters pmin

Table 6 .
Start simplex for the task 1 of stage III Results were achieved with solution accuracy of 1e-6.

Table 7 .
Values of minimized parameters pmin www.intechopen.comIdentification of Thermal Conductivity of Modern MaterialsUsing the Finite Element Method and Nelder -Mead's Optimization Algorithm 297

Table 8 .
Start simplex for the stage III (task 2) Results were achieved with solution accuracy of 1e-6.Required accuracy of solution was obtained after 151 steps with F=2.9336 for the following set of parameters.

Table 9 .
Values of minimized parameters pmin

Table 10 .
Task 1 -Different start simplexes taken for calculation

Table 11 .
. Numerical results are presented in table below.Values of minimized parameters pmin for simplexes A, B, C

Table 12 .
. Numerical results are presented in table below.Values of minimized parameters pmin for simplexes A, B, C

Table 13 .
Task 2 -Different start simplexes taken for calculation

Table 14 .
. Numerical results are presented in table below.Values of minimized parameters pmin for simplexes A, B, C

Table 15 .
. Numerical results are presented in table below.Values of minimized parameters pmin for simplexes A, B, C