diff --git a/bin/qp_gaussian b/bin/qp_gaussian
new file mode 100755
index 00000000..a059a023
--- /dev/null
+++ b/bin/qp_gaussian
@@ -0,0 +1,83 @@
+#!/usr/bin/env python3
+#
+"""
+Runs a Quantum Package calculation using a Gaussian input file.
+
+Usage:
+      qp_gaussian INPUT
+
+"""
+
+# Requires pymatgen (https://pymatgen.org/)
+# pip install pymatgen
+
+
+import os
+import sys
+import os.path
+
+try:
+    import qp_path
+except ImportError:
+    print("source quantum_package.rc")
+
+from docopt import docopt
+import pymatgen
+from pymatgen.io.gaussian import GaussianInput
+
+
+def main(arguments):
+
+    filename = arguments["INPUT"]
+
+    with open(filename,'r') as f:
+        text = f.read()
+
+    in_file = GaussianInput.from_string(text)
+
+    d = in_file.as_dict()
+    charge = ("%d"%(d["charge"])).replace('-','m')
+    basis  = d["basis_set"]
+    mult   = d["spin_multiplicity"]
+    natoms = len(d["molecule"]["sites"])
+    with open("g09.xyz","w") as f:
+      f.write("%d\n"%natoms)
+      f.write("%s\n"%d["title"])
+      f.write("%s\n"%in_file.get_cart_coords())
+
+    if basis is None:
+       print("Basis set not found. Use '/' before basis set")
+       sys.exit(1)
+
+    command = f"rm -rf g09.ezfio"
+    os.system(command)
+
+    command = f"qp_create_ezfio -c {charge} -m {mult} g09.xyz -b {basis} -o g09.ezfio"
+    os.system(command)
+
+    command = f"rm -rf g09.xyz"
+    os.system(command)
+
+    command = f"qp_run scf g09.ezfio"
+    os.system(command)
+
+    command = f"qp_set_frozen_core g09.ezfio"
+    os.system(command)
+
+    if d["functional"] == "FCI":
+      command = f"qp_run fci g09.ezfio"
+    elif d["functional"] == "CIS":
+      command = f"qp_run cis g09.ezfio"
+    elif d["functional"] == "CISD":
+      command = f"qp_run cisd g09.ezfio"
+
+    os.system(command)
+
+    
+
+
+
+
+if __name__ == '__main__':
+    ARGUMENTS = docopt(__doc__)
+    main(ARGUMENTS)
diff --git a/config/gfortran.cfg b/config/gfortran.cfg
index 342acae9..ec72e722 100644
--- a/config/gfortran.cfg
+++ b/config/gfortran.cfg
@@ -13,7 +13,7 @@
 FC           : gfortran -g -ffree-line-length-none -I . -fPIC
 LAPACK_LIB   : -lblas -llapack
 IRPF90       : irpf90
-IRPF90_FLAGS : --ninja --align=32 --assert 
+IRPF90_FLAGS : --ninja --align=32 --assert -DGNU_CHECK_OMP
 
 # Global options
 ################
diff --git a/config/ifort.cfg b/config/ifort.cfg
index 714c4b10..0382360a 100644
--- a/config/ifort.cfg
+++ b/config/ifort.cfg
@@ -9,7 +9,7 @@
 FC           : ifort -fpic
 LAPACK_LIB   : -mkl=parallel -lirc -lsvml -limf -lipps
 IRPF90       : irpf90
-IRPF90_FLAGS : --ninja --align=32  -DINTEL
+IRPF90_FLAGS : --ninja --align=32 -DINTEL
 
 # Global options
 ################
diff --git a/config/ifort_debug.cfg b/config/ifort_debug.cfg
index 9b718380..d70b1465 100644
--- a/config/ifort_debug.cfg
+++ b/config/ifort_debug.cfg
@@ -9,7 +9,7 @@
 FC           : ifort -fpic
 LAPACK_LIB   : -mkl=parallel -lirc -lsvml -limf -lipps
 IRPF90       : irpf90
-IRPF90_FLAGS : --ninja --align=32  --assert -DINTEL
+IRPF90_FLAGS : --ninja --align=32 --assert -DINTEL
 
 # Global options
 ################
