Wf analyzis

ocaml/Basis.ml

@@ -36,9 +36,11 @@ let read_element in_channel at_number element =
 
 
 
-let to_string_general ~fmt ~atom_sep b =
+let to_string_general ~fmt ~atom_sep ?ele_array b =
   let new_nucleus n = 
-    Printf.sprintf "Atom %d" n
+    match ele_array with
+    | None -> Printf.sprintf "Atom %d" n
+    | Some x -> Printf.sprintf "%s" (Element.to_string x.(n-1))
   in
   let rec do_work accu current_nucleus = function
   | [] -> List.rev accu
@@ -56,12 +58,12 @@ let to_string_general ~fmt ~atom_sep b =
   do_work [new_nucleus 1] 1 b
   |> String.concat ~sep:"\n"
 
-let to_string_gamess =
+let to_string_gamess ?ele_array =
-    to_string_general ~fmt:Gto.Gamess ~atom_sep:""
+    to_string_general ?ele_array ~fmt:Gto.Gamess ~atom_sep:""
 
-let to_string_gaussian b =
+let to_string_gaussian ?ele_array b =
   String.concat ~sep:"\n"
-  [ to_string_general ~fmt:Gto.Gaussian ~atom_sep:"****" b ; "****" ]
+  [ to_string_general ?ele_array ~fmt:Gto.Gaussian ~atom_sep:"****" b ; "****" ]
 
 let to_string ?(fmt=Gto.Gamess) =
   match fmt with

ocaml/Basis.mli

@@ -14,7 +14,7 @@ val read_element :
   in_channel -> Nucl_number.t -> Element.t -> (Gto.t * Nucl_number.t) list
 
 (** Convert the basis to a string *)
-val to_string : ?fmt:Gto.fmt -> (Gto.t * Nucl_number.t) list -> string
+val to_string : ?fmt:Gto.fmt -> ?ele_array:Element.t array -> (Gto.t * Nucl_number.t) list -> string
 
 (** Convert the basis to an MD5 hash *)
 val to_md5 :  (Gto.t * Nucl_number.t) list -> MD5.t

plugins/analyze_wf/NEEDED_CHILDREN_MODULES

@@ -0,0 +1 @@
+Determinants

plugins/analyze_wf/README.rst

@@ -0,0 +1,12 @@
+==========
+analyze_wf
+==========
+
+Needed Modules
+==============
+.. Do not edit this section It was auto-generated
+.. by the `update_README.py` script.
+Documentation
+=============
+.. Do not edit this section It was auto-generated
+.. by the `update_README.py` script.

plugins/analyze_wf/analyze_wf.irp.f

@@ -0,0 +1,70 @@
+program analyze_wf
+  implicit none
+  BEGIN_DOC
+! Wave function analyzis
+  END_DOC
+  read_wf = .True.
+  SOFT_TOUCH read_wf
+  call run()
+end
+
+subroutine run
+  implicit none
+  integer                        :: istate, i
+  integer                        :: class(0:mo_tot_num,5)
+  double precision               :: occupation(mo_tot_num)
+
+  write(*,'(A)')  'MO Occupation'
+  write(*,'(A)')  '============='
+  write(*,'(A)')  ''
+  do istate=1,N_states
+    call get_occupation_from_dets(occupation,1)
+    write(*,'(A)')  ''
+    write(*,'(A,I3)'),  'State ', istate
+    write(*,'(A)')  '---------------'
+    write(*,'(A)')  ''
+    write (*,'(A)')  '======== ================'
+    class = 0
+    do i=1,mo_tot_num
+      write (*,'(I8,X,F16.10)')  i, occupation(i)
+      if (occupation(i) > 1.999d0) then
+        class(0,1) += 1
+        class( class(0,1), 1) = i
+      else if (occupation(i) > 1.95d0) then
+        class(0,2) += 1
+        class( class(0,2), 2) = i
+      else if (occupation(i) < 0.001d0) then
+        class(0,5) += 1
+        class( class(0,5), 5) = i
+      else if (occupation(i) < 0.01d0) then
+        class(0,4) += 1
+        class( class(0,4), 4) = i
+      else 
+        class(0,3) += 1
+        class( class(0,3), 3) = i
+      endif
+    enddo
+    write (*,'(A)')  '======== ================'
+    write (*,'(A)')  ''
+
+    write (*,'(A)')  'Suggested classes'
+    write (*,'(A)')  '-----------------'
+    write (*,'(A)')  ''
+    write (*,'(A)')  'Core :'
+    write (*,*)      (class(i,1), ',', i=1,class(0,1))
+    write (*,*)      ''
+    write (*,'(A)')  'Inactive :'
+    write (*,*)      (class(i,2), ',', i=1,class(0,2))
+    write (*,'(A)')  ''
+    write (*,'(A)')  'Active :'
+    write (*,*)      (class(i,3), ',', i=1,class(0,3))
+    write (*,'(A)')  ''
+    write (*,'(A)')  'Virtual :'
+    write (*,*)      (class(i,4), ',', i=1,class(0,4))
+    write (*,'(A)')  ''
+    write (*,'(A)')  'Deleted :'
+    write (*,*)      (class(i,5), ',', i=1,class(0,5))
+    write (*,'(A)')  ''
+  enddo
+
+end

plugins/analyze_wf/occupation.irp.f

@@ -0,0 +1,23 @@
+subroutine get_occupation_from_dets(occupation, istate)
+  implicit none
+  double precision, intent(out)  :: occupation(mo_tot_num)
+  integer, intent(in)            :: istate
+  BEGIN_DOC
+  ! Returns the average occupation of the MOs
+  END_DOC
+  integer                        :: i,j, ispin
+  integer                        :: list(N_int*bit_kind_size,2)
+  integer                        :: n_elements(2)
+  double precision               :: c
+  
+  occupation = 0.d0
+  do i=1,N_det
+    c = psi_coef(i,istate)*psi_coef(i,istate)
+    call bitstring_to_list_ab(psi_det(1,1,i), list, n_elements, N_int)
+    do ispin=1,2
+      do j=1,n_elements(ispin)
+        occupation( list(j,ispin) ) += c
+      enddo
+    enddo
+  enddo
+end

src/Determinants/density_matrix.irp.f

@@ -194,7 +194,6 @@ subroutine set_natural_mos
  double precision, allocatable :: tmp(:,:)
 
  label = "Natural"
-! call mo_as_eigvectors_of_mo_matrix(one_body_dm_mo,size(one_body_dm_mo,1),mo_tot_num,label,-1)
  call mo_as_svd_vectors_of_mo_matrix(one_body_dm_mo,size(one_body_dm_mo,1),mo_tot_num,mo_tot_num,label)
 
 end