Changeset bd2278d for partem_p.f

Timestamp:

09/05/08 11:49:42 (16 years ago)

Author:

baerbaer <baerbaer@…>

Branches:

master

Children:

fafe4d6

Parents:

2ebb8b6

Message:

Reformatting comments and continuation marks.

Fortran 90 and higher use ! to mark comments no matter where they are in the
code. The only valid continuation marker is &.
I also added the SMMP.kdevelop.filelist to the repository to make it easier
to use kdevelop.

git-svn-id: ​svn+ssh://svn.berlios.de/svnroot/repos/smmp/trunk@12 26dc1dd8-5c4e-0410-9ffe-d298b4865968

File:

• partem_p.f (modified) (31 diffs)

partem_p.f

@@ -1 @@
-c**************************************************************
-c     
-c This file contains the subroutines: partem_p 
-C USE WITH main_p, NOT WITH main!!!!!!
-c     
-c Copyright 2003-2005  Frank Eisenmenger, U.H.E. Hansmann,
-c                      Shura Hayryan, Chin-Ku 
-c Copyright 2007       Frank Eisenmenger, U.H.E. Hansmann,
-c                      Jan H. Meinke, Sandipan Mohanty
-c     
-c     **************************************************************
+!**************************************************************
+!     
+! This file contains the subroutines: partem_p 
+! USE WITH main_p, NOT WITH main!!!!!!
+!     
+! Copyright 2003-2005  Frank Eisenmenger, U.H.E. Hansmann,
+!                      Shura Hayryan, Chin-Ku 
+! Copyright 2007       Frank Eisenmenger, U.H.E. Hansmann,
+!                      Jan H. Meinke, Sandipan Mohanty
+!     
+!     **************************************************************
 
       subroutine  partem_p(num_rep, nequi, nswp, nmes, nsave, newsta,
      &                     switch, rep_id, partem_comm)
-C     
-C     PURPOSE: SIMULATION OF PROTEINS BY PARALLEL TEMPERING ALGORITHM 
-C     ON PARALLEL COMPUTERS USING MPI
-C     
-C     switch: Choses the starting configuration:
-C     -1 - stretched configuration
-C     0 - don't change anything
-C     1 - random start configuration
-C     
-c     CALLS:  addang,contacts,energy,hbond,helix,iendst,metropolis,
-c     outvar,(rand),rgyr
-C     
+!     
+!     PURPOSE: SIMULATION OF PROTEINS BY PARALLEL TEMPERING ALGORITHM 
+!     ON PARALLEL COMPUTERS USING MPI
+!     
+!     switch: Choses the starting configuration:
+!     -1 - stretched configuration
+!     0 - don't change anything
+!     1 - random start configuration
+!     
+!     CALLS:  addang,contacts,energy,hbond,helix,iendst,metropolis,
+!     outvar,(rand),rgyr
+!     
       include 'INCL.H'
       include 'INCP.H'
@@ -31 +31 @@
       logical newsta
       integer switch, partem_comm, rep_id, nsave
-c     external rand
+!     external rand
       external can_weight
 
-C     nequi:  number of Monte Carlo sweeps for thermalization
-C     nswp:   number of Monte Carlo sweeps
-C     nmes:   number of Monte Carlo sweeps between measurments
-C     newsta: .true. for new simulations, .false. for re-start
+!     nequi:  number of Monte Carlo sweeps for thermalization
+!     nswp:   number of Monte Carlo sweeps
+!     nmes:   number of Monte Carlo sweeps between measurments
+!     newsta: .true. for new simulations, .false. for re-start
 
       dimension  eavm(MAX_PROC),sph(MAX_PROC),intem(MAX_PROC),
@@ -47 +47 @@
       double precision    e_min, e_minp(MAX_PROC), e_minpt(MAX_PROC)
       integer   h_max, h_maxp(MAX_PROC)
-c     Order of replica exchange
+!     Order of replica exchange
       integer   odd
 !     Counter to keep random number generators in sync
       integer randomCount
       
