/* nag_opt_handle_solve_bounds_foas (e04kfc) Example Program.
 * 
 * Copyright 2019 Numerical Algorithms Group.
 *
 * Mark 27.0, 2019.
 */

/* 
 * NLP example: Nonlinear objective + box constraint
 */

#include <stdio.h>
#include <math.h>
#include <nag.h>

#ifdef __cplusplus
extern "C" {
#endif

    static void NAG_CALL objfun(Integer nvar,
                                const double x[],
                                double *fx,
                                Integer * inform, Nag_Comm * comm);
    static void NAG_CALL objgrd(Integer nvar,
                                const double x[],
                                Integer nnzfd,
                                double fdx[],
                                Integer * inform, Nag_Comm * comm);
    static void NAG_CALL monit(Integer nvar,
                               const double x[],
                               Integer * inform,
                               const double rinfo[],
                               const double stats[], Nag_Comm * comm);

#ifdef __cplusplus
}
#endif
static void NAG_CALL objfun(Integer nvar,
                            const double x[],
                            double *fx, Integer * inform, Nag_Comm * comm)
{
    *fx = pow(1.0 - x[0], 2.0) + 100.0 * pow(x[1] - pow(x[0], 2.0), 2.0);
}

static void NAG_CALL objgrd(Integer nvar,
                            const double x[],
                            Integer nnzfd,
                            double fdx[],
                            Integer * inform, Nag_Comm * comm)
{
    fdx[0] = 2.0 * x[0] - 400.0 * x[0] * (x[1] - pow(x[0], 2.0)) - 2.0;
    fdx[1] = 200.0 * (x[1] - pow(x[0], 2.0));
}

static void NAG_CALL monit(Integer nvar,
                           const double x[],
                           Integer * inform,
                           const double rinfo[],
                           const double stats[], Nag_Comm * comm)
{
    return;
}

int main(void)
{
/* Define problem size */
#define NVAR 2

    /* Scalars */
    Integer nvar = NVAR;
    Integer nu = 2*nvar;

    /* Arrays */
    /* Lower bounds */
    double blx[NVAR] = { -1.0, -2.0 };
    /* Upper bounds */
    double bux[NVAR] = { 0.8, 2.0 };
    /* Indexes for the gradient entries */
    Integer idxfd[NVAR] = { 1, 2 };
    /* Initial guess x0 */
    double x[NVAR] = { -1.5, 1.9 };
    /* space holder for the Lagrange multiplyers */
    double u[2*NVAR];
    double rinfo[100], stats[100];
    void *handle = NULL;
    char opt[80];

    /* Nag Types */
    Nag_Comm comm;
    NagError fail;
    Nag_FileID mon_out;

    /* nag_file_open (x04acc).
     *  Open unit number for reading, writing or appending, and
     *  associate unit with named file  
     */
    nag_file_open("e04kf_cmon.txt", 1, &mon_out, NAGERR_DEFAULT);

    /* nag_file_line_write (x04bac).
     * Write formatted record to external file
     */
    Vprintf("nag_opt_handle_solve_bounds_foas (e04kfc) "
                   "Example Program Results\n");
    fflush(stdout);
    
    /* nag_opt_handle_init (e04rac).
     * Initialize an empty problem handle with NVAR variables. 
     */
    nag_opt_handle_init(&handle, nvar, NAGERR_DEFAULT);

    /* nag_opt_handle_set_simplebounds (e04rhc).
     * Define bounds on the variables
     */
    nag_opt_handle_set_simplebounds(handle, nvar, blx, bux, NAGERR_DEFAULT);

    /* nag_opt_handle_set_nlnobj (e04rgc).
     * Define nonlinear objective 
     */
    nag_opt_handle_set_nlnobj(handle, nvar, idxfd, NAGERR_DEFAULT);

    /* nag_opt_handle_opt_set (e04zmc).
     * Set optional arguments of the solver 
     */
    nag_opt_handle_opt_set(handle, "FOAS Print Frequency = 5", NAGERR_DEFAULT);

    nag_opt_handle_opt_set(handle, "Print Solution = yes", NAGERR_DEFAULT);

    sprintf(opt, "Monitoring File = %" NAG_IFMT, mon_out);
    nag_opt_handle_opt_set(handle, opt, NAGERR_DEFAULT);

    nag_opt_handle_opt_set(handle, "Monitoring Level = 3", NAGERR_DEFAULT);

    nag_opt_handle_opt_set(handle, "Print Level = 1", NAGERR_DEFAULT);

    /* Solver the problem 
     * nag_opt_handle_solve_bounds_foas (e04kfc).
     */
    SET_FAIL(fail);
    nag_opt_handle_solve_bounds_foas(handle, objfun, objgrd, monit, nvar,
                                     x, rinfo, stats, &comm, &fail);

    /* Print objective value at solution */
    if (fail.code==NE_NOERROR || fail.code==NW_NOT_CONVERGED)
    {
      Vprintf("\n\n  Solution found.\n"
        "  Objective function value at solution: %8.1e\n", rinfo[0]);
      nag_opt_handle_set_get_real(handle, "Dual Variables", 1, 
                                  &nu, u, NAGERR_DEFAULT);
      /* Print gradient vector FDX */
      Vprintf("  Gradient at solution:                 %8.1e %8.1e\n\n",
             u[0]-u[1], u[2]-u[3]);
      /* Print Lagrange multipliers */
      Vprintf("  Estimated Lagrange multipliers: blx   %8.1e %8.1e\n",
             u[0], u[2]);
      Vprintf("  Estimated Lagrange multipliers: bux   %8.1e %8.1e\n\n",
             u[1], u[3]);
    }

    nag_file_close(mon_out, NAGERR_DEFAULT);

    /* Clean up */
    if (handle)
        /* nag_opt_handle_free (e04rzc).
         * Destroy the problem handle and deallocate all the memory used
         */
        nag_opt_handle_free(&handle, NAGERR_DEFAULT);

    return(0);

}