Shared2017-09-19-120509.sagewsOpen in CoCalc
def simpsons_rule(fcn,a,b,n):
# output: simpsons rule approx of int_a^b fcn(x)dx using n steps (n must be an even integer)
Deltax = (b-a)*1.0/n
n2=int(n/2)
coeffs = [4,2]*n2
coeffs = [1] +coeffs[:n-1]+[1]
valsf = [f(a+Deltax*i) for i in range(n+1)]
return (Deltax/3)*sum([coeffs[i]*valsf[i] for i in range(n+1)])

#Example
n=10;
f(x)=1/x
print simpsons_rule(f(x),1,2,n).n();

0.693150230688930




def trapezoid_rule(fcn,a,b,n):
# output: trapezoid rule approx of int_a^b fcn(x)dx using n steps
Deltax = (b-a)*1.0/n
coeffs = [2]*(n-1)
coeffs = [1]+coeffs+[1]
valsf = [f(a+Deltax*i) for i in range(n+1)]
return (Deltax/2)*sum([coeffs[i]*valsf[i] for i in range(n+1)])

#Example
n=10;
f(x)=1/x
print trapezoid_rule(f(x),1,2,n).n();

0.693771403175428


def midpoint_rule(fcn,a,b,n):
# output: midpoint rule approx of int_a^b fcn(x)dx using n steps
Deltax = (b-a)*1.0/n
xs=[a+Deltax*i for i in range(n+1)]
ysmid=[f((xs[i]+xs[i+1])/2) for i in range(n)]
return Deltax*sum(ysmid)

#Example
n=10;
f(x)=1/x
print midpoint_rule(f(x),1,2,n).n();

0.692835360409960


def righthand_rs(fcn,a,b,n):
# output: right-hand riemann sum approx of int_a^b fcn(x)dx using n steps
Deltax = (b-a)*1.0/n
return Deltax*sum([f(a+Deltax*(i+1)) for i in range(n)])

#Example
n=10
f(x)=1/x
print righthand_rs(f(x),1,2,n).n();

0.668771403175428

def lefthand_rs(fcn,a,b,n):
# output: left-hand riemann sum approx of int_a^n fcn(x)dx using n steps
Deltax = (b-a)/n
return Deltax*sum([f(a+Deltax*i) for i in range(n)])

#Example
n=6
f(x)=1/x
print lefthand_rs(f(x),1,2,n).n();

0.736544011544012

# the actual integral gives
f(x)=1/x
integral(f , (x,1,2))
N(log(2))

log(2) 0.693147180559945


x = var('x')
f1(x) = 1/x
f = Piecewise([[(1,2),f1]])
trapezoid_sum = f.trapezoid(10)
P = f.plot(rgbcolor=(0,0,1), plot_points=40)
Q = trapezoid_sum.plot(rgbcolor=(1,0,0), plot_points=40)
L = add([line([[a,0],[a,f(x=a)]], rgbcolor=(1,0,0)) for (a,b), f in trapezoid_sum.list()])
M = line([[2,0],[2,f1(2)]], rgbcolor=(1,0,0))
show(P + Q + L + M)



f(x) = x^2-5*x+10
f = Piecewise([[(0,10), f]])
g = f.riemann_sum(6, mode="midpoint")
F = f.plot(color="blue")
R= add([line([[a,0],[a,f(x=a)],[b,f(x=b)],[b,0]], color="red") for (a,b), f in g.list()])
show(F+R)


integral(1/(x^2+1),(x, 0, oo)) # we getting pi/2
plot(1/(x^2+1),x,0,200) # this gives us that the idea that integral of the function will converge from 0 to infinity.

1/2*pi

plot(1/(sqrt(x-2)),2,5) # discontinuous integrands (Improper integral of type 2)
integral(1/(sqrt(x-2)),2,5) # this integral is convergent


2*sqrt(3)


f(x)=1/(sqrt(x^2-x))
g(x)=1/x
plot(f, 3, 7,gridlines='minor', frame=True, color='red')+ plot(g,3 ,7)



f(x)=x/(sqrt(x^2+1))
g(x)=x^4-x
plot(f, 0, 1.4,gridlines='minor', frame=True)+ plot(g,0 ,1.4)+ text('$\int_0^{1.18}(f-g)(x)\,dx$',(0.5,2),fontsize=20)



f(x)=x^2
g(x)=2*x-x^2
plot(f(x), 0, 2)+ plot(g(x),0 ,2)
integral(g(x)-f(x),(x,0,1))
q(x)=f(x)-g(x)
solve(q(x),x)

1/3 [x == 0, x == 1]

