Skip to content

Commit

Permalink
quick update lappy
Browse files Browse the repository at this point in the history
  • Loading branch information
@kfgarrity
kfgarrity committed Mar 19, 2024
1 parent 87c1f6d commit dd29947
Show file tree
Hide file tree
Showing 2 changed files with 44 additions and 35 deletions.
35 changes: 20 additions & 15 deletions src/CalcTB_laguerre.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -528,7 +528,7 @@ function make_coefs(at_list, dim; datH=missing, datS=missing, cutoff=18.01, min_

for a in at_arr
#println("a $a")
println(typeof(a))
# println(typeof(a))
if !haskey(lim, a)
lim[a] = 0.02
end
Expand Down Expand Up @@ -568,9 +568,9 @@ function make_coefs(at_list, dim; datH=missing, datS=missing, cutoff=18.01, min_
end
end

for x in [dim, datH, datS, totH, totS, data_info, inds_int, at_list, orbs, cutoff, min_dist, dist_frontier2, version, lim, repval]
println(typeof(x))
end
# for x in [dim, datH, datS, totH, totS, data_info, inds_int, at_list, orbs, cutoff, min_dist, dist_frontier2, version, lim, repval]
# println(typeof(x))
# end

return coefs(dim, datH, datS, totH, totS, data_info, inds_int, at_list, orbs, cutoff, min_dist, dist_frontier2, version, lim, repval)

Expand Down Expand Up @@ -3938,12 +3938,14 @@ function calc_tb_prepare_fast(reference_tbc::tb_crys; use_threebody=false, use_t
cut = cutoff_fn(dist, cutoff_onX - cutoff_length, cutoff_onX)
end

ad = 2.0*dist
expa=exp.(-0.5*ad)
if dist > 1e-5
ad = 2.0*dist
expa=exp.(-0.5*ad)

rho[a1, 1] += (1.0 * expa) * cut
rho[a1, 2] += (1.0 .- ad) * expa * cut

rho[a1, 1] += (1.0 * expa) * cut
rho[a1, 2] += (1.0 .- ad) * expa * cut
# println("ADD RHO FIT $a1 dist $dist 1 $((1.0 * expa) * cut) 2 $( (1.0 .- ad) * expa * cut) cut $cut")
end

for o1 = orb2ind[a1]
a1a,t1a,s1 = ind2orb[o1]
Expand Down Expand Up @@ -4006,12 +4008,12 @@ function calc_tb_prepare_fast(reference_tbc::tb_crys; use_threebody=false, use_t
for o = orb2ind[a]
aa,t,s = ind2orb[o]
ind = ind_conversion[(o,o,c_zero)]
println("t $t ", typeof(t))
# println("t $t ", typeof(t))
at_set = Set((t,))
println(at_set)
println("at_set ", at_set)
println(keys(eam_arrays))
println("RHO $(rho[a,:]) vals $([rho[a,1]^2, rho[a,2]^2, rho[a,1]*rho[a,2]])")
# println(at_set)
# println("at_set ", at_set)
# println(keys(eam_arrays))
# println("RHO $(rho[a,:]) vals $([rho[a,1]^2, rho[a,2]^2, rho[a,1]*rho[a,2]])")
eam_arrays[at_set][1][ind,:] += [rho[a,1]^2, rho[a,2]^2, rho[a,1]*rho[a,2]]
end
end
Expand Down Expand Up @@ -7892,6 +7894,8 @@ function calc_tb_LV(crys::crystal, database=missing; reference_tbc=missing, verb

rho_th[a1, 1, id] += lag_arr[1] * cut_on
rho_th[a1, 2, id] += lag_arr[2] * cut_on

# println("add rho $a1 dist $dist_a 1 $(lag_arr[1] * cut_on) 2 $(lag_arr[2] * cut_on)")

core!(cham, a1, a2, t1, t2, norb, orbs_arr, DAT_IND_ARR, lag_arr, DAT_ARR, cut_a, H, S, lmn_arr, sym_arr, sym_arrS)
core_onsite!(c_zero, a1, a2, t1, t2, norb, orbs_arr, DAT_IND_ARR_O, lag_arr, DAT_ARR, cut_on, H, lmn_arr, sym_arr, sym_arrS)
Expand All @@ -7906,11 +7910,12 @@ function calc_tb_LV(crys::crystal, database=missing; reference_tbc=missing, verb

if use_eam
rho = sum(rho_th, dims=3)
# println("rho $rho")
for a = 1:crys.nat
t = crys.stypes[a]
d = database[(:eam, t)].datH
temp = d[1]*rho[a,1]^2 + d[2]*rho[a,2]^2 + d[3]*rho[a,1]*rho[a,2]
println("a $a d $d rho[a,:] $(rho[a,:]) temp $temp vals $([rho[a,1]^2, rho[a,2]^2, rho[a,1]*rho[a,2]]) ")
# println("a $a d $d rho[a,:] $(rho[a,:]) temp $temp vals $([rho[a,1]^2, rho[a,2]^2, rho[a,1]*rho[a,2]]) ")
for ox = 1:norb[a]
oxx = orbs_arr[a,ox,1]
H[ oxx, oxx, c_zero] += temp
Expand Down
44 changes: 24 additions & 20 deletions test/test_laguerre_eam.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,31 +21,33 @@ function testEAM()
tbc_list = []
@testset "testing laguerre fitting fake example" begin

atom = :H

@suppress begin

c1 = makecrys([7.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0], [:H], units="Bohr");
c2 = makecrys([8.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0], [:H], units="Bohr");
c3 = makecrys([8.5 0 0; 0 8.5 0; 0 0 7.0], [0 0 0], [:H], units="Bohr");
c4 = makecrys([6.0 0 0; 0 6.0 0; 0 0 6.0], [0 0 0], [:H], units="Bohr");
c1 = makecrys([7.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0], [atom], units="Bohr");
c2 = makecrys([8.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0], [atom], units="Bohr");
c3 = makecrys([8.5 0 0; 0 8.5 0; 0 0 7.0], [0 0 0], [atom], units="Bohr");
c4 = makecrys([6.0 0 0; 0 6.0 0; 0 0 6.0], [0 0 0], [atom], units="Bohr");

c5 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.49 0 0 ], [:H, :H], units="Bohr");
c6 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.45 0 0 ], [:H, :H], units="Bohr");
c7 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.40 0 0 ], [:H, :H], units="Bohr");
c8 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.35 0 0 ], [:H, :H], units="Bohr");
c9 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.47 0 0 ], [:H, :H], units="Bohr");
c5 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.49 0 0 ], [atom, atom], units="Bohr");
c6 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.45 0 0 ], [atom, atom], units="Bohr");
c7 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.40 0 0 ], [atom, atom], units="Bohr");
c8 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.35 0 0 ], [atom, atom], units="Bohr");
c9 = makecrys([14.0 0 0; 0 7.0 0; 0 0 7.0], [0 0 0; 0.47 0 0 ], [atom, atom], units="Bohr");

