大学化学, 2021, 36(2): 2003046-0 doi: 10.3866/PKU.DXHX202003046

自学之友

SymPy在结构化学课程教学中的应用

侯成,

Application of SymPy in the Teaching of Structural Chemistry Course

Hou Cheng,

通讯作者: 侯成, Email: houcheng@gxnu.edu.com

收稿日期: 2020-03-14   接受日期: 2020-04-13  

基金资助: 广西师范大学本科教学改革工程项目.  2019JGB32

Received: 2020-03-14   Accepted: 2020-04-13  

Abstract

Structural chemistry is one of the core basic courses of chemistry-related majors. In the course of teaching structural chemistry, complex formula derivation, mathematical operations and function image drawing have always been the most difficult points of the course. This article introduces the application of the Python scientific computing library SymPy for symbolic calculation and function image drawing in the teaching of structural chemistry through specific examples. The programming language involved in the library is easy to learn, understand, and operate. It can effectively help students overcome their fear to learning, deepen their understanding of structural chemistry, and also help subsequent learning of computational chemistry courses.

Keywords: Structural chemistry ; Teaching reform ; Python ; SymPy

PDF (2535KB) 元数据 多维度评价 相关文章 导出 EndNote| Ris| Bibtex  收藏本文

本文引用格式

侯成. SymPy在结构化学课程教学中的应用. 大学化学[J], 2021, 36(2): 2003046-0 doi:10.3866/PKU.DXHX202003046

Hou Cheng. Application of SymPy in the Teaching of Structural Chemistry Course. University Chemistry[J], 2021, 36(2): 2003046-0 doi:10.3866/PKU.DXHX202003046

结构化学是本科化学课程中最为重要的理论课程,在化学专业课程体系中扮演着承上启下的作用[1]。在过去的教学过程中,学生学习困难的地方主要表现为两个方面,第一个是公式的推导涉及较多复杂运算,而数学一直是部分同学学习过程中的“拦路虎”,另一方面,对于课程中所涉及复杂函数的图像较为抽象,不少同学在理解上还有所欠缺。

针对这两方面的问题,除了积极补充数学基础知识以外,在授课中融合计算机编程语言,可以有效消除学生对数学的恐惧,加深学生对物理图像本质的理解。对于公式的推导所涉及的符号运算和函数图像绘制,常见的科学计算软件如Matlab、Mathmatica、Maple等均可以实现[2]。然而这些软件也有着以下缺点:其均为商业性收费软件;运行对计算机配置有要求;每一种软件都有各自的编程语言,增加学习难度。近年来,Python语言已经成为目前最为流行的编程语言,其优势主要表现在以下几个方面:免费开源;语言简洁,易学易懂; 相关程序运行轻盈,占用内存资源较少;配合科学计算库SymPy,可以集符号计算与图像绘制功能为一体,适合课堂教学演示和计算化学实验的开展[3]

本文介绍了Python开源库SymPy的特点及安装方法,通过具体案例介绍其在结构化学课程授课过程的应用,并给出代码以供参考。

1 SymPy库安装与使用

Sympy库是基于Python语言的开源库(图 1(a)),主要针对于科学计算过程中的符号计算,如微积分运算,微分方程求解,因式分解和矩阵运算等。它安装非常简便,对于已经安装了Python运行环境的计算机,直接在cmd命令行中运行pip install sympy命令即可完成安装。对于没有安装python环境的计算机,笔者推荐直接安装Python环境管理工具Anaconda,在该套装软件中,包含了常用的Python科学计算库,同时提供了Web交互式编程程序Jupyter Notebook以及强大的集成开发环境Spyder,特别是Jupyter Notebook仅仅通过浏览器就可以编写运行Python代码(图 1(b)),非常适合课堂教学演示,接下来我们将通过具体实例介绍其基本功能,实例运行环境为Python 3.7.4版本。

图1

图1   SymPy库介绍页面(a)以及Jupyter Notebook运行界面(b)


