We model the latent heats of crystallization and fusion in phase change materials with a unified latent heat of phase change, ensuring energy conservation by coupling the heat of phase change with amorphous and crystalline specific heats. We demonstrate the model with 2-D finite element simulations of Ge2Sb2Te5 and find that the heat of phase change increases local temperature up to 180 K in 300 nm × 300 nm structures during crystallization, significantly impacting grain distributions. We also show in electrothermal simulations of 45 nm confined and 10 nm mushroom cells that the higher amorphous specific heat predicted by this model increases nucleation probability at the end of reset operations. These nuclei can decrease set time, leading to variability, as demonstrated for the mushroom cell.

Phase change memory (PCM) is a non-volatile memory technology that offers faster write speeds than flash, is CMOS compatible, and can be used for logic operations, allowing significantly faster computations by circumventing the Von Neumann bottleneck.1–5 PCM stores information in the crystalline (high conductivity) or amorphous (low conductivity) phase of materials such as Ge2Sb2Te5 (GST). We have previously reported a phase change modeling framework6–11 which captures stochastic nucleation, growth, amorphization, and grain boundary dynamics coupled with electrothermal physics to simulate PCM device operation. These models use the same specific heat [Cp(T)] for crystalline, amorphous, and molten phase change materials, sharply increasing Cp(T) near the melting temperature (Tmelt to Tmelt + 10 K) to account for the latent heat of fusion (ΔHf). In Woods et al.,7 we additionally incorporated the latent heat of crystallization using the value measured at the glass transition temperature [ΔHcrys(Tglass)]. However, the specific heats of amorphous and crystalline phases are expected to diverge for T > Tglass, and materials can crystallize at any temperature where Tglass < T < Tmelt. Since ΔHcrys(Tglass) ≠ ΔHf,12 phase change models to date which incorporate these latent heats do not conserve energy in amorphization-crystallization cycles. Here, we model both the latent heat of crystallization and the latent heat of fusion as the latent heat of phase change ΔHa,c(T). We ensure that energy is conserved in amorphization-crystallization cycles by coupling ΔHa,c(T) with specific heat. We consider melting to be the crystalline-to-amorphous transition and treat the melt as a continuum of the amorphous phase.

We demonstrate the proposed model with 2-D finite element simulations in COMSOL Multiphysics13 using the framework described in Woods et al.,6,7 which tracks grain orientation and crystallinity using a 3-vector CDCDǁ1 = 1 or 0 corresponding to crystalline or amorphous, respectively) with coupled rate equations

(1)

where CDi is a component of CD and Nucleation, Growth, and Amorphization are the component, phase, and temperature dependent. (1) is coupled with electric current14 

(2)

and heat transfer physics7,14

(3)

Here, J is the current density, σ is the electrical conductivity, V is the electric potential, S is the Seebeck coefficient, d is the mass density, and k is the thermal conductivity. QH accounts for the latent heat of phase change

(4)

where previously7 we set ΔHa,c(T) = ΔHcrys(Tglass) and (4) was non-zero only for positive dǁCDǁ1/dt (crystallization). Here, we update ΔHa,c(T) and use (4) for both crystallization and amorphization.

We use material parameters as in Woods et al.6,7 but update σGST, which affects electronic thermal conductivities (κe = LσT, where L = 2.44 × 10−8 W Ω K−2 is the Lorenz number), and thermal boundary conductances GTh based on our latest metastable amorphous15 and thin film16 measurements (Fig. 1). We model GTh as the sum of electronic and phonon contributions, based on measured values,17–19 as described in Silva et al.20 

(5)

GGSTTiNelec is expected to scale with κe and become dominant at elevated temperatures [Fig. 1(c)]. We also update ΔHa,c(T), Cp,c(T), and Cp,a(T) (c and a denoting crystalline and amorphous phases), as described below (Fig. 2).

FIG. 1.

σGST (a) and GTh [(b) and (c)] phonon (dashed), electron (dotted), and total (solid) components and measured values15–19 (markers).

FIG. 1.

σGST (a) and GTh [(b) and (c)] phonon (dashed), electron (dotted), and total (solid) components and measured values15–19 (markers).

