ralna / RALFit

A non-linear least squares solver
Other
25 stars 6 forks source link

Inconsistent behaviour returning from `eval_HP` #43

Closed talassio closed 3 years ago

talassio commented 5 years ago

In eval_HP return status is not monitored yet HP parameter is of intent=OUT. See https://github.com/ralna/RALFit/blob/2b0dfc8bb9d34e8d1c302842964afa951bcfcacb/libRALFit/src/ral_nlls_internal.f90#L3597

Compiling with nagfor -nan will return NANs in HP illustrating the issue and thus breaking.

Issues are:

  1. even if eval_hp_provided=.True. return status should be monitored and allow calls to dgemv if appropriate, otherwise it should abort with ERROR status -which at me moment is not possible-.

  2. current behaviour is to return HP=0 and this is not correct.

subroutine calculate_sHs( n, m, s, params) ! change to
subroutine calculate_sHs( n, m, s, params, status)

       If (params%eval_hp_provided) Then
          Call params%eval_hp(en_data=en_data,objhprd=objhprd,                 &
            objhprdh=objhprdh,iuser=iuser,ruser=ruser,status=status,n=n,       &
            m=params%m,x=params%x,y=s(1:n),hp=params%tenj%hs,                  &
            params=params%parent_params)
          If (status/=0) Then
!           Error eval_hp returned non-zero status, propagate back!
            Go To 100
          End If
        Else
          Do ii = 1, params%m
            params%tenj%h(1:n,1:n) = params%hi(1:n,1:n,ii)
            params%tenj%hs(1:n,ii) = 0.0_wp
            Call dgemv(trans='N',m=n,n=n,alpha=1.0_wp,a=params%tenj%h(1,1),    &
              lda=n,x=s,incx=1,beta=0.0_wp,y=params%tenj%hs(1,ii),incy=1)
          End Do
        End If
        Call dgemv(trans='T',m=n,n=params%m,alpha=1.0_wp,                      &
          a=params%tenj%hs(1,1),lda=n,x=s,incx=1,beta=0.0_wp,                  &
          y=params%tenj%sths(1),incy=1)
100     Continue