Fixed Makefiles, handle all force fields thru appropriate options
[unres.git] / source / unres / src_CSA / 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