Close modal
FIG. 2.

(a) Cp,a(T) and Cp,c(T) models and measured values.26 The shaded area is the enthalpy difference between phases minus ΔHcrys(Tglass). (b) Enthalpy difference between phases.

FIG. 2.

(a) Cp,a(T) and Cp,c(T) models and measured values.26 The shaded area is the enthalpy difference between phases minus ΔHcrys(Tglass). (b) Enthalpy difference between phases.

Close modal

The heat content of a system, enthalpy (H(T)), cannot be directly measured. However, its temperature derivative Cp(T) can

(6)

Enthalpy depends on phase and temperature. The enthalpy difference between the amorphous and crystalline phases, ΔHa,c(T), is released or absorbed during phase change

(7)

The change in enthalpy between two temperatures T1 and T2 is the integral of (6)

(8)

The change in ΔHa,c(T) from T1 to T2 is

(9)

Rearranging (9) and substituting in (8), we obtain

(10)

ΔHa,c(T) is typically reported as ΔHf at Tmelt and ΔHcrys at Tglass. ΔHf ≈ 68–129 J/g12,21–25 and ΔHcrys ≈ 35 J/g12,21,25 have been reported for GST and similar phase change materials. Here, we use reported GST values ΔHf = 128.9 J/g12 and ΔHcrys = 34.2 J/g12 for ΔHa,c(Tmelt) and ΔHa,c(Tglass). Substituting ΔHf = ΔHa,c(Tmelt) for ΔHa,c(T2) in (10), we get

(11)

We model Cp,c(T) with a quadratic function fit to GST measurements performed from 300 K to 650 K26 [Fig. 2(a)], with Cp,c(T) constant for T > 700 K (maximum of the quadratic). We assume Cp,c(T) = Cp,a(T) for T < Tglass, which is typical for glass formers,27 and model Cp,a(T) for T > Tglass as a constant such that (11) accounts for the difference in ΔHf and ΔHcrys(Tglass)

(12)

We then calculate ΔHa,c(T) using (11) [Fig. 2(b)].

In our simulations, we implement ΔHa,c(T) with (4), leading to local temperature variations during phase change (Fig. 3). We fix the initial and boundary temperatures (Tinit, Tbnd), allowing heat exchange at the boundaries of the simulated area but not out-of-plane. We use the steady state nucleation rate Iss(T) and growth velocity vg(T) from Burr et al.28Iss(T) is the rate at which stable nuclei appear and is distinct from the evolution of sub-critical nuclei (which we do not model) proposed by Loke et al.34 and others35,36 as a means of reducing set times. Iss(T) peaks at ∼600 K, while vg(T) increases until ∼720 K (Fig. 4). Hence, heating at crystal growth fronts increases vg for T < 720 K, increasing dǁCDǁ1/dt and creating positive feedback [Fig. 3(a)]. Heating growth fronts beyond 720 K decreases vg (negative feedback). If cGST is heated above Tmelt, heat absorbed during amorphization cools the crystalline-amorphous interface and slows the amorphization process [Fig. 3(b)].

FIG. 3.

Simulated crystal growth (a) and amorphization (b) with initial and boundary temperatures (Tinit, Tbnd) fixed at 600 K (a) and 900 K (b). The out-of-plane depth is 20 nm.

FIG. 3.

Simulated crystal growth (a) and amorphization (b) with initial and boundary temperatures (Tinit, Tbnd) fixed at 600 K (a) and 900 K (b). The out-of-plane depth is 20 nm.

Close modal
FIG. 4.

Average grain radius (blue markers, left axis) depends on Iss(T) and vg(T)28 (solid and dashed lines, scaled by maximum values, right axis). Radii are from 300 nm × 300 nm simulations using no heat of phase change, see Fig. 5 for simulation details.

FIG. 4.

Average grain radius (blue markers, left axis) depends on Iss(T) and vg(T)28 (solid and dashed lines, scaled by maximum values, right axis). Radii are from 300 nm × 300 nm simulations using no heat of phase change, see Fig. 5 for simulation details.

Close modal