diff --git a/config/ifort_rome.cfg b/config/ifort_rome.cfg
index 5ed01227..1ac78717 100644
--- a/config/ifort_rome.cfg
+++ b/config/ifort_rome.cfg
@@ -9,7 +9,7 @@
 FC           : ifort -fpic
 LAPACK_LIB   : -mkl=parallel -lirc -lsvml -limf -lipps
 IRPF90       : irpf90
-IRPF90_FLAGS : --ninja --align=32  -DINTEL
+IRPF90_FLAGS : --ninja --align=32 -DINTEL
 
 # Global options
 ################
diff --git a/configure b/configure
index 052e5257..65146c22 100755
--- a/configure
+++ b/configure
@@ -3,8 +3,6 @@
 # Quantum Package configuration script
 #
 
-TEMP=$(getopt -o d:c:i:h -l download:,config:,install:,help -n $0 -- "$@") || exit 1
-eval set -- "$TEMP"
 
 export QP_ROOT="$( cd "$(dirname "$0")" ; pwd -P )"
 echo "QP_ROOT="$QP_ROOT
@@ -24,17 +22,17 @@ function help()
 Quantum Package configuration script.
 
 Usage:
-  $(basename $0) -c <file>    | --config=<file>
-  $(basename $0) -h           | --help
-  $(basename $0) -i <package> | --install=<package>
+  $(basename $0) -c <file>    
+  $(basename $0) -h           
+  $(basename $0) -i <package> 
 
 Options:
-  -c, --config=<file>      Define a COMPILATION configuration file,
-                           in "${QP_ROOT}/config/".
-  -h, --help               Print the HELP message
-  -i, --install=<package>  INSTALL <package>. Use at your OWN RISK:
-                           no support will be provided for the installation of
-                           dependencies.
+  -c  <file>    Define a COMPILATION configuration file,
+                in "${QP_ROOT}/config/".
+  -h            Print the HELP message
+  -i <package>  INSTALL <package>. Use at your OWN RISK:
+                no support will be provided for the installation of
+                dependencies.
 
 Example:
   ./$(basename $0) -c config/gfortran.cfg
@@ -68,32 +66,31 @@ function execute () {
 }
 
 PACKAGES=""
+echo $@
 
-while true ; do
-    case "$1" in
-        -c|--config)
-            case "$2" in
+
+while getopts "d:c:i:h" c ; do
+    case "$c" in
+        c)
+            case "$OPTARG" in
                 "") help ; break;;
-                *) if [[ -f $2 ]] ; then
-                      CONFIG="$2"
+                *) if [[ -f $OPTARG ]] ; then
+                      CONFIG="$OPTARG"
                    else
-                      error "error: configuration file $2 not found."
+                      error "error: configuration file $OPTARG not found."
                       exit 1
                    fi
-            esac
-            shift 2;;
-        -i|--install)
-            case "$2" in
+            esac;;
+        i)
+            case "$OPTARG" in
                 "") help ; break;;
-                *) PACKAGES="${PACKAGE} $2"
-            esac
-            shift 2;;
-        -h|-help|--help)
+                *) PACKAGES="${PACKAGE} $OPTARG"
+            esac;;
+        h)
             help
             exit 0;;
-        --) shift ; break ;;
         *)
-            error $(basename $0)": unknown option $1, try --help"
+            error $(basename $0)": unknown option $c, try --help"
             exit 2;;
     esac
 done
