Setting up VeloxChem on a HPC cluster

Objectives

  • Learn how to use the VeloxChem module on Dardel.

Keypoints

  • VeloxChem and its Python dependencies are installed in a virtual environment.

Dardel (PDC)

Dardel is a Cray XE system at PDC.

VeloxChem is available as a module on Dardel:

$ module load PDC/21.11
$ module load veloxchem/1.0dev

Example submission script:

#!/bin/bash

#SBATCH -J myjob
#SBATCH -t 01:00:00
#SBATCH -p main
#SBATCH -A edu22.veloxchem
#SBATCH --reservation velox-lab1

#SBATCH --nodes=1
#SBATCH --ntasks-per-node=8

module load PDC/21.11
module load veloxchem/1.0dev

export OMP_NUM_THREADS=16
export OMP_PLACES=cores

job=porphyrin
srun vlx ${job}.inp ${job}.out

Note

  • There are two ways to log in to Dardel.

  • The time allocation for the workshop is edu22.veloxchem

    • For the first day we can use the reservation velox-lab1

    • For the second day we can use the reservation velox-lab2

Exercise

The purpose of the exercise is to check that your VeloxChem can run a simple calculation. You can test with the cpp.inp file.

  • On Dardel:

    salloc -N 1 --ntasks-per-node 8 -t 00:20:00 -p main -A edu22.veloxchem --reservation velox-lab1
    module load PDC/21.11 veloxchem/1.0dev
    export OMP_NUM_THREADS=16
    export OMP_PLACES=cores
    wget https://gitlab.com/veloxchem/veloxchem/-/raw/master/docs/inputs/cpp.inp
    srun vlx cpp.inp
    exit
    

Script for analyzing singlet excitations from linear response output

from pathlib import Path
import numpy as np
import h5py
import sys

out_file = Path(sys.argv[1])
rsp_file = out_file.with_suffix('.rsp.solutions.h5')

# threshold for priting excitations

coef_thresh = 0.1

# read output file for number of occupied and virtual orbitals

nocc, n_ao = None, None
n_lindep = 0

with out_file.open('r') as f_out:
    for line in f_out:
        if 'Number of alpha electrons   :' in line:
            nocc = int(line.split(':')[1].split()[0])
        if 'Contracted Basis Functions :' in line:
            n_ao = int(line.split(':')[1].split()[0])
        if 'Removed' in line and 'linearly dependent vectors' in line:
            n_lindep = int(line.split('Removed')[1].split()[0])

n_mo = n_ao - n_lindep
nvir = n_mo - nocc
n_ov = nocc * nvir

# read hdf5 file for excitations

hf = h5py.File(str(rsp_file), 'r')
keys = list(hf.keys())

state = 0
while True:
    state += 1
    key = f'S{state}'
    if key not in keys:
        break
    array = np.array(hf.get(key))
    assert array.size == n_ov * 2

    excitations = []
    de_excitations = []

    for i in range(nocc):
        for a in range(nvir):
            ia = i * nvir + a
            exc_coef = array[ia]
            de_exc_coef = array[n_ov + ia]

            homo = f'HOMO-{nocc-1-i}' if i != nocc - 1 else 'HOMO'
            lumo = f'LUMO+{a}' if a != 0 else 'LUMO'

            if abs(exc_coef) > coef_thresh:
                excitations.append((
                    abs(exc_coef),
                    f'{homo:<8s} -> {lumo:<8s} {exc_coef:10.4f}',
                ))

            if abs(de_exc_coef) > coef_thresh:
                de_excitations.append((
                    abs(de_exc_coef),
                    f'{homo:<8s} <- {lumo:<8s} {de_exc_coef:10.4f}',
                ))

    print(f'S{state}:')
    for exc in sorted(excitations, reverse=True):
        print('  ', exc[1])
    for de_exc in sorted(de_excitations, reverse=True):
        print('  ', de_exc[1])

hf.close()