速途網10月10日消息（報道：李楠）為什么無論多高的雨滴落下來都有一個不大的末速度？普通物理課會告訴學生，因為雨滴在下落時會受到與速度正相關的空氣阻力，導致雨滴的下落末速度趨于某個定值。但這個阻力與速度的關系是怎樣的呢？一百七十年前，愛爾蘭物理學家斯托克斯對類似的問題產生了興趣，并提出了在小雷諾數條件下流體阻力與速度成正比的斯托克斯定律。斯托克斯當時又是如何推導的呢？

10月6日14時，《張朝陽的物理課》第二百二十五期開播，搜狐創始人、董事局主席兼首席執行官、物理學博士張朝陽坐鎮搜狐視頻直播間，從流體力學最基本的納維爾-斯托克斯方程（Navier-Stokes equations，下文簡稱NS方程）出發，“馬拉松”式地一步步推導出斯托克斯定律，課堂全程持續近5個小時。

回顧雨滴下落速度的求解

假設有一個密度為ρ_d的球形雨滴，它在密度為ρ_a的空氣中從靜止開始下落，這個過程中雨滴會受到三個力：向下的重力G、向上的空氣浮力F_{浮}和向上的空氣阻力F?，其中空氣阻力滿足斯托克斯定律，將它們依次寫出來是

注意看最后一個式子所描述的空氣阻力，它正是斯托克斯定律所給出的表達式，其中μ是空氣的粘滯系數，R是雨滴半徑，v是雨滴相對空氣的速度，這三個量的一次方相乘再乘以6π就得到了球形物體在不可壓縮定常流體中所受到的阻力。本節課的任務，正是從基本的流體力學最基本的NS方程出發，推導出斯托克斯定律。不過在進入正題前，先看看斯托克斯定律在雨滴下落情景中發揮的作用。

(對雨滴的受力分析)

由剛剛分析出的三個力可以寫下雨滴的運動方程

為了化簡方程，可以將系數用兩個字母代替

這樣就能把雨滴運動方程化簡為

可以猜測這個微分方程的解具有如下形式

代回原方程，得到

再考慮到初始條件，即雨滴在零時刻靜止，可以得到

所以

和自由落體的速度線性增長不同的是，在考慮空氣阻力的情況下，雨滴的下落速度在t→∞時無限逼近一個有限值，它就是雨滴下落的最終速度

為了更加直觀地感受空氣阻力的效果，張朝陽取了一些具體數據代入

算出雨滴下落的最終速度為10m/s，和大家日常生活的經驗十分吻合，說明斯托克斯定律對空氣阻力的描述相當地好。

以上計算只能是針對0.3mm的小雨滴，或者說，它更接近云層中的水霧?，F實中，云層中的水霧會匯聚成毫米量級的大水珠再落下來，水珠下落時受到空氣阻力還會發生形變，所以斯托克斯定律并不能完美地描述雨滴實際受到的阻力。事實上，此時占主導地位的阻力是一個與速度平方成正比的壓降阻力。不過在本堂課，我們更關心斯托克斯定律的導出，至于斯托克斯定律在雨滴下落問題的適用性，可以認為，在雨滴很小的情況下，它總是近似成立的。

關于斯托克斯力并不是大部分情況下雨滴所受空氣阻力的論述，是課后由孫昌璞院士友情指正的，特此致謝。

尋根問底：從流體力學最基礎的NS方程出發

上一節求解了雨滴在受到與速度成正比的空氣阻力下的運動過程，那為什么空氣阻力可以具有與速度的線性關系？斯托克斯定律為什么能在小雨滴在空氣中下落的情形下成立？它又是否是對所有的流體導致的阻力都適用的呢？

帶著這樣的疑問，張朝陽開啟了本場“馬拉松”課程的正式內容。球形雨滴在空氣中下落，這個過程具有以下落方向為軸的旋轉對稱性，而雨滴本身又是球對稱的，所以自然地可以想到要把這些拉普拉斯算子在柱坐標或者球坐標下寫開。張朝陽選擇先以下落方向為極軸，建立了如下所示的球坐標系。

