>>883
別スレの問題を「プログラムで面積を数値計算させていたらエラー発生

体験したければ読みにいコードの下記をどうぞ。

r= 0.66017210648907 ; s=2.89088405538017

# 面積計算にintegrate で エラー
sub <- function(c,Area1=1){ # c 円Cの中心のx座標
b0=(c^2+r^2-s^2)/(2*c) # B(b0,b1) # 円Cと円Oの交点Bの座標
b1=sqrt(r^2-(c^2+r^2-s^2)^2/(4*c^2))

blue = function(x) sqrt(s^2-(x-c)^2) # 円Cの方程式
BLUE= integrate(blue,b0,c+s)$value*2 # 青の面積
green = function(x) sqrt(r^2-x^2) # 円Oの方程式
GREEN=integrate(green,-r,b0)$value*2 # 緑の面積

GREEN+BLUE-Area1 # GREEN+BLUE=1
}
sub=Vectorize(sub)
x=seq(-r-s,r-s,length=100) ; plot(x,sub(x),pch=19) ; abline(h=0,lty=3)


# 面積計算にcubintegrate でエラーなし
library(cubature)
sub <- function(c,Area1=1){ # c 円Cの中心のx座標
b0=(c^2+r^2-s^2)/(2*c) # B(b0,b1) # 円Cと円Oの交点Bの座標
b1=sqrt(r^2-(c^2+r^2-s^2)^2/(4*c^2))

blue = function(x) sqrt(s^2-(x-c)^2) # 円Cの方程式
BLUE= cubintegrate(blue,b0,c+s,method="pcubature")$integral*2 # 青の面積
green = function(x) sqrt(r^2-x^2) # 円Oの方程式
GREEN=cubintegrate(green,-r,b0,method='pcubature')$integral*2 # 緑の面積

GREEN+BLUE-Area1 # GREEN+BLUE=1
}
sub=Vectorize(sub)
x=seq(-r-s,r-s,length=100) ; plot(x,sub(x),pch=19) ; abline(h=0,lty=3)