Merge pull request #16 from scemama/develop

Develop

ocaml/Message.ml

@@ -288,13 +288,14 @@ module Psi : sig
       n_det_selectors : Strictly_positive_int.t option;
       psi_det   :  string                    ;
       psi_coef  :  string                    ;
+      energy    :  string;
   }
   val create : n_state:Strictly_positive_int.t
      -> n_det:Strictly_positive_int.t 
      -> psi_det_size:Strictly_positive_int.t 
      -> n_det_generators:Strictly_positive_int.t option
      -> n_det_selectors:Strictly_positive_int.t option
-     -> psi_det:string -> psi_coef:string -> t
+     -> psi_det:string -> psi_coef:string -> energy:string -> t
 end = struct
   type t = 
   {
@@ -305,14 +306,16 @@ end = struct
       n_det_selectors : Strictly_positive_int.t option;
       psi_det   :  string                    ;
       psi_coef  :  string                    ;
+      energy    :  string                    ;
   }
   let create ~n_state ~n_det ~psi_det_size
-    ~n_det_generators ~n_det_selectors ~psi_det ~psi_coef =
+    ~n_det_generators ~n_det_selectors ~psi_det ~psi_coef 
+    ~energy =
     assert (Strictly_positive_int.to_int n_det <=
             Strictly_positive_int.to_int psi_det_size );
     {  n_state; n_det ; psi_det_size ;
        n_det_generators ; n_det_selectors ;
-       psi_det ; psi_coef }
+       psi_det ; psi_coef ; energy }
 end
 
 (** GetPsiReply_msg : Reply to the GetPsi message *)
@@ -329,19 +332,6 @@ end = struct
     psi       :  Psi.t }
   let create ~client_id ~psi =
     {  client_id ; psi }
-  let to_string_list x =
-    let g, s = 
-      match x.psi.Psi.n_det_generators, x.psi.Psi.n_det_selectors with
-      | Some g, Some s -> Strictly_positive_int.to_int g, Strictly_positive_int.to_int s
-      | _ -> -1, -1
-    in
-    [ Printf.sprintf "get_psi_reply %d %d %d %d %d %d"
-      (Id.Client.to_int x.client_id)
-      (Strictly_positive_int.to_int x.psi.Psi.n_state)
-      (Strictly_positive_int.to_int x.psi.Psi.n_det) 
-      (Strictly_positive_int.to_int x.psi.Psi.psi_det_size)
-      g s ;
-      x.psi.Psi.psi_det ; x.psi.Psi.psi_coef ]
   let to_string x =
     let g, s = 
       match x.psi.Psi.n_det_generators, x.psi.Psi.n_det_selectors with
@@ -354,6 +344,9 @@ end = struct
       (Strictly_positive_int.to_int x.psi.Psi.n_det) 
       (Strictly_positive_int.to_int x.psi.Psi.psi_det_size) 
       g s
+  let to_string_list x =
+    [ to_string x ;
+      x.psi.Psi.psi_det ; x.psi.Psi.psi_coef ; x.psi.Psi.energy ]
 end
 
 
@@ -375,7 +368,8 @@ module PutPsi_msg : sig
      psi_det:string option ->
      psi_coef:string option ->
      n_det_generators: string option -> 
-     n_det_selectors:string option ->  t
+     n_det_selectors:string option ->
+     energy:string option -> t
   val to_string_list : t -> string list
   val to_string : t -> string 
 end = struct
@@ -388,7 +382,7 @@ end = struct
     n_det_selectors  : Strictly_positive_int.t option;
     psi       :  Psi.t option }
   let create ~client_id ~n_state ~n_det ~psi_det_size ~psi_det ~psi_coef 
-    ~n_det_generators ~n_det_selectors  =
+    ~n_det_generators ~n_det_selectors ~energy =
     let n_state, n_det, psi_det_size = 
        Int.of_string n_state 
        |> Strictly_positive_int.of_int ,
@@ -407,45 +401,19 @@ end = struct
       | _ -> None, None
     in
     let psi =
-      match (psi_det, psi_coef) with
+      match (psi_det, psi_coef, energy) with
-      | (Some psi_det, Some psi_coef) ->
+      | (Some psi_det, Some psi_coef, Some energy) ->
         Some (Psi.create ~n_state ~n_det ~psi_det_size ~psi_det
-          ~psi_coef ~n_det_generators ~n_det_selectors)
+          ~psi_coef ~n_det_generators ~n_det_selectors ~energy)
       | _ -> None
     in
     { client_id = Id.Client.of_string client_id ;
       n_state ; n_det ; psi_det_size ; n_det_generators ;
       n_det_selectors ; psi }
