《計算物理學》課件全套 第1-6章 計算物理學簡介-有限元方法

計算物理學

—理論物理與實驗物理的橋梁計算物理是過去二十年來物理學中發(fā)展最迅速的一個領域，在物理學的各個分支學科中發(fā)揮了及其重要的作用，已成為實驗物理和理論物理同等重要并可持續(xù)發(fā)展的二級學科。計算物理是物理學研究的重要方向，已成為解決傳統(tǒng)理論物理所無法解決的問題、替代或減少實驗成本、揭示新的物理效應和規(guī)律的必要手段。計算物理的發(fā)展也對材料、信息、能源、化學、生物學科及其應用的發(fā)展起到推動作用，在國家科技與國防戰(zhàn)略發(fā)展中至關重要。計算物理以計算機為工具，以計算方法和計算軟件為手段，研究和發(fā)現(xiàn)物質(zhì)結(jié)構及其運動規(guī)律，是物理學的重要組成部分。

——《中國學科發(fā)展戰(zhàn)略：計算物理學》第一章計算物理學簡介1.1、計算物理學的誕生1.2、計算物理學的研究內(nèi)容和方法1.3、計算物理學的計算機實現(xiàn)1.4、計算物理學的發(fā)展1.5、計算物理學在相關學科中的應用1.6、本課程的授課內(nèi)容、作業(yè)和考核方式傳統(tǒng)物理學實驗科學—19世紀中葉前：原始積累時期，大部分的物理規(guī)律是基于實驗歸納得出的。實驗+理論科學—1886年Maxwell（麥克斯韋）給出了電磁場麥克斯韋方程組，預言了電磁波的存在，后被實驗證實。初步展示了物理理論思維和歸納演繹方法的強大之處。從這時起，有別于實驗物理，理論物理相對獨立，開始形成了物理學的另一個重要分支。實驗+理論科學—20世紀初，量子力學和相對論的誕生，使物理學進入了一個全新的時代，理論物理發(fā)展成為一門成熟的物理學分支學科。由此物理學正式形成了實驗、理論兩大分支，并相互促進、相互發(fā)展，引發(fā)了20世紀科學技術的重大革命。

一切靠實驗。

理論也能管點用。

理論能頂半邊天。好景不長隨著研究體系復雜性的增加，理論與實驗遇到了難以克服的困難：過于繁瑣的計算：例如，求解H2O電子波函數(shù)大量的實驗數(shù)據(jù)：例如，一秒鐘產(chǎn)生100x100個數(shù)據(jù)點，實驗持續(xù)了一個小時，共計36000000個數(shù)據(jù)點。系統(tǒng)復雜性體現(xiàn)在各個方面：理論：由單體問題轉(zhuǎn)變到多體問題、線性系統(tǒng)發(fā)展到非線性系統(tǒng)、低維體系到高維體系、標量系統(tǒng)擴展到矢量系統(tǒng)、常微分方程轉(zhuǎn)變到偏微分方程、低級微擾轉(zhuǎn)變到高級微擾、理想化模型擴展到實際復雜模型、單一物理學科發(fā)展到綜合學科的研究。實驗：研究對象和范圍的拓廣、研究精度、極限更高、設備更復雜、以及海量的的實驗數(shù)據(jù)。

花費畢生精力？

記錄這些數(shù)據(jù)需要多少紙張？每頁記錄500個數(shù)據(jù)點、1頁厚度0.1mm，需要約7米厚的紙記錄這些數(shù)據(jù)計算物理學省時、省力、省錢、靈活自由、可以模擬極端條件等。。。

數(shù)學不好不能怨天尤人

成本很高、難度很大原始推動力：據(jù)查證：“計算物理（ComputationalPhysics）”一詞首次正式出現(xiàn)是在美國1963年出版的“計算物理方法”叢書中。該套叢書從1963年至1977年相繼出版了17卷，是美國Manhattan（曼哈頓）計劃解密的研究成果。其內(nèi)容涉及統(tǒng)計物理、量子力學、流體力學、核物理、天體物理、固體物理、等離子體物理、受控熱核反應等方面的物理問題，介紹了有關“計算物理”的概貌。因此，可以稱“計算物理學”的原始推動力是（核）武器研制。直至今日，世界上速度最快的計算機仍然擔負著武器研制和國家安全的任務。研制我國第一顆原子彈所采用的手搖計算機核武器研制？第一章、計算物理學簡介1.1、計算物理學的誕生1.2、計算物理學的研究內(nèi)容和方法1.3、計算物理學的計算機實現(xiàn)1.4、計算物理學的發(fā)展1.5、計算物理學在相關學科中的應用1.6、本課程的授課內(nèi)容、作業(yè)和考核方式計算物理學–定義以計算機和計算機科學技術為工具和手段，應用適當?shù)臄?shù)學方法對物理問題進行數(shù)值分析與研究、對物理過程進行數(shù)值模擬的一門新興學科，是物理、數(shù)學、計算機應用三者結(jié)合的產(chǎn)物。是物理學的第三大分支。關鍵詞：計算機為工具，物理問題與過程、數(shù)值分析與模擬其他定義：A.J.Freeman：

復雜物理體系的數(shù)值研究秦元勛：利用現(xiàn)代大型快速計算機對物理過程進行數(shù)值模擬D.Biskamp：物理學中計算機應用的一切方面Computationalphysicscombinesphysics,computerscienceandappliedmathematicsinordertoprovidescientificsolutionstorealisticandoftencomplexproblems.Areasofapplicationincludeenvironmentalmodeling,nuclearcleanup,thedesignofmaterials,groundwatertransport,thenatureofelementaryparticles,medicalimaging,andenergymanagement.Acomputationalphysicistunderstandsnotonlytheworkingsofcomputersandtherelevantscienceandmathematics,butalsohowcomputeralgorithmsandsimulationsconnectthetwo.計算物理學—一些英文書籍中的定義或描述計算物理學—用計算機做實驗互聯(lián)網(wǎng)計算物理學的研究范圍計算理論物理：1、復雜物理體系的數(shù)值計算與模擬

2、復雜物理體系的解析計算與分析（符號運算）計算實驗物理3、物理實驗數(shù)據(jù)的采集、分析、與處理4、物理實驗過程與實驗系統(tǒng)的模擬與控制5、物理圖像的獲得、識別與處理。計算物理學的研究方法--三板斧三板斧：建模(modeling)、模擬(simulation)、計算(computation)Modeling：偏重于物理、數(shù)學模型的建立，這是計算物理的基礎。Simulation：是指對物理過程的描述、對物理現(xiàn)象的表達和對物理規(guī)律的探索。Computation：應用計算機的數(shù)值研究，包括數(shù)值計算和數(shù)值分析。與Simulation一通為計算物理的基本核心內(nèi)容Modeling主要是人腦完成的Simulation&Computation是電腦完成的

極度重要+非常辛苦

強悍的計算機可以彌補人腦的不足

辛苦程度不詳計算物理學的研究方法–四步走第一步：確定物理模型

對物理問題進行分析，抓住主要因素、忽略次要因素，建立相應的物理模擬第二步：選取數(shù)值方法：

在第一步的基礎上，選擇合適的數(shù)值方法在計算機上實現(xiàn)。

算法的重要性：計算能力不足以滿足直接求解的需求，因此需要一個合適的算法，將復雜問題化簡為簡單問題、簡單問題化簡為基本問題，基本問題化簡為計算機可以執(zhí)行的運算。第三步：分析計算結(jié)果

對計算得到的結(jié)果進行分析，得到有價值的物理信息。同時要考慮到計算產(chǎn)生的誤差、收斂性、穩(wěn)定性等等因素第四步：得出物理結(jié)論

整理大量數(shù)據(jù)的基礎上、得出相關物理結(jié)論、反應物理規(guī)律、解釋物理機制、進而給出可能的物理預言。

可以推測算法的選取對人的生理和心理均有很大影響

正確的分析可以避免因為計算誤差而給出錯誤的物理信息、不迷信數(shù)值結(jié)果計算物理學和理論物理、實驗物理的關系實驗物理理論物理計算物理提供原理、解釋結(jié)果

