Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.

TR Solver performance

OP1
New Contributor II
269 Views

Hi:

I am trying to use the TR solver to solve a trivial set of two nonlinear equations. The equations are f1 = cosd(x1)+0.5 and f2=sind(x2)+0.5 (they are not even coupled). I provide an initial guess of x1 = 100 and x2 = -20, which is not that far from the solution x1 = 120 and x2 = -30, as well as a reasonable value of EPS (1.0D-5).
The solver takes 597 iterations to solve this and requires as many jacobian calculations?? I am puzzled. It's probably the case of me looking for too long at this example and not seing the obvious anymore :) .

The code is attached below. Thanks for your help in advance!

Olivier

[fortran]PROGRAM TEST_MKL_TR_SOLVER
IMPLICIT NONE
INTEGER,PARAMETER :: DP = 8

INTEGER TR_SUCCESS	        
PARAMETER (TR_SUCCESS = 1501)
INTEGER TR_INVALID_OPTION
PARAMETER (TR_INVALID_OPTION = 1502)
INTEGER TR_OUT_OF_MEMORY
PARAMETER (TR_OUT_OF_MEMORY = 1503)
EXTERNAL DTRNLSPBC_INIT
INTEGER DTRNLSPBC_INIT
EXTERNAL DTRNLSPBC_SOLVE
INTEGER DTRNLSPBC_SOLVE
EXTERNAL DTRNLSPBC_GET
INTEGER DTRNLSPBC_GET
EXTERNAL DTRNLSPBC_DELETE
INTEGER DTRNLSPBC_DELETE

INTEGER :: STATUS,N,M,ITER1,ITER2,REQUEST,RES,ITER,ST_CR,N_FVEC,N_FJAC
INTEGER(8) :: HANDLE
REAL(KIND=DP) :: X(2),EPS(6),RS,FVEC(2),FJAC(2,2),R1,R2,LB(2),UB(2)
LOGICAL :: KEEP_GOING

N = 2
M = 2
N_FVEC = 0
N_FJAC = 0

LB = -180.0_DP
UB =  180.0_DP

X(1)  = 100.0_DP
X(2)  = -20.0_DP
ITER1 = 1000
ITER2 = 100
RS    = 100.0_DP
EPS   = 1.0E-5_DP

STATUS = DTRNLSPBC_INIT(HANDLE,N,M,X,LB,UB,EPS,ITER1,ITER2,RS)

FVEC = 0.0_DP
FJAC = 0.0_DP

KEEP_GOING = .TRUE.
REQUEST = 0

DO WHILE (KEEP_GOING)

    RES = DTRNLSPBC_SOLVE(HANDLE,FVEC,FJAC,REQUEST)
    
    IF (RES/=TR_SUCCESS) THEN
    
        WRITE(*,*) 'Error in DTRNLSP_SOLVE'
        STOP
    
    END IF
    
    SELECT CASE (REQUEST)

        CASE (-1,-2,-3,-4,-5,-6)
        
            KEEP_GOING = .FALSE.

        CASE (1)

            FVEC(1) = COSD(X(1))+0.5_DP
            FVEC(2) = SIND(X(2))+0.5_DP
            N_FVEC = N_FVEC+1
    
        CASE (2)
        
            FJAC(1,1) = -SIND(X(1))
            FJAC(2,1) = 0.0_DP
            FJAC(1,2) = 0.0_DP
            FJAC(2,2) = COSD(X(2))
            N_FJAC = N_FJAC+1

    END SELECT
   
END DO

RES = DTRNLSPBC_GET(HANDLE,ITER,ST_CR,R1,R2)

WRITE(*,*) 'Last request: ',REQUEST
WRITE(*,*) 'Number of iterations: ',ITER
WRITE(*,*) 'Number of FVEC calculations: ',N_FVEC
WRITE(*,*) 'Number of FJAC calculations: ',N_FJAC
WRITE(*,*) 'Stop criterion: ',ST_CR
WRITE(*,*) 'Last X: ',X

STATUS = DTRNLSPBC_DELETE(HANDLE)

END PROGRAM TEST_MKL_TR_SOLVER[/fortran]
0 Kudos
2 Replies
mecej4
Honored Contributor III
269 Views
The derivative of sind(x) w.r.t. x is not cosd(x), nor is that of cosd(x) equal to -sind(x).

The penalty for providing wrong expressions for derivatives is slow convergence (or even divergence). Without correct derivatives, there is no hope of reaching ultimate second-order convergence.

0 Kudos
OP1
New Contributor II
269 Views

Duh!!! Thanks a lot Mecej4... Can't believe I missed that one... after correction this converges in 3 iterations...

Olivier

0 Kudos
Reply