#include <ruby.h>
#include <math.h>
#include <stdio.h>
#include <narray.h>
#include "ssl2.h"


extern void rqdr_(float*,float*,float*,float*,int*);
extern void cqdr_(scomplex*,scomplex*,scomplex*,scomplex*,int*);
extern void lowp_(float*,int*,float*,int*);
extern void rjetr_(float*,int*,float*,float*,int*);
extern void cjart_(scomplex*,int*,scomplex*,int*);
extern void tsd1_(float*,float*,float(*)(float*),float*,float*,int*);
extern void tsdm_(float*,float(*)(float*),int*,float*,float*,int*,int*);
/* extern void ctsdm_(float*,float*(*)(float*),int*,float*,float*,int*,int*); */
extern void nolbr_(float*,int*,float(*)(float*),float*,float*,float*,int*,float*,float*,int*);

static VALUE
rb_rqdr(self, a0, a1, a2)
     VALUE self, a0, a1, a2;
{
  float ca0, ca1, ca2;
  float *cz;
  int icon;

  VALUE z;
  int shape[2];


  ca0 = (float) NUM2DBL(a0);
  ca1 = (float) NUM2DBL(a1);
  ca2 = (float) NUM2DBL(a2);

  cz = ALLOC_N(float,2*2);
  rqdr_(&ca0,&ca1,&ca2,cz,&icon);

  shape[0] = 2;
  shape[1] = 2;
  z = ssl2_getary(cz,2,shape);
  free(cz);

  ssl2_error(icon);
  return z;
}

static VALUE
rb_cqdr(self,z0,z1,z2)
     VALUE self, z0, z1, z2;
{
  scomplex *cz0, *cz1, *cz2;
  scomplex cz[2];
  int icon;

  VALUE z;
  int shape[1];
  int nn;
    
  cz0 = ssl2_getcaryc(z0,&nn);
  if (nn!=1) rb_raise(rb_eRuntimeError, "z0 != scomplex");
  cz1 = ssl2_getcaryc(z1,&nn);
  if (nn!=1) rb_raise(rb_eRuntimeError, "z1 != scomplex");
  cz2 = ssl2_getcaryc(z2,&nn);
  if (nn!=1) rb_raise(rb_eRuntimeError, "z2 != scomplex");

  cqdr_(cz0,cz1,cz2,cz,&icon);
  free(cz0);
  free(cz1);
  free(cz2);

  shape[0] = 2;
  z = ssl2_getaryc(cz,1,shape);

  ssl2_error(icon);
  return z;
}

static VALUE
rb_lowp(self,a)
     VALUE self, a;
{
  float *ca;
  int n;
  float *cz;
  int icon;

  VALUE z;
  int shape[2];
  int nn;
    
  ca = ssl2_getcary(a,&nn);
  if (ca[0]==0)
    rb_raise(rb_eRuntimeError, "a[0] must be 0");
  n = nn-1;
  if (n<=0||n>5)
    reise(rb_eRuntimeError, "n must be 0<n<=5");

  cz = ALLOC_N(float, n*2);
  lowp_(ca,&n,cz,&icon);
  free(ca);

  shape[0]=2;
  shape[1]=n;
  z = ssl2_getary(cz,2,shape);
  free(cz);

  ssl2_error(icon);
  return z;
}

static VALUE
rb_rjetr(self, a)
     VALUE self, a;
{
  float *ca;
  int n;
  float *cz;
  float *vw;
  int icon;

  VALUE z;
  int shape[2];
  int nn;
    
  ca = ssl2_getcary(a,&nn);
  n = nn-1;
  if (n<1) reise(rb_eRuntimeError, "n must be >= 1");

  cz = ALLOC_N(float, n*2);
  vw = ALLOC_N(float, 6*(n+1));
  rjetr_(ca,&n,cz,vw,&icon);
  free(ca);
  free(vw);
    
  shape[0]=2;
  shape[1]=n;
  z = ssl2_getary(cz,2,shape);
  free(cz);

  ssl2_error(icon);
  return z;
}

static VALUE
rb_cjart(self, za)
     VALUE self, za;
{
  scomplex *cza;
  int n;
  scomplex *cz;
  int icon;

  VALUE z;
  int shape[1];
  int nn;

  cza = ssl2_getcaryc(za,&nn);
  if (cza[0].r==0&&cza[0].i==0) rb_raise(rb_eRuntimeError, "za[0]==[0,0]");
  n = nn-1;
  if (n<1) rb_raise(rb_eRuntimeError, "n<1");

  cz = ALLOC_N(scomplex, n);
  cjart_(cza,&n,cz,&icon);
  free(cza);
    
  shape[0]=n;
  z = ssl2_getaryc(cz,1,shape);
  free(cz);

  ssl2_error(icon);
  return z;
}

