diff --git a/diag.f90 b/diag.f90 index 90041cb..a070a1c 100644 --- a/diag.f90 +++ b/diag.f90 @@ -15,7 +15,7 @@ module diag private verbose ! ! - public diag_tridiag,diag_tridiag_pack,diag_dsyev_i8,diag_dgelss,diag_syev_ilp,diag_syev_i8,diag_dseupd,diag_dseupd_p + public diag_tridiag,diag_tridiag_pack,diag_dsyev_i8,diag_dgelss,diag_syev_ilp,diag_syev_i8 public dseupd_omp_arpack,diag_propack,daprod ! integer,parameter:: verbose = 4 @@ -2558,456 +2558,6 @@ end subroutine matvec_p - subroutine diag_dseupd_p(n,bterm,nroots,factor,maxitr_,tol,h,e) - - integer,intent(in) :: n,bterm(n,2),maxitr_ - integer,intent(inout) :: nroots - double precision,intent(in) :: factor - double precision,intent(in) :: tol - double precision,intent(inout) :: h(:,:) - double precision,intent(out) :: e(nroots) - - double precision,allocatable :: v(:,:),workl(:),workd(:),d(:,:),resid(:) - ! - logical,allocatable :: select(:) - integer(ik) :: iparam(11), ipntr(11),ldv,iter - ! - integer,parameter :: maxnprocs=1024 - ! - integer(ik) :: kstart(0:maxnprocs),iproc,dx,i,jproc - ! - character(len=1) :: bmat - character(len=2) :: which - character(len=cl) :: blacs_or_mpi = 'NONE' - ! - integer(ik) :: ido, nev, ncv, lworkl, info, ierr, j, & - mode, ishfts, alloc, maxitr - ! - logical rvec - double precision :: sigma - - ! - !double precision,external :: dnrm2 - ! - -! -!----------------------------------------------------------------------- -! - integer logfil, ndigit, mgetv0,& - msaupd, msaup2, msaitr, mseigt, msapps, msgets, mseupd,& - mnaupd, mnaup2, mnaitr, mneigh, mnapps, mngets, mneupd,& - mcaupd, mcaup2, mcaitr, mceigh, mcapps, mcgets, mceupd - common /debug/ & - logfil, ndigit, mgetv0,& - msaupd, msaup2, msaitr, mseigt, msapps, msgets, mseupd,& - mnaupd, mnaup2, mnaitr, mneigh, mnapps, mngets, mneupd,& - mcaupd, mcaup2, mcaitr, mceigh, mcapps, mcgets, mceupd - - real t0, t1, t2, t3, t4, t5 - save t0, t1, t2, t3, t4, t5 -! - integer nopx, nbx, nrorth, nitref, nrstrt - real tsaupd, tsaup2, tsaitr, tseigt, tsgets, tsapps, tsconv,& - tnaupd, tnaup2, tnaitr, tneigh, tngets, tnapps, tnconv,& - tcaupd, tcaup2, tcaitr, tceigh, tcgets, tcapps, tcconv,& - tmvopx, tmvbx, tgetv0, titref, trvec - common /timing/ & - nopx, nbx, nrorth, nitref, nrstrt,& - tsaupd, tsaup2, tsaitr, tseigt, tsgets, tsapps, tsconv,& - tnaupd, tnaup2, tnaitr, tneigh, tngets, tnapps, tnconv,& - tcaupd, tcaup2, tcaitr, tceigh, tcgets, tcapps, tcconv,& - tmvopx, tmvbx, tgetv0, titref, trvec - - - -! %-----------------% -! | BLACS INTERFACE | -! %-----------------% -! - integer comm, iam, nprocs, nloc, & - nprow, npcol, myprow, mypcol -! - !external BLACS_PINFO, BLACS_SETUP, BLACS_GET, & - ! BLACS_GRIDINIT, BLACS_GRIDINFO - !integer,parameter :: MPI_COMM_WORLD=0 -! - -! %---------------% -! | MPI INTERFACE | -! %---------------% - - integer myid, rc -! include 'mpif.h' - - -! -! %----------------------------------------------% -! | Local Buffers needed for BLACS communication | -! %----------------------------------------------% -! - Double precision,allocatable :: mv_buf(:) -! -! %------------% -! | Parameters | -! %------------% -! - Double precision zero - parameter (zero = 0.0) -! -! %-----------------------------% -! | BLAS & LAPACK routines used | -! %-----------------------------% -! - Double precision pdnorm2 - external pdnorm2, daxpy -! -! %---------------------% -! | Intrinsic Functions | -! %---------------------% -! - intrinsic abs -! -! %-----------------------% -! | Executable Statements | -! %-----------------------% -! - - write(out,"('Start PARPACK-diagonalization')") - - -#if (blacs_ > 0) - write(out,"('BLAS-PINFO-start')") - call BLACS_PINFO( iam, nprocs ) - print *,nprocs - blacs_or_mpi = 'BLACS' -#endif - ! - write(out,"('BLAS-PINFO-done')") - ! - !call BLACS_PINFO( iam, nprocs ) - !print *,nprocs - -#if (mpi_ > 0) - call MPI_INIT( ierr ) - comm = MPI_COMM_WORLD - call MPI_COMM_RANK( comm, myid, ierr ) - call MPI_COMM_SIZE( comm, nprocs, ierr ) - ! - print *,comm,myid,nprocs - ! - if (trim(blacs_or_mpi)=='BLACS') then - write(out,"('diag_dseupd_p: IT IS ILLEGAL TO USE MPI AND BLACS AT THE SAME TIME')") - stop 'diag_dseupd_p: IT IS ILLEGAL TO USE MPI AND BLACS AT THE SAME TIME' - endif - ! - blacs_or_mpi = 'MPI' - ! -#endif - ! - - -! -! If in PVM, create virtual machine if it doesn't exist -! - if (nprocs .lt. 1) then - nprocs = 1 -#if (blacs_ > 0) - call BLACS_SETUP( iam, nprocs ) -#endif - ! - print *,nprocs - ! - endif - if (nprocs >maxnprocs) stop 'nprocs > maxnprocs' - ! - ! Set up processors in 1D Grid - ! - nprow = nprocs - npcol = 1 - ! - do iproc = 0,nprocs - ! - kstart(iproc) = iproc*(n/nprocs+1)+1 - if (iproc>=mod(n,nprocs)) kstart(iproc) = iproc*(n/nprocs)+1+mod(n,nprocs) - print *,iproc,kstart(iproc) - ! - enddo - ! - if (verbose>=2) call TimerStart('diag_dseupd_p: diagonalization') - ! - nev = nroots - ! - if (nroots==0) nev=max(100,n) - ! - ncv = factor*nev - ! - if (ncv==nev) ncv = int(21*nev/10)+1 - ! - ncv = min(n,ncv) - ! - bmat = 'I' - which = 'SM' - ! - ! The work array WORKL is used in DSAUPD as workspace. - ! The parameter TOL determines the stopping criterion. - ! If TOL<=0, machine precision is used. The variable IDO is used for - ! reverse communication and is initially set to 0. - ! Setting INFO=0 indicates that a random vector is generated in DSAUPD - ! to start the Arnoldi iteration. - ! - - lworkl = ncv*(ncv+10) - ldv = n - ! - info = 0 - ido = 0 - ! - allocate(v(ldv,ncv), workl(lworkl),workd(3*n),d(ncv,2),resid(n),select(ncv),mv_buf(n),stat=alloc) - call ArrayStart('diag_dseupd_p',alloc,size(v),kind(v)) - call ArrayStart('diag_dseupd_p',alloc,size(workl),kind(workl)) - call ArrayStart('diag_dseupd_p',alloc,size(workd),kind(workd)) - call ArrayStart('diag_dseupd_p',alloc,size(d),kind(d)) - call ArrayStart('diag_dseupd_p',alloc,size(resid),kind(resid)) - call ArrayStart('diag_dseupd_p',alloc,size(select),kind(select)) - call ArrayStart('diag_dseupd_p',alloc,size(mv_buf),kind(mv_buf)) - - ! - ! Estimate block overlapping: - ! - dx = 0 - ! - do iproc = 0,nprocs-1 - ! - do i = kstart(iproc),kstart(iproc+1)-1 - ! - do jproc = 0,iproc-1 - ! - if (bterm(i,1)kstart(jproc)) dx = max(dx,jproc-iproc) - ! - enddo - ! - enddo - ! - enddo - ! -! -! Get default system context, and define grid -! - ! - myprow = 1 ; mypcol = 1 ; myid = 1 - ! -#if (blacs_ > 0) - call BLACS_GET( 0, 0, comm ) - call BLACS_GRIDINIT( comm, 'Row', nprow, npcol ) - call BLACS_GRIDINFO( comm, nprow, npcol, myprow, mypcol ) - ! - write(out,"(' myprow, nprow, mypcol, npcol = ',4i8)") myprow, nprow, mypcol, npcol - ! -#endif - ! - if (verbose>=2.and.trim(blacs_or_mpi)=='BLACS') then - write(out,"('myprow, nprow, mypcol, npcol, nprocs = ',5i8)") myprow, nprow, mypcol, npcol, nprocs - endif - ! - if (verbose>=2.and.trim(blacs_or_mpi)=='MPI') write(out,"('myid,nproc = ',2i8)") myid,nprocs -! ! -! If I'm not in grid, go to end of program -! - if ( trim(blacs_or_mpi)=='BLACS'.and.(myprow .ge. nprow) .or. (mypcol .ge. npcol) ) then - write(out,"('not in grid, myprow, nprow, mypcol, npcol = ',4i8)") myprow, nprow, mypcol, npcol - return - endif -! - ndigit = -3 - logfil = 6 - msaupd = 1 - - - !--------------------------------------! - ! Set up distribution of data to nodes ! - !--------------------------------------! - - nloc = (n/nprocs) - if ( mod(n,nprocs)>myprow ) nloc = nloc + 1 - ! - ! - !allocate(mv_buf(n),stat=alloc) - !call ArrayStart('diag_dseupd_p',alloc,size(mv_buf),kind(mv_buf)) - ! - iparam = 0 - ipntr = 0 - ! - ! This routone uses exact shifts with respect to - ! the current Hessenberg matrix (IPARAM(1) = 1). - ! IPARAM(3) specifies the maximum number of Arnoldi - ! iterations allowed. Mode 1 of DSAUPD is used - ! (IPARAM(7) = 1). - ! - ishfts = 1 - maxitr = maxitr_*nev - mode = 1 - ! - iparam(1) = ishfts - iparam(3) = maxitr - iparam(7) = mode - ! - ! - ! M A I N L O O P (Reverse communication) - ! - !if (verbose>=4) write(out,"(/'Arpack: eigenvalues computed:')") - ! - if (verbose>=3) write(out,"(/'Arpack: n = ',i8,' nev = ',i8,' maxitr = ',i8,' tol = ',e12.5)") n,nev,maxitr,tol - ! - iter = 0 - ! - do while(ido<=1) - ! - iter = iter + 1 - ! - if (verbose>=3.and.mod(iter,10)==0) write(out,"(' iter = ',i8)") iter - ! - ! Repeatedly call the routine DSAUPD and take - ! actions indicated by parameter IDO until - ! either convergence is indicated or maxitr - ! has been exceeded. - ! - -#if (blacs_ > 0 || mpi_ > 0) - ! - call pdsaupd ( comm, ido, bmat, nloc, which, nev, tol, resid, & - ncv, v, ldv, iparam, ipntr, workd, workl, & - lworkl, info ) - ! -#else - ! - call dsaupd ( ido, bmat, n, which, nev, tol, resid, & - ncv, v, ldv, iparam, ipntr, workd, workl, & - lworkl, info ) - ! -#endif - ! - !if (verbose>=4.and.iparam(5)>0) then - ! ! - ! write(out,"(i8)") iparam(5) - ! ! - !endif - ! - if (ido==-1 .or. ido==1) then - ! - ! Perform matrix vector multiplication - ! y <--- OP*x - ! - call matvec_p(comm,n,nloc,nprocs,kstart(0:nprocs),dx,bterm,h,mv_buf,workd(ipntr(1)),workd(ipntr(2))) - ! - endif - ! - enddo - ! - ! Either we have convergence or there is an error. - ! - if ( info < 0 ) then - write(out,"(/'Error with _saupd, info = ',i8)") info - stop 'Error with _saupd' - endif - ! - ! - ! No fatal errors occurred. Post-Process using DSEUPD. - ! Computed eigenvalues may be extracted. - ! - ! Eigenvectors may also be computed now if desired. (indicated by rvec = .true.) - ! - rvec = .true. - ! -#if (blacs_ > 0) - ! - if (verbose>=5) write(out,"(/'Arpack: dseupd')") - ! - call pdseupd ( comm, rvec, 'All', select, d, v, ldv, sigma, & - bmat, n, which, nev, tol, resid, ncv, v, ldv, & - iparam, ipntr, workd, workl, lworkl, ierr ) - - - if (verbose>=5) write(out,"(/'Arpack: done!')") - - ! -#else - ! - write(out,"(/'Arpack was not activated yet. Please uncomment dsaupd and dseupd bellow')") - stop 'Arpack was not activated' - ! -#endif - ! - if ( ierr < 0 ) then - write(out,"(/'Error with_seupd, info = ',i8)") ierr - stop 'Error with _saupd' - endif - ! - nroots = iparam(5) - ! - e(1:nroots) = d(1:nroots,1) - ! - !$omp parallel do private(j) shared(h) schedule(dynamic) - do j=1,n - h(j,1:nroots) = v(j,1:nroots) - enddo - !$omp end parallel do - ! - ! Eigenvalues are returned in the first column - ! of the two dimensional array D and the - ! corresponding eigenvectors are returned in - ! the first NEV columns of the two dimensional - ! array V if requested. Otherwise, an - ! orthogonal basis for the invariant subspace - ! corresponding to the eigenvalues in D is - ! returned in V. - ! - !do j=1, nconv - ! - ! Compute the residual norm || A*x - lambda*x || - ! for the NCONV accurately computed eigenvalues and - ! eigenvectors. (iparam(5) indicates how many are - ! accurate to the requested tolerance) - ! - !call matvec(n,bterm,v(1,j), ax) - ! - !call daxpy(n, -d(j,1), v(1,j), 1, ax, 1) - !d(j,2) = dnrm2(n, ax, 1) - !d(j,2) = d(j,2) / abs(d(j,1)) - ! - !enddo - -! -! %---------------------------% -! | Done with program pdsdrv1.| -! %---------------------------% -! - 9000 continue -! -#if (blacs_ > 0) - call BLACS_GRIDEXIT ( comm ) - call BLACS_EXIT(0) -#endif -#if (mpi_ > 0) - call MPI_FINALIZE(rc) -#endif - - - deallocate(v,workl,workd,d,resid,select,mv_buf) - ! - call ArrayStop('diag_dseupd_p') - ! - if (verbose>=2) call TimerStop('diag_dseupd_p: diagonalization') - ! - end subroutine diag_dseupd_p - - - subroutine dseupd_omp_arpack(n,bterm,nroots,factor,maxitr_,tol,h,e) integer,intent(in) :: n,bterm(n,2),maxitr_ @@ -3270,9 +2820,6 @@ subroutine dseupd_omp_arpack(n,bterm,nroots,factor,maxitr_,tol,h,e) end subroutine dseupd_omp_arpack - - - subroutine diag_propack(n,bterm,nroots,factor,maxiter,iverbose,tol,h,e) ! integer,intent(in) :: n,bterm(n,2),maxiter,iverbose diff --git a/lapack.f90 b/lapack.f90 index b52b0e0..7aad724 100644 --- a/lapack.f90 +++ b/lapack.f90 @@ -1657,433 +1657,6 @@ end subroutine matvec_p - subroutine dseupd_p_arpack(n,bterm,nroots,factor,maxitr_,tol,h,e) - - integer,intent(in) :: n,bterm(n,2),maxitr_ - integer,intent(inout) :: nroots - double precision,intent(in) :: factor - double precision,intent(in) :: tol - double precision,intent(inout) :: h(:,:) - double precision,intent(out) :: e(nroots) - - double precision,allocatable :: v(:,:),workl(:),workd(:),d(:,:),resid(:) - ! - logical,allocatable :: select(:) - integer(ik) :: iparam(11), ipntr(11),ldv,iter - ! - integer,parameter :: maxnprocs=256 - ! - integer(ik) :: kstart(0:maxnprocs),iproc,dx,i,jproc - ! - character(len=1) :: bmat - character(len=2) :: which - character(len=cl) :: blacs_or_mpi = 'NONE' - ! - integer(ik) :: ido, nev, ncv, lworkl, info, ierr, j, & - mode, ishfts, alloc, maxitr - ! - logical rvec - double precision :: sigma - - ! - !double precision,external :: dnrm2 - ! - -! -!----------------------------------------------------------------------- -! - integer logfil, ndigit, mgetv0,& - msaupd, msaup2, msaitr, mseigt, msapps, msgets, mseupd,& - mnaupd, mnaup2, mnaitr, mneigh, mnapps, mngets, mneupd,& - mcaupd, mcaup2, mcaitr, mceigh, mcapps, mcgets, mceupd - common /debug/ & - logfil, ndigit, mgetv0,& - msaupd, msaup2, msaitr, mseigt, msapps, msgets, mseupd,& - mnaupd, mnaup2, mnaitr, mneigh, mnapps, mngets, mneupd,& - mcaupd, mcaup2, mcaitr, mceigh, mcapps, mcgets, mceupd - - real t0, t1, t2, t3, t4, t5 - save t0, t1, t2, t3, t4, t5 -! - integer nopx, nbx, nrorth, nitref, nrstrt - real tsaupd, tsaup2, tsaitr, tseigt, tsgets, tsapps, tsconv,& - tnaupd, tnaup2, tnaitr, tneigh, tngets, tnapps, tnconv,& - tcaupd, tcaup2, tcaitr, tceigh, tcgets, tcapps, tcconv,& - tmvopx, tmvbx, tgetv0, titref, trvec - common /timing/ & - nopx, nbx, nrorth, nitref, nrstrt,& - tsaupd, tsaup2, tsaitr, tseigt, tsgets, tsapps, tsconv,& - tnaupd, tnaup2, tnaitr, tneigh, tngets, tnapps, tnconv,& - tcaupd, tcaup2, tcaitr, tceigh, tcgets, tcapps, tcconv,& - tmvopx, tmvbx, tgetv0, titref, trvec - - - -! %-----------------% -! | BLACS INTERFACE | -! %-----------------% -! - integer comm, iam, nprocs, nloc, & - nprow, npcol, myprow, mypcol -! - !external BLACS_PINFO, BLACS_SETUP, BLACS_GET, & - ! BLACS_GRIDINIT, BLACS_GRIDINFO - integer,parameter :: MPI_COMM_WORLD=0 -! - -! %---------------% -! | MPI INTERFACE | -! %---------------% - - integer myid, rc - - -! -! %----------------------------------------------% -! | Local Buffers needed for BLACS communication | -! %----------------------------------------------% -! - Double precision,allocatable :: mv_buf(:) -! -! %------------% -! | Parameters | -! %------------% -! - Double precision zero - parameter (zero = 0.0) -! -! %-----------------------------% -! | BLAS & LAPACK routines used | -! %-----------------------------% -! - Double precision pdnorm2 - external pdnorm2, daxpy -! -! %---------------------% -! | Intrinsic Functions | -! %---------------------% -! - intrinsic abs -! -! %-----------------------% -! | Executable Statements | -! %-----------------------% -! - -#if (blacs_ > 0) - call BLACS_PINFO( iam, nprocs ) - blacs_or_mpi = 'BLACS' -#endif - -#if (mpi_ > 0) - call MPI_INIT( ierr ) - comm = MPI_COMM_WORLD - call MPI_COMM_RANK( comm, myid, ierr ) - call MPI_COMM_SIZE( comm, nprocs, ierr ) - ! - print *,comm,myid,nprocs - ! - if (trim(blacs_or_mpi)=='BLACS') then - write(out,"('dseupd_p_arpack: IT IS ILLEGAL TO USE MPI AND BLACS AT THE SAME TIME')") - stop 'dseupd_p_arpack: IT IS ILLEGAL TO USE MPI AND BLACS AT THE SAME TIME' - endif - ! - blacs_or_mpi = 'MPI' - ! - ! -#endif - - -! -! If in PVM, create virtual machine if it doesn't exist -! - if (nprocs .lt. 1) then - nprocs = 1 -#if (blacs_ > 0) - call BLACS_SETUP( iam, nprocs ) -#endif - endif - if (nprocs >maxnprocs) stop 'nprocs > maxnprocs' - ! - ! Set up processors in 1D Grid - ! - nprow = nprocs - npcol = 1 - ! - do iproc = 0,nprocs - ! - kstart(iproc) = iproc*(n/nprocs+1)+1 - if (iproc>=mod(n,nprocs)) kstart(iproc) = iproc*(n/nprocs)+1+mod(n,nprocs) - print *,iproc,kstart(iproc) - ! - enddo - ! - if (verbose>=2) call TimerStart('dseupd_p_arpack: diagonalization') - ! - nev = nroots - ! - if (nroots==0) nev=max(100,n) - ! - ncv = factor*nev - ! - if (ncv==nev) ncv = int(21*nev/10)+1 - ! - ncv = min(n,ncv) - ! - bmat = 'I' - which = 'SM' - ! - ! The work array WORKL is used in DSAUPD as workspace. - ! The parameter TOL determines the stopping criterion. - ! If TOL<=0, machine precision is used. The variable IDO is used for - ! reverse communication and is initially set to 0. - ! Setting INFO=0 indicates that a random vector is generated in DSAUPD - ! to start the Arnoldi iteration. - ! - - lworkl = ncv*(ncv+10) - ldv = n - ! - info = 0 - ido = 0 - ! - allocate(v(ldv,ncv), workl(lworkl),workd(3*n),d(ncv,2),resid(n),select(ncv),mv_buf(n),stat=alloc) - call ArrayStart('dseupd_p_arpack',alloc,size(v),kind(v)) - call ArrayStart('dseupd_p_arpack',alloc,size(workl),kind(workl)) - call ArrayStart('dseupd_p_arpack',alloc,size(workd),kind(workd)) - call ArrayStart('dseupd_p_arpack',alloc,size(d),kind(d)) - call ArrayStart('dseupd_p_arpack',alloc,size(resid),kind(resid)) - call ArrayStart('dseupd_p_arpack',alloc,size(select),kind(select)) - call ArrayStart('dseupd_p_arpack',alloc,size(mv_buf),kind(mv_buf)) - - ! - ! Estimate block overlapping: - ! - dx = 0 - ! - do iproc = 0,nprocs-1 - ! - do i = kstart(iproc),kstart(iproc+1)-1 - ! - do jproc = 0,iproc-1 - ! - if (bterm(i,1)kstart(jproc)) dx = max(dx,jproc-iproc) - ! - enddo - ! - enddo - ! - enddo - ! -! -! Get default system context, and define grid -! - ! - myprow = 1 ; mypcol = 1 ; myid = 1 -#if (blacs_ > 0) - call BLACS_GET( 0, 0, comm ) - call BLACS_GRIDINIT( comm, 'Row', nprow, npcol ) - call BLACS_GRIDINFO( comm, nprow, npcol, myprow, mypcol ) -#endif - ! - if (verbose>=2.and.trim(blacs_or_mpi)=='BLACS') write(out,"('myprow, nprow, mypcol, npcol, nprocs = ',5i8)") & - myprow, nprow, mypcol, npcol, nprocs - if (verbose>=2.and.trim(blacs_or_mpi)=='MPI') write(out,"('myid,nproc = ',2i8)") myid,nprocs -! ! -! If I'm not in grid, go to end of program -! - if ( trim(blacs_or_mpi)=='BLACS'.and.(myprow .ge. nprow) .or. (mypcol .ge. npcol) ) then - write(out,"('not in grid, myprow, nprow, mypcol, npcol = ',4i8)") myprow, nprow, mypcol, npcol - return - endif -! - ndigit = -3 - logfil = 6 - msaupd = 1 - - - !--------------------------------------! - ! Set up distribution of data to nodes ! - !--------------------------------------! - - nloc = (n/nprocs) - if ( mod(n,nprocs)>myprow ) nloc = nloc + 1 - ! - ! - !allocate(mv_buf(n),stat=alloc) - !call ArrayStart('dseupd_p_arpack',alloc,size(mv_buf),kind(mv_buf)) - ! - iparam = 0 - ipntr = 0 - ! - ! This routone uses exact shifts with respect to - ! the current Hessenberg matrix (IPARAM(1) = 1). - ! IPARAM(3) specifies the maximum number of Arnoldi - ! iterations allowed. Mode 1 of DSAUPD is used - ! (IPARAM(7) = 1). - ! - ishfts = 1 - maxitr = maxitr_*nev - mode = 1 - ! - iparam(1) = ishfts - iparam(3) = maxitr - iparam(7) = mode - ! - ! - ! M A I N L O O P (Reverse communication) - ! - !if (verbose>=4) write(out,"(/'Arpack: eigenvalues computed:')") - ! - if (verbose>=3) write(out,"(/'Arpack: n = ',i8,' nev = ',i8,' maxitr = ',i8,' tol = ',e12.5)") n,nev,maxitr,tol - ! - iter = 0 - ! - do while(ido<=1) - ! - iter = iter + 1 - ! - if (verbose>=3.and.mod(iter,10)==0) write(out,"(' iter = ',i8)") iter - ! - ! Repeatedly call the routine DSAUPD and take - ! actions indicated by parameter IDO until - ! either convergence is indicated or maxitr - ! has been exceeded. - ! - -#if (blacs_ > 0 || mpi_ > 0) - ! - call pdsaupd ( comm, ido, bmat, nloc, which, nev, tol, resid, & - ncv, v, ldv, iparam, ipntr, workd, workl, & - lworkl, info ) - ! -#else - ! - call dsaupd ( ido, bmat, n, which, nev, tol, resid, & - ncv, v, ldv, iparam, ipntr, workd, workl, & - lworkl, info ) - ! -#endif - ! - !if (verbose>=4.and.iparam(5)>0) then - ! ! - ! write(out,"(i8)") iparam(5) - ! ! - !endif - ! - if (ido==-1 .or. ido==1) then - ! - ! Perform matrix vector multiplication - ! y <--- OP*x - ! - call matvec_p(comm,n,nloc,nprocs,kstart(0:nprocs),dx,bterm,h,mv_buf,workd(ipntr(1)),workd(ipntr(2))) - ! - endif - ! - enddo - ! - ! Either we have convergence or there is an error. - ! - if ( info < 0 ) then - write(out,"(/'Error with _saupd, info = ',i8)") info - stop 'Error with _saupd' - endif - ! - ! - ! No fatal errors occurred. Post-Process using DSEUPD. - ! Computed eigenvalues may be extracted. - ! - ! Eigenvectors may also be computed now if desired. (indicated by rvec = .true.) - ! - rvec = .true. - ! -#if (blacs_ > 0) - ! - if (verbose>=5) write(out,"(/'Arpack: dseupd')") - ! - call pdseupd ( comm, rvec, 'All', select, d, v, ldv, sigma, & - bmat, n, which, nev, tol, resid, ncv, v, ldv, & - iparam, ipntr, workd, workl, lworkl, ierr ) - - - if (verbose>=5) write(out,"(/'Arpack: done!')") - - ! -#else - ! - write(out,"(/'Arpack was not activated yet. Please uncomment dsaupd and dseupd bellow')") - stop 'Arpack was not activated' - ! -#endif - ! - if ( ierr < 0 ) then - write(out,"(/'Error with_seupd, info = ',i8)") ierr - stop 'Error with _saupd' - endif - ! - nroots = iparam(5) - ! - e(1:nroots) = d(1:nroots,1) - ! - !$omp parallel do private(j) shared(h) schedule(dynamic) - do j=1,n - h(j,1:nroots) = v(j,1:nroots) - enddo - !$omp end parallel do - ! - ! Eigenvalues are returned in the first column - ! of the two dimensional array D and the - ! corresponding eigenvectors are returned in - ! the first NEV columns of the two dimensional - ! array V if requested. Otherwise, an - ! orthogonal basis for the invariant subspace - ! corresponding to the eigenvalues in D is - ! returned in V. - ! - !do j=1, nconv - ! - ! Compute the residual norm || A*x - lambda*x || - ! for the NCONV accurately computed eigenvalues and - ! eigenvectors. (iparam(5) indicates how many are - ! accurate to the requested tolerance) - ! - !call matvec(n,bterm,v(1,j), ax) - ! - !call daxpy(n, -d(j,1), v(1,j), 1, ax, 1) - !d(j,2) = dnrm2(n, ax, 1) - !d(j,2) = d(j,2) / abs(d(j,1)) - ! - !enddo - -! -! %---------------------------% -! | Done with program pdsdrv1.| -! %---------------------------% -! - 9000 continue -! -#if (blacs_ > 0) - call BLACS_GRIDEXIT ( comm ) - call BLACS_EXIT(0) -#endif -#if (mpi_ > 0) - call MPI_FINALIZE(rc) -#endif - - - deallocate(v,workl,workd,d,resid,select,mv_buf) - ! - call ArrayStop('dseupd_p_arpack') - ! - if (verbose>=2) call TimerStop('dseupd_p_arpack: diagonalization') - ! - end subroutine dseupd_p_arpack !subroutine BLACS_GRIDINFO( comm, nprow, npcol, myprow, mypcol ) diff --git a/makefile b/makefile index c7bb7f3..710f203 100644 --- a/makefile +++ b/makefile @@ -13,17 +13,15 @@ pot_user = pot_ch4 PLAT = #PLAT = _2205_i17 FOR = mpif90 -FFLAGS = -qopenmp -xHost -O3 -ip -g3 -traceback -#FFLAGS = -fopenmp -ffree-line-length-none -march=native -O3 -fcray-pointer -g3 +#FFLAGS = -qopenmp -xHost -O3 -ip -g3 -traceback -DTROVE_USE_MPI_ +FFLAGS = -fopenmp -ffree-line-length-512 -march=native -O0 -fcray-pointer -g -fallow-argument-mismatch -fbacktrace -DTROVE_USE_MPI_ +#LAPACK = -mkl=parallel -lmkl_scalapack_lp64 -lmkl_blacs_intelmpi_lp64 +LAPACK = -L${MKLROOT}/lib/intel64 -Wl,--no-as-needed -lmkl_gf_lp64 -lmkl_gnu_thread -lmkl_blacs_intelmpi_lp64 -lmkl_scalapack_lp64 -lmkl_core -lgomp -lpthread -lm -ldl -#ARPACK = ~/libraries/ARPACK/libarpack_omp_64.a +ARPACK = -larpack -LAPACK = -mkl=parallel -lmkl_scalapack_lp64 -lmkl_blacs_intelmpi_lp64 -#LAPACK = -L${MKLROOT}/lib/intel64 -Wl,--no-as-needed -lmkl_gf_lp64 -lmkl_sequential -lmkl_core -lpthread -lm -ldl - - -LIB = $(LAPACK) +LIB = $(LAPACK) $(ARPACK) %.o : %.f90 $(FOR) -cpp -c $(FFLAGS) $< diff --git a/mpi_aux.f90 b/mpi_aux.f90 index 457ec8d..140d5e2 100644 --- a/mpi_aux.f90 +++ b/mpi_aux.f90 @@ -1,5 +1,7 @@ module mpi_aux +#ifdef TROVE_USE_MPI_ use mpi_f08 +#endif use timer use accuracy implicit none @@ -10,14 +12,43 @@ module mpi_aux public send_or_recv, comm_size, mpi_rank public co_startdim, co_enddim - public blacs_size, blacs_rank, blacs_ctxt + public blacs_size, blacs_rank public nprow,npcol,myprow,mypcol public mpi_real_size, mpi_int_size +#ifndef TROVE_USE_MPI_ + public MPI_Datatype, MPI_File, MPI_Status, MPI_Request +#endif + interface co_sum module procedure :: co_sum_double end interface +#ifndef TROVE_USE_MPI_ + ! + ! When not using MPI, create a dummy type to maintain the interface of this module, + ! specifically the last argument of co_block_type_init. + ! The value passed will not be accessed. + ! + type MPI_Datatype + integer :: dummy = 0 + end type MPI_Datatype + + type MPI_File + integer :: dummy = 0 + end type MPI_File + + type MPI_Status + integer :: dummy = 0 + end type MPI_Status + + type MPI_Request + integer :: dummy = 0 + end type MPI_Request + + parameter MPI_OFFSET_KIND=8 +#endif + integer,dimension(:),allocatable :: proc_sizes, proc_offsets, send_or_recv integer :: comm_size, mpi_rank integer :: co_startdim, co_enddim @@ -38,6 +69,7 @@ subroutine co_init_blacs() if (.not. comms_inited) stop "CO_INIT_BLACS COMMS NOT INITED" +#ifdef TROVE_USE_MPI_ ! Must be initialised to zero - if stack contains garbage here MPI_Dims_create WILL fail blacs_dims = 0 @@ -51,6 +83,14 @@ subroutine co_init_blacs() call blacs_gridinfo(blacs_ctxt, nprow, npcol, myprow, mypcol) !write(*,"('BLACS: [',i2,',',i2'](',i4,i4,i4,i4',)')") mpi_rank,blacs_rank,nprow,npcol,myprow,mypcol +#else + blacs_size = 1 + blacs_rank = 0 + nprow = 1 + npcol = 1 + myprow = 0 + mypcol = 0 +#endif end subroutine co_init_blacs subroutine co_block_type_init(smat, dimx, dimy, descr, allocinfo, mpi_type) @@ -58,12 +98,14 @@ subroutine co_block_type_init(smat, dimx, dimy, descr, allocinfo, mpi_type) real(rk),intent(out),dimension(:,:),allocatable :: smat + ! Note: allocated matrix will have total dimensions (dimy x dimx) integer,intent(in) :: dimx, dimy integer,intent(out),dimension(9) :: descr integer,intent(out) :: allocinfo type(MPI_Datatype),intent(out),optional :: mpi_type +#ifdef TROVE_USE_MPI_ integer,dimension(2) :: global_size, distr, dargs integer :: MB,NB,MLOC,NLOC,ierr @@ -88,6 +130,10 @@ subroutine co_block_type_init(smat, dimx, dimy, descr, allocinfo, mpi_type) MPI_ORDER_FORTRAN, MPI_DOUBLE_PRECISION, mpi_type, ierr) call MPI_Type_commit(mpi_type, ierr) endif +#else + allocate(smat(dimy,dimx), stat=allocinfo) + if(allocinfo.ne.0) return +#endif end subroutine co_block_type_init @@ -97,6 +143,7 @@ subroutine co_sum_double(x, root_process) integer, optional :: root_process if (comm_size.eq.1) return +#ifdef TROVE_USE_MPI_ call TimerStart('co_sum_double') if (present(root_process)) then @@ -111,10 +158,12 @@ subroutine co_sum_double(x, root_process) end if call TimerStop('co_sum_double') +#endif end subroutine subroutine co_init_comms() +#ifdef TROVE_USE_MPI_ integer :: ierr call mpi_init(ierr) @@ -130,11 +179,17 @@ subroutine co_init_comms() if (mpi_rank.ne.0) then open(newunit=out, file='/dev/null', status='replace', iostat=ierr, action="write") endif - +#else + comm_size = 1 + mpi_rank = 0 +#endif comms_inited = .true. call co_init_blacs() + + + end subroutine co_init_comms subroutine co_finalize_comms() @@ -142,18 +197,24 @@ subroutine co_finalize_comms() if (.not. comms_inited) stop "CO_FINALIZE_COMMS COMMS NOT INITED" +#ifdef TROVE_USE_MPI_ call mpi_finalize(ierr) if (ierr .gt. 0) stop "MPI_FINALIZE" +#endif end subroutine co_finalize_comms + subroutine co_init_distr(dimen, startdim, enddim, blocksize) integer,intent(in) :: dimen integer,intent(out) :: startdim, enddim, blocksize + integer :: ierr + +#ifdef TROVE_USE_MPI_ integer,dimension(:),allocatable :: starts, ends integer :: localsize, proc_index, localsize_ - integer :: i, ierr, to_calc, ioslice_width, ioslice_maxwidth + integer :: i, to_calc, ioslice_width, ioslice_maxwidth if (.not. comms_inited) stop "COMMS NOT INITIALISED" !if (distr_inited) stop "DISTRIBUTION ALREADY INITIALISED" @@ -255,15 +316,34 @@ subroutine co_init_distr(dimen, startdim, enddim, blocksize) deallocate(starts,ends) +#else + if (.not. comms_inited) stop "COMMS NOT INITIALISED" + if (.not. distr_inited) then + allocate(send_or_recv(1),stat=ierr) + if (ierr .gt. 0) stop "CO_INIT_DISTR ALLOCATION FAILED" + endif + startdim = 1 + enddim = dimen + co_startdim = 1 + co_enddim = dimen + blocksize = dimen*dimen + send_or_recv(1) = 0 +#endif + distr_inited = .true. end subroutine co_init_distr + ! + ! Distribute the contents of an array among processes. + ! If only using one process or not using MPI, do nothing. + ! subroutine co_distr_data(x, tmp, blocksize, lb, ub) real(rk),dimension(:,lb:),intent(inout) :: x real(rk),dimension(:,:,:),intent(inout) :: tmp integer,intent(in) :: blocksize, lb, ub +#ifdef TROVE_USE_MPI_ integer :: i, icoeff, jcoeff, offset, ierr, k type(MPI_Request) :: reqs(comm_size) @@ -304,6 +384,7 @@ subroutine co_distr_data(x, tmp, blocksize, lb, ub) enddo call TimerStop('MPI_transpose_local') call TimerStop('MPI_transpose') +#endif end subroutine co_distr_data @@ -311,8 +392,9 @@ subroutine co_read_matrix_distr_ordered(x, longdim, lb, ub, infile) real(rk),dimension(:,lb:),intent(out) :: x integer,intent(in) :: longdim, lb, ub - type(MPI_File),intent(inout) :: infile + +#ifdef TROVE_USE_MPI_ type(MPI_Status) :: writestat integer(kind=MPI_Offset_kind) :: offset_start,offset_end integer :: readcount, ierr @@ -324,6 +406,7 @@ subroutine co_read_matrix_distr_ordered(x, longdim, lb, ub, infile) call MPI_File_read_ordered(infile,x,1,mpitype_column,writestat,ierr) call TimerStop('MPI_read_matrix') +#endif end subroutine co_read_matrix_distr_ordered @@ -331,8 +414,9 @@ subroutine co_read_matrix_distr(x, longdim, lb, ub, infile) real(rk),dimension(:,lb:),intent(out) :: x integer,intent(in) :: longdim, lb, ub - type(MPI_File),intent(inout) :: infile + +#ifdef TROVE_USE_MPI_ type(MPI_Status) :: writestat integer(kind=MPI_Offset_kind) :: offset_start,offset_end integer :: readcount, ierr @@ -353,6 +437,7 @@ subroutine co_read_matrix_distr(x, longdim, lb, ub, infile) endif call TimerStop('MPI_read_matrix') +#endif end subroutine co_read_matrix_distr @@ -361,6 +446,8 @@ subroutine co_write_matrix_distr(x, longdim, lb, ub, outfile) real(rk),dimension(:,lb:),intent(in) :: x integer,intent(in) :: longdim, lb, ub type(MPI_File),intent(inout) :: outfile + +#ifdef TROVE_USE_MPI_ integer :: writecount, ierr integer(kind=MPI_Offset_kind) :: offset_start, offset_end type(MPI_Status) :: writestat @@ -372,9 +459,11 @@ subroutine co_write_matrix_distr(x, longdim, lb, ub, outfile) call MPI_File_write_ordered(outfile,x,1,mpitype_column,writestat,ierr) call TimerStop('MPI_write_matrix') +#endif end subroutine co_write_matrix_distr +#ifdef TROVE_USE_MPI_ subroutine co_create_type_column(extent, blocksize, ncols) integer, intent(in) :: extent, blocksize, ncols integer :: ierr,writecount @@ -401,5 +490,6 @@ subroutine co_create_type_subarray(extent, coldim, rowdim, blockid, mpi_newtype) call MPI_Type_commit(mpi_newtype, ierr) end subroutine co_create_type_subarray +#endif end module diff --git a/perturbation.f90 b/perturbation.f90 index c1d5a5f..ba17da5 100644 --- a/perturbation.f90 +++ b/perturbation.f90 @@ -7725,12 +7725,15 @@ end subroutine PThamiltonian_contract !!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine open_chkptfile_mpi(fileh, filename, mode) +#ifdef TROVE_USE_MPI_ use mpi_f08 +#endif use mpi_aux type(MPI_File),intent(inout) :: fileh character(len=*),intent(in) :: filename, mode +#ifdef TROVE_USE_MPI_ integer :: amode, ierr select case(mode) @@ -7748,13 +7751,17 @@ subroutine open_chkptfile_mpi(fileh, filename, mode) !File errors indicate big trouble, so we set errors to be fatal - not likely to be recoverable call MPI_File_set_errhandler(fileh, MPI_ERRORS_ARE_FATAL) +#endif end subroutine open_chkptfile_mpi subroutine close_chkptfile_mpi(fileh) +#ifdef TROVE_USE_MPI_ use mpi_f08 +#endif use mpi_aux type(MPI_File), intent(inout) :: fileh +#ifdef TROVE_USE_MPI_ integer :: ierr call mpi_file_close(fileh, ierr) @@ -7762,27 +7769,25 @@ subroutine close_chkptfile_mpi(fileh) if (mpi_rank .eq. 0) write(*,*) "Error closing MPI-IO-formatted Vib. kinetic checkpoint file." stop "MPI_PTrestore_rot_kinetic_matrix_elements - Error closing MATELEM MPI-IO input file" endif +#endif end subroutine close_chkptfile_mpi subroutine PTrestore_rot_kinetic_matrix_elements_mpi(jrot, task, fileh, dimen, & - !PT, PTvibrational_me_calc,grot,gcor,hvib, & ncontr, maxcontr, icontr) +#ifdef TROVE_USE_MPI_ use mpi_f08 +#endif use mpi_aux integer(ik),intent(in) :: jrot character(len=cl),intent(in) :: task type(MPI_File),intent(inout) :: fileh integer(ik),intent(in) :: dimen - !type(PTelementsT),intent(in) :: PT - !logical,intent(in) :: PTvibrational_me_calc - !type(PTcontrME),pointer :: grot(:,:),gcor(:) ! rot. kinetic part - !type(PTcontrME) :: hvib ! rot. kinetic part - integer(ik),intent(inout),optional :: ncontr integer(ik),intent(inout),optional :: maxcontr integer(ik),intent(in),optional :: icontr +#ifdef TROVE_USE_MPI_ type(MPI_File) :: fileh_slice character(len=cl) :: job_id,filename,readbuf integer(kind=MPI_Offset_kind) :: file_offset @@ -8337,16 +8342,20 @@ subroutine PTrestore_rot_kinetic_matrix_elements_mpi(jrot, task, fileh, dimen, & !! endif !! ! end select +#endif end subroutine subroutine divided_slice_open_mpi(islice,fileh,chkpt_type,suffix) +#ifdef TROVE_USE_MPI_ use mpi_f08 +#endif use mpi_aux implicit none integer(ik),intent(in) :: islice type(MPI_File),intent(inout) :: fileh character(len=*),intent(in) :: chkpt_type,suffix +#ifdef TROVE_USE_MPI_ integer(ik) :: ilen character(len=cl) :: readbuf,filename,jchar integer :: ierr @@ -8370,15 +8379,20 @@ subroutine divided_slice_open_mpi(islice,fileh,chkpt_type,suffix) if (mpi_rank .eq. 0) write (out,"(' [MPI] kinetic checkpoint slice ',a20,': header is missing or wrong',a)") filename,readbuf(1:ilen) stop 'PTrestore_rot_kinetic_matrix_elements_mpi - in slice - header missing or wrong' end if +#endif end subroutine divided_slice_open_mpi subroutine divided_slice_close_mpi(islice,fileh,chkpt_type) +#ifdef TROVE_USE_MPI_ use mpi_f08 +#endif + use mpi_aux integer(ik),intent(in) :: islice type(MPI_File),intent(inout) :: fileh character(len=*),intent(in) :: chkpt_type +#ifdef TROVE_USE_MPI_ integer(ik) :: ilen character(len=cl) :: readbuf integer :: ierr @@ -8394,6 +8408,7 @@ subroutine divided_slice_close_mpi(islice,fileh,chkpt_type) end if ! call close_chkptfile_mpi(fileh) +#endif end subroutine divided_slice_close_mpi !!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -10617,10 +10632,6 @@ subroutine diagonalization_contract(jrot,gamma,dimen_s,mat,zpe,rlevel,total_root stop 'diagonalization_contract, mat - too small' endif ! - !call dseupd_p_arpack(dimen_s,bterm,nroots,job%factor,job%maxiter,job%tolerance,mat,energy) - ! - !call diag_dseupd_p(dimen_s,bterm,nroots,job%factor,job%maxiter,job%tolerance,mat,energy) - ! call dseupd_omp_arpack(dimen_s,bterm,nroots,job%factor,job%maxiter,job%tolerance,mat,energy) ! case('SEUPD') @@ -15413,6 +15424,7 @@ subroutine PTcontracted_matelem_class(jrot) ! job_is ='Vib. matrix elements of the rot. kinetic part' if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_open(mpi_comm_world, job%kinetmat_file, mpi_mode_wronly+mpi_mode_create, mpi_info_null, chkptMPIIO, ierr) call MPI_File_set_errhandler(chkptMPIIO, MPI_ERRORS_ARE_FATAL) @@ -15437,6 +15449,7 @@ subroutine PTcontracted_matelem_class(jrot) endif call MPI_Barrier(mpi_comm_world, ierr) call MPI_File_seek_shared(chkptMPIIO, mpioffset, MPI_SEEK_END, ierr) +#endif else call IOStart(trim(job_is),chkptIO) ! @@ -15625,12 +15638,14 @@ subroutine PTcontracted_matelem_class(jrot) if (trim(job%IOkinet_action)=='SAVE'.and..not.job%IOmatelem_split) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ if(mpi_rank.eq.0) then !call MPI_File_seek(chkptMPIIO, mpioffset, MPI_SEEK_END) call MPI_File_write_shared(chkptMPIIO,'g_rot',5,mpi_character,mpi_status_ignore,ierr) !else ! call MPI_File_write_shared(chkptMPIIO,'',0,mpi_character,mpi_status_ignore,ierr) endif +#endif else write(chkptIO) 'g_rot' endif @@ -15713,12 +15728,14 @@ subroutine PTcontracted_matelem_class(jrot) if (trim(job%IOkinet_action)=='SAVE'.and..not.job%IOmatelem_split) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ ! call MPI_File_seek(chkptMPIIO, mpioffset, MPI_SEEK_END) if(mpi_rank.eq.0) then call MPI_File_write_shared(chkptMPIIO,'g_cor',5,mpi_character,mpi_status_ignore,ierr) !else ! call MPI_File_write_shared(chkptMPIIO,'',0,mpi_character,mpi_status_ignore,ierr) endif +#endif else write(chkptIO) 'g_cor' endif @@ -16184,6 +16201,7 @@ subroutine PTcontracted_matelem_class(jrot) (.not.job%IOmatelem_split.or.job%iswap(1)==0)) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ if(mpi_rank.eq.0) then !call MPI_File_seek(chkptMPIIO, mpioffset, MPI_SEEK_END) call MPI_File_write_shared(chkptMPIIO,'hvib',4,mpi_character,mpi_status_ignore,ierr) @@ -16191,6 +16209,7 @@ subroutine PTcontracted_matelem_class(jrot) ! call MPI_File_write_shared(chkptMPIIO,'',0,mpi_character,mpi_status_ignore,ierr) endif call co_write_matrix_distr(hvib%me,mdimen, startdim, enddim,chkptMPIIO) +#endif else write(chkptIO) 'hvib' write(chkptIO) hvib%me @@ -16207,12 +16226,14 @@ subroutine PTcontracted_matelem_class(jrot) (.not.job%IOmatelem_split.or.job%iswap(1)==0 ) ) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call mpi_barrier(mpi_comm_world, ierr) call MPI_File_seek_shared(chkptMPIIO, mpioffset, MPI_SEEK_END) if(mpi_rank.eq.0) then call MPI_File_write_shared(chkptMPIIO,'End Kinetic part',16,mpi_character,mpi_status_ignore,ierr) endif call MPI_File_close(chkptMPIIO, ierr) +#endif else write(chkptIO) 'End Kinetic part' close(chkptIO,status='keep') @@ -16243,6 +16264,7 @@ subroutine PTcontracted_matelem_class(jrot) ! job_is ='external field contracted matrix elements for J=0' if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_open(mpi_comm_world, job%extFmat_file, mpi_mode_wronly+mpi_mode_create, mpi_info_null, chkptMPIIO, ierr) mpioffset=0 call MPI_File_set_size(chkptMPIIO, mpioffset, ierr) @@ -16267,6 +16289,7 @@ subroutine PTcontracted_matelem_class(jrot) ! ! ! call TimerStop('mpiiosingle') !AT endif +#endif else call IOStart(trim(job_is),chkptIO) ! @@ -16333,12 +16356,14 @@ subroutine PTcontracted_matelem_class(jrot) ! always store the matrix elements of the extF moment ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ if(mpi_rank.eq.0) then call MPI_File_write_shared(chkptMPIIO,imu,1,mpi_integer,mpi_status_ignore,ierr) endif !call mpi_barrier(mpi_comm_world,ierr) ! call co_write_matrix_distr(extF_t,mdimen, startdim, enddim,chkptMPIIO) +#endif else write(chkptIO) imu ! @@ -16386,8 +16411,11 @@ subroutine PTcontracted_matelem_class(jrot) ! endif ! +#ifdef TROVE_USE_MPI_ call mpi_barrier(mpi_comm_world,ierr) +#endif if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ if (mpi_rank.eq.0) then !AT if(.not.job%IOextF_divide) then !call MPI_File_seek(chkptMPIIO, mpioffset, MPI_SEEK_END) @@ -16397,6 +16425,7 @@ subroutine PTcontracted_matelem_class(jrot) endif endif call MPI_File_close(chkptMPIIO, ierr) +#endif else if (.not.job%IOextF_divide) write(chkptIO) 'End external field' endif @@ -16526,6 +16555,8 @@ subroutine write_divided_slice_mpi(islice,name,suffix,N,field) character(len=*),intent(in) :: name,suffix integer(ik),intent(in) :: N real(rk),dimension(:,:),intent(in) :: field + +#ifdef TROVE_USE_MPI_ character(len=4) :: jchar character(len=cl) :: filename character(len=cl) :: job_is @@ -16555,6 +16586,7 @@ subroutine write_divided_slice_mpi(islice,name,suffix,N,field) ! call MPI_File_close(chkptMPIIO, ierr) ! +#endif end subroutine write_divided_slice_mpi @@ -16602,6 +16634,8 @@ subroutine divided_slice_open_mpi(islice,chkptIO,name,suffix) integer(ik),intent(in) :: islice type(MPI_File),intent(inout) :: chkptIO character(len=*),intent(in) :: name,suffix + +#ifdef TROVE_USE_MPI_ character(len=4) :: jchar character(len=cl) :: buf,filename,job_is integer(ik) :: ilen @@ -16633,6 +16667,7 @@ subroutine divided_slice_open_mpi(islice,chkptIO,name,suffix) !stop 'PTrestore_rot_kinetic_matrix_elements - in slice - header missing or wrong' endif endif +#endif end subroutine divided_slice_open_mpi ! subroutine divided_slice_close(islice,chkptIO,name) @@ -16673,15 +16708,19 @@ subroutine divided_slice_close_mpi(islice,chkptIO,name) ilen = LEN_TRIM(name) ! if(mpi_rank .eq. 0) then +#ifdef TROVE_USE_MPI_ call MPI_File_read_shared(chkptIO, buf, ilen, mpi_character, mpi_status_ignore, ierr) if ( trim(buf(1:ilen))/=trim(name) ) then write (out,"(' divided_slice_close, kinetic checkpoint slice ',a,': footer is missing or wrong',a)") trim(name),buf(1:ilen) call MPI_Abort(mpi_comm_world, 1) !stop 'divided_slice_close - in slice - footer missing or wrong' endif +#endif endif ! +#ifdef TROVE_USE_MPI_ call MPI_File_close(chkptIO, ierr) +#endif ! end subroutine divided_slice_close_mpi @@ -34056,6 +34095,7 @@ subroutine PTstoreMPI_icontr_cnu(maxcontracts,iunit,dir) ! case ('SAVE') ! +#ifdef TROVE_USE_MPI_ call mpi_file_write(iunit, maxcontracts, 1,mpi_integer,mpi_status_ignore,ierr) ! call mpi_file_write(iunit, 'icontr_cnu', 10,mpi_character,mpi_status_ignore,ierr) @@ -34068,6 +34108,7 @@ subroutine PTstoreMPI_icontr_cnu(maxcontracts,iunit,dir) call mpi_file_write(iunit, pt%icontr_ideg(0:pt%nclasses,1:maxcontracts), (1+pt%nclasses)*maxcontracts, mpi_integer, & mpi_status_ignore, ierr) +#endif end select end subroutine PTstorempi_icontr_cnu diff --git a/tran.f90 b/tran.f90 index 2e4694a..3d89204 100644 --- a/tran.f90 +++ b/tran.f90 @@ -1369,7 +1369,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) !$omp parallel do private(icoeff,irow,ib,iterm,ielem,blacs_row,blacs_col,i_local,j_local) default(shared) schedule(dynamic) do icoeff = 1,dimen ! +#ifdef TROVE_USE_MPI_ call infog2l(icoeff,iroot,desc_psi,nprow,npcol,myprow,mypcol,i_local,j_local,blacs_row,blacs_col) +#endif if (myprow.eq.blacs_row.and.mypcol.eq.blacs_col) then ! psi(i_local,j_local) = 0 @@ -1469,7 +1471,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) read(iunit, rec = irec) vec ! do i=1,dimen +#ifdef TROVE_USE_MPI_ call infog2l(i,iroot,desc_psi,nprow,npcol,myprow,mypcol,i_local,j_local,blacs_row,blacs_col) +#endif if (myprow.eq.blacs_row.and.mypcol.eq.blacs_col) then psi(i_local,j_local) = vec(i) endif @@ -1486,9 +1490,11 @@ subroutine TRconvert_matel_j0_eigen(jrot) ! ! TODO TEMP intel bugaround - explicitly transpose psi into psi_t if (blacs_size .gt. 1) then +#ifdef TROVE_USE_MPI_ write(out,*) "Explicitly transposing psi into psi_t" call pdtran(Neigenroots, dimen, 1.0d0, psi, 1, 1, desc_psi, 0.0d0, psi_t, & 1, 1, desc_mat_t) +#endif endif ! if (job%verbose>=3) write(out,"(' ...done!')") @@ -1507,6 +1513,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) job_is ='Eigen-vib. matrix elements of the rot. kinetic part' ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_open(mpi_comm_world, job%kineteigen_file, mpi_mode_wronly+mpi_mode_create, mpi_info_null, fileh_w, ierr) call MPI_File_set_errhandler(fileh_w, MPI_ERRORS_ARE_FATAL) mpioffset = 0 @@ -1527,6 +1534,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) mpioffset = 0 treat_vibration = .false. endif +#endif else call IOStart(trim(job_is),chkptIO) ! @@ -1576,7 +1584,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) endif ! +#ifdef TROVE_USE_MPI_ call MPI_File_open(mpi_comm_world, job%kinetmat_file, mpi_mode_rdonly, mpi_info_null, fileh, ierr) +#endif ! The eigen-vibrational (J=0) matrix elements of the rotational and coriolis ! kinetic parts are being computed here. ! @@ -1587,11 +1597,13 @@ subroutine TRconvert_matel_j0_eigen(jrot) task = 'rot' ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_seek_shared(fileh_w, int(0,MPI_OFFSET_KIND), MPI_SEEK_END) if(mpi_rank.eq.0) call MPI_File_write_shared(fileh_w, 'g_rot', 5, mpi_character, mpi_status_ignore, ierr) call mpi_barrier(MPI_COMM_WORLD, ierr) ! call restore_rot_kinetic_matrix_elements_mpi(jrot,treat_vibration,task,fileh) +#endif else write(chkptIO) 'g_rot' ! @@ -1617,12 +1629,14 @@ subroutine TRconvert_matel_j0_eigen(jrot) islice = 0 ! if ((.not.job%IOmatelem_split) .and. blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ call MPI_File_get_position(fileh, read_offset, ierr) call MPI_File_set_view(fileh, read_offset, mpi_byte, gmat_block_type, "native", MPI_INFO_NULL, ierr) !call MPI_File_seek_shared(fileh_w, int(0,MPI_OFFSET_KIND),MPI_SEEK_END) call MPI_File_get_position_shared(fileh_w, write_offset, ierr) call MPI_File_set_view(fileh_w, write_offset, mpi_byte, mat_s_block_type, "native", MPI_INFO_NULL, ierr) +#endif endif ! do k1 = 1,3 @@ -1648,7 +1662,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) else ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_read_all(fileh, gmat, size(gmat), mpi_double_precision, mpi_status_ignore, ierr) +#endif else read(iunit) gmat endif @@ -1657,10 +1673,12 @@ subroutine TRconvert_matel_j0_eigen(jrot) ! ! if (blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ call pdgemm('T','N',Neigenroots,dimen,dimen,alpha,psi,1,1,desc_psi,& gmat,1,1,desc_gmat,beta,mat_t,1,1,desc_mat_t) call pdgemm('N','T',Neigenroots,Neigenroots,dimen,alpha,mat_t,1,1,desc_mat_t,& psi_t,1,1,desc_mat_t,beta,mat_s,1,1,desc_mat_s) +#endif else call dgemm('T','N',Neigenroots,dimen,dimen,alpha,psi,dimen,& gmat,dimen,beta,mat_t,Neigenroots) @@ -1683,7 +1701,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) else ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_write_all(fileh_w, mat_s, size(mat_s), mpi_double_precision, mpi_status_ignore, ierr) +#endif else write (chkptIO) mat_s endif @@ -1696,6 +1716,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) ! ! Reset view to flat file if ((.not.job%IOmatelem_split) .and. blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ read_offset = read_offset + 9*int(dimen,MPI_OFFSET_KIND)*dimen*mpi_real_size call MPI_File_set_view(fileh, int(0,MPI_OFFSET_KIND), mpi_byte, mpi_byte, "native", MPI_INFO_NULL, ierr) call MPI_File_seek(fileh, read_offset, MPI_SEEK_SET) @@ -1706,6 +1727,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) call MPI_File_seek_shared(fileh_w, write_offset, MPI_SEEK_SET) !write_offset = 0 !call MPI_File_seek_shared(fileh_w, write_offset, MPI_SEEK_END) +#endif endif ! @@ -1719,11 +1741,13 @@ subroutine TRconvert_matel_j0_eigen(jrot) task = 'cor' if (trim(job%kinetmat_format).eq.'MPIIO') then ! +#ifdef TROVE_USE_MPI_ call restore_rot_kinetic_matrix_elements_mpi(jrot,treat_vibration,task,fileh) ! if(mpi_rank.eq.0) call MPI_File_write_shared(fileh_w, 'g_cor', 5, mpi_character, mpi_status_ignore, ierr) call MPI_Barrier(MPI_COMM_WORLD, ierr) !call MPI_File_seek_shared(fileh_w, int(0,MPI_OFFSET_KIND), MPI_SEEK_END) +#endif else call restore_rot_kinetic_matrix_elements(jrot,treat_vibration,task,iunit) ! @@ -1739,11 +1763,13 @@ subroutine TRconvert_matel_j0_eigen(jrot) islice = 9 ! if ((.not.job%IOmatelem_split) .and. blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ call MPI_File_get_position(fileh, read_offset, ierr) call MPI_File_set_view(fileh, read_offset, mpi_byte, gmat_block_type, "native", MPI_INFO_NULL, ierr) call MPI_File_get_position_shared(fileh_w, write_offset, ierr) call MPI_File_set_view(fileh_w, write_offset, mpi_byte, mat_s_block_type, "native", MPI_INFO_NULL, ierr) +#endif endif ! !do k1 = 1,FLNmodes @@ -1771,7 +1797,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) else ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_read_all(fileh, gmat, size(gmat), mpi_double_precision, mpi_status_ignore, ierr) +#endif else read(iunit) gmat endif @@ -1779,10 +1807,12 @@ subroutine TRconvert_matel_j0_eigen(jrot) endif ! if (blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ call pdgemm('T','N',Neigenroots,dimen,dimen,alpha,psi,1,1,desc_psi,& gmat,1,1,desc_gmat,beta,mat_t,1,1,desc_mat_t) call pdgemm('N','T',Neigenroots,Neigenroots,dimen,alpha,mat_t,1,1,desc_mat_t,& psi_t,1,1,desc_mat_t,beta,mat_s,1,1,desc_mat_s) +#endif else call dgemm('T','N',Neigenroots,dimen,dimen,alpha,psi,dimen,& gmat,dimen,beta,mat_t,Neigenroots) @@ -1806,7 +1836,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) else ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_write_all(fileh_w, mat_s, size(mat_s), mpi_double_precision, mpi_status_ignore, ierr) +#endif else write (chkptIO) mat_s endif @@ -1819,6 +1851,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) ! ! Reset view to flat file if ((.not.job%IOmatelem_split) .and. blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ read_offset = read_offset + 3*int(dimen,MPI_OFFSET_KIND)*dimen*mpi_real_size call MPI_File_set_view(fileh, int(0,MPI_OFFSET_KIND), mpi_byte, mpi_byte, "native", MPI_INFO_NULL, ierr) call MPI_File_seek(fileh, read_offset, MPI_SEEK_SET) @@ -1828,13 +1861,16 @@ subroutine TRconvert_matel_j0_eigen(jrot) call MPI_File_set_view(fileh_w, int(0,MPI_OFFSET_KIND), mpi_byte, mpi_byte, "native", MPI_INFO_NULL, ierr) call MPI_File_seek_shared(fileh_w, write_offset, MPI_SEEK_END) !write_offset = 0 +#endif endif ! if (job%verbose>=5) call TimerStop('J0-convertion for g_cor') ! if ((.not.job%IOmatelem_split.or.job%iswap(1)==1).and.(mpi_rank.eq.0)) then if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_write_shared(fileh_w, 'End Kinetic part', 16, mpi_character, mpi_status_ignore, ierr) +#endif else write(chkptIO) 'End Kinetic part' endif @@ -1842,7 +1878,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) ! if (.not.job%vib_rot_contr) then if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_close(fileh_w, ierr) +#endif else close(chkptIO,status='keep') endif @@ -1864,7 +1902,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) endif ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_close(fileh, ierr) +#endif else close(chkptIO,status='keep') endif @@ -1898,6 +1938,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) ! if (trim(job%kinetmat_format).eq.'MPIIO') then ! +#ifdef TROVE_USE_MPI_ call MPI_File_open(mpi_comm_world, filename, mpi_mode_rdonly, mpi_info_null, fileh, ierr) ! call MPI_File_read_all(fileh, buf20, 7, mpi_character, mpi_status_ignore, ierr) @@ -1920,6 +1961,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) stop 'restore_Extvib_matrix_elements - in file - illegal ncontracts ' end if ! +#endif else open(iunit,form='unformatted',action='read',position='rewind',status='old',file=filename) ! @@ -1948,6 +1990,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) call IOStart(trim(job_is),chkptIO) ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ ! call mpi_file_open(mpi_comm_world, job%exteigen_file, mpi_mode_wronly+mpi_mode_create, mpi_info_null, fileh_w, ierr) call mpi_file_set_errhandler(fileh_w, mpi_errors_are_fatal) @@ -1966,6 +2009,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) call mpi_barrier(mpi_comm_world, ierr) call mpi_file_seek(fileh_w, int(0,mpi_offset_kind), mpi_seek_end) ! +#endif else ! open(chkptio,form='unformatted',action='write',position='rewind',status='replace',file=job%exteigen_file) @@ -1994,11 +2038,13 @@ subroutine TRconvert_matel_j0_eigen(jrot) endif ! if ((.not.job%IOextF_divide) .and. blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ call MPI_File_get_position(fileh, read_offset, ierr) call MPI_File_set_view(fileh, read_offset, mpi_byte, extF_block_type, "native", MPI_INFO_NULL, ierr) call MPI_File_get_position_shared(fileh_w, write_offset, ierr) call MPI_File_set_view(fileh_w, write_offset, mpi_byte, mat_s_block_type, "native", MPI_INFO_NULL, ierr) +#endif endif ! do imu = fitting%iparam(1),fitting%iparam(2) @@ -2020,9 +2066,11 @@ subroutine TRconvert_matel_j0_eigen(jrot) else ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_read_all(fileh, imu_t, 1, mpi_integer, mpi_status_ignore, ierr) ! call MPI_File_read_all(fileh, extF_me, size(extF_me), mpi_double_precision, mpi_status_ignore, ierr) +#endif else read(iunit) imu_t ! @@ -2032,10 +2080,12 @@ subroutine TRconvert_matel_j0_eigen(jrot) endif ! if(blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ call pdgemm('T','N',Neigenroots,dimen,dimen,alpha,psi,1,1,desc_psi,& extF_me,1,1,desc_gmat,beta,mat_t,1,1,desc_mat_t) call pdgemm('N','T',Neigenroots,Neigenroots,dimen,alpha,mat_t,1,1,desc_mat_t,& psi_t,1,1,desc_mat_t,beta,mat_s,1,1,desc_mat_s) +#endif else call dgemm('T','N',Neigenroots,dimen,dimen,alpha,psi,dimen,& extF_me,dimen,beta,mat_t,Neigenroots) @@ -2072,10 +2122,12 @@ subroutine TRconvert_matel_j0_eigen(jrot) if (.not.job%IOextF_divide.or.job%IOextF_stitch) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ if(mpi_rank.eq.0) call MPI_File_write(fileh_w, imu, 1, mpi_integer, mpi_status_ignore, ierr) call MPI_Barrier(mpi_comm_world, ierr) call MPI_File_seek_shared(fileh_w, int(0,MPI_OFFSET_KIND), MPI_SEEK_END) call MPI_File_write_all(fileh_w, mat_s, size(mat_s), mpi_double_precision, mpi_status_ignore, ierr) +#endif else write(chkptIO) imu write(chkptIO) mat_s @@ -2094,6 +2146,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) enddo ! Reset view to flat file if ((.not.job%IOextF_divide) .and. blacs_size.gt.1) then +#ifdef TROVE_USE_MPI_ read_offset = read_offset + (fitting%iparam(2)-fitting%iparam(1)+1)*int(ncontr_t,MPI_OFFSET_KIND)*ncontr_t*mpi_real_size & + (fitting%iparam(2)-fitting%iparam(1)+1)*mpi_int_size call MPI_File_set_view(fileh, int(0,MPI_OFFSET_KIND), mpi_byte, mpi_byte, "native", MPI_INFO_NULL, ierr) @@ -2104,6 +2157,7 @@ subroutine TRconvert_matel_j0_eigen(jrot) call MPI_File_seek_shared(fileh_w, write_offset, MPI_SEEK_SET) !write_offset = 0 !call MPI_File_seek_shared(fileh_w, write_offset, MPI_SEEK_END) +#endif endif ! if (allocated(extF_me)) deallocate(extF_me) @@ -2112,7 +2166,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) if (.not.job%IOextF_divide) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_read_all(fileh, buf20, 18, mpi_character, mpi_status_ignore, ierr) +#endif else read(iunit) buf20(1:18) endif @@ -2123,7 +2179,9 @@ subroutine TRconvert_matel_j0_eigen(jrot) end if ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ call MPI_File_close(fileh, ierr) +#endif else close(iunit,status='keep') endif @@ -2136,8 +2194,10 @@ subroutine TRconvert_matel_j0_eigen(jrot) if (.not.job%IOextF_divide.or.job%IOextF_stitch) then ! if (trim(job%kinetmat_format).eq.'MPIIO') then +#ifdef TROVE_USE_MPI_ if(mpi_rank.eq.0) call MPI_File_write(fileh_w, 'End external field', 18, mpi_character, mpi_status_ignore, ierr) call MPI_File_close(fileh_w, ierr) +#endif else write(chkptIO) 'End external field' close(chkptIO,status='keep') @@ -2213,6 +2273,8 @@ subroutine divided_slice_write_mpi(islice,name,suffix,N,field,block_type) integer(ik),intent(in) :: N real(rk),dimension(:,:),intent(in) :: field type(MPI_Datatype),intent(in) :: block_type + +#ifdef TROVE_USE_MPI_ character(len=4) :: jchar character(len=cl) :: filename type(MPI_File) :: chkptMPIIO @@ -2245,6 +2307,7 @@ subroutine divided_slice_write_mpi(islice,name,suffix,N,field,block_type) ! call MPI_File_close(chkptMPIIO, ierr) ! +#endif end subroutine divided_slice_write_mpi ! ! @@ -2321,6 +2384,7 @@ subroutine divided_slice_read_mpi(islice,name,suffix,N,field,block_type,ierr) real(rk),dimension(:,:),intent(out) :: field type(MPI_Datatype),intent(in) :: block_type integer(ik),intent(out) :: ierr +#ifdef TROVE_USE_MPI_ character(len=4) :: jchar type(MPI_File) :: chkptMPIIO character(len=cl) :: buf,filename,job_is @@ -2379,6 +2443,7 @@ subroutine divided_slice_read_mpi(islice,name,suffix,N,field,block_type,ierr) ! call MPI_File_close(chkptMPIIO, ierr) ! +#endif end subroutine divided_slice_read_mpi ! ! @@ -2738,6 +2803,7 @@ subroutine restore_rot_kinetic_matrix_elements_mpi(jrot,treat_vibration,kinetic_ logical,intent(in) :: treat_vibration character(len=3),intent(in) :: kinetic_part type(MPI_File),intent(in) :: fileh +#ifdef TROVE_USE_MPI_ ! integer(ik) :: iclasses,alloc,ilevel,ib,max_deg_size,nclasses,islice character(len=cl) :: job_is @@ -2968,6 +3034,7 @@ subroutine find_groundstate_icontr(Nclasses) ! end subroutine find_groundstate_icontr ! +#endif end subroutine restore_rot_kinetic_matrix_elements_mpi ! !