subroutine pmfread(*) include "DIMENSIONS" include "DIMENSIONS.ZSCOPT" #ifdef MPI include 'mpif.h' include 'COMMON.MPI' include 'COMMON.MPISET' integer ierror double precision bufin(maxpoint),bufout(maxpoint) #endif include "COMMON.PMFDERIV" include "COMMON.IOUNITS" include "COMMON.WEIGHTS" c include "COMMON.MARQ" c character*64 fname(8) c common /kuku/ it,jt character*1 typ(-2:2) /'p','a','G','A','P'/ character*64 torfile character*16 pmf_print integer ilen external ilen logical lprn double precision rjunk,sumw(4) write (iout,*) "Calling PMFREAD",torsion_pmf,turn_pmf,eello_pmf call getenv("PMF_DIR",pmf_dir) call getenv("PMF_INFIX",pmf_infix) write (iout,*) "PMF_DIR ",PMF_DIR write (iout,*) "PMF_INFIX",PMF_INFIX call getenv("PMF_PRINT",pmf_print) lprn=index(pmf_print,"PMF_PRINT").gt.0 c logical lder c do i=1,n c read(iparin,'(e20.10,1x,a)') x(i),parname(i) c print *,i,x(i)," ",parname(i) c enddo c maxit=100 c maxmar=10 c rl=1.0d2 c vmarq=1.0d1 c tolx=1.0d-3 c tolf=1.0d-2 c rtolf=1.0d-3 c tollam=1.0d2 c read (iparin,*,end=12,err=12) c & maxit,maxmar,rl,vmarq,tolx,tolf,rtolf,tollam c 12 continue c write (iout,'(//21(1h*)/a/21(1h*))') "Initial parameters:" c do i=1,n c write (iout,'(i3,1x,a16,f10.5)') i,parname(i),x(i) c enddo ifunmode=1 m=3 c write (iout,'(a,i3)') 'number of variables ',m c write (iout,'(a,i3)') 'number of parameters ',n ii=1 i1=1 c close(inp) npoint=0 do it1=0,2 do it2=-2,2 c do it1=1,1 c do it2=1,1 c goto 1111 if (torsion_pmf .or. turn_pmf) then #ifdef MPI write (iout,*) "Master",Master if (me.eq.Master) then #endif torfile = PMF_DIR(:ilen(PMF_DIR))//'/'//typ(it1)//typ(it2)//"_" & //PMF_INFIX(:ilen(PMF_INFIX)) & //".eturn3" write (iout,*) "Reading torsional and turn PMFs" write (iout,'(2a)') "Torsional and turn3 potential file ",torfile open(89,file=torfile,status="old",err=99) iii = 1 i1 = ii do i=1,maxpoint read (89,*,end=13,err=13) zz(1,ii),zz(2,ii),zz(3,ii),y(1,ii), & rjunk,y(2,ii) if (y(1,ii).ge.999.0 .or. y(2,ii).ge.999.0) cycle c if (zz(1,ii).ne.75.0d0 .or. zz(2,ii).ne.130.0d0 .or. c & zz(3,ii).ne.-45.0d0) cycle if (lprn) & write(iout,'(i7,3f7.1,2f10.5)')ii,zz(1,ii),zz(2,ii),zz(3,ii), & y(1,ii),y(2,ii) ii=ii+1 iii = iii+1 enddo 13 np(it1,it2)=iii-1 npoint=npoint+2*np(it1,it2) close (89) write (iout,*)np(it1,it2)," data points found in file" write (iout,*)npoint," data points found so far" call flush(iout) sumw=0.0d0 ymin=y(1,i1) ymax=y(1,i1) do i=i1,ii-1 if (y(1,i).gt.ymax) ymax=y(1,i) if (y(1,i).lt.ymin) ymin=y(1,i) enddo do i=i1,ii-1 if (ifunmode.eq.0) then w(1,i)=dexp(-0.169*(y(1,i)-ymin)) w(2,i)=1.0*dexp(-0.168*y(2,i)) else w(1,i)=1.0d0 w(2,i)=1.0d0 if (y(1,i).ge.999.0) w(1,i)=0.0d0 if (y(2,i).ge.999.0) w(2,i)=0.0d0 endif do j=1,2 sumw(j)=sumw(j)+w(j,i) enddo enddo do i=i1,ii-1 do j=1,2 w(j,i)=w(j,i)/sumw(j) enddo enddo close(89) #ifdef MPI np_=np(it1,it2) call work_partition_pmf(np_,.true.) endif write (iout,*) "After scatter" write (iout,*) "scount_pmf",(scount_pmf(i),i=1,procs) call flush(iout) call MPI_Scatter(scount_pmf,1,MPI_INTEGER,np(it1,it2),1, & MPI_INTEGER,Master,WHAM_COMM,ierror) write (iout,*) "After scatter scount np",np(it1,it2) call flush(iout) do j=1,2 do i=i1,ii-1 bufin(i)=y(j,i) enddo call MPI_Scatterv(bufin(i1),scount_pmf,idispl_pmf, & MPI_DOUBLE_PRECISION, & bufout(1),np(it1,it2),MPI_DOUBLE_PRECISION,Master,WHAM_COMM, & ierror) do i=1,np(it1,it2) y(j,i+i1-1)=bufout(i) enddo write (iout,*) "After scatter y",j call flush(iout) do i=i1,ii-1 bufin(i)=w(j,i) enddo call MPI_Scatterv(bufin(i1),scount_pmf,idispl_pmf, & MPI_DOUBLE_PRECISION, & bufout(1),np(it1,it2),MPI_DOUBLE_PRECISION,Master,WHAM_COMM, & ierror) do i=1,np(it1,it2) w(j,i+i1-1)=bufout(i) enddo write (iout,*) "After scatter w",j call flush(iout) enddo do j=1,4 do i=i1,ii-1 bufin(i)=zz(j,i) enddo call MPI_Scatterv(bufin(i1),scount_pmf,idispl_pmf, & MPI_DOUBLE_PRECISION, & bufout(1),np(it1,it2),MPI_DOUBLE_PRECISION,Master,WHAM_COMM, & ierror) write (iout,*) "After scatter z",j call flush(iout) do i=1,np(it1,it2) zz(j,i+i1-1)=bufout(i) enddo enddo if (lprn) then write (iout,*) "Distributed etor and eturn" do i=i1,i1+np(it1,it2)-1 write(iout,'(i7,3f7.1,2f10.5)')i,(zz(j,i),j=1,3),(y(j,i),j=1,2) enddo endif ii=i1+np(it1,it2) i1=ii #endif 1111 continue write (iout,*)np(it1,it2)," data points found in file" endif ! torsion_pmf .or. turn_pmf if (eello_pmf) then #define PIZDUN #ifdef PIZDUN #ifdef MPI if (me.eq.Master) then #endif torfile = PMF_DIR(:ilen(PMF_DIR))//'/'//typ(it1)//typ(it2)//"_"// & PMF_INFIX(:ilen(PMF_INFIX)) & //".eello3" write (iout,*) "Reading eello3 PMFs" write (iout,'(2a)') "local-electrostatic PMF file ",torfile open(89,file=torfile,status="old",err=99) iii = 1 i1 = ii do i=1,maxpoint read (89,*,end=133,err=133) (zz(j,ii),j=1,8),(y(j,ii),j=1,4) if (y(1,ii).ge.999.0 .or. y(2,ii).ge.999.0 .or. y(3,ii).ge.999.0 & .or. y(4,ii).ge.999.0) cycle if (lprn) write(iout,'(i7,8f7.1,4f10.5)')ii,(zz(j,ii),j=1,8), & (y(j,ii),j=1,4) ii=ii+1 iii = iii+1 enddo 133 npel(it1,it2)=iii-1 close (89) npoint=npoint+npel(it1,it2)*4 sumw=0.0d0 do i=i1,ii do j=1,4 if (ifunmode.eq.0) then w(j,ii)=10*dexp(-0.169*(y(j,ii))) else w(j,ii)=1.0d0 endif enddo do j=1,4 sumw(j)=sumw(j)+w(j,i) enddo enddo do i=i1,ii-1 do j=1,4 w(j,i)=w(j,i)/sumw(j) enddo enddo #ifdef MPI np_=npel(it1,it2) call work_partition_pmf(np_,.true.) endif write (iout,*) "Before scatter npel" call MPI_Scatter(scount_pmf,1,MPI_INTEGER,npel(it1,it2),1, & MPI_INTEGER,Master,WHAM_COMM,ierror) write (iout,*) "After scatter scount npel",npel(it1,it2) call flush(iout) do j=1,4 do i=i1,ii-1 bufin(i)=y(j,i) enddo call MPI_Scatterv(bufin(i1),scount_pmf,idispl_pmf, & MPI_DOUBLE_PRECISION, & bufout(1),npel(it1,it2),MPI_DOUBLE_PRECISION,Master,WHAM_COMM, & ierror) write (iout,*) "After scatter y",j call flush(iout) do i=1,npel(it1,it2) y(j,i+i1-1)=bufout(i) enddo do i=i1,ii-1 bufin(i)=w(j,i) enddo call MPI_Scatterv(bufin(i1),scount_pmf,idispl_pmf, & MPI_DOUBLE_PRECISION, & bufout(1),npel(it1,it2),MPI_DOUBLE_PRECISION,Master,WHAM_COMM, & ierror) write (iout,*) "After scatter w",j call flush(iout) do i=1,npel(it1,it2) w(j,i+i1-1)=bufout(i) enddo enddo do j=1,8 do i=i1,ii-1 bufin(i)=zz(j,i) enddo call MPI_Scatterv(bufin(i1),scount_pmf,idispl_pmf, & MPI_DOUBLE_PRECISION, & bufout(1),npel(it1,it2),MPI_DOUBLE_PRECISION,Master,WHAM_COMM, & ierror) do i=1,npel(it1,it2) zz(j,i+i1-1)=bufout(i) enddo write (iout,*) "After scatter z",j call flush(iout) enddo if (lprn) then write (iout,*) "Distributed eelloc" do i=i1,i1+npel(it1,it2)-1 write(iout,'(i7,8f7.1,4f10.5)')i,(zz(j,i),j=1,8),(y(j,i),j=1,4) enddo endif ii=i1+npel(it1,it2) i1=ii #endif #endif write (iout,*)npel(it1,it2)," data points found in file" write (iout,*)npoint," data points found so far" call flush(iout) endif ! eello_pmf enddo ! it2 enddo ! it1 return 99 write (iout,*) "Error opening PMF fiie(s)" return1 end #ifdef MPI c------------------------------------------------------------------------------- subroutine work_partition_pmf(n,lprint) implicit none include 'DIMENSIONS' include "DIMENSIONS.ZSCOPT" include 'mpif.h' include 'COMMON.MPI' include 'COMMON.MPISET' integer ierror,errcode integer n,ntot,nn,chunk,i,remainder logical lprint c nprocs1=numtasks C C Divide conformations between processors; for each proteins the first and C the last conformation to handle by ith processor is stored in C indstart(i) and indend(i), respectively. C C First try to assign equal number of conformations to each processor. C ntot=n indstart_pmf(0)=1 chunk = N/nprocs1 scount_pmf(0) = chunk do i=1,nprocs1-1 indstart_pmf(i)=chunk+indstart_pmf(i-1) scount_pmf(i)=scount_pmf(i-1) enddo C C Determine how many conformations remained yet unassigned. C remainder=N-(indstart_pmf(nprocs1-1)+scount_pmf(nprocs1-1)-1) c print *,"remainder",remainder C C Assign the remainder conformations to consecutive processors, starting C from the lowest rank; this continues until the list is exhausted. C if (remainder .gt. 0) then do i=1,remainder scount_pmf(i-1) = scount_pmf(i-1) + 1 indstart_pmf(i) = indstart_pmf(i) + i enddo do i=remainder+1,nprocs1-1 indstart_pmf(i) = indstart_pmf(i) + remainder enddo endif indstart_pmf(nprocs1)=N+1 scount_pmf(nprocs1)=0 do i=0,NProcs1 indend_pmf(i)=indstart_pmf(i)+scount_pmf(i)-1 idispl_pmf(i)=indstart_pmf(i)-1 enddo N=0 do i=0,Nprocs1-1 N=N+indend_pmf(i)-indstart_pmf(i)+1 enddo if (N.ne.ntot) then write (2,*) "!!! Checksum error on processor",me call flush(2) call MPI_Abort( MPI_COMM_WORLD, Ierror, Errcode ) endif if (lprint) then write (2,*) "Partition of work between processors" do i=0,nprocs1-1 write (2,'(a,i5,a,i7,a,i7,a,i7)') & "Processor",i," indstart",indstart_pmf(i), & " indend",indend_pmf(i)," count",scount_pmf(i)," idispl", & idispl_pmf(i) enddo endif return end #endif