提供模擬結(jié)果、數(shù)值計算、

提供計算結(jié)果、啟發(fā)新理論

模擬實驗、提供理論數(shù)據(jù)、數(shù)據(jù)分析、采集、實驗自動化產(chǎn)生數(shù)據(jù)、檢驗結(jié)果

提出實驗、解釋實驗

提供數(shù)據(jù)、驗證理論第一章、計算物理學簡介1.1、計算物理學的誕生1.2、計算物理學的研究內(nèi)容和方法1.3、計算物理學的計算機實現(xiàn)1.4、計算物理學的發(fā)展1.5、計算物理學在相關學科中的應用1.6、本課程的授課內(nèi)容、作業(yè)和考核方式計算機實現(xiàn)

中央處理器(CPU)

內(nèi)存HD系統(tǒng)總線(主板)網(wǎng)絡(可選)輸入輸出設備可視化設備操作系統(tǒng)編譯器計算程序數(shù)學庫并行環(huán)境可視化軟件硬件軟件中央處理器(CPU)

–完成計算內(nèi)存–存儲數(shù)據(jù)供CPU調(diào)用HD–長期存儲數(shù)據(jù)系統(tǒng)總線(主板)–連接各個單元網(wǎng)絡(可選)–連接各個計算節(jié)點輸入輸出設備–與人交互可視化設備–分析結(jié)果操作系統(tǒng)--讓硬件聽指揮編譯器–告訴計算機如何工作計算程序–指揮計算機工作數(shù)學庫–避免重復勞動并行環(huán)境–人多力量大可視化軟件計算能力的發(fā)展單個計算機的計算能力，包括大型機、小型機、工作站、個人電腦等。199619972000200320062009CPU(flop/s)10G100G1T10T100T1PRAM10GB100GB1TB10TB100TB1PB（超級）并行計算機19441946~1960~1970~1990~2000~2010CPU(flop/s)3flop/s5000~1M~10M~100M~2G~40GRAM32MB1GB32GB1M=106;1G=109;1T=1012;1P=1015并行計算機并行計算機是由多個處理器組成(在這些處理器之間可以相互通訊和協(xié)調(diào))將一個應用任務分解成多個子任務，分配給不同的處理器，各個處理器之間相互協(xié)同的，同時執(zhí)行子任務的過程并能夠高速、高效地進行復雜計算的計算機系統(tǒng)。并行計算機得名是相對于串行計算機而言的。串行計算機只有單個處理器，順序執(zhí)行計算機程序并行計算的概念是在20世紀70年代中后期提出的，已有30多年的發(fā)展歷史并行計算的必要性1、可以大大加快運算速度，即在更短的時間內(nèi)完成相同的計算量，或者解決原來根本不能計算的非常復雜的問題。2、提高傳統(tǒng)計算機的計算速度受到諸多限制：物理極限、量子效應、加工工藝、散熱、成本，等等3、并行計算較之SMP機器投入較低、靈活性強。SMP=SymmetricMulti-ProcessingSharedMemory：多個處理器公用內(nèi)存

并行計算機的發(fā)展大型機：SMP(如IBMrs6k)PC

Cluster：多節(jié)點+SMP+多核心刀片服務器多節(jié)點+多處理器+多核心+GPU硬件廠商整機：曙光：魔方@SSC(曙光5000A)聯(lián)想：深騰7000@SCCASIBM：Roadrunner(BladeCenter)Cray：Jaguar@ORNL(XT5)SGI：Pleiades@NASA(Altix)HP：TSUBAME@JapanSun：RedSky@NREL(SunBlade)DELL…CPU+GPU：Intel：Xeon/ItaniumAMD：OpteronIBM：POWERNVidia：Tesla…網(wǎng)絡連接：QsNetMyrinetInfiniband…并行機的軟件操作系統(tǒng)：LinuxUnixWindows…數(shù)學庫：MKL(Intel)：MathKernelLibraryACML(AMD)：AMDCoreMathLibraryBLAS：BasicLinearAlgebraSubprogramsLAPACK：LinearAlgebraPACKage…編譯器：Ifort(Fortran)GCC(C/CPP)PGI(Fortran/C/CPP)VTune(調(diào)試器)…并行環(huán)境：MPI(消息傳遞)OpenMP(共享存儲)…編程語言：FortranC/CPPMatlabPython…科學計算軟件：PWSCFVASPANSYSLAMMPS…MPI什么是MPI：MessagePassingInterface：消息傳遞函庫函數(shù)的標準規(guī)范、支持Fortran、C、Matlab等計算機語言。不是計算機語言，具有上百個函數(shù)庫可以被直接調(diào)用，實現(xiàn)消息傳遞功能是一種技術標準、規(guī)范的代表，不是特定指某個具體實現(xiàn)MPI的優(yōu)勢：高可移植性：可以在所有主流計算機上，例如Unix、Liunx、MacOS、Windows等，實現(xiàn)其并行功能，程序無需修改。第一章、計算物理學簡介1.1、計算物理學的誕生1.2、計算物理學的研究內(nèi)容和方法1.3、計算物理學的計算機實現(xiàn)1.4、計算物理學的發(fā)展1.5、計算物理學在相關學科中的應用1.6、本課程的授課內(nèi)容、作業(yè)和考核方式計算方法的發(fā)展19世紀：常微分和偏微分方程

50年代：MonteCarlomethod(MC，蒙特卡羅)

60－70年代：DensityFunctionalTheory(DFT，密度泛函理論)、Moleculardynamics(MD，分子動力學)、Wilsonnumericalrenormalizationgroup(NRG，數(shù)值化重整化群)…

80－90年代：

QuantumMonteCarlomethods(QMC，量子蒙特卡羅)、Car-Parrinellomethod(CP)、Densitymatrixrenormalizationgroupmethod(DMRG，密度矩陣重整化群)…

21世紀：non-equilibriumGreen‘sfunction(非平衡格林函數(shù))，…O/PDE

常/偏微分方程–19世紀微分方程論是數(shù)學的重要分支之一。大致和微積分同時產(chǎn)生，并隨實際需要而發(fā)展。含自變量、未知函數(shù)和它的微商（或偏微商）的方程稱為常（或偏）微分方程。微分方程是研究自然科學和社會科學中的事件、物體和現(xiàn)象運動、演化和變化規(guī)律的最為基本的數(shù)學理論和方法。物理、化學、生物、工程、航空航天、醫(yī)學、經(jīng)濟和金融領域中的許多原理和規(guī)律都可以描述成適當?shù)奈⒎址匠?。然而，我們知道，只有少?shù)十分簡單的微分方程能夠用初等方法求得它們的解，多數(shù)情形只能利用近似方法求解，如漸進法。還有一類近似方法稱為數(shù)值方法，它可以給出解在一些離散點上的近似值。微分方程數(shù)值解就是在計算機上使用的，研究并解決數(shù)學問題的數(shù)值近似解的方法。它既有理論上的抽象性和嚴謹性，又有實用性和實驗性的技術特征。因此，微分方程數(shù)值解已用到科學技術和社會生活的各個領域中。MonteCarlo(蒙特卡羅)

