Job arrays and parameter sweeps

The trouble with a qsub loop

With several job directories, the first idea is to submit them one by one in a loop.

for d in */; do
    (cd "$d" && qsub ../run.pbs)
done

For ten jobs that is fine. For hundreds you end up with hundreds of job IDs to track, extra load on the scheduler, and on many clusters a per-user limit on queued jobs. A job array is a single job that runs the same script many times with a different index. You submit once, get one job ID, and each run (a subjob) receives its index in an environment variable.

PBS Pro job arrays

#PBS -J 1-6 creates six subjobs, each with $PBS_ARRAY_INDEX set to a value from 1 to 6. The resource request (select, walltime) applies to one subjob.

#!/bin/bash
#PBS -N dft-grid
#PBS -l select=1:ncpus=4:mpiprocs=4:mem=8gb
#PBS -l walltime=04:00:00
#PBS -j oe
#PBS -J 1-6

cd "$PBS_O_WORKDIR" || exit 1
module load openmpi/4.1.6
ORCA="$HOME/opt/orca_6_1_0/orca"

NAME=$(sed -n "${PBS_ARRAY_INDEX}p" list.txt)
cd "$NAME" || exit 1
"$ORCA" "$NAME.inp" > "$NAME.out" 2>&1

Array job IDs carry brackets, like 1234[].server. qstat -u $USER shows the whole array as one line; add -t to see the individual subjobs. Each subjob writes its own log.

qstat -t -u $USER          # include subjobs
qdel "1234[]"               # cancel the whole array
qdel "1234[3]"              # cancel subjob 3 only

Brackets are glob characters in the shell, so quote them. PBS Pro requires at least two subjobs in the range, so -J 1-1 is rejected. A step is allowed, as in -J 1-20:2.

Slurm job arrays

In Slurm it is --array=1-20, and the index arrives in $SLURM_ARRAY_TASK_ID. Appending %5 limits the array to five subjobs running at once, which is good manners for large arrays on a shared cluster. In the log name, %A is the array job ID and %a the subjob index.

#!/bin/bash
#SBATCH --job-name=dft-grid
#SBATCH --nodes=1
#SBATCH --ntasks=4
#SBATCH --mem=8G
#SBATCH --time=04:00:00
#SBATCH --array=1-6%3
#SBATCH --output=dft-grid.%A_%a.log

module load openmpi/4.1.6
ORCA="$HOME/opt/orca_6_1_0/orca"

mapfile -t names < list.txt
NAME=${names[$((SLURM_ARRAY_TASK_ID - 1))]}

cd "$NAME" || exit 1
"$ORCA" "$NAME.inp" > "$NAME.out" 2>&1

squeue lists the pending range, such as 1234_[4-6%3], separately from running subjobs such as 1234_1. scancel 1234 cancels everything, scancel 1234_5 just one subjob. Whether PBS Pro offers a similar limit on concurrent subjobs depends on the version and site configuration; check your cluster's documentation.

Turning an index into an input

A subjob only gets a number. The simplest way to turn it into a job name is a list file with one name per line. The two scripts above use two ways of reading it.

$ cat list.txt
B3LYP_def2-SVP
B3LYP_def2-TZVP
PBE0_def2-SVP
PBE0_def2-TZVP
TPSS_def2-SVP
TPSS_def2-TZVP
$ IDX=4; sed -n "${IDX}p" list.txt
PBE0_def2-TZVP
$ mapfile -t names < list.txt; echo "${#names[@]} jobs, #4 is ${names[3]}"
6 jobs, #4 is PBE0_def2-TZVP

sed -n "4p" prints only line 4. mapfile -t reads the file into a bash array, one line per element; arrays count from 0, hence the minus one. The array range must match the number of lines in the list. Passing -J or --array on the command line saves editing the script, since command-line options override the directives inside it.

N=$(wc -l < list.txt)
qsub -J 1-"$N" array.pbs            # PBS Pro
sbatch --array=1-"$N"%5 array.sh     # Slurm

Generating inputs for a parameter sweep

Build water optimisation inputs for every combination of three functionals and two basis sets. Each combination gets its own directory, and the list file is written along the way. The input files are written with a here-doc (<<EOF … EOF).

#!/bin/bash
# make_inputs.sh: one directory and ORCA input per functional/basis pair
set -euo pipefail

functionals=(B3LYP PBE0 TPSS)
basis_sets=(def2-SVP def2-TZVP)
coords=$(tail -n +3 water.xyz)

: > list.txt
for f in "${functionals[@]}"; do
    for b in "${basis_sets[@]}"; do
        name="${f}_${b}"
        mkdir -p "$name"
        cat > "$name/$name.inp" <<EOF
! $f D4 $b Opt
%pal nprocs 4 end
%maxcore 1500
* xyz 0 1
$coords
*
EOF
        echo "$name" >> list.txt
    done
done