(建立球坐標系)

建立好坐標系后，就可以寫出空氣的各點性質關于坐標的函數，比如各個位置上的壓強p(r,θ,?)和速度v(r,θ,?)。物理學家通常用場的概念描述這種隨空間分布的性質，比如壓強就是一個標量場，速度就是一個矢量場。一般的場還會帶有時間坐標t，但這里研究的是達到穩定狀態后的現象，所以可以忽略坐標t。

流體力學中的壓強場和速度場，就好像質點力學中的力和速度。上一節課將雨滴看成一個質點，用牛頓第二定律寫下了它的運動方程。本節要研究的是，空氣作為一個黏滯流體，它會給在其中做相對運動的雨滴一個什么樣的阻力（有時也叫流體曳力，drag）。為了求解空氣阻力，需要先求出空氣自己的運動狀態，而空氣作為一個流體，它的運動方程可以由NS方程給出

NS方程本質上就是牛頓第二定律在流體上的應用，方程的左邊類似于ma，其中加速度包含速度場對時間的偏導項，和速度場對空間偏導后再對時間求偏導，方程的右邊類似于F，其中包含了壓強梯度所導致的正向壓力差，和流體運動的粘性所導致的粘滯力。

現在假設流體已經達到了穩態，所以NS方程左邊第一項，也就是關于時間求偏導的項為0。而方程左邊的第二項是一個速度場關于空間分布變化率的項，可以進一步假設流體微元在隨著流線運動的過程中速度的空間變化率是緩慢的，也就是近似認為NS方程左邊第二項為0。經過穩態和空間緩變的這樣兩個假設，NS方程被簡化為了一個線性的微分方程

類比電動力學，巧妙引入渦度

方程的左邊是一階導，右邊是二階導，有沒有可能將右邊降階為一階導呢？回憶矢量微積分中有這樣一條公式

而流體的質量守恒給出

這里引入第三個假設：在速度比較小時空氣是不可壓縮的流體。這個假設在流速小于0.3倍音速的情形下通常是成立的。在不可壓縮假設下，空氣密度是一個常量，所以質量守恒導出速度場是無散的

將(1)-(3)式聯立，得到

這個等式的右邊看起來還是二階導，但與(1)式不同的是，這里的nabla算子▽是依次以叉乘的形式作用在后面的矢量上的，而(1)式是兩個nabla算子以點乘成拉普拉斯算子的形式作用到速度矢量上，前者的兩次求導操作是容易拆分的，后者要拆分的話比較困難，需要先作用一次導出二階張量再求散度來縮并回一階矢量。受到(4)式的啟發，可以定義一個中間變量

它是速度場的旋度，稱為渦度(vorticity)。這樣就能把(1)式右邊也降階為一階導的形式

這里巧妙地引入了渦度場來代替速度場，從而讓NS方程兩邊的求導階數相等。而巧合的是，電動力學中也曾有類似的操作。張朝陽將電動力學的量類比成三層樓，第一層是電勢Φ和磁矢勢A，第二層是電場E和磁場B，第三層是電荷密度ρ和電流密度j。每上升一層就對應求一次導，例如

比較(5)和(7)，不難看出對于同樣是無散的速度場和磁場，都可以利用它們的旋度，也就是渦度和電流密度，來作為中間變量進而簡化方程。

(類比電動力學的“三層樓”)

經過三個假設后，NS方程簡化為了(6)式，它是一個關于壓強場和渦度場的方程，有沒有可能將這兩個場分開呢？首先，注意到任何一個矢量場的旋度無散，即

所以對(6)式兩邊求散度，得到只關于壓強場的泊松方程

單獨一個(8)式只給出了壓強場的信息，顯然它損失了渦度場的信息。張朝陽打比方道，這就像是有兩個數3和4，將它們乘在一起得到了12，但單獨一個12是還原不出3和4的，它也可能對應2和6。類似的道理，對方程(6)求散度得到(8)式是會損失信息的，還需要再從方程(6)導出渦度場的信息。同樣地，注意到任意一個標量場的梯度無旋，即