–上世紀4、50年代計算機隨機模擬方法，是一種基于“隨機數(shù)”的計算方法。這一方法源于美國在第二次世界大戰(zhàn)進研制原子彈的“曼哈頓計劃”。馮·諾伊曼用摩納哥賭城—MonteCarlo—來命名這種方法。MonteCarlo方法的基本思想很早以前就被人們所發(fā)現(xiàn)和利用。早在17世紀，人們就知道用事件發(fā)生的“頻率”來決定事件的“概率”。19世紀人們用投針試驗的方法來決定圓周率π。本世紀40年代電子計算機的出現(xiàn)，特別是近年來高速電子計算機的出現(xiàn)，使得用數(shù)學方法在計算機上大量、快速地模擬這樣的試驗成為可能?？萍加嬎阒械膯栴}比這要復雜得多。比如金融衍生產(chǎn)品（期權等）的定價及交易風險估算，問題的維數(shù)（即變量的個數(shù)）可能高達數(shù)百甚至數(shù)千。對這類問題，難度隨維數(shù)的增加呈指數(shù)增長，這就是所謂的“維數(shù)的災難”(CourseDimensionality)，傳統(tǒng)的數(shù)值方法難以對付（即使使用速度最快的計算機）。MonteCarlo方法能很好地用來對付維數(shù)的災難，因為該方法的計算復雜性不再依賴于維數(shù)。以前那些本來是無法計算的問題現(xiàn)在也能夠計算量。為提高方法的效率，科學家們提出了許多所謂的“方差縮減”技巧。Moleculardynamics(分子動力學)–上世紀6、70年代分子動力學是一門結(jié)合物理，數(shù)學和化學的綜合技術。分子動力學是一套分子模擬方法，該方法主要是依靠牛頓力學來模擬分子體系的運動，以在由分子體系的不同狀態(tài)構成的系綜中抽取樣本，從而計算體系的構型積分，并以構型積分的結(jié)果為基礎進一步計算體系的熱力學量和其他宏觀性質(zhì)。分子動力學可以用于NPT，NVE，NVT等系綜的計算，是一種基于牛頓力學確定論的熱力學計算方法，與蒙特卡洛法相比在宏觀性質(zhì)計算上具有更高的準確度和有效性，可以廣泛應用于物理，化學，生物，材料，醫(yī)學等各個領域。在實際應用中，經(jīng)常把分子動力學方法和蒙特卡羅法聯(lián)合使用。MD與MC類似，都是“由小見大”的方法：刻畫微觀過程，統(tǒng)計的描述體系宏觀性質(zhì)。MD由牛頓方程支配，描述微觀過程更準確。第一性原理密度泛函(DFT)方法–上世紀6、70年代根據(jù)原子核和電子互相作用的原理及其基本運動規(guī)律，運用量子力學原理，從具體要求出發(fā)，經(jīng)過一些近似處理后直接求解薛定諤方程的算法，習慣上稱為第一原理。第一性原理通常是跟計算聯(lián)系在一起的，是指在進行計算的時候除了告訴程序你所使用的原子和他們的位置外，沒有其他的實驗的，經(jīng)驗的或者半經(jīng)驗的參量，且具有很好的移植性。密度泛函理論是一種研究多電子體系電子結(jié)構的量子力學方法。核心思想是：體系的基態(tài)能量是體系電荷密度的泛函—這將最難處理的多體問題（由于處在一個外部靜電勢中的電子相互作用而產(chǎn)生的）簡化成了一個沒有相互作用的電子在有效勢場中運動的問題。密度泛函理論在物理和化學上都有廣泛的應用，特別是用來研究分子和凝聚態(tài)的性質(zhì)，是凝聚態(tài)物理和計算化學領域最常用的方法之一。在多數(shù)情況下，與其他解決量子力學多體問題的方法相比，采用局域密度近似的密度泛函理論給出了非常令人滿意的結(jié)果，同時計算相比實驗的費用要少。???計算量如何???難以承受直接求解Nanoscale計算物理的研究尺度：時間和空間尺度DFTMDMCO/PDE、FEQCDQMCMacroscaleMesoscaleAtomicscale具體的界限并不嚴格，且可能隨著計算機技術的發(fā)展而改變！計算物理的研究尺度：時間和空間尺度具體的界限并不嚴格，且可能隨著計算機技術的發(fā)展而改變！第一性原理方法(從頭算方法,

Ab-initio)從第一性原理出發(fā)，通過解Schr?dinger(orDirac)量子力學方程，求解材料的各種性質(zhì)。優(yōu)點：1、可以給出原子幾何（力學）和電子結(jié)構信息。2、可以處理化學鍵的斷裂、形成或者電子結(jié)構的重新分布3、知道計算誤差的來源，可以較好的評估計算精度4、原則上可以只輸入原子種類和坐標就基本可以得到準確的物質(zhì)性質(zhì)不足：1、只能處理小的體系，O(102)個原子2、只能處理較快的過程，O(10)ps.3、計算機比較辛苦半經(jīng)驗方法,

Semi-empiricalmethods包含一些從實驗擬合或者嚴格重頭算得到的參數(shù)，通過求解簡化的Schr?dinger(orDirac)量子力學方程，求解材料的各種性質(zhì)。優(yōu)點：1、同樣可以處理化學鍵的斷裂、形成或者電子結(jié)構的重新分布2、可以求解更大的體系，O(103)個原子3、可以求解更長時間尺度的物理問題，O(10)ns不足：1、計算誤差的來源復雜，無法較好的評估計算精度2、需要從實驗擬合或者嚴格的第一性原理計算得到經(jīng)驗參數(shù)。原子尺度模擬方法,

Atomisticsimulationmethods使用第一性原理方法或者半經(jīng)驗方法得到原子之間的相互作用勢，結(jié)合半經(jīng)典的統(tǒng)計力學，采用MC或者MD方法求解體系的熱力學或輸運性質(zhì)等。適用從納米到介觀尺度。優(yōu)點：1、可以求解微觀體系的性質(zhì)甚至更復雜一些的體系，O(104-6)atoms.2、可以求解更、更長時間尺度的物理問題，O(1)μs不足：1、計算精度取決于原子之間相互作用勢的好壞2、無電子相關的任何信息3、無法描述一些物理過程，例如時間尺度更短的：固體中的原子擴散、化學反應；時間尺度更長的：蛋白質(zhì)折疊等。介觀尺度方法,Mesoscalemethods在原子尺度方法的基礎上進行簡化，例如，不考慮時間尺度更短的自由度，將一些原子看成一個原子集合(‘blobsofmatter’)，考慮集合-集合之間的相互作用—有效勢。適用范圍從介觀到宏觀尺度。優(yōu)點：1、可以求解介觀尺度復雜體系的結(jié)構特性，O(108-9)atoms2、可以描述前面幾種方法無法達到的長時間尺度，O(1)s不足：1、通常只能定性的給出趨勢，較難給出定量的結(jié)果2、因為引入的太多的近似，因此對結(jié)果的解釋并不是非常確定宏觀尺度(Macroscale)方法,Continuummethods假設物質(zhì)是連續(xù)的（不考慮原子細節(jié)），且體系物理性質(zhì)可以用場量來描述。數(shù)值求解與該現(xiàn)象（彈性力學、流體力學等等）相關的平衡方程（常、偏微分方程），預測體系的物理性質(zhì)。優(yōu)點：1、可以原則上解決任何宏觀幾何、時間適度的問題不足：1、需要輸入大量參數(shù)(彈性張量、擴散系數(shù)、物態(tài)方程等等)，這些參數(shù)需要從實驗或者更小尺度的模擬中得到。2、無法從電子、原子、分子尺度的細節(jié)解釋計算結(jié)果。不同尺度間的聯(lián)系,

Connectionbetweenthescales從低到高(Up-scaling)：從更低尺度的計算結(jié)果中得到一些經(jīng)驗參數(shù)用來進行較大尺度的模擬。相對簡單，是一個演繹問題，有確定解。從高到低(down-scaling)：用較大尺度的信息（通常是實驗結(jié)果）去擬合較小尺度的參數(shù)。相對困難，因為解并不確定，是一個歸納問題。例如介觀模擬中并沒有原子細節(jié)（僅有原子集合），如何利用這些信息得到原子細節(jié)是較為困難的。第一章、計算物理學簡介1.1、計算物理學的誕生1.2、計算物理學的研究內(nèi)容和方法1.3、計算物理學的計算機實現(xiàn)1.4、計算物理學的發(fā)展1.5、計算物理學在相關學科中的應用1.6、本課程的授課內(nèi)容、作業(yè)和考核方式計算物理的應用領域物理學領域凝聚態(tài)物理：模擬材料的力學、熱學、電學、磁學、光學性質(zhì)高能物理：模擬核反應、粒子撞擊大氣物理：模擬天氣預報、臺風地球物理：模擬地球內(nèi)部運動、火山、板塊活動天體物理：模擬星體的運動、誕生、湮滅理論物理：符號運算，使用計算機“推公式”實驗物理：數(shù)據(jù)采集、處理、擬合材料學領域：多尺度模擬：量子力學

