Thermal conductivity and thermal diffusivity of materials must be known at high accuracy for thermal engineering applications, in order to understand energy dissipation in devices and engines. Thermal conductivity libraries can be unreliable since their reported values may not reflect the complexities of the samples under investigation, and new materials may not yet be listed. Over the past three decades, modulated thermoreflectance (MTR) has emerged and evolved as a reliable, non-contact, and noninvasive technique to measure the thermal properties of a wide range of bulk and thin film materials and their interfaces. This Tutorial discusses the basis of, and recent advances in, the MTR technique, and its applications in the thermal characterization of a variety of materials. MTR experiments use two micrometer-wide laser beams. Experimental data (amplitude and phase vs the offset between the pump and the probe) give information about heat diffusion around the heat source along several tens of micrometers. Heat diffusion equations are used to fit the experimental data and extract the required thermal properties. Importantly, best fit procedures are not always needed because some analytical approximations provide evidence of the required thermal properties. We present many examples (bulk materials, layered sample, etc.) which illustrate this.

## I. INTRODUCTION

Nowadays, fundamental research demands accurate measurement of the thermal properties of complex structures at various scales. For instance, in integrated electronics and optoelectronics, power dissipates at the microscale, limiting the performance of devices and contributing significantly to global warming. In this context, the density of electronic components generates heat that must be evacuated by more sophisticated means than fans. Engineering optimization examines the choice of material and geometry of components. Consequently, a measurement of reliable and accurate local thermal properties is very important for applied and fundamental research. Unfortunately, thermal conductivity values reported in libraries and in the literature are often at odds with the ones observed in the samples being studied. The microstructure and state of the surfaces of samples have a strong influence on heat conduction, affecting the thermal conductivity (W/mK) and thermal diffusivity (m^{2}/s).

Heat can also be a source of valuable energy by harvesting dissipated power, even at the microscale. Research on thermoelectric devices has demonstrated that new geometries and materials may ameliorate the thermoelectric yield. For example, very thin layers exhibit very low thermal conductivities and are good candidates for efficient thermoelectric devices.^{1} Another example is the storage of hydrogen, which requires complex materials with high thermal conductivity. It is increasingly important to develop experimental setups to measure high and low thermal conductivities on very small samples, down to the micrometer scale.

Whatever the geometry, the principle of thermal measurement is always based on applying a heat source that diffuses in and around the sample, monitoring the diffusion, and comparing it with a model based on heat diffusion laws, to obtain the desired thermal properties. The heat source can use contact or contactless methods for massive, thin films, and even for nanowire samples. In the former case, the source can be a bolometer, strips deposited on the surface (3 Ω method)^{2} or the tip of a scanning thermal microscope.^{3} Local measurements of temperature are usually done by the device that creates the heat source. In the latter case, a non-contact method in the form of a photothermal effect usually creates the heat source: the surface of the sample is illuminated by a light beam whose intensity is unstationary (pulse or beam modulated in intensity) and whose geometry may be a point or a line, or may cover the entire surface. The sample may be opaque (surface source) or semi-transparent (volume source). Many physical parameters (refractive index, infrared emission, local deformation of the surface, etc.) are affected by heat diffusion from the modulated source, permitting temperature field measurement.

Modulated heat diffusion is characterized by a typical length, the thermal diffusion length. It must be compared with various other lengths—the size of the sample, thickness of the different layers, size of the source, and so on—in order to determine the depth reached by the heat in the sample. Complete heat diffusion equations, taking into account losses at the interfaces (radiation, convection, and conduction) and interfacial resistances, based on the Fourier law, have to be written into the geometry of the experiment in order to extract the desired thermal properties. Let us underline that, when very high frequencies or very small distances are envisaged, the mean free path of the heat carriers has to be considered. The lifetimes of the heat carriers must also be taken into account as they can be of the order of the pulse widths. In all these cases, the heat conduction equation and Fourier's law become inadequate. In the following discussion, we consider only the case where the local thermodynamic equilibrium is reached and the thermal flux is proportional to the temperature gradient (Fourier's law).

To see how a photothermally modulated excitation associated with local optical reflection measurements allows the thermal properties of samples to be determined, consider Fig. 1.

The total or partial absorption of light energy locally heats the sample's surface and triggers heat diffusion in the sample and in the surrounding media, if the de-excitation is not totally radiative [Fig. 1 (left)]. If a constant light intensity is used, the temperature field is strongly dependent on the geometrical boundary conditions (size support and surroundings). Typically, if the light beam intensity is modulated, a new parameter appears, called the thermal diffusion length, which is the length associated with the modulated heat diffusion regime. It is well known^{4} that in the case of a modulated diffusion regime, the heat propagation behaves as a wave (“thermal wave”) characterized by a wavelength of 2π*μ* and very strong damping. This thermal diffusion length *μ* (m) depends upon the thermal diffusivity *D* (m^{2}/s) and the modulation frequency *f* where $\mu =D/\pi f(m)$. We notice that the higher the frequency, the shorter is the diffusion length. As a consequence, in the modulated heat propagation regime, boundary effects due to the finite dimensions of the sample may vanish when the modulation frequency is increased.

As an illustration, Table I shows the thermal diffusion length in micrometers at different modulation frequencies for a range of materials, from a good thermal conductor (diamond) to a poor one (ZrO_{2}).

Sample . | Diffusivity D(m ^{2}/s) × 10^{−4}
. | Thermal diffusion length (μm)
. | |||
---|---|---|---|---|---|

1 Hz . | 1 kHz . | 1 MHz . | 1 GHz . | ||

Diamond | 10 | 17 840 | 570 | 18 | 0.57 |

Graphite | 1.8 | 7570 | 240 | 7 | 0.24 |

AlN | 0.5 | 3980 | 125 | 4 | 0.125 |

Air | 0.2 | 2500 | 80 | 2 | 0.08 |

Amorphous C | 0.05 | 1260 | 40 | 1 | 0.04 |

ZrO_{2} | 0.008 | 500 | 16 | 0.3 | 0.016 |

Sample . | Diffusivity D(m ^{2}/s) × 10^{−4}
. | Thermal diffusion length (μm)
. | |||
---|---|---|---|---|---|

1 Hz . | 1 kHz . | 1 MHz . | 1 GHz . | ||

Diamond | 10 | 17 840 | 570 | 18 | 0.57 |

Graphite | 1.8 | 7570 | 240 | 7 | 0.24 |

AlN | 0.5 | 3980 | 125 | 4 | 0.125 |