static VALUE
rb_tsd1(self,ai,bi,epst)
     VALUE self,ai,bi,epst;
{
  float cai,cbi;
  float cepst;
  float cx;
  int icon;


  cai = (float)NUM2DBL(ai);
  cbi = (float)NUM2DBL(bi);
  if (ssl2_func(&cai)*ssl2_func(&cbi)>0)
    rb_raise(rb_eRuntimeError, "f(a)f(b) > 0");
  cepst = (float)NUM2DBL(epst);
  if (epst<0)
    rb_raise(rb_eRuntimeError, "epst<0");

  tsd1_(&cai,&cbi,&ssl2_func,&cepst,&cx,&icon);

  ssl2_error(icon);
  return rb_float_new((double)cx);
}

static VALUE
rb_tsdm(self,x,isw,eps,eta,m)
     VALUE self,x,isw,eps,eta,m;
{
  float cx;
  int cisw;
  float ceps,ceta;
  int cm;
  int icon;

  cx = (float)NUM2DBL(x);
  cisw = NUM2INT(isw);
  if (cisw<1||cisw>3)
    rb_raise(rb_eRuntimeError, "isw must be 1,2,or3");
  ceps = (float)NUM2DBL(eps);
  ceta = (float)NUM2DBL(eta);
  if (isw==1&&eps<0)
    rb_raise(rb_eRuntimeError, "isw==1 && eps<0");
  if (isw==2&&eta<0)
    rb_raise(rb_eRuntimeError, "isw==2 && eta<0");
  if (isw==3&&eps<0&&eta<0)
    rb_raise(rb_eRuntimeError, "isw==3 && eps<0 && eta<0");
  cm = NUM2INT(m);
  if (m<1)
    rb_raise(rb_eRuntimeError, "m<1");

  tsdm_(&cx,&ssl2_func,&cisw,&ceps,&ceta,&cm,&icon);

  return rb_ary_new3(3,rb_float_new((double)cx),INT2NUM(cm),INT2NUM(icon));
}

static VALUE
rb_nolbr(self,x,epsz,epst,fc,m)
     VALUE self,x,epsz,epst,fc,m;
{
  float *cx;
  int n;
  float cepsz,cepst;
  float cfc;
  int cm;
  float cfnor;
  float *vw;
  int icon;

  int shape[1];

  cx = ssl2_getcary(x,&n);
  if (n<=0)
    rb_raise(rb_eRuntimeError, "n<=0");
  cepsz = (float)NUM2DBL(epsz);
  if (cepsz<0)
    rb_raise(rb_eRuntimeError, "epsz<0");
  cepst = (float)NUM2DBL(epst);
  if (cepst<0)
    rb_raise(rb_eRuntimeError, "epst<0");
  cfc = (float)NUM2DBL(fc);
  if (cfc<=0)
    rb_raise(rb_eRuntimeError, "fc<0");
  cm = NUM2INT(m);
  if (cm<=0)
    rb_raise(rb_eRuntimeError, "cm<=0");

  vw = ALLOC_N(float,n*(n+3));
  nolbr_(cx,&n,&ssl2_func,&cepsz,&cepst,&cfc,&cm,&cfnor,vw,&icon);
  free(vw);

  shape[0]=n;
  x = ssl2_getary(cx,1,shape);
  free(cx);

  return rb_ary_new3(4,x,rb_float_new((double)cfc),INT2NUM(cm),INT2NUM(icon));
}


void init_nonlinear_equations(mSSL2)
     VALUE mSSL2;
{
    rb_define_module_function(mSSL2, "rqdr", rb_rqdr, 3);
    rb_define_module_function(mSSL2, "cqdr", rb_cqdr, 3);
    rb_define_module_function(mSSL2, "lowp", rb_lowp, 1);
    rb_define_module_function(mSSL2, "rjetr", rb_rjetr, 1);
    rb_define_module_function(mSSL2, "cjart", rb_cjart, 1);
    rb_define_module_function(mSSL2, "tsd1", rb_tsd1, 3);
    rb_define_module_function(mSSL2, "tsdm", rb_tsdm, 5);
    rb_define_module_function(mSSL2, "nolbr", rb_nolbr, 5);
}