15 - Workflows

This script demonstrates using the vasp-ase recipe system for automated calculation workflows, including multi-step relaxations and complex pipelines. python run.py

import numpy as np
from ase.build import bulk

from vasp import Vasp
from vasp.recipes.core import double_relax_flow, relax_job, static_job
from vasp.recipes.decorators import flow, job

print("=" * 60)
print("Automated Workflows with Recipes")
print("=" * 60)
print()

============================================================
Automated Workflows with Recipes
============================================================

Part 1: Recipe basics

print("Part 1: Recipe system overview")
print("-" * 40)
print()

print("The recipe system provides reusable calculation patterns:")
print()
print("  @job: Single VASP calculation")
print("    - static_job: Single-point energy")
print("    - relax_job: Geometry optimization")
print()
print("  @flow: Multi-step workflow")
print("    - double_relax_flow: Coarse → fine relaxation")
print("    - bulk_to_slabs_flow: Bulk → slab generation → calculation")
print()
print("  @subflow: Parallel sub-calculations")
print("    - Process multiple structures in parallel")
print()

Part 1: Recipe system overview
----------------------------------------

The recipe system provides reusable calculation patterns:

  @job: Single VASP calculation
    - static_job: Single-point energy
    - relax_job: Geometry optimization

  @flow: Multi-step workflow
    - double_relax_flow: Coarse → fine relaxation
    - bulk_to_slabs_flow: Bulk → slab generation → calculation

  @subflow: Parallel sub-calculations
    - Process multiple structures in parallel

Part 2: Static job workflow

print("Part 2: Static job workflow")
print("-" * 40)
print()

# Create silicon structure
si = bulk('Si', 'diamond', a=5.43)

# Run static calculation using built-in recipe
print("Running static_job on Si...")
result = static_job(
    si,
    label='results/workflows/si_static',
    kpts=(4, 4, 4),
    encut=300,
)
print(f"  Energy: {result.energy:.4f} eV")
print(f"  Converged: {result.converged}")
print()

Part 2: Static job workflow
----------------------------------------

Running static_job on Si...
  Energy: -10.8354 eV
  Converged: True

Part 3: Relaxation workflow

print("Part 3: Relaxation workflow")
print("-" * 40)
print()

# Slightly perturb Si lattice to show relaxation
si_perturbed = bulk('Si', 'diamond', a=5.50)  # Slightly expanded

print("Running relax_job on expanded Si (a=5.50 Å)...")
relax_result = relax_job(
    si_perturbed,
    label='results/workflows/si_relax',
    relax_cell=True,
    kpts=(4, 4, 4),
    encut=300,
    ediffg=-0.02,
    nsw=20,
)
print(f"  Initial a: 5.50 Å")
print(f"  Final energy: {relax_result.energy:.4f} eV")
a_final = relax_result.atoms.cell[0, 0] * np.sqrt(2)
print(f"  Final a: {a_final:.3f} Å")
print()

print("Running double_relax_flow (coarse → fine)...")
double_result = double_relax_flow(
    si_perturbed.copy(),
    label='results/workflows/si_double_relax',
    kpts=(4, 4, 4),
)
print(f"  Final energy: {double_result.energy:.4f} eV")
print()

Part 3: Relaxation workflow
----------------------------------------

Running relax_job on expanded Si (a=5.50 Å)...
  Initial a: 5.50 Å
  Final energy: -10.8417 eV
  Final a: 0.000 Å

Running double_relax_flow (coarse → fine)...
  Final energy: -10.8273 eV

Part 4: Multi-structure workflow

print("Part 4: Multi-structure workflow")
print("-" * 40)
print()

# Compare energies of different structures
structures = [
    ('Si', bulk('Si', 'diamond', a=5.43)),
    ('Ge', bulk('Ge', 'diamond', a=5.66)),
    ('C', bulk('C', 'diamond', a=3.57)),
]

print("Running static_job on multiple structures...")
print()

energies = {}
for name, atoms in structures:
    result = static_job(
        atoms,
        label=f'results/workflows/{name.lower()}_static',
        kpts=(4, 4, 4),
        encut=300,
    )
    e_per_atom = result.energy / len(atoms)
    energies[name] = e_per_atom
    print(f"  {name}: E = {e_per_atom:.4f} eV/atom")

print()
print("Cohesive energy comparison (relative to Si):")
for name, e in energies.items():
    delta = e - energies['Si']
    print(f"  {name}: ΔE = {delta:+.3f} eV/atom")
print()

Part 4: Multi-structure workflow
----------------------------------------

Running static_job on multiple structures...

  Si: E = -5.4177 eV/atom
  Ge: E = -4.4679 eV/atom
  C: E = -9.1649 eV/atom

Cohesive energy comparison (relative to Si):
  Si: ΔE = +0.000 eV/atom
  Ge: ΔE = +0.950 eV/atom
  C: ΔE = -3.747 eV/atom

Part 5: Custom recipe

print("Part 5: Creating custom recipes")
print("-" * 40)
print()

@job
def convergence_test_job(atoms, label='results/conv_test', encuts=None, **kwargs):
    """Test ENCUT convergence."""
    if encuts is None:
        encuts = [250, 300, 350, 400]

    results = []
    for encut in encuts:
        calc = Vasp(
            atoms=atoms.copy(),
            label=f'{label}_encut{encut}',
            encut=encut,
            **kwargs
        )
        energy = calc.potential_energy
        results.append({'encut': encut, 'energy': energy})

    return results

print("Custom recipe: convergence_test_job")
print()

# Run actual convergence test
al = bulk('Al', 'fcc', a=4.05)
test_results = convergence_test_job(
    al,
    label='results/workflows/al_conv',
    encuts=[250, 300, 350],
    kpts=(6, 6, 6),
    ismear=1,
    sigma=0.1,
)