經(jīng)典力學化學領域：化學反應動力學生命科學領域：生物大分子折疊、藥物篩選、生物信息學環(huán)境科學領域：環(huán)境風險評估、污染檢測工程領域：船舶、飛行器、汽車、建筑設計與評估社會科學領域：金融物理學、最優(yōu)化策略凝聚態(tài)物理–DFT/MD方法晶格常數(shù)的預測計算物理應用舉例(1)結(jié)構、電子（導電性、光學性質(zhì)）、聲子（熱學性質(zhì)）、磁性。。。過渡金屬氧化物中O2s態(tài)的波函數(shù)實空間分布金剛石的能帶結(jié)構Al-Mn-N化合物的自旋極化電子態(tài)密度

半金屬計算物理應用舉例(2)計算物理應用舉例(3)計算物理應用舉例(4)磁性計算物理應用舉例(5)Diffusion-LimitedAggregationRandomaggregateof3600particlesonasquarelattice.PhysicalReviewLetter47,1400(1981).Citedover2000times.計算物理應用舉例(6)計算物理應用舉例(7)材料模擬–分子動力學方法金屬玻璃在形變過程中非仿射形變的分布計算物理應用舉例(8)生命科學–DFT/MD/MC

方法致病機理研究帕金森氏病癥會殺死控制協(xié)調(diào)移動的大腦細胞，患者將遭受肌肉震顫和僵硬?？茖W家最初認為一種叫做阿爾法核素的蛋白質(zhì)分子能夠形成束狀，具有破壞性的纖維物質(zhì)，然而其具體的形成過程尚不清楚。模擬顯示阿爾法核素分子可以在大腦細胞膜上穿孔(圖像中綠色區(qū)域)聚集形成環(huán)狀結(jié)構，從而導致大腦細胞死亡，引發(fā)帕金森氏癥。這些模擬研究有助于研發(fā)新一代預防疾病藥物。藥物篩選：新藥研發(fā)過程費用昂貴、時間冗長、淘汰率高。大約有90%的候選藥物在臨床期間被淘汰，這是新藥研發(fā)過程費用昂貴的主要因素。所謂的藥物虛擬篩選(virtualscreening)，是指對化合物在其合成之前通過計算機模擬預測其藥動學相關的特性而進行篩選。右圖是對H1N1病毒受體的藥物篩選。H蛋白質(zhì)負責寄主識別，N蛋白質(zhì)負責病毒傳播，因此抑制N蛋白就可以有效防止細胞感染。采用特定化學分子使N蛋白(神經(jīng)氨酸)失活計算物理應用舉例(9)生命科學--蛋白質(zhì)折疊

largetimescale(~s)moleculardynamicsmodelingDNA解旋mRNA轉(zhuǎn)錄

核糖體翻譯

多肽鏈

折疊

具有功能的蛋白質(zhì)計算物理應用舉例(10)高能物理這張圖片模擬的是ATLAS里產(chǎn)生的黑洞。這個軌跡是模仿大型強子對撞機(LHC)上的ATLAS探測器得出的模擬數(shù)據(jù)。如果質(zhì)子-質(zhì)子撞擊期間產(chǎn)生了微型黑洞，這些軌跡就會形成。這種小型黑洞會通過霍金輻射(Hawkingradiation)方式，立刻消失不見。大氣物理全球濕度與降水的演化，白色的云代表濕度，色點代表降水。溫度升高對格陵蘭島和南極洲上的巨大冰蓋將產(chǎn)生什么樣的影響？

未來幾十年中，海洋和大氣運動方式對區(qū)域氣候?qū)a(chǎn)生什么樣的影響？氣候變化對包括颶風在內(nèi)的熱帶氣旋的強度和頻率將產(chǎn)生何種影響？計算物理應用舉例(11)地球物理GeoPhysics：多種地球物理現(xiàn)象的耦合——火山內(nèi)部電場、磁場、流體流動耦合分析。天體物理位于地球上低緯度和中緯度地帶的人們可看到，春季太陽西沉，黃昏過后之時，西方地平線有著微末的三角形光錐；而秋季太陽東升，晨曦未現(xiàn)之時，東方自地平線向上伸展出些許“火舌”，這就是黃道光?？茖W家們認為黃道光的起因主要是行星際塵埃對太陽光的散射，類似一個龐大的以太陽為中心的“塵埃餅”，但至今仍對這些塵埃的來源摸不著頭腦。美國行星動力學家用計算機模型破解了神秘的夜空輝光——黃道光形成之謎。結(jié)論認為，黃道光塵埃幾乎全部來自短周期彗星，推翻了長久以來的推測。計算物理應用舉例(12)實驗物理控制系統(tǒng)運行采集實驗數(shù)據(jù)監(jiān)視儀器狀態(tài)數(shù)據(jù)在線分析…理論物理符號運算周期性？使用快速傅立葉變換(FFT)給出周期性計算物理應用舉例(13)環(huán)境評估辦公室內(nèi)的風環(huán)境：風水師必會解復雜體系的流體力學方程據(jù)物理學家組織網(wǎng)報道，研究人員根據(jù)美國宇航局衛(wèi)星監(jiān)測數(shù)據(jù)、地面空氣質(zhì)量平行監(jiān)測數(shù)據(jù)和最近公布的農(nóng)作物產(chǎn)量統(tǒng)計資料，通過一個模擬全球臭氧污染形成和轉(zhuǎn)運的計算機模型顯示，地表臭氧濃度不斷提高導致了農(nóng)作物產(chǎn)量下降，全球農(nóng)作物損失每年約達260億美元甚至更多。計算物理應用舉例(14)工程領域J20紅外信號模擬–紅外隱身解熱傳導方程賽車行駛時氣流旋渦的概況：氣流流過車輪處開始紊亂并形成旋渦--如圖中藍色湍流所示；氣流流過車頂時還保持良好秩序，經(jīng)過尾翼下方時迅速紊亂并形成巨大的旋渦--這也是賽車前進阻力的最大源泉。

解空氣動力學方程船舶水動力特性計算模擬，提高適航性，航速等。解流體力學方程計算物理應用舉例(15)金融物理金融物理學的英文為Econophysics,是由波士頓大學的物理學教授H.E.Stanley在1995年首先提出的,從而解決了“為什么物理專業(yè)的學生可以從事金融學研究并取得物理學位”這一實際問題。金融物理學是用統(tǒng)計物理、理論物理、復雜系統(tǒng)理論、非線性科學、應用數(shù)學等的概念、方法和理論研究金融市場通過自組織而涌現(xiàn)的宏觀規(guī)律及其復雜性的一門新興交叉學科。簡言之,金融物理學家將金融市場看作一個復雜系統(tǒng),把其中的各種數(shù)據(jù)如個股價格、指數(shù)、房價等看作是物理實驗數(shù)據(jù),力圖尋找和闡釋其中的“物理”規(guī)律。計算物理應用舉例(16)第一章、計算物理學簡介1.1、計算物理學的誕生1.2、計算物理學的研究內(nèi)容和方法1.3、計算物理學的計算機實現(xiàn)1.4、計算物理學的發(fā)展1.5、計算物理學在相關學科中的應用1.6、本課程的授課內(nèi)容、作業(yè)和考核方式授課內(nèi)容第二章、常微分方程數(shù)值解OrdinaryDifferentialEquation(ODE)第三章、偏微分方程的數(shù)值解PartialDifferentialEquation(PDE)第四章、分子動力學方法MolecularDynamicsMethod(MD)第五章、蒙特卡洛方法MonteCarloMethod(MC)第六章、有限元方法FiniteElementsMethod(FE)第七章、密度泛函理論Density

Functional

Theory

(DFT)（根據(jù)課程進展情況再定）參考書目《計算物理學》李茂枝，季威，郭茵，盧仲毅編著，中國人民大學出版社，2014《計算方法引論》

徐萃薇，孫繩武編著，高等教育出版社，2007《計算物理學》

馬文淦編著，科學出版社，2005《計算物理學》

