Merge pull request #64 from TApplencourt/master

Pseudo! ❤️

Makefile

@@ -45,5 +45,5 @@ ocaml:
 
 veryclean:
 	rm -rf EZFIO
-	$(MAKE) EZFIO
+	rm -rf resultsFile
 	$(MAKE) -C src veryclean

ocaml/.merlin

@@ -0,0 +1,4 @@
+PKG core ZMQ cryptokit
+B _build/
+
+

ocaml/Makefile

@@ -35,7 +35,7 @@ default: $(ALL_TESTS) $(ALL_EXE) .gitignore
 executables: $(QPACKAGE_ROOT)/data/executables
 
 $(QPACKAGE_ROOT)/data/executables:
-	$(QPACKAGE_ROOT)/scripts/create_executables_list.sh
+	$(QPACKAGE_ROOT)/scripts/create/create_executables_list.sh
 
 external_libs:
 	opam install cryptokit core

ocaml/qp_create_ezfio_from_xyz.ml

@@ -5,18 +5,20 @@ open Core.Std;;
 let spec =
   let open Command.Spec in
   empty 
-  +> flag "o"  (optional string)
+  +> flag "o" (optional string)
      ~doc:"file Name of the created EZFIO file."
   +> flag "b" (required string)
      ~doc:"string Name of basis set."
   +> flag "c" (optional_with_default 0 int)
      ~doc:"int Total charge of the molecule. Default is 0."
-  +> flag "m"  (optional_with_default 1 int)
+  +> flag "m" (optional_with_default 1 int)
      ~doc:"int Spin multiplicity (2S+1) of the molecule. Default is 1."
+  +> flag "p" no_arg
+     ~doc:"Using pseudo."
   +> anon ("xyz_file" %: string)
 ;;
 
-let run ?o b c m xyz_file =
+let run ?o b c m p xyz_file =
 
   (* Read molecule *)
   let molecule =
@@ -60,6 +62,10 @@ let run ?o b c m xyz_file =
       | None -> (* Principal basis *)
         let basis = elem_and_basis_name in
         let command =
+          if (p) then  
+            Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename 
+              ^ "\" \"" ^ basis ^"\" pseudo"
+          else
             Qpackage.root ^ "/scripts/get_basis.sh \"" ^ temp_filename 
               ^ "\" \"" ^ basis ^"\""
         in
@@ -247,6 +253,13 @@ let run ?o b c m xyz_file =
   ~rank:2 ~dim:[| ao_num ; ao_prim_num_max |] ~data:ao_expo) ;
 
 
