interp_method_nnn.c

Go to the documentation of this file.

// Copyright (c) 2024 The YAC Authors
//
// SPDX-License-Identifier: BSD-3-Clause
 
#ifdef HAVE_CONFIG_H
// Get the definition of the 'restrict' keyword.
#include "config.h"
#endif
 
#include <limits.h>
#include <string.h>
 
#include "interpolation/methods/interp_method_internal.h"
#include "interpolation/methods/interp_method_config_internal.h"
#include "interpolation/methods/interp_method_nnn.h"
#include "interpolation/methods/interp_method_utils.h"
#include "parameter/param_enum.h"
#include "parameter/param_int.h"
#include "parameter/param_dble.h"
#include "parameter/param_distance.h"
#include "parameter/param_struct.h"
#include "instrument_internal.h"
#include "yac_mpi_internal.h"
#include "yac_mpi_common.h"
 
static size_t do_search_nnn(struct interp_method * method,
                            struct yac_interp_grid * interp_grid,
                            size_t * tgt_points, size_t count,
                            struct yac_interp_weights * weights,
                            int * interpolation_complete);
static size_t do_search_1nn(struct interp_method * method,
                            struct yac_interp_grid * interp_grid,
                            size_t * tgt_points, size_t count,
                            struct yac_interp_weights * weights,
                            int * interpolation_complete);
static size_t do_search_zero(struct interp_method * method,
                            struct yac_interp_grid * interp_grid,
                            size_t * tgt_points, size_t count,
                            struct yac_interp_weights * weights,
                            int * interpolation_complete);
static void delete_nnn(struct interp_method * method);
 
static struct interp_method_vtable
  interp_method_nnn_vtable = {
    .do_search = do_search_nnn,
    .delete = delete_nnn},
  interp_method_1nn_vtable = {
    .do_search = do_search_1nn,
    .delete = delete_nnn},
  interp_method_zero_vtable = {
    .do_search = do_search_zero,
    .delete = delete_nnn};
 
struct interp_method_nnn {
 
  struct interp_method_vtable * vtable;
  yac_interp_compute_weights_func compute_weights;
  size_t n;
  size_t n_min;
  double max_search_distance;
  double scale;
};
 
struct interp_method_1nn {
 
  struct interp_method_vtable * vtable;
  void (*interp_weights_add)(
    struct yac_interp_weights * weights, struct remote_points * tgts,
    struct remote_point * srcs);
  double max_search_distance;
};
 
struct interp_method_zero {
 
  struct interp_method_vtable * vtable;
};
 
