1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
|
! Simple self-test based on the driver*.f90 examples. Returns 1 in
! the incorrect value is reached. Only other change is the removal of
! long headers and comments.
!
! Copyright Ciyou Zhu, Richard Byrd, Jorge Nocedal, Jose Luis
! Morales. 3-clause BSD license.
program driver
implicit none
!
integer, parameter :: n = 25, m = 5, iprint = 1
integer, parameter :: dp = kind(1.0d0)
real(dp), parameter :: factr = 1.0d+7, pgtol = 1.0d-5
!
character(len=60) :: task, csave
logical :: lsave(4)
integer :: isave(44)
real(dp) :: f
real(dp) :: dsave(29)
integer, allocatable :: nbd(:), iwa(:)
real(dp), allocatable :: x(:), l(:), u(:), g(:), wa(:)
! Declare a few additional variables for this sample problem
real(dp) :: t1, t2
integer :: i
! Allocate dynamic arrays
allocate ( nbd(n), x(n), l(n), u(n), g(n) )
allocate ( iwa(3*n) )
allocate ( wa(2*m*n + 5*n + 11*m*m + 8*m) )
!
do 10 i=1, n, 2
nbd(i) = 2
l(i) = 1.0d0
u(i) = 1.0d2
10 continue
! Next set bounds on the even-numbered variables.
do 12 i=2, n, 2
nbd(i) = 2
l(i) = -1.0d2
u(i) = 1.0d2
12 continue
! We now define the starting point.
do 14 i=1, n
x(i) = 3.0d0
14 continue
write (6,16)
16 format(/,5x, 'Solving sample problem.', &
/,5x, ' (f = 0.0 at the optimal solution.)',/)
! We start the iteration by initializing task.
task = 'START'
! The beginning of the loop
do while(task(1:2).eq.'FG'.or.task.eq.'NEW_X'.or. &
task.eq.'START')
! This is the call to the L-BFGS-B code.
call setulb ( n, m, x, l, u, nbd, f, g, factr, pgtol, &
wa, iwa, task, iprint,&
csave, lsave, isave, dsave )
if (task(1:2) .eq. 'FG') then
f=.25d0*( x(1)-1.d0 )**2
do 20 i=2, n
f = f + ( x(i)-x(i-1 )**2 )**2
20 continue
f = 4.d0*f
! Compute gradient g for the sample problem.
t1 = x(2) - x(1)**2
g(1) = 2.d0*(x(1) - 1.d0) - 1.6d1*x(1)*t1
do 22 i=2, n-1
t2 = t1
t1 = x(i+1) - x(i)**2
g(i) = 8.d0*t2 - 1.6d1*x(i)*t1
22 continue
g(n) = 8.d0*t1
end if
end do
if (abs(f) < 1e-6) then
stop 0
else
stop 1
end if
! end of loop do while
end program driver
|