2 符号计算功能应用示例

2.1 微分方程求解

在结构化学课程中,微分方程的求解是反复出现的内容。SymPy强大的符号计算功能可以轻松完成微分方程通解的计算,我们以一维势箱模型薛定谔方程求解为例[4],在Jupyter Notebook介绍其求解微分方程的功能,一维势箱模型中势箱内的薛定谔方程为:

若设k2 = 2mE/$ {\hslash }^{2} $ = $ {{p}^{2}/\hslash }^{2} $,则方程可简化为:

利用SymPy求解该微分方程代码如下,

1. from sympy import *            #导入SymPy库

2. x, k=symbols('x k')            #定义x, k为符号

3. f=Function('f')            #定义f为函数

4. dsolve(diff(f(x), x, 2)+k**2*f(x), f(x))            #dsolve函数求解微分方程

在Jupyter Notebook中点击运行,即可直接得到微分方程的复数形式通解y(x) = c1exp(−ikx) + c2exp(ikx),运行结果如图 2(a)所示。

图2

图2   SymPy符号计算功能应用示例

(a)微分方程求解应用;(b)积分运算应用;(c)矩阵运算应用


2.2 积分运算

除了微分方程以外,积分运算也是结构化学的公示推导和习题练习中常常出现的内容,在这里我们以求一维势箱粒子位置的期望值为例[5],对SymPy积分运算过程进行介绍。

对于一维势箱,其不含时的波函数为:

其中n为正整数,位置期望值可以利用期望值公式求得:

利用SymPy求解该微分方程代码如下,

1. from sympy import *            #导入SymPy库

2. x, l = symbols('x l')            #定义x, l, 为符号

3. n = Symbol("n", positive=True, integer=True, nonzero=True)            #定义n为符号, 且为非零正整数

4. y = integrate(sin(n*pi*x/l)*x*sin(n*pi*x/l), (x, 0, l))            #进行积分运算

5. result = y*2/l            #乘以归一化系数

6. result

在Jupyter Notebook中点击运行,即可直接得到积分值为l/2,运行结果如图 2(b)所示。

2.3 矩阵运算和行列式展开

结构化学课程中涉及分子轨道理论的课程章节中,涉及行列式以及矩阵的运算较多,和常见的科学计算软件一样,SymPy也可以进行矩阵运算和行列式展开的操作,我们以苯分子的久期行列式为例,对这一功能进行简要说明。对于苯分子,在休克尔分子轨道近似条件下,其久期行列式为:

该六阶行列式其展开过程较为繁琐,在授课过程中我们通常利用对称性进行简化,减少计算量。利用SymPy可以直接将其展开,其代码如下:

1. from sympy import *            #导入SymPy库

2. x = symbols('x')            #定义x为符号

3. M = Matrix([

4. [x, 1, 0, 0, 0, 1],

5. [1, x, 1, 0, 0, 0],

6. [0, 1, x, 1, 0, 0],

7. [0, 0, 1, x, 1, 0],

8. [0, 0, 0, 1, x, 1],

9. [1, 0, 0, 0, 1, x]])            #定义矩阵

10. M.det()            #行列式展开

11. solve(M.det(), x)            #求解方程

点击运行后,可以直接给出其展开式x6 − 6x4 + 9x2 − 4 = 0,进一步利用solve函数以x为变量求解方程x6 − 6x4 + 9x2 − 4 = 0,可得方程的根为[−2, −1, 1, 2],运行结果如图 2(c)所示。

3 函数图像绘制应用示例

3.1 势箱模型波函数绘制

波函数的绘制对于理解粒子的波动性以及能量的量子化等抽象概念有着重要的意义[6]。我们首先以一维、二维势箱的波函数图像绘制,其具体代码如下,

一维势箱波函数绘制代码如下:

1. from sympy import *            #导入SymPy库

2. n=3            #给出量子数

3. x = Symbol('x')            #定义x为符号