static size_t do_search_nnn (struct interp_method * method,
                             struct yac_interp_grid * interp_grid,
                             size_t * tgt_points, size_t count,
                             struct yac_interp_weights * weights,
                             int * interpolation_complete) {
 
  if (*interpolation_complete) return 0;
 
  YAC_INSTRUMENT_START(YAC_INSTR_INTERP_SEARCH_NNN);
 
  struct interp_method_nnn * method_nnn = (struct interp_method_nnn *)method;
 
  CHECK_SRC_FIELD_COUNT_SINGLE(interp_grid)
 
  // get coordinates of target points
  yac_coordinate_pointer tgt_coords = xmalloc(count * sizeof(*tgt_coords));
  yac_interp_grid_get_tgt_coordinates(
    interp_grid, tgt_points, count, tgt_coords);
 
  size_t n = method_nnn->n;
  size_t n_min = method_nnn->n_min;
 
  size_t * size_t_buffer = xmalloc((2 * n + 1) * count * sizeof(*size_t_buffer));
  size_t * src_points = size_t_buffer;
  size_t * src_points_reorder = size_t_buffer + n * count;
  size_t * reorder_idx = size_t_buffer + 2 * n * count;
  int * fail_flag = xmalloc(count * sizeof(*fail_flag));
  size_t * num_src_per_tgt = xmalloc(count * sizeof(*num_src_per_tgt));
 
  // search for the n nearest source points
  yac_interp_grid_do_nnn_search_src(
    interp_grid, tgt_coords, count, n, src_points,
    method_nnn->max_search_distance);
 
  size_t result_count = 0;
  size_t total_num_weights = 0;
 
  // count the number of source points found for each target point
  // (the search returns the found source points at the beginning of the
  //  result buffer of each target point and fills up the remaining entries
  //  with SIZE_MAX)
  // the search is considered to be successful, if at least n_min source
  // points were found
  for (size_t i = 0, k = 0; i < count; ++i) {
    size_t curr_num_src = 0;
    for (size_t j = 0; j < n; ++j, ++k)
      if (src_points[k] != SIZE_MAX) ++curr_num_src;
    num_src_per_tgt[i] = curr_num_src;
    int flag = curr_num_src < n_min;
    fail_flag[i] = flag;
    if (!flag) {
      ++result_count;
      total_num_weights += curr_num_src;
    }
    reorder_idx[i] = i;
  }
 
  // sort target points, for which the search was successful to the beginning of
  // the array
  yac_quicksort_index_int_size_t(fail_flag, count, reorder_idx);
  free(fail_flag);
 
  yac_interp_compute_weights_func compute_weights =
    method_nnn->compute_weights;
  double scale = method_nnn->scale;
  double * w = xmalloc(total_num_weights * sizeof(*w));
  size_t * num_weights_per_tgt =
    xmalloc(result_count * sizeof(*num_weights_per_tgt));
  yac_coordinate_pointer src_coord_buffer = xmalloc(n * sizeof(*src_coord_buffer));
  yac_const_coordinate_pointer src_field_coordinates =
    yac_interp_grid_get_src_field_coords(interp_grid, 0);
 
  // compute the weights
  for (size_t i = 0, weight_offset = 0; i < result_count; ++i) {
    size_t curr_reorder_idx = reorder_idx[i];
    size_t curr_num_src = num_src_per_tgt[curr_reorder_idx];
    for (size_t j = 0; j < curr_num_src; ++j) {
      size_t curr_src_point = src_points[curr_reorder_idx * n + j];
      memcpy(src_coord_buffer[j], src_field_coordinates[curr_src_point],
             3 * sizeof(double));
      src_points_reorder[weight_offset + j] = curr_src_point;
    }
 
    compute_weights(
      tgt_coords[curr_reorder_idx],
      (yac_const_coordinate_pointer)src_coord_buffer, curr_num_src,
      w + weight_offset, scale);
    num_weights_per_tgt[i] = curr_num_src;
    weight_offset += curr_num_src;
  }
 
  free(num_src_per_tgt);
 
  // move the non-interpolated target points to the end
  for (size_t i = 0; i < count; ++i) src_points[reorder_idx[i]] = i;
  yac_quicksort_index_size_t_size_t(src_points, count, tgt_points);
 
  free(src_coord_buffer);
  free(tgt_coords);
 
  struct remote_point * srcs =
    yac_interp_grid_get_src_remote_points(
      interp_grid, 0, src_points_reorder, total_num_weights);
 
  struct remote_points tgts = {
    .data =
      yac_interp_grid_get_tgt_remote_points(
        interp_grid, tgt_points, result_count),
    .count = result_count};
 
  // store weights
  yac_interp_weights_add_wsum(
    weights, &tgts, num_weights_per_tgt, srcs, w);
 
  free(size_t_buffer);
  free(tgts.data);
  free(srcs);
  free(w);
  free(num_weights_per_tgt);
 
  YAC_INSTRUMENT_STOP(YAC_INSTR_INTERP_SEARCH_NNN);
 
  return result_count;
}
 