diff --git a/scripts/verif_omp/check_actual_setup.sh b/scripts/verif_omp/check_actual_setup.sh
new file mode 100755
index 00000000..6eaa4517
--- /dev/null
+++ b/scripts/verif_omp/check_actual_setup.sh
@@ -0,0 +1,12 @@
+#!/bin/sh
+
+# go in qp2/src/fci to run check_omp_actual_setup
+# to see if we can run in parallel an omp section in another one
+echo ""
+echo "Please wait..."
+echo ""
+cd ../../src/fci
+ninja || echo "Please recompile from the root"  
+echo ""
+./check_omp_actual_setup
+cd ../../scripts/verif_omp
diff --git a/scripts/verif_omp/check_omp_v2.f90 b/scripts/verif_omp/check_omp_v2.f90
new file mode 100644
index 00000000..ca6af8bd
--- /dev/null
+++ b/scripts/verif_omp/check_omp_v2.f90
@@ -0,0 +1,175 @@
+program check_omp_v2
+ 
+  use omp_lib
+
+  implicit none
+
+  integer :: accu, accu2
+  integer :: s, n_setting
+  logical :: verbose, test_versions
+  logical, allocatable :: is_working(:)
+
+  verbose = .False.
+  test_versions = .True.
+  n_setting = 4
+
+  allocate(is_working(n_setting))
+
+  is_working = .False.
+
+  ! set the number of threads
+  call omp_set_num_threads(2) 
+
+  do s = 1, n_setting
+
+    accu = 0
+    accu2 = 0
+
+    call omp_set_max_active_levels(1)
+    call omp_set_nested(.False.)
+
+    if (s==1) then
+      !call set_multiple_levels_omp()
+      cycle
+    elseif (s==2) then
+      call omp_set_max_active_levels(5)
+    elseif (s==3) then
+      call omp_set_nested(.True.)
+    else
+      call omp_set_nested(.True.)
+      call omp_set_max_active_levels(5)
+    endif
+
+    ! Level 1
+    !$OMP PARALLEL 
+      if (verbose) then
+        print*,'Num threads level 1:',omp_get_num_threads()
+      endif
+
+      ! Level 2
+      !$OMP PARALLEL 
+        if (verbose) then
+          print*,'Num threads level 2:',omp_get_num_threads()
+        endif
+
+        ! Level 3
+        !$OMP PARALLEL 
+          if (verbose) then
+            print*,'Num threads level 3:',omp_get_num_threads()
+          endif
+
+          call check_omp_in_subroutine(accu2)
+     
+          ! Level 4
+          !$OMP PARALLEL 
+
+            if (verbose) then
+              print*,'Num threads level 4:',omp_get_num_threads()
+            endif
+
+           !$OMP ATOMIC    
+           accu = accu + 1
+           !$OMP END ATOMIC 
+       
+          !$OMP END PARALLEL
+     
+     
+        !$OMP END PARALLEL
+
+
+      !$OMP END PARALLEL
+
+
+    !$OMP END PARALLEL
+
+  if (verbose) then
+    print*,'Setting:',s,'accu=',accu
+    print*,'Setting:',s,'accu2=',accu2
+  endif
+
+  if (accu == 16 .and. accu2 == 16) then
+    is_working(s) = .True.
+  endif
+
+  enddo
+
+  if (verbose) then
+    if (is_working(2)) then
+      print*,'The parallelization works on 4 levels with:'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,''
+      print*,'Please use the irpf90 flags -DSET_MAX_ACT in qp2/config/${compiler_name}.cfg'
+    elseif (is_working(3)) then
+      print*,'The parallelization works on 4 levels with:'
+      print*,'call omp_set_nested(.True.)'
+      print*,''
+      print*,'Please use the irpf90 flag -DSET_NESTED in qp2/config/${compiler_name}.cfg'
+    elseif (is_working(4)) then
+      print*,'The parallelization works on 4 levels with:'
+      print*,'call omp_set_nested(.True.)'
+      print*,'+'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,''
+      print*,'Please use the irpf90 flags -DSET_NESTED -DSET_MAX_ACT in qp2/config/${compiler_name}.cfg'
+    else
+      print*,'The parallelization on multiple levels does not work with:'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,'or'
+      print*,'call omp_set_nested(.True.)'
+      print*,'or'
+      print*,'call omp_set_nested(.True.)'
+      print*,'+'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,''
+      print*,'Try an other compiler and good luck...'
+    endif
+
+   ! if (is_working(1)) then
+   !   print*,''
+   !   print*,'=========================================================='
+   !   print*,'Your actual set up works for parallelization with 4 levels'
+   !   print*,'=========================================================='
+   !   print*,''
+   ! else
+   !   print*,''
+   !   print*,'==================================================================='
+   !   print*,'Your actual set up does not work for parallelization with 4 levels'
+   !   print*,'Please look at the previous messages to understand the requirements'
+   !   print*,'==================================================================='
+   !   print*,''
+   ! endif
+  endif
+
+  ! List of working flags
+  if (test_versions) then
+    print*,'Tests:',is_working(2:4)
+  endif
+
+  ! IRPF90_FLAGS
+  if (is_working(2)) then
+    print*,'-DSET_MAX_ACT'
+  elseif (is_working(3)) then
+    print*,'-DSET_NESTED'
+  elseif (is_working(4)) then
+    print*,'-DSET_MAX_ACT -DSET_NESTED'
+  else
+    print*,'ERROR'
+  endif
+
+end
+
+subroutine check_omp_in_subroutine(accu2)
+
+  implicit none
+ 
+  integer, intent(inout) :: accu2
+
+  !$OMP PARALLEL
+
+  !$OMP ATOMIC
+  accu2 = accu2 + 1
+  !$OMP END ATOMIC
+
+  !$OMP END PARALLEL
+
+end
diff --git a/scripts/verif_omp/check_required_setup.sh b/scripts/verif_omp/check_required_setup.sh
new file mode 100755
index 00000000..facb6cbb
--- /dev/null
+++ b/scripts/verif_omp/check_required_setup.sh
@@ -0,0 +1,19 @@
+#!/bin/sh
+
+# take one argument which is the compiler used
+# return the required IRPF90_FLAGS for the $1 compiler
+
+if [ -z "$1" ]
+then
+    echo "Give the compiler in argument"
+else
+
+$1 --version > /dev/null \
+&& $1 -O0 -fopenmp check_omp_v2.f90 \
+&& ./a.out | tail -n 1
+
+
+# if there is an error or if the compiler is not found
+$1 --version > /dev/null || echo 'compiler not found'
+
+fi
diff --git a/scripts/verif_omp/study_omp.sh b/scripts/verif_omp/study_omp.sh
new file mode 100755
index 00000000..00668d59
--- /dev/null
+++ b/scripts/verif_omp/study_omp.sh
@@ -0,0 +1,30 @@
+#!/bin/sh
+
+# list of compilers
+list_comp="ifort gfortran-7 gfortran-8 gfortran-9"
+
+# file to store the results
+FILE=results.dat
+
+touch $FILE
+rm $FILE
+
+# Comments
+echo "1: omp_set_max_active_levels(5)" >> $FILE
+echo "2: omp_set_nested(.True.)" >> $FILE
+echo "3: 1 + 2" >> $FILE
+echo "" >> $FILE
+echo "1 2 3" >> $FILE
+
+# loop on the comp
+for comp in $list_comp
+do
+	$comp --version > /dev/null \
+        && $comp -O0 -fopenmp check_omp_v2.f90 \
+	&& echo $(./a.out | grep "Tests:" | cut -d ":" -f2- ) $(echo " : ") $($comp --version | head -n 1) >> $FILE
+
+done
+
+# Display
+cat $FILE
+
diff --git a/scripts/verif_omp/update_comp.sh b/scripts/verif_omp/update_comp.sh
new file mode 100755
index 00000000..14b644de
--- /dev/null
+++ b/scripts/verif_omp/update_comp.sh
@@ -0,0 +1,49 @@
+#!/bin/bash
+
+# Compiler
+COMP=$1
+
+# Path to file.cfg
+config_PATH="../../config/"
+END="*.cfg"
+CONFIG="/config/"
+
+#LIST=${config_PATH}${COMP}${END} # without ${QP_ROOT}
+LIST=${QP_ROOT}${CONFIG}${COMP}${END}
+
+if [ -z "$1" ]
+then
+    echo "Give the compiler in argument"
+else
+
+    # List of the config files for the compiler
+    #list_files=$(ls ../../config/$comp*.cfg) #does not give the right list
+    list_files=${LIST}
+    echo "Files that will be modified:"
+    echo $list_files
+   
+    # Flags that must be added 
+    FLAGS=$(./check_required_setup.sh $COMP)
+
+    # Add the flags
+    for file in $list_files
+    do
+        echo $file
+        BASE="IRPF90_FLAGS : --ninja"
+        ACTUAL=$(grep "$BASE" $file)
+        
+        # To have only one time each flag	
+	grep " -DSET_MAX_ACT" $file && ${ACTUAL/" -DSET_MAX"/""}
+	grep " -DSET_NESTED" $file && ${ACTUAL/" -DSET_NESTED"/""}
+        SPACE=" "
+
+	NEW=${ACTUAL}${SPACE}${FLAGS}
+
+	# Debug
+        #echo ${NEW}
+
+        sed "s/${ACTUAL}/${NEW}/" $file
+        # -i # to change the files
+    done
+
+fi
diff --git a/src/ao_one_e_ints/pot_ao_pseudo_ints.irp.f b/src/ao_one_e_ints/pot_ao_pseudo_ints.irp.f
index 24f43311..e75ca056 100644
--- a/src/ao_one_e_ints/pot_ao_pseudo_ints.irp.f
+++ b/src/ao_one_e_ints/pot_ao_pseudo_ints.irp.f
@@ -28,6 +28,7 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integrals, (ao_num,ao_num)]
 END_PROVIDER
 
 BEGIN_PROVIDER [ double precision, ao_pseudo_integrals_local, (ao_num,ao_num)]
