.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/selection/GCH-ROY.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_examples_selection_GCH-ROY.py: Generalized Convex Hull construction for the polymorphs of ROY ============================================================== This notebook analyzes the structures of 264 polymorphs of ROY, from `Beran et Al, Chemical Science (2022) `__, comparing the conventional density-energy convex hull with a Generalized Convex Hull (GCH) analysis (see `Anelli et al., Phys. Rev. Materials (2018) `__). .. GENERATED FROM PYTHON SOURCE LINES 13-23 .. code-block:: Python import matplotlib.tri as mtri import numpy as np from matplotlib import pyplot as plt from sklearn.decomposition import PCA from skmatter.datasets import load_roy_dataset from skmatter.sample_selection import DirectionalConvexHull .. GENERATED FROM PYTHON SOURCE LINES 24-27 Loads the structures (that also contain properties in the ``info`` field) .. GENERATED FROM PYTHON SOURCE LINES 27-37 .. code-block:: Python roy_data = load_roy_dataset() density = roy_data["densities"] energy = roy_data["energies"] structype = roy_data["structure_types"] iknown = np.where(structype == "known")[0] iothers = np.where(structype != "known")[0] .. GENERATED FROM PYTHON SOURCE LINES 38-44 Energy-density hull ------------------- The Directional Convex Hull routines can be used to compute a conventional density-energy hull .. GENERATED FROM PYTHON SOURCE LINES 44-49 .. code-block:: Python dch_builder = DirectionalConvexHull(low_dim_idx=[0]) dch_builder.fit(density.reshape(-1, 1), energy) .. raw:: html
DirectionalConvexHull(low_dim_idx=[0])
In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook.
On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.