static size_t do_search_1nn(struct interp_method * method,
                            struct yac_interp_grid * interp_grid,
                            size_t * tgt_points, size_t count,
                            struct yac_interp_weights * weights,
                            int * interpolation_complete) {
 
  if (*interpolation_complete) return 0;
 
  YAC_INSTRUMENT_START(YAC_INSTR_INTERP_SEARCH_1NN);
 
  CHECK_SRC_FIELD_COUNT_SINGLE(interp_grid)
 
  // get coordinates of target points
  yac_coordinate_pointer tgt_coords = xmalloc(count * sizeof(*tgt_coords));
  yac_interp_grid_get_tgt_coordinates(
    interp_grid, tgt_points, count, tgt_coords);
 
  size_t * src_points = xmalloc(count * sizeof(*src_points));
 
  // search for the nearest source points
  yac_interp_grid_do_nnn_search_src(
    interp_grid, tgt_coords, count, 1, src_points,
    ((struct interp_method_1nn *)method)->max_search_distance);
 
  // move the non-interpolated target points to the end
  // (source points == SIZE_MAX if no source point was found)
  yac_quicksort_index_size_t_size_t(src_points, count, tgt_points);
  free(tgt_coords);
 
  // count the number of target points for which we found a source point
  size_t result_count = 0;
  for (; (result_count < count) && (src_points[result_count] != SIZE_MAX);
       ++result_count);
 
  struct remote_point * srcs =
    yac_interp_grid_get_src_remote_points(
      interp_grid, 0, src_points, result_count);
 
  struct remote_points tgts = {
    .data =
      yac_interp_grid_get_tgt_remote_points(
        interp_grid, tgt_points, result_count),
    .count = result_count};
 
  // store weights
  ((struct interp_method_1nn *)method)->interp_weights_add(
    weights, &tgts, srcs);
 
  free(src_points);
  free(tgts.data);
  free(srcs);
 
  YAC_INSTRUMENT_STOP(YAC_INSTR_INTERP_SEARCH_1NN);
 
  return result_count;
}
 
static void interp_weights_add_zero(
    struct yac_interp_weights * weights, struct remote_points * tgts,
    struct remote_point src) {
 
  size_t * num_src_per_tgt = xmalloc(tgts->count * sizeof(*num_src_per_tgt));
  double * w = xmalloc(tgts->count * sizeof(*w));
  struct remote_point * srcs = xmalloc(tgts->count * sizeof(*srcs));
 
  for (size_t i = 0; i < tgts->count; ++i) {
    num_src_per_tgt[i] = 1;
    w[i] = 0.0;
    srcs[i] = src;
  }
 
  yac_interp_weights_add_wsum(weights, tgts, num_src_per_tgt, srcs, w);
 
  free(srcs);
  free(w);
  free(num_src_per_tgt);
}
 
static size_t do_search_zero(struct interp_method * method,
                             struct yac_interp_grid * interp_grid,
                             size_t * tgt_points, size_t count,
                             struct yac_interp_weights * weights,
                             int * interpolation_complete) {
 
  if (*interpolation_complete) return 0;
 
  YAC_INSTRUMENT_START(YAC_INSTR_INTERP_SEARCH_ZERO);
 
  CHECK_SRC_FIELD_COUNT_SINGLE(interp_grid)
 
  UNUSED(method);
 
  // get global ids of all valid local source points
  size_t * src_points, src_count;
  yac_interp_grid_get_src_points(interp_grid, 0, &src_points, &src_count);
  yac_int * src_global_ids = xmalloc(src_count * sizeof(*src_global_ids));
  if (src_count > 0)
    yac_interp_grid_get_src_global_ids(
      interp_grid, src_points, src_count, 0, src_global_ids);
  free(src_points);
 
  // determine smallest global id of all unmasked source points in the
  // local data
  yac_int min_src_global_id = YAC_INT_MAX;
  for (size_t i = 0; i < src_count; ++i)
    if (src_global_ids[i] < min_src_global_id)
      min_src_global_id = src_global_ids[i];
  free(src_global_ids);
 
  // determine smallest global id of all unmasked source points in the
  // global data
  MPI_Comm comm = yac_interp_grid_get_MPI_Comm(interp_grid);
  yac_mpi_call(
    MPI_Allreduce(
      MPI_IN_PLACE, &min_src_global_id, 1, yac_int_dt, MPI_MIN, comm), comm);
  YAC_ASSERT(
    min_src_global_id != YAC_INT_MAX, "unable to find unmasked source points");
 
  // get the local id of the source point with the smallest global ids
  // (add it to the local data if it does not yet exist there)
  size_t min_src_local_id;
  yac_interp_grid_src_global_to_local(
    interp_grid, 0, &min_src_global_id, 1, &min_src_local_id);
 
  struct remote_point * src =
    yac_interp_grid_get_src_remote_points(
      interp_grid, 0, &min_src_local_id, 1);
 
  struct remote_points tgts = {
    .data =
      yac_interp_grid_get_tgt_remote_points(
        interp_grid, tgt_points, count),
    .count = count};
 
  // store weights
  interp_weights_add_zero(weights, &tgts, *src);
 
  free(tgts.data);
  free(src);
 
  YAC_INSTRUMENT_STOP(YAC_INSTR_INTERP_SEARCH_ZERO);
 
  return count;
}
 
