source/unres/src_MD/newconf.F

   1 #ifdef MPI
   2 ccccccccccccccccccccccccccccccccccccccccccccccccc
   3 ccccccccccccccccccccccccccccccccccccccccccccccccc
   4       subroutine make_var(n,idum,iter_csa)
   5       implicit real*8 (a-h,o-z)
   6       include 'DIMENSIONS'
   7       include 'COMMON.IOUNITS'
   8       include 'COMMON.CSA'
   9       include 'COMMON.BANK'
  10       include 'COMMON.CHAIN'
  11       include 'COMMON.INTERACT'
  12       include 'COMMON.HAIRPIN'
  13       include 'COMMON.VAR'
  14       include 'COMMON.DISTFIT'
  15       include 'COMMON.GEO'
  16       include 'COMMON.CONTROL'
  17       logical nicht_getan,nicht_getan1,fail,lfound
  18       integer nharp,iharp(4,maxres/3),nconf_harp
  19       integer iisucc(mxio)
  20       logical ifused(mxio)
  21       integer nhx_seed(max_seed),ihx_seed(4,maxres/3,max_seed)
  22       integer nhx_use(max_seed),ihx_use(0:4,maxres/3,max_seed)
  23       integer nlx_seed(max_seed),ilx_seed(2,maxres/3,max_seed),
  24      &         nlx_use(max_seed),ilx_use(maxres/3,max_seed)
  25       real ran1,ran2
  26
  27       write (iout,*) 'make_var : nseed=',nseed,'ntry=',n
  28       index=0
  29
  30 c-----------------------------------------
  31       if (n7.gt.0.or.n8.gt.0.or.n9.gt.0.or.n14.gt.0.or.n15.gt.0
  32      &    .or.n16.gt.0.or.n17.gt.0.or.n18.gt.0)
  33      &           call select_frag(n7frag,n8frag,n14frag,
  34      &           n15frag,nbefrag,iter_csa)
  35
  36 c---------------------------------------------------
  37 c N18 - random perturbation of one phi(=gamma) angle in a loop
  38 c
  39       IF (n18.gt.0) THEN
  40       nlx_tot=0
  41       do iters=1,nseed
  42         i1=is(iters)
  43         nlx_seed(iters)=0
  44         do i2=1,n14frag
  45           if (lvar_frag(i2,1).eq.i1) then
  46             nlx_seed(iters)=nlx_seed(iters)+5
  47             ilx_seed(1,nlx_seed(iters),iters)=lvar_frag(i2,2)
  48             ilx_seed(2,nlx_seed(iters),iters)=lvar_frag(i2,3)
  49             ilx_use(nlx_seed(iters),iters)=5
  50           endif
  51         enddo
  52         nlx_use(iters)=nlx_seed(iters)
  53         nlx_tot=nlx_tot+nlx_seed(iters)
  54       enddo
  55
  56       if (nlx_tot .ge. n18*nseed) then
  57         ntot_gen=n18*nseed
  58       else
  59         ntot_gen=(nlx_tot/nseed)*nseed
  60       endif
  61
  62       ngen=0
  63       do while (ngen.lt.ntot_gen)
  64       do iters=1,nseed
  65         iseed=is(iters)
  66         if (nlx_use(iters).gt.0) then
  67           nicht_getan=.true.
  68           do while (nicht_getan)
  69             iih=iran_num(1,nlx_seed(iters))
  70             if (ilx_use(iih,iters).gt.0) then
  71               nicht_getan=.false.
  72               ilx_use(iih,iters)=ilx_use(iih,iters)-1
  73               nlx_use(iters)=nlx_use(iters)-1
  74             endif
  75           enddo
  76           ngen=ngen+1
  77           index=index+1
  78           movenx(index)=18
  79           parent(1,index)=iseed
  80           parent(2,index)=0
  81
  82
  83            if (vdisulf) then
  84              nss_in(index)=bvar_nss(iseed)
  85              do ij=1,nss_in(index)
  86                iss_in(ij,index)=bvar_ss(1,ij,iseed)
  87                jss_in(ij,index)=bvar_ss(2,ij,iseed)
  88              enddo
  89            endif
  90
  91
  92           do k=1,numch
  93             do j=2,nres-1
  94               do i=1,4
  95                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
  96               enddo
  97             enddo
  98           enddo
  99
 100           jr=iran_num(ilx_seed(1,iih,iters),ilx_seed(2,iih,iters))
 101           d=ran_number(-pi,pi)
 102           dihang_in(2,jr-2,1,index)=pinorm(dihang_in(2,jr-2,1,index)+d)
 103
 104
 105           if (ngen.eq.ntot_gen) goto 145
 106         endif
 107       enddo
 108       enddo
 109   145 continue
 110
 111       ENDIF
 112
 113
 114 c-----------------------------------------
 115 c N17 : zip a beta in a seed by forcing one additional p-p contact
 116 c
 117       IF (n17.gt.0) THEN
 118       nhx_tot=0
 119       do iters=1,nseed
 120         i1=is(iters)
 121         nhx_seed(iters)=0
 122         nhx_use(iters)=0
 123         do i2=1,nbefrag
 124           if (avar_frag(i2,1).eq.i1) then
 125             nhx_seed(iters)=nhx_seed(iters)+1
 126             ihx_use(2,nhx_seed(iters),iters)=1
 127             if (avar_frag(i2,5)-avar_frag(i2,3).le.3.and.
 128      &           avar_frag(i2,2).gt.1.and.avar_frag(i2,4).lt.nres) then
 129              ihx_seed(1,nhx_seed(iters),iters)=avar_frag(i2,2)-1
 130              ihx_seed(2,nhx_seed(iters),iters)=avar_frag(i2,4)+1
 131              ihx_use(0,nhx_seed(iters),iters)=1
 132              ihx_use(1,nhx_seed(iters),iters)=0
 133              nhx_use(iters)=nhx_use(iters)+1
 134             else
 135               if (avar_frag(i2,4).gt.avar_frag(i2,5)) then
 136                if (avar_frag(i2,2).gt.1.and.
 137      &                     avar_frag(i2,4).lt.nres) then
 138                 ihx_seed(1,nhx_seed(iters),iters)=avar_frag(i2,2)-1
 139                 ihx_seed(2,nhx_seed(iters),iters)=avar_frag(i2,4)+1
 140                 ihx_use(0,nhx_seed(iters),iters)=1
 141                 ihx_use(1,nhx_seed(iters),iters)=0
 142                 nhx_use(iters)=nhx_use(iters)+1
 143                endif
 144                if (avar_frag(i2,3).lt.nres.and.
 145      &                     avar_frag(i2,5).gt.1) then
 146                 ihx_seed(1,nhx_seed(iters),iters)=avar_frag(i2,3)+1
 147                 ihx_seed(2,nhx_seed(iters),iters)=avar_frag(i2,5)-1
 148                 ihx_use(0,nhx_seed(iters),iters)=
 149      &                    ihx_use(0,nhx_seed(iters),iters)+1
 150                 ihx_use(2,nhx_seed(iters),iters)=0
 151                 nhx_use(iters)=nhx_use(iters)+1
 152                endif
 153               else
 154                if (avar_frag(i2,2).gt.1.and.
 155      &                     avar_frag(i2,4).gt.1) then
 156                 ihx_seed(1,nhx_seed(iters),iters)=avar_frag(i2,2)-1
 157                 ihx_seed(2,nhx_seed(iters),iters)=avar_frag(i2,4)-1
 158                 ihx_use(0,nhx_seed(iters),iters)=1
 159                 ihx_use(1,nhx_seed(iters),iters)=0
 160                 nhx_use(iters)=nhx_use(iters)+1
 161                endif
 162                if (avar_frag(i2,3).lt.nres.and.
 163      &                     avar_frag(i2,5).lt.nres) then
 164                 ihx_seed(1,nhx_seed(iters),iters)=avar_frag(i2,3)+1
 165                 ihx_seed(2,nhx_seed(iters),iters)=avar_frag(i2,5)+1
 166                 ihx_use(0,nhx_seed(iters),iters)=
 167      &                    ihx_use(0,nhx_seed(iters),iters)+1
 168                 ihx_use(2,nhx_seed(iters),iters)=0
 169                 nhx_use(iters)=nhx_use(iters)+1
 170                endif
 171               endif
 172             endif
 173           endif
 174         enddo
 175
 176         nhx_tot=nhx_tot+nhx_use(iters)
 177 cd        write (iout,*) "debug N17",iters,nhx_seed(iters),
 178 cd     &                     nhx_use(iters),nhx_tot
 179       enddo
 180
 181       if (nhx_tot .ge. n17*nseed) then
 182         ntot_gen=n17*nseed
 183       else if (nhx_tot .ge. nseed) then
 184         ntot_gen=(nhx_tot/nseed)*nseed
 185       else
 186         ntot_gen=nhx_tot
 187       endif
 188 cd      write (iout,*) "debug N17==",ntot_gen,nhx_tot,nseed
 189
 190       ngen=0
 191       do while (ngen.lt.ntot_gen)
 192       do iters=1,nseed
 193         iseed=is(iters)
 194         if (nhx_use(iters).gt.0) then
 195 cd          write (iout,*) "debug N17",nhx_use(iters),ngen,ntot_gen
 196 cd          write (iout,*) "debugN17^",
 197 cd     &                    (ihx_use(0,k,iters),k=1,nhx_use(iters))
 198           nicht_getan=.true.
 199           do while (nicht_getan)
 200             iih=iran_num(1,nhx_seed(iters))
 201 cd            write (iout,*) "debugN17^",iih
 202             if (ihx_use(0,iih,iters).gt.0) then
 203               iim=iran_num(1,2)
 204 cd                write (iout,*) "debugN17=",iih,nhx_seed(iters)
 205 cd                write (iout,*) "debugN17-",iim,'##',
 206 cd     &                           (ihx_use(k,iih,iters),k=0,2)
 207 cd                call flush(iout)
 208               do while (ihx_use(iim,iih,iters).eq.1)
 209                 iim=iran_num(1,2)
 210 cd                write (iout,*) "debugN17-",iim,'##',
 211 cd     &                           (ihx_use(k,iih,iters),k=0,2)
 212 cd                call flush(iout)
 213               enddo
 214               nicht_getan=.false.
 215               ihx_use(iim,iih,iters)=1
 216               ihx_use(0,iih,iters)=ihx_use(0,iih,iters)-1
 217               nhx_use(iters)=nhx_use(iters)-1
 218             endif
 219           enddo
 220           ngen=ngen+1
 221           index=index+1
 222           movenx(index)=17
 223           parent(1,index)=iseed
 224           parent(2,index)=0
 225
 226            if (vdisulf) then
 227              nss_in(index)=bvar_nss(iseed)
 228              do ij=1,nss_in(index)
 229                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 230                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 231              enddo
 232            endif
 233
 234           do k=1,numch
 235             do j=2,nres-1
 236               do i=1,4
 237                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
 238               enddo
 239             enddo
 240           enddo
 241
 242           if (iim.eq.1) then
 243             idata(1,index)=ihx_seed(1,iih,iters)
 244             idata(2,index)=ihx_seed(2,iih,iters)
 245           else
 246             idata(1,index)=ihx_seed(3,iih,iters)
 247             idata(2,index)=ihx_seed(4,iih,iters)
 248           endif
 249
 250           if (ngen.eq.ntot_gen) goto 115
 251         endif
 252       enddo
 253       enddo
 254   115 continue
 255       write (iout,*) "N17",n17," ngen/nseed",ngen/nseed,
 256      &                           ngen,nseed
 257
 258
 259       ENDIF
 260 c-----------------------------------------
 261 c N16 : slide non local beta in a seed by +/- 1 or +/- 2
 262 c
 263       IF (n16.gt.0) THEN
 264       nhx_tot=0
 265       do iters=1,nseed
 266         i1=is(iters)
 267         nhx_seed(iters)=0
 268         do i2=1,n7frag
 269           if (bvar_frag(i2,1).eq.i1) then
 270             nhx_seed(iters)=nhx_seed(iters)+1
 271             ihx_seed(1,nhx_seed(iters),iters)=bvar_frag(i2,3)
 272             ihx_seed(2,nhx_seed(iters),iters)=bvar_frag(i2,4)
 273             ihx_seed(3,nhx_seed(iters),iters)=bvar_frag(i2,5)
 274             ihx_seed(4,nhx_seed(iters),iters)=bvar_frag(i2,6)
 275             ihx_use(0,nhx_seed(iters),iters)=4
 276             do i3=1,4
 277               ihx_use(i3,nhx_seed(iters),iters)=0
 278             enddo
 279           endif
 280         enddo
 281         nhx_use(iters)=4*nhx_seed(iters)
 282         nhx_tot=nhx_tot+nhx_seed(iters)
 283 cd        write (iout,*) "debug N16",iters,nhx_seed(iters)
 284       enddo
 285
 286       if (4*nhx_tot .ge. n16*nseed) then
 287         ntot_gen=n16*nseed
 288       else if (4*nhx_tot .ge. nseed) then
 289         ntot_gen=(4*nhx_tot/nseed)*nseed
 290       else
 291         ntot_gen=4*nhx_tot
 292       endif
 293       write (iout,*) "debug N16",ntot_gen,4*nhx_tot,nseed
 294
 295       ngen=0
 296       do while (ngen.lt.ntot_gen)
 297       do iters=1,nseed
 298         iseed=is(iters)
 299         if (nhx_use(iters).gt.0) then
 300           nicht_getan=.true.
 301           do while (nicht_getan)
 302             iih=iran_num(1,nhx_seed(iters))
 303             if (ihx_use(0,iih,iters).gt.0) then
 304               iim=iran_num(1,4)
 305               do while (ihx_use(iim,iih,iters).eq.1)
 306 cd                write (iout,*) iim,
 307 cd     &                ihx_use(0,iih,iters),ihx_use(iim,iih,iters)
 308                 iim=iran_num(1,4)
 309               enddo
 310               nicht_getan=.false.
 311               ihx_use(iim,iih,iters)=1
 312               ihx_use(0,iih,iters)=ihx_use(0,iih,iters)-1
 313               nhx_use(iters)=nhx_use(iters)-1
 314             endif
 315           enddo
 316           ngen=ngen+1
 317           index=index+1
 318           movenx(index)=16
 319           parent(1,index)=iseed
 320           parent(2,index)=0
 321
 322            if (vdisulf) then
 323              nss_in(index)=bvar_nss(iseed)
 324              do ij=1,nss_in(index)
 325                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 326                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 327              enddo
 328            endif
 329
 330           do k=1,numch
 331             do j=2,nres-1
 332               do i=1,4
 333                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
 334               enddo
 335             enddo
 336           enddo
 337
 338           do i=1,4
 339            idata(i,index)=ihx_seed(i,iih,iters)
 340           enddo
 341           idata(5,index)=iim
 342
 343           if (ngen.eq.ntot_gen) goto 116
 344         endif
 345       enddo
 346       enddo
 347   116 continue
 348       write (iout,*) "N16",n16," ngen/nseed",ngen/nseed,
 349      &                           ngen,nseed
 350       ENDIF
 351 c-----------------------------------------
 352 c N15 : copy two 2nd structure elements from 1 or 2 conf. in bank to a seed
 353 c
 354       IF (n15.gt.0) THEN
 355
 356        do iters=1,nseed
 357          iseed=is(iters)
 358          do i=1,mxio
 359            ifused(i)=.false.
 360          enddo
 361
 362          do idummy=1,n15
 363            iter=0
 364    84      continue
 365
 366            iran=0
 367            iif=iran_num(1,n15frag)
 368            do while( (ifused(iif) .or. svar_frag(iif,1).eq.iseed) .and.
 369      &                      iran.le.mxio )
 370             iif=iran_num(1,n15frag)
 371             iran=iran+1
 372            enddo
 373            if(iran.ge.mxio) goto 811
 374
 375            iran=0
 376            iig=iran_num(1,n15frag)
 377            do while( (ifused(iig) .or. svar_frag(iig,1).eq.iseed  .or.
 378      &             .not.(svar_frag(iif,3).lt.svar_frag(iig,2).or.
 379      &                   svar_frag(iig,3).lt.svar_frag(iif,2)) ) .and.
 380      &                      iran.le.mxio )
 381             iig=iran_num(1,n15frag)
 382             iran=iran+1
 383            enddo
 384            if(iran.ge.mxio) goto 811
 385
 386            index=index+1
 387            movenx(index)=15
 388            parent(1,index)=iseed
 389            parent(2,index)=svar_frag(iif,1)
 390            parent(3,index)=svar_frag(iig,1)
 391
 392
 393            if (vdisulf) then
 394              nss_in(index)=bvar_nss(iseed)
 395              do ij=1,nss_in(index)
 396                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 397                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 398              enddo
 399            endif
 400
 401            ifused(iif)=.true.
 402            ifused(iig)=.true.
 403            call newconf_copy(idum,dihang_in(1,1,1,index),
 404      &             svar_frag(iif,1),svar_frag(iif,2),svar_frag(iif,3))
 405
 406            do j=svar_frag(iig,2),svar_frag(iig,3)
 407              do i=1,4
 408               dihang_in(i,j,1,index)=bvar(i,j,1,svar_frag(iig,1))
 409              enddo
 410            enddo
 411
 412
 413            if(iter.lt.10) then
 414             call check_old(icheck,index)
 415             if(icheck.eq.1) then
 416               index=index-1
 417               ifused(iif)=.false.
 418               goto 84
 419             endif
 420            endif
 421
 422   811     continue
 423          enddo
 424        enddo
 425       ENDIF
 426
 427 c-----------------------------------------
 428 c N14 local_move (Maurizio) for loops in a seed
 429 c
 430       IF (n14.gt.0) THEN
 431       nlx_tot=0
 432       do iters=1,nseed
 433         i1=is(iters)
 434         nlx_seed(iters)=0
 435         do i2=1,n14frag
 436           if (lvar_frag(i2,1).eq.i1) then
 437             nlx_seed(iters)=nlx_seed(iters)+3
 438             ilx_seed(1,nlx_seed(iters),iters)=lvar_frag(i2,2)
 439             ilx_seed(2,nlx_seed(iters),iters)=lvar_frag(i2,3)
 440             ilx_use(nlx_seed(iters),iters)=3
 441           endif
 442         enddo
 443         nlx_use(iters)=nlx_seed(iters)
 444         nlx_tot=nlx_tot+nlx_seed(iters)
 445 cd        write (iout,*) "debug N14",iters,nlx_seed(iters)
 446       enddo
 447
 448       if (nlx_tot .ge. n14*nseed) then
 449         ntot_gen=n14*nseed
 450       else
 451         ntot_gen=(nlx_tot/nseed)*nseed
 452       endif
 453 cd      write (iout,*) "debug N14",ntot_gen,n14frag,nseed
 454
 455       ngen=0
 456       do while (ngen.lt.ntot_gen)
 457       do iters=1,nseed
 458         iseed=is(iters)
 459         if (nlx_use(iters).gt.0) then
 460           nicht_getan=.true.
 461           do while (nicht_getan)
 462             iih=iran_num(1,nlx_seed(iters))
 463             if (ilx_use(iih,iters).gt.0) then
 464               nicht_getan=.false.
 465               ilx_use(iih,iters)=ilx_use(iih,iters)-1
 466               nlx_use(iters)=nlx_use(iters)-1
 467             endif
 468           enddo
 469           ngen=ngen+1
 470           index=index+1
 471           movenx(index)=14
 472           parent(1,index)=iseed
 473           parent(2,index)=0
 474
 475           idata(1,index)=ilx_seed(1,iih,iters)
 476           idata(2,index)=ilx_seed(2,iih,iters)
 477
 478
 479            if (vdisulf) then
 480              nss_in(index)=bvar_nss(iseed)
 481              do ij=1,nss_in(index)
 482                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 483                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 484              enddo
 485            endif
 486
 487
 488           do k=1,numch
 489             do j=2,nres-1
 490               do i=1,4
 491                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
 492               enddo
 493             enddo
 494           enddo
 495
 496           if (ngen.eq.ntot_gen) goto 131
 497         endif
 498       enddo
 499       enddo
 500   131 continue
 501 cd      write (iout,*) "N14",n14," ngen/nseed",ngen/nseed,
 502 cd     &                           ngen,nseed
 503
 504       ENDIF
 505 c-----------------------------------------
 506 c N9 : shift a helix in a seed
 507 c
 508       IF (n9.gt.0) THEN
 509       nhx_tot=0
 510       do iters=1,nseed
 511         i1=is(iters)
 512         nhx_seed(iters)=0
 513         do i2=1,n8frag
 514           if (hvar_frag(i2,1).eq.i1) then
 515             nhx_seed(iters)=nhx_seed(iters)+1
 516             ihx_seed(1,nhx_seed(iters),iters)=hvar_frag(i2,2)
 517             ihx_seed(2,nhx_seed(iters),iters)=hvar_frag(i2,3)
 518             ihx_use(0,nhx_seed(iters),iters)=4
 519             do i3=1,4
 520               ihx_use(i3,nhx_seed(iters),iters)=0
 521             enddo
 522           endif
 523         enddo
 524         nhx_use(iters)=4*nhx_seed(iters)
 525         nhx_tot=nhx_tot+nhx_seed(iters)
 526 cd        write (iout,*) "debug N9",iters,nhx_seed(iters)
 527       enddo
 528
 529       if (4*nhx_tot .ge. n9*nseed) then
 530         ntot_gen=n9*nseed
 531       else
 532         ntot_gen=(4*nhx_tot/nseed)*nseed
 533       endif
 534 cd      write (iout,*) "debug N9",ntot_gen,n8frag,nseed
 535
 536       ngen=0
 537       do while (ngen.lt.ntot_gen)
 538       do iters=1,nseed
 539         iseed=is(iters)
 540         if (nhx_use(iters).gt.0) then
 541           nicht_getan=.true.
 542           do while (nicht_getan)
 543             iih=iran_num(1,nhx_seed(iters))
 544             if (ihx_use(0,iih,iters).gt.0) then
 545               iim=iran_num(1,4)
 546               do while (ihx_use(iim,iih,iters).eq.1)
 547 cd                write (iout,*) iim,
 548 cd     &                ihx_use(0,iih,iters),ihx_use(iim,iih,iters)
 549                 iim=iran_num(1,4)
 550               enddo
 551               nicht_getan=.false.
 552               ihx_use(iim,iih,iters)=1
 553               ihx_use(0,iih,iters)=ihx_use(0,iih,iters)-1
 554               nhx_use(iters)=nhx_use(iters)-1
 555             endif
 556           enddo
 557           ngen=ngen+1
 558           index=index+1
 559           movenx(index)=9
 560           parent(1,index)=iseed
 561           parent(2,index)=0
 562
 563            if (vdisulf) then
 564              nss_in(index)=bvar_nss(iseed)
 565              do ij=1,nss_in(index)
 566                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 567                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 568              enddo
 569            endif
 570
 571           do k=1,numch
 572             do j=2,nres-1
 573               do i=1,4
 574                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
 575               enddo
 576             enddo
 577           enddo
 578
 579           jstart=max(nnt,ihx_seed(1,iih,iters)+1)
 580           jend=min(nct,ihx_seed(2,iih,iters))
 581 cd          write (iout,*) "debug N9",iters,iih,jstart,jend
 582           if (iim.eq.1) then
 583             ishift=-2
 584           else if (iim.eq.2) then
 585             ishift=-1
 586           else if (iim.eq.3) then
 587             ishift=1
 588           else if (iim.eq.4) then
 589             ishift=2
 590           else
 591             write (iout,*) 'CHUJ NASTAPIL: iim=',iim
 592             call mpi_abort(mpi_comm_world,ierror,ierrcode)
 593           endif
 594           do j=jstart,jend
 595             if (itype(j).eq.10) then
 596               iang=2
 597             else
 598               iang=4
 599             endif
 600             do i=1,iang
 601               if (j+ishift.ge.nnt.and.j+ishift.le.nct)
 602      &          dihang_in(i,j+ishift,1,index)=bvar(i,j,1,iseed)
 603             enddo
 604           enddo
 605           if (ishift.gt.0) then
 606            do j=0,ishift-1
 607             if (itype(jend+j).eq.10) then
 608               iang=2
 609             else
 610               iang=4
 611             endif
 612             do i=1,iang
 613               if (jend+j.ge.nnt.and.jend+j.le.nct)
 614      &          dihang_in(i,jstart+j,1,index)=bvar(i,jend+j,1,iseed)
 615             enddo
 616            enddo
 617           else
 618            do j=0,-ishift-1
 619             if (itype(jstart+j).eq.10) then
 620               iang=2
 621             else
 622               iang=4
 623             endif
 624             do i=1,iang
 625               if (jend+j.ge.nnt.and.jend+j.le.nct)
 626      &          dihang_in(i,jend+j,1,index)=bvar(i,jstart+j,1,iseed)
 627             enddo
 628            enddo
 629           endif
 630           if (ngen.eq.ntot_gen) goto 133
 631         endif
 632       enddo
 633       enddo
 634   133 continue
 635 cd      write (iout,*) "N9",n9," ngen/nseed",ngen/nseed,
 636 cd     &                           ngen,nseed
 637
 638       ENDIF
 639 c-----------------------------------------
 640 c N8 : copy a helix from bank to seed
 641 c
 642       if (n8.gt.0) then
 643        if (n8frag.lt.n8) then
 644           write (iout,*) "N8: only ",n8frag,'helices'
 645           n8c=n8frag
 646        else
 647           n8c=n8
 648        endif
 649
 650        do iters=1,nseed
 651          iseed=is(iters)
 652          do i=1,mxio
 653            ifused(i)=.false.
 654          enddo
 655
 656
 657          do idummy=1,n8c
 658            iter=0
 659    94      continue
 660            iran=0
 661            iif=iran_num(1,n8frag)
 662            do while( (ifused(iif) .or. hvar_frag(iif,1).eq.iseed) .and.
 663      &                      iran.le.mxio )
 664             iif=iran_num(1,n8frag)
 665             iran=iran+1
 666            enddo
 667
 668            if(iran.ge.mxio) goto 911
 669
 670            index=index+1
 671            movenx(index)=8
 672            parent(1,index)=iseed
 673            parent(2,index)=hvar_frag(iif,1)
 674
 675
 676            if (vdisulf) then
 677              nss_in(index)=bvar_nss(iseed)
 678              do ij=1,nss_in(index)
 679                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 680                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 681              enddo
 682            endif
 683
 684            ifused(iif)=.true.
 685            if (hvar_frag(iif,3)-hvar_frag(iif,2).le.6) then
 686              call newconf_copy(idum,dihang_in(1,1,1,index),
 687      &             hvar_frag(iif,1),hvar_frag(iif,2),hvar_frag(iif,3))
 688            else
 689               ih_start=iran_num(hvar_frag(iif,2),hvar_frag(iif,3)-6)
 690               ih_end=iran_num(ih_start,hvar_frag(iif,3))
 691              call newconf_copy(idum,dihang_in(1,1,1,index),
 692      &             hvar_frag(iif,1),ih_start,ih_end)
 693            endif
 694            iter=iter+1
 695            if(iter.lt.10) then
 696             call check_old(icheck,index)
 697             if(icheck.eq.1) then
 698               index=index-1
 699               ifused(iif)=.false.
 700               goto 94
 701             endif
 702            endif
 703
 704
 705   911     continue
 706
 707          enddo
 708        enddo
 709
 710       endif
 711
 712 c-----------------------------------------
 713 c N7 : copy nonlocal beta fragment from bank to seed
 714 c
 715       if (n7.gt.0) then
 716        if (n7frag.lt.n7) then
 717           write (iout,*) "N7: only ",n7frag,'nonlocal fragments'
 718           n7c=n7frag
 719        else
 720           n7c=n7
 721        endif
 722
 723        do i=1,maxres
 724          do j=1,mxio2
 725             iff_in(i,j)=0
 726          enddo
 727        enddo
 728        index2=0
 729        do i=1,mxio
 730          isend2(i)=0
 731        enddo
 732
 733        do iters=1,nseed
 734          iseed=is(iters)
 735          do i=1,mxio
 736            ifused(i)=.false.
 737          enddo
 738
 739          do idummy=1,n7c
 740            iran=0
 741            iif=iran_num(1,n7frag)
 742            do while( (ifused(iif) .or. bvar_frag(iif,1).eq.iseed) .and.
 743      &                      iran.le.mxio )
 744             iif=iran_num(1,n7frag)
 745             iran=iran+1
 746            enddo
 747
 748 cd       write (*,'(3i5,l,4i5)'),iters,idummy,iif,ifused(iif),
 749 cd     &                    bvar_frag(iif,1),iseed,iran,index2
 750
 751            if(iran.ge.mxio) goto 999
 752            if(index2.ge.mxio2) goto 999
 753
 754            index=index+1
 755            movenx(index)=7
 756            parent(1,index)=iseed
 757            parent(2,index)=bvar_frag(iif,1)
 758            index2=index2+1
 759            isend2(index)=index2
 760            ifused(iif)=.true.
 761
 762            if (vdisulf) then
 763              nss_in(index)=bvar_nss(iseed)
 764              do ij=1,nss_in(index)
 765                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 766                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 767              enddo
 768            endif
 769
 770            do k=1,numch
 771             do j=2,nres-1
 772              do i=1,4
 773               dihang_in2(i,j,k,index2)=bvar(i,j,k,bvar_frag(iif,1))
 774              enddo
 775             enddo
 776            enddo
 777
 778            if (bvar_frag(iif,2).eq.4) then
 779              do i=bvar_frag(iif,3),bvar_frag(iif,4)
 780                iff_in(i,index2)=1
 781              enddo
 782              if (bvar_frag(iif,5).lt.bvar_frag(iif,6)) then
 783 cd               print *,'###',bvar_frag(iif,3),bvar_frag(iif,4),
 784 cd     &                 bvar_frag(iif,5),bvar_frag(iif,6)
 785                do i=bvar_frag(iif,5),bvar_frag(iif,6)
 786                  iff_in(i,index2)=1
 787                enddo
 788              else
 789 cd               print *,'###',bvar_frag(iif,3),bvar_frag(iif,4),
 790 cd     &                 bvar_frag(iif,6),bvar_frag(iif,5)
 791                do i=bvar_frag(iif,6),bvar_frag(iif,5)
 792                  iff_in(i,index2)=1
 793                enddo
 794              endif
 795            endif
 796
 797            do k=1,numch
 798             do j=2,nres-1
 799              do i=1,4
 800               dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
 801              enddo
 802             enddo
 803            enddo
 804
 805
 806   999     continue
 807
 808          enddo
 809        enddo
 810
 811       endif
 812 c-----------------------------------------------
 813 c N6 : copy random continues fragment from bank to seed
 814 c
 815       do iters=1,nseed
 816        iseed=is(iters)
 817        do idummy=1,n6
 818         isize=(is2-is1+1)*ran1(idum)+is1
 819         index=index+1
 820         movenx(index)=6
 821
 822
 823            if (vdisulf) then
 824              nss_in(index)=bvar_nss(iseed)
 825              do ij=1,nss_in(index)
 826                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 827                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 828              enddo
 829            endif
 830
 831         iter=0
 832   104   continue
 833         if(icycle.le.0) then
 834          i1=nconf*  ran1(idum)+1
 835          i1=nbank-nconf+i1
 836         else
 837          i1=nbank*  ran1(idum)+1
 838         endif
 839         if(i1.eq.iseed) goto 104
 840         iter=iter+1
 841         call newconf_residue(idum,dihang_in(1,1,1,index),i1,isize)
 842         parent(1,index)=iseed
 843         parent(2,index)=i1
 844         if(iter.lt.10) then
 845          call check_old(icheck,index)
 846          if(icheck.eq.1) goto 104
 847         endif
 848        enddo
 849       enddo
 850 c-----------------------------------------
 851       if (n3.gt.0.or.n4.gt.0) call gen_hairpin
 852       nconf_harp=0
 853       do iters=1,nseed
 854         if (nharp_seed(iters).gt.0) nconf_harp=nconf_harp+1
 855       enddo
 856 c-----------------------------------------
 857 c N3 : copy hairpin from bank to seed
 858 c
 859       do iters=1,nseed
 860        iseed=is(iters)
 861        nsucc=0
 862        nacc=0
 863        do idummy=1,n3
 864         index=index+1
 865         iter=0
 866   124   continue
 867         if(icycle.le.0) then
 868          i1=nconf*  ran1(idum)+1
 869          i1=nbank-nconf+i1
 870         else
 871          i1=nbank*  ran1(idum)+1
 872         endif
 873         if(i1.eq.iseed) goto 124
 874         do k=1,nsucc
 875           if (i1.eq.iisucc(k).and.nsucc.lt.nconf_harp-1) goto 124
 876         enddo
 877         nsucc=nsucc+1
 878         iisucc(nsucc)=i1
 879         iter=iter+1
 880         call newconf_residue_hairpin(idum,dihang_in(1,1,1,index),
 881      &     i1,fail)
 882         if (fail) then
 883           if (icycle.le.0 .and. nsucc.eq.nconf .or.
 884      &        icycle.gt.0 .and. nsucc.eq.nbank)  then
 885             index=index-1
 886             goto 125
 887           else
 888             goto 124
 889           endif
 890         endif
 891         if(iter.lt.10) then
 892          call check_old(icheck,index)
 893          if(icheck.eq.1) goto 124
 894         endif
 895         movenx(index)=3
 896         parent(1,index)=iseed
 897         parent(2,index)=i1
 898
 899
 900            if (vdisulf) then
 901              nss_in(index)=bvar_nss(iseed)
 902              do ij=1,nss_in(index)
 903                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 904                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 905              enddo
 906            endif
 907
 908         nacc=nacc+1
 909        enddo
 910 c if not enough hairpins, supplement with windows
 911   125  continue
 912 cdd       if (n3.ne.0) write (iout,*) "N3",n3," nsucc",nsucc," nacc",nacc
 913        do idummy=nacc+1,n3
 914         isize=(is2-is1+1)*ran1(idum)+is1
 915         index=index+1
 916         movenx(index)=6
 917         parent(1,index)=iseed
 918         parent(2,index)=i1
 919
 920
 921            if (vdisulf) then
 922              nss_in(index)=bvar_nss(iseed)
 923              do ij=1,nss_in(index)
 924                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 925                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 926              enddo
 927            endif
 928
 929         iter=0
 930   114   continue
 931         if(icycle.le.0) then
 932          i1=nconf*  ran1(idum)+1
 933          i1=nbank-nconf+i1
 934         else
 935          i1=nbank*  ran1(idum)+1
 936         endif
 937         if(i1.eq.iseed) goto 114
 938         iter=iter+1
 939         call newconf_residue(idum,dihang_in(1,1,1,index),i1,isize)
 940         if(iter.lt.10) then
 941          call check_old(icheck,index)
 942          if(icheck.eq.1) goto 114
 943         endif
 944        enddo
 945       enddo
 946 c-----------------------------------------
 947 c N4 : shift a turn in hairpin in seed
 948 c
 949       IF (N4.GT.0) THEN
 950       if (4*nharp_tot .ge. n4*nseed) then
 951         ntot_gen=n4*nseed
 952       else
 953         ntot_gen=(4*nharp_tot/nseed)*nseed
 954       endif
 955       ngen=0
 956       do while (ngen.lt.ntot_gen)
 957       do iters=1,nseed
 958         iseed=is(iters)
 959 c        write (iout,*) 'iters',iters,' iseed',iseed,' nharp_seed',
 960 c     &     nharp_seed(iters),' nharp_use',nharp_use(iters),
 961 c     &     ' ntot_gen',ntot_gen
 962 c        write (iout,*) 'iharp_use(0)',
 963 c     &     (iharp_use(0,k,iters),k=1,nharp_seed(iters))
 964         if (nharp_use(iters).gt.0) then
 965           nicht_getan=.true.
 966           do while (nicht_getan)
 967             iih=iran_num(1,nharp_seed(iters))
 968 c            write (iout,*) 'iih',iih,' iharp_use',
 969 c     &            (iharp_use(k,iih,iters),k=1,4)
 970             if (iharp_use(0,iih,iters).gt.0) then
 971               nicht_getan1=.true.
 972               do while (nicht_getan1)
 973                 iim=iran_num(1,4)
 974                 nicht_getan1=iharp_use(iim,iih,iters).eq.1
 975               enddo
 976               nicht_getan=.false.
 977               iharp_use(iim,iih,iters)=1
 978               iharp_use(0,iih,iters)=iharp_use(0,iih,iters)-1
 979               nharp_use(iters)=nharp_use(iters)-1
 980 cdd              write (iout,'(a16,i3,a5,i2,a10,2i4)')
 981 cdd     &           'N4 selected hairpin',iih,' move',iim,' iharp_seed',
 982 cdd     &            iharp_seed(1,iih,iters),iharp_seed(2,iih,iters)
 983             endif
 984           enddo
 985           ngen=ngen+1
 986           index=index+1
 987           movenx(index)=4
 988           parent(1,index)=iseed
 989           parent(2,index)=0
 990
 991
 992            if (vdisulf) then
 993              nss_in(index)=bvar_nss(iseed)
 994              do ij=1,nss_in(index)
 995                iss_in(ij,index)=bvar_ss(1,ij,iseed)
 996                jss_in(ij,index)=bvar_ss(2,ij,iseed)
 997              enddo
 998            endif
 999