-  let to_string_list x =
+
-    match x.n_det_generators, x.n_det_selectors, x.psi with
-    | Some g, Some s, Some psi ->
-      [ Printf.sprintf "put_psi %d %d %d %d %d %d"
-        (Id.Client.to_int x.client_id)
-        (Strictly_positive_int.to_int x.n_state)
-        (Strictly_positive_int.to_int x.n_det) 
-        (Strictly_positive_int.to_int x.psi_det_size)  
-        (Strictly_positive_int.to_int g)  
-        (Strictly_positive_int.to_int s) ; 
-          psi.Psi.psi_det ; psi.Psi.psi_coef ]
-    | Some g, Some s, None ->
-      [ Printf.sprintf "put_psi %d %d %d %d %d %d"
-        (Id.Client.to_int x.client_id)
-        (Strictly_positive_int.to_int x.n_state)
-        (Strictly_positive_int.to_int x.n_det) 
-        (Strictly_positive_int.to_int x.psi_det_size)  
-        (Strictly_positive_int.to_int g)  
-        (Strictly_positive_int.to_int s) ; 
-          "None" ; "None" ]
-    | _ ->
-      [ Printf.sprintf "put_psi %d %d %d %d -1 -1"
-        (Id.Client.to_int x.client_id)
-        (Strictly_positive_int.to_int x.n_state)
-        (Strictly_positive_int.to_int x.n_det) 
-        (Strictly_positive_int.to_int x.psi_det_size) ;
-          "None" ; "None" ]
   let to_string x =
-    match x.n_det_generators, x.n_det_selectors, x.psi with
+    match x.n_det_generators, x.n_det_selectors with
-    | Some g, Some s, Some psi ->
+    | Some g, Some s ->
       Printf.sprintf "put_psi %d %d %d %d %d %d"
         (Id.Client.to_int x.client_id)
         (Strictly_positive_int.to_int x.n_state)
@@ -453,21 +421,20 @@ end = struct
         (Strictly_positive_int.to_int x.psi_det_size)  
         (Strictly_positive_int.to_int g)  
         (Strictly_positive_int.to_int s) 
-    | Some g, Some s, None ->
+    | _, _ ->
-      Printf.sprintf "put_psi %d %d %d %d %d %d"
-        (Id.Client.to_int x.client_id)
-        (Strictly_positive_int.to_int x.n_state)
-        (Strictly_positive_int.to_int x.n_det) 
-        (Strictly_positive_int.to_int x.psi_det_size)  
-        (Strictly_positive_int.to_int g)  
-        (Strictly_positive_int.to_int s) 
-    | _, _, _ ->
       Printf.sprintf "put_psi %d %d %d %d %d %d"
         (Id.Client.to_int x.client_id)
         (Strictly_positive_int.to_int x.n_state)
         (Strictly_positive_int.to_int x.n_det) 
         (Strictly_positive_int.to_int x.psi_det_size)  
         (-1) (-1)
+
+  let to_string_list x =
+    match x.psi with
+    | Some psi ->
+      [ to_string x ; psi.Psi.psi_det ; psi.Psi.psi_coef ; psi.Psi.energy ]
+    | None ->
+      [ to_string x ; "None" ; "None" ; "None" ]
 end
 
 (** PutPsiReply_msg : Reply to the PutPsi message *)
@@ -606,10 +573,10 @@ let of_string s =
   | "put_psi"    :: client_id :: n_state :: n_det :: psi_det_size :: n_det_generators :: n_det_selectors :: [] ->
        PutPsi   (PutPsi_msg.create ~client_id ~n_state ~n_det ~psi_det_size
                  ~n_det_generators:(Some n_det_generators) ~n_det_selectors:(Some n_det_selectors)
-                 ~psi_det:None ~psi_coef:None )
+                 ~psi_det:None ~psi_coef:None ~energy:None )
   | "put_psi"    :: client_id :: n_state :: n_det :: psi_det_size :: [] ->
        PutPsi   (PutPsi_msg.create ~client_id ~n_state ~n_det ~psi_det_size ~n_det_generators:None
-                ~n_det_selectors:None ~psi_det:None ~psi_coef:None )
+                ~n_det_selectors:None ~psi_det:None ~psi_coef:None ~energy:None)
   | "ok"         :: [] ->
        Ok (Ok_msg.create ())
   | "error"      :: rest ->