$ bash make_inputs.sh
$ ls
B3LYP_def2-SVP
B3LYP_def2-TZVP
PBE0_def2-SVP
PBE0_def2-TZVP
TPSS_def2-SVP
TPSS_def2-TZVP
list.txt
make_inputs.sh
water.xyz
$ cat PBE0_def2-TZVP/PBE0_def2-TZVP.inp
! PBE0 D4 def2-TZVP Opt
%pal nprocs 4 end
%maxcore 1500
* xyz 0 1
O      0.000000    0.000000    0.119262
H      0.000000    0.763239   -0.477047
H      0.000000   -0.763239   -0.477047
*

: > list.txt truncates the list file, so rerunning the script does not pile up duplicate names. The coordinates are read from the xyz file and pasted in, so each input is self-contained and still works when only the input is copied to scratch.

A LAMMPS temperature series

LAMMPS inputs use variables like ${T}, and bash substitutes the same syntax. Quoting the here-doc delimiter, 'EOF', tells bash to leave the body alone.

cat > in.nvt <<'EOF'
units           lj
atom_style      atomic
lattice         fcc 0.8442
region          box block 0 10 0 10 0 10
create_box      1 box
create_atoms    1 box
mass            1 1.0
velocity        all create ${T} 87287 loop geom
pair_style      lj/cut 2.5
pair_coeff      1 1 1.0 1.0 2.5
fix             1 all nvt temp ${T} ${T} 0.5
thermo          100
run             5000
EOF

Without the quotes bash replaces ${T} with its own variable, which is empty if unset:

fix             1 all nvt temp   0.5

Pass the temperature with LAMMPS's -var option. Keep one input file and change only the directory and temperature per subjob.

# temps.txt: 0.8, 1.0, 1.2, 1.4, one per line
T=$(sed -n "${SLURM_ARRAY_TASK_ID}p" temps.txt)
mkdir -p "T$T" && cd "T$T" || exit 1
lmp -var T "$T" -in ../in.nvt -log log.lammps

Collecting results

With one directory per calculation, the collection script stays simple. Looping over directories, checking for normal termination and gathering energies into a CSV is shown in chapter 06.

Resubmitting only the failures

If a few subjobs did not finish because of walltime or convergence problems, pick the failed names out of the list, write them to a new list and resubmit with that. A missing output file (a subjob that never started) also counts as a failure.

#!/bin/bash
# find_failed.sh: names from list.txt whose ORCA output is missing or unfinished
while read -r name; do
    out="$name/$name.out"
    if [[ ! -f $out ]] || ! grep -q 'ORCA TERMINATED NORMALLY' "$out"; then
        echo "$name"
    fi
done < list.txt > failed.txt

$ bash find_failed.sh
$ cat failed.txt
PBE0_def2-SVP
TPSS_def2-SVP
$ grep -L "ORCA TERMINATED NORMALLY" */*.out
PBE0_def2-SVP/PBE0_def2-SVP.out

Here PBE0_def2-SVP has a truncated output and TPSS_def2-SVP has none. grep -L (files without a match) only finds the first, because grep cannot look inside a file that does not exist.

To resubmit, let the array script take the list file name from a variable, falling back to list.txt.

LIST=${LIST:-list.txt}
IDX=${PBS_ARRAY_INDEX:-1}
NAME=$(sed -n "${IDX}p" "$LIST")
[[ -n $NAME ]] || { echo "no line $IDX in $LIST"; exit 1; }

N=$(wc -l < failed.txt)
qsub -J 1-"$N" -v LIST=failed.txt array.pbs              # PBS Pro, N >= 2
sbatch --array=1-"$N" --export=ALL,LIST=failed.txt array.sh  # Slurm

If only one job failed, submit it to PBS Pro as a normal job without -J: qsub -v LIST=failed.txt array.pbs. PBS_ARRAY_INDEX is then unset and IDX falls back to 1. A #PBS -J line left in the script would turn it back into an array, so it is easier to always give the range on the command line. Before rerunning, look at why the jobs failed: raise the walltime if it ran out, fix the input if the SCF did not converge, or the same failures come back.

Many small jobs inside one allocation

Calculations that finish in seconds, like xtb single points, spend longer in the queue than running if each gets its own job. Request several cores in one job and run them side by side with xargs -P or GNU parallel. Start with a script that handles one xyz file.

#!/bin/bash
# run_one.sh: xtb single point for one xyz file, in its own directory
set -euo pipefail
xyz=$(realpath "$1")
name=$(basename "$xyz" .xyz)
mkdir -p "sp/$name"
cd "sp/$name"
xtb "$xyz" > sp.out 2>&1
awk -v n="$name" '/TOTAL ENERGY/ {print n, $4}' sp.out

$ export OMP_NUM_THREADS=1
$ printf '%s\n' mols/*.xyz | xargs -P 4 -n 1 ./run_one.sh
ammonia -4.426093057513
methanol -8.225938944071
methane -4.175074573917
water -5.070222286727
$ parallel -j 4 --keep-order ./run_one.sh ::: mols/*.xyz
ammonia -4.426093057513
methane -4.175074573917
methanol -8.225938944071
water -5.070222286727

xargs -P 4 -n 1 takes one input line per command and runs up to four at a time. Output appears in completion order, which can change from run to run. GNU parallel's --keep-order prints in input order. Each calculation gets its own directory because xtb writes its result files under fixed names; running in one directory at the same time would overwrite them.