4. plot((sqrt(2)*sin(n*pi*x)), (x, 0, 1))            #画出函数图像,x取值范围为0到1

二维势箱波函数绘制代码如下:

1. from sympy import *            #导入SymPy库

2. from sympy.plotting import plot3d            #导入SymPy三维函数图像库

3. nx=3;ny=3            #给出量子数

4. x, y = symbols('x y')            #定义x, y为符号

5. f = 2*sin(nx*pi*x)*sin(ny*pi*y)            #给出波函数表达式

6. plot3d(f, (x, 0, 1), (y, 0, 1), title='2D particle in a box', xlabel='x', ylabel='y')            #画出函数图像

点击运行后,即可绘制对应波函数图像,如图 3所示。

图3

图3   SymPy绘制波函数图像

(a)一维势箱波函数图像绘制;(b)二维势箱波函数图像绘制


3.2 原子轨道波函数图像绘制

对于复杂函数的图像绘制而言,Sympy库中已经内置了物理上常见的特殊函数和氢原子轨道波函数,可以方便调用,调用球谐函数的运行代码如下:

1. from sympy import Ynm, Symbol, simplify            #导入SymPy库中球谐函数

2. from sympy.abc import n, m            #导入参数

3. theta = Symbol("theta")            #定义theta为符号

4. phi = Symbol("phi")            #定义phi为符号

5. Ynm(1, 0, theta, phi).expand(func=True)            #给出n=1, m=0时的球谐函数表达式

点击运行后,即可输出球谐函数Y1, 0表达式,也是Pz轨道波函数的角度部分,见图 4(a)

图4

图4   SymPy特殊函数调用和氢原子轨道波函数等高线绘制

(a)球谐函数调用实例;(b)氢原子轨道波函数调用实例;(c) 2Pz轨道波函数等高线图绘制


调用氢原子轨道波函数代码如下:

1. from sympy import *            #导入SymPy库

2. from sympy.physics.hydrogen import Psi_nlm            #导入SymPy库中氢原子波函数

3. Z=Symbol("Z", positive=True, integer=True, nonzero=True)            #Z为核电荷数

4. r=Symbol("r", real=True, positive=True)            #定义r为正实数

5. theta=Symbol("theta", real=True)            #定义theta为符号

6. phi=Symbol("phi", real=True)            #定义phi为符号

7. Psi_nlm(2, 1, 0, r, phi, theta, 1)            #给出n=2, l=1, m=0时的氢原子波函数表达式

点击运行后,即可输出氢原子2Pz轨道波函数y2, 1, 0表达式,见图 4(b)

不可否认的是,SymPy库和专业科学计算软件在作图上还有一定差距,譬如无法做轨道等高线图和等值面图,但由于其基于Python语言,可以配合其它的Python计算绘图库,对原子轨道的复杂波函数图形进行绘制,我们以Pz轨道等高线图为例,基于matplotlib和numpy库进行绘制,具体代码如下:

1. import matplotlib.pyplot as plt            #导入matplotlib图形库

2. import numpy as np            #导入numpy库

3.

4. n = 256

5. x = np.linspace(-10, 10, n); z = np.linspace(-10, 10, n)

6. X, Z= np.meshgrid(x, z)            #生成二维网格

7. r = np.sqrt(X**2+Z**2)

8. theta = np.arctan(X/Z)

9. f = np.sqrt(2)/(8*np.sqrt(np.pi))*r*np.exp(-r/2)*np.cos(theta)            #给出函数表达式

10. plt.figure(figsize=(5, 5))

11. plt.ylabel('z'); plt.xlabel('x')

12. plt.grid(True)

13. C = plt.contour(X, Z, f, 8)            #做等高线图

14. plt.clabel(C, inline=True, fontsize=10)

15. plt.show()

点击运行后,即可输出Pz轨道波函数的等高线图,见图 4(c)

4 教学设计样例:利用SymPy处理烯烃电子结构

4.1 学情分析

