Examples¶
We look at the examples in the ‘examples’ folder of the Github repository. For all scripts, the user must configure the parameters (mostly directories) at the top of the script under the “User input” section.
Converting DFT data to Zacros input files¶
This script reads VASP data and computes parameters for the input files.
# -*- coding: utf-8 -*-
"""
Created on Fri Apr 07 09:39:17 2017
@author: wittregr
"""
import DFT_to_Thermochemistry as _thermo
import os
import platform
if platform.system() == 'Windows':
Base_path = 'C:\Users\wittregr\Documents\Python Scripts'
else:
Base_path = '/home/1729'
Input = 'Input'
Output = 'Output'
'''
Read reference species and calculate all thermodynamic quantities
'''
#filepath = os.path.join(Base_path, Input, 'Reference_set_info.txt')
filepath = os.path.join(Base_path, Input, 'Zacros_Species_Energy.txt')
[lines, dict] = _thermo.DFTFileRead(filepath)
T_ref = []
Path = os.path.join(Base_path, Input)
for s in lines[3:]:
T_ref.append(_thermo.Reference(s.split('\t'), dict, Path))
'''
Read target species and calculate all thermodynamic quantities
'''
filepath = os.path.join(Base_path, 'Input', 'Tobe_Referenced.txt')
[lines, dict] = _thermo.DFTFileRead(filepath)
T_target = []
for s in lines[3:]:
T_target.append(_thermo.Target(s.split('\t'), dict, Base_path))
'''
Create basis set from reference molecules
'''
Basis = _thermo.Reference.BasisSet(T_ref)
'''
Read surface slab energy from input file
'''
filepath = os.path.join(Base_path, Input, 'Surface_energy_data.txt')
[lines, dict] = _thermo.DFTFileRead(filepath)
T_surface = []
for s in lines[3:]:
T_surface.append(_thermo.Surface(s.split('\t'), dict))
'''
Apply reference molecule correct to thermodynamic data
'''
T_target = _thermo.Target.ReferenceDFT(T_target, T_surface, Basis)
'''
Output thermat file from calculated thermodynamic data
'''
T_target = _thermo.Target.CreateThermdat(T_target, Base_path, Output)
'''
Output Zacros input file from thermodynamic data
'''
pass
'''
Test output (To be removed in final version)
'''
print '\n\nHeat Capacities [kcal/mol K] over a range of Temperatures [k]'
print '-------------------------------------------------------------\n'
print 'Name ',
for y in range(0, len(T_target[0].Cp_Range)):
print ' %4d ' % (T_target[0].Cp_Range[y]),
print '\n',
print '---------------- ',
for y in range(0, len(T_target[0].Cp_Range)):
print ' %s ' % ('----'),
print '\n',
for x in range(0, len(T_target)):
print '%16s ' % (T_target[x].name),
for y in range(0, len(T_target[x].Cp)):
print ' %5.2f' % (T_target[x].Cp[y]),
print '\n',
print '\n\nStandard Heats of Formation and Dispersion'
print '------------------------------------------\n'
print 'Name hf(%6.2f) edisp S(%6.2f)' % (T_target[0].Tstp,
T_target[0].Tstp)
print '---------------- ---------- ------- ---------'
for x in range(0, len(T_target)):
print '%-16s %7.3f %6.3f %6.3f' % (T_target[x].name,
T_target[x].hf_Tstp,
T_target[x].convedisp,
T_target[x].S_Tstp)
print '\n\n'
Analyzing a single trajectory¶
Single_trajectory.py analyzes a single Zacros trajectory and plots some data.
# Read the output for a Zacros trajectory and plot data
import zacros_wrapper as zw
''' ------------ User input section ------------ '''
RunPath = '/home/vlachos/wangyf/Alumina/mechanismI/5'
''' -------------------------------------------- '''
''' Read simulation results '''
my_trajectory = zw.kmc_traj() # Create single trajectory object
my_trajectory.Path = RunPath # Set directory for files
my_trajectory.ReadAllOutput(build_lattice=True) # Read input and output files
''' Plot data '''
my_trajectory.PlotSurfSpecVsTime() # Plot surface species populations versus time
my_trajectory.PlotGasSpecVsTime() # Plot gas species populations versus time
my_trajectory.PlotElemStepFreqs(time_norm = True) # Plot elementary step frequencies
my_trajectory.PlotLattice() # Plot the lattice
my_trajectory.LatticeMovie() # Plot the lattice snapshots
This produces several .png images with data within the folder with the Zacros files.
Running and analyzing multiple trajectories¶
Multiple_trajectories_MPI.py analyzes multiple trajectories, averages them, computes the rate, and performs sensitivity analysis. MPI is used to parallelize
running the trajectories. Therefore it must be launced with the mpiexec command.
# Read multiple trajectories and perform statistical analysis
import zacros_wrapper as zw
from mpi4py import MPI
''' ------------ User input section ------------ '''
zacros_exe = '/home/vlachos/mpnunez/bin/zacros_ZW.x'
KMC_source = '/home/vlachos/mpnunez/Github/Zacros-Wrapper/input_files/AtoB/'
BatchPath = '/home/vlachos/mpnunez/test'
Product = 'B'
n_runs = 16
''' -------------------------------------------- '''
if __name__ == '__main__':
COMM = MPI.COMM_WORLD
COMM.Barrier()
# Prepare template
x = zw.Replicates()
x.runtemplate = zw.kmc_traj()
x.runtemplate.Path = KMC_source
x.runtemplate.exe_file = zacros_exe
x.runtemplate.ReadAllInput()
x.ParentFolder = BatchPath
# Build files and run
x.n_runs = n_runs
if COMM.rank == 0:
x.BuildJobFiles()
'''
Run all trajectories using MPI parallelization
'''
# Collect whatever has to be done in a list. Here we'll just collect a list of
# numbers. Only the first rank has to do this.
if COMM.rank == 0:
jobs = x.run_dirs
jobs = [jobs[_i::COMM.size] for _i in range(COMM.size)] # Split into however many cores are available.
else:
jobs = None
jobs = COMM.scatter(jobs, root=0) # Scatter jobs across cores.
# Now each rank just does its jobs and collects everything in a results list.
# Make sure to not use super big objects in there as they will be pickled to be
# exchanged over MPI.
for job in jobs:
x.runtemplate.Path = job
x.runtemplate.Run_sim()
'''
Read and analyze results
'''
if COMM.rank == 0:
x.ReadMultipleRuns()
x.AverageRuns()
x.runAvg.PlotSurfSpecVsTime()
x.runAvg.PlotGasSpecVsTime()
x.runAvg.PlotElemStepFreqs(time_norm = True)
This produces several .png files with statistics from the average of the trajectories.
Multiple_trajectories.py is usable on squidward.che.udel.edu and uses multiple job submissions to parallelize the jobs, as MPI does not allow for jobs spread over multiple nodes (i.e. > 16 cores).
# Read multiple trajectories and perform statistical analysis
import zacros_wrapper as zw
''' ------------ User input section ------------ '''
zacros_exe = '/home/vlachos/mpnunez/bin/zacros_ZW.x'
KMC_source = '/home/vlachos/mpnunez/Github/Zacros-Wrapper/input_files/AtoB/'
BatchPath = '/home/vlachos/mpnunez/test'
Product = 'B'
n_runs = 16
''' -------------------------------------------- '''
# Prepare template
x = zw.Replicates()
x.runtemplate = zw.kmc_traj()
x.runtemplate.Path = KMC_source
x.runtemplate.exe_file = zacros_exe
x.runtemplate.ReadAllInput()
x.ParentFolder = BatchPath
# Build files and run
x.n_runs = n_runs
x.BuildJobFiles()
x.RunAllJobs_parallel_JobArray()
# Read results
x.ReadMultipleRuns()
# Perform sensitivity analysis
x.AverageRuns()
x.runAvg.PlotSurfSpecVsTime()
x.runAvg.PlotGasSpecVsTime()
x.runAvg.PlotElemStepFreqs(time_norm = True)
Rate constant rescaling and steady state detection¶
Scale_SS_Main.py takes Zacros input files and runs the statistical procedure described in our publication.
Parallelization is used to run multiple trajectories simultaneously. By default, it will use MPI parallelization, but
the flag parallel_mode can be set to use the gridengine submission scripts on Squidward or Farber.
# Read Zacros input files from a folder and then perform rate constant rescaling in a different folder
import zacros_wrapper as zw
''' ------------ User input section ------------ '''
exe_file = '/home/vlachos/mpnunez/bin/zacros_ZW.x'
KMC_source = '/home/vlachos/mpnunez/ZacrosWrapper/sample_systems/AtoB/stiff_input'
RunPath = '/home/vlachos/mpnunez/ZacrosWrapper/sample_systems/AtoB/ScaledownV3'
product_spec = 'B' # product species
number_of_runs = 96
''' -------------------------------------------- '''
if __name__ == '__main__': # Make commands safe for MPI parallelization
x = zw.kmc_traj(path = KMC_source)
x.gas_prod = product_spec
x.exe_file = exe_file
x.ReadAllInput()
final_data = zw.ReachSteadyStateAndRescale(x, RunPath, n_runs = 96,
n_batches = 1000, n_samples = 100, max_iterations = 10, parallel_mode = 'MPI')