Line data Source code
1 : MODULE m_iomatrix_hdf
2 : USE m_judft
3 : USE hdf5
4 : USE m_hdf_tools
5 : USE m_types_mat
6 : USE m_types_mpimat
7 : IMPLICIT NONE
8 : PRIVATE
9 : PUBLIC iomatrix_hdf_close,iomatrix_hdf_open,iomatrix_hdf_write,iomatrix_hdf_read
10 :
11 : CONTAINS
12 0 : SUBROUTINE iomatrix_hdf_read(mat,nrec,did)
13 : CLASS(t_Mat),INTENT(INOUT) :: mat
14 : INTEGER,INTENT(IN) :: nrec
15 : INTEGER(HID_t),INTENT(in) :: did
16 :
17 : INTEGER::mpi_comm,dim(4)
18 :
19 :
20 : INTEGER(HID_t) :: memspace,fspace,trans
21 : INTEGER(hsize_t):: dims(4)
22 : INTEGER :: err
23 0 : REAL,ALLOCATABLE :: dat(:,:,:,:)
24 : SELECT TYPE(mat)
25 : TYPE is (t_mpimat)
26 0 : mpi_comm=mat%blacsdata%mpi_com !Only information used from mat intent(in)
27 : CLASS default
28 0 : mpi_comm=0
29 : END SELECT
30 0 : CALL io_datadim(did,dim)
31 0 : CALL mat%init(DIM(1)==1,DIM(2),DIM(3),mpi_comm,.TRUE.)
32 : SELECT TYPE(mat)
33 : TYPE is (t_mat)
34 0 : IF (mat%l_real) THEN
35 0 : CALL io_read(did,(/1,1,1,nrec/),(/1,mat%matsize1,mat%matsize2,1/),"data_r",mat%data_r)
36 : ELSE
37 0 : CALL io_read(did,(/-1,1,1,nrec/),(/1,mat%matsize1,mat%matsize2,1/),"data_c",mat%data_c)
38 : END IF
39 : TYPE is (t_mpimat)
40 0 : ALLOCATE(dat(MERGE(1,2,mat%l_real),mat%matsize1,mat%matsize2,1))
41 :
42 0 : CALL h5dget_space_f(did,fspace,err)
43 0 : CALL priv_create_hyperslab_from_blacsdesc(mat%l_real,nrec,fspace,mat%blacsdata%blacs_desc)
44 0 : dims=SHAPE(dat)
45 0 : CALL h5screate_simple_f(4,dims,memspace,err)
46 0 : trans=gettransprop()
47 0 : CALL h5dread_f(did,H5T_NATIVE_DOUBLE,dat,dims,err,memspace,fspace,trans)
48 0 : CALL h5sclose_f(memspace,err)
49 0 : CALL h5sclose_f(fspace,err)
50 0 : CALL cleartransprop(trans)
51 0 : IF (mat%l_real) THEN
52 0 : mat%data_r=dat(1,:,:,1)
53 : ELSE
54 0 : mat%data_c=CMPLX(dat(1,:,:,1),dat(2,:,:,1))
55 : ENDIF
56 : END SELECT
57 0 : END SUBROUTINE iomatrix_hdf_read
58 :
59 0 : SUBROUTINE iomatrix_hdf_write(mat,rec,did)
60 : CLASS(t_Mat),INTENT(IN) :: mat
61 : INTEGER,INTENT(IN) :: rec
62 : INTEGER(HID_t),INTENT(in)::did
63 :
64 : INTEGER(HID_t) :: memspace,fspace,trans
65 : INTEGER(HSIZE_t):: dims(4)
66 : INTEGER :: err
67 0 : REAL,ALLOCATABLE :: dat(:,:,:,:)
68 :
69 :
70 : SELECT TYPE(mat)
71 : TYPE is (t_mat)
72 0 : IF (mat%l_real) THEN
73 0 : CALL io_write(did,(/1,1,1,rec/),(/1,mat%matsize1,mat%matsize2,1/),"data_r",mat%data_r)
74 : ELSE
75 0 : CALL io_write(did,(/1,1,1,rec/),(/1,mat%matsize1,mat%matsize2,1/),"data_c_re",REAL(mat%data_c))
76 0 : CALL io_write(did,(/2,1,1,rec/),(/1,mat%matsize1,mat%matsize2,1/),"data_c_im",AIMAG(mat%data_c))
77 : END IF
78 : TYPE is (t_mpimat)
79 0 : ALLOCATE(dat(MERGE(1,2,mat%l_real),mat%matsize1,mat%matsize2,1))
80 0 : IF (mat%l_real) THEN
81 0 : dat(1,:,:,1)=mat%data_r
82 : ELSE
83 0 : dat(1,:,:,1)=REAL(mat%data_c)
84 0 : dat(2,:,:,1)=AIMAG(mat%data_c)
85 : ENDIF
86 0 : CALL h5dget_space_f(did,fspace,err)
87 0 : CALL priv_create_hyperslab_from_blacsdesc(mat%l_real,rec,fspace,mat%blacsdata%blacs_desc)
88 0 : dims=SHAPE(dat)
89 0 : CALL h5screate_simple_f(4,dims,memspace,err)
90 0 : trans=gettransprop()
91 0 : CALL h5dwrite_f(did,H5T_NATIVE_DOUBLE,dat,dims,err,memspace,fspace,trans)
92 0 : CALL h5sclose_f(memspace,err)
93 0 : CALL h5sclose_f(fspace,err)
94 0 : CALL cleartransprop(trans)
95 : END SELECT
96 :
97 0 : END SUBROUTINE iomatrix_hdf_write
98 :
99 0 : SUBROUTINE iomatrix_hdf_close(fid,did)
100 : INTEGER(hid_t),INTENT(inout):: fid,did
101 : INTEGER:: err
102 0 : CALL h5dclose_f(did,err)
103 0 : CALL h5fclose_f(fid,err)
104 0 : END SUBROUTINE iomatrix_hdf_close
105 :
106 0 : SUBROUTINE iomatrix_hdf_open(l_real,matsize,no_rec,filename,fid,did)
107 : #ifdef CPP_HDFMPI
108 : USE mpi
109 : #endif
110 : LOGICAL,INTENT(IN) :: l_real
111 : INTEGER,INTENT(in) :: matsize,no_rec
112 : CHARACTER(len=*),INTENT(in) :: filename
113 : INTEGER(hid_t),INTENT(out) :: fid,did
114 :
115 : INTEGER :: dims(4),err
116 : LOGICAL :: l_exist
117 : INTEGER(HID_T) :: access_prp
118 : #if defined(CPP_HDFMPI) && defined(CPP_MPI)
119 0 : CALL h5pcreate_f(H5P_FILE_ACCESS_F, access_prp, err)
120 0 : CALL h5pset_fapl_mpio_f(access_prp, MPI_COMM_WORLD, MPI_INFO_NULL,err)
121 : #else
122 : access_prp=H5P_DEFAULT_f
123 : #endif
124 :
125 :
126 0 : INQUIRE(file=filename//'.hdf',exist=l_exist)
127 0 : IF (l_exist) THEN
128 0 : CALL h5fopen_f(filename//'.hdf',H5F_ACC_RDWR_F,fid,err,access_prp)
129 : ELSE
130 0 : CALL h5fcreate_f(filename//'.hdf',H5F_ACC_TRUNC_F,fid,err,H5P_DEFAULT_f,access_prp)
131 : ENDIF
132 0 : IF (io_dataexists(fid,'Matrix')) THEN
133 0 : CALL h5dopen_f(fid,"Matrix", did, err)
134 : ELSE
135 : !Create data-space
136 0 : dims(1) = MERGE(1,2,l_real)
137 0 : dims(2:3)= matsize
138 0 : dims(4) = no_rec
139 0 : call io_createvar(fid,"Matrix",H5T_NATIVE_DOUBLE,dims,did)
140 : END IF
141 0 : END SUBROUTINE iomatrix_hdf_open
142 :
143 :
144 0 : SUBROUTINE priv_create_hyperslab_from_blacsdesc(l_real,nrec,sid,blacsdesc)
145 : LOGICAL,INTENT(IN) :: l_real
146 : INTEGER,INTENT(in) :: nrec,blacsdesc(9)
147 : INTEGER(hid_t),INTENT(in):: sid
148 :
149 : INTEGER(hsize_t):: start(4),COUNT(4),stride(4),bloc(4)
150 : INTEGER :: nprow,npcol,myrow,mycol,block_row,block_col,matsize,blacs_ctxt,err
151 : LOGICAL :: ok
152 : !For readability get data from blacsdesc
153 : #ifdef CPP_SCALAPACK
154 0 : blacs_ctxt=blacsdesc(2)
155 0 : block_row=blacsdesc(5)
156 0 : block_col=blacsdesc(6)
157 0 : matsize=blacsdesc(4)
158 0 : CALL blacs_gridinfo(blacs_ctxt,nprow,npcol,myrow,mycol)
159 :
160 0 : CALL h5Sselect_none_f(sid,err) !unselect all elements in dataspace
161 : !Select blocks of blacs-grid
162 0 : start=(/0,myrow*block_row,mycol*block_col,nrec-1/)
163 0 : count=(/1,FLOOR(1.*matsize/block_row)+1,FLOOR(1.*matsize/block_col)+1,1/)
164 0 : stride=(/1,nprow*block_row,npcol*block_col,1/)
165 0 : bloc=(/MERGE(1,2,l_real),block_row,block_col,1/)
166 0 : CALL h5sselect_hyperslab_f(sid,H5S_SELECT_OR_F,start,count,err,stride,bloc)
167 0 : CALL h5sselect_valid_f(sid,ok,err)
168 0 : IF (.NOT.ok) THEN
169 0 : CALL h5sget_simple_extent_dims_f(sid,start,stride,err)
170 : !Cut to actual sizes
171 0 : start=(/0,0,0,0/)
172 0 : count=(/MERGE(1,2,l_real),matsize,matsize,int(stride(4))/)
173 0 : CALL h5sselect_hyperslab_f(sid,H5S_SELECT_AND_F,start,count,err)
174 0 : CALL h5sselect_valid_f(sid,ok,err)
175 0 : IF (.NOT.ok) CALL judft_error("Writing of matrix failed, BUG in parallel HDF-IO")
176 : ENDIF
177 : #endif
178 0 : END SUBROUTINE priv_create_hyperslab_from_blacsdesc
179 0 : END MODULE m_iomatrix_hdf
|