CC BY
13
Discrete & Continuous Models 2021 29 (1) 5—13
. „ . http://journals.rudn.ru/miph
& Applied Computational Science
UDC 519.6
DOI: 10.22363/2658-4670-2021-29-1-5-13
Numerical simulation of thermal processes occurring in materials under the action of femtosecond laser pulses
Ilkizar V. Amirkhanov, Nil R. Sarker, Ibrohim Sarkhadov
Laboratory of Information Technologies Joint Institute for Nuclear Research 6, Joliot-Curie St., Dubna, Moscow Region, 1^1980, Russian Federation
(received: December 26, 2020; accepted: March 12, 2021)
In this work, a numerical study of the solutions of the parabolic and hyperbolic equations of heat conduction with the same physical parameters is carried out and a comparative analysis of the results obtained is carried out. The mathematical formulation of the problem is discussed. The action of the laser is taken into account through the source function, which was chosen as a double femtosecond laser pulse. In the hyperbolic equation, in contrast to the parabolic one, there is an additional parameter that characterizes the relaxation time of the heat flux. In addition, the source of the hyperbolic equation contains an additional term — the derivative of the power density of the source of the parabolic equation. This means that the temperature of the sample is influenced not only by the power density of the source, but also by the rate of its change. The profiles of the sample temperature at different times and its dynamics at different target depths are shown. The calculations were carried out for different time delays between pulses and for different relaxation parameters.
Key words and phrases: parabolic and hyperbolic heat equations, femtosecond laser pulse, numerical simulation
1. Introduction
The study of the interaction of femtosecond laser pulses with matter is important in connection with many fundamental problems (physics of non-equilibrium processes, generation of shock waves, laser acceleration of ions, modification of the properties of the irradiated material, etc.) [1]-[3].
Currently, there is a growing need for the creation and improvement of physical models capable of describing various processes in matter. Moreover, computer modeling now occupies one of the main places in the study of such problems. There are two approaches to the study and creation of physical models — atomistic and continuous.
Atomistic approaches (molecular dynamics method) allow natural consideration of the atomic structure of the crystal lattice, the effect of impurities,
© Amirkhanov I. V., Sarker N. R., Sarkhadov I., 2021
This work is licensed under a Creative Commons Attribution 4.0 International License http://creativecommons.org/licenses/by/4.0/
the presence of dislocations, the kinetics of phase transitions, etc. The continual approach (solving the equations of continuum mechanics) includes the parabolic and hyperbolic heat equation, the two-temperature model of heat conduction, the two-temperature hydrodynamic model, etc. [2].
The molecular dynamics (MD) method [4] can be used to describe the dynamics of fast processes that arise in a substance under the action of a laser pulse. MD is quite effective for microscopic analysis of the mechanisms of melting and evaporation [5], [6]. The appearance and propagation of pressure waves generated by laser radiation [7], [8], as well as the dynamics of laser ablation [9], are well modeled using the MD.
Each approach has its own problems. When studying transport processes within the framework of a parabolic equation, a problem that arises is the infinitely high speed of thermal perturbation propagation (a consequence of the Fourier law). Generalizing the Fourier law, taking into account the relaxation time of the heat flux, we obtain the hyperbolic equation of heat conduction. The relaxation time is a characteristic of nonequilibrium of the heat conduction process. Under exposure to femtosecond pulses, non-equilibrium heating of the material occurs. Therefore, the study of such processes may turn out to be more adequate using the hyperbolic heat equation.
In this work, we carried out a numerical study of the physical processes arising under the action of femtosecond laser pulses within the framework of the parabolic and hyperbolic equations of heat conduction and carried out a comparative analysis of the results obtained.
2. Setting of the problem
When simulating thermal processes arising in materials under the action of femtosecond laser pulses, we use a hyperbolic model of the heat conduction equation:
(dT d2 t\ a2 t A/ N dA(x,t) Cp(m+T^ W )=Xdx2 + A(x, t) + rr — . (1)
Here c, p, X are the specific heat capacity, density, and heat conductivity of the sample material, respectively. T(x,t) is the sample temperature, A(x,t) is the source function, which determines the heat release power density at the point with the coordinate x at the time moment t, rr is the characteristic time of energy flux relaxation.
The second term in the left-hand side of equation (1) reflects the fact that the thermal process is actually hyperbolic rather than parabolic, and this model of heat conduction is widely used in practice [1], [10]—[12].
The relaxation time rr of the heat flux is related to the velocity of heat
propagation by the formula v = ^X/cprr. If v ^ ^ (i.e., rr ^ 0), then we get an equation of the parabolic type. The term rrdA/dt means that the temperature is affected by not only the power density of its sources, but also by the rate of its change. For metals [12] rr = 10-11 s; for steel v = 1800 m/s, for aluminum v = 2830 m/s, for amorphous bodies like glass and polymers the relaxation time attains 10-7 — 10-5 s; in this case v can exceed the velocity of sound propagation vs in these media.
In general, the heat capacity, thermal conductivity, and material density depend on temperature. In this work, the temperature dependence of the parameters of the sample material is disregarded.
Equation (1) is solved with the following initial and boundary conditions:
T(x,0) — To, T(xmax,t) — To,
dT(x,t)
dt
— 0,
dT(x,t)
t=o
dx
— 0.
x=0
(2)
The source function is chosen in the factorized form
A(x,t) = I0[1 - R(TS)]f1 (x)f2(t), Ts = T(0,t).
Here f1 (x), f2(t) are the spatial and temporal shape of the source, respectively, I0 is the source intensity, R(TS) is the coefficient of reflection of the laser pulse from the material surface.
In the present work, f1(x) u f2 (t) are chosen the same as in Ref. [13]:
A (x) —
exp(-x/L}
Î2(t) —
exp
(t-to) 2*2
21
+ exp
(t-to - Td) ■ 2a?
2
Here Lp is the depth of penetration of laser radiation into the substance, t0 is the time moment when the first pulse of the source takes the maximum value, td is the time shift of the second pulse of the source with respect to the first pulse. The radiation dose is
$ —
f? (t)dt — 21 oo at.
When numerically solving equation (1) with initial and boundary conditions (2), it is convenient to replace the dimensional variables and quantities with their dimensionless counterparts. This is carried out as follows:
r-— --— - - + -h - XAt
1—T0 ; X—Ax; t—At; at — At; t(0 — At' ko — cpAx?;
- — A(x,t)At Tr aj_ ; si\X, Z)
cfT ~dt
d 2 T
dt2
d 2 T
dx2
r<7 + rr — ko f-,-2 + A(x, t) + t,
cpTo
_ dÀ(x,t)
dt
dt
dT(0,t) dx
0; T(Xmax 1.
(3)
(4)
The dimensionless source function and the normalization conditions in this case take the form
A(x,t)=AQ h(x)f2(t),
I0[l - R(Ta)]At j . . /r
An = ——;-—; fi(x) = exp(—ax), a = Ax/Lp,
Lp cpT0
f2(t) =
V
exp
(t-to)2 2a+2
+ exp
(t-to-rd )2 2a+
$ = 2I0 Atat.
3. Discussion of numerical results
Numerical experiments were carried out for aluminum irradiated by the double-pulse laser with the following parameters:
w X = 236 —, Km P = 2700 kg 3 , C ' m3 = 920 kgR,
t = 3 •l0-7 •^max ° ±u m, T0 = 300 K, R(TS) = 0,
$ = 4^l05 Ar, m2 = 5 • l0 -14 s, t 0 = 3 • l0-13 s
Ax = 3 • l0 -8 m, At = l0- -12 s.
The total dose $ = 4 • l05 J/m2 for the specified source corresponds to the intensity I0 ^ 1.5957 • l017 W/m2. Dimensionless constants k0, A0, a, t0, at take the following values:
k0 ^ 0.10556; ^ 8404.34l37; a = l; t0 = 0.3; at =0.05.
Below
f(t) = f2 (t)+Tr ^
describes the time dependence of the source. For rr = 0, we get the source for a parabolic equation.
Equation (3) with the initial and boundary conditions (4) was solved using a finite-difference three-layer explicit scheme.
Figures 1 and 3 show the time dependence of the source function, temperature profiles at different times and the dynamics of the sample temperature at different depths. The times ti, i = l,2,... ,10 are selected in such a way that the first five of them correspond to the action times of the source first pulse, and the rest correspond to the action times of the second pulse. The calculations were carried out until the moment the source was turned off at different times of the delay between the pulses rd.
Figures 2 and 4 show the temperature profiles at long times, when the sources are turned off, i.e., f(t) = 0.
0,0 0,2 0,4 0,6 0,8
'(t' nTd=0.4 ps
t=0.1 ps
25 SO 75 100 125 150
T(X,t) 105K
0 25 50 75 100 125 150 T(X,t) 105K
0,0 0,2 0,4 0,6 0,8
0,0 0,2 0,4 0,6 0,8 1,0
T(*„t) 105 K
25 50 75 100 125 150
0,0 0,2 0,4 0,6 0,8 1,0
Figure 1. Time dependence of function f(t) = f2(t) + Trdf2(t)/dt, temperature profiles at different time moments T(x, tj), j = 1,2,..., 10, t1 = 0.25 ps, t2 = 0.3 ps, t3 = 0.35 ps,i4 = 0.45 ps, t5 = 0.55 ps, t6 = 0.65 ps, t7 = 0.7 ps, t8 = 0.75 ps, tg = 0.85 ps, t10 = 1 ps, and dynamics of sample temperature at different depths (T(xt, t), i = 1,2,3, x1 = 0 nm, x2 = 3 nm, x3 =6 nm), obtained in the framework of the hyperbolic heat conduction equation for different values of the parameter Tr (rr = 0 ps, 0.1 ps, 10 ps) h Td = 0.4 ps
Figure 2. Temperature profiles at different time moments T(x, tj), j = 1, 2,..., 5, tj = j ps and the sample temperature dynamics at different depths (T(xi, t), i = 1,2,3, x1 =0 nm, x2 = 3 nm, x3 =6 nm), obtained in the frameworks of the hyperbolic heat conduction equation at different values of parameter Tr (rr = 0 ps, 0.1 ps, 10 ps) and Td = 0.4 ps
1p0 '(') . xd=0.6 ps Tr=0 ps 0,8 0,6 0,4
0,2
0,0
0,0 0,2 0,4 0,6 0,8 1,0 1,2
f(t) Td=0.6 ps x=0.1 ps
tps
0,0 0,2 0,4 0,6 0,8 1,0 v
150 '(') xd=0.6ps x=10ps 100
25 50 75 100 125 150
5 T(X,t) 10 K
43 2
25 50 75 100 125 150
T(X,t) 105K
0,0 0,2 0,4 0,6 0,8 1,0 1,2
t......................*i
0 25 50 75 100 125 150
0,0 0,2 0,4 0,6 0,8 1,0 1,2
T(X„t) 105K
0,0 0,2 0,4 0,6 0,8 1,0 1,2
T(X„t) 105K
0,0 0,2 0,4 0,6 0,8 1,0 1,2
Figure 3. Time dependence of function f(t) = f2(t) + Trdf2(t)/dt, temperature profiles at different time moments T(x, t^), j = 1,2,..., 10, t1 = 0.25 ps, t2 = 0.3 ps, t3 = 0.35 ps,i4 = 0.45 ps, t5 = 0.65 ps, t6 = 0.75 ps, t7 = 0.85 ps, t8 = 0.9 ps, tg = 0.95 ps, t10 = 1.05 ps, and the dynamics of sample temperature at different depths (T(xi, t), i = 1,2,3, x1 =0 nm, x2 =3 nm, x3 = 6 nm), obtained in the frameworks of hyperbolic heat conduction equation at different values of parameter Tr (rr = 0 ps, 0.1 ps, 10 ps) and Td = 0.6 ps
Figure 4. Temperature profiles at different time moments T(x, tj), j = 1, 2,..., 5, tj = j ps and sample temperature dynamics at different depths (T(xt, t), i = 1,2,3, x1 =0 nm, x2 = 3 nm, x3 =6 nm), obtained in the frameworks of the hyperbolic heat conduction equation at different values of parameter Tr (rr = 0 ps, 0.1 ps, 10 ps) and Td = 6 ps
4. Conclusion
In contrast to the parabolic equation, the hyperbolic one includes an additional parameter that characterizes the heat flux relaxation time. A derivative of the power density of the source of the parabolic equation is additionally present in the source of the hyperbolic equation. This fact means that the sample temperature is affected not only by the source power density, but also by the rate of its variation. Due to this dependence, at some time moments the source takes negative values depending on the relaxation time parameter. Nevertheless, the temperature at the sample surface given by the solution of the hyperbolic equation is higher than that given by the solution of the parabolic equation.
Acknowledgments
The work was carried out under the financial support from th Russian Foundation for Basic Research, grants No. 19-01-00645a and No. 20-51-44001 mong-a.
References
[1] S. L. Sobolev, "Local non-equilibrium transport models," Physics Us-pekhi, vol. 40, no. 10, pp. 1043-1053, 1997, in Russian. DOI: 10.1070/ PU1997v040n10ABEH000292.
[2] S. I. Anisimov and B. S. Luk'yanchuk, "Selected problems of laser ablation theory," Usp. Fiz. Nauk, vol. 172, no. 3, pp. 301-333, 2002. DOI: 10.3367/UFNr.0172.200203b.0301.
[3] V. P. Veiko, M. N. Libensonm, G. G. Chervyakov, and E. B. Yakovlev, Interaction of laser radiation with matter. Power optics [Vzaimodeystviye lazernogo izlucheniya s veshchestvom. Silovaya optika], V. I. Konov, Ed. Moscow: Fizmatlit, 2008, in Russian.
[4] M. P. Allen and D. J. Tildesley, Computer simulation of liquids. Clarendon Press, 1991.
[5] Z. H. Jin, P. Gumbsch, K. Lu, and E. Ma, "Melting mechanisms at the limit of superheating," Physical Review Letters, vol. 87, p. 055 703, 5 Jul. 2001. DOI: 10.1103/PhysRevLett.87.055703.
[6] F. F. Abraham and J. Q. Broughton, "Pulsed melting of silicon (111) and (100) surfaces simulated by molecular dynamics," Phys. Rev. Lett., vol. 56, pp. 734-737, 7 Feb. 1986. DOI: 10.1103/PhysRevLett.56.734.
[7] V. Zhigilei and B. J. Garrison, "Pressure Waves in Microscopic Simulations of Laser Ablation," in Materials Research Society (MRS) Proceedings, vol. 538, Cambridge University Press, 1998, pp. 491-496. DOI: 10.1557/ PROC-538-491.
[8] J. I. Etcheverry and M. Mesaros, "Molecular dynamics simulation of the production of acoustic waves by pulsed laser irradiation," Phys. Rev. B, vol. 60, pp. 9430-9434, 13 Oct. 1999. DOI: 10.1103/PhysRevB.60.9430.
[9] L. V. Zhigilei and B. J. Garrison, "Microscopic mechanisms of laser ablation of organic solids in the thermal and stress confinement irradiation regimes," Journal of Applied Physics, vol. 88, no. 3, pp. 1281-1298, 2000. DOI: 10.1063/1.373816.
[10] A. V. Lykov, Heat and Mass Transfer [Teplomassoobmen], 2nd. Moscow: Energiya, 1978, in Russian.
[11] P. Vernott, "Les paradoxes de la théorie continue de l'équation de la chaleur," Comptes rendus de l'Académie des Sciences, vol. 246, no. 22, pp. 3154-3155, 1958.
[12] E. M. Kartashov and V. A. Kudinov, Analytical methods of the theory of heat conduction and its applications [Analiticheskiye metody teorii teploprovodnosti i yeye prilozheniy]. Moscow: LENAND, 2018.
[13] V. B. Fokin, "Continuous-automaton model and its application for numerical calculation of the effect of single and double femtosecond laser pulses on metals [Kontinual'no-atomaticheskaya model' i yeye prime-neniye dlya chislennogo rascheta vozdeystviya odinochnogo i dvoynogo femtosekundnogo lazernogo impul'sa na metally]," in Russian, Candidate of Sci. in Phys. and Math. (PhD) Thesis, Joint Institute for High Temperatures of the Russian Academy of Sciences, Moscow, 2017.
For citation:
I. V. Amirkhanov, N. R. Sarker, I. Sarkhadov, Numerical simulation of thermal processes occurring in materials under the action of femtosecond laser pulses, Discrete and Continuous Models and Applied Computational Science 29 (1) (2021) 5-13. DOI: 10.22363/2658-4670-2021-29-1-5-13.
Information about the authors:
Amirkhanov, Ilkizar V. — Candidate of Physical and Mathematical Sciences, Head of Sector "Scientific Division of Computational Physics". Laboratory of Information Technologies of the Joint Institute for Nuclear Research (e-mail: camir@jinr.ru, ORCID: https://orcid.org/0000-0003-2621-144X, Scopus Author ID: 6507929197)
Sarker, Nil R. — Candidate of Physical and Mathematical Sciences, Senior Researcher "Scientific Division of Computational Physics". Laboratory of Information Technologies of the Joint Institute for Nuclear Research (e-mail: sarker@jinr.ru, ORCID: https://orcid.org/0000-0003-0690-2534, Scopus Author ID: 14829498100)
Sarkhadov, Ibrohim — Candidate of Physical and Mathematical Sciences, Senior Researcher "Scientific Division of Computational Physics". Laboratory of Information Technologies of the Joint Institute for Nuclear Research (e-mail: ibrohim@jinr.ru, ORCID: https://orcid.org/0000-0001-5534-3332, Scopus Author ID: 14829718300)
УДК 519.6
DOI: 10.22363/2658-4670-2021-29-1-5-13
Численное моделирование тепловых процессов, возникающих в материалах при воздействии фемтосекундных лазерных импульсов
И. В. Амирханов, Н. Р. Саркер, И. Сархадов
Лаборатория информационных технологий Объединенный институт ядерных исследований ул. Жолио-Кюри, д. 6, Дубна, Московская область, 1^1980, Россия
В работе проведено численное исследование решений параболического и гиперболического уравнений теплопроводности при одинаковых физических параметрах, а также сравнительный анализ полученных результатов. Обсуждена математическая постановка задачи. Действие лазера учтено через функцию источника, которую выбрали в виде двойного фемтосекундного лазерного импульса. В гиперболическом уравнении, в отличие от параболического, присутствует дополнительный параметр, который характеризует время релаксации потока тепла. Кроме этого, в источнике гиперболического уравнения присутствует дополнительное слагаемое — производная от плотности мощности источника параболического уравнения. Это означает, что на температуру образца оказывает влияние не только плотность мощности источника, но и скорости его изменения. Приведены профили температуры образца в разные моменты времени и её динамика на разных глубинах мишени. Расчёты проводились при различных временах задержки между импульсами и при различных параметрах релаксации.
Ключевые слова: параболическое и гиперболическое уравнения теплопроводности, фемтосекундный лазерный импульс, численное моделирование
