=Paper= {{Paper |id=Vol-2426/paper6 |storemode=property |title=Numerical Simulation of Heat Transfer in Semiconductor Heterostructures |pdfUrl=https://ceur-ws.org/Vol-2426/paper6.pdf |volume=Vol-2426 |authors=Karine K. Abgaryan,Ilya S. Kolbin }} ==Numerical Simulation of Heat Transfer in Semiconductor Heterostructures== https://ceur-ws.org/Vol-2426/paper6.pdf
          Numerical Simulation of Heat Transfer in Semiconductor
                            Heterostructures

                                          Karine K. Abgaryan1,2, Ilya S. Kolbin1
                           1
                          Dorodnitsyn Computing Center FRC CSC RAS, Moscow, Russia
                    2
                     Moscow Aviation Institute (NRU) (MAI), Moscow, Russia, kristal83@mail.ru



                                                            Abstract
                        The paper considers the problems of building-up numerical models of
                        heat transfer in multilayer nanostructure AlAs/GaAs. The calculations
                        relied that the temperature inside each layer is equal through its full
                        width, while the changes proceed in discrete steps at the interfaces. The
                        problem solving was carried out on the basis of hybrid finite-difference
                        meshless method, for time discretization an implicit finite-difference
                        scheme was chosen. In these settings dimensional approximation was
                        performed by applying the radial basis function package. The choice of
                        weight coefficients for meshless solution was carried on through the
                        system of linear equations, and the centers and the “widths” of the bases
                        were fixed. The obtained results were compared for different bases
                        (Gaussian function, multiquadric and inverse multiquadric) together with
                        and without normalization approximating model. It was shown that
                        normalization gives higher accuracy results in some cases. For the quality
                        assessment of the obtained results, they were compared with those,
                        calculated according classical method of finite differences. The offered
                        algorithm gave similar results. Investigation the matter of increasing
                        speed method using concurrent computation revealed that application of
                        multiprocessor system can allow to reach significant growth of efficiency
                        of the developed algorithm.

1        Introduction
   Nowadays there has been a significant amount of interest in multilayer heterostructures with nanoscale layer
thickness. In a point of fact, these structures have become a frequent practice in micro and nanoscale semiconductor
electronics production. The reason for this is high values of thermal conduction, electrical conductivity etc. as
compared to traditional materials. Thus, simulation of procedures of physical processes in corresponding structures is
an urgent science and technology task.
   One of the promising approaches for heat transfer simulation is the usage of meshless methods. Application of
corresponding algorithms allows to penetrate some difficulties commonly found in classical methods of finite
differences, finite elements, and etc. such as exponential increase of complexity considering the problems of high
dimensionality, solving problems with compound range and etc. [3]. The widely accepted meshless approach is
approximation of simulating processes by means of radial basic functions (RBF) [19]. Several equivalent forms of
recording approximating RBF. This paper presents the following:
                                           u ( x ) =  i =1  i  i ( x ) ,
                                                          M
                                                                                                                           (1)

                                            i ( x ) =  ( ri ( x ),  i ) ,                                               (2)

                                              ri ( x ) = x − x ic ,                                                        (3)
whereω – weight coefficients, ε – “widths” of bases, xc – centers coordinates. For example, Hardy multiquadrics are
used as radial functions φ [7, 8, 17]:


    Copyright © 2019 for the individual papers by the papers' authors. Copyright © 2019 for the volume as a collection by its
editors. This volume and its papers are published under the Creative Commons License Attribution 4.0 International (CC BY 4.0).
   In: Sergey I. Smagin, Alexander A. Zatsarinnyy (eds.): V International Conference Information Technologies and High-
Performance Computing (ITHPC-2019), Khabarovsk, Russia, 16-19 Sep, 2019, published at http://ceur-ws.org

                                                                  37
                                  Computer Design of New Materials
    ___________________________________________________________________________________________


                                                  ( x) = 1 +  2 r 2 .                                                (4)
    The algorithm, suggested by Kanza [9, 11, 14, 16], involves selection of weight coefficients of approximating
model at fixed bases location in the given widths. In this setting the choice of the weights is done through solving
simultaneous linear algebraic equations. Considering that multiquadrics are non-local approximators, the matrices of
weight systems comes nonsparse. In case the centers and widths were chosen unsuccessfully, they become ill-
conditioned. Nevertheless, this paper demonstrates [14] that problem solving can be accomplished with a small set of
bases due to the high accuracy of approximating functions. In these conditions the algorithm gave a lower magnitude
of error and performed quicker than the finite differences method.
    A more general meshless approach to problem solving of heat transfer is the method recommended by A. N.
Vasiliev and D. N. Tarkhov. It involves complete adjustment of all parameters of approximating solution. For this
purpose, methods of unconditional nonlinear optimization are used [3, 18]. Discrete quadratic error functional is
formed to achieve this. Then some set of check points is chosen in space, and the approximation model is plugged
into initial equations and boundary conditions. As a rule, the error functional is smooth and in that regard it is rational
to use first-order methods for accelerated search of local optimum. In a point of fact, papers [1, 2] reveal effective
application of conjugate gradients non-linear method CG_DESCENT [4, 15]. An important advantage of this method
is its guaranteed convergence. Alternatively, there may be considered quasi-Newton methods [5, 6], which are also
very effective, but much costlier as they need more memory resources in order to keep approximated Hessian matrix
(or its part, for methods with limited memory, such as L-BFGS [6]). It is necessary to understand that despite high
efficiency of the optimization methods mentioned above, they are local, i.e. resulting accuracy of the solution and the
time spent on its search depend heavily on initial data. Generally, the search of global optimum doesn’t bring any
practical utility due to high computing costs. In a number of cases it may not even exist, particularly while
considering ill-posed problems. Various approaches to global optimization, valid for parameter choice for meshless
simulations, are stated in the monograph [3].
    This paper presents numerical simulation of heat transfer in heterogenic multilayer structure on the macroscopic
level. The paper [10] introduces an approximative model for calculating temperature balance in the layers. A finite-
difference approach was used for problem solving. The computing has been exercised on implicit scheme due to
physical limitation of the observed materials. In this paper we propose a hybrid method of problem solving with
meshless approximation in space and finite-difference in time.

2        Problem Setting
   It is required to build numerical model of heat transfer in multilayer nanostructure AlAs/GaAs. It is assumed that
the temperature is equal inside each layer full width, and the changes proceed in discrete steps at the interfaces. The
model of the heat transfer is as follows [10]:
                                    dTi
                       C i  i hi         =  i −1 (Ti −1 − Ti ) −  i (Ti − Ti +1 ), i = 1K N ,                       (5)
                                    dt
whereC –heatcapacity, ρ –density, h –layerwidth, σ –heatconductivity, N – numberoflayers. The outer boundaries of
the structure are isolated. That was achieved by zeroing of extreme values of heat conductivity:
                                                   0 =  N = 0 , then
                                                        dT1
                                           C 1  1 h1         = − i (T1 − T 2 ) ,                                     (6)
                                                        dt
                                                    dT N
                                     C N  N hN               =  N −1 (T N −1 − T N ) .                               (7)
                                                        dt
Along with that initial temperature distribution by layers was defined:
                                                   Ti ( x , 0 ) = Ti 0 ( x ) .                                         (8)


3        Problem-Solving Method
   A hybrid finite-difference meshless approach in used for solving the problem. The main principle of this method
involves finite-difference splitting by times and continuous approximation [1, 2] at each time layer. Thus, the
algorithm performance will result in a set of continuous models, which give approximate solution at each moment of
simulation. In the general case, the architecture of each approximating model may be individual, but for the
simplicity, homogeneous networks with equal number of elements and identical basis function are used in this paper.
   The implicit finite-difference scheme serves for discretization of the first temperature derivative by time:


                                                                       38
                                  Computer Design of New Materials
    ___________________________________________________________________________________________


                                       dT i k +1 ( x )        T k +1 ( x ) − T k ( x )
                                                          =                              .                         (9)
                                             dt                                 
   For dimensional approximation the radial basis functions network was used. Significantly, it is rather difficult to
gain high accuracy throughout x-coordinate using RBF as infinitely smooth approximation of step functions gives a
bigger error at a point of opens. In order to avoid this limitation, the bases were fixed in the centers of layers.
Supposedly, the value in the center defines the temperature of the whole layer.
   As radial basis functions, beside of multiquadrics, the following were used for comparison:
                                 ( x ) = exp( −  2 r 2 ) , Gaussian function,                                   (10)

                                                      1
                                 ( x) =                         ,inverse multiquadric.                           (11)
                                                1 +  2r 2
   Practically, normalization [13, 20] can be applied while using RBF-aproximation. In a number of cases it allows to
obtain a more accurate result [1, 2, 12, 13], therefore approximative model takes the form:

                                                         ( x) .
                                                                 M
                                                                 i =1       i       i
                                             u ( x) =                                                             (12)
                                                        ( x)
                                                                    M
                                                                     i =1       i

   The choice of parameters of the approximative model is carried out according to Kanza’s method, i.e. the centers
and the widths remain fixed, and the algorithm must calculate only weight values. For this purpose, let us plug (9) in
(5). After elementary transformations it becomes possible to make a tridiagonal matrix of the system of linear
algebraic equations relative to the temperature of the layers T.
                                                  A i , i = 1 + q i −1 + q i ,
                                                  A i , i −1 = − q i −1 ,
                                                  Ai ,i + `1 = − q i ,                                            (13)

                                                               i
                                                  qi =                  .
                                                          Ci ih
   In case of the classical method of finite differences the temperature calculations on the successive time layer will
be conducted according the formula given below:
                                                     T k +1 = A −1T k .                                           (14)
  This obtained solution will be used as a standard for the accuracy assessment of meshless approximations.
Numerical calculation of weight coefficients for approximating model needs the following auxiliary matrices:
                                Bi , j = Ai ,i −1C i −1, j + Ai ,i C i , j + Ai ,i +1C i +1, j ,
                                                                                                                  (15)
                                Ci, j =  i, j .
    In this context, the weight coefficients on the successive time step can be found in the following way:
                                                   W k + 1 = B − 1T k ,
                                                                                                                  (16)
                                                    T k +1 = CW k +1 .


4        Test Problem
   For testing purposes of the suggested problem solving method let us consider a heterogenic structure with N =
1000 alternating layers AlAs/GaAs, equal width for all layers h = 50 nm, thermophysical properties of materials
tested:
                                       CAlAs= 424 J/(kg·К), CGaAs= 327 J/(kg·К),
                                         ρAlAs= 3730 kg/m3, ρGaAs= 5320 kg/m3 ,
heat conductivity on interfaces is preset as:
                              σAlAs/GaAs= 9.6·109 Wt/(m2·К), σGaAs/AlAs= 1.92 109 Wt/(m2·К).




                                                                                39
                                  Computer Design of New Materials
    ___________________________________________________________________________________________

    The problem solving involved normalized and unnormalized radial-basis approximators. The bases centers were
                                                                                                                      1
located in the centers of the nanostructure layers: x i = (i − 0 .5) h . The bases widths were fixed:  i =
                                                           c
                                                                                                                          .
                                                                                                                      h
Timestepwasacceptedasτ = 10-10s.
   Accordingly, all necessary data for setting up the system lineal algebraic equation needed for the search of weight
coefficients for each time layer are given. Here are the results of the computing.
   The quality assessment was conducted by means of comparison with the results, obtained while using classical
method of finite differences. The Figure 1, which is given below, reveals the solutions found at different moments of
simulation. The numeric evaluation of accuracy needed computing of average squared deviation of the obtained
solutions in the centers of layers at the considered time period. Absolute values of errors for various meshless
approximations are presented in Table 1.

                       Table 1: Absolute values of errors for various meshless approximations

                                        Approximator                                       Error, К
                                      Multiquadric (Mq)                                     1.5·10-4
                                    Gaussian function (G)                                   6.2·10-6
                                  Inverse multiquadric (IMq)                               3.22·10-6
                          Multiquadric, with normalization (Mq-N)                          4.53·10-5
                         Gaussian function, with normalization (G-N)                        3.3·10-6
                      Inverse multiquadric, with normalization (IMq-N)                     3.52·10-6




    Figure 1: Model problem solving: FDM – finite difference method; meshless methods with bases: G – Gaussian
         function, Мq –multiquadric, IMq – Inverse multiquadric; normalized approximations with suffix “-N

   It is obvious that the results shown in the diagram, which were obtained using finite differences method stiffly
accurate coincide with meshless approximation. This fact demonstrates the efficiency of the method under the
discussion. Table 1 shows that normalization usage allows, in some cases, to increase the accuracy of the results
significantly. In fact, multiquadric has lowered the error three times, while Gaussian function has decreased its error 2
times.

5         Parallel Implementation
   Nowadays, one of the common ways of increasing the rate of calculation is considered to be a concurrent
computation on a great number of processors. Typical usage includes multicore central processing units and graphic
processing units. The algorithms, considered in this paper, utilize effectively parallelized procedures. That is the
natural reason for estimating possible growth of productivity stimulated by application of concurrent computation.

                                                           40
                                  Computer Design of New Materials
    ___________________________________________________________________________________________

Within the presented paper the tests were conducted in order to estimate the growth of calculation rate of the applied
algorithms at parallel performance on various systems:
     − iMac 2012, Intel Core i5-3470S (2 cores), OSX 10.2, Clang 9;
     − MSI GT70, Intel Core i7-4800MQ (4 cores), Lubuntu 17.04, GCC 5.4;
     − IBM PowerNV, Power8 (8 cores), Ubuntu 16.04, GCC 5.4;
     − IBM PowerNV, Power8 (10 cores), RHEL 4.8.5, XLC_R 13.1.
   In order to compare productivity, synthetic tests were held, involving computation of matrix products with the help
of the dgmm procedure from the standard function set for linear algebra BLAS, as execution of this function
consumes a larger amount of operating hours of the program. The calculations were carried out for square matrices
sized 10000х10000. Only direct calculating time was measured. For instance, time used for memory allocation was
not taken into consideration. Relative results are given in Table 2. The first configuration (iMac, running single
thread), was accepted as a standard.

                                 Table 2: Parallel implementation efficiency growth

                                 System (number of cores)                           Efficiency Growth
                                Intel Core i5-3470S (1 core)                                 1
                                Intel Core i5-3470S (2 cores)                               1.98
                               Intel Core i7-4800MQ (4 cores)                               4.34
                                      Power8 (8 cores)                                      8.94
                                     Power8 (10 cores)                                     14.85

   The review of the quoted results reveals that employing multicore systems allows to increase the algorithm
efficiency considerably. It should be emphasized that only concurrency on CPU was under the research in this paper.
The implementation of the methods mentioned above on the graphical processors will show more significant
efficiency growth.


6        Conclusion
   The paper considered the hybrid finite difference method for heat transfer simulation in heterogenic
nanostructures. The performance of the suggested algorithm was demonstrated on the example of solving a model
problem. The possibility was investigated of computing rate increase by means of parallelization of calculations. It
was demonstrated that usage of multiprocessing system can lead to a significant efficiency growth.
   The calculations were performed by Hybrid high-performance computing cluster of FRC CS RAS) [21, 22] and
Shared Facility Center “Data Center of FEB RAS” (Khabarovsk) [23].

References
1. Vasilyev, A.N., Kolbin, I.S., Reviznikov, D.L.: Meshfree Computational Algorithms Based on Normalized Radial
         Basis Functions. Advances in Neural Networks – ISNN 2016. Lecture Notes in Computer Science, vol.
         9719, pp. 583-591. (2016). DOI: 10.1007/978-3-319-40663-3 67.

2. Kolbin, I.S., Reviznikov, D.L: .Resheniye zadach matematicheskoy fiziki s ispolzovanieym normalizovannikh
        radialno-basisnikh setey. // «Neirokompyuteri»: razrabotka, primineniye. № 2. pp. 12-19. (2012). (In
        Russian)

3. Vasilyev, A.N., Tarkhov, D.A.: Neirosetevoye modelirovaniye. Principi, algorimi, prilozheniya – SPb.:Izd. Politeh.
         (2009). (In Russian)

4. Hager, W.W., Zhang, H.: A new conjugate gradient method with guaranteed descent and an efficient line search.
        SIAM Journal on Optimimization, vol. 16, pp. 170-192. (2005). DOI: 10.1137/030601880.

5. Nocedal, J., Wright, S.J.: Numerical optimization. New York: Springer Verlag, p. 634. (1999).

6. Nocedal, J.: Updating quasi-Newton matrices with limited storage. Mathematics of Computation, vol. 35, pp. 773-
        782. (1980). DOI: 10.1090/S0025-5718-1980-0572855-7.



                                                         41
                                Computer Design of New Materials
  ___________________________________________________________________________________________

7. Hardy, R.L.: Multiquadric equations of topography and other irregular surfaces. Journal of Geophysical Research,
        vol. 76.no. 8, pp. 1905–1915. (1971). DOI: 10.1029/JB076i008p01905.

8. Hardy, R.L.: Theory and applications of the multiquadricbiharmonic method — 20 years of discovery 1968–1988.
        Computers and Mathematics with Applications, vol. 19, no. 8/9, pp. 163–208. (1990). DOI: 10.1016/0898-
        1221(90)90272-L.

9. Kansa, E.J.: Multiquadrics — a scattered data approximation scheme with applications to computational fluid-
        dynamics. 1. Surface approximations and partial derivative estimates. Computers & Mathematics with
        Applications, vol. 19, no. 8/9, pp. 127–45. (1990). DOI: 10.1016/0898-1221(90)90270-T.

10. Vorobyov, D.A, Khvesyuk, V.I.: Мetod rascheta nestacionarnogo nagreva nanostructur. // Nauka i obrazovaniye.
        – pp. 541-550. (2013). DOI: 10.7463/1013.0617255. (In Russian)

11. Chen, W.: New RBF Collocation Methods and Kernel RBF with Applications: Meshfree Methods for Partial
        Differential Equations. Lecture Notes in Computational Science and Engineering, vol. 26, pp. 75–86. (2002).
        DOI: 10.1007/978-3-642-56103-0_6.

12. Benaim, M.: On the functional approximation with normalized Gaussian units. Neural Computation, vol. 6, pp.
        314-333. (1994). DOI: 10.1162/neco.1994.6.2.319.

13. Bugmann, G.: Normalized Gaussian radial basis function networks. Neurocomputing, vol. 20, no. 1/3, pp. 97-110.
       (1998). DOI: 10.1016/S0925-2312(98)00027-7.

14. Kansa, E.J.: Multiquadrics – a scattered data approximation scheme with applications to computational fluid
        dynamics-II. Solutions to hyperbolic, parabolic, and elliptic partial differential equations.Computers &
        Mathematics with Applications, vol. 9, no. 8/9, pp. 147-161. (1990). DOI: 10.1016/0898-1221(90)90271-K.

15. Hager, W.W., Zhang, H.: Algorithm 851: CG_DESCENT, A conjugate gradient method with guaranteed descent.
        ACM Transactions on Mathematical Software, no. 32, pp. 113-137. (2006). DOI:
        10.1145/1132973.1132979.

16. Kansa, E.J. Sarra, S.A.: Multiquadric Radial Basis Function Approximation Methods for the Numerical Solution
        of Partial Differential Equations. (2009). URL: http://www.techscience.com/acm/2009/v2.html. (accessed
        20.12.2016)

17. Franke, R.: Scattered data interpolation: Tests of some methods. Mathematics of Computation, vol38., no. 157,
        pp. 181-192. DOI: 10.1090/S0025-5718-1982-0637296-4. (1982)

18. Vasilyev, A.N., Tarkhov, D.A.: Postroyeniyepriblizhennikhneirosetvikhmodeleyporaznorodnimdannim. // Mat.
        Model. – 19:12. – pp. 43–51. (2007). (In Russian)

19. Haykin, S.S.: Neural Networks: A Comprehensive Foundation, 2nd Edition. Prentice-Hall, p. 936. (1999).

20. Bugmann, G., Koay, K.L., Barlow, N., Phillips M., Rodney, D.: Stable Encoding of Robot Trajectories using
       Normalized Radial Basis Functions: Application to an Autonomous Wheelchair. // Proc. 29th Intl. Symp.
       Robotics ISR'98, Birmingham, UK. – pp. 232-235. DOI: 10.1504/IJVAS.2006.01221. (1998)

21. Federal Research Center ‘Computer Science and Control’ of the Russian Academy of Sciences. Available at:
        http://hhpcc.frccsc.ru (accessed 09/12/2018).

22. Zatsarinny, A.A., Gorshenin, A.K., Kondrashev, V.A., Volovich, K.I., Denisov, S.A.: Toward high performance
         solutions as services of research digital platform. // Procedia Computer Science. Volume 150. p. 622-627.
         (2019)

23. Sorokin, A.A., Makogonov, S.V., Korolev, S.P.: The Information Infrastructure for Collective Scientific Work in
        the Far East of Russia // Scientific and Technical Information Processing. Vol. 44. Number 4. P. 302-304.
        (2017)




                                                        42