Fix #48 add ecp arguments to main function and cli parser (#49)

* Add ECP argument to SPAHM(a,b) CLI arguments
* Includes ecp check to xyz_to_mol() and raises error

Author: YAY-C

qstack/compound.py

@@ -50,10 +50,14 @@ def xyz_to_mol(fin, basis="def2-svp", charge=0, spin=0, ignore=False, unit='ANG'
         mol.charge = -(sum(map(data.elements.charge, elements))%2)
         mol.spin = 0
 
-    if ecp:
+    if ecp is not None:
         mol.ecp = ecp
 
     mol.build()
+    species_charges = [data.elements.charge(z) for z in mol.elements]
+    if mol.basis == 'minao' and ecp is None and (np.array(species_charges) > 36).any():
+        msg = f"{mol.basis} basis set requires the use of effective core potentials for atoms with Z>36"
+        raise RuntimeError(msg)
     return mol
 
 

qstack/spahm/rho/atom.py

@@ -68,12 +68,13 @@ def main():
     parser.add_argument('--charge',    dest='charge',    default=0,                            type=int, help='total charge of the system (default: 0)')
     parser.add_argument('--spin',      dest='spin',      default=None,                         type=int, help='number of unpaired electrons (default: None) (use 0 to treat a closed-shell system in a UHF manner)')
     parser.add_argument('--xc',        dest='xc',        default=defaults.xc,                  type=str, help=f'DFT functional for the SAD guess (default: {defaults.xc})')
+    parser.add_argument('--ecp',        dest='ecp',        default=None,                  type=str, help=f'effective core potential to use (default: None)')
     parser.add_argument('--nameout',   dest='NameOut',   default=None,                         type=str, help='name of the output representations file.')
     parser.add_argument('--omod',      dest='omod',      default=defaults.omod,     nargs='+', type=str, help=f'model(s) for open-shell systems (alpha, beta, sum, diff, default: {defaults.omod})')
     args = parser.parse_args()
     print(vars(args))
 
-    mol = compound.xyz_to_mol(check_file(args.mol), args.basis, charge=args.charge, spin=args.spin, unit=args.units)
+    mol = compound.xyz_to_mol(check_file(args.mol), args.basis, charge=args.charge, spin=args.spin, unit=args.units, ecp=args.ecp)
     dm = None if args.dm is None else np.load(args.dm)
 
     representations = get_repr(mol, args.elements, args.charge, args.spin,

qstack/spahm/rho/bond.py

@@ -57,6 +57,7 @@ def main():
     parser.add_argument('--guess',         type=str,            dest='guess',          default=defaults.guess,           help='initial guess')
     parser.add_argument('--units',         type=str,            dest='units',          default='Angstrom',               help='the units of the input coordinates (default: Angstrom)')
     parser.add_argument('--basis',         type=str,            dest='basis'  ,        default=defaults.basis,           help='AO basis set (default=MINAO)')
+    parser.add_argument('--ecp',        dest='ecp',        default=None,                  type=str, help=f'effective core potential to use (default: None)')
     parser.add_argument('--charge',        type=str,            dest='charge',         default=None,                     help='charge / path to a file with a list of thereof')
     parser.add_argument('--spin',          type=str,            dest='spin',           default=None,                     help='number of unpaired electrons / path to a file with a list of thereof')
     parser.add_argument('--xc',            type=str,            dest='xc',             default=defaults.xc,              help=f'DFT functional for the SAD guess (default={defaults.xc})')
@@ -90,7 +91,7 @@ def main():
         xyzlist = utils.get_xyzlist(xyzlistfile)
         charge  = utils.get_chsp(args.charge, len(xyzlist))
         spin    = utils.get_chsp(args.spin,   len(xyzlist))
-    mols    = utils.load_mols(xyzlist, charge, spin, args.basis, args.print, units=args.units)
+    mols    = utils.load_mols(xyzlist, charge, spin, args.basis, args.print, units=args.units, ecp=args.ecp)
 
     allvec = get_repr(mols, xyzlist, args.guess, xc=args.xc, spin=args.spin, readdm=args.readdm, printlevel=args.print,
                       pairfile=args.pairfile, dump_and_exit=args.dump_and_exit,