+  use omp_lib
   implicit none
   BEGIN_DOC
   ! Local pseudo-potential
@@ -42,7 +43,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integrals_local, (ao_num,ao_num)]
 
   double precision               :: wall_1, wall_2, wall_0
   integer                        :: thread_num
-  integer                        :: omp_get_thread_num
   double precision               :: c
   double precision               :: Z
 
@@ -158,6 +158,7 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integrals_local, (ao_num,ao_num)]
 
 
  BEGIN_PROVIDER [ double precision, ao_pseudo_integrals_non_local, (ao_num,ao_num)]
+  use omp_lib
   implicit none
   BEGIN_DOC
   ! Non-local pseudo-potential
@@ -169,7 +170,6 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integrals_local, (ao_num,ao_num)]
   integer                        :: power_A(3),power_B(3)
   integer                        :: i,j,k,l,m
   double precision               :: Vloc, Vpseudo
-  integer                        :: omp_get_thread_num
 
   double precision               :: wall_1, wall_2, wall_0
   integer                        :: thread_num
diff --git a/src/cis/EZFIO.cfg b/src/cis/EZFIO.cfg
index 7e0eeb03..955d1bef 100644
--- a/src/cis/EZFIO.cfg
+++ b/src/cis/EZFIO.cfg
@@ -5,4 +5,3 @@ interface: ezfio
 size: (determinants.n_states)
 
 
