diff --git a/config.yml b/config.yml index fcd3a261..bb4b54c7 100644 --- a/config.yml +++ b/config.yml @@ -72,3 +72,8 @@ notebooks: - mat3ra-ade - mat3ra-wode - mat3ra-prode + - name: torch + packages_pyodide: + - emfs:/drive/packages/torch-2.1.0a0-cp311-cp311-emscripten_3_1_45_wasm32.whl + # Needed by MACE + - ssl diff --git a/other/experiments/create_interface_with_relaxation_mace.ipynb b/other/experiments/create_interface_with_relaxation_mace.ipynb new file mode 100644 index 00000000..67d71539 --- /dev/null +++ b/other/experiments/create_interface_with_relaxation_mace.ipynb @@ -0,0 +1,654 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Create an Interface with ZSL and Relax it with MACE\n", + "\n", + "Use Zur and McGill superlattices matching [algorithm](https://doi.org/10.1063/1.3330840) to create interfaces between two materials using the `mat3ra-made` [implementation](https://github.com/Exabyte-io/made).\n", + "\n", + "

Usage

\n", + "\n", + "1. Drop the materials files into the \"uploads\" folder in the JupyterLab file browser\n", + "1. Set Input Parameters below or use the default values\n", + "1. Click \"Run\" > \"Run All\" to run all cells\n", + "1. Wait for the run to complete. Scroll down to view cell results.\n", + "1. Review the strain plot and modify parameters as needed\n", + "\n", + "## Methodology\n", + "\n", + "1. Load materials from JSON files and create substrate and film slabs via `mat3ra-made`\n", + "2. Run ZSL strain matching and plot strain vs interface area\n", + "3. Create the selected interface and convert to ASE atoms with `to_ase()`\n", + "4. Relax the interface with MACE-MP-0 and visualize convergence\n", + "5. Compute interface binding energy decomposition" + ], + "id": "9515ed910c085db8" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 1. Set Input Parameters\n", + "\n", + "### 1.1. Substrate and Film Materials" + ], + "id": "4737d145950b1cc8" + }, + { + "metadata": {}, + "cell_type": "code", + "source": [ + "\n", + "SUBSTRATE_NAME = \"Nickel\"\n", + "FILM_NAME = \"Graphene\"\n", + "SUBSTRATE_INDEX = 0\n", + "FILM_INDEX = 1" + ], + "id": "ef377fcae54adf52", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 1.2. Slab Parameters" + ], + "id": "401cb093b302dbb1" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "SUBSTRATE_MILLER_INDICES = (1, 1, 1)\n", + "SUBSTRATE_THICKNESS = 6 # in atomic layers\n", + "SUBSTRATE_TERMINATION_FORMULA = None # if None, the first termination is used\n", + "\n", + "FILM_MILLER_INDICES = (0, 0, 1)\n", + "FILM_THICKNESS = 1 # in atomic layers\n", + "FILM_TERMINATION_FORMULA = None # if None, the first termination is used\n", + "\n", + "USE_CONVENTIONAL_CELL = True" + ], + "id": "edd2f829bc0a72ce", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 1.3. Interface and Relaxation Parameters" + ], + "id": "8fecb819abaa2ca7" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "INTERFACE_DISTANCE = None # gap between substrate and film, in Angstrom, if None, the distance will be set to the sum of the covalent radii of the two materials\n", + "INTERFACE_VACUUM = 10.0 # vacuum over film, in Angstrom\n", + "REDUCE_RESULT_CELL_TO_PRIMITIVE = True\n", + "\n", + "MAX_AREA = 150 # in Angstrom^2\n", + "MAX_AREA_TOLERANCE = 0.09\n", + "MAX_LENGTH_TOLERANCE = 0.05\n", + "MAX_ANGLE_TOLERANCE = 0.02\n", + "\n", + "RELAXATION_PARAMETERS = {\n", + " \"FMAX\": 0.01,\n", + "}\n", + "MACE_MODEL = \"large\" # \"small\", \"medium\", or \"large\" MACE-MP-0 foundation model" + ], + "id": "57dc565952aa2e4e", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 2. Install Packages" + ], + "id": "4e89f2d820acb1ed" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "# !pip install mat3ra-made torch \"mace-torch\" ase \"e3nn==0.4.4\" \"numpy<=1.26.4\" pymatgen" + ], + "id": "810243b01ba923b5", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 3. Load Materials" + ], + "id": "f89d3c98ddce2ab5" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from mat3ra.made.material import Material\n", + "from mat3ra.standata.materials import Materials\n", + "\n", + "substrate = Material.create(Materials.get_by_name_first_match(SUBSTRATE_NAME))\n", + "film = Material.create(Materials.get_by_name_first_match(FILM_NAME))\n", + "\n", + "print(\"Substrate:\", substrate.name)\n", + "print(\"Film: \", film.name)" + ], + "id": "8fd400dace70549e", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 3.1. Visualize Input Materials" + ], + "id": "42f12abf6b65aa2c" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from utils.visualize import visualize_materials as visualize\n", + "\n", + "visualize([substrate, film], repetitions=[3, 3, 1], rotation=\"0x\")" + ], + "id": "88b4ce8e27118174", + "outputs": [], + "execution_count": null + }, + { + "metadata": {}, + "cell_type": "markdown", + "source": "## 3.2. Calculate nearest neighbor distance for each material to inform interface distance choice\n", + "id": "51c20f951b402e48" + }, + { + "metadata": {}, + "cell_type": "code", + "source": [ + "from mat3ra.made.tools.build_components.entities.reusable.three_dimensional.supercell.helpers import create_supercell\n", + "from mat3ra.made.tools.analyze.rdf import RadialDistributionFunction\n", + "### 7.3. Plot Radial Distribution Functions\n", + "from utils.plot import plot_rdf\n", + "\n", + "substrate_supercell = create_supercell(substrate, scaling_factor=[3, 3, 3])\n", + "film_supercell = create_supercell(film, scaling_factor=[3, 3, 3])\n", + "\n", + "rdf_substrate = RadialDistributionFunction.from_material(substrate_supercell, cutoff=5.0)\n", + "rdf_film = RadialDistributionFunction.from_material(film_supercell, cutoff=5.0)\n", + "\n", + "first_peak_substrate = rdf_substrate.first_peak_distance\n", + "first_peak_film = rdf_film.first_peak_distance\n", + "\n", + "print(f\"First RDF peak for substrate ({substrate.name}): {first_peak_substrate:.3f} Å\")\n", + "print(f\"First RDF peak for film ({film.name}): {first_peak_film:.3f} Å\")\n", + "\n", + "if INTERFACE_DISTANCE is None:\n", + " INTERFACE_DISTANCE = (first_peak_substrate + first_peak_film) / 2\n", + " print(f\"Setting interface distance to {INTERFACE_DISTANCE:.3f} Å based on RDF peaks\")\n", + "\n", + "plot_rdf(substrate_supercell, cutoff=5.0)\n", + "plot_rdf(film_supercell, cutoff=5.0)" + ], + "id": "ca5955d87b780158", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 4. Configure Slabs\n", + "\n", + "### 4.1. Get Possible Terminations" + ], + "id": "46b23cc60cd0454e" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from mat3ra.made.tools.helpers import get_slab_terminations\n", + "\n", + "film_slab_terminations = get_slab_terminations(material=film, miller_indices=FILM_MILLER_INDICES)\n", + "substrate_slab_terminations = get_slab_terminations(material=substrate, miller_indices=SUBSTRATE_MILLER_INDICES)\n", + "print(\"Film slab terminations: \", film_slab_terminations)\n", + "print(\"Substrate slab terminations:\", substrate_slab_terminations)" + ], + "id": "be35fa07cf206e76", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 4.2. Visualize Slabs for All Possible Terminations" + ], + "id": "caf19a75fb14f280" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from mat3ra.made.tools.helpers import create_slab, select_slab_termination\n", + "from mat3ra.made.tools.helpers import create_interface_zsl_between_slabs\n", + "\n", + "film_slabs = [\n", + " create_slab(film, miller_indices=FILM_MILLER_INDICES, termination_top=t, vacuum=0)\n", + " for t in film_slab_terminations\n", + "]\n", + "substrate_slabs = [\n", + " create_slab(substrate, miller_indices=SUBSTRATE_MILLER_INDICES, termination_top=t, vacuum=0, number_of_layers=4)\n", + " for t in substrate_slab_terminations\n", + "]\n", + "\n", + "visualize(\n", + " [{\"material\": s, \"title\": str(t)} for s, t in zip(film_slabs, film_slab_terminations)],\n", + " repetitions=[3, 3, 1], rotation=\"-90x\",\n", + ")\n", + "visualize(\n", + " [{\"material\": s, \"title\": str(t)} for s, t in zip(substrate_slabs, substrate_slab_terminations)],\n", + " repetitions=[3, 3, 1], rotation=\"-90x\",\n", + ")" + ], + "id": "7bae67a0b5821aa7", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 4.3. Create Substrate and Film Slabs" + ], + "id": "a692280e024d5883" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from mat3ra.made.tools.build.pristine_structures.two_dimensional.slab import SlabConfiguration, SlabBuilder\n", + "\n", + "substrate_slab_config = SlabConfiguration.from_parameters(\n", + " material_or_dict=substrate,\n", + " miller_indices=SUBSTRATE_MILLER_INDICES,\n", + " number_of_layers=SUBSTRATE_THICKNESS,\n", + " vacuum=0.0,\n", + " termination_top_formula=SUBSTRATE_TERMINATION_FORMULA,\n", + " use_conventional_cell=USE_CONVENTIONAL_CELL,\n", + ")\n", + "film_slab_config = SlabConfiguration.from_parameters(\n", + " material_or_dict=film,\n", + " miller_indices=FILM_MILLER_INDICES,\n", + " number_of_layers=FILM_THICKNESS,\n", + " vacuum=0.0,\n", + " termination_bottom_formula=FILM_TERMINATION_FORMULA,\n", + " use_conventional_cell=USE_CONVENTIONAL_CELL,\n", + ")\n", + "\n", + "substrate_slab = SlabBuilder().get_material(substrate_slab_config)\n", + "film_slab = SlabBuilder().get_material(film_slab_config)" + ], + "id": "89f15267bd2f6f66", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 5. Find Interfaces with ZSL Strain Matching\n", + "\n", + "### 5.1. Initialize ZSL Analyzer" + ], + "id": "4c4b81397d2eb68" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from mat3ra.made.tools.analyze.interface import ZSLInterfaceAnalyzer\n", + "\n", + "zsl_analyzer = ZSLInterfaceAnalyzer(\n", + " substrate_slab_configuration=substrate_slab_config,\n", + " film_slab_configuration=film_slab_config,\n", + " max_area=MAX_AREA,\n", + " max_area_ratio_tol=MAX_AREA_TOLERANCE,\n", + " max_length_tol=MAX_LENGTH_TOLERANCE,\n", + " max_angle_tol=MAX_ANGLE_TOLERANCE,\n", + " reduce_result_cell=False,\n", + ")" + ], + "id": "b04a8543e8fcf455", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 5.2. Generate and Plot Matches" + ], + "id": "a1956430738a9a9a" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from utils.plot import plot_strain_vs_area\n", + "\n", + "PLOT_SETTINGS = {\n", + " \"HEIGHT\": 600,\n", + " \"X_SCALE\": \"log\",\n", + " \"Y_SCALE\": \"log\",\n", + "}\n", + "\n", + "matches = zsl_analyzer.zsl_match_holders\n", + "print(f\"Found {len(matches)} matches\")\n", + "plot_strain_vs_area(matches, PLOT_SETTINGS)" + ], + "id": "8682bbecf48aa30b", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 5.3. Select the Interface\n", + "\n", + "Choose the match index from the plot above (index 0 has the lowest strain)." + ], + "id": "245b45e7bb4a7ac3" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "selected_index = 0" + ], + "id": "a83a9d7a43391187", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 6. Create the Interface" + ], + "id": "5360346b47cf21ed" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "interface = create_interface_zsl_between_slabs(\n", + " substrate_slab=substrate_slab,\n", + " film_slab=film_slab,\n", + " gap=INTERFACE_DISTANCE,\n", + " vacuum=INTERFACE_VACUUM,\n", + " match_id=selected_index,\n", + " max_area=MAX_AREA,\n", + " max_area_ratio_tol=MAX_AREA_TOLERANCE,\n", + " max_length_tol=MAX_LENGTH_TOLERANCE,\n", + " max_angle_tol=MAX_ANGLE_TOLERANCE,\n", + " reduce_result_cell_to_primitive=REDUCE_RESULT_CELL_TO_PRIMITIVE,\n", + ")" + ], + "id": "4d413413dfaa3f6d", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 6.1. Visualize Interface" + ], + "id": "24cc0a761f161676" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from utils.visualize import ViewersEnum\n", + "\n", + "visualize([{\"material\": interface, \"title\": interface.name}], viewer=ViewersEnum.wave)\n", + "visualize(interface, repetitions=[1, 1, 1], rotation=\"-90x\")" + ], + "id": "fcbb4e6c1de21233", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## 7. Apply Relaxation\n", + "### 7.1. Relax with MACE" + ], + "id": "e64688fc18c49bb6" + }, + { + "cell_type": "code", + "metadata": { + "jupyter": { + "is_executing": true + } + }, + "source": [ + "import plotly.graph_objs as go\n", + "from IPython.display import display\n", + "from plotly.subplots import make_subplots\n", + "from src.utils import ase_to_poscar\n", + "from mat3ra.made.tools.convert import to_ase\n", + "from ase.optimize import BFGS\n", + "from mace.calculators import mace_mp\n", + "\n", + "calculator = mace_mp(model=MACE_MODEL, dispersion=True, default_dtype=\"float32\", device=\"cpu\")\n", + "\n", + "ase_interface = to_ase(interface)\n", + "ase_interface.set_calculator(calculator)\n", + "dyn = BFGS(ase_interface)\n", + "\n", + "steps = []\n", + "energies = []\n", + "\n", + "fig = make_subplots(rows=1, cols=1, specs=[[{\"type\": \"scatter\"}]])\n", + "scatter = go.Scatter(x=[], y=[], mode=\"lines+markers\", name=\"Energy\")\n", + "fig.add_trace(scatter)\n", + "fig.update_layout(title_text=\"Real-time Optimization Progress\", xaxis_title=\"Step\", yaxis_title=\"Energy (eV)\")\n", + "\n", + "f = go.FigureWidget(fig)\n", + "display(f)\n", + "\n", + "\n", + "def plotly_callback():\n", + " step = dyn.nsteps\n", + " energy = ase_interface.get_total_energy()\n", + " steps.append(step)\n", + " energies.append(energy)\n", + " print(f\"Step: {step}, Energy: {energy:.4f} eV\")\n", + " with f.batch_update():\n", + " f.data[0].x = steps\n", + " f.data[0].y = energies\n", + "\n", + "\n", + "dyn.attach(plotly_callback, interval=1)\n", + "dyn.run(fmax=RELAXATION_PARAMETERS[\"FMAX\"])\n", + "\n", + "ase_original_interface = to_ase(interface)\n", + "ase_original_interface.set_calculator(calculator)\n", + "ase_final_interface = ase_interface\n", + "\n", + "original_energy = ase_original_interface.get_total_energy()\n", + "relaxed_energy = ase_interface.get_total_energy()\n", + "\n", + "print(\"Original structure:\\n\", ase_to_poscar(ase_original_interface))\n", + "print(\"\\nRelaxed structure:\\n\", ase_to_poscar(ase_final_interface))\n", + "print(f\"The final energy is {float(relaxed_energy):.3f} eV.\")" + ], + "id": "3d8746a77f71bab5", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 7.2. View Structure Before and After Relaxation" + ], + "id": "abfa372909a96bf8" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "from mat3ra.made.tools.convert import from_ase\n", + "\n", + "\n", + "def atoms_to_material(atoms, title):\n", + " material = Material.create(from_ase(atoms))\n", + " material.name = title\n", + " return material\n", + "\n", + "\n", + "material_original = atoms_to_material(ase_original_interface, f\"Original E={original_energy:.3f} eV\")\n", + "material_relaxed = atoms_to_material(ase_final_interface, f\"Relaxed E={relaxed_energy:.3f} eV\")\n", + "\n", + "visualize(\n", + " [\n", + " {\"material\": material_original, \"title\": material_original.name},\n", + " {\"material\": material_relaxed, \"title\": material_relaxed.name},\n", + " ],\n", + " viewer=ViewersEnum.wave,\n", + ")" + ], + "id": "9565d0931b198f63", + "outputs": [], + "execution_count": null + }, + { + "metadata": {}, + "cell_type": "markdown", + "source": "## 7.4. Output interlayer distance before and after relaxation", + "id": "e4b49774283e5517" + }, + { + "metadata": {}, + "cell_type": "code", + "source": [ + "from mat3ra.made.tools.analyze.other import get_average_interlayer_distance\n", + "\n", + "print(f\"Interlayer distance before relaxation: {get_average_interlayer_distance(material_original, 0, 1):.4f} Å\")\n", + "print(f\"Interlayer distance after relaxation: {get_average_interlayer_distance(material_relaxed, 0, 1):.4f} Å\")" + ], + "id": "6dd00402bc2e9d59", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### 7.4. Calculate Interface Energy" + ], + "id": "7ee3d26311a30687" + }, + { + "metadata": {}, + "cell_type": "code", + "source": [ + "def filter_atoms_by_tag(atoms, material_index):\n", + " return atoms[atoms.get_tags() == material_index]\n", + "\n", + "\n", + "def calculate_energy(atoms, calc):\n", + " atoms.set_calculator(calc)\n", + " return atoms.get_total_energy()\n", + "\n", + "\n", + "def calculate_delta_energy(total_energy, *component_energies):\n", + " return total_energy - sum(component_energies)\n", + "\n", + "\n", + "substrate_original = filter_atoms_by_tag(ase_original_interface, SUBSTRATE_INDEX)\n", + "layer_original = filter_atoms_by_tag(ase_original_interface, FILM_INDEX)\n", + "substrate_relaxed = filter_atoms_by_tag(ase_final_interface, SUBSTRATE_INDEX)\n", + "layer_relaxed = filter_atoms_by_tag(ase_final_interface, FILM_INDEX)\n", + "\n", + "original_substrate_energy = calculate_energy(substrate_original, calculator)\n", + "original_layer_energy = calculate_energy(layer_original, calculator)\n", + "relaxed_substrate_energy = calculate_energy(substrate_relaxed, calculator)\n", + "relaxed_layer_energy = calculate_energy(layer_relaxed, calculator)\n", + "\n", + "delta_original = calculate_delta_energy(original_energy, original_substrate_energy, original_layer_energy)\n", + "delta_relaxed = calculate_delta_energy(relaxed_energy, relaxed_substrate_energy, relaxed_layer_energy)\n", + "\n", + "area = ase_original_interface.get_volume() / ase_original_interface.cell[2, 2]\n", + "n_interface = ase_final_interface.get_global_number_of_atoms()\n", + "n_substrate = substrate_relaxed.get_global_number_of_atoms()\n", + "n_layer = layer_relaxed.get_global_number_of_atoms()\n", + "effective_delta_relaxed = (\n", + " relaxed_energy / n_interface\n", + " - (relaxed_substrate_energy / n_substrate + relaxed_layer_energy / n_layer)\n", + " ) / (2 * area)\n", + "\n", + "print(f\"Original Substrate energy: {original_substrate_energy:.4f} eV\")\n", + "print(f\"Relaxed Substrate energy: {relaxed_substrate_energy:.4f} eV\")\n", + "print(f\"Original Layer energy: {original_layer_energy:.4f} eV\")\n", + "print(f\"Relaxed Layer energy: {relaxed_layer_energy:.4f} eV\")\n", + "print(\"\\nDelta between interface energy and sum of component energies\")\n", + "print(f\"Original Delta: {delta_original:.4f} eV\")\n", + "print(f\"Relaxed Delta: {delta_relaxed:.4f} eV\")\n", + "print(f\"Original Delta per area: {delta_original / area:.4f} eV/Ang^2\")\n", + "print(f\"Relaxed Delta per area: {delta_relaxed / area:.4f} eV/Ang^2\")\n", + "print(f\"Relaxed interface energy: {relaxed_energy:.4f} eV\")\n", + "print(\n", + " f\"Effective relaxed Delta per area: {effective_delta_relaxed:.4f} eV/Ang^2 ({effective_delta_relaxed / 0.16:.4f} J/m^2)\")" + ], + "id": "79ea902feda4a8d3", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## References\n", + "\n", + "[1] mat3ra-made interface builder: https://github.com/Exabyte-io/made \n", + "[2] MACE-MP-0 foundation model: https://github.com/ACEsuit/mace?tab=readme-ov-file#foundation-models " + ], + "id": "2f60fdb73e44c09c" + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "name": "python", + "version": "3.10.0" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/other/experiments/jupyterlite/create_interface_with_relaxation_mace.ipynb b/other/experiments/jupyterlite/create_interface_with_relaxation_mace.ipynb new file mode 100644 index 00000000..dd48696f --- /dev/null +++ b/other/experiments/jupyterlite/create_interface_with_relaxation_mace.ipynb @@ -0,0 +1,651 @@ +{ + "metadata": { + "kernelspec": { + "name": "python", + "display_name": "Python (Pyodide)", + "language": "python" + }, + "language_info": { + "codemirror_mode": { + "name": "python", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8" + } + }, + "nbformat_minor": 5, + "nbformat": 4, + "cells": [ + { + "cell_type": "markdown", + "source": "# Create an Interface with ZSL and Relax it with MACE\n\nUse Zur and McGill superlattices matching [algorithm](https://doi.org/10.1063/1.3330840) to create interfaces between two materials using the `mat3ra-made` [implementation](https://github.com/Exabyte-io/made).\n\n