Grain distributions after annealing are a function of both Iss(T) and vg(T): higher Iss leads to more, smaller grains, while higher vg leads to fewer, larger grains (Fig. 4). In 300 nm × 300 nm simulations with a 20 nm out-of-plane depth (Fig. 5) (Multimedia view), we obtain the smallest grains29 when annealing with Tinit = Tbnd ≈ 550 K: the increase in Iss(T) from 550 K to 600 K is outweighed by the increase in vg(T) (Figs. 4 and 5). This is true even when ΔHa,c(T) is not accounted for, but including ΔHa,c(T) results in fewer grains at or above Tinit = Tbnd ≈ 550 K as local temperatures increase up to ∼180 K above Tbnd (Fig. 5). ΔHa,c(Tglass) is released when crystallizing at any temperature in Woods et al.,6,7 which underestimates ΔHa,c(T) above Tglass and thus affects grain distributions to a lesser extent than the model presented here (Fig. 5).

FIG. 5.

Snapshots of simulated grain (left column) and temperature (right column) maps when annealing with Tinit = Tbnd = 650 K, comparing this work (a) and (b), where QH gives rise to ∼180 K increase in local temperature, to Woods et al.6,7 (c) and (d) and with no latent heat of phase change (e) and (f). (g) shows grains from 5 simulations at each temperature and QH model ranked by size. The out-of-plane depth is 20 nm. Grain analysis is performed with ImageJ.29 Multimedia view: https://doi.org/10.1063/1.5025331.1

FIG. 5.

Snapshots of simulated grain (left column) and temperature (right column) maps when annealing with Tinit = Tbnd = 650 K, comparing this work (a) and (b), where QH gives rise to ∼180 K increase in local temperature, to Woods et al.6,7 (c) and (d) and with no latent heat of phase change (e) and (f). (g) shows grains from 5 simulations at each temperature and QH model ranked by size. The out-of-plane depth is 20 nm. Grain analysis is performed with ImageJ.29 Multimedia view: https://doi.org/10.1063/1.5025331.1

Close modal

The relative contribution of the latent heat of phase change is a function of device structure, thermal boundary conditions, and operation speed and voltage. We simulate reset and set operations in nanoscale PCM cells to compare the presented model to our earlier work.6,7 There is no substantial change in the voltage or current necessary to reset a 45 nm × 50 nm confined cell (45 nm out-of-plane depth), and the temperature profiles while melting are similar using either model (Fig. 6). However, the maximum GST temperature hovers near Tmelt for ∼2 ns when cooling at the end of reset using our previous model due to the spike in Cp(T) from Tmelt to Tmelt + 10 K. This spike releases ΔHf, modeling a liquid-to-solid transition, but here we treat the molten phase as a continuum of the amorphous phase. The model presented here predicts a higher Cp,a(T) for Tmelt > T > Tglass, and thus, aGST cools more slowly in this temperature range and remains near the peak nucleation temperature for a longer duration [Fig. 6(b)]. This increases the nucleation probability while quenching at the end of reset (pre-nucleation), which can reduce set time and increase device-to-device and cycle-to-cycle variability30–33 (Fig. 7) (Multimedia view). These nuclei assist set by providing crystallization sites and establishing a preferential current path: current percolates through the conductive nuclei rather than flowing uniformly throughout the amorphous dome. In 30 reset simulations of a 10 nm mushroom cell (Fig. 7), the presented model shows more nuclei in 10 cases and fewer nuclei in 2 compared to the model in Woods et al.6,7 (Fig. 7). Pre-nucleation cannot be leveraged to directly reduce set time and voltage in conventional memory implementations, where a uniform waveform should set all cells, due to its stochastic nature. However, randomness introduced by pre-nucleation makes PCM cells suitable for hardware security applications (e.g., true random number generation, physical uncloneable functions, and physical obfuscated keys37–40)

FIG. 6.

Reset of a confined cell (45 nm out-of-plane depth) (a) with square Vapp(t), amplitude Vreset = Vdd = 6 V (b). Current (b) and temperature (c) in the GST lag slightly when modeling the latent heat of fusion as in this work (yellow spheres) or Woods et al.6,7 (blue spheres) compared to when the latent heat is neglected (solid line). Woods et al.6,7 releases the latent heat of fusion when cooling, giving a relatively flat maximum temperature from 8 to 10 ns (c). This work uses a higher Cp,a(T) for Tmelt > T > Tglass, resulting in slower cooling and more time spent near the peak nucleation temperature [highlighted region in (c)].

