File: help_basis

package info (click to toggle)
openmolcas 25.02-2
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 170,204 kB
  • sloc: f90: 498,088; fortran: 139,779; python: 13,587; ansic: 5,745; sh: 745; javascript: 660; pascal: 460; perl: 325; makefile: 17
file content (285 lines) | stat: -rwxr-xr-x 10,112 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
#!/usr/bin/env python3
# -*- coding: utf-8 -*-

#***********************************************************************
# This file is part of OpenMolcas.                                     *
#                                                                      *
# OpenMolcas is free software; you can redistribute it and/or modify   *
# it under the terms of the GNU Lesser General Public License, v. 2.1. *
# OpenMolcas is distributed in the hope that it will be useful, but it *
# is provided "as is" and without any express or implied warranties.   *
# For more details see the full text of the license in the file        *
# LICENSE or in <http://www.gnu.org/licenses/>.                        *
#                                                                      *
# Copyright (C) 2017, Ignacio Fdez. Galván                             *
#***********************************************************************

import sys
from glob import glob
from os.path import abspath, join, getmtime
from os import environ
import re
from fnmatch import translate
import json

try:
  from colorama import init, Fore, Style
except ImportError:
  def init(*args, **kwargs):
    pass
  class Dummy(object):
    pass
  Fore = Dummy()
  Fore.RED = ''
  Fore.BLUE = ''
  Fore.GREEN = ''
  Style = Dummy()
  Style.RESET_ALL = ''

helptext = '''
{1}################################################################################
#                  Querying the OpenMolcas basis set library                   #
################################################################################{0}

Use with a command-line argument that is a basis label in Molcas format:

       {2}[element].[name].[author].[primitives].[contraction].[aux. labels]*
          (1)     (2)     (3)        (4)           (5)          (6){0}

trailing dots can be omitted.

Examples:

In the simplest form, just write an element symbol to get a list of basis sets
for that element:
  {3}Fe{0}               lists all basis sets for iron
  
Additional fields can be included to restrict the search:
  {3}Fe..jensen{0}       lists all basis sets for iron, with "jensen" in the author
                   field (case insensitive)
  {3}Fe.....ECP{0}       lists all basis sets for silver that use ECP
                   
Wildcards are supported:
  {3}Fe.*cc-p*v?z*{0}    lists all basis sets for iron belonging to the "cc" families
                   
A single space matches empty fields (needs quotes):
  {3}"Fe..... ."{0}      lists all basis sets for silver with an empty 6th field
                   (probably all-electron)
                   
Minimum number of functions for ANO-type basis sets can be specified:
  {3}Fe....8s7p7d4f{0}   lists all basis sets for iron with a [8s7p7d4f] contraction or
                   larger if the basis set supports using fewer functions

Exact matches are printed in green, "fuzzy" matches in red
'''.format(Style.RESET_ALL, Fore.BLUE, Fore.RED, Fore.GREEN)

compat = False
try:
  flag = sys.argv[1]
except IndexError:
  flag = None
  pass
if (flag == '-bw'):
  sys.argv.pop(1)
  compat = True

# First get the template from the command line,
# Use colors in terminal, not if output is redirected to file
if (sys.stdout.isatty() and not compat):
  init()
else:
  init(strip=True)

atom_basis = re.compile(r'^/(\S*)\s*$')
alias = re.compile(r'^\s*([^#]\S*)\s+(\S*)\s*$')
props = ['atom', 'name', 'auth', 'prim', 'func', 'aux1', 'aux2']
angs = ['s', 'p', 'd', 'f', 'g', 'h', 'i', 'k', 'l', 'm', 'n', 'o', 'q', 'r', 't', 'u']

# Parse a basis set label into a dict:
# - The label has fields separated by "." in specific order
# - The fields correspond to the elements in "prop"
# - Missing fields are empty
def parse_basis(label):
  parts = label.rstrip('.').split('.')
  parts.extend((len(props)-len(parts))*[''])
  return {k: v for k, v in zip(props, parts)}

# Decide if a given field in a basis matches the template
# A field matches if:
# - The corresponding template field is empty
# - The field is equal to the template (case-insensitive, wildcards are parsed,
#   a single space in the template matches empty fields only)
# - For the "func" field, and if the "name" field is right, if the numbers of
#   contracted functions is equal or greater than in the template for all
#   angular momenta (i.e., one can select fewer functions for these types of bases)
def match_field(field, basis, template):
  if (template[field] == ''):
    return True
  elif (template[field] == ' '):
    return (basis[field] == '')
  if ((field == 'func') and
      (basis['name'].upper().startswith('ANO') or
       (basis['name'].upper() in ['ECP', 'PSD', 'RYDBERG'])
      )
     ):
    # Parse the "func" field into an array of integers,
    # for basis and template
    funcs_basis = []
    funcs_template = []
    # Catch wrong angular momentum labels (e.g. "j")
    if (any(i not in angs for i in re.sub(r'\d', '', template[field]))):
      print('The [contraction] field ({0}) contains an invalid label.'.format(Fore.RED + template[field] + Style.RESET_ALL))
      print('Valid labels are: {0}'.format(Fore.BLUE + ', '.join(angs) + Style.RESET_ALL))
      print()
      sys.exit(1)
    for i in angs:
      match = re.search(r'\d+(?={0})'.format(i), basis[field].lower())
      if match:
        funcs_basis.append(int(match.group()))
      else:
        funcs_basis.append(0)
      match = re.search(r'\d+(?={0})'.format(i), template[field].lower())
      if match:
        funcs_template.append(int(match.group()))
      else:
        funcs_template.append(0)
    return all(i >= j for i,j in zip(funcs_basis, funcs_template))
  else:
    # Compare with wildcard parsing
    return bool(re.match(translate(template[field].upper()), basis[field].upper()))

