Lots of minor fixes

This commit is contained in:
Anthony Scemama 2015-04-09 21:46:28 +02:00
parent dd1a51af63
commit 102bbb0b4f
16 changed files with 121 additions and 184 deletions

View File

@ -32,6 +32,12 @@ full_ci
do_pt2_end True
var_pt2_ratio 0.75
cas_sd
n_det_max_cas_sd 100000
pt2_max 1.e-4
do_pt2_end True
var_pt2_ratio 0.75
all_singles
n_det_max_fci 50000
pt2_max 1.e-8

View File

@ -176,7 +176,7 @@ def get_dict_config_file(config_file_path, module_lower):
d[pvd][option] = d_default[option]
# If interface is output we need a default value information
if d[pvd]["interface"] == "output":
if d[pvd]["interface"] == "input":
try:
d[pvd]["default"] = config_file.get(section, "default")
except ConfigParser.NoOptionError:
@ -210,8 +210,8 @@ def create_ezfio_provider(dict_ezfio_cfg):
ez_p.set_ezfio_dir(dict_info['ezfio_dir'])
ez_p.set_ezfio_name(dict_info['ezfio_name'])
ez_p.set_default(dict_info['default'])
ez_p.set_output("output_%s" % dict_info['ezfio_dir'])
dict_code_provider[provider_name] = str(ez_p)
return dict_code_provider
@ -239,7 +239,7 @@ def save_ezfio_provider(path_head, dict_code_provider):
"! from file {0}/EZFIO.cfg\n".format(path_head) + \
"\n"
for provider_name, code in dict_code_provider.iteritems():
output += code + "\n"
output += str(code) + "\n"
if output != old_output:
with open(path, "w") as f:

View File

@ -112,16 +112,19 @@ END_PROVIDER
break
v = buffer[1]
name = self.name
true = True
false= False
try:
v_eval = eval(v)
except:
v = "call ezfio_get_%(v)s(%(name)s)"%locals()
else:
if type(v_eval) == bool:
v = '.%s.'%(v)
elif type(v_eval) == float:
v = v.replace('e','d')
v = v.replace('E','D')
v = "%(name)s = %(v)s"%locals()
except:
v = "call ezfio_get_%(v)s(%(name)s)"%locals()
self.default = v

View File

@ -53,7 +53,7 @@ class H_apply(object):
!$OMP N_elec_in_key_hole_2,ia_ja_pairs) &
!$OMP SHARED(key_in,N_int,elec_num_tab,mo_tot_num, &
!$OMP hole_1, particl_1, hole_2, particl_2, &
!$OMP elec_alpha_num,i_generator)"""
!$OMP elec_alpha_num,i_generator) FIRSTPRIVATE(iproc)"""
s["omp_end_parallel"] = "!$OMP END PARALLEL"
s["omp_master"] = "!$OMP MASTER"
s["omp_end_master"] = "!$OMP END MASTER"

View File

