1!dalton_copyright_start 2! 3! 4!dalton_copyright_end 5 6#ifdef VAR_MPI 7module file_io_model 8 9! stefan: - this module provides all necessary functionality 10! to setup a file i/o model in parallel mcscf/ci 11! calculations. 12! 13! written by sknecht, may 2007 for DIRAC MCSCF/KR-CI/LUCITA 14! adapted for DALTON by sknecht, november 2010. 15#ifdef USE_MPI_MOD_F90 16 use mpi 17 implicit none 18#else 19 implicit none 20#include "mpif.h" 21#endif 22 23 public setup_file_io_model 24 public close_file_io_model 25 public set_file_io_offset 26 27 private 28 29 save 30 31 integer(kind=MPI_INTEGER_KIND) :: my_MPI_REAL8 = MPI_REAL8 32 integer(kind=MPI_INTEGER_KIND) :: istat(MPI_STATUS_SIZE) 33 integer(kind=MPI_INTEGER_KIND) :: ierr_mpi 34 integer, parameter :: flabel_length = 14 35 36contains 37 38 subroutine setup_file_io_model(communicator_io_group, & 39 nr_files, & 40 fh_array, & 41 fh_offset, & 42 group_id, & 43 group_io_size, & 44 file_identification, & 45 print_unit) 46!****************************************************************************** 47! 48! purpose: return valid fh for parallel MPI-I/O in CI/MCSCF runs: 49! - provides: 50! a. file handles 51! b. opened files ready for reading and writing 52! - requires: 53! a. global (group) communication handle 54! a. file identification string 55! 56! for a detailed description of the I/O model see references: 57! S. Knecht, H. J. Aa. Jensen, and T. Fleig 58! JCP, 128, 014108 (2008) 59! JCP, 132, 014108 (2010) 60! 61!******************************************************************************* 62 integer, intent(in) :: communicator_io_group 63 integer, intent(in) :: nr_files 64 integer, intent(inout) :: fh_array(nr_files) 65 integer, intent(in) :: fh_offset 66 integer, intent(in) :: group_id 67 integer, intent(in) :: group_io_size 68 integer, intent(in) :: print_unit 69 character (len=5), intent(in) :: file_identification 70!------------------------------------------------------------------------------- 71 integer :: i 72 integer :: i_relative 73 integer :: j 74 integer(kind=MPI_INTEGER_KIND) :: file_info_obj 75 integer(kind=MPI_INTEGER_KIND) :: comm_iogrp_mpi 76 integer(kind=MPI_INTEGER_KIND) :: fh_array_mpi 77 integer(kind=MPI_OFFSET_KIND) :: displacement 78 character (len= 4) :: file_info_groupsz 79 character (len= 4) :: fstring 80 character (len= 4) :: gstring 81 character (len= flabel_length) :: flabel 82!------------------------------------------------------------------------------- 83 84! initial displacement in newly created file 85 displacement = 0 86 87! group id appended to each file 88 call int2char_converter(group_id,gstring) 89 90! file info object - provide hints to the MPI implementation 91 call mpi_info_create(file_info_obj,ierr_mpi) 92 93! 1. number of processes sharing the following MPI-I/O files 94 write(file_info_groupsz,'(i4)') group_io_size 95 call mpi_info_set(file_info_obj,"nb_proc",file_info_groupsz,ierr_mpi) 96! 97#ifdef VAR_PFS 98! 2. extra information on IBMs GPFS to enhance I/O performance 99 call mpi_info_set(file_info_obj,"IBM_largeblock_io","true",ierr_mpi) 100#endif 101 102 do i = 1, nr_files 103 104 i_relative = i + fh_offset 105 106! step a. setup individual file identifier 107 call int2char_converter(i_relative,fstring) 108 109! step b. determine full file name 110 write(flabel,'(a5,a4,a1,a4)') file_identification,fstring,'.',gstring 111 112! step c. open the file 113 comm_iogrp_mpi = communicator_io_group 114 call mpi_file_open(comm_iogrp_mpi,flabel(1:flabel_length), & 115 MPI_MODE_CREATE + MPI_MODE_RDWR + MPI_MODE_DELETE_ON_CLOSE, & 116 file_info_obj,fh_array_mpi,ierr_mpi) 117 fh_array(i) = fh_array_mpi 118! step d. set fileview 119 call mpi_file_set_view(fh_array_mpi,displacement,my_MPI_REAL8,my_MPI_REAL8, & 120 "native",file_info_obj,ierr_mpi) 121 end do 122 123! free info object 124 call mpi_info_free(file_info_obj,ierr_mpi) 125 126 end subroutine setup_file_io_model 127!******************************************************************************* 128 129 subroutine close_file_io_model(number_of_files, & 130 fh_offset, & 131 fh_array) 132!******************************************************************************* 133! 134! purpose: close MPI-I/O files and "nullify" file handles. 135! 136!******************************************************************************* 137 integer, intent(in ) :: number_of_files 138 integer, intent(in ) :: fh_offset 139 integer, intent(inout) :: fh_array(*)!(nr_files+fh_offset) 140!------------------------------------------------------------------------------- 141 integer :: i 142 integer(kind=MPI_INTEGER_KIND) :: fh_array_mpi 143!------------------------------------------------------------------------------- 144 145 do i = 1, number_of_files 146 fh_array_mpi = fh_array(i+fh_offset) 147 call mpi_file_close(fh_array_mpi,ierr_mpi) 148 end do 149 150 end subroutine close_file_io_model 151!******************************************************************************* 152 153 subroutine set_file_io_offset(nr_files, & 154 file_offset_off, & 155 file_offset_array, & 156 file_offset_fac, & 157 gvec_offset_real, & 158 gvec_offset_cplx, & 159 gvec_ablock_real, & 160 gvec_ablock_cplx, & 161 gvec_tblock_gnrl, & 162 complex_algebra, & 163 my_process_id, & 164 nr_gvec_ttss_blocks, & 165 intra_node_size, & 166 group_list, & 167 par_dist_block_list, & 168 block_list) 169!******************************************************************************* 170! 171! purpose: set absolute offset for each process in MPI-I/O files and 172! provide general offset parameters. 173! 174!******************************************************************************* 175 integer, intent(in ) :: nr_files 176 integer, intent(in ) :: file_offset_off 177 integer, intent(in ) :: my_process_id 178 integer, intent(in ) :: nr_gvec_ttss_blocks 179 integer, intent(in ) :: intra_node_size 180 integer, intent(in ) :: group_list(*) 181 integer, intent(in ) :: par_dist_block_list(nr_gvec_ttss_blocks) 182 integer, intent(in ) :: block_list(nr_gvec_ttss_blocks) 183 integer, intent(out) :: gvec_ablock_real 184 integer, intent(out) :: gvec_ablock_cplx 185 integer, intent(out) :: gvec_tblock_gnrl 186 integer(kind=MPI_OFFSET_KIND), intent(out) :: gvec_offset_real 187 integer(kind=MPI_OFFSET_KIND), intent(out) :: gvec_offset_cplx 188 integer(kind=MPI_OFFSET_KIND), intent(out) :: file_offset_array(*)!(nr_files+file_offset_off) 189! integer(kind=MPI_OFFSET_KIND), intent(in) :: file_offset_fac(*)!(nr_files+file_offset_off) 190 integer , intent(in) :: file_offset_fac(*)!(nr_files+file_offset_off) 191 logical, intent(in ) :: complex_algebra 192!------------------------------------------------------------------------------- 193 integer :: i 194 integer :: mult_fac_gvec_blocks 195 integer(kind=MPI_OFFSET_KIND) :: mult_fac_gvec_offset 196 integer(kind=MPI_OFFSET_KIND) :: save_vec_offset 197!------------------------------------------------------------------------------- 198 199 mult_fac_gvec_blocks = 1 200 mult_fac_gvec_offset = 1 201 202 if(complex_algebra)then 203 mult_fac_gvec_blocks = 2 204 mult_fac_gvec_offset = 2 205 end if 206 207 call calc_general_offset_param(gvec_offset_real, & 208 gvec_offset_cplx, & 209 gvec_ablock_real, & 210 gvec_ablock_cplx, & 211 gvec_tblock_gnrl, & 212 save_vec_offset, & 213 mult_fac_gvec_blocks, & 214 mult_fac_gvec_offset, & 215 my_process_id, & 216 nr_gvec_ttss_blocks, & 217 intra_node_size, & 218 group_list, & 219 par_dist_block_list, & 220 block_list) 221 do i = 1, nr_files 222 file_offset_array(i+file_offset_off) = & 223 save_vec_offset * mult_fac_gvec_offset * file_offset_fac(i+file_offset_off) 224 end do 225 226 end subroutine set_file_io_offset 227!******************************************************************************* 228 229 subroutine calc_general_offset_param(gvec_offset_real, & 230 gvec_offset_cplx, & 231 gvec_ablock_real, & 232 gvec_ablock_cplx, & 233 gvec_tblock_gnrl, & 234 save_vec_offset, & 235 mult_fac_gvec_blocks, & 236 mult_fac_gvec_offset, & 237 my_process_id, & 238 nr_gvec_ttss_blocks, & 239 intra_node_size, & 240 group_list, & 241 par_dist_block_list, & 242 block_list) 243!******************************************************************************* 244! 245! purpose: provide general offset parameters for MPI-I/O files. 246! 247!******************************************************************************* 248 integer, intent(in ) :: my_process_id 249 integer, intent(in ) :: nr_gvec_ttss_blocks 250 integer, intent(in ) :: intra_node_size 251 integer, intent(in ) :: group_list(*) 252 integer, intent(in ) :: par_dist_block_list(nr_gvec_ttss_blocks) 253 integer, intent(in ) :: block_list(nr_gvec_ttss_blocks) 254 integer, intent(in ) :: mult_fac_gvec_blocks 255 integer, intent(out) :: gvec_ablock_real 256 integer, intent(out) :: gvec_ablock_cplx 257 integer, intent(out) :: gvec_tblock_gnrl 258 integer(kind=MPI_OFFSET_KIND), intent(out) :: gvec_offset_real 259 integer(kind=MPI_OFFSET_KIND), intent(out) :: gvec_offset_cplx 260 integer(kind=MPI_OFFSET_KIND), intent(out) :: save_vec_offset 261 integer(kind=MPI_OFFSET_KIND), intent(in ) :: mult_fac_gvec_offset 262!------------------------------------------------------------------------------- 263 integer :: current_proc 264 integer :: process_id 265 integer :: current_block 266 integer :: tmp_active_blocks 267 integer(kind=MPI_OFFSET_KIND) :: tmp_gvec_offset 268 integer(kind=MPI_OFFSET_KIND) :: tmp_gvec_offset_2 269!------------------------------------------------------------------------------- 270 271 gvec_ablock_real = 0 272 gvec_ablock_cplx = 0 273 gvec_tblock_gnrl = 0 274 gvec_offset_real = 0 275 gvec_offset_cplx = 0 276 277! total number of ttss-blocks 278 gvec_tblock_gnrl = nr_gvec_ttss_blocks 279 280! calculate general offset parameters for each process in a given group: 281! a. total number of active ttss-blocks (real and complex) 282! b. vector offset (real and complex) 283 284 tmp_active_blocks = 0 285 tmp_gvec_offset = 0 286 tmp_gvec_offset_2 = 0 287 288 do current_proc = 1, intra_node_size 289 290 process_id = group_list(current_proc) 291 292 do current_block = 1, nr_gvec_ttss_blocks 293 294 if(par_dist_block_list(current_block) == process_id)then 295 tmp_gvec_offset = tmp_gvec_offset + block_list(current_block) 296 if(my_process_id == process_id)then 297 gvec_offset_real = gvec_offset_real + block_list(current_block) 298 gvec_offset_cplx = gvec_offset_cplx + block_list(current_block) 299 tmp_active_blocks = tmp_active_blocks + 1 300 end if 301 end if 302 303 end do !nr_gvec_ttss_blocks 304 305 if(my_process_id == process_id)then 306 307 gvec_ablock_real = tmp_active_blocks 308 gvec_ablock_cplx = mult_fac_gvec_blocks * tmp_active_blocks 309 gvec_offset_cplx = mult_fac_gvec_offset * gvec_offset_cplx 310 311 save_vec_offset = tmp_gvec_offset_2 312 313 end if 314 315! keep track of vector and block offsets 316 tmp_active_blocks = 0 317 tmp_gvec_offset_2 = tmp_gvec_offset 318 end do ! intra_node_size 319 320 end subroutine calc_general_offset_param 321!******************************************************************************* 322 323 subroutine int2char_converter(int_number,string_rep) 324!******************************************************************************* 325! 326! purpose: convert positive integer number int_number (< 10 000) 327! into the 4-byte string string_rep. 328! 329! based on the routine num2str originally written 330! by C.V. Larsen in Dirac. 331! 332!******************************************************************************* 333 integer, intent(in ) :: int_number 334 character (len = 4), intent(inout) :: string_rep 335!------------------------------------------------------------------------------- 336 character (len = 1) :: tmp_str(4) 337 integer :: num1 338 integer :: num2 339 integer :: num3 340 integer :: num4 341!------------------------------------------------------------------------------- 342 343 num3 = int_number 344 num2 = 1 345 num1 = 1000 346 347 do 348 num4 = num3/num1 349 tmp_str(num2) = char(num4 + 48) 350 num3 = mod(num3,num1) 351 num2 = num2 + 1 352 num1 = num1/10 353 354 if(num1 < 1) exit 355 end do 356 357 string_rep=tmp_str(1)//tmp_str(2)//tmp_str(3)//tmp_str(4) 358 359 end subroutine int2char_converter 360 361end module 362#else 363subroutine fhio_model 364! dummy routine for non-mpi compilation 365end 366#endif 367