LCPQ/QCaml
10
1
Fork 0
mirror of https://gitlab.com/scemama/QCaml.git synced 2024-11-07 22:53:41 +01:00
Code Releases Activity

Working n F12

Browse Source
This commit is contained in:
Anthony Scemama 2019-03-22 22:42:20 +01:00
parent 29ddf41c0b
commit b40a9b8db8
1 changed files with 23 additions and 19 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

42
CI/F12CI.ml
View File

@ -241,23 +241,29 @@ Util.debug_matrix "psi" psi;
Format.printf "Amplitude (1,1) : %f@." (f12_amplitudes psi).{1,1}; Format.printf "Amplitude (1,1) : %f@." (f12_amplitudes psi).{1,1};
Format.printf "Dressing vector(1,1) : %f@." (Matrix.get delta 1 1); Format.printf "Dressing vector(1,1) : %f@." (Matrix.get delta 1 1);
let f = 1.0 /. psi.{1,1} in
let delta_00 =
Util.list_range 2 (Mat.dim1 psi)
|> List.fold_left (fun accu i -> accu +.
(Matrix.get delta i 1) *. psi.{i,1} *. f) 0.
in
let delta = Matrix.to_mat delta in
delta.{1,1} <- delta.{1,1} -. delta_00;
(*------ (*------
TODO SINGLE STATE HERE TODO SINGLE STATE HERE
*) *)
(*
let m_H_dressed = Matrix.to_mat m_H in let m_H_dressed = Matrix.to_mat m_H in
let f = 1.0 /. psi.{1,1} in let f = 1.0 /. psi.{1,1} in
Util.list_range 1 (Mat.dim1 psi) Util.list_range 1 (Mat.dim1 psi)
|> List.iter (fun i -> |> List.iter (fun i ->
let delta = (Matrix.get delta i 1) *. f in let delta = delta.{i,1} *. f in
m_H_dressed.{i,1} <- m_H_dressed.{i,1} +. delta; m_H_dressed.{i,1} <- m_H_dressed.{i,1} +. delta;
if i <> 1 then if i <> 1 then
begin begin
m_H_dressed.{1,i} <- m_H_dressed.{1,i} +. delta; m_H_dressed.{1,i} <- m_H_dressed.{1,i} +. delta;
m_H_dressed.{1,1} <- m_H_dressed.{1,1} -. delta *. psi.{i,1} *. f;
end end
(* DIAGONAL DRESSING
m_H_dressed.{i,i} <- m_H_dressed.{i,i} +. (Matrix.get delta i 1) /. (psi.{i,1} +. 1.e-5);
*)
); );
let eigenvectors, eigenvalues = let eigenvectors, eigenvalues =
Util.diagonalize_symm m_H_dressed Util.diagonalize_symm m_H_dressed
@ -277,26 +283,25 @@ TODO SINGLE STATE HERE
eigenvectors, eigenvalues eigenvectors, eigenvalues
in in
iteration ci_coef iteration ci_coef
*)
(* (*
------- *) ------- *)
(*
let n_states = ci.CI.n_states in let n_states = ci.CI.n_states in
let f = 1.0 /. (psi.{1,1} ) in
let delta_00 =
Util.list_range 2 (Mat.dim1 psi)
|> List.fold_left (fun accu i -> accu +.
(Matrix.get delta i 1) *. psi.{i,1} *. f) 0.
in
let delta = Matrix.to_mat delta in
delta.{1,1} <- delta.{1,1} -. delta_00;
let diagonal = let diagonal =
Vec.init (Matrix.dim1 m_H) (fun i -> Matrix.get m_H i i +. (if i=1 then delta.{1,1} *. psi.{1,1} else 0.) ) Vec.init (Matrix.dim1 m_H) (fun i -> Matrix.get m_H i i +. (if i=1 then delta.{1,1} *. psi.{1,1} else 0.) )
in in
let delta = Matrix.dense_of_mat delta in
let matrix_prod c = let matrix_prod c =
Matrix.add let w =
(Matrix.mm ~transa:`T m_H c) Matrix.mm ~transa:`T m_H c
delta |> Matrix.to_mat
in
Util.list_range 1 (Mat.dim1 w)
|> List.iter (fun i ->
w.{i,1} <- w.{i,1} +. delta.{i,1} ;
if (i <> 1) then
w.{1,1} <- w.{1,1} +. delta.{i,1} *. f *. (Matrix.get c i 1);
);
Matrix.dense_of_mat w
in in
let eigenvectors, eigenvalues = let eigenvectors, eigenvalues =
Parallel.broadcast (lazy ( Parallel.broadcast (lazy (
@ -318,7 +323,6 @@ TODO SINGLE STATE HERE
eigenvectors, eigenvalues eigenvectors, eigenvalues
in in
iteration ci_coef iteration ci_coef
*)
) )
in in