-c     Collect partial energies. Only the root writes to disk. We have to
-c     collect the information from the different replicas and provide 
-c     arrays to store them.
-c     eyslr    storage array for solvent energy
-c     eyelp     -      "        - coulomb energy
-c     eyvwp     -      "        - van-der-Waals energy
-c     eyhbp     -      "        - hydrogen bonding energy
-c     eysmi    -      "        - intermolecular interaction energy
-c     eyabp     -      "        - Abagyan correction term
+!     Collect partial energies. Only the root writes to disk. We have to
+!     collect the information from the different replicas and provide 
+!     arrays to store them.
+!     eyslr    storage array for solvent energy
+!     eyelp     -      "        - coulomb energy
+!     eyvwp     -      "        - van-der-Waals energy
+!     eyhbp     -      "        - hydrogen bonding energy
+!     eysmi    -      "        - intermolecular interaction energy
+!     eyabp     -      "        - Abagyan correction term
       double precision eyslr(MAX_PROC)
       double precision eyelp(MAX_PROC),eyvwp(MAX_PROC),eyhbp(MAX_PROC), 
      &     eyvrp(MAX_PROC),eysmip(MAX_PROC), eyabp(MAX_PROC)
-c     Collect information about accessible surface and van-der-Waals volume
-c     asap      storage array for solvent accessible surface
-c     vdvolp     storage array for van-der-Waals volume
+!     Collect information about accessible surface and van-der-Waals volume
+!     asap      storage array for solvent accessible surface
+!     vdvolp     storage array for van-der-Waals volume
       double precision asap(MAX_PROC), vdvolp(MAX_PROC)
 
@@ -73 +73 @@
       integer imhbp(MAX_PROC)
       character*80 filebase, fileNameMP, tbase0,tbase1
-c     frame     frame number for writing configurations
-c     trackID   configuration that should be tracked and written out
-c     dir          direction in random walk
-c     -1 - visited highest temperature last
-c     1 - visited lowest temperature last
-c     0 - haven't visited the boundaries yet.
-c     dirp      storage array for directions.
+!     frame     frame number for writing configurations
+!     trackID   configuration that should be tracked and written out
+!     dir          direction in random walk
+!     -1 - visited highest temperature last
+!     1 - visited lowest temperature last
+!     0 - haven't visited the boundaries yet.
+!     dirp      storage array for directions.
       integer frame, trackID, dir
       integer dirp(MAX_PROC)
@@ -90 +90 @@
      &            rep_id, num_rep, partem_comm, myrank
       call flush(6)
-C     
-c     
-C     File with temperatures 
+!     
+!     
+!     File with temperatures 
       open(11,file='temperatures',status='old')
-C     File with reference conformation
+!     File with reference conformation
       tbase0='trj_00000'
       open(18,file=fileNameMP(tbase0,5,9,rep_id),status='unknown')
       if (rep_id.eq.0.and.myrank.eq.0) then
-c     File with time series of simulation
+!     File with time series of simulation
          open(14,file='ts.d',status='unknown')
-c     Track weights
-c      open(16, file='weights.dat', status='unknown')
+!     Track weights
+!      open(16, file='weights.dat', status='unknown')
       endif
       
-C     READ IN TEMPERATURES
+!     READ IN TEMPERATURES
       do i=1,num_rep
          read(11,*) j,temp
@@ -111 +111 @@
       close(11)
 
-c     nresi:  number of residues
+!     nresi:  number of residues
       nresi=irsml2(1)-irsml1(1)+1
-C     
-C     Initialize variables
+!     
+!     Initialize variables
       do i=1,num_rep      
          acx1(i) = 0.0d0
@@ -132 +132 @@
       dir = dirp(rep_id + 1)
 
-c     _________________________________ Initialize Variables
+!     _________________________________ Initialize Variables
       if(newsta) then
          iold=0
@@ -139 +139 @@
             intem(i) = i
          end do
-c     _________________________________ initialize starting configuration
+!     _________________________________ initialize starting configuration
          if (switch.ne.0) then
             do i=1,nvr
@@ -178 +178 @@
          CALL MPI_BCAST(INODE,num_rep,MPI_INTEGER,0,MPI_COMM_WORLD,IERR)
          CALL MPI_BCAST(YOL,num_rep,MPI_DOUBLE_PRECISION,0,
-     #        MPI_COMM_WORLD,IERR)
+     &        MPI_COMM_WORLD,IERR)
          CALL MPI_BCAST(E_MINP, num_rep, MPI_DOUBLE_PRECISION, 0, 
-     #        MPI_COMM_WORLD, IERR)
+     &        MPI_COMM_WORLD, IERR)
          CALL MPI_BCAST(h_maxp,num_rep,MPI_INTEGER,0,MPI_COMM_WORLD,
-     $        IERR)
+     &        IERR)
       end if
       