1000           do k=1,numch
1001             do j=2,nres-1
1002               do i=1,4
1003                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1004               enddo
1005             enddo
1006           enddo
1007           jstart=iharp_seed(1,iih,iters)+1
1008           jend=iharp_seed(2,iih,iters)
1009           if (iim.eq.1) then
1010             ishift=-2
1011           else if (iim.eq.2) then
1012             ishift=-1
1013           else if (iim.eq.3) then
1014             ishift=1
1015           else if (iim.eq.4) then
1016             ishift=2
1017           else
1018             write (iout,*) 'CHUJ NASTAPIL: iim=',iim
1019             call mpi_abort(mpi_comm_world,ierror,ierrcode)
1020           endif
1021 c          write (iout,*) 'jstart',jstart,' jend',jend,' ishift',ishift
1022 c          write (iout,*) 'Before turn shift'
1023 c          do j=2,nres-1
1024 c            theta(j+1)=dihang_in(1,j,1,index)
1025 c            phi(j+2)=dihang_in(2,j,1,index)
1026 c            alph(j)=dihang_in(3,j,1,index)
1027 c            omeg(j)=dihang_in(4,j,1,index)
1028 c          enddo
1029 c          call intout
1030           do j=jstart,jend
1031             if (itype(j).eq.10) then
1032               iang=2
1033             else
1034               iang=4
1035             endif
1036             do i=1,iang
1037               if (j+ishift.ge.nnt.and.j+ishift.le.nct)
1038      &          dihang_in(i,j+ishift,1,index)=bvar(i,j,1,iseed)
1039             enddo
1040           enddo
1041 c          write (iout,*) 'After turn shift'
1042 c          do j=2,nres-1
1043 c            theta(j+1)=dihang_in(1,j,1,index)
1044 c            phi(j+2)=dihang_in(2,j,1,index)
1045 c            alph(j)=dihang_in(3,j,1,index)
1046 c            omeg(j)=dihang_in(4,j,1,index)
1047 c          enddo
1048 c          call intout
1049           if (ngen.eq.ntot_gen) goto 135
1050         endif
1051       enddo
1052       enddo
1053 c if not enough hairpins, supplement with windows
1054 c      write (iout,*) 'end of enddo'
1055   135 continue
1056 cdd      write (iout,*) "N4",n4," ngen/nseed",ngen/nseed,
1057 cdd     &                           ngen,nseed
1058       do iters=1,nseed
1059         iseed=is(iters)
1060         do idummy=ngen/nseed+1,n4
1061           isize=(is2-is1+1)*ran1(idum)+is1
1062           index=index+1
1063           movenx(index)=6
1064
1065            if (vdisulf) then
1066              nss_in(index)=bvar_nss(iseed)
1067              do ij=1,nss_in(index)
1068                iss_in(ij,index)=bvar_ss(1,ij,iseed)
1069                jss_in(ij,index)=bvar_ss(2,ij,iseed)
1070              enddo
1071            endif
1072
1073
1074           iter=0
1075   134     continue
1076           if(icycle.le.0) then
1077           i1=nconf*  ran1(idum)+1
1078           i1=nbank-nconf+i1
1079           else
1080             i1=nbank*  ran1(idum)+1
1081           endif
1082           if(i1.eq.iseed) goto 134
1083             iter=iter+1
1084             call newconf_residue(idum,dihang_in(1,1,1,index),i1,isize)
1085             parent(1,index)=iseed
1086             parent(2,index)=i1
1087             if(iter.lt.10) then
1088               call check_old(icheck,index)
1089             if(icheck.eq.1) goto 134
1090           endif
1091         enddo
1092       enddo
1093       ENDIF
1094 c-----------------------------------------
1095 c N5 : copy one residue from bank to seed (normally switched off - use N1)
1096 c
1097       do iters=1,nseed
1098        iseed=is(iters)
1099        isize=1
1100        do i=1,n5
1101         index=index+1
1102         movenx(index)=5
1103
1104            if (vdisulf) then
1105              nss_in(index)=bvar_nss(iseed)
1106              do ij=1,nss_in(index)
1107                iss_in(ij,index)=bvar_ss(1,ij,iseed)
1108                jss_in(ij,index)=bvar_ss(2,ij,iseed)
1109              enddo
1110            endif
1111
1112
1113         iter=0
1114   105   continue
1115         if(icycle.le.0) then
1116          i1=nconf*  ran1(idum)+1
1117          i1=nbank-nconf+i1
1118         else
1119          i1=nbank*  ran1(idum)+1
1120         endif
1121         if(i1.eq.iseed) goto 105
1122         iter=iter+1
1123         call newconf_residue(idum,dihang_in(1,1,1,index),i1,isize)
1124         parent(1,index)=iseed
1125         parent(2,index)=i1
1126         if(iter.lt.10) then
1127          call check_old(icheck,index)
1128          if(icheck.eq.1) goto 105
1129         endif
1130        enddo
1131       enddo
1132 c-----------------------------------------
1133 c N2 : copy backbone of one residue from bank or first bank to seed
1134 c (normally switched off - use N1)
1135 c
1136       do iters=1,nseed
1137        iseed=is(iters)
1138        do i=n2,1,-1
1139         if(icycle.le.0.and.iuse.gt.nconf-irr) then
1140          iseed=ran1(idum)*nconf+1
1141          iseed=nbank-nconf+iseed
1142         endif
1143         index=index+1
1144         movenx(index)=2
1145
1146         if (vdisulf) then
1147              nss_in(index)=bvar_nss(iseed)
1148              do ij=1,nss_in(index)
1149                iss_in(ij,index)=bvar_ss(1,ij,iseed)
1150                jss_in(ij,index)=bvar_ss(2,ij,iseed)
1151              enddo
1152         endif
1153
1154         iter=0
1155   102   i1=  ran1(idum)*nbank+1
1156         if(i1.eq.iseed) goto 102
1157         iter=iter+1
1158         if(icycle.le.0.and.iuse.gt.nconf-irr) then
1159          nran=mod(i-1,nran0)+3
1160          call newconf1arr(idum,dihang_in(1,1,1,index),nran,i1)
1161          parent(1,index)=-iseed
1162          parent(2,index)=-i1
1163         else if(icycle.le.0.and.iters.le.iuse) then
1164          nran=mod(i-1,nran0)+1
1165          call newconf1abr(idum,dihang_in(1,1,1,index),nran,i1)
1166          parent(1,index)=iseed
1167          parent(2,index)=-i1
1168         else
1169          nran=mod(i-1,nran1)+1
1170          if(ran1(idum).lt.0.5) then
1171           call newconf1abr(idum,dihang_in(1,1,1,index),nran,i1)
1172           parent(1,index)=iseed
1173           parent(2,index)=-i1
1174          else
1175           call newconf1abb(idum,dihang_in(1,1,1,index),nran,i1)
1176           parent(1,index)=iseed
1177           parent(2,index)=i1
1178          endif
1179         endif
1180         if(iter.lt.10) then
1181          call check_old(icheck,index)
1182          if(icheck.eq.1) goto 102
1183         endif
1184        enddo
1185       enddo
1186 c-----------------------------------------
1187 c N1 : copy backbone or sidechain of one residue from bank or
1188 c      first bank to seed
1189 c
1190       do iters=1,nseed
1191        iseed=is(iters)
1192        do i=n1,1,-1
1193         if(icycle.le.0.and.iuse.gt.nconf-irr) then
1194          iseed=ran1(idum)*nconf+1
1195          iseed=nbank-nconf+iseed
1196         endif
1197         index=index+1
1198         movenx(index)=1
1199
1200         if (vdisulf) then
1201              nss_in(index)=bvar_nss(iseed)
1202              do ij=1,nss_in(index)
1203                iss_in(ij,index)=bvar_ss(1,ij,iseed)
1204                jss_in(ij,index)=bvar_ss(2,ij,iseed)
1205              enddo
1206         endif
1207
1208         iter=0
1209   101   i1=  ran1(idum)*nbank+1
1210
1211         if(i1.eq.iseed) goto 101
1212         iter=iter+1
1213         if(icycle.le.0.and.iuse.gt.nconf-irr) then
1214          nran=mod(i-1,nran0)+3
1215          call newconf1rr(idum,dihang_in(1,1,1,index),nran,i1)
1216          parent(1,index)=-iseed
1217          parent(2,index)=-i1
1218         else if(icycle.le.0.and.iters.le.iuse) then
1219          nran=mod(i-1,nran0)+1
1220          call newconf1br(idum,dihang_in(1,1,1,index),nran,i1)
1221          parent(1,index)=iseed
1222          parent(2,index)=-i1
1223         else
1224          nran=mod(i-1,nran1)+1
1225          if(ran1(idum).lt.0.5) then
1226           call newconf1br(idum,dihang_in(1,1,1,index),nran,i1)
1227           parent(1,index)=iseed
1228           parent(2,index)=-i1
1229          else
1230           call newconf1bb(idum,dihang_in(1,1,1,index),nran,i1)
1231           parent(1,index)=iseed
1232           parent(2,index)=i1
1233          endif
1234         endif
1235         if(iter.lt.10) then
1236          call check_old(icheck,index)
1237          if(icheck.eq.1) goto 101
1238         endif
1239        enddo
1240       enddo
1241 c-----------------------------------------
1242 c N0 just all seeds
1243 c
1244       IF (n0.gt.0) THEN
1245       do iters=1,nseed
1246        iseed=is(iters)
1247        index=index+1
1248        movenx(index)=0
1249        parent(1,index)=iseed
1250        parent(2,index)=0
1251
1252        if (vdisulf) then
1253              nss_in(index)=bvar_nss(iseed)
1254              do ij=1,nss_in(index)
1255                iss_in(ij,index)=bvar_ss(1,ij,iseed)
1256                jss_in(ij,index)=bvar_ss(2,ij,iseed)
1257              enddo
1258        endif
1259
1260        do k=1,numch
1261         do j=2,nres-1
1262          do i=1,4
1263           dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1264          enddo
1265         enddo
1266        enddo
1267       enddo
1268       ENDIF
1269 c-----------------------------------------
1270       if (vdisulf) then
1271       do iters=1,nseed
1272        iseed=is(iters)
1273
1274         do k=1,numch
1275           do j=2,nres-1
1276             theta(j+1)=bvar(1,j,k,iseed)
1277             phi(j+2)=bvar(2,j,k,iseed)
1278             alph(j)=bvar(3,j,k,iseed)
1279             omeg(j)=bvar(4,j,k,iseed)
1280           enddo
1281         enddo
1282         call chainbuild
1283
1284 cd       write(iout,*) 'makevar DYNSS',iseed,'#',bvar_ns(iseed),
1285 cd     &     (bvar_s(k,iseed),k=1,bvar_ns(iseed)),
1286 cd     &     bvar_nss(iseed),
1287 cd     &     (bvar_ss(1,k,iseed)-nres,'-',
1288 cd     &      bvar_ss(2,k,iseed)-nres,k=1,bvar_nss(iseed))
1289
1290        do i1=1,bvar_ns(iseed)
1291 c
1292 c N10 fussion of free halfcysteines in seed
1293 c      first select CYS with distance < 7A
1294 c
1295          do j1=i1+1,bvar_ns(iseed)
1296            if (dist(bvar_s(i1,iseed)+nres,bvar_s(j1,iseed)+nres)
1297      &         .lt.7.0.and.
1298      &         iabs(bvar_s(i1,iseed)-bvar_s(j1,iseed)).gt.3) then
1299
1300              index=index+1
1301              movenx(index)=10
1302              parent(1,index)=iseed
1303              parent(2,index)=0
1304              do ij=1,bvar_nss(iseed)
1305                iss_in(ij,index)=bvar_ss(1,ij,iseed)
1306                jss_in(ij,index)=bvar_ss(2,ij,iseed)
1307              enddo
1308              ij=bvar_nss(iseed)+1
1309              nss_in(index)=ij
1310              iss_in(ij,index)=bvar_s(i1,iseed)+nres
1311              jss_in(ij,index)=bvar_s(j1,iseed)+nres
1312
1313 cd             write(iout,*) 'makevar NSS0',index,
1314 cd     &   dist(bvar_s(i1,iseed)+nres,bvar_s(j1,iseed)+nres),
1315 cd     &   nss_in(index),iss_in(ij,index)-nres,'-',jss_in(ij,index)-nres
1316
1317              do k=1,numch
1318               do j=2,nres-1
1319                do i=1,4
1320                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1321                enddo
1322               enddo
1323              enddo
1324
1325            endif
1326          enddo
1327 c
1328 c N11 type I transdisulfidation
1329 c
1330          do j1=1,bvar_nss(iseed)
1331            if (dist(bvar_s(i1,iseed)+nres,bvar_ss(1,j1,iseed))
1332      &         .lt.7.0.and.
1333      &         iabs(bvar_s(i1,iseed)-(bvar_ss(1,j1,iseed)-nres))
1334      &         .gt.3) then
1335
1336              index=index+1
1337              movenx(index)=11
1338              parent(1,index)=iseed
1339              parent(2,index)=0
1340              do ij=1,bvar_nss(iseed)
1341                if (ij.ne.j1) then
1342                 iss_in(ij,index)=bvar_ss(1,ij,iseed)
1343                 jss_in(ij,index)=bvar_ss(2,ij,iseed)
1344                endif
1345              enddo
1346              nss_in(index)=bvar_nss(iseed)
1347              iss_in(j1,index)=bvar_s(i1,iseed)+nres
1348              jss_in(j1,index)=bvar_ss(1,j1,iseed)
1349              if (iss_in(j1,index).gt.jss_in(j1,index)) then
1350                iss_in(j1,index)=bvar_ss(1,j1,iseed)
1351                jss_in(j1,index)=bvar_s(i1,iseed)+nres
1352              endif
1353
1354 cd             write(iout,*) 'makevar NSS1 #1',index,
1355 cd     &       bvar_s(i1,iseed),bvar_ss(1,j1,iseed)-nres,
1356 cd     &       dist(bvar_s(i1,iseed)+nres,bvar_ss(1,j1,iseed)),
1357 cd     &       (iss_in(ij,index)-nres,'-',jss_in(ij,index)-nres,
1358 cd     &       ij=1,nss_in(index))
1359
1360              do k=1,numch
1361               do j=2,nres-1
1362                do i=1,4
1363                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1364                enddo
1365               enddo
1366              enddo
1367            endif
1368            if (dist(bvar_s(i1,iseed)+nres,bvar_ss(2,j1,iseed))
1369      &         .lt.7.0.and.
1370      &         iabs(bvar_s(i1,iseed)-(bvar_ss(2,j1,iseed)-nres))
1371      &         .gt.3) then
1372
1373              index=index+1
1374              movenx(index)=11
1375              parent(1,index)=iseed
1376              parent(2,index)=0
1377              do ij=1,bvar_nss(iseed)
1378                if (ij.ne.j1) then
1379                 iss_in(ij,index)=bvar_ss(1,ij,iseed)
1380                 jss_in(ij,index)=bvar_ss(2,ij,iseed)
1381                endif
1382              enddo
1383              nss_in(index)=bvar_nss(iseed)
1384              iss_in(j1,index)=bvar_s(i1,iseed)+nres
1385              jss_in(j1,index)=bvar_ss(2,j1,iseed)
1386              if (iss_in(j1,index).gt.jss_in(j1,index)) then
1387                iss_in(j1,index)=bvar_ss(2,j1,iseed)
1388                jss_in(j1,index)=bvar_s(i1,iseed)+nres
1389              endif
1390
1391
1392 cd             write(iout,*) 'makevar NSS1 #2',index,
1393 cd     &       bvar_s(i1,iseed),bvar_ss(2,j1,iseed)-nres,
1394 cd     &       dist(bvar_s(i1,iseed)+nres,bvar_ss(2,j1,iseed)),
1395 cd     &       (iss_in(ij,index)-nres,'-',jss_in(ij,index)-nres,
1396 cd     &       ij=1,nss_in(index))
1397
1398              do k=1,numch
1399               do j=2,nres-1
1400                do i=1,4
1401                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1402                enddo
1403               enddo
1404              enddo
1405
1406            endif
1407          enddo
1408        enddo
1409
1410 c
1411 c N12 type II transdisulfidation
1412 c
1413        do i1=1,bvar_nss(iseed)
1414          do j1=i1+1,bvar_nss(iseed)
1415             if (dist(bvar_ss(1,i1,iseed),bvar_ss(1,j1,iseed))
1416      &          .lt.7.0.and.
1417      &          dist(bvar_ss(2,i1,iseed),bvar_ss(2,j1,iseed))
1418      &          .lt.7.0.and.
1419      &          iabs(bvar_ss(1,i1,iseed)-bvar_ss(1,j1,iseed))
1420      &          .gt.3.and.
1421      &          iabs(bvar_ss(2,i1,iseed)-bvar_ss(2,j1,iseed))
1422      &          .gt.3) then
1423              index=index+1
1424              movenx(index)=12
1425              parent(1,index)=iseed
1426              parent(2,index)=0
1427              do ij=1,bvar_nss(iseed)
1428                if (ij.ne.i1 .and. ij.ne.j1) then
1429                 iss_in(ij,index)=bvar_ss(1,ij,iseed)
1430                 jss_in(ij,index)=bvar_ss(2,ij,iseed)
1431                endif
1432              enddo
1433              nss_in(index)=bvar_nss(iseed)
1434              iss_in(i1,index)=bvar_ss(1,i1,iseed)
1435              jss_in(i1,index)=bvar_ss(1,j1,iseed)
1436              if (iss_in(i1,index).gt.jss_in(i1,index)) then
1437                iss_in(i1,index)=bvar_ss(1,j1,iseed)
1438                jss_in(i1,index)=bvar_ss(1,i1,iseed)
1439              endif
1440              iss_in(j1,index)=bvar_ss(2,i1,iseed)
1441              jss_in(j1,index)=bvar_ss(2,j1,iseed)
1442              if (iss_in(j1,index).gt.jss_in(j1,index)) then
1443                iss_in(j1,index)=bvar_ss(2,j1,iseed)
1444                jss_in(j1,index)=bvar_ss(2,i1,iseed)
1445              endif
1446
1447
1448 cd             write(iout,*) 'makevar NSS2 #1',index,
1449 cd     &       bvar_ss(1,i1,iseed)-nres,bvar_ss(1,j1,iseed)-nres,
1450 cd     &       dist(bvar_ss(1,i1,iseed),bvar_ss(1,j1,iseed)),
1451 cd     &       bvar_ss(2,i1,iseed)-nres,bvar_ss(2,j1,iseed)-nres,
1452 cd     &       dist(bvar_ss(2,i1,iseed),bvar_ss(2,j1,iseed)),
1453 cd     &       (iss_in(ij,index)-nres,'-',jss_in(ij,index)-nres,
1454 cd     &       ij=1,nss_in(index))
1455
1456              do k=1,numch
1457               do j=2,nres-1
1458                do i=1,4
1459                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1460                enddo
1461               enddo
1462              enddo
1463
1464            endif
1465
1466             if (dist(bvar_ss(1,i1,iseed),bvar_ss(2,j1,iseed))
1467      &          .lt.7.0.and.
1468      &          dist(bvar_ss(2,i1,iseed),bvar_ss(1,j1,iseed))
1469      &          .lt.7.0.and.
1470      &          iabs(bvar_ss(1,i1,iseed)-bvar_ss(2,j1,iseed))
1471      &          .gt.3.and.
1472      &          iabs(bvar_ss(2,i1,iseed)-bvar_ss(1,j1,iseed))
1473      &          .gt.3) then
1474              index=index+1
1475              movenx(index)=12
1476              parent(1,index)=iseed
1477              parent(2,index)=0
1478              do ij=1,bvar_nss(iseed)
1479                if (ij.ne.i1 .and. ij.ne.j1) then
1480                 iss_in(ij,index)=bvar_ss(1,ij,iseed)
1481                 jss_in(ij,index)=bvar_ss(2,ij,iseed)
1482                endif
1483              enddo
1484              nss_in(index)=bvar_nss(iseed)
1485              iss_in(i1,index)=bvar_ss(1,i1,iseed)
1486              jss_in(i1,index)=bvar_ss(2,j1,iseed)
1487              if (iss_in(i1,index).gt.jss_in(i1,index)) then
1488                iss_in(i1,index)=bvar_ss(2,j1,iseed)
1489                jss_in(i1,index)=bvar_ss(1,i1,iseed)
1490              endif
1491              iss_in(j1,index)=bvar_ss(2,i1,iseed)
1492              jss_in(j1,index)=bvar_ss(1,j1,iseed)
1493              if (iss_in(j1,index).gt.jss_in(j1,index)) then
1494                iss_in(j1,index)=bvar_ss(1,j1,iseed)
1495                jss_in(j1,index)=bvar_ss(2,i1,iseed)
1496              endif
1497
1498
1499 cd             write(iout,*) 'makevar NSS2 #2',index,
1500 cd     &       bvar_ss(1,i1,iseed)-nres,bvar_ss(2,j1,iseed)-nres,
1501 cd     &       dist(bvar_ss(1,i1,iseed),bvar_ss(2,j1,iseed)),
1502 cd     &       bvar_ss(2,i1,iseed)-nres,bvar_ss(1,j1,iseed)-nres,
1503 cd     &       dist(bvar_ss(2,i1,iseed),bvar_ss(1,j1,iseed)),
1504 cd     &       (iss_in(ij,index)-nres,'-',jss_in(ij,index)-nres,
1505 cd     &       ij=1,nss_in(index))
1506
1507              do k=1,numch
1508               do j=2,nres-1
1509                do i=1,4
1510                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1511                enddo
1512               enddo
1513              enddo
1514
1515            endif
1516
1517
1518          enddo
1519        enddo
1520 c
1521 c N13 removal of disulfide bond
1522 c
1523       if (bvar_nss(iseed).gt.0) then
1524         i1=bvar_nss(iseed)*ran1(idum)+1
1525
1526              index=index+1
1527              movenx(index)=13
1528              parent(1,index)=iseed
1529              parent(2,index)=0
1530              ij=0
1531              do j1=1,bvar_nss(iseed)
1532               if (j1.ne.i1) then
1533                ij=ij+1
1534                iss_in(ij,index)=bvar_ss(1,j1,iseed)
1535                jss_in(ij,index)=bvar_ss(2,j1,iseed)
1536               endif
1537              enddo
1538              nss_in(index)=bvar_nss(iseed)-1
1539
1540 cd             write(iout,*) 'NSS3',index,i1,
1541 cd     &   bvar_ss(1,i1,iseed)-nres,'=',bvar_ss(2,i1,iseed)-nres,'#',
1542 cd     &   (iss_in(ij,index)-nres,'-',jss_in(ij,index)-nres,
1543 cd     &   ij=1,nss_in(index))
1544
1545              do k=1,numch
1546               do j=2,nres-1
1547                do i=1,4
1548                 dihang_in(i,j,k,index)=bvar(i,j,k,iseed)
1549                enddo
1550               enddo
1551              enddo
1552
1553       endif
1554
1555       enddo
1556       endif
1557 c-----------------------------------------
1558
1559
1560
1561       if(index.ne.n) write(iout,*)'make_var : ntry=',index
1562
1563       n=index
1564 cd      do ii=1,n
1565 cd        write (istat,*) "======== ii=",ii," the dihang array"
1566 cd        do i=1,nres
1567 cd       write (istat,'(i5,4f15.5)') i,(dihang_in(k,i,1,ii)*rad2deg,k=1,4)
1568 cd        enddo
1569 cd      enddo
1570       return
1571       end
1572 ccccccccccccccccccccccccccccccccccccccccccccccccc
1573 ccccccccccccccccccccccccccccccccccccccccccccccccc
1574       subroutine check_old(icheck,n)
1575       implicit real*8 (a-h,o-z)
1576       include 'DIMENSIONS'
1577       include 'COMMON.CSA'
1578       include 'COMMON.BANK'
1579       include 'COMMON.CHAIN'
1580       include 'COMMON.GEO'
1581       data ctdif /10./
1582       data ctdiff /60./
1583
1584       i1=n
1585       do i2=1,n-1
1586        diff=0.d0
1587        do m=1,numch
1588         do j=2,nres-1
1589          do i=1,4
1590           dif=rad2deg*dabs(dihang_in(i,j,m,i1)-dihang_in(i,j,m,i2))
1591           if(dif.gt.180.0) dif=360.0-dif
1592           if(dif.gt.ctdif) goto 100
1593           diff=diff+dif
1594           if(diff.gt.ctdiff) goto 100
1595          enddo
1596         enddo
1597        enddo
1598        icheck=1
1599        return
1600   100  continue
1601       enddo
1602
1603       icheck=0
1604
1605       return
1606       end
1607 ccccccccccccccccccccccccccccccccccccccccccccccccc
1608 ccccccccccccccccccccccccccccccccccccccccccccccccc
1609       subroutine newconf1rr(idum,vvar,nran,i1)
1610       implicit real*8 (a-h,o-z)
1611       include 'DIMENSIONS'
1612       include 'COMMON.IOUNITS'
1613       include 'COMMON.CSA'
1614       include 'COMMON.BANK'
1615       include 'COMMON.CHAIN'
1616       include 'COMMON.GEO'
1617       real ran1,ran2
1618       dimension vvar(mxang,maxres,mxch),iold(ntotal)
1619       ctdif=10.
1620
1621       do k=1,numch
1622        do j=2,nres-1
1623         do i=1,4
1624          vvar(i,j,k)=rvar(i,j,k,iseed)
1625         enddo
1626        enddo
1627       enddo
1628
1629       do index=1,nran
1630        iold(index) = 0
1631       enddo
1632
1633       number=ntotgr
1634
1635       iter=0
1636       do index=1,nran
1637    10  iran=  ran1(idum)*number+1
1638        if(iter.gt.number) return
1639        iter=iter+1
1640        if(iter.eq.1) goto 11
1641        do ind=1,index-1
1642         if(iran.eq.iold(ind)) goto 10
1643        enddo
1644    11  continue
1645
1646        do ind=1,ngroup(iran)
1647         i=igroup(1,ind,iran)
1648         j=igroup(2,ind,iran)
1649         k=igroup(3,ind,iran)
1650         dif=rad2deg*dabs(vvar(i,j,k)-rvar(i,j,k,i1))
1651         if(dif.gt.180.) dif=360.-dif
1652         if(dif.gt.ctdif) goto 20
1653        enddo
1654        if(iter.gt.number) goto 20
1655        goto 10
1656    20  continue
1657        do ind=1,ngroup(iran)
1658         i=igroup(1,ind,iran)
1659         j=igroup(2,ind,iran)
1660         k=igroup(3,ind,iran)
1661         vvar(i,j,k)=rvar(i,j,k,i1)
1662        enddo
1663        iold(index)=iran
1664       enddo
1665
1666       return
1667       end
1668 ccccccccccccccccccccccccccccccccccccccccccccccccc
1669 ccccccccccccccccccccccccccccccccccccccccccccccccc
1670       subroutine newconf1br(idum,vvar,nran,i1)
1671       implicit real*8 (a-h,o-z)
1672       include 'DIMENSIONS'
1673       include 'COMMON.IOUNITS'
1674       include 'COMMON.CSA'
1675       include 'COMMON.BANK'
1676       include 'COMMON.CHAIN'
1677       include 'COMMON.GEO'
1678       include 'COMMON.TORCNSTR'
1679       include 'COMMON.CONTROL'
1680       real ran1,ran2
1681       dimension vvar(mxang,maxres,mxch),iold(ntotal)
1682       ctdif=10.
1683
1684       do k=1,numch
1685        do j=2,nres-1
1686         do i=1,4
1687          vvar(i,j,k)=bvar(i,j,k,iseed)
1688         enddo
1689        enddo
1690       enddo
1691
1692       do index=1,nran
1693        iold(index) = 0
1694       enddo
1695
1696       number=ntotgr
1697
1698       iter=0
1699       do index=1,nran
1700    10  iran=  ran1(idum)*number+1
1701        if(i2ndstr.gt.0) then
1702         rtmp=ran1(idum)
1703         if(rtmp.le.rdih_bias) then
1704          i=0
1705          do j=1,ndih_nconstr
1706           if(igroup(2,1,iran).eq.idih_nconstr(j))i=j
1707          enddo
1708          if(i.eq.0) then
1709           juhc=0
1710 4321      juhc=juhc+1
1711           iran=  ran1(idum)*number+1
1712           i=0
1713           do j=1,ndih_nconstr
1714            if(igroup(2,1,iran).eq.idih_nconstr(j))i=j
1715           enddo
1716           if(i.eq.0.or.juhc.lt.1000)goto 4321
1717           if(juhc.eq.1000) then
1718            print *, 'move 6 : failed to find unconstrained group'
1719            write(iout,*) 'move 6 : failed to find unconstrained group'
1720           endif
1721          endif
1722         endif
1723        endif
1724        if(iter.gt.number) return
1725        iter=iter+1
1726        if(iter.eq.1) goto 11
1727        do ind=1,index-1
1728         if(iran.eq.iold(ind)) goto 10
1729        enddo
1730    11  continue
1731
1732        do ind=1,ngroup(iran)
1733         i=igroup(1,ind,iran)
1734         j=igroup(2,ind,iran)
1735         k=igroup(3,ind,iran)
1736         dif=rad2deg*dabs(vvar(i,j,k)-rvar(i,j,k,i1))
1737         if(dif.gt.180.) dif=360.-dif
1738         if(dif.gt.ctdif) goto 20
1739        enddo
1740        if(iter.gt.number) goto 20
1741        goto 10
1742    20  continue
1743        do ind=1,ngroup(iran)
1744         i=igroup(1,ind,iran)
1745         j=igroup(2,ind,iran)
1746         k=igroup(3,ind,iran)
1747         vvar(i,j,k)=rvar(i,j,k,i1)
1748        enddo
1749        iold(index)=iran
1750       enddo
1751
1752       return
1753       end
1754 ccccccccccccccccccccccccccccccccccccccccccccccccc
1755 ccccccccccccccccccccccccccccccccccccccccccccccccc
1756       subroutine newconf1bb(idum,vvar,nran,i1)
1757       implicit real*8 (a-h,o-z)
1758       include 'DIMENSIONS'
1759       include 'COMMON.IOUNITS'
1760       include 'COMMON.CSA'
1761       include 'COMMON.BANK'
1762       include 'COMMON.CHAIN'
1763       include 'COMMON.GEO'
1764       real ran1,ran2
1765       dimension vvar(mxang,maxres,mxch),iold(ntotal)
1766       ctdif=10.
1767
1768       do k=1,numch
1769        do j=2,nres-1
1770         do i=1,4
1771          vvar(i,j,k)=bvar(i,j,k,iseed)
1772         enddo
1773        enddo
1774       enddo
1775
1776       do index=1,nran
1777        iold(index) = 0
1778       enddo
1779
1780       number=ntotgr
1781
1782       iter=0
1783       do index=1,nran
1784    10  iran=  ran1(idum)*number+1
1785        if(iter.gt.number) return
1786        iter=iter+1
1787        if(iter.eq.1) goto 11
1788        do ind=1,index-1
1789         if(iran.eq.iold(ind)) goto 10
1790        enddo
1791    11  continue
1792
1793        do ind=1,ngroup(iran)
1794         i=igroup(1,ind,iran)
1795         j=igroup(2,ind,iran)
1796         k=igroup(3,ind,iran)
1797         dif=rad2deg*dabs(vvar(i,j,k)-bvar(i,j,k,i1))
1798         if(dif.gt.180.) dif=360.-dif
1799         if(dif.gt.ctdif) goto 20
1800        enddo
1801        if(iter.gt.number) goto 20
1802        goto 10
1803    20  continue
1804        do ind=1,ngroup(iran)
1805         i=igroup(1,ind,iran)
1806         j=igroup(2,ind,iran)
1807         k=igroup(3,ind,iran)
1808         vvar(i,j,k)=bvar(i,j,k,i1)
1809        enddo
1810        iold(index)=iran
1811       enddo
1812
1813       return
1814       end
1815 ccccccccccccccccccccccccccccccccccccccccccccccccc
1816 ccccccccccccccccccccccccccccccccccccccccccccccccc
1817       subroutine newconf1arr(idum,vvar,nran,i1)
1818       implicit real*8 (a-h,o-z)
1819       include 'DIMENSIONS'
1820       include 'COMMON.IOUNITS'
1821       include 'COMMON.CSA'
1822       include 'COMMON.BANK'
1823       include 'COMMON.CHAIN'
1824       include 'COMMON.GEO'
1825       real ran1,ran2
1826       dimension vvar(mxang,maxres,mxch),iold(ntotal)
1827       ctdif=10.
1828
1829       do k=1,numch
1830        do j=2,nres-1
1831         do i=1,4
1832          vvar(i,j,k)=rvar(i,j,k,iseed)
1833         enddo
1834        enddo
1835       enddo
1836
1837       do index=1,nran
1838        iold(index) = 0
1839       enddo
1840
1841       number=nres-2
1842
1843       iter=0
1844       do index=1,nran
1845    10  iran=  ran1(idum)*number+1
1846        if(iter.gt.number) return
1847        iter=iter+1
1848        if(iter.eq.1) goto 11
1849        do ind=1,index-1
1850         if(iran.eq.iold(ind)) goto 10
1851        enddo
1852    11  continue
1853
1854        do ind=1,ngroup(iran)
1855         i=igroup(1,ind,iran)
1856         j=igroup(2,ind,iran)
1857         k=igroup(3,ind,iran)
1858         dif=rad2deg*dabs(vvar(i,j,k)-rvar(i,j,k,i1))
1859         if(dif.gt.180.) dif=360.-dif
1860         if(dif.gt.ctdif) goto 20
1861        enddo
1862        if(iter.gt.number) goto 20
1863        goto 10
1864    20  continue
1865        do ind=1,ngroup(iran)
1866         i=igroup(1,ind,iran)
1867         j=igroup(2,ind,iran)
1868         k=igroup(3,ind,iran)
1869         vvar(i,j,k)=rvar(i,j,k,i1)
1870        enddo
1871        iold(index)=iran
1872       enddo
1873
1874       return
1875       end
1876 ccccccccccccccccccccccccccccccccccccccccccccccccc
1877 ccccccccccccccccccccccccccccccccccccccccccccccccc
1878       subroutine newconf1abr(idum,vvar,nran,i1)
1879       implicit real*8 (a-h,o-z)
1880       include 'DIMENSIONS'
1881       include 'COMMON.IOUNITS'
1882       include 'COMMON.CSA'
1883       include 'COMMON.BANK'
1884       include 'COMMON.CHAIN'
1885       include 'COMMON.GEO'
1886       include 'COMMON.TORCNSTR'
1887       include 'COMMON.CONTROL'
1888       real ran1,ran2
1889       dimension vvar(mxang,maxres,mxch),iold(ntotal)
1890       ctdif=10.
1891
1892       do k=1,numch
1893        do j=2,nres-1
1894         do i=1,4
1895          vvar(i,j,k)=bvar(i,j,k,iseed)
1896         enddo
1897        enddo
1898       enddo
1899
1900       do index=1,nran
1901        iold(index) = 0
1902       enddo
1903
1904       number=nres-2
1905
1906       iter=0
1907       do index=1,nran
1908    10  iran=  ran1(idum)*number+1
1909       if(i2ndstr.gt.0) then
1910        rtmp=ran1(idum)
1911        if(rtmp.le.rdih_bias) then
1912         iran=ran1(idum)*ndih_nconstr+1
1913         iran=idih_nconstr(iran)
1914        endif
1915       endif
1916        if(iter.gt.number) return
1917        iter=iter+1
1918        if(iter.eq.1) goto 11
1919        do ind=1,index-1
1920         if(iran.eq.iold(ind)) goto 10
1921        enddo
1922    11  continue
1923
1924        do ind=1,ngroup(iran)
1925         i=igroup(1,ind,iran)
1926         j=igroup(2,ind,iran)
1927         k=igroup(3,ind,iran)
1928         dif=rad2deg*dabs(vvar(i,j,k)-rvar(i,j,k,i1))
1929         if(dif.gt.180.) dif=360.-dif
1930         if(dif.gt.ctdif) goto 20
1931        enddo
1932        if(iter.gt.number) goto 20
1933        goto 10
1934    20  continue
1935        do ind=1,ngroup(iran)
1936         i=igroup(1,ind,iran)
1937         j=igroup(2,ind,iran)
1938         k=igroup(3,ind,iran)
1939         vvar(i,j,k)=rvar(i,j,k,i1)
1940        enddo
1941        iold(index)=iran
1942       enddo
1943
1944       return
1945       end
1946 ccccccccccccccccccccccccccccccccccccccccccccccccc
1947 ccccccccccccccccccccccccccccccccccccccccccccccccc
1948       subroutine newconf1abb(idum,vvar,nran,i1)
1949       implicit real*8 (a-h,o-z)
1950       include 'DIMENSIONS'
1951       include 'COMMON.IOUNITS'
1952       include 'COMMON.CSA'
1953       include 'COMMON.BANK'
1954       include 'COMMON.CHAIN'
1955       include 'COMMON.GEO'
1956       include 'COMMON.TORCNSTR'
1957       include 'COMMON.CONTROL'
1958       real ran1,ran2
1959       dimension vvar(mxang,maxres,mxch),iold(ntotal)
1960       ctdif=10.
1961
1962       do k=1,numch
1963        do j=2,nres-1
1964         do i=1,4
1965          vvar(i,j,k)=bvar(i,j,k,iseed)
1966         enddo
1967        enddo
1968       enddo
1969
1970       do index=1,nran
1971        iold(index) = 0
1972       enddo
1973
1974       number=nres-2
1975
1976       iter=0
1977       do index=1,nran
1978    10  iran=  ran1(idum)*number+1
1979       if(i2ndstr.gt.0) then
1980        rtmp=ran1(idum)
1981        if(rtmp.le.rdih_bias) then
1982         iran=ran1(idum)*ndih_nconstr+1
1983         iran=idih_nconstr(iran)
1984        endif
1985       endif
1986        if(iter.gt.number) return
1987        iter=iter+1
1988        if(iter.eq.1) goto 11
1989        do ind=1,index-1
1990         if(iran.eq.iold(ind)) goto 10
1991        enddo
1992    11  continue
1993
1994        do ind=1,ngroup(iran)
1995         i=igroup(1,ind,iran)
1996         j=igroup(2,ind,iran)
1997         k=igroup(3,ind,iran)
1998         dif=rad2deg*dabs(vvar(i,j,k)-bvar(i,j,k,i1))
1999         if(dif.gt.180.) dif=360.-dif
2000         if(dif.gt.ctdif) goto 20
2001        enddo
2002        if(iter.gt.number) goto 20
2003        goto 10
2004    20  continue
2005        do ind=1,ngroup(iran)
2006         i=igroup(1,ind,iran)
2007         j=igroup(2,ind,iran)
2008         k=igroup(3,ind,iran)
2009         vvar(i,j,k)=bvar(i,j,k,i1)
2010        enddo
2011        iold(index)=iran
2012       enddo
2013
2014       return
2015       end
2016 ccccccccccccccccccccccccccccccccccccccccccccccccc
2017 ccccccccccccccccccccccccccccccccccccccccccccccccc
2018       subroutine newconf_residue(idum,vvar,i1,isize)
2019       implicit real*8 (a-h,o-z)
2020       include 'DIMENSIONS'
2021       include 'COMMON.IOUNITS'
2022       include 'COMMON.CSA'
2023       include 'COMMON.BANK'
2024       include 'COMMON.CHAIN'
2025       include 'COMMON.GEO'
2026       include 'COMMON.TORCNSTR'
2027       include 'COMMON.CONTROL'
2028       real ran1,ran2
2029       dimension vvar(mxang,maxres,mxch),iold(ntotal)
2030       ctdif=10.
2031
2032       if (iseed.gt.mxio .or. iseed.lt.1) then
2033         write (iout,*) 'Dimension ERROR in NEWCONF: ISEED',iseed
2034         call mpi_abort(mpi_comm_world,ierror,ierrcode)
2035       endif
2036       do k=1,numch
2037        do j=2,nres-1
2038         do i=1,4
2039          vvar(i,j,k)=bvar(i,j,k,iseed)
2040         enddo
2041        enddo
2042       enddo
2043
2044
2045       k=1
2046       number=nres+isize-2
2047       iter=1
2048    10 iran=  ran1(idum)*number+1
2049       if(i2ndstr.gt.0) then
2050        rtmp=ran1(idum)
2051        if(rtmp.le.rdih_bias) then
2052         iran=ran1(idum)*ndih_nconstr+1
2053         iran=idih_nconstr(iran)
2054        endif
2055       endif
2056       istart=iran-isize+1
2057       iend=iran
2058       if(istart.lt.2) istart=2
2059       if(iend.gt.nres-1) iend=nres-1
2060
2061       if(iter.eq.1) goto 11
2062       do ind=1,iter-1
2063        if(iran.eq.iold(ind)) goto 10
2064       enddo
2065    11 continue
2066
2067       do j=istart,iend
2068        do i=1,4
2069          dif=rad2deg*dabs(vvar(i,j,k)-bvar(i,j,k,i1))
2070          if(dif.gt.180.) dif=360.-dif
2071          if(dif.gt.ctdif) goto 20
2072        enddo
2073       enddo
2074       iold(iter)=iran
2075       iter=iter+1
2076       if(iter.gt.number) goto 20
2077       goto 10
2078
2079    20 continue
2080       do j=istart,iend
2081        do i=1,4
2082         vvar(i,j,k)=bvar(i,j,k,i1)
2083        enddo
2084       enddo
2085
2086       return
2087       end
2088
2089 ccccccccccccccccccccccccccccccccccccccccccccccccc
2090 ccccccccccccccccccccccccccccccccccccccccccccccccc
2091       subroutine newconf_copy(idum,vvar,i1,istart,iend)
2092       implicit real*8 (a-h,o-z)
2093       include 'DIMENSIONS'
2094       include 'COMMON.IOUNITS'
2095       include 'COMMON.CSA'
2096       include 'COMMON.BANK'
2097       include 'COMMON.CHAIN'
2098       include 'COMMON.GEO'
2099       include 'COMMON.TORCNSTR'
2100       include 'COMMON.CONTROL'
2101       real ran1,ran2
2102       dimension vvar(mxang,maxres,mxch),iold(ntotal)
2103       ctdif=10.
2104
2105       if (iseed.gt.mxio .or. iseed.lt.1) then
2106         write (iout,*) 'Dimension ERROR in NEWCONF: ISEED',iseed
2107         call mpi_abort(mpi_comm_world,ierror,ierrcode)
2108       endif
2109       do k=1,numch
2110        do j=2,nres-1
2111         do i=1,4
2112          vvar(i,j,k)=bvar(i,j,k,iseed)
2113         enddo
2114        enddo
2115       enddo
2116
2117
2118       do j=istart,iend
2119        do i=1,4
2120         vvar(i,j,1)=bvar(i,j,1,i1)
2121        enddo
2122       enddo
2123
2124       return
2125       end
2126 ccccccccccccccccccccccccccccccccccccccccccccccccc
2127 ccccccccccccccccccccccccccccccccccccccccccccccccc
2128       subroutine newconf_residue_hairpin(idum,vvar,i1,fail)
2129       implicit real*8 (a-h,o-z)
2130       include 'DIMENSIONS'
2131       include 'COMMON.IOUNITS'
2132       include 'COMMON.CSA'
2133       include 'COMMON.BANK'
2134       include 'COMMON.CHAIN'
2135       include 'COMMON.GEO'
2136       include 'COMMON.VAR'
2137       real ran1,ran2
2138       dimension vvar(mxang,maxres,mxch),iold(ntotal)
2139       integer nharp,iharp(4,maxres/3),icipa(maxres/3)
2140       logical fail,not_done
2141       ctdif=10.
2142
2143       fail=.false.
2144       if (iseed.gt.mxio .or. iseed.lt.1) then
2145         write (iout,*) 'Dimension ERROR in NEWCONF: ISEED',iseed
2146         call mpi_abort(mpi_comm_world,ierror,ierrcode)
2147       endif
2148       do k=1,numch
2149        do j=2,nres-1
2150         do i=1,4
2151          vvar(i,j,k)=bvar(i,j,k,iseed)
2152         enddo
2153        enddo
2154       enddo
2155       do k=1,numch
2156       do j=2,nres-1
2157         theta(j+1)=bvar(1,j,k,i1)
2158         phi(j+2)=bvar(2,j,k,i1)
2159         alph(j)=bvar(3,j,k,i1)
2160         omeg(j)=bvar(4,j,k,i1)
2161       enddo
2162       enddo
2163 c      call intout
2164       call chainbuild
2165       call hairpin(.false.,nharp,iharp)
2166
2167       if (nharp.eq.0) then
2168         fail=.true.
2169         return
2170       endif
2171
2172       n_used=0
2173
2174       DO III=1,NHARP
2175
2176       not_done = .true.
2177       icount=0
2178       do while (not_done)
2179         icount=icount+1
2180         iih=iran_num(1,nharp)
2181         do k=1,n_used
2182           if (iih.eq.icipa(k)) then
2183             iih=0
2184             goto 22
2185           endif
2186         enddo
2187         not_done=.false.
2188         n_used=n_used+1
2189         icipa(n_used)=iih
2190    22   continue
2191         not_done = not_done .and. icount.le.nharp
2192       enddo
2193
2194       if (iih.eq.0) then
2195         write (iout,*) "CHUJ NASTAPIL W NEWCONF_RESIDUE_HAIRPIN!!!!"
2196         fail=.true.
2197         return
2198       endif
2199
2200       istart=iharp(1,iih)+1
2201       iend=iharp(2,iih)
2202
2203 cdd      write (iout,*) "newconf_residue_hairpin: iih",iih,
2204 cdd     &   " istart",istart," iend",iend
2205
2206       do k=1,numch
2207       do j=istart,iend
2208        do i=1,4
2209          dif=rad2deg*dabs(vvar(i,j,k)-bvar(i,j,k,i1))
2210          if(dif.gt.180.) dif=360.-dif
2211          if(dif.gt.ctdif) goto 20
2212        enddo
2213       enddo
2214       enddo
2215       goto 10
2216    20 continue
2217       do k=1,numch
2218       do j=istart,iend
2219        do i=1,4
2220         vvar(i,j,k)=bvar(i,j,k,i1)
2221        enddo
2222       enddo
2223       enddo
2224 c      do j=1,numch
2225 c       do l=2,nres-1
2226 c        write (iout,'(4f8.3)') (rad2deg*vvar(i,l,j),i=1,4)
2227 c       enddo
2228 c      enddo
2229       return
2230    10 continue
2231       ENDDO
2232
2233       fail=.true.
2234
2235       return
2236       end
2237 ccccccccccccccccccccccccccccccccccccccccccccccccc
2238 ccccccccccccccccccccccccccccccccccccccccccccccccc
2239       subroutine gen_hairpin
2240       implicit real*8 (a-h,o-z)
2241       include 'DIMENSIONS'
2242       include 'COMMON.IOUNITS'
2243       include 'COMMON.CSA'
2244       include 'COMMON.BANK'
2245       include 'COMMON.CHAIN'
2246       include 'COMMON.GEO'
2247       include 'COMMON.VAR'
2248       include 'COMMON.HAIRPIN'
2249
2250 c      write (iout,*) 'Entering GEN_HAIRPIN'
2251       do iters=1,nseed
2252         i1=is(iters)
2253         do k=1,numch
2254           do j=2,nres-1
2255             theta(j+1)=bvar(1,j,k,i1)
2256             phi(j+2)=bvar(2,j,k,i1)
2257             alph(j)=bvar(3,j,k,i1)
2258             omeg(j)=bvar(4,j,k,i1)
2259           enddo
2260         enddo
2261         call chainbuild
2262         call hairpin(.false.,nharp_seed(iters),iharp_seed(1,1,iters))
2263       enddo
2264
2265       nharp_tot=0
2266       do iters=1,nseed
2267         nharp_tot=nharp_tot+nharp_seed(iters)
2268         nharp_use(iters)=4*nharp_seed(iters)
2269         do j=1,nharp_seed(iters)
2270           iharp_use(0,j,iters)=4
2271           do k=1,4
2272             iharp_use(k,j,iters)=0
2273           enddo
2274         enddo
2275       enddo
2276
2277       write (iout,*) 'GEN_HAIRPIN: nharp_tot',nharp_tot
2278 cdd      do i=1,nseed
2279 cdd        write (iout,*) 'seed',i
2280 cdd        write (iout,*) 'nharp_seed',nharp_seed(i),
2281 cdd     &     ' nharp_use',nharp_use(i)
2282 cd        write (iout,*) 'iharp_seed, iharp_use'
2283 cd        do j=1,nharp_seed(i)
2284 cd          write (iout,'(7i3)') iharp_seed(1,j,i),iharp_seed(2,j,i),
2285 cd     &      (iharp_use(k,j,i),k=0,4)
2286 cd        enddo
2287 cdd      enddo
2288       return
2289       end
2290
2291 ccccccccccccccccccccccccccccccccccccccccccccccccc
2292 ccccccccccccccccccccccccccccccccccccccccccccccccc
2293       subroutine select_frag(nn,nh,nl,ns,nb,i_csa)
2294       implicit real*8 (a-h,o-z)
2295       include 'DIMENSIONS'
2296       include 'COMMON.IOUNITS'
2297       include 'COMMON.CSA'
2298       include 'COMMON.BANK'
2299       include 'COMMON.CHAIN'
2300       include 'COMMON.GEO'
2301       include 'COMMON.VAR'
2302       include 'COMMON.HAIRPIN'
2303       include 'COMMON.DISTFIT'
2304       character*50 linia
2305       integer isec(maxres)
2306
2307
2308       nn=0
2309       nh=0
2310       nl=0
2311       ns=0
2312       nb=0
2313 cd      write (iout,*) 'Entering select_frag'
2314       do i1=1,nbank
2315         do i=1,nres
2316           isec(i)=0
2317         enddo
2318         do k=1,numch
2319           do j=2,nres-1
2320             theta(j+1)=bvar(1,j,k,i1)
2321             phi(j+2)=bvar(2,j,k,i1)
2322             alph(j)=bvar(3,j,k,i1)
2323             omeg(j)=bvar(4,j,k,i1)
2324           enddo
2325         enddo
2326         call chainbuild
2327 cd        write (iout,*) ' -- ',i1,' -- '
2328         call secondary2(.false.)
2329 c
2330 c bvar_frag nn==pair of nonlocal strands in beta sheet (loop>4)
2331 c               strands > 4 residues; used by N7 and N16
2332 c
2333         do j=1,nbfrag
2334 c
2335 Ctest 09/12/02 bfrag(2,j)-bfrag(1,j).gt.3
2336 c
2337               do i=bfrag(1,j),bfrag(2,j)
2338                 isec(i)=1
2339               enddo
2340               do i=bfrag(4,j),bfrag(3,j),sign(1,bfrag(3,j)-bfrag(4,j))
2341                 isec(i)=1
2342               enddo
2343
2344            if ( (bfrag(3,j).lt.bfrag(4,j) .or.
2345      &          bfrag(4,j)-bfrag(2,j).gt.4) .and.
2346      &          bfrag(2,j)-bfrag(1,j).gt.4 ) then
2347              nn=nn+1
2348
2349
2350              if (bfrag(3,j).lt.bfrag(4,j)) then
2351                write(linia,'(a6,i3,a1,i3,a1,i3,a1,i3)')
2352      &                     "select",bfrag(1,j)-1,"-",bfrag(2,j)-1
2353      &                         ,",",bfrag(3,j)-1,"-",bfrag(4,j)-1
2354              else
2355                write(linia,'(a6,i3,a1,i3,a1,i3,a1,i3)')
2356      &                     "select",bfrag(1,j)-1,"-",bfrag(2,j)-1
2357      &                         ,",",bfrag(4,j)-1,"-",bfrag(3,j)-1
2358
2359              endif
2360 cd             call write_pdb(i_csa*1000+nn+nh,linia,0d0)
2361
2362              bvar_frag(nn,1)=i1
2363              bvar_frag(nn,2)=4
2364              do i=1,4
2365                bvar_frag(nn,i+2)=bfrag(i,j)
2366              enddo
2367            endif
2368         enddo
2369
2370 c
2371 c hvar_frag nh==helices; used by N8 and N9
2372 c
2373         do j=1,nhfrag
2374
2375               do i=hfrag(1,j),hfrag(2,j)
2376                 isec(i)=2
2377               enddo
2378
2379           if ( hfrag(2,j)-hfrag(1,j).gt.4 ) then
2380             nh=nh+1
2381
2382 cd            write(linia,'(a6,i3,a1,i3)')
2383 cd     &                "select",hfrag(1,j)-1,"-",hfrag(2,j)-1
2384 cd            call write_pdb(i_csa*1000+nn+nh,linia,0d0)
2385
2386             hvar_frag(nh,1)=i1
2387             hvar_frag(nh,2)=hfrag(1,j)
2388             hvar_frag(nh,3)=hfrag(2,j)
2389           endif
2390         enddo
2391
2392
2393 cv        write(iout,'(i4,1pe12.4,1x,1000i1)')
2394 cv     &         i1,bene(i1),(isec(i),i=1,nres)
2395 cv            write(linia,'(i4,1x,1000i1)')
2396 cv     &         i1,(isec(i),i=1,nres)
2397 cv            call write_pdb(i_csa*1000+i1,linia,bene(i1))
2398 c
2399 c lvar_frag nl==loops; used by N14
2400 c
2401         i=1
2402         nl1=nl
2403         do while (i.lt.nres)
2404           if (isec(i).eq.0) then
2405            nl=nl+1
2406            lvar_frag(nl,1)=i1
2407            lvar_frag(nl,2)=i
2408            i=i+1
2409            do while (isec(i).eq.0.and.i.le.nres)
2410              i=i+1
2411            enddo
2412            lvar_frag(nl,3)=i-1
2413            if (lvar_frag(nl,3)-lvar_frag(nl,2).lt.1) nl=nl-1
2414           endif
2415           i=i+1
2416         enddo
2417 cd        write(iout,'(4i5)') (i,(lvar_frag(i,ii),ii=1,3),i=nl1+1,nl)
2418
2419 c
2420 c svar_frag ns==an secondary structure element; used by N15
2421 c
2422         i=1
2423         ns1=ns
2424         do while (i.lt.nres)
2425           if (isec(i).gt.0) then
2426            ns=ns+1
2427            svar_frag(ns,1)=i1
2428            svar_frag(ns,2)=i
2429            i=i+1
2430            do while (isec(i).gt.0.and.isec(i-1).eq.isec(i)
2431      &                         .and.i.le.nres)
2432              i=i+1
2433            enddo
2434            svar_frag(ns,3)=i-1
2435            if (svar_frag(ns,3)-svar_frag(ns,2).lt.1) ns=ns-1
2436           endif
2437           if (isec(i).eq.0) i=i+1
2438         enddo
2439 cd        write(iout,'(4i5)') (i,(svar_frag(i,ii),ii=1,3),i=ns1+1,ns)
2440
2441 c
2442 c avar_frag nb==any pair of beta strands; used by N17
2443 c
2444         do j=1,nbfrag
2445              nb=nb+1
2446              avar_frag(nb,1)=i1
2447              do i=1,4
2448                avar_frag(nb,i+1)=bfrag(i,j)
2449              enddo
2450         enddo
2451
2452       enddo
2453
2454       return
2455       end
2456 #endif