Usage

\n\n1. Drop the materials files into the \"uploads\" folder in the JupyterLab file browser\n1. Set Input Parameters below or use the default values\n1. Click \"Run\" > \"Run All\" to run all cells\n1. Wait for the run to complete. Scroll down to view cell results.\n1. Review the strain plot and modify parameters as needed\n\n## Methodology\n\n1. Load materials from JSON files and create substrate and film slabs via `mat3ra-made`\n2. Run ZSL strain matching and plot strain vs interface area\n3. Create the selected interface and convert to ASE atoms with `to_ase()`\n4. Relax the interface with MACE-MP-0 and visualize convergence\n5. Compute interface binding energy decomposition", + "metadata": {}, + "id": "9515ed910c085db8" + }, + { + "cell_type": "markdown", + "source": "## 1. Set Input Parameters\n\n### 1.1. Substrate and Film Materials", + "metadata": {}, + "id": "4737d145950b1cc8" + }, + { + "cell_type": "code", + "source": "\nSUBSTRATE_NAME = \"Nickel\"\nFILM_NAME = \"Graphene\"\nSUBSTRATE_INDEX = 0\nFILM_INDEX = 1", + "metadata": { + "trusted": true + }, + "id": "ef377fcae54adf52", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 1.2. Slab Parameters", + "metadata": {}, + "id": "401cb093b302dbb1" + }, + { + "cell_type": "code", + "source": "SUBSTRATE_MILLER_INDICES = (1, 1, 1)\nSUBSTRATE_THICKNESS = 6 # in atomic layers\nSUBSTRATE_TERMINATION_FORMULA = None # if None, the first termination is used\n\nFILM_MILLER_INDICES = (0, 0, 1)\nFILM_THICKNESS = 1 # in atomic layers\nFILM_TERMINATION_FORMULA = None # if None, the first termination is used\n\nUSE_CONVENTIONAL_CELL = True", + "metadata": { + "trusted": true + }, + "id": "edd2f829bc0a72ce", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 1.3. Interface and Relaxation Parameters", + "metadata": {}, + "id": "8fecb819abaa2ca7" + }, + { + "cell_type": "code", + "source": "INTERFACE_DISTANCE = None # gap between substrate and film, in Angstrom, if None, the distance will be set to the sum of the covalent radii of the two materials\nINTERFACE_VACUUM = 10.0 # vacuum over film, in Angstrom\nREDUCE_RESULT_CELL_TO_PRIMITIVE = True\n\nMAX_AREA = 150 # in Angstrom^2\nMAX_AREA_TOLERANCE = 0.09\nMAX_LENGTH_TOLERANCE = 0.05\nMAX_ANGLE_TOLERANCE = 0.02\n\nRELAXATION_PARAMETERS = {\n \"FMAX\": 0.05,\n}\nMACE_MODEL = \"large\" # \"small\", \"medium\", or \"large\" MACE-MP-0 foundation model", + "metadata": { + "trusted": true + }, + "id": "57dc565952aa2e4e", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 2. Install Packages", + "metadata": {}, + "id": "4e89f2d820acb1ed" + }, + { + "cell_type": "code", + "source": [ + "\n", + "try:\n", + " from pyodide.http import pyfetch as _pyodide_pyfetch\n", + "except ImportError:\n", + " _pyodide_pyfetch = None\n", + "\n", + "# Mirrors mace.calculators.foundations_models.mace_mp_urls (mace-torch 0.3.x).\n", + "MODEL_PATHS_MAP = {\n", + "# \"small\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0/2023-12-10-mace-128-L0_energy_epoch-249.model\",\n", + "# \"small\": \"/drive/packages/models/2023-12-10-mace-128-L0_energy_epoch-249.model\",\n", + "# \"medium\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0/2023-12-03-mace-128-L1_epoch-199.model\",\n", + " \"large\": \"/drive/packages/models/MACE_MPtrj_2022.9.model\",\n", + "# \"large\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0/MACE_MPtrj_2022.9.model\",\n", + "# \"small-0b\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0b/mace_agnesi_small.model\",\n", + "# \"medium-0b\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0b/mace_agnesi_medium.model\",\n", + "# \"small-0b2\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0b2/mace-small-density-agnesi-stress.model\",\n", + "# \"medium-0b2\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0b2/mace-medium-density-agnesi-stress.model\",\n", + "# \"large-0b2\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0b2/mace-large-density-agnesi-stress.model\",\n", + "# \"medium-0b3\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mp_0b3/mace-mp-0b3-medium.model\",\n", + "# \"medium-mpa-0\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_mpa_0/mace-mpa-0-medium.model\",\n", + "# \"small-omat-0\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_omat_0/mace-omat-0-small.model\",\n", + "# \"medium-omat-0\": \"https://github.com/ACEsuit/mace-mp/releases/download/mace_omat_0/mace-omat-0-medium.model\",\n", + "# \"mace-matpes-pbe-0\": \"https://github.com/ACEsuit/mace-foundations/releases/download/mace_matpes_0/MACE-matpes-pbe-omat-ft.model\",\n", + "# \"mace-matpes-r2scan-0\": \"https://github.com/ACEsuit/mace-foundations/releases/download/mace_matpes_0/MACE-matpes-r2scan-omat-ft.model\",\n", + "}\n", + "\n", + "def _matscipy_neighbour_list_compat(quantities, atoms=None, cutoff=None, positions=None, cell=None, pbc=None, **_):\n", + " \"\"\"Translate matscipy-style keyword-arg call into an ASE neighbor_list call.\"\"\"\n", + " from ase import Atoms as _Atoms\n", + " from ase.neighborlist import neighbor_list as _ase_neighbor_list\n", + " if atoms is None:\n", + " atoms = _Atoms(positions=positions, cell=cell, pbc=pbc if pbc is not None else [False, False, False])\n", + " return _ase_neighbor_list(quantities, atoms, cutoff)\n", + "\n", + "def _tensor_array_compat(self, dtype=None):\n", + " \"\"\"Replacement for Tensor.__array__ in Pyodide where tensor.numpy() is unavailable.\"\"\"\n", + " import numpy as _np\n", + " arr = _np.array(self.tolist())\n", + " return arr.astype(dtype) if dtype is not None else arr\n", + "\n", + "\n", + "class _TorchDFTD3CalculatorStub:\n", + " def __init__(self, *args, **kwargs):\n", + " raise RuntimeError(\"torch_dftd is unavailable in Pyodide; use dispersion=False with mace_mp()\")\n", + "\n", + " \n", + "\n", + "class _LoggingTensorModeStub:\n", + " def __enter__(self):\n", + " return self\n", + "\n", + " def __exit__(self, *args):\n", + " return False\n", + "\n", + "\n", + "def _capture_logs_stub(*args, **kwargs):\n", + " return _LoggingTensorModeStub()\n", + "\n", + "\n", + "def patch_mace_for_pyodide():\n", + "\n", + " \"\"\"Stub modules absent in pyodide's stripped torch build and missing C-extension packages.\"\"\"\n", + " _internal = types.ModuleType(\"torch.testing._internal\")\n", + " _internal.__path__ = [] # __path__ marks it as a package so sub-imports resolve\n", + " _internal.__package__ = \"torch.testing._internal\"\n", + "\n", + " _common_utils = types.ModuleType(\"torch.testing._internal.common_utils\")\n", + " _common_utils.dtype_abbrs = {}\n", + "\n", + " _logging_tensor = types.ModuleType(\"torch.testing._internal.logging_tensor\")\n", + " _logging_tensor.LoggingTensorMode = _LoggingTensorModeStub\n", + " _logging_tensor.capture_logs = _capture_logs_stub\n", + "\n", + " _internal.common_utils = _common_utils\n", + " _internal.logging_tensor = _logging_tensor\n", + " sys.modules[\"torch.testing._internal\"] = _internal\n", + " sys.modules[\"torch.testing._internal.common_utils\"] = _common_utils\n", + " sys.modules[\"torch.testing._internal.logging_tensor\"] = _logging_tensor\n", + "\n", + " _matscipy = types.ModuleType(\"matscipy\")\n", + " _matscipy.__path__ = []\n", + " _matscipy.__package__ = \"matscipy\"\n", + " _matscipy_neighbours = types.ModuleType(\"matscipy.neighbours\")\n", + " _matscipy_neighbours.neighbour_list = _matscipy_neighbour_list_compat\n", + " _matscipy.neighbours = _matscipy_neighbours\n", + " sys.modules[\"matscipy\"] = _matscipy\n", + " sys.modules[\"matscipy.neighbours\"] = _matscipy_neighbours\n", + " \n", + " \n", + " # lmdb, orjson, h5py are C-extension packages that cannot be compiled in\n", + " # Emscripten (no subprocess / cc support). They are only used by the\n", + " # training/dataset-loading paths (mace.data.lmdb_dataset,\n", + " # mace.data.hdf5_dataset) which are pulled in by mace.data.__init__ but\n", + " # are never exercised during Pyodide inference. Stubs let the import\n", + " # succeed; any actual call to these APIs would raise at runtime, which is\n", + " # the correct behaviour for an unsupported operation.\n", + " for _pkg in (\"lmdb\", \"h5py\"):\n", + " if _pkg not in sys.modules:\n", + " sys.modules[_pkg] = types.ModuleType(_pkg)\n", + "\n", + " # mace.data.atomic_data accesses torch_geometric.data.Data at class-definition\n", + " # time. In Pyodide the submodule attribute can be unset when mace.tools.__init__\n", + " # is still mid-import (circular via mace.cli.visualise_train) at the point where\n", + " # train.py runs `from . import torch_geometric` – so the subpackage exists in\n", + " # sys.modules but .data has not yet been stamped onto it. Pre-importing here,\n", + " # with the matscipy stub already in place, forces the full init to complete and\n", + " # guarantees the attribute is set before any mace.data import runs.\n", + " try:\n", + " import importlib as _importlib\n", + " _tg = _importlib.import_module(\"mace.tools.torch_geometric\")\n", + " _tg_data = _importlib.import_module(\"mace.tools.torch_geometric.data\")\n", + " _tg.data = _tg_data\n", + " except Exception:\n", + " pass\n", + " \n", + " \n", + " # Pyodide's WASM torch build disables Tensor.numpy() and __array__.\n", + " # MACECalculator.__init__ calls np.array([model.r_max.cpu(), ...]) which\n", + " # triggers __array__ on each tensor element. Route both through tolist().\n", + " try:\n", + " import torch as _torch\n", + " except ImportError:\n", + " _torch = None\n", + " # Pyodide's WASM torch build disables Tensor.numpy() and __array__.\n", + " # MACECalculator.__init__ calls np.array([model.r_max.cpu(), ...]) which\n", + " # triggers __array__ on each tensor element. Route both through tolist().\n", + " if _torch is not None:\n", + " import numpy as _np\n", + " _torch.Tensor.__array__ = _tensor_array_compat\n", + " _torch.Tensor.numpy = lambda self: _np.array(self.detach().tolist())\n", + "\n", + " # torch.compiler exists in the WASM build but is_compiling is absent;\n", + " # mace.modules.utils.prepare_graph calls it unconditionally.\n", + " if not hasattr(_torch, \"compiler\"):\n", + " _torch.compiler = types.ModuleType(\"torch.compiler\")\n", + " if not hasattr(_torch.compiler, \"is_compiling\"):\n", + " _torch.compiler.is_compiling = lambda: False\n", + " \n", + " # torch.linalg.det calls LAPACK, which is absent in the WASM build.\n", + " # mace.modules.utils.compute_forces_virials uses it to compute the\n", + " # unit-cell volume (a 3×3 determinant) when compute_stress=True.\n", + " # Redirect to NumPy for the fallback case.\n", + " _orig_linalg_det = _torch.linalg.det\n", + "\n", + " def _safe_linalg_det(A):\n", + " try:\n", + " return _orig_linalg_det(A)\n", + " except RuntimeError:\n", + " import numpy as _np\n", + " result = _np.linalg.det(_np.array(A.detach().tolist()))\n", + " return _torch.tensor(result, dtype=A.dtype, device=A.device)\n", + "\n", + " _torch.linalg.det = _safe_linalg_det\n", + "\n", + "\n", + "\n", + "async def download_mace_model_pyodide(model = None) -> str:\n", + " \"\"\"Download a MACE foundation model in Pyodide using the browser's fetch API.\n", + "\n", + " urllib.request does not support HTTPS in Pyodide WASM; this function uses\n", + " pyodide.http.pyfetch instead, which delegates to the browser's native fetch.\n", + "\n", + " Args:\n", + " model: MACE model name (e.g. \"medium\", \"medium-mpa-0\") or a direct HTTPS URL.\n", + " None defaults to \"medium-mpa-0\".\n", + " Returns:\n", + " Local filesystem path to the downloaded model file.\n", + " \"\"\"\n", + " if _pyodide_pyfetch is None:\n", + " raise RuntimeError(\"pyodide.http.pyfetch unavailable; call mace_mp() directly in standard Python.\")\n", + " url = MODEL_PATHS_MAP.get(model or \"medium-mpa-0\", model)\n", + " if url is None or not url.startswith(\"http\"):\n", + " return model # already a local path\n", + " cache_dir = os.path.join(os.path.expanduser(\"~\"), \".cache\", \"mace\")\n", + " os.makedirs(cache_dir, exist_ok=True)\n", + " filename = \"\".join(c for c in os.path.basename(url) if c.isalnum() or c in \"_\")\n", + " local_path = os.path.join(cache_dir, filename)\n", + " if not os.path.isfile(local_path):\n", + " print(f\"Downloading MACE model from {url!r}\")\n", + " response = await _pyodide_pyfetch(url)\n", + " with open(local_path, \"wb\") as f:\n", + " f.write(await response.bytes())\n", + " print(f\"Cached MACE model to {local_path}\")\n", + " return local_path\n", + "\n" + ], + "metadata": { + "trusted": true + }, + "id": "e58d503d-7686-4320-b182-1a9be23a3f56", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "source": [ + "import sys\n", + "\n", + "if sys.platform == \"emscripten\":\n", + " import micropip\n", + "\n", + "\n", + " await micropip.install(\"mat3ra-api-examples\", deps=False)\n", + " await micropip.install(\"mat3ra-utils\")\n", + " from mat3ra.utils.jupyterlite.packages import install_packages\n", + " \n", + " import types\n", + "\n", + " await install_packages(\"api_examples|torch\")\n", + " await micropip.install(\"torch-dftd\")\n", + "\n", + " await micropip.install(\"opt_einsum\")\n", + " await micropip.install(\"opt_einsum_fx\", deps=False)\n", + " await micropip.install(\"e3nn==0.4.4\", deps=False)\n", + "\n", + " # MACE training deps\n", + " await micropip.install(\"prettytable\")\n", + " await micropip.install(\"torch_ema\", deps=False)\n", + " await micropip.install(\"lightning-utilities\", deps=False)\n", + " await micropip.install(\"torchmetrics\", deps=False)\n", + " \n", + " #\n", + " await micropip.install(\"ssl\")\n", + " await micropip.install(\"h5py\")\n", + " await micropip.install(\"lmdb\")\n", + "\n", + " await micropip.install(\"orjson\")\n", + " await micropip.install(\"anywidget\")\n", + " \n", + " # MACE\n", + " await micropip.install(\"mace-torch\", deps=False)\n", + " patch_mace_for_pyodide()\n" + ], + "metadata": { + "trusted": true + }, + "id": "3cdce8f3-7508-49bd-9ee1-e9dc99bf1a2c", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "source": "", + "metadata": { + "trusted": true + }, + "id": "a660f174-a0d3-4beb-bc84-8cfc36471281", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "code", + "source": "# micropip.list()", + "metadata": { + "trusted": true + }, + "id": "810243b01ba923b5", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 3. Load Materials", + "metadata": {}, + "id": "f89d3c98ddce2ab5" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.material import Material\nfrom mat3ra.standata.materials import Materials\n\nsubstrate = Material.create(Materials.get_by_name_first_match(SUBSTRATE_NAME))\nfilm = Material.create(Materials.get_by_name_first_match(FILM_NAME))\n\nprint(\"Substrate:\", substrate.name)\nprint(\"Film: \", film.name)", + "metadata": { + "trusted": true + }, + "id": "8fd400dace70549e", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 3.1. Visualize Input Materials", + "metadata": {}, + "id": "42f12abf6b65aa2c" + }, + { + "cell_type": "code", + "source": "from utils.visualize import visualize_materials as visualize\n\nvisualize([substrate, film], repetitions=[3, 3, 1], rotation=\"0x\")", + "metadata": { + "trusted": true + }, + "id": "88b4ce8e27118174", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 3.2. Calculate nearest neighbor distance for each material to inform interface distance choice\n", + "metadata": {}, + "id": "51c20f951b402e48" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.build_components.entities.reusable.three_dimensional.supercell.helpers import create_supercell\nfrom mat3ra.made.tools.analyze.rdf import RadialDistributionFunction\n### 7.3. Plot Radial Distribution Functions\nfrom utils.plot import plot_rdf\n\nsubstrate_supercell = create_supercell(substrate, scaling_factor=[3, 3, 3])\nfilm_supercell = create_supercell(film, scaling_factor=[3, 3, 3])\n\nrdf_substrate = RadialDistributionFunction.from_material(substrate_supercell, cutoff=5.0)\nrdf_film = RadialDistributionFunction.from_material(film_supercell, cutoff=5.0)\n\nfirst_peak_substrate = rdf_substrate.first_peak_distance\nfirst_peak_film = rdf_film.first_peak_distance\n\nprint(f\"First RDF peak for substrate ({substrate.name}): {first_peak_substrate:.3f} Å\")\nprint(f\"First RDF peak for film ({film.name}): {first_peak_film:.3f} Å\")\n\nif INTERFACE_DISTANCE is None:\n INTERFACE_DISTANCE = (first_peak_substrate + first_peak_film) / 2\n print(f\"Setting interface distance to {INTERFACE_DISTANCE:.3f} Å based on RDF peaks\")\n\n# plot_rdf(substrate_supercell, cutoff=5.0)\n# plot_rdf(film_supercell, cutoff=5.0)", + "metadata": { + "trusted": true + }, + "id": "ca5955d87b780158", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 4. Configure Slabs\n\n### 4.1. Get Possible Terminations", + "metadata": {}, + "id": "46b23cc60cd0454e" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.helpers import get_slab_terminations\n\nfilm_slab_terminations = get_slab_terminations(material=film, miller_indices=FILM_MILLER_INDICES)\nsubstrate_slab_terminations = get_slab_terminations(material=substrate, miller_indices=SUBSTRATE_MILLER_INDICES)\nprint(\"Film slab terminations: \", film_slab_terminations)\nprint(\"Substrate slab terminations:\", substrate_slab_terminations)", + "metadata": { + "trusted": true + }, + "id": "be35fa07cf206e76", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 4.2. Visualize Slabs for All Possible Terminations", + "metadata": {}, + "id": "caf19a75fb14f280" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.helpers import create_slab, select_slab_termination\nfrom mat3ra.made.tools.helpers import create_interface_zsl_between_slabs\n\nfilm_slabs = [\n create_slab(film, miller_indices=FILM_MILLER_INDICES, termination_top=t, vacuum=0)\n for t in film_slab_terminations\n]\nsubstrate_slabs = [\n create_slab(substrate, miller_indices=SUBSTRATE_MILLER_INDICES, termination_top=t, vacuum=0, number_of_layers=4)\n for t in substrate_slab_terminations\n]\n\n# visualize(\n# [{\"material\": s, \"title\": str(t)} for s, t in zip(film_slabs, film_slab_terminations)],\n# repetitions=[3, 3, 1], rotation=\"-90x\",\n# )\n# visualize(\n# [{\"material\": s, \"title\": str(t)} for s, t in zip(substrate_slabs, substrate_slab_terminations)],\n# repetitions=[3, 3, 1], rotation=\"-90x\",\n# )", + "metadata": { + "trusted": true + }, + "id": "7bae67a0b5821aa7", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 4.3. Create Substrate and Film Slabs", + "metadata": {}, + "id": "a692280e024d5883" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.build.pristine_structures.two_dimensional.slab import SlabConfiguration, SlabBuilder\n\nsubstrate_slab_config = SlabConfiguration.from_parameters(\n material_or_dict=substrate,\n miller_indices=SUBSTRATE_MILLER_INDICES,\n number_of_layers=SUBSTRATE_THICKNESS,\n vacuum=0.0,\n termination_top_formula=SUBSTRATE_TERMINATION_FORMULA,\n use_conventional_cell=USE_CONVENTIONAL_CELL,\n)\nfilm_slab_config = SlabConfiguration.from_parameters(\n material_or_dict=film,\n miller_indices=FILM_MILLER_INDICES,\n number_of_layers=FILM_THICKNESS,\n vacuum=0.0,\n termination_bottom_formula=FILM_TERMINATION_FORMULA,\n use_conventional_cell=USE_CONVENTIONAL_CELL,\n)\n\nsubstrate_slab = SlabBuilder().get_material(substrate_slab_config)\nfilm_slab = SlabBuilder().get_material(film_slab_config)", + "metadata": { + "trusted": true + }, + "id": "89f15267bd2f6f66", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 5. Find Interfaces with ZSL Strain Matching\n\n### 5.1. Initialize ZSL Analyzer", + "metadata": {}, + "id": "4c4b81397d2eb68" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.analyze.interface import ZSLInterfaceAnalyzer\n\nzsl_analyzer = ZSLInterfaceAnalyzer(\n substrate_slab_configuration=substrate_slab_config,\n film_slab_configuration=film_slab_config,\n max_area=MAX_AREA,\n max_area_ratio_tol=MAX_AREA_TOLERANCE,\n max_length_tol=MAX_LENGTH_TOLERANCE,\n max_angle_tol=MAX_ANGLE_TOLERANCE,\n reduce_result_cell=False,\n)", + "metadata": { + "trusted": true + }, + "id": "b04a8543e8fcf455", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 5.2. Generate and Plot Matches", + "metadata": {}, + "id": "a1956430738a9a9a" + }, + { + "cell_type": "code", + "source": "from utils.plot import plot_strain_vs_area\n\nPLOT_SETTINGS = {\n \"HEIGHT\": 600,\n \"X_SCALE\": \"log\",\n \"Y_SCALE\": \"log\",\n}\n\nmatches = zsl_analyzer.zsl_match_holders\nprint(f\"Found {len(matches)} matches\")\n# plot_strain_vs_area(matches, PLOT_SETTINGS)", + "metadata": { + "trusted": true + }, + "id": "8682bbecf48aa30b", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 5.3. Select the Interface\n\nChoose the match index from the plot above (index 0 has the lowest strain).", + "metadata": {}, + "id": "245b45e7bb4a7ac3" + }, + { + "cell_type": "code", + "source": "selected_index = 0", + "metadata": { + "trusted": true + }, + "id": "a83a9d7a43391187", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 6. Create the Interface", + "metadata": {}, + "id": "5360346b47cf21ed" + }, + { + "cell_type": "code", + "source": "interface = create_interface_zsl_between_slabs(\n substrate_slab=substrate_slab,\n film_slab=film_slab,\n gap=INTERFACE_DISTANCE,\n vacuum=INTERFACE_VACUUM,\n match_id=selected_index,\n max_area=MAX_AREA,\n max_area_ratio_tol=MAX_AREA_TOLERANCE,\n max_length_tol=MAX_LENGTH_TOLERANCE,\n max_angle_tol=MAX_ANGLE_TOLERANCE,\n reduce_result_cell_to_primitive=REDUCE_RESULT_CELL_TO_PRIMITIVE,\n)", + "metadata": { + "trusted": true + }, + "id": "4d413413dfaa3f6d", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 6.1. Visualize Interface", + "metadata": {}, + "id": "24cc0a761f161676" + }, + { + "cell_type": "code", + "source": "from utils.visualize import ViewersEnum\n\n# visualize([{\"material\": interface, \"title\": interface.name}], viewer=ViewersEnum.wave)\n# visualize(interface, repetitions=[1, 1, 1], rotation=\"-90x\")", + "metadata": { + "trusted": true + }, + "id": "fcbb4e6c1de21233", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 7. Apply Relaxation\n### 7.1. Relax with MACE", + "metadata": {}, + "id": "e64688fc18c49bb6" + }, + { + "cell_type": "code", + "source": [ + "\n", + "import os\n", + "import plotly.graph_objs as go\n", + "from IPython.display import display\n", + "from plotly.subplots import make_subplots\n", + "\n", + "from mat3ra.made.tools.convert import to_ase\n", + "from ase.optimize import BFGS\n", + "\n", + "from mace.calculators import MACECalculator\n", + "\n", + "calculator = MACECalculator(model_path=MODEL_PATHS_MAP[\"large\"], dispersion=False, default_dtype=\"float32\", device=\"cpu\")\n", + "\n", + "ase_interface = to_ase(interface)\n", + "ase_interface.set_calculator(calculator)\n", + "dyn = BFGS(ase_interface)\n", + "\n", + "steps = []\n", + "energies = []\n", + "\n", + "fig = make_subplots(rows=1, cols=1, specs=[[{\"type\": \"scatter\"}]])\n", + "scatter = go.Scatter(x=[], y=[], mode=\"lines+markers\", name=\"Energy\")\n", + "fig.add_trace(scatter)\n", + "fig.update_layout(title_text=\"Real-time Optimization Progress\", xaxis_title=\"Step\", yaxis_title=\"Energy (eV)\")\n", + "\n", + "f = go.FigureWidget(fig)\n", + "display(f)\n", + "\n", + "\n", + "def plotly_callback():\n", + " step = dyn.nsteps\n", + " energy = ase_interface.get_total_energy()\n", + " steps.append(step)\n", + " energies.append(energy)\n", + " print(f\"Step: {step}, Energy: {energy:.4f} eV\")\n", + " with f.batch_update():\n", + " f.data[0].x = steps\n", + " f.data[0].y = energies\n", + "\n", + "\n", + "dyn.attach(plotly_callback, interval=1)\n", + "dyn.run(fmax=RELAXATION_PARAMETERS[\"FMAX\"])\n", + "\n", + "ase_original_interface = to_ase(interface)\n", + "ase_original_interface.set_calculator(calculator)\n", + "ase_final_interface = ase_interface\n", + "\n", + "original_energy = ase_original_interface.get_total_energy()\n", + "relaxed_energy = ase_interface.get_total_energy()\n", + "\n", + "# print(\"Original structure:\\n\", ase_to_poscar(ase_original_interface))\n", + "# print(\"\\nRelaxed structure:\\n\", ase_to_poscar(ase_final_interface))\n", + "print(f\"The final energy is {float(relaxed_energy):.3f} eV.\")" + ], + "metadata": { + "jupyter": { + "is_executing": true + }, + "trusted": true + }, + "id": "3d8746a77f71bab5", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 7.2. View Structure Before and After Relaxation", + "metadata": {}, + "id": "abfa372909a96bf8" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.convert import from_ase\n\n\ndef atoms_to_material(atoms, title):\n material = Material.create(from_ase(atoms))\n material.name = title\n return material\n\n\nmaterial_original = atoms_to_material(ase_original_interface, f\"Original E={original_energy:.3f} eV\")\nmaterial_relaxed = atoms_to_material(ase_final_interface, f\"Relaxed E={relaxed_energy:.3f} eV\")\n\nvisualize(\n [\n {\"material\": material_original, \"title\": material_original.name},\n {\"material\": material_relaxed, \"title\": material_relaxed.name},\n ],\n viewer=ViewersEnum.wave,\n)", + "metadata": { + "trusted": true + }, + "id": "9565d0931b198f63", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## 7.4. Output interlayer distance before and after relaxation", + "metadata": {}, + "id": "e4b49774283e5517" + }, + { + "cell_type": "code", + "source": "from mat3ra.made.tools.analyze.other import get_average_interlayer_distance\n\nprint(f\"Interlayer distance before relaxation: {get_average_interlayer_distance(material_original, 0, 1):.4f} Å\")\nprint(f\"Interlayer distance after relaxation: {get_average_interlayer_distance(material_relaxed, 0, 1):.4f} Å\")", + "metadata": { + "trusted": true + }, + "id": "6dd00402bc2e9d59", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "### 7.4. Calculate Interface Energy", + "metadata": {}, + "id": "7ee3d26311a30687" + }, + { + "cell_type": "code", + "source": "def filter_atoms_by_tag(atoms, material_index):\n return atoms[atoms.get_tags() == material_index]\n\n\ndef calculate_energy(atoms, calc):\n atoms.set_calculator(calc)\n return atoms.get_total_energy()\n\n\ndef calculate_delta_energy(total_energy, *component_energies):\n return total_energy - sum(component_energies)\n\n\nsubstrate_original = filter_atoms_by_tag(ase_original_interface, SUBSTRATE_INDEX)\nlayer_original = filter_atoms_by_tag(ase_original_interface, FILM_INDEX)\nsubstrate_relaxed = filter_atoms_by_tag(ase_final_interface, SUBSTRATE_INDEX)\nlayer_relaxed = filter_atoms_by_tag(ase_final_interface, FILM_INDEX)\n\noriginal_substrate_energy = calculate_energy(substrate_original, calculator)\noriginal_layer_energy = calculate_energy(layer_original, calculator)\nrelaxed_substrate_energy = calculate_energy(substrate_relaxed, calculator)\nrelaxed_layer_energy = calculate_energy(layer_relaxed, calculator)\n\ndelta_original = calculate_delta_energy(original_energy, original_substrate_energy, original_layer_energy)\ndelta_relaxed = calculate_delta_energy(relaxed_energy, relaxed_substrate_energy, relaxed_layer_energy)\n\narea = ase_original_interface.get_volume() / ase_original_interface.cell[2, 2]\nn_interface = ase_final_interface.get_global_number_of_atoms()\nn_substrate = substrate_relaxed.get_global_number_of_atoms()\nn_layer = layer_relaxed.get_global_number_of_atoms()\neffective_delta_relaxed = (\n relaxed_energy / n_interface\n - (relaxed_substrate_energy / n_substrate + relaxed_layer_energy / n_layer)\n ) / (2 * area)\n\nprint(f\"Original Substrate energy: {original_substrate_energy:.4f} eV\")\nprint(f\"Relaxed Substrate energy: {relaxed_substrate_energy:.4f} eV\")\nprint(f\"Original Layer energy: {original_layer_energy:.4f} eV\")\nprint(f\"Relaxed Layer energy: {relaxed_layer_energy:.4f} eV\")\nprint(\"\\nDelta between interface energy and sum of component energies\")\nprint(f\"Original Delta: {delta_original:.4f} eV\")\nprint(f\"Relaxed Delta: {delta_relaxed:.4f} eV\")\nprint(f\"Original Delta per area: {delta_original / area:.4f} eV/Ang^2\")\nprint(f\"Relaxed Delta per area: {delta_relaxed / area:.4f} eV/Ang^2\")\nprint(f\"Relaxed interface energy: {relaxed_energy:.4f} eV\")\nprint(\n f\"Effective relaxed Delta per area: {effective_delta_relaxed:.4f} eV/Ang^2 ({effective_delta_relaxed / 0.16:.4f} J/m^2)\")", + "metadata": { + "trusted": true + }, + "id": "79ea902feda4a8d3", + "outputs": [], + "execution_count": null + }, + { + "cell_type": "markdown", + "source": "## References\n\n[1] mat3ra-made interface builder: https://github.com/Exabyte-io/made \n[2] MACE-MP-0 foundation model: https://github.com/ACEsuit/mace?tab=readme-ov-file#foundation-models ", + "metadata": {}, + "id": "2f60fdb73e44c09c" + } + ] +} diff --git a/other/experiments/jupyterlite/uploads/.gitkeep b/other/experiments/jupyterlite/uploads/.gitkeep new file mode 100644 index 00000000..e69de29b