---

gnn_evaluation.c

Go to the documentation of this file.

/***************************************************************************
 *  @file gnn_evaluation.c
 *  @brief Evaluation Implementation.
 *
 *  @date   : 16-09-03 23:47
 *  @author : Pedro Ortega C. <peortega@dcc.uchile.cl>
 *  Copyright  2003  Pedro Ortega C.
 ****************************************************************************/
/*
 *  This program is free software; you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation; either version 2 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU Library General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program; if not, write to the Free Software
 *  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
 */

/**
 * @brief High-Level \ref gnn_node evaluation interface.
 * @defgroup gnn_evaluation_doc Evaluations.
 * @ingroup  libgnn
 *
 * The \ref gnn_node structure, taken the way it is, provides a very low-level
 * implementation of gradient models. Its programming interface it very
 * versatile, but sometimes, for simpler problems, it can be hard to use.
 *
 * The present module provides a high-level programming interface, where
 * the 3 most frequent tasks have been packed into complex, but easy-to-use
 * routines.
 *
 * The API provides routines for the following tasks:
 * - Evaluation of input patterns: Given the model \f$F(x,w)\f$ implemented
 *   by a particular \ref gnn_node, take a set of input vectors
 *   \f$X = \{x^p \in \mathbb{R}^n, p=0,1,\ldots,P-1\}\f$,
 *   and evaluate all (or some) outputs \f$y^p=F(x^p, \overline{w})\f$
 *   for a fixed parameter vector \f$\overline{w}\f$.
 * - Evaluation of the error and its gradients: Given the model \f$F(x,w)\f$
 *   implemented by a particular \ref gnn_node, an error criterion \f$E(y,t)\f$
 *   implemented by a \ref gnn_criterion and a pattern set (or a subset of)
 *   \f${(x^k, t^k, p^k)|k=0,\ldots,P-1}\f$, evaluate the averaged error
 *   \f$<E>\f$, gradient \f$<\frac{\partial E}{\partial x}>\f$ and
 *   gradient \f$<\frac{\partial E}{\partial w}>\f$.
 * - Line evaluation procedures: Consider the one-dimensional function
 *   \f$g(\alpha) = F(x, w_0 + \alpha d)\f$, where \f$w_0\f$ is a vector
 *   called the \em origin, \f$d\f$ is a vector called the directional vector,
 *   and \f$\alpha\f$ is a scalar. Then, compute the averaged error
 *   \f$<E>\f$ and the derivative along the given line
 *   \f$<\frac{d E}{d \alpha}>
 *   = <\frac{\partial E}{\partial w} \frac{\partial w}{\partial \alpha}>\f$.
 *
 * All the previous evaluation tasks need a special structure
 * containing all the needed buffers in order to perform the required
 * evaluations. There are three types of buffers:
 * - \ref gnn_eval : A buffer for computing the output evaluations.
 * - \ref gnn_grad : A buffer for computing the error and gradients.
 * - \ref gnn_line : A buffer for computing the line evaluation procedures.
 *
 * Before calling the functions, you should first create an appropiate
 * buffer. This is done by calling \ref gnn_eval_new for outputs,
 * \ref gnn_grad_new for gradients and \ref gnn_line_new for line evaluations.
 * These buffers are then given to the evaluation functions. After using
 * them, you should free the buffers with \ref gnn_eval_destroy,
 * \ref gnn_grad_destroy or \ref gnn_line_destroy respectively.
 *
 * Usually, you should call the evaluation functions giving them the
 * needed buffers, and after then access the its fields in order to obtain
 * the result using the provided macros and functions.
 *
 * When an evaluation buffer is destroyed, its associated objects, like
 * the \ref gnn_node or \ref gnn_crit used for its allocation, won't be
 * destroyed. These should be deallocated separately.
 * <b> These three type of buffers do not own the associated \ref gnn_node,
 * \ref gnn_dataset, and all required structures to build them. This
 * means that these have to be destroyed separately. </b>
 *
 *
 *
 * \section gnn_eval_use Using gnn_eval.
 *
 * The \ref gnn_eval evaluation interface takes a given model, and a set
 * of inputs, and computes the corresponding outputs.
 *
 * Example use:
   \code
   gnn_node   *node;
   gnn_input  *input;
   gnn_output *output;
   gnn_eval   *eval;

   // Create the needed objects.
   node   = my_new_custom_node ();
   input  = my_new_custom_input ();
   output = my_new_custom_output ();

   // Allocate the evaluation buffer.
   eval = gnn_eval_new (node, input, output);

   // Compute all estimated outputs: the gnn_node evaluates each
   // input example and puts the results into the gnn_output.
   gnn_eval_all (eval);
   
   // Now, "output" contains the results. Do something with them.
   ...
   
   // Destroy buffer and objects.
   gnn_eval_destroy (eval);
   gnn_node_destroy (node);
   gnn_input_destroy (input);
   gnn_output_destroy (output);
   
   \endcode
 *
 *
 *
 * \section gnn_grad_use Using gnn_grad.
 *
 * The \ref gnn_grad evaluation interface takes a given model, a criterion and
 * a dataset, and computes the weighted sum of errors \f$<E>\f$,
 * the weighted sum of the gradient with respect to the inputs
 * \f$<\frac{\partial E}{\partial x}>\f$,
 * or the weighted sum of the gradients with respect to the parameters
 * \f$<\frac{\partial E}{\partial w}>\f$.
 *
 * When calling one of the evaluation routines (that is, either by calling
 * \ref gnn_grad_pats or \ref gnn_grad_all), you have to specify how
 * far the computation should be performed. That's the purpose of the last
 * flag parameter:
 *
 * - gnnGradE : Compute only the weighted sum of the resulting errors \f$<E>\f$.
 * - gnnGradDx : Compute \f$<E>\f$ \b and \f$<\frac{\partial E}{\partial x}>\f$.
 * - gnnGradDw : Compute \f$<E>\f$, \f$<\frac{\partial E}{\partial x}>\f$ and
 *               \f$<\frac{\partial E}{\partial w}>\f$.
 *
 * The values that haven't been computed shouldn't be accessed, because
 * their values are undefined.
 *
 * Example use:
   \code
   double E;
   gsl_vector *Dx;
   
   gnn_node      *node;
   gnn_criterion *crit;
   gnn_dataset   *data;
   gnn_grad      *grad;

   // Create the needed objects.
   node = my_new_custom_node ();
   crit = my_new_custom_criterion ();
   data = my_new_custom_dataset ();

   // Allocate the gradient evaluation buffer.
   grad = gnn_grad_new (node, crit, data);

   // In this example, we need to compute the gradient with respect
   // to the inputs. So we will execute the evaluation routine with the
   // "gnnGradDx" flag.
   gnn_grad_all (grad, gnnGradDx);

   // We can retrieve the results now by using the macros.
   
   // First, get the weighted error which was also computed.
   E = GNN_GRAD_E (grad);
   
   // Now, get the gradient. First allocate a vector were to store the
   // result, and then use the macro to access the buffer's internal
   // Dx vector.
   Dx = gsl_vector_alloc (GNN_GRAD_INPUT_SIZE (grad));
   gsl_vector_memcpy (Dx, GNN_GRAD_DX (grad));
   
   // Do something else here...
   ...

   // Destroy buffer and objects.
   gnn_grad_destroy (grad);
   gnn_node_destroy (node);
   gnn_criterion_destroy (crit);
   gnn_dataset_destroy (data);

   \endcode
 *
 *
 *
 * \section gnn_line_use Using gnn_line.
 *
 * The \ref gnn_line evaluation interface takes an already built \ref gnn_grad
 * buffer plus a direction, and computes the error or the error derivative
 * along the given line.
 *
 * When calling one of the evaluation routines (that is, either by calling
 * \ref gnn_line_pats or \ref gnn_line_all), you have to specify how
 * far the computation should be performed. That's the purpose of the last
 * flag parameter:
 *
 * - gnnLineE  : Compute only the weighted sum of the resulting error \f$<E>\f$.
 * - gnnLineDE : Compute \f$<E>\f$ \b and \f$<\frac{d E}{d \alpha}>\f$.
 *
 * The values that haven't been computed shouldn't be accessed, because
 * their values are undefined.
 *
 * Example use:
   \code
   double E;
   double DE;
   gsl_vector *d;

   gnn_grad *grad;
   gnn_line *line;

   // Create the needed objects.
   grad = my_new_custom_grad ();
   
   // In this example, we want to compute the error derivative along
   // the first parameter. So we have to prepare the corresponding
   // directional vector.
   d = gsl_vector_calloc (GNN_GRAD_PARAMETER_SIZE (grad));
   gsl_vector_set (d, 0, 1.0);
   
   // Allocate the line evaluation buffer.
   line = gnn_line_new (grad, d);

   // We have to execute the evaluation routine with the "gnnLineDE"
   // flag in order to get the derivative computed.
   gnn_line_all (line, 0.0, gnnLineDE);

   // Get the results by using the macros.
   E  = GNN_LINE_E (line);
   DE = GNN_LINE_DE (line);
   
   // Of course, the evaluation can be performed at another point which
   // can be more convenient for the application.
   gnn_line_all (line, 1.0, gnnLineDE);
   ...

   // Do something else here...
   ...

   // Destroy buffer and objects.
   gnn_line_destroy (line);
   gnn_node_destroy (GNN_GRAD_GET_NODE (grad));
   gnn_criterion_destroy (GNN_GRAD_GET_CRITERION (grad));
   gnn_dataset_destroy (GNN_GRAD_GET_DATASET (grad));
   gnn_grad_destroy (grad);

   \endcode
 *
 */

