Skip to content

Commit

Permalink
redundant intercluster ras1 and ras3 pairs added in projection vector
Browse files Browse the repository at this point in the history
  • Loading branch information
@arnab82
arnab82 committed Oct 27, 2023
1 parent d272284 commit c8b3b2e
Show file tree
Hide file tree
Showing 2 changed files with 60 additions and 22 deletions.
12 changes: 6 additions & 6 deletions src/finite_difference.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -754,12 +754,12 @@ function orbital_hessian_numerical(ints, clusters, kappa, fspace, ansatze,d::RDM
function step_numerical!(ints, d1, k)
norb = n_orb(ints)
K = unpack_gradient(k, norb)
if zero_intra_rots
# Remove intracluster rotations
for ci in clusters
K[ci.orb_list, ci.orb_list] .= 0
end
end
# if zero_intra_rots
# # Remove intracluster rotations
# for ci in clusters
# K[ci.orb_list, ci.orb_list] .= 0
# end
# end

Ui = exp(K)

Expand Down
70 changes: 54 additions & 16 deletions src/incore_cmf.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -440,6 +440,7 @@ function cmf_ci(ints, clusters, fspace, ansatze::Vector{<:Ansatz}, in_rdm1::RDM1
use_pyscf = true,
sequential = false)
rdm1 = deepcopy(in_rdm1)
# println(size(rdm1))
energies = []
e_prev = 0

Expand All @@ -463,7 +464,7 @@ function cmf_ci(ints, clusters, fspace, ansatze::Vector{<:Ansatz}, in_rdm1::RDM1
sequential = sequential
)
rdm1_curr = assemble_full_rdm(clusters, rdm1_dict)

# println(size(rdm1_curr))
append!(energies,e_curr)
error = (rdm1_curr.a+rdm1_curr.b) - (rdm1.a+rdm1.b)
d_err = norm(error)
Expand Down Expand Up @@ -1375,10 +1376,12 @@ function cmf_oo_diis(ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fspace
max_ss_size = 8,
diis_start = 1,
alpha = .1,
use_pyscf=true,
step_trust_region=0.95,
use_pyscf = true,
zero_intra_rots = true,
orb_hessian = true,
sequential = false
sequential = false,
trust_region = false
) where T
#={{{=#
println(" Solve OO-CMF with DIIS")
Expand Down Expand Up @@ -1410,8 +1413,9 @@ function cmf_oo_diis(ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fspace
tol_d1 = tol_d1,
tol_ci = tol_ci,
use_pyscf=use_pyscf,
verbose = 0,
sequential = sequential)
verbose = verbose,
sequential = sequential
)
d1_i, d2_i = assemble_full_rdm(clusters, rdm1_dict, rdm2_dict)

