From df370e24e912011721ed16099d9caab199dec15c Mon Sep 17 00:00:00 2001
From: Anton Khirnov <anton@khirnov.net>
Date: Wed, 3 Oct 2018 15:20:03 +0200
Subject: Export the solver functionality as an external API.

To be used for initializing the multigrid solver.
---
 src/maximal_slicing_axi.c | 27 +++++++++++++++------------
 src/ms_solve.c            |  5 +++--
 src/ms_solve.h            |  2 +-
 3 files changed, 19 insertions(+), 15 deletions(-)

diff --git a/src/maximal_slicing_axi.c b/src/maximal_slicing_axi.c
index 5e7f07d..4ed0acc 100644
--- a/src/maximal_slicing_axi.c
+++ b/src/maximal_slicing_axi.c
@@ -198,22 +198,20 @@ static CoordPatch *get_coord_patch(MaximalSlicingContext *ms, int level)
     return cp;
 }
 
-static int context_init(cGH *cctkGH, MaximalSlicingContext **ctx)
+int msa_context_init(cGH *gh, MaximalSlicingContext **ctx,
+                     int basis_order0, int basis_order1, double outer_bound, double filter_power)
 {
     MaximalSlicingContext *ms;
     int ret;
 
-    DECLARE_CCTK_ARGUMENTS;
-    DECLARE_CCTK_PARAMETERS;
-
     ms = calloc(1, sizeof(*ms));
     if (!ms)
         return -ENOMEM;
 
-    ms->gh = cctkGH;
+    ms->gh = gh;
 
-    ret = ms_solver_init(&ms->solver, cctkGH, basis_order_r, basis_order_z,
-                         scale_factor, filter_power, 0.0);
+    ret = ms_solver_init(&ms->solver, gh, basis_order0, basis_order1,
+                         outer_bound, filter_power, 0.0);
     if (ret < 0)
         return ret;
 
@@ -222,7 +220,7 @@ static int context_init(cGH *cctkGH, MaximalSlicingContext **ctx)
     return 0;
 }
 
-static void context_free(MaximalSlicingContext **pms)
+void msa_context_free(MaximalSlicingContext **pms)
 {
     MaximalSlicingContext *ms = *pms;
 
@@ -294,6 +292,11 @@ void msa_lapse_eval(MaximalSlicingContext *ms, double *dst)
     //            1.0, alp, cctk_lsh[0] * cctk_lsh[1]);
 }
 
+void msa_lapse_solve(MaximalSlicingContext *ms)
+{
+    ms_solver_solve(ms->solver);
+}
+
 void maximal_slicing_axi_initlapse(CCTK_ARGUMENTS)
 {
     MaximalSlicingContext *ms;
@@ -303,7 +306,7 @@ void maximal_slicing_axi_initlapse(CCTK_ARGUMENTS)
 
     int i, ret;
 
-    context_init(cctkGH, &ms);
+    msa_context_init(cctkGH, &ms, basis_order_r, basis_order_z, 64.0, filter_power);
 
     ms_solver_solve(ms->solver);
 
@@ -311,7 +314,7 @@ void maximal_slicing_axi_initlapse(CCTK_ARGUMENTS)
 
     memcpy(alp, alpha, sizeof(*alp) * cctk_lsh[0] * cctk_lsh[2]);
 
-    context_free(&ms);
+    msa_context_free(&ms);
 }
 
 void maximal_slicing_axi(CCTK_ARGUMENTS)
@@ -329,10 +332,10 @@ void maximal_slicing_axi(CCTK_ARGUMENTS)
 
     /* on the first run, init the solver */
     if (!ms)
-        context_init(cctkGH, &ms);
+        msa_context_init(cctkGH, &ms, basis_order_r, basis_order_z, 64.0, filter_power);
 
     CCTK_TimerStart("MaximalSlicingAxi_Solve");
-    ms_solver_solve(ms->solver);
+    msa_lapse_solve(ms);
     CCTK_TimerStop("MaximalSlicingAxi_Solve");
 
     if (export_coeffs)
diff --git a/src/ms_solve.c b/src/ms_solve.c
index 90502be..1bf87a2 100644
--- a/src/ms_solve.c
+++ b/src/ms_solve.c
@@ -19,6 +19,7 @@
 #include "common.h"
 
 #include <errno.h>
+#include <limits.h>
 #include <math.h>
 #include <stdint.h>
 #include <stdio.h>
@@ -492,7 +493,7 @@ fail:
 int ms_solver_init(MSSolver **pctx,
                     cGH *cctkGH,
                     int basis_order_r, int basis_order_z,
-                    double sf, double filter_power, double input_filter_power)
+                    double outer_bound, double filter_power, double input_filter_power)
 {
     MSSolver *ctx;
     MSSolverPriv *s;
@@ -535,7 +536,7 @@ int ms_solver_init(MSSolver **pctx,
 
     //FIXME
     scale_factor = 1.0;
-    scale_factor = (64.0 / ctx->basis[0]->colloc_point(s->colloc_grid_order[0], ctx->nb_colloc_points[0] - 1));
+    scale_factor = (outer_bound / ctx->basis[0]->colloc_point(s->colloc_grid_order[0], ctx->nb_colloc_points[0] - 1));
     fprintf(stderr, "scale factor %16.16g\n", scale_factor);
 
     init_opencl(ctx);
diff --git a/src/ms_solve.h b/src/ms_solve.h
index 2ebf1e1..c0c3e65 100644
--- a/src/ms_solve.h
+++ b/src/ms_solve.h
@@ -41,7 +41,7 @@ typedef struct MSSolver {
 int ms_solver_init(MSSolver **ctx,
                    cGH *cctkGH,
                    int basis_order_r, int basis_order_z,
-                   double sf, double filter_power, double input_filter_power);
+                   double outer_bound, double filter_power, double input_filter_power);
 
 void ms_solver_free(MSSolver **ctx);
 
-- 
cgit v1.2.3