#include <oxstd.h>
#import <maximize>

fRosenbrock(const vP, const adFunc, const avScore, const amHessian)
{
    decl ab2 = vP[1] - vP[0] ^ 2, a1 = 1 - vP[0];

    adFunc[0] = -100 * ab2 ^ 2 - a1 ^ 2;
    if (avScore)                   // if !0: compute score
    {
        (avScore[0])[0] = 400 * ab2 * vP[0] + 2 * a1;
        (avScore[0])[1] = -200 * ab2;
    }
    if (amHessian)               // if !0: compute Hessian
    {  // this bit is not needed for numerical derivatives
        (amHessian[0])[0][0] =
            400 * vP[1] - 1200 * vP[0]^2 - 2;
        (amHessian[0])[1][1] = -200;
        (amHessian[0])[1][0] =
            (amHessian[0])[0][1] = 400 * vP[0];
    }
return 1;                           // 1 indicates success
}

main()
{
    decl vp, dfunc, ir, mhess;

    MaxControl(100, 1, 1);

    vp = zeros(2, 1);                   // starting values
    ir = MaxNewton(fRosenbrock, &vp, &dfunc,&mhess, TRUE);

    print("\n", MaxConvergenceMsg(ir),
        " using numerical 2nd derivatives",
        "\nFunction value = ", dfunc, "; parameters:", vp,
        "final Hessian:", mhess);   

    vp = zeros(2, 1);                   // starting values
    ir = MaxNewton(fRosenbrock, &vp, &dfunc,&mhess,FALSE);

    print("\n", MaxConvergenceMsg(ir),
        " using analytical 2nd derivatives",
        "\nFunction value = ", dfunc, "; parameters:", vp,
        "final Hessian:", mhess);   
}