// resolves the configured minimum number of source points
// (a value of 0 is interpreted as "use n")
// (the validity of n_min is checked in yac_interp_method_nnn_new)
static inline int get_nnn_n_min(struct yac_interp_method_nnn_config config) {
 
  return (config.n_min <= 0)?config.n:config.n_min;
}
 
struct interp_method * yac_interp_method_nnn_new(
  struct yac_interp_method_nnn_config config) {
 
  YAC_ASSERT_F(
    (config.max_search_distance >= 0.0) &&
    (config.max_search_distance <= M_PI),
    "unsupported value for max_search_distance (%lf), has to be in [0.0,PI]",
    config.max_search_distance)
 
  YAC_ASSERT_F(
    config.n >= 1, "unsupported value for n (%d), has to be >= 1", config.n)
 
  YAC_ASSERT_F(
    (config.n_min >= 0) && (config.n_min <= config.n),
    "unsupported value for n_min (%d), has to be in [0,n] (n: %d)",
    config.n_min, config.n)
 
  if (config.type == YAC_INTERP_NNN_ZERO) {
 
    YAC_ASSERT(
      config.n == 1, "n != 1 does not make sense for zero")
    YAC_ASSERT(
      config.max_search_distance == 0.0,
      "max_search_distance != 0.0 does not make sense for zero")
 
    struct interp_method_zero * method = xmalloc(1 * sizeof(*method));
    method->vtable = &interp_method_zero_vtable;
    return (struct interp_method*)method;
  }
 
  double max_search_distance =
    (config.max_search_distance == 0.0)?M_PI:config.max_search_distance;
 
  if (config.n == 1) {
 
    struct interp_method_1nn * method = xmalloc(1 * sizeof(*method));
    method->vtable = &interp_method_1nn_vtable;
    method->interp_weights_add = yac_interp_weights_add_direct;
    method->max_search_distance = max_search_distance;
    return (struct interp_method*)method;
  }
 
  struct interp_method_nnn * method = xmalloc(1 * sizeof(*method));
 
  method->vtable = &interp_method_nnn_vtable;
 
  method->max_search_distance = max_search_distance;
  method->n = (size_t)config.n;
  method->n_min = (size_t)get_nnn_n_min(config);
 
  switch (config.type) {
    YAC_UNREACHABLE_DEFAULT("invalid n-nearest-neighbour weight type")
    case(YAC_INTERP_NNN_AVG):
      YAC_ASSERT(
        config.n >= 1, "n needs to be at least 1 for distance weighted average")
      method->scale = 0.0;
      method->compute_weights = yac_interp_compute_weights_avg;
      break;
    case(YAC_INTERP_NNN_DIST):
      YAC_ASSERT(
        config.n >= 2, "n needs to be at least 2 for distance weighted average")
      method->scale = 0.0;
      method->compute_weights = yac_interp_compute_weights_dist;
      break;
    case(YAC_INTERP_NNN_GAUSS):
      YAC_ASSERT(
        config.n >= 2, "n needs to be at least 2 for gauss")
      method->scale = config.data.gauss_scale;
      method->compute_weights = yac_interp_compute_weights_gauss;
      break;
    case(YAC_INTERP_NNN_RBF):
      YAC_ASSERT(
        config.n >= 2, "n needs to be at least 2 for RBF")
      method->scale = config.data.rbf_scale;
      method->compute_weights = yac_interp_compute_weights_rbf;
      break;
  };
 
  return (struct interp_method*)method;
}
 
static void delete_nnn(struct interp_method * method) {
  free(method);
}
 
// --- Modular config struct and vtable for NNN method ---
 
struct yac_interp_method_config_nnn {
  struct yac_interp_method_config base;
  struct yac_interp_method_nnn_config config;
  int rbf_kernel;
};
 
static void config_nnn_delete(
  struct yac_interp_method_config *config) {
  free(config);
}
 
