Calculis3

R Markdown

這是R 簡書, http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

We want to check the statement by R codes. \[\displaystyle\sum_{k=1}^\infty \frac{1}{2^{-2k+1}}=\frac{1}{2} + 2^{-3} +\frac{1}{2^5}+ \cdots =\frac{2}{3}\]

kk=2*(1:100)-1 # kk= seq(1,199, by=2)
kk[1:10]       # 列出前10個元素

##  [1]  1  3  5  7  9 11 13 15 17 19

sum(2^(-kk))   # 理論值=2/3

## [1] 0.6666667

sum=0           # C style program
for(i in 1:100)
  sum=sum+1/2^(2*i-1)
sum            # 理論值=2/3

## [1] 0.6666667

Including Plots

The parameteric equations \[ \begin{cases} x(t)=t^2\\ y(t)=t^3-3t \end{cases} \] When \(t=\sqrt{3}\), the slope at \((3,0)\) is \[ m=\frac{dy}{dx}=\frac{dy/dt}{dx/dt}=\frac{3t^2-3}{2t}|_{t=\sqrt{3}}=\sqrt{3}\]

tt=seq(-3, 3, by=0.01)
xx=tt^2; yy=tt^3-3*tt

#plot(xx, yy)
plot(xx, yy, type="l", col=2, ylim=c(-5,5), xlab="x",ylab="y")
#plot(xx, yy, type="h")
hh=0.2
t1=sqrt(3);
x1=t1^2; y1=t1^3-3*t1
 t2=t1+hh
 x2=t2^2; y2=t2^3-3*t2
points(x1, y1, col="blue")
points(x2, y2, col="blue")
lines(c(x1,x2), c(y1,y2), col=4, lwd = 3)
m=(y2-y1)/(x2-x1)
cat("理論值=", sqrt(3)," \t 而此例兩點斜率為", m, "\n")

## 理論值= 1.732051     而此例兩點斜率為 1.932051

hh=0.02
t1=sqrt(3);
x1=t1^2; y1=t1^3-3*t1
 t2=t1+hh
 x2=t2^2; y2=t2^3-3*t2
points(x1, y1, col="green")
points(x2, y2, col="green")

m=(y2-y1)/(x2-x1)
cat("理論值=", sqrt(3)," \t 而此例兩點斜率為", m, "\n")

## 理論值= 1.732051     而此例兩點斜率為 1.752051

Example

Use R to find the area under the function \(f(x)=x^2\) from \(x=0\) to \(x=1\). \[ A= \sum_{i=1}^n (x_i)^2 \triangle x_i\]

nn=100; 
xx=seq(0, 1, by= 0.005)
yy=xx^2
plot(xx, yy, type="h", col=4,  xlab="x",ylab="y")

A=0; for(i in 1:nn) A=A+i^2/nn^3; 
cat("\n The area is\t", A, "\n")

## 
##  The area is  0.33835

# Method 2
w.x=(1-0)/nn
xx=seq(0,1, by=w.x) #0, 2/nn, 4/nn, ... ,2
cat("\n The area is\t",  sum(xx^2*w.x), "\n")

## 
##  The area is  0.33835

Formula of arc length of \(r=f(\theta)\), see for example LibreTexts
\[ L=\int_{\alpha}^\beta \sqrt{f(\theta)^2+\left(f'(\theta)\right)^2}d\theta \]

• \(f(x)=x^2\), \(f'(x)=2x\) \[\begin{aligned} L&=\int_0^1 \sqrt{1+(2x)^2}\,dx\\ &=\frac{1}{2}\int_0^{\pi/4} \sec^3\theta d\theta \\ &= 1.4789428 \end{aligned} \]

# Call wxMaxima
'integrate(sqrt(1+4*x^2), x,0,1);
integrate(sqrt(1+4*x^2), x,0,1);
float(integrate(sqrt(1+4*x^2), x,0,1););

# Calculating the arc length 
L=0
xx=seq(0, 1, by= 0.005)
yy=xx^2
nn=length(xx)
for ( i in 1:(nn-1)){
  L=L+sqrt((xx[i+1]-xx[i])^2+(yy[i+1]-yy[i])^2)
}
cat("\n Thearc length is\t", L, "\n") # L=1.4789428575445975

## 
##  Thearc length is     1.478941

Example 3:

Find the area under one arch of the cycloid \[ \begin{cases} x=r(\theta-\sin\,\theta)\\ y=r(1-\cos\,\theta) \end{cases}\]

theta=seq(-pi, 3*pi+.1, by=0.05)
r=1;
xx=theta-sin(theta); 
yy=1-cos(theta); 

#plot(xx, yy)
plot(xx, yy, type="h", col=2, ylim=c(-.1, 2), xlab="x",ylab="y")

\[ A=\int_0^{2pr} r(1-\cos\,\theta)\frac{d}{d\theta}r(\theta-\sin\,\theta)=\cdots= d\theta=3\pi r^2 \]

Example 5:

Find the length of one arch of the cycloid \[ \begin{cases} x=r(\theta-\sin\,\theta)\\ y=r(1-\cos\,\theta) \end{cases}\]