@@ -194 +194 @@
          write(*,*) rep_id, yol(rep_id + 1), eol
       endif
-C     Start of simulation
+!     Start of simulation
       write (*,*) '[',rep_id, myrank, beta, partem_comm,
      &            '] Energy before equilibration:', eol
-c     =====================Equilibration by canonical Metropolis
+!     =====================Equilibration by canonical Metropolis
       do nsw=1,nequi
          call metropolis(eol,acz,can_weight)
@@ -204 +204 @@
       write (*,*) '[',rep_id,'] Energy after equilibration:', eol
       call flush(6)
-C     
-C======================Multiple Markov Chains
+!     
+!======================Multiple Markov Chains
       acz = 0
       do nsw=1,nswp
-c------------First ordinary Metropolis 
+!------------First ordinary Metropolis 
          call metropolis(eol,acz,can_weight)
          iold = iold + 1        
@@ -223 +223 @@
             endif
             acz0 = acz
-c     Evaluate RMSD
+!     Evaluate RMSD
             nml = 1
             rmsv = rmsdfun(nml,irsml1(nml),irsml2(nml),ixatp,xatp,yatp, &
      &             zatp,0)
-c            print *,myrank,'received RMSD,energy ',rmsv,eyab,beta
-C     Measure global radius of gyration
+!            print *,myrank,'received RMSD,energy ',rmsv,eyab,beta
+!     Measure global radius of gyration
             call rgyr(0,rgy,ee)  
             rgyp = rgy
-C     Measure Helicity and Sheetness
+!     Measure Helicity and Sheetness
             call helix(nhel,mhel,nbet,mbet)
-C     Measure Number of hydrogen bonds
+!     Measure Number of hydrogen bonds
             mhb = 0
             do i = 1, ntlml
@@ -240 +240 @@
             enddo
             call interhbond(imhb) 
-C     Measure total number of contacts (NCTOT) and number of
-C     native contacts (NCNAT)
+!     Measure total number of contacts (NCTOT) and number of
+!     native contacts (NCNAT)
             call contacts(nctot,ncnat,dham)
-c     Add tracking of lowest energy configuration
+!     Add tracking of lowest energy configuration
             if (eol.lt.e_min) then
-c     Write out configuration
+!     Write out configuration
                i=rep_id+1
                j=inode(i)
@@ -262 +262 @@
                close(15)
             endif
-c     Add tracking of configuration with larges hydrogen contents.
+!     Add tracking of configuration with larges hydrogen contents.
             if ((mhb + imhb).gt.h_max) then
-c     Write out configuration
+!     Write out configuration
                i = rep_id + 1
                j = inode(i)
@@ -282 +282 @@
             endif
 
-C     
-C--------------------Gather measurement data
+!     
+!--------------------Gather measurement data
 ! I only use the master node of each replica for data collection. The
 ! variable partem_comm provides the appropriate communicator.
             if (partem_comm.ne.MPI_COMM_NULL) then
                CALL MPI_GATHER(rmsv,1,MPI_DOUBLE_PRECISION,rmsdp,1,
-     #              MPI_DOUBLE_PRECISION, 0,partem_comm,IERR)
+     &              MPI_DOUBLE_PRECISION, 0,partem_comm,IERR)
                CALL MPI_GATHER(eyab,1,MPI_DOUBLE_PRECISION,eyabp,1,
-     #              MPI_DOUBLE_PRECISION, 0,partem_comm,IERR)
+     &              MPI_DOUBLE_PRECISION, 0,partem_comm,IERR)
                CALL MPI_GATHER(RGYP,1,MPI_DOUBLE_PRECISION,RGYRP,1,
-     #              MPI_DOUBLE_PRECISION, 0,partem_comm,IERR)
+     &              MPI_DOUBLE_PRECISION, 0,partem_comm,IERR)
                CALL MPI_GATHER(NHEL,1,MPI_INTEGER,NHELP,1,MPI_INTEGER,
      &              0,partem_comm,IERR)
@@ -335 +335 @@
 !     &                0,MPI_COMM_WORLD, IERR)                
 
-c     Write trajectory
+!     Write trajectory
                write (18,*) '@@@',iold,inode(rep_id+1)
                call outvbs(0,18)
                write (18,*) '###'
 !                call flush(18)
-c     Write current configuration
+!     Write current configuration
                if ((mod(iold, nsave).eq.0)) then
                   filebase = "conf_0000.var"
