diff --git a/.gitattributes b/.gitattributes new file mode 100644 index 00000000..99956fdc --- /dev/null +++ b/.gitattributes @@ -0,0 +1,5 @@ +*.irp.f linguist-language=IRPF90 +*.f linguist-language=Fortran +*.ml linguist-language=Ocaml +*.sh linguist-language=Bash +*.py linguist-language=Python diff --git a/INSTALL.md b/INSTALL.md index 712e3d8f..47319c9c 100644 --- a/INSTALL.md +++ b/INSTALL.md @@ -6,8 +6,9 @@ * m4 * GNU make * Fortran compiler (ifort or gfortran are tested) -* Python >= 2.7 +* Python >= 2.6 * Bash +* Patch (for opam) ## Standard installation diff --git a/data/ezfio_defaults/q_package.ezfio_default b/data/ezfio_defaults/WILL_BE_DELETED.ezfio_default similarity index 95% rename from data/ezfio_defaults/q_package.ezfio_default rename to data/ezfio_defaults/WILL_BE_DELETED.ezfio_default index 5da4d813..cf54a1dd 100644 --- a/data/ezfio_defaults/q_package.ezfio_default +++ b/data/ezfio_defaults/WILL_BE_DELETED.ezfio_default @@ -32,16 +32,17 @@ full_ci do_pt2_end True var_pt2_ratio 0.75 +cas_sd + n_det_max_cas_sd 100000 + pt2_max 1.e-4 + do_pt2_end True + var_pt2_ratio 0.75 + all_singles n_det_max_fci 50000 pt2_max 1.e-8 do_pt2_end False -cassd - n_det_max_cassd 10000 - pt2_max 1.e-4 - do_pt2_end True - hartree_fock n_it_scf_max 200 thresh_scf 1.e-10 diff --git a/data/ezfio_defaults/determinants.ezfio_default b/data/ezfio_defaults/determinants.ezfio_default index 42f6b384..2cfbe3ea 100644 --- a/data/ezfio_defaults/determinants.ezfio_default +++ b/data/ezfio_defaults/determinants.ezfio_default @@ -4,6 +4,6 @@ determinants n_det_max_jacobi 1000 threshold_generators 0.99 threshold_selectors 0.999 - read_wf False - s2_eig False - only_single_double_dm False + read_wf false + s2_eig false + only_single_double_dm false diff --git a/ocaml/Makefile b/ocaml/Makefile index 25598e1f..eaa8e382 100644 --- a/ocaml/Makefile +++ b/ocaml/Makefile @@ -25,9 +25,12 @@ default: $(ALL_TESTS) $(ALL_EXE) .gitignore .gitignore: $(MLFILES) @for i in .gitignore ezfio.ml Qptypes.ml qptypes_generator.byte _build $(ALL_EXE) $(ALL_TESTS) \ - $(patsubst %.ml,%,$(wildcard test_*.ml)) $(patsubst %.ml,%,$(wildcard qp_*.ml)) ; do \ + $(patsubst %.ml,%,$(wildcard test_*.ml)) $(patsubst %.ml,%,$(wildcard qp_*.ml)) \ + $(shell grep Input Input_auto_generated.ml | awk '{print $$2 ".ml"}') \ + Input_auto_generated.ml;\ + do \ echo $$i ; \ - done >> .gitignore + done > .gitignore executables: $(QPACKAGE_ROOT)/data/executables diff --git a/ocaml/Qpackage.ml b/ocaml/Qpackage.ml index 4d08192a..a3862f11 100644 --- a/ocaml/Qpackage.ml +++ b/ocaml/Qpackage.ml @@ -94,7 +94,7 @@ let get_ezfio_default_in_file ~directory ~data ~filename = match (String.lsplit2 ~on:' ' (String.strip line)) with | Some (l,r) -> if (l = data) then - String.lowercase (String.strip r) + String.strip r else find_data rest | None -> raise Not_found diff --git a/ocaml/qp_run.ml b/ocaml/qp_run.ml index 341a4ce6..2abea107 100644 --- a/ocaml/qp_run.ml +++ b/ocaml/qp_run.ml @@ -20,6 +20,10 @@ let run exe ezfio_file = Printf.printf "===============\nQuantum Package\n===============\n\n"; Printf.printf "Date : %s\n\n%!" (Time.to_string time_start); + match (Sys.command ("qp_edit -c "^ezfio_file)) with + | 0 -> () + | i -> failwith "Error: Input inconsistent\n"; + ; let exe = match (List.find ~f:(fun (x,_) -> x = exe) executables) with | None -> assert false diff --git a/ocaml/qptypes_generator.ml b/ocaml/qptypes_generator.ml index 633727f8..fd1739f1 100644 --- a/ocaml/qptypes_generator.ml +++ b/ocaml/qptypes_generator.ml @@ -150,15 +150,12 @@ end = struct let to_string = function | Huckel -> \"Huckel\" - | HCore -> \"HCore\" + | HCore -> \"Hcore\" let of_string s = - let s = - String.lowercase s - in match s with - | \"huckel\" -> Huckel - | \"hcore\" -> HCore + | \"Huckel\" -> Huckel + | \"Hcore\" -> HCore | _ -> failwith (\"Wrong Guess type : \"^s) end @@ -179,13 +176,10 @@ end = struct | Write -> \"Write\" | None -> \"None\" let of_string s = - let s = - String.lowercase s - in match s with - | \"read\" -> Read - | \"write\" -> Write - | \"none\" -> None + | \"Read\" -> Read + | \"Write\" -> Write + | \"None\" -> None | _ -> failwith (\"Wrong IO type : \"^s) end diff --git a/ocaml/test_input.ml b/ocaml/test_input.ml deleted file mode 100644 index 956fd745..00000000 --- a/ocaml/test_input.ml +++ /dev/null @@ -1,177 +0,0 @@ -open Qptypes;; - -let test_ao () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Ao_basis.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Ao_basis.to_string b); - print_endline (Input.Ao_basis.to_rst b |> Rst_string.to_string); -;; - -let test_bielec_intergals () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Bielec_integrals.read () with - | Some x -> x - | None -> assert false - in - let output = Input.Bielec_integrals.to_string b - in - print_endline output; - let rst = Input.Bielec_integrals.to_rst b in - let b2 = match Input.Bielec_integrals.of_rst rst with - | Some x -> x - | None -> assert false - in - if (b = b2) then - print_endline "OK" - else - print_endline "rst failed"; -;; - -let test_bitmasks () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Bitmasks.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Bitmasks.to_string b); -;; - - -let test_dets () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Determinants.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Determinants.to_rst b |> Rst_string.to_string ) ; - print_endline (Input.Determinants.sexp_of_t b |> Sexplib.Sexp.to_string ) ; - let rst = Input.Determinants.to_rst b in - let b2 = match Input.Determinants.of_rst rst with - | Some x -> x - | None -> assert false - in - if (b2 = b) then - print_endline "OK" - else - print_endline "Failed" -;; - -let test_cisd_sc2 () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Cisd_sc2_selected.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Cisd_sc2_selected.to_string b); - let rst = Input.Cisd_sc2_selected.to_rst b in - let b2 = match Input.Cisd_sc2_selected.of_rst rst with - | Some x -> x - | None -> assert false - in - if (b = b2) then - print_endline "OK" - else - print_endline "rst failed"; - -;; - -let test_electrons () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Electrons.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Electrons.to_string b); - let rst = Input.Electrons.to_rst b in - let b2 = match Input.Electrons.of_rst rst with - | Some x -> x - | None -> assert false - in - if (b = b2) then - print_endline "OK" - else - print_endline "Failed in rst" -;; - -let test_fci () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Full_ci.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Full_ci.to_string b); - let rst = Input.Full_ci.to_rst b in - let b2 = match Input.Full_ci.of_rst rst with - | Some x -> x - | None -> assert false - in - print_endline (Input.Full_ci.to_string b); - if (b = b2) then - print_endline "OK" - else - print_endline "Failed in rst" -;; - -let test_hf () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Hartree_fock.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Hartree_fock.to_string b); - let rst = Input.Hartree_fock.to_rst b in - let b2 = match Input.Hartree_fock.of_rst rst with - | Some x -> x - | None -> assert false - in - print_endline (Input.Hartree_fock.to_string b); - if (b = b2) then - print_endline "OK" - else - print_endline "Failed in rst" -;; - -let test_mo () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Mo_basis.read () with - | Some x -> x - | None -> assert false - in - print_endline (Input.Mo_basis.to_string b); -;; - -let test_nucl () = - Ezfio.set_file "F2.ezfio" ; - let b = match Input.Nuclei.read () with - | Some x -> x - | None -> assert false - in - let rst = Input.Nuclei.to_rst b in - let b2 = match Input.Nuclei.of_rst rst with - | Some x -> x - | None -> assert false - in - print_endline (Input.Nuclei.to_string b); - if (b = b2) then - print_endline "OK" - else - print_endline "Failed in rst" -;; - -(* -test_ao ();; -test_bitmasks (); -test_cis (); -test_cisd_sc2 (); -test_dets (); -test_hf ();; -test_mo ();; -test_nucl (); -test_bielec_intergals ();; -test_electrons(); -*) -test_dets (); - diff --git a/scripts/build_modules.sh b/scripts/build_modules.sh index 5607c72b..c3e0cda5 100755 --- a/scripts/build_modules.sh +++ b/scripts/build_modules.sh @@ -30,6 +30,7 @@ Build failed for module $MODULE " fi fi + ${QPACKAGE_ROOT}/scripts/create_gitignore.sh cd ${OLDPWD} done ${QPACKAGE_ROOT}/scripts/create_executables_list.sh diff --git a/scripts/create_gitignore.sh b/scripts/create_gitignore.sh index 605b0ab8..b42cb24b 100755 --- a/scripts/create_gitignore.sh +++ b/scripts/create_gitignore.sh @@ -35,4 +35,4 @@ then fi find . -type l | sed "s@./@@" >> .gitignore - +find . -type f -executable -print | sed "s@./@@" >> .gitignore \ No newline at end of file diff --git a/scripts/create_rst_templates.sh b/scripts/create_rst_templates.sh index 931701ee..4d85b827 100755 --- a/scripts/create_rst_templates.sh +++ b/scripts/create_rst_templates.sh @@ -12,6 +12,8 @@ fi source ${QPACKAGE_ROOT}/scripts/qp_include.sh check_current_dir_is_module +MODULE=$(basename $PWD) + README="True" if [[ -f README.rst ]] diff --git a/scripts/ezfio_interface/ei_handler.py b/scripts/ezfio_interface/ei_handler.py index 8fad8511..a02b13d5 100755 --- a/scripts/ezfio_interface/ei_handler.py +++ b/scripts/ezfio_interface/ei_handler.py @@ -12,8 +12,8 @@ By default all the option are executed. Options: -h --help --path The path of the `EZFIO.cfg`, by default will look in the ${pwd} - --irpf90 Create the `ezfio_interface.ipf90` - who containt all the provider needed + --irpf90 Create the `ezfio_interface.irpf90` + which contains all the providers needed (aka all with the `interface: input` parameter) in `${pwd}` --ezfio_config Create the `${module_lower}_ezfio_interface_config` in @@ -28,7 +28,7 @@ Options: Format specification : [provider_name] | the name of the provider in irp.f90 doc:{str} | Is the doc - Type:{str} | Is a fancy_type support by the ocaml + Type:{str} | Is a fancy_type supported by the ocaml ezfio_name:{str} | Will be the name of the file for the ezfio (optional by default is the name of the provider) interface:{str} | The provider is a imput or a output @@ -67,13 +67,13 @@ Type = namedtuple('Type', 'fancy ocaml fortran') def is_bool(str_): """ - Take a string, if is a bool return the convert into - fortran and ocaml one. + Take a string, if is a bool return the conversion into + fortran and ocaml. """ - if str_.lower() in ['true', '.true.']: - return Type(None, "True", ".True.") - elif str_.lower() in ['false', '.False.']: - return Type(None, "False", ".False") + if "true" in str_.lower(): + return Type(None, "true", ".True.") + elif "false" in str_.lower(): + return Type(None, "false", ".False") else: raise TypeError @@ -81,7 +81,7 @@ def is_bool(str_): def get_type_dict(): """ This function makes the correspondance between the type of value read in - ezfio.cfg into the f90 and Ocam Type + ezfio.cfg into the f90 and Ocaml Type return fancy_type[fancy_type] = namedtuple('Type', 'ocaml fortran') For example fancy_type['Ndet'].fortran = interger .ocaml = int @@ -154,7 +154,7 @@ def get_type_dict(): # q p _ t y p e s _ g e n e r a t e # # ~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~#~ # - # Read the fancy_type, the ocaml. and convert the ocam to the fortran + # Read the fancy_type, the ocaml. and convert the ocaml to the fortran for i in l_gen + l_un: str_fancy_type = i.split()[1].strip() str_ocaml_type = i.split()[3] @@ -507,6 +507,9 @@ def create_ocaml_input(dict_ezfio_cfg, module_lower): l_type.append(v["type"]) l_doc.append(v["doc"]) + if not l_ezfio_name: + raise ValueError + e_glob = EZFIO_ocaml(l_ezfio_name=l_ezfio_name, l_type=l_type, l_doc=l_doc) @@ -602,6 +605,49 @@ def save_ocaml_input(module_lower, str_ocaml_input): f.write(str_ocaml_input) +def get_l_module_lower(): + """ + Get all module who have EZFIO.cfg with input data + (NB `search` in all the ligne and `match` only in one) + """ + + # ~#~#~#~ # + # I n i t # + # ~#~#~#~ # + + mypath = "{0}/src".format(os.environ['QPACKAGE_ROOT']) + + # ~#~#~#~#~#~#~#~ # + # L _ f o l d e r # + # ~#~#~#~#~#~#~#~ # + + from os import listdir + from os.path import isdir, join, exists + + l_folder = [f for f in listdir(mypath) if isdir(join(mypath, f))] + + # ~#~#~#~#~#~#~#~#~#~#~#~#~#~ # + # L _ m o d u l e _ l o w e r # + # ~#~#~#~#~#~#~#~#~#~#~#~#~#~ # + + l_module_lower = [] + import re + p = re.compile(ur'interface:\s+input') + + for f in l_folder: + path = "{0}/{1}/EZFIO.cfg".format(mypath, f) + if exists(path): + with open(path, 'r') as file_: + if p.search(file_.read()): + l_module_lower.append(f.lower()) + + # ~#~#~#~#~#~ # + # R e t u r n # + # ~#~#~#~#~#~ # + + return l_module_lower + + def create_ocaml_input_global(): """ Check for all the EZFIO.cfg get the module lower @@ -612,15 +658,7 @@ def create_ocaml_input_global(): # I n i t # # ~#~#~#~# # - from os import listdir - from os.path import isdir, join, exists - - mypath = "{0}/src".format(os.environ['QPACKAGE_ROOT']) - - onlyfodler = [f for f in listdir(mypath) if isdir(join(mypath, f))] - - l_module_lower = [f.lower() for f in onlyfodler - if exists("{0}/{1}/EZFIO.cfg".format(mypath, f))] + l_module_lower = get_l_module_lower() # ~#~#~#~#~#~#~#~# # # C r e a t i o n # @@ -761,8 +799,12 @@ if __name__ == "__main__": # O c a m l # # ~#~#~#~#~#~# if do_all or arguments["--ocaml"]: - str_ocaml_input = create_ocaml_input(dict_ezfio_cfg, module_lower) - save_ocaml_input(module_lower, str_ocaml_input) + try: + str_ocaml_input = create_ocaml_input(dict_ezfio_cfg, module_lower) + except ValueError: + pass + else: + save_ocaml_input(module_lower, str_ocaml_input) # ~#~#~#~#~#~#~#~#~#~#~#~#~ # # e z f i o _ d e f a u l t # diff --git a/scripts/ezfio_with_default.py b/scripts/ezfio_with_default.py index 46981132..1b5f01a8 100755 --- a/scripts/ezfio_with_default.py +++ b/scripts/ezfio_with_default.py @@ -90,19 +90,14 @@ END_PROVIDER self.default = t def get_default(self): + filename = '/'.join( [os.environ['QPACKAGE_ROOT'], 'data', + 'ezfio_defaults', + 'WILL_BE_DELETED.ezfio_default'] ) - from os import listdir - from os.path import isfile, join - - mypath = '/'.join( [os.environ['QPACKAGE_ROOT'], 'data', 'ezfio_defaults','/']) - onlyfiles = [ f for f in listdir(mypath) if isfile(join(mypath,f)) ] - - lines= [] - for filename in onlyfiles: - file = open(mypath+filename,'r') - lines.extend(file.readlines()[:]) - file.close() - + file = open(filename,'r') + lines = file.readlines() + file.close() + k=-1 # Search directory for k,line in enumerate(lines): if line[0] != ' ': @@ -120,16 +115,21 @@ END_PROVIDER break v = buffer[1] name = self.name + true = True + false= False try: + true = True + false = False v_eval = eval(v) + except: + v = "call ezfio_get_%(v)s(%(name)s)"%locals() + else: if type(v_eval) == bool: v = '.%s.'%(v) elif type(v_eval) == float: v = v.replace('e','d') v = v.replace('E','D') v = "%(name)s = %(v)s"%locals() - except: - v = "call ezfio_get_%(v)s(%(name)s)"%locals() self.default = v diff --git a/scripts/generate_h_apply.py b/scripts/generate_h_apply.py index ad6a57ee..280c9f72 100755 --- a/scripts/generate_h_apply.py +++ b/scripts/generate_h_apply.py @@ -53,7 +53,7 @@ class H_apply(object): !$OMP N_elec_in_key_hole_2,ia_ja_pairs) & !$OMP SHARED(key_in,N_int,elec_num_tab,mo_tot_num, & !$OMP hole_1, particl_1, hole_2, particl_2, & - !$OMP elec_alpha_num,i_generator)""" + !$OMP elec_alpha_num,i_generator) FIRSTPRIVATE(iproc)""" s["omp_end_parallel"] = "!$OMP END PARALLEL" s["omp_master"] = "!$OMP MASTER" s["omp_end_master"] = "!$OMP END MASTER" diff --git a/scripts/qp_convert.py b/scripts/qp_convert.py index 54a6d2da..ab008e9e 100755 --- a/scripts/qp_convert.py +++ b/scripts/qp_convert.py @@ -53,16 +53,68 @@ def write_ezfioFile(res,filename): basis = res.uncontracted_basis geom = res.geometry + res.clean_contractions() + # AO Basis + import string + at = [] + num_prim = [] + magnetic_number = [] + angular_number = [] + power_x = [] + power_y = [] + power_z = [] + coefficient = [] + exponent = [] + res.convert_to_cartesian() + for b in res.basis: + c = b.center + for i,atom in enumerate(res.geometry): + if atom.coord == c: + at.append(i+1) + num_prim.append(len(b.prim)) + s = b.sym + power_x.append( string.count(s,"x") ) + power_y.append( string.count(s,"y") ) + power_z.append( string.count(s,"z") ) + coefficient.append( b.coef ) + exponent.append( [ p.expo for p in b.prim ] ) + ezfio.set_ao_basis_ao_num(len(res.basis)) + ezfio.set_ao_basis_ao_nucl(at) + ezfio.set_ao_basis_ao_prim_num(num_prim) + ezfio.set_ao_basis_ao_power(power_x+power_y+power_z) + prim_num_max = ezfio.get_ao_basis_ao_prim_num_max() + len_res_basis = len(res.basis) + for i in range(len(res.basis)): + coefficient[i] += [ 0. for j in range(len(coefficient[i]),prim_num_max) ] + exponent[i] += [ 0. for j in range(len(exponent[i]),prim_num_max) ] + coefficient = reduce(lambda x, y: x+y, coefficient, []) + exponent = reduce(lambda x, y: x+y, exponent, []) + coef = [] + expo = [] + for i in range(prim_num_max): + for j in range(i,len(coefficient),prim_num_max): + coef.append ( coefficient[j] ) + expo.append ( exponent[j] ) + ezfio.set_ao_basis_ao_coef(coef) + ezfio.set_ao_basis_ao_expo(expo) + ezfio.set_ao_basis_ao_basis("Read by resultsFile") + + # MO MoTag = res.determinants_mo_type ezfio.set_mo_basis_mo_label('Orthonormalized') MO_type = MoTag - allMOs = res.uncontracted_mo_sets[MO_type] + allMOs = res.mo_sets[MO_type] - closed = [ (allMOs[i].eigenvalue,i) for i in res.closed_mos ] - active = [ (allMOs[i].eigenvalue,i) for i in res.active_mos ] - virtual =[ (allMOs[i].eigenvalue,i) for i in res.virtual_mos ] + try: + closed = [ (allMOs[i].eigenvalue,i) for i in res.closed_mos ] + active = [ (allMOs[i].eigenvalue,i) for i in res.active_mos ] + virtual =[ (allMOs[i].eigenvalue,i) for i in res.virtual_mos ] + except: + closed = [] + virtual = [] + active = [ (allMOs[i].eigenvalue,i) for i in range(len(allMOs)) ] # closed.sort() # active.sort() @@ -111,117 +163,6 @@ def write_ezfioFile(res,filename): while len(MoMatrix) < len(MOs[0].vector)**2: MoMatrix.append(0.) - ezfio.set_mo_basis_mo_tot_num(mo_tot_num) - ezfio.set_mo_basis_mo_occ(OccNum) - - - res.clean_contractions() - # AO Basis - import string - at = [] - num_prim = [] - magnetic_number = [] - angular_number = [] - power_x = [] - power_y = [] - power_z = [] - coefficient = [] - exponent = [] - res.convert_to_cartesian() - for b in res.basis: - c = b.center - for i,atom in enumerate(res.geometry): - if atom.coord == c: - at.append(i+1) - num_prim.append(len(b.prim)) - s = b.sym - power_x.append( string.count(s,"x") ) - power_y.append( string.count(s,"y") ) - power_z.append( string.count(s,"z") ) - coefficient.append( b.coef ) - exponent.append( [ p.expo for p in b.prim ] ) - ezfio.set_ao_basis_ao_num(len(res.basis)) - ezfio.set_ao_basis_ao_nucl(at) - ezfio.set_ao_basis_ao_prim_num(num_prim) - ezfio.set_ao_basis_ao_power(power_x+power_y+power_z) - prim_num_max = ezfio.get_ao_basis_ao_prim_num_max() - len_res_basis = len(res.basis) - for i in range(len(res.basis)): - coefficient[i] += [ 0. for j in range(len(coefficient[i]),prim_num_max) ] - exponent[i] += [ 0. for j in range(len(exponent[i]),prim_num_max) ] - coefficient = reduce(lambda x, y: x+y, coefficient, []) - exponent = reduce(lambda x, y: x+y, exponent, []) - coef = [] - expo = [] - for i in range(prim_num_max): - for j in range(i,len(coefficient),prim_num_max): - coef.append ( coefficient[j] ) - expo.append ( exponent[j] ) - ezfio.set_ao_basis_ao_coef(coef) - ezfio.set_ao_basis_ao_expo(expo) - ezfio.set_ao_basis_ao_basis("Read by resultsFile") - -# Apply threshold to determinants - if len(res.determinants) == 1: - sorted_determinants = [ (-1.,1.,res.determinants[0]) ] - else: - sorted_determinants = [] - for i,j in zip(res.det_coefficients[0],res.determinants): - sorted_determinants.append((-abs(i),i,j)) - sorted_determinants.sort() - norm = 0.0 - for length, (a,b,c) in enumerate(sorted_determinants): - if -a < det_threshold: - length -=1 - break - norm += a**2 - norm = sqrt(norm) - length += 1 - for i in xrange(length): - a = sorted_determinants[i] - sorted_determinants[i] = (a[0],a[1]/norm,a[2]) - sorted_determinants = sorted_determinants[:length] - -# MOs - mo_tot_num = len(res.mo_sets[MoTag]) - closed_mos = res.closed_mos - active_mos = res.active_mos - virtual_mos = res.virtual_mos - to_remove = [] - to_add = [] - for i in active_mos: - found = False - for (a,b,c) in sorted_determinants: - if i in c['alpha']+c['beta']: - found = True - break - if not found: - to_remove.append(i) - to_add.append(i) - virtual_mos = to_add + virtual_mos - for i in active_mos: - always = True - for (a,b,c) in sorted_determinants: - if not (i in c['alpha'] and i in c['beta']): - always = False - break - if always: - to_remove.append(i) - closed_mos.append(i) - for i in to_remove: - active_mos.remove(i) - - - MOindices = closed_mos + active_mos + virtual_mos - while len(MOindices) < mo_tot_num: - MOindices.append(len(MOindices)) - MOmap = list(MOindices) - for i in range(len(MOindices)): - MOmap[i] = MOindices.index(i) - - - ezfio.set_mo_basis_mo_tot_num(mo_tot_num) - mo = [] for i in MOindices: mo.append(res.mo_sets[MoTag][i]) @@ -234,35 +175,11 @@ def write_ezfioFile(res,filename): while len(mo) < mo_tot_num: mo.append(newmo) Energies = [ m.eigenvalue for m in mo ] - - if res.occ_num is not None: - OccNum = [] - for i in MOindices: - OccNum.append(res.occ_num[MoTag][i]) - - while len(OccNum) < mo_tot_num: - OccNum.append(0.) - ezfio.set_mo_basis_mo_occ(OccNum) - - cls = [ "v" for i in mo ] - for i in closed_mos: - cls[MOmap[i]] = 'c' - for i in active_mos: - cls[MOmap[i]] = 'a' - - sym0 = [ i.sym for i in res.mo_sets[MoTag] ] - sym = [ i.sym for i in res.mo_sets[MoTag] ] - for i in xrange(len(sym)): - sym[MOmap[i]] = sym0[i] - - MoMatrix = [] - for m in mo: - for coef in m.vector: - MoMatrix.append(coef) - while len(MoMatrix) < len(mo[0].vector)**2: - MoMatrix.append(0.) + + ezfio.set_mo_basis_mo_tot_num(mo_tot_num) + ezfio.set_mo_basis_mo_occ(OccNum) ezfio.set_mo_basis_mo_coef(MoMatrix) - del MoMatrix + diff --git a/scripts/update_README.py b/scripts/update_README.py index 268c960e..0f32404d 100755 --- a/scripts/update_README.py +++ b/scripts/update_README.py @@ -1,16 +1,16 @@ #!/usr/bin/env python """Updates the README.rst file as the include directive is disabled on GitHub.""" -__date__ = "Thu Apr 3 23:06:18 CEST 2014" +__date__ = "Thu Apr 3 23:06:18 CEST 2014" __author__ = "Anthony Scemama " -README="README.rst" -Assum_key="Assumptions\n===========\n" -Needed_key="Needed Modules\n==============\n" -Doc_key="Documentation\n=============\n" -Sentinel="@@$%&@@" -URL="http://github.com/LCPQ/quantum_package/tree/master/src/" +README = "README.rst" +Assum_key = "Assumptions\n===========\n" +Needed_key = "Needed Modules\n==============\n" +Doc_key = "Documentation\n=============\n" +Sentinel = "@@$%&@@" +URL = "http://github.com/LCPQ/quantum_package/tree/master/src/" import os import subprocess @@ -22,30 +22,30 @@ header = """ """ try: -# subprocess.check_output("git status".split()) -# has_git = True - pass + # subprocess.check_output("git status".split()) + has_git = True except OSError: - has_git = False + has_git = False + def fetch_splitted_data(): """Read the README.rst file and split it in strings: * The description * The assumptions - * The documentation + * The documentation * The needed modules The result is given as a list of strings """ - file = open(README,'r') + file = open(README, 'r') data = file.read() file.close() # Place sentinels - data = data.replace(Assum_key,Sentinel+Assum_key) - data = data.replace(Doc_key,Sentinel+Doc_key) - data = data.replace(Needed_key,Sentinel+Needed_key) - + data = data.replace(Assum_key, Sentinel + Assum_key) + data = data.replace(Doc_key, Sentinel + Doc_key) + data = data.replace(Needed_key, Sentinel + Needed_key) + # Now Split data using the sentinels result = data.split(Sentinel) @@ -56,9 +56,9 @@ def update_assumptions(data): """Read the ASSUMPTIONS.rst file, and replace the data with it.""" try: - file = open('ASSUMPTIONS.rst','r') + file = open('ASSUMPTIONS.rst', 'r') except IOError: - file = open('ASSUMPTIONS.rst','w') + file = open('ASSUMPTIONS.rst', 'w') assumptions = "" else: assumptions = file.read() @@ -74,7 +74,7 @@ def update_assumptions(data): data[i] = assumptions if not has_assumptions: - data.insert(1,assumptions) + data.insert(1, assumptions) return data @@ -83,12 +83,12 @@ 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_MODULES', 'r') modules = file.read() file.close() if modules.strip() != "": - modules = [ '* `%s <%s%s>`_'%(x,URL,x) for x in modules.split() ] + modules = ['* `%s <%s%s>`_' % (x, URL, x) for x in modules.split()] modules = "\n".join(modules) modules = Needed_key + header + modules + '\n\n' @@ -105,112 +105,114 @@ def update_needed(data): def update_documentation(data): - """Reads the BEGIN_DOC ... END_DOC blocks and builds the documentation""" + """Reads the BEGIN_DOC ... END_DOC blocks and builds the documentation""" - # If the file does not exist, don't do anything - try: - file = open('tags','r') - except: - return - tags = file.readlines() - file.close() - - def extract_doc(item): - """Extracts the documentation contained in IRPF90_man file""" - file = open("IRPF90_man/%s.l"%(item),'r') - lines = file.readlines() + # If the file does not exist, don't do anything + try: + file = open('tags', 'r') + except: + return + tags = file.readlines() file.close() - result = [] - inside = False - for line in lines: - if not inside: - inside = line.startswith(".SH Description") - else: - if line.startswith(".SH"): - break - result.append(" "+line.strip()) - if result == []: - result = [" Undocumented"] - return "\n".join(result)+'\n' - - - - items = [] - dirname = os.path.basename(os.getcwd()) - command = "git ls-tree --full-tree --name-only HEAD:src/%s" - command = command%(dirname) - try: - if dirname != 'src': - p = subprocess.Popen(command.split(),stdout=subprocess.PIPE) - tracked_files = p.stdout.read() - else: - tracked_files = "" - tracked_files = tracked_files.splitlines() - except: - tracked_files = [] - for filename in tracked_files: - if filename.endswith('.irp.f'): - # Search for providers, subroutines and functions in each file using - # the tags file - search = "\t"+filename+"\t" - tmp = filter(lambda line: search in line, tags) + def extract_doc(item): + """Extracts the documentation contained in IRPF90_man file""" + file = open("IRPF90_man/%s.l" % (item), 'r') + lines = file.readlines() + file.close() + result = [] + inside = False + for line in lines: + if not inside: + inside = line.startswith(".SH Description") + else: + if line.startswith(".SH"): + break + result.append(" " + line.strip()) - # Search for the documentation in the IRPF90_man directory - for item in tmp : - item, _, line = item.strip().split('\t') - doc = extract_doc(item) - items.append( (item, filename, doc, line) ) + if result == []: + result = [" Undocumented"] + return "\n".join(result) + '\n' - dirname = os.path.basename(os.getcwd()) - # Write the documentation in the README - template = "`%(item)s <%(url)s%(dirname)s/%(filename)s#L%(line)s>`_\n%(doc)s\n" + items = [] + dirname = os.path.basename(os.getcwd()) + command = "git ls-tree --full-tree --name-only HEAD:src/%s" + command = command % (dirname) + try: + if dirname != 'src': + p = subprocess.Popen(command.split(), stdout=subprocess.PIPE) + tracked_files = p.stdout.read() + else: + tracked_files = "" + tracked_files = tracked_files.splitlines() + except: + tracked_files = [] + for filename in tracked_files: + if filename.endswith('.irp.f'): + # Search for providers, subroutines and functions in each file using + # the tags file + search = "\t" + filename + "\t" + tmp = filter(lambda line: search in line, tags) - documentation = Doc_key + header - url = URL - for item, filename, doc, line in items: - documentation += template%locals() - documentation += '\n\n' + # Search for the documentation in the IRPF90_man directory + for item in tmp: + item, _, line = item.strip().split('\t') + doc = extract_doc(item) + items.append((item, filename, doc, line)) - has_doc = False - for i in range(len(data)): - if data[i].startswith(Doc_key): - has_doc = True - data[i] = documentation + dirname = os.path.basename(os.getcwd()) + # Write the documentation in the README + template = "`%(item)s <%(url)s%(dirname)s/%(filename)s#L%(line)s>`_\n%(doc)s\n" - if not has_doc: - data.append(documentation) + documentation = Doc_key + header + url = URL + for item, filename, doc, line in items: + documentation += template % locals() + documentation += '\n\n' - return data + has_doc = False + for i in range(len(data)): + if data[i].startswith(Doc_key): + has_doc = True + data[i] = documentation + + if not has_doc: + data.append(documentation) + + return data def git_add(): """Executes: git add README.rst - if git is present on the machine.""" - command = "git add "+README - os.system(command+" &> /dev/null") + throw an error if git is not precent""" + + import subprocess + + try: + # pipe output to /dev/null for silence + null = open("/dev/null", "w") + subprocess.Popen("git add README.rst", stdout=null, stderr=null) + null.close() + + except OSError: + raise def main(): - if not has_git: - return data = fetch_splitted_data() data = update_assumptions(data) data = update_documentation(data) data = update_needed(data) output = ''.join(data) - - file = open(README,'w') - file.write(output) - file.close() - git_add() + with open(README, 'w') as f: + f.write(output) + try: + git_add() + except OSError: + pass if __name__ == '__main__': - try: - main() - except: - pass - + main() diff --git a/setup_environment.sh b/setup_environment.sh index 59f3af70..12fcf4a5 100755 --- a/setup_environment.sh +++ b/setup_environment.sh @@ -50,16 +50,6 @@ then fi -echo "${BLUE}===== Installing EZFIO ===== ${BLACK}" -${QPACKAGE_ROOT}/scripts/install_ezfio.sh | tee install_ezfio.log - -if [[ ! -d ${QPACKAGE_ROOT}/EZFIO ]] -then - echo $RED "Error in EZFIO installation" $BLACK - exit 1 -fi - - echo "${BLUE}===== Installing Zlib ===== ${BLACK}" ${QPACKAGE_ROOT}/scripts/install_zlib.sh | tee install_zlib.log diff --git a/src/AOs/README.rst b/src/AOs/README.rst index 1a1c1378..f9f81f5f 100644 --- a/src/AOs/README.rst +++ b/src/AOs/README.rst @@ -50,21 +50,65 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`ao_overlap `_ + Overlap between atomic basis functions: + :math:`\int \chi_i(r) \chi_j(r) dr)` + +`ao_overlap_abs `_ + Overlap between absolute value of atomic basis functions: + :math:`\int |\chi_i(r)| |\chi_j(r)| dr)` + +`ao_overlap_x `_ + Overlap between atomic basis functions: + :math:`\int \chi_i(r) \chi_j(r) dr)` + +`ao_overlap_y `_ + Overlap between atomic basis functions: + :math:`\int \chi_i(r) \chi_j(r) dr)` + +`ao_overlap_z `_ + Overlap between atomic basis functions: + :math:`\int \chi_i(r) \chi_j(r) dr)` + `ao_coef `_ Coefficients, exponents and powers of x,y and z - ao_coef(i,j) = coefficient of the jth primitive on the ith ao -`ao_coef_transp `_ +`ao_coef_transp `_ Transposed ao_coef and ao_expo +`ao_coef_unnormalized `_ + 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_expo `_ Coefficients, exponents and powers of x,y and z - ao_coef(i,j) = coefficient of the jth primitive on the ith ao -`ao_expo_transp `_ +`ao_expo_transp `_ Transposed ao_coef and ao_expo -`ao_nucl `_ +`ao_expo_unsorted `_ + 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 `_ + 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 `_ + 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 `_ + Undocumented + +`ao_md5 `_ + MD5 key characteristic of the AO basis + +`ao_nucl `_ Index of the nuclei on which the ao is centered `ao_num `_ @@ -75,16 +119,37 @@ Documentation `ao_power `_ Coefficients, exponents and powers of x,y and z - ao_coef(i,j) = coefficient of the jth primitive on the ith ao -`ao_prim_num `_ +`ao_prim_num `_ Number of primitives per atomic orbital -`ao_prim_num_max `_ +`ao_prim_num_max `_ Undocumented -`ao_prim_num_max_align `_ +`ao_prim_num_max_align `_ Undocumented +`l_to_charater `_ + character corresponding to the "L" value of an AO orbital + +`n_aos_max `_ + Number of AOs per atom + +`nucl_aos `_ + List of AOs attached on each atom + +`nucl_list_shell_aos `_ + 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 `_ + Number of AOs per atom + +`nucl_num_shell_aos `_ + 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 + diff --git a/src/Bielec_integrals/EZFIO.cfg b/src/Bielec_integrals/EZFIO.cfg index bd1d774a..eaada232 100644 --- a/src/Bielec_integrals/EZFIO.cfg +++ b/src/Bielec_integrals/EZFIO.cfg @@ -7,26 +7,26 @@ ezfio_name: direct [disk_access_mo_integrals] type: Disk_access -doc: Write, Read, None for MO integrals from disk (EZFIO folder) +doc: Read/Write MO integrals from/to disk [ Write | Read | None ] interface: input default: None [disk_access_ao_integrals] type: Disk_access -doc: Write or Read or None for AO integrals from disk (EZFIO folder) +doc: Read/Write AO integrals from/to disk [ Write | Read | None ] interface: input default: None [ao_integrals_threshold] type: Threshold -doc: If < ao_integrals_threshold then is null +doc: If || < ao_integrals_threshold then is zero interface: input default: 1.e-15 ezfio_name: threshold_ao [mo_integrals_threshold] type: Threshold -doc: If < ao_integrals_threshold then is null +doc: If || < ao_integrals_threshold then is zero interface: input default: 1.e-15 -ezfio_name: threshold_mo \ No newline at end of file +ezfio_name: threshold_mo diff --git a/src/Bielec_integrals/README.rst b/src/Bielec_integrals/README.rst index 6287ffad..38dc9e96 100644 --- a/src/Bielec_integrals/README.rst +++ b/src/Bielec_integrals/README.rst @@ -32,152 +32,193 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`ao_bielec_integral `_ +`ao_bielec_integral `_ integral of the AO basis or (ij|kl) i(r1) j(r1) 1/r12 k(r2) l(r2) -`ao_bielec_integral_schwartz `_ - Needed to compuet Schwartz inequalities +`ao_bielec_integral_schwartz `_ + Needed to compute Schwartz inequalities -`ao_bielec_integrals_in_map `_ +`ao_bielec_integral_schwartz_accel `_ + integral of the AO basis or (ij|kl) + i(r1) j(r1) 1/r12 k(r2) l(r2) + +`ao_bielec_integrals_in_map `_ Map of Atomic integrals i(r1) j(r2) 1/r12 k(r1) l(r2) -`compute_ao_bielec_integrals `_ +`ao_l4 `_ + Computes the product of l values of i,j,k,and l + +`compute_ao_bielec_integrals `_ Compute AO 1/r12 integrals for all i and fixed j,k,l -`eri `_ +`eri `_ 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 `_ +`general_primitive_integral `_ Computes the integral where p,q,r,s are Gaussian primitives -`give_polynom_mult_center_x `_ +`give_polynom_mult_center_x `_ 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 `_ +`i_x1_new `_ recursive function involved in the bielectronic integral -`i_x1_pol_mult `_ +`i_x1_pol_mult `_ recursive function involved in the bielectronic integral -`i_x1_pol_mult_a1 `_ +`i_x1_pol_mult_a1 `_ recursive function involved in the bielectronic integral -`i_x1_pol_mult_a2 `_ +`i_x1_pol_mult_a2 `_ recursive function involved in the bielectronic integral -`i_x1_pol_mult_recurs `_ +`i_x1_pol_mult_recurs `_ recursive function involved in the bielectronic integral -`i_x2_new `_ +`i_x2_new `_ recursive function involved in the bielectronic integral -`i_x2_pol_mult `_ +`i_x2_pol_mult `_ recursive function involved in the bielectronic integral -`integrale_new `_ +`integrale_new `_ 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 `_ - Returns the upper boundary of the degree of the polynom involved in the +`n_pt_sup `_ + 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) -`gauleg `_ +`gauleg `_ Gauss-Legendre -`gauleg_t2 `_ +`gauleg_t2 `_ t_w(i,1,k) = w(i) t_w(i,2,k) = t(i) -`gauleg_w `_ +`gauleg_w `_ t_w(i,1,k) = w(i) t_w(i,2,k) = t(i) -`ao_integrals_map `_ +`n_pt_max_integrals_16 `_ + Aligned n_pt_max_integrals + +`ao_integrals_map `_ AO integrals -`bielec_integrals_index `_ +`bielec_integrals_index `_ Undocumented -`clear_ao_map `_ +`bielec_integrals_index_reverse `_ + Undocumented + +`clear_ao_map `_ Frees the memory of the AO map -`clear_mo_map `_ +`clear_mo_map `_ Frees the memory of the MO map -`get_ao_bielec_integral `_ +`get_ao_bielec_integral `_ Gets one AO bi-electronic integral from the AO map -`get_ao_bielec_integrals `_ +`get_ao_bielec_integrals `_ Gets multiple AO bi-electronic integral from the AO map . All i are retrieved for j,k,l fixed. -`get_ao_bielec_integrals_non_zero `_ +`get_ao_bielec_integrals_non_zero `_ Gets multiple AO bi-electronic integral from the AO map . All non-zero i are retrieved for j,k,l fixed. -`get_ao_map_size `_ +`get_ao_map_size `_ Returns the number of elements in the AO map -`get_mo_bielec_integral `_ +`get_mo_bielec_integral `_ Returns one integral in the MO basis -`get_mo_bielec_integrals `_ +`get_mo_bielec_integrals `_ Returns multiple integrals in the MO basis, all i for j,k,l fixed. -`get_mo_bielec_integrals_existing_ik `_ +`get_mo_bielec_integrals_existing_ik `_ Returns multiple integrals in the MO basis, all + i(1)j(1) 1/r12 k(2)l(2) i for j,k,l fixed. -`get_mo_map_size `_ +`get_mo_map_size `_ Return the number of elements in the MO map -`insert_into_ao_integrals_map `_ +`insert_into_ao_integrals_map `_ Create new entry into AO map -`insert_into_mo_integrals_map `_ +`insert_into_mo_integrals_map `_ Create new entry into MO map, or accumulate in an existing entry -`mo_bielec_integral `_ +`mo_bielec_integral `_ Returns one integral in the MO basis -`mo_integrals_map `_ +`mo_integrals_map `_ MO integrals -`add_integrals_to_map `_ +`add_integrals_to_map `_ Adds integrals to tha MO map according to some bitmask -`mo_bielec_integral_jj `_ +`mo_bielec_integral_jj `_ mo_bielec_integral_jj(i,j) = J_ij - mo_bielec_integral_jj_exchange(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 `_ +`mo_bielec_integral_jj_anti `_ mo_bielec_integral_jj(i,j) = J_ij - mo_bielec_integral_jj_exchange(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 `_ +`mo_bielec_integral_jj_anti_from_ao `_ + 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 `_ mo_bielec_integral_jj(i,j) = J_ij - mo_bielec_integral_jj_exchange(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_integrals_in_map `_ +`mo_bielec_integral_jj_exchange_from_ao `_ + 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 `_ + 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_integrals_in_map `_ If True, the map of MO bielectronic integrals is provided -`mo_bielec_integrals_index `_ +`mo_bielec_integrals_index `_ Computes an unique index for i,j,k,l integrals +`read_ao_integrals `_ + One level of abstraction for disk_access_ao_integrals and disk_access_mo_integrals + +`read_mo_integrals `_ + One level of abstraction for disk_access_ao_integrals and disk_access_mo_integrals + +`write_ao_integrals `_ + One level of abstraction for disk_access_ao_integrals and disk_access_mo_integrals + +`write_mo_integrals `_ + One level of abstraction for disk_access_ao_integrals and disk_access_mo_integrals + diff --git a/src/Bielec_integrals/map_integrals.irp.f b/src/Bielec_integrals/map_integrals.irp.f index 1426288a..8fdff799 100644 --- a/src/Bielec_integrals/map_integrals.irp.f +++ b/src/Bielec_integrals/map_integrals.irp.f @@ -368,11 +368,11 @@ subroutine get_mo_bielec_integrals_existing_ik(j,l,sze,out_array,map) enddo logical :: integral_is_in_map - if (cache_key_kind == 8) then + if (key_kind == 8) then call i8radix_sort(hash,iorder,kk,-1) - else if (cache_key_kind == 4) then + else if (key_kind == 4) then call iradix_sort(hash,iorder,kk,-1) - else if (cache_key_kind == 2) then + else if (key_kind == 2) then call i2radix_sort(hash,iorder,kk,-1) endif diff --git a/src/Bitmask/README.rst b/src/Bitmask/README.rst index e7e45997..b8e3aa57 100644 --- a/src/Bitmask/README.rst +++ b/src/Bitmask/README.rst @@ -54,6 +54,24 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`is_a_two_holes_two_particles `_ + Undocumented + +`number_of_holes `_ + Undocumented + +`number_of_holes_verbose `_ + Undocumented + +`number_of_particles `_ + Undocumented + +`number_of_particles_verbose `_ + Undocumented + +`cas_bitmask `_ + Bitmasks for CAS reference determinants. (N_int, alpha/beta, CAS reference) + `cis_ijkl_bitmask `_ Bitmask to include all possible single excitations from Hartree-Fock @@ -61,23 +79,35 @@ Documentation Bitmask to include all possible MOs `generators_bitmask `_ - Bitmasks for generator determinants. (N_int, alpha/beta, hole/particle, generator). + Bitmasks for generator determinants. + (N_int, alpha/beta, hole/particle, generator). + .br 3rd index is : + .br * 1 : hole for single exc - * 1 : particle for single exc + .br + * 2 : particle for single exc + .br * 3 : hole for 1st exc of double + .br * 4 : particle for 1st exc of double - * 5 : hole for 2dn exc of double - * 6 : particle for 2dn exc of double + .br + * 5 : hole for 2nd exc of double + .br + * 6 : particle for 2nd exc of double + .br `hf_bitmask `_ Hartree Fock bit mask -`i_bitmask_gen `_ +`i_bitmask_gen `_ Current bitmask for the generators -`i_bitmask_ref `_ - Current bitmask for the reference +`inact_bitmask `_ + Bitmasks for the inactive orbitals that are excited in post CAS method + +`n_cas_bitmask `_ + Number of bitmasks for CAS `n_generators_bitmask `_ Number of bitmasks for generators @@ -85,36 +115,31 @@ Documentation `n_int `_ Number of 64-bit integers needed to represent determinants as binary strings -`n_reference_bitmask `_ - Number of bitmasks for reference - `ref_bitmask `_ Reference bit mask, used in Slater rules, chosen as Hartree-Fock bitmask -`reference_bitmask `_ - Bitmasks for reference determinants. (N_int, alpha/beta, hole/particle, reference) +`virt_bitmask `_ + Bitmasks for the inactive orbitals that are excited in post CAS method -`bitstring_to_hexa `_ +`bitstring_to_hexa `_ Transform a bit string to a string in hexadecimal format for printing `bitstring_to_list `_ Gives the inidices(+1) of the bits set to 1 in the bit string -`bitstring_to_str `_ +`bitstring_to_str `_ Transform a bit string to a string for printing -`debug_det `_ - Undocumented +`debug_det `_ + Subroutine to print the content of a determinant in '+-' notation and + hexadecimal representation. `list_to_bitstring `_ - return the physical string "string(N_int,2)" from the array of occupations "list(N_int*bit_kind_size,2) - list - <== ipos ==> - | - v - string :|------------------------|-------------------------|------------------------| - <==== bit_kind_size ====> <==== bit_kind_size ====> <==== bit_kind_size ====> - { iint } { iint } { iint } + Returns the physical string "string(N_int,2)" from the array of + occupations "list(N_int*bit_kind_size,2) + +`print_det `_ + Subroutine to print the content of a determinant using the '+-' notation diff --git a/src/Bitmask/bitmasks.irp.f b/src/Bitmask/bitmasks.irp.f index e0e04614..2d044ca5 100644 --- a/src/Bitmask/bitmasks.irp.f +++ b/src/Bitmask/bitmasks.irp.f @@ -123,7 +123,6 @@ BEGIN_PROVIDER [ integer(bit_kind), generators_bitmask, (N_int,2,6,N_generators_ call ezfio_has_bitmasks_generators(exists) if (exists) then - print*,'EXIST !!' call ezfio_get_bitmasks_generators(generators_bitmask) else integer :: k, ispin @@ -181,7 +180,6 @@ BEGIN_PROVIDER [ integer(bit_kind), cas_bitmask, (N_int,2,N_cas_bitmask) ] PROVIDE ezfio_filename call ezfio_has_bitmasks_cas(exists) - print*,'exists = ',exists if (exists) then call ezfio_get_bitmasks_cas(cas_bitmask) else diff --git a/src/CAS_SD_selected/ASSUMPTIONS.rst b/src/CAS_SD/ASSUMPTIONS.rst similarity index 100% rename from src/CAS_SD_selected/ASSUMPTIONS.rst rename to src/CAS_SD/ASSUMPTIONS.rst diff --git a/src/CAS_SD/EZFIO.cfg b/src/CAS_SD/EZFIO.cfg new file mode 100644 index 00000000..7f27f95a --- /dev/null +++ b/src/CAS_SD/EZFIO.cfg @@ -0,0 +1,36 @@ +[N_det_max_cas_sd] +type: Det_number_max +doc: Max number of determinants in the wave function +interface: input +default: 10000 + +[do_pt2_end] +type: logical +doc: If true, compute the PT2 at the end of the selection +interface: input +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: input +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: input +default: 0.75 + +[energy] +type: double precision +doc: "Calculated CAS-SD energy" +interface: output + +[energy_pt2] +type: double precision +doc: "Calculated selected CAS-SD energy with PT2 correction" +interface: output + diff --git a/src/CAS_SD_selected/H_apply.irp.f b/src/CAS_SD/H_apply.irp.f similarity index 68% rename from src/CAS_SD_selected/H_apply.irp.f rename to src/CAS_SD/H_apply.irp.f index 76c046ae..85b2f1c6 100644 --- a/src/CAS_SD_selected/H_apply.irp.f +++ b/src/CAS_SD/H_apply.irp.f @@ -2,11 +2,14 @@ use bitmasks BEGIN_SHELL [ /usr/bin/env python ] from generate_h_apply import * -s = H_apply("FCI") +s = H_apply("CAS_SD") +print s + +s = H_apply("CAS_SD_selected") s.set_selection_pt2("epstein_nesbet_2x2") print s -s = H_apply("FCI_PT2") +s = H_apply("CAS_SD_PT2") s.set_perturbation("epstein_nesbet_2x2") print s diff --git a/src/CAS_SD_selected/Makefile b/src/CAS_SD/Makefile similarity index 100% rename from src/CAS_SD_selected/Makefile rename to src/CAS_SD/Makefile diff --git a/src/CAS_SD_selected/NEEDED_MODULES b/src/CAS_SD/NEEDED_MODULES similarity index 100% rename from src/CAS_SD_selected/NEEDED_MODULES rename to src/CAS_SD/NEEDED_MODULES diff --git a/src/CAS_SD/README.rst b/src/CAS_SD/README.rst new file mode 100644 index 00000000..0dc4ea56 --- /dev/null +++ b/src/CAS_SD/README.rst @@ -0,0 +1,44 @@ +====================== +CAS_SD_selected Module +====================== + +Selected CAS + SD module. + +1) Set the different MO classes using the ``qp_set_mo_class`` command +2) Run the selected CAS+SD program + +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +`full_ci `_ + Undocumented + + + +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Bielec_integrals `_ +* `Bitmask `_ +* `Dets `_ +* `Electrons `_ +* `Ezfio_files `_ +* `Generators_CAS `_ +* `Hartree_Fock `_ +* `MOGuess `_ +* `MonoInts `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Perturbation `_ +* `Properties `_ +* `Selectors_full `_ +* `Utils `_ + diff --git a/src/CAS_SD_selected/cas_sd.irp.f b/src/CAS_SD/cas_sd.irp.f similarity index 87% rename from src/CAS_SD_selected/cas_sd.irp.f rename to src/CAS_SD/cas_sd.irp.f index e5e0497a..144eec88 100644 --- a/src/CAS_SD_selected/cas_sd.irp.f +++ b/src/CAS_SD/cas_sd.irp.f @@ -11,12 +11,12 @@ program full_ci pt2 = 1.d0 diag_algorithm = "Lapack" - if (N_det > n_det_max_fci) then + if (N_det > n_det_max_cas_sd) then call diagonalize_CI call save_wavefunction psi_det = psi_det_sorted psi_coef = psi_coef_sorted - N_det = n_det_max_fci + N_det = n_det_max_cas_sd soft_touch N_det psi_det psi_coef call diagonalize_CI call save_wavefunction @@ -28,17 +28,17 @@ program full_ci print *, '-----' endif - do while (N_det < n_det_max_fci.and.maxval(abs(pt2(1:N_st))) > pt2_max) - call H_apply_FCI(pt2, norm_pert, H_pert_diag, N_st) + do while (N_det < n_det_max_cas_sd.and.maxval(abs(pt2(1:N_st))) > pt2_max) + call H_apply_CAS_SD(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 + if (N_det > n_det_max_cas_sd) then psi_det = psi_det_sorted psi_coef = psi_coef_sorted - N_det = n_det_max_fci + N_det = n_det_max_cas_sd soft_touch N_det psi_det psi_coef endif call diagonalize_CI @@ -60,7 +60,7 @@ program full_ci integer :: exc_max, degree_min exc_max = 0 print *, 'CAS determinants : ', N_det_generators - do i=1,N_det_generators + do i=1,min(N_det_generators,10) do k=i,N_det_generators call get_excitation_degree(psi_det_generators(1,1,k),psi_det_generators(1,1,i),degree,N_int) exc_max = max(exc_max,degree) diff --git a/src/CAS_SD/cas_sd_selected.irp.f b/src/CAS_SD/cas_sd_selected.irp.f new file mode 100644 index 00000000..86cac969 --- /dev/null +++ b/src/CAS_SD/cas_sd_selected.irp.f @@ -0,0 +1,87 @@ +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 + PROVIDE N_det_cas + + pt2 = 1.d0 + diag_algorithm = "Lapack" + if (N_det > n_det_max_cas_sd) then + call diagonalize_CI + call save_wavefunction + psi_det = psi_det_sorted + psi_coef = psi_coef_sorted + N_det = n_det_max_cas_sd + 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 + + do while (N_det < n_det_max_cas_sd.and.maxval(abs(pt2(1:N_st))) > pt2_max) + call H_apply_CAS_SD_selected(pt2, norm_pert, H_pert_diag, N_st) + + PROVIDE psi_coef + PROVIDE psi_det + PROVIDE psi_det_sorted + + if (N_det > n_det_max_cas_sd) then + psi_det = psi_det_sorted + psi_coef = psi_coef_sorted + N_det = n_det_max_cas_sd + soft_touch N_det psi_det psi_coef + endif + 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 *, '-----' + call ezfio_set_full_ci_energy(CI_energy) + if (abort_all) 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_cas + do i=1,min(N_det_cas,10) + do k=i,N_det_cas + call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int) + exc_max = max(exc_max,degree) + enddo + call debug_det(psi_cas(1,1,i),N_int) + 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 diff --git a/src/CAS_SD_selected/README.rst b/src/CAS_SD_selected/README.rst deleted file mode 100644 index f66fe229..00000000 --- a/src/CAS_SD_selected/README.rst +++ /dev/null @@ -1,5 +0,0 @@ -====================== -CAS_SD_selected Module -====================== - -Selected CAS + SD module diff --git a/src/CAS_SD_selected/options.irp.f b/src/CAS_SD_selected/options.irp.f deleted file mode 100644 index 553fbe9f..00000000 --- a/src/CAS_SD_selected/options.irp.f +++ /dev/null @@ -1,32 +0,0 @@ -BEGIN_SHELL [ /usr/bin/python ] -from ezfio_with_default import EZFIO_Provider -T = EZFIO_Provider() -T.set_type ( "integer" ) -T.set_name ( "N_det_max_fci" ) -T.set_doc ( "Max number of determinants in the wave function" ) -T.set_ezfio_dir ( "full_ci" ) -T.set_ezfio_name( "N_det_max_fci" ) -T.set_output ( "output_full_ci" ) -print T - -T.set_type ( "logical" ) -T.set_name ( "do_pt2_end" ) -T.set_doc ( "If true, compute the PT2 at the end of the selection" ) -T.set_ezfio_name( "do_pt2_end" ) -print T - -T.set_type ( "double precision" ) -T.set_name ( "pt2_max" ) -T.set_doc ( """The selection process stops when the largest PT2 (for all the states) -is lower than pt2_max in absolute value""" ) -T.set_ezfio_name( "pt2_max" ) -print T - -T.set_type ( "double precision" ) -T.set_name ( "var_pt2_ratio" ) -T.set_doc ( """The selection process stops when the energy ratio variational/(variational+PT2) -is equal to var_pt2_ratio""" ) -T.set_ezfio_name( "var_pt2_ratio" ) -print T -END_SHELL - diff --git a/src/CID/README.rst b/src/CID/README.rst index 5a0c5026..6adc4dcd 100644 --- a/src/CID/README.rst +++ b/src/CID/README.rst @@ -16,19 +16,20 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ +* `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ -* `Selectors_full `_ Documentation ============= diff --git a/src/CID_SC2_selected/README.rst b/src/CID_SC2_selected/README.rst index 8688fd41..720b6385 100644 --- a/src/CID_SC2_selected/README.rst +++ b/src/CID_SC2_selected/README.rst @@ -20,20 +20,21 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `CISD `_ -* `SC2 `_ * `CISD_selected `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ * `Perturbation `_ +* `Properties `_ * `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ diff --git a/src/CID_selected/README.rst b/src/CID_selected/README.rst index 8c64bc25..2ee45ac6 100644 --- a/src/CID_selected/README.rst +++ b/src/CID_selected/README.rst @@ -23,19 +23,21 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `CISD `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ * `Perturbation `_ +* `Properties `_ +* `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ -* `Selectors_full `_ diff --git a/src/CIS/README.rst b/src/CIS/README.rst index bb22445a..e9ba93db 100644 --- a/src/CIS/README.rst +++ b/src/CIS/README.rst @@ -32,17 +32,18 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ +* `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ -* `Selectors_full `_ diff --git a/src/CISD/README.rst b/src/CISD/README.rst index 003b51fa..07528e59 100644 --- a/src/CISD/README.rst +++ b/src/CISD/README.rst @@ -16,19 +16,20 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ +* `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ -* `Selectors_full `_ Documentation ============= diff --git a/src/CISD/cisd.irp.f b/src/CISD/cisd.irp.f index 0b6bc9fd..6d310d95 100644 --- a/src/CISD/cisd.irp.f +++ b/src/CISD/cisd.irp.f @@ -13,6 +13,7 @@ program cisd print *, 'E_corr = ',CI_electronic_energy(i) - ref_bitmask_energy enddo + call save_wavefunction ! call CISD_SC2(psi_det,psi_coef,eigvalues,size(psi_coef,1),N_det,N_states,N_int) ! do i = 1, N_states ! print*,'eigvalues(i) = ',eigvalues(i) diff --git a/src/CISD_SC2_selected/EZFIO.cfg b/src/CISD_SC2_selected/EZFIO.cfg index 660f9fd6..b76eac5d 100644 --- a/src/CISD_SC2_selected/EZFIO.cfg +++ b/src/CISD_SC2_selected/EZFIO.cfg @@ -1,6 +1,6 @@ [N_det_max_cisd_sc2] type: Det_number_max -doc: Get n_det_max_cisd_sc2 from EZFIO file +doc: Max number of determinants interface: input default: 10000 @@ -8,12 +8,12 @@ default: 10000 type: logical doc: If true, compute the PT2 at the end of the selection interface: input -default: true +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 +doc: The selection process stops when the largest PT2 (for all the states) is lower + than abs(pt2_max) interface: input default: 0.0001 diff --git a/src/CISD_SC2_selected/README.rst b/src/CISD_SC2_selected/README.rst index 4c67caaa..25f4368f 100644 --- a/src/CISD_SC2_selected/README.rst +++ b/src/CISD_SC2_selected/README.rst @@ -8,6 +8,9 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`cisd_sc2_selected `_ + Undocumented + Needed Modules @@ -17,20 +20,21 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `CISD `_ -* `SC2 `_ * `CISD_selected `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ * `Perturbation `_ +* `Properties `_ * `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ diff --git a/src/CISD_selected/README.rst b/src/CISD_selected/README.rst index 4fed26e5..1ba5f9c5 100644 --- a/src/CISD_selected/README.rst +++ b/src/CISD_selected/README.rst @@ -14,6 +14,12 @@ Documentation `cisd `_ Undocumented +`n_det_max_cisd `_ + Get n_det_max_cisd from EZFIO file + +`pt2_max `_ + Get pt2_max from EZFIO file + Needed Modules @@ -23,19 +29,21 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `CISD `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ * `Perturbation `_ +* `Properties `_ +* `Selectors_full `_ * `SingleRefMethod `_ * `Utils `_ -* `Selectors_full `_ diff --git a/src/DDCI_selected/README.rst b/src/DDCI_selected/README.rst index 33d21997..db75b101 100644 --- a/src/DDCI_selected/README.rst +++ b/src/DDCI_selected/README.rst @@ -2,3 +2,38 @@ DDCI_selected Module ==================== +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +`full_ci `_ + Undocumented + + + +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Bielec_integrals `_ +* `Bitmask `_ +* `Dets `_ +* `Electrons `_ +* `Ezfio_files `_ +* `Generators_CAS `_ +* `Hartree_Fock `_ +* `MOGuess `_ +* `MonoInts `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Perturbation `_ +* `Properties `_ +* `Selectors_full `_ +* `Utils `_ + diff --git a/src/Dets/H_apply_template.f b/src/Dets/H_apply_template.f index e44562e7..a9a282ae 100644 --- a/src/Dets/H_apply_template.f +++ b/src/Dets/H_apply_template.f @@ -1,4 +1,4 @@ -subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_generator, iproc $parameters ) +subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_generator, iproc_in $parameters ) use omp_lib use bitmasks implicit none @@ -14,7 +14,7 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene integer(bit_kind),allocatable :: keys_out(:,:,:) integer(bit_kind), intent(in) :: hole_1(N_int,2), particl_1(N_int,2) integer(bit_kind), intent(in) :: hole_2(N_int,2), particl_2(N_int,2) - integer, intent(in) :: iproc + integer, intent(in) :: iproc_in integer(bit_kind), allocatable :: hole_save(:,:) integer(bit_kind), allocatable :: key(:,:),hole(:,:), particle(:,:) integer(bit_kind), allocatable :: hole_tmp(:,:), particle_tmp(:,:) @@ -30,6 +30,7 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene integer, allocatable :: ia_ja_pairs(:,:,:) integer, allocatable :: ib_jb_pairs(:,:) double precision :: diag_H_mat_elem + integer :: iproc integer(omp_lock_kind), save :: lck, ifirst=0 if (ifirst == 0) then !$ call omp_init_lock(lck) @@ -38,12 +39,13 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene logical :: check_double_excitation check_double_excitation = .True. - + iproc = iproc_in $initialization $omp_parallel +!$ iproc = omp_get_thread_num() allocate (keys_out(N_int,2,size_max), hole_save(N_int,2), & key(N_int,2),hole(N_int,2), particle(N_int,2), hole_tmp(N_int,2),& particle_tmp(N_int,2), occ_particle(N_int*bit_kind_size,2), & @@ -248,7 +250,7 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene $finalization end -subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $parameters ) +subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc_in $parameters ) use omp_lib use bitmasks implicit none @@ -262,7 +264,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param integer ,intent(in) :: i_generator integer(bit_kind),intent(in) :: key_in(N_int,2) integer(bit_kind),intent(in) :: hole_1(N_int,2), particl_1(N_int,2) - integer, intent(in) :: iproc + integer, intent(in) :: iproc_in integer(bit_kind),allocatable :: keys_out(:,:,:) integer(bit_kind),allocatable :: hole_save(:,:) integer(bit_kind),allocatable :: key(:,:),hole(:,:), particle(:,:) @@ -281,8 +283,11 @@ 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 + integer :: iproc logical :: check_double_excitation + iproc = iproc_in + check_double_excitation = .True. $check_double_excitation @@ -295,6 +300,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param $initialization $omp_parallel +!$ iproc = omp_get_thread_num() allocate (keys_out(N_int,2,size_max), hole_save(N_int,2), & key(N_int,2),hole(N_int,2), particle(N_int,2), hole_tmp(N_int,2),& particle_tmp(N_int,2), occ_particle(N_int*bit_kind_size,2), & @@ -396,7 +402,8 @@ subroutine $subroutine($params_main) integer :: iproc $initialization - PROVIDE H_apply_buffer_allocated mo_bielec_integrals_in_map + PROVIDE H_apply_buffer_allocated mo_bielec_integrals_in_map psi_det_generators psi_coef_generators + nmax = mod( N_det_generators,nproc ) @@ -406,6 +413,7 @@ subroutine $subroutine($params_main) call wall_time(wall_0) + iproc = 0 allocate( mask(N_int,2,6) ) do i_generator=1,nmax @@ -443,12 +451,12 @@ subroutine $subroutine($params_main) call $subroutine_diexc(psi_det_generators(1,1,i_generator), & mask(1,1,d_hole1), mask(1,1,d_part1), & mask(1,1,d_hole2), mask(1,1,d_part2), & - i_generator, 0 $params_post) + i_generator, iproc $params_post) endif if($do_mono_excitations)then call $subroutine_monoexc(psi_det_generators(1,1,i_generator), & mask(1,1,s_hole ), mask(1,1,s_part ), & - i_generator, 0 $params_post) + i_generator, iproc $params_post) endif call wall_time(wall_1) $printout_always @@ -463,7 +471,6 @@ subroutine $subroutine($params_main) !$OMP PARALLEL DEFAULT(SHARED) & !$OMP PRIVATE(i_generator,wall_1,wall_0,ispin,k,mask,iproc) call wall_time(wall_0) - iproc = 0 !$ iproc = omp_get_thread_num() allocate( mask(N_int,2,6) ) !$OMP DO SCHEDULE(dynamic,1) diff --git a/src/Dets/README.rst b/src/Dets/README.rst index de6b5069..f03df8da 100644 --- a/src/Dets/README.rst +++ b/src/Dets/README.rst @@ -33,11 +33,10 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Electrons `_ * `Ezfio_files `_ -* `Hartree_Fock `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ @@ -50,18 +49,23 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`copy_h_apply_buffer_to_wf `_ +`copy_h_apply_buffer_to_wf `_ Copies the H_apply buffer to psi_coef. You need to touch psi_det, psi_coef and N_det 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 `_ +`fill_h_apply_buffer_no_selection `_ Fill the H_apply buffer with determiants for CISD -`h_apply_buffer_allocated `_ +`h_apply_buffer_allocated `_ Buffer of determinants/coefficients/perturbative energy for H_apply. Uninitialized. Filled by H_apply subroutines. -`resize_h_apply_buffer `_ +`h_apply_buffer_lock `_ + Buffer of determinants/coefficients/perturbative energy for H_apply. + Uninitialized. Filled by H_apply subroutines. + +`resize_h_apply_buffer `_ Undocumented `cisd_sc2 `_ @@ -80,26 +84,43 @@ Documentation .br Initial guess vectors are not necessarily orthonormal -`repeat_excitation `_ +`connected_to_ref `_ Undocumented -`connected_to_ref `_ +`connected_to_ref_by_mono `_ Undocumented -`det_is_not_or_may_be_in_ref `_ +`det_is_not_or_may_be_in_ref `_ If true, det is not in ref If false, det may be in ref -`is_in_wavefunction `_ - Undocumented +`det_search_key `_ + Return an integer*8 corresponding to a determinant index for searching -`key_pattern_not_in_ref `_ +`get_index_in_psi_det_sorted_bit `_ + Returns the index of the determinant in the ``psi_det_sorted_bit`` array + +`is_in_wavefunction `_ + True if the determinant ``det`` is in the wave function + +`key_pattern_not_in_ref `_ Min and max values of the integers of the keys of the reference -`davidson_converged `_ +`occ_pattern_search_key `_ + Return an integer*8 corresponding to a determinant index for searching + +`do_mono_excitation `_ + Apply the mono excitation operator : a^{dager}_(i_particle) a_(i_hole) of spin = ispin + on key_in + ispin = 1 == alpha + ispin = 2 == beta + i_ok = 1 == the excitation is possible + i_ok = -1 == the excitation is not possible + +`davidson_converged `_ True if the Davidson algorithm is converged -`davidson_criterion `_ +`davidson_criterion `_ Can be : [ energy | residual | both | wall_time | cpu_time | iterations ] `davidson_diag `_ @@ -146,10 +167,10 @@ Documentation `davidson_sze_max `_ Max number of Davidson sizes -`davidson_threshold `_ +`davidson_threshold `_ Can be : [ energy | residual | both | wall_time | cpu_time | iterations ] -`one_body_dm_mo `_ +`one_body_dm_mo `_ One-body density matrix `one_body_dm_mo_alpha `_ @@ -158,63 +179,162 @@ Documentation `one_body_dm_mo_beta `_ Alpha and beta one-body density matrix for each state -`save_natural_mos `_ +`one_body_single_double_dm_mo_alpha `_ + Alpha and beta one-body density matrix for each state + +`one_body_single_double_dm_mo_beta `_ + Alpha and beta one-body density matrix for each state + +`one_body_spin_density_mo `_ + rho(alpha) - rho(beta) + +`save_natural_mos `_ Save natural orbitals, obtained by diagonalization of the one-body density matrix in the MO basis -`set_natural_mos `_ +`set_natural_mos `_ Set natural orbitals, obtained by diagonalization of the one-body density matrix in the MO basis -`state_average_weight `_ +`state_average_weight `_ Weights in the state-average calculation of the density matrix -`det_search_key `_ - Return an integer*8 corresponding to a determinant index for searching +`det_svd `_ + Computes the SVD of the Alpha x Beta determinant coefficient matrix + +`create_wf_of_psi_svd_matrix `_ + Matrix of wf coefficients. Outer product of alpha and beta determinants + +`filter_3_highest_electrons `_ + Returns a determinant with only the 3 highest electrons + +`generate_all_alpha_beta_det_products `_ + Create a wave function from all possible alpha x beta determinants + +`int_of_3_highest_electrons `_ + Returns an integer*8 as : + .br + |_<--- 21 bits ---><--- 21 bits ---><--- 21 bits --->| + .br + |0<--- i1 ---><--- i2 ---><--- i3 --->| + .br + It encodes the value of the indices of the 3 highest MOs + in descending order + .br + +`max_degree_exc `_ + Maximum degree of excitation in the wf `n_det `_ Number of determinants in the wave function -`psi_average_norm_contrib `_ +`n_det_alpha_unique `_ + Unique alpha determinants + +`n_det_beta_unique `_ + Unique beta determinants + +`psi_average_norm_contrib `_ Contribution of determinants to the state-averaged density -`psi_average_norm_contrib_sorted `_ +`psi_average_norm_contrib_sorted `_ Wave function sorted by determinants contribution to the norm (state-averaged) -`psi_coef `_ +`psi_coef `_ The wave function coefficients. Initialized with Hartree-Fock if the EZFIO file is empty -`psi_coef_sorted `_ +`psi_coef_sorted `_ Wave function sorted by determinants contribution to the norm (state-averaged) -`psi_coef_sorted_bit `_ - Determinants on which we apply for perturbation. - o They are sorted by determinants interpreted as integers. Useful - to accelerate the search of a determinant +`psi_coef_sorted_ab `_ + Determinants on which we apply . + 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 `_ +`psi_coef_sorted_bit `_ + Determinants on which we apply 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 `_ The wave function determinants. Initialized with Hartree-Fock if the EZFIO file is empty -`psi_det_size `_ +`psi_det_alpha `_ + List of alpha determinants of psi_det + +`psi_det_alpha_unique `_ + Unique alpha determinants + +`psi_det_beta `_ + List of beta determinants of psi_det + +`psi_det_beta_unique `_ + Unique beta determinants + +`psi_det_size `_ Size of the psi_det/psi_coef arrays -`psi_det_sorted `_ +`psi_det_sorted `_ Wave function sorted by determinants contribution to the norm (state-averaged) -`psi_det_sorted_bit `_ - Determinants on which we apply for perturbation. - o They are sorted by determinants interpreted as integers. Useful - to accelerate the search of a determinant +`psi_det_sorted_ab `_ + Determinants on which we apply . + 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 `_ +`psi_det_sorted_bit `_ + Determinants on which we apply 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 `_ + Determinants on which we apply . + 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_svd_alpha `_ + SVD wave function + +`psi_svd_beta `_ + SVD wave function + +`psi_svd_coefs `_ + SVD wave function + +`psi_svd_matrix `_ + Matrix of wf coefficients. Outer product of alpha and beta determinants + +`read_dets `_ Reads the determinants from the EZFIO file -`save_wavefunction `_ +`save_wavefunction `_ Save the wave function into the EZFIO file -`save_wavefunction_general `_ +`save_wavefunction_general `_ Save the wave function into the EZFIO file +`save_wavefunction_unsorted `_ + Save the wave function into the EZFIO file + +`sort_dets_by_3_highest_electrons `_ + Determinants on which we apply . + 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 `_ + Determinants are sorted are sorted according to their det_search_key. + Useful to accelerate the search of a random determinant in the wave + function. + +`spin_det_search_key `_ + Return an integer*8 corresponding to a determinant index for searching + `double_exc_bitmask `_ double_exc_bitmask(:,1,i) is the bitmask for holes of excitation 1 double_exc_bitmask(:,2,i) is the bitmask for particles of excitation 1 @@ -233,10 +353,13 @@ Documentation single_exc_bitmask(:,2,i) is the bitmask for particles for a given couple of hole/particle excitations i. -`ci_eigenvectors `_ +`ci_eigenvectors `_ Eigenvectors/values of the CI matrix -`ci_electronic_energy `_ +`ci_eigenvectors_s2 `_ + Eigenvectors/values of the CI matrix + +`ci_electronic_energy `_ Eigenvectors/values of the CI matrix `ci_energy `_ @@ -245,7 +368,7 @@ Documentation `diag_algorithm `_ Diagonalization algorithm (Davidson or Lapack) -`diagonalize_ci `_ +`diagonalize_ci `_ Replace the coefficients of the CI states by the coefficients of the eigenstates of the CI matrix @@ -256,15 +379,31 @@ Documentation Eigenvectors/values of the CI matrix `ci_sc2_energy `_ - N_states lowest eigenvalues of the CI matrix + N_states_diag lowest eigenvalues of the CI matrix `diagonalize_ci_sc2 `_ - Replace the coefficients of the CI states by the coefficients of the + Replace the coefficients of the CI states_diag by the coefficients of the eigenstates of the CI matrix `threshold_convergence_sc2 `_ convergence of the correlation energy of SC2 iterations +`ci_eigenvectors_mono `_ + Eigenvectors/values of the CI matrix + +`ci_eigenvectors_s2_mono `_ + Eigenvectors/values of the CI matrix + +`ci_electronic_energy_mono `_ + Eigenvectors/values of the CI matrix + +`diagonalize_ci_mono `_ + Replace the coefficients of the CI states by the coefficients of the + eigenstates of the CI matrix + +`apply_mono `_ + Undocumented + `filter_connected `_ Filters out the determinants that are not connected by H .br @@ -276,9 +415,16 @@ Documentation .br idx(0) is the number of determinants that interact with key1 -`filter_connected_davidson `_ +`filter_connected_davidson `_ Filters out the determinants that are not connected by H + returns the array idx which contains the index of the + determinants in the array key1 that interact + via the H operator with key2. .br + idx(0) is the number of determinants that interact with key1 + key1 should come from psi_det_sorted_ab. + +`filter_connected_i_h_psi0 `_ returns the array idx which contains the index of the .br determinants in the array key1 that interact @@ -287,16 +433,7 @@ Documentation .br idx(0) is the number of determinants that interact with key1 -`filter_connected_i_h_psi0 `_ - returns the array idx which contains the index of the - .br - determinants in the array key1 that interact - .br - via the H operator with key2. - .br - idx(0) is the number of determinants that interact with key1 - -`filter_connected_i_h_psi0_sc2 `_ +`filter_connected_i_h_psi0_sc2 `_ standard filter_connected_i_H_psi but returns in addition .br the array of the index of the non connected determinants to key1 @@ -307,31 +444,156 @@ Documentation .br to repeat the excitations +`filter_connected_sorted_ab `_ + Filters out the determinants that are not connected by H + returns the array idx which contains the index of the + determinants in the array key1 that interact + via the H operator with key2. + idx(0) is the number of determinants that interact with key1 + .br + Determinants are taken from the psi_det_sorted_ab array + +`put_gess `_ + Undocumented + +`det_to_occ_pattern `_ + Transform a determinant to an occupation pattern + +`make_s2_eigenfunction `_ + Undocumented + +`n_occ_pattern `_ + array of the occ_pattern present in the wf + psi_occ_pattern(:,1,j) = jth occ_pattern of the wave function : represent all the single occupation + psi_occ_pattern(:,2,j) = jth occ_pattern of the wave function : represent all the double occupation + +`occ_pattern_to_dets `_ + Generate all possible determinants for a give occ_pattern + +`occ_pattern_to_dets_size `_ + Number of possible determinants for a given occ_pattern + +`psi_occ_pattern `_ + array of the occ_pattern present in the wf + psi_occ_pattern(:,1,j) = jth occ_pattern of the wave function : represent all the single occupation + psi_occ_pattern(:,2,j) = jth occ_pattern of the wave function : represent all the double occupation + +`rec_occ_pattern_to_dets `_ + Undocumented + +`n_states_diag `_ + Number of states to consider for the diagonalization + +`pouet `_ + Undocumented + +`routine `_ + Undocumented + +`idx_cas `_ + CAS wave function, defined from the application of the CAS bitmask on the + determinants. idx_cas gives the indice of the CAS determinant in psi_det. + +`idx_non_cas `_ + Set of determinants which are not part of the CAS, defined from the application + of the CAS bitmask on the determinants. + idx_non_cas gives the indice of the determinant in psi_det. + +`n_det_cas `_ + CAS wave function, defined from the application of the CAS bitmask on the + determinants. idx_cas gives the indice of the CAS determinant in psi_det. + +`n_det_non_cas `_ + Set of determinants which are not part of the CAS, defined from the application + of the CAS bitmask on the determinants. + idx_non_cas gives the indice of the determinant in psi_det. + +`psi_cas `_ + CAS wave function, defined from the application of the CAS bitmask on the + determinants. idx_cas gives the indice of the CAS determinant in psi_det. + +`psi_cas_coef `_ + CAS wave function, defined from the application of the CAS bitmask on the + determinants. idx_cas gives the indice of the CAS determinant in psi_det. + +`psi_cas_coef_sorted_bit `_ + CAS determinants sorted to accelerate the search of a random determinant in the wave + function. + +`psi_cas_sorted_bit `_ + CAS determinants sorted to accelerate the search of a random determinant in the wave + function. + +`psi_non_cas `_ + Set of determinants which are not part of the CAS, defined from the application + of the CAS bitmask on the determinants. + idx_non_cas gives the indice of the determinant in psi_det. + +`psi_non_cas_coef `_ + Set of determinants which are not part of the CAS, defined from the application + of the CAS bitmask on the determinants. + idx_non_cas gives the indice of the determinant in psi_det. + +`psi_non_cas_coef_sorted_bit `_ + CAS determinants sorted to accelerate the search of a random determinant in the wave + function. + +`psi_non_cas_sorted_bit `_ + CAS determinants sorted to accelerate the search of a random determinant in the wave + function. + +`bi_elec_ref_bitmask_energy `_ + Energy of the reference bitmask used in Slater rules + +`kinetic_ref_bitmask_energy `_ + Energy of the reference bitmask used in Slater rules + +`mono_elec_ref_bitmask_energy `_ + Energy of the reference bitmask used in Slater rules + +`nucl_elec_ref_bitmask_energy `_ + Energy of the reference bitmask used in Slater rules + +`ref_bitmask_energy `_ + Energy of the reference bitmask used in Slater rules + +`expected_s2 `_ + Expected value of S2 : S*(S+1) + `get_s2 `_ Returns -`get_s2_u0 `_ +`get_s2_u0 `_ Undocumented +`s2_values `_ + array of the averaged values of the S^2 operator on the various states + `s_z `_ - Undocumented + z component of the Spin `s_z2_sz `_ + z component of the Spin + +`prog_save_casino `_ + Undocumented + +`save_casino `_ Undocumented `save_dets_qmcchem `_ Undocumented -`save_for_qmc `_ +`save_for_qmc `_ Undocumented `save_natorb `_ Undocumented -`a_operator `_ +`a_operator `_ Needed for diag_H_mat_elem -`ac_operator `_ +`ac_operator `_ Needed for diag_H_mat_elem `decode_exc `_ @@ -341,10 +603,10 @@ Documentation s1,s2 : Spins (1:alpha, 2:beta) degree : Degree of excitation -`det_connections `_ - .br +`det_connections `_ + Build connection proxy between determinants -`diag_h_mat_elem `_ +`diag_h_mat_elem `_ Computes `get_double_excitation `_ @@ -356,16 +618,16 @@ Documentation `get_excitation_degree `_ Returns the excitation degree between two determinants -`get_excitation_degree_vector `_ +`get_excitation_degree_vector `_ Applies get_excitation_degree to an array of determinants `get_mono_excitation `_ Returns the excitation operator between two singly excited determinants and the phase -`get_occ_from_key `_ +`get_occ_from_key `_ Returns a list of occupation numbers from a bitstring -`h_u_0 `_ +`h_u_0 `_ Computes v_0 = H|u_0> .br n : number of determinants @@ -375,10 +637,13 @@ Documentation `i_h_j `_ Returns where i and j are determinants -`i_h_psi `_ +`i_h_j_verbose `_ + Returns where i and j are determinants + +`i_h_psi `_ for the various Nstates -`i_h_psi_sc2 `_ +`i_h_psi_sc2 `_ for the various Nstate .br returns in addition @@ -391,11 +656,33 @@ Documentation .br to repeat the excitations -`n_con_int `_ +`i_h_psi_sc2_verbose `_ + for the various Nstate + .br + returns in addition + .br + the array of the index of the non connected determinants to key1 + .br + in order to know what double excitation can be repeated on key1 + .br + idx_repeat(0) is the number of determinants that can be used + .br + to repeat the excitations + +`i_h_psi_sec_ord `_ + for the various Nstates + +`n_con_int `_ Number of integers to represent the connections between determinants +`write_spindeterminants `_ + Undocumented + +`cisd `_ + Undocumented + `h_matrix_all_dets `_ - H matrix on the basis of the slater deter;inants defined by psi_det + H matrix on the basis of the slater determinants defined by psi_det diff --git a/src/Dets/connected_to_ref.irp.f b/src/Dets/connected_to_ref.irp.f index 07331782..3c7eb581 100644 --- a/src/Dets/connected_to_ref.irp.f +++ b/src/Dets/connected_to_ref.irp.f @@ -95,9 +95,9 @@ integer function get_index_in_psi_det_sorted_bit(key,Nint) enddo i += 1 -! if (i > N_det) then -! return -! endif + if (i > N_det) then + return + endif !DIR$ FORCEINLINE do while (det_search_key(psi_det_sorted_bit(1,1,i),Nint) == det_ref) @@ -116,39 +116,39 @@ integer function get_index_in_psi_det_sorted_bit(key,Nint) enddo if (is_in_wavefunction) then get_index_in_psi_det_sorted_bit = i - exit -! return +! exit + return endif endif i += 1 if (i > N_det) then - exit -! return +! exit + return endif enddo ! DEBUG is_in_wf - if (is_in_wavefunction) then - degree = 1 - do i=1,N_det - integer :: degree - call get_excitation_degree(key,psi_det(1,1,i),degree,N_int) - if (degree == 0) then - exit - endif - enddo - if (degree /=0) then - stop 'pouet 1' - endif - else - do i=1,N_det - call get_excitation_degree(key,psi_det(1,1,i),degree,N_int) - if (degree == 0) then - stop 'pouet 2' - endif - enddo - endif +! if (is_in_wavefunction) then +! degree = 1 +! do i=1,N_det +! integer :: degree +! call get_excitation_degree(key,psi_det(1,1,i),degree,N_int) +! if (degree == 0) then +! exit +! endif +! enddo +! if (degree /=0) then +! stop 'pouet 1' +! endif +! else +! do i=1,N_det +! call get_excitation_degree(key,psi_det(1,1,i),degree,N_int) +! if (degree == 0) then +! stop 'pouet 2' +! endif +! enddo +! endif ! END DEBUG is_in_wf end diff --git a/src/Dets/det_svd.irp.f b/src/Dets/det_svd.irp.f index 8365b42c..0a57acf3 100644 --- a/src/Dets/det_svd.irp.f +++ b/src/Dets/det_svd.irp.f @@ -3,7 +3,7 @@ program det_svd BEGIN_DOC ! Computes the SVD of the Alpha x Beta determinant coefficient matrix END_DOC - integer :: i,j + integer :: i,j,k read_wf = .True. TOUCH read_wf @@ -40,17 +40,22 @@ program det_svd print *, 'N_det_alpha = ', N_det_alpha_unique print *, 'N_det_beta = ', N_det_beta_unique print *, '' - -! do i=1,N_det_alpha_unique -! do j=1,N_det_beta_unique -! print *, i,j,psi_svd_matrix(i,j,:) -! enddo -! enddo - print *, '' + call diagonalize_ci print *, 'Energy = ', ci_energy + + do i=1,N_det_alpha_unique + do j=1,N_det_beta_unique + do k=1,N_states + if (dabs(psi_svd_matrix(i,j,k)) < 1.d-15) then + psi_svd_matrix(i,j,k) = 0.d0 + endif + enddo + enddo + enddo + print *, '' print *, psi_svd_coefs(1:20,1) -! call save_wavefunction + call save_wavefunction end diff --git a/src/Dets/determinants.irp.f b/src/Dets/determinants.irp.f index 0a0ec864..00e683fc 100644 --- a/src/Dets/determinants.irp.f +++ b/src/Dets/determinants.irp.f @@ -467,33 +467,52 @@ END_PROVIDER ! to accelerate the search of a random determinant in the wave ! function. END_DOC + + call sort_dets_by_det_search_key(N_det, psi_det, psi_coef, & + psi_det_sorted_bit, psi_coef_sorted_bit) +END_PROVIDER + +subroutine sort_dets_by_det_search_key(Ndet, det_in, coef_in, det_out, coef_out) + use bitmasks + implicit none + integer, intent(in) :: Ndet + integer(bit_kind), intent(in) :: det_in (N_int,2,psi_det_size) + double precision , intent(in) :: coef_in(psi_det_size,N_states) + integer(bit_kind), intent(out) :: det_out (N_int,2,psi_det_size) + double precision , intent(out) :: coef_out(psi_det_size,N_states) + BEGIN_DOC + ! Determinants are sorted are sorted according to their det_search_key. + ! Useful to accelerate the search of a random determinant in the wave + ! function. + END_DOC integer :: i,j,k integer, allocatable :: iorder(:) integer*8, allocatable :: bit_tmp(:) integer*8, external :: det_search_key - allocate ( iorder(N_det), bit_tmp(N_det) ) + allocate ( iorder(Ndet), bit_tmp(Ndet) ) - do i=1,N_det + do i=1,Ndet iorder(i) = i !$DIR FORCEINLINE - bit_tmp(i) = det_search_key(psi_det(1,1,i),N_int) + bit_tmp(i) = det_search_key(det_in(1,1,i),N_int) enddo - call i8sort(bit_tmp,iorder,N_det) + call i8sort(bit_tmp,iorder,Ndet) !DIR$ IVDEP - do i=1,N_det + do i=1,Ndet do j=1,N_int - psi_det_sorted_bit(j,1,i) = psi_det(j,1,iorder(i)) - psi_det_sorted_bit(j,2,i) = psi_det(j,2,iorder(i)) + det_out(j,1,i) = det_in(j,1,iorder(i)) + det_out(j,2,i) = det_in(j,2,iorder(i)) enddo do k=1,N_states - psi_coef_sorted_bit(i,k) = psi_coef(iorder(i),k) + coef_out(i,k) = coef_in(iorder(i),k) enddo enddo deallocate(iorder, bit_tmp) -END_PROVIDER +end + subroutine int_of_3_highest_electrons( det_in, res, Nint ) implicit none diff --git a/src/Dets/psi_cas.irp.f b/src/Dets/psi_cas.irp.f new file mode 100644 index 00000000..299e5e8f --- /dev/null +++ b/src/Dets/psi_cas.irp.f @@ -0,0 +1,114 @@ +use bitmasks + + BEGIN_PROVIDER [ integer(bit_kind), psi_cas, (N_int,2,psi_det_size) ] +&BEGIN_PROVIDER [ double precision, psi_cas_coef, (psi_det_size,n_states) ] +&BEGIN_PROVIDER [ integer, idx_cas, (psi_det_size) ] +&BEGIN_PROVIDER [ integer, N_det_cas ] + implicit none + BEGIN_DOC + ! CAS wave function, defined from the application of the CAS bitmask on the + ! determinants. idx_cas gives the indice of the CAS determinant in psi_det. + END_DOC + integer :: i, k, l + logical :: good + N_det_cas = 0 + do i=1,N_det + do l=1,n_cas_bitmask + good = .True. + do k=1,N_int + good = good .and. ( & + iand(not(cas_bitmask(k,1,l)), psi_det(k,1,i)) == & + iand(not(cas_bitmask(k,1,l)), psi_det(k,1,1)) ) .and. ( & + iand(not(cas_bitmask(k,2,l)), psi_det(k,2,i)) == & + iand(not(cas_bitmask(k,2,l)), psi_det(k,2,1)) ) + enddo + if (good) then + exit + endif + enddo + if (good) then + N_det_cas = N_det_cas+1 + do k=1,N_int + psi_cas(k,1,N_det_cas) = psi_det(k,1,i) + psi_cas(k,2,N_det_cas) = psi_det(k,2,i) + enddo + idx_cas(N_det_cas) = i + do k=1,N_states + psi_cas_coef(N_det_cas,k) = psi_coef(i,k) + enddo + endif + enddo + call write_int(output_dets,N_det_cas, 'Number of determinants in the CAS') + +END_PROVIDER + + + BEGIN_PROVIDER [ integer(bit_kind), psi_cas_sorted_bit, (N_int,2,psi_det_size) ] +&BEGIN_PROVIDER [ double precision, psi_cas_coef_sorted_bit, (psi_det_size,N_states) ] + implicit none + BEGIN_DOC + ! CAS determinants sorted to accelerate the search of a random determinant in the wave + ! function. + END_DOC + call sort_dets_by_det_search_key(N_det_cas, psi_cas, psi_cas_coef, & + psi_cas_sorted_bit, psi_cas_coef_sorted_bit) + +END_PROVIDER + + + + BEGIN_PROVIDER [ integer(bit_kind), psi_non_cas, (N_int,2,psi_det_size) ] +&BEGIN_PROVIDER [ double precision, psi_non_cas_coef, (psi_det_size,n_states) ] +&BEGIN_PROVIDER [ integer, idx_non_cas, (psi_det_size) ] +&BEGIN_PROVIDER [ integer, N_det_non_cas ] + implicit none + BEGIN_DOC + ! Set of determinants which are not part of the CAS, defined from the application + ! of the CAS bitmask on the determinants. + ! idx_non_cas gives the indice of the determinant in psi_det. + END_DOC + integer :: i_non_cas,j,k + integer :: degree + logical :: in_cas + i_non_cas =0 + do k=1,N_det + in_cas = .False. + do j=1,N_det_cas + call get_excitation_degree(psi_cas(1,1,j), psi_det(1,1,k), degree, N_int) + if (degree == 0) then + in_cas = .True. + exit + endif + enddo + if (.not.in_cas) then + double precision :: hij + i_non_cas += 1 + do j=1,N_int + psi_non_cas(j,1,i_non_cas) = psi_det(j,1,k) + psi_non_cas(j,2,i_non_cas) = psi_det(j,2,k) + enddo + do j=1,N_states + psi_non_cas_coef(i_non_cas,j) = psi_coef(k,j) + enddo + idx_non_cas(i_non_cas) = k + endif + enddo + N_det_non_cas = i_non_cas +END_PROVIDER + + BEGIN_PROVIDER [ integer(bit_kind), psi_non_cas_sorted_bit, (N_int,2,psi_det_size) ] +&BEGIN_PROVIDER [ double precision, psi_non_cas_coef_sorted_bit, (psi_det_size,N_states) ] + implicit none + BEGIN_DOC + ! CAS determinants sorted to accelerate the search of a random determinant in the wave + ! function. + END_DOC + call sort_dets_by_det_search_key(N_det_cas, psi_non_cas, psi_non_cas_coef, & + psi_non_cas_sorted_bit, psi_non_cas_coef_sorted_bit) + +END_PROVIDER + + + + + diff --git a/src/Dets/save_for_qmcchem.irp.f b/src/Dets/save_for_qmcchem.irp.f index 363e0ba2..7dea70c6 100644 --- a/src/Dets/save_for_qmcchem.irp.f +++ b/src/Dets/save_for_qmcchem.irp.f @@ -42,8 +42,10 @@ subroutine save_dets_qmcchem enddo close(31) call system('gzip -f '//trim(ezfio_filename)//'/mo_basis/mo_classif') + end program save_for_qmc call save_dets_qmcchem + call write_spindeterminants end diff --git a/src/Dets/spindeterminants.ezfio_config b/src/Dets/spindeterminants.ezfio_config index e6a0ed9a..1c7d81e3 100644 --- a/src/Dets/spindeterminants.ezfio_config +++ b/src/Dets/spindeterminants.ezfio_config @@ -7,4 +7,8 @@ spindeterminants psi_det_alpha integer*8 (spindeterminants_n_int*spindeterminants_bit_kind/8,spindeterminants_n_det_alpha) psi_det_beta integer*8 (spindeterminants_n_int*spindeterminants_bit_kind/8,spindeterminants_n_det_beta) psi_coef_matrix double precision (spindeterminants_n_det_alpha,spindeterminants_n_det_beta,spindeterminants_n_states) + n_svd_coefs integer + psi_svd_alpha double precision (spindeterminants_n_det_alpha,spindeterminants_n_svd_coefs,spindeterminants_n_states) + psi_svd_beta double precision (spindeterminants_n_det_beta,spindeterminants_n_svd_coefs,spindeterminants_n_states) + psi_svd_coefs double precision (spindeterminants_n_svd_coefs,spindeterminants_n_states) diff --git a/src/Dets/spindeterminants.irp.f b/src/Dets/spindeterminants.irp.f index 236dc1fc..4b426faa 100644 --- a/src/Dets/spindeterminants.irp.f +++ b/src/Dets/spindeterminants.irp.f @@ -26,7 +26,7 @@ subroutine write_spindeterminants enddo call ezfio_set_spindeterminants_psi_det_alpha(psi_det_alpha_unique) deallocate(tmpdet) - + allocate(tmpdet(N_int2,N_det_beta_unique)) do i=1,N_det_beta_unique do k=1,N_int @@ -38,7 +38,54 @@ subroutine write_spindeterminants enddo call ezfio_set_spindeterminants_psi_det_beta(psi_det_beta_unique) deallocate(tmpdet) - + call ezfio_set_spindeterminants_psi_coef_matrix(psi_svd_matrix) + + integer :: n_svd_coefs + double precision :: norm, f + f = 1.d0/dble(N_states) + norm = 1.d0 + do n_svd_coefs=1,N_det_alpha_unique + do k=1,N_states + norm -= psi_svd_coefs(n_svd_coefs,k)*psi_svd_coefs(n_svd_coefs,k) + enddo + if (norm < 1.d-6) then + exit + endif + enddo + n_svd_coefs -= 1 + call ezfio_set_spindeterminants_n_svd_coefs(n_svd_coefs) + + double precision, allocatable :: dtmp(:,:,:) + allocate(dtmp(N_det_alpha_unique,n_svd_coefs,N_states)) + do k=1,N_states + do j=1,n_svd_coefs + do i=1,N_det_alpha_unique + dtmp(i,j,k) = psi_svd_alpha(i,j,k) + enddo + enddo + enddo + call ezfio_set_spindeterminants_psi_svd_alpha(dtmp) + deallocate(dtmp) + + allocate(dtmp(N_det_beta_unique,n_svd_coefs,N_states)) + do k=1,N_states + do j=1,n_svd_coefs + do i=1,N_det_beta_unique + dtmp(i,j,k) = psi_svd_beta(i,j,k) + enddo + enddo + enddo + call ezfio_set_spindeterminants_psi_svd_beta(dtmp) + deallocate(dtmp) + + allocate(dtmp(n_svd_coefs,N_states,1)) + do k=1,N_states + do j=1,n_svd_coefs + dtmp(j,k,1) = psi_svd_coefs(j,k) + enddo + enddo + call ezfio_set_spindeterminants_psi_svd_coefs(dtmp) + deallocate(dtmp) end diff --git a/src/FCIdump/README.rst b/src/FCIdump/README.rst index c53bee01..1fdd9660 100644 --- a/src/FCIdump/README.rst +++ b/src/FCIdump/README.rst @@ -4,3 +4,32 @@ FCIdump Module Interface for the `NECI `_ Full-CI QMC program. +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +`fcidump `_ + Undocumented + + + +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Bielec_integrals `_ +* `Bitmask `_ +* `Dets `_ +* `Electrons `_ +* `Ezfio_files `_ +* `MonoInts `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Utils `_ + diff --git a/src/Full_CI/EZFIO.cfg b/src/Full_CI/EZFIO.cfg index 30528fd8..89679d1a 100644 --- a/src/Full_CI/EZFIO.cfg +++ b/src/Full_CI/EZFIO.cfg @@ -14,7 +14,7 @@ default: 10000 type: logical doc: If true, compute the PT2 at the end of the selection interface: input -default: true +default: True [PT2_max] type: PT2_energy @@ -32,11 +32,11 @@ default: 0.75 [energy] type: double precision -doc: "Calculated Full CI energy" +doc: Calculated Selected FCI energy interface: output [energy_pt2] type: double precision -doc: "Calculated Full CI energy" +doc: Calculated FCI energy + PT2 interface: output diff --git a/src/Full_CI/README.rst b/src/Full_CI/README.rst index 38b30dea..0b37ca69 100644 --- a/src/Full_CI/README.rst +++ b/src/Full_CI/README.rst @@ -10,7 +10,10 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`full_ci `_ +`full_ci `_ + Undocumented + +`var_pt2_ratio_run `_ Undocumented @@ -22,18 +25,20 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Generators_full `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ +* `Perturbation `_ +* `Properties `_ * `Selectors_full `_ * `Utils `_ -* `Perturbation `_ diff --git a/src/Generators_CAS/README.rst b/src/Generators_CAS/README.rst index 06db0139..53e8c5a0 100644 --- a/src/Generators_CAS/README.rst +++ b/src/Generators_CAS/README.rst @@ -3,3 +3,55 @@ Generators_CAS Module ===================== The generator determinants are those filtered by the ``cas_bitmask`` mask. +Assumptions +=========== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +The active space is defined by the ``reference_bitmask``. + + +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +`n_det_generators `_ + Number of generator detetrminants + +`psi_coef_generators `_ + For Single reference wave functions, the generator is the + Hartree-Fock determinant + +`psi_det_generators `_ + For Single reference wave functions, the generator is the + Hartree-Fock determinant + +`select_max `_ + Memo to skip useless selectors + +`size_select_max `_ + Size of the select_max array + + + +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Bielec_integrals `_ +* `Bitmask `_ +* `Dets `_ +* `Electrons `_ +* `Ezfio_files `_ +* `MonoInts `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Utils `_ + diff --git a/src/Generators_CAS/generators.irp.f b/src/Generators_CAS/generators.irp.f index 5dcfc1b0..511877a0 100644 --- a/src/Generators_CAS/generators.irp.f +++ b/src/Generators_CAS/generators.irp.f @@ -62,7 +62,6 @@ END_PROVIDER psi_det_generators(k,2,m) = psi_det(k,2,i) enddo psi_coef_generators(m,:) = psi_coef(m,:) -! call debug_det(psi_det_generators(1,1,m),N_int) endif enddo diff --git a/src/Generators_full/README.rst b/src/Generators_full/README.rst index 3e7f8134..a8492dbc 100644 --- a/src/Generators_full/README.rst +++ b/src/Generators_full/README.rst @@ -11,20 +11,26 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`n_det_generators `_ +`degree_max_generators `_ + Max degree of excitation (respect to HF) of the generators + +`n_det_generators `_ For Single reference wave functions, the number of generators is 1 : the Hartree-Fock determinant -`psi_generators `_ +`psi_coef_generators `_ For Single reference wave functions, the generator is the Hartree-Fock determinant -`select_max `_ +`psi_det_generators `_ + For Single reference wave functions, the generator is the + Hartree-Fock determinant + +`select_max `_ Memo to skip useless selectors -`threshold_generators `_ - Percentage of the norm of the state-averaged wave function to - consider for the generators +`size_select_max `_ + Size of the select_max array @@ -35,12 +41,13 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ diff --git a/src/Hartree_Fock/README.rst b/src/Hartree_Fock/README.rst index 8060f092..a1030e4a 100644 --- a/src/Hartree_Fock/README.rst +++ b/src/Hartree_Fock/README.rst @@ -11,11 +11,12 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Electrons `_ * `Ezfio_files `_ * `MonoInts `_ +* `MOGuess `_ * `MOs `_ * `Nuclei `_ * `Output `_ @@ -36,16 +37,16 @@ Documentation `fock_matrix_alpha_ao `_ Alpha Fock matrix in AO basis set -`fock_matrix_alpha_mo `_ +`fock_matrix_alpha_mo `_ Fock matrix on the MO basis -`fock_matrix_ao `_ +`fock_matrix_ao `_ Fock matrix in AO basis set `fock_matrix_beta_ao `_ Alpha Fock matrix in AO basis set -`fock_matrix_beta_mo `_ +`fock_matrix_beta_mo `_ Fock matrix on the MO basis `fock_matrix_diag_mo `_ @@ -78,24 +79,30 @@ Documentation K = Fb - Fa .br -`fock_mo_to_ao `_ +`fock_mo_to_ao `_ Undocumented -`hf_energy `_ +`hf_energy `_ Hartree-Fock energy `hf_density_matrix_ao `_ - Density matrix in the AO basis + S^-1 Density matrix in the AO basis S^-1 `hf_density_matrix_ao_alpha `_ - Alpha density matrix in the AO basis + S^-1 x Alpha density matrix in the AO basis x S^-1 `hf_density_matrix_ao_beta `_ - Beta density matrix in the AO basis + S^-1 Beta density matrix in the AO basis x S^-1 -`run `_ +`guess `_ Undocumented +`create_guess `_ + Create an MO guess if no MOs are present in the EZFIO directory + +`run `_ + Run SCF calculation + `scf `_ Undocumented @@ -105,7 +112,7 @@ Documentation `diagonal_fock_matrix_mo `_ Diagonal Fock matrix in the MO basis -`diagonal_fock_matrix_mo_sum `_ +`diagonal_fock_matrix_mo_sum `_ diagonal element of the fock matrix calculated as the sum over all the interactions with all the electrons in the RHF determinant diagonal_Fock_matrix_mo_sum(i) = sum_{j=1, N_elec} 2 J_ij -K_ij @@ -113,20 +120,8 @@ Documentation `eigenvectors_fock_matrix_mo `_ Diagonal Fock matrix in the MO basis -`bi_elec_ref_bitmask_energy `_ - Energy of the reference bitmask used in Slater rules - -`kinetic_ref_bitmask_energy `_ - Energy of the reference bitmask used in Slater rules - -`mono_elec_ref_bitmask_energy `_ - Energy of the reference bitmask used in Slater rules - -`nucl_elec_ref_bitmask_energy `_ - Energy of the reference bitmask used in Slater rules - -`ref_bitmask_energy `_ - Energy of the reference bitmask used in Slater rules +`huckel_guess `_ + Build the MOs using the extended Huckel model diff --git a/src/MOs/README.rst b/src/MOs/README.rst index 9e3e7e92..d7a15219 100644 --- a/src/MOs/README.rst +++ b/src/MOs/README.rst @@ -23,11 +23,11 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ +* `Electrons `_ * `Ezfio_files `_ * `Nuclei `_ * `Output `_ * `Utils `_ -* `Electrons `_ Documentation ============= @@ -45,29 +45,58 @@ Documentation `mo_density_matrix_virtual `_ Density matrix in MO basis (virtual MOs) -`mo_coef `_ +`mo_overlap `_ + Undocumented + +`ao_to_mo `_ + Transform A from the AO basis to the MO basis + +`mix_mo_jk `_ + subroutine that rotates the jth MO with the kth MO + to give two new MO's that are + '+' = 1/sqrt(2) (|j> + |k>) + '-' = 1/sqrt(2) (|j> - |k>) + by convention, the '+' MO is in the lower index (min(j,k)) + by convention, the '-' MO is in the greater index (max(j,k)) + +`mo_coef `_ Molecular orbital coefficients on AO basis set mo_coef(i,j) = coefficient of the ith ao on the jth mo mo_label : Label characterizing the MOS (local, canonical, natural, etc) -`mo_coef_transp `_ +`mo_coef_transp `_ Molecular orbital coefficients on AO basis set -`mo_label `_ +`mo_label `_ Molecular orbital coefficients on AO basis set mo_coef(i,j) = coefficient of the ith ao on the jth mo mo_label : Label characterizing the MOS (local, canonical, natural, etc) -`mo_occ `_ +`mo_occ `_ MO occupation numbers +`mo_to_ao `_ + Transform A from the MO basis to the AO basis + +`mo_to_ao_no_overlap `_ + Transform A from the MO basis to the S^-1 AO basis + `mo_tot_num `_ Total number of molecular orbitals and the size of the keys corresponding -`mo_tot_num_align `_ +`mo_tot_num_align `_ Aligned variable for dimensioning of arrays -`mo_as_eigvectors_of_mo_matrix `_ +`s_mo_coef `_ + Product S.C where S is the overlap matrix in the AO basis and C the mo_coef matrix. + +`mo_as_eigvectors_of_mo_matrix `_ + Undocumented + +`mo_as_eigvectors_of_mo_matrix_sort_by_observable `_ + Undocumented + +`mo_sort_by_observable `_ Undocumented `save_mos `_ diff --git a/src/MP2/README.rst b/src/MP2/README.rst index 7a909344..92d915b6 100644 --- a/src/MP2/README.rst +++ b/src/MP2/README.rst @@ -8,7 +8,7 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`cisd `_ +`mp2 `_ Undocumented @@ -20,18 +20,20 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ * `Output `_ * `Perturbation `_ -* `Utils `_ +* `Properties `_ * `Selectors_full `_ * `SingleRefMethod `_ +* `Utils `_ diff --git a/src/MRCC/EZFIO.cfg b/src/MRCC/EZFIO.cfg new file mode 100644 index 00000000..ff586985 --- /dev/null +++ b/src/MRCC/EZFIO.cfg @@ -0,0 +1,4 @@ +[energy] +type: double precision +doc: Calculated MRCC energy +interface: output \ No newline at end of file diff --git a/src/MRCC/H_apply.irp.f b/src/MRCC/H_apply.irp.f index 4f284112..cadaca32 100644 --- a/src/MRCC/H_apply.irp.f +++ b/src/MRCC/H_apply.irp.f @@ -2,13 +2,13 @@ use bitmasks BEGIN_SHELL [ /usr/bin/env python ] from generate_h_apply import * -s = H_apply("mrcc") +s = H_apply("mrcc_simple") s.data["parameters"] = ", delta_ij_sd_, Ndet_sd" s.data["declarations"] += """ integer, intent(in) :: Ndet_sd double precision, intent(in) :: delta_ij_sd_(Ndet_sd,Ndet_sd,*) """ -s.data["keys_work"] = "call mrcc_dress(delta_ij_sd_,Ndet_sd,i_generator,key_idx,keys_out,N_int,iproc)" +s.data["keys_work"] = "call mrcc_dress_simple(delta_ij_sd_,Ndet_sd,i_generator,key_idx,keys_out,N_int,iproc)" s.data["params_post"] += ", delta_ij_sd_, Ndet_sd" s.data["params_main"] += "delta_ij_sd_, Ndet_sd" s.data["decls_main"] += """ @@ -18,8 +18,14 @@ s.data["decls_main"] += """ s.data["finalization"] = "" s.data["copy_buffer"] = "" s.data["generate_psi_guess"] = "" -s.data["size_max"] = "256" +s.data["size_max"] = "3072" +print s + + + +s.data["subroutine"] = "H_apply_mrcc" +s.data["keys_work"] = "call mrcc_dress(delta_ij_sd_,Ndet_sd,i_generator,key_idx,keys_out,N_int,iproc)" print s END_SHELL diff --git a/src/MRCC/NEEDED_MODULES b/src/MRCC/NEEDED_MODULES index 8ae0a127..5e074d3c 100644 --- a/src/MRCC/NEEDED_MODULES +++ b/src/MRCC/NEEDED_MODULES @@ -1,2 +1,2 @@ -AOs Bielec_integrals Bitmask CAS_SD_selected Dets Electrons Ezfio_files Generators_CAS Generators_full Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full Utils +AOs Bielec_integrals Bitmask Dets Electrons Ezfio_files Generators_full Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full Utils diff --git a/src/MRCC/README.rst b/src/MRCC/README.rst index 1fbcccdb..702d19aa 100644 --- a/src/MRCC/README.rst +++ b/src/MRCC/README.rst @@ -1,14 +1,91 @@ --4.142795384334731 +=========== +MRCC Module +=========== -4.695183071437694E-002 - Determinant 64 - --------------------------------------------- - 000000000000002E|000000000000002E - |-+++-+----------------------------------------------------------| - |-+++-+----------------------------------------------------------| +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Bielec_integrals `_ +* `Bitmask `_ +* `Dets `_ +* `Electrons `_ +* `Ezfio_files `_ +* `Generators_full `_ +* `Hartree_Fock `_ +* `MOGuess `_ +* `MonoInts `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Perturbation `_ +* `Properties `_ +* `Selectors_full `_ +* `Utils `_ + +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +`mrcc `_ + Undocumented + +`run `_ + Undocumented + +`run_mrcc `_ + Undocumented + +`run_mrcc_test `_ + Undocumented + +`find_triples_and_quadruples `_ + Undocumented + +`mrcc_dress `_ + Undocumented + +`mrcc_dress_simple `_ + Undocumented + +`psi_cas_lock `_ + Locks on CAS determinants to fill delta_ij + +`ci_eigenvectors_dressed `_ + Eigenvectors/values of the CI matrix + +`ci_eigenvectors_s2_dressed `_ + Eigenvectors/values of the CI matrix + +`ci_electronic_energy_dressed `_ + Eigenvectors/values of the CI matrix + +`ci_energy_dressed `_ + N_states lowest eigenvalues of the dressed CI matrix + +`delta_ij `_ + Dressing matrix in N_det basis + +`delta_ij_non_cas `_ + Dressing matrix in SD basis + +`diagonalize_ci_dressed `_ + Replace the coefficients of the CI states by the coefficients of the + eigenstates of the CI matrix + +`dressing_type `_ + [ Simple | MRCC ] + +`h_matrix_dressed `_ + Dressed H with Delta_ij + +`lambda_mrcc `_ + cm/ -CAS-SD: -4.14214374069306 - : -4.14230904320551 -E0 = -11.5634986758976 diff --git a/src/MRCC/mrcc.irp.f b/src/MRCC/mrcc.irp.f index 663b9f5a..2b879882 100644 --- a/src/MRCC/mrcc.irp.f +++ b/src/MRCC/mrcc.irp.f @@ -3,41 +3,65 @@ program mrcc read_wf = .True. TOUCH read_wf call run + call run_mrcc +! call run_mrcc_test end subroutine run implicit none - print *, N_det - print *, N_det_cas - print *, N_det_sd - -! call update_generators - integer :: i + integer :: i,j print *, 'CAS' print *, '===' do i=1,N_det_cas - print *, psi_cas_coefs(i,:) + print *, psi_cas_coef(i,:) call debug_det(psi_cas(1,1,i),N_int) enddo ! print *, 'SD' ! print *, '==' -! do i=1,N_det_sd -! print *, psi_sd_coefs(i,:) -! call debug_det(psi_sd(1,1,i),N_int) -! enddo -! - print *, 'MRCC' - print *, '====' - print *, ci_energy(:) - print *, h_matrix_all_dets(3,3), delta_ij(3,3,1) - print *, h_matrix_all_dets(3,3), delta_ij(3,3,1) - print *, ci_energy_dressed(:) -! print *, 'max', maxval(delta_ij_sd) -! print *, 'min', minval(delta_ij_sd) -! -! do i=1,N_det -! print '(10(F10.6,X))', delta_ij(i,1:N_det,1) +! do i=1,N_det_non_cas +! print *, psi_non_cas_coef(i,:) +! call debug_det(psi_non_cas(1,1,i),N_int) ! enddo + call write_double(6,ci_energy(1),"Initial CI energy") +end +subroutine run_mrcc_test + implicit none + integer :: i,j + double precision :: pt2 + pt2 = 0.d0 + do j=1,N_det + do i=1,N_det + pt2 += psi_coef(i,1)*psi_coef(j,1) * delta_ij(i,j,1) + enddo + enddo + print *, ci_energy(1) + print *, ci_energy(1)+pt2 +end +subroutine run_mrcc + implicit none + integer :: i,j + + double precision :: E_new, E_old, delta_e + integer :: iteration + E_new = 0.d0 + delta_E = 1.d0 + iteration = 0 + do while (delta_E > 1.d-8) + iteration += 1 + print *, '===========================' + print *, 'MRCC Iteration', iteration + print *, '===========================' + print *, '' + E_old = sum(ci_energy_dressed) + call write_double(6,ci_energy_dressed(1),"MRCC energy") + call diagonalize_ci_dressed + E_new = sum(ci_energy_dressed) + delta_E = dabs(E_new - E_old) + enddo + call write_double(6,ci_energy_dressed(1),"Final MRCC energy") + call ezfio_set_mrcc_energy(ci_energy_dressed(1)) +! call save_wavefunction + end diff --git a/src/MRCC/mrcc_dress.irp.f b/src/MRCC/mrcc_dress.irp.f index e5c8d233..083cea5a 100644 --- a/src/MRCC/mrcc_dress.irp.f +++ b/src/MRCC/mrcc_dress.irp.f @@ -1,15 +1,169 @@ -subroutine mrcc_dress(delta_ij_sd_,Ndet_sd,i_generator,n_selected,det_buffer,Nint,iproc) +use omp_lib + +BEGIN_PROVIDER [ integer(omp_lock_kind), psi_cas_lock, (psi_det_size) ] + implicit none + BEGIN_DOC + ! Locks on CAS determinants to fill delta_ij + END_DOC + integer :: i + do i=1,psi_det_size + call omp_init_lock( psi_cas_lock(i) ) + enddo + +END_PROVIDER + +subroutine mrcc_dress(delta_ij_,Ndet,i_generator,n_selected,det_buffer,Nint,iproc) use bitmasks implicit none integer, intent(in) :: i_generator,n_selected, Nint, iproc - integer, intent(in) :: Ndet_sd - double precision, intent(inout) :: delta_ij_sd_(Ndet_sd,Ndet_sd,*) + integer, intent(in) :: Ndet + double precision, intent(inout) :: delta_ij_(Ndet,Ndet,*) + + integer(bit_kind), intent(in) :: det_buffer(Nint,2,n_selected) + integer :: i,j,k,l + integer :: degree_alpha(psi_det_size) + integer :: idx_alpha(0:psi_det_size) + logical :: good + + integer(bit_kind) :: tq(Nint,2,n_selected) + integer :: N_tq, c_ref ,degree + integer :: connected_to_ref + + double precision :: hIk, hla, hIl, dIk(N_states), dka(N_states), dIa(N_states) + double precision, allocatable :: dIa_hla(:,:) + double precision :: haj, phase, phase2 + double precision :: f(N_states), ci_inv(N_states) + integer :: exc(0:2,2,2) + integer :: h1,h2,p1,p2,s1,s2 + integer(bit_kind) :: tmp_det(Nint,2) + integer :: iint, ipos + integer :: i_state, k_sd, l_sd, i_I, i_alpha + + call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq) + + allocate (dIa_hla(N_states,Ndet)) + + ! |I> + + ! |alpha> + do i_alpha=1,N_tq + call get_excitation_degree_vector(psi_non_cas,tq(1,1,i_alpha),degree_alpha,Nint,N_det_non_cas,idx_alpha) + + ! |I> + do i_I=1,N_det_cas + ! Find triples and quadruple grand parents + call get_excitation_degree(tq(1,1,i_alpha),psi_cas(1,1,i_I),degree,Nint) + if (degree > 4) then + cycle + endif + + do i_state=1,N_states + dIa(i_state) = 0.d0 + enddo + + ! |alpha> + do k_sd=1,idx_alpha(0) + call get_excitation_degree(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(k_sd)),degree,Nint) + if (degree > 2) then + cycle + endif + ! + ! + call i_h_j(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(k_sd)),Nint,hIk) + do i_state=1,N_states + dIk(i_state) = hIk * lambda_mrcc(i_state,idx_alpha(k_sd)) + enddo + ! |l> = Exc(k -> alpha) |I> + call get_excitation(psi_non_cas(1,1,idx_alpha(k_sd)),tq(1,1,i_alpha),exc,degree,phase,Nint) + call decode_exc(exc,degree,h1,p1,h2,p2,s1,s2) + do k=1,N_int + tmp_det(k,1) = psi_cas(k,1,i_I) + tmp_det(k,2) = psi_cas(k,2,i_I) + enddo + ! Hole (see list_to_bitstring) + iint = ishft(h1-1,-bit_kind_shift) + 1 + ipos = h1-ishft((iint-1),bit_kind_shift)-1 + tmp_det(iint,s1) = ibclr(tmp_det(iint,s1),ipos) + + ! Particle + iint = ishft(p1-1,-bit_kind_shift) + 1 + ipos = p1-ishft((iint-1),bit_kind_shift)-1 + tmp_det(iint,s1) = ibset(tmp_det(iint,s1),ipos) + if (degree_alpha(k_sd) == 2) then + ! Hole (see list_to_bitstring) + iint = ishft(h2-1,-bit_kind_shift) + 1 + ipos = h2-ishft((iint-1),bit_kind_shift)-1 + tmp_det(iint,s2) = ibclr(tmp_det(iint,s2),ipos) + + ! Particle + iint = ishft(p2-1,-bit_kind_shift) + 1 + ipos = p2-ishft((iint-1),bit_kind_shift)-1 + tmp_det(iint,s2) = ibset(tmp_det(iint,s2),ipos) + endif + + ! + do i_state=1,N_states + dka(i_state) = 0.d0 + enddo + do l_sd=k_sd+1,idx_alpha(0) + call get_excitation_degree(tmp_det,psi_non_cas(1,1,idx_alpha(l_sd)),degree,Nint) + if (degree == 0) then + call get_excitation(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(l_sd)),exc,degree,phase2,Nint) + call i_h_j(psi_cas(1,1,i_I),psi_non_cas(1,1,idx_alpha(l_sd)),Nint,hIl) + do i_state=1,N_states + dka(i_state) = hIl * lambda_mrcc(i_state,idx_alpha(l_sd)) * phase * phase2 + enddo + exit + endif + enddo + do i_state=1,N_states + dIa(i_state) = dIa(i_state) + dIk(i_state) * dka(i_state) + enddo + enddo + + do i_state=1,N_states + ci_inv(i_state) = 1.d0/psi_cas_coef(i_I,i_state) + enddo + do l_sd=1,idx_alpha(0) + k_sd = idx_alpha(l_sd) + call i_h_j(tq(1,1,i_alpha),psi_non_cas(1,1,idx_alpha(l_sd)),Nint,hla) + do i_state=1,N_states + dIa_hla(i_state,k_sd) = dIa(i_state) * hla + enddo + enddo + call omp_set_lock( psi_cas_lock(i_I) ) + do l_sd=1,idx_alpha(0) + k_sd = idx_alpha(l_sd) + do i_state=1,N_states + delta_ij_(idx_non_cas(k_sd),idx_cas(i_I),i_state) += dIa_hla(i_state,k_sd) + delta_ij_(idx_cas(i_I),idx_non_cas(k_sd),i_state) += dIa_hla(i_state,k_sd) + delta_ij_(idx_cas(i_I),idx_cas(i_I),i_state) -= dIa_hla(i_state,k_sd) * ci_inv(i_state) * psi_non_cas_coef(k_sd,i_state) + enddo + enddo + call omp_unset_lock( psi_cas_lock(i_I) ) + enddo + enddo + deallocate (dIa_hla) +end + + + + + + + +subroutine mrcc_dress_simple(delta_ij_non_cas_,Ndet_non_cas,i_generator,n_selected,det_buffer,Nint,iproc) + use bitmasks + implicit none + + integer, intent(in) :: i_generator,n_selected, Nint, iproc + integer, intent(in) :: Ndet_non_cas + double precision, intent(inout) :: delta_ij_non_cas_(Ndet_non_cas,Ndet_non_cas,*) integer(bit_kind), intent(in) :: det_buffer(Nint,2,n_selected) integer :: i,j,k,m integer :: new_size - logical :: is_in_wavefunction integer :: degree(psi_det_size) integer :: idx(0:psi_det_size) logical :: good @@ -18,6 +172,51 @@ subroutine mrcc_dress(delta_ij_sd_,Ndet_sd,i_generator,n_selected,det_buffer,Nin integer :: N_tq, c_ref integer :: connected_to_ref + call find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq) + + ! Compute / (E0 - Haa) + double precision :: hka, haa + double precision :: haj + double precision :: f(N_states) + + do i=1,N_tq + call get_excitation_degree_vector(psi_non_cas,tq(1,1,i),degree,Nint,Ndet_non_cas,idx) + call i_h_j(tq(1,1,i),tq(1,1,i),Nint,haa) + do m=1,N_states + f(m) = 1.d0/(ci_electronic_energy(m)-haa) + enddo + do k=1,idx(0) + call i_h_j(tq(1,1,i),psi_non_cas(1,1,idx(k)),Nint,hka) + do j=k,idx(0) + call i_h_j(tq(1,1,i),psi_non_cas(1,1,idx(j)),Nint,haj) + do m=1,N_states + delta_ij_non_cas_(idx(k), idx(j),m) += haj*hka* f(m) + delta_ij_non_cas_(idx(j), idx(k),m) += haj*hka* f(m) + enddo + enddo + enddo + enddo +end + + +subroutine find_triples_and_quadruples(i_generator,n_selected,det_buffer,Nint,tq,N_tq) + use bitmasks + implicit none + + integer, intent(in) :: i_generator,n_selected, Nint + + integer(bit_kind), intent(in) :: det_buffer(Nint,2,n_selected) + integer :: i,j,k,m + logical :: is_in_wavefunction + integer :: degree(psi_det_size) + integer :: idx(0:psi_det_size) + logical :: good + + integer(bit_kind), intent(out) :: tq(Nint,2,n_selected) + integer, intent(out) :: N_tq + integer :: c_ref + integer :: connected_to_ref + N_tq = 0 do i=1,N_selected c_ref = connected_to_ref(det_buffer(1,1,i),psi_det_generators,Nint, & @@ -48,143 +247,11 @@ subroutine mrcc_dress(delta_ij_sd_,Ndet_sd,i_generator,n_selected,det_buffer,Nin endif enddo - ! Compute / (E0 - Haa) - double precision :: hka, haa - double precision :: haj - double precision :: f(N_states) - - do i=1,N_tq - call get_excitation_degree_vector(psi_sd,tq(1,1,i),degree,Nint,Ndet_sd,idx) - call i_h_j(tq(1,1,i),tq(1,1,i),Nint,haa) - do m=1,N_states - f(m) = 1.d0/(ci_electronic_energy(m)-haa) - enddo - do k=1,idx(0) - call i_h_j(tq(1,1,i),psi_sd(1,1,idx(k)),Nint,hka) - do j=k,idx(0) - call i_h_j(tq(1,1,i),psi_sd(1,1,idx(j)),Nint,haj) - do m=1,N_states - delta_ij_sd_(idx(k), idx(j),m) += haj*hka* f(m) - delta_ij_sd_(idx(j), idx(k),m) += haj*hka* f(m) - enddo - enddo - enddo - enddo end -BEGIN_PROVIDER [ double precision, delta_ij_sd, (N_det_sd, N_det_sd,N_states) ] - implicit none - BEGIN_DOC - ! Dressing matrix in SD basis - END_DOC - delta_ij_sd = 0.d0 - call H_apply_mrcc(delta_ij_sd,N_det_sd) - -END_PROVIDER - -BEGIN_PROVIDER [ double precision, delta_ij, (N_det,N_det,N_states) ] - implicit none - BEGIN_DOC - ! Dressing matrix in N_det basis - END_DOC - integer :: i,j,m - delta_ij = 0.d0 - do m=1,N_states - do j=1,N_det_sd - do i=1,N_det_sd - delta_ij(idx_sd(i),idx_sd(j),m) = delta_ij_sd(i,j,m) - enddo - enddo - enddo -END_PROVIDER - -BEGIN_PROVIDER [ double precision, h_matrix_dressed, (N_det,N_det) ] - implicit none - BEGIN_DOC - ! Dressed H with Delta_ij - END_DOC - integer :: i, j - do j=1,N_det - do i=1,N_det - h_matrix_dressed(i,j) = h_matrix_all_dets(i,j) + delta_ij(i,j,1) - if (i==j) then - print *, i, delta_ij(i,j,1), h_matrix_all_dets(i,j) - endif - enddo - enddo - -END_PROVIDER - BEGIN_PROVIDER [ double precision, CI_electronic_energy_dressed, (N_states_diag) ] -&BEGIN_PROVIDER [ double precision, CI_eigenvectors_dressed, (N_det,N_states_diag) ] -&BEGIN_PROVIDER [ double precision, CI_eigenvectors_s2_dressed, (N_states_diag) ] - implicit none - BEGIN_DOC - ! Eigenvectors/values of the CI matrix - END_DOC - integer :: i,j - - do j=1,N_states_diag - do i=1,N_det - CI_eigenvectors_dressed(i,j) = psi_coef(i,j) - enddo - enddo - - if (diag_algorithm == "Davidson") then - - stop 'use Lapack' -! call davidson_diag(psi_det,CI_eigenvectors_dressed,CI_electronic_energy_dressed, & -! size(CI_eigenvectors_dressed,1),N_det,N_states_diag,N_int,output_Dets) - - else if (diag_algorithm == "Lapack") then - - double precision, allocatable :: eigenvectors(:,:), eigenvalues(:) - allocate (eigenvectors(size(H_matrix_dressed,1),N_det)) - allocate (eigenvalues(N_det)) - call lapack_diag(eigenvalues,eigenvectors, & - H_matrix_dressed,size(H_matrix_dressed,1),N_det) - CI_electronic_energy_dressed(:) = 0.d0 - do i=1,N_det - CI_eigenvectors_dressed(i,1) = eigenvectors(i,1) - enddo - integer :: i_state - double precision :: s2 - i_state = 0 - do j=1,N_det - call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2) - if(dabs(s2-expected_s2).le.0.3d0)then - i_state += 1 - do i=1,N_det - CI_eigenvectors_dressed(i,i_state) = eigenvectors(i,j) - enddo - CI_electronic_energy_dressed(i_state) = eigenvalues(j) - CI_eigenvectors_s2_dressed(i_state) = s2 - endif - if (i_state.ge.N_states_diag) then - exit - endif - enddo - deallocate(eigenvectors,eigenvalues) - endif - -END_PROVIDER -BEGIN_PROVIDER [ double precision, CI_energy_dressed, (N_states_diag) ] - implicit none - BEGIN_DOC - ! N_states lowest eigenvalues of the dressed CI matrix - END_DOC - - integer :: j - character*(8) :: st - call write_time(output_Dets) - do j=1,N_states_diag - CI_energy_dressed(j) = CI_electronic_energy_dressed(j) + nuclear_repulsion - write(st,'(I4)') j - call write_double(output_Dets,CI_energy(j),'Energy of state '//trim(st)) - call write_double(output_Dets,CI_eigenvectors_s2(j),'S^2 of state '//trim(st)) - enddo -END_PROVIDER + diff --git a/src/MRCC/mrcc_utils.irp.f b/src/MRCC/mrcc_utils.irp.f index 2c5dc811..d33b7902 100644 --- a/src/MRCC/mrcc_utils.irp.f +++ b/src/MRCC/mrcc_utils.irp.f @@ -1,150 +1,161 @@ - use bitmasks -BEGIN_PROVIDER [ integer(bit_kind), cas_bitmask, (N_int,2,N_cas_bitmask) ] - implicit none - BEGIN_DOC - ! Bitmasks for CAS reference determinants. (N_int, alpha/beta, CAS reference) - END_DOC - logical :: exists - integer :: i - PROVIDE ezfio_filename - - call ezfio_has_bitmasks_cas(exists) - if (exists) then - call ezfio_get_bitmasks_cas(cas_bitmask) - else - do i=1,N_cas_bitmask - cas_bitmask(:,:,i) = iand(not(HF_bitmask(:,:)),full_ijkl_bitmask(:,:)) - enddo - endif - -END_PROVIDER - -BEGIN_PROVIDER [ integer, N_det_cas ] - implicit none - BEGIN_DOC - ! Number of generator detetrminants - END_DOC - integer :: i,k,l - logical :: good - call write_time(output_dets) - N_det_cas = 0 - do i=1,N_det - do l=1,n_cas_bitmask - good = .True. - do k=1,N_int - good = good .and. ( & - iand(not(cas_bitmask(k,1,l)), psi_det(k,1,i)) == & - iand(not(cas_bitmask(k,1,l)), psi_det(k,1,1)) ) .and. ( & - iand(not(cas_bitmask(k,2,l)), psi_det(k,2,i)) == & - iand(not(cas_bitmask(k,2,l)), psi_det(k,2,1)) ) - enddo - if (good) then - exit - endif - enddo - if (good) then - N_det_cas += 1 - endif - enddo - N_det_cas = max(N_det_cas, 1) - call write_int(output_dets,N_det_cas, 'Number of determinants in the CAS') -END_PROVIDER - - BEGIN_PROVIDER [ integer(bit_kind), psi_cas, (N_int,2,N_det_cas) ] -&BEGIN_PROVIDER [ double precision, psi_cas_coefs, (N_det_cas,n_states) ] -&BEGIN_PROVIDER [ integer, idx_cas, (N_det_cas) ] - implicit none - BEGIN_DOC - ! For Single reference wave functions, the generator is the - ! Hartree-Fock determinant - END_DOC - integer :: i, k, l, m - logical :: good - m=0 - do i=1,N_det - do l=1,n_cas_bitmask - good = .True. - do k=1,N_int - good = good .and. ( & - iand(not(cas_bitmask(k,1,l)), psi_det(k,1,i)) == & - iand(not(cas_bitmask(k,1,l)), psi_det(k,1,1)) ) .and. ( & - iand(not(cas_bitmask(k,2,l)), psi_det(k,2,i)) == & - iand(not(cas_bitmask(k,2,l)), psi_det(k,2,1)) ) - enddo - if (good) then - exit - endif - enddo - if (good) then - m = m+1 - do k=1,N_int - psi_cas(k,1,m) = psi_det(k,1,i) - psi_cas(k,2,m) = psi_det(k,2,i) - enddo - idx_cas(m) = i - do k=1,N_states - psi_cas_coefs(m,k) = psi_coef(i,k) - enddo - endif - enddo - -END_PROVIDER - - - - - BEGIN_PROVIDER [ integer(bit_kind), psi_sd, (N_int,2,N_det) ] -&BEGIN_PROVIDER [ double precision, psi_sd_coefs, (N_det,n_states) ] -&BEGIN_PROVIDER [ integer, idx_sd, (N_det) ] -&BEGIN_PROVIDER [ integer, N_det_sd] - implicit none - BEGIN_DOC - ! SD - END_DOC - integer :: i_sd,j,k - integer :: degree - logical :: in_cas - i_sd =0 - do k=1,N_det - in_cas = .False. - do j=1,N_det_cas - call get_excitation_degree(psi_cas(1,1,j), psi_det(1,1,k), degree, N_int) - if (degree == 0) then - in_cas = .True. - exit - endif - enddo - if (.not.in_cas) then - double precision :: hij - i_sd += 1 - psi_sd(1:N_int,1:2,i_sd) = psi_det(1:N_int,1:2,k) - psi_sd_coefs(i_sd,1:N_states) = psi_coef(k,1:N_states) - idx_sd(i_sd) = k - endif - enddo - N_det_sd = i_sd -END_PROVIDER - -BEGIN_PROVIDER [ double precision, lambda_mrcc, (psi_det_size,n_states) ] +BEGIN_PROVIDER [ double precision, lambda_mrcc, (N_states,psi_det_size) ] implicit none BEGIN_DOC ! cm/ END_DOC integer :: i,k double precision :: ihpsi(N_states) - do i=1,N_det_sd - call i_h_psi(psi_sd(1,1,i), psi_cas, psi_cas_coefs, N_int, N_det_cas, & - size(psi_cas_coefs,1), n_states, ihpsi) + do i=1,N_det_non_cas + call i_h_psi(psi_non_cas(1,1,i), psi_cas, psi_cas_coef, N_int, N_det_cas, & + size(psi_cas_coef,1), n_states, ihpsi) double precision :: hij do k=1,N_states - if (dabs(ihpsi(k)) < 1.d-6) then - lambda_mrcc(i,k) = 0.d0 + if (dabs(ihpsi(k)) > 1.d-5) then + lambda_mrcc(k,i) = psi_non_cas_coef(i,k)/ihpsi(k) + lambda_mrcc(k,i) = min( lambda_mrcc (k,i),0.d0 ) else - lambda_mrcc(i,k) = psi_sd_coefs(i,k)/ihpsi(k) - lambda_mrcc(i,k) = min( lambda_mrcc (i,k),0.d0 ) + lambda_mrcc(k,i) = 0.d0 endif enddo enddo END_PROVIDER + + +BEGIN_PROVIDER [ character*(32), dressing_type ] + implicit none + BEGIN_DOC + ! [ Simple | MRCC ] + END_DOC + dressing_type = "MRCC" +END_PROVIDER + +BEGIN_PROVIDER [ double precision, delta_ij_non_cas, (N_det_non_cas, N_det_non_cas,N_states) ] + implicit none + BEGIN_DOC + ! Dressing matrix in SD basis + END_DOC + delta_ij_non_cas = 0.d0 + call H_apply_mrcc_simple(delta_ij_non_cas,N_det_non_cas) +END_PROVIDER + +BEGIN_PROVIDER [ double precision, delta_ij, (N_det,N_det,N_states) ] + implicit none + BEGIN_DOC + ! Dressing matrix in N_det basis + END_DOC + integer :: i,j,m + delta_ij = 0.d0 + if (dressing_type == "MRCC") then + call H_apply_mrcc(delta_ij,N_det) + else if (dressing_type == "Simple") then + do m=1,N_states + do j=1,N_det_non_cas + do i=1,N_det_non_cas + delta_ij(idx_non_cas(i),idx_non_cas(j),m) = delta_ij_non_cas(i,j,m) + enddo + enddo + enddo + endif +END_PROVIDER + +BEGIN_PROVIDER [ double precision, h_matrix_dressed, (N_det,N_det) ] + implicit none + BEGIN_DOC + ! Dressed H with Delta_ij + END_DOC + integer :: i, j + do j=1,N_det + do i=1,N_det + h_matrix_dressed(i,j) = h_matrix_all_dets(i,j) + delta_ij(i,j,1) + enddo + enddo + +END_PROVIDER + + + BEGIN_PROVIDER [ double precision, CI_electronic_energy_dressed, (N_states_diag) ] +&BEGIN_PROVIDER [ double precision, CI_eigenvectors_dressed, (N_det,N_states_diag) ] +&BEGIN_PROVIDER [ double precision, CI_eigenvectors_s2_dressed, (N_states_diag) ] + implicit none + BEGIN_DOC + ! Eigenvectors/values of the CI matrix + END_DOC + integer :: i,j + + do j=1,N_states_diag + do i=1,N_det + CI_eigenvectors_dressed(i,j) = psi_coef(i,j) + enddo + enddo + + if (diag_algorithm == "Davidson") then + + stop 'use Lapack' +! call davidson_diag(psi_det,CI_eigenvectors_dressed,CI_electronic_energy_dressed, & +! size(CI_eigenvectors_dressed,1),N_det,N_states_diag,N_int,output_Dets) + + else if (diag_algorithm == "Lapack") then + + double precision, allocatable :: eigenvectors(:,:), eigenvalues(:) + allocate (eigenvectors(size(H_matrix_dressed,1),N_det)) + allocate (eigenvalues(N_det)) + call lapack_diag(eigenvalues,eigenvectors, & + H_matrix_dressed,size(H_matrix_dressed,1),N_det) + CI_electronic_energy_dressed(:) = 0.d0 + do i=1,N_det + CI_eigenvectors_dressed(i,1) = eigenvectors(i,1) + enddo + integer :: i_state + double precision :: s2 + i_state = 0 + do j=1,N_det + call get_s2_u0(psi_det,eigenvectors(1,j),N_det,N_det,s2) + if(dabs(s2-expected_s2).le.0.3d0)then + i_state += 1 + do i=1,N_det + CI_eigenvectors_dressed(i,i_state) = eigenvectors(i,j) + enddo + CI_electronic_energy_dressed(i_state) = eigenvalues(j) + CI_eigenvectors_s2_dressed(i_state) = s2 + endif + if (i_state.ge.N_states_diag) then + exit + endif + enddo + deallocate(eigenvectors,eigenvalues) + endif + +END_PROVIDER + +BEGIN_PROVIDER [ double precision, CI_energy_dressed, (N_states_diag) ] + implicit none + BEGIN_DOC + ! N_states lowest eigenvalues of the dressed CI matrix + END_DOC + + integer :: j + character*(8) :: st + call write_time(output_Dets) + do j=1,N_states_diag + CI_energy_dressed(j) = CI_electronic_energy_dressed(j) + nuclear_repulsion + enddo + +END_PROVIDER + +subroutine diagonalize_CI_dressed + implicit none + BEGIN_DOC +! Replace the coefficients of the CI states by the coefficients of the +! eigenstates of the CI matrix + END_DOC + integer :: i,j + do j=1,N_states_diag + do i=1,N_det + psi_coef(i,j) = CI_eigenvectors_dressed(i,j) + enddo + enddo + SOFT_TOUCH psi_coef + +end diff --git a/src/Makefile b/src/Makefile index 245868a8..237b3ae4 100644 --- a/src/Makefile +++ b/src/Makefile @@ -16,6 +16,10 @@ default: ezfio veryclean: $(QPACKAGE_ROOT)/scripts/clean_modules.sh $(ALL_MODULES) + # Define the dict [type in EZFIO.cfg] = ocaml type , f90 type + # If you change the qptypes_generator.ml, you need to rm this + # For simplicity add this to the veryclean rule + rm $(QPACKAGE_ROOT)/scripts/ezfio_interface/fancy_type.p $(ALL_MODULES): ezfio $(QPACKAGE_ROOT)/scripts/build_modules.sh $@ diff --git a/src/Molden/README.rst b/src/Molden/README.rst index 8660286d..d0e2343d 100644 --- a/src/Molden/README.rst +++ b/src/Molden/README.rst @@ -2,3 +2,40 @@ Molden Module ============= +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +`print_mos `_ + Undocumented + +`write_ao_basis `_ + Undocumented + +`write_geometry `_ + Undocumented + +`write_intro_gamess `_ + Undocumented + +`write_mo_basis `_ + Undocumented + + + +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +* `AOs `_ +* `Electrons `_ +* `Ezfio_files `_ +* `MOs `_ +* `Nuclei `_ +* `Output `_ +* `Utils `_ + diff --git a/src/MonoInts/README.rst b/src/MonoInts/README.rst index d10cee82..fdbb086b 100644 --- a/src/MonoInts/README.rst +++ b/src/MonoInts/README.rst @@ -18,30 +18,10 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`ao_mono_elec_integral `_ +`ao_mono_elec_integral `_ array of the mono electronic hamiltonian on the AOs basis : sum of the kinetic and nuclear electronic potential -`ao_overlap `_ - Overlap between atomic basis functions: - :math:`\int \chi_i(r) \chi_j(r) dr)` - -`ao_overlap_abs `_ - Overlap between absolute value of atomic basis functions: - :math:`\int |\chi_i(r)| |\chi_j(r)| dr)` - -`ao_overlap_x `_ - Overlap between atomic basis functions: - :math:`\int \chi_i(r) \chi_j(r) dr)` - -`ao_overlap_y `_ - Overlap between atomic basis functions: - :math:`\int \chi_i(r) \chi_j(r) dr)` - -`ao_overlap_z `_ - Overlap between atomic basis functions: - :math:`\int \chi_i(r) \chi_j(r) dr)` - `check_ortho `_ Undocumented @@ -72,58 +52,63 @@ Documentation .br {\tt ao_deriv2_x} = \langle \chi_i(x,y,z) \frac{\partial^2}{\partial x^2} |\chi_j (x,y,z) \rangle -`ao_kinetic_integral `_ +`ao_kinetic_integral `_ array of the priminitve basis kinetic integrals \langle \chi_i |\hat{T}| \chi_j \rangle `mo_kinetic_integral `_ Undocumented -`mo_mono_elec_integral `_ +`mo_mono_elec_integral `_ array of the mono electronic hamiltonian on the MOs basis : sum of the kinetic and nuclear electronic potential -`mo_overlap `_ - Undocumented - `orthonormalize_mos `_ Undocumented `ao_nucl_elec_integral `_ interaction nuclear electron -`give_polynom_mult_center_mono_elec `_ +`ao_nucl_elec_integral_per_atom `_ + ao_nucl_elec_integral_per_atom(i,j,k) = - + where Rk is the geometry of the kth atom + +`give_polynom_mult_center_mono_elec `_ Undocumented -`i_x1_pol_mult_mono_elec `_ +`i_x1_pol_mult_mono_elec `_ Undocumented -`i_x2_pol_mult_mono_elec `_ +`i_x2_pol_mult_mono_elec `_ Undocumented -`int_gaus_pol `_ +`int_gaus_pol `_ Undocumented -`nai_pol_mult `_ +`nai_pol_mult `_ Undocumented -`v_e_n `_ +`v_e_n `_ Undocumented -`v_phi `_ +`v_phi `_ Undocumented -`v_r `_ +`v_r `_ Undocumented -`v_theta `_ +`v_theta `_ Undocumented -`wallis `_ +`wallis `_ Undocumented `mo_nucl_elec_integral `_ - Undocumented + interaction nuclear electron on the MO basis + +`mo_nucl_elec_integral_per_atom `_ + mo_nucl_elec_integral_per_atom(i,j,k) = - + where Rk is the geometry of the kth atom `save_ortho_mos `_ Undocumented @@ -173,7 +158,7 @@ Documentation array of the integrals of AO_i * y^2 AO_j array of the integrals of AO_i * z^2 AO_j -`overlap_bourrin_deriv_x `_ +`overlap_bourrin_deriv_x `_ Undocumented `overlap_bourrin_dipole `_ @@ -182,7 +167,7 @@ Documentation `overlap_bourrin_spread `_ Undocumented -`overlap_bourrin_x `_ +`overlap_bourrin_x `_ Undocumented `overlap_bourrin_x_abs `_ diff --git a/src/NEEDED_MODULES b/src/NEEDED_MODULES index a144c42c..efe8b8f8 100644 --- a/src/NEEDED_MODULES +++ b/src/NEEDED_MODULES @@ -1 +1 @@ -AOs Bielec_integrals Bitmask CID CID_SC2_selected CID_selected CIS CISD CISD_selected CISD_SC2_selected Dets Electrons Ezfio_files Full_CI Generators_full Hartree_Fock MOGuess MonoInts MOs MP2 Nuclei Output Selectors_full Utils Molden FCIdump Generators_CAS CAS_SD_selected DDCI_selected MRCC +AOs Bielec_integrals Bitmask CID CID_SC2_selected CID_selected CIS CISD CISD_selected CISD_SC2_selected Dets Electrons Ezfio_files Full_CI Generators_full Hartree_Fock MOGuess MonoInts MOs MP2 Nuclei Output Selectors_full Utils Molden FCIdump Generators_CAS CAS_SD DDCI_selected MRCC diff --git a/src/Nuclei/README.rst b/src/Nuclei/README.rst index cf164529..b21d02ee 100644 --- a/src/Nuclei/README.rst +++ b/src/Nuclei/README.rst @@ -22,6 +22,9 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`element_name `_ + Array of the name of element, sorted by nuclear charge (integer) + `nucl_charge `_ Nuclear charges @@ -65,8 +68,11 @@ Documentation `nucl_num_aligned `_ Number of nuclei -`nuclear_repulsion `_ +`nuclear_repulsion `_ Nuclear repulsion energy +`positive_charge_barycentre `_ + Centroid of the positive charges + diff --git a/src/Output/README.rst b/src/Output/README.rst index edbb7ff4..adcae302 100644 --- a/src/Output/README.rst +++ b/src/Output/README.rst @@ -39,5 +39,23 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`output_cpu_time_0 `_ + Initial CPU and wall times when printing in the output files + +`output_wall_time_0 `_ + Initial CPU and wall times when printing in the output files + +`write_bool `_ + Write an logical value in output + +`write_double `_ + Write a double precision value in output + +`write_int `_ + Write an integer value in output + +`write_time `_ + Write a time stamp in the output for chronological reconstruction + diff --git a/src/Selectors_full/README.rst b/src/Selectors_full/README.rst index 50010c0a..aaa07bbd 100644 --- a/src/Selectors_full/README.rst +++ b/src/Selectors_full/README.rst @@ -8,7 +8,18 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. -`coef_hf_selector `_ +`coef_hf_selector `_ + energy of correlation per determinant respect to the Hartree Fock determinant + .br + for the all the double excitations in the selectors determinants + .br + E_corr_per_selectors(i) = * c(D_i)/c(HF) if |D_i> is a double excitation + .br + E_corr_per_selectors(i) = -1000.d0 if it is not a double excitation + .br + coef_hf_selector = coefficient of the Hartree Fock determinant in the selectors determinants + +`delta_e_per_selector `_ energy of correlation per determinant respect to the Hartree Fock determinant .br for the all the double excitations in the selectors determinants @@ -28,7 +39,29 @@ Documentation .br n_double_selectors = number of double excitations in the selectors determinants -`e_corr_per_selectors `_ +`e_corr_double_only `_ + energy of correlation per determinant respect to the Hartree Fock determinant + .br + for the all the double excitations in the selectors determinants + .br + E_corr_per_selectors(i) = * c(D_i)/c(HF) if |D_i> is a double excitation + .br + E_corr_per_selectors(i) = -1000.d0 if it is not a double excitation + .br + coef_hf_selector = coefficient of the Hartree Fock determinant in the selectors determinants + +`e_corr_per_selectors `_ + energy of correlation per determinant respect to the Hartree Fock determinant + .br + for the all the double excitations in the selectors determinants + .br + E_corr_per_selectors(i) = * c(D_i)/c(HF) if |D_i> is a double excitation + .br + E_corr_per_selectors(i) = -1000.d0 if it is not a double excitation + .br + coef_hf_selector = coefficient of the Hartree Fock determinant in the selectors determinants + +`e_corr_second_order `_ energy of correlation per determinant respect to the Hartree Fock determinant .br for the all the double excitations in the selectors determinants @@ -48,6 +81,39 @@ Documentation .br n_double_selectors = number of double excitations in the selectors determinants +`i_h_hf_per_selectors `_ + energy of correlation per determinant respect to the Hartree Fock determinant + .br + for the all the double excitations in the selectors determinants + .br + E_corr_per_selectors(i) = * c(D_i)/c(HF) if |D_i> is a double excitation + .br + E_corr_per_selectors(i) = -1000.d0 if it is not a double excitation + .br + coef_hf_selector = coefficient of the Hartree Fock determinant in the selectors determinants + +`inv_selectors_coef_hf `_ + energy of correlation per determinant respect to the Hartree Fock determinant + .br + for the all the double excitations in the selectors determinants + .br + E_corr_per_selectors(i) = * c(D_i)/c(HF) if |D_i> is a double excitation + .br + E_corr_per_selectors(i) = -1000.d0 if it is not a double excitation + .br + coef_hf_selector = coefficient of the Hartree Fock determinant in the selectors determinants + +`inv_selectors_coef_hf_squared `_ + energy of correlation per determinant respect to the Hartree Fock determinant + .br + for the all the double excitations in the selectors determinants + .br + E_corr_per_selectors(i) = * c(D_i)/c(HF) if |D_i> is a double excitation + .br + E_corr_per_selectors(i) = -1000.d0 if it is not a double excitation + .br + coef_hf_selector = coefficient of the Hartree Fock determinant in the selectors determinants + `n_double_selectors `_ degree of excitation respect to Hartree Fock for the wave function .br @@ -57,23 +123,40 @@ Documentation .br n_double_selectors = number of double excitations in the selectors determinants -`n_det_selectors `_ +`n_det_selectors `_ For Single reference wave functions, the number of selectors is 1 : the Hartree-Fock determinant -`psi_selectors `_ +`psi_selectors `_ Determinants on which we apply for perturbation. -`psi_selectors_coef `_ +`psi_selectors_ab `_ + Determinants on which we apply . + 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_selectors_coef `_ Determinants on which we apply for perturbation. -`psi_selectors_size `_ +`psi_selectors_coef_ab `_ + Determinants on which we apply . + 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_selectors_diag_h_mat `_ + Diagonal elements of the H matrix for each selectors + +`psi_selectors_next_ab `_ + Determinants on which we apply . + 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_selectors_size `_ Undocumented -`threshold_selectors `_ - Percentage of the norm of the state-averaged wave function to - consider for the selectors - Needed Modules @@ -83,12 +166,13 @@ Needed Modules .. NEEDED_MODULES file. * `AOs `_ -* `BiInts `_ +* `Bielec_integrals `_ * `Bitmask `_ * `Dets `_ * `Electrons `_ * `Ezfio_files `_ * `Hartree_Fock `_ +* `MOGuess `_ * `MonoInts `_ * `MOs `_ * `Nuclei `_ diff --git a/src/Utils/LinearAlgebra.irp.f b/src/Utils/LinearAlgebra.irp.f index 7d173246..a7462f94 100644 --- a/src/Utils/LinearAlgebra.irp.f +++ b/src/Utils/LinearAlgebra.irp.f @@ -451,3 +451,6 @@ subroutine set_zero_extra_diag(i1,i2,matrix,lda,m) end + + + diff --git a/src/Utils/README.rst b/src/Utils/README.rst index 283a3779..68a4f080 100644 --- a/src/Utils/README.rst +++ b/src/Utils/README.rst @@ -10,35 +10,104 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`apply_rotation `_ + Apply the rotation found by find_rotation + +`find_rotation `_ + Find A.C = B + +`get_pseudo_inverse `_ + Find C = A^-1 + +`lapack_diag `_ + Diagonalize matrix H + .br + H is untouched between input and ouptut + .br + eigevalues(i) = ith lowest eigenvalue of the H matrix + .br + eigvectors(i,j) = where i is the basis function and psi_j is the j th eigenvector + .br + +`lapack_diag_s2 `_ + Diagonalize matrix H + .br + H is untouched between input and ouptut + .br + eigevalues(i) = ith lowest eigenvalue of the H matrix + .br + eigvectors(i,j) = where i is the basis function and psi_j is the j th eigenvector + .br + +`lapack_diagd `_ + Diagonalize matrix H + .br + H is untouched between input and ouptut + .br + eigevalues(i) = ith lowest eigenvalue of the H matrix + .br + eigvectors(i,j) = where i is the basis function and psi_j is the j th eigenvector + .br + +`lapack_partial_diag `_ + Diagonalize matrix H + .br + H is untouched between input and ouptut + .br + eigevalues(i) = ith lowest eigenvalue of the H matrix + .br + eigvectors(i,j) = where i is the basis function and psi_j is the j th eigenvector + .br + +`ortho_lowdin `_ + Compute C_new=C_old.S^-1/2 canonical orthogonalization. + .br + overlap : overlap matrix + .br + LDA : leftmost dimension of overlap array + .br + N : Overlap matrix is NxN (array is (LDA,N) ) + .br + C : Coefficients of the vectors to orthogonalize. On exit, + orthogonal vectors + .br + LDC : leftmost dimension of C + .br + m : Coefficients matrix is MxN, ( array is (LDC,N) ) + .br + +`set_zero_extra_diag `_ + Undocumented + `abort_all `_ If True, all the calculation is aborted -`abort_here `_ +`abort_here `_ If True, all the calculation is aborted -`catch_signal `_ +`catch_signal `_ What to do on Ctrl-C. If two Ctrl-C are pressed within 1 sec, the calculation if aborted. -`trap_signals `_ +`trap_signals `_ What to do when a signal is caught. Here, trap Ctrl-C and call the control_C subroutine. -`add_poly `_ +`add_poly `_ Add two polynomials D(t) =! D(t) +( B(t)+C(t)) -`add_poly_multiply `_ +`add_poly_multiply `_ Add a polynomial multiplied by a constant D(t) =! D(t) +( cst * B(t)) -`f_integral `_ +`f_integral `_ function that calculates the following integral \int_{\-infty}^{+\infty} x^n \exp(-p x^2) dx -`gaussian_product `_ +`gaussian_product `_ Gaussian product in 1D. e^{-a (x-x_A)^2} e^{-b (x-x_B)^2} = K_{ab}^x e^{-p (x-x_P)^2} -`gaussian_product_x `_ +`gaussian_product_x `_ Gaussian product in 1D. e^{-a (x-x_A)^2} e^{-b (x-x_B)^2} = K_{ab}^x e^{-p (x-x_P)^2} @@ -50,108 +119,148 @@ Documentation * [ sum (l_y = 0,i_order(2)) P_new(l_y,2) * (y-P_center(2))^l_y ] exp (- p (y-P_center(2))^2 ) * [ sum (l_z = 0,i_order(3)) P_new(l_z,3) * (z-P_center(3))^l_z ] exp (- p (z-P_center(3))^2 ) +`give_explicit_poly_and_gaussian_double `_ + Transforms the product of + (x-x_A)^a(1) (x-x_B)^b(1) (x-x_A)^a(2) (y-y_B)^b(2) (z-z_A)^a(3) (z-z_B)^b(3) + exp(-(r-A)^2 alpha) exp(-(r-B)^2 beta) exp(-(r-Nucl_center)^2 gama + .br + into + fact_k * [ sum (l_x = 0,i_order(1)) P_new(l_x,1) * (x-P_center(1))^l_x ] exp (- p (x-P_center(1))^2 ) + * [ sum (l_y = 0,i_order(2)) P_new(l_y,2) * (y-P_center(2))^l_y ] exp (- p (y-P_center(2))^2 ) + * [ sum (l_z = 0,i_order(3)) P_new(l_z,3) * (z-P_center(3))^l_z ] exp (- p (z-P_center(3))^2 ) + `give_explicit_poly_and_gaussian_x `_ Transform the product of (x-x_A)^a(1) (x-x_B)^b(1) (x-x_A)^a(2) (y-y_B)^b(2) (z-z_A)^a(3) (z-z_B)^b(3) exp(-(r-A)^2 alpha) exp(-(r-B)^2 beta) into fact_k (x-x_P)^iorder(1) (y-y_P)^iorder(2) (z-z_P)^iorder(3) exp(-p(r-P)^2) -`hermite `_ +`hermite `_ Hermite polynomial -`multiply_poly `_ +`multiply_poly `_ Multiply two polynomials D(t) =! D(t) +( B(t)*C(t)) -`recentered_poly2 `_ +`recentered_poly2 `_ Recenter two polynomials -`rint `_ +`rint `_ .. math:: .br \int_0^1 dx \exp(-p x^2) x^n .br -`rint1 `_ +`rint1 `_ Standard version of rint -`rint_large_n `_ +`rint_large_n `_ Version of rint for large values of n -`rint_sum `_ +`rint_sum `_ Needed for the calculation of two-electron integrals. +`overlap_a_b_c `_ + Undocumented + `overlap_gaussian_x `_ .. math:: .br \sum_{-infty}^{+infty} (x-A_x)^ax (x-B_x)^bx exp(-alpha(x-A_x)^2) exp(-beta(x-B_X)^2) dx .br -`overlap_gaussian_xyz `_ +`overlap_gaussian_xyz `_ .. math:: .br S_x = \int (x-A_x)^{a_x} exp(-\alpha(x-A_x)^2) (x-B_x)^{b_x} exp(-beta(x-B_x)^2) dx \\ S = S_x S_y S_z .br -`overlap_x_abs `_ +`overlap_x_abs `_ .. math :: .br \int_{-infty}^{+infty} (x-A_center)^(power_A) * (x-B_center)^power_B * exp(-alpha(x-A_center)^2) * exp(-beta(x-B_center)^2) dx .br -`align_double `_ +`progress_active `_ + Current status for displaying progress bars. Global variable. + +`progress_bar `_ + Current status for displaying progress bars. Global variable. + +`progress_timeout `_ + Current status for displaying progress bars. Global variable. + +`progress_title `_ + Current status for displaying progress bars. Global variable. + +`progress_value `_ + Current status for displaying progress bars. Global variable. + +`run_progress `_ + Display a progress bar with documentation of what is happening + +`start_progress `_ + Starts the progress bar + +`stop_progress `_ + Stop the progress bar + +`align_double `_ Compute 1st dimension such that it is aligned for vectorization. -`all_utils `_ - Dummy provider to provide all utils - -`approx_dble `_ +`approx_dble `_ Undocumented -`binom `_ +`binom `_ Binomial coefficients -`binom_func `_ +`binom_func `_ .. math :: .br \frac{i!}{j!(i-j)!} .br -`binom_transp `_ +`binom_transp `_ Binomial coefficients -`dble_fact `_ +`dble_fact `_ n!! -`fact `_ +`dble_logfact `_ + n!! + +`fact `_ n! -`fact_inv `_ +`fact_inv `_ 1/n! -`inv_int `_ +`inv_int `_ 1/i -`normalize `_ +`logfact `_ + n! + +`normalize `_ Normalizes vector u u is expected to be aligned in memory. -`nproc `_ +`nproc `_ Number of current OpenMP threads -`u_dot_u `_ +`u_dot_u `_ Compute -`u_dot_v `_ +`u_dot_v `_ Compute -`wall_time `_ +`wall_time `_ The equivalent of cpu_time, but for the wall time. -`write_git_log `_ +`write_git_log `_ Write the last git commit in file iunit. - +