仕様上はできそうなので試してみたらできたと言うだけの話なのですが、Rで呼び出したFortranのプロシージャからRの関数を呼ぶ方法を記しておきます。
技術的にどうこうと言うのは無いのですが、流れはややこしくて以下のようになります。
- .CallでRのコールバック関数のある環境を引数にCの関数（invoker）を呼ぶ
- Cの関数（invoker）で環境のポインターを引数に、Fortranのプロシージャ（mainroutine）を呼ぶ
- iso_c_bindingsしたFortranのプロシージャ（mainroutine）から、環境のポインターを引数に（別の）Cの関数（call_r）を呼ぶ
- Cの関数（call_r）から、Rの関数（r_callback）を呼ぶ
Rの環境（名前空間）の変数を連携するところがポイントで、そうしないとRのローカル変数を呼ぶことが難しくなります*1。
Fortranのコード
Fortranのプロシージャ（mainroutine）は、呼び出し前のデータと、呼び出し後のデータを表示するだけです。
注意すべきは、use iso_c_bindingをしてCの変数を誤動作のリスクを減らして扱えるようにすること、Cの関数の定義をinterfaceで書くこと、Rの環境（名前空間）はFortranから見たらCのポインターであることぐらいです。
Subroutine mainroutine(envir_ptr) bind(C, name = "mainroutine") use iso_c_binding implicit None type(c_ptr) :: envir_ptr interface INTEGER(c_int) function call_r(nr, r, ptr) bind(C, name = "call_r") use iso_c_binding implicit none INTEGER(c_int), intent(in) :: nr REAL(c_double), dimension(:), intent(inout) :: r type(c_ptr), intent(in) :: ptr end function end interface INTEGER(c_int), parameter :: nr = 3 REAL(c_double), dimension(nr) :: r = (/1., 2., 3./) write(*, *) "r = ", r write(*, *) "callback: ", call_r(nr, r, envir_ptr) write(*, *) "r = ", r End Subroutine
ファイル名をRCF_F.f90として保存します。
Cのコード
Rから呼ばれるinvoker関数と、Fortranから呼ばれるcall_r関数を用意しています。RとCの連携用の手続きが長めにありますが、そこは仕様どおりです。invokerで環境のポインターをFortranに渡していることに注意してください。
#include <R.h> #include <Rinternals.h> /* *nrv: rvの大きさ **rv: 呼び出し元（Fortran）で確保した領域 **envir_ptr: Rの環境（名前空間） */ int call_r(int *nrv, double **rv, void **envir_ptr) { SEXP R_fcall, ans, s, t, r; unsigned i; /* r_callback(theta)のAbstract Syntax Treeを構築 */ /* 命令用の領域を確保; tとsには同じポインタが入る */ t = s = PROTECT(allocList(2)); /* LANGSXPを指定 */ SET_TYPEOF(s, LANGSXP); /* 最初の領域には "r_callback" を入れる; 後の環境と同様に、Rから渡された変数へのポインターのセットもできるが、C言語でRの関数名を指定した。 */ SETCAR(t, install("r_callback")); /* tを次の領域に移動させる; sは最初の領域のまま */ t = CDR(t); /* 次の領域にはベクターを入れる */ PROTECT(r = allocVector(REALSXP, *nrv)); for(i = 0; i<*nrv; i++) REAL(r)[i] = (*rv)[i]; /* R用の変数を引数にセット */ SETCAR(t, r); /* thetaと名前をつけておく */ SET_TAG(t, install("theta")); /* 環境*envir_ptrでevalを行う。代わりに.R_GlobalEnvを指定しても呼び出せるが、呼ばれた関数がローカル変数を読めなくなる。*/ PROTECT(ans = eval(s, *envir_ptr)); /* 戻り値の型を確認 */ if(!isReal(ans)){ UNPROTECT(3); error("'ans' should be of numeric."); } /* Rの仕様で戻り値はベクトル */ for(i = 0; i<length(ans) && i<*nrv; i++) (*rv)[i] = REAL(ans)[i]; /* PROTECTを外して、ガーベッジコレクターが領域を開放できるようにする */ UNPROTECT(3); return 0.0; } /* Fortranのプロシージャを関数として宣言 */ void mainroutine(void *addr); /* Rから呼び出されるCの関数 */ SEXP invoker(SEXP rho){ /* 引数が環境か確認 */ if(!isEnvironment(rho)) error("'rho' should be an environment"); /* Rの環境を引数に、Fortranのプロシージャを呼ぶ */ mainroutine(&rho); return R_NilValue; }
RCF_C.cと名前をつけて保存します。なお、RFC.f90とRFC.cにすると、オブジェクトファイルがRCF.oとRCF.oで被ってDLLの作成ができません。
Rのコード
できたRCF.so（もしくはRCF.dll）をRから呼び出します。
dll <- paste("RCF", .Platform$dynlib.ext, sep="") dyn.load(dll) # ローカル環境（名前空間）をつくってCのコードを呼び出す example <- function(){ # ローカル変数 local_variables <- sqrt(c(2:4)) # C経由でFortranにコールバックされる関数 r_callback <- function(theta){ local_variables + theta } # r_callbackの書かれた環境を引き数に、invokerを呼び出す .Call("invoker", environment(r_callback)) } example() dyn.unload(dll)
引数とローカル変数を合計した値が、表示されましたね。
まとめ
実際に書いてみるとオーバーヘッドが大きく、環境の受け渡しが煩雑なので安易に使いたくはない感じですが、準備しているFortranでHMC法によるベイズ推定（と言うかサンプリング）をするコードで、Rで定義した目的関数を呼び出したいので試してみました。
*1:逆に言うと、Rの関数が定義されたスコープの変数を参照する必要がなければ、環境の受け渡しの部分は無くてすっきりします。
