第三单元用python学习微积分（二十二）功、平均值、概率 ...

酹月誰邊 · 发表于 2022-4-11 15:07:04

本文内容来自于学习麻省理工学院公开课：单变量微积分-数值积分-网易公开课
用python解方程组_星-耀的博客-CSDN博客_python解方程组
做功_百度百科
首先老师阐释了这个：热量 = 能量 = 做功
热量的单位有千卡(kcal) 能量的单位有千瓦时( kwh ) 做功的单位有焦耳( j )

一、概率的经典问题，投飞标问题，问如果一个人向靶心投飞镖，击中的次数与  (c是常数,假设为1)成正比，
有多大的概率会射到靶子旁边的人。至于这个计算有没有实际意义呢？这里老师提到，二战时有人研究德国的V-2导弹瞄准伦敦发射，会击中哪里，用的就是这个函数（结果相近！）。
这里击中的概率与距离靶心的长度 r 相关，类正态分布
from sympy import *
import numpy as np
from scipy.optimize import fsolve
import matplotlib.pyplot as plt

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
ax.spines[&#39;left&#39;].set_position(&#39;zero&#39;)
ax.spines[&#39;bottom&#39;].set_position(&#39;zero&#39;)
ax.spines[&#39;right&#39;].set_color(&#39;none&#39;)
ax.spines[&#39;top&#39;].set_color(&#39;none&#39;)
ax.xaxis.set_ticks_position(&#39;bottom&#39;)
ax.yaxis.set_ticks_position(&#39;left&#39;)
ax.set_aspect(1)

def DrawXY(xFrom,xTo,steps,expr,color,label,plt):
yarr = []
xarr = np.linspace(xFrom ,xTo, steps)
for xval in xarr:
      yval = expr.subs(x,xval)
      yarr.append(yval)
y_nparr = np.array(yarr)
plt.plot(xarr, y_nparr, c=color, label=label)

def DrawInt(xFrom,xTo,steps,expr,color,plt, label=&#39;&#39;):
if(xFrom < 0 and xTo < 0):
      DrawIntNegative(xFrom,xTo,steps,expr,color,plt, label)
else:
      if(xFrom > 0 and xTo > 0):
         DrawIntPostive(xFrom,xTo,steps,expr,color,plt, label)
      else:
         DrawIntNegative(xFrom,0,steps,expr,color,plt, label)
         DrawIntPostive(0,xTo,steps,expr,color,plt, label)

def DrawIntNegative(xFrom1,xTo1,steps,expr,color,plt, label=&#39;&#39;):
xFrom = 0 - xTo1
xTo = 0 - xFrom1
width = (xTo - xFrom)/steps
xarr = []
yarr = []
area = 0
xprev = xFrom
yvalAll =  0
xarr.append(0)
yarr.append(0)
for step in range(steps):
      yval = expr.subs(x,xprev)
      area += width * yval
      xarr.append(0-xprev)
      yarr.append(0-area)
      xprev= xprev + width
plt.plot(xarr, yarr, c=color, label =label)

def DrawIntPostive(xFrom,xTo,steps,expr,color,plt, label=&#39;&#39;):
width = (xTo - xFrom)/steps
xarr = []
yarr = []
area = 0
xprev = xFrom
yvalAll =  0
xarr.append(0)
yarr.append(0)
for step in range(steps):
      yval = expr.subs(x,xprev)
      area += width * yval
      xarr.append(xprev)
      yarr.append(area)
      xprev= xprev + width
plt.plot(xarr, yarr, c=color, label =label)

def DrawParabola(xFrom,xTo,steps,expr,color,plt, label=&#39;&#39;):
width = (xTo - xFrom)/steps
xarrRect = []
yarrRect = []
area = 0
xprev = xFrom
yvalAll =  0
for step in range(int(steps/2)):
      x1 = xprev
      x2 = xprev+width
      x3 = xprev+width+width
      y1 = expr.subs(x,x1)
      y2 = expr.subs(x,x2)
      y3 = expr.subs(x,x3)

      def solveFuncsInDrawParabola(paramlist):
         a,b,c =paramlist[0],paramlist[1],paramlist[2]
         return [a*x1**2+b*x1+c-y1,
                  a*x2**2+b*x2+c-y2,
                  a*x3**2+b*x3+c-y3]
      s=fsolve(solveFuncsInDrawParabola,[0,0,0])
      a,b,c = s[0],s[1],s[2]
      expr1 = a*x**2+b*x+c
      DrawXY(x1,x3,20,expr1,color,&#39;&#39;,plt)
      xarrRect.append(x1)
      xarrRect.append(x1)
      xarrRect.append(x3)
      xarrRect.append(x3)
      yarrRect.append(y1)
      yarrRect.append(0)
      yarrRect.append(0)
      yarrRect.append(y3)
      area += 2*width * (y1+4*y2 + y3)/6
      plt.plot(xarrRect, yarrRect, c=color, label = &#39;&#39;)
      plt.plot([x2,x2], [0,y2], &#39;gray&#39;, label = &#39;&#39;)

      xprev = x3
print(&#39;============================&#39;)
if len(label)!=0:
      print(label)
print(&#39;============================&#39;)
print(&#39;width = &#39;, width)
print(&#39;area = &#39;,area)
areaFinal = (integrate(expr, (x,xFrom,xTo)))
print (&#39;area final = &#39;,areaFinal)
print (&#39;ave final = &#39;, areaFinal / (xTo - xFrom))




def DrawTrapezoidal(xFrom,xTo,steps,expr,color,plt, label=&#39;&#39;):
width = (xTo - xFrom)/steps
xarrRect = []
yarrRect = []
area = 0
xprev = xFrom
yvalAll =  0
for step in range(steps):
      yval1 =expr.subs(x,xprev)
      yval2 = expr.subs(x,xprev + width)
      xarrRect.append(xprev)
      xarrRect.append(xprev)
      xarrRect.append(xprev + width)
      xarrRect.append(xprev + width)
      xarrRect.append(xprev)
      yarrRect.append(0)
      yarrRect.append(yval1)
      yarrRect.append(yval2)
      yarrRect.append(0)
      yarrRect.append(0)
      area += width * (yval1+yval2)/2
      plt.plot(xarrRect, yarrRect, c=color)
      xprev= xprev + width
print(&#39;============================&#39;)
if len(label)!=0:
      print(label)
print(&#39;============================&#39;)
print(&#39;width = &#39;, width)
print(&#39;area = &#39;,area)
areaFinal = (integrate(expr, (x,xFrom,xTo)))
print (&#39;area final = &#39;,areaFinal)
print (&#39;ave final = &#39;, areaFinal / (xTo - xFrom))

def DrawRects(xFrom,xTo,steps,expr,color,plt, label=&#39;&#39;):
width = (xTo - xFrom)/steps
xarrRect = []
yarrRect = []
area = 0
xprev = xFrom
yvalAll =  0
for step in range(steps):
      yval = expr.subs(x,xprev + width)
      xarrRect.append(xprev)
      xarrRect.append(xprev)
      xarrRect.append(xprev + width)
      xarrRect.append(xprev + width)
      xarrRect.append(xprev)
      yarrRect.append(0)
      yarrRect.append(yval)
      yarrRect.append(yval)
      yarrRect.append(0)
      yarrRect.append(0)
      area += width * yval
      plt.plot(xarrRect, yarrRect, c=color)
      xprev= xprev + width
      yvalAll += yval
print(&#39;============================&#39;)
if len(label)!=0:
      print(label)
print(&#39;============================&#39;)
print(&#39;width = &#39;, width)
print(&#39;ave = &#39;, yvalAll / steps)
print(&#39;area = &#39;,area)
areaFinal = (integrate(expr, (x,xFrom,xTo)))
print (&#39;area final = &#39;,areaFinal)
print (&#39;ave final = &#39;, areaFinal / (xTo - xFrom))

def TangentLine(exprY,x0Val,xVal):
diffExpr = diff(exprY)
x1,y1,xo,yo = symbols(&#39;x1 y1 xo yo&#39;)
expr = (y1-yo)/(x1-xo) - diffExpr.subs(x,x0Val)
eq = expr.subs(xo,x0Val).subs(x1,xVal).subs(yo,exprY.subs(x,x0Val))
eq1 = Eq(eq,0)
solveY = solve(eq1)
return xVal,solveY

def DrawTangentLine(exprY, x0Val,xVal1, xVal2, clr, txt):
x1,y1 = TangentLine(exprY, x0Val, xVal1)
x2,y2 = TangentLine(exprY, x0Val, xVal2)
plt.plot([x1,x2],[y1,y2], color = clr, label=txt)

def Newton(expr, x0):
ret = x0 - expr.subs(x, x0)/ expr.diff().subs(x,x0)
return ret

def DrawText(x,y, text, fsize=9, clr = &#39;black&#39;):
plt.text(x,y, text, fontsize=fsize, color=clr)

c = 1
x = symbols(&#39;x&#39;)
expr = c*((np.e)**(-x**2))

DrawXY(-2,2,100,expr,&#39;blue&#39;,&#39;c*((np.e)**(-x**2))&#39;,plt)
plt.plot([1,1],[0,expr.subs(x, 1)], &#39;red&#39;, label=&#39;r1&#39;)
plt.plot([1.5,1.5],[0,expr.subs(x, 1.5)], &#39;green&#39;, label=&#39;r2&#39;)

plt.legend(loc=&#39;lower left&#39;)
plt.show()

上图蓝色曲线围绕 y 轴旋转一周可以形成击中靶子的概率模型，而我们需要考虑的则是r1到r2之间的部分。
看看顶视图（从上向下看）
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
ax.spines[&#39;left&#39;].set_position(&#39;zero&#39;)
ax.spines[&#39;bottom&#39;].set_position(&#39;zero&#39;)
ax.spines[&#39;right&#39;].set_color(&#39;none&#39;)
ax.spines[&#39;top&#39;].set_color(&#39;none&#39;)
ax.xaxis.set_ticks_position(&#39;bottom&#39;)
ax.yaxis.set_ticks_position(&#39;left&#39;)
ax.set_aspect(1)
x = symbols(&#39;x&#39;)

expr = (1 - x**2)**0.5

DrawXY(-1,1,100,expr,&#39;red&#39;,&#39;r1&#39;,plt)

expr = -(1 - x**2)**0.5
DrawXY(-1,1,100,expr,&#39;red&#39;,&#39;&#39;,plt)

expr1 = (1.5*1.5 - x**2)**0.5

DrawXY(-1.5,1.5,100,expr1,&#39;green&#39;,&#39;r2&#39;,plt)

expr1 = -(1.5*1.5 - x**2)**0.5
DrawXY(-1.5,1.5,100,expr1,&#39;green&#39;,&#39;&#39;,plt)

plt.legend(loc=&#39;lower left&#39;)
plt.show()

这里我给的c=1，老师说这个常量等于几并不重要，因为会被消掉....

想象蓝色曲线绕y轴旋转，左侧前视图，右侧顶视图

使用壳层法计算r1到r2之间的体积：第三单元用python学习微积分（二十二）功、平均值、概率 ... 第4张图片

这里老师一眼看出答案：第三单元用python学习微积分（二十二）功、平均值、概率 ... 第5张图片

显然这个c确实被削掉了。
同时，把所有的半径都考虑进去，那击中的概率为：第三单元用python学习微积分（二十二）功、平均值、概率 ... 第9张图片

这时老师添加了新的假设：
这个投掷标箭的人在7岁时，靶子的半径为 a , 他中靶的概率为 0.5

fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
ax.spines[&#39;left&#39;].set_position(&#39;zero&#39;)
ax.spines[&#39;bottom&#39;].set_position(&#39;zero&#39;)
ax.spines[&#39;right&#39;].set_color(&#39;none&#39;)
ax.spines[&#39;top&#39;].set_color(&#39;none&#39;)
ax.xaxis.set_ticks_position(&#39;bottom&#39;)
ax.yaxis.set_ticks_position(&#39;left&#39;)
ax.set_aspect(1)
x = symbols(&#39;x&#39;)

a = 1
expr = (a - x**2)**0.5

DrawXY(-a,a,100,expr,&#39;red&#39;,&#39;r = a&#39;,plt)

expr = -(a - x**2)**0.5
DrawXY(-a,a,100,expr,&#39;red&#39;,&#39;&#39;,plt)

r = symbols(&#39;r&#39;)
xminExpr = r*cos(-30*2/180*np.pi)

expr2 = -((2*a)**2 - x**2)**0.5
x1Min = float(xminExpr.subs(r,2*a))
DrawXY(x1Min,2*a,100,expr2,&#39;gray&#39;,&#39;area where the people stand&#39;,plt)

xminExpr = r*cos(-30*2/180*np.pi)
expr3 = -((3*a)**2 - x**2)**0.5
x2Min=float(xminExpr.subs(r,3*a))
DrawXY(x2Min,3*a,100,expr3,&#39;gray&#39;,&#39;&#39;,plt)

plt.plot([2*a, 3*a],[0,0], &#39;gray&#39;, &#39;&#39;)
plt.plot([x1Min, x2Min],[expr2.subs(x,x1Min),expr3.subs(x,x2Min)], &#39;gray&#39;, &#39;&#39;)

DrawText(3.2*a, 0, &#39;3 o&#34;clock&#39;,12)
DrawText(x1Min, expr3.subs(x,x2Min) - 0.2, &#39;5 o&#34;clock&#39;,12)

plt.legend(loc=&#39;upper left&#39;)
plt.show()

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第10张图片

由于上图灰色这段环占整个完整的圆环第三单元用python学习微积分（二十二）功、平均值、概率 ... 第11张图片

(....类似2个小时比12个小时...)，所以有
本问题就是求第三单元用python学习微积分（二十二）功、平均值、概率 ... 第12张图片

由已知：

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第13张图片

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第14张图片

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第15张图片

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第16张图片

再添加了各种预设后，这个题目的解，站在标靶旁的小人被击中的概率约等于 1%
这个问题是讨论加权的平均值或带权重的积分
由第三单元用python学习微积分（二十二）功、平均值、概率 ... 第17张图片

可知

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第18张图片

，而当 r=0时，权重为0。
二、数值积分
1、黎曼集合(Riemann sums) 。。。。。（老师语：效率低, 而且结果不令人满意）
这个就是把函数曲线下面到x轴按分为n个矩形相加得到面积的那个

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第19张图片

目标：平均值和累加y 来得到第三单元用python学习微积分（二十二）功、平均值、概率 ... 第20张图片

的近似值
左侧和：第三单元用python学习微积分（二十二）功、平均值、概率 ... 第21张图片

右侧和：

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第22张图片

用刚刚的权重函数做个例子
x = symbols(&#39;x&#39;)
expr = 2 *np.pi * x * np.e**(-x**2)
DrawXY(0,3,100,expr,&#39;blue&#39;,&#39;2 *np.pi * x * e**(-x**2)&#39;,plt)
DrawRects(0,3,20,expr,&#39;gray&#39;,plt, label=&#39;rect&#39;)
plt.legend(loc=&#39;upper left&#39;)
plt.show()

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第23张图片

2、梯形法。。。。。。（老师语：效率低, 而且结果不令人满意）
这个方法和黎曼和的方法类似，就是不简单的切分矩形，而是换成根据平局的间隔切分出不同的梯形，再累加。
x = symbols(&#39;x&#39;)
expr = 2 *np.pi * x * np.e**(-x**2)
DrawXY(0,3,100,expr,&#39;blue&#39;,&#39;2 *np.pi * x * e**(-x**2)&#39;,plt)
DrawTrapezoidal(0,3,20,expr,&#39;gray&#39;,plt, label=&#39;Trapezoidal&#39;)
plt.legend(loc=&#39;upper left&#39;)
plt.show()

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第24张图片

从结果看明显比黎曼和要好些....

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第25张图片

第三单元用python学习微积分（二十二）功、平均值、概率 ... 第26张图片

(黎曼左侧和＋黎曼右侧和) / 2
页面不能太大，剩下的请看Bullseye：第三单元用python学习微积分（二十二）功、平均值、概率（下）和数值积分 (2)

上一篇：完整版：全球著名咨询公司管理常用模型（值得收藏！）
下一篇：MIT 6.858: Computer Systems Security 计算机系统安全实验1

白金会员	积分	兔币	帖子
白金会员, 积分 3464, 距离下一级还需 1536 积分	3464	1897	1549
在线时间：0 小时	最后登录：2024-1-25

[问答] 第三单元用python学习微积分（二十二）功、平均值、概率 ...

关联主题

三一重工董监高“大撤退”：减持1亿砸出75

政府收购商品房来了！泰安征集已建成存量商

前英皇男星儿子生病入院，被老婆炮轰曝光冷

未来10天北方迎三轮降雨南方多地高温或破

市值蒸发3万亿！英伟达遭美国反垄断调查，

奥迪A4L销量不如奔驰C级和宝马3系，压力依

30万级别, 买运动计时选它准没错！

丰田与宝马的野心：氢燃料汽车将取代电动车

房价何时见底

暴跌后的油价何去何从?

何小鹏：中国新能源车企组球队 10年内说不

星二代，被扎堆送出国

恒大地产，又添诉讼！

存量房贷，调整的脚步声近了？

北京买房：理清思路，购房建议1294

美欧说“遗憾”，乌克兰发出“威胁”，蒙古

山东东平一公交车失控，家长逃生后发视频时

“断崖式衰老”发生在这两个岁数！Nature子

三一重工董监高“大撤退”：减持1亿砸出75

荣耀宣布即将发布首款 AI Agent 手机，Magi

悟空出世爆火，深圳老板入账15个亿

与操盘手深聊4小时，我们找到2024下半年经

成就非凡：家长该如何放手，让孩子书写自己

大S恋爱脑上头，消费汪小菲给具俊晔赚流量

新机：Mate70全系列亮相；红米K80 Pro长这

“ 到公司第一件事就是…”小米王腾晒游戏

iOS18新测试版推送达1.2GB ，新功能不多，

何小鹏：中国新能源车企组球队 10年内说不

乌军越境俄罗斯激战持续一周：沉默中行动，

小鹏汽车官宣！完成收购

财神驾到

一抹伤

哇哇的哭

绿林道的

willlin

[问答] 第三单元 用python学习微积分（二十二）功、平均值、概率 ...

关联主题

[问答] 第三单元用python学习微积分（二十二）功、平均值、概率 ...