Air | 0.2 | 2500 | 80 | 2 | 0.08 |

Amorphous C | 0.05 | 1260 | 40 | 1 | 0.04 |

ZrO_{2} | 0.008 | 500 | 16 | 0.3 | 0.016 |

To explore the temperature field at the sample's surface, the temperature dependence of numerous physical parameters can be considered: emissivity, refractive index, local deformation, and so on.^{5–15} Refractive index variation and, more precisely, optical reflection coefficient variation, are considered in this paper. The optical reflection that is dependent upon only the local temperature will be taken into consideration in the following. This method is called Modulated Thermo Reflectance (MTR). Figure 1 (right) shows that the surface temperature field can be explored with a second beam reflected off the heated area.^{16–17} It is important to note that these fully optical MTR experiments are performed by non-contact, *in situ*, and non-destructive means by using two diode lasers.

Then, in addition to the frequency modulation, which fixes the characteristic diffusion length, MTR includes another experimental length: the pump-probe spacing. The collected data are amplitude and phase vs the offset between the two 1 *μ*m diameter spots, and modulation frequency is in the range of 1 Hz to 10 MHz.

Two lengths are thus experimentally available in MTR experiments: thermal diffusion length and probe beam spacing, which can be easily adjusted. Sweeps spacing around localized heating are fundamental for MTR. In this Tutorial, we will show that, most of the time, the thermal properties are deduced from the slope of the experimental data plots, without any absolute measurements of the heat flow nor of the probe diameter. In practical situations, scans of the spacing of the probe beams are recorded for different frequencies.

Moreover, by varying the thermal diffusion length on three orders of magnitude, we will show that a depth profiling analysis of multilayered samples can be obtained.

By playing with the spacing between the spots that ranges from a few micrometers to a few tens of micrometers and the frequency, it is then possible to determine thermal parameters at various spatial scales, which is fundamental to study heterogeneous microstructures, that is, ceramics and films with granular microstructures having different phases. The macroscopic thermal properties of these materials are strongly dependent upon the microscopic properties of their grains and grain boundaries.^{61}

It is important to report that a similar experimental scheme can be used in the time domain; that is, the so-called Time Domain Thermo Reflectance setup (TDTR).^{18} Indeed, a more complex setup is needed since the TDTR uses ps or fs lasers that deliver short pulses with a repetition rate of around ten ns (Ti-sapphire laser for instance). Each pulse of the Ti-sapphire laser is split into pump and probe pulses. The pump beam passes through a crystal that doubles the pump light frequency and the probe pulses are delayed relative to the pump. The pump and the probe overlap on the sample surface. The collected data are recorded vs time. Acoustic information is superimposed on a thermal profile between two excitation pulses. Let us underline that to increase the sensitivity of the method in a TDTR setup, the pump beam passes through an electro-optic modulator (EOM) and a lock-in is set to detect the TDTR signal at the EOM frequency. Let us underline that due to the high pulse repetition of the Ti-sapphire laser, there is insufficient time for the sample to return to thermal equilibrium between two laser pulses, and a thermal wave at EOM modulation frequency is created. So, the TDTR setup can also be used to measure thermal parameters at EOM frequency and is then called FDTR (Frequency Domain Thermo Reflectance). The experimental FDTR data are generally phase results vs frequency (10–100 MHz range).^{19} This frequency range prevents thermal waves from exploring thick samples in depth (see Table I). In FDTR, both beams are generally set on the same location on the surface. The diameters of the beams are about 10 *μ*m large. Recently, some experiments were reported where a spacing has been set between the two beams.^{26}

It is worth underlining that the same complete heat diffusion model based on thermal waves with the same thermal parameters was used in MTR and FDTR experiments.^{20,21} In FDTR experiments, a sophisticated best fit was used to extract all thermally isotropic or anisotropic parameters of bare substrates, films, and layered samples.^{22–27}Numerical sensitivity studies allow to comfort the thermal parameters determination.

To summarize the comparison between MTR and FDTR, in both experiments, the temperature field at the surface of the sample is associated with the propagation of a thermal wave created with the pump beam. The same heat diffusion equations are taken into consideration with the same thermal parameters. The main experimental differences are the size of the diameter beams, frequency ranges, and the complexity (and the cost) of the experimental setups.

The aim of this Tutorial is to describe the temperature field resulting from localized heating on the surface of the sample in order to discuss the sensitivity of MTR to extract different thermal properties from experiments. Whenever possible, we present analytical formula useful for these analyses and which allow a simpler qualitative and quantitative description of the temperature field.

Section II is devoted to an analysis of the heat equation in cylindrical symmetry for bulk samples, to enlighten the MTR experimental setup presented in Sec. III. Sections IV and V show theoretical and experimental MTR results on bulk and layered samples. In Sec. IV, experiments are presented on good and poor thermal semi-infinite samples (semiconductors Si_{3}N_{4} and Bi_{2}Se_{3}). We show that the scan of the heated area on tens of micrometers from the heat source directly gives the thermal diffusion length and the lateral thermal diffusivity. We also illustrate the thermal wave probing on heterogenous materials, that is, AlN ceramics, where a knowledge of thermal parameters at various spatial scales is necessary. Section V is devoted to layered samples. The layer can be a metallic film or not be one. We focus on examples with a high thermal contrast between the layer and the substrate. When the layer is much more conductive than the substrate, analytical solutions for the heat diffusion equations exist, whatever be the thickness of the good conductive layer. Moreover, we demonstrate that the relevant information is found when the probe beam scans at sufficiently large distances from the heat source. In the opposite case (i.e., a poor conductive layer on a good conductive substrate), the thermal properties of the layer and of the substrate can be measured by scanning the probe beam sufficiently close or far away from the heat source, respectively. Finally, we demonstrate that through the MTR experiments, we were able to achieve depth profiling analysis in multilayered samples.

## II. HEAT DIFFUSION EQUATION IN CYLINDRICAL GEOMETRY (BULK SAMPLE)

Prior to discussing the experimental setup of the MTR experiments, it is interesting to determine the temperature field at the surface of a semi-infinite bulk sample illuminated by a point source located at the surface, that is, *r = 0* and *z = 0*, as shown in Fig. 2, which describes the geometry of the experiment. The local temperature written *T(r, z, t)* is the solution of the heat diffusion equation written by adopting a cylindrical geometry^{7–14}

where *k* is thermal conductivity (W/m^{2 }K), and *D* is thermal diffusivity (m^{2}/s). In a homogenous sample, *k* and *D* are related by $D=k/\rho C$, where *ρ* (kg/m^{3}) is the mass density and *C* (J/kg K) is the specific heat.