/******************************************/
/* Include Files                          */
/******************************************/

#include "gnn_evaluation.h"



/******************************************/
/* Public Interface Implementation        */
/******************************************/

/*************************/
/* a) Pattern Evaluation */
/*************************/

/**
 * @brief Builds a new buffer for output evaluation.
 * @ingroup gnn_evaluation_doc
 *
 * This function builds a new \ref gnn_eval buffer structure for computing
 * outputs.
 *
 * @param node A pointer to a \ref gnn_node.
 * @param in   A pointer to a \ref gnn_input input sampler.
 * @param out  A pointer to a \ref gnn_output output device.
 * @return Returns a pointer to a new buffer structure of NULL if failed.
 */
gnn_eval *
gnn_eval_new (gnn_node *node, gnn_input *in, gnn_output *out)
{
    size_t n;
    size_t m;
    size_t l;
    size_t P;
    gnn_eval *eval;
    
    assert (node != NULL);
    assert (in   != NULL);
    assert (out  != NULL);
    
    /* get sizes */
    n = gnn_node_input_get_size (node);
    m = gnn_node_output_get_size (node);
    l = gnn_node_param_get_size (node);
    P = gnn_input_get_size (in);
    
    /* check sizes */
    if (n != gnn_input_sample_get_size (in))
    {
        GSL_ERROR_VAL ("the node's input and the input device's sample size "
                       "should match", GSL_EINVAL, NULL);
    }
    if (gnn_output_is_random_access (out) && m != gnn_output_sample_get_size (out))
    {
        GSL_ERROR_VAL ("the node's and the dataset's output size "
                       "should match", GSL_EINVAL, NULL);
    }
    if (gnn_output_is_random_access (out) && P != gnn_output_get_size (out))
    {
        GSL_ERROR_VAL ("a random access output device should handle the same "
                       "amount of outputs as available inputs",
                       GSL_EINVAL, NULL);
    }

    /* allocate buffer structure */
    eval = (gnn_eval *) malloc (sizeof (*eval));
    if (eval == NULL)
    {
        GSL_ERROR_VAL ("could not allocate memory for gnn_eval",
                       GSL_ENOMEM, NULL);
    }

    /* allocate buffer output vector */
    eval->y = gsl_vector_alloc (m);
    if (eval->y == NULL)
    {
        GSL_ERROR_VAL ("could not allocate memory for gnn_eval buffers",
                       GSL_ENOMEM, NULL);
    }

    /* set fields */
    eval->n      = n;
    eval->m      = m;
    eval->l      = l;
    eval->P      = P;
    eval->node   = node;
    eval->input  = in;
    eval->output = out;
    
    return eval;
}

