From 39366f7ec1df4635cd3fe4bffc6a232a0149d75b Mon Sep 17 00:00:00 2001 From: Anthony Scemama Date: Tue, 14 Jun 2016 13:59:05 +0200 Subject: [PATCH] mmap works for integrals --- plugins/Mmap/fortran_mmap.c | 72 ++++++++++ plugins/Mmap/mmap.f90 | 69 ++++++++++ src/Integrals_Bielec/NEEDED_CHILDREN_MODULES | 3 +- src/Integrals_Bielec/ao_bi_integrals.irp.f | 7 +- src/Integrals_Bielec/mo_bi_integrals.irp.f | 10 +- src/Utils/NEEDED_CHILDREN_MODULES | 2 +- src/Utils/map_module.f90 | 133 ++++++++++++++++++- 7 files changed, 279 insertions(+), 17 deletions(-) create mode 100644 plugins/Mmap/fortran_mmap.c create mode 100644 plugins/Mmap/mmap.f90 diff --git a/plugins/Mmap/fortran_mmap.c b/plugins/Mmap/fortran_mmap.c new file mode 100644 index 00000000..2748dcba --- /dev/null +++ b/plugins/Mmap/fortran_mmap.c @@ -0,0 +1,72 @@ +#include +#include +#include +#include +#include +#include +#include + + +void* mmap_fortran(char* filename, size_t bytes, int* file_descr, int read_only) +{ + int i; + int fd; + int result; + void* map; + + if (read_only == 1) + { + fd = open(filename, O_RDONLY, (mode_t)0600); + if (fd == -1) { + printf("%s:\n", filename); + perror("Error opening mmap file for reading"); + exit(EXIT_FAILURE); + } + map = mmap(0, bytes, PROT_READ, MAP_SHARED, fd, 0); + } + else + { + fd = open(filename, O_RDWR | O_CREAT, (mode_t)0600); + if (fd == -1) { + printf("%s:\n", filename); + perror("Error opening mmap file for writing"); + exit(EXIT_FAILURE); + } + + result = lseek(fd, bytes, SEEK_SET); + if (result == -1) { + close(fd); + printf("%s:\n", filename); + perror("Error calling lseek() to stretch the file"); + exit(EXIT_FAILURE); + } + + result = write(fd, "", 1); + if (result != 1) { + close(fd); + printf("%s:\n", filename); + perror("Error writing last byte of the file"); + exit(EXIT_FAILURE); + } + + map = mmap(0, bytes, PROT_READ | PROT_WRITE, MAP_SHARED, fd, 0); + } + + if (map == MAP_FAILED) { + close(fd); + printf("%s:\n", filename); + perror("Error mmapping the file"); + exit(EXIT_FAILURE); + } + + *file_descr = fd; + return map; +} + +void munmap_fortran(size_t bytes, int fd, void* map) +{ + if (munmap(map, bytes) == -1) { + perror("Error un-mmapping the file"); + } + close(fd); +} diff --git a/plugins/Mmap/mmap.f90 b/plugins/Mmap/mmap.f90 new file mode 100644 index 00000000..ce33e301 --- /dev/null +++ b/plugins/Mmap/mmap.f90 @@ -0,0 +1,69 @@ +module mmap_module + + use iso_c_binding + + interface + + ! File descriptors + ! ---------------- + + type(c_ptr) function c_mmap_fortran(filename, length, fd, read_only) bind(c,name='mmap_fortran') + use iso_c_binding + character(c_char), intent(in) :: filename(*) + integer(c_size_t), intent(in), value :: length + integer(c_int), intent(out) :: fd + integer(c_int), intent(in), value :: read_only + end function + + subroutine c_munmap(length, fd, map) bind(c,name='munmap_fortran') + use iso_c_binding + integer(c_size_t), intent(in), value :: length + integer(c_int), intent(in), value :: fd + type(c_ptr), intent(in), value :: map + end subroutine + + end interface + + contains + + subroutine mmap(filename, shape, bytes, fd, read_only, map) + use iso_c_binding + implicit none + character*(*), intent(in) :: filename ! Name of the mapped file + integer*8, intent(in) :: shape(:) ! Shape of the array to map + integer, intent(in) :: bytes ! Number of bytes per element + logical, intent(in) :: read_only ! If true, mmap is read-only + integer, intent(out) :: fd ! File descriptor + type(c_ptr), intent(out) :: map ! C Pointer + + integer(c_long) :: length + integer(c_int) :: fd_ + + length = PRODUCT( shape(:) ) * bytes + if (read_only) then + map = c_mmap_fortran( trim(filename)//char(0), length, fd_, 1) + else + map = c_mmap_fortran( trim(filename)//char(0), length, fd_, 0) + endif + fd = fd_ + end subroutine + + subroutine munmap(shape, bytes, fd, map) + use iso_c_binding + implicit none + integer*8, intent(in) :: shape(:) ! Shape of the array to map + integer, intent(in) :: bytes ! Number of bytes per element + integer, intent(in) :: fd ! File descriptor + type(c_ptr), intent(in) :: map ! C pointer + + integer(c_long) :: length + integer(c_int) :: fd_ + + length = PRODUCT( shape(:) ) * bytes + fd_ = fd + call c_munmap( length, fd_, map) + end + +end module mmap_module + + diff --git a/src/Integrals_Bielec/NEEDED_CHILDREN_MODULES b/src/Integrals_Bielec/NEEDED_CHILDREN_MODULES index 152711f3..c17850c6 100644 --- a/src/Integrals_Bielec/NEEDED_CHILDREN_MODULES +++ b/src/Integrals_Bielec/NEEDED_CHILDREN_MODULES @@ -1 +1,2 @@ -Pseudo Bitmask ZMQ +Pseudo Bitmask ZMQ Mmap + diff --git a/src/Integrals_Bielec/ao_bi_integrals.irp.f b/src/Integrals_Bielec/ao_bi_integrals.irp.f index b7c75fb8..e9700c76 100644 --- a/src/Integrals_Bielec/ao_bi_integrals.irp.f +++ b/src/Integrals_Bielec/ao_bi_integrals.irp.f @@ -353,11 +353,10 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_in_map ] if (read_ao_integrals) then integer :: load_ao_integrals print*,'Reading the AO integrals' - if (load_ao_integrals(trim(ezfio_filename)//'/work/ao_integrals.bin') == 0) then + call map_load_from_disk(trim(ezfio_filename)//'/work/ao_ints',ao_integrals_map) print*, 'AO integrals provided' ao_bielec_integrals_in_map = .True. return - endif endif print*, 'Providing the AO integrals' @@ -402,8 +401,10 @@ BEGIN_PROVIDER [ logical, ao_bielec_integrals_in_map ] print*, ' wall time :',wall_2 - wall_1, 's ( x ', (cpu_2-cpu_1)/(wall_2-wall_1+tiny(1.d0)), ' )' ao_bielec_integrals_in_map = .True. + if (write_ao_integrals) then - call dump_ao_integrals(trim(ezfio_filename)//'/work/ao_integrals.bin') + call ezfio_set_work_empty(.False.) + call map_save_to_disk(trim(ezfio_filename)//'/work/ao_ints',ao_integrals_map) call ezfio_set_integrals_bielec_disk_access_ao_integrals("Read") endif diff --git a/src/Integrals_Bielec/mo_bi_integrals.irp.f b/src/Integrals_Bielec/mo_bi_integrals.irp.f index 4d471545..868a18d7 100644 --- a/src/Integrals_Bielec/mo_bi_integrals.irp.f +++ b/src/Integrals_Bielec/mo_bi_integrals.irp.f @@ -30,10 +30,9 @@ BEGIN_PROVIDER [ logical, mo_bielec_integrals_in_map ] if (read_mo_integrals) then integer :: load_mo_integrals print*,'Reading the MO integrals' - if (load_mo_integrals(trim(ezfio_filename)//'/work/mo_integrals.bin') == 0) then - print*, 'MO integrals provided' - return - endif + call map_load_from_disk(trim(ezfio_filename)//'/work/mo_ints',mo_integrals_map) + print*, 'MO integrals provided' + return endif call add_integrals_to_map(full_ijkl_bitmask_4) @@ -299,7 +298,8 @@ subroutine add_integrals_to_map(mask_ijkl) print*,' wall time :',wall_2 - wall_1, 's ( x ', (cpu_2-cpu_1)/(wall_2-wall_1), ')' if (write_mo_integrals) then - call dump_mo_integrals(trim(ezfio_filename)//'/work/mo_integrals.bin') + call ezfio_set_work_empty(.False.) + call map_save_to_disk(trim(ezfio_filename)//'/work/mo_ints',mo_integrals_map) call ezfio_set_integrals_bielec_disk_access_mo_integrals("Read") endif diff --git a/src/Utils/NEEDED_CHILDREN_MODULES b/src/Utils/NEEDED_CHILDREN_MODULES index 8b137891..993ebdae 100644 --- a/src/Utils/NEEDED_CHILDREN_MODULES +++ b/src/Utils/NEEDED_CHILDREN_MODULES @@ -1 +1 @@ - +Mmap diff --git a/src/Utils/map_module.f90 b/src/Utils/map_module.f90 index 47adc83e..38af0830 100644 --- a/src/Utils/map_module.f90 +++ b/src/Utils/map_module.f90 @@ -30,8 +30,8 @@ module map_module integer*8, parameter :: map_mask = ibset(0_8,15)-1_8 type cache_map_type - integer(cache_key_kind), pointer :: key(:) real(integral_kind), pointer :: value(:) + integer(cache_key_kind), pointer :: key(:) logical :: sorted integer(cache_map_size_kind) :: map_size integer(cache_map_size_kind) :: n_elements @@ -39,10 +39,14 @@ module map_module end type cache_map_type type map_type - type(cache_map_type), allocatable :: map(:) + type(cache_map_type), pointer :: map(:) + real(integral_kind), pointer :: consolidated_value(:) + integer(cache_key_kind), pointer :: consolidated_key(:) + integer*8, pointer :: consolidated_idx(:) + logical :: sorted + logical :: consolidated integer(map_size_kind) :: map_size integer(map_size_kind) :: n_elements - logical :: sorted integer(omp_lock_kind) :: lock end type map_type @@ -92,6 +96,7 @@ subroutine map_init(map,keymax) map%n_elements = 0_8 map%map_size = ishft(keymax,map_shift) + map%consolidated = .False. allocate(map%map(0_8:map%map_size),stat=err) if (err /= 0) then @@ -655,12 +660,10 @@ subroutine search_key_big_interval(key,X,sze,idx,ibegin_in,iend_in) idx = ibegin if (min(iend_in,sze) > ibegin+16) then iend = ibegin+16 - !DIR$ VECTOR ALIGNED do while (cache_key > X(idx)) idx = idx+1 end do else - !DIR$ VECTOR ALIGNED do while (cache_key > X(idx)) idx = idx+1 if (idx /= iend) then @@ -768,13 +771,11 @@ subroutine search_key_value_big_interval(key,value,X,Y,sze,idx,ibegin_in,iend_in value = Y(idx) if (min(iend_in,sze) > ibegin+16) then iend = ibegin+16 - !DIR$ VECTOR ALIGNED do while (cache_key > X(idx)) idx = idx+1 value = Y(idx) end do else - !DIR$ VECTOR ALIGNED do while (cache_key > X(idx)) idx = idx+1 value = Y(idx) @@ -853,3 +854,121 @@ subroutine get_cache_map(map,map_idx,keys,values,n_elements) enddo end + + +subroutine map_save_to_disk(filename,map) + use map_module + use mmap_module + implicit none + character*(*), intent(in) :: filename + type(map_type), intent(inout) :: map + type(c_ptr) :: c_pointer(3) + integer :: fd(3) + integer*8 :: i,k + integer :: j + + + + if (map % consolidated) then + stop 'map already consolidated' + endif + + call mmap(trim(filename)//'_consolidated_idx', (/ map % map_size + 2_8 /), 8, fd(1), .False., c_pointer(1)) + call c_f_pointer(c_pointer(1),map % consolidated_idx, (/ map % map_size +2_8/)) + + call mmap(trim(filename)//'_consolidated_key', (/ map % n_elements /), cache_key_kind, fd(2), .False., c_pointer(2)) + call c_f_pointer(c_pointer(2),map % consolidated_key, (/ map % n_elements /)) + + call mmap(trim(filename)//'_consolidated_value', (/ map % n_elements /), integral_kind, fd(3), .False., c_pointer(3)) + call c_f_pointer(c_pointer(3),map % consolidated_value, (/ map % n_elements /)) + + if (.not.associated(map%consolidated_key)) then + stop 'cannot consolidate map : consolidated_key not associated' + endif + + if (.not.associated(map%consolidated_value)) then + stop 'cannot consolidate map : consolidated_value not associated' + endif + + if (.not.associated(map%consolidated_idx)) then + stop 'cannot consolidate map : consolidated_idx not associated' + endif + + call map_sort(map) + k = 1_8 + do i=0_8, map % map_size + map % consolidated_idx (i+1) = k + do j=1, map % map(i) % n_elements + map % consolidated_value(k) = map % map(i) % value(j) + map % consolidated_key (k) = map % map(i) % key(j) + k = k+1_8 + enddo + deallocate(map % map(i) % value) + deallocate(map % map(i) % key) + map % map(i) % value => map % consolidated_value ( map % consolidated_idx (i+1) :) + map % map(i) % key => map % consolidated_key ( map % consolidated_idx (i+1) :) + enddo + map % consolidated_idx (map % map_size + 2_8) = k + map % consolidated = .True. + + + call munmap( (/ map % map_size + 2_8 /), 8, fd(1), c_pointer(1)) + call mmap(trim(filename)//'_consolidated_idx', (/ map % map_size + 2_8 /), 8, fd(1), .True., c_pointer(1)) + call c_f_pointer(c_pointer(1),map % consolidated_idx, (/ map % map_size +2_8/)) + + call munmap( (/ map % n_elements /), cache_key_kind, fd(2), c_pointer(2)) + call mmap(trim(filename)//'_consolidated_key', (/ map % n_elements /), cache_key_kind, fd(2), .True., c_pointer(2)) + call c_f_pointer(c_pointer(2),map % consolidated_key, (/ map % n_elements /)) + + call munmap( (/ map % n_elements /), integral_kind, fd(3), c_pointer(3)) + call mmap(trim(filename)//'_consolidated_value', (/ map % n_elements /), integral_kind, fd(3), .True., c_pointer(3)) + call c_f_pointer(c_pointer(3),map % consolidated_value, (/ map % n_elements /)) + +end + +subroutine map_load_from_disk(filename,map) + use map_module + use mmap_module + implicit none + character*(*), intent(in) :: filename + type(map_type), intent(inout) :: map + type(c_ptr) :: c_pointer(3) + integer :: fd(3) + integer*8 :: i,k + integer :: n_elements + + + + if (map % consolidated) then + stop 'map already consolidated' + endif + + call mmap(trim(filename)//'_consolidated_idx', (/ map % map_size + 2_8 /), 8, fd(1), .True., c_pointer(1)) + call c_f_pointer(c_pointer(1),map % consolidated_idx, (/ map % map_size + 2_8/)) + + map% n_elements = map % consolidated_idx (map % map_size+2_8)-1 + + call mmap(trim(filename)//'_consolidated_key', (/ map % n_elements /), cache_key_kind, fd(2), .True., c_pointer(2)) + call c_f_pointer(c_pointer(2),map % consolidated_key, (/ map % n_elements /)) + + call mmap(trim(filename)//'_consolidated_value', (/ map % n_elements /), integral_kind, fd(3), .True., c_pointer(3)) + call c_f_pointer(c_pointer(3),map % consolidated_value, (/ map % n_elements /)) + + k = 1_8 + do i=0_8, map % map_size + deallocate(map % map(i) % value) + deallocate(map % map(i) % key) + map % map(i) % value => map % consolidated_value ( map % consolidated_idx (i+1) :) + map % map(i) % key => map % consolidated_key ( map % consolidated_idx (i+1) :) + map % map(i) % sorted = .True. + n_elements = map % consolidated_idx (i+2) - k + k = map % consolidated_idx (i+2) + map % map(i) % map_size = n_elements + map % map(i) % n_elements = n_elements + enddo + map % n_elements = k-1 + map % sorted = .True. + map % consolidated = .True. + +end +