diff --git a/.gitignore b/.gitignore index e69de29bb2d1d6434b8b29ae775ad8c2e48c5391..b0541cc72e077713a561e3e1bf1be74ac19dc450 100644 --- a/.gitignore +++ b/.gitignore @@ -0,0 +1,2 @@ +*~ +#~ diff --git a/cdn/cdnval.F90 b/cdn/cdnval.F90 index a9cf083868adad4483b3b5b0e9825ed9cd3a586b..ccca2c241a1e2e09bd3a05526edebb4d7bc4b7be 100644 --- a/cdn/cdnval.F90 +++ b/cdn/cdnval.F90 @@ -44,6 +44,7 @@ CONTAINS ! and bands !*********************************************************************** ! + USE m_eig66_io,ONLY: write_dos USE m_radfun USE m_radflo USE m_rhomt @@ -158,8 +159,8 @@ CONTAINS INTEGER, ALLOCATABLE :: nmtsl(:,:),nslat(:,:) REAL, ALLOCATABLE :: zsl(:,:),volsl(:) REAL, ALLOCATABLE :: volintsl(:) - REAL, ALLOCATABLE :: qintsl(:,:,:,:),qmtsl(:,:,:,:) - REAL, ALLOCATABLE :: orbcomp(:,:,:,:,:),qmtp(:,:,:,:) + REAL, ALLOCATABLE :: qintsl(:,:),qmtsl(:,:) + REAL, ALLOCATABLE :: orbcomp(:,:,:),qmtp(:,:) !-new_sl !-dw INTEGER, ALLOCATABLE :: gvac1d(:),gvac2d(:) ,kveclo(:) @@ -429,10 +430,10 @@ CONTAINS atoms,cell,nsld,& nsl,zsl,nmtsl,nslat,volsl,volintsl) ! - ALLOCATE ( qintsl(nsld,dimension%neigd,kpts%nkptd,dimension%jspd)) - ALLOCATE ( qmtsl(nsld,dimension%neigd,kpts%nkptd,jsp_start:jsp_end)) - ALLOCATE ( orbcomp(dimension%neigd,23,atoms%natd,kpts%nkptd,jsp_start:jsp_end) ) - ALLOCATE ( qmtp(dimension%neigd,atoms%natd,kpts%nkptd,jsp_start:jsp_end) ) + ALLOCATE ( qintsl(nsld,dimension%neigd)) + ALLOCATE ( qmtsl(nsld,dimension%neigd)) + ALLOCATE ( orbcomp(dimension%neigd,23,atoms%natd) ) + ALLOCATE ( qmtp(dimension%neigd,atoms%natd) ) IF (.not.input%film) qvac(:,:,:,jspin) = 0.0 ENDIF !-q_sl @@ -626,7 +627,7 @@ CONTAINS cell,& z,noccbd,lapw,& nsl,zsl,nmtsl,oneD,& - qintsl(:,:,ikpt,jspin)) + qintsl(:,:)) ! ENDIF ENDIF @@ -727,7 +728,7 @@ CONTAINS IF (banddos%dos.AND.(banddos%ndir.EQ.-3)) THEN CALL q_mt_sl(ispin, atoms,noccbd,nsld, ikpt,noccbd,ccof(-atoms%llod,1,1,1,ispin),& skip_t,noccbd, acof(:,0:,:,ispin),bcof(:,0:,:,ispin),usdus,& - nmtsl,nsl, qmtsl(:,:,ikpt,ispin)) + nmtsl,nsl, qmtsl(:,:)) ! INQUIRE (file='orbcomprot',exist=l_orbcomprot) @@ -738,7 +739,7 @@ CONTAINS ENDIF CALL orb_comp(ispin,noccbd,atoms,noccbd,usdus,acof(1:,0:,1:,ispin),bcof(1:,0:,1:,ispin),& - ccof(-atoms%llod:,1:,1:,1:,ispin), orbcomp(:,:,:,ikpt,ispin), qmtp(:,:,ikpt,ispin)) + ccof(-atoms%llod:,1:,1:,1:,ispin), orbcomp, qmtp) ! ENDIF !-new @@ -829,56 +830,10 @@ CONTAINS ! !--dw now write k-point data to tmp_dos ! - CALL judft_error("IO in cdnval not implemented") - IF (mpi%irank==0) THEN - IF ( .not.l_mcd ) THEN - DO ispin = jsp_start,jsp_end - WRITE (84,rec=kpts%nkpt*(ispin-1)+ikpt) bkpt,wk,nbands,eig,& - & qal(:,:,:,ispin),qvac(:,:,ikpt,ispin),& - & qis(:,ikpt,ispin),qvlay(:,:,:,ikpt,ispin),& - & qstars,ksym,jsym - ENDDO - ELSE - WRITE (84,rec=kpts%nkpt*(jspin-1)+ikpt) bkpt,wk,nbands,eig,& - & qal(:,:,:,jspin),qvac(:,:,ikpt,jspin),& - & qis(:,ikpt,jspin),qvlay(:,:,:,ikpt,jspin),& - & qstars,ksym,jsym,mcd(:,1:ncored,:) - ENDIF ! not.l_mcd - !+new_sl - IF (banddos%ndir.EQ.-3) THEN - IF ( .not.l_mcd ) THEN - DO ispin = jsp_start,jsp_end - WRITE (129,rec=kpts%nkpt*(ispin-1)+ikpt) bkpt,wk,& - nbands,eig,qvac(:,:,ikpt,ispin),& - qintsl(:,:,ikpt,ispin),& - qmtsl(:,:,ikpt,ispin),& - orbcomp(:,:,:,ikpt,ispin),& - qmtp(:,:,ikpt,ispin) - ENDDO - ELSE - WRITE (6,FMT=500) - WRITE (16,FMT=500) -500 FORMAT('A calculation of an orbital composition& - &of input%film states (s,px,py,pz,....)& - & for the l_mcd case is not implemented.') - ENDIF ! not.l_mcd - ENDIF ! banddos%ndir == -3 - ENDIF ! mpi%irank==0 -#ifdef CPP_MPI - DO ispin = jsp_start,jsp_end - CALL mpi_col_dos(& - mpi,l_evp,noccbd,noccbd_l,& - dimension,atoms,atoms,& - nsld,ncored,& - vacuum,vacuum,l_mcd,banddos,kpts,ispin,ikpt,& - bkpt,wk,nbands,eig,qal(0,1,1,ispin),& - qvac(1,1,ikpt,ispin),qis(1,nkpt,ispin),& - qvlay(:,:,:,ikpt,ispin),qstars,ksym,& - jsym,m_mcd,qintsl(:,:,ikpt,ispin),& - qmtsl(:,:,ikpt,ispin),qmtp(:,:,ikpt,ispin),& - orbcomp(:,:,:,ikpt,ispin)) - ENDDO -#endif + call write_dos(eig_id,ikpt,ispin,qal(:,:,:,ispin),qvac(:,:,ikpt,ispin),qis(:,ikpt,ispin),& + qvlay(:,:,:,ikpt,ispin),qstars,ksym,jsym,mcd,qintsl,& + qmtsl(:,:),qmtp(:,:),orbcomp) + CALL timestop("cdnval: write_info") !-new_sl ENDIF @@ -917,7 +872,7 @@ enddo call timestart("cdnval: dos") IF (mpi%irank==0) THEN CALL doswrite(& - dimension,kpts,atoms,vacuum,& + eig_id,dimension,kpts,atoms,vacuum,& input,banddos,& sliceplot,noco,sym,& cell,& @@ -925,7 +880,7 @@ enddo results%ef,nsld,oneD) IF (banddos%dos.AND.(banddos%ndir.EQ.-3)) THEN CALL Ek_write_sl(& - dimension,kpts,atoms,vacuum,& + eig_id,dimension,kpts,atoms,vacuum,& nsld,input,jspin,& sym,cell,& nsl,nslat) diff --git a/dos/Ek_write_sl.f90 b/dos/Ek_write_sl.f90 index 2904d286f9c650c292e8ca262839f74a15342e8f..00d53f47502b78d088be838c5e038da1dde4d79f 100644 --- a/dos/Ek_write_sl.f90 +++ b/dos/Ek_write_sl.f90 @@ -1,16 +1,14 @@ MODULE m_Ekwritesl use m_juDFT CONTAINS - SUBROUTINE Ek_write_sl(& - & dimension,kpts,atoms,vacuum,nsld,& - & input,jspin,& - & sym,cell,& - & nsl,nslat) + SUBROUTINE Ek_write_sl(eig_id,dimension,kpts,atoms,vacuum,nsld,& + input,jspin, sym,cell, nsl,nslat) !----------------------------------------------------------------- !-- now write E(k) for all kpts if on T3E !-- now read data from tmp_dos and write of E(k) in ek_orbcomp !----------------------------------------------------------------- USE m_types + USE m_eig66_io IMPLICIT NONE TYPE(t_dimension),INTENT(IN) :: dimension TYPE(t_input),INTENT(IN) :: input @@ -21,7 +19,7 @@ CONTAINS TYPE(t_atoms),INTENT(IN) :: atoms ! .. ! .. Scalar Arguments .. - INTEGER, INTENT (IN) :: nsld + INTEGER, INTENT (IN) :: nsld,eig_id INTEGER, INTENT (IN) :: nsl ,jspin ! .. ! .. Array Arguments .. @@ -36,19 +34,22 @@ CONTAINS INTEGER norb(23),iqsl(nsld),iqvacpc(2) REAL bkpt(3),qvact(2) REAL, ALLOCATABLE :: eig(:),qvac(:,:,:,:),orbcomp(:,:,:,:,:) + REAL, ALLOCATABLE :: qal(:,:,:),qis(:),qvlay(:,:,:) + COMPLEX,ALLOCATABLE::qstars(:,:,:,:) + INTEGER,ALLOCATABLE::ksym(:),jsym(:) REAL, ALLOCATABLE :: qintsl(:,:,:,:),qmtsl(:,:,:,:),qmtp(:,:,:,:) CHARACTER (len=2) :: chntype CHARACTER (len=99) :: chform ! .. - ! .. Intrinsic Functions - INTRINSIC nint - ! IF (nsl.GT.nsld) THEN CALL juDFT_error("nsl.GT.nsld",calledby="Ek_write_sl") ENDIF ALLOCATE(eig(dimension%neigd),orbcomp(dimension%neigd,23,atoms%natd,kpts%nkptd,dimension%jspd)) ALLOCATE(qvac(dimension%neigd,2,kpts%nkptd,dimension%jspd),qintsl(nsld,dimension%neigd,kpts%nkptd,dimension%jspd)) ALLOCATE(qmtsl(nsld,dimension%neigd,kpts%nkptd,dimension%jspd),qmtp(dimension%neigd,atoms%natd,kpts%nkptd,dimension%jspd)) + ALLOCATE(qal(4,atoms%ntype,dimension%neigd),qis(dimension%neigd),qvlay(dimension%neigd,vacuum%layerd,2)) + ALLOCATE(qstars(vacuum%nstars,dimension%neigd,vacuum%layerd,2)) + ALLOCATE(ksym(dimension%neigd),jsym(dimension%neigd)) ! ! ---> open files for a bandstucture with an orbital composition ! ---> in the case of the film geometry @@ -83,10 +84,9 @@ CONTAINS !============================================================== DO ikpt=1,kpts%nkpt ! - READ (129,rec=kpts%nkpt*(kspin-1)+ikpt) bkpt,wk,nbands,eig,& - & qvac(:,:,ikpt,kspin),qintsl(:,:,ikpt,kspin),& - & qmtsl(:,:,ikpt,kspin),& - & orbcomp(:,:,:,ikpt,kspin),qmtp(:,:,ikpt,kspin) + call read_eig(eig_id,ikpt,kspin,bk=bkpt,neig=nbands,eig=eig) + call read_dos(eig_id,ikpt,kspin,qal,qvac(:,:,ikpt,kspin),qis,qvlay,qstars,ksym,jsym,& + qintsl=qintsl(:,:,ikpt,kspin),qmtsl= qmtsl(:,:,ikpt,kspin),qmtp=qmtp(:,:,ikpt,kspin),orbcomp=orbcomp(:,:,:,ikpt,kspin)) ! write(*,*) kspin,nkpt,qmtp(1,:,ikpt,kspin) ! WRITE (130,FMT=8000) (bkpt(i),i=1,3) diff --git a/dos/doswrite.f90 b/dos/doswrite.f90 index 37e5c4fc37d256df5bb8f91e49757d489eb5bc5c..9a367199d26621d7dbea379d89ca51099b67a30b 100644 --- a/dos/doswrite.f90 +++ b/dos/doswrite.f90 @@ -6,17 +6,18 @@ MODULE m_doswrite ! CONTAINS SUBROUTINE doswrite(& - & DIMENSION,kpts,atoms,vacuum,& + & eig_id,DIMENSION,kpts,atoms,vacuum,& & input,banddos,& & sliceplot,noco,sym,& & cell,& & l_mcd,ncored,ncore,e_mcd,& & efermi,nsld,oneD) - + USE m_eig66_io,ONLY:read_dos USE m_evaldos USE m_cdninf USE m_types IMPLICIT NONE + TYPE(t_dimension),INTENT(IN) :: DIMENSION TYPE(t_oneD),INTENT(IN) :: oneD TYPE(t_banddos),INTENT(IN) :: banddos @@ -31,7 +32,7 @@ CONTAINS ! .. Scalar Arguments .. INTEGER,PARAMETER :: n2max=13 - INTEGER, INTENT (IN) :: nsld + INTEGER, INTENT (IN) :: nsld,eig_id INTEGER, INTENT (IN) :: ncored REAL, INTENT (IN) :: efermi LOGICAL, INTENT (IN) :: l_mcd @@ -98,10 +99,10 @@ CONTAINS ENDIF DO ikpt=1,kpts%nkpt - READ (84,rec=kpts%nkpt*(kspin-1)+ikpt) bkpt,wk,ne,eig,& + call read_dos(eig_id,ikpt,kspin,& & qal(:,:,:,kspin),qvac(:,:,ikpt,kspin),& & qis(:,ikpt,kspin),& - & qvlay(:,:,:),qstars,ksym,jsym + & qvlay(:,:,:),qstars,ksym,jsym) CALL cdninf(& & input,sym,noco,kspin,atoms,& @@ -127,7 +128,7 @@ CONTAINS ! IF (banddos%dos.AND.(banddos%ndir.LT.0)) THEN CALL evaldos(& - & input,banddos,vacuum,kpts,atoms,sym,noco,oneD,cell,& + & eig_id,input,banddos,vacuum,kpts,atoms,sym,noco,oneD,cell,& & DIMENSION,efermi,& & l_mcd,ncored,ncore,e_mcd,nsld) ENDIF diff --git a/dos/evaldos.f90 b/dos/evaldos.f90 index accc5829519e4cb0277247602b8a45d91c5ef711..49ad7555b1ca7b9c54f8d21eacd43c2a0621c415 100644 --- a/dos/evaldos.f90 +++ b/dos/evaldos.f90 @@ -1,6 +1,6 @@ MODULE m_evaldos CONTAINS - SUBROUTINE evaldos(input,banddos,vacuum,kpts,atoms,sym,noco,oneD,cell,& + SUBROUTINE evaldos(eig_id,input,banddos,vacuum,kpts,atoms,sym,noco,oneD,cell,& dimension,efermiarg, l_mcd,ncored,ncore,e_mcd,nsld) !---------------------------------------------------------------------- ! @@ -19,6 +19,7 @@ ! ntb=max(nevk) ! !---------------------------------------------------------------------- + USE m_eig66_io,ONLY:read_dos USE m_triang USE m_maketetra USE m_tetrados @@ -27,7 +28,7 @@ USE m_smooth USE m_types IMPLICIT NONE - + INTEGER,INTENT(IN) :: eig_id TYPE(t_dimension),INTENT(IN) :: dimension TYPE(t_oneD),INTENT(IN) :: oneD TYPE(t_banddos),INTENT(IN) :: banddos @@ -62,7 +63,7 @@ REAL vk(3,kpts%nkptd),wt(kpts%nkptd),voltet(6*kpts%nkpt),kx(kpts%nkpt),vkr(3,kpts%nkpt) REAL ev(dimension%neigd,kpts%nkptd),e(ned),gpart(ned,atoms%ntype),atr(2*kpts%nkpt) REAL e_grid(ned+1),spect(ned,3*atoms%ntypd),ferwe(dimension%neigd,kpts%nkptd) - REAL, ALLOCATABLE :: qal(:,:,:),qval(:,:,:),qlay(:,:,:),g(:,:) + REAL, ALLOCATABLE :: qal(:,:,:),qval(:,:,:),qlay(:,:,:),g(:,:),qal_tmp(:,:,:),qis(:),qvlay(:,:,:) REAL, ALLOCATABLE :: mcd(:,:,:),orbcomp(:,:,:),qmtp(:,:) REAL, ALLOCATABLE :: qintsl(:,:),qmtsl(:,:),qvac(:,:) COMPLEX, ALLOCATABLE :: qstars(:,:,:,:) @@ -155,44 +156,35 @@ ENDDO ENDDO ! -! read data from dos_tmp file! +! read data from file! ! - IF (( .not.l_mcd ).AND.(banddos%ndir.NE.-3)) THEN - READ(84,rec=kpts%nkpt*(jspin-1)+k) vk(:,k),wt(k),nevk(k),ev(:,k),& - qal(1:lmax*atoms%ntype,:,k), & ! atoms 1...atoms%ntype& - qal(lmax*atoms%ntype+2,:,k), & ! vacuum 1& - qal(lmax*atoms%ntype+3,:,k), & ! vacuum 2& - qal(lmax*atoms%ntype+1,:,k), & ! interstitial& - qlay,qstars - ELSEIF (banddos%ndir.NE.-3) THEN ALLOCATE( ksym(dimension%neigd),jsym(dimension%neigd) ) - READ(84,rec=kpts%nkpt*(jspin-1)+k) vk(:,k),wt(k),nevk(k),ev(:,k),& - qal(1:lmax*atoms%ntype,:,k), & ! atoms 1...atoms%ntype& - qal(lmax*atoms%ntype+2,:,k), & ! vacuum 1& - qal(lmax*atoms%ntype+3,:,k), & ! vacuum 2& - qal(lmax*atoms%ntype+1,:,k), & ! interstitial& - qlay,qstars,ksym,jsym,& - mcd(:,:,k) + ALLOCATE( qal_tmp(1:lmax,atoms%ntype,dimension%neigd)) + ALLOCATE( orbcomp(dimension%neigd,23,atoms%natd),qintsl(nsld,dimension%neigd)) + ALLOCATE( qmtsl(nsld,dimension%neigd),qmtp(dimension%neigd,atoms%natd),qvac(dimension%neigd,2)) + ALLOCATE( qis(dimension%neigd),qvlay(dimension%neigd,vacuum%layerd,2)) + CALL read_dos(eig_id,k,jspin,qal_tmp,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + qal(1:lmax*atoms%ntype,:,k)=reshape(qal_tmp,(/lmax*atoms%ntype,size(qal_tmp,3)/)) + qal(lmax*atoms%ntype+2,:,k)=qvac(:,1) ! vacuum 1 + qal(lmax*atoms%ntype+3,:,k)=qvac(:,2) ! vacuum 2 + qal(lmax*atoms%ntype+1,:,k)=qis ! interstitial DEALLOCATE( ksym,jsym ) - ELSE - ALLOCATE( orbcomp(dimension%neigd,23,atoms%natd),qintsl(nsld,dimension%neigd)) - ALLOCATE( qmtsl(nsld,dimension%neigd),qmtp(dimension%neigd,atoms%natd),qvac(dimension%neigd,2)) - READ (129,rec=kpts%nkpt*(jspin-1)+k) vk(:,k),wt(k),nevk(k),ev(:,k),qvac,qintsl,qmtsl,orbcomp,qmtp - IF (n_orb == 0) THEN - qal(1:nsld,:,k) = qintsl(:,:) - qal(nsld+1:2*nsld,:,k) = qmtsl(:,:) - ELSE - DO i = 1, 23 - DO l = 1, nevk(k) - qal(i,l,k) = orbcomp(l,i,n_orb)*qmtp(l,n_orb)/10000. - ENDDO - DO l = nevk(k)+1, dimension%neigd - qal(i,l,k) = 0.0 - ENDDO - ENDDO - ENDIF - DEALLOCATE( orbcomp,qintsl,qmtsl,qmtp,qvac) - ENDIF + IF (l_orbcomp)THEN + IF (n_orb == 0) THEN + qal(1:nsld,:,k) = qintsl(:,:) + qal(nsld+1:2*nsld,:,k) = qmtsl(:,:) + ELSE + DO i = 1, 23 + DO l = 1, nevk(k) + qal(i,l,k) = orbcomp(l,i,n_orb)*qmtp(l,n_orb)/10000. + ENDDO + DO l = nevk(k)+1, dimension%neigd + qal(i,l,k) = 0.0 + ENDDO + ENDDO + ENDIF + ENDIF + DEALLOCATE( orbcomp,qintsl,qmtsl,qmtp,qvac,qis,qal_tmp,qvlay) ntb = max(ntb,nevk(k)) ! ! set vacuum partial charge zero, if bulk calculation diff --git a/io/eig66_da.F90 b/io/eig66_da.F90 index e47e1306b35383833f330f88c27f76e8308c6455..6382616ed968db0b603d168182ff112d1dbb1473 100644 --- a/io/eig66_da.F90 +++ b/io/eig66_da.F90 @@ -1,281 +1,364 @@ MODULE m_eig66_da #include "juDFT_env.h" - ! Do the IO of the eig-file in fortran direct-access - ! The eig-file is split into two parts: - ! eig.bas contains the basis-set information - ! eig.vec contains the eigenvalues and the eigenvectors - ! The record number is given by nrec=nk+(jspin-1)*nkpts - ! each record contains: - ! eig.bas: el,evac,ello,bkpt,wtkpt,nv,nmat,k1,k2,k3,kveclo - ! eig.vec: ne,eig,z** - !**: real or complex depending on calculation type - USE m_eig66_data - IMPLICIT NONE + ! Do the IO of the eig-file in fortran direct-access + ! The eig-file is split into two parts: + ! eig.bas contains the basis-set information + ! eig.vec contains the eigenvalues and the eigenvectors + ! The record number is given by nrec=nk+(jspin-1)*nkpts + ! each record contains: + ! eig.bas: el,evac,ello,bkpt,wtkpt,nv,nmat,k1,k2,k3,kveclo + ! eig.vec: ne,eig,z** + !**: real or complex depending on calculation type + USE m_eig66_data + IMPLICIT NONE CONTAINS - subroutine priv_find_data(id,d) - INTEGER,INTENT(IN) :: id - TYPE(t_data_DA),POINTER,INTENT(out) :: d - - class(t_data),pointer ::dp - call eig66_find_data(dp,id) - select type(dp) - type is (t_data_da) - d=>dp - class default - call judft_error("BUG: wrong datatype in eig66_da") - END SELECT - END subroutine - - SUBROUTINE open_eig(id,nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,create,filename) - INTEGER, INTENT(IN) :: id,nmat,neig,nkpts,jspins,nlo,ntype,lmax,nlotot - LOGICAL, INTENT(IN) :: create - CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename - !locals - LOGICAL :: l_file - INTEGER :: i1,recl_z,recl_eig - REAL :: r1,r3(3) - COMPLEX :: c1 - TYPE(t_data_DA),POINTER:: d - - call priv_find_data(id,d) - - IF (present(filename)) d%fname=filename - call eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype) - - !Allocate the storage for the DATA always read/write - ALLOCATE(d%el_s(0:lmax,ntype),d%ello_s(nlo,ntype),d%evac_s(2)) - ALLOCATE(d%kvec_s(nmat,3),d%kveclo_s(nlotot)) - !Calculate the record length - INQUIRE(IOLENGTH=recl_eig) d%el_s,d%evac_s,d%ello_s,r3,r1,i1,i1,d%kvec_s,d%kveclo_s - d%recl_bas=recl_eig - INQUIRE(IOLENGTH=recl_eig) r1 - recl_eig=recl_eig*(neig+2) ! add a 2 for integer 'neig' + SUBROUTINE priv_find_data(id,d) + INTEGER,INTENT(IN) :: id + TYPE(t_data_DA),POINTER,INTENT(out) :: d + + CLASS(t_data),POINTER ::dp + CALL eig66_find_data(dp,id) + SELECT TYPE(dp) + TYPE is (t_data_da) + d=>dp + CLASS default + CALL judft_error("BUG: wrong datatype in eig66_da") + END SELECT + END SUBROUTINE priv_find_data + + SUBROUTINE open_eig(id,nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,create,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + INTEGER, INTENT(IN) :: id,nmat,neig,nkpts,jspins,nlo,ntype,lmax,nlotot + LOGICAL, INTENT(IN) :: create + LOGICAL,INTENT(IN) :: l_dos,l_mcd,l_orb + CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename + INTEGER,INTENT(IN),OPTIONAL :: layers,nstars,ncored,nsld,nat + !locals + LOGICAL :: l_file + INTEGER :: i1,recl_z,recl_eig,recl_dos + REAL :: r1,r3(3) + COMPLEX :: c1 + TYPE(t_data_DA),POINTER:: d + + CALL priv_find_data(id,d) + + IF (PRESENT(filename)) d%fname=filename + CALL eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype,l_dos,l_mcd,l_orb) + + !Allocate the storage for the DATA always read/write + ALLOCATE(d%el_s(0:lmax,ntype),d%ello_s(nlo,ntype),d%evac_s(2)) + ALLOCATE(d%kvec_s(nmat,3),d%kveclo_s(nlotot)) + !Calculate the record length + INQUIRE(IOLENGTH=recl_eig) d%el_s,d%evac_s,d%ello_s,r3,r1,i1,i1,d%kvec_s,d%kveclo_s + d%recl_bas=recl_eig + INQUIRE(IOLENGTH=recl_eig) r1 + recl_eig=recl_eig*(neig+2) ! add a 2 for integer 'neig' #if ( defined(CPP_INVERSION) && !defined(CPP_SOC) ) - INQUIRE(IOLENGTH=recl_z) r1 + INQUIRE(IOLENGTH=recl_z) r1 #else - INQUIRE(IOLENGTH=recl_z) c1 + INQUIRE(IOLENGTH=recl_z) c1 #endif - recl_z=recl_z*nmat*neig - - d%recl_vec=recl_eig+recl_z - - - IF (create) THEN - d%file_io_id_bas=priv_free_uid() - INQUIRE(file=trim(d%fname)//".bas",opened=l_file) - DO while(l_file) - print *,"eig66_open_da:",d%fname," in use" - d%fname=trim(d%fname)//"6" - INQUIRE(file=trim(d%fname)//".bas",opened=l_file) - ENDDO - OPEN(d%file_io_id_bas,FILE=trim(d%fname)//".bas",ACCESS='direct',FORM='unformatted',RECL=d%recl_bas,STATUS='unknown') - d%file_io_id_vec=priv_free_uid() - OPEN(d%file_io_id_vec,FILE=trim(d%fname)//".vec",ACCESS='direct',FORM='unformatted',RECL=d%recl_vec,STATUS='unknown') - ELSE - d%file_io_id_bas=priv_free_uid() - OPEN(d%file_io_id_bas,FILE=trim(d%fname)//".bas",ACCESS='direct',FORM='unformatted',RECL=d%recl_bas,STATUS='old') - d%file_io_id_vec=priv_free_uid() - OPEN(d%file_io_id_vec,FILE=trim(d%fname)//".vec",ACCESS='direct',FORM='unformatted',RECL=d%recl_vec,STATUS='old') - ENDIF - CONTAINS - INTEGER FUNCTION priv_free_uid() RESULT(uid) - IMPLICIT NONE - LOGICAL::used - used=.TRUE. - uid=665 - DO WHILE(used) - uid=uid+1 - INQUIRE(UNIT=uid,OPENED=used) - END DO - END FUNCTION - END SUBROUTINE open_eig - SUBROUTINE close_eig(id,filename) - INTEGER,INTENT(IN)::id - CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename - type(t_data_DA),pointer:: d - - call priv_find_data(id,d) - - DEALLOCATE(d%el_s,d%ello_s,d%evac_s,d%kvec_s,d%kveclo_s) - CLOSE(d%file_io_id_bas) - CLOSE(d%file_io_id_vec) - d%recl_vec=0 - d%recl_bas=0 - - !If a filename was given and the name is not the current filename then rename - IF (present(filename)) THEN - IF (filename.NE.d%fname) THEN - CALL system("mv "//trim(d%fname)//".bas "//trim(filename)//".bas") - CALL system("mv "//trim(d%fname)//".vec "//trim(filename)//".vec") - ENDIF - ENDIF - d%fname="eig" - CALL eig66_remove_data(id) - END SUBROUTINE close_eig - SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) - IMPLICIT NONE - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat - INTEGER, INTENT(OUT),OPTIONAL :: neig - REAL, INTENT(OUT),OPTIONAL :: eig(:) - INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) - REAL, INTENT(OUT),OPTIONAL :: bk(:),wk - INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end - CLASS(*),TARGET,INTENT(OUT),OPTIONAL :: z(:,:) - - !Local variables - INTEGER:: nv_s,nmat_s,n,nrec,neig_s - REAL :: bkpt(3),wtkpt - REAL,ALLOCATABLE::eig_s(:),zr_s(:,:) - COMPLEX,ALLOCATABLE::zc_s(:,:) - REAL,POINTER :: zr(:,:) - COMPLEX,POINTER:: zc(:,:) - type(t_data_DA),POINTER:: d - - - - call priv_find_data(id,d) - ! check if io is performed correctly - IF (present(n_start)) THEN - IF (n_start/=1) & - CALL juDFT_error("In direct access mode only all eigenstates can be read") - ENDIF - - nrec=nk+(jspin-1)*d%nkpts - IF (present(el).OR.present(ello).OR.present(evac).OR.present(bk).OR.present(wk).OR.& - present(nv).OR.present(nmat).OR.present(k1).OR.present(k2).OR.present(k3).OR.& - present(kveclo)) THEN - !IO of basis-set information - READ(d%file_io_id_bas,REC=nrec) nmat_s,d%el_s,d%evac_s,d%ello_s,bkpt,wtkpt,nv_s,d%kvec_s,d%kveclo_s - IF (present(el)) el=d%el_s - IF (present(evac)) evac=d%evac_s - IF (present(ello)) ello=d%ello_s - IF (present(bk)) bk=bkpt - IF (present(wk)) wk=wtkpt - IF (present(nv)) nv=nv_s - IF (present(nmat)) nmat=nmat_s - IF (present(k1)) k1=d%kvec_s(:,1) - IF (present(k2)) k2=d%kvec_s(:,2) - IF (present(k3)) k3=d%kvec_s(:,3) - IF (present(kveclo)) kveclo=d%kveclo_s - ENDIF - - IF (.NOT.(present(eig).OR.present(neig).OR.present(z))) RETURN - - IF (.NOT.(present(eig).OR.present(z))) THEN - READ(d%file_io_id_vec,REC=nrec) neig - RETURN - ENDIF - IF (present(eig)) THEN - ALLOCATE(eig_s(size(eig))) - ENDIF - IF (present(z).and..not.present(eig)) THEN - ALLOCATE(eig_s(size(z,2))) - ENDIF - IF (present(z)) THEN - SELECT TYPE(z) - TYPE IS (real) - INQUIRE(IOLENGTH=n) neig_s,eig_s,real(z) - IF (n>d%recl_vec) CALL juDFT_error("BUG: Too long record") - zr=>z - READ(d%file_io_id_vec,REC=nrec) neig_s,eig_s,zr - TYPE IS (complex) - INQUIRE(IOLENGTH=n) neig_s,eig_s,cmplx(z) - IF (n>d%recl_vec) CALL juDFT_error("BUG: Too long record") - zc=>z - READ(d%file_io_id_vec,REC=nrec) neig_s,eig_s,zc - END SELECT - ELSE - INQUIRE(IOLENGTH=n) neig_s,eig_s - IF (n>d%recl_vec) CALL juDFT_error("BUG: Too long record") - READ(d%file_io_id_vec,REC=nrec) neig_s,eig_s - ENDIF - IF (present(neig)) neig=neig_s - IF (present(eig)) eig=eig_s - - END SUBROUTINE read_eig - - SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk, & - eig,el,ello,evac,nlotot,kveclo,n_size,n_rank,z) - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank - REAL, INTENT(IN),OPTIONAL :: wk - INTEGER, INTENT(IN),OPTIONAL :: neig,nv,nmat,nlotot,neig_total - INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:) - REAL, INTENT(IN),OPTIONAL :: evac(:),ello(:,:) - CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) - - INTEGER:: nrec,r_len - INTEGER:: nv_s,nmat_s - REAL :: bkpt(3),wtkpt - type(t_data_DA),POINTER:: d - - call priv_find_data(id,d) - !This mode requires all data to be written at once!! - - IF (present(n_size).AND.present(n_rank)) THEN - IF (n_size/=1.OR.n_rank/=0) & - CALL juDFT_error("Direct Access IO not possible in eigenvalue parallel code") - ENDIF - !check record length - !INQUIRE(iolength=r_len) nmat,el,evac,ello,bk,wk,nv,d%kvec_s,kveclo - !if (r_len>recl_bas) call juDFT_error("BUG: too long record") - - !Now it is time for the IO :-) - nrec=nk+(jspin-1)*d%nkpts - IF (present(nmat).AND..NOT.present(el)) THEN - !IO of basis-set information - READ(d%file_io_id_bas,REC=nrec,ERR=88) nmat_s,d%el_s,d%evac_s,d%ello_s,bkpt,wtkpt,nv_s,d%kvec_s,d%kveclo_s -88 WRITE(d%file_io_id_bas,REC=nrec) nmat ,d%el_s,d%evac_s,d%ello_s,bkpt,wtkpt,nv_s,d%kvec_s,d%kveclo_s - IF (present(wk).OR.present(nv).OR.present(nlotot) & - .OR.present(k1).OR.present(k2).OR.present(k3).OR.present(kveclo).OR.& - present(bk).OR.present(ello).OR.present(evac)) THEN - CALL juDFT_error("BUG:Direct access IO of eig-file only with all scalar data") - ENDIF - ELSE - IF (.NOT.(present(wk).AND.present(nv).AND.present(nmat).AND.present(nlotot) & - .AND.present(k1).AND.present(k2).AND.present(k3).AND.present(kveclo).AND.& - present(bk).AND.present(el).AND.present(ello).AND.present(evac))) THEN - CALL juDFT_error("BUG:Direct access IO of eig-file only with all data") - ENDIF - d%kvec_s(:,1)=k1 - d%kvec_s(:,2)=k2 - d%kvec_s(:,3)=k3 - if ((size(el).ne.size(d%el_s)).or.(size(ello).ne.size(d%ello_s).or.(size(evac).ne.size(d%evac_s)))) THEN - write(*,*) shape(el),shape(d%el_s) - write(*,*) shape(ello),shape(d%ello_s) - write(*,*) shape(evac),shape(d%evac_s) - call juDFT_error("Mismatch of sizes") - endif - WRITE(d%file_io_id_bas,REC=nrec) nmat,el,evac,ello,bk,wk,nv,d%kvec_s,kveclo - - ENDIF - IF (present(neig).AND.present(neig_total)) THEN - IF (neig.NE.neig_total) THEN - CALL juDFT_error("Neig and neig_total have to be equal in DA mode",calledby="eig66_da") - ENDIF - ENDIF - IF (.NOT.present(eig).OR..NOT.present(neig)) RETURN - !Now the IO of the eigenvalues/vectors - IF (present(z)) THEN - SELECT TYPE(z) - TYPE IS (real) - INQUIRE(IOLENGTH=r_len) neig,eig,real(z) - IF (r_len>d%recl_vec) CALL juDFT_error("BUG: too long record") - WRITE(d%file_io_id_vec,REC=nrec) neig,eig,real(z) - TYPE IS (complex) - INQUIRE(IOLENGTH=r_len) neig,eig,cmplx(z) - IF (r_len>d%recl_vec) CALL juDFT_error("BUG: too long record") - WRITE(d%file_io_id_vec,REC=nrec) neig,eig,cmplx(z) - END SELECT - ELSE - INQUIRE(IOLENGTH=r_len) neig,eig - IF (r_len>d%recl_vec) CALL juDFT_error("BUG: too long record") - WRITE(d%file_io_id_vec,REC=nrec) neig,eig - ENDIF - - END SUBROUTINE write_eig + recl_z=recl_z*nmat*neig + + d%recl_vec=recl_eig+recl_z + + IF (l_dos) THEN + IF (.NOT.(PRESENT(layers).AND.PRESENT(nstars).AND.PRESENT(ncored).AND.PRESENT(nsld).AND.PRESENT(nat))) & + CALL judft_error("BUG:Could not open file for DOS-data",calledby="eig66_da") + INQUIRE(IOLENGTH=i1) i1 + recl_dos=i1*2*neig !ksym&jsym + INQUIRE(IOLENGTH=i1) r1 + recl_dos=recl_dos+i1*3*neig !qvac&qis + recl_dos=recl_dos+i1*4*ntype !qal + recl_dos=recl_dos+i1*neig*2*layers !qvlay + IF (l_orb) THEN + recl_dos=recl_dos+i1*2*nsld*neig !qintsl,qmtsl + recl_dos=recl_dos+i1*24*neig*nat !qmtp,orbcomp + ENDIF + INQUIRE(IOLENGTH=i1) c1 + recl_dos=recl_dos+i1*nstars*neig*layers*2 !qstars + IF (l_mcd) recl_dos=recl_dos+i1*3*ntype*ncored*neig !mcd + ELSE + recl_dos=-1 + ENDIF + d%recl_dos=recl_dos + + + IF (create) THEN + d%file_io_id_bas=priv_free_uid() + INQUIRE(file=TRIM(d%fname)//".bas",opened=l_file) + DO WHILE(l_file) + PRINT *,"eig66_open_da:",d%fname," in use" + d%fname=TRIM(d%fname)//"6" + INQUIRE(file=TRIM(d%fname)//".bas",opened=l_file) + ENDDO + OPEN(d%file_io_id_bas,FILE=TRIM(d%fname)//".bas",ACCESS='direct',FORM='unformatted',RECL=d%recl_bas,STATUS='unknown') + d%file_io_id_vec=priv_free_uid() + OPEN(d%file_io_id_vec,FILE=TRIM(d%fname)//".vec",ACCESS='direct',FORM='unformatted',RECL=d%recl_vec,STATUS='unknown') + IF(d%recl_dos>0) THEN + d%file_io_id_dos=priv_free_uid() + OPEN(d%file_io_id_dos,FILE=TRIM(d%fname)//".dos",ACCESS='direct',FORM='unformatted',RECL=d%recl_dos,STATUS='unknown') + ENDIF + + ELSE + d%file_io_id_bas=priv_free_uid() + OPEN(d%file_io_id_bas,FILE=TRIM(d%fname)//".bas",ACCESS='direct',FORM='unformatted',RECL=d%recl_bas,STATUS='old') + d%file_io_id_vec=priv_free_uid() + OPEN(d%file_io_id_vec,FILE=TRIM(d%fname)//".vec",ACCESS='direct',FORM='unformatted',RECL=d%recl_vec,STATUS='old') + IF(d%recl_dos>0) THEN + d%file_io_id_dos=priv_free_uid() + OPEN(d%file_io_id_dos,FILE=TRIM(d%fname)//".dos",ACCESS='direct',FORM='unformatted',RECL=d%recl_dos,STATUS='old') + ENDIF + ENDIF + CONTAINS + INTEGER FUNCTION priv_free_uid() RESULT(uid) + IMPLICIT NONE + LOGICAL::used + used=.TRUE. + uid=665 + DO WHILE(used) + uid=uid+1 + INQUIRE(UNIT=uid,OPENED=used) + END DO + END FUNCTION priv_free_uid + END SUBROUTINE open_eig + SUBROUTINE close_eig(id,filename) + INTEGER,INTENT(IN)::id + CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename + TYPE(t_data_DA),POINTER:: d + + CALL priv_find_data(id,d) + + DEALLOCATE(d%el_s,d%ello_s,d%evac_s,d%kvec_s,d%kveclo_s) + CLOSE(d%file_io_id_bas) + CLOSE(d%file_io_id_vec) + d%recl_vec=0 + d%recl_bas=0 + + !If a filename was given and the name is not the current filename then rename + IF (PRESENT(filename)) THEN + IF (filename.NE.d%fname) THEN + CALL system("mv "//TRIM(d%fname)//".bas "//TRIM(filename)//".bas") + CALL system("mv "//TRIM(d%fname)//".vec "//TRIM(filename)//".vec") + ENDIF + ENDIF + d%fname="eig" + CALL eig66_remove_data(id) + END SUBROUTINE close_eig + SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat + INTEGER, INTENT(OUT),OPTIONAL :: neig + REAL, INTENT(OUT),OPTIONAL :: eig(:) + INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) + REAL, INTENT(OUT),OPTIONAL :: bk(:),wk + INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end + CLASS(*),TARGET,INTENT(OUT),OPTIONAL :: z(:,:) + + !Local variables + INTEGER:: nv_s,nmat_s,n,nrec,neig_s + REAL :: bkpt(3),wtkpt + REAL,ALLOCATABLE::eig_s(:),zr_s(:,:) + COMPLEX,ALLOCATABLE::zc_s(:,:) + REAL,POINTER :: zr(:,:) + COMPLEX,POINTER:: zc(:,:) + TYPE(t_data_DA),POINTER:: d + + + + CALL priv_find_data(id,d) + ! check if io is performed correctly + IF (PRESENT(n_start)) THEN + IF (n_start/=1) & + CALL juDFT_error("In direct access mode only all eigenstates can be read") + ENDIF + + nrec=nk+(jspin-1)*d%nkpts + IF (PRESENT(el).OR.PRESENT(ello).OR.PRESENT(evac).OR.PRESENT(bk).OR.PRESENT(wk).OR.& + PRESENT(nv).OR.PRESENT(nmat).OR.PRESENT(k1).OR.PRESENT(k2).OR.PRESENT(k3).OR.& + PRESENT(kveclo)) THEN + !IO of basis-set information + READ(d%file_io_id_bas,REC=nrec) nmat_s,d%el_s,d%evac_s,d%ello_s,bkpt,wtkpt,nv_s,d%kvec_s,d%kveclo_s + IF (PRESENT(el)) el=d%el_s + IF (PRESENT(evac)) evac=d%evac_s + IF (PRESENT(ello)) ello=d%ello_s + IF (PRESENT(bk)) bk=bkpt + IF (PRESENT(wk)) wk=wtkpt + IF (PRESENT(nv)) nv=nv_s + IF (PRESENT(nmat)) nmat=nmat_s + IF (PRESENT(k1)) k1=d%kvec_s(:,1) + IF (PRESENT(k2)) k2=d%kvec_s(:,2) + IF (PRESENT(k3)) k3=d%kvec_s(:,3) + IF (PRESENT(kveclo)) kveclo=d%kveclo_s + ENDIF + + IF (.NOT.(PRESENT(eig).OR.PRESENT(neig).OR.PRESENT(z))) RETURN + + IF (.NOT.(PRESENT(eig).OR.PRESENT(z))) THEN + READ(d%file_io_id_vec,REC=nrec) neig + RETURN + ENDIF + IF (PRESENT(eig)) THEN + ALLOCATE(eig_s(SIZE(eig))) + ENDIF + IF (PRESENT(z).AND..NOT.PRESENT(eig)) THEN + ALLOCATE(eig_s(SIZE(z,2))) + ENDIF + IF (PRESENT(z)) THEN + SELECT TYPE(z) + TYPE IS (REAL) + INQUIRE(IOLENGTH=n) neig_s,eig_s,REAL(z) + IF (n>d%recl_vec) CALL juDFT_error("BUG: Too long record") + zr=>z + READ(d%file_io_id_vec,REC=nrec) neig_s,eig_s,zr + TYPE IS (COMPLEX) + INQUIRE(IOLENGTH=n) neig_s,eig_s,CMPLX(z) + IF (n>d%recl_vec) CALL juDFT_error("BUG: Too long record") + zc=>z + READ(d%file_io_id_vec,REC=nrec) neig_s,eig_s,zc + END SELECT + ELSE + INQUIRE(IOLENGTH=n) neig_s,eig_s + IF (n>d%recl_vec) CALL juDFT_error("BUG: Too long record") + READ(d%file_io_id_vec,REC=nrec) neig_s,eig_s + ENDIF + IF (PRESENT(neig)) neig=neig_s + IF (PRESENT(eig)) eig=eig_s + + END SUBROUTINE read_eig + + SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk, & + eig,el,ello,evac,nlotot,kveclo,n_size,n_rank,z) + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank + REAL, INTENT(IN),OPTIONAL :: wk + INTEGER, INTENT(IN),OPTIONAL :: neig,nv,nmat,nlotot,neig_total + INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:) + REAL, INTENT(IN),OPTIONAL :: evac(:),ello(:,:) + CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) + + INTEGER:: nrec,r_len + INTEGER:: nv_s,nmat_s + REAL :: bkpt(3),wtkpt + TYPE(t_data_DA),POINTER:: d + + CALL priv_find_data(id,d) + !This mode requires all data to be written at once!! + + IF (PRESENT(n_size).AND.PRESENT(n_rank)) THEN + IF (n_size/=1.OR.n_rank/=0) & + CALL juDFT_error("Direct Access IO not possible in eigenvalue parallel code") + ENDIF + !check record length + !INQUIRE(iolength=r_len) nmat,el,evac,ello,bk,wk,nv,d%kvec_s,kveclo + !if (r_len>recl_bas) call juDFT_error("BUG: too long record") + + !Now it is time for the IO :-) + nrec=nk+(jspin-1)*d%nkpts + IF (PRESENT(nmat).AND..NOT.PRESENT(el)) THEN + !IO of basis-set information + READ(d%file_io_id_bas,REC=nrec,ERR=88) nmat_s,d%el_s,d%evac_s,d%ello_s,bkpt,wtkpt,nv_s,d%kvec_s,d%kveclo_s +88 WRITE(d%file_io_id_bas,REC=nrec) nmat ,d%el_s,d%evac_s,d%ello_s,bkpt,wtkpt,nv_s,d%kvec_s,d%kveclo_s + IF (PRESENT(wk).OR.PRESENT(nv).OR.PRESENT(nlotot) & + .OR.PRESENT(k1).OR.PRESENT(k2).OR.PRESENT(k3).OR.PRESENT(kveclo).OR.& + PRESENT(bk).OR.PRESENT(ello).OR.PRESENT(evac)) THEN + CALL juDFT_error("BUG:Direct access IO of eig-file only with all scalar data") + ENDIF + ELSE + IF (.NOT.(PRESENT(wk).AND.PRESENT(nv).AND.PRESENT(nmat).AND.PRESENT(nlotot) & + .AND.PRESENT(k1).AND.PRESENT(k2).AND.PRESENT(k3).AND.PRESENT(kveclo).AND.& + PRESENT(bk).AND.PRESENT(el).AND.PRESENT(ello).AND.PRESENT(evac))) THEN + CALL juDFT_error("BUG:Direct access IO of eig-file only with all data") + ENDIF + d%kvec_s(:,1)=k1 + d%kvec_s(:,2)=k2 + d%kvec_s(:,3)=k3 + IF ((SIZE(el).NE.SIZE(d%el_s)).OR.(SIZE(ello).NE.SIZE(d%ello_s).OR.(SIZE(evac).NE.SIZE(d%evac_s)))) THEN + WRITE(*,*) SHAPE(el),SHAPE(d%el_s) + WRITE(*,*) SHAPE(ello),SHAPE(d%ello_s) + WRITE(*,*) SHAPE(evac),SHAPE(d%evac_s) + CALL juDFT_error("Mismatch of sizes") + ENDIF + WRITE(d%file_io_id_bas,REC=nrec) nmat,el,evac,ello,bk,wk,nv,d%kvec_s,kveclo + + ENDIF + IF (PRESENT(neig).AND.PRESENT(neig_total)) THEN + IF (neig.NE.neig_total) THEN + CALL juDFT_error("Neig and neig_total have to be equal in DA mode",calledby="eig66_da") + ENDIF + ENDIF + IF (.NOT.PRESENT(eig).OR..NOT.PRESENT(neig)) RETURN + !Now the IO of the eigenvalues/vectors + IF (PRESENT(z)) THEN + SELECT TYPE(z) + TYPE IS (REAL) + INQUIRE(IOLENGTH=r_len) neig,eig,REAL(z) + IF (r_len>d%recl_vec) CALL juDFT_error("BUG: too long record") + WRITE(d%file_io_id_vec,REC=nrec) neig,eig,REAL(z) + TYPE IS (COMPLEX) + INQUIRE(IOLENGTH=r_len) neig,eig,CMPLX(z) + IF (r_len>d%recl_vec) CALL juDFT_error("BUG: too long record") + WRITE(d%file_io_id_vec,REC=nrec) neig,eig,CMPLX(z) + END SELECT + ELSE + INQUIRE(IOLENGTH=r_len) neig,eig + IF (r_len>d%recl_vec) CALL juDFT_error("BUG: too long record") + WRITE(d%file_io_id_vec,REC=nrec) neig,eig + ENDIF + + END SUBROUTINE write_eig + + + SUBROUTINE write_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(IN) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(IN) :: qstars(:,:,:,:) + INTEGER,INTENT(IN) :: ksym(:),jsym(:) + REAL,INTENT(IN),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(IN),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + TYPE(t_data_DA),POINTER:: d + INTEGER:: nrec + CALL priv_find_data(id,d) + + nrec=nk+(jspin-1)*d%nkpts + + IF (d%l_orb.AND.PRESENT(qmtsl)) THEN + IF (d%l_mcd) CPP_error("mcd & orbital decomposition not implemented in IO") + WRITE(d%file_io_id_dos,REC=nrec) qal,qvac,qis,qvlay,qstars,ksym,jsym,qintsl,qmtsl,qmtp,orbcomp + ELSEIF(d%l_mcd.AND.PRESENT(mcd)) THEN + WRITE(d%file_io_id_dos,REC=nrec) qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd + ELSE + WRITE(d%file_io_id_dos,REC=nrec) qal,qvac,qis,qvlay,qstars,ksym,jsym + END IF + END SUBROUTINE write_dos + + SUBROUTINE read_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(OUT) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(OUT) :: qstars(:,:,:,:) + INTEGER,INTENT(OUT) :: ksym(:),jsym(:) + REAL,INTENT(OUT),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(OUT),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + TYPE(t_data_DA),POINTER:: d + INTEGER:: nrec + CALL priv_find_data(id,d) + + nrec=nk+(jspin-1)*d%nkpts + + IF (d%l_orb.AND.PRESENT(qmtsl)) THEN + IF (d%l_mcd) CPP_error("mcd & orbital decomposition not implemented in IO") + READ(d%file_io_id_dos,REC=nrec) qal,qvac,qis,qvlay,qstars,ksym,jsym,qintsl,qmtsl,qmtp,orbcomp + ELSEIF(d%l_mcd.AND.PRESENT(mcd)) THEN + READ(d%file_io_id_dos,REC=nrec) qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd + ELSE + READ(d%file_io_id_dos,REC=nrec) qal,qvac,qis,qvlay,qstars,ksym,jsym + END IF + END SUBROUTINE read_dos + + END MODULE m_eig66_da diff --git a/io/eig66_data.F90 b/io/eig66_data.F90 index d758aa32557f4304c350168e69a104a415642444..31eed0a07dd30626c90e094dd61d9a372fe209b5 100644 --- a/io/eig66_data.F90 +++ b/io/eig66_data.F90 @@ -8,25 +8,32 @@ module m_eig66_data TYPE :: t_data INTEGER:: io_mode INTEGER:: jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype + LOGICAL:: l_dos,l_mcd,l_orb END TYPE TYPE,EXTENDS(t_data):: t_data_DA REAL, ALLOCATABLE:: el_s(:,:),ello_s(:,:),evac_s(:) INTEGER,ALLOCATABLE:: kvec_s(:,:),kveclo_s(:) - INTEGER :: recl_bas=0,recl_vec=0 + INTEGER :: recl_bas=0,recl_vec=0,recl_dos CHARACTER(LEN=20) :: fname="eig" - INTEGER :: file_io_id_bas,file_io_id_vec + INTEGER :: file_io_id_bas,file_io_id_vec,file_io_id_dos END TYPE TYPE,extends(t_data):: t_data_MPI INTEGER :: n_size=1 INTEGER :: size_k,size_el,size_ello,size_eig INTEGER :: int_handle,real_handle,eig_handle,zr_handle,zc_handle,neig_handle + INTEGER :: qal_handle,qvac_handle,qis_handle,qvlay_handle,qintsl_handle,qmtsl_handle + INTEGER :: qmtp_handle,orbcomp_handle,qstars_handle,mcd_handle,jsym_handle,ksym_handle INTEGER,ALLOCATABLE :: pe_basis(:,:),slot_basis(:,:) INTEGER,ALLOCATABLE :: pe_ev(:,:,:),slot_ev(:,:,:) INTEGER :: irank INTEGER,POINTER :: neig_data(:),int_data(:) REAL,POINTER :: eig_data(:),zr_data(:),real_data(:) + REAL,POINTER :: qal_data(:),qvac_data(:),qis_data(:),qvlay_data(:) + REAL,POINTER :: qintsl_data(:),qmtsl_data(:),qmtp_data(:),orbcomp_data(:),mcd_data(:) + COMPLEX,POINTER :: qstars_data(:) + INTEGER,POINTER :: jsym_data(:),ksym_data(:) COMPLEX,POINTER :: zc_data(:) END TYPE TYPE,EXTENDS(t_data):: t_data_hdf @@ -36,6 +43,9 @@ module m_eig66_data INTEGER(HID_T) :: bksetid,wksetid,ksetid INTEGER(HID_T) :: neigsetid,nvsetid,nmatsetid INTEGER(HID_T) :: energysetid,evsetid + INTEGER(HID_T) :: qalsetid,qvacsetid,qissetid,qvlaysetid + INTEGER(HID_T) :: qstarssetid,ksymsetid,jsymsetid,mcdsetid + INTEGER(HID_T) :: qintslsetid,qmtslsetid,qmtpstid,orbcompsetid CHARACTER(LEN=20) :: fname="eig" #endif END TYPE @@ -46,6 +56,18 @@ module m_eig66_data REAL,ALLOCATABLE :: eig_eig(:,:) REAL,ALLOCATABLE :: eig_vecr(:,:) COMPLEX,ALLOCATABLE :: eig_vecc(:,:) + REAL,ALLOCATABLE :: qal(:,:,:,:) + REAL,ALLOCATABLE :: qvac(:,:,:) + REAL,ALLOCATABLE :: qis(:,:) + REAL,ALLOCATABLE :: qvlay(:,:,:,:) + COMPLEX,ALLOCATABLE :: qstars(:,:,:,:,:) + INTEGER,ALLOCATABLE :: ksym(:,:) + INTEGER,ALLOCATABLE :: jsym(:,:) + REAL,ALLOCATABLE :: mcd(:,:,:,:) + REAL,ALLOCATABLE :: qintsl(:,:,:) + REAL,ALLOCATABLE :: qmtsl(:,:,:) + REAL,ALLOCATABLE :: qmtp(:,:,:) + REAL,ALLOCATABLE :: orbcomp(:,:,:,:) END TYPE TYPE t_list @@ -62,9 +84,10 @@ module m_eig66_data contains - subroutine eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype) + subroutine eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype,l_dos,l_mcd,l_orb) CLASS(t_data)::d INTEGER,INTENT(IN)::jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype + LOGICAL,INTENT(IN)::l_dos,l_mcd,l_orb d%jspins=jspins d%nkpts=nkpts d%nmat=nmat @@ -73,6 +96,9 @@ module m_eig66_data d%nlotot=nlotot d%nlo=nlo d%ntype=ntype + d%l_dos=l_dos + d%l_mcd=l_mcd + d%l_orb=l_orb END SUBROUTINE subroutine eig66_find_data(d,id,io_mode) diff --git a/io/eig66_hdf.F90 b/io/eig66_hdf.F90 index 6be19efbba9447591e753d520e4979901dc1c5ce..e25ad18cc3967b8a037871fc442418c94184569d 100644 --- a/io/eig66_hdf.F90 +++ b/io/eig66_hdf.F90 @@ -1,530 +1,675 @@ MODULE m_eig66_hdf #include "juDFT_env.h" - !***************************************************************** - ! DESC:Module for hdf-io of eig-file - ! To be compatible with f90 interface of HDF, use kind for vars - ! - ! !ATTENTION before calling openeig and after calling closeeig! - ! !the hdf library has to be initialized or finalized, respectively - ! - ! CONTAINS the following subroutines: - ! openeig opens file - ! closeeig closes file - ! read_keb reads kpt, enpara and basis data - ! read_neig read no of eigenvalues (and eigenvalues itself) - ! read_eig reads eigenvectors - ! writeeig saves all data for kpt - ! writesingleeig saves data for one kpt and energy - ! - ! - ! Daniel Wortmann, Tue Nov 512:07:522002 - !***************************************************************** - USE m_eig66_data + !***************************************************************** + ! DESC:Module for hdf-io of eig-file + ! To be compatible with f90 interface of HDF, use kind for vars + ! + ! !ATTENTION before calling openeig and after calling closeeig! + ! !the hdf library has to be initialized or finalized, respectively + ! + ! CONTAINS the following subroutines: + ! openeig opens file + ! closeeig closes file + ! read_keb reads kpt, enpara and basis data + ! read_neig read no of eigenvalues (and eigenvalues itself) + ! read_eig reads eigenvectors + ! writeeig saves all data for kpt + ! writesingleeig saves data for one kpt and energy + ! + ! + ! Daniel Wortmann, Tue Nov 512:07:522002 + !***************************************************************** + USE m_eig66_data #ifdef CPP_HDF - USE hdf5 - USE m_hdf_tools - IMPLICIT NONE + USE hdf5 + USE m_hdf_tools + IMPLICIT NONE - PRIVATE - INTEGER, PARAMETER :: one=1,two=2,three=3,zero=0 - !to have the correct - !type for array constructors + PRIVATE + INTEGER, PARAMETER :: one=1,two=2,three=3,zero=0 + !to have the correct + !type for array constructors #endif - PUBLIC open_eig,close_eig - PUBLIC read_eig - PUBLIC write_eig!,writesingleeig,writeeigc,writebas + PUBLIC open_eig,close_eig + PUBLIC read_eig + PUBLIC write_eig!,writesingleeig,writeeigc,writebas CONTAINS - subroutine priv_find_data(id,d) - INTEGER,INTENT(IN)::id - TYPE(t_data_hdf),pointer:: d - - class(t_data),pointer ::dp - call eig66_find_data(dp,id) - select type(dp) - type is (t_data_hdf) - d=>dp - class default - call judft_error("BUG: wrong datatype in eig66_hdf") - END SELECT - END subroutine - !---------------------------------------------------------------------- - SUBROUTINE open_eig(id,mpi_comm,nmat,neig,nkpts,jspins,lmax,nlo,ntype,& - & create,readonly,filename) - - !***************************************************************** - ! opens hdf-file for eigenvectors+values - !***************************************************************** - IMPLICIT NONE - - INTEGER, INTENT(IN) :: id,mpi_comm - INTEGER, INTENT(IN) :: nmat,neig,nkpts,jspins,nlo,ntype,lmax - LOGICAL, INTENT(IN) :: create,readonly - CHARACTER(LEN=*),OPTIONAL :: filename + SUBROUTINE priv_find_data(id,d) + INTEGER,INTENT(IN)::id + TYPE(t_data_hdf),POINTER:: d + + CLASS(t_data),POINTER ::dp + CALL eig66_find_data(dp,id) + SELECT TYPE(dp) + TYPE is (t_data_hdf) + d=>dp + CLASS default + CALL judft_error("BUG: wrong datatype in eig66_hdf") + END SELECT + END SUBROUTINE priv_find_data + !---------------------------------------------------------------------- + SUBROUTINE open_eig(id,mpi_comm,nmat,neig,nkpts,jspins,lmax,nlo,ntype,create,readonly,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + + !***************************************************************** + ! opens hdf-file for eigenvectors+values + !***************************************************************** + IMPLICIT NONE + + INTEGER, INTENT(IN) :: id,mpi_comm + INTEGER, INTENT(IN) :: nmat,neig,nkpts,jspins,nlo,ntype,lmax + LOGICAL, INTENT(IN) :: create,readonly,l_dos,l_mcd,l_orb + CHARACTER(LEN=*),OPTIONAL :: filename + INTEGER,INTENT(IN),OPTIONAL :: layers,nstars,ncored,nsld,nat #ifdef CPP_HDF - - INTEGER :: hdferr,access_mode - INTEGER(HID_T) :: creation_prp,access_prp,spaceid - LOGICAL :: l_exist - INTEGER(HSIZE_T):: dims(5) - TYPE(t_data_HDF),pointer::d - !Set creation and access properties + + INTEGER :: hdferr,access_mode + INTEGER(HID_T) :: creation_prp,access_prp,spaceid + LOGICAL :: l_exist + INTEGER(HSIZE_T):: dims(5) + TYPE(t_data_HDF),POINTER::d + !Set creation and access properties #ifdef CPP_MPI - INCLUDE 'mpif.h' - IF (readonly) THEN - access_prp=H5P_DEFAULT_f - creation_prp=H5P_DEFAULT_f - ELSE - CALL h5pcreate_f(H5P_FILE_ACCESS_F, access_prp, hdferr) - ! CALL h5pset_fapl_mpiposix_f(access_prp,MPI_COMM, - ! +.false.,hdferr) - CALL h5pset_fapl_mpio_f(access_prp, MPI_COMM, MPI_INFO_NULL,hdferr) - creation_prp=H5P_DEFAULT_f !no special creation property - ENDIF + INCLUDE 'mpif.h' + IF (readonly) THEN + access_prp=H5P_DEFAULT_f + creation_prp=H5P_DEFAULT_f + ELSE + CALL h5pcreate_f(H5P_FILE_ACCESS_F, access_prp, hdferr) + ! CALL h5pset_fapl_mpiposix_f(access_prp,MPI_COMM, + ! +.false.,hdferr) + CALL h5pset_fapl_mpio_f(access_prp, MPI_COMM, MPI_INFO_NULL,hdferr) + creation_prp=H5P_DEFAULT_f !no special creation property + ENDIF #else - access_prp=H5P_DEFAULT_f - creation_prp=H5P_DEFAULT_f + access_prp=H5P_DEFAULT_f + creation_prp=H5P_DEFAULT_f #endif - IF (present(filename)) d%fname=filename - call priv_find_data(id,d) - !set access_flags according - IF (readonly) THEN - access_mode=H5F_ACC_RDONLY_F - ELSE - access_mode=H5F_ACC_RDWR_F - ENDIF - ! OPEN FILE and get D%FID's - IF (create) THEN - INQUIRE(FILE=trim(d%fname)//'.hdf',EXIST=l_exist) - access_mode=H5F_ACC_TRUNC_F - ! IF (l_exist) WRITE (*,*)'Warning: eig.hdf was overwritten' - CALL h5fcreate_f(trim(d%fname)//'.hdf',access_Mode, d%fid, hdferr& - & ,creation_prp,access_prp) - ! create dataspaces and datasets - ! scalars - dims(:2)=(/nkpts,jspins/) - CALL h5screate_simple_f(2,dims(:2),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "neig", H5T_NATIVE_INTEGER, spaceid,& - & d%neigsetid, hdferr) - CALL h5dcreate_f(d%fid, "wk", H5T_NATIVE_DOUBLE,spaceid,& - & d%wksetid, hdferr) - CALL h5dcreate_f(d%fid, "nv", H5T_NATIVE_INTEGER, spaceid,& - & d%nvsetid, hdferr) - CALL h5dcreate_f(d%fid, "nmat", H5T_NATIVE_INTEGER, spaceid,& - & d%nmatsetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - ! vectors - dims(1:3)=(/two,nkpts,jspins/) - CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "evac", H5T_NATIVE_DOUBLE, spaceid,& - & d%evacsetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - dims(:3)=(/three,nkpts,jspins/) - CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "bk", H5T_NATIVE_DOUBLE, spaceid,& - & d%bksetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - dims(:3)=(/neig,nkpts,jspins/) - CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) - ! ew - CALL h5dcreate_f(d%fid, "energy", H5T_NATIVE_DOUBLE, spaceid,& - & d%energysetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - ! enparas - dims(1:4)=(/lmax+1,ntype,nkpts,jspins/) - CALL h5screate_simple_f(4,dims(1:4),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "el", H5T_NATIVE_DOUBLE, spaceid,& - & d%esetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - - dims(:4)=(/nlo,ntype,nkpts,jspins/) - CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "ello", H5T_NATIVE_DOUBLE, spaceid,& - & d%ellosetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - ! ev + CALL priv_find_data(id,d) + IF (PRESENT(filename)) d%fname=filename + CALL eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype,l_dos,l_mcd,l_orb) + !set access_flags according + IF (readonly) THEN + access_mode=H5F_ACC_RDONLY_F + ELSE + access_mode=H5F_ACC_RDWR_F + ENDIF + ! OPEN FILE and get D%FID's + IF (create) THEN + INQUIRE(FILE=TRIM(d%fname)//'.hdf',EXIST=l_exist) + access_mode=H5F_ACC_TRUNC_F + ! IF (l_exist) WRITE (*,*)'Warning: eig.hdf was overwritten' + CALL h5fcreate_f(TRIM(d%fname)//'.hdf',access_Mode, d%fid, hdferr ,creation_prp,access_prp) + ! create dataspaces and datasets + ! scalars + dims(:2)=(/nkpts,jspins/) + CALL h5screate_simple_f(2,dims(:2),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "neig", H5T_NATIVE_INTEGER, spaceid, d%neigsetid, hdferr) + CALL h5dcreate_f(d%fid, "wk", H5T_NATIVE_DOUBLE,spaceid, d%wksetid, hdferr) + CALL h5dcreate_f(d%fid, "nv", H5T_NATIVE_INTEGER, spaceid, d%nvsetid, hdferr) + CALL h5dcreate_f(d%fid, "nmat", H5T_NATIVE_INTEGER, spaceid, d%nmatsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + ! vectors + dims(1:3)=(/two,nkpts,jspins/) + CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "evac", H5T_NATIVE_DOUBLE, spaceid, d%evacsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:3)=(/three,nkpts,jspins/) + CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "bk", H5T_NATIVE_DOUBLE, spaceid, d%bksetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:3)=(/neig,nkpts,jspins/) + CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) + ! ew + CALL h5dcreate_f(d%fid, "energy", H5T_NATIVE_DOUBLE, spaceid, d%energysetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + ! enparas + dims(1:4)=(/lmax+1,ntype,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(1:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "el", H5T_NATIVE_DOUBLE, spaceid, d%esetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + + dims(:4)=(/nlo,ntype,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "ello", H5T_NATIVE_DOUBLE, spaceid, d%ellosetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + ! ev #if ( defined(CPP_INVERSION) && !defined(CPP_SOC) ) - dims(:5)=(/one,nmat,neig,nkpts,jspins/) + dims(:5)=(/one,nmat,neig,nkpts,jspins/) #else - dims(:5)=(/two,nmat,neig,nkpts,jspins/) + dims(:5)=(/two,nmat,neig,nkpts,jspins/) #endif - CALL h5screate_simple_f(5,dims(:5),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "ev", H5T_NATIVE_DOUBLE, spaceid,& - & d%evsetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - ! basis - dims(:4)=(/nmat,three,nkpts,jspins/) - CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) - CALL h5dcreate_f(d%fid, "k", H5T_NATIVE_INTEGER, spaceid,& - & d%ksetid, hdferr) - CALL h5sclose_f(spaceid,hdferr) - ELSE - CALL h5fopen_f (trim(d%fname)//'.hdf', access_Mode, d%fid,& - & hdferr,access_prp) - !get dataset-ids - CALL h5dopen_f(d%fid, 'el', d%esetid, hdferr) - CALL h5dopen_f(d%fid, 'evac', d%evacsetid, hdferr) - CALL h5dopen_f(d%fid, 'ello', d%ellosetid, hdferr) - CALL h5dopen_f(d%fid, 'bk', d%bksetid, hdferr) - CALL h5dopen_f(d%fid, 'wk', d%wksetid, hdferr) - CALL h5dopen_f(d%fid, 'energy', d%energysetid, hdferr) - CALL h5dopen_f(d%fid, 'k', d%ksetid, hdferr) - CALL h5dopen_f(d%fid, 'neig', d%neigsetid, hdferr) - CALL h5dopen_f(d%fid, 'ev', d%evsetid, hdferr) - CALL h5dopen_f(d%fid, 'nv', d%nvsetid, hdferr) - CALL h5dopen_f(d%fid, 'nmat', d%nmatsetid, hdferr) - ENDIF - IF (.NOT.access_prp==H5P_DEFAULT_f) CALL H5Pclose_f(access_prp& + CALL h5screate_simple_f(5,dims(:5),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "ev", H5T_NATIVE_DOUBLE, spaceid, d%evsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + ! basis + dims(:4)=(/nmat,three,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "k", H5T_NATIVE_INTEGER, spaceid, d%ksetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + !stuff for dos etc + IF (l_dos) THEN + dims(:5)=(/4,ntype,neig,nkpts,jspins/) + CALL h5screate_simple_f(5,dims(:5),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qal", H5T_NATIVE_DOUBLE, spaceid, d%qalsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:4)=(/neig,2,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qvac", H5T_NATIVE_DOUBLE, spaceid, d%qvacsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:3)=(/neig,nkpts,jspins/) + CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qis", H5T_NATIVE_DOUBLE, spaceid, d%qissetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:5)=(/neig,layers,2,nkpts,jspins/) + CALL h5screate_simple_f(5,dims(:5),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qvlay", H5T_NATIVE_DOUBLE, spaceid, d%qvlaysetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:7)=(/2,nstars,neig,layers,2,nkpts,jspins/) + CALL h5screate_simple_f(7,dims(:7),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qstars", H5T_NATIVE_DOUBLE, spaceid, d%qstarssetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:3)=(/neig,nkpts,jspins/) + CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "ksym", H5T_NATIVE_DOUBLE, spaceid, d%ksymsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:3)=(/neig,nkpts,jspins/) + CALL h5screate_simple_f(3,dims(:3),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "jsym", H5T_NATIVE_DOUBLE, spaceid, d%jsymsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + IF (l_mcd) THEN + dims(:5)=(/3*ntype,ncored,neig,nkpts,jspins/) + CALL h5screate_simple_f(5,dims(:5),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "mcd", H5T_NATIVE_DOUBLE, spaceid, d%mcdsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + ENDIF + IF (l_orb) THEN + dims(:4)=(/nsld,neig,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qintsl", H5T_NATIVE_DOUBLE, spaceid, d%qintslsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:4)=(/nsld,neig,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qmtsl", H5T_NATIVE_DOUBLE, spaceid, d%qmtslsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:4)=(/neig,natd,nkpts,jspins/) + CALL h5screate_simple_f(4,dims(:4),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "qmtp", H5T_NATIVE_DOUBLE, spaceid, d%qmtpsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + dims(:5)=(/neig,23,natd,nkpts,jspins/) + CALL h5screate_simple_f(5,dims(:5),spaceid,hdferr) + CALL h5dcreate_f(d%fid, "orbcomp", H5T_NATIVE_DOUBLE, spaceid, d%orbcompsetid, hdferr) + CALL h5sclose_f(spaceid,hdferr) + ENDIF + ENDIF + ELSE + CALL h5fopen_f (TRIM(d%fname)//'.hdf', access_Mode, d%fid, hdferr,access_prp) + !get dataset-ids + CALL h5dopen_f(d%fid, 'el', d%esetid, hdferr) + CALL h5dopen_f(d%fid, 'evac', d%evacsetid, hdferr) + CALL h5dopen_f(d%fid, 'ello', d%ellosetid, hdferr) + CALL h5dopen_f(d%fid, 'bk', d%bksetid, hdferr) + CALL h5dopen_f(d%fid, 'wk', d%wksetid, hdferr) + CALL h5dopen_f(d%fid, 'energy', d%energysetid, hdferr) + CALL h5dopen_f(d%fid, 'k', d%ksetid, hdferr) + CALL h5dopen_f(d%fid, 'neig', d%neigsetid, hdferr) + CALL h5dopen_f(d%fid, 'ev', d%evsetid, hdferr) + CALL h5dopen_f(d%fid, 'nv', d%nvsetid, hdferr) + CALL h5dopen_f(d%fid, 'nmat', d%nmatsetid, hdferr) + IF (l_dos) THEN + CALL h5dopen_f(d%fid, 'qal', d%qalsetid, hdferr) + CALL h5dopen_f(d%fid, 'qvac', d%qvacsetid, hdferr) + CALL h5dopen_f(d%fid, 'qis', d%qissetid, hdferr) + CALL h5dopen_f(d%fid, 'qvlay', d%qvlaysetid, hdferr) + CALL h5dopen_f(d%fid, 'qstars', d%qstarssetid, hdferr) + CALL h5dopen_f(d%fid, 'ksym', d%ksymsetid, hdferr) + CALL h5dopen_f(d%fid, 'jsym', d%jsymsetid, hdferr) + IF (l_mcd) THEN + CALL h5dopen_f(d%fid, 'mcd', d%mcdsetid, hdferr) + ENDIF + IF (l_orb) THEN + CALL h5dopen_f(d%fid, 'qintsl', d%qintslsetid, hdferr) + CALL h5dopen_f(d%fid, 'qmtsl', d%qmtslsetid, hdferr) + CALL h5dopen_f(d%fid, 'qmtp', d%qmtpsetid, hdferr) + CALL h5dopen_f(d%fid, 'orbcomp', d%orbcompsetid, hdferr) + ENDIF + ENDIF + IF (.NOT.access_prp==H5P_DEFAULT_f) CALL H5Pclose_f(access_prp& & ,hdferr) #else - CALL juDFT_error("Could not use HDF5 for IO, please recompile") + CALL juDFT_error("Could not use HDF5 for IO, please recompile") #endif - END SUBROUTINE open_eig - !---------------------------------------------------------------------- - SUBROUTINE close_eig(id,filename) - !***************************************************************** - ! closes hdf-file for eigenvectors+values - !***************************************************************** - IMPLICIT NONE - INTEGER,INTENT(IN) :: id - CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: filename - - INTEGER::hdferr - TYPE(t_data_HDF),pointer::d - - !close datasets + END SUBROUTINE open_eig + !---------------------------------------------------------------------- + SUBROUTINE close_eig(id,filename) + !***************************************************************** + ! closes hdf-file for eigenvectors+values + !***************************************************************** + IMPLICIT NONE + INTEGER,INTENT(IN) :: id + CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: filename + + INTEGER::hdferr + TYPE(t_data_HDF),POINTER::d + + !close datasets #ifdef CPP_HDF - call priv_find_data(id,d) - - CALL h5dclose_f(d%esetid,hdferr) - CALL h5dclose_f(d%evacsetid,hdferr) - CALL h5dclose_f(d%ellosetid,hdferr) - CALL h5dclose_f(d%bksetid,hdferr) - CALL h5dclose_f(d%wksetid,hdferr) - CALL h5dclose_f(d%energysetid,hdferr) - CALL h5dclose_f(d%ksetid,hdferr) - CALL h5dclose_f(d%neigsetid,hdferr) - CALL h5dclose_f(d%evsetid,hdferr) - CALL h5dclose_f(d%nvsetid,hdferr) - CALL h5dclose_f(d%nmatsetid,hdferr) - !close file - CALL h5fclose_f(d%fid,hdferr) - !If a filename was given and the name is not the current filename - IF (present(filename)) THEN - IF (filename.NE.d%fname) THEN - CALL system("mv "//trim(d%fname)//".hdf "//trim(filename)//".hdf") - ENDIF - ENDIF - d%fname="eig" - call eig66_remove_data(id) + CALL priv_find_data(id,d) + + CALL h5dclose_f(d%esetid,hdferr) + CALL h5dclose_f(d%evacsetid,hdferr) + CALL h5dclose_f(d%ellosetid,hdferr) + CALL h5dclose_f(d%bksetid,hdferr) + CALL h5dclose_f(d%wksetid,hdferr) + CALL h5dclose_f(d%energysetid,hdferr) + CALL h5dclose_f(d%ksetid,hdferr) + CALL h5dclose_f(d%neigsetid,hdferr) + CALL h5dclose_f(d%evsetid,hdferr) + CALL h5dclose_f(d%nvsetid,hdferr) + CALL h5dclose_f(d%nmatsetid,hdferr) + IF (l_dos) THEN + CALL h5dclose_f(d%qalsetid, hdferr) + CALL h5dclose_f(d%qvacsetid, hdferr) + CALL h5dclose_f(d%qissetid, hdferr) + CALL h5dclose_f(d%qvlaysetid, hdferr) + CALL h5dclose_f(d%qstarssetid, hdferr) + CALL h5dclose_f(d%ksymsetid, hdferr) + CALL h5dclose_f(d%jsymsetid, hdferr) + IF (l_mcd) THEN + CALL h5dclose_f(d%mcdsetid, hdferr) + ENDIF + IF (l_orb) THEN + CALL h5dclose_f(d%qintslsetid, hdferr) + CALL h5dclose_f(d%qmtslsetid, hdferr) + CALL h5dclose_f(d%qmtpsetid, hdferr) + CALL h5dclose_f(d%orbcompsetid, hdferr) + ENDIF + ENDIF + !close file + CALL h5fclose_f(d%fid,hdferr) + !If a filename was given and the name is not the current filename + IF (PRESENT(filename)) THEN + IF (filename.NE.d%fname) THEN + CALL system("mv "//TRIM(d%fname)//".hdf "//TRIM(filename)//".hdf") + ENDIF + ENDIF + d%fname="eig" + CALL eig66_remove_data(id) #endif - END SUBROUTINE + END SUBROUTINE close_eig #ifdef CPP_HDF - !---------------------------------------------------------------------- - SUBROUTINE priv_r_vec(d,nk,jspin,n_start,n_end,nmat,z) + !---------------------------------------------------------------------- + SUBROUTINE priv_r_vec(d,nk,jspin,n_start,n_end,nmat,z) - USE m_hdf_tools - IMPLICIT NONE - TYPE(t_data_HDF),INTENT(IN)::d - INTEGER, INTENT(IN) :: nk,jspin - INTEGER, INTENT(IN) :: n_start,n_end - INTEGER, INTENT(OUT) :: nmat - REAL, INTENT(OUT) :: z(:,:) + USE m_hdf_tools + IMPLICIT NONE + TYPE(t_data_HDF),INTENT(IN)::d + INTEGER, INTENT(IN) :: nk,jspin + INTEGER, INTENT(IN) :: n_start,n_end + INTEGER, INTENT(OUT) :: nmat + REAL, INTENT(OUT) :: z(:,:) - INTEGER i,j,neig_l + INTEGER i,j,neig_l - neig_l = n_end - n_start + 1 + neig_l = n_end - n_start + 1 - ! read matrix size - CALL io_read_integer0(d%nmatsetid,(/nk,jspin/),(/1,1/),nmat) + ! read matrix size + CALL io_read_integer0(d%nmatsetid,(/nk,jspin/),(/1,1/),nmat) - IF ( nmat > size(z,1) .OR. neig_l > size(z,2) ) THEN - WRITE (6,*) nmat,size(z,1),size(z,2) - CALL juDFT_error("eig66_hdf$read_vec",calledby ="eig66_hdf") - ENDIF + IF ( nmat > SIZE(z,1) .OR. neig_l > SIZE(z,2) ) THEN + WRITE (6,*) nmat,SIZE(z,1),SIZE(z,2) + CALL juDFT_error("eig66_hdf$read_vec",calledby ="eig66_hdf") + ENDIF - !read eigenvectors - CALL io_read_real2(d%evsetid,(/1,1,n_start,nk,jspin/),& + !read eigenvectors + CALL io_read_real2(d%evsetid,(/1,1,n_start,nk,jspin/),& & (/1,nmat,neig_l,1,1/),& & z(:nmat,:neig_l) ) - END SUBROUTINE priv_r_vec + END SUBROUTINE priv_r_vec #endif + SUBROUTINE read_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(OUT) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(OUT) :: qstars(:,:,:,:) + INTEGER,INTENT(OUT) :: ksym(:),jsym(:) + REAL,INTENT(OUT),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(OUT),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + TYPE(t_data_HDF),POINTER :: d + REAL,ALLOCATABLE :: r_tmp5(:,:,:,:,:) + CALL priv_find_data(id,d) +#ifdef CPP_HDF + CALL io_read_real3(d%qalsetid,(/1,1,1,nk,jspin/),(/SIZE(qal,1),SIZE(qal,2),SIZE(qal,3),1,1/),qal) + CALL io_read_real2(d%qvacsetid,(/1,1,nk,jspin/),(/SIZE(qvac,1),SIZE(qvac,2),1,1/),qvac) + CALL io_read_real1(d%qissetid,(/1,nk,jspin/),(/SIZE(qis,1),1,1/),qis) + CALL io_read_real3(d%qvlaysetid,(/1,1,1,nk,jspin/),(/SIZE(qvlay,1),SIZE(qvlay,2),SIZE(qvlay,3),1,1/),qvlay) + ALLOCATE(r_tmp5(2,SIZE(qstars,1),SIZE(qstars,2),SIZE(qstars,3),SIZE(qstars,4))) + CALL io_read_real5(d%qstarssetid,(/1,1,1,1,1,nk,jspin/),(/2,SIZE(qstars,1),SIZE(qstars,2),SIZE(qstars,3),SIZE(qstars,4),1,1/),r_tmp5(:,:,:,:,:)) + qstars=CMPLX(r_tmp5(1,:,:,:,:),r_tmp5(2,:,:,:,:)) + DEALLOCATE(r_tmp5) + CALL io_read_integer1(d%ksymsetid,(/1,nk,jspin/),(/SIZE(ksym,1),1,1/),ksym) + CALL io_read_integer1(d%jsymsetid,(/1,nk,jspin/),(/SIZE(jsym,1),1,1/),jsym) + IF (d%l_mcd.AND.PRESENT(mcd)) THEN + CALL io_read_real3(d%mcdsetid,(/1,1,1,nk,jspin/),(/SIZE(mcd,1),SIZE(mcd,2),SIZE(mcd,3),1,1/),mcd) + ENDIF + IF (d%l_orb.AND.PRESENT(qintsl)) THEN + CALL io_read_real2(d%qintslsetid,(/1,1,nk,jspin/),(/SIZE(qintsl,1),SIZE(qintsl,2),1,1/),qintsl) + CALL io_read_real2(d%qmtslsetid,(/1,1,nk,jspin/),(/SIZE(qmtsl,1),SIZE(qmtsl,2),1,1/),qmtsl) + CALL io_read_real2(d%qmtpsetid,(/1,1,nk,jspin/),(/SIZE(qmtp,1),SIZE(qmtp,2),1,1/),qmtp) + CALL io_read_real2(d%orbcompsetid,(/1,1,1,nk,jspin/),(/SIZE(orbcomp,1),23,SIZE(orbcomp,3),1,1/),orbcomp) + ENDIF +#endif + END SUBROUTINE read_dos + + + SUBROUTINE write_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(IN) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(IN) :: qstars(:,:,:,:) + INTEGER,INTENT(IN) :: ksym(:),jsym(:) + REAL,INTENT(IN),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(IN),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + TYPE(t_data_HDF),POINTER ::d + CALL priv_find_data(id,d) +#ifdef CPP_HDF + CALL io_write_real3(d%qalsetid,(/1,1,1,nk,jspin/),(/SIZE(qal,1),SIZE(qal,2),SIZE(qal,3),1,1/),qal) + CALL io_write_real2(d%qvacsetid,(/1,1,nk,jspin/),(/SIZE(qvac,1),SIZE(qvac,2),1,1/),qvac) + CALL io_write_real1(d%qissetid,(/1,nk,jspin/),(/SIZE(qis,1),1,1/),qis) + CALL io_write_real3(d%qvlaysetid,(/1,1,1,nk,jspin/),(/SIZE(qvlay,1),SIZE(qvlay,2),SIZE(qvlay,3),1,1/),qvlay) + CALL io_write_real5(d%qstarssetid,(/1,1,1,1,1,nk,jspin/),(/1,SIZE(qstars,1),SIZE(qstars,2),SIZE(qstars,3),SIZE(qstars,4),1,1/),REAL(qstars)) + CALL io_write_real5(d%qstarssetid,(/2,1,1,1,1,nk,jspin/),(/1,SIZE(qstars,1),SIZE(qstars,2),SIZE(qstars,3),SIZE(qstars,4),1,1/),AIMAG(qstars)) + + CALL io_write_integer1(d%ksymsetid,(/1,nk,jspin/),(/SIZE(ksym,1),1,1/),ksym) + CALL io_write_integer1(d%jsymsetid,(/1,nk,jspin/),(/SIZE(jsym,1),1,1/),jsym) + IF (d%l_mcd.AND.PRESENT(mcd)) THEN + CALL io_write_real3(d%mcdsetid,(/1,1,1,nk,jspin/),(/SIZE(mcd,1),SIZE(mcd,2),SIZE(mcd,3),1,1/),mcd) + ENDIF + IF (d%l_orb.AND.PRESENT(qintsl)) THEN + CALL io_write_real2(d%qintslsetid,(/1,1,nk,jspin/),(/SIZE(qintsl,1),SIZE(qintsl,2),1,1/),qintsl) + CALL io_write_real2(d%qmtslsetid,(/1,1,nk,jspin/),(/SIZE(qmtsl,1),SIZE(qmtsl,2),1,1/),qmtsl) + CALL io_write_real2(d%qmtpsetid,(/1,1,nk,jspin/),(/SIZE(qmtp,1),SIZE(qmtp,2),1,1/),qmtp) + CALL io_write_real2(d%orbcompsetid,(/1,1,1,nk,jspin/),(/SIZE(orbcomp,1),23,SIZE(orbcomp,3),1,1/),orbcomp) + ENDIF +#endif + END SUBROUTINE write_dos + - SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,& - & eig,el,ello,evac,& - & nlotot,kveclo,n_size,n_rank,z) + SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,& + & eig,el,ello,evac,& + & nlotot,kveclo,n_size,n_rank,z) - !***************************************************************** - ! writes all eignevecs for the nk-th kpoint - !***************************************************************** - IMPLICIT NONE + !***************************************************************** + ! writes all eignevecs for the nk-th kpoint + !***************************************************************** + IMPLICIT NONE - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank - REAL, INTENT(IN),OPTIONAL :: wk - INTEGER, INTENT(IN),OPTIONAL :: neig,nv,nmat,nlotot,neig_total - INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:) - REAL, INTENT(IN),OPTIONAL :: evac(2),ello(:,:) - CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank + REAL, INTENT(IN),OPTIONAL :: wk + INTEGER, INTENT(IN),OPTIONAL :: neig,nv,nmat,nlotot,neig_total + INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:) + REAL, INTENT(IN),OPTIONAL :: evac(2),ello(:,:) + CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) - INTEGER i,j,k,nv_local,n1,n2,ne - TYPE(t_data_HDF),pointer::d - call priv_find_data(id,d) + INTEGER i,j,k,nv_local,n1,n2,ne + TYPE(t_data_HDF),POINTER::d + CALL priv_find_data(id,d) #ifdef CPP_HDF - ! - !write enparas - ! - nv_local=huge(1) + ! + !write enparas + ! + nv_local=HUGE(1) - IF (present(el))& + IF (PRESENT(el))& & CALL io_write_real2(& & d%esetid,(/1,1,nk,jspin/),& - & (/size(el,1),size(el,2),1,1/),el) + & (/SIZE(el,1),SIZE(el,2),1,1/),el) - IF (present(ello))& + IF (PRESENT(ello))& & CALL io_write_real2(& & d%ellosetid,(/1,1,nk,jspin/),& - & (/size(ello,1),size(ello,2),1,1/),ello) + & (/SIZE(ello,1),SIZE(ello,2),1,1/),ello) - IF (present(evac)) CALL io_write_real1(& + IF (PRESENT(evac)) CALL io_write_real1(& & d%evacsetid,(/1,nk,jspin/),(/2,1,1/),evac) - ! - !write kpts - ! + ! + !write kpts + ! - IF (present(bk)) CALL io_write_real1(& + IF (PRESENT(bk)) CALL io_write_real1(& & d%bksetid,(/1,nk,jspin/),(/3,1,1/),bk) - IF (present(wk)) CALL io_write_real0(& + IF (PRESENT(wk)) CALL io_write_real0(& & d%wksetid,(/nk,jspin/),(/1,1/),wk) - ! - !write basis - ! + ! + !write basis + ! - IF (present(nv)) THEN - nv_local=nv - CALL io_write_integer0(d%nvsetid,(/nk,jspin/),(/1,1/),nv) - ENDIF + IF (PRESENT(nv)) THEN + nv_local=nv + CALL io_write_integer0(d%nvsetid,(/nk,jspin/),(/1,1/),nv) + ENDIF - IF (present(nmat)) CALL io_write_integer0(& + IF (PRESENT(nmat)) CALL io_write_integer0(& & d%nmatsetid,(/nk,jspin/),(/1,1/),nmat) - IF (present(k1)) CALL io_write_integer1(& + IF (PRESENT(k1)) CALL io_write_integer1(& & d%ksetid,(/1,1,nk,jspin/),& - & (/min(nv_local,size(k1)),1,1,1/),k1(:min(nv_local,size(k1)))) + & (/MIN(nv_local,SIZE(k1)),1,1,1/),k1(:MIN(nv_local,SIZE(k1)))) - IF (present(k2)) CALL io_write_integer1(& + IF (PRESENT(k2)) CALL io_write_integer1(& & d%ksetid,(/1,2,nk,jspin/),& - & (/min(nv_local,size(k2)),1,1,1/),k2(:min(nv_local,size(k2)))) + & (/MIN(nv_local,SIZE(k2)),1,1,1/),k2(:MIN(nv_local,SIZE(k2)))) - IF (present(k3)) CALL io_write_integer1(& + IF (PRESENT(k3)) CALL io_write_integer1(& & d%ksetid,(/1,3,nk,jspin/),& - & (/min(nv_local,size(k3)),1,1,1/),k3(:min(nv_local,size(k3)))) - - IF (present(kveclo).AND.present(nlotot).AND.& - & present(k1).AND.present(k2).AND.present(k3)) THEN - - WRITE(*,*) kveclo,nlotot - DO k = 1, nlotot - CALL io_write_integer0(& - & d%ksetid,(/nv+k,1,nk,jspin/),(/1,1,1,1/),k1(kveclo(k))) - CALL io_write_integer0(& - & d%ksetid,(/nv+k,2,nk,jspin/),(/1,1,1,1/),k2(kveclo(k))) - CALL io_write_integer0(& - & d%ksetid,(/nv+k,3,nk,jspin/),(/1,1,1,1/),k3(kveclo(k))) - ENDDO - ENDIF - ! - !write eigenvalues - ! - - IF (present(neig_total)) THEN - CALL io_write_integer0(d%neigsetid,(/nk,jspin/),(/1,1/),neig_total) - ENDIF - - IF (present(n_rank).AND.present(n_size).AND.& - & present(eig).AND.present(neig)) THEN - CALL io_write_real1s(& - & d%energysetid,(/n_rank+1,nk,jspin/), & - & (/neig,1,1/),eig(:neig),(/n_size,1,1/)) - !write eigenvectors - ! - ELSEIF (present(eig).AND.present(neig)) THEN - CALL io_write_real1s(& - & d%energysetid,(/1,nk,jspin/),& - & (/neig,1,1/),eig(:neig),(/1,1,1/)) - ELSE - IF (present(eig)) CALL juDFT_error("BUG in calling write_eig") - ENDIF - IF (present(z).AND..NOT.present(neig))& + & (/MIN(nv_local,SIZE(k3)),1,1,1/),k3(:MIN(nv_local,SIZE(k3)))) + + IF (PRESENT(kveclo).AND.PRESENT(nlotot).AND.& + & PRESENT(k1).AND.PRESENT(k2).AND.PRESENT(k3)) THEN + + WRITE(*,*) kveclo,nlotot + DO k = 1, nlotot + CALL io_write_integer0(& + & d%ksetid,(/nv+k,1,nk,jspin/),(/1,1,1,1/),k1(kveclo(k))) + CALL io_write_integer0(& + & d%ksetid,(/nv+k,2,nk,jspin/),(/1,1,1,1/),k2(kveclo(k))) + CALL io_write_integer0(& + & d%ksetid,(/nv+k,3,nk,jspin/),(/1,1,1,1/),k3(kveclo(k))) + ENDDO + ENDIF + ! + !write eigenvalues + ! + + IF (PRESENT(neig_total)) THEN + CALL io_write_integer0(d%neigsetid,(/nk,jspin/),(/1,1/),neig_total) + ENDIF + + IF (PRESENT(n_rank).AND.PRESENT(n_size).AND.& + & PRESENT(eig).AND.PRESENT(neig)) THEN + CALL io_write_real1s(& + & d%energysetid,(/n_rank+1,nk,jspin/), & + & (/neig,1,1/),eig(:neig),(/n_size,1,1/)) + !write eigenvectors + ! + ELSEIF (PRESENT(eig).AND.PRESENT(neig)) THEN + CALL io_write_real1s(& + & d%energysetid,(/1,nk,jspin/),& + & (/neig,1,1/),eig(:neig),(/1,1,1/)) + ELSE + IF (PRESENT(eig)) CALL juDFT_error("BUG in calling write_eig") + ENDIF + IF (PRESENT(z).AND..NOT.PRESENT(neig))& & CALL juDFT_error("BUG in calling write_eig with eigenvector") - n1=1;n2=0 - IF (present(n_size)) n1=n_size - IF (present(n_rank)) n2=n_rank - IF (present(z)) THEN - - SELECT TYPE(z) - TYPE IS (real) - CALL io_write_real2s(& - & d%evsetid,(/1,1,n2+1,nk,jspin/),& - & (/1,nmat,neig,1,1/),real(z(:nmat,:neig)),(/1,1,n1,1,1/)) - TYPE IS (complex) - CALL io_write_real2s(& - & d%evsetid,(/1,1,n2+1,nk,jspin/),& - & (/1,nmat,neig,1,1/),real(z(:nmat,:neig)),(/1,1,n1,1,1/)) - CALL io_write_real2s(& - & d%evsetid,(/2,1,n2+1,nk,jspin/),& - & (/1,nmat,neig,1,1/),aimag(z(:nmat,:neig)),& - & (/1,1,n1,1,1/)) - END SELECT - ENDIF + n1=1;n2=0 + IF (PRESENT(n_size)) n1=n_size + IF (PRESENT(n_rank)) n2=n_rank + IF (PRESENT(z)) THEN + + SELECT TYPE(z) + TYPE IS (REAL) + CALL io_write_real2s(& + & d%evsetid,(/1,1,n2+1,nk,jspin/),& + & (/1,nmat,neig,1,1/),REAL(z(:nmat,:neig)),(/1,1,n1,1,1/)) + TYPE IS (COMPLEX) + CALL io_write_real2s(& + & d%evsetid,(/1,1,n2+1,nk,jspin/),& + & (/1,nmat,neig,1,1/),REAL(z(:nmat,:neig)),(/1,1,n1,1,1/)) + CALL io_write_real2s(& + & d%evsetid,(/2,1,n2+1,nk,jspin/),& + & (/1,nmat,neig,1,1/),AIMAG(z(:nmat,:neig)),& + & (/1,1,n1,1,1/)) + END SELECT + ENDIF #endif - END SUBROUTINE write_eig + END SUBROUTINE write_eig #ifdef CPP_HDF - !---------------------------------------------------------------------- - SUBROUTINE priv_r_vecc(& - & d,nk,jspin,n_start,n_end,& - & nmat,z) + !---------------------------------------------------------------------- + SUBROUTINE priv_r_vecc(& + & d,nk,jspin,n_start,n_end,& + & nmat,z) - USE m_hdf_tools - IMPLICIT NONE - TYPE(t_data_HDF),INTENT(IN)::d - INTEGER, INTENT(IN) :: nk,jspin - INTEGER, INTENT(IN) :: n_start,n_end - INTEGER, INTENT(OUT) :: nmat - COMPLEX, INTENT(OUT) :: z(:,:) + USE m_hdf_tools + IMPLICIT NONE + TYPE(t_data_HDF),INTENT(IN)::d + INTEGER, INTENT(IN) :: nk,jspin + INTEGER, INTENT(IN) :: n_start,n_end + INTEGER, INTENT(OUT) :: nmat + COMPLEX, INTENT(OUT) :: z(:,:) - REAL, ALLOCATABLE :: z1(:,:,:) - INTEGER i,j,neig_l + REAL, ALLOCATABLE :: z1(:,:,:) + INTEGER i,j,neig_l - neig_l = n_end - n_start + 1 + neig_l = n_end - n_start + 1 - ! read matrix size - CALL io_read_integer0(& + ! read matrix size + CALL io_read_integer0(& & d%nmatsetid,(/nk,jspin/),(/1,1/),& & nmat) - IF ( nmat > size(z,1) .OR. neig_l > size(z,2) ) THEN - WRITE (6,*) nmat,size(z,1),size(z,2) - CALL juDFT_error("eig66_hdf$read_vec",calledby ="eig66_hdf") - ENDIF + IF ( nmat > SIZE(z,1) .OR. neig_l > SIZE(z,2) ) THEN + WRITE (6,*) nmat,SIZE(z,1),SIZE(z,2) + CALL juDFT_error("eig66_hdf$read_vec",calledby ="eig66_hdf") + ENDIF - ! read eigenvectors - ALLOCATE (z1(2,nmat,neig_l)) - CALL io_read_real3(d%evsetid,(/1,1,n_start,nk,jspin/),& + ! read eigenvectors + ALLOCATE (z1(2,nmat,neig_l)) + CALL io_read_real3(d%evsetid,(/1,1,n_start,nk,jspin/),& & (/2,nmat,neig_l,1,1/),z1) - - DO i=1,neig_l - DO j=1,nmat - z(j,i) = cmplx( z1(1,j,i) ,z1(2,j,i) ) - ENDDO - ENDDO - DEALLOCATE (z1) + DO i=1,neig_l + DO j=1,nmat + z(j,i) = CMPLX( z1(1,j,i) ,z1(2,j,i) ) + ENDDO + ENDDO - END SUBROUTINE priv_r_vecc - !----------------------------------------------------------------------- + DEALLOCATE (z1) + + END SUBROUTINE priv_r_vecc + !----------------------------------------------------------------------- #endif - SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,& - & ello,evac,kveclo,n_start,n_end,z) - IMPLICIT NONE - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat - INTEGER, INTENT(OUT),OPTIONAL :: neig - REAL, INTENT(OUT),OPTIONAL :: eig(:) - INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) - REAL, INTENT(OUT),OPTIONAL :: bk(:),wk - INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end - CLASS(*),TARGET,INTENT(OUT),OPTIONAL :: z(:,:) + SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,& + & ello,evac,kveclo,n_start,n_end,z) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat + INTEGER, INTENT(OUT),OPTIONAL :: neig + REAL, INTENT(OUT),OPTIONAL :: eig(:) + INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) + REAL, INTENT(OUT),OPTIONAL :: bk(:),wk + INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end + CLASS(*),TARGET,INTENT(OUT),OPTIONAL :: z(:,:) #ifdef CPP_HDF - INTEGER:: n1,n,k,k1_t,k2_t,k3_t - TYPE(t_data_HDF),pointer::d - call priv_find_data(id,d) - - - IF (present(neig)) THEN - CALL io_read_integer0(d%neigsetid,(/nk,jspin/),(/1,1/),neig) - - IF ( present(eig) ) THEN ! read eigenv - IF ( neig > size(eig) ) THEN - WRITE(*,*) neig,size(eig) - CALL juDFT_error("eig66_hdf$readeig",calledby ="eig66_hdf") + INTEGER:: n1,n,k,k1_t,k2_t,k3_t + TYPE(t_data_HDF),POINTER::d + CALL priv_find_data(id,d) + + + IF (PRESENT(neig)) THEN + CALL io_read_integer0(d%neigsetid,(/nk,jspin/),(/1,1/),neig) + + IF ( PRESENT(eig) ) THEN ! read eigenv + IF ( neig > SIZE(eig) ) THEN + WRITE(*,*) neig,SIZE(eig) + CALL juDFT_error("eig66_hdf$readeig",calledby ="eig66_hdf") + ENDIF + CALL io_read_real1(d%energysetid,(/1,nk,jspin/),(/neig,1,1/),& + & eig(:neig)) + ENDIF + ENDIF + + IF (PRESENT(k1)) THEN + IF (.NOT.(PRESENT(k2).AND.PRESENT(k3).AND.PRESENT(kveclo)))& + & CALL juDFT_error("BUG1 in calling read_eig") + + CALL io_read_integer0(d%nvsetid,(/nk,jspin/),(/1,1/),n1) + IF (n1>SIZE(k1)) CALL juDFT_error("eig66_hdf$read_basis",& + & calledby="eig66_hdf") + !read basis + CALL io_read_integer1(d%ksetid,(/1,1,nk,jspin/),(/n1,1,1,1/),k1(:n1)) + CALL io_read_integer1(d%ksetid,(/1,2,nk,jspin/),(/n1,1,1,1/),k2(:n1)) + CALL io_read_integer1(d%ksetid,(/1,3,nk,jspin/),(/n1,1,1,1/),k3(:n1)) + DO k = 1, SIZE(kveclo) + CALL io_read_integer0(d%ksetid,(/n1+k,1,nk,jspin/),(/1,1,1,1/),k1_t) + CALL io_read_integer0(d%ksetid,(/n1+k,2,nk,jspin/),(/1,1,1,1/),k2_t) + CALL io_read_integer0(d%ksetid,(/n1+k,3,nk,jspin/),(/1,1,1,1/),k3_t) + DO n = 1, n1 + IF (( (k1_t == k1(n)).AND.(k2_t == k2(n)) ).AND.(k3_t == k3(n)) ) THEN + kveclo(k) = n + CYCLE ENDIF - CALL io_read_real1(d%energysetid,(/1,nk,jspin/),(/neig,1,1/),& - & eig(:neig)) - ENDIF - ENDIF - - IF (present(k1)) THEN - IF (.NOT.(present(k2).AND.present(k3).AND.present(kveclo)))& - & CALL juDFT_error("BUG1 in calling read_eig") - - CALL io_read_integer0(d%nvsetid,(/nk,jspin/),(/1,1/),n1) - IF (n1>size(k1)) CALL juDFT_error("eig66_hdf$read_basis",& - & calledby="eig66_hdf") - !read basis - CALL io_read_integer1(d%ksetid,(/1,1,nk,jspin/),(/n1,1,1,1/),k1(:n1)) - CALL io_read_integer1(d%ksetid,(/1,2,nk,jspin/),(/n1,1,1,1/),k2(:n1)) - CALL io_read_integer1(d%ksetid,(/1,3,nk,jspin/),(/n1,1,1,1/),k3(:n1)) - DO k = 1, size(kveclo) - CALL io_read_integer0(d%ksetid,(/n1+k,1,nk,jspin/),(/1,1,1,1/),k1_t) - CALL io_read_integer0(d%ksetid,(/n1+k,2,nk,jspin/),(/1,1,1,1/),k2_t) - CALL io_read_integer0(d%ksetid,(/n1+k,3,nk,jspin/),(/1,1,1,1/),k3_t) - DO n = 1, n1 - IF (( (k1_t == k1(n)).AND.(k2_t == k2(n)) ).AND.(k3_t == k3(n)) ) THEN - kveclo(k) = n - CYCLE - ENDIF - ENDDO - ENDDO - IF (present(nv)) nv=n1 - ELSE - IF (present(nv)) CALL io_read_integer0(d%nvsetid,(/nk,jspin/),(/1,1/),nv) - - ENDIF - IF (present(nmat)) & + ENDDO + ENDDO + IF (PRESENT(nv)) nv=n1 + ELSE + IF (PRESENT(nv)) CALL io_read_integer0(d%nvsetid,(/nk,jspin/),(/1,1/),nv) + + ENDIF + IF (PRESENT(nmat)) & & CALL io_read_integer0(d%nmatsetid,(/nk,jspin/),(/1,1/),nmat) - IF (present(el)) CALL io_read_real2(d%esetid,(/1,1,nk,jspin/),& - & (/size(el,1),size(el,2),1,1/),el(:,:)) - IF (present(ello)) CALL io_read_real2(d%ellosetid,(/1,1,nk,jspin/),& - & (/size(ello,1),size(ello,2),1,1/),ello(:,:)) - IF (present(evac)) CALL io_read_real1(d%evacsetid,(/1,nk,jspin/),& + IF (PRESENT(el)) CALL io_read_real2(d%esetid,(/1,1,nk,jspin/),& + & (/SIZE(el,1),SIZE(el,2),1,1/),el(:,:)) + IF (PRESENT(ello)) CALL io_read_real2(d%ellosetid,(/1,1,nk,jspin/),& + & (/SIZE(ello,1),SIZE(ello,2),1,1/),ello(:,:)) + IF (PRESENT(evac)) CALL io_read_real1(d%evacsetid,(/1,nk,jspin/),& & (/2,1,1/),evac) - IF (present(bk)) CALL& - io_read_real1(d%bksetid,(/1,nk,jspin/),(/3,1,1/),bk) - IF (present(wk)) CALL& - io_read_real0(d%wksetid,(/nk,jspin/),(/1,1/),wk) - - IF (present(n_start)) THEN - IF (.NOT.present(n_end)) CALL juDFT_error("BUG3 in read_eig") - IF (present(z)) THEN - SELECT TYPE(z) - TYPE IS (real) - CALL priv_r_vec(d,nk,jspin,n_start,n_end,n1,z) - TYPE is (complex) - CALL priv_r_vecc(d,nk,jspin,n_start,n_end,n1,z) - END SELECT - ENDIF - IF (present(nmat)) nmat=n1 - ENDIF + IF (PRESENT(bk)) CALL& + io_read_real1(d%bksetid,(/1,nk,jspin/),(/3,1,1/),bk) + IF (PRESENT(wk)) CALL& + io_read_real0(d%wksetid,(/nk,jspin/),(/1,1/),wk) + + IF (PRESENT(n_start)) THEN + IF (.NOT.PRESENT(n_end)) CALL juDFT_error("BUG3 in read_eig") + IF (PRESENT(z)) THEN + SELECT TYPE(z) + TYPE IS (REAL) + CALL priv_r_vec(d,nk,jspin,n_start,n_end,n1,z) + TYPE is (COMPLEX) + CALL priv_r_vecc(d,nk,jspin,n_start,n_end,n1,z) + END SELECT + ENDIF + IF (PRESENT(nmat)) nmat=n1 + ENDIF #endif - END SUBROUTINE + END SUBROUTINE read_eig -END MODULE + END MODULE diff --git a/io/eig66_io.F90 b/io/eig66_io.F90 index a45b0ba9c3ec6b3c1301bb7e39f2a16f760a5530..3c3cda56138fe97769b7ecc9c0ba5fb2b5c53735 100644 --- a/io/eig66_io.F90 +++ b/io/eig66_io.F90 @@ -1,14 +1,15 @@ MODULE m_eig66_io #include "juDFT_env.h" - use m_eig66_data - IMPLICIT NONE - PRIVATE + USE m_eig66_data + IMPLICIT NONE + PRIVATE - PUBLIC open_eig,close_eig - PUBLIC read_eig, write_eig - CONTAINS + PUBLIC open_eig,close_eig + PUBLIC read_eig, write_eig + PUBLIC read_dos,write_dos +CONTAINS - FUNCTION open_eig(mpi_comm,nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,l_noco,l_create,l_readonly,n_size,mode_in,filename)result(id) + FUNCTION open_eig(mpi_comm,nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,l_noco,l_create,l_readonly,n_size,mode_in,filename,layers,nstars,ncored,nsld,nat,l_dos,l_mcd,l_orb)RESULT(id) USE m_eig66_hdf,ONLY:open_eig_hdf=>open_eig USE m_eig66_DA ,ONLY:open_eig_DA=>open_eig USE m_eig66_mem,ONLY:open_eig_mem=>open_eig @@ -17,7 +18,9 @@ MODULE m_eig66_io INTEGER,INTENT(IN) :: nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,mpi_comm LOGICAL,INTENT(IN) :: l_noco,l_readonly,l_create INTEGER,INTENT(IN),OPTIONAL :: n_size,mode_in + LOGICAL,INTENT(IN),OPTIONAL :: l_dos,l_mcd,l_orb CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename + INTEGER,INTENT(IN),OPTIONAL :: layers,nstars,ncored,nsld,nat INTEGER:: id,mode INTEGER:: neig_local,isize,err @@ -27,7 +30,7 @@ MODULE m_eig66_io neig_local=neig #endif mode=-1 - if (present(mode_in)) mode=mode_in + IF (PRESENT(mode_in)) mode=mode_in IF (mode<0) THEN !Use default mode @@ -37,44 +40,44 @@ MODULE m_eig66_io mode=MEM_mode #endif !check if default was given on command-line - if (juDFT_was_argument("-mpi")) mode=MPI_mode - if (juDFT_was_argument("-mem")) mode=MEM_mode - if (juDFT_was_argument("-da")) mode=DA_mode - if (juDFT_was_argument("-hdf")) mode=HDF_mode + IF (juDFT_was_argument("-mpi")) mode=MPI_mode + IF (juDFT_was_argument("-mem")) mode=MEM_mode + IF (juDFT_was_argument("-da")) mode=DA_mode + IF (juDFT_was_argument("-hdf")) mode=HDF_mode ENDIF !Check if mode is available #ifndef CPP_MPI - if (mode==MPI_mode) call juDFT_error("MPI-mode not available. Recompile with CPP_MPI",calledby="eig66_io") + IF (mode==MPI_mode) CALL juDFT_error("MPI-mode not available. Recompile with CPP_MPI",calledby="eig66_io") #else - CALL MPI_COMM_SIZE(mpi_comm,isize,err) - if (isize>1.and.((mode==DA_mode.or.mode==mem_mode))) & - call juDFT_error("In a parallel calculation MEM/DA-mode are not available",calledby="eig66_io") + CALL MPI_COMM_SIZE(mpi_comm,isize,err) + IF (isize>1.AND.((mode==DA_mode.OR.mode==mem_mode))) & + CALL juDFT_error("In a parallel calculation MEM/DA-mode are not available",calledby="eig66_io") #endif #ifndef CPP_HDF - if (mode==HDF_mode) call juDFT_error("HDF-mode not available. Recompile with CPP_HDF",calledby="eig66_io") + IF (mode==HDF_mode) CALL juDFT_error("HDF-mode not available. Recompile with CPP_HDF",calledby="eig66_io") #endif - id=eig66_data_newid(mode) + id=eig66_data_newid(mode) - print *,"open_eig:",id,mode + PRINT *,"open_eig:",id,mode - call timestart("Open file/memory for IO of eig66") + CALL timestart("Open file/memory for IO of eig66") SELECT CASE (eig66_data_mode(id)) - CASE (DA_mode) - CALL open_eig_DA(id,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,nlotot,l_create,filename) - CASE (hdf_mode) - CALL open_eig_HDF(id,mpi_comm,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,l_create,l_readonly,filename) - CASE (mem_mode) - CALL open_eig_MEM(id,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,l_create,nlotot,l_noco,filename) - CASE (mpi_mode) - CALL open_eig_MPI(id,mpi_comm,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,l_create,nlotot,l_noco,n_size,filename) - CASE DEFAULT - CALL juDFT_error("Invalid IO-mode in eig66_io") + CASE (DA_mode) + CALL open_eig_DA(id,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,nlotot,l_create,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + CASE (hdf_mode) + CALL open_eig_HDF(id,mpi_comm,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,l_create,l_readonly,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + CASE (mem_mode) + CALL open_eig_MEM(id,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,l_create,nlotot,l_noco,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + CASE (mpi_mode) + CALL open_eig_MPI(id,mpi_comm,nmat,neig_local,nkpts,jspins,lmax,nlo,ntype,l_create,nlotot,l_noco,n_size,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + CASE DEFAULT + CALL juDFT_error("Invalid IO-mode in eig66_io") END SELECT - call timestop("Open file/memory for IO of eig66") - END FUNCTION + CALL timestop("Open file/memory for IO of eig66") + END FUNCTION open_eig - SUBROUTINE close_eig(id,filename) + SUBROUTINE close_eig(id,filename) USE m_eig66_hdf,ONLY:close_eig_hdf=>close_eig USE m_eig66_DA ,ONLY:close_eig_DA=>close_eig USE m_eig66_mem,ONLY:close_eig_MEM=>close_eig @@ -84,23 +87,23 @@ MODULE m_eig66_io CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename INTEGER :: mode mode=eig66_data_mode(id) - print*,"close_eig:",id,mode + PRINT*,"close_eig:",id,mode SELECT CASE (mode) - CASE (DA_mode) - CALL close_eig_DA(id,filename) - CASE (hdf_mode) - CALL close_eig_HDF(id,filename) - CASE (mem_mode) - CALL close_eig_Mem(id,filename=filename) - CASE (MPI_mode) - CALL close_eig_MPI(id,filename=filename) - CASE (-1) - CALL juDFT_error("ID not assigned in close_eig",calledby="eig66_io") + CASE (DA_mode) + CALL close_eig_DA(id,filename) + CASE (hdf_mode) + CALL close_eig_HDF(id,filename) + CASE (mem_mode) + CALL close_eig_Mem(id,filename=filename) + CASE (MPI_mode) + CALL close_eig_MPI(id,filename=filename) + CASE (-1) + CALL juDFT_error("ID not assigned in close_eig",calledby="eig66_io") END SELECT - END SUBROUTINE + END SUBROUTINE close_eig - SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) + SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) USE m_eig66_hdf,ONLY:read_eig_hdf=>read_eig USE m_eig66_DA ,ONLY:read_eig_DA=>read_eig USE m_eig66_mem,ONLY:read_eig_mem=>read_eig @@ -118,21 +121,21 @@ MODULE m_eig66_io INTEGER::n CALL timestart("IO (read)") SELECT CASE (eig66_data_mode(id)) - CASE (DA_mode) - CALL read_eig_DA(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) - CASE (hdf_mode) - CALL read_eig_hdf(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) - CASE (mem_mode) - CALL read_eig_mem(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) - CASE (mpi_mode) - CALL read_eig_mpi(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) - CASE (-1) - CALL juDFT_error("Could not read eig-file before opening") + CASE (DA_mode) + CALL read_eig_DA(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) + CASE (hdf_mode) + CALL read_eig_hdf(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) + CASE (mem_mode) + CALL read_eig_mem(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) + CASE (mpi_mode) + CALL read_eig_mpi(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,ello,evac,kveclo,n_start,n_end,z) + CASE (-1) + CALL juDFT_error("Could not read eig-file before opening") END SELECT - call timestop("IO (read)") - END SUBROUTINE + CALL timestop("IO (read)") + END SUBROUTINE read_eig - SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) + SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) USE m_eig66_hdf,ONLY:write_eig_hdf=>write_eig USE m_eig66_DA ,ONLY:write_eig_DA=>write_eig USE m_eig66_mem,ONLY:write_eig_MEM=>write_eig @@ -144,20 +147,77 @@ MODULE m_eig66_io INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:),evac(2),ello(:,:) CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) - call timestart("IO (write)") + CALL timestart("IO (write)") + SELECT CASE (eig66_data_mode(id)) + CASE (da_mode) + CALL write_eig_DA(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) + CASE (hdf_mode) + CALL write_eig_HDF(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) + CASE (mem_mode) + CALL write_eig_Mem(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) + CASE (MPI_mode) + CALL write_eig_MPI(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) + CASE (-1) + CALL juDFT_error("Could not write eig-file before opening") + END SELECT + CALL timestop("IO (write)") + END SUBROUTINE write_eig + + SUBROUTINE write_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + USE m_eig66_hdf,ONLY:write_dos_hdf=>write_dos + USE m_eig66_DA ,ONLY:write_dos_DA=>write_dos + USE m_eig66_mem,ONLY:write_dos_MEM=>write_dos + USE m_eig66_MPI,ONLY:write_dos_MPI=>write_dos + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(IN) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(IN) :: qstars(:,:,:,:) + INTEGER,INTENT(IN) :: ksym(:),jsym(:) + REAL,INTENT(IN),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(IN),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + CALL timestart("IO (dos-write)") + SELECT CASE (eig66_data_mode(id)) + CASE (da_mode) + CALL write_dos_DA(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (hdf_mode) + CALL write_dos_HDF(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (mem_mode) + CALL write_dos_Mem(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (MPI_mode) + CALL write_dos_MPI(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (-1) + CALL juDFT_error("Could not write eig-file before opening") + END SELECT + CALL timestop("IO (dos-write)") + END SUBROUTINE write_dos + + + SUBROUTINE read_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + USE m_eig66_hdf,ONLY:read_dos_hdf=>read_dos + USE m_eig66_DA ,ONLY:read_dos_DA=>read_dos + USE m_eig66_mem,ONLY:read_dos_MEM=>read_dos + USE m_eig66_MPI,ONLY:read_dos_MPI=>read_dos + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(OUT) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(OUT) :: qstars(:,:,:,:) + INTEGER,INTENT(OUT) :: ksym(:),jsym(:) + REAL,INTENT(OUT),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(OUT),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + CALL timestart("IO (dos-read)") SELECT CASE (eig66_data_mode(id)) - CASE (da_mode) - CALL write_eig_DA(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) - CASE (hdf_mode) - CALL write_eig_HDF(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) - CASE (mem_mode) - CALL write_eig_Mem(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) - CASE (MPI_mode) - CALL write_eig_MPI(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk,eig,el,ello,evac,nlotot,kveclo,n_start,n_end,z) - CASE (-1) - CALL juDFT_error("Could not write eig-file before opening") + CASE (da_mode) + CALL read_dos_DA(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (hdf_mode) + CALL read_dos_HDF(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (mem_mode) + CALL read_dos_Mem(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (MPI_mode) + CALL read_dos_MPI(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + CASE (-1) + CALL juDFT_error("Could not read eig-file before opening") END SELECT - call timestop("IO (write)") - END SUBROUTINE + CALL timestop("IO (dos-read)") + END SUBROUTINE read_dos END MODULE m_eig66_io diff --git a/io/eig66_mem.F90 b/io/eig66_mem.F90 index 75dccb160914112d4818a90c0cf5c59ceb2617db..bf4eaf35a98db6de89c87db5b8cc7d716c9eceb5 100644 --- a/io/eig66_mem.F90 +++ b/io/eig66_mem.F90 @@ -1,56 +1,57 @@ -module m_eig66_mem +MODULE m_eig66_mem #include "juDFT_env.h" -! Do the IO of the eig-file into memory -! The eig-file is split into four arrays: -! eig_int contains the basis-set information/integers (nv,nmat,ne,k1,k2,k3,kveclo) -! eig_real contains the basis-set information/real (el,evac,ello,bkpt,wtkpt) -! eig_eig contains the eigenvalues -! eig_vec contains the eigenvectors -! The record number is given by nrec=nk+(jspin-1)*nkpts - use m_eig66_data - implicit none - contains - - subroutine priv_find_data(id,d) - INTEGER,INTENT(IN)::id - TYPE(t_data_mem),pointer,intent(out):: d - - class(t_data),pointer ::dp - call eig66_find_data(dp,id) - select type(dp) - type is (t_data_mem) - d=>dp - class default - call judft_error("BUG: wrong datatype in eig66_mem") - END SELECT - END subroutine - - subroutine open_eig(id,nmat,neig,nkpts,jspins,lmax,nlo,ntype,l_create,nlotot,l_noco,filename) + ! Do the IO of the eig-file into memory + ! The eig-file is split into four arrays: + ! eig_int contains the basis-set information/integers (nv,nmat,ne,k1,k2,k3,kveclo) + ! eig_real contains the basis-set information/real (el,evac,ello,bkpt,wtkpt) + ! eig_eig contains the eigenvalues + ! eig_vec contains the eigenvectors + ! The record number is given by nrec=nk+(jspin-1)*nkpts + USE m_eig66_data + IMPLICIT NONE +CONTAINS + + SUBROUTINE priv_find_data(id,d) + INTEGER,INTENT(IN)::id + TYPE(t_data_mem),POINTER,INTENT(out):: d + + CLASS(t_data),POINTER ::dp + CALL eig66_find_data(dp,id) + SELECT TYPE(dp) + TYPE is (t_data_mem) + d=>dp + CLASS default + CALL judft_error("BUG: wrong datatype in eig66_mem") + END SELECT + END SUBROUTINE priv_find_data + + SUBROUTINE open_eig(id,nmat,neig,nkpts,jspins,lmax,nlo,ntype,l_create,nlotot,l_noco,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) INTEGER, INTENT(IN) :: id,nmat,neig,nkpts,jspins,nlo,ntype,lmax,nlotot - LOGICAL, INTENT(IN) :: l_noco,l_create + LOGICAL, INTENT(IN) :: l_noco,l_create,l_dos,l_mcd,l_orb CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename + INTEGER,INTENT(IN),OPTIONAL :: layers,nstars,ncored,nsld,nat + !locals + INTEGER:: length + TYPE(t_data_mem),POINTER:: d + CALL priv_find_data(id,d) - integer:: length - TYPE(t_data_mem),pointer:: d - call priv_find_data(id,d) + IF (ALLOCATED(d%eig_int)) THEN + IF (.NOT.l_create) THEN + IF (PRESENT(filename)) CALL priv_readfromfile() + RETURN + ENDIF + CALL close_eig(id,.TRUE.) - if (allocated(d%eig_int)) then - if (.not.l_create) THEN - if (present(filename)) call priv_readfromfile() - return - endif - call close_eig(id,.true.) + ENDIF - endif - - call eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype) + CALL eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype,l_dos,l_mcd,l_orb) !d%eig_int length=3 !nv+nmat+ne length=length+nmat*3 !k1,k2,k3 length=length+nlotot !kveclo - allocate(d%eig_int(length,jspins*nkpts)) + ALLOCATE(d%eig_int(length,jspins*nkpts)) !d%eig_real length=3+1+2 !bk,wk,evac @@ -59,215 +60,297 @@ module m_eig66_mem ALLOCATE(d%eig_real(length,jspins*nkpts)) !d%eig_eig length=jspins - if (l_noco) length=1 - allocate(d%eig_eig(neig,jspins*nkpts)) + IF (l_noco) length=1 + ALLOCATE(d%eig_eig(neig,jspins*nkpts)) !d%eig_vec #ifdef CPP_INVERSION - allocate(d%eig_vecr(nmat*neig,length*nkpts)) + ALLOCATE(d%eig_vecr(nmat*neig,length*nkpts)) #ifdef CPP_SOC - call judft_error("SOC+INVERSION can not be used with eigenvalues stored in memory") + CALL judft_error("SOC+INVERSION can not be used with eigenvalues stored in memory") #endif #else - allocate(d%eig_vecc(nmat*neig,length*nkpts)) + ALLOCATE(d%eig_vecc(nmat*neig,length*nkpts)) #endif - if (present(filename)) call priv_readfromfile() - contains - subroutine priv_readfromfile() - use m_eig66_da,ONLY:open_eig_IO=>open_eig,read_eig_IO=>read_eig,close_eig_IO=>close_eig - integer:: jspin,nk,i,ii,iii,nv,tmp_id - real :: wk,bk3(3),evac(2) - integer :: k1(nmat),k2(nmat),k3(nmat),kveclo(nlotot) - real :: eig(neig),ello(d%nlo,d%ntype),el(d%lmax,d%ntype) + length=length*nkpts + IF (l_dos) THEN + ALLOCATE(d%qal(0:3,ntype,neig,length)) + ALLOCATE(d%qvac(neig,2,length)) + ALLOCATE(d%qis(neig,length)) + ALLOCATE(d%qvlay(neig,layers,2,length)) + ALLOCATE(d%qstars(nstars,neig,layers,2,length)) + ALLOCATE(d%ksym(neig,length)) + ALLOCATE(d%jsym(neig,length)) + IF (l_mcd) ALLOCATE(d%mcd(3*ntype,ncored,neig,length)) + IF (l_orb) THEN + ALLOCATE(d%qintsl(nsld,neig,length)) + ALLOCATE(d%qmtsl(nsld,neig,length)) + ALLOCATE(d%qmtp(neig,nat,length)) + ALLOCATE(d%orbcomp(neig,23,nat,length)) + ENDIF + ENDIF + IF (PRESENT(filename)) CALL priv_readfromfile() + CONTAINS + SUBROUTINE priv_readfromfile() + USE m_eig66_da,ONLY:open_eig_IO=>open_eig,read_eig_IO=>read_eig,close_eig_IO=>close_eig + INTEGER:: jspin,nk,i,ii,iii,nv,tmp_id + REAL :: wk,bk3(3),evac(2) + INTEGER :: k1(nmat),k2(nmat),k3(nmat),kveclo(nlotot) + REAL :: eig(neig),ello(d%nlo,d%ntype),el(d%lmax,d%ntype) #ifdef CPP_INVERSION - real :: z(nmat,neig) + REAL :: z(nmat,neig) #else - complex :: z(nmat,neig) + COMPLEX :: z(nmat,neig) #endif - tmp_id=eig66_data_newid(DA_mode) - call open_eig_IO(tmp_id,nmat,neig,nkpts,jspins,d%lmax,d%nlo,d%ntype,nlotot,.false.,filename) - DO jspin=1,jspins - DO nk=1,nkpts - call read_eig_IO(tmp_id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) - call write_eig(id,nk,jspin,ii,ii,nv,i,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) - ENDDO - ENDDO - call close_eig_IO(tmp_id) - end subroutine - - end subroutine open_eig - - subroutine close_eig(id,delete,filename) - integer,intent(in) :: id - logical,intent(in),optional::delete - character(len=*),optional,intent(in)::filename - TYPE(t_data_mem),pointer:: d - call priv_find_data(id,d) - - if (present(filename)) call priv_writetofile() - - if (present(delete)) THEN - if (delete) THEN - if (allocated(d%eig_int)) deallocate(d%eig_int) - if (allocated(d%eig_real)) deallocate(d%eig_real) - if (allocated(d%eig_eig)) deallocate(d%eig_eig) - if (allocated(d%eig_vecr)) deallocate(d%eig_vecr) - if (allocated(d%eig_vecc)) deallocate(d%eig_vecc) - endif - endif - contains - subroutine priv_writetofile() - use m_eig66_DA,ONLY:open_eig_DA=>open_eig,write_eig_DA=>write_eig,close_eig_DA=>close_eig - implicit none - - integer:: nlotot,nk,jspin,nv,i,ii,tmp_id - real :: wk,bk3(3),evac(2) - integer :: k1(d%nmat),k2(d%nmat),k3(d%nmat),kveclo(size(d%eig_int,1)-3-3*d%nmat) - real :: eig(size(d%eig_eig,1)),ello(d%nlo,d%ntype),el(d%lmax,d%ntype) + tmp_id=eig66_data_newid(DA_mode) + IF (l_dos) CPP_error("Can not read DOS-data") + CALL open_eig_IO(tmp_id,nmat,neig,nkpts,jspins,d%lmax,d%nlo,d%ntype,nlotot,.FALSE.,.FALSE.,.FALSE.,.FALSE.,filename) + DO jspin=1,jspins + DO nk=1,nkpts + CALL read_eig_IO(tmp_id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) + CALL write_eig(id,nk,jspin,ii,ii,nv,i,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) + ENDDO + ENDDO + CALL close_eig_IO(tmp_id) + END SUBROUTINE priv_readfromfile + + END SUBROUTINE open_eig + + SUBROUTINE close_eig(id,delete,filename) + INTEGER,INTENT(in) :: id + LOGICAL,INTENT(in),OPTIONAL::delete + CHARACTER(len=*),OPTIONAL,INTENT(in)::filename + TYPE(t_data_mem),POINTER:: d + CALL priv_find_data(id,d) + + IF (PRESENT(filename)) CALL priv_writetofile() + + IF (PRESENT(delete)) THEN + IF (delete) THEN + IF (ALLOCATED(d%eig_int)) DEALLOCATE(d%eig_int) + IF (ALLOCATED(d%eig_real)) DEALLOCATE(d%eig_real) + IF (ALLOCATED(d%eig_eig)) DEALLOCATE(d%eig_eig) + IF (ALLOCATED(d%eig_vecr)) DEALLOCATE(d%eig_vecr) + IF (ALLOCATED(d%eig_vecc)) DEALLOCATE(d%eig_vecc) + ENDIF + ENDIF + CONTAINS + SUBROUTINE priv_writetofile() + USE m_eig66_DA,ONLY:open_eig_DA=>open_eig,write_eig_DA=>write_eig,close_eig_DA=>close_eig + IMPLICIT NONE + + INTEGER:: nlotot,nk,jspin,nv,i,ii,tmp_id + REAL :: wk,bk3(3),evac(2) + INTEGER :: k1(d%nmat),k2(d%nmat),k3(d%nmat),kveclo(SIZE(d%eig_int,1)-3-3*d%nmat) + REAL :: eig(SIZE(d%eig_eig,1)),ello(d%nlo,d%ntype),el(d%lmax,d%ntype) #ifdef CPP_INVERSION - real :: z(d%nmat,size(d%eig_eig,1)) + REAL :: z(d%nmat,SIZE(d%eig_eig,1)) #else - complex :: z(d%nmat,size(d%eig_eig,1)) + COMPLEX :: z(d%nmat,SIZE(d%eig_eig,1)) #endif - tmp_id=eig66_data_newid(DA_mode) - call open_eig_DA(tmp_id,d%nmat,d%neig,d%nkpts,d%jspins,d%lmax,d%nlo,d%ntype,d%nlotot,.false.,filename) - DO jspin=1,d%jspins - DO nk=1,d%nkpts - call read_eig(id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) - call write_eig_DA(tmp_id,nk,jspin,ii,ii,nv,i,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) - ENDDO - ENDDO - call close_eig_DA(tmp_id) - call eig66_remove_data(id) - end subroutine - end subroutine close_eig - - subroutine read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,& - ello,evac,kveclo,n_start,n_end,z) - implicit none - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat - INTEGER, INTENT(OUT),OPTIONAL :: neig - REAL, INTENT(OUT),OPTIONAL :: eig(:) - INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) - REAL, INTENT(OUT),OPTIONAL :: bk(:),wk - INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end - CLASS(*),INTENT(OUT),OPTIONAL :: z(:,:) - - INTEGER::nrec - TYPE(t_data_mem),pointer:: d - call priv_find_data(id,d) - - nrec=nk+(jspin-1)*d%nkpts - ! data from d%eig_int - if (present(nv)) nv=d%eig_int(1,nrec) - if (present(nmat)) nmat=d%eig_int(2,nrec) - if (present(neig)) then - neig=d%eig_int(3,nrec) - endif - if (present(k1)) then - if (.not.present(k2).or..not.present(k3)) call juDFT_error("BUG: always read k1,k2,k3") - k1=d%eig_int(3+1:3+d%nmat,nrec) - k2=d%eig_int(3+d%nmat+1:3+2*d%nmat,nrec) - k3=d%eig_int(3+2*d%nmat+1:3+3*d%nmat,nrec) - endif - if (present(kveclo)) kveclo=d%eig_int(4+3*d%nmat:3+3*d%nmat+size(kveclo),nrec) - - !data from d%eig_real - if (present(bk)) bk=d%eig_real(1:3,nrec) - if (present(wk)) wk=d%eig_real(4,nrec) - if (present(evac)) evac=d%eig_real(5:6,nrec) - if (present(el)) el=reshape(d%eig_real(7:7+size(el)-1,nrec),shape(el)) - if (present(ello)) ello=reshape(d%eig_real(size(d%eig_real,1)-size(ello)+1:,nrec),shape(ello)) - - !data from d%eig_eig - if (present(eig)) THEN - eig=0.0 - eig=d%eig_eig(:size(eig),nrec) - !print *,"R-eig:",nrec,shape(eig) - !print*,"R-eig(data):",shape(d%eig_eig) - !print*,"R:",eig - ENDIF - !data from d%eig_vec - - if (present(z)) then - write(*,*) "R-Z:",nrec,shape(z) - SELECT TYPE(z) - TYPE is (real) - if (.not.allocated(d%eig_vecr)) call juDFT_error("BUG: can not read real vectors from memory") - z=reshape(d%eig_vecr(:size(z),nrec),shape(z)) - TYPE is (complex) - if (.not.allocated(d%eig_vecc)) call juDFT_error("BUG: can not read complex vectors from memory") - z=reshape(d%eig_vecc(:size(z),nrec),shape(z)) - END SELECT - endif - end subroutine read_eig - - - subroutine write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk, & - eig,el,ello,evac, & - nlotot,kveclo,n_size,n_rank,z) - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank - REAL, INTENT(IN),OPTIONAL :: wk - INTEGER, INTENT(IN),OPTIONAL :: neig,neig_total,nv,nmat,nlotot - INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:) - REAL, INTENT(IN),OPTIONAL :: evac(:),ello(:,:) - CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) - INTEGER::nrec - TYPE(t_data_mem),pointer:: d - call priv_find_data(id,d) - - nrec=nk+(jspin-1)*d%nkpts - ! data from d%eig_int - if (present(nv)) d%eig_int(1,nrec)=nv - if (present(nmat)) d%eig_int(2,nrec)=nmat - if (present(neig)) THEN - if (present(neig_total)) THEN - if (neig.ne.neig_total) STOP "BUG in eig_mem" - d%eig_int(3,nrec)=neig_total - else - STOP "BUG2 in eig_mem" - endif - endif - - if (present(k1)) then - if (.not.present(k2).or..not.present(k3)) call juDFT_error("BUG: always write k1,k2,k3") - d%eig_int(3+1:3+d%nmat,nrec)=k1 - d%eig_int(3+d%nmat+1:3+2*d%nmat,nrec)=k2 - d%eig_int(3+2*d%nmat+1:3+3*d%nmat,nrec)=k3 - endif - if (present(kveclo)) d%eig_int(4+3*d%nmat:3+3*d%nmat+size(kveclo),nrec)=kveclo - - !data from d%eig_real - if (present(bk)) d%eig_real(1:3,nrec)=bk - if (present(wk)) d%eig_real(4,nrec)=wk - if (present(evac)) d%eig_real(5:6,nrec)=evac - if (present(el)) d%eig_real(7:7+size(el)-1,nrec)=reshape(el,(/size(el)/)) - if (present(ello)) d%eig_real(size(d%eig_real,1)-size(ello)+1:,nrec)=reshape(ello,(/size(ello)/)) - !data from d%eig_eig - if (present(eig)) then - d%eig_eig(:size(eig),nrec)=eig - !print*,"W-eig:",nrec,shape(eig) - !print*,"W:",eig - endif - !data from d%eig_vec - if (present(z)) then - write(*,*) "W-Z:",nrec,shape(z) - - SELECT TYPE(z) - TYPE IS (real) - if (.not.allocated(d%eig_vecr)) call juDFT_error("BUG: can not write real vectors to memory") - d%eig_vecr(:size(z),nrec)=reshape(real(z),(/size(z)/)) - TYPE IS(complex) - if (.not.allocated(d%eig_vecc)) call juDFT_error("BUG: can not write complex vectors to memory") - d%eig_vecc(:size(z),nrec)=reshape(cmplx(z),(/size(z)/)) - END SELECT - endif - - - end subroutine write_eig - - -end module m_eig66_mem + tmp_id=eig66_data_newid(DA_mode) + IF (d%l_dos) CPP_error("Could not write DOS data") + CALL open_eig_DA(tmp_id,d%nmat,d%neig,d%nkpts,d%jspins,d%lmax,d%nlo,d%ntype,d%nlotot,.FALSE.,.FALSE.,.FALSE.,.FALSE.,filename) + DO jspin=1,d%jspins + DO nk=1,d%nkpts + CALL read_eig(id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) + CALL write_eig_DA(tmp_id,nk,jspin,ii,ii,nv,i,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) + ENDDO + ENDDO + CALL close_eig_DA(tmp_id) + CALL eig66_remove_data(id) + END SUBROUTINE priv_writetofile + END SUBROUTINE close_eig + + SUBROUTINE write_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(IN) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(IN) :: qstars(:,:,:,:) + INTEGER,INTENT(IN) :: ksym(:),jsym(:) + REAL,INTENT(IN),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(IN),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + + INTEGER::nrec + TYPE(t_data_mem),POINTER:: d + CALL priv_find_data(id,d) + + nrec=nk+(jspin-1)*d%nkpts + + d%qal(:,:,:,nrec)=qal + d%qvac(:,:,nrec)=qvac + d%qis(:,nrec)=qis + d%qvlay(:,:,:,nrec)=qvlay + d%qstars(:,:,:,:,nrec)=qstars + d%ksym(:,nrec)=ksym + d%jsym(:,nrec)=jsym + IF (d%l_mcd.AND.PRESENT(mcd)) d%mcd(:,:,:,nrec)=mcd + IF (d%l_orb.AND.PRESENT(qintsl)) THEN + d%qintsl(:,:,nrec)=qintsl + d%qmtsl(:,:,nrec)=qmtsl + d%qmtp(:,:,nrec)=qmtp + d%orbcomp(:,:,:,nrec)=orbcomp + ENDIF + END SUBROUTINE write_dos + + SUBROUTINE read_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(OUT) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(OUT) :: qstars(:,:,:,:) + INTEGER,INTENT(OUT) :: ksym(:),jsym(:) + REAL,INTENT(OUT),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(OUT),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + + INTEGER::nrec + TYPE(t_data_mem),POINTER:: d + CALL priv_find_data(id,d) + + nrec=nk+(jspin-1)*d%nkpts + + qal=d%qal(:,:,:,nrec) + qvac=d%qvac(:,:,nrec) + qis=d%qis(:,nrec) + qvlay=d%qvlay(:,:,:,nrec) + qstars=d%qstars(:,:,:,:,nrec) + ksym=d%ksym(:,nrec) + jsym=d%jsym(:,nrec) + IF (d%l_mcd.AND.PRESENT(mcd)) mcd=d%mcd(:,:,:,nrec) + IF (d%l_orb.AND.PRESENT(qintsl)) THEN + qintsl=d%qintsl(:,:,nrec) + qmtsl=d%qmtsl(:,:,nrec) + qmtp=d%qmtp(:,:,nrec) + orbcomp=d%orbcomp(:,:,:,nrec) + ENDIF + END SUBROUTINE read_dos + + + SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk,wk,neig,eig,el,& + ello,evac,kveclo,n_start,n_end,z) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat + INTEGER, INTENT(OUT),OPTIONAL :: neig + REAL, INTENT(OUT),OPTIONAL :: eig(:) + INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) + REAL, INTENT(OUT),OPTIONAL :: bk(:),wk + INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end + CLASS(*),INTENT(OUT),OPTIONAL :: z(:,:) + + INTEGER::nrec + TYPE(t_data_mem),POINTER:: d + CALL priv_find_data(id,d) + + nrec=nk+(jspin-1)*d%nkpts + ! data from d%eig_int + IF (PRESENT(nv)) nv=d%eig_int(1,nrec) + IF (PRESENT(nmat)) nmat=d%eig_int(2,nrec) + IF (PRESENT(neig)) THEN + neig=d%eig_int(3,nrec) + ENDIF + IF (PRESENT(k1)) THEN + IF (.NOT.PRESENT(k2).OR..NOT.PRESENT(k3)) CALL juDFT_error("BUG: always read k1,k2,k3") + k1=d%eig_int(3+1:3+d%nmat,nrec) + k2=d%eig_int(3+d%nmat+1:3+2*d%nmat,nrec) + k3=d%eig_int(3+2*d%nmat+1:3+3*d%nmat,nrec) + ENDIF + IF (PRESENT(kveclo)) kveclo=d%eig_int(4+3*d%nmat:3+3*d%nmat+SIZE(kveclo),nrec) + + !data from d%eig_real + IF (PRESENT(bk)) bk=d%eig_real(1:3,nrec) + IF (PRESENT(wk)) wk=d%eig_real(4,nrec) + IF (PRESENT(evac)) evac=d%eig_real(5:6,nrec) + IF (PRESENT(el)) el=RESHAPE(d%eig_real(7:7+SIZE(el)-1,nrec),SHAPE(el)) + IF (PRESENT(ello)) ello=RESHAPE(d%eig_real(SIZE(d%eig_real,1)-SIZE(ello)+1:,nrec),SHAPE(ello)) + + !data from d%eig_eig + IF (PRESENT(eig)) THEN + eig=0.0 + eig=d%eig_eig(:SIZE(eig),nrec) + !print *,"R-eig:",nrec,shape(eig) + !print*,"R-eig(data):",shape(d%eig_eig) + !print*,"R:",eig + ENDIF + !data from d%eig_vec + + IF (PRESENT(z)) THEN + WRITE(*,*) "R-Z:",nrec,SHAPE(z) + SELECT TYPE(z) + TYPE is (REAL) + IF (.NOT.ALLOCATED(d%eig_vecr)) CALL juDFT_error("BUG: can not read real vectors from memory") + z=RESHAPE(d%eig_vecr(:SIZE(z),nrec),SHAPE(z)) + TYPE is (COMPLEX) + IF (.NOT.ALLOCATED(d%eig_vecc)) CALL juDFT_error("BUG: can not read complex vectors from memory") + z=RESHAPE(d%eig_vecc(:SIZE(z),nrec),SHAPE(z)) + END SELECT + ENDIF + END SUBROUTINE read_eig + + + SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk,wk, & + eig,el,ello,evac, & + nlotot,kveclo,n_size,n_rank,z) + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank + REAL, INTENT(IN),OPTIONAL :: wk + INTEGER, INTENT(IN),OPTIONAL :: neig,neig_total,nv,nmat,nlotot + INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(IN),OPTIONAL :: bk(3),eig(:),el(:,:) + REAL, INTENT(IN),OPTIONAL :: evac(:),ello(:,:) + CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) + INTEGER::nrec + TYPE(t_data_mem),POINTER:: d + CALL priv_find_data(id,d) + + nrec=nk+(jspin-1)*d%nkpts + ! data from d%eig_int + IF (PRESENT(nv)) d%eig_int(1,nrec)=nv + IF (PRESENT(nmat)) d%eig_int(2,nrec)=nmat + IF (PRESENT(neig)) THEN + IF (PRESENT(neig_total)) THEN + IF (neig.NE.neig_total) STOP "BUG in eig_mem" + d%eig_int(3,nrec)=neig_total + ELSE + STOP "BUG2 in eig_mem" + ENDIF + ENDIF + + IF (PRESENT(k1)) THEN + IF (.NOT.PRESENT(k2).OR..NOT.PRESENT(k3)) CALL juDFT_error("BUG: always write k1,k2,k3") + d%eig_int(3+1:3+d%nmat,nrec)=k1 + d%eig_int(3+d%nmat+1:3+2*d%nmat,nrec)=k2 + d%eig_int(3+2*d%nmat+1:3+3*d%nmat,nrec)=k3 + ENDIF + IF (PRESENT(kveclo)) d%eig_int(4+3*d%nmat:3+3*d%nmat+SIZE(kveclo),nrec)=kveclo + + !data from d%eig_real + IF (PRESENT(bk)) d%eig_real(1:3,nrec)=bk + IF (PRESENT(wk)) d%eig_real(4,nrec)=wk + IF (PRESENT(evac)) d%eig_real(5:6,nrec)=evac + IF (PRESENT(el)) d%eig_real(7:7+SIZE(el)-1,nrec)=RESHAPE(el,(/SIZE(el)/)) + IF (PRESENT(ello)) d%eig_real(SIZE(d%eig_real,1)-SIZE(ello)+1:,nrec)=RESHAPE(ello,(/SIZE(ello)/)) + !data from d%eig_eig + IF (PRESENT(eig)) THEN + d%eig_eig(:SIZE(eig),nrec)=eig + !print*,"W-eig:",nrec,shape(eig) + !print*,"W:",eig + ENDIF + !data from d%eig_vec + IF (PRESENT(z)) THEN + WRITE(*,*) "W-Z:",nrec,SHAPE(z) + + SELECT TYPE(z) + TYPE IS (REAL) + IF (.NOT.ALLOCATED(d%eig_vecr)) CALL juDFT_error("BUG: can not write real vectors to memory") + d%eig_vecr(:SIZE(z),nrec)=RESHAPE(REAL(z),(/SIZE(z)/)) + TYPE IS(COMPLEX) + IF (.NOT.ALLOCATED(d%eig_vecc)) CALL juDFT_error("BUG: can not write complex vectors to memory") + d%eig_vecc(:SIZE(z),nrec)=RESHAPE(CMPLX(z),(/SIZE(z)/)) + END SELECT + ENDIF + + + END SUBROUTINE write_eig + + +END MODULE m_eig66_mem diff --git a/io/eig66_mpi.F90 b/io/eig66_mpi.F90 index e47e30bc2dcaa122bb83ca0679e5158be9609267..4bdc66704f22a8b9315a8e0f86a7e3fa95c6fffc 100644 --- a/io/eig66_mpi.F90 +++ b/io/eig66_mpi.F90 @@ -1,517 +1,652 @@ MODULE m_eig66_mpi #include "juDFT_env.h" - use m_eig66_data - IMPLICIT NONE + USE m_eig66_data + IMPLICIT NONE #ifdef CPP_MPI - include 'mpif.h' - PRIVATE + INCLUDE 'mpif.h' + PRIVATE #endif - PUBLIC open_eig,read_eig,write_eig,close_eig + PUBLIC open_eig,read_eig,write_eig,close_eig CONTAINS - subroutine priv_find_data(id,d) - INTEGER,INTENT(IN)::id - TYPE(t_data_mpi),POINTER:: d - - class(t_data),pointer ::dp - call eig66_find_data(dp,id) - select type(dp) - type is (t_data_mpi) - d=>dp - class default - call judft_error("BUG: wrong datatype in eig66_mpi") - END SELECT - END subroutine - - - SUBROUTINE open_eig(id,mpi_comm,nmat,neig,nkpts,jspins,lmax,nlo,ntype,create,nlotot,l_noco,n_size_opt,filename) - USE,INTRINSIC::iso_c_binding - IMPLICIT NONE - INTEGER, INTENT(IN) :: id,mpi_comm,nmat,neig,nkpts,jspins,nlo,ntype,lmax,nlotot - LOGICAL, INTENT(IN) :: l_noco,create - INTEGER,INTENT(IN),OPTIONAL:: n_size_opt - CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename -#ifdef CPP_MPI - INTEGER:: isize,e,slot_size,type_size,local_slots - INTEGER(MPI_ADDRESS_KIND):: length - TYPE(c_ptr)::ptr - TYPE(t_data_MPI),pointer :: d - - call priv_find_data(id,d) - - IF (present(n_size_opt)) d%n_size=n_size_opt - IF (allocated(d%pe_ev)) THEN - IF (create) THEN - d%neig_data=0 - d%eig_data=1E99 - d%int_data=9999999 - d%real_data=1E99 -#ifdef CPP_INVERSION - d%zr_data=0.0 -#else - d%zc_data=0.0 -#endif - ENDIF - if (present(filename)) call priv_readfromfile() - RETURN !everything already done! - ENDIF - - CALL timestart("preperation") - CALL MPI_COMM_RANK(MPI_COMM,d%irank,e) - CALL MPI_COMM_SIZE(MPI_COMM,isize,e) - d%jspins=jspins - d%nkpts=nkpts - d%nmat=nmat - d%neig=neig - d%nlotot=nlotot - d%lmax=lmax - d%nlo=nlo - d%ntype=ntype - - - CALL timestart("map creation") - CALL create_maps(d,isize,nkpts,jspins,neig,d%n_size) - CALL timestop("map creation") - local_slots=count(d%pe_basis==d%irank) - !Now create the windows - - !Window for neig - CALL MPI_TYPE_SIZE(MPI_INTEGER,type_size,e) - slot_size=1*type_size - length=slot_size*local_slots - CALL timestop("preperation") - !Allocate memory - CALL timestart("Allocate") - CALL MPI_ALLOC_MEM(length,MPI_INFO_NULL,ptr,e) - CALL c_f_pointer(ptr,d%neig_data,(/length/type_size/)) - IF (e.NE.0) CPP_error("Could not allocated MPI-Data (neig) in eig66_mpi") - !create window - d%neig_data=0 - CALL timestop("Allocate") - CALL timestart("win create") - CALL MPI_WIN_CREATE(d%neig_data,length,slot_size,Mpi_INFO_NULL,MPI_COMM, d%neig_handle, e) - CALL timestop("win create") - - - !The integer values - d%size_k=nmat - slot_size=(5+3*d%size_k+1+nlotot)*type_size - length=slot_size*local_slots - !Allocate memory - CALL timestart("Allocate") - CALL MPI_ALLOC_MEM(length,MPI_INFO_NULL,ptr,e) - CALL c_f_pointer(ptr,d%int_data,(/length/type_size/)) - IF (e.NE.0) CPP_error("Could not allocated MPI-Data (int) in eig66_mpi") - !create window - d%int_data=9999999 - CALL timestop("Allocate") - CALL timestart("win create") - CALL MPI_WIN_CREATE(d%int_data, length,slot_size,Mpi_INFO_NULL, MPI_COMM, d%int_handle, e) - CALL timestop("win create") - CALL timestart("Allocate") - - !The real values - CALL MPI_TYPE_SIZE(MPI_DOUBLE_PRECISION,type_size,e) - d%size_el=(1+lmax)*ntype - d%size_ello=nlo*ntype - slot_size=(6+d%size_el+d%size_ello)*type_size - length=slot_size*local_slots - !Allocate memory - CALL MPI_ALLOC_MEM(length,MPI_INFO_NULL,ptr,e) - CALL c_f_pointer(ptr,d%real_data,(/length/type_size/)) - IF (e.NE.0) CPP_error("Could not allocated MPI-Data (real) in eig66_mpi") - !create window - d%real_data=1E99 - CALL timestop("Allocate") - CALL timestart("win create") - CALL MPI_WIN_CREATE(d%real_data, length,slot_size,Mpi_INFO_NULL, MPI_COMM, d%real_handle, e) - CALL timestop("win create") - CALL timestart("Allocate") - - !The eigenvalues - d%size_eig=neig - slot_size=d%size_eig*type_size - length=slot_size*local_slots - !Allocate memory - CALL MPI_ALLOC_MEM(length,MPI_INFO_NULL,ptr,e) - CALL c_f_pointer(ptr,d%eig_data,(/length/type_size/)) - IF (e.NE.0) CPP_error("Could not allocated MPI-Data (real) in eig66_mpi") - !create window - d%eig_data=1E99 - CALL timestop("Allocate") - CALL timestart("win create") - CALL MPI_WIN_CREATE(d%eiG_data, length,slot_size,Mpi_INFO_NULL, MPI_COMM, d%eig_handle, e) - CALL timestop("win create") - !The eigenvectors -#ifndef CPP_INVERSION - CALL MPI_TYPE_SIZE(MPI_DOUBLE_COMPLEX,type_size,e) -#endif - local_slots=count(d%pe_ev==d%irank) - slot_size=nmat*type_size - length=slot_size*local_slots - !Allocate memory - CALL timestart("Allocate") - CALL MPI_ALLOC_MEM(length,MPI_INFO_NULL,ptr,e) - IF (e.NE.0) CPP_error("Could not allocated MPI-Data (real) in eig66_mpi") + SUBROUTINE priv_find_data(id,d) + INTEGER,INTENT(IN)::id + TYPE(t_data_mpi),POINTER:: d -#ifdef CPP_INVERSION - !create window + CLASS(t_data),POINTER ::dp + CALL eig66_find_data(dp,id) + SELECT TYPE(dp) + TYPE is (t_data_mpi) + d=>dp + CLASS default + CALL judft_error("BUG: wrong datatype in eig66_mpi") + END SELECT + END SUBROUTINE priv_find_data - CALL c_f_pointer(ptr,d%zr_data,(/length/type_size/)) - CALL timestop("Allocate") - CALL timestart("win create") - CALL MPI_WIN_CREATE(d%zr_data, length,slot_size,Mpi_INFO_NULL, MPI_COMM,d%zr_handle, e) - CALL timestop("win create") -#ifdef CPP_SOC - CALL judft_error("SOC+INVERSION can not be used with eigenvalues stored in memory") -#endif + SUBROUTINE open_eig(id,mpi_comm,nmat,neig,nkpts,jspins,lmax,nlo,ntype,create,nlotot,l_noco,n_size_opt,l_dos,l_mcd,l_orb,filename,layers,nstars,ncored,nsld,nat) + USE,INTRINSIC::iso_c_binding + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,mpi_comm,nmat,neig,nkpts,jspins,nlo,ntype,lmax,nlotot + LOGICAL, INTENT(IN) :: l_noco,create + INTEGER,INTENT(IN),OPTIONAL:: n_size_opt + LOGICAL,INTENT(IN) ::l_dos,l_mcd,l_orb + CHARACTER(LEN=*),INTENT(IN),OPTIONAL :: filename + INTEGER,INTENT(IN),OPTIONAL :: layers,nstars,ncored,nsld,nat +#ifdef CPP_MPI + INTEGER:: isize,e,slot_size,local_slots + TYPE(t_data_MPI),POINTER :: d + + CALL priv_find_data(id,d) + CALL eig66_data_storedefault(d,jspins,nkpts,nmat,neig,lmax,nlotot,nlo,ntype,l_dos,l_mcd,l_orb) + + IF (PRESENT(n_size_opt)) d%n_size=n_size_opt + IF (ALLOCATED(d%pe_ev)) THEN + IF (create) THEN + d%neig_data=0 + d%eig_data=1E99 + d%int_data=9999999 + d%real_data=1E99 +#if defined(CPP_INVERSION)&&!defined(CPP_SOC) + d%zr_data=0.0 #else - !create window - - CALL c_f_pointer(ptr,d%zc_data,(/length/type_size/)) - CALL timestop("Allocate") - CALL timestart("win create") + d%zc_data=0.0 +#endif + d%qal_data=0.0 + d%qvac_data=0.0 + d%qvlay_data=0.0 + d%qstars_data=0.0 + d%ksym_data=0.0 + d%jsym_data=0.0 + d%mcd_data=0.0 + d%qintsl_data=0.0 + d%qmtsl_data=0.0 + d%qmtp_data=0.0 + d%orbcomp_data=0.0 + ENDIF + IF (PRESENT(filename)) CALL priv_readfromfile() + RETURN !everything already done! + ENDIF - CALL MPI_WIN_CREATE(d%zc_data, length,slot_size,Mpi_INFO_NULL, MPI_COMM, d%zc_handle, e) - CALL timestop("win create") + CALL timestart("create data spaces in ei66_mpi") + CALL MPI_COMM_RANK(MPI_COMM,d%irank,e) + CALL MPI_COMM_SIZE(MPI_COMM,isize,e) + + CALL create_maps(d,isize,nkpts,jspins,neig,d%n_size) + local_slots=COUNT(d%pe_basis==d%irank) + !Now create the windows + + !Window for neig + slot_size=1 + CALL priv_create_memory(1,local_slots,d%neig_data,d%neig_handle) + d%neig_data=0 + + !The integer values + d%size_k=nmat + slot_size=(5+3*d%size_k+1+nlotot) + CALL priv_create_memory(slot_size,local_slots,d%int_data,d%int_handle) + d%int_data=9999999 + + !The real values + d%size_el=(1+lmax)*ntype + d%size_ello=nlo*ntype + slot_size=(6+d%size_el+d%size_ello) + CALL priv_create_memory(slot_size,local_slots,d%real_data,d%real_handle) + d%real_data=1E99 + + !The eigenvalues + d%size_eig=neig + CALL priv_create_memory(d%size_eig,local_slots,d%eig_data,d%eig_handle) + d%eig_data=1E99 + !The eigenvectors + local_slots=COUNT(d%pe_ev==d%irank) + slot_size=nmat + +#if !defined(CPP_INVERSION)||defined(CPP_SOC) + CALL priv_create_memory(slot_size,local_slots,d%zr_data,d%zr_handle) +#else + CALL priv_create_memory(slot_size,local_slots,d%zc_data,d%zc_handle) #endif - if (present(filename).and..not.create) call priv_readfromfile() - contains - subroutine priv_readfromfile() - use m_eig66_DA,ONLY:open_eig_DA=>open_eig,read_eig_DA=>read_eig,close_eig_da=>close_eig - integer:: jspin,nk,i,ii,iii,nv,tmp_id - real :: wk,bk3(3),evac(2) - integer :: k1(nmat),k2(nmat),k3(nmat),kveclo(nlotot) - real :: eig(neig),ello(nlo,ntype),el(lmax,ntype) + !Data for DOS etc + IF (l_dos) THEN + local_slots=COUNT(d%pe_basis==d%irank) + CALL priv_create_memory(4*ntype*neig,local_slots,d%qal_data,d%qal_handle) + CALL priv_create_memory(neig*2,local_slots,d%qvac_data,d%qvac_handle) + CALL priv_create_memory(neig,local_slots,d%qis_data,d%qis_handle) + CALL priv_create_memory(neig*layers*2,local_slots,d%qvlay_data,d%qvlay_handle) + CALL priv_create_memory(nstars,neigd*layers*2,local_slots,d%qstars_data,d%qstars_handle) + CALL priv_create_memory(neig,local_slots,d%jsym_data,d%jsym_handle) + CALL priv_create_memory(neig,local_slots,d%ksym_data,d%ksym_handle) + IF (l_mcd) CALL priv_create_memory(3*ntype*mcored,neig,local_slots,d%mcd_data,d%mcd_handle) + IF (l_orb) THEN + CALL priv_create_memory(nsld*neig,local_slots,d%qintsl_data,d%qintsl_handle) + CALL priv_create_memory(nsld*neig,local_slots,d%qmtsl_data,d%qmtsl_handle) + CALL priv_create_memory(nat*neig,local_slots,d%qmtp_data,d%qmtp_handle) + CALL priv_create_memory(23*nat*neig,local_slots,d%orbcomp_data,d%orbcomp_handle) + ENDIF + ELSE + ALLOCATE(d%qal_data(1),d%qvac_data(1),d%qis_data(1),d%qvlay_data(1),d%qstars_data(1),& + d%jsym_data(1),d%ksym_data(1),d%mcd_data(1),d%qintsl_data(1),d%qmtsl_data(1),& + d%qmtp_data(1),d%orbcomp_data(1)) + ENDIF + IF (PRESENT(filename).AND..NOT.create) CALL priv_readfromfile() + CALL timestop("create data spaces in ei66_mpi") + CONTAINS + SUBROUTINE priv_create_memory(slot_size,local_slots,data_ptr,handle) + IMPLICIT NONE + INTEGER,INTENT(IN) :: slot_size,local_slots + CLASS(*),POINTER,INTENT(OUT) :: data_ptr + INTEGER,INTENT(OUT) :: handle + + TYPE(c_ptr)::ptr + INTEGER:: e + INTEGER(MPI_ADDRESS_KIND) :: length + INTEGER :: type_size,slot_size + + length=slot_size*local_slots + + SELECT TYPE(data_ptr) + TYPE IS (REAL) + CALL MPI_TYPE_SIZE(MPI_DOUBLE_PRECISION,type_size,e) + TYPE IS (COMPLEX) + CALL MPI_TYPE_SIZE(MPI_DOUBLE_COMPLEX,type_size,e) + TYPE IS (INTEGER) + CALL MPI_TYPE_SIZE(MPI_INTEGER,type_size,e) + END SELECT + length=length*type_size + + CALL MPI_ALLOC_MEM(length,MPI_INFO_NULL,ptr,e) + IF (e.NE.0) CPP_error("Could not allocated MPI-Data in eig66_mpi") + + CALL C_F_POINTER(ptr,data_ptr,(/length/type_size/)) + CALL MPI_WIN_CREATE(data_ptr, length,slot_size,Mpi_INFO_NULL, MPI_COMM,handle, e) + END SUBROUTINE priv_create_memory + SUBROUTINE priv_readfromfile() + USE m_eig66_DA,ONLY:open_eig_DA=>open_eig,read_eig_DA=>read_eig,close_eig_da=>close_eig + INTEGER:: jspin,nk,i,ii,iii,nv,tmp_id + REAL :: wk,bk3(3),evac(2) + INTEGER :: k1(nmat),k2(nmat),k3(nmat),kveclo(nlotot) + REAL :: eig(neig),ello(nlo,ntype),el(lmax,ntype) #ifdef CPP_INVERSION - real :: z(nmat,neig) + REAL :: z(nmat,neig) #else - complex :: z(nmat,neig) + COMPLEX :: z(nmat,neig) #endif - !only do this with PE=0 - if (d%irank==0) THEN - tmp_id=eig66_data_newid(DA_mode) - call open_eig_DA(tmp_id,nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,.false.,filename) - DO jspin=1,jspins - DO nk=1,nkpts - call read_eig_DA(tmp_id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) - call write_eig(id,nk,jspin,ii,ii,nv,nmat,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) - ENDDO + !only do this with PE=0 + IF (d%irank==0) THEN + tmp_id=eig66_data_newid(DA_mode) + IF (l_dos) CPP_error("Could not read DOS data") + CALL open_eig_DA(tmp_id,nmat,neig,nkpts,jspins,lmax,nlo,ntype,nlotot,.FALSE.,.FALSE.,.FALSE.,.FALSE.,filename) + DO jspin=1,jspins + DO nk=1,nkpts + CALL read_eig_DA(tmp_id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) + CALL write_eig(id,nk,jspin,ii,ii,nv,nmat,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) ENDDO - call close_eig_DA(tmp_id) - ENDIF - end subroutine + ENDDO + CALL close_eig_DA(tmp_id) + ENDIF + END SUBROUTINE priv_readfromfile #endif - END SUBROUTINE open_eig - SUBROUTINE close_eig(id,delete,filename) + END SUBROUTINE open_eig + SUBROUTINE close_eig(id,delete,filename) INTEGER,INTENT(IN) :: id LOGICAL,INTENT(IN),OPTIONAL:: delete CHARACTER(LEN=*),INTENT(IN),OPTIONAL::filename - TYPE(t_data_MPI),pointer :: d - call priv_find_data(id,d) + TYPE(t_data_MPI),POINTER :: d + CALL priv_find_data(id,d) - IF (present(delete)) THEN - IF (delete) WRITE(*,*) "No deallocation of memory implemented in eig66_mpi" + IF (PRESENT(delete)) THEN + IF (delete) WRITE(*,*) "No deallocation of memory implemented in eig66_mpi" ENDIF - if (present(filename)) call priv_writetofile() - contains - subroutine priv_writetofile() - use m_eig66_DA,ONLY:open_eig_DA=>open_eig,write_eig_DA=>write_eig,close_eig_DA=>close_eig - implicit none - - integer:: nlotot,nk,jspin,nv,i,ii,tmp_id - real :: wk,bk3(3),evac(2) - integer :: k1(d%nmat),k2(d%nmat),k3(d%nmat),kveclo(d%nlotot) - real :: eig(d%neig),ello(d%nlo,d%ntype),el(d%lmax,d%ntype) + IF (PRESENT(filename)) CALL priv_writetofile() + CONTAINS + SUBROUTINE priv_writetofile() + USE m_eig66_DA,ONLY:open_eig_DA=>open_eig,write_eig_DA=>write_eig,close_eig_DA=>close_eig + IMPLICIT NONE + + INTEGER:: nlotot,nk,jspin,nv,i,ii,tmp_id + REAL :: wk,bk3(3),evac(2) + INTEGER :: k1(d%nmat),k2(d%nmat),k3(d%nmat),kveclo(d%nlotot) + REAL :: eig(d%neig),ello(d%nlo,d%ntype),el(d%lmax,d%ntype) #ifdef CPP_INVERSION - real :: z(d%nmat,d%neig) + REAL :: z(d%nmat,d%neig) #else - complex :: z(d%nmat,d%neig) + COMPLEX :: z(d%nmat,d%neig) #endif - nlotot=d%nlotot + nlotot=d%nlotot - if (d%irank==0) THEN - tmp_id=eig66_data_newid(DA_mode) - call open_eig_DA(tmp_id,d%nmat,d%neig,d%nkpts,d%jspins,d%lmax,d%nlo,d%ntype,d%nlotot,.false.,filename) - DO jspin=1,d%jspins - DO nk=1,d%nkpts - call read_eig(id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) - call write_eig_DA(tmp_id,nk,jspin,ii,ii,nv,i,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) + IF (d%irank==0) THEN + tmp_id=eig66_data_newid(DA_mode) + IF (d%l_dos) CPP_error("Could not write DOS data") + CALL open_eig_DA(tmp_id,d%nmat,d%neig,d%nkpts,d%jspins,d%lmax,d%nlo,d%ntype,d%nlotot,.FALSE.,.FALSE.,.FALSE.,.FALSE.,filename) + DO jspin=1,d%jspins + DO nk=1,d%nkpts + CALL read_eig(id,nk,jspin,nv,i,k1,k2,k3,bk3,wk,ii,eig,el,ello,evac,kveclo,z=z) + CALL write_eig_DA(tmp_id,nk,jspin,ii,ii,nv,i,k1,k2,k3,bk3,wk,eig,el,ello,evac,nlotot,kveclo,z=z) + ENDDO ENDDO - ENDDO - CALL close_eig_DA(tmp_id) - ENDIF - call eig66_remove_data(id) - end subroutine - - END SUBROUTINE close_eig - - SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk3,wk,neig,eig,el,& - ello,evac,kveclo,n_start,n_end,z) - IMPLICIT NONE - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat - INTEGER, INTENT(OUT),OPTIONAL :: neig - REAL, INTENT(OUT),OPTIONAL :: eig(:) - INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) - REAL, INTENT(OUT),OPTIONAL :: bk3(:),wk - INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end - CLASS(*),TARGET,INTENT(OUT),OPTIONAL :: z(:,:) + CALL close_eig_DA(tmp_id) + ENDIF + CALL eig66_remove_data(id) + END SUBROUTINE priv_writetofile + + END SUBROUTINE close_eig + + SUBROUTINE read_eig(id,nk,jspin,nv,nmat,k1,k2,k3,bk3,wk,neig,eig,el,& + ello,evac,kveclo,n_start,n_end,z) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(OUT),OPTIONAL :: nv,nmat + INTEGER, INTENT(OUT),OPTIONAL :: neig + REAL, INTENT(OUT),OPTIONAL :: eig(:) + INTEGER, INTENT(OUT),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(OUT),OPTIONAL :: evac(:),ello(:,:),el(:,:) + REAL, INTENT(OUT),OPTIONAL :: bk3(:),wk + INTEGER, INTENT(IN),OPTIONAL :: n_start,n_end + CLASS(*),TARGET,INTENT(OUT),OPTIONAL :: z(:,:) #ifdef CPP_MPI - INTEGER :: pe,tmp_size,e - INTEGER(MPI_ADDRESS_KIND) :: slot - INTEGER :: n1,n2,n3,n - INTEGER,ALLOCATABLE :: tmp_int(:) - REAL,ALLOCATABLE :: tmp_real(:) - COMPLEX,ALLOCATABLE :: tmp_cmplx(:) - TYPE(t_data_MPI),pointer :: d - call priv_find_data(id,d) - pe=d%pe_basis(nk,jspin) - slot=d%slot_basis(nk,jspin) - IF (present(neig))THEN - CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%neig_handle,e) - ! Get current values - CALL MPI_GET(neig,1,MPI_INTEGER,pe,slot,1,MPI_INTEGER,d%neig_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%neig_handle,e) - - ENDIF - !read the integer values - IF (any((/present(nv),present(nmat),present(k1),present(k2),present(k3),present(kveclo)/))) THEN - tmp_size=4+3*d%size_k - IF (present(kveclo)) tmp_size=tmp_size+size(kveclo) - ALLOCATE(tmp_int(tmp_size)) - CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%int_handle,e) - ! Get current values - CALL MPI_GET(tmp_int,tmp_size,MPI_INTEGER,pe,slot,tmp_size,MPI_INTEGER,d%int_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%int_handle,e) - !IF (present(neig)) neig=tmp_int(1) - IF (present(nv)) nv=tmp_int(2) - IF (present(nmat)) nmat=tmp_int(3) - IF (present(k1)) k1=tmp_int(4+1:4+size(k1)) - IF (present(k2)) k2=tmp_int(4+d%size_k+1:4+d%size_k+size(k2)) - IF (present(k3)) k3=tmp_int(4+2*d%size_k+1:4+2*d%size_k+size(k3)) - IF (present(kveclo)) kveclo=tmp_int(4+3*d%size_k+1:4+3*d%size_k+size(kveclo)) - - ENDIF - !read the real-values - IF (any((/present(wk),present(bk3),present(el),present(ello),present(evac)/))) THEN - tmp_size=6+d%size_el+d%size_ello - ALLOCATE(tmp_real(tmp_size)) - CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%real_handle,e) - ! Get current values - CALL MPI_GET(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,d%real_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%real_handle,e) - IF (present(wk)) wk=tmp_real(1) - IF (present(bk3)) bk3=tmp_real(2:4) - IF (present(evac)) evac=tmp_real(5:6) - IF (present(el)) el=reshape(tmp_real(6+1:6+size(el)),shape(el)) - IF (present(ello)) ello=reshape(tmp_real(6+d%size_el+1:6+d%size_el+size(ello)),shape(ello)) - DEALLOCATE(tmp_real) - ENDIF - IF (present(eig)) THEN - CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%eig_handle,e) - ALLOCATE(tmp_real(d%size_eig)) - CALL MPI_GET(tmp_real,d%size_eig,MPI_DOUBLE_PRECISION,pe,slot,d%size_eig,MPI_DOUBLE_PRECISION,d%eig_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%eig_handle,e) - n1=1;n3=1;n2=size(eig) - IF (present(n_start)) n1=n_start - IF (present(n_end)) n2=n_end - eig(:n2-n1+1)=tmp_real(n1:n2) - DEALLOCATE(tmp_real) - ENDIF - - IF (present(z)) THEN - tmp_size=size(z,1) - ALLOCATE(tmp_real(tmp_size)) - ALLOCATE(tmp_cmplx(tmp_size)) - DO n=1,size(z,2) - n1=n - IF (present(n_start)) n1=n_start+n-1 - IF (present(n_end)) THEN - IF (n1>n_end) CYCLE - ENDIF - slot=d%slot_ev(nk,jspin,n1) - pe=d%pe_ev(nk,jspin,n1) - SELECT TYPE(z) - TYPE IS(real) - CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%zr_handle,e) - CALL MPI_GET(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,d%zr_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%zr_handle,e) - !print *, nk,jspin,n1,"r PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_real(1) - z(:,n)=tmp_real - TYPE IS (complex) - CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%zc_handle,e) - CALL MPI_GET(tmp_cmplx,tmp_size,MPI_DOUBLE_COMPLEX,pe,slot,tmp_size,MPI_DOUBLE_COMPLEX,d%zc_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%zc_handle,e) - !print *, nk,jspin,n1,"r PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_cmplx(1) - z(:,n)=tmp_cmplx - END SELECT - ENDDO - ENDIF + INTEGER :: pe,tmp_size,e + INTEGER(MPI_ADDRESS_KIND) :: slot + INTEGER :: n1,n2,n3,n + INTEGER,ALLOCATABLE :: tmp_int(:) + REAL,ALLOCATABLE :: tmp_real(:) + COMPLEX,ALLOCATABLE :: tmp_cmplx(:) + TYPE(t_data_MPI),POINTER :: d + CALL priv_find_data(id,d) + pe=d%pe_basis(nk,jspin) + slot=d%slot_basis(nk,jspin) + IF (PRESENT(neig))THEN + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%neig_handle,e) + ! Get current values + CALL MPI_GET(neig,1,MPI_INTEGER,pe,slot,1,MPI_INTEGER,d%neig_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%neig_handle,e) + + ENDIF + !read the integer values + IF (ANY((/PRESENT(nv),PRESENT(nmat),PRESENT(k1),PRESENT(k2),PRESENT(k3),PRESENT(kveclo)/))) THEN + tmp_size=4+3*d%size_k + IF (PRESENT(kveclo)) tmp_size=tmp_size+SIZE(kveclo) + ALLOCATE(tmp_int(tmp_size)) + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%int_handle,e) + ! Get current values + CALL MPI_GET(tmp_int,tmp_size,MPI_INTEGER,pe,slot,tmp_size,MPI_INTEGER,d%int_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%int_handle,e) + !IF (present(neig)) neig=tmp_int(1) + IF (PRESENT(nv)) nv=tmp_int(2) + IF (PRESENT(nmat)) nmat=tmp_int(3) + IF (PRESENT(k1)) k1=tmp_int(4+1:4+SIZE(k1)) + IF (PRESENT(k2)) k2=tmp_int(4+d%size_k+1:4+d%size_k+SIZE(k2)) + IF (PRESENT(k3)) k3=tmp_int(4+2*d%size_k+1:4+2*d%size_k+SIZE(k3)) + IF (PRESENT(kveclo)) kveclo=tmp_int(4+3*d%size_k+1:4+3*d%size_k+SIZE(kveclo)) + + ENDIF + !read the real-values + IF (ANY((/PRESENT(wk),PRESENT(bk3),PRESENT(el),PRESENT(ello),PRESENT(evac)/))) THEN + tmp_size=6+d%size_el+d%size_ello + ALLOCATE(tmp_real(tmp_size)) + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%real_handle,e) + ! Get current values + CALL MPI_GET(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,d%real_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%real_handle,e) + IF (PRESENT(wk)) wk=tmp_real(1) + IF (PRESENT(bk3)) bk3=tmp_real(2:4) + IF (PRESENT(evac)) evac=tmp_real(5:6) + IF (PRESENT(el)) el=RESHAPE(tmp_real(6+1:6+SIZE(el)),SHAPE(el)) + IF (PRESENT(ello)) ello=RESHAPE(tmp_real(6+d%size_el+1:6+d%size_el+SIZE(ello)),SHAPE(ello)) + DEALLOCATE(tmp_real) + ENDIF + IF (PRESENT(eig)) THEN + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%eig_handle,e) + ALLOCATE(tmp_real(d%size_eig)) + CALL MPI_GET(tmp_real,d%size_eig,MPI_DOUBLE_PRECISION,pe,slot,d%size_eig,MPI_DOUBLE_PRECISION,d%eig_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%eig_handle,e) + n1=1;n3=1;n2=SIZE(eig) + IF (PRESENT(n_start)) n1=n_start + IF (PRESENT(n_end)) n2=n_end + eig(:n2-n1+1)=tmp_real(n1:n2) + DEALLOCATE(tmp_real) + ENDIF + + IF (PRESENT(z)) THEN + tmp_size=SIZE(z,1) + ALLOCATE(tmp_real(tmp_size)) + ALLOCATE(tmp_cmplx(tmp_size)) + DO n=1,SIZE(z,2) + n1=n + IF (PRESENT(n_start)) n1=n_start+n-1 + IF (PRESENT(n_end)) THEN + IF (n1>n_end) CYCLE + ENDIF + slot=d%slot_ev(nk,jspin,n1) + pe=d%pe_ev(nk,jspin,n1) + SELECT TYPE(z) + TYPE IS(REAL) +#ifdef CPP_SOC + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%zc_handle,e) + CALL MPI_GET(tmp_cmplx,tmp_size,MPI_DOUBLE_COMPLEX,pe,slot,tmp_size,MPI_DOUBLE_COMPLEX,d%zc_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%zc_handle,e) + !print *, nk,jspin,n1,"r PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_cmplx(1) + z(:,n)=REAL(tmp_cmplx) +#else + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%zr_handle,e) + CALL MPI_GET(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,d%zr_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%zr_handle,e) + !print *, nk,jspin,n1,"r PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_real(1) + z(:,n)=tmp_real +#endif + TYPE IS (COMPLEX) + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,d%zc_handle,e) + CALL MPI_GET(tmp_cmplx,tmp_size,MPI_DOUBLE_COMPLEX,pe,slot,tmp_size,MPI_DOUBLE_COMPLEX,d%zc_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%zc_handle,e) + !print *, nk,jspin,n1,"r PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_cmplx(1) + z(:,n)=tmp_cmplx + END SELECT + ENDDO + ENDIF #endif - END SUBROUTINE read_eig - - SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk3,wk, & - eig,el,ello,evac, & - nlotot,kveclo,n_size,n_rank,z) - INTEGER, INTENT(IN) :: id,nk,jspin - INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank - REAL, INTENT(IN),OPTIONAL :: wk - INTEGER, INTENT(IN),OPTIONAL :: neig,nv,nmat,nlotot,neig_total - INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) - REAL, INTENT(IN),OPTIONAL :: bk3(3),eig(:),el(:,:) - REAL, INTENT(IN),OPTIONAL :: evac(:),ello(:,:) - CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) + END SUBROUTINE read_eig + + SUBROUTINE write_eig(id,nk,jspin,neig,neig_total,nv,nmat,k1,k2,k3,bk3,wk, & + eig,el,ello,evac, & + nlotot,kveclo,n_size,n_rank,z) + INTEGER, INTENT(IN) :: id,nk,jspin + INTEGER, INTENT(IN),OPTIONAL :: n_size,n_rank + REAL, INTENT(IN),OPTIONAL :: wk + INTEGER, INTENT(IN),OPTIONAL :: neig,nv,nmat,nlotot,neig_total + INTEGER, INTENT(IN),OPTIONAL :: k1(:),k2(:),k3(:),kveclo(:) + REAL, INTENT(IN),OPTIONAL :: bk3(3),eig(:),el(:,:) + REAL, INTENT(IN),OPTIONAL :: evac(:),ello(:,:) + CLASS(*),INTENT(IN),OPTIONAL :: z(:,:) #ifdef CPP_MPI - INTEGER :: pe,tmp_size,e - INTEGER(MPI_ADDRESS_KIND) :: slot - INTEGER :: n1,n2,n3,n - INTEGER,ALLOCATABLE :: tmp_int(:) - REAL,ALLOCATABLE :: tmp_real(:) - COMPLEX,ALLOCATABLE :: tmp_cmplx(:) - LOGICAL :: acc - TYPE(t_data_MPI),pointer :: d - - call priv_find_data(id,d) - - pe=d%pe_basis(nk,jspin) - slot=d%slot_basis(nk,jspin) - !write the number of eigenvalues values - IF (present(neig_total)) THEN - CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%neig_handle,e) - ALLOCATE(tmp_int(1)) - tmp_int(1)=neig_total - CALL MPI_PUT(tmp_int,1,MPI_INTEGER,pe,slot,1,MPI_INTEGER,d%neig_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%neig_handle,e) - DEALLOCATE(tmp_int) - ENDIF - - IF (any((/present(nv),present(nmat),present(nlotot),present(k1),present(k2),present(k3),present(kveclo)/))) THEN - tmp_size=5+3*d%size_k - IF (present(kveclo)) tmp_size=tmp_size+size(kveclo) - ALLOCATE(tmp_int(tmp_size)) - tmp_int=9999999 - tmp_int(1)=0 - IF (present(nv)) tmp_int(2)=nv - IF (present(nmat)) tmp_int(3)=nmat - IF (present(nlotot)) tmp_int(4)=nlotot - IF (present(k1)) tmp_int(4+1:4+size(k1))=k1 - IF (present(k2)) tmp_int(4+d%size_k+1:4+d%size_k+size(k2))=k2 - IF (present(k3)) tmp_int(4+2*d%size_k+1:4+2*d%size_k+size(k3))=k3 - IF (present(kveclo)) tmp_int(4+3*d%size_k+1:4+3*d%size_k+size(kveclo))=kveclo - CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%int_handle,e) - CALL MPI_ACCUMULATE(tmp_int,tmp_size,MPI_INTEGER,pe,slot,tmp_size,MPI_INTEGER,MPI_MIN,d%int_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%int_handle,e) - ENDIF - !write the real-values - IF (any((/present(wk),present(bk3),present(el),present(ello),present(evac)/))) THEN - tmp_size=6+d%size_el+d%size_ello - ALLOCATE(tmp_real(tmp_size)) - tmp_real=1E99 - IF (present(wk)) tmp_real(1)=wk - IF (present(bk3)) tmp_real(2:4)=bk3 - IF (present(evac)) tmp_real(5:6)=evac - IF (present(el)) tmp_real(6+1:6+size(el))=reshape(el,(/size(el)/)) - IF (present(ello)) tmp_real(6+d%size_el+1:6+d%size_el+size(ello))=reshape(ello,(/size(ello)/)) - - CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%real_handle,e) - CALL MPI_ACCUMULATE(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,MPI_MIN,d%real_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%real_handle,e) - DEALLOCATE(tmp_real) - ENDIF - IF (present(eig)) THEN - ALLOCATE(tmp_real(d%size_eig)) - tmp_real=1E99 - n1=1;n3=1 - IF (present(n_rank)) n1=n_rank+1 - IF (present(n_size)) n3=n_size - n2=size(eig)*n3+n1-1 - tmp_real(n1:n2:n3)=eig - CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%eig_handle,e) - CALL MPI_ACCUMULATE(tmp_real,d%size_eig,MPI_DOUBLE_PRECISION,pe,slot,d%size_eig,MPI_DOUBLE_PRECISION,MPI_MIN,d%eig_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%eig_handle,e) - DEALLOCATE(tmp_real) - ENDIF - IF (present(z)) THEN - tmp_size=size(z,1) - ALLOCATE(tmp_real(tmp_size)) - ALLOCATE(tmp_cmplx(tmp_size)) - DO n=1,size(z,2) - n1=n-1 - IF (present(n_size)) n1=n_size*n1 - IF (present(n_rank)) n1=n1+n_rank - slot=d%slot_ev(nk,jspin,n1+1) - pe=d%pe_ev(nk,jspin,n1+1) - !print *, "PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_real(1) - SELECT TYPE(z) - TYPE IS(REAL) - tmp_real=z(:,n) - CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%zr_handle,e) - CALL MPI_PUT(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,d%zr_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%zr_handle,e) - TYPE IS(complex) - tmp_cmplx=z(:,n) - CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%zc_handle,e) - CALL MPI_PUT(tmp_cmplx,tmp_size,MPI_DOUBLE_COMPLEX,pe,slot,tmp_size,MPI_DOUBLE_COMPLEX,d%zc_handle,e) - CALL MPI_WIN_UNLOCK(pe,d%zc_handle,e) - END SELECT - ENDDO - ENDIF + INTEGER :: pe,tmp_size,e + INTEGER(MPI_ADDRESS_KIND) :: slot + INTEGER :: n1,n2,n3,n + INTEGER,ALLOCATABLE :: tmp_int(:) + REAL,ALLOCATABLE :: tmp_real(:) + COMPLEX,ALLOCATABLE :: tmp_cmplx(:) + LOGICAL :: acc + TYPE(t_data_MPI),POINTER :: d + + CALL priv_find_data(id,d) + + pe=d%pe_basis(nk,jspin) + slot=d%slot_basis(nk,jspin) + !write the number of eigenvalues values + IF (PRESENT(neig_total)) THEN + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%neig_handle,e) + ALLOCATE(tmp_int(1)) + tmp_int(1)=neig_total + CALL MPI_PUT(tmp_int,1,MPI_INTEGER,pe,slot,1,MPI_INTEGER,d%neig_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%neig_handle,e) + DEALLOCATE(tmp_int) + ENDIF + + IF (ANY((/PRESENT(nv),PRESENT(nmat),PRESENT(nlotot),PRESENT(k1),PRESENT(k2),PRESENT(k3),PRESENT(kveclo)/))) THEN + tmp_size=5+3*d%size_k + IF (PRESENT(kveclo)) tmp_size=tmp_size+SIZE(kveclo) + ALLOCATE(tmp_int(tmp_size)) + tmp_int=9999999 + tmp_int(1)=0 + IF (PRESENT(nv)) tmp_int(2)=nv + IF (PRESENT(nmat)) tmp_int(3)=nmat + IF (PRESENT(nlotot)) tmp_int(4)=nlotot + IF (PRESENT(k1)) tmp_int(4+1:4+SIZE(k1))=k1 + IF (PRESENT(k2)) tmp_int(4+d%size_k+1:4+d%size_k+SIZE(k2))=k2 + IF (PRESENT(k3)) tmp_int(4+2*d%size_k+1:4+2*d%size_k+SIZE(k3))=k3 + IF (PRESENT(kveclo)) tmp_int(4+3*d%size_k+1:4+3*d%size_k+SIZE(kveclo))=kveclo + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%int_handle,e) + CALL MPI_ACCUMULATE(tmp_int,tmp_size,MPI_INTEGER,pe,slot,tmp_size,MPI_INTEGER,MPI_MIN,d%int_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%int_handle,e) + ENDIF + !write the real-values + IF (ANY((/PRESENT(wk),PRESENT(bk3),PRESENT(el),PRESENT(ello),PRESENT(evac)/))) THEN + tmp_size=6+d%size_el+d%size_ello + ALLOCATE(tmp_real(tmp_size)) + tmp_real=1E99 + IF (PRESENT(wk)) tmp_real(1)=wk + IF (PRESENT(bk3)) tmp_real(2:4)=bk3 + IF (PRESENT(evac)) tmp_real(5:6)=evac + IF (PRESENT(el)) tmp_real(6+1:6+SIZE(el))=RESHAPE(el,(/SIZE(el)/)) + IF (PRESENT(ello)) tmp_real(6+d%size_el+1:6+d%size_el+SIZE(ello))=RESHAPE(ello,(/SIZE(ello)/)) + + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%real_handle,e) + CALL MPI_ACCUMULATE(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,MPI_MIN,d%real_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%real_handle,e) + DEALLOCATE(tmp_real) + ENDIF + IF (PRESENT(eig)) THEN + ALLOCATE(tmp_real(d%size_eig)) + tmp_real=1E99 + n1=1;n3=1 + IF (PRESENT(n_rank)) n1=n_rank+1 + IF (PRESENT(n_size)) n3=n_size + n2=SIZE(eig)*n3+n1-1 + tmp_real(n1:n2:n3)=eig + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%eig_handle,e) + CALL MPI_ACCUMULATE(tmp_real,d%size_eig,MPI_DOUBLE_PRECISION,pe,slot,d%size_eig,MPI_DOUBLE_PRECISION,MPI_MIN,d%eig_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%eig_handle,e) + DEALLOCATE(tmp_real) + ENDIF + IF (PRESENT(z)) THEN + tmp_size=SIZE(z,1) + ALLOCATE(tmp_real(tmp_size)) + ALLOCATE(tmp_cmplx(tmp_size)) + DO n=1,SIZE(z,2) + n1=n-1 + IF (PRESENT(n_size)) n1=n_size*n1 + IF (PRESENT(n_rank)) n1=n1+n_rank + slot=d%slot_ev(nk,jspin,n1+1) + pe=d%pe_ev(nk,jspin,n1+1) + !print *, "PE:",pe," Slot: ",slot," Size:",tmp_size,tmp_real(1) + SELECT TYPE(z) + TYPE IS(REAL) +#ifdef CPP_SOC + tmp_cmplx=z(:,n) + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%zc_handle,e) + CALL MPI_PUT(tmp_cmplx,tmp_size,MPI_DOUBLE_COMPLEX,pe,slot,tmp_size,MPI_DOUBLE_COMPLEX,d%zc_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%zc_handle,e) +#else + tmp_real=z(:,n) + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%zr_handle,e) + CALL MPI_PUT(tmp_real,tmp_size,MPI_DOUBLE_PRECISION,pe,slot,tmp_size,MPI_DOUBLE_PRECISION,d%zr_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%zr_handle,e) +#endif + TYPE IS(COMPLEX) + tmp_cmplx=z(:,n) + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,d%zc_handle,e) + CALL MPI_PUT(tmp_cmplx,tmp_size,MPI_DOUBLE_COMPLEX,pe,slot,tmp_size,MPI_DOUBLE_COMPLEX,d%zc_handle,e) + CALL MPI_WIN_UNLOCK(pe,d%zc_handle,e) + END SELECT + ENDDO + ENDIF #endif - END SUBROUTINE write_eig + END SUBROUTINE write_eig + + + SUBROUTINE priv_put_data(pe,slot,DATA,handle) + IMPLICIT NONE + INTEGER,INTENT(IN) :: pe,slot + CLASS(*),INTENT(IN) :: DATA(:) + INTEGER,INTENT(IN) :: handle + + INTEGER :: len,e + INTEGER,ALLOCATABLE :: int_tmp(:) + REAL,ALLOCATABLE :: real_tmp(:) + COMPLEX,ALLOCATABLE:: cmplx_tmp(:) + INCLUDE 'mpif.h' + len=SIZE(DATA) + SELECT TYPE(DATA) + TYPE IS (INTEGER) + ALLOCATE(int_tmp(len)) + int_tmp=DATA + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,handle,e) + CALL MPI_PUT(int_tmp,len,MPI_INTEGER,pe,slot,len,MPI_INTEGER,handle,e) + CALL MPI_WIN_UNLOCK(pe,handle,e) + TYPE is (REAL) + ALLOCATE(real_tmp(len)) + real_tmp=DATA + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,handle,e) + CALL MPI_PUT(real_tmp,len,MPI_DOUBLE_PRECISION,pe,slot,len,MPI_DOUBLE_PRECISION,handle,e) + CALL MPI_WIN_UNLOCK(pe,handle,e) + TYPE is (COMPLEX) + ALLOCATE(cmplx_tmp(len)) + cmplx_tmp=DATA + CALL MPI_WIN_LOCK(MPI_LOCK_EXCLUSIVE,pe,0,handle,e) + CALL MPI_PUT(cmplx_tmp,len,MPI_DOUBLE_COMPLEX,pe,slot,len,MPI_DOUBLE_COMPLEX,handle,e) + CALL MPI_WIN_UNLOCK(pe,handle,e) + END SELECT + END SUBROUTINE priv_put_data + + SUBROUTINE priv_get_data(pe,slot,DATA,len,handle) + IMPLICIT NONE + INTEGER,INTENT(IN) :: pe,slot,len + CLASS(*),INTENT(OUT) :: DATA(len) + INTEGER,INTENT(IN) :: handle + + INTEGER :: e + INTEGER,ALLOCATABLE :: int_tmp(:) + REAL,ALLOCATABLE :: real_tmp(:) + COMPLEX,ALLOCATABLE:: cmplx_tmp(:) + INCLUDE 'mpif.h' + + SELECT TYPE(DATA) + TYPE IS (INTEGER) + ALLOCATE(int_tmp(len)) + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,handle,e) + CALL MPI_GET(int_tmp,len,MPI_INTEGER,pe,slot,len,MPI_INTEGER,handle,e) + CALL MPI_WIN_UNLOCK(pe,handle,e) + DATA=int_tmp + TYPE is (REAL) + ALLOCATE(real_tmp(len)) + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,handle,e) + CALL MPI_GET(real_tmp,len,MPI_DOUBLE_PRECISION,pe,slot,len,MPI_DOUBLE_PRECISION,handle,e) + CALL MPI_WIN_UNLOCK(pe,handle,e) + DATA=real_tmp + TYPE is (COMPLEX) + ALLOCATE(cmplx_tmp(len)) + CALL MPI_WIN_LOCK(MPI_LOCK_SHARED,pe,0,handle,e) + CALL MPI_GET(cmplx_tmp,len,MPI_DOUBLE_COMPLEX,pe,slot,len,MPI_DOUBLE_COMPLEX,handle,e) + CALL MPI_WIN_UNLOCK(pe,handle,e) + DATA=cmplx_tmp + END SELECT + END SUBROUTINE priv_get_data + + + SUBROUTINE write_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(IN) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(IN) :: qstars(:,:,:,:) + INTEGER,INTENT(IN) :: ksym(:),jsym(:) + REAL,INTENT(IN),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(IN),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + + TYPE(t_data_MPI),POINTER :: d + INTEGER:: pe,slot + + CALL priv_find_data(id,d) + pe=d%pe_basis(nk,jspin) + slot=d%slot_basis(nk,jspin) + + CALL priv_put_data(pe,slot,RESHAPE(qal,(/SIZE(qal)/)),d%qal_handle) + CALL priv_put_data(pe,slot,RESHAPE(qvac,(/SIZE(qvac)/)),d%qvac_handle) + CALL priv_put_data(pe,slot,RESHAPE(qis,(/SIZE(qis)/)),d%qis_handle) + CALL priv_put_data(pe,slot,RESHAPE(qvlay,(/SIZE(qvlay)/)),d%qvlay_handle) + CALL priv_put_data(pe,slot,RESHAPE(qstars,(/SIZE(qstars)/)),d%qstars_handle) + CALL priv_put_data(pe,slot,RESHAPE(ksym,(/SIZE(ksym)/)),d%ksym_handle) + CALL priv_put_data(pe,slot,RESHAPE(jsym,(/SIZE(jsym)/)),d%jsym_handle) + IF (d%l_mcd.AND.PRESENT(mcd)) CALL priv_put_data(pe,slot,RESHAPE(mcd,(/SIZE(mcd)/)),d%mcd_handle) + IF (d%l_orb.AND.PRESENT(qintsl)) THEN + CALL priv_put_data(pe,slot,RESHAPE(qintsl,(/SIZE(qintsl)/)),d%qintsl_handle) + CALL priv_put_data(pe,slot,RESHAPE(qmtsl,(/SIZE(qmtsl)/)),d%qmtsl_handle) + CALL priv_put_data(pe,slot,RESHAPE(qmtp,(/SIZE(qmtp)/)),d%qmtp_handle) + CALL priv_put_data(pe,slot,RESHAPE(orbcomp,(/SIZE(orbcomp)/)),d%orbcomp_handle) + ENDIF + END SUBROUTINE write_dos + + SUBROUTINE read_dos(id,nk,jspin,qal,qvac,qis,qvlay,qstars,ksym,jsym,mcd,qintsl,qmtsl,qmtp,orbcomp) + IMPLICIT NONE + INTEGER, INTENT(IN) :: id,nk,jspin + REAL,INTENT(out) :: qal(:,:,:),qvac(:,:),qis(:),qvlay(:,:,:) + COMPLEX,INTENT(out) :: qstars(:,:,:,:) + INTEGER,INTENT(out) :: ksym(:),jsym(:) + REAL,INTENT(out),OPTIONAL :: mcd(:,:,:) + REAL,INTENT(out),OPTIONAL :: qintsl(:,:),qmtsl(:,:),qmtp(:,:),orbcomp(:,:,:) + + TYPE(t_data_MPI),POINTER :: d + INTEGER:: pe,slot + + CALL priv_find_data(id,d) + pe=d%pe_basis(nk,jspin) + slot=d%slot_basis(nk,jspin) + + CALL priv_get_data(pe,slot,qal,SIZE(qal),d%qal_handle) + CALL priv_get_data(pe,slot,qvac,SIZE(qvac),d%qvac_handle) + CALL priv_get_data(pe,slot,qis,SIZE(qis),d%qis_handle) + CALL priv_get_data(pe,slot,qvlay,SIZE(qvlay),d%qvlay_handle) + CALL priv_get_data(pe,slot,qstars,SIZE(qstars),d%qstars_handle) + CALL priv_get_data(pe,slot,ksym,SIZE(ksym),d%ksym_handle) + CALL priv_get_data(pe,slot,jsym,SIZE(jsym),d%jsym_handle) + IF (d%l_mcd.AND.PRESENT(mcd)) CALL priv_get_data(pe,slot,mcd,SIZE(mcd),d%mcd_handle) + IF (d%l_orb.AND.PRESENT(qintsl)) THEN + CALL priv_get_data(pe,slot,qintsl,SIZE(qintsl),d%qintsl_handle) + CALL priv_get_data(pe,slot,qmtsl,SIZE(qmtsl),d%qmtsl_handle) + CALL priv_get_data(pe,slot,qmtp,SIZE(qmtp),d%qmtp_handle) + CALL priv_get_data(pe,slot,orbcomp,SIZE(orbcomp),d%orbcomp_handle) + ENDIF + END SUBROUTINE read_dos + #ifdef CPP_MPI - SUBROUTINE create_maps(d,isize,nkpts,jspins,neig,n_size) - IMPLICIT NONE - TYPE(t_data_MPI),INTENT(INOUT):: d - INTEGER,INTENT(IN):: isize,nkpts,jspins,neig,n_size - - INTEGER:: nk,j,n1,n2,n,pe,n_members - INTEGER::used(0:isize) - - ALLOCATE(d%pe_basis(nkpts,jspins),d%slot_basis(nkpts,jspins)) - ALLOCATE(d%pe_ev(nkpts,jspins,neig),d%slot_ev(nkpts,jspins,neig)) - - !basis contains a total of nkpts*jspins entries - d%pe_basis=-1 - d%pe_ev=-1 - used=0 - n_members=isize/n_size !no of k-points in parallel - DO j=1,jspins - DO nk=1,nkpts - n1=nk+(j-1)*nkpts-1 - pe=MOD(n1,n_members)*n_size - d%pe_basis(nk,j)=pe - d%slot_basis(nk,j)=used(pe) - used(pe)=used(pe)+1 - ENDDO - ENDDO - used=0 - DO n=1,neig - DO j=1,jspins - DO nk=1,nkpts - n1=nk+(j-1)*nkpts-1 - !eigenvectors have more entries - pe=MOD(n1,n_members)*n_size+mod(n,n_size) - d%pe_ev(nk,j,n)=pe - d%slot_ev(nk,j,n)=used(pe) - used(pe)=used(pe)+1 - ENDDO - ENDDO - ENDDO + SUBROUTINE create_maps(d,isize,nkpts,jspins,neig,n_size) + IMPLICIT NONE + TYPE(t_data_MPI),INTENT(INOUT):: d + INTEGER,INTENT(IN):: isize,nkpts,jspins,neig,n_size + + INTEGER:: nk,j,n1,n2,n,pe,n_members + INTEGER::used(0:isize) + + ALLOCATE(d%pe_basis(nkpts,jspins),d%slot_basis(nkpts,jspins)) + ALLOCATE(d%pe_ev(nkpts,jspins,neig),d%slot_ev(nkpts,jspins,neig)) + + !basis contains a total of nkpts*jspins entries + d%pe_basis=-1 + d%pe_ev=-1 + used=0 + n_members=isize/n_size !no of k-points in parallel + DO j=1,jspins + DO nk=1,nkpts + n1=nk+(j-1)*nkpts-1 + pe=MOD(n1,n_members)*n_size + d%pe_basis(nk,j)=pe + d%slot_basis(nk,j)=used(pe) + used(pe)=used(pe)+1 + ENDDO + ENDDO + used=0 + DO n=1,neig + DO j=1,jspins + DO nk=1,nkpts + n1=nk+(j-1)*nkpts-1 + !eigenvectors have more entries + pe=MOD(n1,n_members)*n_size+MOD(n,n_size) + d%pe_ev(nk,j,n)=pe + d%slot_ev(nk,j,n)=used(pe) + used(pe)=used(pe)+1 + ENDDO + ENDDO + ENDDO - END SUBROUTINE + END SUBROUTINE create_maps #endif