c4a = makecrys([5.0 0 0; 0 5.0 0; 0 0 5.0], [0 0 0], [:H], units="Bohr");
c4b = makecrys([4.0 0 0; 0 4.0 0; 0 0 4.0], [0 0 0], [:H], units="Bohr");
c4c = makecrys([3.0 0 0; 0 3.0 0; 0 0 3.0], [0 0 0], [:H], units="Bohr");
c4a = makecrys([5.0 0 0; 0 5.0 0; 0 0 5.0], [0 0 0], [atom], units="Bohr");
c4b = makecrys([4.0 0 0; 0 4.0 0; 0 0 4.0], [0 0 0], [atom], units="Bohr");
c4c = makecrys([3.0 0 0; 0 3.0 0; 0 0 3.0], [0 0 0], [atom], units="Bohr");


begin
# if true

database = Dict()
database[(:H, :H)] = ThreeBodyTB.CalcTB.make_coefs(Set([:H, :H]), 2)
database[(:eam, :H)] = ThreeBodyTB.CalcTB.make_coefs(Set([:H, :H]), 0)
database[(:H, :H)].datH .= 0.0
database[(atom, atom)] = ThreeBodyTB.CalcTB.make_coefs(Set([atom, atom]), 2)
database[(:eam, atom)] = ThreeBodyTB.CalcTB.make_coefs(Set([atom, atom]), 0)
# database[(atom, atom)].datH .= 0.0

tbc1 = ThreeBodyTB.CalcTB.calc_tb_LV(c1, database, use_threebody=false, use_threebody_onsite=false, use_eam=true);
tbc2 = ThreeBodyTB.CalcTB.calc_tb_LV(c2, database, use_threebody=false, use_threebody_onsite=false, use_eam=true);
Expand All @@ -62,16 +64,18 @@ function testEAM()

tbc_list = [tbc1, tbc2, tbc3, tbc4, tbc5, tbc6, tbc7, tbc8, tbc9, tbc4a, tbc4b, tbc4c]

newdatabase = ThreeBodyTB.FitTB.do_fitting(tbc_list, fit_threebody=false, fit_threebody_onsite=false,do_plot=false)
newdatabase = ThreeBodyTB.FitTB.do_fitting(tbc_list, fit_threebody=false, fit_threebody_onsite=false,do_plot=false, fit_eam=true)
# newdatabase = ThreeBodyTB.FitTB.do_fitting_recursive(tbc_list, fit_threebody=false, do_plot=false)


# println("newdartabase")
# println(newdatabase[("Li", "Li")])
println(newdatabase[(:H, :H)].datH)
println(database[(:H, :H)].datH)
println(sum(abs.(newdatabase[(:H, :H)].datH .- database[(:H, :H)].datH)))
@test sum(abs.(newdatabase[(:H, :H)].datH .- database[(:H, :H)].datH)) ≤ 1.5e-3
# println(newdatabase[(atom, atom)].datH)
# println(database[(atom, atom)].datH)
# println(sum(abs.(newdatabase[(atom, atom)].datH .- database[(atom, atom)].datH)))
@test sum(abs.(newdatabase[(atom, atom)].datH .- database[(atom, atom)].datH)) ≤ 1.5e-3

@test sum(abs.(newdatabase[(:eam, atom)].datH .- database[(:eam, atom)].datH)) ≤ 1.5e-3

end

Expand Down

0 comments on commit dd29947

Please sign in to comment.