diff --git a/README.md b/README.md index d5496f5e..b8141b88 100644 --- a/README.md +++ b/README.md @@ -36,8 +36,8 @@ https://arxiv.org/abs/1902.08154 # Build status -* Master [![master build status](https://travis-ci.org/QuantumPackage/qp2.svg?branch=master)](https://travis-ci.org/QuantumPackage/qp2) -* Development [![dev build status](https://travis-ci.org/QuantumPackage/qp2.svg?branch=dev)](https://travis-ci.org/QuantumPackage/qp2) +* Master [![master build status](https://travis-ci.com/QuantumPackage/qp2.svg?branch=master)](https://travis-ci.org/QuantumPackage/qp2) +* Development [![dev build status](https://travis-ci.com/QuantumPackage/qp2.svg?branch=dev)](https://travis-ci.org/QuantumPackage/qp2) * Documentation [![Documentation Status](https://readthedocs.org/projects/quantum-package/badge/?version=master)](https://quantum-package.readthedocs.io/en/master/?badge=master) diff --git a/RELEASE_NOTES.org b/RELEASE_NOTES.org index f9dbcbb0..758618a9 100644 --- a/RELEASE_NOTES.org +++ b/RELEASE_NOTES.org @@ -29,6 +29,7 @@ - Disk-based Davidson when too much memory is required - Fixed bug in DIIS - Fixed bug in molden (Au -> Angs) + - Fixed bug with non-contiguous MOs in active space and deleter MOs *** User interface @@ -52,6 +53,10 @@ - Added ~print_hamiltonian~ - Added input for two body RDM - Added keyword ~save_wf_after_selection~ + - Added a ~restore_symm~ flag to enforce the restoration of + symmetry in matrices + - qp_export_as_tgz exports also plugin codes + - Added a basis module containing basis set information *** Code @@ -76,6 +81,8 @@ - Added ~h_core_guess~ routine - Fixed Laplacians in real space (indices) - Added LIB file to add extra libs in plugin + - Using Intel IPP for sorting when using Intel compiler + - Removed parallelism in sorting ao_one_e_integral_zero banned_excitations diff --git a/VERSION b/VERSION new file mode 100644 index 00000000..c043eea7 --- /dev/null +++ b/VERSION @@ -0,0 +1 @@ +2.2.1 diff --git a/bin/qp_export_as_tgz b/bin/qp_export_as_tgz index 7332d1ed..20624dba 100755 --- a/bin/qp_export_as_tgz +++ b/bin/qp_export_as_tgz @@ -99,7 +99,9 @@ function find_libs () { } function find_exec () { - find ${QP_ROOT}/$1 -perm /u+x -type f + for i in $@ ; do + find ${QP_ROOT}/$i -perm /u+x -type f + done } @@ -119,7 +121,7 @@ fi echo "Copying binary files" # -------------------- -FORTRAN_EXEC=$(find_exec src) +FORTRAN_EXEC=$(find_exec src/*/) if [[ -z $FORTRAN_EXEC ]] ; then error 'No Fortran binaries found.' exit 1 diff --git a/config/ifort.cfg b/config/ifort.cfg index 866aae3d..714c4b10 100644 --- a/config/ifort.cfg +++ b/config/ifort.cfg @@ -7,9 +7,9 @@ # [COMMON] FC : ifort -fpic -LAPACK_LIB : -mkl=parallel +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=32 +IRPF90_FLAGS : --ninja --align=32 -DINTEL # Global options ################ diff --git a/config/ifort_avx.cfg b/config/ifort_avx.cfg index d689050e..a2cb4c8a 100644 --- a/config/ifort_avx.cfg +++ b/config/ifort_avx.cfg @@ -7,9 +7,9 @@ # [COMMON] FC : ifort -fpic -LAPACK_LIB : -mkl=parallel +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=32 +IRPF90_FLAGS : --ninja --align=32 -DINTEL # Global options ################ diff --git a/config/ifort_avx_mpi.cfg b/config/ifort_avx_mpi.cfg index 9a839e66..f2bb8889 100644 --- a/config/ifort_avx_mpi.cfg +++ b/config/ifort_avx_mpi.cfg @@ -7,9 +7,9 @@ # [COMMON] FC : mpiifort -fpic -LAPACK_LIB : -mkl=parallel +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=32 -DMPI +IRPF90_FLAGS : --ninja --align=32 -DMPI -DINTEL # Global options ################ diff --git a/config/ifort_debug.cfg b/config/ifort_debug.cfg index f2fbf8ce..9b718380 100644 --- a/config/ifort_debug.cfg +++ b/config/ifort_debug.cfg @@ -7,9 +7,9 @@ # [COMMON] FC : ifort -fpic -LAPACK_LIB : -mkl=parallel +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=32 --assert +IRPF90_FLAGS : --ninja --align=32 --assert -DINTEL # Global options ################ diff --git a/config/ifort_mpi.cfg b/config/ifort_mpi.cfg index 57087847..e0d489a0 100644 --- a/config/ifort_mpi.cfg +++ b/config/ifort_mpi.cfg @@ -7,9 +7,9 @@ # [COMMON] FC : mpiifort -fpic -LAPACK_LIB : -mkl=parallel +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=32 -DMPI +IRPF90_FLAGS : --ninja --align=32 -DMPI -DINTEL # Global options ################ diff --git a/config/ifort_rome.cfg b/config/ifort_rome.cfg index a4bce680..5ed01227 100644 --- a/config/ifort_rome.cfg +++ b/config/ifort_rome.cfg @@ -7,9 +7,9 @@ # [COMMON] FC : ifort -fpic -LAPACK_LIB : -mkl=parallel +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=32 +IRPF90_FLAGS : --ninja --align=32 -DINTEL # Global options ################ @@ -31,8 +31,8 @@ OPENMP : 1 ; Append OpenMP flags # -ftz : Flushes denormal results to zero # [OPT] -FC : -traceback -FCFLAGS : -O2 -ip -g -march=core-avx2 -align array64byte -fma -ftz -fomit-frame-pointer +FC : -traceback -shared-intel +FCFLAGS : -O2 -ip -g -march=core-avx2 -align array64byte -fma -ftz -fomit-frame-pointer # Profiling flags ################# diff --git a/config/ifort_xHost.cfg b/config/ifort_xHost.cfg index 5d952e54..ddb4aa2d 100644 --- a/config/ifort_xHost.cfg +++ b/config/ifort_xHost.cfg @@ -6,10 +6,10 @@ # --align=32 : Align all provided arrays on a 32-byte boundary # [COMMON] -FC : ifort -fpic -LAPACK_LIB : -mkl=parallel +FC : ifort -fpic +LAPACK_LIB : -mkl=parallel -lirc -lsvml -limf -lipps IRPF90 : irpf90 -IRPF90_FLAGS : --ninja --align=64 +IRPF90_FLAGS : --ninja --align=64 -DINTEL # Global options ################ diff --git a/data/basis/cc-pvdz b/data/basis/cc-pvdz index 3d1a77e6..156e46ac 100644 --- a/data/basis/cc-pvdz +++ b/data/basis/cc-pvdz @@ -3606,4 +3606,4 @@ D 5 5 1.507524E+00 2.667560E-01 D 1 1 5.030000E-01 1.000000E+00 -$END \ No newline at end of file +$END diff --git a/etc/local.rc b/etc/local.rc index 954b315b..182599c6 100644 --- a/etc/local.rc +++ b/etc/local.rc @@ -19,4 +19,3 @@ # export QP_NIC=lo # export QP_NIC=ib0 - diff --git a/external/ezfio b/external/ezfio index ccee52d0..ed1df9f3 160000 --- a/external/ezfio +++ b/external/ezfio @@ -1 +1 @@ -Subproject commit ccee52d00c2cde1d628b0d34f4a247143747bf36 +Subproject commit ed1df9f3c1f51752656ca98da5693a4119add05c diff --git a/external/irpf90 b/external/irpf90 index 132a4a16..33ca5e10 160000 --- a/external/irpf90 +++ b/external/irpf90 @@ -1 +1 @@ -Subproject commit 132a4a1661c9878d21dcbf0ac14f7fe9a3b110d0 +Subproject commit 33ca5e1018f3bbb5e695e6ee558f5dac0753b271 diff --git a/ocaml/Symmetry.ml b/ocaml/Angmom.ml similarity index 95% rename from ocaml/Symmetry.ml rename to ocaml/Angmom.ml index eb4b637b..ed13e8dc 100644 --- a/ocaml/Symmetry.ml +++ b/ocaml/Angmom.ml @@ -26,7 +26,7 @@ let of_string = function | "J" | "j" -> J | "K" | "k" -> K | "L" | "l" -> L - | x -> raise (Failure ("Symmetry should be S|P|D|F|G|H|I|J|K|L, + | x -> raise (Failure ("Angmom should be S|P|D|F|G|H|I|J|K|L, not "^x^".")) let of_char = function @@ -40,7 +40,7 @@ let of_char = function | 'J' | 'j' -> J | 'K' | 'k' -> K | 'L' | 'l' -> L - | x -> raise (Failure ("Symmetry should be S|P|D|F|G|H|I|J|K|L")) + | x -> raise (Failure ("Angmom should be S|P|D|F|G|H|I|J|K|L")) let to_l = function | S -> Positive_int.of_int 0 @@ -68,7 +68,7 @@ let of_l i = | 7 -> J | 8 -> K | 9 -> L - | x -> raise (Failure ("Symmetry should be S|P|D|F|G|H|I|J|K|L")) + | x -> raise (Failure ("Angmom should be S|P|D|F|G|H|I|J|K|L")) type st = t diff --git a/ocaml/Symmetry.mli b/ocaml/Angmom.mli similarity index 100% rename from ocaml/Symmetry.mli rename to ocaml/Angmom.mli diff --git a/ocaml/GaussianPrimitive.ml b/ocaml/GaussianPrimitive.ml index d2144dec..c59a68de 100644 --- a/ocaml/GaussianPrimitive.ml +++ b/ocaml/GaussianPrimitive.ml @@ -2,14 +2,14 @@ open Qptypes open Sexplib.Std type t = - { sym : Symmetry.t ; + { sym : Angmom.t ; expo : AO_expo.t ; } [@@deriving sexp] let to_string p = let { sym = s ; expo = e } = p in Printf.sprintf "(%s, %22e)" - (Symmetry.to_string s) + (Angmom.to_string s) (AO_expo.to_float e) diff --git a/ocaml/Gto.ml b/ocaml/Gto.ml index 465e5627..200b4fa6 100644 --- a/ocaml/Gto.ml +++ b/ocaml/Gto.ml @@ -9,7 +9,7 @@ type fmt = | Gaussian type t = -{ sym : Symmetry.t ; +{ sym : Angmom.t ; lc : ((GaussianPrimitive.t * AO_coef.t) list) } [@@deriving sexp] @@ -47,7 +47,7 @@ let read_one in_channel = in let sym_str = String.sub buffer 0 2 in let n_str = String.sub buffer 2 ((String.length buffer)-2) in - let sym = Symmetry.of_string (String_ext.strip sym_str) in + let sym = Angmom.of_string (String_ext.strip sym_str) in let n = int_of_string (String_ext.strip n_str) in (* Read all the primitives *) let rec read_lines result = function @@ -82,7 +82,7 @@ let read_one in_channel = (** Write the GTO in Gamess format *) let to_string_gamess { sym = sym ; lc = lc } = let result = - Printf.sprintf "%s %3d" (Symmetry.to_string sym) (List.length lc) + Printf.sprintf "%s %3d" (Angmom.to_string sym) (List.length lc) in let rec do_work accu i = function | [] -> List.rev accu @@ -102,7 +102,7 @@ let to_string_gamess { sym = sym ; lc = lc } = (** Write the GTO in Gaussian format *) let to_string_gaussian { sym = sym ; lc = lc } = let result = - Printf.sprintf "%s %3d 1.00" (Symmetry.to_string sym) (List.length lc) + Printf.sprintf "%s %3d 1.00" (Angmom.to_string sym) (List.length lc) in let rec do_work accu i = function | [] -> List.rev accu diff --git a/ocaml/Gto.mli b/ocaml/Gto.mli index 91534ebe..ccbac7fe 100644 --- a/ocaml/Gto.mli +++ b/ocaml/Gto.mli @@ -5,7 +5,7 @@ type fmt = | Gaussian type t = - { sym : Symmetry.t ; + { sym : Angmom.t ; lc : (GaussianPrimitive.t * Qptypes.AO_coef.t) list; } [@@deriving sexp] diff --git a/ocaml/Input_ao_basis.ml b/ocaml/Input_ao_basis.ml index a73d641c..95d37a7a 100644 --- a/ocaml/Input_ao_basis.ml +++ b/ocaml/Input_ao_basis.ml @@ -9,7 +9,7 @@ module Ao_basis : sig ao_prim_num : AO_prim_number.t array; ao_prim_num_max : AO_prim_number.t; ao_nucl : Nucl_number.t array; - ao_power : Symmetry.Xyz.t array; + ao_power : Angmom.Xyz.t array; ao_coef : AO_coef.t array; ao_expo : AO_expo.t array; ao_cartesian : bool; @@ -32,7 +32,7 @@ end = struct ao_prim_num : AO_prim_number.t array; ao_prim_num_max : AO_prim_number.t; ao_nucl : Nucl_number.t array; - ao_power : Symmetry.Xyz.t array; + ao_power : Angmom.Xyz.t array; ao_coef : AO_coef.t array; ao_expo : AO_expo.t array; ao_cartesian : bool; @@ -87,7 +87,7 @@ end = struct if (data.(2*dim+i-1) > 0) then result.(i-1) <- result.(i-1)^"z"^(string_of_int data.(2*dim+i-1)); done; - Array.map Symmetry.Xyz.of_string result + Array.map Angmom.Xyz.of_string result ;; let read_ao_coef () = @@ -133,7 +133,7 @@ end = struct let ao_num = AO_number.to_int b.ao_num in let gto_array = Array.init (AO_number.to_int b.ao_num) (fun i -> - let s = Symmetry.Xyz.to_symmetry b.ao_power.(i) in + let s = Angmom.Xyz.to_symmetry b.ao_power.(i) in let ao_prim_num = AO_prim_number.to_int b.ao_prim_num.(i) in let prims = List.init ao_prim_num (fun j -> let prim = { GaussianPrimitive.sym = s ; @@ -217,9 +217,9 @@ end = struct let ao_power = let l = Array.to_list ao_power in List.concat [ - (list_map (fun a -> Positive_int.to_int a.Symmetry.Xyz.x) l) ; - (list_map (fun a -> Positive_int.to_int a.Symmetry.Xyz.y) l) ; - (list_map (fun a -> Positive_int.to_int a.Symmetry.Xyz.z) l) ] + (list_map (fun a -> Positive_int.to_int a.Angmom.Xyz.x) l) ; + (list_map (fun a -> Positive_int.to_int a.Angmom.Xyz.y) l) ; + (list_map (fun a -> Positive_int.to_int a.Angmom.Xyz.z) l) ] in Ezfio.set_ao_basis_ao_power(Ezfio.ezfio_array_of_list ~rank:2 ~dim:[| ao_num ; 3 |] ~data:ao_power) ; @@ -409,7 +409,7 @@ end = struct | [] -> [] | (i,n,x)::tail -> (Printf.sprintf " %5d %6d %-8s\n" i (Nucl_number.to_int n) - (Symmetry.Xyz.to_string x) + (Angmom.Xyz.to_string x) )::(do_work tail) in do_work l |> String.concat "" @@ -496,7 +496,7 @@ md5 = %s (b.ao_nucl |> Array.to_list |> list_map Nucl_number.to_string |> String.concat ", ") (b.ao_power |> Array.to_list |> list_map (fun x-> - "("^(Symmetry.Xyz.to_string x)^")" )|> String.concat ", ") + "("^(Angmom.Xyz.to_string x)^")" )|> String.concat ", ") (b.ao_coef |> Array.to_list |> list_map AO_coef.to_string |> String.concat ", ") (b.ao_expo |> Array.to_list |> list_map AO_expo.to_string diff --git a/ocaml/Long_basis.ml b/ocaml/Long_basis.ml index 6c88954a..a8ea3c66 100644 --- a/ocaml/Long_basis.ml +++ b/ocaml/Long_basis.ml @@ -2,7 +2,7 @@ open Qptypes open Qputils open Sexplib.Std -type t = (Symmetry.Xyz.t * Gto.t * Nucl_number.t ) list [@@deriving sexp] +type t = (Angmom.Xyz.t * Gto.t * Nucl_number.t ) list [@@deriving sexp] let of_basis b = let rec do_work accu = function @@ -10,7 +10,7 @@ let of_basis b = | (g,n)::tail -> begin let new_accu = - Symmetry.Xyz.of_symmetry g.Gto.sym + Angmom.Xyz.of_symmetry g.Gto.sym |> List.rev_map (fun x-> (x,g,n)) in do_work (new_accu@accu) tail @@ -25,7 +25,7 @@ let to_basis b = | [] -> List.rev accu | (s,g,n)::tail -> let first_sym = - Symmetry.Xyz.of_symmetry g.Gto.sym + Angmom.Xyz.of_symmetry g.Gto.sym |> List.hd in let new_accu = @@ -42,7 +42,7 @@ let to_basis b = let to_string b = let middle = list_map (fun (x,y,z) -> "( "^((string_of_int (Nucl_number.to_int z)))^", "^ - (Symmetry.Xyz.to_string x)^", "^(Gto.to_string y) + (Angmom.Xyz.to_string x)^", "^(Gto.to_string y) ^" )" ) b |> String.concat ",\n" diff --git a/ocaml/Long_basis.mli b/ocaml/Long_basis.mli index 26009c10..2c41ad74 100644 --- a/ocaml/Long_basis.mli +++ b/ocaml/Long_basis.mli @@ -5,16 +5,16 @@ open Qptypes;; * all the D orbitals are converted to xx, xy, xz, yy, yx * etc *) -type t = (Symmetry.Xyz.t * Gto.t * Nucl_number.t) list [@@deriving sexp] +type t = (Angmom.Xyz.t * Gto.t * Nucl_number.t) list [@@deriving sexp] (** Transform a basis to a long basis *) val of_basis : - (Gto.t * Nucl_number.t) list -> (Symmetry.Xyz.t * Gto.t * Nucl_number.t) list + (Gto.t * Nucl_number.t) list -> (Angmom.Xyz.t * Gto.t * Nucl_number.t) list (** Transform a long basis to a basis *) val to_basis : - (Symmetry.Xyz.t * Gto.t * Nucl_number.t) list -> (Gto.t * Nucl_number.t) list + (Angmom.Xyz.t * Gto.t * Nucl_number.t) list -> (Gto.t * Nucl_number.t) list (** Convert the basis into its string representation *) val to_string : - (Symmetry.Xyz.t * Gto.t * Nucl_number.t) list -> string + (Angmom.Xyz.t * Gto.t * Nucl_number.t) list -> string diff --git a/ocaml/Primitive.mli b/ocaml/Primitive.mli index f7d8809d..f63ecd8c 100644 --- a/ocaml/Primitive.mli +++ b/ocaml/Primitive.mli @@ -1,5 +1,5 @@ type t = -{ sym : Symmetry.t; +{ sym : Angmom.t; expo : Qptypes.AO_expo.t; } [@@deriving sexp] @@ -7,5 +7,5 @@ type t = val to_string : t -> string (** Creation *) -val of_sym_expo : Symmetry.t -> Qptypes.AO_expo.t -> t +val of_sym_expo : Angmom.t -> Qptypes.AO_expo.t -> t diff --git a/ocaml/qp_create_ezfio.ml b/ocaml/qp_create_ezfio.ml index f4c7f3b6..a4865e2b 100644 --- a/ocaml/qp_create_ezfio.ml +++ b/ocaml/qp_create_ezfio.ml @@ -512,13 +512,14 @@ let run ?o b au c d m p cart xyz_file = let ao_num = List.length long_basis in Ezfio.set_ao_basis_ao_num ao_num; Ezfio.set_ao_basis_ao_basis b; + Ezfio.set_basis_basis b; let ao_prim_num = list_map (fun (_,g,_) -> List.length g.Gto.lc) long_basis and ao_nucl = list_map (fun (_,_,n) -> Nucl_number.to_int n) long_basis and ao_power= let l = list_map (fun (x,_,_) -> x) long_basis in - (list_map (fun t -> Positive_int.to_int Symmetry.Xyz.(t.x)) l)@ - (list_map (fun t -> Positive_int.to_int Symmetry.Xyz.(t.y)) l)@ - (list_map (fun t -> Positive_int.to_int Symmetry.Xyz.(t.z)) l) + (list_map (fun t -> Positive_int.to_int Angmom.Xyz.(t.x)) l)@ + (list_map (fun t -> Positive_int.to_int Angmom.Xyz.(t.y)) l)@ + (list_map (fun t -> Positive_int.to_int Angmom.Xyz.(t.z)) l) in let ao_prim_num_max = List.fold_left (fun s x -> if x > s then x @@ -561,6 +562,78 @@ let run ?o b au c d m p cart xyz_file = and ao_expo = create_expo_coef `Expos in let () = + let shell_num = List.length basis in + let lc : (GaussianPrimitive.t * Qptypes.AO_coef.t) list list = + list_map ( fun (g,_) -> g.Gto.lc ) basis + in + let ang_mom = + list_map (fun (l : (GaussianPrimitive.t * Qptypes.AO_coef.t) list) -> + let x, _ = List.hd l in + Angmom.to_l x.GaussianPrimitive.sym |> Qptypes.Positive_int.to_int + ) lc + in + let expo = + list_map (fun l -> list_map (fun (x,_) -> Qptypes.AO_expo.to_float x.GaussianPrimitive.expo) l ) lc + |> List.concat + in + let coef = + list_map (fun l -> + list_map (fun (_,x) -> Qptypes.AO_coef.to_float x) l + ) lc + |> List.concat + in + let shell_prim_num = + list_map List.length lc + in + let shell_prim_idx = + let rec aux count accu = function + | [] -> List.rev accu + | l::rest -> + let newcount = count+(List.length l) in + aux newcount (count::accu) rest + in + aux 1 [] lc + in + let prim_num = List.length coef in + Ezfio.set_basis_typ "Gaussian"; + Ezfio.set_basis_shell_num shell_num; + Ezfio.set_basis_prim_num prim_num ; + Ezfio.set_basis_shell_prim_num (Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| shell_num |] ~data:shell_prim_num); + Ezfio.set_basis_shell_ang_mom (Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| shell_num |] ~data:ang_mom ) ; + Ezfio.set_basis_shell_prim_index (Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| shell_num |] ~data:shell_prim_idx) ; + Ezfio.set_basis_basis_nucleus_index (Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| nucl_num |] + ~data:( + list_map (fun (_,n) -> Nucl_number.to_int n) basis + |> List.fold_left (fun accu i -> + match accu with + | [] -> [] + | (h,j) :: rest -> if j == i then ((h+1,j)::rest) else ((h+1,i)::(h+1,j)::rest) + ) [(0,0)] + |> List.rev + |> List.map fst + )) ; + Ezfio.set_basis_nucleus_shell_num(Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| nucl_num |] + ~data:( + list_map (fun (_,n) -> Nucl_number.to_int n) basis + |> List.fold_left (fun accu i -> + match accu with + | [] -> [(1,i)] + | (h,j) :: rest -> if j == i then ((h+1,j)::rest) else ((1,i)::(h,j)::rest) + ) [] + |> List.rev + |> List.map fst + )) ; + Ezfio.set_basis_prim_coef (Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| prim_num |] ~data:coef) ; + Ezfio.set_basis_prim_expo (Ezfio.ezfio_array_of_list + ~rank:1 ~dim:[| prim_num |] ~data:expo) ; + + Ezfio.set_ao_basis_ao_prim_num (Ezfio.ezfio_array_of_list ~rank:1 ~dim:[| ao_num |] ~data:ao_prim_num) ; Ezfio.set_ao_basis_ao_nucl(Ezfio.ezfio_array_of_list diff --git a/src/ao_basis/NEED b/src/ao_basis/NEED index d2066b18..a69fe768 100644 --- a/src/ao_basis/NEED +++ b/src/ao_basis/NEED @@ -1 +1,2 @@ nuclei +basis diff --git a/src/ao_basis/aos.irp.f b/src/ao_basis/aos.irp.f index e95ac711..17aa784d 100644 --- a/src/ao_basis/aos.irp.f +++ b/src/ao_basis/aos.irp.f @@ -6,6 +6,23 @@ BEGIN_PROVIDER [ integer, ao_prim_num_max ] ao_prim_num_max = maxval(ao_prim_num) END_PROVIDER +BEGIN_PROVIDER [ integer, ao_shell, (ao_num) ] + implicit none + BEGIN_DOC + ! Index of the shell to which the AO corresponds + END_DOC + integer :: i, j, k, n + k=0 + do i=1,shell_num + n = shell_ang_mom(i)+1 + do j=1,(n*(n+1))/2 + k = k+1 + ao_shell(k) = i + enddo + enddo + +END_PROVIDER + BEGIN_PROVIDER [ double precision, ao_coef_normalized, (ao_num,ao_prim_num_max) ] &BEGIN_PROVIDER [ double precision, ao_coef_normalization_factor, (ao_num) ] implicit none @@ -23,14 +40,15 @@ END_PROVIDER do i=1,ao_num - powA(1) = ao_power(i,1) - powA(2) = ao_power(i,2) - powA(3) = ao_power(i,3) + powA(1) = ao_power(i,1) + ao_power(i,2) + ao_power(i,3) + powA(2) = 0 + powA(3) = 0 ! Normalization of the primitives if (primitives_normalized) then do j=1,ao_prim_num(i) - call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j),powA,powA,overlap_x,overlap_y,overlap_z,norm,nz) + call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,j), & + powA,powA,overlap_x,overlap_y,overlap_z,norm,nz) ao_coef_normalized(i,j) = ao_coef(i,j)/sqrt(norm) enddo else @@ -39,6 +57,10 @@ END_PROVIDER enddo endif + powA(1) = ao_power(i,1) + powA(2) = ao_power(i,2) + powA(3) = ao_power(i,3) + ! Normalization of the contracted basis functions if (ao_normalized) then norm = 0.d0 @@ -56,39 +78,6 @@ END_PROVIDER END_PROVIDER -BEGIN_PROVIDER [ double precision, ao_coef_normalization_libint_factor, (ao_num) ] - implicit none - BEGIN_DOC - ! |AO| normalization for interfacing with libint - END_DOC - double precision :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c - integer :: l, powA(3), nz - integer :: i,j,k - nz=100 - C_A(1) = 0.d0 - C_A(2) = 0.d0 - C_A(3) = 0.d0 - - do i=1,ao_num - powA(1) = ao_l(i) - powA(2) = 0 - powA(3) = 0 - - ! Normalization of the contracted basis functions - norm = 0.d0 - do j=1,ao_prim_num(i) - do k=1,ao_prim_num(i) - call overlap_gaussian_xyz(C_A,C_A,ao_expo(i,j),ao_expo(i,k),powA,powA,overlap_x,overlap_y,overlap_z,c,nz) - norm = norm+c*ao_coef_normalized(i,j)*ao_coef_normalized(i,k) - enddo - enddo - ao_coef_normalization_libint_factor(i) = ao_coef_normalization_factor(i) * sqrt(norm) - - enddo - -END_PROVIDER - - BEGIN_PROVIDER [ double precision, ao_coef_normalized_ordered, (ao_num,ao_prim_num_max) ] &BEGIN_PROVIDER [ double precision, ao_expo_ordered, (ao_num,ao_prim_num_max) ] implicit none @@ -231,38 +220,11 @@ END_PROVIDER END_DOC do i = 1, nucl_num Nucl_num_shell_Aos(i) = 0 - do j = 1, Nucl_N_Aos(i) - if(ao_l(Nucl_Aos(i,j))==0)then - ! S type function - Nucl_num_shell_Aos(i)+=1 - Nucl_list_shell_Aos(i,Nucl_num_shell_Aos(i))=Nucl_Aos(i,j) - elseif(ao_l(Nucl_Aos(i,j))==1)then - ! P type function - if(ao_power(Nucl_Aos(i,j),1)==1)then + if (ao_power(Nucl_Aos(i,j),1) == ao_l(Nucl_Aos(i,j))) then Nucl_num_shell_Aos(i)+=1 Nucl_list_shell_Aos(i,Nucl_num_shell_Aos(i))=Nucl_Aos(i,j) endif - elseif(ao_l(Nucl_Aos(i,j))==2)then - ! D type function - if(ao_power(Nucl_Aos(i,j),1)==2)then - Nucl_num_shell_Aos(i)+=1 - Nucl_list_shell_Aos(i,Nucl_num_shell_Aos(i))=Nucl_Aos(i,j) - endif - elseif(ao_l(Nucl_Aos(i,j))==3)then - ! F type function - if(ao_power(Nucl_Aos(i,j),1)==3)then - Nucl_num_shell_Aos(i)+=1 - Nucl_list_shell_Aos(i,Nucl_num_shell_Aos(i))=Nucl_Aos(i,j) - endif - elseif(ao_l(Nucl_Aos(i,j))==4)then - ! G type function - if(ao_power(Nucl_Aos(i,j),1)==4)then - Nucl_num_shell_Aos(i)+=1 - Nucl_list_shell_Aos(i,Nucl_num_shell_Aos(i))=Nucl_Aos(i,j) - endif - endif - enddo enddo diff --git a/src/ao_one_e_ints/pseudopot.f90 b/src/ao_one_e_ints/pseudopot.f90 index 6c8e5c83..48e3803e 100644 --- a/src/ao_one_e_ints/pseudopot.f90 +++ b/src/ao_one_e_ints/pseudopot.f90 @@ -96,8 +96,12 @@ end ! x=cos(theta) double precision function ylm_real(l,m,x,phi) - implicit double precision (a-h,o-z) - DIMENSION PM(0:100,0:100) + implicit none + integer :: MM, iabs_m, m, l + double precision :: pi, fourpi, factor, x, phi, coef + double precision :: xchap, ychap, zchap + double precision, external :: fact + double precision :: PM(0:100,0:100), plm MM=100 pi=dacos(-1.d0) fourpi=4.d0*pi @@ -1150,8 +1154,10 @@ end ! Output: PM(m,n) --- Pmn(x) ! ===================================================== ! - IMPLICIT DOUBLE PRECISION (P,X) - DIMENSION PM(0:MM,0:(N+1)) + implicit none +! IMPLICIT DOUBLE PRECISION (P,X) + integer :: MM, N, I, J, M + double precision :: PM(0:MM,0:(N+1)), X, XQ, XS DOUBLE PRECISION, SAVE :: INVERSE(100) = 0.D0 DOUBLE PRECISION :: LS, II, JJ IF (INVERSE(1) == 0.d0) THEN @@ -1202,8 +1208,9 @@ end ! P_l^|m|(cos(theta)) exp(i m phi) subroutine erreur(x,n,rmoy,error) - implicit double precision(a-h,o-z) - dimension x(n) + implicit none + integer :: i, n + double precision :: x(n), rn, rn1, error, rmoy ! calcul de la moyenne rmoy=0.d0 do i=1,n diff --git a/src/basis/EZFIO.cfg b/src/basis/EZFIO.cfg new file mode 100644 index 00000000..7f2ede4c --- /dev/null +++ b/src/basis/EZFIO.cfg @@ -0,0 +1,75 @@ +[basis] +type: character*(256) +doc: Name of the |AO| basis set +interface: ezfio, provider + +[typ] +type: character*(32) +doc: Type of basis set. Only 'Gaussian' is supported +interface: ezfio, provider + +[shell_num] +type: integer +doc: Number of shells +interface: ezfio, provider + +[nucleus_shell_num] +type: integer +doc: Number of shells per nucleus +size: (nuclei.nucl_num) +interface: ezfio, provider + +[shell_normalization_factor] +type: double precision +doc: Normalization factor applied to the whole shell, ex $1/\sqrt{ }$ +size: (basis.shell_num) +interface: ezfio + +[shell_ang_mom] +type: integer +doc: Angular momentum of each shell +size: (basis.shell_num) +interface: ezfio, provider + +[shell_prim_num] +type: integer +doc: Number of primitives in a shell +size: (basis.shell_num) +interface: ezfio, provider + +[shell_prim_index] +type: integer +doc: Max number of primitives in a shell +size: (basis.shell_num) +interface: ezfio, provider + +[basis_nucleus_index] +type: integer +doc: Index of the nucleus on which the shell is centered +size: (nuclei.nucl_num) +interface: ezfio, provider + +[prim_normalization_factor] +type: double precision +doc: Normalization factor applied to each primitive +size: (basis.prim_num) +interface: ezfio + +[prim_num] +type: integer +doc: Total number of primitives +interface: ezfio, provider + +[prim_coef] +type: double precision +doc: Primitive coefficients +size: (basis.prim_num) +interface: ezfio, provider + +[prim_expo] +type: double precision +doc: Exponents in the shell +size: (basis.prim_num) +interface: ezfio, provider + + diff --git a/src/basis/NEED b/src/basis/NEED new file mode 100644 index 00000000..d2066b18 --- /dev/null +++ b/src/basis/NEED @@ -0,0 +1 @@ +nuclei diff --git a/src/basis/README.rst b/src/basis/README.rst new file mode 100644 index 00000000..17493b4f --- /dev/null +++ b/src/basis/README.rst @@ -0,0 +1,8 @@ +====== +basis +====== + +This module contains the basis set information, which will then be used to build the atomic orbitals. + + + diff --git a/src/basis/basis.irp.f b/src/basis/basis.irp.f new file mode 100644 index 00000000..6a406e28 --- /dev/null +++ b/src/basis/basis.irp.f @@ -0,0 +1,119 @@ +BEGIN_PROVIDER [ double precision, shell_normalization_factor , (shell_num) ] + implicit none + BEGIN_DOC + ! Number of primitives per |AO| + END_DOC + + logical :: has + PROVIDE ezfio_filename + if (mpi_master) then + if (size(shell_normalization_factor) == 0) return + + call ezfio_has_basis_shell_normalization_factor(has) + if (has) then + write(6,'(A)') '.. >>>>> [ IO READ: shell_normalization_factor ] <<<<< ..' + call ezfio_get_basis_shell_normalization_factor(shell_normalization_factor) + else + + double precision :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c + integer :: l, powA(3), nz + integer :: i,j,k + nz=100 + C_A(1) = 0.d0 + C_A(2) = 0.d0 + C_A(3) = 0.d0 + + do i=1,shell_num + + powA(1) = shell_ang_mom(i) + powA(2) = 0 + powA(3) = 0 + + norm = 0.d0 + do k=shell_prim_index(i),shell_prim_index(i)+shell_prim_num(i)-1 + do j=shell_prim_index(i),shell_prim_index(i)+shell_prim_num(i)-1 + call overlap_gaussian_xyz(C_A,C_A,prim_expo(j),prim_expo(k), & + powA,powA,overlap_x,overlap_y,overlap_z,c,nz) + norm = norm+c*prim_coef(j)*prim_coef(k) * prim_normalization_factor(j) * prim_normalization_factor(k) + enddo + enddo + shell_normalization_factor(i) = 1.d0/dsqrt(norm) + enddo + + endif + endif + IRP_IF MPI_DEBUG + print *, irp_here, mpi_rank + call MPI_BARRIER(MPI_COMM_WORLD, ierr) + IRP_ENDIF + IRP_IF MPI + include 'mpif.h' + integer :: ierr + call MPI_BCAST( shell_normalization_factor, (shell_num), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr) + if (ierr /= MPI_SUCCESS) then + stop 'Unable to read shell_normalization_factor with MPI' + endif + IRP_ENDIF + + call write_time(6) + +END_PROVIDER + + +BEGIN_PROVIDER [ double precision, prim_normalization_factor , (prim_num) ] + implicit none + BEGIN_DOC + ! Number of primitives per |AO| + END_DOC + + logical :: has + PROVIDE ezfio_filename + if (mpi_master) then + if (size(prim_normalization_factor) == 0) return + + call ezfio_has_basis_prim_normalization_factor(has) + if (has) then + write(6,'(A)') '.. >>>>> [ IO READ: prim_normalization_factor ] <<<<< ..' + call ezfio_get_basis_prim_normalization_factor(prim_normalization_factor) + else + + double precision :: norm,overlap_x,overlap_y,overlap_z,C_A(3), c + integer :: l, powA(3), nz + integer :: i,j,k + nz=100 + C_A(1) = 0.d0 + C_A(2) = 0.d0 + C_A(3) = 0.d0 + + do i=1,shell_num + + powA(1) = shell_ang_mom(i) + powA(2) = 0 + powA(3) = 0 + + do k=shell_prim_index(i),shell_prim_index(i)+shell_prim_num(i)-1 + call overlap_gaussian_xyz(C_A,C_A,prim_expo(k),prim_expo(k), & + powA,powA,overlap_x,overlap_y,overlap_z,norm,nz) + prim_normalization_factor(k) = 1.d0/dsqrt(norm) + enddo + enddo + + + endif + endif + IRP_IF MPI_DEBUG + print *, irp_here, mpi_rank + call MPI_BARRIER(MPI_COMM_WORLD, ierr) + IRP_ENDIF + IRP_IF MPI + include 'mpif.h' + integer :: ierr + call MPI_BCAST( prim_normalization_factor, (prim_num), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr) + if (ierr /= MPI_SUCCESS) then + stop 'Unable to read prim_normalization_factor with MPI' + endif + IRP_ENDIF + + call write_time(6) + +END_PROVIDER diff --git a/src/bitmask/bitmasks.irp.f b/src/bitmask/bitmasks.irp.f index 91617397..27c7dc84 100644 --- a/src/bitmask/bitmasks.irp.f +++ b/src/bitmask/bitmasks.irp.f @@ -27,9 +27,7 @@ BEGIN_PROVIDER [ integer(bit_kind), full_ijkl_bitmask, (N_int) ] full_ijkl_bitmask(j) = 0_bit_kind do i=0,bit_kind_size-1 k=k+1 - if (mo_class(k) /= 'Deleted') then - full_ijkl_bitmask(j) = ibset(full_ijkl_bitmask(j),i) - endif + full_ijkl_bitmask(j) = ibset(full_ijkl_bitmask(j),i) if (k == mo_num) exit enddo enddo diff --git a/src/cipsi/pt2_stoch_routines.irp.f b/src/cipsi/pt2_stoch_routines.irp.f index 7554c39e..6f0d6683 100644 --- a/src/cipsi/pt2_stoch_routines.irp.f +++ b/src/cipsi/pt2_stoch_routines.irp.f @@ -286,7 +286,7 @@ subroutine ZMQ_pt2(E, pt2_data, pt2_data_err, relative_error, N_in) call write_int(6,nproc_target,'Number of threads for PT2') call write_double(6,mem,'Memory (Gb)') - call omp_set_nested(.false.) + call omp_set_max_active_levels(1) print '(A)', '========== ======================= ===================== ===================== ===========' @@ -313,6 +313,7 @@ subroutine ZMQ_pt2(E, pt2_data, pt2_data_err, relative_error, N_in) endif !$OMP END PARALLEL call end_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'pt2') + call omp_set_max_active_levels(8) print '(A)', '========== ======================= ===================== ===================== ===========' diff --git a/src/cipsi/selection.irp.f b/src/cipsi/selection.irp.f index 59792ea9..eda9642c 100644 --- a/src/cipsi/selection.irp.f +++ b/src/cipsi/selection.irp.f @@ -253,12 +253,7 @@ subroutine select_singles_and_doubles(i_generator,hole_mask,particle_mask,fock_d deallocate(exc_degree) nmax=k-1 - allocate(iorder(nmax)) - do i=1,nmax - iorder(i) = i - enddo - call isort(indices,iorder,nmax) - deallocate(iorder) + call isort_noidx(indices,nmax) ! Start with 32 elements. Size will double along with the filtering. allocate(preinteresting(0:32), prefullinteresting(0:32), & @@ -749,7 +744,7 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d double precision :: eigvalues(N_states+1) double precision :: work(1+6*(N_states+1)+2*(N_states+1)**2) - integer :: iwork(3+5*(N_states+1)), info, k + integer :: info, k , iwork(N_states+1) if (do_diag) then double precision :: pt2_matrix(N_states+1,N_states+1) @@ -761,8 +756,8 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d pt2_matrix(N_states+1,istate) = mat(istate,p1,p2) enddo - call DSYEVD( 'V', 'U', N_states+1, pt2_matrix, N_states+1, eigvalues, & - work, size(work), iwork, size(iwork), info ) + call DSYEV( 'V', 'U', N_states+1, pt2_matrix, N_states+1, eigvalues, & + work, size(work), info ) if (info /= 0) then print *, 'error in '//irp_here stop -1 @@ -770,7 +765,7 @@ subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_d pt2_matrix = dabs(pt2_matrix) iwork(1:N_states+1) = maxloc(pt2_matrix,DIM=1) do k=1,N_states - e_pert(iwork(k)) = eigvalues(k) - E0(iwork(k)) + e_pert(k) = eigvalues(iwork(k)) - E0(k) enddo endif diff --git a/src/csf/sigma_vector.irp.f b/src/csf/sigma_vector.irp.f index 265d2384..85ed5f84 100644 --- a/src/csf/sigma_vector.irp.f +++ b/src/csf/sigma_vector.irp.f @@ -1,3 +1,9 @@ + real*8 function logabsgamma(x) + implicit none + real*8, intent(in) :: x + logabsgamma = log(abs(gamma(x))) + end function logabsgamma + BEGIN_PROVIDER [ integer, NSOMOMax] &BEGIN_PROVIDER [ integer, NCSFMax] &BEGIN_PROVIDER [ integer*8, NMO] @@ -22,33 +28,65 @@ integer NSOMO integer dimcsfpercfg integer detDimperBF - real*8 :: coeff + real*8 :: coeff, binom1, binom2 integer MS integer ncfgpersomo + real*8, external :: logabsgamma detDimperBF = 0 MS = elec_alpha_num-elec_beta_num ! number of cfgs = number of dets for 0 somos - n_CSF = cfg_seniority_index(0)-1 + n_CSF = 0 ncfgprev = cfg_seniority_index(0) - do i = 0-iand(MS,1)+2, NSOMOMax,2 - if(cfg_seniority_index(i) .EQ. -1)then - ncfgpersomo = N_configuration + 1 - else - ncfgpersomo = cfg_seniority_index(i) - endif - ncfg = ncfgpersomo - ncfgprev - !detDimperBF = max(1,nint((binom(i,(i+1)/2)))) - if (i > 2) then - dimcsfpercfg = max(1,nint((binom(i-2,(i-2+1)/2)-binom(i-2,((i-2+1)/2)+1)))) - else - dimcsfpercfg = 1 - endif - n_CSF += ncfg * dimcsfpercfg - !if(cfg_seniority_index(i+2) == -1) EXIT - !if(detDimperBF > maxDetDimPerBF) maxDetDimPerBF = detDimperBF - ncfgprev = cfg_seniority_index(i) + ncfgpersomo = ncfgprev + do i = iand(MS,1), NSOMOMax-2,2 + if(cfg_seniority_index(i) .EQ. -1) then + cycle + endif + if(cfg_seniority_index(i+2) .EQ. -1) then + ncfgpersomo = N_configuration + 1 + else + if(cfg_seniority_index(i+2) > ncfgpersomo) then + ncfgpersomo = cfg_seniority_index(i+2) + else + k = 0 + do while(cfg_seniority_index(i+2+k) < ncfgpersomo) + k = k + 2 + ncfgpersomo = cfg_seniority_index(i+2+k) + enddo + endif + endif + ncfg = ncfgpersomo - ncfgprev + if(iand(MS,1) .EQ. 0) then + !dimcsfpercfg = max(1,nint((binom(i,i/2)-binom(i,i/2+1)))) + binom1 = dexp(logabsgamma(1.0d0*(i+1)) & + - logabsgamma(1.0d0*((i/2)+1)) & + - logabsgamma(1.0d0*(i-((i/2))+1))); + binom2 = dexp(logabsgamma(1.0d0*(i+1)) & + - logabsgamma(1.0d0*(((i/2)+1)+1)) & + - logabsgamma(1.0d0*(i-((i/2)+1)+1))); + dimcsfpercfg = max(1,nint(binom1 - binom2)) + else + !dimcsfpercfg = max(1,nint((binom(i,(i+1)/2)-binom(i,(i+3)/2)))) + binom1 = dexp(logabsgamma(1.0d0*(i+1)) & + - logabsgamma(1.0d0*(((i+1)/2)+1)) & + - logabsgamma(1.0d0*(i-(((i+1)/2))+1))); + binom2 = dexp(logabsgamma(1.0d0*(i+1)) & + - logabsgamma(1.0d0*((((i+3)/2)+1)+1)) & + - logabsgamma(1.0d0*(i-(((i+3)/2)+1)+1))); + dimcsfpercfg = max(1,nint(binom1 - binom2)) + endif + n_CSF += ncfg * dimcsfpercfg + if(cfg_seniority_index(i+2) > ncfgprev) then + ncfgprev = cfg_seniority_index(i+2) + else + k = 0 + do while(cfg_seniority_index(i+2+k) < ncfgprev) + k = k + 2 + ncfgprev = cfg_seniority_index(i+2+k) + enddo + endif enddo - END_PROVIDER +END_PROVIDER subroutine get_phase_qp_to_cfg(Ialpha, Ibeta, phaseout) diff --git a/src/davidson/diagonalization_hcsf_dressed.irp.f b/src/davidson/diagonalization_hcsf_dressed.irp.f index 2a83cc28..da23b919 100644 --- a/src/davidson/diagonalization_hcsf_dressed.irp.f +++ b/src/davidson/diagonalization_hcsf_dressed.irp.f @@ -197,6 +197,7 @@ subroutine davidson_diag_csf_hjj(dets_in,u_in,H_jj,energies,dim_in,sze,sze_csf,N call write_int(6,N_st,'Number of states') call write_int(6,N_st_diag,'Number of states in diagonalization') call write_int(6,sze,'Number of determinants') + call write_int(6,sze_csf,'Number of CSFs') call write_int(6,nproc_target,'Number of threads for diagonalization') call write_double(6, r1, 'Memory(Gb)') if (disk_based) then diff --git a/src/dressing/run_dress_slave.irp.f b/src/dressing/run_dress_slave.irp.f index 8a92962c..a33fb1dd 100644 --- a/src/dressing/run_dress_slave.irp.f +++ b/src/dressing/run_dress_slave.irp.f @@ -72,7 +72,7 @@ subroutine run_dress_slave(thread,iproce,energy) provide psi_energy ending = dress_N_cp+1 ntask_tbd = 0 - call omp_set_nested(.true.) + call omp_set_max_active_levels(8) !$OMP PARALLEL DEFAULT(SHARED) & !$OMP PRIVATE(interesting, breve_delta_m, task_id) & @@ -84,7 +84,7 @@ subroutine run_dress_slave(thread,iproce,energy) zmq_socket_push = new_zmq_push_socket(thread) integer, external :: connect_to_taskserver !$OMP CRITICAL - call omp_set_nested(.false.) + call omp_set_max_active_levels(1) if (connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread) == -1) then print *, irp_here, ': Unable to connect to task server' stop -1 @@ -296,7 +296,7 @@ subroutine run_dress_slave(thread,iproce,energy) !$OMP END CRITICAL !$OMP END PARALLEL - call omp_set_nested(.false.) + call omp_set_max_active_levels(1) ! do i=0,dress_N_cp+1 ! call omp_destroy_lock(lck_sto(i)) ! end do diff --git a/src/mo_one_e_ints/EZFIO.cfg b/src/mo_one_e_ints/EZFIO.cfg index ca4dabf4..d58b3da1 100644 --- a/src/mo_one_e_ints/EZFIO.cfg +++ b/src/mo_one_e_ints/EZFIO.cfg @@ -47,10 +47,3 @@ type: Disk_access doc: Read/Write |MO| one-electron integrals from/to disk [ Write | Read | None ] interface: ezfio,provider,ocaml default: None - - -[restore_symm] -type: logical -doc: If true, try to find symmetry in the MO coefficient matrices -interface: ezfio,provider,ocaml -default: True diff --git a/src/two_body_rdm/two_e_dm_mo.irp.f b/src/two_body_rdm/two_e_dm_mo.irp.f index 19e8d75e..4dadd2e6 100644 --- a/src/two_body_rdm/two_e_dm_mo.irp.f +++ b/src/two_body_rdm/two_e_dm_mo.irp.f @@ -5,18 +5,14 @@ BEGIN_PROVIDER [double precision, two_e_dm_mo, (mo_num,mo_num,mo_num,mo_num)] ! ! ! - ! WHERE ALL ORBITALS (i,j,k,l) BELONGS TO ALL OCCUPIED ORBITALS : core, inactive and active + ! where the indices (i,j,k,l) belong to all MOs. ! - ! THE NORMALIZATION (i.e. sum of diagonal elements) IS SET TO N_{elec} * (N_{elec} - 1)/2 + ! The normalization (i.e. sum of diagonal elements) is set to $N_{elec} * (N_{elec} - 1)/2$ ! - ! !!!!! WARNING !!!!! ALL SLATER DETERMINANTS IN PSI_DET MUST BELONG TO AN ACTIVE SPACE DEFINED BY "list_act" + ! !!!!! WARNING !!!!! IF "no_core_density" then all elements involving at least one CORE MO are set to zero + ! The state-averaged two-electron energy : ! - ! !!!!! WARNING !!!!! IF "no_core_density" then all elements involving at least one CORE MO is set to zero - ! The two-electron energy of each state can be computed as: - ! - ! \sum_{i,j,k,l = 1, n_core_inact_act_orb} two_e_dm_mo(i,j,k,l,istate) * < ii jj | kk ll > - ! - ! with ii = list_core_inact_act(i), jj = list_core_inact_act(j), kk = list_core_inact_act(k), ll = list_core_inact_act(l) + ! \sum_{i,j,k,l = 1, mo_num} two_e_dm_mo(i,j,k,l) * < ii jj | kk ll > END_DOC two_e_dm_mo = 0.d0 integer :: i,j,k,l,iorb,jorb,korb,lorb,istate diff --git a/src/utils/EZFIO.cfg b/src/utils/EZFIO.cfg new file mode 100644 index 00000000..7d367f0c --- /dev/null +++ b/src/utils/EZFIO.cfg @@ -0,0 +1,5 @@ +[restore_symm] +type: logical +doc: If true, try to find symmetry in the MO coefficient matrices +interface: ezfio,provider,ocaml +default: False diff --git a/src/utils/intel.f90 b/src/utils/intel.f90 new file mode 100644 index 00000000..2b61cb38 --- /dev/null +++ b/src/utils/intel.f90 @@ -0,0 +1,173 @@ +module intel + use, intrinsic :: iso_c_binding + interface + subroutine ippsSortAscend_32s_I(pSrc, len) bind(C, name='ippsSortAscend_32s_I') + use iso_c_binding + integer, intent(in), value :: len + integer, intent(inout) :: pSrc(len) + end + end interface + interface + subroutine ippsSortAscend_32f_I(pSrc, len) bind(C, name='ippsSortAscend_32f_I') + use iso_c_binding + integer, intent(in), value :: len + real, intent(inout) :: pSrc(len) + end + end interface + interface + subroutine ippsSortAscend_64s_I(pSrc, len) bind(C, name='ippsSortAscend_64s_I') + use iso_c_binding + integer, intent(in), value :: len + integer*8, intent(inout) :: pSrc(len) + end + end interface + interface + subroutine ippsSortAscend_64f_I(pSrc, len) bind(C, name='ippsSortAscend_64f_I') + use iso_c_binding + integer, intent(in), value :: len + double precision, intent(inout) :: pSrc(len) + end + end interface + + interface + subroutine ippsSortRadixIndexGetBufferSize(len, dataType, pBufSize) bind(C, name='ippsSortRadixIndexGetBufferSize') + use iso_c_binding + integer, intent(in), value :: len + integer, intent(in), value :: dataType + integer, intent(out) :: pBufSize + end + end interface + + interface + subroutine ippsSortRadixAscend_16s_I(pSrc, len, pTmp) bind(C, name='ippsSortRadixAscend_16s_I') + use iso_c_binding + integer, intent(in), value :: len + integer*2, intent(inout) :: pSrc(len) + character, intent(inout) :: pTmp(len) + end + end interface + interface + subroutine ippsSortRadixAscend_32s_I(pSrc, len, pTmp) bind(C, name='ippsSortRadixAscend_32s_I') + use iso_c_binding + integer, intent(in), value :: len + integer, intent(inout) :: pSrc(len) + character, intent(inout) :: pTmp(len) + end + end interface + interface + subroutine ippsSortRadixAscend_32f_I(pSrc, len, pTmp) bind(C, name='ippsSortRadixAscend_32f_I') + use iso_c_binding + integer, intent(in), value :: len + real, intent(inout) :: pSrc(len) + character, intent(inout) :: pTmp(len) + end + end interface + interface + subroutine ippsSortRadixAscend_64s_I(pSrc, len, pTmp) bind(C, name='ippsSortRadixAscend_64s_I') + use iso_c_binding + integer, intent(in), value :: len + integer*8, intent(inout) :: pSrc(len) + character, intent(inout) :: pTmp(len) + end + end interface + interface + subroutine ippsSortRadixAscend_64f_I(pSrc, len, pTmp) bind(C, name='ippsSortRadixAscend_64f_I') + use iso_c_binding + integer, intent(in), value :: len + double precision, intent(inout) :: pSrc(len) + character, intent(inout) :: pTmp(len) + end + end interface + + interface + subroutine ippsSortRadixIndexAscend_16s(pSrc, srcStrideBytes, pDstIndx, len, pTmpIndx) bind(C, name='ippsSortRadixIndexAscend_16s') + use iso_c_binding + integer, intent(in), value :: len + integer*2, intent(inout) :: pSrc(len) + integer, intent(in), value :: srcStrideBytes + integer, intent(inout) :: pDstIndx(len) + character, intent(inout) :: pTmpIndx(len) + end + end interface + interface + subroutine ippsSortRadixIndexAscend_32s(pSrc, srcStrideBytes, pDstIndx, len, pTmpIndx) bind(C, name='ippsSortRadixIndexAscend_32s') + use iso_c_binding + integer, intent(in), value :: len + integer, intent(inout) :: pSrc(len) + integer, intent(in), value :: srcStrideBytes + integer, intent(inout) :: pDstIndx(len) + character, intent(inout) :: pTmpIndx(len) + end + end interface + interface + subroutine ippsSortRadixIndexAscend_32f(pSrc, srcStrideBytes, pDstIndx, len, pTmpIndx) bind(C,name='ippsSortRadixIndexAscend_32f') + use iso_c_binding + integer, intent(in), value :: len + real , intent(inout) :: pSrc(len) + integer, intent(in), value :: srcStrideBytes + integer, intent(inout) :: pDstIndx(len) + character, intent(inout) :: pTmpIndx(len) + end + end interface + interface + subroutine ippsSortRadixIndexAscend_64s(pSrc, srcStrideBytes, pDstIndx, len, pTmpIndx) bind(C, name='ippsSortRadixIndexAscend_64s') + use iso_c_binding + integer, intent(in), value :: len + integer*8, intent(inout) :: pSrc(len) + integer, intent(in), value :: srcStrideBytes + integer, intent(inout) :: pDstIndx(len) + character, intent(inout) :: pTmpIndx(len) + end + end interface + interface + subroutine ippsSortRadixIndexAscend_64f(pSrc, srcStrideBytes, pDstIndx, len, pTmpIndx) bind(C,name='ippsSortRadixIndexAscend_64f') + use iso_c_binding + integer, intent(in), value :: len + real*8 , intent(inout) :: pSrc(len) + integer, intent(in), value :: srcStrideBytes + integer, intent(inout) :: pDstIndx(len) + character, intent(inout) :: pTmpIndx(len) + end + end interface + + interface + subroutine ippsSortIndexAscend_32f_I(pSrcDst, pDstIndx, len) bind(C,name='ippsSortIndexAscend_32f_I') + use iso_c_binding + real(4), intent(in) :: pSrcDst(*) + integer(4), intent(inout) :: pDstIndx(*) + integer(4), intent(in), value :: len + end + end interface + interface + subroutine ippsSortIndexAscend_32s_I(pSrcDst, pDstIndx, len) bind(C,name='ippsSortIndexAscend_32s_I') + use iso_c_binding + integer(4), intent(in) :: pSrcDst(*) + integer(4), intent(inout) :: pDstIndx(*) + integer(4), intent(in), value :: len + end + end interface + interface + subroutine ippsSortIndexAscend_64f_I(pSrcDst, pDstIndx, len) bind(C,name='ippsSortIndexAscend_64f_I') + use iso_c_binding + real(8), intent(in) :: pSrcDst(*) + integer(4), intent(inout) :: pDstIndx(*) + integer(4), intent(in), value :: len + end + end interface + interface + subroutine ippsSortIndexAscend_64s_I(pSrcDst, pDstIndx, len) bind(C,name='ippsSortIndexAscend_64s_I') + use iso_c_binding + integer(8), intent(in) :: pSrcDst(*) + integer(4), intent(inout) :: pDstIndx(*) + integer(4), intent(in), value :: len + end + end interface + interface + subroutine ippsSortIndexAscend_16s_I(pSrcDst, pDstIndx, len) bind(C,name='ippsSortIndexAscend_16s_I') + use iso_c_binding + integer(2), intent(in) :: pSrcDst(*) + integer(4), intent(inout) :: pDstIndx(*) + integer(4), intent(in), value :: len + end + end interface +end module diff --git a/src/utils/sort.irp.f b/src/utils/sort.irp.f index 2a655eed..a63eb4a3 100644 --- a/src/utils/sort.irp.f +++ b/src/utils/sort.irp.f @@ -38,15 +38,7 @@ BEGIN_TEMPLATE $type,intent(inout) :: x(isize) integer,intent(inout) :: iorder(isize) integer, external :: omp_get_num_threads - if (omp_get_num_threads() == 1) then - !$OMP PARALLEL DEFAULT(SHARED) - !$OMP SINGLE - call rec_$X_quicksort(x,iorder,isize,1,isize,nproc) - !$OMP END SINGLE - !$OMP END PARALLEL - else - call rec_$X_quicksort(x,iorder,isize,1,isize,nproc) - endif + call rec_$X_quicksort(x,iorder,isize,1,isize,nproc) end recursive subroutine rec_$X_quicksort(x, iorder, isize, first, last, level) @@ -57,7 +49,7 @@ BEGIN_TEMPLATE $type :: c, tmp integer :: itmp integer :: i, j - + if(isize<2)return c = x( shiftr(first+last,1) ) @@ -89,16 +81,11 @@ BEGIN_TEMPLATE endif else if (first < i-1) then - !$OMP TASK DEFAULT(SHARED) FIRSTPRIVATE(isize,first,i,level) call rec_$X_quicksort(x, iorder, isize, first, i-1,level/2) - !$OMP END TASK endif if (j+1 < last) then - !$OMP TASK DEFAULT(SHARED) FIRSTPRIVATE(isize,last,j,level) call rec_$X_quicksort(x, iorder, isize, j+1, last,level/2) - !$OMP END TASK endif - !$OMP TASKWAIT endif end @@ -262,7 +249,60 @@ SUBST [ X, type ] i2 ; integer*2 ;; END_TEMPLATE + +!---------------------- INTEL +IRP_IF INTEL + BEGIN_TEMPLATE + subroutine $Xsort(x,iorder,isize) + use intel + implicit none + BEGIN_DOC + ! Sort array x(isize). + ! iorder in input should be (1,2,3,...,isize), and in output + ! contains the new order of the elements. + END_DOC + integer,intent(in) :: isize + $type,intent(inout) :: x(isize) + integer,intent(inout) :: iorder(isize) + integer :: n + character, allocatable :: tmp(:) + if (isize < 2) return + call ippsSortRadixIndexGetBufferSize(isize, $ippsz, n) + allocate(tmp(n)) + call ippsSortRadixIndexAscend_$ityp(x, $n, iorder, isize, tmp) + deallocate(tmp) + iorder(1:isize) = iorder(1:isize)+1 + call $Xset_order(x,iorder,isize) + end + + subroutine $Xsort_noidx(x,isize) + use intel + implicit none + BEGIN_DOC + ! Sort array x(isize). + ! iorder in input should be (1,2,3,...,isize), and in output + ! contains the new order of the elements. + END_DOC + integer,intent(in) :: isize + $type,intent(inout) :: x(isize) + integer :: n + character, allocatable :: tmp(:) + if (isize < 2) return + call ippsSortRadixIndexGetBufferSize(isize, $ippsz, n) + allocate(tmp(n)) + call ippsSortRadixAscend_$ityp_I(x, isize, tmp) + deallocate(tmp) + end + +SUBST [ X, type, ityp, n, ippsz ] + ; real ; 32f ; 4 ; 13 ;; + i ; integer ; 32s ; 4 ; 11 ;; + i2 ; integer*2 ; 16s ; 2 ; 7 ;; +END_TEMPLATE + +BEGIN_TEMPLATE + subroutine $Xsort(x,iorder,isize) implicit none BEGIN_DOC @@ -289,12 +329,12 @@ BEGIN_TEMPLATE endif end subroutine $Xsort -SUBST [ X, type, Y ] - ; real ; i ;; - d ; double precision ; i8 ;; +SUBST [ X, type ] + d ; double precision ;; END_TEMPLATE BEGIN_TEMPLATE + subroutine $Xsort(x,iorder,isize) implicit none BEGIN_DOC @@ -306,8 +346,112 @@ BEGIN_TEMPLATE $type,intent(inout) :: x(isize) integer,intent(inout) :: iorder(isize) integer :: n -! call $Xradix_sort(x,iorder,isize,-1) - call quick_$Xsort(x,iorder,isize) + if (isize < 2) then + return + endif + call sorted_$Xnumber(x,isize,n) + if (isize == n) then + return + endif + if ( isize < 32) then + call insertion_$Xsort(x,iorder,isize) + else + call $Xradix_sort(x,iorder,isize,-1) + endif + end subroutine $Xsort + +SUBST [ X, type ] + i8 ; integer*8 ;; +END_TEMPLATE + +!---------------------- END INTEL +IRP_ELSE +!---------------------- NON-INTEL +BEGIN_TEMPLATE + + subroutine $Xsort_noidx(x,isize) + implicit none + BEGIN_DOC + ! Sort array x(isize). + END_DOC + integer,intent(in) :: isize + $type,intent(inout) :: x(isize) + integer, allocatable :: iorder(:) + integer :: i + allocate(iorder(isize)) + do i=1,isize + iorder(i)=i + enddo + call $Xsort(x,iorder,isize) + deallocate(iorder) + end subroutine $Xsort_noidx + +SUBST [ X, type ] + ; real ;; + d ; double precision ;; + i ; integer ;; + i8 ; integer*8 ;; + i2 ; integer*2 ;; +END_TEMPLATE + +BEGIN_TEMPLATE + + subroutine $Xsort(x,iorder,isize) + implicit none + BEGIN_DOC + ! Sort array x(isize). + ! iorder in input should be (1,2,3,...,isize), and in output + ! contains the new order of the elements. + END_DOC + integer,intent(in) :: isize + $type,intent(inout) :: x(isize) + integer,intent(inout) :: iorder(isize) + integer :: n + if (isize < 2) then + return + endif +! call sorted_$Xnumber(x,isize,n) +! if (isize == n) then +! return +! endif + if ( isize < 32) then + call insertion_$Xsort(x,iorder,isize) + else +! call heap_$Xsort(x,iorder,isize) + call quick_$Xsort(x,iorder,isize) + endif + end subroutine $Xsort + +SUBST [ X, type ] + ; real ;; + d ; double precision ;; +END_TEMPLATE + +BEGIN_TEMPLATE + + subroutine $Xsort(x,iorder,isize) + implicit none + BEGIN_DOC + ! Sort array x(isize). + ! iorder in input should be (1,2,3,...,isize), and in output + ! contains the new order of the elements. + END_DOC + integer,intent(in) :: isize + $type,intent(inout) :: x(isize) + integer,intent(inout) :: iorder(isize) + integer :: n + if (isize < 2) then + return + endif + call sorted_$Xnumber(x,isize,n) + if (isize == n) then + return + endif + if ( isize < 32) then + call insertion_$Xsort(x,iorder,isize) + else + call $Xradix_sort(x,iorder,isize,-1) + endif end subroutine $Xsort SUBST [ X, type ] @@ -316,6 +460,11 @@ SUBST [ X, type ] i2 ; integer*2 ;; END_TEMPLATE +IRP_ENDIF +!---------------------- END NON-INTEL + + + BEGIN_TEMPLATE subroutine $Xset_order(x,iorder,isize) implicit none @@ -413,10 +562,12 @@ SUBST [ X, type ] i2; integer*2 ;; END_TEMPLATE + BEGIN_TEMPLATE - recursive subroutine $Xradix_sort$big(x,iorder,isize,iradix) +recursive subroutine $Xradix_sort$big(x,iorder,isize,iradix) implicit none + BEGIN_DOC ! Sort integer array x(isize) using the radix sort algorithm. ! iorder in input should be (1,2,3,...,isize), and in output @@ -552,24 +703,14 @@ BEGIN_TEMPLATE endif -! !$OMP PARALLEL DEFAULT(SHARED) if (isize > 1000000) -! !$OMP SINGLE if (i3>1_$int_type) then -! !$OMP TASK FIRSTPRIVATE(iradix_new,i3) SHARED(x,iorder) if(i3 > 1000000) call $Xradix_sort$big(x,iorder,i3,iradix_new-1) -! !$OMP END TASK endif if (isize-i3>1_$int_type) then -! !$OMP TASK FIRSTPRIVATE(iradix_new,i3) SHARED(x,iorder) if(isize-i3 > 1000000) call $Xradix_sort$big(x(i3+1_$int_type),iorder(i3+1_$int_type),isize-i3,iradix_new-1) -! !$OMP END TASK endif -! !$OMP TASKWAIT -! !$OMP END SINGLE -! !$OMP END PARALLEL - return endif @@ -624,16 +765,11 @@ BEGIN_TEMPLATE if (i1>1_$int_type) then - !$OMP TASK FIRSTPRIVATE(i0,iradix,i1) SHARED(x,iorder) if(i1 >1000000) call $Xradix_sort$big(x(i0+1_$int_type),iorder(i0+1_$int_type),i1,iradix-1) - !$OMP END TASK endif if (i0>1) then - !$OMP TASK FIRSTPRIVATE(i0,iradix) SHARED(x,iorder) if(i0 >1000000) call $Xradix_sort$big(x,iorder,i0,iradix-1) - !$OMP END TASK endif - !$OMP TASKWAIT end @@ -646,3 +782,4 @@ SUBST [ X, type, integer_size, is_big, big, int_type ] END_TEMPLATE + diff --git a/src/zmq/utils.irp.f b/src/zmq/utils.irp.f index 93dbd16a..7cb6c896 100644 --- a/src/zmq/utils.irp.f +++ b/src/zmq/utils.irp.f @@ -127,7 +127,9 @@ function zmq_port(ishift) END_DOC integer, intent(in) :: ishift character*(8) :: zmq_port + !$OMP CRITICAL(write) write(zmq_port,'(I8)') zmq_port_start+ishift + !$OMP END CRITICAL(write) zmq_port = adjustl(trim(zmq_port)) end @@ -518,7 +520,9 @@ subroutine new_parallel_job(zmq_to_qp_run_socket,zmq_socket_pull,name_in) zmq_to_qp_run_socket = new_zmq_to_qp_run_socket() zmq_socket_pull = new_zmq_pull_socket () + !$OMP CRITICAL(write) write(name,'(A,I8.8)') trim(name_in)//'.', icount + !$OMP END CRITICAL(write) sze = len(trim(name)) zmq_state = trim(name) call lowercase(name,sze) @@ -582,7 +586,9 @@ subroutine end_parallel_job(zmq_to_qp_run_socket,zmq_socket_pull,name_in) integer, save :: icount=0 icount = icount+1 + !$OMP CRITICAL(write) write(name,'(A,I8.8)') trim(name_in)//'.', icount + !$OMP END CRITICAL(write) sze = len(trim(name)) call lowercase(name,sze) if (name /= zmq_state) then @@ -704,7 +710,9 @@ integer function disconnect_from_taskserver_state(zmq_to_qp_run_socket, worker_i disconnect_from_taskserver_state = -1 + !$OMP CRITICAL(write) write(message,*) 'disconnect '//trim(state), worker_id + !$OMP END CRITICAL(write) sze = min(510,len(trim(message))) rc = f77_zmq_send(zmq_to_qp_run_socket, trim(message), sze, 0) @@ -781,7 +789,9 @@ integer function zmq_abort(zmq_to_qp_run_socket) character*(512) :: message zmq_abort = 0 + !$OMP CRITICAL(write) write(message,*) 'abort ' + !$OMP END CRITICAL(write) sze = len(trim(message)) @@ -823,7 +833,9 @@ integer function task_done_to_taskserver(zmq_to_qp_run_socket, worker_id, task_i task_done_to_taskserver = 0 + !$OMP CRITICAL(write) write(message,*) 'task_done '//trim(zmq_state), worker_id, task_id + !$OMP END CRITICAL(write) sze = len(trim(message)) rc = f77_zmq_send(zmq_to_qp_run_socket, trim(message), sze, 0) @@ -856,9 +868,11 @@ integer function tasks_done_to_taskserver(zmq_to_qp_run_socket, worker_id, task_ tasks_done_to_taskserver = 0 + !$OMP CRITICAL(write) allocate(character(LEN=64+n_tasks*12) :: message) write(fmt,*) '(A,X,A,I10,X,', n_tasks, '(I11,1X))' write(message,*) 'task_done '//trim(zmq_state), worker_id, (task_id(k), k=1,n_tasks) + !$OMP END CRITICAL(write) sze = len(trim(message)) rc = f77_zmq_send(zmq_to_qp_run_socket, trim(message), sze, 0) @@ -900,7 +914,9 @@ integer function get_task_from_taskserver(zmq_to_qp_run_socket,worker_id,task_id get_task_from_taskserver = 0 + !$OMP CRITICAL(write) write(message,*) 'get_task '//trim(zmq_state), worker_id + !$OMP END CRITICAL(write) sze = len(trim(message)) rc = f77_zmq_send(zmq_to_qp_run_socket, message, sze, 0) @@ -961,7 +977,9 @@ integer function get_tasks_from_taskserver(zmq_to_qp_run_socket,worker_id,task_i get_tasks_from_taskserver = 0 + !$OMP CRITICAL(write) write(message,'(A,A,X,I10,I10)') 'get_tasks ', trim(zmq_state), worker_id, n_tasks + !$OMP END CRITICAL(write) sze = len(trim(message)) rc = f77_zmq_send(zmq_to_qp_run_socket, message, sze, 0) @@ -1061,7 +1079,9 @@ integer function zmq_delete_task(zmq_to_qp_run_socket,zmq_socket_pull,task_id,mo zmq_delete_task = 0 + !$OMP CRITICAL(write) write(message,*) 'del_task ', zmq_state, task_id + !$OMP END CRITICAL(write) rc = f77_zmq_send(zmq_to_qp_run_socket,trim(message),len(trim(message)),0) if (rc /= len(trim(message))) then zmq_delete_task = -1 @@ -1101,7 +1121,9 @@ integer function zmq_delete_task_async_send(zmq_to_qp_run_socket,task_id,sending endif zmq_delete_task_async_send = 0 + !$OMP CRITICAL(write) write(message,*) 'del_task ', zmq_state, task_id + !$OMP END CRITICAL(write) rc = f77_zmq_send(zmq_to_qp_run_socket,trim(message),len(trim(message)),0) if (rc /= len(trim(message))) then zmq_delete_task_async_send = -1 @@ -1159,8 +1181,10 @@ integer function zmq_delete_tasks(zmq_to_qp_run_socket,zmq_socket_pull,task_id,n allocate(character(LEN=64+n_tasks*12) :: message) + !$OMP CRITICAL(write) write(fmt,*) '(A,1X,A,1X,', n_tasks, '(I11,1X))' write(message,*) 'del_task '//trim(zmq_state), (task_id(k), k=1,n_tasks) + !$OMP END CRITICAL(write) rc = f77_zmq_send(zmq_to_qp_run_socket,trim(message),len(trim(message)),0) @@ -1206,8 +1230,10 @@ integer function zmq_delete_tasks_async_send(zmq_to_qp_run_socket,task_id,n_task allocate(character(LEN=64+n_tasks*12) :: message) + !$OMP CRITICAL(write) write(fmt,*) '(A,1X,A,1X,', n_tasks, '(I11,1X))' write(message,*) 'del_task '//trim(zmq_state), (task_id(k), k=1,n_tasks) + !$OMP END CRITICAL(write) rc = f77_zmq_send(zmq_to_qp_run_socket,trim(message),len(trim(message)),0)