劉金遠，段萍，鄂鵬，科學出版社，2012《An

Introduction

to

Computational

Physics》Tao

Pang,

Cambridge

University

press,

2006作業(yè)和考核平時：60%出勤、課堂情況、作業(yè)(應用能力)期末：40%閉卷考試(概念+建模能力)計算物理

計算物理是指使用現(xiàn)代計算技術（計算機、軟件和硬件）來——探索、研究和驗證新的物理現(xiàn)象或物理特性；一方面，它作為理論的一部分被用來驗證和解釋實驗發(fā)現(xiàn)；另一方面，它本身就是一種實驗，被用來檢驗理論模型的正確性；在許多情況下，它被用來取代實驗，降低科研成本。對于計算物理而言，計算方法或算法的發(fā)展和應用是其兩個重要的內(nèi)涵。計算物理涉及到物理學的各個分支（高能與粒子物理、原子與分子物理、統(tǒng)計物理和凝聚態(tài)物理等），已經(jīng)成為一個研究領域，而且計算物理方法在工程、流體力學、材料科學、金融與經(jīng)濟和生命科學等交叉學科研究領域的有著廣泛應用潛力。

缺點：不能獲得物理定律和理論公式，且計算結(jié)果缺乏嚴格的論證，其結(jié)果仍需實驗驗證。小結(jié)和展望計算物理學是一門蓬勃發(fā)展的“新興”學科歷史悠久、發(fā)展迅猛、前景廣闊廣闊的應用領域：材料、化學、生命科學、清潔能源、生物醫(yī)藥等等任何一個領域均與人們的生活息息相關還面臨著很多挑戰(zhàn)算準、算大、算快希望更多的有志青年投身到這一領域之中尤其是在座的各位作業(yè)1、計算物理的含義是什么？2、計算物理的研究范疇是什么？3、計算物理的研究方法是什么？(1、2、3三選二)4、闡述計算機數(shù)值模擬方法與理論、實驗方法相比有什么特殊的優(yōu)點和局限性。5、闡述計算物理學和實驗物理及理論物理的關系(4、5二選一)6、闡述在計算機上進行計算物理研究都需要什么樣的硬件和軟件7、闡述并行計算的必要性、優(yōu)點和所需軟、硬件(6、7二選一)8、舉例說明除凝聚態(tài)物理、材料科學以外計算物理的應用領域（至少三個，(必做))9、簡述多尺度模擬中幾個尺度的特點及其所使用方法(必做)第二章常微分方程的數(shù)值求解方法1、引言2、數(shù)值求解的基本思想3、歐拉方法4、R-K方法5、線性多步法6、預估校正方法7、步長選擇8、常微分方程組和高階微分方程9、邊值問題10、有限差分方法內(nèi)容提綱1、引言微分方程：定義：包含自變量、未知函數(shù)及未知函數(shù)的導數(shù)或微分的方程稱為微分方程?！俺！焙汀捌保涸谖⒎址匠讨?自變量的個數(shù)只有一個,

稱為常微分方程。自變量的個數(shù)為兩個或兩個以上的微分方程叫偏微分方程。微分方程中出現(xiàn)的未知函數(shù)最高階導數(shù)的階數(shù)稱為微分方程的階數(shù)。如果未知函數(shù)y及其各階導數(shù)。都是一次的，則稱它是線性的，否則稱為非線性的。解析求解?在高等數(shù)學中，對于常微分方程的求解，給出了一些典型方程求解析解的基本方法，如可分離變量法、常系數(shù)齊次線性方程的解法、常系數(shù)非齊次線性方程的解法等。

但能求解的常微分方程仍然是有限的（1%?)。

大多數(shù)的常微分方程是不可能給出解析解。譬如：這個一階微分方程就不能用初等函數(shù)及其積分來表達它的解。一個具體的例子

解析求解無能為力例：多體運動問題多體運動問題是物理研究中的一個普遍問題。設有n個物體（粒子）在靜電場或引力場中運動，因而有6n個未知數(shù)xi(坐標)，vi（速度），（i=1,2,…,n）。根據(jù)牛頓定律或庫倫定律，應該滿足6n個一階常微分方程組：其中mi是第i個物體（或粒子）的質(zhì)量，ei在靜電場中是第i個粒子的電荷，在引力場中事實上，當n≥3時，無法找到解析解。數(shù)值求解根據(jù)常微分方程的定解條件，常微分方程可分為初值問題和邊值問題。如果定解條件是描述函數(shù)在初始點的狀態(tài)，則稱為初值問題。一個典型的常微分方程初值問題的形式如下：

具有以下形式的常微分方程及定解條件，如

則為常微分方程的邊值問題。常微分方程初值問題的數(shù)值求解

由此可見，常微分方程數(shù)值求解的基本出發(fā)點就是將求解區(qū)間離散化。

3、歐拉（Euler）方法這是一種最為基礎、簡單的方法，它的精確度不高，因此在實際計算中并不常使用。由于其能反映出很多復雜方法的特征，被廣泛用于微分方程解法的講解。根據(jù)計算方法不同，分以下三種：歐拉（Euler）方法：根據(jù)前一個點的函數(shù)值，可以計算后一個點的函數(shù)值。

三種簡單的數(shù)值解法1、化導數(shù)為差商；2、數(shù)值積分法；3、Taylor展開

數(shù)值積分法

Talyor級數(shù)展開

歐拉公式及其幾何意義

歐拉方法的穩(wěn)定性分析由歐拉公式可以看出，在常微分方程初值問題的數(shù)值求解過程中，每一步計算都包含了前一步的計算結(jié)果，因此前面的計算誤差將會影響后面的數(shù)值計算結(jié)果。需要對數(shù)值求解方法的誤差進行系統(tǒng)分析，考察計算過程的誤差積累會不會掩蓋真解，這就是常微分方程數(shù)值解的穩(wěn)定性問題。需要對所采用的數(shù)值方法進行穩(wěn)定性分析，確定其穩(wěn)定性區(qū)域，保證該方法得到的數(shù)值解是合理的。

總體誤差和局部誤差的關系

由此可得：

收斂性和穩(wěn)定性是兩個不同的概念。收斂性是指數(shù)值方法的截斷誤差對計算結(jié)果的影響穩(wěn)定性指的是某一步的計算誤差對計算結(jié)果的影響此外，穩(wěn)定性與步長密切相關。

數(shù)值方法的穩(wěn)定性

0-1-2ReImg歐拉方法的絕對穩(wěn)定區(qū)域絕對穩(wěn)定區(qū)域越大，這個方法的適應性就越強。若絕對穩(wěn)定區(qū)域包含復平面的整個左半平面，就稱這個數(shù)值方法是A穩(wěn)定的。210ReImg后退歐拉方法的絕對穩(wěn)定區(qū)域：向后歐拉方法是A穩(wěn)定的，但是仍受到迭代的限制，即0<hL<1。

總體截斷誤差局部截斷誤差穩(wěn)定性歐拉方法向后歐拉方法A穩(wěn)定改進歐拉方法A穩(wěn)定

以上三種方法的誤差和穩(wěn)定性總結(jié)如下：RUNGE-KUTTA方法

該方程組存在無窮多個解。所有滿足該方程組的解都為二階Runge-Kutta公式。

這就是改進歐拉公式。因此，改進歐拉公式也是二階Runge-Kutta公式。

具體推導參考《計算物理學》李茂枝等編著（2024）

四階Runge-Kutta公式的經(jīng)典形式

一步法在計算時只用到了前面一步的近似值，但要提高精度，需要增加中間函數(shù)值的計算，這就加大了計算量。R-K方法的優(yōu)缺點是否可以尋找一種精度較高、計算量適中的方法呢？注意到：計算yn+1時，yn，yn-1，yn-2等等均已知，那么可否利用這些值提高精度而不顯著加大計算量呢？優(yōu)點：一步法，在給定初值后可以逐步計算下去；精度高；在計算過程中便于改變步長缺點：計算量大線性多步法

四階Adams方法

Adams內(nèi)插公式

107內(nèi)插值迭代求解108Adams預估--矯正公式