-
diff --git a/src/cis/cis.irp.f b/src/cis/cis.irp.f
index acec29c2..ab2294ad 100644
--- a/src/cis/cis.irp.f
+++ b/src/cis/cis.irp.f
@@ -79,6 +79,6 @@ subroutine run
   call ezfio_set_cis_energy(CI_energy)
   psi_coef = ci_eigenvectors
   SOFT_TOUCH psi_coef
-  call save_wavefunction_truncated(1.d-12)
+  call save_wavefunction_truncated(save_threshold)
 
 end
diff --git a/src/cisd/cisd.irp.f b/src/cisd/cisd.irp.f
index 6c55e2ff..fca3b10e 100644
--- a/src/cisd/cisd.irp.f
+++ b/src/cisd/cisd.irp.f
@@ -63,7 +63,7 @@ subroutine run
   endif
   psi_coef = ci_eigenvectors
   SOFT_TOUCH psi_coef
-  call save_wavefunction
+  call save_wavefunction_truncated(save_threshold)
   call ezfio_set_cisd_energy(CI_energy)
 
   do i = 1,N_states
diff --git a/src/csf/sigma_vector.irp.f b/src/csf/sigma_vector.irp.f
index 77f6190b..0d24ae57 100644
--- a/src/csf/sigma_vector.irp.f
+++ b/src/csf/sigma_vector.irp.f
@@ -59,24 +59,30 @@
       endif
     endif
     ncfg = ncfgpersomo - ncfgprev
-    if(iand(MS,1) .EQ. 0) then
-      !dimcsfpercfg = max(1,nint((binom(i,i/2)-binom(i,i/2+1))))
-      binom1 = dexp(logabsgamma(1.0d0*(i+1))                            &
-                  - logabsgamma(1.0d0*((i/2)+1))                        &
-                  - logabsgamma(1.0d0*(i-((i/2))+1)));
-      binom2 = dexp(logabsgamma(1.0d0*(i+1))                            &
-                  - logabsgamma(1.0d0*(((i/2)+1)+1))                    &
-                  - logabsgamma(1.0d0*(i-((i/2)+1)+1)));
-      dimcsfpercfg = max(1,nint(binom1 - binom2))
+    if(i .EQ. 0 .OR. i .EQ. 1) then
+      dimcsfpercfg = 1
+    elseif( i .EQ. 3) then
+      dimcsfpercfg = 2
     else
