add funcs

src/init.c

@@ -100,8 +100,45 @@ double chi(double b){
   return find_root_5ord(y_previous[0], y[0], y_previous[1], y[1], phi_previous, phi);
 }
 
-static void refine_interval(double xl, double xr, double yl, double yr, SampleData *sample){}
+void sample_push(SampleData *sample, double x, double y){
+  if (sample->size + 1 >= sample->capacity) {
+    sample->x = realloc(sample->x, sizeof(double) * (sample->capacity)*2);
+    sample->y = realloc(sample->y, sizeof(double) * (sample->capacity)*2);
+  }
+  sample->x[sample->size] = x;
+  sample->y[sample->size] = y;
+  sample->size += 1;
+}
+
+static void refine_interval(double xl, double xr, double yl, double yr, double rela_err_limit, SampleData *sample){
+  double xm = (xl + yl)/2;
+  double ym = chi(xm);
+  sample_push(sample, xm, ym);
+  if (fabs((ym-(yl+yr)/2)/ym) >= rela_err_limit) {
+    refine_interval(xl, xm, yl, ym, rela_err_limit, sample);
+    refine_interval(xm, xr, ym, yr, rela_err_limit, sample);
+  }
+  return ;
+}
 
 static int init(double bmax, double rela_err_limit){
+  SampleData sample;
+  sample.capacity = 2;
+  sample.size=2;
+  sample.x = malloc(sizeof(double) * 2);
+  sample.y = malloc(sizeof(double) * 2);
+  sample.x[0] = bmin;
+  sample.x[1] = bmax;
+  sample.y[0] = chi(bmin);
+  sample.y[1] = chi(bmax);
+
+  refine_interval(bmin, bmax, sample.y[0], sample.y[1], rela_err_limit, &sample);
+
+  //sort
+
+  //output
+
+  free(sample.x);
+  free(sample.y);
   return 0;
 }