.. GENERATED FROM PYTHON SOURCE LINES 50-53 We can get the indices of the selection, and compute the distance from the hull .. GENERATED FROM PYTHON SOURCE LINES 53-58 .. code-block:: Python sel = dch_builder.selected_idx_ dch_dist = dch_builder.score_samples(density.reshape(-1, 1), energy) .. GENERATED FROM PYTHON SOURCE LINES 59-68 Structures on the hull are stable with respect to synthesis at constant molar volume. Any other structure would lower the energy by decomposing into a mixture of the two nearest structures along the hull. Given that the lattice energy is an imperfect proxy for the free energy, and that synthesis can be performed in other ways than by fixing the density, structures that are not exactly on the hull might also be stable. One can compute a “hull energy” as an indication of how close these structures are to being stable . .. GENERATED FROM PYTHON SOURCE LINES 68-82 .. code-block:: Python fig, ax = plt.subplots(1, 1, figsize=(6, 4)) ax.scatter(density, energy, c=dch_dist, marker=".") ssel = sel[np.argsort(density[sel])] ax.plot(density[ssel], energy[ssel], "k--") ax.set_xlabel("density / g/cm$^3$") ax.set_ylabel("energy / kJ/mol") print( f"Mean hull energy for 'known' stable structures {dch_dist[iknown].mean()} kJ/mol" ) print(f"Mean hull energy for 'other' structures {dch_dist[iothers].mean()} kJ/mol") .. image-sg:: /examples/selection/images/sphx_glr_GCH-ROY_001.png :alt: GCH ROY :srcset: /examples/selection/images/sphx_glr_GCH-ROY_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none Mean hull energy for 'known' stable structures 1.816657381014075 kJ/mol Mean hull energy for 'other' structures 6.312730486304906 kJ/mol .. GENERATED FROM PYTHON SOURCE LINES 83-87 You can also visualize the hull with ``chemiscope``. This runs only in a notebook, and requires having the ``chemiscope`` package installed. .. GENERATED FROM PYTHON SOURCE LINES 87-110 .. code-block:: Python """ import chemiscope chemiscope.show( structures, dict( energy=energy, density=density, hull_energy=dch_dist, structure_type=structype ), settings={ "map": { "x": {"property": "density"}, "y": {"property": "energy"}, "color": {"property": "hull_energy"}, "symbol": "structure_type", "size": {"factor": 35}, }, "structure": [{"unitCell": True, "supercell": {"0": 2, "1": 2, "2": 2}}], }, ) """ .. rst-class:: sphx-glr-script-out .. code-block:: none '\nimport chemiscope\nchemiscope.show(\n structures,\n dict(\n energy=energy, density=density,\n hull_energy=dch_dist, structure_type=structype\n ),\n settings={\n "map": {\n "x": {"property": "density"},\n "y": {"property": "energy"},\n "color": {"property": "hull_energy"},\n "symbol": "structure_type",\n "size": {"factor": 35},\n },\n "structure": [{"unitCell": True, "supercell": {"0": 2, "1": 2, "2": 2}}],\n },\n)\n' .. GENERATED FROM PYTHON SOURCE LINES 111-123 Generalized Convex Hull ----------------------- A GCH is a similar construction, in which generic structural descriptors are used in lieu of composition, density or other thermodynamic constraints. The idea is that configurations that are found close to the GCH are locally stable with respect to structurally-similar configurations. In other terms, one can hope to find a thermodynamic constraint (i.e. synthesis conditions) that act differently on these structures in comparison with the others, and may potentially stabilize them. .. GENERATED FROM PYTHON SOURCE LINES 126-132 A first step is to computes suitable ML descriptors. Here we have used ``rascaline`` to evaluate average SOAP features for the structures. We will load pre-computed features to reduce the dependencies of these examples on external packages, but you can use this as a stub to apply this analysis to other chemical systems .. GENERATED FROM PYTHON SOURCE LINES 132-158 .. code-block:: Python """ from rascaline import SoapPowerSpectrum from equistore import mean_over_samples hypers = { "cutoff": 4, "max_radial": 6, "max_angular": 4, "atomic_gaussian_width": 0.7, "cutoff_function": {"ShiftedCosine": {"width": 0.5}}, "radial_basis": {"Gto": {"accuracy": 1e-6}}, "center_atom_weight": 1.0 } calculator = SoapPowerSpectrum(**hypers) rho2i = calculator.compute(structures) rho2i=rho2i.keys_to_samples(['species_center']).keys_to_properties( ['species_neighbor_1', 'species_neighbor_2']) rho2i_structure = mean_over_samples(rho2i, sample_names=["center", "species_center"]) np.savez("roy_features.npz", feats=rho2i_structure.block(0).values) """ features = roy_data["features"] .. GENERATED FROM PYTHON SOURCE LINES 159-163 Computes PCA projection to generate low-dimensional descriptors that reflect structural diversity. Any other dimensionality reduction scheme could be used in a similar fashion. .. GENERATED FROM PYTHON SOURCE LINES 163-175 .. code-block:: Python pca = PCA(n_components=4) pca_features = pca.fit_transform(features) fig, ax = plt.subplots(1, 1, figsize=(6, 4)) scatter = ax.scatter(pca_features[:, 0], pca_features[:, 1], c=energy) ax.set_xlabel("PCA[1]") ax.set_ylabel("PCA[2]") cbar = fig.colorbar(scatter, ax=ax) cbar.set_label("energy / kJ/mol") .. image-sg:: /examples/selection/images/sphx_glr_GCH-ROY_002.png :alt: GCH ROY :srcset: /examples/selection/images/sphx_glr_GCH-ROY_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 176-178 Builds a convex hull on the first two PCA features .. GENERATED FROM PYTHON SOURCE LINES 178-185 .. code-block:: Python dch_builder = DirectionalConvexHull(low_dim_idx=[0, 1]) dch_builder.fit(pca_features, energy) sel = dch_builder.selected_idx_ dch_dist = dch_builder.score_samples(pca_features, energy) .. GENERATED FROM PYTHON SOURCE LINES 186-188 Generates a 3D Plot .. GENERATED FROM PYTHON SOURCE LINES 188-200 .. code-block:: Python triang = mtri.Triangulation(pca_features[sel, 0], pca_features[sel, 1]) fig = plt.figure(figsize=(7, 5), tight_layout=True) ax = fig.add_subplot(projection="3d") ax.plot_trisurf(triang, energy[sel], color="gray") ax.scatter(pca_features[:, 0], pca_features[:, 1], energy, c=dch_dist) ax.set_xlabel("PCA[1]") ax.set_ylabel("PCA[2]") ax.set_zlabel("energy / kJ/mol\n \n", labelpad=11) ax.view_init(25, 110) .. image-sg:: /examples/selection/images/sphx_glr_GCH-ROY_003.png :alt: GCH ROY :srcset: /examples/selection/images/sphx_glr_GCH-ROY_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 201-205 The GCH construction improves the separation between the hull energies of “known” and hypothetical polymorphs (compare with the density-energy values above) .. GENERATED FROM PYTHON SOURCE LINES 205-212 .. code-block:: Python print( f"Mean hull energy for 'known' stable structures {dch_dist[iknown].mean()} kJ/mol" ) print(f"Mean hull energy for 'other' structures {dch_dist[iothers].mean()} kJ/mol") .. rst-class:: sphx-glr-script-out .. code-block:: none Mean hull energy for 'known' stable structures 0.8539156846959327 kJ/mol Mean hull energy for 'other' structures 5.19873860454363 kJ/mol .. GENERATED FROM PYTHON SOURCE LINES 213-216 Visualize in ``chemiscope``. This runs only in a notebook, and requires having the ``chemiscope`` package installed. .. GENERATED FROM PYTHON SOURCE LINES 216-249 .. code-block:: Python """ import chemiscope for i, f in enumerate(structures): for j in range(len(pca_features[i])): f.info["pca_" + str(j + 1)] = pca_features[i, j] structure_properties = chemiscope.extract_properties(structures) structure_properties.update({"per_atom_energy": energy, "hull_energy": dch_dist}) chemiscope.show( frames=structures, properties=structure_properties, settings={ "map": { "x": {"property": "pca_1"}, "y": {"property": "pca_2"}, "z": {"property": "energy"}, "symbol": "type", "symbol": "type", "color": {"property": "hull_energy"}, "size": {"factor": 35, "mode": "linear", "property": "", "reverse": True}, }, "structure": [ { "bonds": True, "unitCell": True, "keepOrientation": True, } ], }, ) """ .. rst-class:: sphx-glr-script-out .. code-block:: none '\nimport chemiscope\nfor i, f in enumerate(structures):\n for j in range(len(pca_features[i])):\n f.info["pca_" + str(j + 1)] = pca_features[i, j]\nstructure_properties = chemiscope.extract_properties(structures)\nstructure_properties.update({"per_atom_energy": energy, "hull_energy": dch_dist})\nchemiscope.show(\n frames=structures,\n properties=structure_properties,\n settings={\n "map": {\n "x": {"property": "pca_1"},\n "y": {"property": "pca_2"},\n "z": {"property": "energy"},\n "symbol": "type",\n "symbol": "type",\n "color": {"property": "hull_energy"},\n "size": {"factor": 35, "mode": "linear",\n "property": "", "reverse": True},\n },\n "structure": [\n {\n "bonds": True,\n "unitCell": True,\n "keepOrientation": True,\n }\n ],\n },\n)\n' .. _sphx_glr_download_examples_selection_GCH-ROY.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: GCH-ROY.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: GCH-ROY.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: GCH-ROY.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_