simulations.heat
#
Classes:
|
Heat simulations including laser excitation and N-temperature model. |
- class udkm1Dsim.simulations.heat.Heat(S, force_recalc, **kwargs)[source]#
Bases:
Simulation
Heat simulations including laser excitation and N-temperature model.
- Parameters:
S (Structure) – sample to do simulations with.
force_recalc (boolean) – force recalculation of results.
- Keyword Arguments:
save_data (boolean) – true to save simulation results.
cache_dir (str) – path to cached data.
disp_messages (boolean) – true to display messages from within the simulations.
progress_bar (boolean) – enable tqdm progress bar.
heat_diffusion (boolean) – true when including heat diffusion in the calculations.
intp_at_interface (int) – number of additional spacial points at the interface of each layer.
backend (str) – pde solver backend - either default scipy or matlab.
- Attributes:
S (Structure) – sample structure to calculate simulations on.
force_recalc (boolean) – force recalculation of results.
save_data (boolean) – true to save simulation results.
cache_dir (str) – path to cached data.
disp_messages (boolean) – true to display messages from within the simulations.
progress_bar (boolean) – enable tqdm progress bar.
heat_diffusion (boolean) – true when including heat diffusion in the calculations.
intp_at_interface (int) – number of additional spacial points at the interface of each layer.
backend (str) – pde solver backend - either default scipy or matlab.
excitation (dict{ndarray[float, Quantity]}) – excitation parameters fluence, delay_pump, pulse_width, wavelength, theta, polarization, multilayer_absorption, backside
boundary_conditions (dict{str, float, Quantity}) – boundary conditions of the top and bottom boundary for the heat diffusion calculation.
top_type
orbottom_type
must be one ofboundary_types
. For the last two types the corresponding value,top_value
andbottom_value
have to be set as \(K \times 1\) array, where \(K\) is the number of sub-systems.boundary_types (list[str]) – description of boundary types.
isolator
temperature
flux
distances (ndarray[float, Quantity]) – spatial grid for heat diffusion [m]
ode_options (dict) – options for scipy solve_ivp ode solver
ode_options_matlab (dict) – dict with options for the MATLAB pdepe solver.
matlab_engine (module) – MATLAB to Python API engine required for calculating heat diffusion.
Methods:
get_hash
(delays, init_temp, **kwargs)Returns a unique hash given by the
delays
andinit_temp
as well as the sample structure hash for relevant thermal parameters.check_initial_temperature
(init_temp[, distances])Initial temperature for heat simulations can be either a scalar temperature which is assumed to be valid for all layers in the structure or a temperature profile is given with one temperature for each layer in the structure and for each subsystem.
check_excitation
(delays)The optical excitation is a dictionary with fluence \(F\) [J/m²], delays \(t\) [s] of the pump events, and pulse width \(\tau\) [s].
get_absorption_profile
([distances, backside])Returns the differential absorption profile \(\mbox{d}A/\mbox{d}z\).
The transmission is given by:
Calculates the intensity, differential absorption and temperature increase profiles in each layer of a multilayers structure for \(p\) -polarized light.
Calculate the final temperature and temperature change for each layer of the sample structure after an optical excitation with a fluence \(F\) [J/m^2] and an initial temperature \(T_1\) [K]:
get_temp_map
(delays, init_temp)Returns a temperature profile for the sample structure after optical excitation.
calc_temp_map
(delays, input_init_temp)Calculates the temperature profile for the sample structure after optical excitation.
calc_heat_diffusion
(init_temp, distances, ...)Calculates a temperature profile including heat diffusion for a given delay array and initial temperature profile.
odefunc
(t, u, N, K, d_x_grid, x, ...)Ordinary differential equation that is solved for 1D heat diffusion.
conv_with_function
(y, x, handle)Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).
disp_message
(message)Wrapper to display messages for that class.
save
(full_filename, data, *args)Save data to file.
- get_hash(delays, init_temp, **kwargs)[source]#
Returns a unique hash given by the
delays
andinit_temp
as well as the sample structure hash for relevant thermal parameters.- Parameters:
delays (ndarray[float]) – delay grid for the simulation.
init_temp (ndarray[float]) – initial spatial temperature profile.
**kwargs (float, optional) – optional parameters.
- Returns:
hash (str) – unique hash.
- check_initial_temperature(init_temp, distances=[])[source]#
Initial temperature for heat simulations can be either a scalar temperature which is assumed to be valid for all layers in the structure or a temperature profile is given with one temperature for each layer in the structure and for each subsystem. Alternatively a spatial grid can be provided.
- Parameters:
init_temp (float, Quantity, ndarray[float, Quantity]) – initial temperature scalar or array [K].
distances (ndarray[float, Quantity], optional) – spatial grid of the initial temperatrue.
- Returns:
init_temp (ndarray[float]) – checked initial temperature as array on the according spatial grid.
- check_excitation(delays)[source]#
The optical excitation is a dictionary with fluence \(F\) [J/m²], delays \(t\) [s] of the pump events, and pulse width \(\tau\) [s]. \(N\) is the number of pump events.
Traverse excitation vector to update the
delay_pump
\(t_p\) vector for finite pulse durations \(w(i)\) as follows\[[t_p(i)-\mbox{window}\cdot w(i)]:[t_p(i)+\mbox{window}\cdot w(i)]: [w(i)/\mbox{intp}]\]and to combine excitations which have overlapping intervals.
Moreover the incidence angle \(\vartheta\) is taken into account for the user-defined incidence fluence in order to project the laser footprint onto the sample surface.
- Parameters:
delays (ndarray[Quantity]) – delays range of simulation [s].
- Returns:
(tuple) –
res (list[ndarray[float]]) - resulting list of excitations with interpolated delay density around excitations.
fluence (ndarray[float]) - projected excitation fluences.
delay_pump (ndarray[float]) - delays of the excitations.
pulse_width (ndarray[float]) - pulse widths of the excitations.
- get_absorption_profile(distances=[], backside=False)[source]#
Returns the differential absorption profile \(\mbox{d}A/\mbox{d}z\).
- Parameters:
distances (ndarray[float], optional) – spatial grid for calculation.
backside (boolean, optional) – backside or frontside excitation.
- Returns:
dAdz (ndarray[float]) – differential absorption within each layer calculated either by Lambert-Beers law or by a multilayer absorption formalism.
- get_Lambert_Beer_absorption_profile(distances=[], backside=False)[source]#
The transmission is given by:
\[\tau = \frac{I}{I_0} = \exp(-z/ \zeta)\]and the absorption by:
\[A = 1 - \tau = 1 - \exp(-z/\zeta)\]The absorption profile can be derived from the spatial derivative:
\[\frac{\mbox{d}A(z)}{\mbox{d}z} = \frac{1}{\zeta} \exp(-z/\zeta)\]- Parameters:
distances (ndarray[float], optional) – spatial grid for calculation.
backside (boolean, optional) – backside or frontside excitation.
- Returns:
dAdz (ndarray[float]) – differential absorption within each layer calculated by Lambert-Beers law.
- get_multilayers_absorption_profile(distances=[], backside=False)[source]#
Calculates the intensity, differential absorption and temperature increase profiles in each layer of a multilayers structure for \(p\) -polarized light.
Calculation based on the method in Ref [5] and code developed Matlab by L. Le Guyader, see Ref [6].
Copyright (2012-2014) Loïc Le Guyader <loic.le_guyader@helmholtz-berlin.de>
- Parameters:
distances (ndarray[float], optional) – spatial grid for calculation.
backside (boolean, optional) – backside or frontside excitation.
- Returns:
(tuple) –
dAdz (ndarray[float]) - differential absorption within each layer.
Ints (ndarray[float]) - intensity profiles within each layer.
R_total (float) - total amount of reflection from the multilayer.
T_total (float) - total transmission in the last layer of the multilayer.
References:
- get_temperature_after_delta_excitation(fluence, init_temp, distances=[])[source]#
Calculate the final temperature and temperature change for each layer of the sample structure after an optical excitation with a fluence \(F\) [J/m^2] and an initial temperature \(T_1\) [K]:
\[\Delta E = \int_{T_1}^{T_2} m \, c(T)\, \mbox{d}T\]where \(\Delta E\) is the absorbed energy in each layer and \(c(T)\) is the temperature-dependent heat capacity [J/kg K] and \(m\) is the mass [kg].
The absorbed energy per layer can be linearized from the absorption profile \(\mbox{d} \alpha / \mbox{d}z\) as
\[\Delta E = \frac{\mbox{d} \alpha}{\mbox{d}z} E_0 \Delta z\]where \(E_0\) is the initial energy impinging on the first layer given by the fluence \(F = E / A\). \(\Delta z\) is equal to the thickness of each layer.
Finally, one has to minimize the following modulus to obtain the final temperature \(T_2\) of each layer:
\[\left| \int_{T_1}^{T_2} m c(T)\, \mbox{d}T - \frac{\mbox{d}\alpha} {\mbox{d}z} E_0 \Delta z \right| \stackrel{!}{=} 0\]- Parameters:
fluence (float, Quantity) – incident fluence [J/m²].
init_temp (float, Quantity, ndarray[float, Quantity]) – initial temperature scalar or array [K].
distances (ndarray[float], optional) – spatial grid for calculation.
- Returns:
(tuple) –
final_temp (ndarray[float]) - final temperature after delta excitation.
delta_T (ndarray[float]) - temperature change.
- get_temp_map(delays, init_temp)[source]#
Returns a temperature profile for the sample structure after optical excitation. The result can be saved using an unique hash of the sample and the simulation parameters in order to reuse it.
- Parameters:
delays (ndarray[Quantity]) – delays range of simulation [s].
init_temp (float, Quantity, ndarray[float, Quantity]) – initial temperature scalar or array [K].
- Returns:
(tuple) –
temp_map (ndarray[float]) - spatio-temporal temperature map.
delta_temp_map (ndarray[float]) - spatio-temporal temperature change map.
- calc_temp_map(delays, input_init_temp)[source]#
Calculates the temperature profile for the sample structure after optical excitation. Heat diffusion can be included if
heat_diffusion = true
.- Parameters:
delays (ndarray[Quantity]) – delays range of simulation [s].
input_init_temp (float, Quantity, ndarray[float, Quantity]) – initial temperature scalar or array [K].
- Returns:
(tuple) –
temp_map (ndarray[float]) - spatio-temporal temperature map.
delta_temp_map (ndarray[float]) - spatio-temporal temperature change map.
checked_excitations (list[ndarray[float]]) - resulting list of checked excitations.
- calc_heat_diffusion(init_temp, distances, delays, delay_pump, pulse_width, fluence)[source]#
Calculates a temperature profile including heat diffusion for a given delay array and initial temperature profile. Here we have to solve the 1D heat diffusion equation for \(N\) subsystems:
\[\begin{split}c_1(T_1) \rho \frac{\partial T_1}{\partial t} & = & \frac{\partial}{\partial z} \left( k_1(T_1) \frac{\partial T_1}{\partial z} \right) + G_1(T_1,...,T_N) + S(z,t) \\ & \vdots & \\ c_N(T_N) \rho \frac{\partial T_N}{\partial t} & = & \frac{\partial}{\partial z} \left( k_N(T_N) \frac{\partial T_N}{\partial z} \right) + G_N(T_1,...,T_N)\end{split}\]where \(T_i\) is the temperature [K], \(z\) the distance [m], \(t\) the delay [s], \(c_i(T)\) the temperature dependent heat capacity [J/kg K], \(\rho\) the density [kg/m³] and \(k_i(T)\) is the temperature-dependent thermal conductivity [W/m K] and \(S(z,t)\) is a source term [W/m³]. The energy flow between the subsystems is given by the
sub_system_coupling
parameter \(G_i(T_1,...,T_N)\) of the individual layers. The index \(i\) refers to the \(i\)-th subsystem.The 1D heat diffusion equation can be either solved with SciPy or Matlab as backend.
- Parameters:
init_temp (float, Quantity, ndarray[float, Quantity]) – initial temperature scalar or array [K].
distances (ndarray[float]) – spatial grid for calculation.
delays (ndarray[Quantity]) – delays range of simulation [s].
delay_pump (ndarray[float]) – delays of the excitations.
pulse_width (ndarray[float]) – pulse widths of the excitations.
fluence (ndarray[float]) – excitation fluences.
- Returns:
temp_map (ndarray[float]) – spatio-temporal temperature map.
- static odefunc(t, u, N, K, d_x_grid, x, thermal_conds, heat_capacities, sub_system_coupling, densities, indices, dAdz, fluence, delay_pump, pulse_length, bc_top_type, bc_top_value, bc_bottom_type, bc_bottom_value, pbar, state)[source]#
Ordinary differential equation that is solved for 1D heat diffusion.
- Parameters:
t (ndarray[float]) – internal time steps of the ode solver.
u (ndarray[float]) – internal variable of the ode solver.
N (int) – number of spatial grid points.
K (int) – number of sub-systems.
d_x_grid (ndarray[float]) – derivative of spatial grid.
x (ndarray[float]) – start point of actual layers.
thermal_conds (ndarray[@lambda]) – T-dependent thermal conductivity function handles.
heat_capacities (ndarray[@lambda]) – T-dependent heat capacity function handles.
sub_system_coupling (ndarray[@lambda]) – T-dependent sub-system coupling.
densities (ndarray[float]) – density of layers.
indices (ndarray[int]) – indices of actual layers in respect to interpolated spatial grid.
dAdz (ndarray[float]) – differential absorption profile.
fluence (ndarray[float]) – excitation fluences.
delay_pump (ndarray[float]) – delay of excitations.
pulse_length (ndarray[float]) – pulse widths of excitations.
bc_top_type (int) – top boundary type.
bc_top_value (ndarray[float]) – top boundary value.
bc_bottom_type (int) – bottom boundary type.
bc_bottom_value (ndarray[float]) – bottom boundary value.
pbar (tqdm) – tqdm progressbar.
state (list[float]) – state variables for progress bar.
- Returns:
dudt (ndarray[float]) – temporal derivative of internal variable.
- static conv_with_function(y, x, handle)#
Convolutes the array \(y(x)\) with a function given by the handle on the argument array \(x\).
- Parameters:
y (ndarray[float]) – y data.
x (ndarray[float]) – x data.
handle (@lamdba) – convolution function.
- Returns:
y_conv (ndarray[float]) – convoluted data.
- disp_message(message)#
Wrapper to display messages for that class.
- Parameters:
message (str) – message to display.
- save(full_filename, data, *args)#
Save data to file. The variable name can be handed as variable argument.
- Parameters:
full_filename (str) – full file name to data file.
data (ndarray) – actual data to save.
*args (str, optional) – variable name within the data file.