f(x)=e^(x)
g(x)=x
plot(f(x), 0, 1,color='green', legend_label='f(x)=e^(x)')+ plot(g(x),0 ,1,color='red', legend_label='g(x)=x')
integral(f(x)-g(x),(x,0,1))

e - 3/2


f(x)=x/(sqrt(x^2+1))
g(x)=x^4-x
plot(f(x), 0, 1.4,gridlines=True, frame=False)+ plot(g(x),0 ,1.4)
integral(f(x)-g(x),(x,0,1.18))
q(x)=f(x)-g(x)
solve(q(x),x)
find_root(q(x), 1, 1.2)

0.7853869527188246 [x == 0, x^3 == (sqrt(x^2 + 1) + 1)/sqrt(x^2 + 1)] 1.1807757031062647

solve(x^3-7*x^2+14*x-8, x)

[x == 2, x == 4, x == 1]


find_root(cos(x)-x^2+1, 1,2)

1.1765019399018324

f(x)=cos(x)
g(x)=x^2-1
plot(f, -2, 2,gridlines='minor', frame=True, color='red')+ plot(g, -2 ,2)
integral(2*(f-g),0,1.18)

3.113857358149374

plot(atan(x)/(2+e^x),x)


s = 'sum (5*(-2/3)^(n-1),n,1,inf), simpsum'
SR(sage.calculus.calculus.maxima(s))

3

s = 'sum (1/n^2,n,1,inf), simpsum'
SR(sage.calculus.calculus.maxima(s))

1/6*pi^2

s = 'sum (1/2^n,n,1,inf), simpsum'
SR(sage.calculus.calculus.maxima(s))

1


s = 'sum (1/2^n,n,1,5), simpsum'
SR(sage.calculus.calculus.maxima(s))
numerical_approx(maxima(s), digits=10)

31/32 0.9687500000


print " Plot Slope Field "
t , v = var ("t,v")
pt = plot_slope_field ( 9.8-0.2*v , (t , 0 , 10) ,(v , 40 , 55) , plot_points = \
20)
show ( pt , aspect_ratio ='0.5')

Plot Slope Field


t = var('t')
x = function('x', t)
de = lambda v: diff(v,t) == 9.8-0.2*v
desolve(de(x),[x,t]);

(_C + 49*e^(1/5*t))*e^(-1/5*t)

print " Plot Slope Field "
t , p = var ("t,p")
pt = plot_slope_field ( 0.5*p-450 , (t , 0 , 10) ,(p , 800 , 1000) , plot_points = \
20)
show ( pt , aspect_ratio ='0.02')

Plot Slope Field


t = var('t')
x = function('x', t)
de = lambda p: diff(p,t) == 0.5*p-450
desolve(de(x),[x,t]);


integrate(e^x/(3-e^(2*x)),x)

-1/6*sqrt(3)*log(-(sqrt(3) - e^x)/(sqrt(3) + e^x))

plot(e^x/(3-e^(2*x)),(x,0.54,0.56))


integrate(1/(x**2),(x,1,oo))

1

k(x)=cos(x)
l(x)=(x^2)-1
solve(k(x)==l(x),x)

[x == -sqrt(cos(x) + 1), x == sqrt(cos(x) + 1)]
find_root(k(x)-l(x),-1.18,1.18)

1.1765019399018324

f(x)=3*x
g(x)=(x^3)-x
plot((f(x),g(x)),x,0,3,figsize=5,ymin=-2,ymax=10,fill=True)


solve(f(x)==g(x),x)

[x == -2, x == 2, x == 0]
integral(f(x)-g(x),x,0,2)

4

u = var('u')
f =u
g = u^2
a=revolution_plot3d(f, (u,0,2), show_curve=True, opacity=0.7).show(aspect_ratio=(1,1,1))
b=revolution_plot3d(g, (u,0,2), show_curve=True, opacity=0.7).show(aspect_ratio=(1,1,1))
(a+b).show()

3D rendering not yet implemented
3D rendering not yet implemented
Error in lines 6-6 Traceback (most recent call last): File "/cocalc/lib/python2.7/site-packages/smc_sagews/sage_server.py", line 1013, in execute exec compile(block+'\n', '', 'single') in namespace, locals File "", line 1, in <module> TypeError: unsupported operand type(s) for +: 'NoneType' and 'NoneType'
var('u,f')
f(u)=u^2+1
plot(f(u),u,0,3,figsize=3,fill=True)

(u, f)

var('u')
f(u)=sqrt(u)
revolution_plot3d(f(u),(u,0,1),show_curve=True,opacity=7)

u
3D rendering not yet implemented
2*pi*integral(x*(1-x^2),x,0,1)

1/2*pi