1 |
|
2 |
subroutine hessupd( nn, mupd, dd, jmin, jmax, xdiff, lphprint ) |
3 |
|
4 |
c ================================================================== |
5 |
c SUBROUTINE hessupd |
6 |
c ================================================================== |
7 |
c |
8 |
c o controls update of descent vector using available |
9 |
c approximation of Hessian Matrix based on gradients of |
10 |
c previous iterations |
11 |
c |
12 |
c o Reference: J.C. Gilbert & C. Lemarechal |
13 |
c Some numerical experiments with variable-storage |
14 |
c quasi-Newton algorithms |
15 |
c Mathematical Programming 45 (1989), pp. 407-435 |
16 |
c |
17 |
c o started: ??? not reproducible |
18 |
c |
19 |
c o changed: Patrick Heimbach, MIT/EAPS |
20 |
c 24-Feb-2000: |
21 |
c - changed some variable names to be consistent |
22 |
c with routine lsoptv, lsline; |
23 |
c |
24 |
c o Version: 2.1.0, 02-Mar-2000: Patrick Heimbach, MIT/EAPS |
25 |
c |
26 |
c ================================================================== |
27 |
c SUBROUTINE hessupd |
28 |
c ================================================================== |
29 |
|
30 |
implicit none |
31 |
|
32 |
#include "blas1.h" |
33 |
include 'mpif.h' |
34 |
c------------------------------------ |
35 |
c declare arguments |
36 |
c------------------------------------ |
37 |
integer nn, mupd, jmin, jmax |
38 |
double precision dd(nn), alpha(100), xdiff(nn) |
39 |
logical lphprint |
40 |
integer status(MPI_STATUS_SIZE),ierr |
41 |
integer myid, nprocs,mystart,myend |
42 |
common /mpi_parm/ myid,nprocs,mystart,myend |
43 |
|
44 |
c------------------------------------ |
45 |
c declare local variables |
46 |
c------------------------------------ |
47 |
external DDOT |
48 |
double precision DDOT |
49 |
|
50 |
integer jfin, i, j, jp |
51 |
double precision r |
52 |
double precision temp |
53 |
|
54 |
c------------------------------------ |
55 |
c initialization |
56 |
c------------------------------------ |
57 |
jfin = jmax |
58 |
|
59 |
if (lphprint .and. myid .eq. 0) |
60 |
& print *, 'pathei-lsopt: in hessupd; ', |
61 |
& 'jmin, jmax, mupd:', jmin, jmax, mupd |
62 |
|
63 |
|
64 |
if (jfin.lt.jmin) jfin = jmax+mupd |
65 |
|
66 |
c------------------------------------ |
67 |
c compute right hand side |
68 |
c------------------------------------ |
69 |
do j = jfin,jmin,-1 |
70 |
|
71 |
if (lphprint .and. myid .eq. 0) |
72 |
& print *, 'pathei-lsopt: in hessupd; loop ', |
73 |
& 'j,jfin,jmin = ', j,jfin,jmin |
74 |
|
75 |
jp = j |
76 |
if (jp.gt.mupd) jp = jp-mupd |
77 |
call dostore( nn, xdiff, .false., 2*jp+3 ) |
78 |
|
79 |
temp = DDOT( nn, dd, 1, xdiff, 1 ) |
80 |
call MPI_ALLREDUCE(temp,r,1,MPI_DOUBLE_PRECISION,MPI_SUM, |
81 |
& MPI_COMM_WORLD,ierr) |
82 |
! r = DDOT( nn, dd, 1, xdiff,1 ) |
83 |
|
84 |
call dostore( nn, xdiff, .false., 2*jp+2 ) |
85 |
alpha(jp) = r |
86 |
do i = 1, nn |
87 |
dd(i) = dd(i) - r*xdiff(i) |
88 |
end do |
89 |
end do |
90 |
|
91 |
c------------------------------------ |
92 |
c multiply precondition matrix |
93 |
c------------------------------------ |
94 |
if (mupd .ne. 0) then |
95 |
call dostore( nn, xdiff, .false., 3 ) |
96 |
do i = 1, nn |
97 |
dd(i) = dd(i)*xdiff(i) |
98 |
end do |
99 |
end if |
100 |
|
101 |
c------------------------------------ |
102 |
c compute left hand side |
103 |
c------------------------------------ |
104 |
do j = jmin,jfin |
105 |
jp = j |
106 |
if (jp .gt. mupd) jp = jp-mupd |
107 |
call dostore( nn, xdiff, .false., 2*jp+2 ) |
108 |
|
109 |
temp = DDOT( nn, dd, 1, xdiff, 1 ) |
110 |
call MPI_ALLREDUCE(temp,r,1,MPI_DOUBLE_PRECISION,MPI_SUM, |
111 |
& MPI_COMM_WORLD,ierr) |
112 |
r = alpha(jp) - r |
113 |
! r = alpha(jp) - DDOT( nn, dd,1 , xdiff, 1 ) |
114 |
|
115 |
call dostore( nn, xdiff, .false., 2*jp+3 ) |
116 |
do i = 1, nn |
117 |
dd(i) = dd(i) + r*xdiff(i) |
118 |
end do |
119 |
end do |
120 |
|
121 |
return |
122 |
|
123 |
end |