NAG Library Manual, Mark 28.3
Interfaces:  FL   CL   CPP   AD 

NAG AD Library Introduction
Example description
!   C05AY_A1W_F Example Program Text
!   Mark 28.3 Release. NAG Copyright 2022.
    Module c05ay_a1w_fe_mod

!     C05AY_A1W_F Example Program Module:
!            Parameters and User-defined Routines

!     .. Use Statements ..
      Use iso_c_binding, Only: c_funloc, c_int, c_ptr
      Use nagad_library, Only: exp, nagad_a1w_get_derivative,                  &
                               nagad_a1w_inc_derivative,                       &
                               nagad_a1w_ir_register_variable,                 &
                               nagad_a1w_w_rtype, nagad_primal, x10ad_a1w_f,   &
                               x10ba_a1w_f, x10bb_a1w_f, x10bd_a1w_f,          &
                               x10be_a1w_f, x10bj_a1w_f, x10ce_a1w_f,          &
                               x10cj_a1w_f, Assignment (=), Operator (-),      &
                               Operator (*)
      Use nag_library, Only: nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Accessibility Statements ..
      Private
      Public                           :: cb_sym, f_a1w
!     .. Parameters ..
      Integer, Parameter, Public       :: nin = 5, nout = 6
    Contains
      Subroutine f_a1w(ad_handle,x,z,iuser,ruser)

!       .. Use Statements ..
        Use nagad_library, Only: nagad_algorithmic, nagad_symbolic
!       .. Scalar Arguments ..
        Type (c_ptr), Intent (Inout)   :: ad_handle
        Type (nagad_a1w_w_rtype), Intent (In) :: x
        Type (nagad_a1w_w_rtype), Intent (Out) :: z
!       .. Array Arguments ..
        Type (nagad_a1w_w_rtype), Intent (Inout) :: ruser(*)
        Integer, Intent (Inout)        :: iuser(*)
!       .. Local Scalars ..
        Type (c_ptr)                   :: cb_handle
        Integer                        :: cb_mode, ifail, mode
!       .. Executable Statements ..
        ifail = 0
        Call x10ad_a1w_f(ad_handle,mode,ifail)

        z = 0.0_nag_wp
        If (mode==nagad_algorithmic) Then
!         Use module extensions to operators and intrinsics.
          z = exp(-x) - x*ruser(1)
        Else If (mode==nagad_symbolic) Then
          Call x10bd_a1w_f(ad_handle,cb_mode,ifail)
          If (cb_mode==nagad_primal) Then
!           Evaluate primal only
            z%value = exp(-x%value) - x%value*ruser(1)%value
            z%tapeindex = 0
          Else
!           Adjoint w.r.t to x or ruser

!           Create callback data object
            ifail = 0
            Call x10ba_a1w_f(cb_handle,ifail)
!           Write cb_mode to object
            Call x10be_a1w_f(cb_handle,cb_mode,ifail)
!           Write inputs to object
            Call x10bj_a1w_f(cb_handle,x,ifail)
            Call x10bj_a1w_f(cb_handle,ruser(1),ifail)

!           Evaluate primal
            z%value = exp(-x%value) - x%value*ruser(1)%value
            z%tapeindex = 0

!           Register output
            Call nagad_a1w_ir_register_variable(z)
!           Write output to object
            Call x10bj_a1w_f(cb_handle,z,ifail)

!           Insert callback
            Call x10bb_a1w_f(cb_handle,c_funloc(cb_sym),ifail)

          End If
        End If
        Return

      End Subroutine f_a1w
      Subroutine cb_sym(callmode,cb_handle) Bind (C)

!       .. Use Statements ..
        Use nagad_library, Only: nagad_dall, nagad_dparam, nagad_dstate
!       .. Scalar Arguments ..
        Type (c_ptr), Value            :: cb_handle
        Integer (Kind=c_int), Value    :: callmode
!       .. Local Scalars ..
        Type (nagad_a1w_w_rtype)       :: ruser, x, z
        Real (Kind=nag_wp)             :: ra, xa, za
        Integer                        :: cb_mode, ifail
!       .. Executable Statements ..

        If (callmode==1) Then
          Return
        End If

!       Extract data from callback data object
        ifail = 0
        Call x10ce_a1w_f(cb_handle,cb_mode,ifail)
        Call x10cj_a1w_f(cb_handle,x,ifail)
        Call x10cj_a1w_f(cb_handle,ruser,ifail)
        Call x10cj_a1w_f(cb_handle,z,ifail)

        za = nagad_a1w_get_derivative(z)

        If (cb_mode==nagad_dstate .Or. cb_mode==nagad_dall) Then