In the case of a surface modulated point heat source $(q(r,0,t)=Q(2\pi r)\u22121\delta (r)\delta (z)exp(i\omega t))$, it is possible to include the source term contribution as a heat flow boundary condition at the interface. In the steady state regime, the thermal field can be decomposed into a zero frequency term and a harmonic one at the angular frequency *ω* = 2*πf*. While the former is principally related to the thermal losses of the sample (radiation, convection…), the latter can be related to the thermal properties of the sample. The amplitude of this term is analyzed below. To solve the partial differential equation in a modulated regime, a Hankel transform is performed over the space parameter (r → λ), leading to the following equation: $\u2202zzT~(\lambda ,z,\omega )=(\lambda 2+i\omega /D)T~(\lambda ,z,\omega )$. The *z* ordinary differential equation can then be solved: $T~(\lambda ,z,\omega )=Aexp(\u2212\sigma z)+Bexp(+\sigma z),$

where the quantity $\sigma =\lambda 2+i\omega /D$ corresponds to a wave vector. The parameters *A* and *B* are determined by considering the boundary conditions (equality of temperature and equality of the flux) at the interface, with the medium being the semi-infinite *B *= 0. Moreover, we suppose that there is no heat transfer in the upper medium, and the boundary condition on the flux gives $A=Q2\pi \sigma $, which is the point source contribution ($A=Q2\pi $ is the Hankel transform of the heat source).

The solution of the surface temperature field is then given by the inverse Hankel transform

where $J0(\lambda r)$ is the Bessel function of the first kind and the thermal diffusion length $\mu =2D/\omega $.

Then, the modulated temperature field *T(r, 0, t)* is obtained by the following formula:

This expression, corresponding to a point source, is called the Green function for the semi-infinite bulk sample. In the following, we ignore the trivial time harmonic dependence of the solution. It is important to note that the temperature at a distance *r* from the source is $r\mu $ out of phase with respect to the heat source.

We observe that the Green function is characterized by

a $exp(\u2212r/\mu )r$ amplitude dependency (3 D regime) and

the same spatial dependence, that is,

*r/μ*, on the phase and the logarithm of the amplitude

However, in real experiments, a pump beam presents a finite diameter size (Gaussian or Airy spot). This can be included by recalling the convolution theorem that permits us to replace the point source Hankel transform $Q2\pi $ in the integral of Eq. (2) with the Hankel transform of the actual heat flow multiplied by $2\pi $. In the case of a Gaussian beam of diameter *d* (beam diameter at 1/e^{2})*,* with a total incoming heat power *Q* (W), the heat flow at the surface of the sample reads $8Q\pi d2exp\u22128r2d2$, corresponding to the Hankel transform $Q2\pi exp\u2212d2\lambda 232$. Then taking into account the finite source size, the thermal field can be written as

There is no analytical expression of the surface temperature field. But, if the beam diameter is much smaller than the thermal diffusion length, the numerical solution is similar to the Green function above. In Fig. 3, the thermal field corresponding to a Diraclike beam (*d* = 10 nm) is well approximated by the Green function solution, while, for larger *d* values, convolution effects are clearly observed.