@ -53,16 +53,68 @@ def write_ezfioFile(res,filename):
basis = res.uncontracted_basis
geom = res.geometry
res.clean_contractions()
# AO Basis
import string
at = []
num_prim = []
magnetic_number = []
angular_number = []
power_x = []
power_y = []
power_z = []
coefficient = []
exponent = []
res.convert_to_cartesian()
for b in res.basis:
c = b.center
for i,atom in enumerate(res.geometry):
if atom.coord == c:
at.append(i+1)
num_prim.append(len(b.prim))
s = b.sym
power_x.append( string.count(s,"x") )
power_y.append( string.count(s,"y") )
power_z.append( string.count(s,"z") )
coefficient.append( b.coef )
exponent.append( [ p.expo for p in b.prim ] )
ezfio.set_ao_basis_ao_num(len(res.basis))
ezfio.set_ao_basis_ao_nucl(at)
ezfio.set_ao_basis_ao_prim_num(num_prim)
ezfio.set_ao_basis_ao_power(power_x+power_y+power_z)
prim_num_max = ezfio.get_ao_basis_ao_prim_num_max()
len_res_basis = len(res.basis)
for i in range(len(res.basis)):
coefficient[i] += [ 0. for j in range(len(coefficient[i]),prim_num_max) ]
exponent[i] += [ 0. for j in range(len(exponent[i]),prim_num_max) ]
coefficient = reduce(lambda x, y: x+y, coefficient, [])
exponent = reduce(lambda x, y: x+y, exponent, [])
coef = []
expo = []
for i in range(prim_num_max):
for j in range(i,len(coefficient),prim_num_max):
coef.append ( coefficient[j] )
expo.append ( exponent[j] )
ezfio.set_ao_basis_ao_coef(coef)
ezfio.set_ao_basis_ao_expo(expo)
ezfio.set_ao_basis_ao_basis("Read by resultsFile")
# MO
MoTag = res.determinants_mo_type
ezfio.set_mo_basis_mo_label('Orthonormalized')
MO_type = MoTag
allMOs = res.uncontracted_mo_sets[MO_type]
allMOs = res.mo_sets[MO_type]
closed = [ (allMOs[i].eigenvalue,i) for i in res.closed_mos ]
active = [ (allMOs[i].eigenvalue,i) for i in res.active_mos ]
virtual =[ (allMOs[i].eigenvalue,i) for i in res.virtual_mos ]
try:
closed = [ (allMOs[i].eigenvalue,i) for i in res.closed_mos ]
active = [ (allMOs[i].eigenvalue,i) for i in res.active_mos ]
virtual =[ (allMOs[i].eigenvalue,i) for i in res.virtual_mos ]
except:
closed = []
virtual = []
active = [ (allMOs[i].eigenvalue,i) for i in range(len(allMOs)) ]
# closed.sort()
# active.sort()
@ -111,117 +163,6 @@ def write_ezfioFile(res,filename):
while len(MoMatrix) < len(MOs[0].vector)**2:
MoMatrix.append(0.)
ezfio.set_mo_basis_mo_tot_num(mo_tot_num)
ezfio.set_mo_basis_mo_occ(OccNum)
res.clean_contractions()
# AO Basis
import string
at = []
num_prim = []
magnetic_number = []
angular_number = []
power_x = []
power_y = []
power_z = []
coefficient = []
exponent = []
res.convert_to_cartesian()
for b in res.basis:
c = b.center
for i,atom in enumerate(res.geometry):
if atom.coord == c:
at.append(i+1)
num_prim.append(len(b.prim))
s = b.sym
power_x.append( string.count(s,"x") )
power_y.append( string.count(s,"y") )
power_z.append( string.count(s,"z") )
coefficient.append( b.coef )
exponent.append( [ p.expo for p in b.prim ] )
ezfio.set_ao_basis_ao_num(len(res.basis))
ezfio.set_ao_basis_ao_nucl(at)
ezfio.set_ao_basis_ao_prim_num(num_prim)
ezfio.set_ao_basis_ao_power(power_x+power_y+power_z)
prim_num_max = ezfio.get_ao_basis_ao_prim_num_max()
len_res_basis = len(res.basis)
for i in range(len(res.basis)):
coefficient[i] += [ 0. for j in range(len(coefficient[i]),prim_num_max) ]
exponent[i] += [ 0. for j in range(len(exponent[i]),prim_num_max) ]
coefficient = reduce(lambda x, y: x+y, coefficient, [])
exponent = reduce(lambda x, y: x+y, exponent, [])
coef = []
expo = []
for i in range(prim_num_max):
for j in range(i,len(coefficient),prim_num_max):
coef.append ( coefficient[j] )
expo.append ( exponent[j] )
ezfio.set_ao_basis_ao_coef(coef)
ezfio.set_ao_basis_ao_expo(expo)
ezfio.set_ao_basis_ao_basis("Read by resultsFile")
# Apply threshold to determinants
if len(res.determinants) == 1:
sorted_determinants = [ (-1.,1.,res.determinants[0]) ]
else:
sorted_determinants = []
for i,j in zip(res.det_coefficients[0],res.determinants):
sorted_determinants.append((-abs(i),i,j))
sorted_determinants.sort()
norm = 0.0
for length, (a,b,c) in enumerate(sorted_determinants):
if -a < det_threshold:
length -=1
break
norm += a**2
norm = sqrt(norm)
length += 1
for i in xrange(length):
a = sorted_determinants[i]
sorted_determinants[i] = (a[0],a[1]/norm,a[2])
sorted_determinants = sorted_determinants[:length]
# MOs
mo_tot_num = len(res.mo_sets[MoTag])
closed_mos = res.closed_mos
active_mos = res.active_mos
virtual_mos = res.virtual_mos
to_remove = []
to_add = []
for i in active_mos:
found = False
for (a,b,c) in sorted_determinants:
if i in c['alpha']+c['beta']:
found = True
break
if not found:
to_remove.append(i)
to_add.append(i)
virtual_mos = to_add + virtual_mos
for i in active_mos:
always = True
for (a,b,c) in sorted_determinants:
if not (i in c['alpha'] and i in c['beta']):
always = False
break
if always:
to_remove.append(i)
closed_mos.append(i)
for i in to_remove:
active_mos.remove(i)
MOindices = closed_mos + active_mos + virtual_mos
while len(MOindices) < mo_tot_num:
MOindices.append(len(MOindices))
MOmap = list(MOindices)
for i in range(len(MOindices)):
MOmap[i] = MOindices.index(i)
ezfio.set_mo_basis_mo_tot_num(mo_tot_num)
mo = []
for i in MOindices:
mo.append(res.mo_sets[MoTag][i])
@ -234,35 +175,11 @@ def write_ezfioFile(res,filename):
while len(mo) < mo_tot_num:
mo.append(newmo)
Energies = [ m.eigenvalue for m in mo ]
if res.occ_num is not None:
OccNum = []
for i in MOindices:
OccNum.append(res.occ_num[MoTag][i])
while len(OccNum) < mo_tot_num:
OccNum.append(0.)
ezfio.set_mo_basis_mo_occ(OccNum)
cls = [ "v" for i in mo ]
for i in closed_mos:
cls[MOmap[i]] = 'c'
for i in active_mos:
cls[MOmap[i]] = 'a'
sym0 = [ i.sym for i in res.mo_sets[MoTag] ]
sym = [ i.sym for i in res.mo_sets[MoTag] ]
for i in xrange(len(sym)):
sym[MOmap[i]] = sym0[i]
MoMatrix = []
for m in mo:
for coef in m.vector:
MoMatrix.append(coef)
while len(MoMatrix) < len(mo[0].vector)**2:
MoMatrix.append(0.)
ezfio.set_mo_basis_mo_tot_num(mo_tot_num)
ezfio.set_mo_basis_mo_occ(OccNum)
ezfio.set_mo_basis_mo_coef(MoMatrix)
del MoMatrix

