I am new to QE and particularly to ESM-RISM. I have read the de facto paper by S. Nishihara and M. Otani, Phys. Rev. B 96, 115429 (2017), as well as reviewed the pw.x documentation. The calculation seems to indeed work for several geometry steps (13, with a final gradient error of 1.6E-03), with the Laue-RISM in the SCF usually converging in 2-10 iterations. Then suddenly, it will use all 5000 rism3d_maxsteps, and the calculation will finish.

The system I am using is extremely simple; it is a lithium slab of 11 layers (11 atoms), which are centred at 0, as the paper suggests. Below, I will provide my ASE script to generate the said input file:

from ase.build import bcc100
from ase.io import write
from pathlib import Path
import numpy as np
from ase.units import Bohr
from ase.visualize import view

pseudopotentials = {'Li': 'li_pbe_v1.4.uspp.F.UPF'}
slab = bcc100('Li', size=(1,1,11), a=3.44, vacuum=7.5, orthogonal=True)
average_z_top_two_layers = (np.mean(np.sort(np.unique(slab.get_positions()[:, 2]))[::-1][:2])) / Bohr
average_z_bottom_two_layers = (np.mean(np.sort(np.unique(slab.get_positions()[:, 2]))[:2])) / Bohr
input_data = {
    'control': {
        'calculation': 'relax',
        'verbosity': 'high',
        'restart_mode': 'from_scratch',
        'nstep': 999,
        'tstress': False,
        'tprnfor': True,
        'outdir': './',
        'prefix': 'pw.dir',
        'max_seconds': 86100,
        'etot_conv_thr': 1.0e-5,
        'forc_conv_thr': 1.0E-4,
        'disk_io': 'low',
        'pseudo_dir': './',
        'trism': True,
    'system': {
        'ibrav': 0,
        'nbnd': int(len(slab) * 3 + (len(slab) * 0.3)),
        'tot_charge': 0.0,
        'tot_magnetization': -10000,
        'ecutwfc': 40.0,
        'ecutrho': 40 * 5,
        'occupations': 'smearing',
        'degauss': 0.00735,
        'smearing': 'cold',
        'input_dft': 'pbe',
        'nspin': 1,
        'assume_isolated': 'esm',
        'esm_bc': 'bc1',
        'vdw_corr': 'none',
    'electrons': {
        'electron_maxstep': 999,
        'scf_must_converge': True,
        'conv_thr': 1.0e-9,
        'mixing_mode': 'local-TF',
        'mixing_beta': 0.20,
        'mixing_ndim': 8,
        'diagonalization': 'rmm-davidson',
        'diago_david_ndim': 2,
        'diago_rmm_ndim': 4,
        'diago_rmm_conv': False,
        'diago_full_acc': False,
    'ions': {
        'ion_dynamics': 'bfgs',
        'upscale': 100,
        'bfgs_ndim': 6,
    'RISM': {
        'nsolv': 3,
        'closure': 'kh',
        'tempv': 298.15,
        'ecutsolv': 40 * 5,
        'solute_lj(1)': 'none',
        'solute_epsilon(1)': 0.01828,
        'solute_sigma(1)': 2.70,
        'starting1d': 'zero',
        'starting3d': 'zero',
        'rism1d_maxstep': 50000,
        'rism3d_maxstep': 5000,
        'rism1d_conv_thr': 1e-08,
        'rism3d_conv_thr': 1e-05,
        'mdiis1d_size': 20,
        'mdiis3d_size': 10,
        'mdiis1d_step': 0.1,
        'mdiis3d_step': 0.8,
        'rism3d_conv_level': 0.5,
        'rism3d_planar_average': True,
        'laue_expand_right': slab.cell[2, 2] * 2 / Bohr,
        'laue_starting_right': average_z_top_two_layers,
write('Li-tester.in', slab, format='espresso-in', input_data=input_data, kpts=[12,12,1], pseudopotentials=pseudopotentials)
view(slab, viewer='x3d')

Solvents (mol/L):

Li+: 1.0 Li+.oplsaa.MOL
PF6-: 1.0 PF6-.oplsaa.MOL
EC: -1.0 EC.oplsua.MOL

These must be manually added. Note that the files for the solvent and electrolyte were obtained from Nishihara's GitHub (https://github.com/nisihara1/MOLs), and the Lennard-Jones parameters for the solvent are obtained from: Jun Haruyama, Tamio Ikeshoji, and Minoru Otani, The Journal of Physical Chemistry C 2018, 122 (18), 9804-9810, DOI: 10.1021/acs.jpcc.8b01979 - supplementary information.

Simply put, I wish to solvate a lithium slab in EC, with Li+ & PF6- electrolyte, then explore GC-DFT in QE from this.

Maybe I misunderstood some keywords, theory, or even obtaining the LJ parameters from the supplementary information. Any and all help would be most appreciated!

Note that the mdiis1d_step was changed from default as this improved 1D-RISM convergence.

Kind regards, Brad



You must log in to answer this question.

Browse other questions tagged .