#1. Repaso Código Integrales (aproximaciones de Área)
\[\int_{0}^{2e}ln(x)dx = lim_{b \to 0}\int_{b}^{2e}ln(x)dx \]
inte <- function(x) x^3
x1 <- seq(0,1,by = 0.01)
plot(x1, inte(x1),type='l', col = 'blue',lwd=2)
polygon(c(0,x1,1),c(0,inte(x1),0), col ='green')
#2. Suma trapezoidal
n <- 20
a <- 0
b <- 1
h <- (b-a)/n
xi <- a+(0:n)*h;xi
## [1] 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.40 0.45 0.50 0.55 0.60 0.65 0.70
## [16] 0.75 0.80 0.85 0.90 0.95 1.00
plot(x1, inte(x1),type='l', col = 'black',lwd=2, main = 'suma trapezoidal')
for(j in 1:n){
segments(xi[j],0,xi[j+1],0,col = 'blue',lwd = 2)
segments(xi[j+1],0,xi[j+1],inte(xi[j+1]),col = 'blue')
segments(xi[j+1],inte(xi[j+1]),xi[j],inte(xi[j]),col = 'blue')
}
pot1 <- (h/2)*(inte(a)+inte(b))
pot2 <- (h/2)*(sum(2*inte(a+(1:(n-1)*h))))
protra <- pot1+pot2
protra
## [1] 0.250625
#3. Sumas de Riemann por la derecha
plot(x1, inte(x1),type='l', col = 'black',lwd=2, main = 'Riemann por la derecha')
for (l in 1:n) {
segments(xi[l],0,xi[l+1],0,col = 'red',lwd=2)
segments(xi[l+1],0,xi[l+1],inte(xi[l+1]),col = 'red',lwd=2)
segments(xi[l+1],inte(xi[l+1]),xi[l],inte(xi[l+1]),col = 'red',lwd=2)
segments(xi[l],inte(xi[l+1]),xi[l],0,col = 'red',lwd=2)
}
prode<-sum(inte(xi[-1]))*h
prode
## [1] 0.275625
Por la izquierda
plot(x1, inte(x1),type='l', col = 'black',lwd=2, main = 'Riemann por la izquierda')
for (l in 1:n) {
segments(xi[l],0,xi[l+1],0,col = 'green',lwd=2)
segments(xi[l],0,xi[l],inte(xi[l]),col = 'green',lwd=2)
segments(xi[l+1],0,xi[l+1],inte(xi[l]),col = 'green',lwd=2)
segments(xi[l+1],inte(xi[l]),xi[l],inte(xi[l]),col = 'green',lwd=2)
}
proiz<-sum(inte(xi[-(n+1)]))*h
proiz
## [1] 0.225625
plot(x1, inte(x1),type='l', col = 'black',lwd=2, main = 'Montecarlo')
mc.in <- function(n1,a,b){
x <- runif(n1,a,b)
y <- runif(n1,0,inte(b))
co <- ifelse(y<=inte(x),'red','black')
points(x,y,col = ifelse(y<=inte(x),'red','black'))
}
n1=5000
mc.in(n1,a,b)
x2 <- runif(n1,a,b)
y2 <- runif(n1,0,inte(b))
co <- ifelse(y2<=inte(x2),'red','black')
table(co)/n1
## co
## black red
## 0.7568 0.2432
mc.3 <- function(n3,a,b){
x <- runif(n3,a,b)
f <- (b-a)*mean(inte(x))
return(f)
}
n3=1000
mc.3(n3,a,b)
## [1] 0.2528943
Resolver por simulación la siguiente integral
\[\int_{0}^{2e}ln(x)dx = lim_{b \to 0}\int_{b}^{2e}ln(x)dx \]
integralT <- function(x) log(x)
x4 <- seq(0,2*exp(1),by = 0.01)
plot(x4, integralT(x4),type='l', col = 'blue',lwd=2,ylab = 'y',xlab = 'x')
abline(v=1,h=0)
plot(x4, integralT(x4),type='l',ylim = c(0,2), col = 'black',lwd=2,ylab = 'y',xlab = 'x')
a <- 0.0001
b <- 2*exp(1)
puntos <- 10000
SimuInte <- function(puntos,a,b){
x <- runif(puntos,a,b)
y <- runif(puntos,0,integralT(b))
co <- ifelse(y<=integralT(x),'yellow','black')
points(x,y,col = ifelse(y<=integralT(x),'yellow','black'))
}
SimuInte(puntos,a,b)
rect(1,0,b,integralT(b), border = 'red')
x5 <- runif(puntos,a,b)
y5 <- runif(puntos,0,integralT(b))
co <- ifelse(y5<=integralT(x5),'yellow','black')
table(co)/puntos
## co
## black yellow
## 0.4843 0.5157
(table(co)/puntos)*(b-1)*(integralT(b))
## co
## black yellow
## 3.637943 3.873812
plot(x4, integralT(x4),type='l',ylim = c(-10,2),xlim = c(0,1), col = 'black',lwd=2,ylab = 'y',xlab = 'x')
a <- 0.0001
b <- 1
puntos <- 50000
SimuInte <- function(puntos,a,b){
x <- runif(puntos,a,b)
y <- runif(puntos,-10,0)
co <- ifelse(y<=integralT(x),'black','yellow')
points(x,y,col = ifelse(y<=integralT(x),'black','yellow'))
}
SimuInte(puntos,a,b)
rect(0,-10,1,0,border = 'red')
x5 <- runif(puntos,a,b)
y5 <- runif(puntos,-10,0)
co <- ifelse(y5<=integralT(x5),'black','yellow')
(table(co)/puntos)*-4
## co
## black yellow
## -3.60288 -0.39712