---

時隔兩天，前一篇文章閱覽量竟然還未破1000（555），但是為了能夠最快速度完成本周微分幾何作業(yè)，我加工1h整出來了本次文章的計算代碼。

上次文章我們解決的問題是給定曲線下曲率與撓率的計算，這次我們反其道而行之。

根據(jù)曲線論基本定理一節(jié)課內容的知識我們知道了在曲率與撓率確定的條件下，空間曲線是唯一的（在剛體變換等價類的意義下）。那么，當我們知道空間曲線的曲率與撓率，我們理應可以按照以下常微分方程求出對應的曲線的參數(shù)表示：
$$?\begin{array}{rl}{r}^{\prime }(s)=& t(s)\\ t^{\prime }(s)=& \kappa (s)?n(s)\\ n^{\prime }(s)=& ?\kappa (s)?t(s)+\tau (s)?b(s)\\ b^{\prime }(s)=& ?\tau (s)?n(s)\end{array}$$
用矩陣形式表示即為：
$$\begin{array}{r}\frac{d}{ds}?\begin{array}{c}r(s)\\ t(s)\\ n(s)\\ b(s)\end{array}?=?\begin{array}{cccc}0& 1& 0& 0\\ 0& 0& \kappa (s)& 0\\ 0& \kappa (s)& 0& \tau (s)\\ 0& 0& ?\tau (s)& 0\end{array}??\begin{array}{c}r(s)\\ t(s)\\ n(s)\\ b(s)\end{array}?\end{array}$$
對于非常系數(shù)線性微分方程，我們并沒有較好的解決方法，但是我們可以通過計算機進行求解。

我們舉下面一個習題作為例子：

---

(1)求曲率$\kappa (s)=a/(a^2+s^2)$($s$為弧長參數(shù))的平面曲線;
(2)求曲率$\kappa (s)=1/\sqrt{a^2?s^2}$($s$為弧長參數(shù))的平面曲線.

---

我們發(fā)現(xiàn)這里只需要求平面上的曲線，而不需要求空間上的曲線，因而我們只需要假定撓率$\tau =0$(撓率為0時，曲線為平面曲線)。

我們直接上代碼：

from sympy import * from sympy.plotting import plot3d_parametric_line,plot_parametric import re

我們這次要用的庫有上次使用過的符號計算庫sympy,以及其對應的繪制三維和二維圖像的兩個函數(shù)，以及進行數(shù)據(jù)清洗的re庫。

#自變量 s = symbols('s', real=True) # r的三個分量 r1 = Function('r1') r2 = Function('r2') r3 = Function('r3') #r=Matrix([r1,r2,r3]) # t的三個分量 t1 = Function('t1') t2 = Function('t2') t3 = Function('t3') #t=Matrix([t1,t2,t3]) # n的三個分量 n1 = Function('n1') n2 = Function('n2') n3 = Function('n3') #n=Matrix([n1,n2,n3]) # b的三個分量 b1 = Function('b1') b2 = Function('b2') b3 = Function('b3') #b=Matrix([b1,b2,b3])

我們用sympy定義出相應的自變量和函數(shù)。需要注意的是，我們上面列出的微分方程中$r,t,n,b$都是三維的列向量，因而總方程數(shù)要有更多，在這里我們也是用分量的形式來進行定義的。

# 曲率與撓率 a=1 k = a/(a**2+s**2)# 曲率 tao = 0 # 撓率

我們輸入題目中給定的曲率與撓率。(因為不會帶參數(shù)的微分方程求解，因而暫時取參數(shù)a為1)（會的朋友歡迎加我好友進行交流）

eq = [r1(s).diff(s) - t1(s),r2(s).diff(s) - t2(s),r3(s).diff(s) - t3(s),t1(s).diff(s) - k * n1(s),t2(s).diff(s) - k * n2(s),t3(s).diff(s) - k * n3(s),n1(s).diff(s) + k * t1(s) - tao * b1(s),n2(s).diff(s) + k * t2(s) - tao * b2(s),n3(s).diff(s) + k * t3(s) - tao * b3(s),b1(s).diff(s) + tao * n1(s),b2(s).diff(s) + tao * n2(s),b3(s).diff(s) + tao * n3(s)] con = {r1(0): 0, r2(0): 0, r3(0): 0, t1(0): 1, t2(0): 0, t3(0): 0, n1(0): 0, n2(0): 1, n3(0): 0, b1(0): 0, b2(0): 0, b3(0): 1}

我們將微分方程組和初值條件分別用列表和字典屬性進行儲存，

a = dsolve(eq, ics=con) b=a[0:3] print(str(b[0]),str(b[1]),str(b[2]),sep='\n')

用dsolve函數(shù)解出相應的方程組，并取出$r$的成分$r1,r2,r3$.
得到的結果為：
$$r_1(s)=asinh(s);r_2(s)=\sqrt{s^2+1}?1;r_3(s)=0.$$

我們還能對其曲線進行繪制：

R1=eval(re.findall(',.*\)',str(b[0]))[0][1:-1]) R2=eval(re.findall(',.*\)',str(b[1]))[0][1:-1]) R3=eval(re.findall(',.*\)',str(b[2]))[0][1:-1]) if R1==0:plot_parametric(R2,R3,(s,-100,100),title='R2-R3') elif R2==0:plot_parametric(R1, R3, (s, -100, 100),title='R1-R3') elif R3==0:plot_parametric(R1, R2, (s, -100, 100),title='R1-R2') else:plot3d_parametric_line(R1,R2,R3,(s,-100,100))

得到了解的積分曲線：

同樣的方法我們可以解答第二問的答案：
$$r_1(s)=s?\sqrt{1?s^2}/2+asin(s)/2;r_2(s)=s^2/2;r_3(s)=0.$$
其積分曲線的圖像為：

---

原創(chuàng)不易，歡迎大家一鍵三連！

---

文案：錦帆遠航
代碼：錦帆遠航
圖片：錦帆遠航 百度百科

總結

以上是生活随笔為你收集整理的用曲率，挠率反求曲线方程！（作业捷径篇 续集）的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內容還不錯，歡迎將生活随笔推薦給好友。