所以對(6)式兩邊求旋度，得到只關于渦度場的方程

這里再次利用了矢量微積分的公式(2)。另外，渦度是速度場的旋度，所以它是無散的。因此渦度場也滿足泊松方程

靈活切換坐標系，轉化渦度場的矢量泊松方程

上一節將NS方程簡化成了兩個泊松方程，一個是關于壓強標量的式(8)，另一個是關于渦度矢量的式(9)。標量的泊松方程在之前的課上有相當多的處理經驗，但矢量的泊松方程是之前從沒遇到過的，為此需要先研究一下矢量被拉普拉斯算符作用后得到的是什么。

根據流體運動的軸向對稱性，速度矢量應該只有徑向分量v_r和極角方向的分量v_θ，而沒有方位角方向上的分量v_?，并且各方向的速度分量對?的求導都為0?；趯λ俣葓龅妮S對稱性的認知，可以在球坐標下寫出渦度

可見在流體速度有軸對稱性的情況下，渦度是一個只有方位角方向分量的矢量。由此(9)式變成了

這是一個矢量場的泊松方程。首先來計算第一次nabla算符作用后的結果，它將被作用的矢量沿不同方向求導，但對求導方向的基矢和被作用后的矢量的基矢這兩個基矢而言做了張量積，張量積既不是點乘也不是叉乘，而是把兩個基矢直接放在一起作為二階張量的基底，以三維空間來看，它包含了3×3=9個系數和基底。用?代表矢量的張量積，可以寫成

(12)式的兩項分別將nabla算符作用在了系數和基矢上，它的第一項再被nabla算符點乘后是

上式第二項中被大括號標出的部分為0，因為球坐標的散度公式為

而基矢\vec{e}_?就相當于g_r=g_θ=0，g_?=1的一個矢量，代入散度公式可知它等于0。

(12)式的第二項涉及到直接對一個矢量求“梯度”得到二階張量，展開來寫是

第二個等式新定義了一個矢量\vec{e}_ρ，它是從\vec{e}_?對?求導出來的。方位角基矢與r和θ無關，只會隨著?變化，且變化的矢量是指向極軸的一個向內的矢量，類似于柱坐標下的徑向矢量。

(柱坐標的切面)

為了方便對這個二階張量求散度，接下來把它切換到直角坐標系下，將\vec{e}_?和\vec{e}_ρ用\vec{i}和\vec{j}展開

直角坐標系下的坐標和球坐標之間有這樣的關系

因此，

或者更顯式地，把它寫成一個2×2的矩陣形式（因為ω沒有z方向分量，可以把問題簡化在xy平面上）

至此，(12)式的第二項化為了直角坐標系下的形式，無論從張量積的形式看還是從矩陣的形式看，它都是一個二階張量，可以用字母F上加兩個箭頭來表示這是一個二階張量。對這個二階張量再求一次散度

最后一步把第一項又寫回到球坐標的形式，第二項則引入了矢量

它是從極軸出發指向空間中某個位置的位矢。注意到f是一個和?無關的函數，雖然它的分子上有ω_?，但從(10)式可以看出ω_?是和?無關的，這一切都來源于流速場具有繞軸旋轉的對稱性，因此▽f在xy平面上的分量只能是平行于ρ的，所以(14)式中被大括號包住的部分為0。

將(13)和(14)式的結果代入(11)式的矢量泊松方程，可以將基底\vec{e}_?提出到拉普拉斯算符外面，從而得到只和系數ω_?有關的泊松方程

巧借氫原子球諧函數，求解球坐標下的壓強場和渦度場

將壓強場的泊松方程(8)式和渦度系數的方程(15)寫在球坐標下，并且注意到它們在?方向上的對稱性，可以得到

壓強場的方程比較簡單，先求解它。假設壓強場方程的解能分離成與徑向有關的部分以及與極角方向有關的部分

代入(16.a)，得到

