2df2658e34918619b8b9bfbd15c27db893047623
[unres.git] / unres.F
1 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
2 C                                                                              C
3 C                                U N R E S                                     C
4 C                                                                              C
5 C Program to carry out conformational search of proteins in an united-residue  C
6 C approximation.                                                               C
7 C                                                                              C
8 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
9       implicit real*8 (a-h,o-z)
10       include 'DIMENSIONS'
11
12
13 #ifdef MPI
14       include 'mpif.h'
15       include 'COMMON.SETUP'
16 #endif
17       include 'COMMON.TIME1'
18       include 'COMMON.INTERACT'
19       include 'COMMON.NAMES'
20       include 'COMMON.GEO'
21       include 'COMMON.HEADER'
22       include 'COMMON.CONTROL'
23       include 'COMMON.CONTACTS'
24       include 'COMMON.CHAIN'
25       include 'COMMON.VAR'
26       include 'COMMON.IOUNITS'
27       include 'COMMON.FFIELD'
28       include 'COMMON.REMD'
29       include 'COMMON.MD'
30       include 'COMMON.SBRIDGE'
31       double precision hrtime,mintime,sectime
32       character*64 text_mode_calc(-2:14) /'test',
33      & 'SC rotamer distribution',
34      & 'Energy evaluation or minimization',
35      & 'Regularization of PDB structure',
36      & 'Threading of a sequence on PDB structures',
37      & 'Monte Carlo (with minimization) ',
38      & 'Energy minimization of multiple conformations',
39      & 'Checking energy gradient',
40      & 'Entropic sampling Monte Carlo (with minimization)',
41      & 'Energy map',
42      & 'CSA calculations',
43      & 'Not used 9',
44      & 'Not used 10',
45      & 'Soft regularization of PDB structure',
46      & 'Mesoscopic molecular dynamics (MD) ',
47      & 'Not used 13',
48      & 'Replica exchange molecular dynamics (REMD)'/
49       external ilen
50
51 c      call memmon_print_usage()
52
53       call init_task
54       if (me.eq.king)
55      & write(iout,*)'### LAST MODIFIED  4/25/08 7:29PM by adam'  
56       if (me.eq.king) call cinfo
57 C Read force field parameters and job setup data
58       call readrtns
59 C
60       if (me.eq.king .or. .not. out1file) then
61        write (iout,'(2a/)') 
62      & text_mode_calc(modecalc)(:ilen(text_mode_calc(modecalc))),
63      & ' calculation.' 
64        if (minim) write (iout,'(a)') 
65      &  'Conformations will be energy-minimized.'
66        write (iout,'(80(1h*)/)') 
67       endif
68       call flush(iout)
69 C
70       if (modecalc.eq.-2) then
71         call test
72         stop
73       else if (modecalc.eq.-1) then
74         write(iout,*) "call check_sc_map next"
75         call check_bond
76         stop
77       endif
78 #ifdef MPI
79       if (fg_rank.gt.0) then
80 C Fine-grain slaves just do energy and gradient components.
81         call ergastulum ! slave workhouse in Latin
82       else
83 #endif
84       if (modecalc.eq.0) then
85         call exec_eeval_or_minim
86       else if (modecalc.eq.1) then
87         call exec_regularize
88       else if (modecalc.eq.2) then
89         call exec_thread
90       else if (modecalc.eq.3 .or. modecalc .eq.6) then
91         call exec_MC
92       else if (modecalc.eq.4) then
93         call exec_mult_eeval_or_minim
94       else if (modecalc.eq.5) then
95         call exec_checkgrad
96       else if (ModeCalc.eq.7) then
97         call exec_map
98       else if (ModeCalc.eq.8) then
99         call exec_CSA
100       else if (modecalc.eq.11) then
101         call exec_softreg
102       else if (modecalc.eq.12) then
103         call exec_MD
104       else if (modecalc.eq.14) then
105 #ifdef MPI
106         call exec_MREMD
107 #else
108         write (iout,*) "Need a parallel version to run MREMD."
109         stop
110 #endif
111       else
112         write (iout,'(a)') 'This calculation type is not supported',
113      &   ModeCalc
114       endif
115 #ifdef MPI
116       endif
117 C Finish task.
118       if (fg_rank.eq.0) call finish_task
119 c      call memmon_print_usage()
120 #ifdef TIMING
121        call print_detailed_timing
122 #endif
123       call MPI_Finalize(ierr)
124       stop 'Bye Bye...'
125 #else
126       call dajczas(tcpu(),hrtime,mintime,sectime)
127       stop '********** Program terminated normally.'
128 #endif
129       end
130 c--------------------------------------------------------------------------
131       subroutine exec_MD
132       include 'DIMENSIONS'
133 #ifdef MPI
134       include "mpif.h"
135 #endif
136       include 'COMMON.SETUP'
137       include 'COMMON.CONTROL'
138       include 'COMMON.IOUNITS'
139 c      if (me.eq.king .or. .not. out1file) then
140 c        write (iout,*) "Calling chainbuild"
141 c        call flush(iout)
142 c      endif
143       call chainbuild
144 c      if (me.eq.king .or. .not. out1file) then
145 c        write (iout,*) "Calling MD"
146 c        call flush(iout)
147 c      endif
148       call MD
149       return
150       end
151 c---------------------------------------------------------------------------
152 #ifdef MPI
153       subroutine exec_MREMD
154       include 'DIMENSIONS'
155 #ifdef MPI
156       include "mpif.h"
157 #endif
158       include 'COMMON.SETUP'
159       include 'COMMON.CONTROL'
160       include 'COMMON.IOUNITS'
161       include 'COMMON.REMD'
162       if (me.eq.king .or. .not. out1file)
163      &   write (iout,*) "Calling chainbuild"
164       call chainbuild
165       if (me.eq.king .or. .not. out1file)
166      &   write (iout,*) "Calling REMD"
167       if (remd_mlist) then 
168         call MREMD
169       else
170         do i=1,nrep
171           remd_m(i)=1
172         enddo
173         call MREMD
174       endif
175       return
176       end
177 #endif
178 c---------------------------------------------------------------------------
179       subroutine exec_eeval_or_minim
180       implicit real*8 (a-h,o-z)
181       include 'DIMENSIONS'
182 #ifdef MPI
183       include 'mpif.h'
184 #endif
185       include 'COMMON.SETUP'
186       include 'COMMON.TIME1'
187       include 'COMMON.INTERACT'
188       include 'COMMON.NAMES'
189       include 'COMMON.GEO'
190       include 'COMMON.HEADER'
191       include 'COMMON.CONTROL'
192       include 'COMMON.CONTACTS'
193       include 'COMMON.CHAIN'
194       include 'COMMON.VAR'
195       include 'COMMON.IOUNITS'
196       include 'COMMON.FFIELD'
197       include 'COMMON.REMD'
198       include 'COMMON.MD'
199       include 'COMMON.SBRIDGE'
200       common /srutu/ icall
201       double precision energy(0:n_ene)
202       double precision energy_long(0:n_ene),energy_short(0:n_ene)
203       double precision varia(maxvar)
204       if (indpdb.eq.0) call chainbuild
205 #ifdef MPI
206       time00=MPI_Wtime()
207 #else
208       time00=tcpu()
209 #endif
210       call chainbuild_cart
211       if (split_ene) then
212        print *,"Processor",myrank," after chainbuild"
213        icall=1
214        call etotal_long(energy_long(0))
215        write (iout,*) "Printing long range energy"
216        call enerprint(energy_long(0))
217        call etotal_short(energy_short(0))
218        write (iout,*) "Printing short range energy"
219        call enerprint(energy_short(0))
220        do i=0,n_ene
221          energy(i)=energy_long(i)+energy_short(i)
222          write (iout,*) i,energy_long(i),energy_short(i),energy(i)
223        enddo
224        write (iout,*) "Printing long+short range energy"
225        call enerprint(energy(0))
226       endif
227       call etotal(energy(0))
228 #ifdef MPI
229       time_ene=MPI_Wtime()-time00
230 #else 
231       time_ene=tcpu()-time00
232 #endif
233       write (iout,*) "Time for energy evaluation",time_ene
234       print *,"after etotal"
235       etota = energy(0)
236       etot =etota
237       call enerprint(energy(0))
238       call hairpin(.true.,nharp,iharp)
239       call secondary2(.true.)
240       if (minim) then
241 crc overlap test
242         if (overlapsc) then 
243           print *, 'Calling OVERLAP_SC'
244           call overlap_sc(fail)
245         endif 
246
247         if (searchsc) then 
248           call sc_move(2,nres-1,10,1d10,nft_sc,etot)
249           print *,'SC_move',nft_sc,etot
250           write(iout,*) 'SC_move',nft_sc,etot
251         endif 
252
253         if (dccart) then
254           print *, 'Calling MINIM_DC'
255 #ifdef MPI
256           time1=MPI_WTIME()
257 #else
258           time1=tcpu()
259 #endif
260           call minim_dc(etot,iretcode,nfun)
261           if(iretcode.eq.8) call check_ecartint
262         else
263           if (indpdb.ne.0) then 
264             call bond_regular
265             call chainbuild
266           endif
267           call geom_to_var(nvar,varia)
268           print *,'Calling MINIMIZE.'
269 #ifdef MPI
270           time1=MPI_WTIME()
271 #else
272           time1=tcpu()
273 #endif
274           call minimize(etot,varia,iretcode,nfun)
275         endif
276         print *,'SUMSL return code is',iretcode,' eval ',nfun
277 #ifdef MPI
278         evals=nfun/(MPI_WTIME()-time1)
279 #else
280         evals=nfun/(tcpu()-time1)
281 #endif
282         print *,'# eval/s',evals
283         print *,'refstr=',refstr
284         call hairpin(.true.,nharp,iharp)
285         call secondary2(.true.)
286         call etotal(energy(0))
287         etot = energy(0)
288         call enerprint(energy(0))
289
290         call intout
291         call briefout(0,etot)
292         if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
293           write (iout,'(a,i3)') 'SUMSL return code:',iretcode
294           write (iout,'(a,i20)') '# of energy evaluations:',nfun+1
295           write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals
296       else
297         print *,'refstr=',refstr
298         if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
299         call briefout(0,etot)
300       endif
301       if (outpdb) call pdbout(etot,titel(:32),ipdb)
302       if (outmol2) call mol2out(etot,titel(:32))
303       return
304       end
305 c---------------------------------------------------------------------------
306       subroutine exec_regularize
307       implicit real*8 (a-h,o-z)
308       include 'DIMENSIONS'
309 #ifdef MPI
310       include 'mpif.h'
311 #endif
312       include 'COMMON.SETUP'
313       include 'COMMON.TIME1'
314       include 'COMMON.INTERACT'
315       include 'COMMON.NAMES'
316       include 'COMMON.GEO'
317       include 'COMMON.HEADER'
318       include 'COMMON.CONTROL'
319       include 'COMMON.CONTACTS'
320       include 'COMMON.CHAIN'
321       include 'COMMON.VAR'
322       include 'COMMON.IOUNITS'
323       include 'COMMON.FFIELD'
324       include 'COMMON.REMD'
325       include 'COMMON.MD'
326       include 'COMMON.SBRIDGE'
327       double precision energy(0:n_ene)
328
329       call gen_dist_constr
330       call sc_conf
331       call intout
332       call regularize(nct-nnt+1,etot,rms,cref(1,nnt,1),iretcode)
333       call etotal(energy(0))
334       energy(0)=energy(0)-energy(14)
335       etot=energy(0)
336       call enerprint(energy(0))
337       call intout
338       call briefout(0,etot)
339       if (outpdb) call pdbout(etot,titel(:32),ipdb)
340       if (outmol2) call mol2out(etot,titel(:32))
341       if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
342       write (iout,'(a,i3)') 'SUMSL return code:',iretcode
343       return
344       end
345 c---------------------------------------------------------------------------
346       subroutine exec_thread
347       include 'DIMENSIONS'
348 #ifdef MP
349       include "mpif.h"
350 #endif
351       include "COMMON.SETUP"
352       call thread_seq
353       return
354       end
355 c---------------------------------------------------------------------------
356       subroutine exec_MC
357       implicit real*8 (a-h,o-z)
358       include 'DIMENSIONS'
359       character*10 nodeinfo
360       double precision varia(maxvar)
361 #ifdef MPI
362       include "mpif.h"
363 #endif
364       include "COMMON.SETUP"
365       include 'COMMON.CONTROL'
366       call mcm_setup
367       if (minim) then
368 #ifdef MPI
369         if (modecalc.eq.3) then
370           call do_mcm(ipar)
371         else
372           call entmcm
373         endif
374 #else
375         if (modecalc.eq.3) then
376           call do_mcm(ipar)
377         else
378           call entmcm
379         endif
380 #endif
381       else
382         call monte_carlo
383       endif
384       return
385       end
386 c---------------------------------------------------------------------------
387       subroutine exec_mult_eeval_or_minim
388       implicit real*8 (a-h,o-z)
389       include 'DIMENSIONS'
390 #ifdef MPI
391       include 'mpif.h'
392       dimension muster(mpi_status_size)
393 #endif
394       include 'COMMON.SETUP'
395       include 'COMMON.TIME1'
396       include 'COMMON.INTERACT'
397       include 'COMMON.NAMES'
398       include 'COMMON.GEO'
399       include 'COMMON.HEADER'
400       include 'COMMON.CONTROL'
401       include 'COMMON.CONTACTS'
402       include 'COMMON.CHAIN'
403       include 'COMMON.VAR'
404       include 'COMMON.IOUNITS'
405       include 'COMMON.FFIELD'
406       include 'COMMON.REMD'
407       include 'COMMON.MD'
408       include 'COMMON.SBRIDGE'
409       double precision varia(maxvar)
410       dimension ind(6)
411       double precision energy(0:max_ene)
412       logical eof
413       eof=.false.
414 #ifdef MPI
415       if(me.ne.king) then
416         call minim_mcmf
417         return
418       endif
419
420       close (intin)
421       open(intin,file=intinname,status='old')
422       write (istat,'(a5,20a12)')"#    ",
423      &  (wname(print_order(i)),i=1,nprint_ene)
424       if (refstr) then
425         write (istat,'(a5,20a12)')"#    ",
426      &   (ename(print_order(i)),i=1,nprint_ene),
427      &   "ETOT total","RMSD","nat.contact","nnt.contact"        
428       else
429         write (istat,'(a5,20a12)')"#    ",
430      &    (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
431       endif
432
433       if (.not.minim) then
434         do while (.not. eof)
435           if (read_cart) then
436             read (intin,'(e15.10,e15.5)',end=1100,err=1100) time,ene
437             call read_x(intin,*11)
438 #ifdef MPI
439 c Broadcast the order to compute internal coordinates to the slaves.
440             if (nfgtasks.gt.1)
441      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
442 #endif
443             call int_from_cart1(.false.)
444           else
445             read (intin,'(i5)',end=1100,err=1100) iconf
446             call read_angles(intin,*11)
447             call geom_to_var(nvar,varia)
448             call chainbuild
449           endif
450           write (iout,'(a,i7)') 'Conformation #',iconf
451           call etotal(energy(0))
452           call briefout(iconf,energy(0))
453           call enerprint(energy(0))
454           etot=energy(0)
455           if (refstr) then 
456             call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
457             write (istat,'(i5,20(f12.3))') iconf,
458      &      (energy(print_order(i)),i=1,nprint_ene),etot,
459      &       rms,frac,frac_nn,co
460 cjlee end
461           else
462             write (istat,'(i5,16(f12.3))') iconf,
463      &     (energy(print_order(i)),i=1,nprint_ene),etot
464           endif
465         enddo
466 1100    continue
467         goto 1101
468       endif
469
470       mm=0
471       imm=0
472       nft=0
473       ene0=0.0d0
474       n=0
475       iconf=0
476 c      do n=1,nzsc
477       do while (.not. eof)
478         mm=mm+1
479         if (mm.lt.nodes) then
480           if (read_cart) then
481             read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
482             call read_x(intin,*11)
483 #ifdef MPI
484 c Broadcast the order to compute internal coordinates to the slaves.
485             if (nfgtasks.gt.1) 
486      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
487 #endif
488             call int_from_cart1(.false.)
489           else
490             read (intin,'(i5)',end=11,err=11) iconf
491             call read_angles(intin,*11)
492             call geom_to_var(nvar,varia)
493             call chainbuild
494           endif
495           write (iout,'(a,i7)') 'Conformation #',iconf
496           n=n+1
497          imm=imm+1
498          ind(1)=1
499          ind(2)=n
500          ind(3)=0
501          ind(4)=0
502          ind(5)=0
503          ind(6)=0
504          ene0=0.0d0
505          call mpi_send(ind,6,mpi_integer,mm,idint,CG_COMM,
506      *                  ierr)
507          call mpi_send(varia,nvar,mpi_double_precision,mm,
508      *                  idreal,CG_COMM,ierr)
509          call mpi_send(ene0,1,mpi_double_precision,mm,
510      *                  idreal,CG_COMM,ierr)
511 c         print *,'task ',n,' sent to worker ',mm,nvar
512         else
513          call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
514      *                 CG_COMM,muster,ierr)
515          man=muster(mpi_source)
516 c         print *,'receiving result from worker ',man,' (',iii1,iii,')'
517          call mpi_recv(varia,nvar,mpi_double_precision, 
518      *               man,idreal,CG_COMM,muster,ierr)
519          call mpi_recv(ene,1,
520      *               mpi_double_precision,man,idreal,
521      *               CG_COMM,muster,ierr)
522          call mpi_recv(ene0,1,
523      *               mpi_double_precision,man,idreal,
524      *               CG_COMM,muster,ierr)
525 c         print *,'result received from worker ',man,' sending now'
526
527           call var_to_geom(nvar,varia)
528           call chainbuild
529           call etotal(energy(0))
530           iconf=ind(2)
531           write (iout,*)
532           write (iout,*)
533           write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
534
535           etot=energy(0)
536           call enerprint(energy(0))
537           call briefout(it,etot)
538 c          if (minim) call briefout(it,etot)
539           if (refstr) then 
540             call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
541             write (istat,'(i5,19(f12.3))') iconf,
542      &     (energy(print_order(i)),i=1,nprint_ene),etot,
543      &     rms,frac,frac_nn,co
544           else
545             write (istat,'(i5,15(f12.3))') iconf,
546      &     (energy(print_order(i)),i=1,nprint_ene),etot
547           endif
548
549           imm=imm-1
550           if (read_cart) then
551             read (intin,'(e15.10,e15.5)',end=1101,err=1101) time,ene
552             call read_x(intin,*11)
553 #ifdef MPI
554 c Broadcast the order to compute internal coordinates to the slaves.
555             if (nfgtasks.gt.1)
556      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
557 #endif
558             call int_from_cart1(.false.)
559           else
560             read (intin,'(i5)',end=1101,err=1101) iconf
561             call read_angles(intin,*11)
562             call geom_to_var(nvar,varia)
563             call chainbuild
564           endif
565           n=n+1
566           imm=imm+1
567           ind(1)=1
568           ind(2)=n
569           ind(3)=0
570           ind(4)=0
571           ind(5)=0
572           ind(6)=0
573           call mpi_send(ind,6,mpi_integer,man,idint,CG_COMM,
574      *                  ierr)
575           call mpi_send(varia,nvar,mpi_double_precision,man, 
576      *                  idreal,CG_COMM,ierr)
577           call mpi_send(ene0,1,mpi_double_precision,man,
578      *                  idreal,CG_COMM,ierr)
579           nf_mcmf=nf_mcmf+ind(4)
580           nmin=nmin+1
581         endif
582       enddo
583 11    continue
584       do j=1,imm
585         call mpi_recv(ind,6,mpi_integer,mpi_any_source,idint,
586      *               CG_COMM,muster,ierr)
587         man=muster(mpi_source)
588         call mpi_recv(varia,nvar,mpi_double_precision, 
589      *               man,idreal,CG_COMM,muster,ierr)
590         call mpi_recv(ene,1,
591      *               mpi_double_precision,man,idreal,
592      *               CG_COMM,muster,ierr)
593         call mpi_recv(ene0,1,
594      *               mpi_double_precision,man,idreal,
595      *               CG_COMM,muster,ierr)
596
597         call var_to_geom(nvar,varia)
598         call chainbuild
599         call etotal(energy(0))
600         iconf=ind(2)
601         write (iout,*)
602         write (iout,*)
603         write (iout,'(a,2i7)') 'Conformation #',iconf,ind(5)
604
605         etot=energy(0)
606         call enerprint(energy(0))
607         call briefout(it,etot)
608         if (refstr) then 
609           call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
610           write (istat,'(i5,19(f12.3))') iconf,
611      &   (energy(print_order(i)),i=1,nprint_ene),etot,
612      &   rms,frac,frac_nn,co
613         else
614           write (istat,'(i5,15(f12.3))') iconf,
615      &    (energy(print_order(i)),i=1,nprint_ene),etot
616         endif
617         nmin=nmin+1
618       enddo
619 1101  continue
620       do i=1, nodes-1
621          ind(1)=0
622          ind(2)=0
623          ind(3)=0
624          ind(4)=0
625          ind(5)=0
626          ind(6)=0
627          call mpi_send(ind,6,mpi_integer,i,idint,CG_COMM,
628      *                  ierr)
629       enddo
630 #else
631       close (intin)
632       open(intin,file=intinname,status='old')
633       write (istat,'(a5,20a12)')"#    ",
634      &   (wname(print_order(i)),i=1,nprint_ene)
635       write (istat,'("#    ",20(1pe12.4))')
636      &   (weights(print_order(i)),i=1,nprint_ene)
637       if (refstr) then
638         write (istat,'(a5,20a12)')"#    ",
639      &   (ename(print_order(i)),i=1,nprint_ene),
640      &   "ETOT total","RMSD","nat.contact","nnt.contact"
641       else
642         write (istat,'(a5,14a12)')"#    ",
643      &   (ename(print_order(i)),i=1,nprint_ene),"ETOT total"
644       endif
645       do while (.not. eof)
646           if (read_cart) then
647             read (intin,'(e15.10,e15.5)',end=11,err=11) time,ene
648             call read_x(intin,*11)
649 #ifdef MPI
650 c Broadcast the order to compute internal coordinates to the slaves.
651             if (nfgtasks.gt.1)
652      &        call MPI_Bcast(6,1,MPI_INTEGER,king,FG_COMM,IERROR)
653 #endif
654             call int_from_cart1(.false.)
655           else
656             read (intin,'(i5)',end=11,err=11) iconf
657             call read_angles(intin,*11)
658             call geom_to_var(nvar,varia)
659             call chainbuild
660           endif
661         write (iout,'(a,i7)') 'Conformation #',iconf
662         if (minim) call minimize(etot,varia,iretcode,nfun)
663         call etotal(energy(0))
664
665         etot=energy(0)
666         call enerprint(energy(0))
667         if (minim) call briefout(it,etot) 
668         if (refstr) then 
669           call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
670           write (istat,'(i5,18(f12.3))') iconf,
671      &   (energy(print_order(i)),i=1,nprint_ene),
672      &   etot,rms,frac,frac_nn,co
673 cjlee end
674         else
675           write (istat,'(i5,14(f12.3))') iconf,
676      &   (energy(print_order(i)),i=1,nprint_ene),etot
677         endif
678       enddo
679    11 continue
680 #endif
681       return
682       end
683 c---------------------------------------------------------------------------
684       subroutine exec_checkgrad
685       implicit real*8 (a-h,o-z)
686       include 'DIMENSIONS'
687 #ifdef MPI
688       include 'mpif.h'
689 #endif
690       include 'COMMON.SETUP'
691       include 'COMMON.TIME1'
692       include 'COMMON.INTERACT'
693       include 'COMMON.NAMES'
694       include 'COMMON.GEO'
695       include 'COMMON.HEADER'
696       include 'COMMON.CONTROL'
697       include 'COMMON.CONTACTS'
698       include 'COMMON.CHAIN'
699       include 'COMMON.VAR'
700       include 'COMMON.IOUNITS'
701       include 'COMMON.FFIELD'
702       include 'COMMON.REMD'
703       include 'COMMON.MD'
704       include 'COMMON.SBRIDGE'
705       common /srutu/ icall
706       double precision energy(0:max_ene)
707 c      do i=2,nres
708 c        vbld(i)=vbld(i)+ran_number(-0.1d0,0.1d0)
709 c        if (itype(i).ne.10) 
710 c     &      vbld(i+nres)=vbld(i+nres)+ran_number(-0.001d0,0.001d0)
711 c      enddo
712       if (indpdb.eq.0) call chainbuild
713 c      do i=0,nres
714 c        do j=1,3
715 c          dc(j,i)=dc(j,i)+ran_number(-0.2d0,0.2d0)
716 c        enddo
717 c      enddo
718 c      do i=1,nres-1
719 c        if (itype(i).ne.10) then
720 c          do j=1,3
721 c            dc(j,i+nres)=dc(j,i+nres)+ran_number(-0.2d0,0.2d0)
722 c          enddo
723 c        endif
724 c      enddo
725 c      do j=1,3
726 c        dc(j,0)=ran_number(-0.2d0,0.2d0)
727 c      enddo
728       usampl=.true.
729       totT=1.d0
730       eq_time=0.0d0
731       call read_fragments
732       call chainbuild_cart
733       call cartprint
734       call intout
735       icall=1
736       call etotal(energy(0))
737       etot = energy(0)
738       call enerprint(energy(0))
739       write (iout,*) "Uconst",Uconst," Uconst_back",uconst_back
740       print *,'icheckgrad=',icheckgrad
741       goto (10,20,30) icheckgrad
742   10  call check_ecartint
743       write(iout,*) "kupadupa"
744       call check_ecartint
745       return
746   20  call check_cartgrad
747       return
748   30  call check_eint
749       return
750       end
751 c---------------------------------------------------------------------------
752       subroutine exec_map
753 C Energy maps
754       call map_read
755       call map
756       return
757       end
758 c---------------------------------------------------------------------------
759       subroutine exec_CSA
760 #ifdef MPI
761       include "mpif.h"
762 #endif
763       include 'DIMENSIONS'
764       include 'COMMON.IOUNITS'
765 C Conformational Space Annealling programmed by Jooyoung Lee.
766 C This method works only with parallel machines!
767 #ifdef MPI
768       call together
769 #else
770       write (iout,*) "CSA works on parallel machines only"
771 #endif
772       return
773       end
774 c---------------------------------------------------------------------------
775       subroutine exec_softreg
776       include 'DIMENSIONS'
777       include 'COMMON.IOUNITS'
778       include 'COMMON.CONTROL'
779       double precision energy(0:max_ene)
780       call chainbuild
781       call etotal(energy(0))
782       call enerprint(energy(0))
783       if (.not.lsecondary) then
784         write(iout,*) 'Calling secondary structure recognition'
785         call secondary2(debug)
786       else
787         write(iout,*) 'Using secondary structure supplied in pdb'
788       endif
789
790       call softreg
791
792       call etotal(energy(0))
793       etot=energy(0)
794       call enerprint(energy(0))
795       call intout
796       call briefout(0,etot)
797       call secondary2(.true.)
798       if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.)
799       return
800       end