source: minim.f@ e40e335

Last change on this file since e40e335 was e40e335, checked in by baerbaer <baerbaer@…>, 16 years ago

Initial import to BerliOS corresponding to 3.0.4

git-svn-id: ​svn+ssh://svn.berlios.de/svnroot/repos/smmp/trunk@1 26dc1dd8-5c4e-0410-9ffe-d298b4865968
  • Property mode set to 100644
File size: 7.3 KB
Line 
1c **************************************************************
2c
3c This file contains the subroutines: minim,move
4c
5c Copyright 2003-2005 Frank Eisenmenger, U.H.E. Hansmann,
6c Shura Hayryan, Chin-Ku
7c Copyright 2007 Frank Eisenmenger, U.H.E. Hansmann,
8c Jan H. Meinke, Sandipan Mohanty
9c
10c **************************************************************
11
12
13 subroutine minim(imin, maxit, eps)
14
15c ......................................................................
16c PURPOSE: Use minimizers
17c
18c imin = 1: use Quasi-Newton
19c imin = 2: use Conjugated Gradients
20c
21c @param maxit maximum number of iterations
22c @param eps acceptance criterium
23cInstitute
24c CALLS: difang,energy,gradient, mincjg,minqsn, nursvr
25c ......................................................................
26
27 include 'INCL.H'
28cf2py intent(in) imin
29cf2py intent(in) maxit
30cf2py intent(in) eps
31 parameter (msvmx=mxvr*(mxvr+5)/(2*(2*mxvr+1)),
32 # msv = 50 )
33
34 dimension w(mxvr*(mxvr+13)/2)
35
36 dimension vlvrn(mxvr),vlvro(mxvr),gdvr(mxvr),scl(mxvr)
37
38c --------------------------- new
39 dimension gbpro(6,mxml)
40
41 mxop=maxit
42 acc=eps
43
44 if (msv.gt.msvmx) then
45 write (*,'(a,i5)') ' minreg> parameter MSV > ',msvmx
46 stop
47 endif
48
49c ----------------------- energy & gradient
50
51 call gradient()
52
53 write(*,'(/,a,/)') ' Energy BEFORE minimization:'
54
55 if (ireg.eq.0) then
56
57 write (*,'(a,e12.5,/,3(a,e11.4),/,2(a,e11.4),/)') ' Total: ',
58 # eysm,
59 # ' Coulomb: ',eyel,' Lennard-Jones: ',eyvw,' HB: ',eyhb,
60 # ' Variables: ',eyvr,' Solvatation: ',eysl
61
62 else
63
64
65 write (*,'(a,e12.5,/,3(a,e11.4),/,3(a,e11.4),/)') ' Total: ',
66 # wtey*eysm + wtrg*eyrg,
67 # ' Coulomb: ',eyel,' Lennard-Jones: ',eyvw,' HB: ',eyhb,
68 # ' Variables: ',eyvr,' Solvatation: ',eysl,
69 # ' Regularization: ',eyrg
70
71 endif
72
73c --------------------------------------- variables
74
75 ntlvr=ivrml1(ntlml)+nvrml(ntlml)-1
76 nv=0
77
78 gdey2=0.d0
79 gdrg2=0.d0
80
81 do i=1,ntlvr
82 if (.not.fxvr(i)) then
83
84 nv=nv+1
85
86 vlvrn(nv) = vlvr(i)
87
88 scl(nv) = 1.d0 ! scale grad_i
89
90 if (ireg.eq.0) then
91 gdvr(nv) = gdeyvr(i)
92 else
93 gdvr(nv) = wtey*gdeyvr(i)+wtrg*gdeyrg(i)
94 gdrg2 = gdrg2+gdeyrg(i)**2
95 endif
96
97 gdey2=gdey2+gdeyvr(i)**2
98
99 endif
100
101 vlvro(i)=vlvr(i)
102 enddo
103
104 n=nv
105
106 if (ireg.eq.0) then
107
108 esm=eysm
109
110 else
111
112 ngbvr=0
113
114 do i=1,ntlml
115 do j=1,6
116
117 n=n+1
118 ngbvr=ngbvr+1
119
120 vlvrn(n) = gbpr(j,i)
121 scl(n) = 1.d0 ! Scale
122
123 gbpro(j,i) = gbpr(j,i) ! save
124
125 gdvr(n) = wtrg*gdeygb(ngbvr)
126
127 gdrg2=gdrg2+gdeygb(ngbvr)**2
128 enddo
129 enddo
130
131c if (abs(wtrg-1.d0).gt.1.d-4.and.abs(wtey-1.d0).gt.1.d-4) then
132c gdey2 = max(acc,gdey2)
133c gdrg2 = max(acc,gdrg2)
134c wtrg = wtrg * sqrt(gdey2/gdrg2)
135c write(*,*) ' --> Wt_energy = ',wtey,' Wt_regul. = ',wtrg
136c write(*,*) ' '
137c endif
138
139 esm=wtey*eysm+wtrg*eyrg
140
141 endif
142
143
144 if (imin.eq.1) then ! Quasi-Newton
145
146 n1=1+(n*(n+1))/2
147 n2=n1+n
148 n3=n2+n
149 n4=n3+n
150 n5=n4+n
151 n6=n5+n
152
153 call minqsn(n,mxvr,vlvrn,esm,gdvr,scl,acc,w,w(n1),w(n2),
154 # w(n3),w(n4),w(n5),w(n6),mxop,nop)
155
156 elseif (imin.eq.2) then ! Conjugated Gradients
157
158 n1=1+n
159 n2=n1+n
160 n3=n2+n
161 n4=n3+n
162 n5=n4+n
163
164 call mincjg(n,mxvr,vlvrn,esm,gdvr,acc,w,w(n1),w(n2), ! no 'scl'
165 # w(n3),w(n4),w(n5),mxop,nop)
166
167 endif
168
169
170 if (nop.lt.mxop) then
171 write (*,'(a)') ' ---- CONVERGENCE ----'
172 else
173 write (*,'(a)') '---- STEP LIMIT ----'
174 endif
175
176
177 write (*,'(/,2a,/)') ' Final energies ',
178 # '__________________________________________________'
179
180 eysm = energy()
181
182 if (ireg.eq.0) then
183
184 write (*,'(a,e12.5,/,3(a,e11.4),/,2(a,e11.4))') ' Total: ',eysm,
185 # ' Coulomb: ',eyel,' Lennard-Jones: ',eyvw,' HB: ',eyhb,
186 # ' Variables: ',eyvr,' Solvatation: ',eysl
187
188 else
189
190 write (*,'(a,e12.5,/,3(a,e11.4),/,3(a,e11.4))') ' Total: ',
191 # wtey*eysm + wtrg*eyrg,
192 # ' Coulomb: ',eyel,' Lennard-Jones: ',eyvw,' HB: ',eyhb,
193 # ' Variables: ',eyvr,' Solvatation: ',eysl,
194 # ' Regularization: ',eyrg
195
196 endif
197
198 write (*,'(/,a,/)') ' Variables _________________'
199
200 nv = 0
201 do i=1,ntlvr !! nvr
202 if (.not.fxvr(i)) then
203
204 nv=nv+1
205 vr=vlvrn(nv)
206 if (abs(vr).gt.pi) vr=vr-sign(pi2,vr)
207
208 write (*,'(1x,a,1x,i4,f8.1,a,f5.1,a)') nmvr(i),nursvr(i),
209 # vr*crd,' (',abs(difang(vr,vlvro(i)))*crd,')'
210
211 vlvr(i) = vr
212 endif
213 enddo
214
215
216 if (ireg.ne.0) then
217
218 write (*,'(/,a,/)') ' Global Variables ___________'
219
220 do i=1,ntlml
221 write(*,*) ' Molecule #',i,' old new'
222 do j=1,3
223 write(*,*) gbpro(j,i),' ',gbpr(j,i)
224 enddo
225 do j=4,6
226 write(*,*) gbpro(j,i)*crd,' ',gbpr(j,i)*crd
227 enddo
228 enddo
229
230 endif
231
232 write (*,'(/,2a)') ' Gradient ',
233 # '______________________________________________________________'
234
235 write (*,'(8(1x,f8.3))') (gdvr(i),i=1,nv)
236
237 if (ireg.ne.0) then
238
239 write (*,*) ' -------------- global variables ------------'
240 write (*,'(6(1x,f8.3))') (gdvr(i+nv),i=1,ngbvr)
241
242 endif
243
244 return
245 end
246c ********************************************
247 subroutine move(nop,nvr1,esm,vlvrn,gdvr)
248c
249c CALLS: gradient
250c
251 include 'INCL.H'
252
253 dimension vlvrn(mxvr),gdvr(mxvr)
254
255
256c ------------------------ compile & new variables
257
258 ntlvr=ivrml1(ntlml)+nvrml(ntlml)-1
259 n=0
260
261 do i=1,ntlvr
262 if (.not.fxvr(i)) then
263 n=n+1
264 vlvr(i)=vlvrn(n)
265 endif
266 enddo
267
268 if (ireg.ne.0) then
269
270 ii=0
271 do i=1,ntlml
272 do j=1,6 ! global vars.
273 ii=ii+1
274 gbpr(j,i)=vlvrn(ii+n)
275 enddo
276 enddo
277
278 endif
279
280c -------------------------- new minimz. gradient
281
282 call gradient()
283
284 gdsmey=0.d0
285 gdsmrg=0.d0
286
287 n=0
288
289 do i=1,ntlvr
290 if (.not.fxvr(i)) then
291 n=n+1
292
293 if (ireg.eq.0) then
294 gdvr(n) = gdeyvr(i)
295 else
296 gdvr(n) = wtey*gdeyvr(i) + wtrg*gdeyrg(i)
297 gdsmrg = gdsmrg + gdeyrg(i)**2
298 endif
299
300 gdsmey = gdsmey + gdeyvr(i)**2
301
302 endif
303 enddo
304
305
306 if (ireg.eq.0) then
307
308 esm=eysm
309
310 write (*,'(a,i5,a,2(e13.6,a))') ' Step ',nop,': energy ',esm
311 # ,' (',gdsmey,' )'
312
313 else
314
315 ii=0
316 do i=1,ntlml ! global vars.
317 do j=1,6
318 n=n+1
319 ii=ii+1
320
321 gdvr(n) = wtrg*gdeygb(ii)
322 gdsmrg = gdsmrg+gdeygb(ii)**2
323
324 enddo
325 enddo
326
327 esm=wtey*eysm+wtrg*eyrg
328
329 write (*,'(a,i5,a,3(e13.6,a))') ' Step ',nop,': energy ',esm
330 # ,' (',gdsmey,',',gdsmrg,' )'
331
332 endif
333
334 return
335 end
336
Note: See TracBrowser for help on using the repository browser.