FIG. 6.

Reset of a confined cell (45 nm out-of-plane depth) (a) with square Vapp(t), amplitude Vreset = Vdd = 6 V (b). Current (b) and temperature (c) in the GST lag slightly when modeling the latent heat of fusion as in this work (yellow spheres) or Woods et al.6,7 (blue spheres) compared to when the latent heat is neglected (solid line). Woods et al.6,7 releases the latent heat of fusion when cooling, giving a relatively flat maximum temperature from 8 to 10 ns (c). This work uses a higher Cp,a(T) for Tmelt > T > Tglass, resulting in slower cooling and more time spent near the peak nucleation temperature [highlighted region in (c)].

Close modal
FIG. 7.

Reading and writing a PCM mushroom cell (10 nm out-of-plane depth) using the models in this work (a)–(d) and Woods et al.6,7 (e)–(h). This work uses a higher Cp,a(T) for Tmelt > T > Tglass, resulting in slower cooling and more probable nucleation after reset (c) and (g). This “pre-nucleation” reduces set time, as evidenced by the third read operation: a 40 ns pulse sets only the pre-nucleated device (d), (h), and (i). Pre-nucleation decreases set time; however, it is a stochastic process. Read, reset, and set are performed with square Vapp(T) pulses with magnitudes Vread = 0.1 V, Vreset= Vdd = 2 V, and Vset = 1.3 V. Multimedia view: https://doi.org/10.1063/1.5025331.2

FIG. 7.

Reading and writing a PCM mushroom cell (10 nm out-of-plane depth) using the models in this work (a)–(d) and Woods et al.6,7 (e)–(h). This work uses a higher Cp,a(T) for Tmelt > T > Tglass, resulting in slower cooling and more probable nucleation after reset (c) and (g). This “pre-nucleation” reduces set time, as evidenced by the third read operation: a 40 ns pulse sets only the pre-nucleated device (d), (h), and (i). Pre-nucleation decreases set time; however, it is a stochastic process. Read, reset, and set are performed with square Vapp(T) pulses with magnitudes Vread = 0.1 V, Vreset= Vdd = 2 V, and Vset = 1.3 V. Multimedia view: https://doi.org/10.1063/1.5025331.2

Close modal

In summary, we have applied thermodynamic principles to model temperature dependent latent heat of phase change and specific heats such that energy is conserved in crystallization-amorphization cycles. There is no distinction between the melt and amorphous phases in this model. Compared to our previous work, this model predicts a higher latent heat of crystallization above the glass transition temperature, and thus, the latent heat of crystallization has a larger impact on grain distributions. Electrothermal device simulations show that the amorphous phase remains longer near the peak nucleation temperature while cooling due to a higher specific heat, resulting in pre-nucleation at the end of reset pulses and thus stochastically decreased minimum set times. Variability introduced by pre-nucleation can be leveraged for hardware security.

This work was supported by AFOSR MURI under Award No. FA9550-14-1-0351.