/**
 * @brief Destroys an evaluation buffer.
 * @ingroup gnn_evaluation_doc
 *
 * @param eval A pointer a \ref gnn_eval buffer.
 */
void
gnn_eval_destroy (gnn_eval *eval)
{
    assert (eval != NULL);
    
    if (eval->y != NULL)
        gsl_vector_free (eval->y);
    free (eval);
}

/**
 * @brief Computes the outputs for a given input set.
 * @ingroup gnn_evaluation_doc
 *
 * This function computes, given a set of input samples, all associated
 * outputs, storing them into the output device provided by the \ref gnn_eval
 * structure. The patterns that are evaluated are those between \f$k=s\f$
 * and \f$k=s+n\f$.
 *
 * @param eval   A pointer to a \ref gnn_eval.
 * @param s      The index of the first input sample to be evaluated.
 * @param n      The amount of samples to be evaluated.
 * @return Returns 0 if suceeded.
 */
int
gnn_eval_pats (gnn_eval *eval, size_t s, size_t n)
{
    size_t k;
    const gsl_vector *x;

    assert (eval != NULL);

    /* compute outputs */
    for (k=s; k<s+n; ++k)
    {
        /* get k-th input sample */
        x = gnn_input_get (eval->input, k);

        /* evaluate */
        gnn_node_evaluate_init (eval->node);
        gnn_node_evaluate_f (eval->node, x, eval->y);
        
        /* store output */
        gnn_output_put (eval->output, k, eval->y);
    }

    return 0;
}

/**
 * @brief Computes the outputs for a given input set.
 * @ingroup gnn_evaluation_doc
 *
 * This function computes, given a set of input samples, all associated
 * outputs, storing them into the output device provided by the \ref gnn_eval
 * structure. The patterns that are evaluated are those between \f$k=s\f$
 * and \f$k=s+n\f$.
 *
 * @param eval   A pointer to a \ref gnn_eval.
 * @param s      The index of the first input sample to be evaluated.
 * @param n      The amount of samples to be evaluated.
 * @return Returns 0 if suceeded.
 */
int
gnn_eval_all (gnn_eval *eval)
{
    return gnn_eval_pats (eval, 0, eval->P);
}



/************************************/
/* b) Error and Gradient Evaluation */
/************************************/