static struct yac_interp_method_config * config_nnn_copy(
  const struct yac_interp_method_config *config) {
  DEF_INTERP_METHOD_CONFIG_TYPE(nnn)
  struct yac_interp_method_config_nnn *copy =
    xmalloc(sizeof(*copy));
  copy->base.vtable = config_nnn->base.vtable;
  copy->base.root_param_cache = NULL;
  copy->config = config_nnn->config;
  copy->rbf_kernel = config_nnn->rbf_kernel;
  return (struct yac_interp_method_config *)copy;
}
 
static int config_nnn_compare(
  void const *a_, void const *b_) {
  DEF_INTERP_METHOD_CONFIG_COMPARE_TYPE(nnn)
  // Compare number of neighbors
  CHECK_INT(config.n);
  // Compare minimum number of neighbors (n_min == 0 is equivalent to n_min == n)
  {
    int a_n_min = get_nnn_n_min(a->config);
    int b_n_min = get_nnn_n_min(b->config);
    if (a_n_min != b_n_min) return (a_n_min > b_n_min) - (a_n_min < b_n_min);
  }
  // Compare max_search_distance with tolerance
  CHECK_ANGLE(config.max_search_distance);
  // Compare weighting type
  CHECK_ENUM(config.type);
  // Compare type-specific scale parameter
  switch (a->config.type) {
    YAC_UNREACHABLE_DEFAULT("invalid n-nearest-neighbour weight type");
    case YAC_INTERP_NNN_GAUSS:
      CHECK_DOUBLE(config.data.gauss_scale);
      break;
    case YAC_INTERP_NNN_RBF:
      CHECK_DOUBLE(config.data.rbf_scale);
      CHECK_INT(rbf_kernel);
      break;
    case YAC_INTERP_NNN_AVG:
    case YAC_INTERP_NNN_DIST:
    case YAC_INTERP_NNN_ZERO:
      // For other types, no further parameter to compare
      break;
  }
  return 0;
}
 
// Gets interpolation method
static enum yac_interpolation_list config_nnn_get_type() {
  return YAC_N_NEAREST_NEIGHBOR;
}
 
static enum yac_interpolation_list config_rbf_get_type() {
  return YAC_RADIAL_BASIS_FUNCTION;
}
 
static struct interp_method * config_nnn_generate(
  struct yac_interp_method_config const *config) {
  const struct yac_interp_method_config_nnn *config_nnn =
    (const struct yac_interp_method_config_nnn *)config;
  return yac_interp_method_nnn_new(config_nnn->config);
}
 
static struct yac_param *config_nnn_get_n_param(
  struct yac_interp_method_config_nnn *config_nnn, int n_default) {
 
  enum {
    N_HAS_DEFAULT = 1,
    N_IS_DEFINED = 1,
  };
 
  int const n_min = 1;
  int const n_max = INT_MAX;
 
  return yac_param_int_new(
    "n",
    &config_nnn->config.n, n_default,
    n_min, n_max, N_HAS_DEFAULT, N_IS_DEFINED);
}
 
static struct yac_param *config_nnn_get_n_min_param(
  struct yac_interp_method_config_nnn *config_nnn, int n_min_default) {
 
  enum {
    N_MIN_HAS_DEFAULT = 1,
    N_MIN_IS_DEFINED = 1,
  };
 
  // a value of 0 is interpreted as "use the value of n"
  int const n_min_min = 0;
  int const n_min_max = INT_MAX;
 
  return yac_param_int_new(
    "n_min",
    &config_nnn->config.n_min, n_min_default,
    n_min_min, n_min_max, N_MIN_HAS_DEFAULT, N_MIN_IS_DEFINED);
}
 
static struct yac_param *config_nnn_get_max_search_distance_param(
  struct yac_interp_method_config_nnn *config_nnn,
  double max_search_distance_default) {
 
  enum {
    MAX_SEARCH_DISTANCE_HAS_DEFAULT = 1,
    MAX_SEARCH_DISTANCE_IS_DEFINED = 1,
  };
 
  double const max_search_distance_min = 0.0;
  double const max_search_distance_max = M_PI - 0.000001;
 
  return yac_param_distance_new(
    "max_search_distance",
    &config_nnn->config.max_search_distance,
    max_search_distance_default,
    max_search_distance_min, max_search_distance_max,
    MAX_SEARCH_DISTANCE_HAS_DEFAULT, MAX_SEARCH_DISTANCE_IS_DEFINED);
}
 
