subroutine lsline( & simul & , nn, ifail, lphprint & , ifunc, nfunc & , ff, dotdg & , sflag, tflag, tact, epsx & , dd, gg, xx, xdiff & ) c ================================================================== c SUBROUTINE lsline c ================================================================== c c o line search algorithm for determining control vector update; c After computing updated control vector from given gradient, c a forward and adjoint model run are performed (simul.F) c using the updated control vector. c Tests are then applied to see whether solution has improved. c c o Reference: J.C. Gilbert & C. Lemarechal c Some numerical experiments with variable-storage c quasi-Newton algorithms c Mathematical Programming 45 (1989), pp. 407-435 c c o started: ??? not reproducible c c o changed: Benny Cheng, NASA/JPL c c ================================================================== c SUBROUTINE lsline c ================================================================== #include "blas1.h" implicit none include 'mpif.h' c---------------------------------- c declare arguments c---------------------------------- integer nmax parameter( nmax = MAX_INDEPEND ) integer nn, ifail, ifunc, nfunc double precision ff, dotdg, tact, epsx double precision xx(nn), gg(nn), xdiff(nn),dd(nn) logical tflag,sflag,isforward logical lphprint external simul integer pidlen,myindx(2) integer status(MPI_STATUS_SIZE),ierr common /mpi_parm/myid,nprocs,mystart,myend integer myid,nprocs,mystart,myend c---------------------------------- c declare local variables c---------------------------------- double precision xpara1, xpara2 parameter( xpara1 = 0.001D0, xpara2=0.9D0 ) double precision factor parameter( factor = 0.2 ) double precision barmax parameter( barmax = 0.3 ) double precision barmul parameter( barmul = 5.0 ) double precision barmin parameter( barmin = 0.01 ) double precision tmax parameter( tmax = 5. ) integer i double precision tg, fg, td, ta double precision fa, fpa, fp double precision fnew, fdiff,ftemp double precision z, test, barr double precision left, right, told double precision t,a,b double precision xtemp(nn),gtemp(nn),gold(nn) external DDOT,DNRM2 double precision DDOT, DNRM2 double precision temp,gnorm0 real*8,allocatable:: vv(:) integer, allocatable:: displs(:) integer, allocatable:: counts(:) c---------------------------------- c check parameters c---------------------------------- if ( (nn.le.0) & .or. (xpara1.le.0.0) .or. (xpara1.ge.0.5) & .or. (xpara2.le.xpara1) .or. (xpara2.ge.1.0) ) then ifail = 9 go to 999 endif c---------------------------------- c initialization c---------------------------------- c======================================================================= c begin of simulation iter. c======================================================================= if (lphprint .and. myid .eq. 0) & print *, 'lsopt: simul.' c------------------------------------ c compute cost function and gradient c------------------------------------ gold = gg temp = DNRM2(nn,gg,1) call MPI_ALLREDUCE(temp*temp,gnorm0,1,MPI_DOUBLE_PRECISION, & MPI_SUM,MPI_COMM_WORLD,ierr) gnorm0 = dsqrt(gnorm0) if (sflag ) then !get estimated cost from hessian update control, this is the trial step isforward = .true. call simul( 0, isforward,nn, xdiff, fnew, gg ) xtemp = xdiff ! new direction vector dd temp = DDOT( nn, xdiff - xx, 1, gold, 1 ) call MPI_ALLREDUCE(temp,b,1,MPI_DOUBLE_PRECISION,MPI_SUM, & MPI_COMM_WORLD,ierr) dotdg = b c----------------------------------------- c quadratic line search c----------------------------------------- a = fnew - ff - b t = -b/(2.*a) t = dmax1(t,0.D0) if (t <= 1.d-2 ) then tact = 1. xdiff = xx - 2.*(.01*ff)*gold*tact/(gnorm0*gnorm0) tflag = .false. if (lphprint .and. myid .eq. 0) & print *, 'lsopt: using steepest descent, t value ',t,'< .01' else xdiff = (1. - t)*xx + t*xdiff tact = t tflag = .true. endif if (lphprint .and. myid .eq. 0) & print *, 'lsopt: trial stepsize =',t ! no hessian update , get new adjoint run goto 999 else !no line search !read new ctrl vector and new cost gradient and new cost isforward = .false. call simul( 0, isforward, nn, xdiff, fnew, gg ) if (tflag) then !read in trial step adjoint output call simul( ifunc+1, isforward, nn, xtemp, ftemp, gtemp) dd = xtemp - xx c----------------------------------------- c apply 1st Wolfe test c----------------------------------------- fdiff = fnew - ff temp = DDOT( nn, dd, 1, gold, 1 ) call MPI_ALLREDUCE(temp,dotdg,1,MPI_DOUBLE_PRECISION,MPI_SUM, & MPI_COMM_WORLD,ierr) if (fdiff .gt. tact*xpara1*dotdg) then ifail = 7 if (lphprint .and. myid .eq. 0) & print *, 'lsopt: fails wolfe test 1' endif c----------------------------------------- c 1st Wolfe test 1 ok, apply 2nd Wolf test c----------------------------------------- temp = DDOT( nn, dd, 1, gg, 1 ) call MPI_ALLREDUCE(temp,fp,1,MPI_DOUBLE_PRECISION,MPI_SUM, & MPI_COMM_WORLD,ierr) if (fp .lt. xpara2*dotdg) then ifail = 7 if (lphprint .and. myid .eq. 0) & print *, 'lsopt: fails wolfe test 2' endif c----------------------------------------- c 2nd Wolfe test 2 ok, donc pas serieux, on sort c----------------------------------------- if (ifail .eq. 7 ) then if (ifunc .lt. nfunc) then t = dmin1(tact,tmax) if (t < 5.D0 ) t = .5*t xdiff = (1. - t)*xx + t*xtemp tact = t ifunc = ifunc + 1 goto 999 endif endif !end ifail = 7 else if ( fnew > .999D0*ff ) then if (ifunc .lt. nfunc) then tact = tact*.5 xdiff = xx - 2.*(.01*ff)*gold*tact/(gnorm0*gnorm0) ifunc = ifunc + 1 goto 999 endif endif endif ! end of tflag dd = xdiff - xx ff = fnew xx =xdiff ifunc = 0 tact = 1. endif ! end of sflag c----------------------------------------- c end of routine c----------------------------------------- 999 continue return end