Merge branch 'master' of github.com:scemama/quantum_package

Conflicts:
	scripts/ezfio_interface.py
	scripts/install_ocaml.sh

INSTALL.md

@@ -3,23 +3,40 @@
 ## Requirements
 
 * curl
-* wget
 * m4
 * GNU make
-* Intel Fortran compiler
+* Fortran compiler (ifort or gfortran are tested)
-* Python
+* Python >= 2.7
 * Bash
 
 
-By default, the Ocaml compiler will be installed in `$HOME/ocamlbrew`.
+## Standard installation
-To install it somewhere else, set the `$OCAMLBREW_BASE` environment
-variable to the required destination, for example:
 
-    export OCAMLBREW_BASE=/usr/local/ocamlbrew
+1) `./setup_environment.sh`
 
-For more info about the Ocaml installation, check the ocamlbrew
+This command will download and install all the requirements.
-website : https://github.com/hcarty/ocamlbrew
+Installing OCaml and the Core library may take somme time
+(up to 20min on an old machine).
 
+2) `source quantum_package.rc`
+
+This file contains all the environement variables neeeded by the quantum package
+both to compile and run. This should also be done before running calculations.
+
+3) `cp ./src/Makefile.config.gfortran  ./src/Makefile.config`
+
+Create the ``Makefile.config`` which contains all the flags needed by the compiler.
+The is also an example for the Intel Compiler (`Makefile.config.ifort`).
+Edit this file and tune the flags as you want.
+
+4) `make build`
+
+It will compile all the executables and tools. 
+
+5) `make binary`
+
+Optional. It will build a `tar.gz` file containing everything needed to run the quantum package on a
+machine where you can't compile.
 
 
 ## Installing behind a firewall
@@ -61,4 +78,3 @@ website : https://github.com/hcarty/ocamlbrew
 
     ssh -fN -D 10000 user@external-server.com
 
-

README.md

@@ -1,8 +1,10 @@
 Quantum package
 ===============
