/* D01UA_T1W_F C++ Header Example Program.
*
* Copyright 2024 Numerical Algorithms Group.
* Mark 30.1, 2024.
*/
#include <dco.hpp>
#include <iostream>
#include <math.h>
#include <nag.h>
#include <nagad.h>
#include <stdio.h>
using namespace std;
int main()
{
// Scalars
int exit_status = 0;
Integer n = 64;
// Arrays
nagad_t1w_w_rtype ruser[1];
Integer iuser[1];
ruser[0] = 1.0;
cout << "D01UA_T1W_F C++ Header Example Program Results\n\n";
cout << "\n Derivatives calculated: First order tangents\n";
cout << " Computational mode : algorithmic\n";
// Create AD configuration Data object
Integer ifail = 0;
nag::ad::handle_t ad_handle;
nagad_t1w_w_rtype a, b;
Integer key;
for (int funid = 1; funid <= 6; funid++)
{
if (funid == 1)
{
a = 0.0;
b = 1.0;
key = 0;
}
if (funid == 2)
{
a = 2.0;
b = 0.0;
key = -5;
}
if (funid == 3)
{
a = 2.0;
b = 1.0;
key = -3;
}
if (funid == 4)
{
a = 2.0;
b = 1.0;
key = 3;
}
if (funid == 5)
{
a = -1.0;
b = 3.0;
key = -4;
}
if (funid == 6)
{
a = -1.0;
b = 3.0;
key = 4;
}
iuser[0] = funid;
auto f = [&](nag::ad::handle_t & ad_handle,
const nagad_t1w_w_rtype *x,
const Integer & nx,
nagad_t1w_w_rtype *fv,
Integer & iflag)
{
// dco/c++ overloading used here to perform AD
Integer i;
nagad_t1w_w_rtype xx, fx, rux;
for (i = 0; i < nx && iflag >= 0; i++)
{
xx = x[i];
rux = ruser[0] * xx;
if (iuser[0] == 1)
{
fx = 4.0 / (1.0 + ruser[0] * xx * xx);
}
if (iuser[0] == 2)
{
fx = 1.0 / (xx * xx * log(ruser[0] * xx));
}
if (iuser[0] == 3)
{
fx = exp(-ruser[0] * xx) / xx;
}
if (iuser[0] == 4)
{
fx = ruser[0] / xx;
}
if (iuser[0] == 5)
{
fx = exp(-3.0 * ruser[0] * xx * xx - 4.0 * xx - 1.0);
}
if (iuser[0] == 6)
{
fx = exp(2.0 * ruser[0] * xx + 2.0);
}
if (iuser[0] < 0 || iuser[0] > 6)
{
iflag = -1;
}
fv[i] = fx;
}
};
// Call the AD routine
ifail = 0;
nagad_t1w_w_rtype dinest;
double inc = 1.0, zero = 0.0;
dco::derivative(ruser[0]) = inc;
nag::ad::d01ua(ad_handle, key, a, b, n, f, dinest, ifail);
dco::derivative(ruser[0]) = zero;
// Get derivatives
double dr;
dr = dco::derivative(dinest);
cout.setf(ios::scientific, ios::floatfield);
cout.setf(ios::right);
cout.precision(4);
cout << "\n dinest = ";
double dinest_value = dco::value(dinest);
cout.width(12);
cout << dinest_value << " ";
cout << " d/druser[0] = ";
cout.width(12);
cout << dr << endl;
}
return exit_status;
}