上式第一個部分只與f(r)有關，第二個部分只與g(θ)有關，所以它們應該都等于一個常數，注意到這個方程和《張朝陽的物理課（第一卷）》中氫原子薛定諤方程的形式一模一樣，不妨就借鑒那里的思路，令這個常數等于角量子數l(l+1)。對于徑向部分，取(17)式的第一部分和第三部分做比奈（Binet）變換

令ψ=rf，上式可化簡為

這個方程有兩個解系，一個是ψ~r^(l+1)，一個是ψ~r^(-l)。

對于極角方向部分，取(17)式的第二部分和第三部分

令x=cos(θ),h(x)=g(θ)

這個方程的解正是勒讓德多項式

結合徑向和角向的結論，可以寫下壓強場的通解

剩下的問題是確定展開系數A和B。為此，張朝陽先從量綱的角度分析了它與哪些物理量有關，對于流體產生的阻力，不難想象它和在其中運動的球狀物體的半徑R、運動速度v?、以及流體粘滯系數μ有關

一個樸素的想法是，先將這些量的一次方相乘起來，看看量綱是什么

和壓強的量綱正好差L2，因此可以大膽地假設，(19)式中只包含r^(-2)的項。當然這里會有一個小問題，R/r是一個無量綱量，完全可以在表達式中再乘上這一無量綱量來維持量綱的平衡。針對這個不確定性，數學家們已經證明了，在給定足夠的邊界條件下，微分方程的解具有唯一性，因此只要求出了其中一個滿足邊界條件的解，就可以認為得到了完整的解。

經過這樣的量綱分析，不難發現只有l=1的情況才滿足假設

再來看(16.b)的渦度場方程，徑向部分的分析和壓強場是一樣的，而角向部分相比(18)式會多出一項

多出來的第三項，正好是球諧函數Y在m=1時方位角分量e^(im?)被拉普拉斯算子作用后出現的項

由此可以確定，(21)式的解是連帶勒讓德多項式

所以渦度場的通解為

渦度是速度的散度，它具有量綱

如果和壓強場一樣，也樸素地假設它含有R的一次項和v?的一次項，那么可以推斷(22)式只含有r^(-2)的項，也就是l=1。由于

因此

降一層臺階：用斯托克斯流函數代入邊界條件

上一節利用之前的物理課解氫原子薛定諤方程的經驗，成功得出了壓強場(20)和渦度場(23)，但它們各自都有一個待定系數，需要用邊界條件來定出。然而渦度的邊界條件是難以給定的，因為根據(10)式，它和速度的兩個分量vr、vθ的導數都有關系。回憶引入渦度時，出發點是為了將(1)式右邊的二階導變成一階導，相當于視角升了一個臺階，從而降低微分的階數。那么有沒有可能再轉化一次視角，把vr和vθ統一成另一個中間量，以便于引入邊界條件呢？

答案是可行的。在引入渦度時，曾將速度場類比成磁場，因為它們都具有無散的性質，可以通過(2)式將拉普拉斯算子轉化為兩個散度的依次作用。無散的矢量場還有一個性質，就是它們總可以寫成另一個矢量場的旋度，正是這一點允許從磁場定義一個磁矢勢。同理，可以定義一個矢量場\vec{Ψ}，它的旋度等于速度場，從而將視角降一個臺階

(張朝陽類比電動力學引入斯托克斯流函數)

在球坐標下寫開各分量

和一般的在球坐標下求散度的情形不同的是，最后一步重新定義了Ψ?，把rsin(θ)吸收了進去，并簡記為Ψ，稱為斯托克斯流函數。再次利用軸對稱的性質，可以知道v沒有?方向的分量，且Ψ是一個與?無關的函數，所以

也就是

到這一步，就把兩個速度分量統一成了一個斯托克斯流函數。將它們代入(10)式

再代入(23)式，就可以從渦度的場方程轉換到斯托克斯流函數的方程

這個偏微分方程相比(16)式要復雜得多，它很難直接地分離出徑向部分和角向部分，但可以先假設它是可分量變量的，從而猜測它應該滿足什么樣的形式，一個最簡單的猜解是

代回檢驗，不難發現這個形式能非常方便地配平(25)式中有關角向的部分的冪次

