diff --git a/ocaml/Input.ml b/ocaml/Input.ml index 5e012996..b1bb1060 100644 --- a/ocaml/Input.ml +++ b/ocaml/Input.ml @@ -7,5 +7,5 @@ include Input_bitmasks;; include Input_determinants_by_hand;; include Input_electrons;; include Input_mo_basis;; -include Input_nuclei;; +include Input_nuclei_by_hand;; include Input_auto_generated;; diff --git a/ocaml/Input_nuclei.ml b/ocaml/Input_nuclei_by_hand.ml similarity index 99% rename from ocaml/Input_nuclei.ml rename to ocaml/Input_nuclei_by_hand.ml index ca81629e..f36b6b82 100644 --- a/ocaml/Input_nuclei.ml +++ b/ocaml/Input_nuclei_by_hand.ml @@ -2,7 +2,7 @@ open Qptypes;; open Qputils;; open Core.Std;; -module Nuclei : sig +module Nuclei_by_hand : sig type t = { nucl_num : Nucl_number.t ; nucl_label : Element.t array; diff --git a/scripts/ezfio_interface/qp_edit_template b/scripts/ezfio_interface/qp_edit_template index af9b295c..9d327124 100644 --- a/scripts/ezfio_interface/qp_edit_template +++ b/scripts/ezfio_interface/qp_edit_template @@ -20,7 +20,7 @@ type keyword = | Determinants_by_hand | Electrons | Mo_basis -| Nuclei +| Nuclei_by_hand {keywords} @@ -30,7 +30,7 @@ let keyword_to_string = function | Determinants_by_hand -> "Determinants_by_hand" | Electrons -> "Electrons" | Mo_basis -> "MO basis" -| Nuclei -> "Molecule" +| Nuclei_by_hand -> "Molecule" {keywords_to_string} @@ -74,8 +74,8 @@ let get s = f Mo_basis.(read, to_rst) | Electrons -> f Electrons.(read, to_rst) - | Nuclei -> - f Nuclei.(read, to_rst) + | Nuclei_by_hand -> + f Nuclei_by_hand.(read, to_rst) | Ao_basis -> f Ao_basis.(read, to_rst) | Determinants_by_hand -> @@ -121,7 +121,7 @@ let set str s = {write} | Electrons -> write Electrons.(of_rst, write) s | Determinants_by_hand -> write Determinants_by_hand.(of_rst, write) s - | Nuclei -> write Nuclei.(of_rst, write) s + | Nuclei_by_hand -> write Nuclei_by_hand.(of_rst, write) s | Ao_basis -> () (* TODO *) | Mo_basis -> () (* TODO *) end @@ -184,7 +184,7 @@ let run check_only ?ndet ?state ezfio_filename = *) let tasks = [ - Nuclei ; + Nuclei_by_hand ; Ao_basis; Electrons ; {tasks} diff --git a/src/Nuclei/EZFIO.cfg b/src/Nuclei/EZFIO.cfg index 81c168ff..dab6b211 100644 --- a/src/Nuclei/EZFIO.cfg +++ b/src/Nuclei/EZFIO.cfg @@ -19,4 +19,15 @@ interface: ezfio, provider doc: Nuclear coordinates in the format (:, {x,y,z}) type: double precision size: (nuclei.nucl_num,3) -interface: ezfio \ No newline at end of file +interface: ezfio + +[disk_access_nuclear_repulsion] +doc: Read/Write Nuclear Repulsion from/to disk [ Write | Read | None ] +type: Disk_access +interface: ezfio,provider,ocaml +default: None + +[nuclear_repulsion] +doc: Nuclear repulsion (Computed automaticaly or Read in the EZFIO) +type:double precision +interface: ezfio diff --git a/src/Nuclei/nuclei.irp.f b/src/Nuclei/nuclei.irp.f index 34fae989..577b8b92 100644 --- a/src/Nuclei/nuclei.irp.f +++ b/src/Nuclei/nuclei.irp.f @@ -141,34 +141,53 @@ BEGIN_PROVIDER [ double precision, positive_charge_barycentre,(3)] enddo END_PROVIDER - BEGIN_PROVIDER [ double precision, nuclear_repulsion ] +BEGIN_PROVIDER [ double precision, nuclear_repulsion ] implicit none BEGIN_DOC ! Nuclear repulsion energy END_DOC - integer :: k,l - double precision :: Z12, r2, x(3) - nuclear_repulsion = 0.d0 - do l = 1, nucl_num - do k = 1, nucl_num - if(k == l) then - cycle - endif - Z12 = nucl_charge(k)*nucl_charge(l) - x(1) = nucl_coord(k,1) - nucl_coord(l,1) - x(2) = nucl_coord(k,2) - nucl_coord(l,2) - x(3) = nucl_coord(k,3) - nucl_coord(l,3) - r2 = x(1)*x(1) + x(2)*x(2) + x(3)*x(3) - nuclear_repulsion += Z12/dsqrt(r2) - enddo - enddo - nuclear_repulsion *= 0.5d0 - + + IF (disk_access_nuclear_repulsion.EQ.'Read') THEN + print*, 'nuclear_repulsion read from disk' + LOGICAL :: has + call ezfio_has_nuclei_nuclear_repulsion(has) + if (has) then + call ezfio_get_nuclei_nuclear_repulsion(nuclear_repulsion) + else + print *, 'nuclei/nuclear_repulsion not found in EZFIO file' + stop 1 + endif + + ELSE + + integer :: k,l + double precision :: Z12, r2, x(3) + nuclear_repulsion = 0.d0 + do l = 1, nucl_num + do k = 1, nucl_num + if(k == l) then + cycle + endif + Z12 = nucl_charge(k)*nucl_charge(l) + x(1) = nucl_coord(k,1) - nucl_coord(l,1) + x(2) = nucl_coord(k,2) - nucl_coord(l,2) + x(3) = nucl_coord(k,3) - nucl_coord(l,3) + r2 = x(1)*x(1) + x(2)*x(2) + x(3)*x(3) + nuclear_repulsion += Z12/dsqrt(r2) + enddo + enddo + nuclear_repulsion *= 0.5d0 + END IF + call write_time(output_Nuclei) call write_double(output_Nuclei,nuclear_repulsion, & 'Nuclear repulsion energy') + + IF (disk_access_nuclear_repulsion.EQ.'Write') THEN + call ezfio_set_nuclei_nuclear_repulsion(nuclear_repulsion) + END IF END_PROVIDER - + BEGIN_PROVIDER [ character*(128), element_name, (78)] BEGIN_DOC ! Array of the name of element, sorted by nuclear charge (integer)