-
Notifications
You must be signed in to change notification settings - Fork 0
/
read_h5.irp.f
133 lines (107 loc) · 4.51 KB
/
read_h5.irp.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
program main
use iso_c_binding
use map_module
implicit none
character (len = 100) :: h5path
character (len = 1) :: basis_id_str
integer :: err
integer :: n_bielec_int, n_hcore_int, n_to_read, o_num_h5
integer :: i_int, chuck_size_cur
real(integral_kind), allocatable :: buffer_values(:)
integer(key_kind), allocatable :: buffer_i(:)
integer(key_kind), allocatable :: buffer_i_hcore(:,:)
double precision, allocatable :: buffer_values_hcore_h5(:)
real(integral_kind), allocatable :: buffer_values_hcore(:,:)
character (len = 20) :: chunk_size_str
integer :: chunk_size
integer :: basis_id = 1
integer :: i,j, a,b
PROVIDE ezfio_filename
! _
! |_) _. ._ _. ._ _
! | (_| | (_| | | |
!
CALL GET_ENVIRONMENT_VARIABLE("h5_int_path", h5path, status = err)
if ( err .ne. 0) then
h5path = "qp_in.qp.h5"
endif
! basis_id == 0 -> MO ; basis_id != 0 -> AO. Default basis_id == 0
CALL GET_ENVIRONMENT_VARIABLE("h5_basis_id", basis_id_str, status = err)
if ( err .ne. 0) then
basis_id = 0
else
read(basis_id_str, '(i)') basis_id
endif
CALL GET_ENVIRONMENT_VARIABLE("h5_chunk_size", chunk_size_str, status = err)
if ( err .ne. 0) then
chunk_size = 100*1000
else
read(chunk_size_str, '(i)') chunk_size
endif
call get_param_h5_int(trim(h5path)//c_null_char, n_bielec_int, n_hcore_int, basis_id)
if ( basis_id.eq. 0 ) then
o_num_h5 = mo_num
else
o_num_h5 = ao_num
endif
! _
! |\/| _ ._ _ _ | _ _ / |_| / _ ._ _ \
! | | (_) | | (_) (/_ | (/_ (_ | | | \_ (_) | (/_ |
! \ /
allocate(buffer_i_hcore(2,n_hcore_int), buffer_values_hcore_h5(n_hcore_int), buffer_values_hcore(o_num_h5,o_num_h5))
call read_hcore_h5(trim(h5path)//c_null_char,buffer_i_hcore, buffer_values_hcore_h5, basis_id)
buffer_values_hcore(:,:) = 0.d0 ! Assume that we may read sparse
! reprensation, but we store dense one
!$OMP SIMD PRIVATE(a,b)
do j = 1, n_hcore_int
a = buffer_i_hcore(1,j)
b = buffer_i_hcore(2,j)
buffer_values_hcore(a+1,b+1) =buffer_values_hcore_h5(j)
end do
if (basis_id == 0) then
call ezfio_set_mo_one_e_ints_mo_one_e_integrals(buffer_values_hcore)
call ezfio_set_mo_one_e_ints_io_mo_one_e_integrals("Read")
else
call ezfio_set_ao_one_e_ints_ao_one_e_integrals(buffer_values_hcore)
call ezfio_set_ao_one_e_ints_io_ao_one_e_integrals("Read")
endif
deallocate(buffer_i_hcore, buffer_values_hcore_h5, buffer_values_hcore)
! _ _
! |_) o |_ | _ _
! |_) | |_ | (/_ (_
!
allocate(buffer_values(chunk_size), buffer_i(chunk_size))
print*, 'Integral two-e: number of to read: ', n_bielec_int
n_to_read = n_bielec_int
i_int = 0
do while ( (n_to_read .ne. 0) .and. (i_int .ne. n_bielec_int))
chuck_size_cur = min(chunk_size, n_to_read)
call read_bielec_h5(trim(h5path)//c_null_char, buffer_i,buffer_values,i_int, chuck_size_cur, basis_id)
if (basis_id == 0) then
call insert_into_mo_integrals_map(chuck_size_cur, buffer_i,buffer_values)
else
call insert_into_ao_integrals_map(chuck_size_cur, buffer_i,buffer_values)
endif
i_int = i_int + chuck_size_cur
n_to_read = n_to_read - chuck_size_cur
!print *, 'Bielect integral: Commited ', i_int, 'integral ', n_to_read, 'to go'
end do
print*, 'Sorting the integrals_map'
if (basis_id == 0) then
call map_sort(mo_integrals_map)
!call map_unique(mo_integrals_map) ! Assume key are already unique
else
call map_sort(ao_integrals_map)
endif
if (basis_id == 0) then
call map_sort(mo_integrals_map)
call map_save_to_disk(trim(ezfio_filename)//'/work/mo_ints',mo_integrals_map)
call ezfio_set_mo_two_e_ints_io_mo_two_e_integrals('Read')
else
call map_sort(ao_integrals_map)
call map_save_to_disk(trim(ezfio_filename)//'/work/ao_ints',ao_integrals_map)
call ezfio_set_ao_two_e_ints_io_ao_two_e_integrals('Read')
endif
!Now n_to_read == 0 and i_int = n_bielec_int)
deallocate(buffer_values,buffer_i)
end program main