function [logmax,dL]=nolin1d_grad_optag_schur_stable_c(x,M,y)

% NOLIN1D_GRAD_OPTAG_SCHUR_STABLE_C     This optimisation script uses only 
% GRADIENT
% information to perform Gaussian Process prior model optimisation for
% time-series data, i.e. one-dimensional data, [M,y] with M as the input
% and y as the output. Hessian is approximated by finite-differencing using
% MATLAB operation.
% The hyperparamters, g represents the length-scale, and V, the noise
% variance.
%
% This is a fast algorithm using Generalised Schur algorithm to perform
% computation task, impossible by standard MATLAB command. The Schur
% algorithm is written in Mex-C coding.
%
% [logmax,dL] = NOLIN1D_GRAD_OPTAG_SCHUR_STABLE_C(x,M,y)
%
% Inputs:   M   : one-dimensional input column vector
%           y   : one-dimensional output column vector
%           x   : hyperparameter for GP, x = [g V] for training
%
% Outputs:  logmax  = log likelihood function value
%           dL      = gradient vector value
%
%
% (C) Gaussian Process Schur Toolbox 1.0
% (C) Copyright 2005-2007, Keith Neo <keith.neo@nuim.ie>
%                          http://www.hamilton.ie

N0=size(M,1);
g=exp(x(1));
V=exp(x(2));
hyp=[g V];
G=get_eigen_block_trace_a(M,hyp);
[trdQ,dQb]=schur_invbfv3b(G,y);
clear G;
G=get_eigen_block_a(M,hyp);
[logdet,yQy,trQ,yQQy,invb]=schur_invbfv3(G,y);
clear G;
logmax=logdet+N0*log(yQy);
Ny=N0/yQy;
dL=zeros(2,1);
dL(1,1)=trdQ-Ny*(invb'*dQb);
dL(2,1)=V*(trQ-Ny*yQQy);