static struct yac_param *config_nnn_get_scale_param(
  char const *name, double *scale_value, double scale_default) {
 
  enum {
    SCALE_HAS_DEFAULT = 1,
    SCALE_IS_DEFINED = 1,
  };
 
  double const scale_min = -DBL_MAX;
  double const scale_max = DBL_MAX;
  enum {
    SCALE_ALLOW_NAN = 0,
    SCALE_ALLOW_INF = 0,
  };
 
  return yac_param_dble_new(
    name, scale_value, scale_default, scale_min, scale_max,
    SCALE_ALLOW_NAN, SCALE_ALLOW_INF, SCALE_HAS_DEFAULT, SCALE_IS_DEFINED);
}
 
static struct yac_param *config_nnn_get_param(
  struct yac_interp_method_config const *config) {
  struct yac_interp_method_config_nnn *config_nnn =
    (struct yac_interp_method_config_nnn *)config;
 
  enum {
    WEIGHTED_HAS_DEFAULT = 1,
    WEIGHTED_IS_DEFINED = 1,
  };
 
  DEF_NAME_TYPE_PAIRS(
    weighted_enum_table,
    DEF_NAME_TYPE_PAIR(distance_weighted, YAC_INTERP_NNN_DIST),
    DEF_NAME_TYPE_PAIR(gauss_weighted, YAC_INTERP_NNN_GAUSS),
    DEF_NAME_TYPE_PAIR(arithmetic_average, YAC_INTERP_NNN_AVG),
    DEF_NAME_TYPE_PAIR(zero, YAC_INTERP_NNN_ZERO))
 
  struct yac_param *weighted_param =
    yac_param_enum_new(
      "weighted",
      (int *)&config_nnn->config.type,
      YAC_INTERP_NNN_WEIGHTED_DEFAULT,
      weighted_enum_table, weighted_enum_table_size,
      WEIGHTED_HAS_DEFAULT, WEIGHTED_IS_DEFINED);
 
  struct yac_param *n_param =
    config_nnn_get_n_param(config_nnn, YAC_INTERP_NNN_N_DEFAULT);
 
  struct yac_param *n_min_param =
    config_nnn_get_n_min_param(config_nnn, YAC_INTERP_NNN_N_MIN_DEFAULT);
 
  struct yac_param *max_search_distance_param =
    config_nnn_get_max_search_distance_param(
      config_nnn, YAC_INTERP_NNN_MAX_SEARCH_DISTANCE_DEFAULT);
 
  struct yac_param *gauss_scale_param =
    config_nnn_get_scale_param(
      "gauss_scale",
      &config_nnn->config.data.gauss_scale, YAC_INTERP_NNN_GAUSS_SCALE_DEFAULT);
 
  struct yac_param *root_param_array[] = {
    weighted_param,
    n_param,
    n_min_param,
    max_search_distance_param,
    gauss_scale_param
  };
  enum {
    ROOT_PARAM_ARRAY_SIZE =
      sizeof(root_param_array) / sizeof(root_param_array[0])
  };
 
  return yac_param_struct_new(
    "nnn", root_param_array, ROOT_PARAM_ARRAY_SIZE);
}
 
static struct yac_interp_method_config_vtable yac_interp_method_config_vtable_nnn = {
  .delete = config_nnn_delete,
  .copy = config_nnn_copy,
  .compare = config_nnn_compare,
  .get_type = config_nnn_get_type,
  .generate = config_nnn_generate,
  .get_param = config_nnn_get_param
};
 
