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