View File

@ -368,11 +368,11 @@ subroutine get_mo_bielec_integrals_existing_ik(j,l,sze,out_array,map)
enddo
logical :: integral_is_in_map
if (cache_key_kind == 8) then
if (key_kind == 8) then
call i8radix_sort(hash,iorder,kk,-1)
else if (cache_key_kind == 4) then
else if (key_kind == 4) then
call iradix_sort(hash,iorder,kk,-1)
else if (cache_key_kind == 2) then
else if (key_kind == 2) then
call i2radix_sort(hash,iorder,kk,-1)
endif

View File

@ -123,7 +123,6 @@ BEGIN_PROVIDER [ integer(bit_kind), generators_bitmask, (N_int,2,6,N_generators_
call ezfio_has_bitmasks_generators(exists)
if (exists) then
print*,'EXIST !!'
call ezfio_get_bitmasks_generators(generators_bitmask)
else
integer :: k, ispin
@ -181,7 +180,6 @@ BEGIN_PROVIDER [ integer(bit_kind), cas_bitmask, (N_int,2,N_cas_bitmask) ]
PROVIDE ezfio_filename
call ezfio_has_bitmasks_cas(exists)
print*,'exists = ',exists
if (exists) then
call ezfio_get_bitmasks_cas(cas_bitmask)
else

View File

@ -2,4 +2,8 @@
CAS_SD_selected Module
======================
Selected CAS + SD module
Selected CAS + SD module.
1) Set the different MO classes using the ``qp_set_mo_class`` command
2) Run the selected CAS+SD program

View File

@ -11,12 +11,12 @@ program full_ci
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > n_det_max_fci) then
if (N_det > n_det_max_cas_sd) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = n_det_max_fci
N_det = n_det_max_cas_sd
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
@ -28,17 +28,17 @@ program full_ci
print *, '-----'
endif
do while (N_det < n_det_max_fci.and.maxval(abs(pt2(1:N_st))) > pt2_max)
do while (N_det < n_det_max_cas_sd.and.maxval(abs(pt2(1:N_st))) > pt2_max)
call H_apply_CAS_SD(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > n_det_max_fci) then
if (N_det > n_det_max_cas_sd) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = n_det_max_fci
N_det = n_det_max_cas_sd
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI

View File

@ -11,12 +11,12 @@ program full_ci
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > n_det_max_fci) then
if (N_det > n_det_max_cas_sd) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = n_det_max_fci
N_det = n_det_max_cas_sd
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
@ -28,17 +28,17 @@ program full_ci
print *, '-----'
endif
do while (N_det < n_det_max_fci.and.maxval(abs(pt2(1:N_st))) > pt2_max)
do while (N_det < n_det_max_cas_sd.and.maxval(abs(pt2(1:N_st))) > pt2_max)
call H_apply_CAS_SD_selected(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > n_det_max_fci) then
if (N_det > n_det_max_cas_sd) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = n_det_max_fci
N_det = n_det_max_cas_sd
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI

View File

@ -1,4 +1,4 @@
subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_generator, iproc $parameters )
subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_generator, iproc_in $parameters )
use omp_lib
use bitmasks
implicit none
@ -14,7 +14,7 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene
integer(bit_kind),allocatable :: keys_out(:,:,:)
integer(bit_kind), intent(in) :: hole_1(N_int,2), particl_1(N_int,2)
integer(bit_kind), intent(in) :: hole_2(N_int,2), particl_2(N_int,2)
integer, intent(in) :: iproc
integer, intent(in) :: iproc_in
integer(bit_kind), allocatable :: hole_save(:,:)
integer(bit_kind), allocatable :: key(:,:),hole(:,:), particle(:,:)
integer(bit_kind), allocatable :: hole_tmp(:,:), particle_tmp(:,:)
@ -30,6 +30,7 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene
integer, allocatable :: ia_ja_pairs(:,:,:)
integer, allocatable :: ib_jb_pairs(:,:)
double precision :: diag_H_mat_elem
integer :: iproc
integer(omp_lock_kind), save :: lck, ifirst=0
if (ifirst == 0) then
!$ call omp_init_lock(lck)
@ -38,12 +39,13 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene
logical :: check_double_excitation
check_double_excitation = .True.
iproc = iproc_in
$initialization
$omp_parallel
!$ iproc = omp_get_thread_num()
allocate (keys_out(N_int,2,size_max), hole_save(N_int,2), &
key(N_int,2),hole(N_int,2), particle(N_int,2), hole_tmp(N_int,2),&
particle_tmp(N_int,2), occ_particle(N_int*bit_kind_size,2), &
@ -248,7 +250,7 @@ subroutine $subroutine_diexc(key_in, hole_1,particl_1, hole_2, particl_2, i_gene
$finalization
end
subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $parameters )
subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc_in $parameters )
use omp_lib
use bitmasks
implicit none
@ -262,7 +264,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param
integer ,intent(in) :: i_generator
integer(bit_kind),intent(in) :: key_in(N_int,2)
integer(bit_kind),intent(in) :: hole_1(N_int,2), particl_1(N_int,2)
integer, intent(in) :: iproc
integer, intent(in) :: iproc_in
integer(bit_kind),allocatable :: keys_out(:,:,:)
integer(bit_kind),allocatable :: hole_save(:,:)
integer(bit_kind),allocatable :: key(:,:),hole(:,:), particle(:,:)
@ -281,8 +283,11 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param
logical, allocatable :: array_pairs(:,:)
double precision :: diag_H_mat_elem
integer(omp_lock_kind), save :: lck, ifirst=0
integer :: iproc
logical :: check_double_excitation
iproc = iproc_in
check_double_excitation = .True.
$check_double_excitation
@ -295,6 +300,7 @@ subroutine $subroutine_monoexc(key_in, hole_1,particl_1,i_generator,iproc $param
$initialization
$omp_parallel
!$ iproc = omp_get_thread_num()
allocate (keys_out(N_int,2,size_max), hole_save(N_int,2), &
key(N_int,2),hole(N_int,2), particle(N_int,2), hole_tmp(N_int,2),&
particle_tmp(N_int,2), occ_particle(N_int*bit_kind_size,2), &
@ -396,7 +402,8 @@ subroutine $subroutine($params_main)
integer :: iproc
$initialization
PROVIDE H_apply_buffer_allocated mo_bielec_integrals_in_map
PROVIDE H_apply_buffer_allocated mo_bielec_integrals_in_map psi_det_generators psi_coef_generators
nmax = mod( N_det_generators,nproc )
@ -406,6 +413,7 @@ subroutine $subroutine($params_main)
call wall_time(wall_0)
iproc = 0
allocate( mask(N_int,2,6) )
do i_generator=1,nmax
@ -443,12 +451,12 @@ subroutine $subroutine($params_main)
call $subroutine_diexc(psi_det_generators(1,1,i_generator), &
mask(1,1,d_hole1), mask(1,1,d_part1), &
mask(1,1,d_hole2), mask(1,1,d_part2), &
i_generator, 0 $params_post)
i_generator, iproc $params_post)
endif
if($do_mono_excitations)then
call $subroutine_monoexc(psi_det_generators(1,1,i_generator), &
mask(1,1,s_hole ), mask(1,1,s_part ), &
i_generator, 0 $params_post)
i_generator, iproc $params_post)
endif
call wall_time(wall_1)
$printout_always
@ -463,7 +471,6 @@ subroutine $subroutine($params_main)
!$OMP PARALLEL DEFAULT(SHARED) &
!$OMP PRIVATE(i_generator,wall_1,wall_0,ispin,k,mask,iproc)
call wall_time(wall_0)
iproc = 0
!$ iproc = omp_get_thread_num()
allocate( mask(N_int,2,6) )
!$OMP DO SCHEDULE(dynamic,1)