# Decide if a parsed basis set matches a template
# All fields have to match, so we can short-circuit
def match_basis(basis, template):
  for k in props:
    if (not match_field(k, basis, template)):
      return False
  return True

# Filter the list of basis aliases to keep only those that match the template
def filter_basis_alias(template):
  return [b for b in basis_alias if match_basis(parse_basis(b[0]), template)]

# Filter the list of basis sets to keep only those that match the template
def filter_basis_sets(template):
  return [b for b in basis_sets if match_basis(b, template)]

# Color a string according to whether or not it matches a template
# - If the string of the template is empty, don't color
# - If it matches (case insensitive), color green
# - If it doesn't, color red
def color(string, template):
  if ((template == '') or (string == '')):
    return string
  elif (template.upper() == string.upper()):
    return Fore.GREEN + string + Style.RESET_ALL
  else:
    return Fore.RED + string + Style.RESET_ALL

# Convert a parsed basis set into a basis label
# - Join the fields given by "props" with ".", and remove any trailing "."
# - Color according to the template, if given
def name_basis(basis, template=None):
  if (template is None):
    return '.'.join([basis[k] for k in props]).rstrip('.')
  else:
    return '.'.join([color(basis[k], template[k]) for k in props]).rstrip('.')

# so we can spare reading the library
try:
  template = parse_basis(sys.argv[1])
except IndexError:
  print(helptext)
  sys.exit(0)

try:
  MOLCAS = environ['MOLCAS']
except KeyError:
  MOLCAS = '.'

# Find out the latest modification date in the basis library
# in order to read/write the database file
basis_lib = abspath(join(MOLCAS, 'basis_library'))
basis_files = glob(join(basis_lib, '*'))
mtime = getmtime(basis_lib)
for f in basis_files:
  mtime = max(mtime, getmtime(f))

# Get a list of all basis sets and all aliases
# "basis_sets": list of parsed basis sets (dicts)
# "basis_alias": list of 2-tuples: (alias, basis set template)
try:
  # If the files are older than the database, read this
  db_file = join(MOLCAS, 'data', 'basis.db')
  if (mtime > getmtime(db_file)):
    raise
  with open(db_file, 'r') as f:
    basis_sets = json.loads(f.readline())
    basis_alias = json.loads(f.readline())
except:
  # Otherwise rebuild the database from the basis library files
  basis_sets = []
  basis_alias = []
  for basis in sorted(basis_files):
    try:
      with open(basis, 'r') as f:
        for line in f:
          match = atom_basis.match(line)
          if (match):
            basis_sets.append(parse_basis(match.group(1)))
    except IOError:
      pass
  with open(join(basis_lib, 'basis.tbl'), 'r') as f:
    for line in f:
      match = alias.match(line)
      if (match):
        basis_alias.append(match.group(1, 2))
  # Attempt to save the database for later use
  try:
    with open(db_file, 'w') as f:
      json.dump(basis_sets, f)
      f.write('\n')
      json.dump(basis_alias, f)
      f.write('\n')
  except:
    pass

# Print the location of the basis library
if (not compat):
  print(Fore.BLUE + 'Basis library:' + Style.RESET_ALL)
  print(basis_lib)

# Print the aliases and basis sets that match the template given as command-line argument
if (compat):
  fmt = '{0}              {1}'
else:
  fmt = '{0}   ->   {1}'
first = True
for alias, basis in filter_basis_alias(template):
  name1 = name_basis(parse_basis(alias), template)
  name2 = name_basis(parse_basis(basis), template)
  if (first):
    if (compat):
      print('#RECOMMENDED')
    else:
      print('')
      print(Fore.BLUE + 'Basis aliases for ' + Fore.GREEN + name_basis(template) + Style.RESET_ALL + '\n')
    first = False
  # Can't use width in format, because it doesn't work with ANSI codes
  name1 += (27 - len(name_basis(parse_basis(alias))))*' '
  print(fmt.format(name1, name2))
if (not first):
  first = True
for basis in filter_basis_sets(template):
  if (first):
    if (compat):
      print('#OTHER')
    else:
      print('')
      print(Fore.BLUE + 'Basis sets for ' + Fore.GREEN + name_basis(template) + Style.RESET_ALL + '\n')
    first = False
  print(name_basis(basis, template))
if (first):
  print('')
  print('No basis sets found matching the query: {0}'.format(Fore.RED + sys.argv[1] + Style.RESET_ALL))
print('')