#!/usr/bin/env python
################################################################################
#
# CDDL HEADER START
#
# The contents of this file are subject to the terms of the Common Development
# and Distribution License Version 1.0 (the "License").
#
# You can obtain a copy of the license at
# http:# www.opensource.org/licenses/CDDL-1.0. See the License for the
# specific language governing permissions and limitations under the License.
#
# When distributing Covered Code, include this CDDL HEADER in each file and
# include the License file in a prominent location with the name LICENSE.CDDL.
# If applicable, add the following below this CDDL HEADER, with the fields
# enclosed by brackets "[]" replaced with your own identifying information:
#
# Portions Copyright (c) [yyyy] [name of copyright owner]. All rights reserved.
#
# CDDL HEADER END
#
# Copyright (c) 2017, Regents of the University of Minnesota.
# All rights reserved.
#
# Contributor(s):
# Ellad B. Tadmor
#
################################################################################
# The docstring below is vc_description
'''This verfication check ensures that a model supports all
of the species that it claims. The verification check attempts to compute the
energy for all clusters up to a preset size for all combinations of supported
elements. Several random samples are generated for each cluster.
The check is successful if all elements appear in at least one
successful calculation. This algorithm accounts for the fact that some models
may only support species in certain combinations. For example a model
supporting species A and B may be for an AB alloy and may only support A-B
interactions and not A-A or B-B. The output provides full information on
which interactions are supported and which are not.'''
# Python 2-3 compatible code issues
from __future__ import print_function
try:
input = raw_input
except NameError:
pass
from kimcalculator import KIMCalculator, KIM_get_supported_species_list
import kimvc
import itertools
from ase import Atoms
import numpy as np
import sys
import math
import random
__version__ = "000"
__author__ = "Ellad Tadmor"
vc_name = "vc-species-supported-as-stated"
vc_description = kimvc.vc_stripall(__doc__)
vc_category = "mandatory"
vc_grade_basis = "passfail"
vc_files = []
################################################################################
#
# FUNCTIONS
#
################################################################################
def pos_has_distance_less_than_dmin(pos,dmin):
for i in range(len(pos)-1):
for j in range(i+1,len(pos)):
dist = np.linalg.norm(pos[j]-pos[i])
if dist<dmin:
return True
return False
################################################################################
def get_cluster_positions(N,rball,dmin):
'''
Select 'N' positions within a sphere of radius 'rball' where no two points
are closer than 'rmin'.
'''
pos = np.zeros((N,3))
points_too_close = True
onethird = 1./3.
while points_too_close:
# create N points uniformly distributed within a sphere of radius 'rball'
for i in range(N):
phi = random.uniform(0,2*math.pi)
costheta = random.uniform(-1.,1.)
theta = np.arccos( costheta )
x = np.sin( theta) * np.cos( phi )
y = np.sin( theta) * np.sin( phi )
z = np.cos( theta )
rad = random.uniform(0,rball)
pos[i] = (x*rad, y*rad, z*rad)
points_too_close = pos_has_distance_less_than_dmin(pos,dmin)
return pos
################################################################################
def do_vc(model, vc):
'''
Perform Species Supported as Stated VC
'''
# Get supported species
species = KIM_get_supported_species_list(model)
species = kimvc.remove_species_not_supported_by_ASE(species)
species.sort()
# Initialize
Nmax = 4 # clusters ranging from 2 to Nmax atoms will be
# tested (as long as 'max_combinations' is not
# exeeded.
max_combinations = 10000 # maximum number of allowed combinations
rball = 1.5 # radius of ball containing clusters
dmin = 1.0 # minimum allowed distance between atoms
Nsamp = 10 # number of samples of each cluster
seed = 13
random.seed(seed)
species_is_supported = {}
for spec in species:
species_is_supported[spec] = False # Assume none of the species
# are supported and now test
yesno = {
True: "yes",
False: "** NO **"}
# Print VC Header
dashwidth = 80
vc.rwrite('')
vc.rwrite('-'*dashwidth)
vc.rwrite('Results for KIM Model : %s' % model.strip())
vc.rwrite('Supported species : %s' % ' '.join(species))
vc.rwrite('')
vc.rwrite('Maximum cluster size considered (#atoms) = %d' % Nmax)
vc.rwrite('Maximum number of allowed combinations = %d' % max_combinations)
vc.rwrite('Radius of ball containing clusters: = %0.3f' % rball)
vc.rwrite('Minimum allowed distance between atoms: = %0.3f' % dmin)
vc.rwrite('Number of random samples for each cluster = %d' % Nsamp)
vc.rwrite('Random number seed: = %d' % seed)
vc.rwrite('-'*dashwidth)
vc.rwrite('')
# Loop over all clusters from 2 to Nmax considering all possible
# species combinations without considering order.
# For example for a model supporting 3 species (A,B,C) the loop would
# be over 2-atom and 3-atom clusters with the following combinations:
# A-A, B-B, C-C, A-B, A-C, B-C
# A-A-A, B-B-B, C-C-C, A-A-B, A-A-C, B-B-A, B-B-C, C-C-A, C-C-B, A-B-C
# Test whether the model supports the configuration. The model passes
# the check if all species appear in at least one supported configuration.
# Report this combination
vc.rwrite("%-*s %7s %s" % (2*Nmax, 'Cluster', 'Success', 'Supported'))
vc.rwrite('-'*(2*Nmax + 24))
for N in range(2,Nmax+1):
# Create all combintations of species for a cluster of size N
species_cluster = list(itertools.combinations_with_replacement(species,N))
# Compute number of combinations and skip if too many
num_combinations = len(species_cluster)
if num_combinations > max_combinations:
vc.rwrite("...... skipping %d-cluster (%d combinations)" % \
(N, num_combinations))
continue
# Loop over combinations
for i in range(len(species_cluster)):
# Loop over samples
num_success = 0
for n in range(Nsamp):
# Assign positions to atoms in a cluster of N atoms in a ball of
# radius 'rball' with no atoms closer than 'dmin'.
pos = get_cluster_positions(N,rball,dmin)
# Define Atoms object for cluster
atoms = Atoms(''.join(species_cluster[i]), positions=pos,
cell=(2*rball, 2*rball, 2*rball),pbc=(0,0,0))
calc = KIMCalculator(model)
# Try to compute energy
try:
atoms.set_calculator(calc)
energy = atoms.get_potential_energy()
num_success += 1
for spec in species_cluster[i]:
species_is_supported[spec] = True
except:
pass
# Report this combination
vc.rwrite("%-*s %2d / %2d %s" % \
(2*Nmax, ''.join(species_cluster[i]), num_success, Nsamp, yesno[num_success>0]))
vc.rwrite('-'*(2*Nmax + 24))
vc.rwrite('')
vc.rwrite('='*dashwidth)
vc.rwrite('To pass this verification check each model element must be supported')
vc.rwrite('in at least one cluster combination.')
vc.rwrite('')
if all(value for value in species_is_supported.values()):
vc_grade = 'P'
vc_comment = 'All elements claimed by model are supported in at least one of the tested configrations.'
else:
vc_grade = 'F'
vc_comment = 'At least one of the elements claimed by the model is not supported in any tested configuration.'
vc.rwrite('Grade: {}'.format(vc_grade))
vc.rwrite('')
vc.rwrite('Comment: '+vc_comment)
vc.rwrite('')
return vc_grade, vc_comment
################################################################################
#
# MAIN PROGRAM
#
###############################################################################
if __name__ == '__main__':
# Get the model extended KIM ID:
model = input("Model Extended KIM ID = ")
# Define VC object and do verification check
vc = kimvc.vc_object(vc_name, vc_description, __author__)
with vc:
# Perform verification check and get grade
try:
vc_grade, vc_comment = do_vc(model, vc)
except:
vc_grade = "N/A"
vc_comment = "Unable to perform verification check due to an error."
sys.stderr.write('ERROR: Unable to perform verification check.\n')
# Pack results in a dictionary and write VC property instance
results = {"vc_name" : vc_name,
"vc_description" : vc_description,
"vc_category" : vc_category,
"vc_grade_basis" : vc_grade_basis,
"vc_grade" : vc_grade,
"vc_comment" : vc_comment,
"vc_files" : vc_files}
vc.write_results(results)