-      !dimcsfpercfg = max(1,nint((binom(i,(i+1)/2)-binom(i,(i+3)/2))))
-      binom1 = dexp(logabsgamma(1.0d0*(i+1))                            &
-                  - logabsgamma(1.0d0*(((i+1)/2)+1))                    &
-                  - logabsgamma(1.0d0*(i-(((i+1)/2))+1)));
-      binom2 = dexp(logabsgamma(1.0d0*(i+1))                            &
-                  - logabsgamma(1.0d0*((((i+3)/2)+1)+1))                &
-                  - logabsgamma(1.0d0*(i-(((i+3)/2)+1)+1)));
-      dimcsfpercfg = max(1,nint(binom1 - binom2))
+      if(iand(MS,1) .EQ. 0) then
+        !dimcsfpercfg = max(1,nint((binom(i,i/2)-binom(i,i/2+1))))
+        binom1 = dexp(logabsgamma(1.0d0*(i+1))                            &
+                    - logabsgamma(1.0d0*((i/2)+1))                        &
+                    - logabsgamma(1.0d0*(i-((i/2))+1)));
+        binom2 = dexp(logabsgamma(1.0d0*(i+1))                            &
+                    - logabsgamma(1.0d0*(((i/2)+1)+1))                    &
+                    - logabsgamma(1.0d0*(i-((i/2)+1)+1)));
+        dimcsfpercfg = max(1,nint(binom1 - binom2))
+      else
+        !dimcsfpercfg = max(1,nint((binom(i,(i+1)/2)-binom(i,(i+3)/2))))
+        binom1 = dexp(logabsgamma(1.0d0*(i+1))                            &
+                    - logabsgamma(1.0d0*(((i+1)/2)+1))                    &
+                    - logabsgamma(1.0d0*(i-(((i+1)/2))+1)));
+        binom2 = dexp(logabsgamma(1.0d0*(i+1))                            &
+                    - logabsgamma(1.0d0*((((i+3)/2)+1)+1))                &
+                    - logabsgamma(1.0d0*(i-(((i+3)/2)+1)+1)));
+        dimcsfpercfg = max(1,nint(binom1 - binom2))
+      endif
     endif
     n_CSF += ncfg * dimcsfpercfg
     if(cfg_seniority_index(i+2) > ncfgprev) then
diff --git a/src/davidson/davidson_parallel.irp.f b/src/davidson/davidson_parallel.irp.f
index 8fd023da..fcee16bc 100644
--- a/src/davidson/davidson_parallel.irp.f
+++ b/src/davidson/davidson_parallel.irp.f
@@ -508,7 +508,8 @@ subroutine H_S2_u_0_nstates_zmq(v_0,s_0,u_0,N_st,sze)
   endif
 
 
-  call omp_set_max_active_levels(5)
+  !call omp_set_max_active_levels(5)
+  call set_multiple_levels_omp()
 
   !$OMP PARALLEL DEFAULT(shared) NUM_THREADS(2) PRIVATE(ithread)
   ithread = omp_get_thread_num()
diff --git a/src/davidson/davidson_parallel_csf.irp.f b/src/davidson/davidson_parallel_csf.irp.f
index fe651b1d..90e4303e 100644
--- a/src/davidson/davidson_parallel_csf.irp.f
+++ b/src/davidson/davidson_parallel_csf.irp.f
@@ -464,7 +464,9 @@ subroutine H_u_0_nstates_zmq(v_0,u_0,N_st,sze)
     print *,  irp_here, ': Failed in zmq_set_running'
   endif
 
-  call omp_set_max_active_levels(4)
+  !call omp_set_max_active_levels(4)
+  call set_multiple_levels_omp()
+
   !$OMP PARALLEL DEFAULT(shared) NUM_THREADS(2) PRIVATE(ithread)
   ithread = omp_get_thread_num()
   if (ithread == 0 ) then
diff --git a/src/davidson/davidson_parallel_nos2.irp.f b/src/davidson/davidson_parallel_nos2.irp.f
index 84cbe3af..091b8666 100644
--- a/src/davidson/davidson_parallel_nos2.irp.f
+++ b/src/davidson/davidson_parallel_nos2.irp.f
@@ -464,7 +464,9 @@ subroutine H_u_0_nstates_zmq(v_0,u_0,N_st,sze)
     print *,  irp_here, ': Failed in zmq_set_running'
   endif
 