1.
A.
Flocke
and
T. G.
Noll
, in
2007 Proceedings of 33rd European Solid-State Circuits Conference, ESSCIRC
(
2007
), p.
328
.
2.
D.
Kau
,
S.
Tang
,
I. V.
Karpov
,
R.
Dodge
,
B.
Klehn
,
J. A.
Kalb
,
J.
Strand
,
A.
Diaz
,
N.
Leung
,
J.
Wu
,
S.
Lee
,
T.
Langtry
,
K. W.
Chang
,
C.
Papagianni
,
J.
Lee
,
J.
Hirst
,
S.
Erra
,
E.
Flores
,
N.
Righos
,
H.
Castro
, and
G.
Spadini
, in
International Electron Devices Meeting (IEDM) Technical Digest
(
IEEE
,
2009
), pp.
1
4
.
3.
E.
Linn
,
R.
Rosezin
,
S.
Tappertzhofen
,
U.
Böttger
, and
R.
Waser
,
Nanotechnology
23
,
305205
(
2012
).
4.
C. D.
Wright
,
P.
Hosseini
, and
J. A. V.
Diosdado
,
Adv. Funct. Mater.
23
,
2248
(
2013
).
5.
N.
Kan'an
,
H.
Silva
, and
A.
Gokirmak
,
IEEE J. Electron Devices Soc.
4
,
72
(
2016
).
6.
Z.
Woods
and
A.
Gokirmak
,
IEEE Trans. Electron Devices
64
,
4466
(
2017
).
7.
Z.
Woods
,
J.
Scoggin
,
A.
Cywar
,
L.
Adnane
, and
A.
Gokirmak
,
IEEE Trans. Electron Devices
64
,
4472
(
2017
).
8.
A.
Faraclas
,
N.
Williams
,
A.
Gokirmak
, and
H.
Silva
,
IEEE Electron Device Lett.
32
,
1737
(
2011
).
9.
A.
Faraclas
,
G.
Bakan
,
L.
Adnane
,
F.
Dirisaglik
,
N. E.
Williams
,
A.
Gokirmak
, and
H.
Silva
,
IEEE Trans. Electron Devices
61
,
372
378
(
2014
).
10.
F.
Dirisaglik
,
G.
Bakan
,
A.
Faraclas
,
A.
Gokirmak
, and
H.
Silva
,
Int. J. High Speed Electron. Syst.
23
,
1450004
(
2014
).
11.
A.
Cywar
,
J.
Li
,
C.
Lam
, and
H.
Silva
,
Nanotechnology
23
,
225201
(
2012
).
12.
J. A.
Kalb
, “
Crystallization kinetics in antimony and tellurium alloys used for phase change recording
,” Ph.D. dissertation (Physikalisches Institut der Rheinisch-Westfälische Technische Hochschule (RWTH), Aachen, Germany,
2006
).
13.
“COMSOL Multiphysics Reference Manual, version 5.3”, COMSOL, Inc, www.comsol.com.
14.
G.
Bakan
,
N.
Khan
,
H.
Silva
, and
A.
Gokirmak
,
Sci. Rep.
3
,
2724
(
2013
).
15.
F.
Dirisaglik
,
G.
Bakan
,
Z.
Jurado
,
S.
Muneer
,
M.
Akbulut
,
J.
Rarey
,
L.
Sullivan
,
M.
Wennberg
,
A.
King
,
L.
Zhang
,
R.
Nowak
,
C.
Lam
,
H.
Silva
, and
A.
Gokirmak
,
Nanoscale
7
,
16625
(
2015
).
16.
L.
Adnane
,
F.
Dirisaglik
,
A.
Cywar
,
K.
Cil
,
Y.
Zhu
,
C.
Lam
,
A. F. M.
Anwar
,
A.
Gokirmak
, and
H.
Silva
,
J. Appl. Phys.
122
,
125104
(
2017
).
17.
J.
Lee
,
J. P.
Reifenberg
,
E.
Bozorg-Grayeli
,
L.
Hom
,
Z.
Li
,
S.
Kim
,
M.
Asheghi
,
H. S. P.
Wong
, and
K. E.
Goodson
, in
2010 12th IEEE Intersociety Conference on Thermal and Thermomechanical Phenomena in Electronic Systems (ITherm)
(IEEE,
2010
), pp.
1
6
.
18.
J.
Lee
,
E.
Bozorg-Grayeli
,
S.
Kim
,
M.
Asheghi
,
H.-S.
Philip Wong
, and
K. E.
Goodson
,
Appl. Phys. Lett.
102
,
191911
(
2013
).
19.
J. P.
Reifenberg
,
K.-W. C. K.-W.
Chang
,
M. a.
Panzer
,
S. K. S.
Kim
,
J. a.
Rowlette
,
M.
Asheghi
,
H.-S. P.
Wong
, and
K. E.
Goodson
,
IEEE Electron Device Lett.
31
,
56
(
2010
).
20.
H.
Silva
,
A.
Faraclas
, and
A.
Gokirmak
, in
Nanoscale Semicond. Memories Technol. Appl.
, edited by
S.
Kurinec
and
K.
Iniewski
, 1st ed. (
CRC Press
,
Boca Raton
,
2014
), pp.
307
332
.
21.
C.
Peng
,
L.
Cheng
, and
M.
Mansuripur
,
J. Appl. Phys.
82
,
4183
(
1997
).
22.
S.
Senkader
and
C. D.
Wright
,
J. Appl. Phys.
95
,
504
(
2004
).
23.
K.
Sonoda
,
A.
Sakai
,
M.
Moniwa
,
K.
Ishikawa
,
O.
Tsuchiya
, and
Y.
Inoue
,
IEEE Trans. Electron Devices
55
,
1672
(
2008
).
24.
J.
Kalb
,
F.
Spaepen
, and
M.
Wuttig
,
J. Appl. Phys.
93
,
2389
(
2003
).
25.
T.
Ohta
,
J. Optoelectron. Adv. Mater.
3
,
609
(
2001
).
26.
J.
Kalb
,
Stresses
,
Viscous Flow and Crystallization Kinetics in Thin Films of Amorphous Chalcogenides Used for Optical Data Storage
(
der Rheinisch-Westfalischen Technischen Hochschule Aachen
,
2002
).
27.
W.
Kauzmann
,
Chem. Rev.
43
,
219
(
1948
).
28.
G. W.
Burr
,
P.
Tchoulfian
,
T.
Topuria
,
C.
Nyffeler
,
K.
Virwani
,
A.
Padilla
,
R. M.
Shelby
,
M.
Eskandari
,
B.
Jackson
, and
B. S.
Lee
,
J. Appl. Phys.
111
,
104308
(
2012
).
29.
C. A.
Schneider
,
W. S.
Rasband
, and
K. W.
Eliceiri
,
Nat. Methods
9
,
671
(
2012
).
30.
U.
Russo
,
D.
Ielmini
,
A.
Redaelli
, and
A.
Lacaita
,
IEEE Trans. Electron Devices
53
,
3032
(
2006
).
31.
A.
Redaelli
,
D.
Ielmini
,
U.
Russo
, and
A. L.
Lacaita
,
IEEE Trans. Electron Devices
53
,
3040
(
2006
).
32.
U.
Russo
,
D.
Ielmini
, and
A. L.
Lacaita
,
IEEE Trans. Electron Devices
54
,
2769
(
2007
).
33.
M.
Rizzi
,
N.
Ciocchini
,
A.
Montefiori
,
M.
Ferro
,
P.
Fantini
,
A. L.
Lacaita
, and
D.
Ielmini
,
IEEE Trans. Electron Devices
62
,
2205
(
2015
).
34.
D.
Loke
,
T. H.
Lee
,
W. J.
Wang
,
L. P.
Shi
,
R.
Zhao
,
Y. C.
Yeo
,
T. C.
Chong
, and
S. R.
Elliott
,
Science
336
,
1566
(
2012
).
35.
T. H.
Lee
,
D.
Loke
,
K. J.
Huang
,
W. J.
Wang
, and
S. R.
Elliott
,
Adv. Mater.
26
,
7493
(
2014
).
36.
J.
Orava
and
A. L.
Greer
,
Acta Mater.
139
,
226
(
2017
).
37.
Security Opportunities in Nano Devices and Emerging Technologies
, edited by
M.
Tehranipoor
,
D.
Forte
,
G. S.
Rose
, and
S.
Bhunia
(
CRC Press
,
2017
).
38.
K.
Kursawe
,
A. R.
Sadeghi
,
D.
Schellekens
,
B.
Škorić
, and
P.
Tuyls
, in
2009 IEEE International Work. Hardware-Oriented Security Trust, HOST 2009, 22
Abstract (
2009
).
39.
L.
Zhang
,
Z. H.
Kong
, and
C. H.
Chang
, in
Proceedings of IEEE International Symposium Circuits Systems
(
2013
), Vol. 1, p.
1444
.
40.
L.
Zhang
,
Z. H.
Kong
,
C. H.
Chang
,
A.
Cabrini
, and
G.
Torelli
,
IEEE Trans. Inf. Forensics Secur.
9
,
921
(
2014
).