#include <oxstd.h>
#include <oxfloat.h>  // defines M_PI, M_2PI
#import <maximize>

static decl s_mY;     // the data sample (T x 1)

fLoglik(const vP, const adFunc, const avScore,const amHess)
{
    decl dsum, dsig2;

    dsum = sumsqrc(s_mY - vP[0]) / rows(s_mY);
    dsig2 = vP[1];
    
    adFunc[0] = -0.5 * (log(M_2PI * dsig2) + dsum / dsig2);

return 1;             // 1 indicates success
}
main()
{
    decl cn, dmu, dsigma2, vtheta, dfunc;

    cn = 50;          // sample size
    dmu = 21;         // distribution parameter: mean
    dsigma2 = 49;     // and variance
                      // generate a sample
    s_mY = dmu + sqrt(dsigma2) * rann(cn, 1);

    vtheta = <20;49>;
    fLoglik(vtheta, &dfunc, 0, 0);     // evaluate

    print("function value is ", cn * dfunc, 
        " at parameter value:", vtheta');

    MaxControl(-1,5);
    MaxBFGS(fLoglik, &vtheta, &dfunc, 0, TRUE);

    print("Function value is ", cn * dfunc,
        " at parameter value:", vtheta');
}