-  call omp_set_max_active_levels(4)
+  !call omp_set_max_active_levels(4)
+  call set_multiple_levels_omp()
+
   !$OMP PARALLEL DEFAULT(shared) NUM_THREADS(2) PRIVATE(ithread)
   ithread = omp_get_thread_num()
   if (ithread == 0 ) then
diff --git a/src/determinants/EZFIO.cfg b/src/determinants/EZFIO.cfg
index 662c6fbb..5e109de8 100644
--- a/src/determinants/EZFIO.cfg
+++ b/src/determinants/EZFIO.cfg
@@ -48,7 +48,7 @@ default: 1
 
 [threshold_generators]
 type: Threshold
-doc: Thresholds on generators (fraction of the square of the norm) 
+doc: Thresholds on generators (fraction of the square of the norm)
 interface: ezfio,provider,ocaml
 default: 0.999
 
@@ -80,7 +80,7 @@ type: integer
 [psi_coef]
 interface: ezfio
 doc: Coefficients of the wave function
-type: double precision 
+type: double precision
 size: (determinants.n_det,determinants.n_states)
 
 [psi_det]
@@ -92,7 +92,7 @@ size: (determinants.n_int*determinants.bit_kind/8,2,determinants.n_det)
 [psi_coef_qp_edit]
 interface: ezfio
 doc: Coefficients of the wave function
-type: double precision 
+type: double precision
 size: (determinants.n_det_qp_edit,determinants.n_states)
 
 [psi_det_qp_edit]
@@ -126,13 +126,19 @@ default: 1.
 
 [thresh_sym]
 type: Threshold
-doc: Thresholds to check if a determinant is connected with HF 
+doc: Thresholds to check if a determinant is connected with HF
 interface: ezfio,provider,ocaml
 default: 1.e-15
 
 [pseudo_sym]
 type: logical
-doc: If |true|, discard any Slater determinants with an interaction smaller than thresh_sym with HF. 
+doc: If |true|, discard any Slater determinants with an interaction smaller than thresh_sym with HF.
 interface: ezfio,provider,ocaml
 default: False
 
