{ "cells": [ { "cell_type": "markdown", "source": [ "# Tutorial using recipes\n", "\n", "Recipes are scripted functions that cover probably 90% of cases. They read and save all information to disk using the HDF5 format and have sensible defaults. If you need more flexibility for your calculation, you can run in full control mode instead. \n", "\n", "In a production run, each one of these cells would most likely best be executed as a job. In that case, you would want to separate the computationally intensive runs from the analysis.\n", "\n", "You will need to install `pyscf`, `ase`, `matplotlib`, `h5py`, `pandas` and `pyqmc` to use these routines. These all can be installed using `pip`. " ], "metadata": {} }, { "cell_type": "code", "execution_count": 1, "source": [ "import os\n", "os.environ[\"MKL_NUM_THREADS\"] = \"1\"\n", "os.environ[\"NUMEXPR_NUM_THREADS\"] = \"1\"\n", "os.environ[\"OMP_NUM_THREADS\"] = \"1\"\n", "import numpy as np\n", "import pandas as pd\n", "from pyscf import gto, scf\n", "import pyqmc.api as pyq\n", "import h5py\n", "import matplotlib.pyplot as plt\n", "from ase.data.pubchem import pubchem_atoms_search\n", "import pyscf.pbc.tools.pyscf_ase as pyscf_ase\n", "# Here we reset the files \n", "for fname in ['mf.hdf5','optimized_wf.hdf5','vmc_data.hdf5','dmc.hdf5']:\n", " if os.path.isfile(fname):\n", " os.remove(fname) " ], "outputs": [ { "output_type": "stream", "name": "stderr", "text": [ "/opt/homebrew/Caskroom/miniforge/base/envs/python3.9/lib/python3.9/site-packages/pyscf/lib/misc.py:46: H5pyDeprecationWarning: Using default_file_mode other than 'r' is deprecated. Pass the mode to h5py.File() instead.\n", " h5py.get_config().default_file_mode = 'a'\n" ] } ], "metadata": {} }, { "cell_type": "markdown", "source": [ "This function computes the mean-field solution and saves the results to the file specified. We recommend using the ccecp pseudopotentials for high accuracy and efficiency.\n", "\n", "Here we search up the atomic positions of water to reduce errors. " ], "metadata": {} }, { "cell_type": "code", "execution_count": 2, "source": [ "def mean_field(chkfile):\n", " atoms = pubchem_atoms_search(name='water')\n", " mol = gto.M(\n", " atom=pyscf_ase.ase_atoms_to_pyscf(atoms),\n", " basis=\"ccECP_cc-pVDZ\", ecp=\"ccecp\",\n", " )\n", " mf = scf.RHF(mol)\n", " mf.chkfile = chkfile\n", " mf.kernel()\n", "mean_field(\"mf.hdf5\")" ], "outputs": [ { "output_type": "stream", "name": "stdout", "text": [ "converged SCF energy = -16.9321257267691\n" ] } ], "metadata": { "tags": [] } }, { "cell_type": "markdown", "source": [ "Now we wish to construct a Slater-Jastrow wave function and optimize its energy. This is done using the OPTIMIZE function in pyq. It's often helpful to do the first optimization with only a few configurations, to get close to the minimum cheaply.\n", "\n", "By default, `OPTIMIZE` uses 10 VMC blocks of 10 steps each, so this is 10x10x100=10,000 VMC samples per iteration. \n", "\n", "To increase precision, it is usually best to increase the number of configurations. " ], "metadata": {} }, { "cell_type": "code", "execution_count": 3, "source": [ "pyq.OPTIMIZE(\"mf.hdf5\", #Construct a Slater-Jastrow wave function from the pyscf output\n", " \"optimized_wf.hdf5\", #Store optimized parameters in this file.\n", " nconfig=100, #Optimize using this many Monte Carlo samples/configurations\n", " max_iterations=10, #10 optimization steps\n", " verbose=False)" ], "outputs": [], "metadata": { "tags": [] } }, { "cell_type": "markdown", "source": [ "The API contains functions such as `read_opt` which assist in reading the output files. Here we print out the energy and uncertainty as the optimization proceeds. You should see that the energy is decreasing (within uncertainty). " ], "metadata": {} }, { "cell_type": "code", "execution_count": 4, "source": [ "df = pyq.read_opt(\"optimized_wf.hdf5\")\n", "df" ], "outputs": [ { "output_type": "execute_result", "data": { "text/plain": [ " energy iteration error fname\n", "0 -16.991545 0 0.017707 optimized_wf.hdf5\n", "1 -17.080173 1 0.018693 optimized_wf.hdf5\n", "2 -17.126727 2 0.010658 optimized_wf.hdf5\n", "3 -17.116194 3 0.009721 optimized_wf.hdf5\n", "4 -17.109356 4 0.011984 optimized_wf.hdf5\n", "5 -17.123817 5 0.010352 optimized_wf.hdf5\n", "6 -17.150079 6 0.011616 optimized_wf.hdf5\n", "7 -17.148353 7 0.010380 optimized_wf.hdf5\n", "8 -17.158911 8 0.009709 optimized_wf.hdf5\n", "9 -17.170465 9 0.008196 optimized_wf.hdf5" ], "text/html": [ "
| \n", " | energy | \n", "iteration | \n", "error | \n", "fname | \n", "
|---|---|---|---|---|
| 0 | \n", "-16.991545 | \n", "0 | \n", "0.017707 | \n", "optimized_wf.hdf5 | \n", "
| 1 | \n", "-17.080173 | \n", "1 | \n", "0.018693 | \n", "optimized_wf.hdf5 | \n", "
| 2 | \n", "-17.126727 | \n", "2 | \n", "0.010658 | \n", "optimized_wf.hdf5 | \n", "
| 3 | \n", "-17.116194 | \n", "3 | \n", "0.009721 | \n", "optimized_wf.hdf5 | \n", "
| 4 | \n", "-17.109356 | \n", "4 | \n", "0.011984 | \n", "optimized_wf.hdf5 | \n", "
| 5 | \n", "-17.123817 | \n", "5 | \n", "0.010352 | \n", "optimized_wf.hdf5 | \n", "
| 6 | \n", "-17.150079 | \n", "6 | \n", "0.011616 | \n", "optimized_wf.hdf5 | \n", "
| 7 | \n", "-17.148353 | \n", "7 | \n", "0.010380 | \n", "optimized_wf.hdf5 | \n", "
| 8 | \n", "-17.158911 | \n", "8 | \n", "0.009709 | \n", "optimized_wf.hdf5 | \n", "
| 9 | \n", "-17.170465 | \n", "9 | \n", "0.008196 | \n", "optimized_wf.hdf5 | \n", "