debuging WHAM for constrains DEBUG ON
[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       if (constr_dist.gt.0) call read_dist_constr
62 #ifdef MPI 
63       write (iout,*) "Calling proc_groups"
64       call proc_groups
65       write (iout,*) "proc_groups exited"
66       call flush(iout)
67 #endif
68       do ipar=1,nParmSet
69         write (iout,*) "Calling parmread",ipar
70         call parmread(ipar,*10)
71         if (.not.separate_parset) then
72           call store_parm(ipar)
73           write (iout,*) "Finished storing parameters",ipar
74         else if (ipar.eq.myparm) then
75           call store_parm(1)
76           write (iout,*) "Finished storing parameters",ipar
77         endif
78         call flush(iout)
79       enddo
80       call read_efree(*10)
81       write (iout,*) "Finished READ_EFREE"
82       call flush(iout)
83       call read_protein_data(*10)
84       write (iout,*) "Finished READ_PROTEIN_DATA"
85       call flush(iout)
86       if (indpdb.gt.0) then
87         call promienie
88         call read_compar
89         call read_ref_structure(*10)
90         call proc_cont
91         call fragment_list
92       endif
93 C      if (constr_dist.gt.0) call read_dist_constr
94       write (iout,*) "Begin read_database"
95       call flush(iout)
96       call read_database(*10)
97       write (iout,*) "Finished read_database"
98       call flush(iout)
99       if (separate_parset) nparmset=1
100       do islice=1,nslice
101         if (ntot(islice).gt.0) then
102 #ifdef MPI 
103           call work_partition(islice,.true.)
104           write (iout,*) "work_partition OK"
105           call flush(iout)
106 #endif
107           call enecalc(islice,*10)
108           write (iout,*) "enecalc OK"
109           call flush(iout)
110           call WHAM_CALC(islice,*10)
111           write (iout,*) "wham_calc OK"
112           call flush(iout)
113           call write_dbase(islice,*10)
114           write (iout,*) "write_dbase OK"
115           call flush(iout)
116           if (ensembles.gt.0) then
117             call make_ensembles(islice,*10)
118             write (iout,*) "make_ensembles OK"
119             call flush(iout)
120           endif
121         endif
122       enddo
123 #ifdef MPI
124       call MPI_Finalize( IERROR )
125 #endif
126       stop
127    10 write (iout,*) "Error termination of the program"
128       call MPI_Finalize( IERROR )
129       stop
130       end
131 c------------------------------------------------------------------------------
132 #ifdef MPI
133       subroutine proc_groups
134 C Split the processors into the Master and Workers group, if needed.
135       implicit none
136       include "DIMENSIONS"
137       include "DIMENSIONS.ZSCOPT"
138       include "DIMENSIONS.FREE"
139       include "mpif.h"
140       include "COMMON.IOUNITS"
141       include "COMMON.MPI"
142       include "COMMON.FREE"
143       integer n,chunk,i,j,ii,remainder
144       integer kolor,key,ierror,errcode
145       logical lprint
146       lprint=.true.
147
148 C Split the communicator if independent runs for different parameter
149 C sets will be performed.
150 C
151       if (nparmset.eq.1 .or. .not.separate_parset) then
152         WHAM_COMM = MPI_COMM_WORLD
153       else if (separate_parset) then
154         if (nprocs.lt.nparmset) then
155           write (iout,*) 
156      & "*** Cannot split parameter sets for fewer processors than sets",
157      &  nprocs,nparmset
158           call MPI_Finalize(ierror)
159           stop
160         endif 
161         write (iout,*) "nparmset",nparmset
162         nprocs = nprocs/nparmset
163         kolor = me/nprocs
164         key = mod(me,nprocs)
165         write (iout,*) "My old rank",me," kolor",kolor," key",key
166         call MPI_Comm_split(MPI_COMM_WORLD,kolor,key,WHAM_COMM,ierror)
167         call MPI_Comm_size(WHAM_COMM,nprocs,ierror)
168         call MPI_Comm_rank(WHAM_COMM,me,ierror)
169         write (iout,*) "My new rank",me," comm size",nprocs
170         write (iout,*) "MPI_COMM_WORLD",MPI_COMM_WORLD,
171      &   " WHAM_COMM",WHAM_COMM
172         myparm=kolor+1
173         write (iout,*) "My parameter set is",myparm
174         call flush(iout)
175       else
176         myparm=nparmset
177       endif
178       Me1 = Me
179       Nprocs1 = Nprocs
180       return
181       end
182 c------------------------------------------------------------------------------
183       subroutine work_partition(islice,lprint)
184 c Split the conformations between processors
185       implicit none
186       include "DIMENSIONS"
187       include "DIMENSIONS.ZSCOPT"
188       include "DIMENSIONS.FREE"
189       include "mpif.h"
190       include "COMMON.IOUNITS"
191       include "COMMON.MPI"
192       include "COMMON.PROT"
193       integer islice
194       integer n,chunk,i,j,ii,remainder
195       integer kolor,key,ierror,errcode
196       logical lprint
197 C
198 C Divide conformations between processors; the first and
199 C the last conformation to handle by ith processor is stored in 
200 C indstart(i) and indend(i), respectively.
201 C
202 C First try to assign equal number of conformations to each processor.
203 C
204         n=ntot(islice)
205         write (iout,*) "n=",n
206         indstart(0)=1
207         chunk = N/nprocs1
208         scount(0) = chunk
209 c        print *,"i",0," indstart",indstart(0)," scount",
210 c     &     scount(0)
211         do i=1,nprocs1-1
212           indstart(i)=chunk+indstart(i-1) 
213           scount(i)=scount(i-1)
214 c          print *,"i",i," indstart",indstart(i)," scount",
215 c     &     scount(i)
216         enddo 
217 C
218 C Determine how many conformations remained yet unassigned.
219 C
220         remainder=N-(indstart(nprocs1-1)
221      &    +scount(nprocs1-1)-1)
222 c        print *,"remainder",remainder
223 C
224 C Assign the remainder conformations to consecutive processors, starting
225 C from the lowest rank; this continues until the list is exhausted.
226 C
227         if (remainder .gt. 0) then 
228           do i=1,remainder
229             scount(i-1) = scount(i-1) + 1
230             indstart(i) = indstart(i) + i
231           enddo
232           do i=remainder+1,nprocs1-1
233             indstart(i) = indstart(i) + remainder
234           enddo
235         endif
236
237         indstart(nprocs1)=N+1
238         scount(nprocs1)=0
239
240         do i=0,NProcs1
241           indend(i)=indstart(i)+scount(i)-1
242           idispl(i)=indstart(i)-1
243         enddo
244
245         N=0
246         do i=0,Nprocs1-1
247           N=N+indend(i)-indstart(i)+1
248         enddo
249
250 c        print *,"N",n," NTOT",ntot(islice)
251         if (N.ne.ntot(islice)) then
252           write (iout,*) "!!! Checksum error on processor",me,
253      &     " slice",islice
254           call flush(iout)
255           call MPI_Abort( MPI_COMM_WORLD, Ierror, Errcode )
256         endif
257
258       if (lprint) then
259         write (iout,*) "Partition of work between processors"
260           do i=0,nprocs1-1
261             write (iout,'(a,i5,a,i7,a,i7,a,i7)')
262      &        "Processor",i," indstart",indstart(i),
263      &        " indend",indend(i)," count",scount(i)
264           enddo
265       endif
266       return
267       end
268 #endif
269 #ifdef AIX
270       subroutine flush(iu)
271       call flush_(iu)
272       return
273       end
274 #endif