2013-02-08 5 views
4

Я пытаюсь написать функцию, которая запускает цикл в C++ из R, используя Rcpp.Итерация R и C++

У меня есть матрица Z, которая на одну строку короче, чем матрица OUT, которую функция должна возвращать, потому что каждое положение первой строки OUT будет задано скалярным sigma_0.

Функция должна реализовывать дифференциальное уравнение. Каждая итерация зависит от значения от матрицы Z, а также от ранее созданного значения матрицы OUT.

То, что я получил это:

cppFunction(' 
    NumericMatrix sim(NumericMatrix Z, long double sigma_0, long double delta, long double omega, long double gamma) { 
     int nrow = Z.nrow() + 1, ncol = Z.ncol(); 
     NumericMatrix out(nrow, ncol); 

     for(int q = 0; q < ncol; q++) { 
      out(0, q) = sigma_0; 
     } 



     for(int i = 0; i < ncol; i++) { 

      for(int j = 1; j < nrow; j++) { 
       long double z = Z(j - 1, i); 
       long double sigma = out(j - 1, i); 
       out(j, i) = pow(abs(z * sigma) - gamma * z * sigma, delta); 

      } 
     } 


     return out; 
    } 
') 

К сожалению, я довольно уверен, что это не работает. Функция выполняется, но рассчитанные значения неверны - я проверил с помощью простых примеров в Excel и простое R-кодирование. Я разделил основное дифференциальное уравнение, пытаясь его постепенно наращивать, чтобы увидеть, когда реализация i Excel и R с использованием C++ начинает различаться. Кажется, что когда я начинаю использовать функцию abs() и функцию power(), но я просто не могу сузить проблему. Любая помощь будет принята с благодарностью - также я могу упомянуть, что это первый раз, когда я использую C++ и C++ вместе с R.

ответ

4

Я думаю, вы хотите fabs, а не abs. abs работает на int с, а fabs работает на парных/двухместных.

+0

Да, это делает трюк. Большое вам спасибо ... просто потратьте полдня на то, чтобы его очень ценили. – user2055639

+0

+1 - Одна из тех вещей, которые являются болезненными, потому что оба они определены в C ... Было бы хорошо, если бы вы (@ user2055639) могли теперь «перевыполнять» и «принимать» ответ, как это принято в StackOverflow. –

+0

Я бы с удовольствием @DirkEddelbuettel. Просто у приёма не хватило репутации для голосования. – user2055639