1 |
heimbach |
1.2 |
|
2 |
|
|
subroutine lsopt_top( nn, xx, ff, gg, simul, optline |
3 |
|
|
$ , epsx, fmin, epsg |
4 |
|
|
$ , iprint, itmax, nfunc, nupdate |
5 |
|
|
$ , dd, gold, xdiff |
6 |
|
|
$ , loffline |
7 |
|
|
$ , ifail ) |
8 |
|
|
|
9 |
|
|
c ================================================================== |
10 |
|
|
c SUBROUTINE lsopt_top |
11 |
|
|
c ================================================================== |
12 |
|
|
c |
13 |
|
|
c o uses a set of control variables, their adjoint variables, |
14 |
|
|
c and a cost function value |
15 |
|
|
c to compute an improved set of controls with respect to the |
16 |
|
|
c cost function via a |
17 |
|
|
c variable-storage Quasi-Newton method |
18 |
|
|
c |
19 |
|
|
c o Reference: J.C. Gilbert & C. Lemarechal |
20 |
|
|
c Some numerical experiments with variable-storage |
21 |
|
|
c quasi-Newton algorithms |
22 |
|
|
c Mathematical Programming 45 (1989), pp. 407-435 |
23 |
|
|
c |
24 |
|
|
c o started: ??? not reproducible |
25 |
|
|
c |
26 |
|
|
c o changed: Patrick Heimbach, MIT/EAPS |
27 |
|
|
c |
28 |
|
|
c o Version: 2.0, 24-Feb-2000: |
29 |
|
|
c (Version 1.0 is considered as version |
30 |
|
|
c starting from which changes were made). |
31 |
|
|
c - severe changes in structure including various |
32 |
|
|
c bug-fixes to enable multi-optimization runs; |
33 |
|
|
c - routine lsoptw incorporated into lsoptv |
34 |
|
|
c - optimization iteration loop restructured |
35 |
|
|
c - complete restructuring of handling |
36 |
|
|
c cold start cases |
37 |
|
|
c - number of 3 control flags for error handling |
38 |
|
|
c (indic, moderl, ifail) reduced to 1 (ifail) |
39 |
|
|
c and homogenized with routine lsline |
40 |
|
|
c |
41 |
|
|
c o Version: 2.1, 29-Feb-2000: |
42 |
|
|
c - handling of case ifail = 6 changed; |
43 |
|
|
c leads to stop of optimization |
44 |
|
|
c (similar to case ifail = 8) |
45 |
|
|
c - logical lphprint included |
46 |
|
|
c |
47 |
|
|
c ================================================================== |
48 |
|
|
c SUBROUTINE lsopt_top |
49 |
|
|
c ================================================================== |
50 |
|
|
|
51 |
|
|
implicit none |
52 |
|
|
|
53 |
|
|
#include <blas1.h> |
54 |
|
|
|
55 |
|
|
c----------------------------------------- |
56 |
|
|
c declare arguments |
57 |
|
|
c----------------------------------------- |
58 |
|
|
integer nn, iprint, itmax, nfunc, nupdate, ifail |
59 |
|
|
|
60 |
|
|
double precision xx(nn), ff, gg(nn), epsx, fmin, epsg |
61 |
|
|
double precision dd(nn), gold(nn), xdiff(nn) |
62 |
|
|
|
63 |
|
|
cph( |
64 |
|
|
integer phniter0, phniterold |
65 |
|
|
double precision phff |
66 |
|
|
COMMON /PH_OPTI/ phniter0, phniterold, phff |
67 |
|
|
cph) |
68 |
|
|
|
69 |
|
|
external simul, optline |
70 |
|
|
|
71 |
|
|
c----------------------------------------- |
72 |
|
|
C declare local variables |
73 |
|
|
c----------------------------------------- |
74 |
|
|
logical cold, lphprint, loffline |
75 |
|
|
parameter (lphprint = .true.) |
76 |
|
|
|
77 |
|
|
integer mm, mupd, jmin, jmax, indic, isize, REAL_BYTE |
78 |
|
|
integer i, iiter, ifunc |
79 |
|
|
|
80 |
|
|
double precision fupd |
81 |
|
|
double precision r1, tmin, tmax, tact, gnorm, gnorm0, eps1 |
82 |
|
|
double precision fold, ys |
83 |
|
|
double precision dotdg |
84 |
|
|
|
85 |
|
|
external SDOT, SNRM2, SSCAL |
86 |
|
|
double precision SDOT, SNRM2 |
87 |
|
|
|
88 |
|
|
c----------------------------------------- |
89 |
|
|
c parameters |
90 |
|
|
c----------------------------------------- |
91 |
|
|
|
92 |
|
|
double precision rmin |
93 |
|
|
parameter( rmin = 1.e-20 ) |
94 |
|
|
|
95 |
|
|
character*(*) iform |
96 |
|
|
parameter(iform='(i5,2x,1pe8.1,1x,i5,4x,1pe11.4,3(2x,1pe8.1))' ) |
97 |
|
|
|
98 |
|
|
c ================================================================== |
99 |
|
|
c |
100 |
|
|
c----------------------------------------- |
101 |
|
|
c initialization |
102 |
|
|
c----------------------------------------- |
103 |
|
|
cold = .true. |
104 |
|
|
fupd = 1.0e+10 |
105 |
|
|
indic = 0 |
106 |
|
|
tmin = 0. |
107 |
|
|
tmax = 1.0e+10 |
108 |
|
|
tact = 1.0 |
109 |
|
|
cph( |
110 |
|
|
phniterold = 0 |
111 |
|
|
cph) |
112 |
|
|
iiter = 0 |
113 |
|
|
eps1 = 1.0 |
114 |
|
|
ifunc = 0 |
115 |
|
|
ifail = 0 |
116 |
|
|
gnorm0 = 1. |
117 |
|
|
|
118 |
|
|
c----------------------------------------- |
119 |
|
|
c initialization for dd and dds |
120 |
|
|
c----------------------------------------- |
121 |
|
|
|
122 |
|
|
jmin = 1 |
123 |
|
|
jmax = 0 |
124 |
|
|
|
125 |
|
|
mm = nn |
126 |
|
|
mupd = nupdate |
127 |
|
|
|
128 |
|
|
REAL_BYTE = 8 |
129 |
|
|
isize = REAL_BYTE |
130 |
|
|
|
131 |
|
|
c----------------------------------------- |
132 |
|
|
c print information |
133 |
|
|
c----------------------------------------- |
134 |
|
|
if (iprint .ge. 1) then |
135 |
|
|
print '(2x,a)', |
136 |
|
|
$ '===============================================' |
137 |
|
|
print '(2x,a)', |
138 |
|
|
$ '=== O P T I M I Z A T I O N ===' |
139 |
|
|
print '(2x,a)', |
140 |
|
|
$ '===============================================' |
141 |
|
|
print '(a,i9)' |
142 |
|
|
$ , ' number of control variables.......', nn |
143 |
|
|
print '(a,e9.2)' |
144 |
|
|
$ , ' precision on x....................', epsx |
145 |
|
|
print '(a,e9.2)' |
146 |
|
|
$ , ' precision on g....................', epsg |
147 |
|
|
print '(a,e9.2)' |
148 |
|
|
$ , ' expected optimal function value...', fmin |
149 |
|
|
print '(a,i9)' |
150 |
|
|
$ , ' maximal number of iterations......', itmax |
151 |
|
|
print '(a,i9)' |
152 |
|
|
$ , ' maximal number of simulations.....', nfunc |
153 |
|
|
print '(a,i9)' |
154 |
|
|
$ , ' information level.................', iprint |
155 |
|
|
print '(a,i9)' |
156 |
|
|
$ , ' number of updates.................', nupdate |
157 |
|
|
print '(a,i9)' |
158 |
|
|
$ , ' size of used memory...............', 3*nn |
159 |
|
|
endif |
160 |
|
|
|
161 |
|
|
c----------------------------------------- |
162 |
|
|
c check arguments |
163 |
|
|
c----------------------------------------- |
164 |
|
|
|
165 |
|
|
if (nn .le. 0) then |
166 |
|
|
if (iprint.ge.1) then |
167 |
|
|
print '(a,i6)' , ' ERROR : n = ', nn |
168 |
|
|
endif |
169 |
|
|
ifail = 1 |
170 |
|
|
goto 999 |
171 |
|
|
endif |
172 |
|
|
|
173 |
|
|
if (itmax .lt. 0) then |
174 |
|
|
if (iprint.ge.1) then |
175 |
|
|
print '(a,i6)' , ' ERROR : itmax = ', itmax |
176 |
|
|
endif |
177 |
|
|
ifail = 1 |
178 |
|
|
goto 999 |
179 |
|
|
endif |
180 |
|
|
|
181 |
|
|
if (nfunc .le. 0) then |
182 |
|
|
if (iprint.ge.10) then |
183 |
|
|
print '(a,i6)' , ' ERROR : nfunc = ', nfunc |
184 |
|
|
endif |
185 |
|
|
ifail = 1 |
186 |
|
|
goto 999 |
187 |
|
|
endif |
188 |
|
|
|
189 |
|
|
if (epsx .le. 0.) then |
190 |
|
|
if (iprint.ge.1) then |
191 |
|
|
print '(a,e9.2)', ' ERROR : epsx = ', epsx |
192 |
|
|
endif |
193 |
|
|
ifail = 1 |
194 |
|
|
goto 999 |
195 |
|
|
endif |
196 |
|
|
|
197 |
|
|
if (epsg .le. 0.) then |
198 |
|
|
if (iprint.ge.1) then |
199 |
|
|
print '(a,e9.2)', ' ERROR : epsg = ', epsg |
200 |
|
|
endif |
201 |
|
|
ifail = 1 |
202 |
|
|
goto 999 |
203 |
|
|
endif |
204 |
|
|
|
205 |
|
|
if (epsg .gt. 1.) then |
206 |
|
|
if (iprint.ge.1) then |
207 |
|
|
print '(a,e9.2)', ' ERROR : epsg = ', epsg |
208 |
|
|
endif |
209 |
|
|
ifail = 1 |
210 |
|
|
goto 999 |
211 |
|
|
endif |
212 |
|
|
|
213 |
|
|
cph( |
214 |
|
|
print *, 'pathei: vor instore ' |
215 |
|
|
cph) |
216 |
|
|
call instore( mm, fupd, gnorm0, isize, mupd, jmin, jmax, cold, |
217 |
|
|
& ifail ) |
218 |
|
|
|
219 |
|
|
cph( |
220 |
|
|
phff = fupd |
221 |
|
|
cph) |
222 |
|
|
|
223 |
|
|
c----------------------------------------- |
224 |
|
|
c check warm start parameters |
225 |
|
|
c----------------------------------------- |
226 |
|
|
if (ifail .ne. 0) then |
227 |
|
|
if (iprint.ge.1) then |
228 |
|
|
print '(a)', ' ERROR : reading restart file' |
229 |
|
|
end if |
230 |
|
|
ifail = 2 |
231 |
|
|
goto 999 |
232 |
|
|
end if |
233 |
|
|
|
234 |
|
|
if (isize .ne. REAL_BYTE) then |
235 |
|
|
if (iprint.ge.1) then |
236 |
|
|
print '(a)', ' ERROR : uncompatible floating point format' |
237 |
|
|
end if |
238 |
|
|
ifail = 2 |
239 |
|
|
goto 999 |
240 |
|
|
end if |
241 |
|
|
|
242 |
|
|
if (mupd .lt. 1) then |
243 |
|
|
if (iprint .ge. 1) then |
244 |
|
|
print '(a)', ' ERROR : m is set too small in instore' |
245 |
|
|
endif |
246 |
|
|
ifail = 2 |
247 |
|
|
goto 999 |
248 |
|
|
endif |
249 |
|
|
|
250 |
|
|
c----------------------------------------- |
251 |
|
|
c cold start or warm restart ? |
252 |
|
|
c----------------------------------------- |
253 |
|
|
|
254 |
|
|
if (cold) then |
255 |
|
|
c--- start if cold start --- |
256 |
|
|
if (lphprint) then |
257 |
|
|
print '(a)', 'pathei-lsopt: cold start' |
258 |
|
|
print * |
259 |
|
|
end if |
260 |
|
|
|
261 |
|
|
call simul( indic, nn, xx, ff, gg ) |
262 |
|
|
cph( |
263 |
|
|
print *, 'pathei: nach simul: nn, ff = ', nn, ff |
264 |
|
|
print *, 'pathei: nach simul: xx, gg = ', xx(1), gg(1) |
265 |
|
|
cph) |
266 |
|
|
|
267 |
|
|
do i = 1, nn |
268 |
|
|
xdiff(i) = 1. |
269 |
|
|
end do |
270 |
|
|
|
271 |
|
|
cph( |
272 |
|
|
print *, 'pathei: vor dostore ' |
273 |
|
|
cph) |
274 |
|
|
call dostore( nn, xx, .true., 1 ) |
275 |
|
|
call dostore( nn, gg, .true., 2 ) |
276 |
|
|
call dostore( nn, xdiff, .true., 3 ) |
277 |
|
|
|
278 |
|
|
cph( |
279 |
|
|
print *, 'pathei: vor lswri ' |
280 |
|
|
cph) |
281 |
|
|
call lswri( iiter, nn, xx, gg, lphprint ) |
282 |
|
|
|
283 |
|
|
cph( |
284 |
|
|
print *, 'pathei: vor gnorm0 ' |
285 |
|
|
cph) |
286 |
|
|
gnorm0 = SNRM2( nn, gg, 1 ) |
287 |
|
|
cph( |
288 |
|
|
print *, 'pathei: gnorm0 = ', gnorm0 |
289 |
|
|
cph) |
290 |
|
|
if (gnorm0 .lt. rmin) then |
291 |
|
|
ifail = 3 |
292 |
|
|
goto 1000 |
293 |
|
|
endif |
294 |
|
|
|
295 |
|
|
cph( |
296 |
|
|
phff = ff |
297 |
|
|
cph) |
298 |
|
|
|
299 |
|
|
if (lphprint) |
300 |
|
|
& print *, 'pathei-lsopt: cold; first call simul: ff = ', |
301 |
|
|
& phff |
302 |
|
|
|
303 |
|
|
c--- end if cold start --- |
304 |
|
|
else |
305 |
|
|
c--- start if warm start --- |
306 |
|
|
if (mm .ne. nn) then |
307 |
|
|
if (iprint.ge.1) then |
308 |
|
|
print '(a,i6)' |
309 |
|
|
$ , ' ERROR : inconsistent nn = ', mm |
310 |
|
|
endif |
311 |
|
|
ifail = 2 |
312 |
|
|
goto 999 |
313 |
|
|
endif |
314 |
|
|
if (mupd .ne. nupdate) then |
315 |
|
|
if (iprint.ge.1) then |
316 |
|
|
print '(a,i6)' |
317 |
|
|
$ , ' ERROR : inconsistent nupdate = ', mupd |
318 |
|
|
endif |
319 |
|
|
ifail = 2 |
320 |
|
|
goto 999 |
321 |
|
|
endif |
322 |
|
|
if (jmin .lt. 1) then |
323 |
|
|
if (iprint.ge.1) then |
324 |
|
|
print '(a,i6)' |
325 |
|
|
$ , ' ERROR : inconsistent jmin = ', jmin |
326 |
|
|
endif |
327 |
|
|
ifail = 2 |
328 |
|
|
goto 999 |
329 |
|
|
endif |
330 |
|
|
if (jmin .gt. mupd) then |
331 |
|
|
if (iprint.ge.1) then |
332 |
|
|
print '(a,i6)' |
333 |
|
|
$ , ' ERROR : inconsistent jmin = ', jmin |
334 |
|
|
endif |
335 |
|
|
ifail = 2 |
336 |
|
|
goto 999 |
337 |
|
|
endif |
338 |
|
|
if (jmax .gt. mupd) then |
339 |
|
|
if (iprint.ge.1) then |
340 |
|
|
print '(a,i6)' |
341 |
|
|
$ , ' ERROR : inconsistent jmax = ', jmax |
342 |
|
|
endif |
343 |
|
|
ifail = 2 |
344 |
|
|
goto 999 |
345 |
|
|
endif |
346 |
|
|
|
347 |
|
|
if (lphprint) then |
348 |
|
|
print *, 'pathei-lsopt: warm start; read via dostore' |
349 |
|
|
print * |
350 |
|
|
endif |
351 |
|
|
|
352 |
|
|
call dostore( nn, xx, .false., 1 ) |
353 |
|
|
call dostore( nn, gg, .false., 2 ) |
354 |
|
|
ff = fupd |
355 |
|
|
|
356 |
|
|
cph( |
357 |
|
|
phff = ff |
358 |
|
|
cph) |
359 |
|
|
|
360 |
|
|
if (lphprint) |
361 |
|
|
& print *, 'pathei-lsopt: warm; first dostore read: ff = ', |
362 |
|
|
& ff |
363 |
|
|
|
364 |
|
|
c--- end if warm start --- |
365 |
|
|
endif |
366 |
|
|
|
367 |
|
|
if (iprint .ge. 1) then |
368 |
|
|
print '(2a)', ' Itn Step Nfun Objective ' |
369 |
|
|
$ , 'Norm G Norm X Norm (X(k-1)-X(k))' |
370 |
|
|
end if |
371 |
|
|
|
372 |
|
|
c----------------------------------------- |
373 |
|
|
c print information line |
374 |
|
|
c----------------------------------------- |
375 |
|
|
if (cold) then |
376 |
|
|
print iform, iiter, tact, ifunc, ff, gnorm0 |
377 |
|
|
$ , SNRM2( nn, xx, 1 ), 0. |
378 |
|
|
|
379 |
|
|
write(94,'(i5,2x,1pe11.4,4(2x,1pe8.1))') |
380 |
|
|
& iiter, ff, gnorm0, tact, |
381 |
|
|
& SNRM2( nn, xx, 1 ), 0. |
382 |
|
|
|
383 |
|
|
if ( itmax .EQ. 0 ) then |
384 |
|
|
ifail = 10 |
385 |
|
|
goto 1000 |
386 |
|
|
end if |
387 |
|
|
end if |
388 |
|
|
|
389 |
|
|
c======================================================================= |
390 |
|
|
c begin of loop |
391 |
|
|
c compute x(k+1) out of x(k) + t*d, where t > 0. |
392 |
|
|
c======================================================================= |
393 |
|
|
|
394 |
|
|
do iiter = 1, itmax |
395 |
|
|
|
396 |
|
|
if (lphprint) then |
397 |
|
|
print *, 'pathei-lsopt: ++++++++++++++++++++++++' |
398 |
|
|
print *, 'pathei-lsopt: entering iiter =', iiter |
399 |
|
|
end if |
400 |
|
|
|
401 |
|
|
c----------------------------------------- |
402 |
|
|
c store old values |
403 |
|
|
c----------------------------------------- |
404 |
|
|
do i = 1, nn |
405 |
|
|
gold(i) = gg(i) |
406 |
|
|
end do |
407 |
|
|
|
408 |
|
|
fold = ff |
409 |
|
|
cph( |
410 |
|
|
phniter0 = iiter |
411 |
|
|
phff = ff |
412 |
|
|
cph) |
413 |
|
|
|
414 |
|
|
c----------------------------------------- |
415 |
|
|
c compute new dd and xx |
416 |
|
|
c----------------------------------------- |
417 |
|
|
c |
418 |
|
|
call lsupdxx( |
419 |
|
|
& nn, ifail, lphprint |
420 |
|
|
& , jmin, jmax, nupdate |
421 |
|
|
& , ff, fmin, fold, gnorm0, dotdg |
422 |
|
|
& , gg, dd, xx, xdiff |
423 |
|
|
& , tmin, tmax, tact, epsx |
424 |
|
|
& ) |
425 |
|
|
|
426 |
|
|
c----------------------------------------- |
427 |
|
|
c check whether new direction is a descent one |
428 |
|
|
c----------------------------------------- |
429 |
|
|
if ( ifail .eq. 4) goto 1000 |
430 |
|
|
|
431 |
|
|
c----------------------------------------- |
432 |
|
|
c optline returns new values of x, f, and g |
433 |
|
|
c----------------------------------------- |
434 |
|
|
c |
435 |
|
|
call optline( |
436 |
|
|
& simul |
437 |
|
|
& , nn, ifail, lphprint |
438 |
|
|
& , ifunc, nfunc |
439 |
|
|
& , ff, dotdg |
440 |
|
|
& , tmin, tmax, tact, epsx |
441 |
|
|
& , dd, gg, xx, xdiff |
442 |
|
|
& ) |
443 |
|
|
|
444 |
|
|
if (lphprint) print *, 'pathei-lsopt: ', |
445 |
|
|
& 'back from optline; ifail = ', ifail |
446 |
|
|
if (lphprint) print *, 'pathei-lsopt: ', |
447 |
|
|
& 'dostore 1,2 at iiter ', iiter |
448 |
|
|
|
449 |
|
|
call dostore( nn, xx, .true., 1 ) |
450 |
|
|
call dostore( nn, gg, .true., 2 ) |
451 |
|
|
cph( |
452 |
|
|
call lswri( iiter, nn, xx, gg, lphprint ) |
453 |
|
|
cph) |
454 |
|
|
|
455 |
|
|
gnorm = SNRM2( nn, gg, 1 ) |
456 |
|
|
|
457 |
|
|
c----------------------------------------- |
458 |
|
|
c print information line |
459 |
|
|
c----------------------------------------- |
460 |
|
|
print iform, iiter, tact, ifunc, ff, gnorm |
461 |
|
|
$ , SNRM2( nn, xx, 1 ), tact*SNRM2( nn, dd, 1 ) |
462 |
|
|
|
463 |
|
|
write(94,'(i5,2x,1pe11.4,4(2x,1pe8.1))') |
464 |
|
|
& iiter, ff, gnorm, tact, |
465 |
|
|
& SNRM2( nn, xx, 1 ), tact*SNRM2( nn, dd, 1 ) |
466 |
|
|
|
467 |
|
|
c----------------------------------------- |
468 |
|
|
c check output mode of ifail |
469 |
|
|
c----------------------------------------- |
470 |
|
|
if (ifail .eq. 7 |
471 |
|
|
& .or. ifail .eq. 8 |
472 |
|
|
& .or. ifail .eq. 9) goto 1000 |
473 |
|
|
|
474 |
|
|
c----------------------------------------- |
475 |
|
|
c maximal number of simulation reached |
476 |
|
|
c no goto in order to update Hessian |
477 |
|
|
c----------------------------------------- |
478 |
|
|
if (ifail .eq. 6) ifail = 0 |
479 |
|
|
|
480 |
|
|
c----------------------------------------- |
481 |
|
|
c NOTE: stopping tests are now done |
482 |
|
|
c after having updated the matrix, |
483 |
|
|
c so that update information can be stored |
484 |
|
|
c in case of a later warm restart |
485 |
|
|
c----------------------------------------- |
486 |
|
|
c----------------------------------------- |
487 |
|
|
c compute new s, y |
488 |
|
|
c----------------------------------------- |
489 |
|
|
do i = 1, nn |
490 |
|
|
xdiff(i) = tact*dd(i) |
491 |
|
|
end do |
492 |
|
|
|
493 |
|
|
c----------------------------------------- |
494 |
|
|
c compute ys |
495 |
|
|
c----------------------------------------- |
496 |
|
|
do i = 1, nn |
497 |
|
|
gold(i) = gg(i)-gold(i) |
498 |
|
|
end do |
499 |
|
|
|
500 |
|
|
ys = SDOT( nn, gold, 1, xdiff, 1 ) |
501 |
|
|
if (ys .le. 0.) then |
502 |
|
|
ifail = 4 |
503 |
|
|
print *, 'pathei-lsopt: ys < 0; ifail = ', ifail |
504 |
|
|
goto 1000 |
505 |
|
|
endif |
506 |
|
|
cph( |
507 |
|
|
cph----------------------------------------- |
508 |
|
|
cph at this point it is clear that xdiff |
509 |
|
|
cph provides a true optimized solution; |
510 |
|
|
cph i.e. take new gradient gg to compute new dd |
511 |
|
|
cph----------------------------------------- |
512 |
|
|
cph) |
513 |
|
|
|
514 |
|
|
c----------------------------------------- |
515 |
|
|
c update pointers for hessupd |
516 |
|
|
c----------------------------------------- |
517 |
|
|
if (nupdate .gt. 0) then |
518 |
|
|
|
519 |
|
|
if (jmax .eq. 0) then |
520 |
|
|
jmax = jmax+1 |
521 |
|
|
if (lphprint) |
522 |
|
|
& print *, 'pathei-lsopt: ', |
523 |
|
|
& 'first pointer update after cold start; ', |
524 |
|
|
& 'iiter, jmin, jmax = ', iiter, jmin, jmax |
525 |
|
|
else |
526 |
|
|
jmax = jmax+1 |
527 |
|
|
if (jmax .gt. nupdate) jmax = jmax-nupdate |
528 |
|
|
|
529 |
|
|
if (jmin .eq. jmax) then |
530 |
|
|
if (lphprint) |
531 |
|
|
& print *, 'pathei-lsopt: pointers updated for ', |
532 |
|
|
& ' iiter, jmin, jmax = ', iiter, jmin, jmax |
533 |
|
|
jmin = jmin+1 |
534 |
|
|
if (jmin .gt. nupdate) jmin = jmin-nupdate |
535 |
|
|
end if |
536 |
|
|
end if |
537 |
|
|
|
538 |
|
|
c----------------------------------------- |
539 |
|
|
c compute sbar, ybar store rec = min 4,5 |
540 |
|
|
c----------------------------------------- |
541 |
|
|
r1 = sqrt(1./ys) |
542 |
|
|
call SSCAL( nn, r1, xdiff, 1 ) |
543 |
|
|
call SSCAL( nn, r1, gold, 1 ) |
544 |
|
|
|
545 |
|
|
if (lphprint) |
546 |
|
|
& print *, 'pathei-lsopt: dostore at iiter, jmin, jmax ', |
547 |
|
|
& iiter, jmin, jmax |
548 |
|
|
|
549 |
|
|
call dostore( nn, gold, .true., 2*jmax+2 ) |
550 |
|
|
call dostore( nn, xdiff, .true., 2*jmax+3 ) |
551 |
|
|
|
552 |
|
|
c----------------------------------------- |
553 |
|
|
c compute the diagonal preconditioner |
554 |
|
|
c use dd as temporary array |
555 |
|
|
c----------------------------------------- |
556 |
|
|
call dgscale( nn, gold, xdiff, dd, rmin ) |
557 |
|
|
|
558 |
|
|
endif |
559 |
|
|
|
560 |
|
|
c----------------------------------------- |
561 |
|
|
c test convergence and stopping criteria |
562 |
|
|
c----------------------------------------- |
563 |
|
|
eps1 = gnorm / gnorm0 |
564 |
|
|
|
565 |
|
|
if (eps1 .lt. epsg) then |
566 |
|
|
|
567 |
|
|
ifail = 11 |
568 |
|
|
goto 1000 |
569 |
|
|
endif |
570 |
|
|
|
571 |
|
|
c======================================================================= |
572 |
|
|
c return of loop |
573 |
|
|
c======================================================================= |
574 |
|
|
|
575 |
|
|
end do |
576 |
|
|
|
577 |
|
|
iiter = iiter - 1 |
578 |
|
|
ifail = 5 |
579 |
|
|
|
580 |
|
|
c----------------------------------------- |
581 |
|
|
c loop exit |
582 |
|
|
c----------------------------------------- |
583 |
|
|
1000 continue |
584 |
|
|
nfunc = ifunc |
585 |
|
|
epsg = eps1 |
586 |
|
|
|
587 |
|
|
c----------------------------------------------------------------------- |
588 |
|
|
c save data for warm start |
589 |
|
|
c----------------------------------------------------------------------- |
590 |
|
|
|
591 |
|
|
call outstore( nn, ff, gnorm0, nupdate, jmin, jmax ) |
592 |
|
|
|
593 |
|
|
c----------------------------------------------------------------------- |
594 |
|
|
c compute dd(i+1), xx(i+1) based on latest available gg(i), xx(i) |
595 |
|
|
c for offline version |
596 |
|
|
c----------------------------------------------------------------------- |
597 |
|
|
|
598 |
|
|
if (loffline) then |
599 |
|
|
|
600 |
|
|
call lsupdxx( |
601 |
|
|
& nn, ifail, lphprint |
602 |
|
|
& , jmin, jmax, nupdate |
603 |
|
|
& , ff, fmin, fold, gnorm0, dotdg |
604 |
|
|
& , gg, dd, xx, xdiff |
605 |
|
|
& , tmin, tmax, tact, epsx |
606 |
|
|
& ) |
607 |
|
|
|
608 |
|
|
c Save xx(i+1) to file for offline version. |
609 |
|
|
call optim_write_control( nn, xdiff ) |
610 |
|
|
|
611 |
|
|
end if |
612 |
|
|
|
613 |
|
|
c----------------------------------------- |
614 |
|
|
c print final information |
615 |
|
|
c----------------------------------------- |
616 |
|
|
if (iprint .ge. 5) then |
617 |
|
|
print * |
618 |
|
|
print '(a,i9)' |
619 |
|
|
$ , ' number of iterations..............', iiter |
620 |
|
|
print '(a,i9)' |
621 |
|
|
$ , ' number of simultations............', nfunc |
622 |
|
|
print '(a,e9.2)' |
623 |
|
|
$ , ' relative precision on g...........', epsg |
624 |
|
|
end if |
625 |
|
|
|
626 |
|
|
if (iprint.ge.10) then |
627 |
|
|
print * |
628 |
|
|
print '(a,e15.8)' |
629 |
|
|
$ , ' cost function...............', ff |
630 |
|
|
print '(a,e15.8)' |
631 |
|
|
$ , ' norm of x...................', SNRM2( nn, xx, 1 ) |
632 |
|
|
print '(a,e15.8)' |
633 |
|
|
$ , ' norm of g...................', SNRM2( nn, gg, 1 ) |
634 |
|
|
end if |
635 |
|
|
|
636 |
|
|
c----------------------------------------- |
637 |
|
|
c print error message |
638 |
|
|
c----------------------------------------- |
639 |
|
|
999 continue |
640 |
|
|
|
641 |
|
|
if (ifail .ne. 0) then |
642 |
|
|
if (iprint .ge. 5) then |
643 |
|
|
print * |
644 |
|
|
print '(a)', ' optimization stopped because :' |
645 |
|
|
if (ifail .lt. 0) then |
646 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
647 |
|
|
$ , ' user request during simul' |
648 |
|
|
else if (ifail .eq. 0) then |
649 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
650 |
|
|
$ , ' optimal solution found' |
651 |
|
|
else if (ifail .eq. 1) then |
652 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
653 |
|
|
$ , ' an input argument is wrong' |
654 |
|
|
else if (ifail .eq. 2) then |
655 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
656 |
|
|
$ , ' warm start file is corrupted' |
657 |
|
|
else if (ifail .eq. 3) then |
658 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
659 |
|
|
$ , ' the initial gradient is too small' |
660 |
|
|
else if (ifail .eq. 4) then |
661 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
662 |
|
|
$ , ' the search direction is not a descent one' |
663 |
|
|
else if (ifail .eq. 5) then |
664 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
665 |
|
|
$ , ' maximal number of iterations reached' |
666 |
|
|
else if (ifail .eq. 6) then |
667 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
668 |
|
|
$ , ' maximal number of simulations reached' |
669 |
|
|
else if (ifail .eq. 7) then |
670 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
671 |
|
|
$ , ' the linesearch failed' |
672 |
|
|
else if (ifail .eq. 8) then |
673 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
674 |
|
|
$ , ' the function could not be improved' |
675 |
|
|
else if (ifail .eq. 9) then |
676 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
677 |
|
|
$ , ' optline parameters wrong' |
678 |
|
|
else if (ifail .eq. 10) then |
679 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
680 |
|
|
$ , ' cold start, no optimization done' |
681 |
|
|
else if (ifail .eq. 11) then |
682 |
|
|
print '(2x,a8,I3,a)', 'ifail = ', ifail |
683 |
|
|
$ , ' convergence achieved within precision' |
684 |
|
|
end if |
685 |
|
|
print * |
686 |
|
|
else if (iprint .ge. 1) then |
687 |
|
|
print * |
688 |
|
|
print '(a,i9)' |
689 |
|
|
$ , ' after optimization ifail..........', ifail |
690 |
|
|
print * |
691 |
|
|
end if |
692 |
|
|
end if |
693 |
|
|
|
694 |
|
|
c----------------------------------------- |
695 |
|
|
c end of subroutine |
696 |
|
|
c----------------------------------------- |
697 |
|
|
return |
698 |
|
|
|
699 |
|
|
end |