本节段为结构化学多原子分子电子结构章节,属于量子化学部分最为重要的内容,主要介绍基于休克尔分子轨道理论近似处理共轭分子体系电子结构的具体方法。该部分知识点主要为理论的基本内容和具体计算方法介绍,是学生进一步学习高阶课程如计算化学的基础。

本节课授课对象为大学三年级学生,开设时间为大三下学期,前节内容已经讲授了休克尔分子轨道理论的主要内容,以及对丁二烯电子结构处理的具体步骤,本节课将对上次内容进行简要复习,并介绍通过SymPy求解类似问题的具体步骤,帮助学生理清思路,消除疑惑,开拓视野,教学过程见表 1

表1   利用SymPy处理丙烯自由基分子的教学设计

索引师生活动教学内容设计意图
导入课程介绍休克尔分子轨道理论内容以及处理丁二烯电子结构的主要步骤。    休克尔分子轨道理论基本假设
(1) σπ键分离;(2)积分简化
    具体计算步骤:以丁二烯为例
(1)根据分子结构构建久期方程组
(2)方程组系数行列式展开
(3)求解行列式展开式中的根并得到能量的表达式
(4)将不同的根依次带入求解该能级分子轨道系数
总结复习上节所学知识,为本节课内容进行铺垫。
提出问题让学生上台黑板书写丙烯自由基的久期行列式,并尝试展开。丙烯自由基的久期行列式为:
$\left| {\begin{array}{*{20}{c}} x&1&0 \\ 1&x&1 \\ 0&1&x \end{array}} \right| = 0$
展开后表达式为:
${x^3} - 2x = 0$
点评学生书写结果,指出随着所处理体系的增大,行列数阶数上升,行列式展开计算量增大。
提出解决方案介绍SymPy矩阵操作函数和求解方程,线性方程组的solve函数。以二阶行列式利用计算机演示行列式展开,代码如下:
1. from sympy import *
2. a, b, c, d = symbols('a b c d')
3. M = Matrix([[a, b], [c, d]])
4. M.det()
运行后得到adbc
以二元一次方程组为例,演示求解线性方程组方法,代码如下:
1. from sympy import *
2. x, y= symbols("x y")
3. f1 = x+y-35
4. f2 = 2*x+4*y-94
5. solve([f1, f2], [x, y])
运行后得到{x: 23, y: 12},即为方程的根
以较为简单的例子,由浅入深讲述程序实现方法,重点帮助学生了解一些基本命令。
具体实例介绍基于休克尔分子轨道理论,利用SymPy处理丙烯自由基电子结构的具体步骤。    构建久期行列式并展开,代码如下:
1. from sympy import *
2. x= symbols('x')
3. M = Matrix([[x, 1, 0], [1, x, 1], [0, 1, x]])
4. M.det()
    利用solve函数求解方程,代码如下:
solve(M.det(), x)
求得方程的根为0,$ \sqrt{2},-\sqrt{2} $,由于x = (aE)/b,分别带入即可解得三个能级的能量表达式

    将每个根分别带入到之前得方程组可以对应能级的波函数系数,在这里我们取x = $ -\sqrt{2} $,同时考虑归一化条件,构建线性方程组,代码如下:
1. c1, c2, c3 = symbols("c1 c2 c3")
2. x = -sqrt(2)
3. f1 = c1*x+c2
4. f2 = c1+c2*x+c3
5. f3 = c2+c3*x
6. f4 = c1**2+c2**2+c3**2-1
7. solve([f1, f2, f3, f4], [c1, c2, c3])
可以得到成键轨道波函数的系数为:
c1 = 1/2,$c2 = \sqrt 2 {\rm{/}}2$c3 = 1/2,波函数表达式为:
${\psi _1} = \frac{1}{2}{\phi _1} + \frac{{\sqrt 2 }}{2}{\phi _2} + \frac{1}{2}{\phi _2}$
同理可得其它轨道波函数表达式:
${\psi _2} = \frac{{\sqrt 2 }}{2}{\phi _1} - \frac{{\sqrt 2 }}{2}{\phi _2}$
${\psi _3} = \frac{1}{2}{\phi _1} - \frac{{\sqrt 2 }}{2}{\phi _2} + \frac{1}{2}{\phi _2}$
通过具体实例讲述,帮助学生进一步熟悉休克尔分子轨道理论处理共轭分子的基本步骤;帮助学生掌握利用计算机处理数学问题的基本方法。
课后练习布置作业,让学生自己编写处理丁二烯电子结构的代码,举一反三。    作业:
1)基于上述代码,基于休克尔分子轨道理论,编写处理1, 3-丁二烯的小程序,并和之前书本例题结果进行核对。
2)进一步编写求1,3-丁二烯原子电荷,键级,自由价的小程序。
让学生对所学到的知识进行巩固。

新窗口打开| 下载CSV


4.2 教学目标设计

a) 掌握休克尔分子轨道内容实质,具体运算步骤;

b) 掌握SymPy在矩阵运算,求解线性方程组等数学问题的具体方法;

c) 提高学生学习兴趣,培养举一反三的能力。

4.3 教学内容设计

4.3.1 内容纲要

2) 复习休克尔分子轨道基本内容以及处理丁二烯电子结构的方法;

2) 介绍SymPy中矩阵运算,线性方程组求解方法;

2) 基于休克尔分子轨道理论,利用SymPy处理丙烯自由基分子的电子结构。

4.3.2 教学重点

a) 休克尔分子轨道理论;

b) SymPy处理该问题的具体方法。

4.3.3 教学难点

数学运算以及代码编写。

4.4 教学策略方法

a) 板书标示教学法;

b) 迁移教学法。

4.5 教学过程设计

4.6 教学反思和学生反馈

本节课从已学知识出发,复习休克尔分子轨道的主要内容,理清计算流程,在此基础上介绍利用SymPy处理丙烯电子结构的具体操作方法,该方法操作方便,学生通过实际动手操作,加深了对于休克尔分子轨道理论的理解,也对计算机在处理量子化学问题有了初步的了解,为后续开展计算化学等相关课程打下了基础。

在学生反馈中,主要存在的问题主要是以下几个方面,一个是代码运行环境的安装,另外是python数学运算符号的意义和作用先后顺序,在进一步了解了这些知识后,由于python语言易学易用的特点,以及所见即所得的程序运行环境,大部分同学都可以完成布置的作业,有了程序设计语言作为有力工具,不少同学也由此对计算化学以及分子模拟等产生了浓厚的兴趣。

5 结语

本文通过SymPy对结构化学课程中进行符号运算和函数图像绘制进行了介绍和实例分析。在具体的教学过程中,可以有效降低学生的畏难情绪,提高学生的学习兴趣,让学生把注意力放到具体的科学问题和物理图像上,不仅对抽象的公式推导和函数图像有了进一步的认识,也接触了编程语言,为以后学习更加复杂的课程内容教学打下了基础。笔者将在未来的教学过程中继续探索Python在结构化学和物理化学课程中复杂问题的应用。

参考文献

万坚. 大学化学, 2017, 32 (4), 11.

URL     [本文引用: 1]

刘平; 张大顺. 计算机与应用化学, 2003, (4), 533.

DOI:10.3969/j.issn.1001-4160.2003.04.041      [本文引用: 1]

苗纯. 合肥师范学院学报, 2018, 36 (3), 47.

URL     [本文引用: 1]

陈兰; 孙宏伟; 赖城明. 大学化学, 2017, 32 (5), 77.

URL     [本文引用: 1]

陈兰; 孙宏伟; 赖城明. 大学化学, 2016, 31 (4), 70.

URL     [本文引用: 1]

孙宏伟; 陈兰. 大学化学, 2019, 34 (7), 100.

URL     [本文引用: 1]

/