Rozgrzebany SCCOR dla wham-M
[unres.git] / source / wham / src-M / wham_multparm.F
1       program WHAM_multparm
2 c Creation/update of the database of conformations
3       implicit none
4 #ifndef ISNAN
5       external proc_proc
6 #endif
7 #ifdef WINPGI
8 cMS$ATTRIBUTES C ::  proc_proc
9 #endif
10       include "DIMENSIONS"
11       include "DIMENSIONS.ZSCOPT"
12       include "DIMENSIONS.FREE"
13 #ifdef MPI
14       include "mpif.h"
15       integer IERROR,ERRCODE
16       include "COMMON.MPI"
17 #endif
18       include "COMMON.IOUNITS"
19       include "COMMON.FREE"
20       include "COMMON.CONTROL"
21       include "COMMON.ALLPARM"
22       include "COMMON.PROT"
23       double precision rr,x(max_paropt)
24       integer idumm
25       integer i,ipar,islice
26 #ifdef MPI
27       call MPI_Init( IERROR )
28       call MPI_Comm_rank( MPI_COMM_WORLD, me, IERROR )
29       call MPI_Comm_size( MPI_COMM_WORLD, nprocs, IERROR )
30       Master = 0
31       if (ierror.gt.0) then
32         write(iout,*) "SEVERE ERROR - Can't initialize MPI."
33         call mpi_finalize(ierror)
34         stop
35       endif
36       if (nprocs.gt.MaxProcs+1) then
37         write (2,*) "Error - too many processors",
38      &   nprocs,MaxProcs+1
39         write (2,*) "Increase MaxProcs and recompile"
40         call MPI_Finalize(IERROR)
41         stop
42       endif
43 #endif
44 c NaNQ initialization
45 #ifndef ISNAN
46       i=-1
47       rr=dacos(100.0d0)
48 #ifdef WINPGI
49       idumm=proc_proc(rr,i)
50 #else
51       call proc_proc(rr,i)
52 #endif
53 #endif
54       call initialize
55       call openunits
56       call cinfo
57       call read_general_data(*10)
58       call flush(iout)
59       call molread(*10)
60       call flush(iout)
61 #ifdef MPI 
62       write (iout,*) "Calling proc_groups"
63       call proc_groups
64       write (iout,*) "proc_groups exited"
65       call flush(iout)
66 #endif
67       do ipar=1,nParmSet
68         write (iout,*) "Calling parmread",ipar
69         call parmread(ipar,*10)
70         if (.not.separate_parset) then
71           call store_parm(ipar)
72           write (iout,*) "Finished storing parameters",ipar
73         else if (ipar.eq.myparm) then
74           call store_parm(1)
75           write (iout,*) "Finished storing parameters",ipar
76         endif
77         call flush(iout)
78       enddo
79       call read_efree(*10)
80       write (iout,*) "Finished READ_EFREE"
81       call flush(iout)
82       call read_protein_data(*10)
83       write (iout,*) "Finished READ_PROTEIN_DATA"
84       call flush(iout)
85       if (indpdb.gt.0) then
86         call promienie
87         call read_compar
88         call read_ref_structure(*10)
89         call proc_cont
90         call fragment_list
91         if (constr_dist.gt.0) call read_dist_constr
92       endif
93       write (iout,*) "Begin read_database"
94       call flush(iout)
95       call read_database(*10)
96       write (iout,*) "Finished read_database"
97       call flush(iout)
98       if (separate_parset) nparmset=1
99       do islice=1,nslice
100         if (ntot(islice).gt.0) then
101 #ifdef MPI 
102           call work_partition(islice,.true.)
103           write (iout,*) "work_partition OK"
104           call flush(iout)
105 #endif
106           call enecalc(islice,*10)
107           write (iout,*) "enecalc OK"
108           call flush(iout)
109           call WHAM_CALC(islice,*10)
110           write (iout,*) "wham_calc OK"
111           call flush(iout)
112           call write_dbase(islice,*10)
113           write (iout,*) "write_dbase OK"
114           call flush(iout)
115           if (ensembles.gt.0) then
116             call make_ensembles(islice,*10)
117             write (iout,*) "make_ensembles OK"
118             call flush(iout)
119           endif
120         endif
121       enddo
122 #ifdef MPI
123       call MPI_Finalize( IERROR )
124 #endif
125       stop
126    10 write (iout,*) "Error termination of the program"
127       call MPI_Finalize( IERROR )
128       stop
129       end
130 c------------------------------------------------------------------------------
131 #ifdef MPI
132       subroutine proc_groups
133 C Split the processors into the Master and Workers group, if needed.
134       implicit none
135       include "DIMENSIONS"
136       include "DIMENSIONS.ZSCOPT"
137       include "DIMENSIONS.FREE"
138       include "mpif.h"
139       include "COMMON.IOUNITS"
140       include "COMMON.MPI"
141       include "COMMON.FREE"
142       integer n,chunk,i,j,ii,remainder
143       integer kolor,key,ierror,errcode
144       logical lprint
145       lprint=.true.
146
147 C Split the communicator if independent runs for different parameter
148 C sets will be performed.
149 C
150       if (nparmset.eq.1 .or. .not.separate_parset) then
151         WHAM_COMM = MPI_COMM_WORLD
152       else if (separate_parset) then
153         if (nprocs.lt.nparmset) then
154           write (iout,*) 
155      & "*** Cannot split parameter sets for fewer processors than sets",
156      &  nprocs,nparmset
157           call MPI_Finalize(ierror)
158           stop
159         endif 
160         write (iout,*) "nparmset",nparmset
161         nprocs = nprocs/nparmset
162         kolor = me/nprocs
163         key = mod(me,nprocs)
164         write (iout,*) "My old rank",me," kolor",kolor," key",key
165         call MPI_Comm_split(MPI_COMM_WORLD,kolor,key,WHAM_COMM,ierror)
166         call MPI_Comm_size(WHAM_COMM,nprocs,ierror)
167         call MPI_Comm_rank(WHAM_COMM,me,ierror)
168         write (iout,*) "My new rank",me," comm size",nprocs
169         write (iout,*) "MPI_COMM_WORLD",MPI_COMM_WORLD,
170      &   " WHAM_COMM",WHAM_COMM
171         myparm=kolor+1
172         write (iout,*) "My parameter set is",myparm
173         call flush(iout)
174       else
175         myparm=nparmset
176       endif
177       Me1 = Me
178       Nprocs1 = Nprocs
179       return
180       end
181 c------------------------------------------------------------------------------
182       subroutine work_partition(islice,lprint)
183 c Split the conformations between processors
184       implicit none
185       include "DIMENSIONS"
186       include "DIMENSIONS.ZSCOPT"
187       include "DIMENSIONS.FREE"
188       include "mpif.h"
189       include "COMMON.IOUNITS"
190       include "COMMON.MPI"
191       include "COMMON.PROT"
192       integer islice
193       integer n,chunk,i,j,ii,remainder
194       integer kolor,key,ierror,errcode
195       logical lprint
196 C
197 C Divide conformations between processors; the first and
198 C the last conformation to handle by ith processor is stored in 
199 C indstart(i) and indend(i), respectively.
200 C
201 C First try to assign equal number of conformations to each processor.
202 C
203         n=ntot(islice)
204         write (iout,*) "n=",n
205         indstart(0)=1
206         chunk = N/nprocs1
207         scount(0) = chunk
208 c        print *,"i",0," indstart",indstart(0)," scount",
209 c     &     scount(0)
210         do i=1,nprocs1-1
211           indstart(i)=chunk+indstart(i-1) 
212           scount(i)=scount(i-1)
213 c          print *,"i",i," indstart",indstart(i)," scount",
214 c     &     scount(i)
215         enddo 
216 C
217 C Determine how many conformations remained yet unassigned.
218 C
219         remainder=N-(indstart(nprocs1-1)
220      &    +scount(nprocs1-1)-1)
221 c        print *,"remainder",remainder
222 C
223 C Assign the remainder conformations to consecutive processors, starting
224 C from the lowest rank; this continues until the list is exhausted.
225 C
226         if (remainder .gt. 0) then 
227           do i=1,remainder
228             scount(i-1) = scount(i-1) + 1
229             indstart(i) = indstart(i) + i
230           enddo
231           do i=remainder+1,nprocs1-1
232             indstart(i) = indstart(i) + remainder
233           enddo
234         endif
235
236         indstart(nprocs1)=N+1
237         scount(nprocs1)=0
238
239         do i=0,NProcs1
240           indend(i)=indstart(i)+scount(i)-1
241           idispl(i)=indstart(i)-1
242         enddo
243
244         N=0
245         do i=0,Nprocs1-1
246           N=N+indend(i)-indstart(i)+1
247         enddo
248
249 c        print *,"N",n," NTOT",ntot(islice)
250         if (N.ne.ntot(islice)) then
251           write (iout,*) "!!! Checksum error on processor",me,
252      &     " slice",islice
253           call flush(iout)
254           call MPI_Abort( MPI_COMM_WORLD, Ierror, Errcode )
255         endif
256
257       if (lprint) then
258         write (iout,*) "Partition of work between processors"
259           do i=0,nprocs1-1
260             write (iout,'(a,i5,a,i7,a,i7,a,i7)')
261      &        "Processor",i," indstart",indstart(i),
262      &        " indend",indend(i)," count",scount(i)
263           enddo
264       endif
265       return
266       end
267 #endif
268 #ifdef AIX
269       subroutine flush(iu)
270       call flush_(iu)
271       return
272       end
273 #endif