科研课堂(三)——Monte Carlo 算法

随着计算机技术的发展,一种可以模拟随机变量与随机过程的方法——蒙特卡洛算法应运而生。与解析解或者近似解相比,Monte Carlo(MC)算法无需复杂的数学推导,实现起来也较为简单。 蒙特卡洛法其实是比较宽泛的一系列算法的统称,它的特点是假设概率分布已知,通过重复的随机采样来获得数值结果。比如根据大数定理,我们可以用采样得到的样本计算得到的样本均值来估计总体期望。又比如,积分的运算往往可以表示为随机变量在一个概率密度函数分布上的期望。 针对Ising模型,往往采用为马尔可夫链蒙特卡洛算法(MCMC)。此种方法具体为随机变量 x 的状态空间 S 上定义一个满足遍历定理的马尔科夫链 X = X1, X2, X3… ,使其平稳分布就是抽样的目标分布 p(x)。然后在这个马尔科夫链上进行随机游走,每个时刻得到一个样本。根据遍历定理,当时间趋向于无穷时,样本的分布趋近于平稳分布,样本的函数均值趋近函数的数学期望。关键就是设计这样一个马尔科夫链。 对于Ising模型,我们每一次随机选取一个格点,计算发生反转前后的能量差,如果能量降低,则使其翻转,如果能量升高,则依概率(\(P = \exp^(-\beta \Delta E)\))翻转。我们可以证明,按照此概率进行演进,系统达到平衡后的状态按照玻尔兹曼函数依能量分布,此过程也符合可逆过程的时间反演对称,即细致平衡原理。 我们可以通过简单的python代码实现用MCMC算法模拟Ising模型中,格点随温度、外加磁场的变化、并且做出相图。 以上蒙特卡洛算法是一种对离散问题的有效解决方案。但是物理上许多常见的相变问题,需采用连续介质假设。否则相对于模拟,问题本身的尺度过大,计算量过于庞大。相场模型(或者time-dependent Ginzburg-Landau模型)是解决连续介质假设下相变问题的有效途径。而该方程的阶数较高,数值实现起来需要一定的技巧。目前常见的方法是FFT-FDM,以及FEM。这将在下一讲中进行简要介绍。 扩展阅读 python skilearn的monte-carlo实现求定积分的库: https://scikit-monaco.readthedocs.io/en/latest/ 一种简单的马尔可夫链——随机行走 https://en.wikipedia.org/wiki/Random_walk 蒙特卡罗方法与机器学习的一个例子: https://www.deeplearningbook.org/contents/monte_carlo.html 阅读数: 1,267

科研课堂(二)——平均场近似和Ising模型

第一节课主要是叙述了“相”的概念以及两种不同的相边界模型(diffusive interface vs sharp interface)。第二节课的主题是相变(phase transition)。 相变是一个非平衡过程,描述起来比平衡过程要复杂。我们从最简单的Ising模型入手,对相变进行刻画。虽然数学模型较为简单,但是好处是可以定性地分析系统的特征温度、相变时对应的宏观物理量的变化。是研究复杂相变过程的基础。 Ising 模型最主要的假设有两个: 状态只有两种 +1 -1 (更为复杂的情况,状态可以描述为矢量等连续变化的量、甚至可以描述为波函数) 假设粒子只与最近的粒子发生相互作用(即只存在近场相互作用,无长程作用) 基于这两点出发,可以给出自旋粒子的哈密顿量。从哈密顿量出发、结合Gibbs概率密度分布(概率依能量指数分配),可推倒出粒子自旋量的概率测度(tanh分布,slide-12)。得到自旋的概率分布,就可以推出一系列的宏观物理量。如net magnetization、能量、磁滞回线等。 因此问题的核心就在于如何解出自旋的概率分布?平均场近似假设:每个粒子所受到的外界作用(课堂上的例子是外界磁场作用)是周围所有粒子的平均。因此在公式(5)中,可以将\(q_i\)用\(<q>\)代替。经过这种方法处理,方程大为简化(不需要考虑周围粒子的反转等瞬态过程)由平均场近似可以得出一系列很好的结果: 相变温度Tc 相变磁场Hc 磁滞回线 区分一级相变二级相变 \(<q> ~ (T_c-T)^{0.5}\) 平均场近似是不精确的。Onsager给出了二维Ising模型的解析接,序参量的critical exponent值不是0.5,而是0.125 。[1] 给同学们留下的作业是 1、回去动手算一下\(\tanh (Ax+b) = x\)的解。我用的是直接迭代。如果用牛顿迭代是不是收敛性更好?可否用程序证明? 2、算一下\(T^* = 0\) 和\(T^* > 0)时候的磁滞回线,看一下与低温时候有何区别? 资料: python环境管理包:anaconda: https://www.anaconda.com/ python的IDE:spyder 可以在anaconda中下载 有用的包: 自行搜索手册、可用anaconda安装或者用pip安装 基本数学工具和数组:numpy 科学计算、求解方程: scipy (部分替代mathematica) 数据表: pandas (可替代excel) 绘图:matplotlib (可替代origin) 参考文献 …

科研课堂(一)——初识相场模型

相场模型(Phase field model)主要是为了解决界面问题。思想来源最早可以追溯到Van der Waals在处理液体表面张力时候的非均匀界面假设,经过Landau&Lifshitz、Allen-Cahn、Cahn-Hillard发展,日趋完善。对于相场模型的表达式,在简单、低维、平衡的条件下,人们已经给出了解析解。但是对于复杂的边界条件以及非平衡状态下序参量的演化过程,计算量十分庞大。 随着计算机科学的不断发展,求解复杂边界条件下的相场方程成为可能。典型的例子是R. Kayabushi(1993)的晶枝生长、以及LQ Chen等的铁电畴变的算例。目前相场方法应用于众多领域,如 铁电 铁磁 合金 马氏体相变 多相合金 粉末冶金 等等。 那么,究竟什么是相场模型?对于一个界面,我们可以有两种处理方式: 对于气、液界面,Tolman [1] 和Ono Kondo [2] 证明,采用如下所示的自由能表示方式,最后得到的界面的序参量是间断变化的。 \(\psi(z) = \psi(z) [\rho(z), T]\) (1) 自由能取最低的条件为表面张力为0。 Van der Waals 提出,在界面处的自由能分布不是均匀的,而是一个跟组份比例有关的参数: \(\Psi(z) = \psi(z) [\rho(z), T]\) – \psi(\rho^{A,B}, T) (2) 与公式(1)相比增加的一部分自由能被称为surface excess。注意这里自由能的分布是位置坐标的函数。 Landau采用相同的思想,将铁磁体自由能的分布写为局部序参量(如magnetization)梯度的函数。Cahn等将梯度引入任意非均匀系统。在此假设下,自由能密度可以写为 \(\psi(c,\nabla c) = \psi _0(c) + G|\nabla c|^{2} \) (3) …