This page was generated from
docs/examples/enthalpy_sorption.ipynb.
To start an interactive version:
Sorption enthalpy calculations¶
Sorption enthalpy, \(\Delta H_{ads}\), is an indication of the strength of the adsorbate-material interaction and can be estimated through several methods.
Isosteric enthalpy (Clausius-Clapeyron)¶
In order to calculate \(\Delta H_{ads}\), at least two isotherms which were taken at different temperatures are required.
First, make sure the data is imported.
[1]:
# import isotherms
%run import.ipynb
%matplotlib inline
# import the characterisation module
import pygaps.characterisation as pgc
Selected 5 isotherms with nitrogen at 77K
Selected 2 room temperature calorimetry isotherms
Selected 2 isotherms for IAST calculation
Selected 3 isotherms for isosteric enthalpy calculation
The file version is None while the parser uses version d546195. Strange things might happen, so double check your data.
Let's quickly plot the isotherms to see how they look. We put the temperature of the experiment in the legend by using the lgd_keys keyword.
[2]:
# import the graphing module
import pygaps.graphing as pgg
pgg.plot_iso(
isotherms_isosteric,
lgd_keys=['adsorbate', 'temperature'],
)
[2]:
<Axes: xlabel='Pressure [$bar$]', ylabel='Loading [$mmol\\/g^{-1}$]'>
The isotherms look good, except perhaps a bit of measurement error in the low pressure region.
The isosteric enthalpy calculation takes the list of the isotherms and returns the results as a dictionary. Using the verbose keyword, we also generate a graph that includes an error bar.
[3]:
result_dict = pgc.isosteric_enthalpy(isotherms_isosteric, verbose=True)
The inaccuracy in the low pressure region has contributed to the odd enthalpy curve. One other option would be to first fit a model to each isotherm, then use it for the enthalpy determination.
Let's try a Double Site Langmuir model and then re-run the isosteric calculation.
[4]:
from pygaps.modelling import model_iso
models_isosteric = [model_iso(iso, model="dslangmuir") for iso in isotherms_isosteric]
result_dict = pgc.isosteric_enthalpy(models_isosteric, verbose=True)
More information about the functions and their use can be found in the manual.
Whittaker method¶
The Whittaker method uses a Tóth’s modification of Polanyi's potential theory to determine enthalpy of adsorption from a single isotherm fitted with a suitable model. We loaded the example isotherm before, so we can now plot it.
[5]:
pgg.plot_iso(isotherms_enth_whittaker, lgd_keys=['adsorbate', 'temperature'], pressure_unit='bar')
[5]:
<Axes: xlabel='Pressure [$bar$]', ylabel='Loading [$mmol\\/g^{-1}$]'>
The code can accept a ModelIsotherm of the required type. Note that this isotherm must be in correct pressure mode (absolute) and units: (Pa).
[6]:
import pygaps.modelling as pgm
for iso in isotherms_enth_whittaker:
iso.convert_pressure(mode_to="absolute", unit_to="Pa")
models_whittaker = pgm.model_iso(
isotherms_enth_whittaker[0],
model=pgm._WHITTAKER_MODELS,
verbose=True
)
result_dict_whitt = pgc.enthalpy_sorption_whittaker(models_whittaker, verbose=True)
Attempting to model using Langmuir.
Model Langmuir success, RMSE is 0.00748
Attempting to model using DSLangmuir.
Model DSLangmuir success, RMSE is 0.00368
Attempting to model using TSLangmuir.
Model TSLangmuir success, RMSE is 0.027
Attempting to model using Toth.
Model Toth success, RMSE is 0.00493
Attempting to model using DSToth.
Model DSToth success, RMSE is 0.0167
Attempting to model using ChemiPhysisorption.
Model ChemiPhysisorption success, RMSE is 0.205
Best model fit is DSLangmuir.
Alternatively, the PointIsotherm and the desired model can be passed as parameters, and fitting is performed automatically before the method is applied.
[7]:
result_dict_whitt = pgc.enthalpy_sorption_whittaker(
isotherms_enth_whittaker[0],
model="guess",
verbose=True
)
Attempting to model using Langmuir.
Model Langmuir success, RMSE is 0.00748
Attempting to model using DSLangmuir.
Model DSLangmuir success, RMSE is 0.00368
Attempting to model using TSLangmuir.
Model TSLangmuir success, RMSE is 0.027
Attempting to model using Toth.
Model Toth success, RMSE is 0.00493
Attempting to model using DSToth.
Model DSToth success, RMSE is 0.0167
Attempting to model using ChemiPhysisorption.
Model ChemiPhysisorption success, RMSE is 0.205
Best model fit is DSLangmuir.
Comparing Whittaker and Clausius-Clapeyron¶
How do these methods compare? We can use the Whittaker method on each of the isotherms in isotherms_isosteric, and compare with the Clausius-Clapeyron results in result_dict.
[42]:
import pygaps.graphing as pgg
import matplotlib.pyplot as plt
_, ax = plt.subplots(1, 1) # set up an ax to plot on
for i in isotherms_isosteric:
# Calculate heat of adsorpation with whittaker method
whittaker = pgc.enthalpy_sorption_whittaker(
i,
model='dslangmuir', # using DSLangmuir as previously showed good fit
verbose=False
)
print(
f"{i.temperature}: {whittaker['model_isotherm'].model.rmse}"
) # check that fit is reasonable for all models
ax.errorbar( # and plotting
whittaker['loading'],
whittaker['enthalpy_sorption'],
yerr=whittaker['std_errs'], # with error bars
label=str(i.temperature),
)
# Now plotting Clausis-Clapeyron enthalpies
pgg.calc_graphs.isosteric_enthalpy_plot(
result_dict['loading'],
result_dict['isosteric_enthalpy'],
result_dict['std_errs'],
units=isotherms_isosteric[0].units,
ax=ax
)
ax.legend(frameon=False)
298.15: 0.009565136818274825
323.15: 0.005256038450545493
348.15: 0.009516772628472047
[42]:
<matplotlib.legend.Legend at 0x7fb6e7a35cd0>
A reasonable match in terms of magnitude of \(\Delta H_{st}\), but not in terms of the overall trend. The Clausius-Clapeyron method predicts increasing enthalpy with loading up to 6 mmol g-1, while the Whittaker method predicts a continuous decrease. Whittaker approximation will always give this decreasing trend.
What happens when we use a different model, say Toth?
[46]:
import pygaps.graphing as pgg
import matplotlib.pyplot as plt
_, ax = plt.subplots(1, 1) # set up an ax to plot on
for i in isotherms_isosteric:
# Calculate heat of adsorpation with whittaker method
whittaker = pgc.enthalpy_sorption_whittaker(
i,
model='toth', # using DSToth
verbose=False
)
print(
f"{i.temperature}: {whittaker['model_isotherm'].model.rmse}"
) # check that fit is reasonable for all models
ax.errorbar( # and plotting
whittaker['loading'],
whittaker['enthalpy_sorption'],
yerr=whittaker['std_errs'], # with error bars
label=str(i.temperature),
)
# Now plotting Clausis-Clapeyron enthalpies
pgg.calc_graphs.isosteric_enthalpy_plot(
result_dict['loading'],
result_dict['isosteric_enthalpy'],
result_dict['std_errs'],
units=isotherms_isosteric[0].units,
ax=ax
)
ax.legend(frameon=False)
---------------------------------------------------------------------------
CalculationError Traceback (most recent call last)
Cell In[46], line 8
4 _, ax = plt.subplots(1, 1) # set up an ax to plot on
6 for i in isotherms_isosteric:
7 # Calculate heat of adsorpation with whittaker method
----> 8 whittaker = pgc.enthalpy_sorption_whittaker(
9 i,
10 model='toth', # using DSToth
11 verbose=False
12 )
13 print(
14 f"{i.temperature}: {whittaker['model_isotherm'].model.rmse}"
15 ) # check that fit is reasonable for all models
16 ax.errorbar( # and plotting
17 whittaker['loading'],
18 whittaker['enthalpy_sorption'],
19 yerr=whittaker['std_errs'], # with error bars
20 label=str(i.temperature),
21 )
File ~/anaconda3/lib/python3.9/site-packages/pygaps/characterisation/enth_sorp_whittaker.py:129, in enthalpy_sorption_whittaker(isotherm, model, loading, verbose, dographs, **kwargs)
126 model = pgm._WHITTAKER_MODELS
128 max_nfev = kwargs.get('max_nfev', None)
--> 129 isotherm = pgm.model_iso(
130 isotherm,
131 branch='ads',
132 model=model,
133 verbose=verbose,
134 optimization_params=dict(max_nfev=max_nfev)
135 )
137 model = isotherm.model.name
139 if not pgm.is_model_whittaker(model):
File ~/anaconda3/lib/python3.9/site-packages/pygaps/modelling/__init__.py:287, in model_iso(isotherm, branch, model, param_guess, param_bounds, optimization_params, verbose)
261 """
262 Fits a PointIsotherm with a model.
263
(...)
284 Prints out extra information about steps taken.
285 """
286 from pygaps.core.modelisotherm import ModelIsotherm
--> 287 return ModelIsotherm.from_pointisotherm(
288 isotherm,
289 branch=branch,
290 model=model,
291 param_guess=param_guess,
292 param_bounds=param_bounds,
293 optimization_params=optimization_params,
294 verbose=verbose,
295 )
File ~/anaconda3/lib/python3.9/site-packages/pygaps/core/modelisotherm.py:343, in ModelIsotherm.from_pointisotherm(cls, isotherm, branch, model, param_guess, param_bounds, optimization_params, verbose)
341 if isinstance(model, str):
342 if model != 'guess':
--> 343 return cls(
344 branch=branch,
345 model=model,
346 param_guess=param_guess,
347 param_bounds=param_bounds,
348 optimization_params=optimization_params,
349 verbose=verbose,
350 **iso_params
351 )
353 return ModelIsotherm.guess(
354 branch=branch,
355 models=model,
(...)
358 **iso_params
359 )
File ~/anaconda3/lib/python3.9/site-packages/pygaps/core/modelisotherm.py:212, in ModelIsotherm.__init__(self, pressure, loading, isotherm_data, pressure_key, loading_key, branch, model, param_guess, param_bounds, optimization_params, verbose, **other_properties)
209 param_guess = self.model.initial_guess(pressure, loading)
211 # fit model to isotherm data
--> 212 self.model.fit(
213 pressure,
214 loading,
215 param_guess,
216 optimization_params,
217 verbose,
218 )
220 # Plot fit if verbose
221 if verbose and other_properties.pop('plot_fit', True):
File ~/anaconda3/lib/python3.9/site-packages/pygaps/modelling/base_model.py:259, in IsothermBaseModel.fit(self, pressure, loading, param_guess, optimization_params, verbose)
256 if optimization_params:
257 fit_args.update(optimization_params)
--> 259 opt_res = self.fit_leastsq(fit_args)
261 # assign params
262 for i, _ in enumerate(param_names):
File ~/anaconda3/lib/python3.9/site-packages/pygaps/modelling/base_model.py:280, in IsothermBaseModel.fit_leastsq(self, leastsq_args)
276 raise CalculationError(
277 f"Fitting routine for {self.name} failed with error:\n\t{err}"
278 ) from err
279 if not opt_res.success:
--> 280 raise CalculationError(
281 f"Fitting routine for {self.name} failed with error:"
282 f"\n\t{opt_res.message}"
283 f"\nTry a different starting point in the nonlinear optimization"
284 f"\nby passing a dictionary of parameter guesses, param_guess, to the constructor."
285 f"\nDefault starting guess for parameters:"
286 f"\n{leastsq_args['x0']}\n"
287 )
288 return opt_res
CalculationError: Fitting routine for Toth failed with error:
The maximum number of function evaluations is exceeded.
Try a different starting point in the nonlinear optimization
by passing a dictionary of parameter guesses, param_guess, to the constructor.
Default starting guess for parameters:
[1.78673000e+01 4.93900934e-05 1.00000000e+00]
This gives an error as pgm.model_iso() can't fit the model. We can try to increase the number of function evaluatins max_nfev to see if that helps.
[51]:
import pygaps.graphing as pgg
import matplotlib.pyplot as plt
_, ax = plt.subplots(1, 1) # set up an ax to plot on
for i in isotherms_isosteric:
# Calculate heat of adsorpation with whittaker method
whittaker = pgc.enthalpy_sorption_whittaker(
i,
model='toth', # using Toth
verbose=False,
**{'max_nfev': 1e5}, # Increasing number of function evaluations
)
print(
f"{i.temperature}: {whittaker['model_isotherm'].model.rmse}"
) # check that fit is reasonable for all models
ax.errorbar( # and plotting
whittaker['loading'],
whittaker['enthalpy_sorption'],
yerr=whittaker['std_errs'], # with error bars
label=str(i.temperature),
)
# Now plotting Clausis-Clapeyron enthalpies
pgg.calc_graphs.isosteric_enthalpy_plot(
result_dict['loading'],
result_dict['isosteric_enthalpy'],
result_dict['std_errs'],
units=isotherms_isosteric[0].units,
ax=ax
)
ax.legend(frameon=False)
298.15: 0.020847257578099514
323.15: 0.01520385363705269
348.15: 0.008170243232597895
[51]:
<matplotlib.legend.Legend at 0x7fb6e7d25190>