+[![Build Status](https://travis-ci.org/LCPQ/quantum_package.svg?branch=master)](https://travis-ci.org/LCPQ/quantum_package)
 
 [![Gitter](https://badges.gitter.im/Join Chat.svg)](https://gitter.im/LCPQ/quantum_package?utm_source=badge&utm_medium=badge&utm_campaign=pr-badge&utm_content=badge)
 
+
 Set of quantum chemistry programs and libraries.
 
 For more information, you can visit the [wiki of the project](http://github.com/LCPQ/quantum_package/wiki>), or the [Install](INSTALL.md) file.

scripts/ezfio_interface.py

@@ -62,13 +62,13 @@ def get_type_dict():
     # ~#~#~#~#~#~#~#~ #
 
     fancy_type['integer'] = Type("int", "integer")
-    fancy_type['int'    ] = Type("int", "integer")
+    fancy_type['int'] = Type("int", "integer")
 
-    fancy_type['float'           ] = Type("float", "double precision")
+    fancy_type['float'] = Type("float", "double precision")
     fancy_type['double precision'] = Type("float", "double precision")
 
     fancy_type['logical'] = Type("bool", "logical")
-    fancy_type['bool'   ] = Type("bool", "logical")
+    fancy_type['bool'] = Type("bool", "logical")
 
     # ~#~#~#~#~#~#~#~ #
     # q p _ t y p e s #
@@ -92,18 +92,39 @@ def get_type_dict():
 type_dict = get_type_dict()
 
 
-def get_dict_config_file(config_file_path,folder):
+def get_dict_config_file(config_file_path, module_lower):
     """
-    Read a ezfio.cfg
+    Input:
-    Return a dict d[provider_name] = {type, default, ezfio_name,ezfio_dir,doc}
+        config_file_path is the config file path
+            (for example FULL_PATH/EZFIO.cfg)
+        module_lower is the MODULE name lowered
+            (Ex fullci)
+
+    Return a dict d[provider_name] = {type,
+                                      doc,
+                                      ezfio_name,
+                                      ezfio_dir,
+                                      interface,
+                                      default}
+
+    Type       : Is a fancy_type named typle who containt fortran and ocaml type
+    doc        : Is the doc
+    ezfio_name : Will be the name of the file
+    ezfio_dir  : Will be the folder who containt the ezfio_name
+        * /ezfio_dir/ezfio_name
+        * equal to MODULE_lower name for the moment.
+    interface  : The provider is a imput or a output
+    if is a output:
+        default    : The default value
+
     """
     # ~#~#~#~ #
     # I n i t #
     # ~#~#~#~ #
 
     d = defaultdict(dict)
-    list_option_required = ["default", "doc"]
+    l_info_required = ["doc", "interface"]
-    list_option          = ["ezfio_name", "output"]
+    l_info_optional = ["ezfio_name"]
 
     # ~#~#~#~#~#~#~#~#~#~#~ #
     # L o a d _ C o n f i g #
@@ -116,43 +137,63 @@ def get_dict_config_file(config_file_path,folder):
     # F i l l _ d i c t #
     # ~#~#~#~#~#~#~#~#~ #
 
-    provider_names = config_file.sections()
+    def error(o, p, c):
-    for p in provider_names:
+        "o option ; p provider_name ;c config_file_path"
-        provider_name = p.lower()
+        print "You need a {0} for {1} in {2}".format(o, p, c)
-        default_d = {"ezfio_name": provider_name, "output": "false" }
+
+    for section in config_file.sections():
+        # pvd = provider
+        pvd = section.lower()
+
+        # Create the dictionary who containt the value per default
+        d_default = {"ezfio_name": pvd}
+
+        # Set the ezfio_dir
+        d[pvd]["ezfio_dir"] = module_lower
+
         # Check if type if avalaible
-        type_ = config_file.get(p, "type")
+        type_ = config_file.get(section, "type")
         if type_ not in type_dict:
             print "{0} not avalaible. Choose in:".format(type_)
             print ", ".join([i for i in type_dict])
             sys.exit(1)
         else:
-            d[provider_name]["type"] = type_dict[type_]
+            d[pvd]["type"] = type_dict[type_]
 
-        # Fill the dict with allother the information
+        # Fill the dict with REQUIRED information
-        for k in list_option_required:
+        for option in l_info_required:
             try:
-                d[provider_name][k] = config_file.get(p, k)
+                d[pvd][option] = config_file.get(section, option)
             except ConfigParser.NoOptionError:
-                print "You need a {0} for {1} in {2}".format(k,
+                error(option, pvd, config_file_path)
-                                                             provider_name,
+                sys.exit(1)
-                                                             config_file_path)
-        d[provider_name]["ezfio_dir"] = folder
-        for k in list_option:
-            try:
-                d[provider_name][k] = config_file.get(p, k).lower()
-            except ConfigParser.NoOptionError:
-                d[provider_name][k] = default_d[k]
 
-        # Convert string to bool
+        # Fill the dict with OPTIONAL information
-        d[provider_name]["output"] = bool_convertor(d[provider_name]["output"])
+        for option in l_info_optional:
+            try:
+                d[pvd][option] = config_file.get(section, option).lower()
+            except ConfigParser.NoOptionError:
+                d[pvd][option] = d_default[option]
+
+        # If interface is output we need a default value information
+        if d[pvd]["interface"] == "output":
+            try:
+                d[pvd]["default"] = config_file.get(section, "default")
+            except ConfigParser.NoOptionError:
+                error("default", pvd, config_file_path)
+                sys.exit(1)
 
     return dict(d)
 
 
 def create_ezfio_provider(dict_ezfio_cfg):
     """
-    From dict d[provider_name] = {type, default, ezfio_name,ezfio_dir,doc}
+    From dict d[provider_name] = {type,
+                                  doc,
+                                  ezfio_name,
+                                  ezfio_dir,
+                                  interface,
+                                  default}
     create the a list who containt all the code for the provider
     return [code, ...]
     """
@@ -162,14 +203,14 @@ def create_ezfio_provider(dict_ezfio_cfg):
 
     ez_p = EZFIO_Provider()
     for provider_name, dict_info in dict_ezfio_cfg.iteritems():
-        if not dict_info["output"]:
+        if "default" in dict_info:
             ez_p.set_type(dict_info['type'].fortran)
             ez_p.set_name(provider_name)
             ez_p.set_doc(dict_info['doc'])
             ez_p.set_ezfio_dir(dict_info['ezfio_dir'])
             ez_p.set_ezfio_name(dict_info['ezfio_name'])
             ez_p.set_default(dict_info['default'])
-     
+
             ez_p.set_output("output_%s" % dict_info['ezfio_dir'])
             dict_code_provider[provider_name] = str(ez_p)
 
@@ -178,8 +219,7 @@ def create_ezfio_provider(dict_ezfio_cfg):
 
 def save_ezfio_provider(path_head, dict_code_provider):
     """
-    Write in "ezfio_interface.irp.f" the
+    Write in path_head/"ezfio_interface.irp.f" the value of dict_code_provider
-    value of dict_code_provider
     """
 
     path = "{0}/ezfio_interface.irp.f".format(path_head)
@@ -206,30 +246,33 @@ def save_ezfio_provider(path_head, dict_code_provider):
             f.write(output)
 
 
-def create_ezfio_config(dict_ezfio_cfg, opt, folder):
+def create_ezfio_config(dict_ezfio_cfg, opt, module_lower):
     """
     From dict_ezfio_cfg[provider_name] = {type, default, ezfio_name,ezfio_dir,doc}
     Return the string ezfio_interface_config
     """
 
-    result = [ folder ]
+    result = [module_lower]
-    lenmax = max( [ len(i) for i in dict_ezfio_cfg ]  )+2
+    lenmax = max([len(i) for i in dict_ezfio_cfg]) + 2
     l = sorted(dict_ezfio_cfg.keys())
     for provider_name in l:
-      provider_info = dict_ezfio_cfg[provider_name]
+        provider_info = dict_ezfio_cfg[provider_name]
-      s = "  {0} {1}".format( provider_name.lower().ljust(lenmax), provider_info["type"].fortran )
+        s = "  {0} {1}".format(
-      result.append(s)
+            provider_name.lower().ljust(lenmax),
+            provider_info["type"].fortran)
+        result.append(s)
     return "\n".join(result)
 
-def save_ezfio_config(folder, str_ezfio_config):
+
+def save_ezfio_config(module_lower, str_ezfio_config):
     """
     Write the str_ezfio_config in
-    $QPACKAGE_ROOT/EZFIO/{0}.ezfio_interface_config".format(folder)
+    $QPACKAGE_ROOT/EZFIO/{0}.ezfio_interface_config".format(module_lower)
     """
 
     ezfio_dir = "{0}/EZFIO".format(os.environ['QPACKAGE_ROOT'])
     path = "{0}/config/{1}.ezfio_interface_config".format(ezfio_dir,
-                                                          folder)
+                                                          module_lower)
 
 #    print "Path = {}".format(path)
 
@@ -245,27 +288,33 @@ def save_ezfio_config(folder, str_ezfio_config):
         with open(path, "w") as f:
             f.write(str_ezfio_config)
 
+
 def main():
-    """Take in argument a EZFIO.cfg"""
+    """
+    Two condition:
+        -Take the EZFIO.cfg path in arg
+        or
+        -Look if EZFIO.cfg is present in the pwd
+    """
 
     try:
-      path = sys.argv[1]
+        config_file_path = sys.argv[1]
     except:
-      path = "EZFIO.cfg"
+        config_file_path = "EZFIO.cfg"
-      if "EZFIO.cfg" not in os.listdir(os.getcwd()):
+        if "EZFIO.cfg" not in os.listdir(os.getcwd()):
-        sys.exit(0)
+            sys.exit(0)
 
-    path = os.path.expanduser(path)
+    config_file_path = os.path.expanduser(config_file_path)
-    path = os.path.expandvars(path)
+    config_file_path = os.path.expandvars(config_file_path)
-    path = os.path.abspath(path)
+    config_file_path = os.path.abspath(config_file_path)
-#    print path
+#    print config_file_path
 
-    path_dirname = os.path.dirname(path)
+    path_dirname = os.path.dirname(config_file_path)
-    folder = [i for i in path_dirname.split("/") if i][-1]
+    module = [i for i in path_dirname.split("/") if i][-1]
-    folder = folder.lower()
+    module_lower = module.lower()
 
-#    print "Found a EZFIO.cfg in {}".format(path)
+#    print "Read {0}".format(config_file_path)
-    dict_info_provider = get_dict_config_file(path,folder)
+    dict_info_provider = get_dict_config_file(config_file_path, module_lower)
 
 #    print "Generating the ezfio_interface.irp.f: \n"
     d_config = create_ezfio_provider(dict_info_provider)
@@ -274,12 +323,11 @@ def main():
     save_ezfio_provider(path_dirname, d_config)
 
 #    print "Generating the ezfio_config"
-    config_ezfio = create_ezfio_config(dict_info_provider, "config", folder)
+    config_ezfio = create_ezfio_config(dict_info_provider, "config", module_lower)
 
 #    print "Saving ezfio_config"
-    save_ezfio_config(folder, config_ezfio)
+    save_ezfio_config(module_lower, config_ezfio)
 
 
 if __name__ == "__main__":
-  main()
+    main()
-

scripts/generate_h_apply.py

@@ -24,6 +24,10 @@ skip
 init_main
 filter_integrals
 filter2h2p
+filterhole
+filterparticle
+do_double_excitations
+check_double_excitation
 """.split()
 
 class H_apply(object):
@@ -116,6 +120,24 @@ class H_apply(object):
       buffer = buffer.replace('$'+key, value)
     return buffer
 
+  def unset_double_excitations(self):
+    self["do_double_excitations"] = ".False."
+    self["check_double_excitation"] = """
+     check_double_excitation = .False.
+    """
+  def set_filter_holes(self):
+    self["filterhole"] = """
+     if(iand(ibset(0_bit_kind,j),hole(k,other_spin)).eq.0_bit_kind )cycle
+    """
+  def set_filter_particl(self):
+    self["filterparticle"] = """
+     if(iand(ibset(0_bit_kind,j_a),hole(k_a,other_spin)).eq.0_bit_kind )cycle
+    """
+  def unset_skip(self):
+    self["skip"] = """
+    """
+
+
   def set_filter_2h_2p(self):
     self["filter2h2p"] = """
 !    ! DIR$ FORCEINLINE

scripts/install_ocaml.sh

@@ -50,11 +50,12 @@ ${QPACKAGE_ROOT}/scripts/fetch_from_web.py \
 chmod +x opam_installer.sh
 echo N | ./opam_installer.sh ${QPACKAGE_ROOT}/bin 
 
-if [[ ! -f ${QPACKAGE_ROOT}/bin/opam ]]; then
+if [[ ! -f ${QPACKAGE_ROOT}/bin/opam ]]
+then
    echo "Installation of OPAM failed"
    exit 2
 fi
-rm opam_installer.sh
+rm -f opam_installer.sh
 ${QPACKAGE_ROOT}/bin/opam config setup -a --dot-profile ${QPACKAGE_ROOT}/quantum_package.rc 
 
 source ${QPACKAGE_ROOT}/quantum_package.rc

scripts/install_zlib.sh

@@ -23,7 +23,7 @@ EOF
 gcc /tmp/main.c -lz
 if [[ $? -eq 0 ]]
 then
-  rm /tmp/main.c
+  rm /tmp/main.c "a.out"
   exit 0
 fi
 rm /tmp/main.c

src/Dets/H_apply_template.f

@@ -36,6 +36,11 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene
     ifirst=1
   endif
   
+  logical :: check_double_excitation 
+  check_double_excitation = .True.
+
+
+
   $initialization
   
   $omp_parallel
@@ -276,6 +281,12 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param
   logical, allocatable           :: array_pairs(:,:)
   double precision               :: diag_H_mat_elem
   integer(omp_lock_kind), save   :: lck, ifirst=0
+
+  logical :: check_double_excitation 
+  check_double_excitation = .True.
+  $check_double_excitation
+
+
   if (ifirst == 0) then
     ifirst=1
 !$    call omp_init_lock(lck)
@@ -333,9 +344,11 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param
       hole = key_in
       k = ishft(i_a-1,-bit_kind_shift)+1
       j = i_a-ishft(k-1,bit_kind_shift)-1
+  $filterhole
       hole(k,ispin) = ibclr(hole(k,ispin),j)
       k_a = ishft(j_a-1,-bit_kind_shift)+1
       l_a = j_a-ishft(k_a-1,bit_kind_shift)-1
+  $filterparticle
       hole(k_a,ispin) = ibset(hole(k_a,ispin),l_a)
       $filter2h2p
       key_idx += 1

src/Full_CI/EZFIO.cfg

@@ -1,37 +1,42 @@
 [N_det_max_fci]
-doc: Max number of determinants in the wave function
 type: Det_number_max
+doc: Max number of determinants in the wave function
+interface: output
 default: 10000
 
 [N_det_max_fci_property]
-doc: Max number of determinants in the wave function when you select for a given property
 type: Det_number_max
+doc: Max number of determinants in the wave function when you select for a given property
+interface: output
 default: 10000
 
 [do_pt2_end]
 type: logical
 doc: If true, compute the PT2 at the end of the selection
+interface: output
 default: true
 
 [PT2_max]
 type: PT2_energy
 doc: The selection process stops when the largest PT2 (for all the state is lower
      than pt2_max in absolute value
+interface: output
 default: 0.0001
 
 [var_pt2_ratio]
 type: Normalized_float
 doc: The selection process stops when the energy ratio variational/(variational+PT2)
      is equal to var_pt2_ratio
+interface: output
 default: 0.75
 
 [energy]
 type: double precision
 doc: "Calculated Full CI energy"
-output: true
+interface: input
 
 [energy_pt2]
 type: double precision
 doc: "Calculated Full CI energy"
-output: true
+interface: input
 

src/Full_CI/H_apply.irp.f

@@ -10,23 +10,26 @@ s = H_apply("FCI_PT2")
 s.set_perturbation("epstein_nesbet_2x2")
 print s
 
+s = H_apply("FCI_no_skip")
+s.set_selection_pt2("epstein_nesbet_2x2")
+s.unset_skip()
+print s
 
-if False:
+s = H_apply("FCI_mono")
-  s = H_apply("FCI_mono")
+s.set_selection_pt2("epstein_nesbet_2x2")
-  s.set_selection_pt2("epstein_nesbet_2x2")
+s.unset_double_excitations()
-  s.unset_double_excitations()
+print s
-  print s
 
 
-  s = H_apply("select_mono_delta_rho")
+s = H_apply("select_mono_delta_rho")
-  s.unset_double_excitations()
+s.unset_double_excitations()
-  s.set_selection_pt2("delta_rho_one_point")
+s.set_selection_pt2("delta_rho_one_point")
-  print s
+print s
 
-  s = H_apply("pt2_mono_delta_rho")
+s = H_apply("pt2_mono_delta_rho")
-  s.unset_double_excitations()
+s.unset_double_excitations()
-  s.set_perturbation("delta_rho_one_point")
+s.set_perturbation("delta_rho_one_point")
-  print s
+print s
 
 s = H_apply("select_mono_di_delta_rho")
 s.set_selection_pt2("delta_rho_one_point")

src/Full_CI/full_ci_no_skip.irp.f

@@ -0,0 +1,91 @@
+program full_ci
+  implicit none
+  integer                        :: i,k
+
+  
+  double precision, allocatable  :: pt2(:), norm_pert(:), H_pert_diag(:)
+  integer                        :: N_st, degree
+  N_st = N_states
+  allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
+  character*(64)                 :: perturbation
+  
+  pt2 = 1.d0
+  diag_algorithm = "Lapack"
+  if (N_det > n_det_max_fci) then
+    call diagonalize_CI
+    call save_wavefunction
+    psi_det = psi_det_sorted
+    psi_coef = psi_coef_sorted
+    N_det = n_det_max_fci
+    soft_touch N_det psi_det psi_coef
+    call diagonalize_CI
+    call save_wavefunction
+    print *,  'N_det    = ', N_det
+    print *,  'N_states = ', N_states
+    print *,  'PT2      = ', pt2
+    print *,  'E        = ', CI_energy
+    print *,  'E+PT2    = ', CI_energy+pt2
+    print *,  '-----'
+  endif
+  double precision :: i_H_psi_array(N_states),diag_H_mat_elem,h,i_O1_psi_array(N_states)
+  if(read_wf)then
+   call i_H_psi(psi_det(1,1,N_det),psi_det,psi_coef,N_int,N_det,psi_det_size,N_states,i_H_psi_array)
+   h = diag_H_mat_elem(psi_det(1,1,N_det),N_int)
+   selection_criterion = dabs(psi_coef(N_det,1) *  (i_H_psi_array(1) - h * psi_coef(N_det,1))) * 0.1d0
+   soft_touch selection_criterion
+  endif
+
+
+  integer :: n_det_before
+  print*,'Beginning the selection ...'
+  do while (N_det < n_det_max_fci.and.maxval(abs(pt2(1:N_st))) > pt2_max)
+    n_det_before = N_det
+    call H_apply_FCI_no_skip(pt2, norm_pert, H_pert_diag,  N_st)
+
+    PROVIDE  psi_coef
+    PROVIDE  psi_det
+    PROVIDE  psi_det_sorted
+
+    if (N_det > n_det_max_fci) then
+       psi_det = psi_det_sorted
+       psi_coef = psi_coef_sorted
+       N_det = n_det_max_fci
+       soft_touch N_det psi_det psi_coef
+    endif
+    call diagonalize_CI
+    call save_wavefunction
+    if(n_det_before == N_det)then
+     selection_criterion = selection_criterion * 0.5d0
+    endif
+    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_full_ci_energy(CI_energy)
+    if (abort_all) then
+      exit
+    endif
+  enddo
+   N_det = min(n_det_max_fci,N_det)
+   touch N_det psi_det psi_coef
+   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_FCI_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_full_ci_energy_pt2(CI_energy+pt2)
+   endif
+   call save_wavefunction
+  deallocate(pt2,norm_pert)
+end

src/Utils/map_module.f90

@@ -11,6 +11,10 @@ module map_module
 ! as integer*2 and is found by applying the map_mask
 ! to the initial key. The element are found in the
 ! cache_map using a binary search
+!
+! When using the map_update subroutine to build the map,
+! the map_unique subroutine
+! should be called before getting data from the map.
 
  use omp_lib
  
@@ -433,6 +437,97 @@ call omp_unset_lock(map%lock)
 
 end
 
+subroutine map_update_verbose(map, key, value, sze, thr)
+  use map_module
+  implicit none
+  type (map_type), intent(inout) :: map
+  integer, intent(in)            :: sze
+  integer(key_kind), intent(inout) :: key(sze)
+  real(integral_kind), intent(inout) :: value(sze)
+  real(integral_kind), intent(in) :: thr
+  
+  integer                        :: i
+  integer(map_size_kind)         :: idx_cache, idx_cache_new
+  integer(cache_map_size_kind)   :: idx
+  integer                        :: sze2
+  integer(cache_key_kind)        :: cache_key
+  integer(map_size_kind)         :: n_elements_temp
+  type (cache_map_type)          :: local_map
+  logical                        :: map_sorted
+! do i = 1, sze
+!  print*,'value in map = ',value(i)
+! enddo
+  
+  sze2 = sze
+  map_sorted = .True.
+  
+  n_elements_temp = 0_8
+  n_elements_temp = n_elements_temp + 1_8
+  do while (sze2>0)
+    i=1
+    do while (i<=sze)
+      if (key(i) /= 0_8) then
+        idx_cache = ishft(key(i),map_shift)
+        if (omp_test_lock(map%map(idx_cache)%lock)) then
+          local_map%key => map%map(idx_cache)%key
+          local_map%value => map%map(idx_cache)%value
+          local_map%sorted = map%map(idx_cache)%sorted
+          local_map%map_size = map%map(idx_cache)%map_size
+          local_map%n_elements = map%map(idx_cache)%n_elements
+          do
+          !DIR$ FORCEINLINE
+          call search_key_big_interval(key(i),local_map%key, local_map%n_elements, idx, 1, local_map%n_elements)
+          if (idx > 0_8) then
+!           print*,'AHAAH'
+!           print*,'local_map%value(idx) = ',local_map%value(idx)
+            local_map%value(idx) = local_map%value(idx) + value(i)
+!           print*,'not a new value !'
+!           print*,'local_map%value(idx) = ',local_map%value(idx)
+          else
+            ! Assert that the map has a proper size
+            if (local_map%n_elements == local_map%map_size) then
+              call cache_map_unique(local_map)
+              call cache_map_reallocate(local_map, local_map%n_elements + local_map%n_elements)
+              call cache_map_shrink(local_map,thr)
+            endif
+            cache_key = iand(key(i),map_mask)
+            local_map%n_elements = local_map%n_elements + 1_8
+            local_map%value(local_map%n_elements) = value(i)
+!           print*,'new value !'
+            local_map%key(local_map%n_elements) = cache_key
+            local_map%sorted = .False.
+            n_elements_temp = n_elements_temp + 1_8
+          endif  ! idx > 0
+          key(i) = 0_8
+          i = i+1
+          sze2 = sze2-1
+          if (i>sze) then
+            i=1
+          endif
+          if ( (ishft(key(i),map_shift) /= idx_cache).or.(key(i)==0_8)) then
+            exit
+          endif
+        enddo
+        map%map(idx_cache)%key => local_map%key
+        map%map(idx_cache)%value => local_map%value
+        map%map(idx_cache)%sorted = local_map%sorted
+        map%map(idx_cache)%n_elements = local_map%n_elements
+        map%map(idx_cache)%map_size = local_map%map_size
+        map_sorted = map_sorted .and. local_map%sorted
+        call omp_unset_lock(map%map(idx_cache)%lock)
+      endif  ! omp_test_lock
+    else
+      i=i+1
+    endif  ! key = 0
+  enddo  ! i
+enddo  ! sze2 > 0
+call omp_set_lock(map%lock)
+map%n_elements = map%n_elements + n_elements_temp
+map%sorted = map%sorted .and. map_sorted
+call omp_unset_lock(map%lock)
+
+end
+
 subroutine map_append(map, key, value, sze)
   use map_module
   implicit none