d1_i = orbital_rotation(d1_i, Ui')
Expand Down Expand Up @@ -1474,8 +1478,13 @@ function cmf_oo_diis(ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fspace
else
step_i = alpha*g_i
end

if trust_region==true
if norm(step_i)> step_trust_region
step_i=step_i*step_trust_region/norm(step_i)
end
end
k_i = k_i - step_i


if nss < max_ss_size
nss += 1
Expand Down Expand Up @@ -1567,7 +1576,7 @@ function cmf_oo_diis(ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fspace


if norm(g_i) < tol_oo
@printf("*ooCMF Iter: %4i Total= %16.15f G= %12.3e #SS: %4s\n", i, e_i, norm(g_i), nss)
@printf("*ooCMF Iter: %4i Total= %16.15f G= %12.3e step_size= %12.3e #SS: %4s\n", i, e_i, norm(g_i),norm(k_i), nss)
break
end

Expand All @@ -1577,7 +1586,7 @@ function cmf_oo_diis(ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fspace



@printf(" ooCMF Iter: %4i Total= %16.15f G= %12.3e #SS: %4s\n", i, e_i, norm(g_i), nss)
@printf("*ooCMF Iter: %4i Total= %16.15f G= %12.3e step_size= %12.3e #SS: %4s\n", i, e_i, norm(g_i),norm(k_i), nss)
end

U = exp(unpack_gradient(k_i,norb))
Expand Down Expand Up @@ -1633,10 +1642,12 @@ function cmf_oo_newton( ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fsp
tol_d1 = 1e-7,
tol_ci = 1e-8,
verbose = 0,
step_trust_region=0.95,
use_pyscf=true,
zero_intra_rots =true,
sequential = false,
step_trust_region=0.95
trust_region=false,

) where T
#={{{=#
println(" Solve OO-CMF with newton")
Expand Down Expand Up @@ -1668,15 +1679,15 @@ function cmf_oo_newton( ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fsp
maxiter_ci = maxiter_ci,
tol_d1 = tol_d1,
tol_ci = tol_ci,
verbose = 0,
verbose = verbose,
use_pyscf=use_pyscf,
sequential = sequential)

gd1, gd2 = assemble_full_rdm(clusters, rdm1_dict, rdm2_dict)
g_i = build_orbital_gradient(ints, gd1, gd2)
# energy_computed=compute_energy(ints, gd1,gd2)
packed_hessian=RDM.build_orbital_hessian(ints,gd1,gd2)

# display(packed_hessian)
# if zero_intra_rots
# g_i = unpack_gradient(g_i, norb)
# for ci in clusters
Expand Down Expand Up @@ -1704,9 +1715,10 @@ function cmf_oo_newton( ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fsp
else
step_i=-(pinv(h_i)*(g_i))#;atol=1e-8
end

if norm(step_i)> step_trust_region
step_i=step_i*0.95/norm(step_i)
if trust_region==true
if norm(step_i)> step_trust_region
step_i=step_i*step_trust_region/norm(step_i)
end
end
e = ei
U = U*Ui
Expand All @@ -1724,7 +1736,7 @@ function cmf_oo_newton( ints_in::InCoreInts{T}, clusters::Vector{MOCluster}, fsp
end

function convert_pairs(original_list, rearranged_pairs)
# Create a mapping from the rearranged index to the original index
# Create a mapping from the rearranged i6dex to the original index
index_mapping = Dict{Int,Int}()
for (original_idx, orbital) in enumerate(original_list)
index_mapping[original_idx] = orbital
Expand All @@ -1746,11 +1758,37 @@ function projection_vector(ansatze::Vector{<:Ansatz}, clusters, norb)
# println(cluster)
count += 1
tmp = ActiveSpaceSolvers.invariant_orbital_rotations(cluster)
# println(tmp)
tmp_global = convert_pairs(clusters_new[count], tmp)
# println(tmp_global)
append!(invar, tmp_global)
# println(invar)
end

for i in ansatze
if typeof(i) == RASCIAnsatz
ras1i, ras2i, ras3i = ActiveSpaceSolvers.RASCI.make_rasorbs(i.ras_spaces[1], i.ras_spaces[2], i.ras_spaces[3], i.no)
for j in ansatze
ras1j, ras2j, ras3j = ActiveSpaceSolvers.RASCI.make_rasorbs(j.ras_spaces[1], j.ras_spaces[2], j.ras_spaces[3], j.no)
pairs = []
for a in 1:length(ras1i)
for b in a:length(ras1j)
# println(ras1i[a])
push!(pairs, (ras1i[a],i.no+ras1j[b]))
end
end

for e in 1:length(ras3i)
for f in e:length(ras3j)
push!(pairs, (ras3i[e],i.no+ras3j[f]))
end
end
# println(pairs)
append!(invar,pairs)
end
end
end
# append!(invar,ras13_ras13_pairs)
# println(invar)
fci = ActiveSpaceSolvers.FCIAnsatz(norb, 0, 0) #dummie FCI anstaz to generate all pairs
full_list = ActiveSpaceSolvers.invariant_orbital_rotations(fci)

Expand Down

0 comments on commit c8b3b2e

Please sign in to comment.