/**
 * @brief Builds a new buffer for gradient evaluation.
 * @ingroup gnn_evaluation_doc
 *
 * This function builds a new \ref gnn_grad buffer structure for computing
 * gradients.
 *
 * @param node A pointer to a \ref gnn_node.
 * @param crit A pointer to a \ref gnn_criterion.
 * @param data A pointer to a \ref gnn_dataset.
 * @return Returns a pointer to a new buffer structure of NULL if failed.
 */
gnn_grad *
gnn_grad_new (gnn_node *node, gnn_criterion *crit, gnn_dataset *data)
{
    size_t n;
    size_t m;
    size_t l;
    size_t P;
    gnn_grad *grad;

    assert (node != NULL);

    /* get sizes */
    n = gnn_node_input_get_size (node);
    m = gnn_node_output_get_size (node);
    l = gnn_node_param_get_size (node);
    P = gnn_dataset_get_size (data);

    /* check sizes */
    if ( n != gnn_dataset_input_get_size (data) )
    {
        GSL_ERROR_VAL ("gnn_node's and gnn_dataset's input sizes should match",
                       GSL_EINVAL, NULL);
    }

    if ( m != gnn_dataset_output_get_size (data) )
    {
        GSL_ERROR_VAL ("gnn_node's and gnn_dataset's output sizes should match",
                       GSL_EINVAL, NULL);
    }

    if ( (crit != NULL) && (m != gnn_criterion_get_size (crit)) )
    {
        GSL_ERROR_VAL ("gnn_node's output size and gnn_criterion's size "
                       "should match", GSL_EINVAL, NULL);
    }

    /* allocate structure */
    grad = (gnn_grad *) malloc (sizeof (*grad));
    if (grad == NULL)
    {
        GSL_ERROR_VAL ("could allocate memory for gnn_grad",
                       GSL_ENOMEM, NULL);
    }

    /* build buffers */
    grad->y   = gsl_vector_alloc (m);
    grad->w   = gsl_vector_alloc (l);
    grad->dy  = gsl_vector_alloc (m);
    grad->dx  = gsl_vector_alloc (n);
    grad->mdx = gsl_vector_alloc (n);
    grad->dw  = gsl_vector_alloc (l);
    grad->mdw = gsl_vector_alloc (l);
    grad->h   = gsl_vector_alloc (l);
    if (    (grad->y   == NULL)
         || (grad->w   == NULL)
         || (grad->dy  == NULL)
         || (grad->dx  == NULL)
         || (grad->mdx == NULL)
         || (grad->dw  == NULL)
         || (grad->mdw == NULL)
         || (grad->h   == NULL) )
    {
        gnn_grad_destroy (grad);
        GSL_ERROR_VAL ("could allocate memory for gnn_grad's buffers",
                       GSL_ENOMEM, NULL);
    }

    /* set pointers and values */
    grad->node = node;
    grad->crit = crit;
    grad->data = data;
    grad->n    = n;
    grad->m    = m;
    grad->l    = l;
    grad->P    = P;

    return grad;
}

/**
 * @brief Destroys a gradient buffer.
 * @ingroup gnn_evaluation_doc
 *
 * @param grad A pointer a \ref gnn_grad buffer.
 */
void
gnn_grad_destroy (gnn_grad *grad)
{
    assert (grad != NULL);

    if (grad->y   != NULL) gsl_vector_free (grad->y);
    if (grad->w   != NULL) gsl_vector_free (grad->w);
    if (grad->dy  != NULL) gsl_vector_free (grad->dy);
    if (grad->dx  != NULL) gsl_vector_free (grad->dx);
    if (grad->mdx != NULL) gsl_vector_free (grad->mdx);
    if (grad->dw  != NULL) gsl_vector_free (grad->dw);
    if (grad->mdw != NULL) gsl_vector_free (grad->mdw);
    if (grad->h   != NULL) gsl_vector_free (grad->h);
    if (grad != NULL) free (grad);
}