由此得到徑向部分的方程

這是一個變系數非齊次線性微分方程，可以確定它的通解是這些冪函數的和

代回原式可以定出a=D?/2。

至此，終于可以考慮速度的邊界條件了。在本問題中，有兩個邊界條件需要考慮，一個是無窮遠處的流體應當沒有受到球狀雨滴的影響，所以保持靜止；另一個是雨滴表面的流體相對雨滴靜止，稱為no-slip condtion，這在小雷諾數情形下是普遍成立的。把(26)代入(24)，并要求無窮遠處速度趨于0，可以定出b=0。再要求球面處有\vec{v}(R)=\vec{v}?，對應r分量和θ分量分別是

由此可以定出

所有的待定系數都已經得到了確定，最終結果是

在用邊界條件定下速度場和渦度場后，就可以聯立(6)和(20)來確定壓強場的待定系數B?了。

比較上下兩式，可以定出

所以壓強場為

當然，空氣應該還有個大氣壓強p?作為背景壓強，但加上一個常數并不影響壓強場的梯度，因此這里可以忽略大氣壓。

計算流體作用在球狀物體上的斯托克斯阻力

在本文最開始介紹完整的NS方程時，就說道等號右側代表流體受力的兩項分別是壓強梯度導致的正向壓力差和流體粘性導致的粘滯力。流體對雨滴的作用力也一樣有這兩個來源，所以對于流體的應力張量，它的表達式是

最終作用在球體上的力，就是這個二階的應力張量與球面的法向量\vec{e}?點乘得到的一階矢量在整個球面上的積分。正是由于二階張量有非對角項，點乘出的力并不一定垂直于法向，而是會有一定的切向分量。

(應力張量與面法向量的點乘)

對于(28)式的第一項，將張量基底寫開并與球面法向點乘是

根據軸對稱性，壓強在球面上的積分只剩下沿極軸z方向的分量，其他方向的分量會互相抵消，所以可以只對壓強的z向分量做積分

同理，(28)式的第二項與球面法向點乘是

上式最后一步代入r=R時，r方向上的速度梯度剛好是一個極值點，所以不貢獻粘滯力，只有θ方向貢獻出來一個切向的壓強。將這個壓強的z分量在球面上積分

(28)式的第三項與球面法向點乘是

這說明這一項在球面不會貢獻壓強。其中第二個等號運用了球坐標單位基矢的求導關系

將以上三項加起來，就得到了斯托克斯定律

回顧整堂課的推導過程，雖然要解決的是流體力學問題，但其中貫穿了各個領域的知識，比如從電動力學的矢量微積分引入了渦度、速度、斯托克斯流函數的“三層樓”，借用量子力學中的球諧函數解出了渦度角向分量的泊松方程。當然，從歷史進程來看，數學家們先是在研究流體力學時發展了這些工具，再到后來發展電動力學和量子力學時，人們就可以直接使用這些已經發展好的工具。而《張朝陽的物理課》研習路徑則剛好反了過來，是先學習了電動力學和量子力學，再來研究這個大部分普通物理課上一筆帶過的斯托克斯定律。張朝陽以此提示大家，學科的發展本就沒有領域的限制，只有將一個領域的內容學扎實了，才能在面對新問題時靈活運用先前的知識。

本堂課也是在向斯托克斯致敬。斯托克斯在劍橋執教期間對流體力學領域做出很多奠基性的貢獻，并早在1851年就得出了今天所講的斯托克斯定律。他的研究對后世產生了深遠的影響，比如矢量分析的工具幫助了麥克斯韋在1865年建立起電磁方程組，密立根油滴實驗用到了斯托克斯定律來確定油滴的質量，風洞試驗、汽車風阻模擬、機翼設計等工程問題也處處有NS方程的身影。NS方程的解的存在性與光滑性至今仍是數學界的研究熱點，并在21世紀初被列為千禧年七大數學難題之一。盡管現在能借助計算機數值求解NS方程，但尋找解析解依舊是幫助人們理解物理圖像的重要手段之一。