步長的自動選擇109小步長高精度110111一個例子112步數(shù)和步長步數(shù)步長10.0735892100.0397480200.00475879300.00169146400.0123375500.356060600.0280562700.00326283790.000106687800.000120953900.005208841000.1201341060.0694882113階數(shù)和步長

二分法選取步長115常微分方程組和高階微分方程116

然后采用合適的數(shù)值方法求解

邊值問題117在實際問題中也常常碰到所謂的兩點邊值問題

要得到方程（1）的解的存在唯一性非常困難。打靶法（shootingmethod）是數(shù)值求解常微分方程兩點邊值問題的常用方法。

打靶法（shootingmethod）118

打靶法的基本思想是將微分方程的邊值問題轉(zhuǎn)化為一系列初值問題。先考慮一個初值問題，即

打靶法就是將微分方程的邊值問題轉(zhuǎn)化為初值問題。通過適當選取和調(diào)整初值條件，求解一系列初值問題，使之逼近給定的邊界條件。如果將描述的曲線視做彈道，那么求解過程就是不斷調(diào)整試射條件使之打到預定的靶子。打靶法的主要思想炮擊中的例子：試射兩發(fā)，調(diào)整后齊射αβ有限差分方法有限差分法是一種求解微分方程數(shù)值解的近似方法，其主要思想是將微分方程中的微分項直接進行差分近似，從而將微分方程轉(zhuǎn)化為代數(shù)方程組求解。隨著計算機技術和數(shù)值算法的快速發(fā)展，有限差分法已成為數(shù)值求解微分方程的重要方法之一，在眾多領域得到了廣泛的應用和發(fā)展。在有限差分方法中，我們放棄了微分方程中獨立變量可以取連續(xù)值的特征，而關注獨立變量離散取值后對應的函數(shù)值。原則上，這種方法仍然可以達到任意滿意的計算精度，因為方程的連續(xù)數(shù)值解可以通過減小獨立變量離散取值的間隔，或者通過離散點上的函數(shù)值插值來近似得到。126有限差分法包含兩部分：1、用差分代替微分方程中的微分，將連續(xù)變化的變量離散化，從而得到差分方程組的數(shù)學形式。2、求解差分方程組。在第一步中，通過所謂的網(wǎng)絡分割法，將函數(shù)定義域分成大量相鄰而不重合的子區(qū)域，通常采用的是規(guī)則的分割方式。網(wǎng)絡線劃分的交點稱為節(jié)點。

具體方法127

將函數(shù)定義域離散化后，需要求得特定問題在所有這些節(jié)點上的近似值，因此數(shù)值求解的關鍵就是要找到適當?shù)臄?shù)值計算方法。設一個函數(shù)在x點上的一階和二階微商，可以近似地用它的近鄰兩點的函數(shù)值的差分來表示，即差分形式128一、二階微商129

利用上面幾個公式，可以構造出微分方程的差分形式。

小結(jié)1、線性多步法(避免大量迭代）

2、預估校正方法（合理選擇步長）

3、步長事后估計4、兩點邊值問題

（打靶法，求解大量實際問題）5、常微分方程組和高階微分方程

高階方程總可以化作一階方程組

一階方程組的求解與一階方程的求解類似6、有限差分方法—更實用的方法131132求解一個實際問題133R-K公式的經(jīng)典形式

134Matlab形式K1=h*feval(@fx,x(i),y(i));K2=h*feval(@fx,x(i)+h/2,y(i)+K1/2);K3=h*feval(@fx,x(i)+h/2,y(i)+K2/2);K4=h*feval(@fx,x(i)+h,y(i)+K3);y(i+1)=y(i)+(K1+2*K2+2*K3+K4)/6;135初值和相應參數(shù)區(qū)間[a,b]

[0,1]

a=0;b=1;步長h=0.2(可調(diào))總共的計算次數(shù)N

N=(b-a)/h初始值：

x0=0;y0=1;136N次求解fori=1:N;x(i)=a+(i-1)*h;K1=h*feval(@fx,x(i),y(i));K2=h*feval(@fx,x(i)+h/2,y(i)+K1/2);K3=h*feval(@fx,x(i)+h/2,y(i)+K2/2);K4=h*feval(@fx,x(i)+h,y(i)+K3);y(i+1)=y(i)+(K1+2*K2+2*K3+K4)/6;end137基本大功告成！???Errorusing==>plotVectorsmustbethesamelengths.Errorin==>RKat23plot(x,y);x(N+1)=b;精確解：y=sqrt(1+2.*x);138結(jié)果令人滿意139增加求解區(qū)間[0,1][0,8]140改變h到0.02141Adams方法的修正方案y(i+1)=c(i+1)-19.*(c(i+1)-p(i+1))./270;p(i+1)=y(i)+h.*(55.*fn-59.*fn_1+37.*fn_2-9.*fn_3)./24;m(i+1)=p(i+1)+251.*(c(i)-p(i))./270;q(i+1)=feval(@fx,x(i+1),m(i+1));c(i+1)=y(i)+h.*(9.*q(i+1)+19.*fn-5.*fn_1+fn_2)./24;142Fn,fn-1,fn-2,fn-3?

根據(jù)定義：fn=feval(@fx,x(i),y(i));fn_1=feval(@fx,x(i-1),y(i-1));fn_2=feval(@fx,x(i-2),y(i-2));fn_3=feval(@fx,x(i-3),y(i-3));如何求解前三個點？R-K方法！143R-K方法求解前三個點fori=1:3;x(i)=a+(i-1).*h;K1=h.*feval(@fx,x(i),y(i));K2=h.*feval(@fx,x(i)+h/2,y(i)+K1/2);K3=h.*feval(@fx,x(i)+h/2,y(i)+K2/2);K4=h.*feval(@fx,x(i)+h,y(i)+K3);y(i+1)=y(i)+(K1+2*K2+2*K3+K4)./6;end144就這么簡單么？???Undefinedfunctionormethod'c'forinputargumentsoftype'double'.Errorin==>Adamsat36m(i+1)=p(i+1)+251.*(c(i)-p(i))./270;i=4,p(4)沒有定義

p(4)=0同樣，c(4)也沒有定義c(4)=0145146在區(qū)間[0,8]上結(jié)果如何？147Physicistscandoanything!h=0.01例：BicycleRacingWhatfactorsdeterminetheultimatespeedofabicycle?Howtoestimatethespeedforarealisticcase?Asimplercase:frictionisignoredEquationofmotion:Constantpowermodel:《ComputationalPhysics》N.J.GiordanoandH.Nakanishi,(2007)選?。簃=70Kg,v0=4m/s,deltat=0.1s.ConsidertheairresistanceC:dragcoefficientA:frontalareaoftheobjectρ:densityofairIntheabovecalculations,theinitialvelocityisnotzero.Whathappensifv0=0?參數(shù)選?。篊=0.5,A=0.33m2,m=70Kg,v0=4m/s,deltat=0.1s例：ThetrajectoryofacannonshellInavacuum,howtodescribethetrajectoryofacannonshell?Equationofmotion:《ComputationalPhysics》N.J.GiordanoandH.Nakanishi,(2007)Initialspeed=700m/sWithairresistanceWhatistheforcedirection?FxFyEquationofmotion:300500B/m=4x10-5/mConsider

Airdensityvs.heightIsothermalmodel:Y0=kBT/mg~1.0x104mAdiabaticapproximation:forairDragforcebecomes159作業(yè)1.2.量子力學中描述一維無限深勢阱的微分方程為：

采用打靶法求解本征值l和本征函數(shù)f161第三章偏微分方程的數(shù)值求解方法偏微分方程

162

雙曲型方程163xt

即對流方程，這里a為常數(shù)。顧名思義，對流方程描述流體中某種物質(zhì)的物理量的變化規(guī)律，例如傳熱過程中溫度的演化規(guī)律，流體中溶質(zhì)的濃度等的變化規(guī)律等。

有限差分方法在有限差分方法中，我們放棄了微分方程中獨立變量可以取連續(xù)值的特征，而關注獨立變量離散取值后對應的函數(shù)值。原則上，這種方法仍然可以達到任意滿意的計算精度，因為方程的連續(xù)數(shù)值解可以通過減小獨立變量離散取值的間隔，或者通過離散點上的函數(shù)值插值來近似得到。有限差分法具體操作分兩部分：1.用差分（差商）代替微分方程中的微分（微商），將連續(xù)變化的變量離散化，從而得到差分方程組的數(shù)學形式。2.求解差分方程組。164在第一步中，通過所謂的網(wǎng)絡分割法，將函數(shù)定義域分成大量相鄰而不重合的子區(qū)域，通常采用的是規(guī)則的分割方式。網(wǎng)絡線劃分的交點稱為節(jié)點。數(shù)值求解的關鍵就是要應用適當?shù)挠嬎惴椒?，求得特定問題在所有這些節(jié)點上的近似值。一個函數(shù)在x點上的一階和二階微商，可以近似地用它所臨近的兩點上的函數(shù)值的差商來表示。具體方法165差分格式的建立hτtxO差分法是處理偏微分方程的一種常用數(shù)值解法166