/**
 * @brief Compute the mean cost and gradients.
 * @ingroup gnn_evaluation_doc
 *
 * This function computes many things simultaneously in one pass. It builds
 * estimations of the gradients \f$\frac{\partial E}{\partial x}\f$
 * and \f$\frac{\partial E}{\partial w}\f$ and the cost \f$E\f$ by
 * averaging over the n patterns starting at s. That is, this function
 * computes:
 *
 * \f[ <\frac{\partial E}{\partial x}>
 *     = \frac{1}{\sum_{k=s}^{s+n} p_k}
 *       \sum_{k=s}^{s+n} p_k \frac{\partial E}{\partial x}(x^k, w),
 * \f]
 * \f[ <\frac{\partial E}{\partial w}>
 *     = \frac{1}{\sum_{k=s}^{s+n} p_k}
 *       \sum_{k=s}^{s+n} p_k \frac{\partial E}{\partial w}(x^k, w)
 * \f]
 * and
 * \f[ <E>
 *     = \frac{1}{\sum_{k=s}^{s+n} p_k}
 *       \sum_{k=s}^{s+n} p_k E(f(x^k, w), t^k)
 * \f]
 * where \f$p_k\f$ is the \f$k\f$-th pattern's weight,
 * \f$ \frac{\partial E}{\partial x}(x^k, w) \f$ is the gradient with
 * respect to the inputs evaluated at the \f$k\f$-th pattern,
 * \f$ \frac{\partial E}{\partial w}(x^k, w) \f$ is the gradient with
 * respect to the parameters evaluated at the \f$k\f$-th pattern, and
 * \f$ E(f(x^k, w), t^k) \f$ is the cost/error for the \f$k\f$-th output.
 *
 * In order to obtain the results, you should tell the function how far
 * it should process the patterns, using a special flag. The evaluation flag
 * can be:
 * - gnnGradE : only compute mean error \f$<E>\f$.
 * - gnnGradDx : compute \f$<E>\f$ and \f$<\frac{\partial E}{\partial x}>\f$.
 * - gnnGradDw : compute \f$<E>\f$ and \f$<\frac{\partial E}{\partial x}>\f$
 *               and \f$<\frac{\partial E}{\partial w}>\f$
 *
 * The computed values can be retrieved directly from the \ref gnn_grad buffer:
 * \code
 * // compute all
 * gnn_grad_pats (g, gnnGradDw, 10, 20);
 *
 * // now, get them
 * e  = GNN_GRAD_E (g);       // e is a double
 * dx = GNN_GRAD_DX (g);      // dx is a gsl_vector *
 * dw = GNN_GRAD_DW (g);      // dw is a gsl_vector *
 *
 * sp = GNN_GRAD_SUMP (g);    // sp is a double
 * se = GNN_GRAD_SUME (g);    // se is a double
 * \endcode
 *
 * The starting index s should be within the valid bounds. The size n of the
 * batch will be adapted if s + n exceeds the valid bounds. That is, e.g. if
 * there are 10 patterns (and so the index of the last pattern to be considered
 * is 9), and s=8 and n=4, then the last pattern in the averaging will be
 * 9, although 8+4=12.
 *
 * @param grad A pointer to a \ref gnn_grad.
 * @param flag The evaluation flag.
 * @param s    The index of the first pattern in the dataset to be considered
 *             in the averaging.
 * @param n    The amount of patterns to be considered.
 * @return Returns 0 if suceeded.
 */
int
gnn_grad_pats (gnn_grad *grad, gnn_grad_eval flag, size_t s, size_t n)
{
    size_t k;
    size_t P;
    double p;
    gsl_vector *x;
    gsl_vector *t;

    assert (grad != NULL);

    /* check batch bounds */
    if (s < 0 || s >= grad->P)
    {
        GSL_ERROR ("invalid starting index", GSL_EINVAL);
    }
    if (n < 1)
    {
        GSL_ERROR ("the number of patterns in the batch "
                   "should be stricly positive",
                   GSL_EINVAL);
    }
    
    /* adjust batch size */
    n = (n+s <= grad->P)? n : grad->P - s;

    /* clear gradients */
    if (flag >= gnnGradDx) gsl_vector_set_zero (grad->mdx);
    if (flag >= gnnGradDw) gsl_vector_set_zero (grad->mdw);

    /* compute gradient */
    grad->mp = 0.0;
    grad->me = 0.0;
    for (k=s; k<s+n; ++k)
    {
        /* get pattern */
        gnn_dataset_get (grad->data, k, &x, &t, &p);

        /* compute gradient */
        gnn_node_evaluate_init (grad->node);
        gnn_node_evaluate_f (grad->node, x, grad->y);
        grad->e = gnn_criterion_evaluate_e (grad->crit, grad->y, t);
        gnn_criterion_evaluate_dy (grad->crit, grad->dy);
        gnn_node_evaluate_dx (grad->node, grad->dy, grad->dx);
        gnn_node_evaluate_dw (grad->node, grad->dw);

        /* scale and sum gradients */
        if (flag >= gnnGradDx)
        {
            gsl_vector_scale (grad->dx, p);
            gsl_vector_add (grad->mdx, grad->dx);
        }
        if (flag >= gnnGradDw)
        {
            gsl_vector_scale (grad->dw, p);
            gsl_vector_add (grad->mdw, grad->dw);
        }

        /* keep track of the weight and cost sum */
        grad->mp += p;
        grad->me += p * grad->e;
    }

    /* scale gradients and cost */
    if (flag >= gnnGradDx)
        gsl_vector_scale (grad->mdx, 1 / grad->mp);

    if (flag >= gnnGradDw)
        gsl_vector_scale (grad->mdw, 1 / grad->mp);

    grad->me = grad->me / grad->mp;
    
    return 0;
}

