void function_cx_1_func(const real_1d_array &c, const real_1d_array &x, double &func, void *ptr)
{
    func = c[0] + x[0]*c[1];
}

int main(int argc, char **argv)
{

    real_2d_array x = "[[0.],[1.],[2.],[3.]]";
    real_1d_array y = "[1., 2., 4., 8.]";
    real_1d_array c = "[0., 0.]";
    double epsf = 0;
    double epsx = 0.000001;
    ae_int_t maxits = 0;
    ae_int_t info;
    lsfitstate state;
    lsfitreport rep;
    double diffstep = 0.0001;

    lsfitcreatef(x, y, c, diffstep, state);
    lsfitsetcond(state, epsf, epsx, maxits);
    lsfitfit(state, function_cx_1_func);
    lsfitresults(state, info, c, rep);

    for(int i=0; i<2; i++)
   for(int j=0; j<2; j++)
       printf("%f\n", state.c_ptr()->optstate.quadraticmodel.ptr.pp_double[i][j]);

    return 0;

}

real_1d_array y = "[1., 2., 4., 7";