nepluno / lbfgsb-gpu

An open source library for the GPU-implementation of L-BFGS-B algorithm
http://dx.doi.org/10.1016/j.cag.2014.01.002
Mozilla Public License 2.0
116 stars 17 forks source link

Logical difference between fortran version and active.cu: upper bound not enforced in active #6

Closed fzimmermann89 closed 1 year ago

fzimmermann89 commented 2 years ago

In https://github.com/nepluno/lbfgsb-gpu/blob/f3522f6ef2903cdb0628a37f346b0e19d5970f70/culbfgsb/active.cu#L34-L40

Whereas in Nocedal et als fortran code (http://users.iems.northwestern.edu/~nocedal/lbfgsb.html) it is

         if (nbd(i) .gt. 0) then
            if (nbd(i) .le. 2 .and. x(i) .le. l(i)) then
               if (x(i) .lt. l(i)) then
                  prjctd = .true.
                  x(i) = l(i)
               endif
               nbdd = nbdd + 1
            else if (nbd(i) .ge. 2 .and. x(i) .ge. u(i)) then
               if (x(i) .gt. u(i)) then
                  prjctd = .true.
                  x(i) = u(i)
               endif
               nbdd = nbdd + 1
            endif
         endif

In the current active.cu code, for nbdi==2 the upper bounds is not enforced.

In the fortran code, the else if (nbd(i) .ge. 2 .and. x(i) .ge. u(i)) will be hit, but in the cuda code, the first if condition already matches, the max is a noop if x>ub>lb and the else will not be hit, so the min is never called and the upper bound not enforced.

Not sure what this means for the rest of the algorithm though..

fzimmermann89 commented 2 years ago

and somewhat related: is the meaning of nbd the same as in the original code?

c     nbd is an INTEGER array of dimension n that must be set by the
c       user to the type of bounds imposed on the variables:
c       nbd(i)=0 if x(i) is unbounded,
c              1 if x(i) has only a lower bound,
c              2 if x(i) has both lower and upper bounds, 
c              3 if x(i) has only an upper bound.

In the readme, the nbd[i]==3 case is not mentioned, but I assume you kept the logic around nbd the same?

nepluno commented 1 year ago

Thanks for reporting the issue! Fixed the bug according to your commit.