/**
 * @brief Compute the mean cost and gradients.
 * @ingroup gnn_evaluation_doc
 *
 * This function behaves like the \ref gnn_grad_pats function, but
 * it computes the gradients and the error over the whole dataset. That is,
 * it computes:
 *
 * \f[ <\frac{\partial E}{\partial x}>
 *     = \frac{1}{\sum_{k=0}^{P-1} p_k}
 *       \sum_{k=0}^{P-1} p_k \frac{\partial E}{\partial x}(x^k, w),
 * \f]
 * \f[ <\frac{\partial E}{\partial w}>
 *     = \frac{1}{\sum_{k=0}^{P-1} p_k}
 *       \sum_{k=0}^{P-1} p_k \frac{\partial E}{\partial w}(x^k, w)
 * \f]
 * and
 * \f[ <E>
 *     = \frac{1}{\sum_{k=0}^{P-1} p_k}
 *       \sum_{k=0}^{P-1} p_k E(f(x^k, w), t^k)
 * \f]
 * where \f$p_k\f$ is the \f$k\f$-th pattern's weight,
 * \f$ \frac{\partial E}{\partial x}(x^k, w) \f$ is the gradient with
 * respect to the inputs evaluated at the \f$k\f$-th pattern,
 * \f$ \frac{\partial E}{\partial w}(x^k, w) \f$ is the gradient with
 * respect to the parameters evaluated at the \f$k\f$-th pattern, and
 * \f$ E(f(x^k, w), t^k) \f$ is the cost/error for the \f$k\f$-th output.
 *
 * Recall that \f$P\f$ is the number of patterns contained in the dataset.
 *
 * The evaluation flag can be:
 * - gnnGradE : only compute mean error \f$<E>\f$.
 * - gnnGradDx : compute \f$<E>\f$ and \f$<\frac{\partial E}{\partial x}>\f$.
 * - gnnGradDw : compute \f$<E>\f$ and \f$<\frac{\partial E}{\partial x}>\f$
 *               and \f$<\frac{\partial E}{\partial w}>\f$
 *
 * The computed values can be retrieved directly from the \ref gnn_grad buffer:
 * \code
 * // compute all
 * gnn_grad_all (g, gnnGradDw);
 *
 * // now, get them
 * e  = GNN_GRAD_E (g);       // e is a double
 * dx = GNN_GRAD_DX (g);      // dx is a gsl_vector *
 * dw = GNN_GRAD_DW (g);      // dw is a gsl_vector *
 *
 * sp = GNN_GRAD_SUMP (g);    // sp is a double
 * se = GNN_GRAD_SUME (g);    // se is a double
 * \endcode
 *
 * @param grad A pointer to a \ref gnn_grad.
 * @param flag The evaluation flag.
 * @return Returns 0 if suceeded.
 */
int
gnn_grad_all (gnn_grad *grad, gnn_grad_eval flag)
{
    assert (grad != NULL);

    return gnn_grad_pats (grad, flag, 0, grad->P);
}



/****************************/
/* c) Error Line Evaluation */
/****************************/

/**
 * @brief Returns a new buffer structure for line evaluations.
 * @ingroup gnn_evaluation_doc
 *
 * This function builds a new buffer structure for line evaluations. The
 * direction is given by the \a direction vector, and the origin is set at
 * the \ref gnn_node's current parameter vector. If the direction is omitted
 * (NULL is given), then it will be a null vector.
 *
 * @param grad A pointer to a \ref gnn_grad.
 * @param direction A pointer to a gsl_vector.
 * @return Returns a pointer to the new \ref gnn_line buffer or NULL if failed.
 */
gnn_line *
gnn_line_new (gnn_grad *grad, gsl_vector *direction)
{
    gnn_line *line;
    
    assert (grad != NULL);

    /* check sizes */
    if (direction != NULL && grad->l != direction->size)
    {
        GSL_ERROR_VAL ("the directional vector's size doesn't match",
                       GSL_EINVAL, NULL);
    }
    
    /* allocate line buffer */
    line = (gnn_line *) malloc (sizeof (*line));
    if (line == NULL)
    {
        GSL_ERROR_VAL ("couldn't allocate memory for gnn_line",
                       GSL_ENOMEM, NULL);
    }

    /* allocate internal buffers and vectors */
    line->w    = gsl_vector_alloc (grad->l);
    line->d    = gsl_vector_alloc (grad->l);
    line->dbuf = gsl_vector_alloc (grad->l);
    if (line->w == NULL || line->d == NULL || line->dbuf == NULL)
    {
        gnn_line_destroy (line);
        GSL_ERROR_VAL ("couldn't allocate memory for gnn_line's buffers",
                       GSL_ENOMEM, NULL);
    }

    /* set remaining fields */
    line->grad       = grad;
    line->derivative = 0.0;
    
    /* set origin */
    gnn_node_param_get (grad->node, line->w);
    
    /* set directional vector */
    if (direction != NULL)
        gsl_vector_memcpy (line->d, direction);
    else
        gsl_vector_set_zero (line->d);

    return line;
}

