NonEquilibriumHeatExchange¶

class NonEquilibriumHeatExchange(configuration, heating_power, heat_source, heat_sink, call_interval=None, update_profile_interval=None, profile_resolution=None, exchange_interval=1)¶

A class that implements a heat flow by constant heating power exchange technique via a hook function.

Parameters:

configuration (BulkConfiguration) – The initial configuration on which the heat flow simulation will be performed.
heating_power (PhysicalQuantity of type energy/time) – The thermal energy that is transferred from the cold to the hot region per time.
heat_source (str | list of ints | tuple of type (float, PhysicalQuantity of type length)) – The heat source atoms can be defined by giving a tag that specifies which atoms are in the heat source, a list of atomic indices, or by defining a spatial region. The spatial region is defined as 2-tuple of the origin and length. The origin is the fractional position along the c-axis and the length is the thickness along the c-axis.
heat_sink (str | list of ints | tuple of type (float, PhysicalQuantity of type length)) – The heat sink atoms can be defined by giving a tag that specifies which atoms are in the heat sink, a list of atomic indices, or by defining a spatial region. The spatial region is defined as 2-tuple of the origin and length. The origin is the fractional position along the c-axis and the length is the thickness along the c-axis.
call_interval (int) – The interval to perform energy exchange by rescaling the velocities.
Default: 1 (every step).
update_profile_interval (int) – [description].
profile_resolution (PhysicalQuantity of type length) – The spatial resolution for the on-the-fly calculation of the temperature profile.
Default: 2.0*Angstrom.
exchange_interval (int) – The interval to perform energy exchange by rescaling the velocities.
Default: 1 (every step).

callInterval()¶

Returns:: The call interval of this hook function.
Return type:: int

dz()¶

Returns:: Float version of the profile resolution.
Return type:: float

exchangeInterval()¶

Returns:: The hot and cold region velocity exchange interval.
Return type:: int

heatCurrent()¶

Returns:: The average heat current.
Return type:: PhysicalQuantity of type energy/time

heatSink()¶

Returns:: The heat sink.
Return type:: SpatialHeatRegion

heatSource()¶

Returns:: The head source.
Return type:: SpatialHeatRegion

masses()¶

Returns:: The masses for the entire system and for the two regions.
Return type:: PhysicalQuantity of type mass

nlprint(stream=<_io.TextIOWrapper name='<stdout>' mode='w' encoding='utf-8'>)¶

Print a formatted string with the average temperature profile, as well as the average heat current.

Parameters:: stream (A stream that supports strings being written to using 'write'.) – The stream the temperature profile is written to.

nslices()¶

Returns:: The number of bins (slices) to store the temperature profile.
Return type:: int

profileResolution()¶

Returns:: The temperature profile resolution.
Return type:: PhysicalQuantity of type length

temperatureProfile()¶

Returns:: The temperature profile as a tuple containing bin-centers and corresponding temperature values.
Return type:: (PhysicalQuantity of type length, PhysicalQuantity of type temperature)

uniqueString()¶: Return a unique string representing the state of the object.

updateProfileInterval()¶

Returns:: The temperature profile update interval.
Return type:: int

Usage Example¶

Perform a non-equilibrium molecular dynamics simulation with a constant heating power on a silicon crystal with a single germanium layer, to obtain the thermal conductance through this impurity layer.

potentialSet = Tersoff_SiGe_1989()
calculator = TremoloXCalculator(parameters=potentialSet)
bulk_configuration.setCalculator(calculator)

# -------------------------------------------------------------
# Molecular Dynamics
# -------------------------------------------------------------

initial_velocity = MaxwellBoltzmannDistribution(
    temperature=300.0*Kelvin,
    remove_center_of_mass_momentum=True
)

method = NVEVelocityVerlet(
    time_step=1*femtoSecond,
    initial_velocity=initial_velocity,
)

# Use a heating power of 50 nW.
heating_power = 50.0e-9*Joule/Second
momentum_exchange_hook = NonEquilibriumHeatExchange(
    configuration=bulk_configuration,
    heating_power=heating_power,
    heat_source='heat_source',
    heat_sink='heat_sink',
    update_profile_interval=10,
    profile_resolution=2.0*Ang
)

md_trajectory = MolecularDynamics(
    bulk_configuration,
    constraints=[],
    trajectory_filename='Si_w_Ge_layer_NEMD.hdf5',
    steps=100000,
    log_interval=500,
    post_step_hook=momentum_exchange_hook,
    method=method
)

Here, only the molecular dynamics block of the script is shown. The full script can be found found in the file nonequilibriumheatexchange.py.

Notes¶

The NonEquilibriumHeatExchange method can be used to run MolecularDynamics simulations with a constant heat flux between the two selected reservoirs. The technique can therefore be used to simulate the thermal conductance of a given configuration (see e.g, [1]).
The reverse non-equilibrium heat exchange method is implemented as a class, which can be used as a post_step_hook in MolecularDynamics.
The velocities in the two reservoir regions are rescaled in such a way that a well-defined amount kinetic energy is added to the atoms in the heat_source region whereas the same amount is taken away from the kinetic energy in the heat_sink region. The magnitude of the kinetic energy that is transferred per time from the heat_sink to the heat_source can be specified by the heating_power parameter. By default the rescaling is performed at every MD step, if desired, the interval can be given by the user via the parameter exchange_interval.
If a non-zero update_profile_interval is specified, a temperature profile is stored on-the-fly during the simulation. After the MD run, the profiles may be accessed using the method temperatureProfile.
The thermal bulk conductivity can be obtained by plotting the temperature profile of the system using the MD-Analyzer tool and fitting its gradient \(\nabla T\) along the transport direction. The thermal conductivity is then obtained via Fourier’s law

\[\kappa = - \frac{\dot{Q}}{A \nabla T}\]

where \(\dot{Q}\) is the specified heating_power and A is the cross-sectional area perpendicular to the transport direction.

The thermal conductance across an interface G (Kapitza conductance) can be obtained as

\[G = -\frac{\dot{Q}}{\Delta T}\]

where \(\Delta T\) is the temperature drop across the interface.