Использование функции gsl_function с использованием поставляемой функции R

0

В настоящее время я пытаюсь использовать gsl_function из библиотеки GSL через RcppGSL, используя файл.cpp и вызывая его с помощью sourceCpp(). Идея состоит в том, чтобы выполнить численное интегрирование с gsl_integration_qags, также из GSL. Мой C-код вызывает определенную пользователем функцию R (SomeRFunction в моем коде ниже), сохраненную в глобальной среде. Код:

#include <RcppGSL.h>
#include <stdio.h>
#include <math.h>
#include <gsl/gsl_integration.h>
#include <gsl/gsl_vector.h>

// [[Rcpp::depends(RcppGSL)]]
// [[Rcpp::export]]
double integfn1(double ub, double lb, double abserr, double Rsq, int pA){ 

  double result, abserror;
  Rcpp::Environment global = Rcpp::Environment::global_env();
  Rcpp::Function f1 = global["SomeRFunction"];

  struct my_f_params { double Rsq; int pA; };
  struct my_f_params params = { Rsq, pA };
  gsl_function F;
  F.function = & f1;
  F.params = &params;
  double lb1 =  lb;
  double ub1 = ub;

  gsl_integration_workspace * w = gsl_integration_workspace_alloc (1000);

  gsl_integration_qags (&F, lb1, ub1, 1e-8, 1e-8, 1000, w, &result, &abserror);

  printf ("result          = % .18f\n", result);
  printf ("estimated error = % .18f\n", abserror);

  gsl_integration_workspace_free (w);

  return result;
}

И появляется следующая ошибка:

"cannot convert 'Rcpp::Function* {aka Rcpp::Function_Impl<Rcpp::PreserveStorage>*}' to 'double (*)(double, void*)' "

Проблема заключается в строке, где я объявляю, что функция для интеграции (т.е. " F.function = & f1; ").

Я искал похожие проблемы, но не смог найти ничего перечисленного... Любой намек будет очень благодарен!

Большое спасибо

Теги:
rcpp
gsl

3 ответа

1
Лучший ответ

Я создал рабочий пример (и приурочен к C), в котором вы можете передать произвольную пользовательскую функцию R функции GSL QAWF. Вы также сможете обобщить его на другие функции gsl.

См. Пример здесь: https://sites.google.com/site/andrassali/computing/user-supplied-functions-in-rcppgsl

Как отмечалось выше, конкурирующая реализация C происходит намного быстрее.

Поделиться
  • 0
    Андрас: Спасибо, что поделился, это отличная помощь!
  • 0
    Рад, я надеюсь, что это помогло. Это дало вам решение в конце концов?
Показать ещё 1 комментарий
1

F.function ожидает функцию этой double (*) (double x, void * params) сигнатуры, поэтому функция принимает double затем void*. Вам нужно помочь Rcpp::Function если вы хотите, чтобы это Rcpp::Function.

Это типичные типы отбрасывания C apis. Итак, что вам нужно, чтобы дать как function это то, что понимает, что вы имеете в качестве params и может вызвать функцию R, с правильными параметрами и преобразовать в двойную. Простым способом сделать это является рассмотрение функции R как данных, поэтому добавьте struct чтобы сделать ее примерно такой:

struct my_f_params { double Rsq; int pA; Function& rfun };

Таким образом, функция, которую вы передаете как function может быть примерно такой:

double proxy_fun( double x, void* params ){
  my_f_params* data = reinterpret_cast<my_f_params*>(params) ;
  return as<double>( data->rfun(x, data->Rsq, data->pA ) ) ;
} ;

Таким образом, вы можете построить gsl_function примерно так:

my_f_params params = { Rsq, pA, f1 };
gsl_function F ; 
F.function = &proxy_fun ;
F.params = &params;

С небольшим количеством вариативных шаблонов и кортежей С++ 11 вы можете обобщить это на что угодно, вместо пары (double, int), но это больше работает.

Поделиться
  • 0
    Это действительно хорошее решение, большое спасибо, Ромен .... Моя R-функция не очень затратна в вычислительном отношении, в этом случае вычислительные затраты будут значительно снижены при переписывании моей функции на C?
  • 1
    Трудно сказать. не верьте тому, что кто-то говорит, измерьте накладные расходы для себя.
Показать ещё 2 комментария
1

Быстрые комментарии:

  1. Вы пишете: мой код C вызывает определенную пользователем функцию R (SomeRFunction в моем коде ниже), сохраненную в глобальной среде. Таким образом, ваш код C все равно будет замедляться при каждой оценке функции из-за стоимости вызова R и более медленной скорости R.

  2. Мой пакет RcppDE (который также предназначен для optmisation) как пример использования внешних указателей (Rcpp::XPtr) для передачи пользовательской функции C до оптимизатора. Та же гибкость. более высокая скорость.

  3. Ваша ошибка компиляции находится именно на этом пересечении - вы не можете просто передать объект Rcpp указателю на пустоту. Rcpp::XPtr помогает вам там, но вам также нужно знать, что вы делаете, поэтому рабочий пример может оказаться полезным.

  4. Это действительно вопрос Rcpp, но вы не добавили тег, поэтому я просто сделал это.

Поделиться
  • 0
    Я перепишу свою "пользовательскую" функцию на C и посмотрю на пакет RcppDE ... Большое спасибо !!!!

Ещё вопросы

Сообщество Overcoder
Наверх
Меню