subroutine read_dist_constr
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
#ifdef MPI
      include 'mpif.h'
#endif
      include 'COMMON.CONTROL'
      include 'COMMON.CHAIN'
      include 'COMMON.IOUNITS'
      include 'COMMON.SBRIDGE'
      integer ifrag_(2,100),ipair_(2,100)
      double precision wfrag_(100),wpair_(100)
      character*500 controlcard
      logical lprn /.true./
      write (iout,*) "Calling read_dist_constr"
      write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
      call flush(iout)
      call card_concat(controlcard)
      call readi(controlcard,"NFRAG",nfrag_,0)
      call readi(controlcard,"NPAIR",npair_,0)
      call readi(controlcard,"NDIST",ndist_,0)
      call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
      call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
      call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
      call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
      call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
      write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
      write (iout,*) "IFRAG"
      do i=1,nfrag_
        write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
      enddo
      write (iout,*) "IPAIR"
      do i=1,npair_
        write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
      enddo
      call flush(iout)
      do i=1,nfrag_
        if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
        if (ifrag_(2,i).gt.nstart_sup+nsup-1)
     &    ifrag_(2,i)=nstart_sup+nsup-1
c        write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
        call flush(iout)
        if (wfrag_(i).gt.0.0d0) then
        do j=ifrag_(1,i),ifrag_(2,i)-1
          do k=j+1,ifrag_(2,i)
            write (iout,*) "j",j," k",k
            ddjk=dist(j,k)
            if (constr_dist.eq.1) then
              nhpb=nhpb+1
              ihpb(nhpb)=j
              jhpb(nhpb)=k
              dhpb(nhpb)=ddjk
              forcon(nhpb)=wfrag_(i) 
            else if (constr_dist.eq.2) then
              if (ddjk.le.dist_cut) then
                nhpb=nhpb+1
                ihpb(nhpb)=j
                jhpb(nhpb)=k
                dhpb(nhpb)=ddjk
                forcon(nhpb)=wfrag_(i) 
              endif
            else
              nhpb=nhpb+1
              ihpb(nhpb)=j
              jhpb(nhpb)=k
              dhpb(nhpb)=ddjk
              forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
            endif
            if (lprn)
     &      write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
     &       nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
          enddo
        enddo
        endif
      enddo
      do i=1,npair_
        if (wpair_(i).gt.0.0d0) then
        ii = ipair_(1,i)
        jj = ipair_(2,i)
        if (ii.gt.jj) then
          itemp=ii
          ii=jj
          jj=itemp
        endif
        do j=ifrag_(1,ii),ifrag_(2,ii)
          do k=ifrag_(1,jj),ifrag_(2,jj)
            nhpb=nhpb+1
            ihpb(nhpb)=j
            jhpb(nhpb)=k
            forcon(nhpb)=wpair_(i)
            dhpb(nhpb)=dist(j,k)
            write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
     &       nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
          enddo
        enddo
        endif
      enddo 
      do i=1,ndist_
        read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
        if (forcon(nhpb+1).gt.0.0d0) then
          nhpb=nhpb+1
          dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
          write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
     &     nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
        endif
      enddo
      call flush(iout)
      return
      end