167其次，分別利用一階微商的向前和向后差商公式，以及中心差商公式構造偏微分方程的差分格式。

將這些差商公式分別代入方程，即可構造出不同的差分格式，并建立相應的差分方程。168

k,jk+1,jk,j+1

169

k-1,jk,jk,j+1k+1,j

k-1,jk,jk,j+1170

171

172差分格式的相容性和收斂性從上面建立的差分格式和他們相應的截斷誤差可以看出，當步長h,τ0時，各差分格式的截斷誤差也0。這表明：差分格式的極限形式就是相應的微分方程這稱為：差分格式和相應的微分方程是相容的通常是希望用差分格式計算出來的解在h,τ0時，能逼近于相應微分方程的解，這就是收斂性問題。相容性和收斂性的概念是不同的，相容是必備條件，收斂才是目的。一般來說相容并不能保證收斂，但在一定條件下，兩者等價。差分格式的收斂性–依賴區(qū)域t=0t=1t=2t=3BCP若用左圖差分格式計算P點的值，依賴于初始線段BC上的網(wǎng)絡節(jié)點的值；對于右圖，則依賴于初始線段AC上節(jié)點的值。BC就是格式左的依賴區(qū)域；AC就是格式右的依賴區(qū)域t=0t=1t=2t=3BCPA173依賴區(qū)域跟：差分格式、時間層相關差分格式的收斂性與偏微分方程的依賴區(qū)域有關。Courant條件差分格式收斂的一個必要條件是，差分格式的依賴區(qū)域應包含微分方程的依賴區(qū)域?！狢ourant條件174

可見：可通過控制t的步長τ和x的步長h來保證收斂。175

差分格式的穩(wěn)定性數(shù)值不穩(wěn)定性：用差分格式求初值問題數(shù)值解，每步計算都會引入舍入誤差，而初值問題還是逐層計算的，所以誤差還會逐層傳播。若誤差傳播越來越大，就可能會淹沒真解。176

差分格式的穩(wěn)定性不但與差分格式有關，還與網(wǎng)格步長的比值有關。如果在一定步長下，差分格式是穩(wěn)定的，則稱該差分格式為條件穩(wěn)定。如果對任何步長，差分格式都是穩(wěn)定的，則稱該差分格式為無條件穩(wěn)定。如果差分格式對任何步長都不穩(wěn)定，則稱該差分格式為完全不穩(wěn)定。177一種討論穩(wěn)定性的常用方法—Fourier方法常用Fourier方法(也稱VonNeumann方法)研究穩(wěn)定性。

178

179

180

精度穩(wěn)定性條件逆風差分格式FTCS差分格式不穩(wěn)定蛙跳差分格式Lax差分格式Lax-Wendroff差分格式Lax等價定理(相容性、穩(wěn)定性、收斂性的關系)：對于一個適定的初值問題，一個與它相容的差分格式收斂的充分必要條件是這個格式的穩(wěn)定性。181對策：182小結(jié)有限差分求解偏微分方程：離散化

建立差分方程討論差分方程的穩(wěn)定性和收斂性具體地，離散化、構造網(wǎng)格、用網(wǎng)格函數(shù)做近似解、用差商代替微商、同Talyor展開建立差分方程、或用特征線構造差分方程截斷誤差的引入，相容的概念差分方程的穩(wěn)定性就是舍入誤差是否會逐層傳播無限增加，判斷其穩(wěn)定性可以用圖表，也可以用Fourier方法Lax定理和Courant條件例：波動方程183略去小量，并代入方程可得，184185

1862、利用二階泰勒公式構造第二行的值。

這樣就得到了一個改進的第二行的近似差分公式。練習：求解兩端固定的弦上的振動形式187188拋物型方程的差分解法一維熱傳導方程：

根據(jù)給定的初值和邊值，該方程有以下兩類定解條件：1892.初邊值混合問題：也就是求解函數(shù)在一個有限區(qū)域內(nèi)既滿足初始條件也滿足邊界條件的定解。

190

191

如果第j層的近似值已經(jīng)求得，通過以上公式可以計算得到網(wǎng)格中第j+1層的近似值。該差分格式稱為顯式差分格式。

192

193

這就是一維熱傳導方程的Crank-Nicholson法，也稱為六點差分格式。194邊界條件的差分公式：

對于第二類邊值問題，含有微商，需要利用差商公式。

195

196

197

198

精度穩(wěn)定性條件顯式差分格式隱式差分格式無條件穩(wěn)定Richardson差分格式完全不穩(wěn)定六點差分格式無條件穩(wěn)定幾種差分格式的精度和穩(wěn)定性條件總結(jié)如下：

199二維熱傳導方程：

原則上可以用前面介紹的幾種差分格式

200方程中有兩個空間變量，一個時間變量，選取網(wǎng)格如下：顯式差分格式：

201

202對于n維熱傳導問題，其顯式差分格式的穩(wěn)定性條件為可以看出，隨著維度的增大，穩(wěn)定性條件要求減小時間步長，這樣會導致計算量越來越大。下面簡要介紹熱傳導方程的交替方向法。

203交替方向法的基本思想是將一個二維問題，逐漸轉(zhuǎn)化為兩個一維問題來處理。

Peaceman-Rachford差分格式

該方法將兩個二階差商分開，在每一個時間層上只計算其中一個，另一個用已知值代替，這樣就得到僅在一個方向上是隱式的差分格式，比較容易求解。對于第二個時間層，在另外一個方向上重復以上過程。這樣相繼的兩個時間層構成了完整的一步。

204

205

206

207橢圓型方程的差分解法橢圓方程最簡單的形式是LapLace方程：

和Poisson方程：

主要定解條件是邊值條件，求函數(shù)u(x,y)在區(qū)域內(nèi)滿足微分方程和邊界條件的解。選取適當網(wǎng)絡，將微分方程離散成差分方程當網(wǎng)絡步長h—>0時，差分方程的準確解是否收斂到微分方程的解解相應的代數(shù)方程組208以Poisson方程第一類邊值問題為例：

差分方程的建立：xOy平面用兩組平行線的網(wǎng)格離散

h,

t為步長，交點為節(jié)點209屬于W的點稱為內(nèi)點，若一個節(jié)點的所有四個相鄰節(jié)點均屬于W+G,

則稱為正則內(nèi)點若一個節(jié)點，其相鄰的四個節(jié)點中至少有一個不屬于W+G,

則稱為非正則內(nèi)點。正則內(nèi)點的差分與前面一致：

210另一種五點格式

當h=t,五點格式

211對于非正則內(nèi)點，需要另行處理：1.直接轉(zhuǎn)移2.作線性差值：可以利用T和R作線性插值來估算S的值，也就是

RSThd圖中S為非正則內(nèi)點，可以將其緊鄰的邊界與網(wǎng)格的交點T的值直接賦給S，即

這里d為S和T之間的距離Q212

分子動力學引言分子動力學基礎知識分子動力學的基本步驟分子間相互作用勢平衡態(tài)分子動力學模擬應用舉例