print("  ENCUT convergence test for Al:")
for r in test_results:
    print(f"    ENCUT={r['encut']}: E={r['energy']:.4f} eV")

# Check convergence
if len(test_results) >= 2:
    delta = abs(test_results[-1]['energy'] - test_results[-2]['energy'])
    print(f"  Energy difference (last two): {delta*1000:.1f} meV")
print()

Part 5: Creating custom recipes
----------------------------------------

Custom recipe: convergence_test_job

  ENCUT convergence test for Al:
    ENCUT=250: E=-3.7285 eV
    ENCUT=300: E=-3.7336 eV
    ENCUT=350: E=-3.7363 eV
  Energy difference (last two): 2.7 meV

Part 6: Workflow composition

print("Part 6: Composing workflows")
print("-" * 40)
print()

@flow
def relax_then_static_flow(atoms, label, **kwargs):
    """Relax structure then run high-quality static calculation."""
    print("  Step 1: Relaxation...")
    relax = relax_job(
        atoms,
        label=f'{label}_relax',
        encut=300,
        ediffg=-0.02,
        nsw=20,
        **kwargs
    )

    print("  Step 2: High-quality static...")
    static = static_job(
        relax.atoms,
        label=f'{label}_static',
        encut=400,  # Higher cutoff for final energy
        **kwargs
    )

    return {
        'relaxed': relax,
        'static': static,
    }

print("Custom flow: relax_then_static_flow")
print()

# Run on copper
cu = bulk('Cu', 'fcc', a=3.65)  # Slightly off equilibrium
results = relax_then_static_flow(
    cu,
    label='results/workflows/cu_flow',
    kpts=(6, 6, 6),
    ismear=1,
    sigma=0.1,
)
print()
print(f"  Relaxed energy: {results['relaxed'].energy:.4f} eV")
print(f"  Static energy (high quality): {results['static'].energy:.4f} eV")
print(f"  Difference: {(results['static'].energy - results['relaxed'].energy)*1000:.1f} meV")
print()

Part 6: Composing workflows
----------------------------------------

Custom flow: relax_then_static_flow

  Step 1: Relaxation...
  Step 2: High-quality static...

  Relaxed energy: -3.7147 eV
  Static energy (high quality): -3.7159 eV
  Difference: -1.2 meV

Part 7: Workflow engine integration

print("Part 7: Workflow engine integration")
print("-" * 40)
print()

print("The recipe decorators support workflow engines:")
print()
print("  Environment variable: VASP_WORKFLOW_ENGINE")
print()
print("  Supported engines:")
print("    - prefect: Prefect workflows")
print("    - dask: Dask distributed")
print("    - parsl: Parsl parallel")
print("    - covalent: Covalent cloud")
print("    - jobflow: Materials Project jobflow")
print()
print("  When an engine is set:")
print("    - @job becomes the engine's task decorator")
print("    - @flow becomes the engine's flow decorator")
print("    - Enables parallel execution, retries, logging")
print()

print("Example with Prefect:")
print('''
  export VASP_WORKFLOW_ENGINE=prefect

  from vasp.recipes import static_job, relax_job

  # Now decorated with @prefect.task / @prefect.flow
  # Can use Prefect features:
  # - Automatic retries
  # - Task caching
  # - Distributed execution
  # - Web UI monitoring
''')
print()

Part 7: Workflow engine integration
----------------------------------------

The recipe decorators support workflow engines:

  Environment variable: VASP_WORKFLOW_ENGINE

  Supported engines:
    - prefect: Prefect workflows
    - dask: Dask distributed
    - parsl: Parsl parallel
    - covalent: Covalent cloud
    - jobflow: Materials Project jobflow

  When an engine is set:
    - @job becomes the engine's task decorator
    - @flow becomes the engine's flow decorator
    - Enables parallel execution, retries, logging

Example with Prefect:

  export VASP_WORKFLOW_ENGINE=prefect

  from vasp.recipes import static_job, relax_job

  # Now decorated with @prefect.task / @prefect.flow
  # Can use Prefect features:
  # - Automatic retries
  # - Task caching
  # - Distributed execution
  # - Web UI monitoring

Summary

print("=" * 60)
print("Summary")
print("=" * 60)
print()

print("Recipe system features:")
print("  - Reusable calculation patterns")
print("  - Consistent parameter handling")
print("  - Result dataclasses for type safety")
print("  - Workflow engine integration")
print("  - Easy composition of complex workflows")
print()

print("Available recipes:")
print("  Core: static_job, relax_job, double_relax_flow")
print("  Slabs: slab_static_job, slab_relax_job, bulk_to_slabs_flow")
print("  Phonons: phonon_job, phonon_flow (requires phonopy)")
print()

print("Best practices:")
print("  - Use recipes for reproducibility")
print("  - Set workflow engine for HPC")
print("  - Create custom @job and @flow for repeated tasks")
print("  - Compose simple jobs into complex workflows")
print()

print("Next: Try 16_neb/ for transition state calculations.")

============================================================
Summary
============================================================

Recipe system features:
  - Reusable calculation patterns
  - Consistent parameter handling
  - Result dataclasses for type safety
  - Workflow engine integration
  - Easy composition of complex workflows

Available recipes:
  Core: static_job, relax_job, double_relax_flow
  Slabs: slab_static_job, slab_relax_job, bulk_to_slabs_flow
  Phonons: phonon_job, phonon_flow (requires phonopy)

Best practices:
  - Use recipes for reproducibility
  - Set workflow engine for HPC
  - Create custom @job and @flow for repeated tasks
  - Compose simple jobs into complex workflows

Next: Try 16_neb/ for transition state calculations.