2 implicit real*8 (a-h,o-z)
9 include 'COMMON.INTERACT'
10 include 'COMMON.IOUNITS'
11 include 'COMMON.DISTFIT'
12 include 'COMMON.SBRIDGE'
13 include 'COMMON.CONTROL'
14 include 'COMMON.FFIELD'
15 include 'COMMON.MINIM'
16 include 'COMMON.CHAIN'
17 double precision time0,time1
18 double precision energy(0:n_ene),ee
19 double precision var(maxvar),var1(maxvar)
21 logical debug,accepted
25 call geom_to_var(nvar,var1)
27 call etotal(energy(0))
30 rms = rmscalc(c,cref,ipermmin)
31 write(iout,*) 'etot=',0,etot,rms
32 call secondary2(.false.)
34 call write_pdb(0,'first structure',etot)
43 betbol=1.0D0/(1.9858D-3*temp)
46 c phi(jr)=pinorm(phi(jr)+d)
48 call etotal(energy(0))
51 rms = rmscalc(c,cref,ipermmin)
52 write(iout,*) 'etot=',1,etot0,rms
53 call write_pdb(1,'perturb structure',etot0)
59 phi(jr)=pinorm(phi(jr)+d)
61 call etotal(energy(0))
64 if (etot.lt.etot0) then
68 xxr=ran_number(0.0D0,1.0D0)
69 xxh=betbol*(etot-etot0)
70 if (xxh.lt.50.0D0) then
72 if (xxh.gt.xxr) accepted=.true.
76 c print *,etot0,etot,accepted
80 rms = rmscalc(c,cref,ipermmin)
81 write(iout,*) 'etot=',i,etot,rms
82 call write_pdb(i,'MC structure',etot)
84 c call geom_to_var(nvar,var1)
85 call sc_move(2,nres-1,1,10d0,nft_sc,etot)
86 call geom_to_var(nvar,var)
87 call minimize(etot,var,iretcode,nfun)
88 write(iout,*)'SUMSL return code is',iretcode,' eval ',nfun
89 call var_to_geom(nvar,var)
92 rms = rmscalc(c,cref,ipermmin)
93 write(iout,*) 'etot mcm=',i,etot,rms
94 call write_pdb(i+1,'MCM structure',etot)
95 call var_to_geom(nvar,var1)
103 c call sc_move(2,nres-1,1,10d0,nft_sc,etot)
104 c call geom_to_var(nvar,var)
107 c call write_pdb(998 ,'sc min',etot)
109 c call minimize(etot,var,iretcode,nfun)
110 c write(iout,*)'------------------------------------------------'
111 c write(iout,*)'SUMSL return code is',iretcode,' eval ',nfun
113 c call var_to_geom(nvar,var)
115 c call write_pdb(999,'full min',etot)
123 implicit real*8 (a-h,o-z)
130 include 'COMMON.INTERACT'
131 include 'COMMON.IOUNITS'
132 include 'COMMON.DISTFIT'
133 include 'COMMON.SBRIDGE'
134 include 'COMMON.CONTROL'
135 include 'COMMON.FFIELD'
136 include 'COMMON.MINIM'
137 include 'COMMON.CHAIN'
138 double precision time0,time1
139 double precision energy(0:n_ene),ee
140 double precision var(maxvar),var1(maxvar)
146 call geom_to_var(nvar,var1)
148 call etotal(energy(0))
150 write(iout,*) nnt,nct,etot
151 call write_pdb(1,'first structure',etot)
152 call secondary2(.true.)
161 call var_to_geom(nvar,var1)
162 write(iout,*) 'N16 test',(jdata(i),i=1,5)
163 call beta_slide(jdata(1),jdata(2),jdata(3),jdata(4),jdata(5)
165 call geom_to_var(nvar,var)
173 call minimize(etot,var,iretcode,nfun)
174 write(iout,*)'------------------------------------------------'
175 write(iout,*)'SUMSL return code is',iretcode,' eval ',nfun,
176 & '+ DIST eval',ieval
183 write (iout,'(a,f6.2,f8.2,a)')' Time for full min.',time1-time0,
184 & nfun/(time1-time0),' eval/s'
186 call var_to_geom(nvar,var)
188 call write_pdb(ij*100+99,'full min',etot)
198 subroutine test_local
199 implicit real*8 (a-h,o-z)
203 include 'COMMON.INTERACT'
204 include 'COMMON.IOUNITS'
205 double precision time0,time1
206 double precision energy(0:n_ene),ee
207 double precision varia(maxvar)
210 c call geom_to_var(nvar,varia)
211 call write_pdb(1,'first structure',0d0)
213 call etotal(energy(0))
215 write(iout,*) nnt,nct,etot
217 write(iout,*) 'calling sc_move'
218 call sc_move(nnt,nct,5,10d0,nft_sc,etot)
219 write(iout,*) nft_sc,etot
220 call write_pdb(2,'second structure',etot)
222 write(iout,*) 'calling local_move'
223 call local_move_init(.false.)
224 call local_move(24,29,20d0,50d0)
226 call write_pdb(3,'third structure',etot)
228 write(iout,*) 'calling sc_move'
229 call sc_move(24,29,5,10d0,nft_sc,etot)
230 write(iout,*) nft_sc,etot
231 call write_pdb(2,'last structure',etot)
238 implicit real*8 (a-h,o-z)
242 include 'COMMON.INTERACT'
243 include 'COMMON.IOUNITS'
244 double precision time0,time1
245 double precision energy(0:n_ene),ee
246 double precision varia(maxvar)
249 c call geom_to_var(nvar,varia)
250 call write_pdb(1,'first structure',0d0)
252 call etotal(energy(0))
254 write(iout,*) nnt,nct,etot
256 write(iout,*) 'calling sc_move'
258 call sc_move(nnt,nct,5,10d0,nft_sc,etot)
259 write(iout,*) nft_sc,etot
260 call write_pdb(2,'second structure',etot)
262 write(iout,*) 'calling sc_move 2nd time'
264 call sc_move(nnt,nct,5,1d0,nft_sc,etot)
265 write(iout,*) nft_sc,etot
266 call write_pdb(3,'last structure',etot)
269 c--------------------------------------------------------
270 subroutine bgrow(bstrand,nbstrand,in,ind,new)
271 implicit real*8 (a-h,o-z)
273 include 'COMMON.CHAIN'
274 integer bstrand(maxres/3,6)
276 ishift=iabs(bstrand(in,ind+4)-new)
278 print *,'bgrow',bstrand(in,ind+4),new,ishift
283 bstrand(nbstrand,5)=bstrand(nbstrand,1)
285 IF (bstrand(nbstrand,3).eq.bstrand(i,3)) THEN
286 if (bstrand(i,5).lt.bstrand(i,6)) then
287 bstrand(i,5)=bstrand(i,5)-ishift
289 bstrand(i,5)=bstrand(i,5)+ishift
294 bstrand(nbstrand,6)=bstrand(nbstrand,2)
296 IF (bstrand(nbstrand,3).eq.bstrand(i,3)) THEN
297 if (bstrand(i,6).lt.bstrand(i,5)) then
298 bstrand(i,6)=bstrand(i,6)-ishift
300 bstrand(i,6)=bstrand(i,6)+ishift
311 c------------------------------------------
313 implicit real*8 (a-h,o-z)
319 include 'COMMON.CHAIN'
320 include 'COMMON.IOUNITS'
322 include 'COMMON.CONTROL'
323 include 'COMMON.SBRIDGE'
324 include 'COMMON.FFIELD'
325 include 'COMMON.MINIM'
327 include 'COMMON.DISTFIT'
328 integer if(20,maxres),nif,ifa(20)
329 integer ibc(0:maxres,0:maxres),istrand(20)
330 integer ibd(maxres),ifb(10,2),nifb,lifb(10),lifb0
331 integer itmp(20,maxres)
332 double precision time0,time1
333 double precision energy(0:n_ene),ee
334 double precision varia(maxvar),vorg(maxvar)
336 logical debug,ltest,usedbfrag(maxres/3)
339 integer betasheet(maxres),ibetasheet(maxres),nbetasheet
340 integer bstrand(maxres/3,6),nbstrand
342 c------------------------
345 c------------------------
352 call geom_to_var(nvar,vorg)
353 call secondary2(debug)
355 if (nbfrag.le.1) return
362 nbetasheet=nbetasheet+1
364 bstrand(1,1)=bfrag(1,1)
365 bstrand(1,2)=bfrag(2,1)
366 bstrand(1,3)=nbetasheet
368 bstrand(1,5)=bfrag(1,1)
369 bstrand(1,6)=bfrag(2,1)
370 do i=bfrag(1,1),bfrag(2,1)
371 betasheet(i)=nbetasheet
375 bstrand(2,1)=bfrag(3,1)
376 bstrand(2,2)=bfrag(4,1)
377 bstrand(2,3)=nbetasheet
378 bstrand(2,5)=bfrag(3,1)
379 bstrand(2,6)=bfrag(4,1)
381 if (bfrag(3,1).le.bfrag(4,1)) then
383 do i=bfrag(3,1),bfrag(4,1)
384 betasheet(i)=nbetasheet
389 do i=bfrag(4,1),bfrag(3,1)
390 betasheet(i)=nbetasheet
397 do while (iused_nbfrag.ne.nbfrag)
401 IF (.not.usedbfrag(j)) THEN
403 write (*,*) j,(bfrag(i,j),i=1,4)
405 write (*,'(i4,a3,10i4)') jk,'B',(bstrand(i,jk),i=1,nbstrand)
407 write (*,*) '------------------'
410 if (bfrag(3,j).le.bfrag(4,j)) then
411 do i=bfrag(3,j),bfrag(4,j)
412 if(betasheet(i).eq.nbetasheet) then
414 do k=bfrag(3,j),bfrag(4,j)
415 betasheet(k)=nbetasheet
420 iused_nbfrag=iused_nbfrag+1
421 do k=bfrag(1,j),bfrag(2,j)
422 betasheet(k)=nbetasheet
423 ibetasheet(k)=nbstrand
425 if (bstrand(in,4).lt.0) then
426 bstrand(nbstrand,1)=bfrag(2,j)
427 bstrand(nbstrand,2)=bfrag(1,j)
428 bstrand(nbstrand,3)=nbetasheet
429 bstrand(nbstrand,4)=-nbstrand
430 bstrand(nbstrand,5)=bstrand(nbstrand,1)
431 bstrand(nbstrand,6)=bstrand(nbstrand,2)
432 if(bstrand(in,1).lt.bfrag(4,j)) then
433 call bgrow(bstrand,nbstrand,in,1,bfrag(4,j))
435 bstrand(nbstrand,5)=bstrand(nbstrand,5)+
436 & (bstrand(in,5)-bfrag(4,j))
438 if(bstrand(in,2).gt.bfrag(3,j)) then
439 call bgrow(bstrand,nbstrand,in,2,bfrag(3,j))
441 bstrand(nbstrand,6)=bstrand(nbstrand,6)-
442 & (-bstrand(in,6)+bfrag(3,j))
445 bstrand(nbstrand,1)=bfrag(1,j)
446 bstrand(nbstrand,2)=bfrag(2,j)
447 bstrand(nbstrand,3)=nbetasheet
448 bstrand(nbstrand,4)=nbstrand
449 bstrand(nbstrand,5)=bstrand(nbstrand,1)
450 bstrand(nbstrand,6)=bstrand(nbstrand,2)
451 if(bstrand(in,1).gt.bfrag(3,j)) then
452 call bgrow(bstrand,nbstrand,in,1,bfrag(3,j))
454 bstrand(nbstrand,5)=bstrand(nbstrand,5)-
455 & (-bstrand(in,5)+bfrag(3,j))
457 if(bstrand(in,2).lt.bfrag(4,j)) then
458 call bgrow(bstrand,nbstrand,in,2,bfrag(4,j))
460 bstrand(nbstrand,6)=bstrand(nbstrand,6)+
461 & (bstrand(in,6)-bfrag(4,j))
466 if(betasheet(bfrag(1,j)+i-bfrag(3,j)).eq.nbetasheet) then
467 in=ibetasheet(bfrag(1,j)+i-bfrag(3,j))
468 do k=bfrag(1,j),bfrag(2,j)
469 betasheet(k)=nbetasheet
474 iused_nbfrag=iused_nbfrag+1
475 do k=bfrag(3,1),bfrag(4,1)
476 betasheet(k)=nbetasheet
477 ibetasheet(k)=nbstrand
479 if (bstrand(in,4).lt.0) then
480 bstrand(nbstrand,1)=bfrag(4,j)
481 bstrand(nbstrand,2)=bfrag(3,j)
482 bstrand(nbstrand,3)=nbetasheet
483 bstrand(nbstrand,4)=-nbstrand
484 bstrand(nbstrand,5)=bstrand(nbstrand,1)
485 bstrand(nbstrand,6)=bstrand(nbstrand,2)
486 if(bstrand(in,1).lt.bfrag(2,j)) then
487 call bgrow(bstrand,nbstrand,in,1,bfrag(2,j))
489 bstrand(nbstrand,5)=bstrand(nbstrand,5)+
490 & (bstrand(in,5)-bfrag(2,j))
492 if(bstrand(in,2).gt.bfrag(1,j)) then
493 call bgrow(bstrand,nbstrand,in,2,bfrag(1,j))
495 bstrand(nbstrand,6)=bstrand(nbstrand,6)-
496 & (-bstrand(in,6)+bfrag(1,j))
499 bstrand(nbstrand,1)=bfrag(3,j)
500 bstrand(nbstrand,2)=bfrag(4,j)
501 bstrand(nbstrand,3)=nbetasheet
502 bstrand(nbstrand,4)=nbstrand
503 bstrand(nbstrand,5)=bstrand(nbstrand,1)
504 bstrand(nbstrand,6)=bstrand(nbstrand,2)
505 if(bstrand(in,1).gt.bfrag(1,j)) then
506 call bgrow(bstrand,nbstrand,in,1,bfrag(1,j))
508 bstrand(nbstrand,5)=bstrand(nbstrand,5)-
509 & (-bstrand(in,5)+bfrag(1,j))
511 if(bstrand(in,2).lt.bfrag(2,j)) then
512 call bgrow(bstrand,nbstrand,in,2,bfrag(2,j))
514 bstrand(nbstrand,6)=bstrand(nbstrand,6)+
515 & (bstrand(in,6)-bfrag(2,j))
522 do i=bfrag(4,j),bfrag(3,j)
523 if(betasheet(i).eq.nbetasheet) then
525 do k=bfrag(4,j),bfrag(3,j)
526 betasheet(k)=nbetasheet
531 iused_nbfrag=iused_nbfrag+1
532 do k=bfrag(1,j),bfrag(2,j)
533 betasheet(k)=nbetasheet
534 ibetasheet(k)=nbstrand
536 if (bstrand(in,4).lt.0) then
537 bstrand(nbstrand,1)=bfrag(1,j)
538 bstrand(nbstrand,2)=bfrag(2,j)
539 bstrand(nbstrand,3)=nbetasheet
540 bstrand(nbstrand,4)=nbstrand
541 bstrand(nbstrand,5)=bstrand(nbstrand,1)
542 bstrand(nbstrand,6)=bstrand(nbstrand,2)
543 if(bstrand(in,1).lt.bfrag(3,j)) then
544 call bgrow(bstrand,nbstrand,in,1,bfrag(3,j))
546 bstrand(nbstrand,5)=bstrand(nbstrand,5)-
547 & (bstrand(in,5)-bfrag(3,j))
549 if(bstrand(in,2).gt.bfrag(4,j)) then
550 call bgrow(bstrand,nbstrand,in,2,bfrag(4,j))
552 bstrand(nbstrand,6)=bstrand(nbstrand,6)+
553 & (-bstrand(in,6)+bfrag(4,j))
556 bstrand(nbstrand,1)=bfrag(2,j)
557 bstrand(nbstrand,2)=bfrag(1,j)
558 bstrand(nbstrand,3)=nbetasheet
559 bstrand(nbstrand,4)=-nbstrand
560 bstrand(nbstrand,5)=bstrand(nbstrand,1)
561 bstrand(nbstrand,6)=bstrand(nbstrand,2)
562 if(bstrand(in,1).gt.bfrag(4,j)) then
563 call bgrow(bstrand,nbstrand,in,1,bfrag(4,j))
565 bstrand(nbstrand,5)=bstrand(nbstrand,5)+
566 & (-bstrand(in,5)+bfrag(4,j))
568 if(bstrand(in,2).lt.bfrag(3,j)) then
569 call bgrow(bstrand,nbstrand,in,2,bfrag(3,j))
571 bstrand(nbstrand,6)=bstrand(nbstrand,6)-
572 & (bstrand(in,6)-bfrag(3,j))
577 if(betasheet(bfrag(2,j)-i+bfrag(4,j)).eq.nbetasheet) then
578 in=ibetasheet(bfrag(2,j)-i+bfrag(4,j))
579 do k=bfrag(1,j),bfrag(2,j)
580 betasheet(k)=nbetasheet
585 iused_nbfrag=iused_nbfrag+1
586 do k=bfrag(4,j),bfrag(3,j)
587 betasheet(k)=nbetasheet
588 ibetasheet(k)=nbstrand
590 if (bstrand(in,4).lt.0) then
591 bstrand(nbstrand,1)=bfrag(4,j)
592 bstrand(nbstrand,2)=bfrag(3,j)
593 bstrand(nbstrand,3)=nbetasheet
594 bstrand(nbstrand,4)=nbstrand
595 bstrand(nbstrand,5)=bstrand(nbstrand,1)
596 bstrand(nbstrand,6)=bstrand(nbstrand,2)
597 if(bstrand(in,1).lt.bfrag(2,j)) then
598 call bgrow(bstrand,nbstrand,in,1,bfrag(2,j))
600 bstrand(nbstrand,5)=bstrand(nbstrand,5)-
601 & (bstrand(in,5)-bfrag(2,j))
603 if(bstrand(in,2).gt.bfrag(1,j)) then
604 call bgrow(bstrand,nbstrand,in,2,bfrag(1,j))
606 bstrand(nbstrand,6)=bstrand(nbstrand,6)+
607 & (-bstrand(in,6)+bfrag(1,j))
610 bstrand(nbstrand,1)=bfrag(3,j)
611 bstrand(nbstrand,2)=bfrag(4,j)
612 bstrand(nbstrand,3)=nbetasheet
613 bstrand(nbstrand,4)=-nbstrand
614 bstrand(nbstrand,5)=bstrand(nbstrand,1)
615 bstrand(nbstrand,6)=bstrand(nbstrand,2)
616 if(bstrand(in,1).gt.bfrag(1,j)) then
617 call bgrow(bstrand,nbstrand,in,1,bfrag(1,j))
619 bstrand(nbstrand,5)=bstrand(nbstrand,5)+
620 & (-bstrand(in,5)+bfrag(1,j))
622 if(bstrand(in,2).lt.bfrag(2,j)) then
623 call bgrow(bstrand,nbstrand,in,2,bfrag(2,j))
625 bstrand(nbstrand,6)=bstrand(nbstrand,6)-
626 & (bstrand(in,6)-bfrag(2,j))
640 do while (usedbfrag(j))
645 nbetasheet=nbetasheet+1
646 bstrand(nbstrand,1)=bfrag(1,j)
647 bstrand(nbstrand,2)=bfrag(2,j)
648 bstrand(nbstrand,3)=nbetasheet
649 bstrand(nbstrand,5)=bfrag(1,j)
650 bstrand(nbstrand,6)=bfrag(2,j)
652 bstrand(nbstrand,4)=nbstrand
653 do i=bfrag(1,j),bfrag(2,j)
654 betasheet(i)=nbetasheet
655 ibetasheet(i)=nbstrand
659 bstrand(nbstrand,1)=bfrag(3,j)
660 bstrand(nbstrand,2)=bfrag(4,j)
661 bstrand(nbstrand,3)=nbetasheet
662 bstrand(nbstrand,5)=bfrag(3,j)
663 bstrand(nbstrand,6)=bfrag(4,j)
665 if (bfrag(3,j).le.bfrag(4,j)) then
666 bstrand(nbstrand,4)=nbstrand
667 do i=bfrag(3,j),bfrag(4,j)
668 betasheet(i)=nbetasheet
669 ibetasheet(i)=nbstrand
672 bstrand(nbstrand,4)=-nbstrand
673 do i=bfrag(4,j),bfrag(3,j)
674 betasheet(i)=nbetasheet
675 ibetasheet(i)=nbstrand
679 iused_nbfrag=iused_nbfrag+1
685 write (*,'(i4,a3,10i4)') jk,'A',(bstrand(i,jk),i=1,nbstrand)
692 if (betasheet(i).ne.0) write(*,*) i,betasheet(i),ibetasheet(i)
696 write (*,'(i4,a3,10i4)') j,':',(bstrand(i,j),i=1,nbstrand)
699 c------------------------
703 if(iabs(bstrand(i,5)-bstrand(j,5)).le.5 .or.
704 & iabs(bstrand(i,6)-bstrand(j,6)).le.5 ) then
706 ifb(nifb,1)=bstrand(i,4)
707 ifb(nifb,2)=bstrand(j,4)
714 write (*,'(a3,20i4)') "ifb",i,ifb(i,1),ifb(i,2)
720 write (*,'(a3,20i4)') "ifa",(ifa(i),i=1,nbstrand)
722 nif=iabs(bstrand(1,6)-bstrand(1,5))+1
724 if (iabs(bstrand(j,6)-bstrand(j,5))+1.gt.nif)
725 & nif=iabs(bstrand(j,6)-bstrand(j,5))+1
731 if(j,i)=bstrand(j,6)+(i-1)*sign(1,bstrand(j,5)-bstrand(j,6))
732 if (if(j,i).gt.0) then
733 if(betasheet(if(j,i)).eq.0 .or.
734 & ibetasheet(if(j,i)).ne.iabs(bstrand(j,4))) if(j,i)=0
739 write(*,'(a3,10i4)') 'if ',(if(j,i),j=1,nbstrand)
742 c read (inp,*) (ifa(i),i=1,4)
744 c read (inp,*,err=20,end=20) (if(j,i),j=1,4)
748 c------------------------
753 cccccccccccccccccccccccccccccccccc
755 cccccccccccccccccccccccccccccccccc
759 istrand(is-j+1)=int(ii/is**(is-j))
760 ii=ii-istrand(is-j+1)*is**(is-j)
764 istrand(k)=istrand(k)+1
765 if(istrand(k).gt.isa) istrand(k)=istrand(k)-2*isa-1
769 if(istrand(k).eq.istrand(l).and.k.ne.l.or.
770 & istrand(k).eq.-istrand(l).and.k.ne.l) ltest=.false.
779 & ifb(m,1).eq.istrand(k).and.ifb(m,2).eq.istrand(k+1).or.
780 & ifb(m,2).eq.istrand(k).and.ifb(m,1).eq.istrand(k+1).or.
781 & -ifb(m,1).eq.istrand(k).and.-ifb(m,2).eq.istrand(k+1).or.
782 & -ifb(m,2).eq.istrand(k).and.-ifb(m,1).eq.istrand(k+1))
788 if (mod(isa,2).eq.0) then
790 if (istrand(k).eq.1) ltest=.false.
794 if (istrand(k).eq.1) ltest=.false.
798 IF (ltest.and.lifb0.eq.1) THEN
801 call var_to_geom(nvar,vorg)
803 write (*,'(i5,i10,10i3)') iconf,ig,(istrand(k),k=1,isa)
804 write (iout,'(i5,i10,10i3)') iconf,ig,(istrand(k),k=1,isa)
805 write (linia,'(10i3)') (istrand(k),k=1,isa)
815 if ( sign(1,istrand(i)).eq.sign(1,ifa(iabs(istrand(i)))) ) then
817 itmp(iabs(istrand(i)),j)=if(iabs(ifa(iabs(istrand(i)))),j)
821 itmp(iabs(istrand(i)),j)=if(iabs(ifa(iabs(istrand(i)))),nif-j+1)
827 write(*,*) (itmp(j,i),j=1,4)
831 c ifa(1),ifa(2),ifa(3),ifa(4)
832 c if(1,i),if(2,i),if(3,i),if(4,i)
837 & ifb(m,1).eq.istrand(k).and.ifb(m,2).eq.istrand(k+1).or.
838 & ifb(m,2).eq.istrand(k).and.ifb(m,1).eq.istrand(k+1).or.
839 & -ifb(m,1).eq.istrand(k).and.-ifb(m,2).eq.istrand(k+1).or.
840 & -ifb(m,2).eq.istrand(k).and.-ifb(m,1).eq.istrand(k+1))
848 ibc(itmp(iabs(istrand(k)),i),itmp(iabs(istrand(k+1)),i))=-1
850 ibc(itmp(iabs(istrand(k)),i),itmp(iabs(istrand(k+1)),i))=-2
854 & ibc(itmp(iabs(istrand(k)),i),itmp(iabs(istrand(k+2)),i))=-3
856 & ibc(itmp(iabs(istrand(k)),i),itmp(iabs(istrand(k+3)),i))=-4
859 c------------------------
862 c freeze sec.elements
872 do i=bfrag(1,j),bfrag(2,j)
877 if (bfrag(3,j).le.bfrag(4,j)) then
878 do i=bfrag(3,j),bfrag(4,j)
884 do i=bfrag(4,j),bfrag(3,j)
892 do i=hfrag(1,j),hfrag(2,j)
900 c------------------------
901 c generate constrains
909 if ( ibc(i,j).eq.-1 .or. ibc(j,i).eq.-1) then
917 else if ( ibc(i,j).eq.-2 .or. ibc(j,i).eq.-2) then
925 else if ( ibc(i,j).eq.-3 .or. ibc(j,i).eq.-3) then
933 else if ( ibc(i,j).eq.-4 .or. ibc(j,i).eq.-4) then
941 else if ( ibc(i,j).gt.0 ) then
942 d0(ind)=DIST(i,ibc(i,j))
949 else if ( ibc(j,i).gt.0 ) then
950 d0(ind)=DIST(ibc(j,i),j)
964 cd--------------------------
966 write(iout,'(i3,2i4,a3,2i4,f7.2)') (i,ibc(ihpb(i),jhpb(i)),
967 & ibc(jhpb(i),ihpb(i)),' --',
968 & ihpb(i),jhpb(i),dhpb(i),i=1,nhpb)
974 call contact_cp_min(varia,ifun,iconf,linia,debug)
981 call minimize(etot,varia,iretcode,nfun)
982 write(iout,*)'------------------------------------------------'
983 write(iout,*)'SUMSL return code is',iretcode,' eval ',nfun,
991 write (iout,'(a,f6.2,f8.2,a)')' Time for full min.',time1-time0,
992 & nfun/(time1-time0),' eval/s'
994 write (linia,'(a10,10i3)') 'full_min',(istrand(k),k=1,isa)
995 call var_to_geom(nvar,varia)
997 call write_pdb(900+iconf,linia,etot)
1000 call etotal(energy(0))
1002 call enerprint(energy(0))
1004 cd call briefout(0,etot)
1005 cd call secondary2(.true.)
1009 cccccccccccccccccccccccccccccccccccc
1012 cccccccccccccccccccccccccccccccccccc
1015 10 write (iout,'(a)') 'Error reading test structure.'
1018 c--------------------------------------------------------
1021 implicit real*8 (a-h,o-z)
1022 include 'DIMENSIONS'
1026 include 'COMMON.GEO'
1027 include 'COMMON.CHAIN'
1028 include 'COMMON.IOUNITS'
1029 include 'COMMON.VAR'
1030 include 'COMMON.CONTROL'
1031 include 'COMMON.SBRIDGE'
1032 include 'COMMON.FFIELD'
1033 include 'COMMON.MINIM'
1035 include 'COMMON.DISTFIT'
1036 integer if(3,maxres),nif
1037 integer ibc(maxres,maxres),istrand(20)
1038 integer ibd(maxres),ifb(10,2),nifb,lifb(10),lifb0
1039 double precision time0,time1
1040 double precision energy(0:n_ene),ee
1041 double precision varia(maxvar)
1047 read (inp,*,err=20,end=20) if(1,i),if(2,i),if(3,i)
1050 write (*,'(a4,3i5)') ('if =',if(1,i),if(2,i),if(3,i),
1054 c------------------------
1055 call secondary2(debug)
1056 c------------------------
1064 c freeze sec.elements and store indexes for beta constrains
1074 do i=bfrag(1,j),bfrag(2,j)
1079 if (bfrag(3,j).le.bfrag(4,j)) then
1080 do i=bfrag(3,j),bfrag(4,j)
1084 ibc(bfrag(1,j)+i-bfrag(3,j),i)=-1
1087 do i=bfrag(4,j),bfrag(3,j)
1091 ibc(bfrag(2,j)-i+bfrag(4,j),i)=-1
1096 do i=hfrag(1,j),hfrag(2,j)
1105 c ---------------- test --------------
1107 if (ibc(if(1,i),if(2,i)).eq.-1) then
1108 ibc(if(1,i),if(2,i))=if(3,i)
1109 ibc(if(1,i),if(3,i))=if(2,i)
1110 else if (ibc(if(2,i),if(1,i)).eq.-1) then
1111 ibc(if(2,i),if(1,i))=0
1112 ibc(if(1,i),if(2,i))=if(3,i)
1113 ibc(if(1,i),if(3,i))=if(2,i)
1115 ibc(if(1,i),if(2,i))=if(3,i)
1116 ibc(if(1,i),if(3,i))=if(2,i)
1122 if (ibc(i,j).ne.0) write(*,'(3i5)') i,j,ibc(i,j)
1125 c------------------------
1131 if ( ibc(i,j).eq.-1 ) then
1139 else if ( ibc(i,j).gt.0 ) then
1140 d0(ind)=DIST(i,ibc(i,j))
1147 else if ( ibc(j,i).gt.0 ) then
1148 d0(ind)=DIST(ibc(j,i),j)
1162 cd--------------------------
1163 write(*,'(i3,2i4,a3,2i4,f7.2)') (i,ibc(ihpb(i),jhpb(i)),
1164 & ibc(jhpb(i),ihpb(i)),' --',
1165 & ihpb(i),jhpb(i),dhpb(i),i=1,nhpb)
1172 call contact_cp_min(varia,ieval,in_pdb,linia,debug)
1179 call minimize(etot,varia,iretcode,nfun)
1180 write(iout,*)'------------------------------------------------'
1181 write(iout,*)'SUMSL return code is',iretcode,' eval ',nfun,
1182 & '+ DIST eval',ieval
1189 write (iout,'(a,f6.2,f8.2,a)')' Time for full min.',time1-time0,
1190 & nfun/(time1-time0),' eval/s'
1193 call var_to_geom(nvar,varia)
1195 call write_pdb(999,'full min',etot)
1198 call etotal(energy(0))
1200 call enerprint(energy(0))
1202 call briefout(0,etot)
1203 call secondary2(.true.)
1206 10 write (iout,'(a)') 'Error reading test structure.'
1214 implicit real*8 (a-h,o-z)
1215 include 'DIMENSIONS'
1219 include 'COMMON.GEO'
1220 include 'COMMON.CHAIN'
1221 include 'COMMON.IOUNITS'
1222 include 'COMMON.VAR'
1223 include 'COMMON.CONTROL'
1224 include 'COMMON.SBRIDGE'
1225 include 'COMMON.FFIELD'
1226 include 'COMMON.MINIM'
1228 include 'COMMON.DISTFIT'
1231 double precision time0,time1
1232 double precision energy(0:n_ene),ee
1233 double precision theta2(maxres),phi2(maxres),alph2(maxres),
1235 & theta1(maxres),phi1(maxres),alph1(maxres),
1237 double precision varia(maxvar),varia2(maxvar)
1241 read (inp,*,err=10,end=10) if(1,1),if(1,2),if(2,1),if(2,2)
1242 write (iout,'(a4,4i5)') 'if =',if(1,1),if(1,2),if(2,1),if(2,2)
1243 read (inp,*,err=10,end=10) (theta2(i),i=3,nres)
1244 read (inp,*,err=10,end=10) (phi2(i),i=4,nres)
1245 read (inp,*,err=10,end=10) (alph2(i),i=2,nres-1)
1246 read (inp,*,err=10,end=10) (omeg2(i),i=2,nres-1)
1248 theta2(i)=deg2rad*theta2(i)
1249 phi2(i)=deg2rad*phi2(i)
1250 alph2(i)=deg2rad*alph2(i)
1251 omeg2(i)=deg2rad*omeg2(i)
1265 c------------------------
1270 do i=if(j,1),if(j,2)
1276 call geom_to_var(nvar,varia)
1277 call write_pdb(1,'first structure',0d0)
1279 call secondary(.true.)
1281 call secondary2(.true.)
1284 if ( (bfrag(3,j).lt.bfrag(4,j) .or.
1285 & bfrag(4,j)-bfrag(2,j).gt.4) .and.
1286 & bfrag(2,j)-bfrag(1,j).gt.3 ) then
1289 if (bfrag(3,j).lt.bfrag(4,j)) then
1290 write(iout,'(a6,i3,a1,i3,a1,i3,a1,i3)')
1291 & "select",bfrag(1,j)-1,"-",bfrag(2,j)-1
1292 & ,",",bfrag(3,j)-1,"-",bfrag(4,j)-1
1294 write(iout,'(a6,i3,a1,i3,a1,i3,a1,i3)')
1295 & "select",bfrag(1,j)-1,"-",bfrag(2,j)-1
1296 & ,",",bfrag(4,j)-1,"-",bfrag(3,j)-1
1309 call geom_to_var(nvar,varia2)
1310 call write_pdb(2,'second structure',0d0)
1314 c-------------------------------------------------------
1317 call contact_cp(varia,varia2,iff,ifun,7)
1324 call minimize(etot,varia,iretcode,nfun)
1325 write(iout,*)'------------------------------------------------'
1326 write(iout,*)'SUMSL return code is',iretcode,' eval ',nfun,
1327 & '+ DIST eval',ifun
1334 write (iout,'(a,f6.2,f8.2,a)')' Time for full min.',time1-time0,
1335 & nfun/(time1-time0),' eval/s'
1338 call var_to_geom(nvar,varia)
1340 call write_pdb(999,'full min',etot)
1343 call etotal(energy(0))
1345 call enerprint(energy(0))
1347 call briefout(0,etot)
1350 10 write (iout,'(a)') 'Error reading test structure.'
1354 c-------------------------------------------------
1355 c-------------------------------------------------
1357 subroutine secondary(lprint)
1358 implicit real*8 (a-h,o-z)
1359 include 'DIMENSIONS'
1360 include 'COMMON.CHAIN'
1361 include 'COMMON.IOUNITS'
1362 include 'COMMON.DISTFIT'
1364 integer ncont,icont(2,maxres*maxres/2),isec(maxres,3)
1365 logical lprint,not_done
1366 real dcont(maxres*maxres/2),d
1371 double precision xpi(3),xpj(3)
1376 cd call write_pdb(99,'sec structure',0d0)
1388 xpi(k)=0.5d0*(c(k,i-1)+c(k,i))
1392 xpj(k)=0.5d0*(c(k,j-1)+c(k,j))
1394 cd d = (c(1,i)-c(1,j))*(c(1,i)-c(1,j)) +
1395 cd & (c(2,i)-c(2,j))*(c(2,i)-c(2,j)) +
1396 cd & (c(3,i)-c(3,j))*(c(3,i)-c(3,j))
1397 cd print *,'CA',i,j,d
1398 d = (xpi(1)-xpj(1))*(xpi(1)-xpj(1)) +
1399 & (xpi(2)-xpj(2))*(xpi(2)-xpj(2)) +
1400 & (xpi(3)-xpj(3))*(xpi(3)-xpj(3))
1401 if ( d.lt.rcomp*rcomp) then
1405 dcont(ncont)=sqrt(d)
1411 write (iout,'(a)') '#PP contact map distances:'
1413 write (iout,'(3i4,f10.5)')
1414 & i,icont(1,i),icont(2,i),dcont(i)
1418 c finding parallel beta
1419 cd write (iout,*) '------- looking for parallel beta -----------'
1425 if(dcont(i).le.rbeta .and. j1-i1.gt.4 .and.
1426 & isec(i1,1).le.1.and.isec(j1,1).le.1.and.
1427 & (isec(i1,2).ne.isec(j1,2).or.isec(i1,2)*isec(j1,2).eq.0).and.
1428 & (isec(i1,3).ne.isec(j1,3).or.isec(i1,3)*isec(j1,3).eq.0).and.
1429 & (isec(i1,2).ne.isec(j1,3).or.isec(i1,2)*isec(j1,3).eq.0).and.
1430 & (isec(i1,3).ne.isec(j1,2).or.isec(i1,3)*isec(j1,2).eq.0)
1434 cd write (iout,*) i1,j1,dcont(i)
1440 if (i1.eq.icont(1,j) .and. j1.eq.icont(2,j)
1441 & .and. dcont(j).le.rbeta .and.
1442 & isec(i1,1).le.1.and.isec(j1,1).le.1.and.
1443 & (isec(i1,2).ne.isec(j1,2).or.isec(i1,2)*isec(j1,2).eq.0).and.
1444 & (isec(i1,3).ne.isec(j1,3).or.isec(i1,3)*isec(j1,3).eq.0).and.
1445 & (isec(i1,2).ne.isec(j1,3).or.isec(i1,2)*isec(j1,3).eq.0).and.
1446 & (isec(i1,3).ne.isec(j1,2).or.isec(i1,3)*isec(j1,2).eq.0)
1451 cd write (iout,*) i1,j1,dcont(j),not_done
1455 if (i1-ii1.gt.1) then
1459 if(lprint)write(iout,*)'parallel beta',nbeta,ii1,i1,jj1,j1
1468 isec(ij,1)=isec(ij,1)+1
1469 isec(ij,1+isec(ij,1))=nbeta
1472 isec(ij,1)=isec(ij,1)+1
1473 isec(ij,1+isec(ij,1))=nbeta
1478 if (nbeta.le.9) then
1479 write(12,'(a18,i1,a9,i3,a2,i3,a1)')
1480 & "DefPropRes 'strand",nstrand,
1481 & "' 'num = ",ii1-1,"..",i1-1,"'"
1483 write(12,'(a18,i2,a9,i3,a2,i3,a1)')
1484 & "DefPropRes 'strand",nstrand,
1485 & "' 'num = ",ii1-1,"..",i1-1,"'"
1488 if (nbeta.le.9) then
1489 write(12,'(a18,i1,a9,i3,a2,i3,a1)')
1490 & "DefPropRes 'strand",nstrand,
1491 & "' 'num = ",jj1-1,"..",j1-1,"'"
1493 write(12,'(a18,i2,a9,i3,a2,i3,a1)')
1494 & "DefPropRes 'strand",nstrand,
1495 & "' 'num = ",jj1-1,"..",j1-1,"'"
1497 write(12,'(a8,4i4)')
1498 & "SetNeigh",ii1-1,i1-1,jj1-1,j1-1
1504 c finding antiparallel beta
1505 cd write (iout,*) '--------- looking for antiparallel beta ---------'
1510 if (dcont(i).le.rbeta.and.
1511 & isec(i1,1).le.1.and.isec(j1,1).le.1.and.
1512 & (isec(i1,2).ne.isec(j1,2).or.isec(i1,2)*isec(j1,2).eq.0).and.
1513 & (isec(i1,3).ne.isec(j1,3).or.isec(i1,3)*isec(j1,3).eq.0).and.
1514 & (isec(i1,2).ne.isec(j1,3).or.isec(i1,2)*isec(j1,3).eq.0).and.
1515 & (isec(i1,3).ne.isec(j1,2).or.isec(i1,3)*isec(j1,2).eq.0)
1519 cd write (iout,*) i1,j1,dcont(i)
1526 if (i1.eq.icont(1,j).and.j1.eq.icont(2,j) .and.
1527 & isec(i1,1).le.1.and.isec(j1,1).le.1.and.
1528 & (isec(i1,2).ne.isec(j1,2).or.isec(i1,2)*isec(j1,2).eq.0).and.
1529 & (isec(i1,3).ne.isec(j1,3).or.isec(i1,3)*isec(j1,3).eq.0).and.
1530 & (isec(i1,2).ne.isec(j1,3).or.isec(i1,2)*isec(j1,3).eq.0).and.
1531 & (isec(i1,3).ne.isec(j1,2).or.isec(i1,3)*isec(j1,2).eq.0)
1532 & .and. dcont(j).le.rbeta ) goto 6
1536 cd write (iout,*) i1,j1,dcont(j),not_done
1540 if (i1-ii1.gt.1) then
1541 if(lprint)write (iout,*)'antiparallel beta',
1542 & nbeta,ii1-1,i1,jj1,j1-1
1545 bfrag(1,nbfrag)=max0(ii1-1,1)
1548 bfrag(4,nbfrag)=max0(j1-1,1)
1553 isec(ij,1)=isec(ij,1)+1
1554 isec(ij,1+isec(ij,1))=nbeta
1558 isec(ij,1)=isec(ij,1)+1
1559 isec(ij,1+isec(ij,1))=nbeta
1565 if (nstrand.le.9) then
1566 write(12,'(a18,i1,a9,i3,a2,i3,a1)')
1567 & "DefPropRes 'strand",nstrand,
1568 & "' 'num = ",ii1-2,"..",i1-1,"'"
1570 write(12,'(a18,i2,a9,i3,a2,i3,a1)')
1571 & "DefPropRes 'strand",nstrand,
1572 & "' 'num = ",ii1-2,"..",i1-1,"'"
1575 if (nstrand.le.9) then
1576 write(12,'(a18,i1,a9,i3,a2,i3,a1)')
1577 & "DefPropRes 'strand",nstrand,
1578 & "' 'num = ",j1-2,"..",jj1-1,"'"
1580 write(12,'(a18,i2,a9,i3,a2,i3,a1)')
1581 & "DefPropRes 'strand",nstrand,
1582 & "' 'num = ",j1-2,"..",jj1-1,"'"
1584 write(12,'(a8,4i4)')
1585 & "SetNeigh",ii1-2,i1-1,jj1-1,j1-2
1591 if (nstrand.gt.0.and.lprint) then
1592 write(12,'(a27,$)') "DefPropRes 'sheet' 'strand1"
1595 write(12,'(a9,i1,$)') " | strand",i
1597 write(12,'(a9,i2,$)') " | strand",i
1600 write(12,'(a1)') "'"
1604 c finding alpha or 310 helix
1610 if (j1.eq.i1+3.and.dcont(i).le.r310
1611 & .or.j1.eq.i1+4.and.dcont(i).le.ralfa ) then
1612 cd if (j1.eq.i1+3) write (iout,*) "found 1-4 ",i1,j1,dcont(i)
1613 cd if (j1.eq.i1+4) write (iout,*) "found 1-5 ",i1,j1,dcont(i)
1616 if (isec(ii1,1).eq.0) then
1625 if (i1.eq.icont(1,j) .and. j1.eq.icont(2,j)) goto 10
1629 cd write (iout,*) i1,j1,not_done
1632 if (j1-ii1.gt.4) then
1634 cd write (iout,*)'helix',nhelix,ii1,j1
1638 hfrag(2,nhfrag)=max0(j1-1,1)
1644 write (iout,'(a6,i3,2i4)') "Helix",nhelix,ii1-1,j1-2
1645 if (nhelix.le.9) then
1646 write(12,'(a17,i1,a9,i3,a2,i3,a1)')
1647 & "DefPropRes 'helix",nhelix,
1648 & "' 'num = ",ii1-1,"..",j1-2,"'"
1650 write(12,'(a17,i2,a9,i3,a2,i3,a1)')
1651 & "DefPropRes 'helix",nhelix,
1652 & "' 'num = ",ii1-1,"..",j1-2,"'"
1659 if (nhelix.gt.0.and.lprint) then
1660 write(12,'(a26,$)') "DefPropRes 'helix' 'helix1"
1662 if (nhelix.le.9) then
1663 write(12,'(a8,i1,$)') " | helix",i
1665 write(12,'(a8,i2,$)') " | helix",i
1668 write(12,'(a1)') "'"
1672 write(12,'(a37)') "DefPropRes 'coil' '! (helix | sheet)'"
1673 write(12,'(a20)') "XMacStand ribbon.mac"
1679 c----------------------------------------------------------------------------
1681 subroutine write_pdb(npdb,titelloc,ee)
1682 implicit real*8 (a-h,o-z)
1683 include 'DIMENSIONS'
1684 include 'COMMON.IOUNITS'
1685 character*50 titelloc1
1686 character*(*) titelloc
1695 if (npdb.lt.1000) then
1696 call numstr(npdb,zahl)
1697 open(ipdb,file=prefix(:lenpre)//'@@'//zahl//'.pdb')
1699 if (npdb.lt.10000) then
1700 write(liczba5,'(i1,i4)') 0,npdb
1702 write(liczba5,'(i5)') npdb
1704 open(ipdb,file=prefix(:lenpre)//'@@'//liczba5//'.pdb')
1706 call pdbout(ee,titelloc1,ipdb)
1711 c-----------------------------------------------------------
1712 subroutine contact_cp2(var,var2,iff,ieval,in_pdb)
1713 implicit real*8 (a-h,o-z)
1714 include 'DIMENSIONS'
1718 include 'COMMON.SBRIDGE'
1719 include 'COMMON.FFIELD'
1720 include 'COMMON.IOUNITS'
1721 include 'COMMON.DISTFIT'
1722 include 'COMMON.VAR'
1723 include 'COMMON.CHAIN'
1724 include 'COMMON.MINIM'
1728 double precision var(maxvar),var2(maxvar)
1729 double precision time0,time1
1730 integer iff(maxres),ieval
1731 double precision theta1(maxres),phi1(maxres),alph1(maxres),
1735 call var_to_geom(nvar,var)
1742 if ( iff(i).eq.1.and.iff(j).eq.1 ) then
1764 call var_to_geom(nvar,var2)
1767 if ( iff(i).eq.1 ) then
1776 cd call write_pdb(3,'combined structure',0d0)
1777 cd time0=MPI_WTIME()
1780 NY=((NRES-4)*(NRES-5))/2
1781 call distfit(.true.,200)
1783 cd time1=MPI_WTIME()
1784 cd write (iout,'(a,f6.2,a)') ' Time for distfit ',time1-time0,' sec'
1794 call geom_to_var(nvar,var)
1795 cd time0=MPI_WTIME()
1796 call minimize(etot,var,iretcode,nfun)
1797 write(iout,*)'SUMSL return code is',iretcode,' eval SOFT',nfun
1799 cd time1=MPI_WTIME()
1800 cd write (iout,'(a,f6.2,f8.2,a)')' Time for soft min.',time1-time0,
1801 cd & nfun/(time1-time0),' SOFT eval/s'
1802 call var_to_geom(nvar,var)
1808 if (iff(1).eq.1) then
1814 if ( iwsk.eq.0.and.iff(i-1).eq.0.and.iff(i).eq.1 ) then
1819 if ( iwsk.eq.1.and.iff(i-1).eq.1.and.iff(i).eq.0 ) then
1825 if (iff(nres).eq.1) then
1831 cd write(linia,'(a6,i3,a1,i3,a1,i3,a1,i3)')
1832 cd & "select",ij(1),"-",ij(2),
1833 cd & ",",ij(3),"-",ij(4)
1834 cd call write_pdb(in_pdb,linia,etot)
1840 cd time0=MPI_WTIME()
1841 call minimize(etot,var,iretcode,nfun)
1842 cd write(iout,*)'SUMSL DIST return code is',iretcode,' eval ',nfun
1845 cd time1=MPI_WTIME()
1846 cd write (iout,'(a,f6.2,f8.2,a)')' Time for DIST min.',time1-time0,
1847 cd & nfun/(time1-time0),' eval/s'
1848 cd call var_to_geom(nvar,var)
1850 cd call write_pdb(6,'dist structure',etot)
1860 c-----------------------------------------------------------