python微元法計(jì)算函數(shù)曲線長度的方法

'''
計(jì)算曲線長度，根據(jù)線積分公式：
\int_A^Bf(x,y,z)dl，令積分函數(shù)為1，即計(jì)算曲線的長度
'''
import numpy as np
from mpl_toolkits.mplot3d import *
import matplotlib.pyplot as plt

## 求二維圓周長，半徑為1，采用參數(shù)形式
def circle_2d(dt=0.001,plot=True):
 dt = dt # 變化率
 t = np.arange(0,2*np.pi, dt)
 x = np.cos(t)
 y = np.sin(t)

 # print(len(t))
 area_list = [] # 存儲(chǔ)每一微小步長的曲線長度

 for i in range(1,len(t)):
  # 計(jì)算每一微小步長的曲線長度，dx = x_{i}-x{i-1}，索引從1開始
  dl_i = np.sqrt( (x[i]-x[i-1])**2 + (y[i]-y[i-1])**2 ) 
  # 將計(jì)算結(jié)果存儲(chǔ)起來
  area_list.append(dl_i)

 area = sum(area_list)# 求和計(jì)算曲線在t:[0,2*pi]的長度

 print("二維圓周長：{:.4f}".format(area))
 if plot:
  fig = plt.figure()
  ax = fig.add_subplot(111)
  ax.plot(x,y)
  plt.title("circle")
  plt.show()


## 二維空間曲線，采用參數(shù)形式
def curve_param_2d(dt=0.0001,plot=True):
 dt = dt # 變化率
 t = np.arange(0,2*np.pi, dt)
 x = t*np.cos(t)
 y = t*np.sin(t)

 # print(len(t))
 area_list = [] # 存儲(chǔ)每一微小步長的曲線長度

 # 下面的方式是循環(huán)實(shí)現(xiàn)
 # for i in range(1,len(t)):
 #  # 計(jì)算每一微小步長的曲線長度，dx = x_{i}-x{i-1}，索引從1開始
 #  dl_i = np.sqrt( (x[i]-x[i-1])**2 + (y[i]-y[i-1])**2 ) 
 #  # 將計(jì)算結(jié)果存儲(chǔ)起來
 #  area_list.append(dl_i)

 # 更加pythonic的寫法
 area_list = [np.sqrt( (x[i]-x[i-1])**2 + (y[i]-y[i-1])**2 ) for i in range(1,len(t))]

 area = sum(area_list)# 求和計(jì)算曲線在t:[0,2*pi]的長度

 print("二維參數(shù)曲線長度：{:.4f}".format(area))

 if plot:

  fig = plt.figure()
  ax = fig.add_subplot(111)
  ax.plot(x,y)
  plt.title("2-D Parameter Curve")
  plt.show()

## 二維空間曲線
def curve_2d(dt=0.0001,plot=True):
 dt = dt # 變化率
 t = np.arange(-6,10, dt)
 x = t
 y = x**3/8 - 4*x + np.sin(3*x)

 # print(len(t))
 area_list = [] # 存儲(chǔ)每一微小步長的曲線長度

 # for i in range(1,len(t)):
 #  # 計(jì)算每一微小步長的曲線長度，dx = x_{i}-x{i-1}，索引從1開始
 #  dl_i = np.sqrt( (x[i]-x[i-1])**2 + (y[i]-y[i-1])**2 ) 
 #  # 將計(jì)算結(jié)果存儲(chǔ)起來
 #  area_list.append(dl_i)

 area_list = [np.sqrt( (x[i]-x[i-1])**2 + (y[i]-y[i-1])**2 ) for i in range(1,len(t))]

 area = sum(area_list)# 求和計(jì)算曲線在t:[0,2*pi]的長度

 print("二維曲線長度：{:.4f}".format(area))

 if plot:
  fig = plt.figure()
  ax = fig.add_subplot(111)
  ax.plot(x,y)
  plt.title("2-D Curve")
  plt.show()

## 三維空間曲線，采用參數(shù)形式
def curve_3d(dt=0.001,plot=True):
 dt = dt # 變化率
 t = np.arange(0,2*np.pi, dt)
 x = t*np.cos(t)
 y = t*np.sin(t)
 z = 2*t

 # print(len(t))
 area_list = [] # 存儲(chǔ)每一微小步長的曲線長度

 for i in range(1,len(t)):
  # 計(jì)算每一微小步長的曲線長度，dx = x_{i}-x{i-1}，索引從1開始
  dl_i = np.sqrt( (x[i]-x[i-1])**2 + (y[i]-y[i-1])**2 + (z[i]-z[i-1])**2 ) 
  # 將計(jì)算結(jié)果存儲(chǔ)起來
  area_list.append(dl_i)

 area = sum(area_list)# 求和計(jì)算曲線在t:[0,2*pi]的長度

 print("三維空間曲線長度：{:.4f}".format(area))

 if plot:
  fig = plt.figure()
  ax = fig.add_subplot(111,projection='3d')
  ax.plot(x,y,z)
  plt.title("3-D Curve")
  plt.show()

if __name__ == '__main__':

 circle_2d(plot=True)
 curve_param_2d(plot=True)
 curve_2d(plot=True)
 curve_3d(plot=True)

二維圓周長：6.2830
二維參數(shù)曲線長度：21.2558
二維曲線長度：128.2037
三維空間曲線長度：25.3421