SCF files the Python Way

Hello All,

Week two of my personal post-a-week challenge. I have been practicing automating my workflow and I have another quick script.

Using my QE module that I talked about last week, I made a quick script that will develop the SCF, NSCF, DOS, and PDOS files. These are the files needed for a density of states calculation. What is nice about this work flow is it allows you, in one run script, to run all these calculations. This is important for running on Super Computers as wait times can sometimes outweigh the actual calculations.


#!/usr/bin/python3.4

import QE as qe
import sys
import numpy as np

class Files:
  def __init__(self):
    self.prefix = ''
    self.pdos = """&projwfc
    outdir='.'
    prefix='PREFIX'
    io_choice='both'
    Emin=-100, Emax=100.0, DeltaE=0.1
    ngauss=1, degauss=0.02
 /
"""
    self.dos = """&dos
  outdir='.'
  prefix='PREFIX'
  fildos='PREFIX.dos'
  Emin=-100.0, Emax=100.0, DeltaE=0.1
/
"""
    self.scf = """"""
    self.nscf = """"""
    self.round = 5
  def Process(self,filestring,QE):
    try:
      f = open(filestring,'r')
    except:
      print("Cannot open %s" % filestring)
      sys.exit(1)
    self.scf = ''
    self.nscf = ''
    for i in f:
      if 'calculation' in i:
        self.scf += i.replace('vc-relax','scf').replace('relax','scf')
        self.nscf += i.replace('vc-relax','nscf').replace('relax','nscf')
      elif '&system' in i:
        tmpstr = '''&system\n  ibrav = 14\n'''
        for j in ['a','b','c']:
          tmpstr += "  " + j + " = " + str(round(QE.norms[j],self.round)) + "\n"
        tmpstr += "  cosBC = " + str(round(np.cos(QE.angles['alpha'] * np.pi/180.),self.round)) + "\n"
        tmpstr += "  cosAC = " + str(round(np.cos(QE.angles['beta'] * np.pi/180.),self.round)) + "\n" 
        tmpstr += "  cosAB = " + str(round(np.cos(QE.angles['gamma'] * np.pi/180.),self.round)) + "\n"
        self.scf += tmpstr
        self.nscf += tmpstr
        for j in range(0,7):
          next(f)
      elif 'ATOMIC_POSITIONS' in i.upper():
        tmpstr = 'ATOMIC_POSITIONS {angstrom}\n'
        for j in QE.atoms:
          tmpstr += j[0] + "\t" + str(j[1]) + "\t" + str(j[2]) + "\t" + str(j[3]) + "\n"
          next(f)
        self.scf += tmpstr
        self.nscf += tmpstr
      elif 'K_POINTS' in i.upper():
        self.scf += i 
        self.nscf += i 
        line = next(f)
        self.scf += line
        ks = line.split()
        for j in range(0,6):
          if j < 3:
            self.nscf += str(2*int(ks[j])) + " " 
          else:
            self.nscf += ks[j] + " "
        self.nscf += "\n"
      elif "prefix" in i.lower():
        self.prefix=i.replace('prefix','').replace('=','').replace("'","").replace('"','').strip()
      else:
        self.scf += i
        self.nscf += i
      self.pdos = self.pdos.replace('PREFIX',self.prefix)
      self.dos = self.dos.replace('PREFIX',self.prefix)
        
def main(command):
  QEStruct= qe.Struct()
  QEStruct.File_Process(command[1])
  outputs = Files()
  outputs.Process(command[0],QEStruct)
  f = open('scf.in','w')
  f.write(outputs.scf)
  f.close()
  f = open('nscf.in','w')
  f.write(outputs.nscf)
  f.close()
  f = open('dos.in','w')
  f.write(outputs.dos)
  f.close()
  f = open('pdos.in','w')
  f.write(outputs.pdos)
  f.close()

if __name__ == "__main__":
  if len(sys.argv) != 3:
    print("Incorrect number of arguments, run as ./Scf.py QEINPUT QEOUTPUT")
    sys.exit(6)
  command = [sys.argv[1],sys.argv[2]]
  main(command)

This generates a scf.in, nscf.in, dos.in, pdos.in, which allows you to run one input script for all files as follows:

pw.x < relax.in > relax.out
./ScfCreator.py relax.in relax.out
pw.x < scf.in > scf.out
pw.x < nscf.in > nscf.out
dos.x < dos.in > dos.out
projwfc.x < pdos.in > pdos.out

Eventually I hope to develop this into a python wrapper for the QE program.

Happy computing,

Levi

6 comments Add yours
  1. Hello Levi,
    I am very new to Quantum Espresso. Recently, i have explored in INTERNET and found your blog. Its really a good place to learn few new tricks with QE. I tried to run the SCFcreator.py file to generate different input files for the calculation. But it shows some error as I have not added the QE to my module. Could you please let me know how to add QE to the module to import it as qe. Looking forward to hear from you.

    1. Hello Amrit,

      Sorry for the delay on responding to this, research gets crazy.

      You need to have my QE.py in your $PYTHONPATH. Once that is in there, you can use the ScfCreator.py file.

      Let me know if you have any other questions.

      Levi

Leave a Reply

Your email address will not be published. Required fields are marked *