分子動力學方法是一種計算機模擬實驗方法。通過模擬實驗或微觀粒子運動過程，研究系統(tǒng)宏觀性質(zhì)的微觀機理等。

什么是分子動力學方法？

MolecularDynamics1957年至1959年間，Alder和Wainwright首次采用分子動力學模擬研究了硬球模型系統(tǒng)中無序-有序轉(zhuǎn)變問題，證實了硬球模型系統(tǒng)的一級相變行為。1964年，Rahman采用分子動力學方法研究了Lennard-Jones相互作用下液態(tài)Ar的結(jié)構和動力學性質(zhì)。1967~1971年，Verlet模擬了更復雜的液體——水、溶鹽、聚合物和蛋白質(zhì)等。1974年，Rahman和Stillinger使用分子動力學模擬研究真實系統(tǒng)，研究對象為液體水。1977年，McCammon第一個蛋白質(zhì)的分子動力學研究，研究對象是牛胰島素抑制劑。分子動力學的簡要發(fā)展歷程分子動力學實例蛋白質(zhì)的折疊通過分子動力學模擬可以研究蛋白質(zhì)的折疊路徑。分子動力學實例水化（溶劑化）描述生物大分子周圍的溶劑（水）分子分布。抹香鯨肌紅蛋白的水化水化位置

（藍，實驗）水化數(shù)密度最大值

（黃，MD模擬）分子動力學實例脆性斷裂分子動力學實例：晶化的分子動力學模擬FCC（綠）HCP（紫）BCC（紅）分子動力學的引入對于大量粒子體系，例如一滴水包含1021個水分子，幾乎不可能計算所有涉及到的狀態(tài)的物理量，再求其平均。而觀測量確是平均值，如何比較計算值和觀測值？近似：假設感興趣的物理量在體系到達一定尺度后，就基本不隨尺度變化而變化了。因此：對一些列有限的態(tài)做統(tǒng)計平均。隨機模擬方法（MC）和確定性模擬方法（MD）MD：由運動方程描述，利用統(tǒng)計方法得到體系的靜態(tài)和動態(tài)性質(zhì)，從而得到體系的宏觀性質(zhì)?？梢钥醋魇求w系在一段時間內(nèi)發(fā)展過程的模擬。不存在任何的隨機因素

確定性模擬方法分子動力學簡介分子動力學是在原子、分子水平上求解多體問題的重要的計算機模擬方法。通過求解所有粒子的運動方程，分子動力學方法可以用于模擬與粒子運動路徑相關的基本過程。在分子動力學中，粒子的運動行為是通過經(jīng)典的Newton運動方程所描述。系統(tǒng)的所有粒子服從經(jīng)典力學的運動規(guī)律，它的動力學方程就是從經(jīng)典力學的運動方程——拉格朗日方程和哈密頓方程導出。分子動力學方法的優(yōu)勢分子動力學得天獨厚的優(yōu)勢是能夠模擬分子的運動軌道，即能夠提供分子運動和變化的最為微觀的、細致的信息。這對傳統(tǒng)的力學、熱力學、光譜學等實驗方法都是難以辦到的。現(xiàn)代分析理論和實驗測量手段要準確地描述物質(zhì)的非平衡態(tài)性質(zhì)是非常困難的。而分子動力學方法能夠計算得到運行軌道，通過非常直接的計算公式給出感興趣的非平衡態(tài)物理量。這種特點對現(xiàn)今所有的理論、實驗、計算手段來講都是獨一無二的。正是這種不可替代的作用，使得分子動力學在當今自然科學研究中獲得了廣泛和深入的應用。實際使用的限制在實際應用中，分子動力學模擬方法和隨機模擬方法都面臨著兩個基本限制：有限觀測時間的限制有限系統(tǒng)大小的限制通常人們感興趣的是體系在熱力學極限下（即粒子數(shù)目趨于無窮時）的性質(zhì)。但是計算機模擬允許的體系大小要比熱力學極限小得多，因此可能會出現(xiàn)有限尺寸效應。為了減小有限尺寸效應，人們往往引入周期性、全反射、漫反射等邊界條件。當然邊界條件的引入顯然會影響體系的某些性質(zhì)。數(shù)值求解中的問題對體系的分子運動方程組采用計算機進行數(shù)值求解時，需要將運動方程離散化為有限差分方程。常用的求解方法有歐拉法、龍格-庫塔法等。數(shù)值計算的誤差階數(shù)顯然取決于所采用的數(shù)值求解方法的近似階數(shù)。原則上，只要計算機計算速度足夠大，內(nèi)存足夠多，我們可以使計算誤差足夠小。分子動力學的基礎知識運動方程的導出求解運動的數(shù)值算法：Verlet算法等相互作用勢周期性邊界條件截斷距離最小像力約定分子動力學運動方程數(shù)值求解系統(tǒng)的動力學機制決定運動方程的形式。在分子動力學方法處理過程中，方程組的建立是通過對物理體系的微觀數(shù)學描述給出的。在這個微觀的物理體系中，每個分子都各自服從經(jīng)典的牛頓力學。每個分子運動的內(nèi)稟動力學是用理論力學上的哈密頓量或者拉格朗日量來描述，也可以直接用牛頓運動方程來描述。這種方法可以處理與時間有關的過程，因而可以處理非平衡態(tài)問題。運動方程采用分子動力學方法時，必須對一組分子運動微分方程做數(shù)值求解。從計算數(shù)學的角度來看，這是個求一個初值問題的微分方程的解。實際上計算數(shù)學為了求解這種問題已經(jīng)發(fā)展了許多的算法。但是并不是所有的這些算法都可以用來解決分子動力學問題?？臻g描述在空間描述如何物體的運動，如果其本身的大小可以忽略時，就可以將其看作是粒子（或質(zhì)點）。粒子描述：空間位置：r速度：v=dr/dt加速度：a=dv/dt=d2r/dt2若一個系統(tǒng)由N個粒子組成，則粒子描述：空間位置：r1，r2，…

，rN笛卡爾坐標系，粒子有3N個自由度設系統(tǒng)有s個自由度廣義坐標：q1，q2，…

，qs廣義速度：

拉格朗日(Lagrange)方程

哈密頓方程粒子運動方程的數(shù)值解法

采用分子動力學方法時，必須對一組分子運動微分方程做數(shù)值求解，這是一個微分方程的初值問題。

速度Verlet算法

速度Verlet算法與Verlet算法是等價的。將改寫為代入可得方程（1）。

利用方程

Verlet的速度形式Verlet速度形式的算法比前一種算法好些。它不僅可以在計算中得到同一時間步上的空間位置和速度，而且數(shù)值計算的穩(wěn)定性也提高了。一般情況下，對于能量確定的系統(tǒng)不可能給出精確的初始條件。這時需要先給出一個合理的初始條件，然后在模擬過程中逐漸調(diào)節(jié)系統(tǒng)能量達到給定值。

如果消去速度項，該算法本質(zhì)上等價于Verlet算法。但比Verlet算法優(yōu)越，其好處在于不需要處理兩個大數(shù)的相減，從而盡可能地減小精度丟失。但相對于速度Verlet算法而言，該算法不能得到同一時刻的坐標和速度。三種Verlet算法對比Verlet算法速度Verlet算法Leap-Frog算法《Computersimulationofliquids》（1986）Gear預測校正法預測校正法是分子動力學模擬中的常用算法之一，其基本思想是Taylor展開：

根據(jù)預測值和誤差可以得到各個量的校正值

方程組

可以表示成矩陣形式方程組

可以表示成矩陣形式Gear預測校正因子

1/65/611/3

19/1203/411/21/12

3/20251/360111/181/61/60保留不同階數(shù)的泰勒展開所得到的的修正系數(shù)如下Gear預測校正因子

1/65/611/3

19/903/411/21/12

3/16251/360111/181/61/60得到的的修正系數(shù)如下原子間相互作用勢

要進行分子動力學模擬就必須知道原子間的相互作用勢。在分子動力學模擬中，我們一般采用經(jīng)驗勢來代替原子間的相互作用勢，如Lennard-Jones勢、Morse勢、嵌入原子勢（Embedded-atommethodpotential）等。