/**
 * @brief Destroys a given line evaluation buffer.
 * @ingroup gnn_evaluation_doc
 *
 * This function frees the \ref gnn_line's internal buffers and deallocates
 * the structure. The installed \ref gnn_grad error and gradient evaluation
 * buffer won't be destroyed.
 *
 * @param line A pointer to a \ref gnn_line.
 */
void
gnn_line_destroy (gnn_line *line)
{
    assert (line != NULL);
    
    if (line->w != NULL)
        gsl_vector_free (line->w);
    if (line->d != NULL)
        gsl_vector_free (line->d);
    if (line->dbuf != NULL)
        gsl_vector_free (line->dbuf);
    free (line);
}

/**
 * @brief Sets a new direction \f$d\f$ for the \ref gnn_line structure.
 * @ingroup gnn_evaluation_doc
 *
 * This function sets a new directional vector \f$d\f$ for the line evaluation
 * buffer. The vector should be of the correct size.
 *
 * @param line A pointer to a \ref gnn_line.
 * @param dir  A pointer to a gsl_vector.
 * @return Returns 0 if succeeded.
 */
int
gnn_line_set_direction (gnn_line *line, const gsl_vector *dir)
{
    assert (line != NULL);
    assert (dir != NULL);

    /* check size */
    if (line->grad->l != dir->size)
    {
        GSL_ERROR ("the size fo the directional vector is incorrect",
                   GSL_EINVAL);
    }

    /* copy vector */
    gsl_vector_memcpy (line->d, dir);

    return 0;
}

/**
 * @brief Gets the installed directional vector \f$d\f$.
 * @ingroup gnn_evaluation_doc
 *
 * This function returns a pointer to the currently used directional buffer
 * \f$d\f$ used by the \ref gnn_line buffer for its line evaluations.
 *
 * @param line A pointer to a \ref gnn_line.
 * @return Returns a pointer to the directional vector.
 */
const gsl_vector *
gnn_line_get_direction (gnn_line *line)
{
    assert (line != NULL);
    return line->d;
}

/**
 * @brief Sets a new origin \f$w_0\f$ for the \ref gnn_line structure.
 * @ingroup gnn_evaluation_doc
 *
 * This function sets a new origin vector \f$w_0\f$ for the line evaluation
 * buffer. The vector should be of the correct size.
 *
 * @param line A pointer to a \ref gnn_line.
 * @param origin A pointer to a gsl_vector.
 * @return Returns 0 if succeeded.
 */
int
gnn_line_set_origin (gnn_line *line, const gsl_vector *origin)
{
    assert (line != NULL);
    assert (origin != NULL);

    /* check size */
    if (line->grad->l != origin->size)
    {
        GSL_ERROR ("the size fo the origin vector is incorrect",
                   GSL_EINVAL);
    }

    /* copy vector */
    gsl_vector_memcpy (line->w, origin);

    return 0;
}

/**
 * @brief Gets the installed origin vector \f$w_0\f$.
 * @ingroup gnn_evaluation_doc
 *
 * This function returns a pointer to the currently used origin vector
 * \f$w_0\f$ by the \ref gnn_line for its line evaluations.
 *
 * @param line A pointer to a \ref gnn_line.
 * @return Returns a pointer to the origin vector.
 */
const gsl_vector *
gnn_line_get_origin (gnn_line *line)
{
    assert (line != NULL);
    return line->w;
}

/**
 * @brief Gets the installed \ref gnn_grad evaluation buffer.
 * @ingroup gnn_evaluation_doc
 *
 * This function returns a buffer for the installed \ref gnn_grad error and
 * gradients evaluation buffer used by the current \ref gnn_line structure
 * for its line evaluation.
 *
 * @param line A pointer to a \ref gnn_line.
 * @return Returns a pointer to the internal \ref gnn_grad structure.
 */
gnn_grad *
gnn_line_get_grad (gnn_line *line)
{
    assert (line != NULL);
    return line->grad;
}