@@ -349 +349 @@
             if(rep_id.eq.0.and.myrank.eq.0) then
                randomCount = 0
-c  Update acceptance, temperature wise average of E and E^2 used to calculate
-c  specific heat. 
+!  Update acceptance, temperature wise average of E and E^2 used to calculate
+!  specific heat. 
                do i=1,num_rep
                   j=intem(i)
                   acy(i)=0.0
-c  Above: contents of acy1 are added to acy(i) a few lines down. 
-c  acy1(intem(i)) contains information received from the node at temperature
-c  i, on how many updates have been accepted in node intem(i). Since acz
-c  is not reset to 0 every cycle, acy(i) must be set to 0 here. Else, there 
-c  will be serious double counting and the values of acceptance printed 
-c  will be simply wrong.
+!  Above: contents of acy1 are added to acy(i) a few lines down. 
+!  acy1(intem(i)) contains information received from the node at temperature
+!  i, on how many updates have been accepted in node intem(i). Since acz
+!  is not reset to 0 every cycle, acy(i) must be set to 0 here. Else, there 
+!  will be serious double counting and the values of acceptance printed 
+!  will be simply wrong.
                   e_minpt(i)=e_minp(intem(i))
                end do
@@ -371 +371 @@
 
 
-C     Write measurements to the time series file ts.d
+!     Write measurements to the time series file ts.d
                do i=1,num_rep
                   j=intem(i)
@@ -382 +382 @@
 !                      call flush(14)
                end do