+[save_threshold]
+type: Threshold
+doc: Cut-off to apply to the CI coefficients when the wave function is stored
+interface: ezfio,provider,ocaml
+default: 0.
+
diff --git a/src/fci/check_omp_actual_setup.irp.f b/src/fci/check_omp_actual_setup.irp.f
new file mode 100644
index 00000000..70514bd3
--- /dev/null
+++ b/src/fci/check_omp_actual_setup.irp.f
@@ -0,0 +1,174 @@
+program check_omp_actual_setup
+ 
+  use omp_lib
+
+  implicit none
+
+  integer :: accu, accu2
+  integer :: s, n_setting
+  logical :: verbose, test_versions
+  logical, allocatable :: is_working(:)
+
+  verbose = .True.
+  test_versions = .False.
+  n_setting = 4
+
+  allocate(is_working(n_setting))
+
+  is_working = .False.
+
+  ! set the number of threads
+  call omp_set_num_threads(2) 
+
+  do s = 1, n_setting
+
+    accu = 0
+    accu2 = 0
+
+    call omp_set_max_active_levels(1)
+    call omp_set_nested(.False.)
+
+    if (s==1) then
+      call set_multiple_levels_omp()
+    elseif (s==2) then
+      call omp_set_max_active_levels(5)
+    elseif (s==3) then
+      call omp_set_nested(.True.)
+    else
+      call omp_set_nested(.True.)
+      call omp_set_max_active_levels(5)
+    endif
+
+    ! Level 1
+    !$OMP PARALLEL 
+      if (verbose) then
+        print*,'Num threads level 1:',omp_get_num_threads()
+      endif
+
+      ! Level 2
+      !$OMP PARALLEL 
+        if (verbose) then
+          print*,'Num threads level 2:',omp_get_num_threads()
+        endif
+
+        ! Level 3
+        !$OMP PARALLEL 
+          if (verbose) then
+            print*,'Num threads level 3:',omp_get_num_threads()
+          endif
+
+          call check_omp_in_subroutine(accu2)
+     
+          ! Level 4
+          !$OMP PARALLEL 
+
+            if (verbose) then
+              print*,'Num threads level 4:',omp_get_num_threads()
+            endif
+
+           !$OMP ATOMIC    
+           accu = accu + 1
+           !$OMP END ATOMIC 
+       
+          !$OMP END PARALLEL
+     
+     
+        !$OMP END PARALLEL
+
+
+      !$OMP END PARALLEL
+
+
+    !$OMP END PARALLEL
+
+  if (verbose) then
+    print*,'Setting:',s,'accu=',accu
+    print*,'Setting:',s,'accu2=',accu2
+  endif
+
+  if (accu == 16 .and. accu2 == 16) then
+    is_working(s) = .True.
+  endif
+
+  enddo
+
+  if (verbose) then
+    if (is_working(2)) then
+      print*,'The parallelization works on 4 levels with:'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,''
+      print*,'Please use the irpf90 flags -DSET_MAX_ACT in qp2/config/${compiler_name}.cfg'
+    elseif (is_working(3)) then
+      print*,'The parallelization works on 4 levels with:'
+      print*,'call omp_set_nested(.True.)'
+      print*,''
+      print*,'Please use the irpf90 flag -DSET_NESTED in qp2/config/${compiler_name}.cfg'
+    elseif (is_working(4)) then
+      print*,'The parallelization works on 4 levels with:'
+      print*,'call omp_set_nested(.True.)'
+      print*,'+'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,''
+      print*,'Please use the irpf90 flags -DSET_NESTED -DSET_MAX_ACT in qp2/config/${compiler_name}.cfg'
+    else
+      print*,'The parallelization on multiple levels does not work with:'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,'or'
+      print*,'call omp_set_nested(.True.)'
+      print*,'or'
+      print*,'call omp_set_nested(.True.)'
+      print*,'+'
+      print*,'call omp_set_max_active_levels(5)'
+      print*,''
+      print*,'Try an other compiler and good luck...'
+    endif
+
+    if (is_working(1)) then
+      print*,''
+      print*,'=========================================================='
+      print*,'Your actual set up works for parallelization with 4 levels'
+      print*,'=========================================================='
+      print*,''
+    else
+      print*,''
+      print*,'==================================================================='
+      print*,'Your actual set up does not work for parallelization with 4 levels'
+      print*,'Please look at the previous messages to understand the requirements'
+      print*,'==================================================================='
+      print*,''
+    endif
+  endif
+
+  ! List of working flags
+  if (test_versions) then
+    print*,is_working(2:4)
+  endif
+
+  ! IRPF90_FLAGS
+  if (is_working(2)) then
+    print*,'-DSET_MAX_ACT'
+  elseif (is_working(3)) then
+    print*,'-DSET_NESTED'
+  elseif (is_working(4)) then
+    print*,'-DSET_MAX_ACT -DSET_NESTED'
+  else
+    print*,'ERROR'
+  endif
+
+end
+
+subroutine check_omp_in_subroutine(accu2)
+
+  implicit none
+ 
+  integer, intent(inout) :: accu2
+
+  !$OMP PARALLEL
+
+  !$OMP ATOMIC
+  accu2 = accu2 + 1
+  !$OMP END ATOMIC
+
+  !$OMP END PARALLEL
+
+end
diff --git a/src/utils/set_multiple_levels_omp.irp.f b/src/utils/set_multiple_levels_omp.irp.f
new file mode 100644
index 00000000..4be3af5b
--- /dev/null
+++ b/src/utils/set_multiple_levels_omp.irp.f
@@ -0,0 +1,16 @@
+subroutine set_multiple_levels_omp()
+
+! Doc : idk
+
+  implicit none
+
+  IRP_IF SET_MAX_ACT
+    !print*,'SET_MAX_ACT: True, call omp_set_max_active_levels(5)'
+    call omp_set_max_active_levels(5)
+  IRP_ENDIF
+  IRP_IF SET_NESTED
+    !print*,'SET_NESTED: True, call omp_set_nested(.True.)'
+    call omp_set_nested(.True.)
+  IRP_ENDIF
+
+end
diff --git a/src/utils/util.irp.f b/src/utils/util.irp.f
index cfb42fd1..ef846bdb 100644
--- a/src/utils/util.irp.f
+++ b/src/utils/util.irp.f
@@ -300,12 +300,12 @@ subroutine wall_time(t)
 end
 
 BEGIN_PROVIDER [ integer, nproc ]
+  use omp_lib
   implicit none
   BEGIN_DOC
   ! Number of current OpenMP threads
   END_DOC
 
-  integer                        :: omp_get_num_threads
   nproc = 1
   !$OMP PARALLEL
   !$OMP MASTER