!         d/dx
          xa = (-ruser%value-exp(-x%value))*za
          Call nagad_a1w_inc_derivative(x,xa)
        End If
        If (cb_mode==nagad_dparam .Or. cb_mode==nagad_dall) Then
!         d/druser
          ra = -x%value*za
          Call nagad_a1w_inc_derivative(ruser,ra)
        End If

      End Subroutine cb_sym

    End Module c05ay_a1w_fe_mod

    Program c05ay_a1w_fe
!     C05AY_A1W_F Example Main Program

!     .. Use Statements ..
      Use c05ay_a1w_fe_mod, Only: f_a1w, nin, nout
      Use iso_c_binding, Only: c_ptr
      Use nagad_library, Only: c05ay_a1w_f, nagad_a1w_get_derivative,          &
                               nagad_a1w_ir_create => x10za_a1w_f,             &
                               nagad_a1w_ir_interpret_adjoint_sparse,          &
                               nagad_a1w_ir_register_variable,                 &
                               nagad_a1w_ir_remove, nagad_a1w_set_derivative,  &
                               nagad_a1w_w_rtype, nagad_algorithmic,           &
                               x10aa_a1w_f, x10ab_a1w_f, x10ac_a1w_f,          &
                               x10ae_a1w_f, x10af_a1w_f, Assignment (=)
      Use nag_library, Only: nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Local Scalars ..
      Type (nagad_a1w_w_rtype)         :: a, b, eps, eta, x
      Type (c_ptr)                     :: ad_handle
      Real (Kind=nag_wp)               :: ar, br, da, db, druser, epsr, etar,  &
                                          rr
      Integer                          :: ifail, mode
!     .. Local Arrays ..
      Type (nagad_a1w_w_rtype)         :: ruser(1)
      Integer                          :: iuser(1)
!     .. Executable Statements ..

      Write (nout,*) 'C05AY_A1W_F Example Program Results'

!     Skip first line of data file
      Read (nin,*)
!     Read problem parameters
      Read (nin,*) ar
      Read (nin,*) br
      Read (nin,*) epsr
      Read (nin,*) etar
      Read (nin,*) rr

      a = ar
      b = br
      eps = epsr
      eta = etar
      ruser(1) = rr

!     Create AD tape
      Call nagad_a1w_ir_create

!     Create AD configuration data object
      ifail = 0
      Call x10aa_a1w_f(ad_handle,ifail)

!     Read AD computational mode
      Read (nin,*) mode
!     Set this mode.
      Call x10ac_a1w_f(ad_handle,mode,ifail)

!     Register variables to differentiate w.r.t.
      Call nagad_a1w_ir_register_variable(a)
      Call nagad_a1w_ir_register_variable(b)
      Call nagad_a1w_ir_register_variable(ruser)

!     Call AD routine
      ifail = 0
      Call x10ae_a1w_f(ad_handle,1,ifail)
      ifail = 0
      Call x10af_a1w_f(ad_handle,-1,ifail)
      Call c05ay_a1w_f(ad_handle,a,b,eps,eta,f_a1w,x,iuser,ruser,ifail)

!     Setup evaluation of derivatives via adjoints
      Call nagad_a1w_set_derivative(x,1.0_nag_wp)
      ifail = 0
      Call nagad_a1w_ir_interpret_adjoint_sparse(ifail)

      Write (nout,*)
      Write (nout,*) ' Derivatives calculated: First order adjoints'
      If (mode==nagad_algorithmic) Then
        Write (nout,*) ' Computational mode    : algorithmic'
      Else
        Write (nout,*) ' Computational mode    : symbolic'
      End If

!     Get derivatives
      da = nagad_a1w_get_derivative(a)
      db = nagad_a1w_get_derivative(b)
      druser = nagad_a1w_get_derivative(ruser(1))

      Write (nout,*)
      Write (nout,99999) ' Solution, x =', x%value
      Write (nout,*) ' Derivatives:'
      Write (nout,99999) '     d/da(x) =', da
      Write (nout,99999) '     d/db(x) =', db
      Write (nout,99999) ' d/druser(x) =', druser
99999 Format (1X,A,1X,E12.5)

!     Remove computational data object and tape
      Call x10ab_a1w_f(ad_handle,ifail)
      Call nagad_a1w_ir_remove

    End Program c05ay_a1w_fe