-c     Write the current parallel tempering information into par_R.in
+!     Write the current parallel tempering information into par_R.in
 !               timeLeft = llwrem(2) ! Time left till hard limit
 !               if ((mod(iold, nsave).eq.0).or.(timeLeft.lt.minTimeLeft)
@@ -393 +393 @@
      &                    h_maxp(i)
                   end do
-C     -------------------------- Various statistics of current run
-c               swp=nswp-nequi
+!     -------------------------- Various statistics of current run
+!               swp=nswp-nequi
                   swp=nsw
                   write(13,*) 'Acceptance rate for change of chains:'
@@ -400 +400 @@
                      temp=1.0d0/pbe(k1)/0.00198773
                      write(13,*) temp, acx1(k1)*2.0d0*nmes/swp 
-c  Above: it's the acceptance rate of exchange of replicas. Since a 
-c  replica exchange is attempted only once every nmes sweeps, the 
-c  rate should be normalized with (nmes/swp).
+!  Above: it's the acceptance rate of exchange of replicas. Since a 
+!  replica exchange is attempted only once every nmes sweeps, the 
+!  rate should be normalized with (nmes/swp).
                   end do 
                   write(13,*)
@@ -411 +411 @@
                      geavm(k1) = nmes*eavm(k1)/swp
                      gsph(k1)  = (nmes*sph(k1)/swp-geavm(k1)**2)
-     #                    *beta*beta/nresi
+     &                    *beta*beta/nresi
                      write(13,'(a,2f9.2,i4,f12.3)') 
      &                    'Temperature, Node,local acceptance rate:',
      &                    beta,temp,k,acy(k1)/dble(nsw*nvr)
-c  Above: Changed (nswp-nequi) in the denominator of acceptance as 
-c  acceptance values are initialized to 0 after equilibration cycles are
-c  finished. Note also that since this is being written in the middle of
-c  the simulation, it is normalized to nsw instead of nswp.
+!  Above: Changed (nswp-nequi) in the denominator of acceptance as 
+!  acceptance values are initialized to 0 after equilibration cycles are
+!  finished. Note also that since this is being written in the middle of
+!  the simulation, it is normalized to nsw instead of nswp.
                      write(13,'(a,3f12.2)') 
      &                    'Last Energy, Average Energy, Spec. Heat:', 
@@ -431 +431 @@
                end if
 
-C--------------------Parallel Tempering  update
-c     Swap with right neighbor (odd, even)            
+!--------------------Parallel Tempering  update
+!     Swap with right neighbor (odd, even)            
                if(odd.eq.1) then
                   nu=1
                   no1 = num_rep-1
-c     Swap with left neighbor (even, odd)
+!     Swap with left neighbor (even, odd)
                else
                   nu = 2
@@ -443 +443 @@
                do i=nu,no1,2
                   j=i+1
-c     Periodic bc for swaps
+!     Periodic bc for swaps
                   if(i.eq.num_rep) j=1
                   in=intem(i)
@@ -449 +449 @@
                   wij=exp(-pbe(i)*yol(jn)-pbe(j)*yol(in)
      &                 +pbe(i)*yol(in)+pbe(j)*yol(jn))
-c The random number generator is getting out of sync here, because
-c        the swap is only done on node 0!
+! The random number generator is getting out of sync here, because
+!        the swap is only done on node 0!
 ! Keep track of number of random numbers used.
                   rd=grnd()
                   randomCount = randomCount + 1
-c                  write (16,*) '>', iold, i,j
-c     &            ,pbe(i),yol(in), pbe(j), yol(jn), wij, rd
+!                  write (16,*) '>', iold, i,j
+!     &            ,pbe(i),yol(in), pbe(j), yol(jn), wij, rd
                   if(wij.ge.rd) then
-c Next line: Replica exchange only happens after equilibration, 
-c which takes place outside this loop over nsw. So, I think nsw.gt.nequi
-c is irrelevant for the calculation of acceptance of replica exchanges.
-c /Sandipan
-c                     if(nsw.gt.nequi) 
+! Next line: Replica exchange only happens after equilibration, 
+! which takes place outside this loop over nsw. So, I think nsw.gt.nequi
+! is irrelevant for the calculation of acceptance of replica exchanges.
+! /Sandipan
+!                     if(nsw.gt.nequi) 
                      acx1(i) = acx1(i)+1
                      intem(i) = jn
@@ -469 +469 @@
                   end if
                end do
-c     ---------------- End Loop over nodes which creates a new temperature
-c     map for all nodes, at the node with rank 0. 
-c     
+!     ---------------- End Loop over nodes which creates a new temperature
+!     map for all nodes, at the node with rank 0. 
+!     
                odd = 1 - odd
             end if
-c     End of "if (myrank.eq.0) ...". The block above includes PT update and 
-c     writing of observables into the time series file etc. 
+!     End of "if (myrank.eq.0) ...". The block above includes PT update and 
+!     writing of observables into the time series file etc. 
            
-c     Below: Communicate new temperature-node map to all nodes
+!     Below: Communicate new temperature-node map to all nodes
             CALL MPI_BCAST(INTEM,num_rep,MPI_INTEGER,0,MPI_COMM_WORLD,
      &           IERR)
@@ -486 +486 @@
             CALL MPI_BCAST(H_MAXP,num_rep,MPI_INTEGER,0,MPI_COMM_WORLD,
      &           IERR)
-c Synchronize random number generators for replica 0
+! Synchronize random number generators for replica 0
             if (rep_id.eq.0) then
                CALL MPI_BCAST(randomCount,1,MPI_INTEGER,0,my_mpi_comm,
@@ -507 +507 @@
 
          endif
-c        End of "if (mod(iold,nmes).eq.0) ..."
+!        End of "if (mod(iold,nmes).eq.0) ..."
       end do
-c-----------End Loop over sweeps
-c     
-C     OUTPUT:
-C--------------------For Re-starts:
+!-----------End Loop over sweeps
+!     
+!     OUTPUT:
+!--------------------For Re-starts:
       nu = rep_id + 1
       filebase = "conf_0000.var"
@@ -524 +524 @@
       if (partem_comm.ne.MPI_COMM_NULL) then
          CALL MPI_GATHER(EOL0,1,MPI_DOUBLE_PRECISION,YOL,1,
-     #        MPI_DOUBLE_PRECISION,0,partem_comm,IERR)
+     &        MPI_DOUBLE_PRECISION,0,partem_comm,IERR)
          CALL MPI_GATHER(acz0,1,MPI_DOUBLE_PRECISION,acy1,1,
-     #     MPI_DOUBLE_PRECISION,0,partem_comm,IERR)
+     &     MPI_DOUBLE_PRECISION,0,partem_comm,IERR)
       endif
       
@@ -536 +536 @@
             write(13,*) i,inode(i),intem(i),yol(i),e_minp(i),h_maxp(i)
          end do
-C     -------------------------- Various statistics of current run
+!     -------------------------- Various statistics of current run
          swp=nswp
          write(13,*) 'Acceptance rate for change of chains:'
@@ -559 +559 @@
          end do
          close(13)
-c         close(16)
+!         close(16)
       end if
       close(18)
 
-c     =====================
+!     =====================
       CALL MPI_BARRIER(MPI_COMM_WORLD,IERR)