Source code for pcmpreprocess

import pcmgetkw as getkw
from copy import deepcopy

isAngstrom = False

[docs]def setup_keywords(): top=getkw.Section('toplevel', callback=verify_top) top.set_status(True) top.add_kw('Units', 'STR', 'AU') cavity=getkw.Section('Cavity', callback=verify_cavity) cavity.add_kw('Type','STR') cavity.add_kw('PatchLevel', 'INT', 2) cavity.add_kw('Coarsity', 'DBL', 0.5) cavity.add_kw('Area','DBL', 0.3) cavity.add_kw('Scaling', 'STR', 'Yes') cavity.add_kw('AddSpheres', 'STR', 'Yes') cavity.add_kw('Mode','STR','Explicit') cavity.add_kw('Atoms','INT_ARRAY') cavity.add_kw('Radii','DBL_ARRAY') cavity.add_kw('RadiiSet', 'STR', 'Bondi') cavity.add_kw('Spheres','DBL_ARRAY', callback=verify_spheres) top.add_sect(cavity) medium=getkw.Section('Medium', callback=verify_medium) medium.add_kw('Solvent', 'STR', 'Water') medium.add_kw('SolverType', 'STR', 'IEFPCM') medium.add_kw('EquationType','STR', 'SecondKind') medium.add_kw('Correction', 'DBL', 0.0) medium.add_kw('ProbeRadius', 'DBL', 1.0) top.add_sect(medium) green=getkw.Section('Green', callback=verify_green) green.add_kw('Type', 'STR', 'Vacuum') green.add_kw('Der', 'STR', 'Derivative') green.add_kw('Eps', 'DBL', 1.0) green.add_kw('EpsRe', 'DBL', 1.0) green.add_kw('EpsImg', 'DBL', 1.0) green.add_kw('SphereRadius', 'DBL', 1.0) green.add_kw('SpherePosition', 'DBL_ARRAY') medium.add_sect(green) green_part = deepcopy(green) green.add_sect(green_part) return top
[docs]def verify_top(section): global isAngstrom allowed_units = ('AU', 'Angstrom') key = section.get('Units') val = key.get() if (val not in allowed_units): print "Allowed units are: ", allowed_units sys.exit(1) if (val == 'Angstrom'): isAngstrom = True
[docs]def verify_cavity(section): allowed = ('GePol', 'Wavelet') type = section.get('Type') if (type.get() not in allowed): print "Allowed types are: ", allowed sys.exit(1) if section['Area'].is_set(): convert_area_scalar(section['Area']) if (type.get() == 'GePol'): area=section.get('Area') a=area.get() if (a < 0.01): print "Area value is too small" print "Minimum value: 0.01" sys.exit(1) elif (type.get() == 'Wavelet'): key = section.get('PatchLevel') if (key.get() < 1): print "Patch level must be > 0" sys.exit(1) key = section.get('Coarsity') if (key.get() < 0.0 or key.get() >= 1.0): print "Coarsity has to be within ]0,1[" sys.exit(1) allowed_modes = ("Explicit", "Atoms", "Implicit") mode = section.get('Mode') if (mode.get() not in allowed_modes): print "Allowed modes are: ", allowed_modes sys.exit(1) atoms=section.get('Atoms') at=atoms.get() radii = section.get('Radii') convert_length_array(radii); r = radii.get() if (mode.get() == 'Atoms'): if (len(r) != len(at) or len(at) == 0): print "Incoherent input for Atoms keyword." print "Check that Atoms and Radii are consistent." sys.exit(1) else: for i, v in enumerate(at): if (at.count(v) > 1): print "Incoherent input for Atoms keyword." print "Too many spheres on the same atom(s)." sys.exit(1)
[docs]def verify_medium(section): allowedSolvents = {'Water': ('Water', 'water', 'H2O'), 'Methanol': ('Methanol', 'methanol', 'CH3OH'), 'Ethanol': ('Ethanol', 'ethanol', 'CH3CH2OH'), 'Chloroform': ('Chloroform', 'chloroform', 'CHCl3'), 'Methylenechloride': ('Methylenechloride', 'methylenechloride', 'CH2Cl2'), '1,2-Dichloroethane': ('1,2-Dichloroethane', '1,2-dichloroethane', 'C2H4Cl2'), 'Carbon Tetrachloride': ('Carbon Tetrachloride', 'carbon tetrachloride', 'CCl4'), 'Benzene': ('Benzene', 'benzene', 'C6H6'), 'Toluene': ('Toluene', 'toluene', 'C6H5CH3'), 'Chlorobenzene': ('Chlorobenzene', 'chlorobenzene', 'C6H5Cl'), 'Nitromethane': ('Nitromethane', 'nitromethane', 'CH3NO2'), 'N-heptane': ('N-heptane', 'n-heptane', 'C7H16'), 'Cyclohexane': ('Cyclohexane', 'cyclohexane', 'C6H12'), 'Aniline': ('Aniline', 'aniline', 'C6H5NH2'), 'Acetone': ('Acetone', 'acetone', 'C2H6CO'), 'Tetrahydrofurane': ('Tetrahydrofurane', 'tetrahydrofurane', 'THF'), 'Dimethylsulfoxide': ('Dimethylsulfoxide', 'dimethylsulfoxide', 'DMSO'), 'Acetonitrile': ('Acetonitrile', 'acetonitrile', 'CH3CN'), 'Explicit': ('Explicit', 'explicit')} solvent = section.get('Solvent') explicitSolvent = solvent.get() in allowedSolvents['Explicit'] if(explicitSolvent): PRF = section.is_set('ProbeRadius') GIF = section.is_set('Green<inside>') GOF = section.is_set('Green<outside>') if (not PRF): print "Error: Explicit solvent chosen but ProbeRadius not specified" if (not GIF): print "Error: Explicit solvent chosen but Green<inside> not specified" if (not GOF): print "Error: Explicit solvent chosen but Green<outside> not specified" if (not GIF or not GOF or not PRF): sys.exit(1) solventFound = False for i, v in allowedSolvents.iteritems(): if (solvent.get() in v): solventName = i solventFound = True break if (not solventFound): print "Unknown solvent" print "Choose a solvent from the following list: " print allowedSolvents.keys() print "or specify the solvent data explicitly." sys.exit(1) correction = section.get('Correction') if (correction.get() < 0.0): print "Correction for CPCM solver must be greater than 0.0" sys.exit(1) convert_length_scalar(section.get('ProbeRadius')) radius = section.get('ProbeRadius') if (radius.get() < 0.1 or radius.get() > 100): print "Probe radius has to be within [0.1,100] Atomic Units" sys.exit(1) allowed_types = ('IEFPCM', 'CPCM', 'Wavelet', 'Linear') key = section.get('SolverType') val = key.get() if (val not in allowed_types): print "Allowed types are: ", allowed_types sys.exit(1) allowed_equations = ('FirstKind', 'SecondKind', 'Full') key = section.get('EquationType') val = key.get() if (val not in allowed_equations): print "Allowed equations are: ", allowed_equations sys.exit(1)
[docs]def verify_green(section): required = ('Type',) allowed = ('Vacuum', 'UniformDielectric', 'MetalSphere', 'GreensFunctionSum') allowed_der = ('Numerical', 'Derivative', 'Gradient', 'Hessian') green1 = section.fetch_sect('Green<one>') green2 = section.fetch_sect('Green<two>') eps = section.get('Eps') epsimg = section.get('EpsImg') epsre = section.get('EpsRe') convert_length_array(section.get('SpherePosition')) position = section.get('SpherePosition') convert_length_scalar(section.get('SphereRadius')) radius = section.get('SphereRadius') type=section.get('Type') if (type.get() not in allowed): print "Allowed Green's functions are:", allowed sys.exit(1) type=section.get('Der') if (type.get() not in allowed_der): print "Allowed Derivatives are:", allowed sys.exit(1) if (type.get() == 'UniformDielectric'): if not eps.is_set(): print "Eps not defined for UniformDielectric" sys.exit(1) if (type.get() == 'MetalSphere'): if not (eps.is_set() and epsre.is_set and epsimg.is_set()): print "Eps and/or EpsImg not defined for MetalSphere" sys.exit(1) if not (position.is_set() and radius.is_set()): print "SpherePosition and/or SphereRadius not defined for MetalSphere" sys.exit(1) if (len(position.get()) != 3): print "SpherePosition error" sys.exit(1) if (radius.get() < 0.1): print "Minimum value allowed for Radius is 0.1" sys.exit(1) if (type.get() == 'GreensFunctionSum'): if not (green1.is_set() and green2.is_set()): print "One or both components not defined for GreensFunctionSum" sys.exit(1)
[docs]def verify_spheres(keyword): length=len(keyword.get()) if (length % 4 != 0): print "Empty or incoherent Spheres list." sys.exit(1) convert_length_array(keyword)
[docs]def convert_length_array(keyword): length=len(keyword.get()) if (isAngstrom): for i in range(length): keyword[i] *= toAtomicUnits
[docs]def convert_length_scalar(keyword): if (isAngstrom): keyword[0] *= toAtomicUnits
[docs]def convert_area_scalar(keyword): if (isAngstrom): keyword[0] *= toAtomicUnits * toAtomicUnits
[docs]def preprocess(): """ Takes the PCM input file in @pcmsolver.inp, and preprocesses to make it machine-readable.""" valid_keywords = setup_keywords() input=getkw.GetkwParser() inkw=input.parseFile('@pcmsolver.inp') inkw.sanitize(valid_keywords) topsect=inkw.get_topsect() inkw.run_callbacks(valid_keywords) xfile='@pcmsolver.inp' fd=open(xfile,'w') fd.write(str(inkw.top)) fd.close()