1 subroutine minimize(etot,x,iretcode,nfun)
2 implicit real*8 (a-h,o-z)
4 parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2))
5 *********************************************************************
6 * OPTIMIZE sets up SUMSL or DFP and provides a simple interface for *
7 * the calling subprogram. *
8 * when d(i)=1.0, then v(35) is the length of the initial step, *
9 * calculated in the usual pythagorean way. *
10 * absolute convergence occurs when the function is within v(31) of *
11 * zero. unless you know the minimum value in advance, abs convg *
12 * is probably not useful. *
13 * relative convergence is when the model predicts that the function *
14 * will decrease by less than v(32)*abs(fun). *
15 *********************************************************************
16 include 'COMMON.IOUNITS'
19 include 'COMMON.MINIM'
20 include 'COMMON.CONTROL'
23 double precision minval,x(maxvar),d(maxvar),v(1:lv),xx(maxvar)
24 double precision energia(0:n_ene)
25 external func,gradient,fdum
26 external func_restr,grad_restr
27 logical not_done,change,reduce
28 c common /przechowalnia/ v
36 call deflt(2,iv,liv,lv,v)
37 * 12 means fresh start, dont call deflt
39 * max num of fun calls
40 if (maxfun.eq.0) maxfun=500
42 * max num of iterations
43 if (maxmin.eq.0) maxmin=1000
49 if (print_min_ini+print_min_stat+print_min_res.gt.0) iv(21)=iout
50 * 1 means to print out result
52 * 1 means to print out summary stats
54 * 1 means to print initial x and d
56 * min val for v(radfac) default is 0.1
58 * max val for v(radfac) default is 4.0
61 * check false conv if (act fnctn decrease) .lt. v(26)*(exp decrease)
62 * the sumsl default is 0.1
64 * false conv if (act fnctn decrease) .lt. v(34)
65 * the sumsl default is 100*machep
67 * absolute convergence
68 if (tolf.eq.0.0D0) tolf=1.0D-4
70 * relative convergence
71 if (rtolf.eq.0.0D0) rtolf=1.0D-4
73 * controls initial step size
75 * large vals of d correspond to small components of step
82 cd print *,'Calling SUMSL'
83 c call var_to_geom(nvar,x)
85 c call etotal(energia(0))
90 call x2xx(x,xx,nvar_restr)
91 call sumsl(nvar_restr,d,xx,func_restr,grad_restr,
92 & iv,liv,lv,v,idum,rdum,fdum)
95 call sumsl(nvar,d,x,func,gradient,iv,liv,lv,v,idum,rdum,fdum)
101 cd print *,'Exit SUMSL; return code:',iretcode,' energy:',etot
102 cd write (iout,'(/a,i4/)') 'SUMSL return code:',iv(1)
105 call var_to_geom(nvar,x)
107 c write (iout,'(a)') 'Reduction worked, minimizing again...'
112 c call etotal(energia(0))
114 c call enerprint(energia(0))
117 c write (*,*) 'Processor',MyID,' leaves MINIMIZE.'
124 c----------------------------------------------------------------------------
125 subroutine ergastulum
126 implicit real*8 (a-h,o-z)
131 include 'COMMON.SETUP'
132 include 'COMMON.DERIV'
134 include 'COMMON.IOUNITS'
135 include 'COMMON.FFIELD'
136 include 'COMMON.INTERACT'
138 include 'COMMON.TIME1'
139 double precision z(maxres6),d_a_tmp(maxres6)
140 double precision edum(0:n_ene),time_order(0:10)
141 double precision Gcopy(maxres2,maxres2)
142 common /przechowalnia/ Gcopy
144 C Workers wait for variables and NF, and NFL from the boss
146 do while (iorder.ge.0)
147 c write (*,*) 'Processor',fg_rank,' CG group',kolor,
148 c & ' receives order from Master'
150 call MPI_Bcast(iorder,1,MPI_INTEGER,king,FG_COMM,IERR)
151 time_Bcast=time_Bcast+MPI_Wtime()-time00
152 if (icall.gt.4 .and. iorder.ge.0)
153 & time_order(iorder)=time_order(iorder)+MPI_Wtime()-time00
156 c & 'Processor',fg_rank,' completed receive MPI_BCAST order',iorder
157 if (iorder.eq.0) then
160 c write (2,*) "After etotal"
161 c write (2,*) "dimen",dimen," dimen3",dimen3
163 else if (iorder.eq.2) then
165 call etotal_short(edum)
166 c write (2,*) "After etotal_short"
167 c write (2,*) "dimen",dimen," dimen3",dimen3
169 else if (iorder.eq.3) then
171 call etotal_long(edum)
172 c write (2,*) "After etotal_long"
173 c write (2,*) "dimen",dimen," dimen3",dimen3
175 else if (iorder.eq.1) then
177 c write (2,*) "After sum_gradient"
178 c write (2,*) "dimen",dimen," dimen3",dimen3
180 else if (iorder.eq.4) then
181 call ginv_mult(z,d_a_tmp)
182 else if (iorder.eq.5) then
183 c Setup MD things for a slave
184 dimen=(nct-nnt+1)+nside
185 dimen1=(nct-nnt)+(nct-nnt+1)
187 c write (2,*) "dimen",dimen," dimen3",dimen3
189 call int_bounds(dimen,igmult_start,igmult_end)
190 igmult_start=igmult_start-1
191 call MPI_Allgather(3*igmult_start,1,MPI_INTEGER,
192 & ng_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
193 my_ng_count=igmult_end-igmult_start
194 call MPI_Allgather(3*my_ng_count,1,MPI_INTEGER,ng_counts(0),1,
195 & MPI_INTEGER,FG_COMM,IERROR)
196 c write (2,*) "ng_start",(ng_start(i),i=0,nfgtasks-1)
197 c write (2,*) "ng_counts",(ng_counts(i),i=0,nfgtasks-1)
198 myginv_ng_count=maxres2*my_ng_count
199 c write (2,*) "igmult_start",igmult_start," igmult_end",
200 c & igmult_end," my_ng_count",my_ng_count
202 call MPI_Allgather(maxres2*igmult_start,1,MPI_INTEGER,
203 & nginv_start(0),1,MPI_INTEGER,FG_COMM,IERROR)
204 call MPI_Allgather(myginv_ng_count,1,MPI_INTEGER,
205 & nginv_counts(0),1,MPI_INTEGER,FG_COMM,IERROR)
206 c write (2,*) "nginv_start",(nginv_start(i),i=0,nfgtasks-1)
207 c write (2,*) "nginv_counts",(nginv_counts(i),i=0,nfgtasks-1)
209 c call MPI_Barrier(FG_COMM,IERROR)
211 call MPI_Scatterv(ginv(1,1),nginv_counts(0),
212 & nginv_start(0),MPI_DOUBLE_PRECISION,gcopy(1,1),
213 & myginv_ng_count,MPI_DOUBLE_PRECISION,king,FG_COMM,IERR)
215 time_scatter_ginv=time_scatter_ginv+MPI_Wtime()-time00
222 c write (2,*) "dimen",dimen," dimen3",dimen3
223 c write (2,*) "End MD setup"
225 c write (iout,*) "My chunk of ginv_block"
226 c call MATOUT2(my_ng_count,dimen3,maxres2,maxers2,ginv_block)
227 else if (iorder.eq.6) then
228 call int_from_cart1(.false.)
229 else if (iorder.eq.7) then
231 else if (iorder.eq.8) then
233 else if (iorder.eq.9) then
234 call fricmat_mult(z,d_a_tmp)
235 else if (iorder.eq.10) then
239 write (*,*) 'Processor',fg_rank,' CG group',kolor,
240 & ' absolute rank',myrank,' leves ERGASTULUM.'
241 write(*,*)'Processor',fg_rank,' wait times for respective orders',
242 & (' order[',i,']',time_order(i),i=0,10)
246 ************************************************************************
247 subroutine func(n,x,nf,f,uiparm,urparm,ufparm)
248 implicit real*8 (a-h,o-z)
250 include 'COMMON.DERIV'
251 include 'COMMON.IOUNITS'
254 double precision energia(0:n_ene)
256 double precision ufparm
262 c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n)
266 cd print *,'func',nf,nfl,icg
267 call var_to_geom(n,x)
270 cd write (iout,*) 'ETOTAL called from FUNC'
271 call etotal(energia(0))
275 c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n)
276 c write (iout,*) 'f=',etot
281 ************************************************************************
282 subroutine func_restr(n,x,nf,f,uiparm,urparm,ufparm)
283 implicit real*8 (a-h,o-z)
285 include 'COMMON.DERIV'
286 include 'COMMON.IOUNITS'
289 double precision energia(0:n_ene)
291 double precision ufparm
297 c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n)
301 call var_to_geom_restr(n,x)
304 cd write (iout,*) 'ETOTAL called from FUNC'
305 call etotal(energia(0))
309 c write (iout,'(10f8.3)') (rad2deg*x(i),i=1,n)
310 c write (iout,*) 'f=',etot
315 c-------------------------------------------------------
316 subroutine x2xx(x,xx,n)
317 implicit real*8 (a-h,o-z)
320 include 'COMMON.CHAIN'
321 include 'COMMON.INTERACT'
322 double precision xx(maxvar),x(maxvar)
332 if (mask_phi(i).eq.1) then
340 if (mask_theta(i).eq.1) then
348 if (itype(i).ne.10) then
350 if (mask_side(i).eq.1) then
363 subroutine xx2x(x,xx)
364 implicit real*8 (a-h,o-z)
367 include 'COMMON.CHAIN'
368 include 'COMMON.INTERACT'
369 double precision xx(maxvar),x(maxvar)
379 if (mask_phi(i).eq.1) then
387 if (mask_theta(i).eq.1) then
395 if (itype(i).ne.10) then
397 if (mask_side(i).eq.1) then
408 c----------------------------------------------------------
409 subroutine minim_dc(etot,iretcode,nfun)
410 implicit real*8 (a-h,o-z)
412 parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2))
416 include 'COMMON.SETUP'
417 include 'COMMON.IOUNITS'
420 include 'COMMON.MINIM'
421 include 'COMMON.CHAIN'
423 double precision minval,x(maxvar),d(maxvar),v(1:lv),xx(maxvar)
424 c common /przechowalnia/ v
426 double precision energia(0:n_ene)
427 external func_dc,grad_dc,fdum
428 logical not_done,change,reduce
429 double precision g(maxvar),f1
431 call deflt(2,iv,liv,lv,v)
432 * 12 means fresh start, dont call deflt
434 * max num of fun calls
435 if (maxfun.eq.0) maxfun=500
437 * max num of iterations
438 if (maxmin.eq.0) maxmin=1000
442 * selects output unit
444 if (print_min_ini+print_min_stat+print_min_res.gt.0) iv(21)=iout
445 * 1 means to print out result
447 * 1 means to print out summary stats
448 iv(23)=print_min_stat
449 * 1 means to print initial x and d
451 * min val for v(radfac) default is 0.1
453 * max val for v(radfac) default is 4.0
456 * check false conv if (act fnctn decrease) .lt. v(26)*(exp decrease)
457 * the sumsl default is 0.1
459 * false conv if (act fnctn decrease) .lt. v(34)
460 * the sumsl default is 100*machep
462 * absolute convergence
463 if (tolf.eq.0.0D0) tolf=1.0D-4
465 * relative convergence
466 if (rtolf.eq.0.0D0) rtolf=1.0D-4
468 * controls initial step size
470 * large vals of d correspond to small components of step
483 if (ialph(i,1).gt.0) then
490 print *,"check_ecart before sumsl"
492 c call exec_checkgrad
494 call sumsl(k,d,x,func_dc,grad_dc,iv,liv,lv,v,idum,rdum,fdum)
504 if (ialph(i,1).gt.0) then
515 cd call func_dc(k,x,nf,f,idum,rdum,fdum)
516 cd call grad_dc(k,x,nf,g,idum,rdum,fdum)
520 cd call func_dc(k,x,nf,f1,idum,rdum,fdum)
522 cd print '(i5,2f15.5)',i,g(i),(f1-f)/1.0D-5
528 print *,"check_ecart"
530 c call exec_checkgrad
534 subroutine func_dc(n,x,nf,f,uiparm,urparm,ufparm)
535 implicit real*8 (a-h,o-z)
540 include 'COMMON.SETUP'
541 include 'COMMON.DERIV'
542 include 'COMMON.IOUNITS'
544 include 'COMMON.CHAIN'
546 double precision energia(0:n_ene)
547 double precision ufparm
564 if (ialph(i,1).gt.0) then
574 call etotal(energia(0))
577 cd print *,'func_dc ',nf,nfl,f
582 subroutine grad_dc(n,x,nf,g,uiparm,urparm,ufparm)
583 implicit real*8 (a-h,o-z)
588 include 'COMMON.SETUP'
589 include 'COMMON.CHAIN'
590 include 'COMMON.DERIV'
592 include 'COMMON.INTERACT'
593 include 'COMMON.FFIELD'
595 include 'COMMON.IOUNITS'
598 double precision urparm(1)
599 dimension x(maxvar),g(maxvar)
605 cd print *,'grad_dc ',nf,nfl,nf-nfl+1,icg
606 if (nf-nfl+1) 20,30,40
607 20 call func_dc(n,x,nf,f,uiparm,urparm,ufparm)
621 if (ialph(i,1).gt.0) then
631 C Evaluate the derivatives of virtual bond lengths and SC vectors in variables.
643 if (ialph(i,1).gt.0) then