View File

@ -1,42 +1,42 @@
[N_det_max_fci]
type: Det_number_max
doc: Max number of determinants in the wave function
interface: output
interface: input
default: 10000
[N_det_max_fci_property]
type: Det_number_max
doc: Max number of determinants in the wave function when you select for a given property
interface: output
interface: input
default: 10000
[do_pt2_end]
type: logical
doc: If true, compute the PT2 at the end of the selection
interface: output
interface: input
default: true
[PT2_max]
type: PT2_energy
doc: The selection process stops when the largest PT2 (for all the state is lower
than pt2_max in absolute value
interface: output
interface: input
default: 0.0001
[var_pt2_ratio]
type: Normalized_float
doc: The selection process stops when the energy ratio variational/(variational+PT2)
is equal to var_pt2_ratio
interface: output
interface: input
default: 0.75
[energy]
type: double precision
doc: "Calculated Full CI energy"
interface: input
interface: output
[energy_pt2]
type: double precision
doc: "Calculated Full CI energy"
interface: input
interface: output

View File

@ -62,7 +62,6 @@ END_PROVIDER
psi_det_generators(k,2,m) = psi_det(k,2,i)
enddo
psi_coef_generators(m,:) = psi_coef(m,:)
! call debug_det(psi_det_generators(1,1,m),N_int)
endif
enddo

View File

@ -1,2 +1,2 @@
AOs BiInts Bitmask CAS_SD_selected Dets Electrons Ezfio_files Generators_CAS Generators_full Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full Utils
AOs BiInts Bitmask CAS_SD Dets Electrons Ezfio_files Generators_CAS Generators_full Hartree_Fock MOGuess MonoInts MOs Nuclei Output Perturbation Properties Selectors_full Utils

View File

@ -1 +1 @@
AOs BiInts Bitmask CID CID_SC2_selected CID_selected CIS CISD CISD_selected CISD_SC2_selected Dets Electrons Ezfio_files Full_CI Generators_full Hartree_Fock MOGuess MonoInts MOs MP2 Nuclei Output Selectors_full Utils Molden FCIdump Generators_CAS CAS_SD_selected DDCI_selected MRCC
AOs BiInts Bitmask CID CID_SC2_selected CID_selected CIS CISD CISD_selected CISD_SC2_selected Dets Electrons Ezfio_files Full_CI Generators_full Hartree_Fock MOGuess MonoInts MOs MP2 Nuclei Output Selectors_full Utils Molden FCIdump Generators_CAS CAS_SD DDCI_selected MRCC

View File

@ -451,3 +451,6 @@ subroutine set_zero_extra_diag(i1,i2,matrix,lda,m)
end