Observe that the beam diameter size effect is relevant at small *r* values (±2 *μ*m) (Fig. 3). Far from this region, the logarithm of the amplitude [log (Ampl)] and the phase shift exhibit the same slope ∼1/*μ* [see Eq. (2) and the insets to Fig. 3]. Note that the phase measurements can be done just outside the diameter of the beam. In regard to log (Ampl), the 1/*μ* slope appears far away from the center, at distances that represent nearly twice the thermal diffusion length in this example. Consequently, the thermal diffusivity can be determined for this semi-infinite sample, provided that the measurement is done sufficiently far away from the heat source. We stress here that no absolute measurement is needed.

Below, we focus on the region that is close to the center. Of course, the diameter of the beam is an important parameter, as shown in Fig. 3: the smaller the diameter of the beam the higher the maximum temperature. Equation (4) permits us to determine the following analytical expression of *T (0,0,ω)*:

In Fig. 4, the amplitude and the phase of the temperature at *r *= 0 are plotted. Interestingly, the amplitude and the phase are strongly dependent upon the beam diameter and the modulation frequency [see curves in black (*d *= 2 *μ*m) and green (*d *= 5 *μ*m) for silica].

We can notice that the low-frequency behavior that corresponds to a 3D regime (*d < μ*) is given by

Equation (6) is particularly useful for estimating the steady state temperature rise at *r = 0*

At high frequency, Eq. (5) can be written as

The beam diameter becoming larger than the thermal diffusion length [Eq. (8)] is the solution for one-dimensional heat flow perpendicular to the surface with a uniform heat flux $8Q/\pi d2$. The phase is then equal to −45° corresponding to the 1D regime.^{4}

To conclude, in order to determine the thermal diffusivity of a semi-infinite sample [Eqs. (2) and (4)] without any absolute measurement, scans around the heat source should be carried out as far as the signal-to-noise ratio allows. The slopes of these scans are direct measurements of the thermal diffusivity.

The thermal properties *k* and *D* can also be determined by measuring the point source temperature *T (0, 0, ω)* with absolute measurements of *Q* and *d* at different modulation frequencies {allowing the determination of *k* at low frequency [Eq. (6)] and of the effusivity $K/D=k\rho C$ at high frequency [Eq. (8)]}.

It is also important to note that lateral heat diffusion governs low frequency experiments (3D regime *μ *> > d), while perpendicular diffusion governs high frequency experiments [1D regime (*μ *< < d) with phase −45°, amplitude as Eq. (8)].

## III. MTR MICROSCOPE

In 1985, Allan Rosencwaig described the first MTR microscope in order to measure the quality of silicon wafer (dopant concentration) with two superimposed beams impinging on the surface.^{16,17} In Sec. II, we showed that scans far from the heat source are needed in order to measure the thermal diffusivity of the semi-infinite samples. As a consequence, the configuration reported in Refs. 16 and 17 (two superimposed beams) is not appropriate for determining the lateral thermal diffusivity of the solid samples, and lateral scans are needed.

Figure 5 illustrates the configuration of an MTR microscope^{28} allowing scanning of the heated area. The setup allows the sample and beams to be observed (Olympus BX-RFA). The pump beam passes through an acousto-optic modulator (MT80-A1 AA Opto-Electronic) and is then reflected on the surface sample with an oscillating mirror, allowing scans to be made. The probe beam is directed onto the surface sample and then reflected onto the photodiode (1801 FS 125 MHz New Focus). An AC signal is recovered with lock-in detection (Signal Recovery 7280). In the case of a semi-infinite sample, fixing the probe and scanning the pump, and the reverse, are fully equivalent. Nevertheless, the fixed probe and scanning pump configuration is more relevant when the surface sample does not present good optical quality everywhere.^{28}

Now, we discuss the choice of the two lasers. First, the pump laser has to be absorbed ($Rpump$ small). Second, the probe beam has to be reflected ($Rprobe$ large). Moreover, in this modulated experiment, the signal is proportional to the optical reflection variation $\Delta Rprobe$,

Only $\u2202Rprobe/\u2202T$ has to be taken into consideration in an MTR setup. To get a photon noise limited setup, the reflected probe beam on the photodiode has to be correctly adjusted. The powers of the probe laser and the optical reflection coefficient have to be matched. Actually, the light absorption is proportional to $(1\u2212Rpump)$, but the signal is also strongly dependent on the substrate.

To achieve experiments on a large range of thermal conductivity samples, with a given laser, it is not reasonable to fix $Rpump$ around 1. To avoid the variation in the optical properties of different samples (transparency, small absorption for the pump wavelength or diffusing sample), a thin layer (100 nm) of gold is often deposited on the sample surface. So, in the setup shown in Fig. 5, the pump beam is a 532 nm Cobolt laser ($Rpump$ = 0, 6) and the probe beam is a 488 nm Oxxius laser diode ($Rprobe$ = 0, 5), $\Delta R\Delta T$ at 488 nm corresponding to a maximum of (5 × 10^{−4}) for gold (Fig. 6). The power applied to the sample is typically in the range of 10–100 *μ*W for the probe beam, and a few mW for the pump beam, both with a spot size of about 1 *μ*m^{2}. In each experimental record, the linearity of the signal is checked. Usually, a ×50 (NA = 0.5) objective is used and the spots on the surface are diffraction limited. The photon noise is about 10^{−7}, so the minimum temperature variation that can be measured is 10^{−3} if $dRdT$ is 10^{−4} K^{−1}.

In Sec. II, we saw that the heated area has to be scanned in order to measure the slopes of the amplitude and phase vs the distance of the two beams. Typically, ± 40 *μ*m scans can be recorded. The time constant of the lock-in detection can be adapted to the signal-to-noise ratio: averaging can be made higher in the wings (small signal) than in the center (high signal).

As compared to the previous expressions of the temperature given in Eqs. (4) and (5), the beam diameter of the probe beam implies a second convolution. If both beams are supposed to be Gaussian, the convolution of two Gaussians is still a Gaussian, with an effective diameter $d=dpump2+dprobe2$. Finally, Eqs. (6) and (8) or Refs. 20 and 29 give an approximation of the DC temperature in r=0, that is, $T(0,0,0)=2Qkd12\pi $for a semi-infinite sample [Eq. (7)].

The working frequency ranges from 1 Hz to 10 MHz. An important consequence of the frequency dependence of the thermal diffusion length is that the larger the modulation frequency, the smaller is the thermal diffusion length. This implies that MTR allows extracting thermal properties close to the surface (high frequencies) or in-depth (low frequencies). This behavior is very useful for heterogeneous samples and is an important benefit of the MTR method.

## IV. BULK SAMPLE EXAMPLES

As introduced in Sec. II, the thermal diffusivity of bulk samples can be obtained by lateral/radial thermal diffusivity measurements (3D regime). In this section, we consider [Eq. (4)] only opaque samples where heat is applied to the surface (i.e., z = 0). Light penetration into the sample affects the response close to the pump light spot (enlargement of the heat source due to optical penetration).^{23} Fortunately, the thermal diffusivity values can be easily extracted at sufficient distances.

Figure 7 shows the MTR experimental results obtained on a GaAs sample at 100 kHz. It is important to recall that the optical reflection coefficient depends on the local temperature and also on the photo-created carrier density when semiconductors are under study. In fact, thermal diffusion can be blurred by carrier diffusion when studying semiconductor samples. In the case of GaAs, photo-created carriers do not interfere with the thermal measurement, since their lifetime is very small (20–300 ns) compared with the modulation frequency (1 MHz). Superimposed onto the circles that represent experimental data are lines showing calculations obtained with Eq. (3). The calculation parameters are as follows: 100 kHz, 1 *μ*m diameter beam, *k *= 45 W/mK *D* = 2, 6 × 10^{−5} m^{2}/s. To illustrate the accuracy of the method, three numerical simulations are superimposed onto the experimental results (=2, 2.6, 3 × 10^{−5} m^{2}/s). The precision is around 5% on the thermal diffusivity measurement.

Note that the thermal diffusivity can also be measured on the slope of the phase (between 4 *μ*m and 20 *μ*m) or the log (Ampl) far from the center (10–20 *μ*m) without any best fit. The slope being $1/\mu $, the thermal diffusivity is given by $D=\pi f\mu 2$.

By varying the modulation frequency, it is possible to corroborate the fit obtained. Indeed, we have already mentioned that heat diffusion can explore the sample in depth by reducing the frequency. This enables the homogeneity of the sample to be checked by running experiments at different frequencies to confirm the sample microstructure. Figure 8 shows MTR amplitude and phase data from a Bi_{2}Se_{3} sample at 50–100–200 kHz.^{30} The thermal diffusivity *D* is obtained from the slope of the phase data. Because the *ρC* product has been previously determined, k is obtained by *k = DρC.* MTR amplitude and phase data can be recorded on more than 3 orders of magnitude (Figs. 7 and 8) with r scans larger than 20 *μ*m. Moreover, in the case of a uniaxial anisotropic bare substrate under investigation, only the radial thermal diffusivity is measured (see Fig. 20 and Sec. V C). The data obtained on Bi_{2}Se_{3}, presented in Fig. 8, are, in reality, measurements of radial diffusivity, since Bi_{2}Se_{3} has a strong thermal anisotropy.^{30}

In the following, we will discuss thermoreflectance measurements performed in heterogeneous materials, particularly in ceramics where MTR has been demonstrated as a valuable tool to extract thermal properties at various scales.^{61} These materials are puzzling from the point of view of thermal properties. The thermal properties are spatially range-dependent. Measurements at the millimeter scale do not give the same values as micrometer range experiments. Ceramics include many defects (such as secondary phase grains and thermal barriers at grain boundaries), so it is legitimate to measure the thermal properties at different spatial ranges in order to determine the limiting factors of heat diffusion. A good example is industrial AlN ceramics. Figure 9 (right) shows the microstructure of this type of material obtained by sintering. The secondary phase Y_{2}O_{3} grains (white) and the grains of AlN (gray) appear with a well-defined geometrical distribution. The MTR microscope allows the measurement of the thermal conductivity of the AlN grains and that of the thermal conductivity of the second phase grains at the micrometer range and also evidences the presence of thermal barriers at grain boundaries.^{31–33} Explorations at the macroscopic scale require another experimental scheme. For instance, a mirage experiment^{34,35} can give heat diffusion information averaged over many grains (over a few hundred micrometers) at low modulation frequency, whereas MTR can provide information at the grain scale and at high modulation frequency. As the grains are around 10 *μ*m in the AlN ceramics under investigation, it is easy to obtain information on the thermal properties at the grain scale (Fig. 9, left). The data at the grain scale associated with measurements at lower modulation frequencies (where a lot of grains are involved) allow us to study the influence of the microstructure and grain distribution on heat propagation. The role of a thermal barrier at the grain boundaries can be a limiting factor. Figures 10 and 11 show the detection of thermal barriers associated with grain boundaries or a secondary phase grain with the MTRM setup. In this case, the two beams have to be superimposed and then the sample scanned.

All these data provide evidence that grain boundaries are, for this sample, the ultimate limitation of heat diffusion. Moreover, it is possible to predict the effective thermal conductivity, thanks to a mesoscopic modeling of the intergranular structure of Y_{2}O_{3}-doped AlN.^{31}

Finally, we show MTR measurements of silicon nitride ceramics, which present excellent mechanical properties associated with a rather high thermal conductivity.^{36–38} Silicon nitride fabricated by tape casting exhibits a highly anisotropic microstructure with large and elongated grains, uniaxially oriented parallel to the casting directions shown in Fig. 12. Experiments done on two samples (parallel to *a*-axis and parallel to *c*-axis) have confirmed the expected anisotropy.^{38}

To conclude, MTR on bulk samples permits lateral thermal diffusivity without any absolute temperature measurement (using slope measurements). Results can be obtained easily on very small areas (such as ceramics whose grains are no larger than a few micrometers, rough samples, and so forth). Whatever be the spatial scale, thermal conductivity can be deduced only after an absolute measurement of the temperature or by a precise evaluation of the *ρC* product.

## V. LAYERED SAMPLES: HEAT DIFFUSION IN CYLINDRICAL GEOMETRY

In this section, we consider samples consisting of a layer and a substrate characterized by different thermal properties. There are two reasons to focus on these systems. The former concerns the already-mentioned thin layer (10–100 nm) of gold that is often deposited on the surface of a bulk sample to prevent light diffusion in the sample. Nevertheless, we will show that it is still possible to get thermal information about the underlying substrate. The latter motivation is that thin layers deposited on a substrate are ubiquitous in modern technology. Consequently, it is extremely important to evaluate the thermal properties of very thin layers, for example, 10 nm–10 *μ*m, in optoelectronics, electronics, and optics and for thermoelectric applications.^{1} It is therefore necessary to solve the heat diffusion equations for a layered structure in cylindrical symmetry (Fig. 13). The layer is characterized by its thickness *h*, its thermal diffusivity *D _{1}*

_{,}and its thermal conductivity

*k*, and the substrate by

_{1}*k*and

_{0}*D*. The diameter of the Gaussian beam is

_{0}*d*at 1/e

^{2}. The heat source is assumed to be applied to the surface, and losses in air and interface resistance between the layer and the substrate are ignored. The surface temperature field is given by Eq. (7),

^{39}

Let us underline that no analytical solutions associated with the surface temperature field can be obtained from Eq. (10). Only a numerical integration can determine the surface temperature field, with a best-fit procedure to determine the unknown parameters. While a best fit can give answers even if all parameters are unknown, two cases have to be discussed to refine the solutions: either the layer properties or the substrate ones are to be determined. In both cases, the beam diameter has to be considered as an unknown parameter, as it is usually difficult to determine accurately. In the former case (layer properties unknown, substrate properties known), the thermal properties of the layer (*k _{1}* and

*D*) and the beam diameter

_{1}*d*are the best fit parameters, while in the latter case (layer properties known, substrate properties unknown), the thermal properties of the substrate are the best fit parameters together with the beam diameter. This discussion highlights that the experimental data have to be recorded with very high accuracy near the center and along the wings (far from the center). In the next sections, we demonstrate that the relevant information could be distributed either in the center or along the wings, depending on whether the thermal conductor is above or below the thermal insulator.

Obviously, in the case of a layered structure, heat propagation is strongly dependent on the thermal contrast α = $k0/k1$. Figure 14 shows what happens in two extreme cases (α being very large or very small): insulator on conductor (IOC) and conductor on insulator (COI). In Fig. 14 (left), the heat propagates through the insulating layer with negligible lateral diffusion. Once it reaches the conducting substrate, the heat propagates through it in all directions. The substrate then heats the surface through the insulating layer over a broad area, not just where the heat was applied. In Fig. 14 (right), the heat disperses rapidly through the conductive layer and can be described as trapped in this layer. Heat barely propagates in the substrate.

A comparison of different lengths (*h, μ _{1},* and

*μ*) and the thermal contrast is the main key to understanding heat propagation in the layered structure. Sections V A and V B are devoted to the study of these two cases, with analytical solutions where possible. Of course, when $h\u226b\mu 1$, the COI and IOC solutions merge into the bulk solution given in the previous section.

_{0}### A. Conductive layer on insulator substrate (COI)

Most conductive layers are metallic. Many papers are devoted to understanding heat propagation in a layered metallic sample (see, for example, Refs. 39–47). In this section, two issues have to be distinguished. In the first case, the substrate's thermal properties must be determined despite the presence of the metallic layer, while in the second case, the layer's thermal properties are the subject of the study. Moreover, it is important to recall that in metals, the measured thermal conductivity (*k*), in an MTR experiment, is the sum of two distinct contributions: the electronic and the phonon thermal conductivities, *k _{el}* and

*k*, respectively. The Wiedemann–Franz law permits to estimate the electronic contribution from electrical conductivity leading to

_{ph}*k*evaluation.

_{ph}In order to understand the thermal behavior of thin metallic layers (gold for instance) deposited on insulator substrates (silica for instance), let us start with the amplitude of surface temperature simulation [Eq. (4)] of a black silica bulk sample and a gold bulk sample at 1 MHz (Fig. 15). As expected, the two maxima are dependent upon the beam diameter, the modulation frequency, and the thermal diffusion length, and are inversely proportional to the thermal conductivity [Eq. (5)].

When a gold layer is deposited on silica, by varying the conductive layer thickness with respect to the thermal diffusion lengths of the substrate and of the layer itself, we can point out three cases (Fig. 16):

The layer is larger than the thermal diffusion length

*μ*_{1}of the layer: $h>\mu 1$ $(h=7.5\mu m>\mu 1=6.3\mu m)$. As already mentioned, this thermally thick regime corresponds to the semi-infinite behavior (red line in Figs. 15 and 16).The layer is very thin: $h\u226a\mu 1(h=0.5nm\u226a\mu 1=6.3\mu m)$: silica behavior is recovered (black line in Figs. 15 and 16).

For intermediate thicknesses (10–100–500–1000 nm): ($h<\mu 1$ thermally thin regime): the simulation reflects the schematic heat diffusion regime reported in Fig. 14 right (heat trapped in the upper layer).

Let us concentrate now on the intermediate region (10 nm < *h*< 1000 nm) in Fig. 16. In Ref. 41, we derived an analytical expression of the asymptotic (larger) surface temperature. In a thermally thin regime $h<\mu 1$ with a large contrast$\alpha <0.1$, the analytical expression for the asymptotic surface temperature is

with

and

Large *r* bulk behaviors are distinguished by

a $1/r$ amplitude dependency that corresponds to a 2D diffusion regime.

The phase and the logarithm of the amplitude present different slopes because the real and the imaginary part of the complex factor in the exponential are different.

Equation (11) expresses that, at a given modulation frequency, by performing a single measurement scan of the phase and of the amplitude, we can determine

*k*and_{1}*D*independently (without knowledge of the density and the specific heat) if the thickness h of the layer and the effusivity [$k0/D0$] of the substrate are known._{1}the thermal effusivity of the substrate ($k0/D0$) if the thickness and the thermal properties of the layer are known (h, k

_{1,}and D_{1}).

Interestingly, by collecting three or four scans at different frequencies, it is possible to determine *k _{1}* and

*D*independently and without the help of numerical simulations:

_{1}^{41}

Combinations of the slopes of the plots of the logarithm of the amplitude and of the phase vs the pump probe distance can be achieved (Fig. 17). Two linear graphs are obtained, one vs $\omega $ and the other vs$\omega $. The slope of the first is inversely proportional to

*D*, while the slope of the second is a combination of the thermal conductivity_{1}*k*, and the thermal effusivity of the substrate ($k0/D0$). The thermal diffusivity and the thermal conductivity of the layer are obtained without any absolute measurements and independently._{1}, hThe thermal diffusivity

*D*of the layer can be obtained even if the substrate and layer thickness, and the thermal conductivity of the sample, are unknown (Fig. 17). The slope q of the combination $b2\u2212a2+2ab=\omega D=q\omega $ permits_{1}*D*to be determined._{1}

It can be shown that the limit of the asymptotic behavior [Eq. (11)] is given by $\alpha \mu 0<h2<min\mu 0\alpha ,\mu 1$. In the case under investigation, that is, gold layer on silica, the thickness range is bracketed between; then, 2 nm *< h* < 6.3 *μ*m for a 1 MHz experiment. The condition for the thickness is then $\alpha \mu 0<h2<\mu 1$.

These experiments using a metallic layer deposited on an insulator substrate clearly evidence the presence of interface resistance.^{46–50} For instance, for four different thicknesses of Al on silica and for four different modulation frequencies for each layer, the above combinations of slopes are given in Fig. 18. The thermal diffusivity obtained for the 55 nm Al layer is larger than the thermal diffusivity of the bulk Al, which is clear evidence of the presence of an interface resistance that prevents heat diffusing in the substrate.^{47} The model is thus incomplete and a two-layer model has to be calculated to take into account the interface. Note that, when an interface resistance is suspected, we have just to replace $1/k0\sigma 0by(1+Rk0\sigma 0)/k0\sigma 0$ in Eq. (10), to take this discontinuity into account.

Next, let us consider the deposition of very thin layers (nanometer range), as are regularly used to protect samples from oxidation (Au capping layer). For applications, it is useful to evaluate the role played by this thin layer on the overall thermal properties. It can be shown that if the thickness is characterized by$h2\u226a\alpha \mu 0$, as in the case of gold on silica (h < < 2 nm), a very good approximation of the thermal behavior can be written as

Here, H_{0} is a modified Struve function and Y_{0} a Bessel function of the second kind. The temperature field (Fig. 19) is the sum of a pure substrate term and of a cap layer perturbation depending on the thermal contrast and the layer thickness. This is a situation where the cap layer acts as a small perturbation.

To conclude this part, when the 1/√r regime is observed on the experimental amplitude, it means that the layer is more conductive than the substrate, and we use an analytical expression of the temperature field to evaluate the thermal properties of the thin layer.

On the other hand, when the layer thermal properties are known, only the thermal effusivity ($k/D$) of the substrate can be extracted from experimental data. Importantly, as the most useful information is contained in the temperature asymptotic behavior, the size of the beam diameter plays no role. It is always possible to achieve a best fit with five parameters in the experimental data to obtain the thermal properties of the layer, substrate, and beam diameter. The presence of an interface resistance is easy to detect in the case of a conductive layer deposited on an insulator.

If the layer thickness is very small, we obtain an analytical expression that shows that the layer acts as a small perturbation for the thermal temperature field. In this case, scanning the temperature field allows us to obtain the thermal diffusivity of the substrate. For very thick metallic layers, the substrate plays no role, and we obtain the thermal diffusivity of the layer. Moreover, the lateral diffusion in the metallic film means that the temperature can be decreased efficiently. A crude estimation of the DC temperature (r = 0) rise is then given by $Qh/\pi kd2$.^{51}

This manuscript presents only the temperature field associated with surface heat source. The optical penetration of the pump beam has to be taken into consideration, and this leads to volumetric heating. Nanometer range optical penetration up to a few tens of nanometers in a metallic layer affects the central part of the probe exploration (near *r *= 0) but does not significantly modify the slopes of the lateral heat diffusion (large *r*). This is a clear advantage of a pump-probe separation scanning when the thermal properties that are under investigation control the long distance temperature field (asymptotic behavior in the case of a good conductor over a layer deposited on a poor conductor substrate, for example). So, we can avoid more complicated diffusion equations.

Finally, when a uniaxial anisotropic substrate is covered with an isotropic metallic layer (Fig. 20, right), the perpendicular thermal properties of the substrate are measured, while on the bare substrate [Fig. 20, left], the parallel properties are measured (see Sec. V C). This particular case has been discussed in Ref. 30 where oriented Bi_{2}Se_{3} was covered by a very thin film of gold.

### B. Insulator layer on conductive substrates (IOC)

As shown in Fig. 14, we expect that the central part of the data, near the pump beam, reflects the thermal layer properties, while bulk properties could be extracted from the wings. Nevertheless, an important point has to be considered. The insulating layer can be semitransparent: the optical penetration of the pump beam has to be taken very carefully into consideration, since the profile of the temperature field is determined by the penetration of the pump beam. Consequently, Eq. (15) has to be corrected.^{23,24}

Below, we limit our presentation to experiments where the optical penetration is much smaller than the thicknesses of the layer so that we can ignore the beam optical penetration.

As in Sec. V A, we start by comparing the amplitude and phase of the bare substrate samples: amorphous silicon (red circles) and silicon (black circles) (Fig. 21), which are poor and good thermal conductors, respectively. When a thin a-Si layer (10 nm–100 nm) is deposited on a silicon substrate, the temperature field is strongly affected near the heat source. As expected, sufficiently far away from the point source, the substrate contribution is recovered, as shown in Fig. 21. We note that the amplitude is strictly the amplitude of the bulk, while the phase is slightly shifted but is parallel to the bulk phase. On the other hand, close to the source point, the amplitude and the phase are strongly affected by the presence of the thermally insulating layer: the amplitude is strongly increased, reflecting the Gaussian beam shape and attesting that there is only weak diffusion in the layer outside the beam diameter.

To explain this behavior, Eq. (10) has to be numerically solved. Then, five parameters can be obtained by a best-fit procedure, that is, *k _{0}, k_{1}, D_{0}, D_{1}*, and

*d*. This is simplified when the thermal properties of the substrate are known.

with $\alpha \u2032=k1/k0=1/\u221d$ and $Z=\sigma 1/\sigma 0$. Note that $\alpha \u2032=1/\u221d$ is small. It is evident in Eq. (15) that the surface temperature field is the superposition of two terms, the first for the substrate and the second for the layer.

If α’ is very small (insulator on conductor) for a thin layer $h<\mu 1$ and $h<d$, the approximation $\alpha \u2032|Ztanh(\sigma 1h)|$ < < 1 allows us to approximate Eq. (15) as

If $h\u226ad$ and $h\u226a\mu 1$ are up to the first order, we can write

It is very important to notice in Eq. (17) that first, the thin film contribution simply adds to the substrate's contribution and second, the integration corresponds to the Gaussian beam,

Consequently, the thermal conductivity of the layer is straightforward to deduce from the real part of the experimental data by a best-fit evaluation of *d* and *k _{1}.* Since the real part is very large compared with the imaginary part, the thermal conductivity can be deduced from the amplitude with high accuracy (Fig. 22).

At the second order, the solution of Eq. (16) could be approximated as

Figure 23 shows the result of three simulations, Eqs. (16), (17), and (19). As clearly demonstrated in Fig. 23 (left), the approximation at the first order is very robust. We recall that on the amplitude/real part, only the beam diameter and thermal conductivity of the layer are pertinent parameters. In regard to phase measurement (Fig. 23, right), we point out that the fitting procedure is more delicate since the approximations [Eqs. (17) and (19)] are relevant only for very thin layers (*h* < 50 nm for this case).

A very simple model [Eq. (17)] can be applied when *h < < μ*1 and *h < < d*, leading to a straightforward evaluation of the layer's thermal conductivity from the real part of the experimental data, just by fitting two parameters, that is, *d* and *k _{1}*. The evaluation of thermal diffusivity is trickier, and experimental data must be fitted using Eq. (15).

For thicker layers only, Eq. (15) has to be used with a best fit with beam diameter, thermal diffusivity of the layer, and thermal conductivity as parameters, if the substrate is known. Experiments at various frequencies and/or various thicknesses are needed to confirm the results. When the amplitude/real part does not present a Gaussian shape (lateral diffusion), only a numerical calculus [Eq. (15)] can be fitted to the data.

Note that the presence of an interface resistance is more difficult to detect than for the COI case. In Eq. (17), the interface resistance R appears clearly by replacing $1/k0\sigma 0$ by $(1+Rk0\sigma 0)/k0\sigma 0$,

The insulator thin layer acts as a resistance and we measure the addition of these two resistances.^{44,45} If the interface resistance value is in the same range as the layer equivalent resistance, it is wise to measure several samples with different layer thicknesses. When the interface resistance is much smaller than the equivalent layer resistance, its evaluation is much more difficult.

An example of an insulating layer on a good conductor substrate is a titanium layer on silicon.^{52} Titanium layers are often used as electrodes in electronic or optoelectronic devices. From a technical point of view, it is very important to know the thermal properties of each part of the device, particularly the electrodes. We have performed MTR measurements for 180 nm titanium deposited on silicon at different modulation frequencies: 50 kHz, 150 kHz, 450 kHz, and 850 kHz. Experimental data are given in Fig. 24 (points). The shape of the amplitude suggests the scenario of a poor conductor deposited on a good conductor substrate.

First, we fitted the real part of the data [Eq. (17)] in order to obtain the beam diameter and the thermal conductivity of the Ti layer.^{52} On the 50 kHz run, we deduced a value of the thermal conductivity of titanium k_{1} = 18.7 W/mK with d = 2.75 *μ*m. Then, we achieved a best fit (lines on Fig. 24) on the experimental data with only one parameter, that is, the thermal diffusivity D_{1}.

The fit of the phase contribution is not perfect (Fig. 24). It can be improved by a two-layer model considering the oxidation of the Ti surface.^{53} This is a good demonstration of the sensitivity of the method.

Another example concerns thermoelectrical materials. Nowadays, new materials are investigated to improve thermoelectrical efficiency. As we have already mentioned, thinner layers in a structure are good candidates since *k* is expected to reduce, increasing the so-called ZT factor.^{1} To explore this, we studied very thin layers (18–191 nm) of Bi_{2}Se_{3} epitaxially deposited on a GaAs (111) substrate.^{54} In this case, we know the ratio between the thermal conductivity and the thermal diffusivity. The best fit was done only on *k _{1}* and the beam diameter as parameters.

Figure 25 shows experimental data for thin layers of Bi_{2}Se_{3} epitaxied on GaAs (111) substrates by molecular beam epitaxy (MBE) for two types of thicknesses, 18 nm and 191 nm. It is interesting to measure the thermal conductivity of different thicknesses of layer in order to confirm their ability to enhance thermoelectric power. We found that in thin films, thermal conductivity measured at room temperature monotonically decreases when reducing film thickness (Fig. 26). We can attribute this reduction to incoherent scattering of phonons whose momentum is out of plane with the film's upper and lower surfaces. This work outlines the crucial role of sample thinning in reducing the out-of-plane thermal conductivity.^{54}

Note that the shape of the MTR experimental data means we can classify the samples under investigation. As already mentioned, the $1/r$ regime is the fingerprint of a good conductive layer, while a very large signal in the center characterizes poor-conductor thin films on good-conductor substrates. Figure 27 concerns a 300 nm YBaCuO layer deposited on a LaAlO_{3} substrate.^{57} Considering the amplitude experimental data, it turns out that the YBaCuO layer is a poor conductor with respect to the substrate.^{55–57} A sophisticated best fit allows us to evidence the thermal properties of the layer and the thermal resistance between the layer and the substrate.^{57} It turns out that the thermal contrast is about 5 (YBaCuO *k *= 2 W/mK– LaAL_{2}O_{3} *k *= 10 W/mK).

To conclude this IOC section, we emphasize that the shape of the experimental amplitude data is a good indicator of thermal contrast. A strong value in the center reveals that the thermal conductivity of the layer is smaller than the thermal conductivity of the substrate, even with a relatively weak thermal contrast (larger than around 2). To evaluate the thermal diffusivity of the substrate, the slopes of the log (Ampl) and the phase data are investigated. To obtain the thermal properties of the layer, a careful investigation of the central part is needed. The beam diameter is an important parameter of the best fit. Moreover, best fits have to be done to extract carefully all parameters with Eq. (15). It is straightforward to see the role of a thermal resistance R in Eq. (19) (which is shown again here for convenience),

### C. Recipe for anisotropic materials

We made it clear in Ref. 30 that when uniaxial anisotropic samples are under study (*k _{r} D_{r}* radial and

*k*perpendicular properties), isotropic heat diffusion equations can be used with effective values for thermal conductivity and thermal diffusivity (Fig. 28, top). When a semi-infinite isotropic sample is covered with an anisotropic layer, an additional effective thickness must be introduced (Fig. 28, bottom).

_{z}D_{z}### D. Depth profiling: Stacks

Here, we discuss a unique and important advantage of MTR: depth profiling. By varying the modulation frequency, since the thermal diffusion length changes, heat can explore the sample in depth. In the case of a multilayered structure, the upper layers are explored at high frequencies and the deep layers make a significant contribution at lower frequencies. It is even possible to explore layer by layer when low thermal conductors are present.^{58,59} Attention must be paid at low frequencies, since the layer/substrate interface has to be considered.

To illustrate this depth profiling, a solid-state electrolyte, Lipon, used in microbatteries, is presented.^{53,58} Experimental work was done on simplified samples very close to the envisioned device. Figure 29 illustrates this interesting topic: on the left part, the unknown layer is sandwiched between a thin conducting gold layer and a substrate whose thermal properties are well known. In the case we present, Lipon is a poor thermal conductor whose thickness can be considered as infinite for high-frequency investigation. Its effusivity is then deduced for high-frequency investigation (see Sec. V A on COI). For the low-frequency best fit, we impose a constraint (effusivity) between the thermal diffusivity and the thermal conductivity since the product density specific heat is not known. In regard to low frequency experiments, we consider four layers: 200 nm Gold, 20 nm Ti, 2 *μ*m Lipon, 250 nm TI, and a silicon substrate. Note that Eq. (10) is valid for only one layer deposited on a substrate and has to be generalized when several layers are superimposed. Reference 60 can help with the calculus. Figure 30 summarizes the results obtained from the microbattery.

## VI. CONCLUSIONS

In this article, we have shown that the modulated thermoreflectance (MTR) technique is the tool of choice for local thermal property evaluation, since it provides a reliable experimental and theoretical arsenal for measurements, data fitting, and numerical simulation. A judicious choice of the experimental parameters enables a measurement of the thermal properties of bulk materials and layered thin structures. The MTR technique is also well adapted to inhomogeneous bulk materials, such as ceramics, permitting the measurement of thermal properties at different spatial scales (for example, grains and grain boundaries).^{53} Two configurations of a highly thermally contrasted layered structure were analyzed. The former is a layer of a good conductor deposited on an insulator substrate. In this case, analytical expressions allow a straightforward evaluation of the thermal properties of the layer. The opposite situation is more challenging. When a thin insulating layer is present, its thermal conductivity makes a significant contribution to the amplitude/real part of the experiments. In most cases, thermal properties can be evaluated by a numerical fit. Finally, we discussed the possibility of depth profiling the thermal properties of stacked structures by MTR.

## ACKNOWLEDGMENTS

This research has received funding from the French National Research Agency under the project HiPerTherMag (No. ANR-18-CE05-0019). The authors gratefully acknowledge the efficient work of many colleagues and students who collaborated in achieving the results presented in this tutorial. LPEM ESPCI/PSL: BinCheng Li, Lionel Pottier, Jean Paul Roger, Karsten Plamann, Benoit C Forget, Claude Boccara, Catherine Pelissonnier-Grosjean, Dominique Jeulin, Alain Thorel, K. Watari, L. Fabbri, Valerie Reita, Tsetuo Ikari, Isabelle Barbereau, Antonio Mansanares, and Lionel Aigouy. INSP UPMC/Sorbonne Université: Mahmoud Eddrief, Paola Atkinson, Jean-Yves Duquesne, Feng Xu, Laurent Belliard, Loïc Beccera, Serge Vincent, Agnes Huynh, Maxime Vabre, and Eric Charron. The authors also want to assure the whole ICPPP community of their recognition; without this, these works would certainly not exist. Many thanks to G. Busse, D. Cahen, G. Diebold, C. Glorieux, V. Gusev, R. Li Voti, A. Mandelis, J. Pelzl, A. Rosencwaig, A. Salazar, J. Thoen, H. Vargas, U. Zamit, etc.

## DATA AVAILABILITY

Data sharing is not applicable to this article as no new data were created or analyzed in this study.