Zero variance almost OK

ocaml/Qmcchem_dataserver.ml

@@ -379,6 +379,7 @@ let run ?(daemon=true) ezfio_filename =
       same host name, but different PIDs. The result is merging all the CPU cores of
       the compute nodes. Happens when [max_file_size] is reached.
   *)
+  (*
   let compress_block_file filename =
     let t0 = 
       Time.now ()
@@ -389,6 +390,7 @@ let run ?(daemon=true) ezfio_filename =
     send_log "status" 0 t0 "Compressed block file";
     block_channel := Out_channel.create ~append:true block_channel_filename_locked
   in
+  *)
 
 
 
@@ -817,7 +819,9 @@ let run ?(daemon=true) ezfio_filename =
   let finalize ~t0 =
     print_string "Finalizing...";
     change_status Status.Stopped;
+    (*
     compress_block_file ();
+    *)
     send_log "status" 0 t0 "Done";
     close_debug_socket ();
     ZMQ.Context.terminate zmq_context;

ocaml/Qmcchem_result.ml

@@ -3,8 +3,20 @@ open Qptypes
 
 (** Display a table that can be plotted by gnuplot *)
 let display_table ~range property =
-  let p = Property.of_string property
-  |> Random_variable.of_raw_data ~range
+  let p =
+  (*---*)
+    if (property = "xxx") then
+    let open Random_variable in
+       (of_raw_data ~range Property.E_ref_zv) 
+                /! (of_raw_data ~range Property.Psi_no_cusp_weight) 
+              -! ( (of_raw_data ~range Property.E_loc_zv) 
+                /! (of_raw_data ~range Property.Weight_zv) )   
+                /! (of_raw_data ~range Property.Psi_no_cusp_weight2) 
+                +@ Input.Trial_wf_energy.to_float (Input.Trial_wf_energy. read () )
+    else
+  (*---*)
+    Property.of_string property
+    |> Random_variable.of_raw_data ~range
   in
   let  conv = Random_variable.convergence p
   and rconv = Random_variable.rev_convergence p
@@ -145,14 +157,19 @@ let display_summary ~range =
   in
   List.iter properties ~f:print_property ;
 
+  (*---*)
   let open Random_variable in
-    let p =   (of_raw_data ~range Property.E_ref_zv)
+    let p =  (of_raw_data ~range Property.E_ref_zv) 
+              /! (of_raw_data ~range Property.Psi_no_cusp_weight) 
             -! ( (of_raw_data ~range Property.E_loc_zv) 
-              /! (of_raw_data ~range Property.Weight_zv) ) 
+              /! (of_raw_data ~range Property.Weight_zv) )   
+              /! (of_raw_data ~range Property.Psi_no_cusp_weight2) 
+              +@ Input.Trial_wf_energy.to_float (Input.Trial_wf_energy. read () )
     in
     Printf.printf "%20s : %s\n"
      ("E_loc_zv(+)")
      (Random_variable.to_string p);
+  (*---*)
 
 (*---
   let open Random_variable in

ocaml/Random_variable.ml

@@ -452,6 +452,40 @@ let merge_per_compute_node_and_block_id =
        (Compute_node.to_string block.Block.compute_node)
        (Block_id.to_int block.Block.block_id) )
 
+let error_x_over_y = function
+| []    -> (Average.of_float 0., None)
+| (x,_)::[] -> (Average.of_float x, None)
+| (x,w)::tail ->
+    begin
+      let avcu0 = ref 0.
+      and ansum = ref w
+      and avsum = ref x
+      and avbl  = ref (x /. w)
+      and avsq  = ref 0.
+      and n     = ref 1.
+      in
+      let avcu  = ref !avbl
+      in
+      List.iter tail ~f:(fun (x,w) ->
+        avcu0 := !avsum /. !ansum;
+        ansum := !ansum +. w;
+        avsum := !avsum +. x;
+        avbl  := x /. w ;
+        if (!ansum <> 0.) then
+          avcu := !avsum /. !ansum
+        else ();
+        avsq := !avsq +. (1. -. w /. !ansum) *. (!avbl -. !avcu0) *.  (!avbl -. !avcu0) *. w;
+        n := !n +. 1.
+      );
+      let arg =
+        Float.abs (!avsq /.(!ansum *. (!n -. 1.)))
+      in
+      let error =
+        sqrt arg
+      in
+      (Average.of_float !avcu, Some (Error.of_float error) )
+    end
+
 
 (** Create float, variable operators *)
 let one_variable_operator ~update_value p f = 

src/PROPERTIES/properties_energy.irp.f

@@ -272,13 +272,42 @@ BEGIN_PROVIDER  [ double precision, E_loc ]
 END_PROVIDER
 
 
+BEGIN_PROVIDER [ double precision, psi_no_cusp_weight ]
+ psi_no_cusp_weight = 1.d0
+END_PROVIDER
+
+BEGIN_PROVIDER [ double precision, psi_no_cusp_weight2 ]
+ psi_no_cusp_weight2 = 1.d0
+END_PROVIDER
+
+BEGIN_PROVIDER [ double precision, E_loc_nocusp ]
+ implicit none
+ BEGIN_DOC
+ ! Local energy without fitcusp
+ END_DOC
+ double precision :: psi_cusp_inv
+ if (do_nucl_fitcusp) then
+   psi_cusp_inv = psi_value_inv
+   do_nucl_fitcusp = .False.
+   TOUCH do_nucl_fitcusp
+   psi_no_cusp_weight = psi_value * psi_cusp_inv
+   psi_no_cusp_weight2 = psi_no_cusp_weight * psi_no_cusp_weight
+   E_loc_nocusp = psi_no_cusp_weight * E_loc
+   do_nucl_fitcusp = .True.
+   TOUCH do_nucl_fitcusp
+ else
+   E_loc_nocusp = E_loc
+   psi_no_cusp_weight = 1.d0
+ endif
+
+END_PROVIDER
 
 BEGIN_PROVIDER [ double precision, E_ref_zv ]
  implicit none
  BEGIN_DOC
  ! E_loc - E_trial
  END_DOC
-  E_ref_zv = E_loc + E_trial
+  E_ref_zv = E_loc_nocusp 
 
   E_ref_zv_min = min(E_ref_zv,E_ref_zv_min)
   E_ref_zv_max = max(E_ref_zv,E_ref_zv_max)
@@ -291,7 +320,7 @@ BEGIN_PROVIDER [ double precision, E_loc_zv ]
  BEGIN_DOC
  ! Estimate of the VMC energy in PDMC 
  END_DOC
-  E_loc_zv = E_loc * dmc_zv_weight 
+  E_loc_zv = E_loc_nocusp * psi_no_cusp_weight * dmc_zv_weight 
 
   E_loc_zv_min = min(E_loc_zv,E_loc_zv_min)
   E_loc_zv_max = max(E_loc_zv,E_loc_zv_max)