ocaml/TaskServer.ml

@@ -483,9 +483,9 @@ let put_psi msg rest_of_msg program_state rep_socket =
         | Some x -> x
         | None ->
           begin
-            let psi_det, psi_coef =
+            let psi_det, psi_coef, energy =
               match rest_of_msg with
-              | [ x ; y ] -> x, y
+              | [ x ; y ; e ] -> x, y, e
               | _ -> failwith "Badly formed put_psi message"
             in
             Message.Psi.create
@@ -496,6 +496,7 @@ let put_psi msg rest_of_msg program_state rep_socket =
               ~n_det_selectors:msg.Message.PutPsi_msg.n_det_selectors
               ~psi_det
               ~psi_coef
+              ~energy
           end
     in
     let new_program_state =

plugins/Full_CI/micro_pt2.irp.f

@@ -22,6 +22,7 @@ subroutine run_wf
 
   integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
   integer(ZMQ_PTR) :: zmq_to_qp_run_socket
+  double precision :: energy(N_states_diag)
 
   print *,  'Getting wave function'
   zmq_context = f77_zmq_ctx_new ()
@@ -30,7 +31,21 @@ subroutine run_wf
 
   ! TODO : do loop here
   ! TODO : wait_state
-  call zmq_get_psi(zmq_to_qp_run_socket, 1)
+  call zmq_get_psi(zmq_to_qp_run_socket,1,energy,size(energy))
+  integer :: j,k
+  do j=1,N_states_diag
+    do k=1,N_det
+     CI_eigenvectors(k,j) = psi_coef(k,j)
+    enddo
+    call get_s2_u0(psi_det,CI_eigenvectors(1,j),N_det,size(CI_eigenvectors,1),CI_eigenvectors_s2(j))
+  enddo
+  if (.True.) then
+    do k=1,size(ci_electronic_energy)
+      ci_electronic_energy(k) = energy(k)
+    enddo
+    SOFT_TOUCH ci_electronic_energy CI_eigenvectors_s2 CI_eigenvectors
+    print *,  energy(:)
+  endif
   call write_double(6,ci_energy,'Energy')
   zmq_state = 'h_apply_fci_pt2'
 

plugins/Selectors_full/zmq.irp.f

@@ -1,4 +1,4 @@
-subroutine zmq_put_psi(zmq_to_qp_run_socket,worker_id)
+subroutine zmq_put_psi(zmq_to_qp_run_socket,worker_id, energy, size_energy)
   use f77_zmq
   implicit none
   BEGIN_DOC
@@ -6,6 +6,8 @@ subroutine zmq_put_psi(zmq_to_qp_run_socket,worker_id)
   END_DOC
   integer(ZMQ_PTR), intent(in)   :: zmq_to_qp_run_socket
   integer, intent(in)            :: worker_id
+  integer, intent(in)            :: size_energy
+  double precision, intent(out)  :: energy(size_energy)
   integer                        :: rc
   character*(256)                :: msg
 
@@ -23,9 +25,15 @@ subroutine zmq_put_psi(zmq_to_qp_run_socket,worker_id)
     stop 'error'
   endif
 
-  rc = f77_zmq_send(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,0)
+  rc = f77_zmq_send(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,ZMQ_SNDMORE)
   if (rc /= psi_det_size*N_states*8) then
-    print *, 'f77_zmq_send(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,0)'
+    print *, 'f77_zmq_send(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,ZMQ_SNDMORE)'
+    stop 'error'
+  endif
+
+  rc = f77_zmq_send(zmq_to_qp_run_socket,energy,size_energy*8,0)
+  if (rc /= size_energy*8) then
+    print *, 'f77_zmq_send(zmq_to_qp_run_socket,energy,size_energy*8,0)'
     stop 'error'
   endif
 
@@ -40,7 +48,7 @@ end
 
 
 
-subroutine zmq_get_psi(zmq_to_qp_run_socket, worker_id)
+subroutine zmq_get_psi(zmq_to_qp_run_socket, worker_id, energy, size_energy)
   use f77_zmq
   implicit none
   BEGIN_DOC
@@ -48,6 +56,8 @@ subroutine zmq_get_psi(zmq_to_qp_run_socket, worker_id)
   END_DOC
   integer(ZMQ_PTR), intent(in)   :: zmq_to_qp_run_socket
   integer, intent(in)            :: worker_id