static struct yac_param *config_rbf_get_param(
  struct yac_interp_method_config const *config) {
  struct yac_interp_method_config_nnn *config_nnn =
    (struct yac_interp_method_config_nnn *)config;
 
  enum {
    RBF_KERNEL_HAS_DEFAULT = 1,
    RBF_KERNEL_IS_DEFINED = 1
  };
 
  DEF_NAME_TYPE_PAIRS(
    rbf_kernel_enum_table,
    DEF_NAME_TYPE_PAIR(gauss_kernel, 0),
    DEF_NAME_TYPE_PAIR(gauss, 0))
 
  struct yac_param *n_param =
    config_nnn_get_n_param(config_nnn, YAC_INTERP_RBF_N_DEFAULT);
 
  struct yac_param *n_min_param =
    config_nnn_get_n_min_param(config_nnn, YAC_INTERP_RBF_N_MIN_DEFAULT);
 
  struct yac_param *max_search_distance_param =
    config_nnn_get_max_search_distance_param(
      config_nnn, YAC_INTERP_RBF_MAX_SEARCH_DISTANCE_DEFAULT);
 
  struct yac_param *rbf_scale_param =
    config_nnn_get_scale_param(
      "rbf_scale",
      &config_nnn->config.data.rbf_scale, YAC_INTERP_RBF_SCALE_DEFAULT);
 
  struct yac_param *rbf_kernel_param =
    yac_param_enum_new(
      "rbf_kernel",
      &config_nnn->rbf_kernel,
      YAC_INTERP_RBF_KERNEL_DEFAULT,
      rbf_kernel_enum_table, rbf_kernel_enum_table_size,
      RBF_KERNEL_HAS_DEFAULT, RBF_KERNEL_IS_DEFINED);
 
  struct yac_param *root_param_array[] = {
    n_param,
    n_min_param,
    max_search_distance_param,
    rbf_scale_param,
    rbf_kernel_param,
  };
  enum {
    ROOT_PARAM_ARRAY_SIZE =
      sizeof(root_param_array) / sizeof(root_param_array[0])
  };
 
  return yac_param_struct_new(
    "rbf", root_param_array, ROOT_PARAM_ARRAY_SIZE);
}
 
static struct yac_interp_method_config_vtable yac_interp_method_config_vtable_rbf = {
  .delete = config_nnn_delete,
  .copy = config_nnn_copy,
  .compare = config_nnn_compare,
  .get_type = config_rbf_get_type,
  .generate = config_nnn_generate,
  .get_param = config_rbf_get_param
};
 
struct yac_interp_method_config *yac_interp_method_config_default_nnn_new(void) {
  struct yac_interp_method_config_nnn *config_nnn = xmalloc(sizeof(*config_nnn));
  config_nnn->base.vtable = &yac_interp_method_config_vtable_nnn;
  config_nnn->base.root_param_cache = NULL;
  config_nnn->rbf_kernel = YAC_INTERP_RBF_KERNEL_DEFAULT;
  config_nnn->config = (struct yac_interp_method_nnn_config){
    .type = YAC_INTERP_NNN_WEIGHTED_DEFAULT,
    .n = YAC_INTERP_NNN_N_DEFAULT,
    .n_min = YAC_INTERP_NNN_N_MIN_DEFAULT,
    .max_search_distance = YAC_INTERP_NNN_MAX_SEARCH_DISTANCE_DEFAULT,
    .data = {.gauss_scale = YAC_INTERP_NNN_GAUSS_SCALE_DEFAULT},
  };
  return (struct yac_interp_method_config *)config_nnn;
}
 
struct yac_interp_method_config *yac_interp_method_config_default_rbf_new(void) {
  struct yac_interp_method_config_nnn *config_nnn = xmalloc(sizeof(*config_nnn));
  config_nnn->base.vtable = &yac_interp_method_config_vtable_rbf;
  config_nnn->base.root_param_cache = NULL;
  config_nnn->rbf_kernel = YAC_INTERP_RBF_KERNEL_DEFAULT;
  config_nnn->config = (struct yac_interp_method_nnn_config){
    .type = YAC_INTERP_NNN_RBF,
    .n = YAC_INTERP_RBF_N_DEFAULT,
    .n_min = YAC_INTERP_RBF_N_MIN_DEFAULT,
    .max_search_distance = YAC_INTERP_RBF_MAX_SEARCH_DISTANCE_DEFAULT,
    .data = {.rbf_scale = YAC_INTERP_RBF_SCALE_DEFAULT},
  };
  return (struct yac_interp_method_config *)config_nnn;
}