/**
 * @brief Compute mean cost and gradients along a direction for a minibatch.
 * @ingroup gnn_evaluation_doc
 *
 * This function behaves like the \ref gnn_grad_pats function, but
 * it computes the error and its derivative along a given line. That is,
 * it performs all evaluations at
 *
 * \f[ \overline{w} = w_0 + \alpha d \f]
 *
 * where \f$w_0\f$ are the node's current parameters, \f$\alpha\f$ is a scalar
 * and \f$d\f$ is a given direction, treating the function as it where just
 * a one-dimensional function. The patterns starting at \a s and ending
 * at \a s+n are used for the averaging.
 *
 * The flag \a flag tells the function if it should compute only the error
 * (\c gnnLineE) or the error with its derivative along the line (\c gnnLineDE).
 * The approximated derivative is apropiately scaled in order to match
 * the real derivative:
 * \f[ \frac{1}{||d||}\frac{\partial E}{\partial \alpha}
 *     = \frac{1}{||d||}\frac{\partial E}{\partial w}
 *                    \frac{\partial w}{\partial \alpha}
 *     = \frac{1}{||d||}\frac{\partial E}{\partial w}
 * \f]
 *
 * The results can be recovered by the macros \ref GNN_LINE_E and
 * \ref GNN_LINE_DE:
 *
 * \code
 * gnn_grad *grad;
 * gnn_line *line;
 * double    error;
 * double    derivative;
 *
 * // build grad and line buffers
 * ...
 *
 * // evaluate the patterns 11-30 at 0.5
 * gnn_line_pats (line, 0.5, gnnLineDE, 11, 30);
 *
 * // get results
 * error      = GNN_LINE_E (line);
 * derivative = GNN_LINE_DE (line);
 * \endcode
 *
 * \warning A line evaluation modifies the \ref gnn_node's internal parameters.
 *          In order to recover the original vector, call this function with
 *          \f$\alpha = 0\f$.
 *
 * @param line  A pointer to a \ref gnn_line.
 * @param alpha The scalar \f$\alpha\f$.
 * @param flag  The evaluation flag.
 * @param s     The first pattern index \f$s\f$.
 * @param n     The size of the minibatch \f$n\f$.
 * @return Returns 0 if suceeded.
 */
int
gnn_line_pats (gnn_line *line,
                 double alpha, gnn_line_eval flag, size_t s, size_t n)
{
    assert (line != NULL);

    /* compute new parameter vector (using the correct directional vector */
    gsl_vector_memcpy (line->dbuf, line->d);
    gsl_vector_scale (line->dbuf, alpha);
    gsl_vector_add (line->dbuf, line->w);

    /* set the new parameter vector */
    gnn_node_param_set (line->grad->node, line->dbuf);

    /* evaluate */
    if (flag == gnnLineE)
    {
        /* compute error */
        gnn_grad_pats (line->grad, gnnGradE, s, n);
        line->error = GNN_GRAD_E (line->grad);
    }
    else if (flag == gnnLineDE)
    {
        double dE;
        double nrmd;
        
        /* compute error */
        gnn_grad_pats (line->grad, gnnGradDw, s, n);
        line->error = GNN_GRAD_E (line->grad);
        
        /* compute derivative */
        gsl_blas_ddot (line->d, GNN_GRAD_DW (line->grad), &dE);
        nrmd = gsl_blas_dnrm2 (line->d);
        line->derivative = dE / nrmd;
    }
    else
        GSL_ERROR ("undefined flag for line operation", GSL_EINVAL);

    return 0;
}

/**
 * @brief Compute mean cost and gradients along a direction for a minibatch.
 * @ingroup gnn_evaluation_doc
 *
 * This function behaves like the \ref gnn_grad_pats function, but
 * it computes the error and its derivative along a given line. That is,
 * it performs all evaluations at
 *
 * \f[ \overline{w} = w_0 + \alpha d \f]
 *
 * where \f$w_0\f$ are the node's current parameters, \f$\alpha\f$ is a scalar
 * and \f$d\f$ is a given direction, treating the function as it where just
 * a one-dimensional function. All patterns are used for this evaluation.
 *
 * The flag \a flag tells the function if it should compute only the error
 * (\c gnnLineE) or the error with its derivative along the line (\c gnnLineDE).
 * The approximated derivative is apropiately scaled in order to match
 * the real derivative:
 * \f[ \frac{1}{||d||}\frac{\partial E}{\partial \alpha}
 *     = \frac{1}{||d||}\frac{\partial E}{\partial w}
 *                      \frac{\partial w}{\partial \alpha}
 *     = \frac{1}{||d||}\frac{\partial E}{\partial w}
 * \f]
 *
 * The results can be recovered by the macros \ref GNN_LINE_E and
 * \ref GNN_LINE_DE:
 *
 * \code
 * gnn_grad *grad;
 * gnn_line *line;
 * double    error;
 *
 * // build grad and line buffers
 * ...
 *
 * // evaluate all patterns at 1.2
 * gnn_line_all (line, 1.2, gnnLineE);
 *
 * // get results
 * error      = GNN_LINE_E (line);
 * \endcode
 *
 * \warning A line evaluation modifies the \ref gnn_node's internal parameters.
 *          In order to recover the original vector, call this function with
 *          \f$\alpha = 0\f$.
 *
 * @param line  A pointer to a \ref gnn_line.
 * @param alpha The scalar \f$\alpha\f$.
 * @param flag  The evaluation flag.
 * @return Returns 0 if suceeded.
 */
int
gnn_line_all (gnn_line *line, double alpha, gnn_line_eval flag)
{
    return gnn_line_pats (line, alpha, flag, 0, line->grad->P);
}

---