+  integer, intent(in)            :: size_energy
+  double precision, intent(out)  :: energy(size_energy)
   integer                        :: rc
   character*(64)                 :: msg
 
@@ -86,9 +96,15 @@ subroutine zmq_get_psi(zmq_to_qp_run_socket, worker_id)
     stop 'error'
   endif
 
-  rc = f77_zmq_recv(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,0)
+  rc = f77_zmq_recv(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,ZMQ_SNDMORE)
   if (rc /= psi_det_size*N_states*8) then
-    print *, '77_zmq_recv(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,0)'
+    print *, '77_zmq_recv(zmq_to_qp_run_socket,psi_coef,psi_det_size*N_states*8,ZMQ_SNDMORE)'
+    stop 'error'
+  endif
+
+  rc = f77_zmq_recv(zmq_to_qp_run_socket,energy,size_energy*8,0)
+  if (rc /= size_energy*8) then
+    print *, '77_zmq_recv(zmq_to_qp_run_socket,energy,size_energy*8,0)'
     stop 'error'
   endif
 

scripts/generate_h_apply.py

@@ -435,8 +435,9 @@ class H_apply_zmq(H_apply):
     norm_pert(k) = 0.d0
     H_pert_diag(k) = 0.d0
     norm_psi(k) = 0.d0
+    energy(k) = %s(k)
   enddo
-     """ 
+     """ % (self.energy)
      self.data["copy_buffer"] = """
   do i=1,N_det_generators
     do k=1,N_st

src/Determinants/H_apply_zmq.template.f

@@ -28,11 +28,12 @@ subroutine $subroutine($params_main)
   integer(ZMQ_PTR) :: zmq_to_qp_run_socket
   double precision, allocatable :: pt2_generators(:,:), norm_pert_generators(:,:)
   double precision, allocatable :: H_pert_diag_generators(:,:)
+  double precision              :: energy(N_st)
 
   call new_parallel_job(zmq_to_qp_run_socket,'$subroutine')
   zmq_socket_pair = new_zmq_pair_socket(.True.)
 
-  call zmq_put_psi(zmq_to_qp_run_socket,1)
+  call zmq_put_psi(zmq_to_qp_run_socket,1,energy,size(energy))
 
   do i_generator=N_det_generators,1,-1
     $skip

src/Integrals_Monoelec/pseudopot.f90

@@ -1068,7 +1068,7 @@ integer n
 double precision g,dble_fact,expo
 double precision, parameter :: sq_pi_ov_2=dsqrt(dacos(-1.d0)*0.5d0)
 expo=0.5d0*dfloat(n+1)
-crochet=dble_fact(n-1)/(2.d0*g)**expo
+crochet=dble_fact(n-1)/(g+g)**expo
 if(mod(n,2).eq.0)crochet=crochet*sq_pi_ov_2
 end
 
@@ -1840,8 +1840,8 @@ double precision function int_prod_bessel(l,gam,n,m,a,b,arg)
     int=0.d0
     done=.false.
 
-    n_1 = 2*(n)+1
+    n_1 = n+n+1
-    m_1 = 2*m+1
+    m_1 = m+m+1
     nlm = n+m+l
     pi=dacos(-1.d0)
     a_over_b_square = (a/b)**2
@@ -1892,7 +1892,7 @@ double precision function int_prod_bessel(l,gam,n,m,a,b,arg)
 
         !Compute the s_q+1_0
 !        s_q_0=s_q_0*(2.d0*q+nlm+1)*b**2/((2.d0*(m+q)+3)*4.d0*(q+1)*gam)
-        s_q_0=s_q_0*(2.d0*q+nlm+1)*b*b/((8.d0*(m+q)+12.d0)*(q+1)*gam)
+        s_q_0=s_q_0*(q+q+nlm+1)*b*b/(dble(8*(m+q)+12)*(q+1)*gam)
         
         if(mod(n+m+l,2).eq.1)s_q_0=s_q_0*dsqrt(pi*.5d0)
         ! Increment q
@@ -1933,7 +1933,7 @@ double precision function int_prod_bessel_large(l,gam,n,m,a,b,arg)
       double precision xq(100),wq(100)
 
       u=(a+b)/(2.d0*dsqrt(gam))
-      factor=dexp(u**2-arg)/dsqrt(gam)
+      factor=dexp(u*u-arg)/dsqrt(gam)
 
 xq(1)= 5.38748089001123    
 xq(2)= 4.60368244955074