+  (* Doesn't work... *)
+  if (p) then 
+    begin
+      Qpackage.root ^ "/scripts/pseudo/put_pseudo_in_ezfio.py " ^ ezfio_file
+      |> Sys.command_exn
+    end;
+
   match Input.Ao_basis.read () with
   | None -> failwith "Error in basis"
   | Some x -> Input.Ao_basis.write x
@@ -266,8 +279,8 @@ elements can be defined as follows:
 
       ")
     spec
-    (fun o b c m xyz_file () ->
+    (fun o b c m p xyz_file () ->
-       run ?o b c m xyz_file )
+       run ?o b c m p xyz_file )
 ;;
 
 let () =

scripts/.gitignore

@@ -1,2 +1,3 @@
 *.pyc
 *.pyo
+docopt.py

scripts/check_dependencies.sh

@@ -1,63 +0,0 @@
-#!/bin/bash
-# 
-# usage:
-#  check_dependencies.sh MOs AOs Electrons
-#
-# Checks that the list of dependencies given in
-# argument is consistent. If the dependencies
-# are OK the exit code is 0, otherwise it is 1. 
-# If no argument is given, the dependencies are
-# read in the Makefile.
-# Thu Apr  3 01:44:23 CEST 2014
-
-if [[ -z ${QPACKAGE_ROOT} ]]
-then
-  print "The QPACKAGE_ROOT environment variable is not set."
-  print "Please reload the quantum_package.rc file."
-  exit -1
-fi
-source ${QPACKAGE_ROOT}/scripts/qp_include.sh
-
-if [[ -z $1 ]]
-then
-  exit 0
-fi
-
-if [[ $1 == "-" ]]
-then
-  COMMAND_LINE=$(cat NEEDED_MODULES)
-else
-  COMMAND_LINE=$(unique_list $@)
-fi
-
-for d in $COMMAND_LINE
-do
-  if [[ ! -d ${QPACKAGE_ROOT}/src/$d ]]
-  then
-    echo Error: Directory $d does not exist 
-    exit 2
-  fi
-
-done
-
-DEPS_LONG=""
-for i in $COMMAND_LINE
-do
-  DEPS_LONG+=" $i "
-  DEPS_LONG+=$(cat "${QPACKAGE_ROOT}/src/${i}/NEEDED_MODULES")
-done
-
-DEPS=$(unique_list $DEPS_LONG)
-
-if [[ ! "$COMMAND_LINE" == "$DEPS" ]]
-then
-  DEPS=$(${QPACKAGE_ROOT}/scripts/check_dependencies.sh ${DEPS})
-fi
-echo "$DEPS"
-
-if [[ "$COMMAND_LINE" == "$DEPS" ]]
-then
-  exit 0
-else
-  exit 1
-fi

scripts/docopt.py

@@ -1,590 +0,0 @@
-"""Pythonic command-line interface parser that will make you smile.
-
- * http://docopt.org
- * Repository and issue-tracker: https://github.com/docopt/docopt
- * Licensed under terms of MIT license (see LICENSE-MIT)
- * Copyright (c) 2013 Vladimir Keleshev, vladimir@keleshev.com
-
-"""
-import sys
-import re
-
-
-__all__ = ['docopt']
-__version__ = '0.6.1'
-
-
-class DocoptLanguageError(Exception):
-
-    """Error in construction of usage-message by developer."""
-
-
-class DocoptExit(SystemExit):
-
-    """Exit in case user invoked program with incorrect arguments."""
-
-    usage = ''
-
-    def __init__(self, message=''):
-        SystemExit.__init__(self, (message + '\n' + self.usage).strip())
-
-
-class Pattern(object):
-
-    def __eq__(self, other):
-        return repr(self) == repr(other)
-
-    def __hash__(self):
-        return hash(repr(self))
-
-    def fix(self):
-        self.fix_identities()
-        self.fix_repeating_arguments()
-        return self
-
-    def fix_identities(self, uniq=None):
-        """Make pattern-tree tips point to same object if they are equal."""
-        if not hasattr(self, 'children'):
-            return self
-        uniq = list(set(self.flat())) if uniq is None else uniq
-        for i, child in enumerate(self.children):
-            if not hasattr(child, 'children'):
-                assert child in uniq
-                self.children[i] = uniq[uniq.index(child)]
-            else:
-                child.fix_identities(uniq)
-
-    def fix_repeating_arguments(self):
-        """Fix elements that should accumulate/increment values."""
-        either = [list(child.children) for child in transform(self).children]
-        for case in either:
-            for e in [child for child in case if case.count(child) > 1]:
-                if isinstance(
-                        e,
-                        Argument) or isinstance(
-                        e,
-                        Option) and e.argcount:
-                    if e.value is None:
-                        e.value = []
-                    elif not isinstance(e.value, list):
-                        e.value = e.value.split()
-                if isinstance(
-                        e,
-                        Command) or isinstance(
-                        e,
-                        Option) and e.argcount == 0:
-                    e.value = 0
-        return self
-
-
-def transform(pattern):
-    """Expand pattern into an (almost) equivalent one, but with single Either.
-
-    Example: ((-a | -b) (-c | -d)) => (-a -c | -a -d | -b -c | -b -d)
-    Quirks: [-a] => (-a), (-a...) => (-a -a)
-
-    """
-    result = []
-    groups = [[pattern]]
-    while groups:
-        children = groups.pop(0)
-        parents = [Required, Optional, OptionsShortcut, Either, OneOrMore]
-        if any(t in map(type, children) for t in parents):
-            child = [c for c in children if type(c) in parents][0]
-            children.remove(child)
-            if isinstance(child, Either):
-                for c in child.children:
-                    groups.append([c] + children)
-            elif isinstance(child, OneOrMore):
-                groups.append(child.children * 2 + children)
-            else:
-                groups.append(child.children + children)
-        else:
-            result.append(children)
-    return Either(*[Required(*e) for e in result])
-
-
-class LeafPattern(Pattern):
-
-    """Leaf/terminal node of a pattern tree."""
-
-    def __init__(self, name, value=None):
-        self.name, self.value = name, value
-
-    def __repr__(self):
-        return '%s(%r, %r)' % (self.__class__.__name__, self.name, self.value)
-
-    def flat(self, *types):
-        return [self] if not types or type(self) in types else []
-
-    def match(self, left, collected=None):
-        collected = [] if collected is None else collected
-        pos, match = self.single_match(left)
-        if match is None:
-            return False, left, collected
-        left_ = left[:pos] + left[pos + 1:]
-        same_name = [a for a in collected if a.name == self.name]
-        if type(self.value) in (int, list):
-            if isinstance(self.value, int):
-                increment = 1
-            else:
-                increment = ([match.value] if isinstance(match.value, str)
-                             else match.value)
-            if not same_name:
-                match.value = increment
-                return True, left_, collected + [match]
-            same_name[0].value += increment
-            return True, left_, collected
-        return True, left_, collected + [match]
-
-
-class BranchPattern(Pattern):
-
-    """Branch/inner node of a pattern tree."""
-
-    def __init__(self, *children):
-        self.children = list(children)
-
-    def __repr__(self):
-        return '%s(%s)' % (self.__class__.__name__,
-                           ', '.join(repr(a) for a in self.children))
-
-    def flat(self, *types):
-        if type(self) in types:
-            return [self]
-        return sum([child.flat(*types) for child in self.children], [])
-
-
-class Argument(LeafPattern):
-
-    def single_match(self, left):
-        for n, pattern in enumerate(left):
-            if isinstance(pattern, Argument):
-                return n, Argument(self.name, pattern.value)
-        return None, None
-
-    @classmethod
-    def parse(class_, source):
-        name = re.findall('(<\S*?>)', source)[0]
-        value = re.findall('\[default: (.*)\]', source, flags=re.I)
-        return class_(name, value[0] if value else None)
-
-
-class Command(Argument):
-
-    def __init__(self, name, value=False):
-        self.name, self.value = name, value
-
-    def single_match(self, left):
-        for n, pattern in enumerate(left):
-            if isinstance(pattern, Argument):
-                if pattern.value == self.name:
-                    return n, Command(self.name, True)
-                else:
-                    break
-        return None, None
-
-
-class Option(LeafPattern):
-
-    def __init__(self, short=None, long=None, argcount=0, value=False):
-        assert argcount in (0, 1)
-        self.short, self.long, self.argcount = short, long, argcount
-        self.value = None if value is False and argcount else value
-
-    @classmethod
-    def parse(class_, option_description):
-        short, long, argcount, value = None, None, 0, False
-        options, _, description = option_description.strip().partition('  ')
-        options = options.replace(',', ' ').replace('=', ' ')
-        for s in options.split():
-            if s.startswith('--'):
-                long = s
-            elif s.startswith('-'):
-                short = s
-            else:
-                argcount = 1
-        if argcount:
-            matched = re.findall('\[default: (.*)\]', description, flags=re.I)
-            value = matched[0] if matched else None
-        return class_(short, long, argcount, value)
-
-    def single_match(self, left):
-        for n, pattern in enumerate(left):
-            if self.name == pattern.name:
-                return n, pattern
-        return None, None
-
-    @property
-    def name(self):
-        return self.long or self.short
-
-    def __repr__(self):
-        return 'Option(%r, %r, %r, %r)' % (self.short, self.long,
-                                           self.argcount, self.value)
-
-
-class Required(BranchPattern):
-
-    def match(self, left, collected=None):
-        collected = [] if collected is None else collected
-        l = left
-        c = collected
-        for pattern in self.children:
-            matched, l, c = pattern.match(l, c)
-            if not matched:
-                return False, left, collected
-        return True, l, c
-
-
-class Optional(BranchPattern):
-
-    def match(self, left, collected=None):
-        collected = [] if collected is None else collected
-        for pattern in self.children:
-            m, left, collected = pattern.match(left, collected)
-        return True, left, collected
-
-
-class OptionsShortcut(Optional):
-
-    """Marker/placeholder for [options] shortcut."""
-
-
-class OneOrMore(BranchPattern):
-
-    def match(self, left, collected=None):
-        assert len(self.children) == 1
-        collected = [] if collected is None else collected
-        l = left
-        c = collected
-        l_ = None
-        matched = True
-        times = 0
-        while matched:
-            # could it be that something didn't match but changed l or c?
-            matched, l, c = self.children[0].match(l, c)
-            times += 1 if matched else 0
-            if l_ == l:
-                break
-            l_ = l
-        if times >= 1:
-            return True, l, c
-        return False, left, collected
-
-
-class Either(BranchPattern):
-
-    def match(self, left, collected=None):
-        collected = [] if collected is None else collected
-        outcomes = []
-        for pattern in self.children:
-            matched, _, _ = outcome = pattern.match(left, collected)
-            if matched:
-                outcomes.append(outcome)
-        if outcomes:
-            return min(outcomes, key=lambda outcome: len(outcome[1]))
-        return False, left, collected
-
-
-class Tokens(list):
-
-    def __init__(self, source, error=DocoptExit):
-        self += source.split() if hasattr(source, 'split') else source
-        self.error = error
-
-    @staticmethod
-    def from_pattern(source):
-        source = re.sub(r'([\[\]\(\)\|]|\.\.\.)', r' \1 ', source)
-        source = [s for s in re.split('\s+|(\S*<.*?>)', source) if s]
-        return Tokens(source, error=DocoptLanguageError)
-
-    def move(self):
-        return self.pop(0) if len(self) else None
-
-    def current(self):
-        return self[0] if len(self) else None
-
-
-def parse_long(tokens, options):
-    """long ::= '--' chars [ ( ' ' | '=' ) chars ] ;"""
-    long, eq, value = tokens.move().partition('=')
-    assert long.startswith('--')
-    value = None if eq == value == '' else value
-    similar = [o for o in options if o.long == long]
-    if tokens.error is DocoptExit and similar == []:  # if no exact match
-        similar = [o for o in options if o.long and o.long.startswith(long)]
-    if len(similar) > 1:  # might be simply specified ambiguously 2+ times?
-        raise tokens.error('%s is not a unique prefix: %s?' %
-                           (long, ', '.join(o.long for o in similar)))
-    elif len(similar) < 1:
-        argcount = 1 if eq == '=' else 0
-        o = Option(None, long, argcount)
-        options.append(o)
-        if tokens.error is DocoptExit:
-            o = Option(None, long, argcount, value if argcount else True)
-    else:
-        o = Option(similar[0].short, similar[0].long,
-                   similar[0].argcount, similar[0].value)
-        if o.argcount == 0:
-            if value is not None:
-                raise tokens.error('%s must not have an argument' % o.long)
-        else:
-            if value is None:
-                if tokens.current() in [None, '--']:
-                    raise tokens.error('%s requires argument' % o.long)
-                value = tokens.move()
-        if tokens.error is DocoptExit:
-            o.value = value if value is not None else True
-    return [o]
-
-
-def parse_shorts(tokens, options):
-    """shorts ::= '-' ( chars )* [ [ ' ' ] chars ] ;"""
-    token = tokens.move()
-    assert token.startswith('-') and not token.startswith('--')
-    left = token.lstrip('-')
-    parsed = []
-    while left != '':
-        short, left = '-' + left[0], left[1:]
-        similar = [o for o in options if o.short == short]
-        if len(similar) > 1:
-            raise tokens.error('%s is specified ambiguously %d times' %
-                               (short, len(similar)))
-        elif len(similar) < 1:
-            o = Option(short, None, 0)
-            options.append(o)
-            if tokens.error is DocoptExit:
-                o = Option(short, None, 0, True)
-        else:  # why copying is necessary here?
-            o = Option(short, similar[0].long,
-                       similar[0].argcount, similar[0].value)
-            value = None
-            if o.argcount != 0:
-                if left == '':
-                    if tokens.current() in [None, '--']:
-                        raise tokens.error('%s requires argument' % short)
-                    value = tokens.move()
-                else:
-                    value = left
-                    left = ''
-            if tokens.error is DocoptExit:
-                o.value = value if value is not None else True
-        parsed.append(o)
-    return parsed
-
-
-def parse_pattern(source, options):
-    tokens = Tokens.from_pattern(source)
-    result = parse_expr(tokens, options)
-    if tokens.current() is not None:
-        raise tokens.error('unexpected ending: %r' % ' '.join(tokens))
-    return Required(*result)
-
-
-def parse_expr(tokens, options):
-    """expr ::= seq ( '|' seq )* ;"""
-    seq = parse_seq(tokens, options)
-    if tokens.current() != '|':
-        return seq
-    result = [Required(*seq)] if len(seq) > 1 else seq
-    while tokens.current() == '|':
-        tokens.move()
-        seq = parse_seq(tokens, options)
-        result += [Required(*seq)] if len(seq) > 1 else seq
-    return [Either(*result)] if len(result) > 1 else result
-
-
-def parse_seq(tokens, options):
-    """seq ::= ( atom [ '...' ] )* ;"""
-    result = []
-    while tokens.current() not in [None, ']', ')', '|']:
-        atom = parse_atom(tokens, options)
-        if tokens.current() == '...':
-            atom = [OneOrMore(*atom)]
-            tokens.move()
-        result += atom
-    return result
-
-
-def parse_atom(tokens, options):
-    """atom ::= '(' expr ')' | '[' expr ']' | 'options'
-             | long | shorts | argument | command ;
-    """
-    token = tokens.current()
-    result = []
-    if token in '([':
-        tokens.move()
-        matching, pattern = {'(': [')', Required], '[': [']', Optional]}[token]
-        result = pattern(*parse_expr(tokens, options))
-        if tokens.move() != matching:
-            raise tokens.error("unmatched '%s'" % token)
-        return [result]
-    elif token == 'options':
-        tokens.move()
-        return [OptionsShortcut()]
-    elif token.startswith('--') and token != '--':
-        return parse_long(tokens, options)
-    elif token.startswith('-') and token not in ('-', '--'):
-        return parse_shorts(tokens, options)
-    elif token.startswith('<') and token.endswith('>') or token.isupper():
-        return [Argument(tokens.move())]
-    else:
-        return [Command(tokens.move())]
-
-
-def parse_argv(tokens, options, options_first=False):
-    """Parse command-line argument vector.
-
-    If options_first:
-        argv ::= [ long | shorts ]* [ argument ]* [ '--' [ argument ]* ] ;
-    else:
-        argv ::= [ long | shorts | argument ]* [ '--' [ argument ]* ] ;
-
-    """
-    parsed = []
-    while tokens.current() is not None:
-        if tokens.current() == '--':
-            return parsed + [Argument(None, v) for v in tokens]
-        elif tokens.current().startswith('--'):
-            parsed += parse_long(tokens, options)
-        elif tokens.current().startswith('-') and tokens.current() != '-':
-            parsed += parse_shorts(tokens, options)
-        elif options_first:
-            return parsed + [Argument(None, v) for v in tokens]
-        else:
-            parsed.append(Argument(None, tokens.move()))
-    return parsed
-
-
-def parse_defaults(doc):
-    defaults = []
-    for s in parse_section('options:', doc):
-        # FIXME corner case "bla: options: --foo"
-        _, _, s = s.partition(':')  # get rid of "options:"
-        split = re.split('\n[ \t]*(-\S+?)', '\n' + s)[1:]
-        split = [s1 + s2 for s1, s2 in zip(split[::2], split[1::2])]
-        options = [Option.parse(s) for s in split if s.startswith('-')]
-        defaults += options
-    return defaults
-
-
-def parse_section(name, source):
-    pattern = re.compile('^([^\n]*' + name + '[^\n]*\n?(?:[ \t].*?(?:\n|$))*)',
-                         re.IGNORECASE | re.MULTILINE)
-    return [s.strip() for s in pattern.findall(source)]
-
-
-def formal_usage(section):
-    _, _, section = section.partition(':')  # drop "usage:"
-    pu = section.split()
-    return '( ' + ' '.join(') | (' if s == pu[0] else s for s in pu[1:]) + ' )'
-
-
-def extras(help, version, options, doc):
-    if help and any((o.name in ('-h', '--help')) and o.value for o in options):
-        print(doc.strip("\n"))
-        sys.exit()
-    if version and any(o.name == '--version' and o.value for o in options):
-        print(version)
-        sys.exit()
-
-
-class Dict(dict):
-
-    def __repr__(self):
-        return '{%s}' % ',\n '.join('%r: %r' % i for i in sorted(self.items()))
-
-
-def docopt(doc, argv=None, help=True, version=None, options_first=False):
-    """Parse `argv` based on command-line interface described in `doc`.
-
-    `docopt` creates your command-line interface based on its
-    description that you pass as `doc`. Such description can contain
-    --options, <positional-argument>, commands, which could be
-    [optional], (required), (mutually | exclusive) or repeated...
-
-    Parameters
-    ----------
-    doc : str
-        Description of your command-line interface.
-    argv : list of str, optional
-        Argument vector to be parsed. sys.argv[1:] is used if not
-        provided.
-    help : bool (default: True)
-        Set to False to disable automatic help on -h or --help
-        options.
-    version : any object
-        If passed, the object will be printed if --version is in
-        `argv`.
-    options_first : bool (default: False)
-        Set to True to require options precede positional arguments,
-        i.e. to forbid options and positional arguments intermix.
-
-    Returns
-    -------
-    args : dict
-        A dictionary, where keys are names of command-line elements
-        such as e.g. "--verbose" and "<path>", and values are the
-        parsed values of those elements.
-
-    Example
-    -------
-    >>> from docopt import docopt
-    >>> doc = '''
-    ... Usage:
-    ...     my_program tcp <host> <port> [--timeout=<seconds>]
-    ...     my_program serial <port> [--baud=<n>] [--timeout=<seconds>]
-    ...     my_program (-h | --help | --version)
-    ...
-    ... Options:
-    ...     -h, --help  Show this screen and exit.
-    ...     --baud=<n>  Baudrate [default: 9600]
-    ... '''
-    >>> argv = ['tcp', '127.0.0.1', '80', '--timeout', '30']
-    >>> docopt(doc, argv)
-    {'--baud': '9600',
-     '--help': False,
-     '--timeout': '30',
-     '--version': False,
-     '<host>': '127.0.0.1',
-     '<port>': '80',
-     'serial': False,
-     'tcp': True}
-
-    See also
-    --------
-    * For video introduction see http://docopt.org
-    * Full documentation is available in README.rst as well as online
-      at https://github.com/docopt/docopt#readme
-
-    """
-    argv = sys.argv[1:] if argv is None else argv
-
-    usage_sections = parse_section('usage:', doc)
-    if len(usage_sections) == 0:
-        raise DocoptLanguageError('"usage:" (case-insensitive) not found.')
-    if len(usage_sections) > 1:
-        raise DocoptLanguageError('More than one "usage:" (case-insensitive).')
-    DocoptExit.usage = usage_sections[0]
-
-    options = parse_defaults(doc)
-    pattern = parse_pattern(formal_usage(DocoptExit.usage), options)
-    # [default] syntax for argument is disabled
-    # for a in pattern.flat(Argument):
-    #    same_name = [d for d in arguments if d.name == a.name]
-    #    if same_name:
-    #        a.value = same_name[0].value
-    argv = parse_argv(Tokens(argv), list(options), options_first)
-    pattern_options = set(pattern.flat(Option))
-    for options_shortcut in pattern.flat(OptionsShortcut):
-        doc_options = parse_defaults(doc)
-        options_shortcut.children = list(set(doc_options) - pattern_options)
-        # if any_options:
-        #    options_shortcut.children += [Option(o.short, o.long, o.argcount)
-        #                    for o in argv if type(o) is Option]
-    extras(help, version, argv, doc)
-    matched, left, collected = pattern.fix().match(argv)
-    if matched and left == []:  # better error message if left?
-        return Dict((a.name, a.value) for a in (pattern.flat() + collected))
-    raise DocoptExit()

scripts/ezfio_interface/ei_handler.py

@@ -62,6 +62,9 @@ import ConfigParser
 from collections import defaultdict
 from collections import namedtuple
 
+
+from qp_utils import cache
+
 Type = namedtuple('Type', 'fancy ocaml fortran')
 
 
@@ -78,6 +81,7 @@ def is_bool(str_):
         raise TypeError
 
 
+@cache
 def get_type_dict():
     """
     This function makes the correspondance between the type of value read in
@@ -89,17 +93,7 @@ def get_type_dict():
     # ~#~#~#~#~ #
     # P i c l e #
     # ~#~#~#~#~ #
-
-    import cPickle as pickle
-
-    from os import listdir
-
     qpackage_root = os.environ['QPACKAGE_ROOT']
-    fancy_type_pickle = qpackage_root + "/scripts/ezfio_interface/fancy_type.p"
-
-    if fancy_type_pickle in listdir(os.getcwd()):
-        fancy_type = pickle.load(open(fancy_type_pickle, "rb"))
-        return fancy_type
 
     # ~#~#~#~ #
     # I n i t #
@@ -148,9 +142,7 @@ def get_type_dict():
         b = r.find('let untouched = "')
         e = r.find(';;', b)
 
-        l_un = [
+        l_un = [i for i in r[b:e].splitlines() if i.strip().startswith("module")]
-            i for i in r[
-                b:e].splitlines() if i.strip().startswith("module")]
 
     # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
     # q p _ t y p e s _ g e n e r a t e #
@@ -174,9 +166,6 @@ def get_type_dict():
     # ~#~#~#~#~#~#~#~ #
     # F i n a l i z e #
     # ~#~#~#~#~#~#~#~ #
-
-    pickle.dump(dict(fancy_type), open(fancy_type_pickle, "wb"))
-
     return dict(fancy_type)
 
 

scripts/ezfio_interface/qp_convert_output_to_ezfio.py

@@ -105,7 +105,7 @@ def write_ezfio(res, filename):
     # Transformt H1 into H
     import re
     p = re.compile(ur'(\d*)$')
-    label = [p.sub("", x.name) for x in res.geometry]
+    label = [p.sub("", x.name).capitalize() for x in res.geometry]
     ezfio.set_nuclei_nucl_label(label)
 
     ezfio.set_nuclei_nucl_coord(coord_x + coord_y + coord_z)
@@ -265,6 +265,16 @@ def write_ezfio(res, filename):
     ezfio.set_mo_basis_mo_occ(OccNum)
     ezfio.set_mo_basis_mo_coef(MoMatrix)
 
+    # ______                   _
+    # | ___ \                 | |
+    # | |_/ /__  ___ _   _  __| | ___
+    # |  __/ __|/ _ \ | | |/ _` |/ _ \
+    # | |  \__ \  __/ |_| | (_| | (_) |
+    # \_|  |___/\___|\__,_|\__,_|\___/
+    #
+
+    ezfio.set_pseudo_integrals_do_pseudo(False)
+
 
 def get_full_path(file_path):
     file_path = os.path.expanduser(file_path)

scripts/ezfio_interface/upgrade_1.0_2.0.sh

@@ -0,0 +1,27 @@
+#!/bin/bash
+#  Convert a old ezfio file (with option.irp.f ezfio_default)
+# into a new EZFIO.cfg type
+
+# Hartree Fock
+# Changin the case, don't know if is needed or not 
+mv $1/Hartree_Fock $1/hartree_fock 2> /dev/null
+
+mv $1/hartree_Fock/thresh_SCF $1/hartree_fock/thresh_scf 2> /dev/null
+
+# BiInts
+mv $1/bi_integrals $1/bielect_integrals 2> /dev/null
+
+if [ -f  $1/bielect_integrals/read_ao_integrals ]; then
+    if [ `cat $1/bielect_integrals/read_ao_integrals` -eq "True" ]
+    then
+            echo "Read" > $1/bielect_integrals/disk_access_ao_integrals
+    
+    elif [ `cat bielect_integrals/write_ao_integrals` -eq "True" ]
+    then
+            echo "Write" > $1/bielect_integrals/disk_access_ao_integrals
+    
+    else
+            echo "None" > $1/bielect_integrals/disk_access_ao_integrals
+    
+    fi
+fi

scripts/generate_h_apply.py

@@ -141,7 +141,7 @@ class H_apply(object):
   def set_filter_2h_2p(self):
     self["filter2h2p"] = """
 !    ! DIR$ FORCEINLINE
-     if(is_a_two_holes_two_particles(key))cycle
+     if (is_a_two_holes_two_particles(key)) cycle
     """
 
 

scripts/get_basis.sh

@@ -42,9 +42,14 @@ then
   echo "ERROR"
   exit 1
 fi
-${EMSL_API_ROOT}/EMSL_api.py get_basis_data --treat_l --save --path="${tmpfile}" --basis="${basis}" $atoms
-
-
 
+pseudo="$1"
+shift
 
+if [[ -z $pseudo ]]
+then
+    ${EMSL_API_ROOT}/EMSL_api.py get_basis_data --treat_l --save --path="${tmpfile}" --basis="${basis}"
+else
+    ${EMSL_API_ROOT}/EMSL_api.py get_basis_data --save --path="${tmpfile}" --basis="${basis}" --db_path="${EMSL_API_ROOT}/db/Pseudo.db"
+fi
 

scripts/install/install_curl.sh

@@ -8,8 +8,8 @@ CURL_URL="http://qmcchem.ups-tlse.fr/files/scemama/${CURL}.tar.bz2"
 
 if [[ -z ${QPACKAGE_ROOT} ]]
 then
-  print "The QPACKAGE_ROOT environment variable is not set."
+  echo "The QPACKAGE_ROOT environment variable is not set."
-  print "Please reload the quantum_package.rc file."
+  echo "Please reload the quantum_package.rc file."
   exit -1
 fi
 

scripts/install/install_docopt.sh

@@ -0,0 +1,21 @@
+#!/bin/bash
+#
+# Installs docopt
+# lundi 27 avril 2015, 16:51:34 (UTC+0200)
+
+DOCOPT="docopt.py"
+DOCOPT_URL="https://raw.githubusercontent.com/docopt/docopt/master/${DOCOPT}"
+
+if [[ -z ${QPACKAGE_ROOT} ]]
+then
+  echo "The QPACKAGE_ROOT environment variable is not set."
+  echo "Please reload the quantum_package.rc file."
+  exit -1
+fi
+
+cd ${QPACKAGE_ROOT}
+
+rm -f -- scripts/${DOCOPT}{,c}
+${QPACKAGE_ROOT}/scripts/fetch_from_web.py ${DOCOPT_URL} ${DOCOPT}
+
+mv ${DOCOPT} scripts/${DOCOPT}

scripts/install/install_emsl.sh

@@ -8,8 +8,8 @@ URL="https://github.com/LCPQ/${BASE}/archive/master.tar.gz"
 
 if [[ -z ${QPACKAGE_ROOT} ]]
 then
-  print "The QPACKAGE_ROOT environment variable is not set."
+  echo "The QPACKAGE_ROOT environment variable is not set."
-  print "Please reload the quantum_package.rc file."
+  echo "Please reload the quantum_package.rc file."
   exit -1
 fi
 

scripts/install/install_ezfio.sh

@@ -8,8 +8,8 @@ URL="https://github.com/LCPQ/${BASE}/archive/master.tar.gz"
 
 if [[ -z ${QPACKAGE_ROOT} ]]
 then
-  print "The QPACKAGE_ROOT environment variable is not set."
+  echo "The QPACKAGE_ROOT environment variable is not set."
-  print "Please reload the quantum_package.rc file."
+  echo "Please reload the quantum_package.rc file."
   exit -1
 fi
 

scripts/install/install_resultsFile.sh

@@ -16,6 +16,7 @@ fi
 cd ${QPACKAGE_ROOT}
 
 
+rm -rf resultsFile-master
 ${QPACKAGE_ROOT}/scripts/fetch_from_web.py ${URL} ${QPACKAGE_ROOT}/resultsFile.tar.gz
 tar -zxf resultsFile.tar.gz && rm resultsFile.tar.gz ||exit 1
 mv resultsFile-master resultsFile

scripts/module/build_modules.sh

@@ -30,7 +30,7 @@ Build failed for module $MODULE
 "
     fi
   fi
-  ${QPACKAGE_ROOT}/scripts/create_gitignore.sh
+  ${QPACKAGE_ROOT}/scripts/module/create_gitignore.sh
   cd ${OLDPWD}
 done
-${QPACKAGE_ROOT}/scripts/create_executables_list.sh
+${QPACKAGE_ROOT}/scripts/module/create_executables_list.sh

scripts/module/clean_modules.sh

@@ -13,8 +13,11 @@ source ${QPACKAGE_ROOT}/scripts/qp_include.sh
 function do_clean()
 {
   rm -rf -- \
-    IRPF90_temp IRPF90_man Makefile.depend $(cat NEEDED_MODULES) include \
+    IRPF90_temp IRPF90_man Makefile.depend \
+    $(module_handler.py print_genealogy) include \
     ezfio_interface.irp.f irpf90.make irpf90_entities tags $(ls_exe) *.mod
+
+  touch -c EZFIO.cfg *.ezfio_config
 }
 
 if [[ -z $1 ]]

scripts/module/create_Needed_modules.sh

@@ -14,7 +14,9 @@ source ${QPACKAGE_ROOT}/scripts/qp_include.sh
 
 check_current_dir_is_module
 
-OUTPUT=$(${QPACKAGE_ROOT}/scripts/check_dependencies.sh $@)
+OUTPUT=$(module_handler.py check_dependencies $@)
-echo ${OUTPUT} > NEEDED_MODULES
-
 
+if [[ $? -eq 0 ]]
+then
+  echo $@ > NEEDED_CHILDREN_MODULES
+fi

scripts/module/module_handler.py

@@ -0,0 +1,185 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+"""
+Create the NEEDED_MODULE
+    aka the genealogy (children module, subchildren module and so on),
+of a NEEDED_CHILDREN_MODULES file
+
+Usage:
+    module_handler.py print_genealogy    [<NEEDED_CHILDREN_MODULES>]
+    module_handler.py check_dependencies   [<module_name>...]
+    module_handler.py create_png         [<NEEDED_CHILDREN_MODULES>]
+
+Options:
+    print_genealogy         Print the genealogy of the NEEDED_CHILDREN_MODULES
+                                 aka (children, subchildren, etc)
+    create_png              Create a png of the file
+    NEEDED_CHILDREN_MODULES The path of NEEDED_CHILDREN_MODULES
+                                by default try to open the file in the current path
+"""
+
+from docopt import docopt
+
+import os
+import sys
+import os.path
+from qp_utils import cache
+
+
+@cache
+def get_dict_genealogy():
+    """Loop over MODULE in  QPACKAGE_ROOT/src, open all the NEEDED_CHILDREN_MODULES
+    and create a dict[MODULE] = [sub module needed, ...]
+    """
+    d_ref = dict()
+
+    qpackage_root = os.environ['QPACKAGE_ROOT']
+    dir_ = os.path.join(qpackage_root, 'src')
+
+    for o in os.listdir(dir_):
+
+        try:
+            with open(os.path.join(dir_, o, "NEEDED_CHILDREN_MODULES"), "r") as f:
+                l_children = f.read().split()
+        except IOError:
+            pass
+        else:
+            d_ref[o] = l_children
+
+    return d_ref
+
+
+def module_genealogy(path):
+    """
+    Take a name of a NEEDED_CHILDREN_MODULES
+    and return a list of all the {sub, subsub, ...}children
+    """
+    try:
+        with open(path, "r") as f:
+            l_children = f.read().split()
+    except IOError as e:
+        print >> sys.stderr, e
+        sys.exit(1)
+    else:
+
+        needed_module = get_it_and_children(l_children)
+
+        return needed_module
+
+
+def get_it_and_children(l_module):
+    """
+    From a list of module return the module and all of the genealogy
+    """
+    d_ref = get_dict_genealogy()
+
+    l = []
+    for module in l_module:
+        if module not in l:
+            l.append(module)
+            try:
+                l.extend(get_it_and_children(d_ref[module]))
+            except KeyError:
+                print >> sys.stderr, "`{0}` in not a good submodule name".format(module)
+                print >> sys.stderr, "Check the corresponding NEEDED_CHILDREN_MODULES"
+                sys.exit(1)
+
+    return list(set(l))
+
+
+def get_all_children(l_module):
+    """
+    From a list of module return all the genealogy
+    """
+
+    it_and_all = get_it_and_children(l_module)
+    return [children for children in it_and_all if children not in l_module]
+
+
+def reduce_(l_module):
+    """
+    Take a l_module and try to find the lower combinaitions
+    of module with the same genealogy
+    """
+    import itertools
+    d_ref = get_dict_genealogy()
+
+    target_genealogy = sorted(get_all_children(l_module))
+
+    for i in xrange(len(d_ref)):
+        for c in itertools.combinations(d_ref, i):
+
+                guess_genealogy = sorted(get_it_and_children(d_ref, c))
+
+                if target_genealogy == guess_genealogy:
+                    return c
+
+
+def create_png_from_path(path):
+    " Change a path like this into a module list"
+    "path = /home/razoa/quantum_package/src/Molden/NEEDED_CHILDREN_MODULES"
+
+    l_module = os.path.split(path)[0].split("/")[-1]
+    create_png([l_module])
+
+
+def create_png(l_module):
+    """Create the png of the dependancy tree for a l_module"""
+
+    # Init
+    import pydot
+    all_ready_done = []
+
+    def draw_module_edge(module, l_children):
+        "Draw all the module recursifly"
+
+        if module not in all_ready_done:
+            for children in l_children:
+                        # Add Edge
+                        edge = pydot.Edge(module, children)
+                        graph.add_edge(edge)
+                        # Recurs
+                        draw_module_edge(children, d_ref[children])
+            all_ready_done.append(module)
+
+    # Init
+    graph = pydot.Dot(graph_type='digraph')
+    d_ref = get_dict_genealogy()
+
+    # Create all the edge
+    for module in l_module:
+        node_a = pydot.Node(module, fontcolor="red")
+        graph.add_node(node_a)
+        draw_module_edge(module, d_ref[module])
+
+    # Save
+    path = '{0}.png'.format("_".join(l_module))
+    print "png saved in {0}".format(path)
+    graph.write_png(path)
+
+if __name__ == '__main__':
+
+    arguments = docopt(__doc__)
+
+    if not arguments['<NEEDED_CHILDREN_MODULES>']:
+        dir_ = os.getcwd()
+        path = os.path.join(dir_, "NEEDED_CHILDREN_MODULES")
+    else:
+        path = os.path.abspath(arguments['<NEEDED_CHILDREN_MODULES>'])
+        path = os.path.expanduser(path)
+        path = os.path.expandvars(path)
+
+    if arguments['print_genealogy']:
+        l_all_needed_molule = module_genealogy(path)
+        print " ".join(sorted(l_all_needed_molule))
+
+    elif arguments["check_dependencies"]:
+        l_module = arguments['<module_name>']
+        if l_module:
+            l_all_needed_molule = get_it_and_children(l_module)
+        else:
+            l_all_needed_molule = module_genealogy(path)
+
+    elif arguments["create_png"]:
+        create_png_from_path(path)

scripts/module/qp_create_module.sh

@@ -71,7 +71,7 @@ debug "Module does not already exist: OK"
 
 
 # Set up dependencies
-ALL_MODULES="${NEEDED_MODULES}"
+ALL_MODULES="$(cat NEEDED_MODULES)"
 echo "Select which modules you are sure you will need:  (press q to quit)"  
 NEEDED_MODULES=""
 select M in ${ALL_MODULES}
@@ -117,7 +117,7 @@ debug "Module directory is created."
 
 
 # Create the Makefile
-"${QPACKAGE_ROOT}/scripts/create_Makefile.sh" || fail "Unable to create Makefile"
+"${QPACKAGE_ROOT}/scripts/module/create_Makefile.sh" || fail "Unable to create Makefile"
 if [[ ! -f Makefile ]]
 then
   fail "Makefile was not created"
@@ -125,8 +125,8 @@ fi
 debug "Makefile created"
 
 # Create the NEEDED_MODULES file
-"${QPACKAGE_ROOT}/scripts/create_Needed_modules.sh" ${NEEDED_MODULES} || fail "Unable to create the NEEDED_MODULES file"
+"${QPACKAGE_ROOT}/scripts/module/create_Needed_modules.sh" ${NEEDED_MODULES} || fail "Unable to create the NEEDED_MODULES file"
-if [[ ! -f NEEDED_MODULES ]]
+if [[ ! -f NEEDED_CHILDREN_MODULES ]]
 then
   fail "NEEDED_MODULES was not created"
 fi
@@ -135,7 +135,7 @@ debug "NEEDED_MODULES created"
 
 
 # Create rst templates
-"${QPACKAGE_ROOT}/scripts/create_rst_templates.sh" || fail "Unable to create rst templates"
+"${QPACKAGE_ROOT}/scripts/module/create_rst_templates.sh" || fail "Unable to create rst templates"
 
 
 # Update module list in main NEEDED_MODULES

scripts/pseudo/elts_num_ele.py

@@ -0,0 +1,118 @@
+name_to_elec = {"H": 1,
+                "He": 2,
+                "Li": 3,
+                "Be": 4,
+                "B": 5,
+                "C": 6,
+                "N": 7,
+                "O": 8,
+                "F": 9,
+                "Ne": 10,
+                "Na": 11,
+                "Mg": 12,
+                "Al": 13,
+                "Si": 14,
+                "P": 15,
+                "S": 16,
+                "Cl": 17,
+                "Ar": 18,
+                "K": 19,
+                "Ca": 20,
+                "Sc": 21,
+                "Ti": 22,
+                "V": 23,
+                "Cr": 24,
+                "Mn": 25,
+                "Fe": 26,
+                "Co": 27,
+                "Ni": 28,
+                "Cu": 29,
+                "Zn": 30,
+                "Ga": 31,
+                "Ge": 32,
+                "As": 33,
+                "Se": 34,
+                "Br": 35,
+                "Kr": 36,
+                "Rb": 37,
+                "Sr": 38,
+                "Y": 39,
+                "Zr": 40,
+                "Nb": 41,
+                "Mo": 42,
+                "Tc": 43,
+                "Ru": 44,
+                "Rh": 45,
+                "Pd": 46,
+                "Ag": 47,
+                "Cd": 48,
+                "In": 49,
+                "Sn": 50,
+                "Sb": 51,
+                "Te": 52,
+                "I": 53,
+                "Xe": 54,
+                "Cs": 55,
+                "Ba": 56,
+                "La": 57,
+                "Ce": 58,
+                "Pr": 59,
+                "Nd": 60,
+                "Pm": 61,
+                "Sm": 62,
+                "Eu": 63,
+                "Gd": 64,
+                "Tb": 65,
+                "Dy": 66,
+                "Ho": 67,
+                "Er": 68,
+                "Tm": 69,
+                "Yb": 70,
+                "Lu": 71,
+                "Hf": 72,
+                "Ta": 73,
+                "W": 74,
+                "Re": 75,
+                "Os": 76,
+                "Ir": 77,
+                "Pt": 78,
+                "Au": 79,
+                "Hg": 80,
+                "Tl": 81,
+                "Pb": 82,
+                "Bi": 83,
+                "Po": 84,
+                "At": 85,
+                "Rn": 86,
+                "Fr": 87,
+                "Ra": 88,
+                "Ac": 89,
+                "Th": 90,
+                "Pa": 91,
+                "U": 92,
+                "Np": 93,
+                "Pu": 94,
+                "Am": 95,
+                "Cm": 96,
+                "Bk": 97,
+                "Cf": 98,
+                "Es": 99,
+                "Fm": 100,
+                "Md": 101,
+                "No": 102,
+                "Lr": 103,
+                "Rf": 104,
+                "Db": 105,
+                "Sg": 106,
+                "Bh": 107,
+                "Hs": 108,
+                "Mt": 109,
+                "Ds": 110,
+                "Rg": 111,
+                "Cn": 112,
+                "Uut": 113,
+                "Fl": 114,
+                "Uup": 115,
+                "Lv": 116,
+                "Uus": 117,
+                "Uuo": 118}

scripts/pseudo/put_pseudo_in_ezfio.py

@@ -0,0 +1,338 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+"""
+Create the pseudo potential for a given atom
+
+Usage:
+    put_pseudo_in_ezfio.py <ezfio_path>
+
+Help:
+    atom is the Abreviation of the atom
+"""
+
+
+import os
+import sys
+from docopt import docopt
+
+from subprocess import Popen, PIPE
+
+qpackage_root = os.environ['QPACKAGE_ROOT']
+
+EZFIO = "{0}/EZFIO".format(qpackage_root)
+sys.path = [EZFIO + "/Python"] + sys.path
+
+from ezfio import ezfio
+
+import re
+p = re.compile(ur'\|(\d+)><\d+\|')
+
+
+def get_pseudo_str(l_atom):
+    """
+    Run EMSL_local for geting the str of the speudo potential
+
+    str_ele :
+        Element Symbol: Na
+        Number of replaced protons: 10
+        Number of projectors: 2
+
+        Pseudopotential data:
+
+        Local component:
+        Coeff.      r^n Exp.
+        1.00000000  -1  5.35838717
+        5.35838717  1   3.67918975
+        -2.07764789 0   1.60507673
+
+        Non-local component:
+        Coeff.      r^n Exp.        Proj.
+        10.69640234 0   1.32389367  |0><0|
+        10.11238853 0   1.14052020  |1><1|
+    """
+
+    EMSL_root = "{0}/EMSL_Basis/".format(qpackage_root)
+    EMSL_path = "{0}/EMSL_api.py".format(EMSL_root)
+    db_path = "{0}/db/Pseudo.db".format(EMSL_root)
+
+    str_ = ""
+
+    for a in l_atom:
+        l_cmd_atom = ["--atom", a]
+
+        l_cmd_head = [EMSL_path, "get_basis_data",
+                      "--db_path", db_path,
+                      "--basis", "BFD-Pseudo"]
+
+        process = Popen(l_cmd_head + l_cmd_atom, stdout=PIPE, stderr=PIPE)
+
+        stdout, _ = process.communicate()
+        str_ += stdout.strip() + "\n"
+
+    return str_
+
+
+def get_v_n_dz_local(str_ele):
+    """
+    From a str_ele of the pseudo (aka only one ele in the str)
+    get the list ussefull for the Local potential : v_k n_k and dz_k
+    """
+    l_v_k = []
+    l_n_k = []
+    l_dz_k = []
+
+    for l in str_ele.splitlines():
+        try:
+            v, n, dz = l.split()
+            v = float(v)
+            n = int(n)
+            dz = float(dz)
+        except ValueError:
+            pass
+        else:
+            l_v_k.append(v)
+            l_n_k.append(n)
+            l_dz_k.append(dz)
+
+    return l_v_k, l_n_k, l_dz_k
+
+
+def get_v_n_dz_l_nonlocal(str_ele):
+    """
+    From a str_ele of the pseudo (aka only one ele in the str)
+    get the list ussefull for the non Local potential
+         v_kl (v, l)
+         n_k (v, l)
+        dz_k (dz ,l)
+    """
+    l_v_kl = []
+    l_n_kl = []
+    l_dz_kl = []
+
+    for l in str_ele.splitlines():
+        try:
+            v, n, dz, proj = l.split()
+            v = float(v)
+            n = int(n)
+            dz = float(dz)
+            l = int(p.match(proj).group(1))
+
+        except ValueError:
+            pass
+        else:
+            l_v_kl.append([v])
+            l_n_kl.append([n])
+            l_dz_kl.append([dz])
+
+    if not l_v_kl:
+        l_v_kl.append([0.])
+        l_n_kl.append([0])
+        l_dz_kl.append([0.])
+
+    return l_v_kl, l_n_kl, l_dz_kl
+
+
+def get_zeff_alpha_beta(str_ele):
+    """
+    Return the the zeff, alpha num elec and beta num elec
+        Assert ezfio_set_file alredy defined
+    """
+
+    import re
+
+    # ___
+    #  |  ._  o _|_
+    # _|_ | | |  |_
+    #
+
+    # ~#~#~#~#~#~#~ #
+    # s t r _ e l e #
+    # ~#~#~#~#~#~#~ #
+
+#    m = re.search('Element Symbol: ([a-zA-Z]+)', str_ele)
+#    name = m.group(1).capitalize()
+    name = str_ele.split("\n")[0].strip().capitalize()
+
+    m = re.search('Number of replaced protons: (\d+)', str_ele)
+    z_remove = int(m.group(1))
+
+    #  _
+    # |_) _. ._ _  _
+    # |  (_| | _> (/_
+    #
+
+    from elts_num_ele import name_to_elec
+    z = name_to_elec[name]
+
+    z_eff = z - z_remove
+
+    alpha = (z_remove / 2)
+    beta = (z_remove / 2)
+
+    #  _
+    # |_)  _ _|_     ._ ._
+    # | \ (/_ |_ |_| |  | |
+    #
+
+    return [z_eff, alpha, beta]
+
+
+def add_zero(array, size, type):
+    for add in xrange(len(array), size):
+        array.append([type(0)])
+
+    return array
+
+
+def make_it_square(matrix, dim, type=float):
+    """
+    matix the matrix to squate
+    dim array  [lmax, kmax]
+    type the null value you want
+        [[[28.59107316], [19.37583724]], [[50.25646328]]]
+            =>
+        [[[28.59107316], [19.37583724]], [[50.25646328], [0.0]]]
+    """
+
+    lmax = dim[0]
+    kmax = dim[1]
+
+    for l_list in matrix:
+
+        l_list = add_zero(l_list, lmax, type)
+
+        for k_list in list_:
+            k_list = add_zero(k_list, kmax, type)
+
+    return matrix
+
+if __name__ == "__main__":
+    arguments = docopt(__doc__)
+    # ___
+    #  |  ._  o _|_
+    # _|_ | | |  |_
+    #
+
+    # ~#~#~#~#~ #
+    # E Z F I O #
+    # ~#~#~#~#~ #
+
+    ezfio_path = arguments["<ezfio_path>"]
+    ezfio_path = os.path.expanduser(ezfio_path)
+    ezfio_path = os.path.expandvars(ezfio_path)
+    ezfio_path = os.path.abspath(ezfio_path)
+
+    ezfio.set_file("{0}".format(ezfio_path))
+
+    # ~#~#~#~#~#~#~#~#~#~#~ #
+    # P s e u d o _ d a t a #
+    # ~#~#~#~#~#~#~#~#~#~#~ #
+
+    l_ele = ezfio.get_nuclei_nucl_label()
+    str_ = get_pseudo_str(l_ele)
+
+    #  _
+    # |_) _. ._ _  _
+    # |  (_| | _> (/_
+    #
+
+    l_str_ele = [str_ele for str_ele in str_.split("Element Symbol: ")
+                 if str_ele]
+
+    for i in "l_zeff v_k n_k dz_k v_kl n_kl dz_kl".split():
+        exec("{0} = []".format(i))
+
+    alpha_tot = 0
+    beta_tot = 0
+
+    for str_ele in l_str_ele:
+
+        # ~#~#~#~#~ #
+        # S p l i t #
+        # ~#~#~#~#~ #
+
+        l = str_ele.find("Local component:")
+        nl = str_ele.find("Non-local component")
+
+        # ~#~#~#~#~ #
+        # L o c a l #
+        # ~#~#~#~#~ #
+
+        l_v, l_n, l_dz = get_v_n_dz_local(str_ele[l:nl])
+
+        v_k.append(l_v)
+        n_k.append(l_n)
+        dz_k.append(l_dz)
+
+        # ~#~#~#~#~#~#~#~#~ #
+        # N o n _ L o c a l #
+        # ~#~#~#~#~#~#~#~#~ #
+
+        l_v_kl, l_n_kl, l_dz_kl = get_v_n_dz_l_nonlocal(str_ele[nl:])
+
+        v_kl.append(l_v_kl)
+        n_kl.append(l_n_kl)
+        dz_kl.append(l_dz_kl)
+
+        # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
+        # Z _ e f f , a l p h a / b e t a _ e l e c #
+        # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
+
+        zeff, alpha, beta = get_zeff_alpha_beta(str_ele)
+
+        alpha_tot += alpha
+        beta_tot += beta
+        l_zeff.append(zeff)
+    #                                  _
+    #  /\   _|  _|   _|_  _     _  _ _|_ o  _
+    # /--\ (_| (_|    |_ (_)   (/_ /_ |  | (_)
+    #
+
+    # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
+    # Z _ e f f , a l p h a / b e t a _ e l e c #
+    # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ #
+
+    ezfio.nuclei_nucl_charge = l_zeff
+
+    alpha_tot = ezfio.get_electrons_elec_alpha_num() - alpha_tot
+    beta_tot = ezfio.get_electrons_elec_beta_num() - beta_tot
+
+    ezfio.electrons_elec_alpha_num = alpha_tot
+    ezfio.electrons_elec_beta_num = beta_tot
+
+    # Change all the array 'cause EZFIO
+    #   v_kl (v, l) => v_kl(l,v)
+    #    v_kl => zip(*_v_kl)
+    # [[7.0, 79.74474797, -49.45159098], [1.0, 5.41040609, -4.60151975]]
+    # [(7.0, 1.0), (79.74474797, 5.41040609), (-49.45159098, -4.60151975)]
+
+    # ~#~#~#~#~ #
+    # L o c a l #
+    # ~#~#~#~#~ #
+
+    klocmax = max([len(i) for i in v_k])
+    ezfio.pseudo_integrals_klocmax = klocmax
+
+    ezfio.pseudo_integrals_v_k = zip(*v_k)
+    ezfio.pseudo_integrals_n_k = zip(*n_k)
+    ezfio.pseudo_integrals_dz_k = zip(*dz_k)
+
+    # ~#~#~#~#~#~#~#~#~ #
+    # N o n _ L o c a l #
+    # ~#~#~#~#~#~#~#~#~ #
+
+    lmax = max([len(i) for i in v_kl])
+    kmax = max([len(sublist) for list_ in v_kl for sublist in list_])
+
+    ezfio.pseudo_integrals_lmaxpo = lmax
+    ezfio.pseudo_integrals_kmax = kmax
+
+    v_kl = make_it_square(v_kl, [lmax, kmax])
+    n_kl = make_it_square(n_kl, [lmax, kmax], int)
+    dz_kl = make_it_square(dz_kl, [lmax, kmax])
+
+    ezfio.pseudo_integrals_v_kl = zip(*v_kl)
+    ezfio.pseudo_integrals_n_kl = zip(*n_kl)
+    ezfio.pseudo_integrals_dz_kl = zip(*dz_kl)
+
+    ezfio.pseudo_integrals_do_pseudo = True

scripts/qp_include.sh

@@ -35,9 +35,10 @@ function check_current_dir_is_module()
      exit -1
   fi
 }
-if [[ -f NEEDED_MODULES ]]
+
+if [[ -f NEEDED_CHILDREN_MODULES ]]
 then
-  NEEDED_MODULES=$(cat NEEDED_MODULES)
+  NEEDED_MODULES=$(module_handler.py print_genealogy NEEDED_CHILDREN_MODULES)
 fi
 
 # List of executables in the current directory

scripts/qp_utils.py

@@ -0,0 +1,17 @@
+from functools import wraps
+
+def cache(func):
+    """
+    A decorator for lazy evaluation off true function
+    """
+    saved = {}
+
+    @wraps(func)
+    def newfunc(*args):
+        if args in saved:
+            return saved[args]
+
+        result = func(*args)
+        saved[args] = result
+        return result
+    return newfunc

scripts/run_Makefile_common.sh

@@ -14,7 +14,7 @@ check_current_dir_is_module
 
 # Check if the NEEDED_MODULES file is consistent
 INCLUDE_DIRS="${NEEDED_MODULES} include"
-NEEDED_MODULES_OK=$( ${QPACKAGE_ROOT}/scripts/check_dependencies.sh ${NEEDED_MODULES} )
+NEEDED_MODULES_OK=$(module_handler.py check_dependencies ${NEEDED_MODULES} )
 if [[ $? -ne 0 ]] 
 then
   error "
@@ -28,7 +28,7 @@ fi
 # Check if README.rst exists
 if [[ ! -f README.rst ]] 
 then
-  ${QPACKAGE_ROOT}/scripts/create_rst_templates.sh
+  ${QPACKAGE_ROOT}/scripts/module/create_rst_templates.sh
   error "
 README.rst was not present, so I created a
 default one for you.
@@ -62,7 +62,7 @@ then
 fi
 
 # Update Makefile.depend
-${QPACKAGE_ROOT}/scripts/create_Makefile_depend.sh
+${QPACKAGE_ROOT}/scripts/module/create_Makefile_depend.sh
 
 # Update EZFIO interface
-  ${QPACKAGE_ROOT}/scripts/ezfio_interface/ei_handler.py
+${QPACKAGE_ROOT}/scripts/ezfio_interface/ei_handler.py

scripts/run_Makefile_global.sh

@@ -52,7 +52,7 @@ ${IRPF90_VERSION}.
 
 IRPF90 version >= ${IRPF90_REQUIRED_VERSION} is required.
 To upgrade IRPF90, run :                       
-  ${QPACKAGE_ROOT}/scripts/upgrade_irpf90.sh   
+  ${QPACKAGE_ROOT}/scripts/upgrade/upgrade_irpf90.sh   
 "
 else
   info "irpf90 version is OK"
@@ -75,7 +75,7 @@ then
 Current EZFIO version : ${EZFIO_VERSION}
 EZFIO version >= ${EZFIO_REQUIRED_VERSION} is required.
 To upgrade EZFIO, run :                       
-  ${QPACKAGE_ROOT}/scripts/upgrade_ezfio.sh   
+  ${QPACKAGE_ROOT}/scripts/upgrade/upgrade_ezfio.sh   
 "
 else
   info "EZFIO version is OK"

scripts/update_README.py

@@ -83,7 +83,7 @@ def update_needed(data):
     """Read the NEEDED_MODULES file, and replace the data with it.
     Create the links to the GitHub pages."""
 
-    file = open('NEEDED_MODULES', 'r')
+    file = open('NEEDED_CHILDREN_MODULES', 'r')
     modules = file.read()
     file.close()
 

scripts/upgrade/upgrade_ezfio.sh

@@ -12,7 +12,7 @@ fi
 cd -- ${QPACKAGE_ROOT}
 mv -- ${QPACKAGE_ROOT}/EZFIO ${QPACKAGE_ROOT}/EZFIO.old
 
-${QPACKAGE_ROOT}/scripts/install_ezfio.sh
+${QPACKAGE_ROOT}/scripts/install/install_ezfio.sh
 
 if [[ $? -eq 0 ]]
 then

scripts/upgrade/upgrade_irpf90.sh

@@ -12,7 +12,7 @@ fi
 cd -- ${QPACKAGE_ROOT}
 mv -f -- ${QPACKAGE_ROOT}/irpf90 ${QPACKAGE_ROOT}/irpf90.old
 
-${QPACKAGE_ROOT}/scripts/install_irpf90.sh
+${QPACKAGE_ROOT}/scripts/install/install_irpf90.sh
 
 if [[ $? -eq 0 ]]
 then

setup_environment.sh

@@ -12,13 +12,15 @@ QPACKAGE_ROOT="$( cd "$(dirname "$BASH_SOURCE")"  ; pwd -P )"
 
 
 cat << EOF > quantum_package.rc
-export IRPF90="${IRPF90}"
 export QPACKAGE_ROOT=\$( cd \$(dirname "\${BASH_SOURCE}")  ; pwd -P )
+export IRPF90="\${QPACKAGE_ROOT}/bin/irpf90"
 export LD_LIBRARY_PATH="\${QPACKAGE_ROOT}"/lib:\${LD_LIBRARY_PATH}
 export LIBRARY_PATH="\${QPACKAGE_ROOT}"/lib:\${LIBRARY_PATH}
 export C_INCLUDE_PATH="\${QPACKAGE_ROOT}"/include:\${C_INCLUDE_PATH}
-export PYTHONPATH=\${PYTHONPATH}:"\${QPACKAGE_ROOT}"/scripts:"\${QPACKAGE_ROOT}"/scripts/ezfio_interface
+
-export PATH=\${PATH}:"\${QPACKAGE_ROOT}"/scripts:"\${QPACKAGE_ROOT}"/scripts/ezfio_interface
+export PYTHONPATH=\${PYTHONPATH}\$(find "\${QPACKAGE_ROOT}"/scripts -type d -printf ":%p")
+
+export PATH=\${PATH}\$(find "\${QPACKAGE_ROOT}"/scripts -type d -printf ":%p")
 export PATH=\${PATH}:"\${QPACKAGE_ROOT}"/bin
 export PATH=\${PATH}:"\${QPACKAGE_ROOT}"/ocaml
 source "\${QPACKAGE_ROOT}"/bin/irpman &> /dev/null
@@ -28,37 +30,29 @@ EOF
 source quantum_package.rc
 
 echo "${BLUE}===== Installing IRPF90 ===== ${BLACK}"
-${QPACKAGE_ROOT}/scripts/install_irpf90.sh     | tee install_irpf90.log
+${QPACKAGE_ROOT}/scripts/install/install_irpf90.sh     | tee ${QPACKAGE_ROOT}/install_logs/install_irpf90.log
-if [[ ! -d ${QPACKAGE_ROOT}/irpf90 ]]
+if [[ ! -d ${QPACKAGE_ROOT}/irpf90 ]]  ||  [[ ! -x ${QPACKAGE_ROOT}/bin/irpf90 ]] ||  [[ ! -x ${QPACKAGE_ROOT}/bin/irpman ]]
-then
-  echo $RED "Error in IRPF90 installation" $BLACK
-  exit 1
-fi
-
-if [[ ! -x ${QPACKAGE_ROOT}/bin/irpf90 ]]
-then
-  echo $RED "Error in IRPF90 installation" $BLACK
-  exit 1
-fi
-
-if [[ ! -x ${QPACKAGE_ROOT}/bin/irpman ]]
 then
   echo $RED "Error in IRPF90 installation" $BLACK
   exit 1
 fi
 
+mkdir -p ${QPACKAGE_ROOT}/install_logs
 
 echo "${BLUE}===== Installing Zlib ===== ${BLACK}"
-${QPACKAGE_ROOT}/scripts/install_zlib.sh       | tee install_zlib.log
+${QPACKAGE_ROOT}/scripts/install/install_zlib.sh       | tee ${QPACKAGE_ROOT}/install_logs/install_zlib.log
 
 echo "${BLUE}===== Installing Curl ===== ${BLACK}"
-${QPACKAGE_ROOT}/scripts/install_curl.sh       | tee install_curl.log
+${QPACKAGE_ROOT}/scripts/install/install_curl.sh       | tee ${QPACKAGE_ROOT}/install_logs/install_curl.log
 
 echo "${BLUE}===== Installing M4 ===== ${BLACK}"
-${QPACKAGE_ROOT}/scripts/install_m4.sh         | tee install_m4.log
+${QPACKAGE_ROOT}/scripts/install/install_m4.sh         | tee ${QPACKAGE_ROOT}/install_logs/install_m4.log
+
+echo "${BLUE}===== Installing Docopt ===== ${BLACK}"
+${QPACKAGE_ROOT}/scripts/install/install_docopt.sh     | tee ${QPACKAGE_ROOT}/install_logs/install_docopt.log
 
 echo "${BLUE}===== Installing EMSL Basis set library ===== ${BLACK}"
-${QPACKAGE_ROOT}/scripts/install_emsl.sh       | tee install_emsl.log
+${QPACKAGE_ROOT}/scripts/install/install_emsl.sh       | tee ${QPACKAGE_ROOT}/install_logs/install_emsl.log
 
 if [[ ! -d ${QPACKAGE_ROOT}/EMSL_Basis ]]
 then
@@ -68,7 +62,7 @@ fi
 
 echo "${BLUE}===== Installing EZFIO ===== ${BLACK}"
 
-${QPACKAGE_ROOT}/scripts/install_ezfio.sh | tee install_ezfio.log
+${QPACKAGE_ROOT}/scripts/install/install_ezfio.sh | tee ${QPACKAGE_ROOT}/install_logs/install_ezfio.log
 if [[ ! -d ${QPACKAGE_ROOT}/EZFIO ]]
 then
   echo $RED "Error in EZFIO installation" $BLACK
@@ -78,7 +72,7 @@ fi
 
 echo "${BLUE}===== Installing Ocaml compiler and libraries ===== ${BLACK}"
 rm -f -- ocaml/Qptypes.ml
-${QPACKAGE_ROOT}/scripts/install_ocaml.sh | tee install_ocaml.log
+${QPACKAGE_ROOT}/scripts/install/install_ocaml.sh | tee ${QPACKAGE_ROOT}/install_logs/install_ocaml.log
 
 if [[ ! -f ${QPACKAGE_ROOT}/ocaml/Qptypes.ml ]]
 then
@@ -87,7 +81,7 @@ then
 fi
 
 echo "${BLUE}===== Installing resultsFile Python library ===== ${BLACK}"
-${QPACKAGE_ROOT}/scripts/install_resultsFile.sh
+${QPACKAGE_ROOT}/scripts/install/install_resultsFile.sh
 if [[ ! -d ${QPACKAGE_ROOT}/resultsFile ]]
 then
   echo $RED "Error in resultsFile installation" $BLACK
@@ -106,8 +100,7 @@ source ${QPACKAGE_ROOT}/quantum_package.rc
 " $BLACK
 
 
-mkdir -p ${QPACKAGE_ROOT}/install_logs
+
-mv ${QPACKAGE_ROOT}/*.log ${QPACKAGE_ROOT}/install_logs/
 
 if [[ $1 == "--robot" ]] ; 
 then

src/AOs/ASSUMPTIONS.rst

@@ -1,8 +1 @@
-* The atomic orbitals are normalized:
-
-  .. math::
-
-   \int \left(\chi_i({\bf r}) \right)^2 dr = 1
-
 * The AO coefficients in the EZFIO files are not necessarily normalized and are normalized after reading
-* The AO coefficients and exponents are ordered in increasing order of exponents

src/AOs/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Nuclei

src/AOs/NEEDED_MODULES

@@ -1,2 +0,0 @@
-Ezfio_files Nuclei Output Utils
-

src/AOs/README.rst

@@ -17,20 +17,20 @@ The AO coefficients are normalized as:
 
   {\tilde c}_{ki} = \frac{c_{ki}}{ \int \left( (x-X_A)^a (y-Y_A)^b (z-Z_A)^c  e^{-\gamma_{ki} |{\bf r} - {\bf R}_A|^2} \right)^2} dr
 
+Warning: ``ao_coef`` contains the AO coefficients given in input. These do not
+include the normalization constant of the AO. The ``ao_coef_normalized`` includes
+this normalization factor.
+
+The AOs are also sorted by increasing exponent to accelerate the calculation of
+the two electron integrals.
+
 Assumptions
 ===========
 
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* The atomic orbitals are normalized:
-
-  .. math::
-
-   \int \left(\chi_i({\bf r}) \right)^2 dr = 1
-
 * The AO coefficients in the EZFIO files are not necessarily normalized and are normalized after reading
-* The AO coefficients and exponents are ordered in increasing order of exponents
 
 
 Needed Modules
@@ -39,10 +39,7 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
 * `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 
 Documentation
 =============
@@ -70,45 +67,41 @@ Documentation
   Overlap between atomic basis functions:
   :math:`\int \chi_i(r) \chi_j(r) dr)`
 
-`ao_coef <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L21>`_
+`ao_coef <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L62>`_
-  Coefficients, exponents and powers of x,y and z
+  AO Coefficients, read from input. Those should not be used directly, as
+  the MOs are expressed on the basis of **normalized** AOs.
 
-`ao_coef_transp <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L157>`_
+`ao_coef_normalized <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L84>`_
-  Transposed ao_coef and ao_expo
+  Coefficients including the AO normalization
 
-`ao_coef_unnormalized <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L116>`_
+`ao_coef_normalized_ordered <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L107>`_
-  Coefficients, exponents and powers of x,y and z as in the EZFIO file
+  Sorted primitives to accelerate 4 index MO transformation
-  ao_coef(i,j) = coefficient of the jth primitive on the ith ao
+
+`ao_coef_normalized_ordered_transp <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L133>`_
+  Transposed ao_coef_normalized_ordered
+
+`ao_expo <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L41>`_
+  AO Exponents read from input
+
+`ao_expo_ordered <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L108>`_
+  Sorted primitives to accelerate 4 index MO transformation
+
+`ao_expo_ordered_transp <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L147>`_
+  Transposed ao_expo_ordered
+
+`ao_l <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L162>`_
   ao_l = l value of the AO: a+b+c in x^a y^b z^c
 
-`ao_expo <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L20>`_
+`ao_l_char <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L163>`_
-  Coefficients, exponents and powers of x,y and z
-
-`ao_expo_transp <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L158>`_
-  Transposed ao_coef and ao_expo
-
-`ao_expo_unsorted <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L117>`_
-  Coefficients, exponents and powers of x,y and z as in the EZFIO file
-  ao_coef(i,j) = coefficient of the jth primitive on the ith ao
   ao_l = l value of the AO: a+b+c in x^a y^b z^c
 
-`ao_l <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L118>`_
+`ao_l_char_space <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L311>`_
-  Coefficients, exponents and powers of x,y and z as in the EZFIO file
-  ao_coef(i,j) = coefficient of the jth primitive on the ith ao
-  ao_l = l value of the AO: a+b+c in x^a y^b z^c
-
-`ao_l_char <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L119>`_
-  Coefficients, exponents and powers of x,y and z as in the EZFIO file
-  ao_coef(i,j) = coefficient of the jth primitive on the ith ao
-  ao_l = l value of the AO: a+b+c in x^a y^b z^c
-
-`ao_l_char_space <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L309>`_
   Undocumented
 
-`ao_md5 <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L400>`_
+`ao_md5 <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L403>`_
   MD5 key characteristic of the AO basis
 
-`ao_nucl <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L207>`_
+`ao_nucl <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L209>`_
   Index of the nuclei on which the ao is centered
 
 `ao_num <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L1>`_
@@ -118,35 +111,35 @@ Documentation
   Number of atomic orbitals
 
 `ao_power <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L19>`_
-  Coefficients, exponents and powers of x,y and z
+  Powers of x,y and z read from input
 
-`ao_prim_num <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L175>`_
+`ao_prim_num <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L177>`_
   Number of primitives per atomic orbital
 
-`ao_prim_num_max <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L197>`_
+`ao_prim_num_max <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L199>`_
   Undocumented
 
-`ao_prim_num_max_align <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L198>`_
+`ao_prim_num_max_align <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L200>`_
   Undocumented
 
-`l_to_charater <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L216>`_
+`l_to_charater <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L218>`_
   character corresponding to the "L" value of an AO orbital
 
-`n_aos_max <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L229>`_
+`n_aos_max <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L231>`_
   Number of AOs per atom
 
-`nucl_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L242>`_
+`nucl_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L244>`_
   List of AOs attached on each atom
 
-`nucl_list_shell_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L260>`_
+`nucl_list_shell_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L262>`_
   Index of the shell type Aos and of the corresponding Aos
   Per convention, for P,D,F and G AOs, we take the index
   of the AO with the the corresponding power in the "X" axis
 
-`nucl_n_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L228>`_
+`nucl_n_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L230>`_
   Number of AOs per atom
 
-`nucl_num_shell_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L261>`_
+`nucl_num_shell_aos <http://github.com/LCPQ/quantum_package/tree/master/src/AOs/aos.irp.f#L263>`_
   Index of the shell type Aos and of the corresponding Aos
   Per convention, for P,D,F and G AOs, we take the index
   of the AO with the the corresponding power in the "X" axis

src/AOs/ao_overlap.irp.f

@@ -21,8 +21,8 @@
   !$OMP  overlap_x,overlap_y, overlap_z, overlap, &
   !$OMP  alpha, beta,i,j,c) &
   !$OMP SHARED(nucl_coord,ao_power,ao_prim_num, &
-  !$OMP  ao_overlap_x,ao_overlap_y,ao_overlap_z,ao_overlap,ao_num,ao_coef_transp,ao_nucl, &
+  !$OMP  ao_overlap_x,ao_overlap_y,ao_overlap_z,ao_overlap,ao_num,ao_coef_normalized_ordered_transp,ao_nucl, &
-  !$OMP  ao_expo_transp,dim1)
+  !$OMP  ao_expo_ordered_transp,dim1)
   do j=1,ao_num
    A_center(1) = nucl_coord( ao_nucl(j), 1 )
    A_center(2) = nucl_coord( ao_nucl(j), 2 )
@@ -44,12 +44,12 @@
     power_B(2)  = ao_power( i, 2 )
     power_B(3)  = ao_power( i, 3 )
     do n = 1,ao_prim_num(j)
-     alpha = ao_expo_transp(n,j)
+     alpha = ao_expo_ordered_transp(n,j)
      !DEC$ VECTOR ALIGNED
      do l = 1, ao_prim_num(i)
-      beta = ao_expo_transp(l,i)
+      beta = ao_expo_ordered_transp(l,i)
       call overlap_gaussian_xyz(A_center,B_center,alpha,beta,power_A,power_B,overlap_x,overlap_y,overlap_z,overlap,dim1)
-      c = ao_coef_transp(n,j) * ao_coef_transp(l,i)
+      c = ao_coef_normalized_ordered_transp(n,j) * ao_coef_normalized_ordered_transp(l,i)
       ao_overlap(i,j) += c * overlap
       ao_overlap_x(i,j) += c * overlap_x
       ao_overlap_y(i,j) += c * overlap_y
@@ -84,8 +84,8 @@ BEGIN_PROVIDER [ double precision, ao_overlap_abs,(ao_num_align,ao_num) ]
   !$OMP  overlap_x,overlap_y, overlap_z, overlap, &
   !$OMP  alpha, beta,i,j,dx) &
   !$OMP SHARED(nucl_coord,ao_power,ao_prim_num, &
-  !$OMP  ao_overlap_abs,ao_num,ao_coef_transp,ao_nucl, &
+  !$OMP  ao_overlap_abs,ao_num,ao_coef_normalized_ordered_transp,ao_nucl, &
-  !$OMP  ao_expo_transp,dim1,lower_exp_val)
+  !$OMP  ao_expo_ordered_transp,dim1,lower_exp_val)
   do j=1,ao_num
    A_center(1) = nucl_coord( ao_nucl(j), 1 )
    A_center(2) = nucl_coord( ao_nucl(j), 2 )
@@ -104,14 +104,14 @@ BEGIN_PROVIDER [ double precision, ao_overlap_abs,(ao_num_align,ao_num) ]
     power_B(2)  = ao_power( i, 2 )
     power_B(3)  = ao_power( i, 3 )
     do n = 1,ao_prim_num(j)
-     alpha = ao_expo_transp(n,j)
+     alpha = ao_expo_ordered_transp(n,j)
      !DEC$ VECTOR ALIGNED
      do l = 1, ao_prim_num(i)
-      beta = ao_expo_transp(l,i)
+      beta = ao_expo_ordered_transp(l,i)
       call overlap_x_abs(A_center(1),B_center(1),alpha,beta,power_A(1),power_B(1),overlap_x,lower_exp_val,dx,dim1)
       call overlap_x_abs(A_center(2),B_center(2),alpha,beta,power_A(2),power_B(2),overlap_y,lower_exp_val,dx,dim1)
       call overlap_x_abs(A_center(3),B_center(3),alpha,beta,power_A(3),power_B(3),overlap_z,lower_exp_val,dx,dim1)
-      ao_overlap_abs(i,j) += abs(ao_coef_transp(n,j) * ao_coef_transp(l,i)) * overlap_x * overlap_y * overlap_z
+      ao_overlap_abs(i,j) += abs(ao_coef_normalized_ordered_transp(n,j) * ao_coef_normalized_ordered_transp(l,i)) * overlap_x * overlap_y * overlap_z
      enddo
     enddo
    enddo

src/AOs/aos.irp.f

@@ -1,152 +1,171 @@
  BEGIN_PROVIDER [ integer, ao_num ]
 &BEGIN_PROVIDER [ integer, ao_num_align ]
- implicit none
+   implicit none
    
- BEGIN_DOC
+   BEGIN_DOC
-! Number of atomic orbitals
+   ! Number of atomic orbitals
- END_DOC
+   END_DOC
    
- ao_num = -1
+   ao_num = -1
- PROVIDE ezfio_filename
+   PROVIDE ezfio_filename
- call ezfio_get_ao_basis_ao_num(ao_num)
+   call ezfio_get_ao_basis_ao_num(ao_num)
- if (ao_num <= 0) then
+   if (ao_num <= 0) then
-   stop 'Number of contracted gaussians should be > 0'
+     stop 'Number of contracted gaussians should be > 0'
- endif
+   endif
- integer :: align_double
+   integer                        :: align_double
- ao_num_align = align_double(ao_num)
+   ao_num_align = align_double(ao_num)
 END_PROVIDER
  
- BEGIN_PROVIDER [ integer, ao_power, (ao_num_align,3) ]
+BEGIN_PROVIDER [ integer, ao_power, (ao_num_align,3) ]
-&BEGIN_PROVIDER [ double precision, ao_expo, (ao_num_align,ao_prim_num_max) ]
+   implicit none
-&BEGIN_PROVIDER [ double precision, ao_coef, (ao_num_align,ao_prim_num_max) ]
+   BEGIN_DOC
- implicit none
+   ! Powers of x,y and z read from input
-
+   END_DOC
- BEGIN_DOC
+   PROVIDE ezfio_filename
-! Coefficients, exponents and powers of x,y and z
- END_DOC
- PROVIDE ezfio_filename
-
- double precision, allocatable :: buffer(:,:)
- allocate ( buffer(ao_num,ao_prim_num_max) )
- integer :: ibuffer(ao_num,3)
- integer :: i,j,k
- character*(128) :: give_ao_character_space
- ibuffer = 0
- call ezfio_get_ao_basis_ao_power(ibuffer)
- ao_power = 0
- do j = 1, 3
-  do i = 1, ao_num
-   ao_power(i,j) = ibuffer(i,j)
-  enddo
- enddo
- ao_expo  = 0.d0
- buffer = 0.d0
- call ezfio_get_ao_basis_ao_expo(buffer)
- do j = 1, ao_prim_num_max
-  do i = 1, ao_num
-   ao_expo(i,j) = buffer(i,j)
-  enddo
- enddo
-
- ao_coef  = 0.d0
- buffer = 0.d0
- call ezfio_get_ao_basis_ao_coef(buffer)
- do j = 1, ao_prim_num_max
-  do i = 1, ao_num
-   ao_coef(i,j) = buffer(i,j)
-  enddo
- enddo
-
- deallocate(buffer)
-
-! Normalization of the AO coefficients
-! ------------------------------------
- double precision :: norm, norm2,overlap_x,overlap_y,overlap_z,C_A(3)
- integer :: l, powA(3), nz
- nz=100
- C_A(1) = 0.d0
- C_A(2) = 0.d0
- C_A(3) = 0.d0
- do i=1,ao_num
-  powA(1) = ao_power(i,1)
-  powA(2) = ao_power(i,2)
-  powA(3) = ao_power(i,3)
-  do j=1,ao_prim_num(i)
-   call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,norm,nz)
-   ao_coef(i,j) = ao_coef(i,j)/sqrt(norm)
-  enddo
- enddo
-
- ! Sorting of the exponents for efficient integral calculations
- integer :: iorder(ao_prim_num_max)
- double precision :: d(ao_prim_num_max,2)
- do i=1,ao_num
-  do j=1,ao_prim_num(i)
-   iorder(j) = j
-   d(j,1) = ao_expo(i,j)
-   d(j,2) = ao_coef(i,j)
-  enddo
-  call dsort(d(1,1),iorder,ao_prim_num(i))
-  call dset_order(d(1,2),iorder,ao_prim_num(i))
-  do j=1,ao_prim_num(i)
-   ao_expo(i,j) = d(j,1)
-   ao_coef(i,j) = d(j,2)
-  enddo
- enddo
-END_PROVIDER
-
-
- BEGIN_PROVIDER [ double precision, ao_coef_transp, (ao_prim_num_max_align,ao_num) ]
-&BEGIN_PROVIDER [ double precision, ao_expo_transp, (ao_prim_num_max_align,ao_num) ]
- implicit none
- BEGIN_DOC
-! Transposed ao_coef and ao_expo
- END_DOC
- integer :: i,j
- do j=1, ao_num
-  do i=1, ao_prim_num_max
-   ao_coef_transp(i,j) = ao_coef(j,i)
-   ao_expo_transp(i,j) = ao_expo(j,i)
-  enddo
- enddo
    
+   integer                        :: i,j,k
+   integer, allocatable           :: ibuffer(:,:)
+   allocate ( ibuffer(ao_num,3) )
+   ibuffer = 0
+   call ezfio_get_ao_basis_ao_power(ibuffer)
+   ao_power = 0
+   do j = 1, 3
+     do i = 1, ao_num
+       ao_power(i,j) = ibuffer(i,j)
+     enddo
+   enddo
+   deallocate(ibuffer)
    
 END_PROVIDER
 
- BEGIN_PROVIDER [ double precision, ao_coef_unnormalized, (ao_num_align,ao_prim_num_max) ]
+BEGIN_PROVIDER [ double precision, ao_expo, (ao_num_align,ao_prim_num_max) ]
-&BEGIN_PROVIDER [ double precision, ao_expo_unsorted, (ao_num_align,ao_prim_num_max) ]
+  implicit none
-&BEGIN_PROVIDER [ integer, ao_l, (ao_num) ]
+  BEGIN_DOC
+  ! AO Exponents read from input
+  END_DOC
+  PROVIDE ezfio_filename
+  
+  double precision, allocatable  :: buffer(:,:)
+  allocate ( buffer(ao_num,ao_prim_num_max) )
+  integer                        :: i,j,k
+  ao_expo  = 0.d0
+  buffer = 0.d0
+  call ezfio_get_ao_basis_ao_expo(buffer)
+  do j = 1, ao_prim_num_max
+    do i = 1, ao_num
+      ao_expo(i,j) = buffer(i,j)
+    enddo
+  enddo
+  deallocate(buffer)
+END_PROVIDER
+
+BEGIN_PROVIDER [ double precision, ao_coef, (ao_num_align,ao_prim_num_max) ]
+  implicit none
+  BEGIN_DOC
+  ! AO Coefficients, read from input. Those should not be used directly, as
+  ! the MOs are expressed on the basis of **normalized** AOs.
+  END_DOC
+  PROVIDE ezfio_filename
+  
+  double precision, allocatable  :: buffer(:,:)
+  allocate ( buffer(ao_num,ao_prim_num_max) )
+  integer                        :: i,j,k
+  ao_coef  = 0.d0
+  buffer = 0.d0
+  call ezfio_get_ao_basis_ao_coef(buffer)
+  do j = 1, ao_prim_num_max
+    do i = 1, ao_num
+      ao_coef(i,j) = buffer(i,j)
+    enddo
+  enddo
+  deallocate(buffer)
+END_PROVIDER
+
+BEGIN_PROVIDER [ double precision, ao_coef_normalized, (ao_num_align,ao_prim_num_max) ]
+  implicit none
+  BEGIN_DOC
+  ! Coefficients including the AO normalization
+  END_DOC
+  double precision               :: norm, norm2,overlap_x,overlap_y,overlap_z,C_A(3)
+  integer                        :: l, powA(3), nz
+  integer                        :: i,j
+  nz=100
+  C_A(1) = 0.d0
+  C_A(2) = 0.d0
+  C_A(3) = 0.d0
+  do i=1,ao_num
+    powA(1) = ao_power(i,1)
+    powA(2) = ao_power(i,2)
+    powA(3) = ao_power(i,3)
+    do j=1,ao_prim_num(i)
+      call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,norm,nz)
+      ao_coef_normalized(i,j) = ao_coef(i,j)/sqrt(norm)
+    enddo
+  enddo
+END_PROVIDER
+
+ BEGIN_PROVIDER [ double precision, ao_coef_normalized_ordered, (ao_num_align,ao_prim_num_max) ]
+&BEGIN_PROVIDER [ double precision, ao_expo_ordered, (ao_num_align,ao_prim_num_max) ]
+   implicit none
+   BEGIN_DOC
+   ! Sorted primitives to accelerate 4 index MO transformation
+   END_DOC
+   
+   integer                        :: iorder(ao_prim_num_max)
+   double precision               :: d(ao_prim_num_max,2)
+   integer                        :: i,j
+   do i=1,ao_num
+     do j=1,ao_prim_num(i)
+       iorder(j) = j
+       d(j,1) = ao_expo(i,j)
+       d(j,2) = ao_coef_normalized(i,j)
+     enddo
+     call dsort(d(1,1),iorder,ao_prim_num(i))
+     call dset_order(d(1,2),iorder,ao_prim_num(i))
+     do j=1,ao_prim_num(i)
+       ao_expo_ordered(i,j) = d(j,1)
+       ao_coef_normalized_ordered(i,j) = d(j,2)
+     enddo
+   enddo
+END_PROVIDER
+
+
+BEGIN_PROVIDER [ double precision, ao_coef_normalized_ordered_transp, (ao_prim_num_max_align,ao_num) ]
+  implicit none
+  BEGIN_DOC
+  ! Transposed ao_coef_normalized_ordered
+  END_DOC
+  integer                        :: i,j
+  do j=1, ao_num
+    do i=1, ao_prim_num_max
+      ao_coef_normalized_ordered_transp(i,j) = ao_coef_normalized_ordered(j,i)
+    enddo
+  enddo
+  
+END_PROVIDER
+
+BEGIN_PROVIDER [ double precision, ao_expo_ordered_transp, (ao_prim_num_max_align,ao_num) ]
+  implicit none
+  BEGIN_DOC
+  ! Transposed ao_expo_ordered
+  END_DOC
+  integer                        :: i,j
+  do j=1, ao_num
+    do i=1, ao_prim_num_max
+      ao_expo_ordered_transp(i,j) = ao_expo_ordered(j,i)
+    enddo
+  enddo
+  
+END_PROVIDER
+
+
+ BEGIN_PROVIDER [ integer, ao_l, (ao_num) ]
 &BEGIN_PROVIDER [ character*(128), ao_l_char, (ao_num) ]
  implicit none
-
  BEGIN_DOC
-! Coefficients, exponents and powers of x,y and z as in the EZFIO file
-! ao_coef(i,j) = coefficient of the jth primitive on the ith ao
 ! ao_l = l value of the AO: a+b+c in x^a y^b z^c
  END_DOC
- PROVIDE ezfio_filename
+ integer :: i
-
- double precision, allocatable :: buffer(:,:)
- allocate ( buffer(ao_num,ao_prim_num_max) )
- integer :: i,j,k
- character*(128) :: give_ao_character_space
- buffer = 0.d0
- call ezfio_get_ao_basis_ao_expo(buffer)
- do j = 1, ao_prim_num_max
-  do i = 1, ao_num
-   ao_expo_unsorted(i,j) = buffer(i,j)
-  enddo
- enddo
-
- buffer = 0.d0
- call ezfio_get_ao_basis_ao_coef(buffer)
- do j = 1, ao_prim_num_max
-  do i = 1, ao_num
-   ao_coef_unnormalized(i,j) = buffer(i,j)
-  enddo
- enddo
- deallocate(buffer)
-
  do i=1,ao_num
    ao_l(i) = ao_power(i,1) + ao_power(i,2) + ao_power(i,3) 
    ao_l_char(i) = l_to_charater(ao_l(i))
@@ -154,23 +173,6 @@ END_PROVIDER
 END_PROVIDER
 
 
- BEGIN_PROVIDER [ double precision, ao_coef_transp, (ao_prim_num_max_align,ao_num) ]
-&BEGIN_PROVIDER [ double precision, ao_expo_transp, (ao_prim_num_max_align,ao_num) ]
- implicit none
- BEGIN_DOC
-! Transposed ao_coef and ao_expo
- END_DOC
- integer :: i,j
- do j=1, ao_num
-  do i=1, ao_prim_num_max
-   ao_coef_transp(i,j) = ao_coef(j,i)
-   ao_expo_transp(i,j) = ao_expo(j,i)
-  enddo
- enddo
-
-
-END_PROVIDER
-
 
 BEGIN_PROVIDER [ integer, ao_prim_num, (ao_num_align) ]
  implicit none
@@ -303,10 +305,10 @@ END_PROVIDER
   enddo
  enddo
 
- END_PROVIDER
+END_PROVIDER
 
 
- BEGIN_PROVIDER [ character*(4), ao_l_char_space, (ao_num) ]
+BEGIN_PROVIDER [ character*(4), ao_l_char_space, (ao_num) ]
  implicit none
  integer :: i
  character*(4) :: give_ao_character_space
@@ -397,6 +399,7 @@ END_PROVIDER
    ao_l_char_space(i) = give_ao_character_space
  enddo
 END_PROVIDER
+
 BEGIN_PROVIDER [ character*(32), ao_md5 ]
  BEGIN_DOC
 ! MD5 key characteristic of the AO basis

src/Bielec_integrals/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+MonoInts Bitmask

src/Bielec_integrals/NEEDED_MODULES

@@ -1 +0,0 @@
-AOs Bitmask Electrons Ezfio_files MOs Nuclei Output Utils MonoInts

src/Bielec_integrals/README.rst

@@ -16,15 +16,8 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 * `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
+* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
 
 Documentation
 =============
@@ -36,7 +29,7 @@ Documentation
   integral of the AO basis <ik|jl> or (ij|kl)
   i(r1) j(r1) 1/r12 k(r2) l(r2)
 
-`ao_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L490>`_
+`ao_bielec_integral_schwartz <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L491>`_
   Needed to compute Schwartz inequalities
 
 `ao_bielec_integral_schwartz_accel <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L107>`_
@@ -53,48 +46,48 @@ Documentation
 `compute_ao_bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L290>`_
   Compute AO 1/r12 integrals for all i and fixed j,k,l
 
-`eri <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L653>`_
+`eri <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L654>`_
   ATOMIC PRIMTIVE bielectronic integral between the 4 primitives ::
   primitive_1 = x1**(a_x) y1**(a_y) z1**(a_z) exp(-alpha * r1**2)
   primitive_2 = x1**(b_x) y1**(b_y) z1**(b_z) exp(- beta * r1**2)
   primitive_3 = x2**(c_x) y2**(c_y) z2**(c_z) exp(-delta * r2**2)
   primitive_4 = x2**(d_x) y2**(d_y) z2**(d_z) exp(- gama * r2**2)
 
-`general_primitive_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L515>`_
+`general_primitive_integral <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L516>`_
   Computes the integral <pq|rs> where p,q,r,s are Gaussian primitives
 
-`give_polynom_mult_center_x <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L851>`_
+`give_polynom_mult_center_x <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L852>`_
   subroutine that returns the explicit polynom in term of the "t"
   variable of the following polynomw :
   I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q)
 
-`i_x1_new <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L772>`_
+`i_x1_new <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L773>`_
   recursive function involved in the bielectronic integral
 
-`i_x1_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L914>`_
+`i_x1_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L915>`_
   recursive function involved in the bielectronic integral
 
-`i_x1_pol_mult_a1 <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L1034>`_
+`i_x1_pol_mult_a1 <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L1035>`_
   recursive function involved in the bielectronic integral
 
-`i_x1_pol_mult_a2 <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L1088>`_
+`i_x1_pol_mult_a2 <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L1089>`_
   recursive function involved in the bielectronic integral
 
-`i_x1_pol_mult_recurs <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L948>`_
+`i_x1_pol_mult_recurs <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L949>`_
   recursive function involved in the bielectronic integral
 
-`i_x2_new <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L807>`_
+`i_x2_new <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L808>`_
   recursive function involved in the bielectronic integral
 
-`i_x2_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L1150>`_
+`i_x2_pol_mult <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L1151>`_
   recursive function involved in the bielectronic integral
 
-`integrale_new <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L698>`_
+`integrale_new <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L699>`_
   calculate the integral of the polynom ::
   I_x1(a_x+b_x, c_x+d_x,p,q) * I_x1(a_y+b_y, c_y+d_y,p,q) * I_x1(a_z+b_z, c_z+d_z,p,q)
   between ( 0 ; 1)
 
-`n_pt_sup <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L837>`_
+`n_pt_sup <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/ao_bi_integrals.irp.f#L838>`_
   Returns the upper boundary of the degree of the polynomial involved in the
   bielctronic integral :
   Ix(a_x,b_x,c_x,d_x) * Iy(a_y,b_y,c_y,d_y) * Iz(a_z,b_z,c_z,d_z)
@@ -172,32 +165,32 @@ Documentation
 `add_integrals_to_map <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L42>`_
   Adds integrals to tha MO map according to some bitmask
 
-`mo_bielec_integral_jj <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L464>`_
+`mo_bielec_integral_jj <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L465>`_
   mo_bielec_integral_jj(i,j) = J_ij
   mo_bielec_integral_jj_exchange(i,j) = K_ij
   mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij
 
-`mo_bielec_integral_jj_anti <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L466>`_
+`mo_bielec_integral_jj_anti <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L467>`_
   mo_bielec_integral_jj(i,j) = J_ij
   mo_bielec_integral_jj_exchange(i,j) = K_ij
   mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij
 
-`mo_bielec_integral_jj_anti_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L326>`_
+`mo_bielec_integral_jj_anti_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L327>`_
   mo_bielec_integral_jj_from_ao(i,j) = J_ij
   mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
   mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij
 
-`mo_bielec_integral_jj_exchange <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L465>`_
+`mo_bielec_integral_jj_exchange <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L466>`_
   mo_bielec_integral_jj(i,j) = J_ij
   mo_bielec_integral_jj_exchange(i,j) = K_ij
   mo_bielec_integral_jj_anti(i,j) = J_ij - K_ij
 
-`mo_bielec_integral_jj_exchange_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L325>`_
+`mo_bielec_integral_jj_exchange_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L326>`_
   mo_bielec_integral_jj_from_ao(i,j) = J_ij
   mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
   mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij
 
-`mo_bielec_integral_jj_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L324>`_
+`mo_bielec_integral_jj_from_ao <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals/mo_bi_integrals.irp.f#L325>`_
   mo_bielec_integral_jj_from_ao(i,j) = J_ij
   mo_bielec_integral_jj_exchange_from_ao(i,j) = J_ij
   mo_bielec_integral_jj_anti_from_ao(i,j) = J_ij - K_ij

src/Bielec_integrals/ao_bi_integrals.irp.f

@@ -42,24 +42,24 @@ double precision function ao_bielec_integral(i,j,k,l)
     
     do p = 1, ao_prim_num(i)
       double precision               :: coef1
-      coef1 = ao_coef_transp(p,i)
+      coef1 = ao_coef_normalized_ordered_transp(p,i)
       do q = 1, ao_prim_num(j)
         double precision               :: coef2
-        coef2 = coef1*ao_coef_transp(q,j)
+        coef2 = coef1*ao_coef_normalized_ordered_transp(q,j)
         double precision               :: p_inv,q_inv
         call give_explicit_poly_and_gaussian(P_new,P_center,pp,fact_p,iorder_p,&
-            ao_expo_transp(p,i),ao_expo_transp(q,j),                 &
+            ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j),                 &
             I_power,J_power,I_center,J_center,dim1)
         p_inv = 1.d0/pp
         do r = 1, ao_prim_num(k)
           double precision               :: coef3
-          coef3 = coef2*ao_coef_transp(r,k)
+          coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
           do s = 1, ao_prim_num(l)
             double precision               :: coef4
-            coef4 = coef3*ao_coef_transp(s,l)
+            coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
             double precision               :: general_primitive_integral
             call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q,&
-                ao_expo_transp(r,k),ao_expo_transp(s,l),             &
+                ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),             &
                 K_power,L_power,K_center,L_center,dim1)
             q_inv = 1.d0/qq
             integral = general_primitive_integral(dim1,              &
@@ -82,15 +82,15 @@ double precision function ao_bielec_integral(i,j,k,l)
     double  precision              :: ERI
 
     do p = 1, ao_prim_num(i)
-      coef1 = ao_coef_transp(p,i)
+      coef1 = ao_coef_normalized_ordered_transp(p,i)
       do q = 1, ao_prim_num(j)
-        coef2 = coef1*ao_coef_transp(q,j)
+        coef2 = coef1*ao_coef_normalized_ordered_transp(q,j)
         do r = 1, ao_prim_num(k)
-          coef3 = coef2*ao_coef_transp(r,k)
+          coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
           do s = 1, ao_prim_num(l)
-            coef4 = coef3*ao_coef_transp(s,l)
+            coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
             integral = ERI(                                          &
-                ao_expo_transp(p,i),ao_expo_transp(q,j),ao_expo_transp(r,k),ao_expo_transp(s,l),&
+                ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j),ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),&
                 I_power(1),J_power(1),K_power(1),L_power(1),         &
                 I_power(2),J_power(2),K_power(2),L_power(2),         &
                 I_power(3),J_power(3),K_power(3),L_power(3))
@@ -149,12 +149,12 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
     
     schwartz_kl(0,0) = 0.d0
     do r = 1, ao_prim_num(k)
-      coef1 = ao_coef_transp(r,k)*ao_coef_transp(r,k)
+      coef1 = ao_coef_normalized_ordered_transp(r,k)*ao_coef_normalized_ordered_transp(r,k)
       schwartz_kl(0,r) = 0.d0
       do s = 1, ao_prim_num(l)
-        coef2 = coef1 * ao_coef_transp(s,l) * ao_coef_transp(s,l)
+        coef2 = coef1 * ao_coef_normalized_ordered_transp(s,l) * ao_coef_normalized_ordered_transp(s,l)
         call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q,&
-            ao_expo_transp(r,k),ao_expo_transp(s,l),                 &
+            ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),                 &
             K_power,L_power,K_center,L_center,dim1)
         q_inv = 1.d0/qq
         schwartz_kl(s,r) = general_primitive_integral(dim1,          &
@@ -168,13 +168,13 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
 
     do p = 1, ao_prim_num(i)
       double precision               :: coef1
-      coef1 = ao_coef_transp(p,i)
+      coef1 = ao_coef_normalized_ordered_transp(p,i)
       do q = 1, ao_prim_num(j)
         double precision               :: coef2
-        coef2 = coef1*ao_coef_transp(q,j)
+        coef2 = coef1*ao_coef_normalized_ordered_transp(q,j)
         double precision               :: p_inv,q_inv
         call give_explicit_poly_and_gaussian(P_new,P_center,pp,fact_p,iorder_p,&
-            ao_expo_transp(p,i),ao_expo_transp(q,j),                 &
+            ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j),                 &
             I_power,J_power,I_center,J_center,dim1)
         p_inv = 1.d0/pp
         schwartz_ij = general_primitive_integral(dim1,               &
@@ -189,16 +189,16 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
              cycle
           endif
           double precision               :: coef3
-          coef3 = coef2*ao_coef_transp(r,k)
+          coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
           do s = 1, ao_prim_num(l)
             double precision               :: coef4
             if (schwartz_kl(s,r)*schwartz_ij < thresh) then
                cycle
             endif
-            coef4 = coef3*ao_coef_transp(s,l)
+            coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
             double precision               :: general_primitive_integral
             call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q,&
-                ao_expo_transp(r,k),ao_expo_transp(s,l),             &
+                ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),             &
                 K_power,L_power,K_center,L_center,dim1)
             q_inv = 1.d0/qq
             integral = general_primitive_integral(dim1,              &
@@ -222,12 +222,12 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
 
     schwartz_kl(0,0) = 0.d0
     do r = 1, ao_prim_num(k)
-      coef1 = ao_coef_transp(r,k)*ao_coef_transp(r,k)
+      coef1 = ao_coef_normalized_ordered_transp(r,k)*ao_coef_normalized_ordered_transp(r,k)
       schwartz_kl(0,r) = 0.d0
       do s = 1, ao_prim_num(l)
-        coef2 = coef1*ao_coef_transp(s,l)*ao_coef_transp(s,l)
+        coef2 = coef1*ao_coef_normalized_ordered_transp(s,l)*ao_coef_normalized_ordered_transp(s,l)
         schwartz_kl(s,r) = ERI(                                      &
-            ao_expo_transp(r,k),ao_expo_transp(s,l),ao_expo_transp(r,k),ao_expo_transp(s,l),&
+            ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),&
             K_power(1),L_power(1),K_power(1),L_power(1),             &
             K_power(2),L_power(2),K_power(2),L_power(2),             &
             K_power(3),L_power(3),K_power(3),L_power(3)) * &
@@ -238,11 +238,11 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
     enddo
 
     do p = 1, ao_prim_num(i)
-      coef1 = ao_coef_transp(p,i)
+      coef1 = ao_coef_normalized_ordered_transp(p,i)
       do q = 1, ao_prim_num(j)
-        coef2 = coef1*ao_coef_transp(q,j)
+        coef2 = coef1*ao_coef_normalized_ordered_transp(q,j)
         schwartz_ij = ERI(                                          &
-                ao_expo_transp(p,i),ao_expo_transp(q,j),ao_expo_transp(p,i),ao_expo_transp(q,j),&
+                ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j),ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j),&
                 I_power(1),J_power(1),I_power(1),J_power(1),         &
                 I_power(2),J_power(2),I_power(2),J_power(2),         &
                 I_power(3),J_power(3),I_power(3),J_power(3))*coef2*coef2
@@ -253,14 +253,14 @@ double precision function ao_bielec_integral_schwartz_accel(i,j,k,l)
           if (schwartz_kl(0,r)*schwartz_ij < thresh) then
              cycle
           endif
-          coef3 = coef2*ao_coef_transp(r,k)
+          coef3 = coef2*ao_coef_normalized_ordered_transp(r,k)
           do s = 1, ao_prim_num(l)
             if (schwartz_kl(s,r)*schwartz_ij < thresh) then
                cycle
             endif
-            coef4 = coef3*ao_coef_transp(s,l)
+            coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)
             integral = ERI(                                          &
-                ao_expo_transp(p,i),ao_expo_transp(q,j),ao_expo_transp(r,k),ao_expo_transp(s,l),&
+                ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j),ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l),&
                 I_power(1),J_power(1),K_power(1),L_power(1),         &
                 I_power(2),J_power(2),K_power(2),L_power(2),         &
                 I_power(3),J_power(3),K_power(3),L_power(3))
@@ -481,6 +481,7 @@ IRP_ENDIF COARRAY
   ao_bielec_integrals_in_map = .True.
   if (write_ao_integrals) then
     call dump_ao_integrals(trim(ezfio_filename)//'/work/ao_integrals.bin')
+    call ezfio_set_bielec_integrals_disk_access_ao_integrals(.True.)
   endif
   
 END_PROVIDER

src/Bielec_integrals/mo_bi_integrals.irp.f

@@ -312,6 +312,7 @@ IRP_ENDIF
   
   if (write_mo_integrals) then
     call dump_mo_integrals(trim(ezfio_filename)//'/work/mo_integrals.bin')
+    call ezfio_set_bielec_integrals_disk_access_mo_integrals(.True.)
   endif
   
   

src/Bitmask/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+MOs

src/Bitmask/NEEDED_MODULES

@@ -1 +0,0 @@
-AOs Electrons Ezfio_files MOs Nuclei Output Utils

src/Bitmask/README.rst

@@ -40,13 +40,7 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
 * `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 
 Documentation
 =============
@@ -134,6 +128,10 @@ Documentation
   Subroutine to print the content of a determinant in '+-' notation and
   hexadecimal representation.
 
+`debug_spindet <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks_routines.irp.f#L155>`_
+  Subroutine to print the content of a determinant in '+-' notation and
+  hexadecimal representation.
+
 `list_to_bitstring <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks_routines.irp.f#L29>`_
   Returns the physical string "string(N_int,2)" from the array of
   occupations "list(N_int*bit_kind_size,2)
@@ -141,5 +139,8 @@ Documentation
 `print_det <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks_routines.irp.f#L138>`_
   Subroutine to print the content of a determinant using the '+-' notation
 
+`print_spindet <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask/bitmasks_routines.irp.f#L171>`_
+  Subroutine to print the content of a determinant using the '+-' notation
+
 
 

src/Bitmask/bitmasks_routines.irp.f

@@ -126,7 +126,7 @@ subroutine debug_det(string,Nint)
   END_DOC
   integer, intent(in)            :: Nint
   integer(bit_kind), intent(in)  :: string(Nint,2)
-  character*(512)                :: output(2)
+  character*(2048)                :: output(2)
   call bitstring_to_hexa( output(1), string(1,1), Nint )
   call bitstring_to_hexa( output(2), string(1,2), Nint )
   print *,  trim(output(1)) , '|', trim(output(2))
@@ -143,7 +143,7 @@ subroutine print_det(string,Nint)
   END_DOC
   integer, intent(in)            :: Nint
   integer(bit_kind), intent(in)  :: string(Nint,2)
-  character*(512)                :: output(2)
+  character*(2048)                :: output(2)
 
   call bitstring_to_str( output(1), string(1,1), Nint )
   call bitstring_to_str( output(2), string(1,2), Nint )
@@ -151,3 +151,34 @@ subroutine print_det(string,Nint)
   print *,  trim(output(2))
 
 end
+
+subroutine debug_spindet(string,Nint)
+  use bitmasks
+  implicit none
+  BEGIN_DOC
+  ! Subroutine to print the content of a determinant in '+-' notation and
+  ! hexadecimal representation.
+  END_DOC
+  integer, intent(in)            :: Nint
+  integer(bit_kind), intent(in)  :: string(Nint,2)
+  character*(2048)                :: output(1)
+  call bitstring_to_hexa( output(1), string(1,1), Nint )
+  print *,  trim(output(1)) 
+  call print_spindet(string,Nint)
+
+end
+
+subroutine print_spindet(string,Nint)
+  use bitmasks
+  implicit none
+  BEGIN_DOC
+  ! Subroutine to print the content of a determinant using the '+-' notation
+  END_DOC
+  integer, intent(in)            :: Nint
+  integer(bit_kind), intent(in)  :: string(Nint,2)
+  character*(2048)                :: output(1)
+
+  call bitstring_to_str( output(1), string(1,1), Nint )
+  print *,  trim(output(1))
+
+end

src/CAS_SD/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Perturbation Selectors_full Generators_CAS

src/CAS_SD/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask Determinants Electrons Ezfio_files Generators_CAS Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full Utils
-

src/CAS_SD/README.rst

@@ -24,21 +24,7 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
-* `Properties <http://github.com/LCPQ/quantum_package/tree/master/src/Properties>`_
 * `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
+* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
 

src/CAS_SD/cas_sd.irp.f

@@ -1,6 +1,7 @@
 program full_ci
   implicit none
   integer                        :: i,k
+  integer                        :: N_det_old
 
   
   double precision, allocatable  :: pt2(:), norm_pert(:), H_pert_diag(:)
@@ -9,6 +10,7 @@ program full_ci
   allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
   character*(64)                 :: perturbation
   
+  N_det_old = 0
   pt2 = 1.d0
   diag_algorithm = "Lapack"
   if (N_det > n_det_max_cas_sd) then
@@ -29,6 +31,7 @@ program full_ci
   endif
 
   do while (N_det < n_det_max_cas_sd.and.maxval(abs(pt2(1:N_st))) > pt2_max)
+    N_det_old = N_det
     call H_apply_CAS_SD(pt2, norm_pert, H_pert_diag,  N_st)
 
     PROVIDE  psi_coef
@@ -53,10 +56,11 @@ program full_ci
     if (abort_all) then
       exit
     endif
+    if (N_det == N_det_old) then
+      exit
+    endif
   enddo
 
-  ! Check that it is a CAS-SD
-  logical :: in_cas
   integer :: exc_max, degree_min
   exc_max = 0
   print *,  'CAS determinants : ', N_det_generators
@@ -69,18 +73,4 @@ program full_ci
     print *,  ''
   enddo
   print *,  'Max excitation degree in the CAS :', exc_max
-  do i=1,N_det
-    in_cas = .False.
-    degree_min = 1000
-    do k=1,N_det_generators
-      call get_excitation_degree(psi_det_generators(1,1,k),psi_det(1,1,i),degree,N_int)
-      degree_min = min(degree_min,degree)
-    enddo
-    if (degree_min > 2) then
-        print *,  'Error : This is not a CAS-SD : '
-        print *,  'Excited determinant:', degree_min
-        call debug_det(psi_det(1,1,k),N_int)
-        stop
-    endif
-  enddo
 end

src/CAS_SD/cas_sd_selected.irp.f

@@ -50,14 +50,30 @@ program full_ci
     print *,  'E        = ', CI_energy
     print *,  'E+PT2    = ', CI_energy+pt2
     print *,  '-----'
-    call ezfio_set_full_ci_energy(CI_energy)
+    call ezfio_set_cas_sd_energy(CI_energy(1))
     if (abort_all) then
       exit
     endif
   enddo
+  call diagonalize_CI
+
+   if(do_pt2_end)then
+    print*,'Last iteration only to compute the PT2'
+    threshold_selectors = 1.d0
+    threshold_generators = 0.999d0
+    call H_apply_CAS_SD_PT2(pt2, norm_pert, H_pert_diag,  N_st)
+
+    print *,  'Final step'
+    print *,  'N_det    = ', N_det
+    print *,  'N_states = ', N_states
+    print *,  'PT2      = ', pt2
+    print *,  'E        = ', CI_energy
+    print *,  'E+PT2    = ', CI_energy+pt2
+    print *,  '-----'
+    call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
+   endif
+
 
-  ! Check that it is a CAS-SD
-  logical :: in_cas
   integer :: exc_max, degree_min
   exc_max = 0
   print *,  'CAS determinants : ', N_det_cas
@@ -70,18 +86,4 @@ program full_ci
     print *,  ''
   enddo
   print *,  'Max excitation degree in the CAS :', exc_max
-  do i=1,N_det
-    in_cas = .False.
-    degree_min = 1000
-    do k=1,N_det_cas
-      call get_excitation_degree(psi_cas(1,1,k),psi_det(1,1,i),degree,N_int)
-      degree_min = min(degree_min,degree)
-    enddo
-    if (degree_min > 2) then
-        print *,  'Error : This is not a CAS-SD : '
-        print *,  'Excited determinant:', degree_min
-        call debug_det(psi_det(1,1,k),N_int)
-        stop
-    endif
-  enddo
 end

src/CID/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Selectors_full SingleRefMethod

src/CID/NEEDED_MODULES

@@ -1,3 +0,0 @@
-AOs Bielec_integrals Bitmask Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Selectors_full SingleRefMethod Utils
-
-

src/CID/README.rst

@@ -15,21 +15,8 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
 * `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 
 Documentation
 =============

src/CID_SC2_selected/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+CID_selected

src/CID_SC2_selected/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask CISD CISD_selected Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full SingleRefMethod Utils
-

src/CID_SC2_selected/README.rst

@@ -19,23 +19,5 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
+* `CID_selected <http://github.com/LCPQ/quantum_package/tree/master/src/CID_selected>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `CISD <http://github.com/LCPQ/quantum_package/tree/master/src/CISD>`_
-* `CISD_selected <http://github.com/LCPQ/quantum_package/tree/master/src/CISD_selected>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
-* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
-* `Properties <http://github.com/LCPQ/quantum_package/tree/master/src/Properties>`_
-* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
-* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 

src/CID_selected/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Perturbation CID

src/CID_selected/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask CISD Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full SingleRefMethod Utils
-

src/CID_selected/README.rst

@@ -22,22 +22,6 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `CISD <http://github.com/LCPQ/quantum_package/tree/master/src/CISD>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
-* `Properties <http://github.com/LCPQ/quantum_package/tree/master/src/Properties>`_
+* `CID <http://github.com/LCPQ/quantum_package/tree/master/src/CID>`_
-* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
-* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 

src/CIS/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Selectors_full SingleRefMethod

src/CIS/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Selectors_full SingleRefMethod Utils
-

src/CIS/README.rst

@@ -31,19 +31,6 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
 * `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 

src/CISD/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Selectors_full SingleRefMethod

src/CISD/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Selectors_full SingleRefMethod Utils
-

src/CISD/README.rst

@@ -15,21 +15,8 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
 * `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 
 Documentation
 =============

src/CISD_SC2_selected/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+CISD_selected

src/CISD_SC2_selected/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask CISD CISD_selected Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full SingleRefMethod Utils
-

src/CISD_SC2_selected/README.rst

@@ -19,23 +19,5 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `CISD <http://github.com/LCPQ/quantum_package/tree/master/src/CISD>`_
 * `CISD_selected <http://github.com/LCPQ/quantum_package/tree/master/src/CISD_selected>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
-* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
-* `Properties <http://github.com/LCPQ/quantum_package/tree/master/src/Properties>`_
-* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
-* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 

src/CISD_selected/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Perturbation CISD

src/CISD_selected/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask CISD Determinants Electrons Ezfio_files Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full SingleRefMethod Utils
-

src/CISD_selected/README.rst

@@ -28,22 +28,6 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `CISD <http://github.com/LCPQ/quantum_package/tree/master/src/CISD>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
-* `Properties <http://github.com/LCPQ/quantum_package/tree/master/src/Properties>`_
+* `CISD <http://github.com/LCPQ/quantum_package/tree/master/src/CISD>`_
-* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
-* `SingleRefMethod <http://github.com/LCPQ/quantum_package/tree/master/src/SingleRefMethod>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 

src/DDCI_selected/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Perturbation Selectors_full Generators_CAS

src/DDCI_selected/NEEDED_MODULES

@@ -1,2 +0,0 @@
-AOs Bielec_integrals Bitmask Determinants Electrons Ezfio_files Generators_CAS Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full Utils
-

src/DDCI_selected/README.rst

@@ -19,21 +19,7 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
-* `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Determinants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
-* `Hartree_Fock <http://github.com/LCPQ/quantum_package/tree/master/src/Hartree_Fock>`_
-* `MOGuess <http://github.com/LCPQ/quantum_package/tree/master/src/MOGuess>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
 * `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
-* `Properties <http://github.com/LCPQ/quantum_package/tree/master/src/Properties>`_
 * `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
+* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
 

src/Determinants/H_apply.irp.f

@@ -53,13 +53,17 @@ subroutine resize_H_apply_buffer(new_size,iproc)
   double precision,  pointer     :: buffer_e2(:,:)
   integer                        :: i,j,k
   integer                        :: Ndet
+
+  BEGIN_DOC
+! Resizes the H_apply buffer of proc iproc. The buffer lock should
+! be set before calling this function.
+  END_DOC
   PROVIDE H_apply_buffer_allocated
   
   ASSERT (new_size > 0)
   ASSERT (iproc >= 0)
   ASSERT (iproc < nproc)
   
-  call omp_set_lock(H_apply_buffer_lock(1,iproc))
   allocate ( buffer_det(N_int,2,new_size),                           &
       buffer_coef(new_size,N_states),                                &
       buffer_e2(new_size,N_states) )
@@ -93,7 +97,6 @@ subroutine resize_H_apply_buffer(new_size,iproc)
   
   H_apply_buffer(iproc)%sze = new_size
   H_apply_buffer(iproc)%N_det = min(new_size,H_apply_buffer(iproc)%N_det)
-  call omp_unset_lock(H_apply_buffer_lock(1,iproc))
   
 end
 
@@ -101,8 +104,7 @@ subroutine copy_H_apply_buffer_to_wf
   use omp_lib
   implicit none
   BEGIN_DOC
-! Copies the H_apply buffer to psi_coef. You need to touch psi_det, psi_coef and N_det
+! Copies the H_apply buffer to psi_coef.
-! after calling this function.
 ! After calling this subroutine, N_det, psi_det and psi_coef need to be touched
   END_DOC
   integer(bit_kind), allocatable :: buffer_det(:,:,:)
@@ -181,8 +183,77 @@ subroutine copy_H_apply_buffer_to_wf
   call normalize(psi_coef,N_det)
   SOFT_TOUCH N_det psi_det psi_coef
   
+  logical :: found_duplicates
+  call remove_duplicates_in_psi_det(found_duplicates)
 end
 
+subroutine remove_duplicates_in_psi_det(found_duplicates)
+  implicit none
+  logical, intent(out) :: found_duplicates
+  BEGIN_DOC
+! Removes duplicate determinants in the wave function.
+  END_DOC
+  integer                        :: i,j,k
+  integer(bit_kind), allocatable :: bit_tmp(:)
+  logical,allocatable            :: duplicate(:)
+
+  allocate (duplicate(N_det), bit_tmp(N_det))
+
+  do i=1,N_det
+    integer, external            :: det_search_key
+    !$DIR FORCEINLINE
+    bit_tmp(i) = det_search_key(psi_det_sorted_bit(1,1,i),N_int)
+    duplicate(i) = .False.
+  enddo
+
+  do i=1,N_det-1
+    if (duplicate(i)) then
+      cycle
+    endif
+    j = i+1
+    do while (bit_tmp(j)==bit_tmp(i))
+      if (duplicate(j)) then
+        j += 1
+        cycle
+      endif
+      duplicate(j) = .True.
+      do k=1,N_int
+        if ( (psi_det_sorted_bit(k,1,i) /= psi_det_sorted_bit(k,1,j) ) &
+        .or. (psi_det_sorted_bit(k,2,i) /= psi_det_sorted_bit(k,2,j) ) ) then
+          duplicate(j) = .False.
+          exit
+        endif
+      enddo
+      j += 1
+      if (j > N_det) then
+        exit
+      endif
+    enddo
+  enddo
+
+  found_duplicates = .False.
+  do i=1,N_det
+    if (duplicate(i)) then
+      found_duplicates = .True.
+      exit
+    endif
+  enddo
+
+  if (found_duplicates) then
+    call write_bool(output_determinants,found_duplicates,'Found duplicate determinants')
+    k=0
+    do i=1,N_det
+      if (.not.duplicate(i)) then
+        k += 1
+        psi_det(:,:,k) = psi_det_sorted_bit (:,:,i)
+        psi_coef(k,:)  = psi_coef_sorted_bit(i,:)
+      endif
+    enddo
+    N_det = k
+    TOUCH N_det psi_det psi_coef
+  endif
+  deallocate (duplicate,bit_tmp)
+end
 
 subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)
   use bitmasks
@@ -196,11 +267,11 @@ subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)
   integer                        :: i,j,k
   integer                        :: new_size
   PROVIDE H_apply_buffer_allocated
+  call omp_set_lock(H_apply_buffer_lock(1,iproc))
   new_size = H_apply_buffer(iproc)%N_det + n_selected
   if (new_size > H_apply_buffer(iproc)%sze) then
     call resize_h_apply_buffer(max(2*H_apply_buffer(iproc)%sze,new_size),iproc)
   endif
-  call omp_set_lock(H_apply_buffer_lock(1,iproc))
   do i=1,H_apply_buffer(iproc)%N_det
     ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,1,i)) )== elec_alpha_num)
     ASSERT (sum(popcnt(H_apply_buffer(iproc)%det(:,2,i))) == elec_beta_num)
@@ -215,7 +286,7 @@ subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)
   enddo
   do j=1,N_states
     do i=1,N_selected
-      H_apply_buffer(iproc)%coef(i,j) = 0.d0
+      H_apply_buffer(iproc)%coef(i+H_apply_buffer(iproc)%N_det,j) = 0.d0
     enddo
   enddo
   H_apply_buffer(iproc)%N_det = new_size

src/Determinants/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Bielec_integrals

src/Determinants/NEEDED_MODULES

@@ -1 +0,0 @@
-AOs Bielec_integrals Bitmask Electrons Ezfio_files MonoInts MOs Nuclei Output Utils

src/Determinants/README.rst

@@ -32,16 +32,7 @@ Needed Modules
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-* `AOs <http://github.com/LCPQ/quantum_package/tree/master/src/AOs>`_
 * `Bielec_integrals <http://github.com/LCPQ/quantum_package/tree/master/src/Bielec_integrals>`_
-* `Bitmask <http://github.com/LCPQ/quantum_package/tree/master/src/Bitmask>`_
-* `Electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Electrons>`_
-* `Ezfio_files <http://github.com/LCPQ/quantum_package/tree/master/src/Ezfio_files>`_
-* `MonoInts <http://github.com/LCPQ/quantum_package/tree/master/src/MonoInts>`_
-* `MOs <http://github.com/LCPQ/quantum_package/tree/master/src/MOs>`_
-* `Nuclei <http://github.com/LCPQ/quantum_package/tree/master/src/Nuclei>`_
-* `Output <http://github.com/LCPQ/quantum_package/tree/master/src/Output>`_
-* `Utils <http://github.com/LCPQ/quantum_package/tree/master/src/Utils>`_
 
 Documentation
 =============
@@ -49,12 +40,11 @@ Documentation
 .. Do not edit this section. It was auto-generated from the
 .. NEEDED_MODULES file.
 
-`copy_h_apply_buffer_to_wf <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L100>`_
+`copy_h_apply_buffer_to_wf <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L103>`_
-  Copies the H_apply buffer to psi_coef. You need to touch psi_det, psi_coef and N_det
+  Copies the H_apply buffer to psi_coef.
-  after calling this function.
   After calling this subroutine, N_det, psi_det and psi_coef need to be touched
 
-`fill_h_apply_buffer_no_selection <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L187>`_
+`fill_h_apply_buffer_no_selection <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L258>`_
   Fill the H_apply buffer with determiants for CISD
 
 `h_apply_buffer_allocated <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L15>`_
@@ -65,8 +55,12 @@ Documentation
   Buffer of determinants/coefficients/perturbative energy for H_apply.
   Uninitialized. Filled by H_apply subroutines.
 
+`remove_duplicates_in_psi_det <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L190>`_
+  Removes duplicate determinants in the wave function.
+
 `resize_h_apply_buffer <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/H_apply.irp.f#L48>`_
-  Undocumented
+  Resizes the H_apply buffer of proc iproc. The buffer lock should
+  be set before calling this function.
 
 `cisd_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/SC2.irp.f#L1>`_
   CISD+SC2 method              :: take off all the disconnected terms of a CISD (selected or not)
@@ -193,10 +187,10 @@ Documentation
 `det_svd <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/det_svd.irp.f#L1>`_
   Computes the SVD of the Alpha x Beta determinant coefficient matrix
 
-`filter_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L426>`_
+`filter_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L423>`_
   Returns a determinant with only the 3 highest electrons
 
-`int_of_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L391>`_
+`int_of_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L388>`_
   Returns an integer*8 as :
   .br
   |_<--- 21 bits ---><--- 21 bits ---><--- 21 bits --->|
@@ -213,26 +207,26 @@ Documentation
 `n_det <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L3>`_
   Number of determinants in the wave function
 
-`psi_average_norm_contrib <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L276>`_
+`psi_average_norm_contrib <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L273>`_
   Contribution of determinants to the state-averaged density
 
-`psi_average_norm_contrib_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L306>`_
+`psi_average_norm_contrib_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L303>`_
   Wave function sorted by determinants contribution to the norm (state-averaged)
 
-`psi_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L230>`_
+`psi_coef <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L227>`_
   The wave function coefficients. Initialized with Hartree-Fock if the EZFIO file
   is empty
 
-`psi_coef_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L305>`_
+`psi_coef_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L302>`_
   Wave function sorted by determinants contribution to the norm (state-averaged)
 
-`psi_coef_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L453>`_
+`psi_coef_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L450>`_
   Determinants on which we apply <i|H|j>.
   They are sorted by the 3 highest electrons in the alpha part,
   then by the 3 highest electrons in the beta part to accelerate
   the research of connected determinants.
 
-`psi_coef_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L336>`_
+`psi_coef_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L333>`_
   Determinants on which we apply <i|H|psi> for perturbation.
   They are sorted by determinants interpreted as integers. Useful
   to accelerate the search of a random determinant in the wave
@@ -245,46 +239,46 @@ Documentation
 `psi_det_size <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L47>`_
   Size of the psi_det/psi_coef arrays
 
-`psi_det_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L304>`_
+`psi_det_sorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L301>`_
   Wave function sorted by determinants contribution to the norm (state-averaged)
 
-`psi_det_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L452>`_
+`psi_det_sorted_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L449>`_
   Determinants on which we apply <i|H|j>.
   They are sorted by the 3 highest electrons in the alpha part,
   then by the 3 highest electrons in the beta part to accelerate
   the research of connected determinants.
 
-`psi_det_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L335>`_
+`psi_det_sorted_bit <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L332>`_
   Determinants on which we apply <i|H|psi> for perturbation.
   They are sorted by determinants interpreted as integers. Useful
   to accelerate the search of a random determinant in the wave
   function.
 
-`psi_det_sorted_next_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L454>`_
+`psi_det_sorted_next_ab <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L451>`_
   Determinants on which we apply <i|H|j>.
   They are sorted by the 3 highest electrons in the alpha part,
   then by the 3 highest electrons in the beta part to accelerate
   the research of connected determinants.
 
-`read_dets <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L583>`_
+`read_dets <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L580>`_
   Reads the determinants from the EZFIO file
 
-`save_wavefunction <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L630>`_
+`save_wavefunction <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L627>`_
   Save the wave function into the EZFIO file
 
-`save_wavefunction_general <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L649>`_
+`save_wavefunction_general <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L646>`_
   Save the wave function into the EZFIO file
 
-`save_wavefunction_unsorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L640>`_
+`save_wavefunction_unsorted <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L637>`_
   Save the wave function into the EZFIO file
 
-`sort_dets_by_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L474>`_
+`sort_dets_by_3_highest_electrons <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L471>`_
   Determinants on which we apply <i|H|j>.
   They are sorted by the 3 highest electrons in the alpha part,
   then by the 3 highest electrons in the beta part to accelerate
   the research of connected determinants.
 
-`sort_dets_by_det_search_key <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L349>`_
+`sort_dets_by_det_search_key <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/determinants.irp.f#L346>`_
   Determinants are sorted are sorted according to their det_search_key.
   Useful to accelerate the search of a random determinant in the wave
   function.
@@ -322,7 +316,7 @@ Documentation
 `diag_algorithm <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI.irp.f#L1>`_
   Diagonalization algorithm (Davidson or Lapack)
 
-`diagonalize_ci <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI.irp.f#L96>`_
+`diagonalize_ci <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI.irp.f#L100>`_
   Replace the coefficients of the CI states by the coefficients of the
   eigenstates of the CI matrix
 
@@ -335,7 +329,7 @@ Documentation
 `ci_sc2_energy <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_SC2.irp.f#L1>`_
   N_states_diag lowest eigenvalues of the CI matrix
 
-`diagonalize_ci_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_SC2.irp.f#L46>`_
+`diagonalize_ci_sc2 <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_SC2.irp.f#L45>`_
   Replace the coefficients of the CI states_diag by the coefficients of the
   eigenstates of the CI matrix
 
@@ -351,7 +345,7 @@ Documentation
 `ci_electronic_energy_mono <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_mono.irp.f#L1>`_
   Eigenvectors/values of the CI matrix
 
-`diagonalize_ci_mono <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_mono.irp.f#L59>`_
+`diagonalize_ci_mono <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/diagonalize_CI_mono.irp.f#L73>`_
   Replace the coefficients of the CI states by the coefficients of the
   eigenstates of the CI matrix
 
@@ -538,7 +532,7 @@ Documentation
 `save_dets_qmcchem <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/save_for_qmcchem.irp.f#L1>`_
   Undocumented
 
-`save_for_qmc <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/save_for_qmcchem.irp.f#L48>`_
+`save_for_qmc <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/save_for_qmcchem.irp.f#L46>`_
   Undocumented
 
 `save_natorb <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/save_natorb.irp.f#L1>`_
@@ -629,61 +623,49 @@ Documentation
 `n_con_int <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/slater_rules.irp.f#L1131>`_
   Number of integers to represent the connections between determinants
 
-`create_wf_of_psi_svd_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L473>`_
+`create_wf_of_psi_svd_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L483>`_
   Matrix of wf coefficients. Outer product of alpha and beta determinants
 
-`generate_all_alpha_beta_det_products <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L528>`_
+`generate_all_alpha_beta_det_products <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L538>`_
   Create a wave function from all possible alpha x beta determinants
 
-`get_index_in_psi_det_alpha_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L131>`_
+`get_index_in_psi_det_alpha_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L139>`_
   Returns the index of the determinant in the ``psi_det_alpha_unique`` array
 
-`get_index_in_psi_det_beta_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L212>`_
+`get_index_in_psi_det_beta_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L221>`_
   Returns the index of the determinant in the ``psi_det_beta_unique`` array
 
-`n_det_alpha_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L54>`_
-  Unique alpha determinants
-
-`n_det_beta_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L91>`_
-  Unique beta determinants
-
 `psi_det_alpha <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L25>`_
   List of alpha determinants of psi_det
 
-`psi_det_alpha_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L53>`_
-  Unique alpha determinants
-
 `psi_det_beta <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L39>`_
   List of beta determinants of psi_det
 
-`psi_det_beta_unique <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L90>`_
+`psi_svd_alpha <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L578>`_
-  Unique beta determinants
-
-`psi_svd_alpha <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L568>`_
   SVD wave function
 
-`psi_svd_beta <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L569>`_
+`psi_svd_beta <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L579>`_
   SVD wave function
 
-`psi_svd_coefs <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L570>`_
+`psi_svd_coefs <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L580>`_
   SVD wave function
 
-`psi_svd_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L457>`_
+`psi_svd_matrix <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L467>`_
   Matrix of wf coefficients. Outer product of alpha and beta determinants
 
-`psi_svd_matrix_columns <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L398>`_
+`psi_svd_matrix_columns <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L407>`_
   Matrix of wf coefficients. Outer product of alpha and beta determinants
 
-`psi_svd_matrix_rows <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L397>`_
+`psi_svd_matrix_rows <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L406>`_
   Matrix of wf coefficients. Outer product of alpha and beta determinants
 
-`psi_svd_matrix_values <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L396>`_
+`psi_svd_matrix_values <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L405>`_
   Matrix of wf coefficients. Outer product of alpha and beta determinants
 
 `spin_det_search_key <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L9>`_
   Return an integer*8 corresponding to a determinant index for searching
 
-`write_spindeterminants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L294>`_
+`write_spindeterminants <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/spindeterminants.irp.f#L303>`_
   Undocumented
 
 `cisd <http://github.com/LCPQ/quantum_package/tree/master/src/Determinants/truncate_wf.irp.f#L1>`_

src/Determinants/determinants.irp.f

@@ -147,9 +147,7 @@ END_PROVIDER
    !$DIR FORCEINLINE
    bit_tmp(i) = occ_pattern_search_key(psi_occ_pattern(1,1,i),N_int)
  enddo
- print*,'passed 1'
  call i8sort(bit_tmp,iorder,N_det)
- print*,'passed 2'
  !DIR$ IVDEP
  do i=1,N_det
   do k=1,N_int
@@ -189,7 +187,6 @@ END_PROVIDER
     endif
   enddo
  enddo
-  print*,'passed 3'
 
  N_occ_pattern=0
  do i=1,N_det

src/Determinants/diagonalize_CI.irp.f

@@ -66,28 +66,32 @@ END_PROVIDER
     enddo
     integer :: i_state
     double precision :: s2
-    i_state = 0
+    if (s2_eig) then
-    do j=1,N_det
+      i_state = 0
-      call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
+      do j=1,N_det
-      if(dabs(s2-expected_s2).le.0.3d0)then
+        call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
-       i_state += 1
+        if(dabs(s2-expected_s2).le.0.3d0)then
-       do i=1,N_det
+        i_state += 1
-         CI_eigenvectors(i,i_state) = eigenvectors(i,j)
+        do i=1,N_det
-       enddo
+          CI_eigenvectors(i,i_state) = eigenvectors(i,j)
-       CI_electronic_energy(i_state) = eigenvalues(j)
+        enddo
-       CI_eigenvectors_s2(i_state) = s2
+        CI_electronic_energy(i_state) = eigenvalues(j)
-      endif
+        CI_eigenvectors_s2(i_state) = s2
-      if (i_state.ge.N_states_diag) then
+        endif
-        exit
+        if (i_state.ge.N_states_diag) then
-      endif
+          exit
-    enddo
+        endif
-!    if(i_state < min(N_states_diag,N_det))then
+      enddo
-!     print *, 'pb with the number of states'
+    else
-!     print *, 'i_state = ',i_state
+      do j=1,N_states_diag
-!     print *, 'N_states_diag ',N_states_diag
+        call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
-!     print *,'stopping ...'
+        do i=1,N_det
-!     stop
+          CI_eigenvectors(i,j) = eigenvectors(i,j)
-!    endif
+        enddo
+        CI_electronic_energy(j) = eigenvalues(j)
+        CI_eigenvectors_s2(j) = s2
+      enddo
+    endif
     deallocate(eigenvectors,eigenvalues)
   endif
   

src/Determinants/diagonalize_CI_SC2.irp.f

@@ -35,8 +35,7 @@ END_PROVIDER
     do i=1,N_det
       CI_SC2_eigenvectors(i,j) = psi_coef(i,j)
     enddo
-! TODO : check comment
+    CI_SC2_electronic_energy(j) = CI_electronic_energy(j) 
-!   CI_SC2_electronic_energy(j) = CI_electronic_energy(j) 
   enddo
   
   call CISD_SC2(psi_det,CI_SC2_eigenvectors,CI_SC2_electronic_energy, &

src/Determinants/diagonalize_CI_mono.irp.f

@@ -32,25 +32,39 @@
     integer :: i_state
     double precision :: s2
     i_state = 0
-    do j=1,N_det
+    if (s2_eig) then
-      call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
+      do j=1,N_det
-      if(dabs(s2-expected_s2).le.0.3d0)then
+        call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
-       print*,'j = ',j
+        if(dabs(s2-expected_s2).le.0.3d0)then
-       print*,'e = ',eigenvalues(j)
+          print*,'j = ',j
-       print*,'c = ',dabs(eigenvectors(1,j))
+          print*,'e = ',eigenvalues(j)
-       if(dabs(eigenvectors(1,j)).gt.0.9d0)then
+          print*,'c = ',dabs(eigenvectors(1,j))
-        i_state += 1
+          if(dabs(eigenvectors(1,j)).gt.0.9d0)then
-        do i=1,N_det
+            i_state += 1
-          CI_eigenvectors_mono(i,i_state) = eigenvectors(i,j)
+            do i=1,N_det
-        enddo
+              CI_eigenvectors_mono(i,i_state) = eigenvectors(i,j)
-        CI_electronic_energy_mono(i_state) = eigenvalues(j)
+            enddo
-        CI_eigenvectors_s2_mono(i_state) = s2
+            CI_electronic_energy_mono(i_state) = eigenvalues(j)
-       endif
+            CI_eigenvectors_s2_mono(i_state) = s2
-      endif
+          endif
-      if (i_state.ge.N_states_diag) then
+        endif
-        exit
+        if (i_state.ge.N_states_diag) then
-      endif
+          exit
-    enddo
+        endif
+      enddo
+    else
+      do j=1,N_states_diag
+        call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2)
+        if(dabs(eigenvectors(1,j)).gt.0.9d0)then
+          i_state += 1
+          do i=1,N_det
+            CI_eigenvectors_mono(i,i_state) = eigenvectors(i,j)
+          enddo
+          CI_electronic_energy_mono(i_state) = eigenvalues(j)
+          CI_eigenvectors_s2_mono(i_state) = s2
+        endif
+      enddo
+    endif
     deallocate(eigenvectors,eigenvalues)
   endif
   

src/Determinants/filter_connected.irp.f

@@ -235,8 +235,8 @@ subroutine filter_connected_davidson(key1,key2,Nint,sze,idx)
     
   else if (Nint==3) then
     
-    !DIR$ LOOP COUNT (1000)
     i = idx(0)
+    !DIR$ LOOP COUNT (1000)
     do j_int=1,N_con_int
       itmp = det_connections(j_int,i)
       do while (itmp /= 0_8)
@@ -261,8 +261,8 @@ subroutine filter_connected_davidson(key1,key2,Nint,sze,idx)
     
   else
     
-    !DIR$ LOOP COUNT (1000)
     i = idx(0)
+    !DIR$ LOOP COUNT (1000)
     do j_int=1,N_con_int
       itmp = det_connections(j_int,i)
       do while (itmp /= 0_8)

src/Determinants/s2.irp.f

@@ -88,19 +88,19 @@ subroutine get_s2_u0(psi_keys_tmp,psi_coefs_tmp,n,nmax,s2)
  double precision, intent(out) :: s2
  integer :: i,j,l
  double precision :: s2_tmp
- s2 = S_z2_Sz
+ s2 = 0.d0
  !$OMP PARALLEL DO DEFAULT(NONE) &
  !$OMP PRIVATE(i,j,s2_tmp) SHARED(n,psi_coefs_tmp,psi_keys_tmp,N_int) &
  !$OMP REDUCTION(+:s2) SCHEDULE(dynamic)
- do i = 1, n
+ do i=1,n
    call get_s2(psi_keys_tmp(1,1,i),psi_keys_tmp(1,1,i),s2_tmp,N_int)
-!  print*,'s2_tmp = ',s2_tmp
+   s2 += psi_coefs_tmp(i)*psi_coefs_tmp(i)*s2_tmp
-  do j = 1, n
+   do j=i+1,n
-    call get_s2(psi_keys_tmp(1,1,i),psi_keys_tmp(1,1,j),s2_tmp,N_int)
+     call get_s2(psi_keys_tmp(1,1,i),psi_keys_tmp(1,1,j),s2_tmp,N_int)
-    if (s2_tmp == 0.d0) cycle
+     s2 += (psi_coefs_tmp(i)+psi_coefs_tmp(i))*psi_coefs_tmp(j)*s2_tmp
-    s2 += psi_coefs_tmp(i)*psi_coefs_tmp(j)*s2_tmp
+   enddo
-  enddo
  enddo
  !$OMP END PARALLEL DO
+ s2 +=  S_z2_Sz
 end
 

src/Determinants/save_for_casino.irp.f

@@ -161,7 +161,7 @@ subroutine save_casino
    if (ao_l(i) == ao_power(i,1)) then
      do j=1,ao_prim_num(i)
        l+=1
-       rtmp(l) = ao_coef(i,ao_prim_num(i)-j+1)
+       rtmp(l) = ao_coef_normalized(i,ao_prim_num(i))
      enddo
    endif
  enddo

src/Determinants/save_for_qmcchem.irp.f

@@ -7,8 +7,6 @@ subroutine save_dets_qmcchem
  integer, allocatable :: occ(:,:,:), occ_tmp(:,:)
  !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: occ, occ_tmp
 
- read_wf = .True.
- TOUCH read_wf
  call ezfio_set_determinants_det_num(N_det)
  call ezfio_set_determinants_det_coef(psi_coef_sorted(1,1))
 
@@ -46,6 +44,8 @@ subroutine save_dets_qmcchem
 end
 
 program save_for_qmc
-  call save_dets_qmcchem
+ read_wf = .True.
+ TOUCH read_wf
+!  call save_dets_qmcchem
  call write_spindeterminants
 end

src/Determinants/spindeterminants.irp.f

@@ -50,80 +50,88 @@ BEGIN_PROVIDER [ integer(bit_kind), psi_det_beta, (N_int,psi_det_size) ]
  enddo
 END_PROVIDER
 
- BEGIN_PROVIDER [ integer(bit_kind), psi_det_alpha_unique, (N_int,psi_det_size) ]
+
-&BEGIN_PROVIDER [ integer, N_det_alpha_unique ]
+BEGIN_TEMPLATE
+
+ BEGIN_PROVIDER [ integer(bit_kind), psi_det_$alpha_unique, (N_int,psi_det_size) ]
+&BEGIN_PROVIDER [ integer, N_det_$alpha_unique ]
  implicit none
  BEGIN_DOC
- ! Unique alpha determinants
+ ! Unique $alpha determinants
  END_DOC
 
- integer                        :: i,k
+ integer                        :: i,j,k
  integer, allocatable           :: iorder(:)
  integer*8, allocatable         :: bit_tmp(:)
  integer*8                      :: last_key
  integer*8, external            :: spin_det_search_key
+ logical,allocatable            :: duplicate(:)
 
- allocate ( iorder(N_det), bit_tmp(N_det))
+ allocate ( iorder(N_det), bit_tmp(N_det), duplicate(N_det) )
 
  do i=1,N_det
    iorder(i) = i
-   bit_tmp(i) = spin_det_search_key(psi_det_alpha(1,i),N_int)
+   bit_tmp(i) = spin_det_search_key(psi_det_$alpha(1,i),N_int)
  enddo
 
  call i8sort(bit_tmp,iorder,N_det)
 
- N_det_alpha_unique = 0
+ N_det_$alpha_unique = 0
  last_key = 0_8
  do i=1,N_det
-  if (bit_tmp(i) /= last_key) then
+   last_key = bit_tmp(i)
-    last_key = bit_tmp(i)
+   N_det_$alpha_unique += 1
-    N_det_alpha_unique += 1
+   do k=1,N_int
-    do k=1,N_int
+     psi_det_$alpha_unique(k,N_det_$alpha_unique) = psi_det_$alpha(k,iorder(i))
-      psi_det_alpha_unique(k,N_det_alpha_unique) = psi_det_alpha(k,iorder(i))
+   enddo
-    enddo
+   duplicate(i) = .False.
+ enddo
+
+ j=1
+ do i=1,N_det_$alpha_unique-1
+   if (duplicate(i)) then
+     cycle
+   endif
+   j = i+1
+   do while (bit_tmp(j)==bit_tmp(i))
+     if (duplicate(j)) then
+       j += 1
+       cycle
+     endif
+     duplicate(j) = .True.
+     do k=1,N_int
+       if (psi_det_$alpha_unique(k,i) /= psi_det_$alpha_unique(k,j)) then
+         duplicate(j) = .False.
+         exit
+       endif
+     enddo
+     j+=1
+     if (j > N_det_$alpha_unique) then
+       exit
+     endif
+   enddo
+ enddo
+
+ j=1
+ do i=2,N_det_$alpha_unique
+   if (duplicate(i)) then
+     cycle
+  else
+    j += 1
+    psi_det_$alpha_unique(:,j) = psi_det_$alpha_unique(:,i)
   endif
  enddo
+ N_det_$alpha_unique = j
 
- deallocate (iorder, bit_tmp)
+ deallocate (iorder, bit_tmp, duplicate)
 END_PROVIDER
 
- BEGIN_PROVIDER [ integer(bit_kind), psi_det_beta_unique, (N_int,psi_det_size) ]
+SUBST [ alpha ]
-&BEGIN_PROVIDER [ integer, N_det_beta_unique ]
- implicit none
- BEGIN_DOC
- ! Unique beta determinants
- END_DOC
 
- integer                        :: i,k
+alpha ;;
- integer, allocatable           :: iorder(:)
+beta ;;
- integer*8, allocatable         :: bit_tmp(:)
- integer*8                      :: last_key
- integer*8, external            :: spin_det_search_key
-
- allocate ( iorder(N_det), bit_tmp(N_det))
-
- do i=1,N_det
-   iorder(i) = i
-   bit_tmp(i) = spin_det_search_key(psi_det_beta(1,i),N_int)
- enddo
-
- call i8sort(bit_tmp,iorder,N_det)
-
- N_det_beta_unique = 0
- last_key = 0_8
- do i=1,N_det
-  if (bit_tmp(i) /= last_key) then
-    last_key = bit_tmp(i)
-    N_det_beta_unique += 1
-    do k=1,N_int
-      psi_det_beta_unique(k,N_det_beta_unique) = psi_det_beta(k,iorder(i))
-    enddo
-  endif
- enddo
-
- deallocate (iorder, bit_tmp)
-END_PROVIDER
 
+END_TEMPLATE
 
 
 
@@ -150,6 +158,7 @@ integer function get_index_in_psi_det_alpha_unique(key,Nint)
 
   !DIR$ FORCEINLINE
   det_ref = spin_det_search_key(key,Nint)
+
   !DIR$ FORCEINLINE
   det_search = spin_det_search_key(psi_det_alpha_unique(1,1),Nint)
 
@@ -337,52 +346,52 @@ subroutine write_spindeterminants
   call ezfio_set_spindeterminants_psi_coef_matrix_rows(psi_svd_matrix_rows)
   call ezfio_set_spindeterminants_psi_coef_matrix_columns(psi_svd_matrix_columns)
   
-  integer                        :: n_svd_coefs
+!  integer                        :: n_svd_coefs
-  double precision               :: norm, f
+!  double precision               :: norm, f
-  f = 1.d0/dble(N_states)
+!  f = 1.d0/dble(N_states)
-  norm = 1.d0
+!  norm = 1.d0
-  do n_svd_coefs=1,N_det_alpha_unique
+!  do n_svd_coefs=1,N_det_alpha_unique
-    do k=1,N_states
+!    do k=1,N_states
-      norm -= psi_svd_coefs(n_svd_coefs,k)*psi_svd_coefs(n_svd_coefs,k)
+!      norm -= psi_svd_coefs(n_svd_coefs,k)*psi_svd_coefs(n_svd_coefs,k)
-    enddo
+!    enddo
-    if (norm < 1.d-4) then
+!    if (norm < 1.d-4) then
-      exit
+!      exit
-    endif
+!    endif
-  enddo
+!  enddo
-  n_svd_coefs -= 1
+!  n_svd_coefs -= 1
-  call ezfio_set_spindeterminants_n_svd_coefs(n_svd_coefs)
+!  call ezfio_set_spindeterminants_n_svd_coefs(n_svd_coefs)
-  
+!  
-  double precision, allocatable  :: dtmp(:,:,:)
+!  double precision, allocatable  :: dtmp(:,:,:)
-  allocate(dtmp(N_det_alpha_unique,n_svd_coefs,N_states))
+!  allocate(dtmp(N_det_alpha_unique,n_svd_coefs,N_states))
-  do k=1,N_states
+!  do k=1,N_states
-    do j=1,n_svd_coefs
+!    do j=1,n_svd_coefs
-      do i=1,N_det_alpha_unique
+!      do i=1,N_det_alpha_unique
-        dtmp(i,j,k) = psi_svd_alpha(i,j,k)
+!        dtmp(i,j,k) = psi_svd_alpha(i,j,k)
-      enddo
+!      enddo
-    enddo
+!    enddo
-  enddo
+!  enddo
-  call ezfio_set_spindeterminants_psi_svd_alpha(dtmp)
+!  call ezfio_set_spindeterminants_psi_svd_alpha(dtmp)
-  deallocate(dtmp)
+!  deallocate(dtmp)
-
+!
-  allocate(dtmp(N_det_beta_unique,n_svd_coefs,N_states))
+!  allocate(dtmp(N_det_beta_unique,n_svd_coefs,N_states))
-  do k=1,N_states
+!  do k=1,N_states
-    do j=1,n_svd_coefs
+!    do j=1,n_svd_coefs
-      do i=1,N_det_beta_unique
+!      do i=1,N_det_beta_unique
-        dtmp(i,j,k) = psi_svd_beta(i,j,k)
+!        dtmp(i,j,k) = psi_svd_beta(i,j,k)
-      enddo
+!      enddo
-    enddo
+!    enddo
-  enddo
+!  enddo
-  call ezfio_set_spindeterminants_psi_svd_beta(dtmp)
+!  call ezfio_set_spindeterminants_psi_svd_beta(dtmp)
-  deallocate(dtmp)
+!  deallocate(dtmp)
-
+!
-  allocate(dtmp(n_svd_coefs,N_states,1))
+!  allocate(dtmp(n_svd_coefs,N_states,1))
-  do k=1,N_states
+!  do k=1,N_states
-    do j=1,n_svd_coefs
+!    do j=1,n_svd_coefs
-        dtmp(j,k,1) = psi_svd_coefs(j,k)
+!        dtmp(j,k,1) = psi_svd_coefs(j,k)
-    enddo
+!    enddo
-  enddo
+!  enddo
-  call ezfio_set_spindeterminants_psi_svd_coefs(dtmp)
+!  call ezfio_set_spindeterminants_psi_svd_coefs(dtmp)
-  deallocate(dtmp)
+!  deallocate(dtmp)
 
 end
 
@@ -439,6 +448,7 @@ BEGIN_PROVIDER  [ double precision, psi_svd_matrix_values, (N_det,N_states) ]
   do k=1,N_det
     i = get_index_in_psi_det_alpha_unique(psi_det(1,1,k),N_int)
     j = get_index_in_psi_det_beta_unique (psi_det(1,2,k),N_int)
+
     do l=1,N_states
       psi_svd_matrix_values(k,l) = psi_coef(k,l)
     enddo