2016-11-03 3 views
0

При вычислении двукратной интеграции я столкнулся с несогласованными результатами интеграции() и suave() из пакета «R2Cuba».Непоследовательные результаты между R2Cuba и integrate()

где

Аналогичный вопрос, связанный с этой теме можно найти here, давайте установим п (х) = 6 * Sin (х) и г (х) = 1 и верхнее время T = 3.

Следующее использование кода интеграция() и обходительный(), чтобы вычислить объективный 2 раза интеграла:

library(R2Cuba) 

integrand = function(x){6*sin(x)} 

phi = function(x){integrate(integrand,lower=x,upper = 3)[["value"]]^2} 

NDIM=1 
NCOMP=1 

phicuba= function(x){suave(NDIM,NCOMP,integrand,lower=x,upper=3)$value^2} 

foldintegral = integrate(Vectorize(phi),lower = 0,upper = 3) 

foldintegralcuba = suave(NDIM,NCOMP,phicuba,lower = 0,upper = 3) 

Результатов:

> foldintegral 
167.3934 with absolute error < 1.9e-12 
> foldintegralcuba 
integral: 6.365749 (+-0.0057) 
nregions: 8; number of evaluations: 10000; probability: 1 

, который не является последовательным. Однако, если мы только сравнить phi и phicuba

> phi(2) 
[1] 11.85476 
> phicuba(2) 
Iteration 1: 1000 integrand evaluations so far 
[1] 3.44367 +- 0.0429822 chisq 0 (0 df) 
Iteration 2: 2000 integrand evaluations so far 
[1] 3.4436 +- 0.00866171 chisq 0.00513973 (2 df) 
Iteration 3: 3000 integrand evaluations so far 
[1] 3.44181 +- 0.00386046 chisq 5.38913 (5 df) 
Iteration 4: 4000 integrand evaluations so far 
[1] 3.44283 +- 0.00206345 chisq 23.0816 (8 df) 
[1] 11.85309 

Мы получили результат, который можно рассматривать как последовательное. Кроме того, если мы будем использовать заменить подынтегральное внутри suave()

> suave(NDIM,NCOMP,phi,lower = 0,upper = 3) 
Iteration 1: 1000 integrand evaluations so far 
[1] 167.426 +- 4.91983  chisq 0 (0 df) 
Iteration 2: 2000 integrand evaluations so far 
[1] 167.315 +- 0.771248  chisq 0.00208764 (2 df) 
Iteration 3: 3000 integrand evaluations so far 
[1] 167.357 +- 0.432941  chisq 5.50239 (5 df) 
Iteration 4: 4000 integrand evaluations so far 
[1] 167.362 +- 0.129012  chisq 5.51661 (8 df) 
integral: 167.3621 (+-0.13) 
nregions: 4; number of evaluations: 4000; probability: 0.2988006 

У нас еще есть последовательный результат.

Мне жаль, что я не могу использовать suave(), поскольку он намного быстрее при работе со сложным интегралом, но почему это противоречивое существует?

===================== ОБНОВЛЕНИЕ ============================ ===========

Похоже, что алгоритм, используемый в suave(), From the Cuba Documentation, позволяет использовать глобально адаптивное подразделение + выбор важности. Если изменить suave() на cuhre(), которые используют только глобально адаптивное подразделение, все должно быть последовательным.

+0

@ZheyuanLi Спасибо, что напомнили! – lz10086lz

ответ

1

Похоже, что алгоритм, используемый в suave(), From the Cuba Documentation, позволяет использовать глобально адаптивное подразделение + выбор важности. Если изменить suave() на cuhre(), которые используют только глобально адаптивное подразделение, все должно быть последовательным.