1、第28卷 第10期2011年lO月28同计算机与应用化学Computers and Applied ChemistryV0128，No10October,201l基于MATLAB的球谐函数M函数库郭仕恒，许庆陵，王东耀，苏育志，张建华(广州大学化学化工学院，广东，广州，510006)摘要：以MATLAB为设计语言，开发了系列用于原子轨道球谐函数Ylnl(e砷)如S、P、d、f原子轨道及sp3和d2sp3杂化轨道角度函数可视化的精确、快速、简便的m函数库。根据球坐标的自变量定义域：币【0，2兀】、0E【0：r】，设定离散的自变量值T，P，对球谐函数Yl。(0，p)及杂化轨道角度蛹数的解析式进行计
2、算机程序化，根据解析式计算对应的R值。然后由MATLAB矩阵变换命令IX，Y,Z=spb-2car“T,ER)和绘图命令mesh(X，Y，z)实现球谐函数YI。(0，神及杂化轨道角度函数的可视化。详细介绍了程序设计思路、编写方法，给出了球谐函数Y“(O砷)及杂化轨道角度函数的解析式数组算法表达式、mfile绘图程序和运行结果。关键词：MATLAB；球谐函数Ylm(0，Q)：nr文件中图分类号：06451 引言文献标识码：B 文章编号：1001-4160(201 1)10-1339-1342MATLAB作为强大的科学计算语言已经在许多领域得到了广泛应用fIl。许多领域的技术专家为了解决专业应用方
3、面的问题，已经开发了各领域的扩展工具箱4-5。MATLAB用于化学化工实验数据处理、图形绘制及化工过程模拟等在国内外已经有许多报道睁14】，但用于原子轨道球谐函数报道尚不多见【15m。因此，开发基于MATLAB的用于原子轨道球谐函数Yl。(9，p)可视化的m函数库，形成一个相对独立的结构化学图形应用软件系统，精确、简便、快速地实现原子轨道、波函数和电子云的可视化对结构化学和量子化学教学和研究具有重要意义。2 方法由文献【l 8】可得原子轨道球谐函数Ylm(O，审)及杂化轨道角度函数【I 8】的解析式。根据球极坐标的自变量定义域：0【O，兀】，甲0,2x，设定离散的自变量值，对球谐函数Yh(0，
4、Q)及杂化轨道角度函数的解析式进行计算机程序化，根据解析式计算对应的R值。为便于用MATLAB的绘图命令绘制三维图形，应将0，Q，r的值由球极坐标变换为直角坐标，此后，便可用MATLAB绘图命令绘制出各个球谐函数YllII(0，币)及的杂化轨道角度函数图形了。基于M朋rLAB的原子轨道球谐函数的可视化是基于MATLAB语言的数组算法，先将目标绘图函数转化为MATLAB可读化形式是关键的步骤，表l给出了基于MATLAB数组算法的原子轨道球谐函数可读化形式，表2给出了sp、dZsp3杂化轨道球谐函数MATLAB可读化形式。表l基于MATLAB数组算法的原子轨道球谐函数可读化形式Table 1 Re
5、adable expressions ofatomic orbital spherical harmonic鱼翌!竺!些i!：生2 11丛垒!生：!璺：函数类型lI。1 基于MATLAB的可读化形式Ys(Yoo)Ypz(V,o)Ypx(Y1+1)sqrt(14pi)sqrt(34pi)abs(cos(theta)gl峨34pi)+ab如叫nleta)cos(phi)Ydz2(Y2,o)Ydxz(Y2,+I)Ydyz(Y2，1)YdxyfY2+2)Ydx2矿fYz-2)Yfx3(Y3o)Yfxz2(Y3+1)Y瞄(Y3,-I)sqnll54#)。ab“sinoheta)*era(thin)cos
6、(phi)sqrttls4pi)+abs(sin(theta)*cos(the)。sin(phi)sq喇1516#)al，s(sin(thm)*sin(thm)sia(2“phi)sq婀lS16pi)abs(sin(dem)*si(ttEta)cos(2。phi)14+sq咂7倒)abs(S*cos(tima)*c06(tl蛐)c05(theta)-3*cot4theta)18+sqrt(42pi)alrsinoheta)(5*cos(theta)*e06(theta)-I)*咖帖)y=l8+sqr吖42pi)abstsinOt*ta)(5*cos(them)*cos(thm)-1)sin删3
7、)Yfz(x：y2)(Y3+2) l4sq咂105倒)abs(s咄岫)*re(thin)co武thcta)*co“2*phi)Yfzxy(Y3-2) !4sqrtfl05：pi)abs(sin(thin)*sin(them)cos(them)*sin(2*phi)Yfx(z23y2)fY，-+O 18*sqrt(70pi)*abs(sin(theta)*sin(tt把ta)sin(theta)co3*phi)YfH3x2-y2)fy3。) I8*sqrt(70pi)*abs(sin(theta)豳(岫回sin(theta)*sin(3*ph0)收稿ri期：2011-0324；修回H期：2011
8、0702作者简介：郭什恒(1967一)，男。广东人，实验师。Email：guoshihengI126com联系人：张建华(1959一)，男，湖北人教授硕士生导师。E-mail：蝴sinaeom万方数据1340 计算机与应用化学 2011，28(1 0)表2基于MATLAB数组算法的杂化轨道球谐函数可读化形式Table 2 Readable expressions of sp3，d2sp3 hybrid orbital!E垒!巫!皇!婴坐!垒唑!皇!垒!虹【!：里2 11坚垒!生垒星：杂化轨道 基于MATLAB的可读化形式”abs(sqrt(18pi)+sqrt(38pi)(sin(them)c
9、os(phi)+sqrt(38pi)+(sm(them)+sm(phi)+sqrt(38pi)(cos(them)abs(sqrt(I8pi卜sqrt(3Spi)+(sin(them)cos(phi)-sqrt(38pi)+(sin(tIIem)+sinlhi)-sqrt(38pi)(cos(them)abs(sqrt(18pi)-sqrt(38pi)+(sin(them)cos(phi)卜sqrt(38pi)+(sin(them)+sinIhi)-sqrt(38pi)+(cos(them)abs(sqrt(18pi)一sqrt(38pi)+(sin(them)cos(phi)sqrt(38p
10、i)+(sin(them)+sinlhi)+sqrt(38pi)+(cos(them)abs(sqrt(16)+sqrt(32)cos(theta)+sqrt(512)f3cos(the哪+cos(them)一1)abs(sqrt(16)-sqrt(32)cos(them)+sqrt(512)f3cos(theta)+cos(them)1)abs(sqrt(16卜sqrt(32)sin(them)。cos(phi)+sqrt(548)(DSP63 3cos(them)+cos(them)一1)+sqrt(1516)sin(them)+sin(theta)*cos(2phi)abs(sqrt(1
11、6)一sqrt(32)sin(them)cos(phi卜sqrt(548)(DSP64 3+cos(them)1lcos(them)-1)+sqrt(1516)+sin(them)sin(thera)*cos(2phi)abs(sqrt(16卜sqrt(32)sin(them)+sin(phi问rt(548)+(DSP65 3+cos(them)+cos(theta)-11sqrt(1516)+sin(them)+sin(them)+cos(2+phi)abs(sqrt(16)-sqrt(32)+sin(them)sin(phi)+sqrt(548)。f3DSP66 +cos(theta)co
12、s(them)-1)-sqn(1516)sin(them)+sin(theta)cos(2+phi)、3原子轨道球谐函数MATLAB程序及其运行结果31原子轨道球谐函数m励e程序编写31一Ys(Yo0)function v：函数声明figure；建立新的图形文件p=一pi：pi100：pi；t-0：2+pi150：pi；设置p,t变化范围和步长P,T-rneshgrid(p，t)；theta=pi2一P；phi=TR=Ys(theta，phi)；用定义的函数计算R值【X，Y,Z=sph2cart(T,P,R)；建立一个矩阵，并以变换后的T只R值为其赋值mesh(X，YZ)；对矩阵X，Y Z】做
13、网格图axis equal；xlabel(X)，ylabel(一Y)，zlabel(Z)；title(Ys)function y=Ys(theta，phi)；y=sqrt(14pi)；函数定义31 2 Vpz(Y,o)function y；figure；p=一pi：pi100：pi；t=0：2牛pi150：pi；【ET】-meshgrid(p，t)；theta=pi2-P；phi-T；R=Ypz(theta，phi)；【X，距】_sph2carKT,P,R)；mesh(X，YZ)；axis equal；xlabel(X)；ylabel(Y)；zlabel(Z)；title(Ypz)functi
14、on y=Ypz(theta，phi)；y=sqrt(34pi)+abs(cos(theta)；313其他原子轨道球谐函数可视化的mfile程序编写将上述Ys(Yo,o)和Vpz(V1o)，2个m-file进行比较发现其结构完全相同，在编写Ypz(V1o)的m-file时只需(1)将Ys(Yo)的 中所有的“Ys”换成“ ”；(2)将)的,omfi中le最后一句 VpzYs(Yoom-file y=sqrt(14pi)；”换成y=sqrt(34pi)+abs(cos(theta)；”，y=sqrt(14pi)和y=sqrt(34pi)abs(cos(theta)分别是原子轨道球谐函数Yoo和Y
15、1o的基于MATLAB的可读化形式。因此，在编写Ypx(Y,+1)、Ypy(Y1,-1)、Yclz：(Y2，o)、Ydxz(Y2+1)、Ydyz(Y2，1)、Ydxy(Y2,+2)、Yd)【2-y2(Y2,2)、Yfj【(Y3o)、YfxzZ(Y3+1)、Yfyz2(Y3_1)、Yfz(x2-y2)(Y3,+2)、Yfzxy(Y3。-2)、Yfx(z2-3y2)(Y3+3)、Yfy(3x2-y2)(Y33)的m-file程序时，仅需根据前述解析函数找出这些原子轨道球谐函数基于MATLAB的可读化形式，然后再做同样的替换即可。表l是这些球谐函数基于MATLAB的可读化形式。32杂化轨道角度函数可
16、视化的mfile程序编写分别以sp3杂化和d2sp3杂化为例说明杂化轨道角度函数可视化的mfile程序编写。sp3fF【ld2sp3杂化轨道角度函数分别由4个函数和6个函数组成见表2。使用hold on指令可以实现在同绘图窗口绘制分别有4个函数和6个函数组成的s口3和d2sp3杂化轨道角度函数的图形。用同样的方法可以实现其他杂化轨道角度函数的可视化。321 spj杂化轨道角度函数figure(1)p=pi2：pi30：pi2：t-0：2pi50：24pit【P，T=meshgrid(p，t)；the协=p沈一Pphi=T： R=SP3 l(theta，phi)； 【X，Y，Z】=sph2car
17、t(T，P，R)；mesh(X，Y，Z)，hold onRl=SP32(theta，phi)【X，Y，Z1=sph2cart(T，P，R1)；mesh(X，Y，Z1)，hold onR2=SP33(theta，phi)： X，Y，Z2=sph2cart(T，P，R2)；mesh(X，Y，Z2)，hold onR3=SP34(theta，phi)：【XY，Z3l=sph2cart(T，P，R3)；mesh(X，Y，Z3)，hold onfunction y=SP3 1(them，phi)y=abs(sqn(18pi)+sqrt(38pi)+(sin(theta)。cos(phi)+sqrt(38
18、pi)4(sin(them)+sin(phi)+sqrt(38pi)(cos(them)function y=SP32(theta，phi)y=abs(sqrt(18pi)+sqrt(38pi)+(sin(them)*cos(phi)sqrt(38pi)+(sin(theta)+sin(phi)一sqrt(38pi)+(cos(them)、function y=SP33(theta，phi)y=abs(sqrt(!8pi)一sqrt(38pi)(sin(theta)*eos(phi)+sqrt(38pi)+(sin(theta)+sin(phi)一sqrt(38pi)+(cos(them)fu
19、nction y=SP34(theta，phi)y=abs(sqrt(18pi)-sqrt(38pi)(sin(theta)cos(phi)-sqrt(38pi)+(sin(them)+sin(phi)+sqrt(38pi)+(cos(them)322 dZsp杂化轨道角度函数figure(1)p=-pi2：pi30：pi2：t-=o：2+pi50：2+pi：【P，T-meshgrid(p，t)；theta=pi2一P；phi=T：R=DSP6 1(theta，phi)：【X，Yz】_sph2cart(T，P，R)；mesh(X，Y，Z)，hold onR1=DSP62(theta,phi)：
20、 【X，Y，zq=sph2cart(T，P，R1)；mesh(X，Y，z 1)，hold onR2=DSP63(theta,phi)；【X，Y,Z2=sph2,car“T，P，R2)：mesh(X，Y，Z2)，hold onR3=DSP64(theta，phi)： 【X，Yz3】sph2cart(T，P，R3)；mesh(X，Yz3)，hold onn驼”拍晰溅娥饿饿兮l宝万方数据呈壁!堕!坐 塑竺堡：苎!薹王丛塑兰些塑蔓堂堕墼坚里墼塞 !兰!R4=DSP65(theta,phi)：X,Y,Z4=sph2e删(T，P，R4)；mesh(X,Y,Z4)holdR5=DSP66(theta,phi
21、)；x,Y,Z5】=sph2cart(T,P,RS)；mesh(X,Y,ZS)holdhillctiOIl y=DSP61(them,phi)y=abs(sq麒l6)+sqrt(32)印he蜘+哪咂5，12)(3d蜊a)+cos(theta)-1)function y=DSP62(theta,phi)，硝s【sqrt(16)-s呻(3国+sob出卜s【叫5，12)0*cos(e埘a)碱血抽)-1)fimcfion y=DSP63(theta,phi)y-=abs(sqrt(16)+sqrt(32)sinoheta)cos(phi)-hqrt(54S)*r3cos(theta)00s(n僦咖1卜
22、sq坦15，16)sin(the)sin(thet曲cos(2phi)function y=DSP64(theta,phi)y-=abs(sqrt(16)-sqrt(32)sin(thin)cosCphi)+sqrt(548)Ocos(them)+comma)-I)+sqrt(151们shi(thcta)sin(the)co“2*phi)fimction y=DSP65(thet&phi)y=abs(sqrt(16)+sqrt(32)+sin(thin)+sin(phi)+sqrt(54S)+(3cos(them)cos(them)1)-sqrt(1516)+sin(theta)sin(the
23、ta)COS(2hi)fimcfion y=DSP66(theta,phi)yffiabs(sqrt(I6)-sqrt(32)sin(rams)sin(phi)+sqrt(548)*仃cce(theta)+co“theta)-D-sqrt05l们sin(thma)+sin(them)+COSr2phi)33运行方法打开MATLAB程序，在命令窗)7(Command window)输入编写的m文件名回车即可，如：运行名为Ysm的程序，即在命令窗口输入Ys，回车。34图像编辑对图形进行编辑，方法是选择运行产生的图像窗口的编辑(Ed盍。下拉菜单中的图像属性(figure properties)命令。
24、使用Rotate 3D可以对图像进行各个角度的旋转，编辑线型线宽等。Ys、Ypx、Ydz2、Y出一y2、Y丘留妒)、Yfx3、yd2咿、Ysp3函数可视化程序运行结果如下图l(a卜图1(0。一I，?蓠：“ 鬟 7、：，雾爹淤 ，f墨，毫蟛7 、一一。，74讨论Refbren：l Nal硼mura SN咖懈ical A加蚴|is and Graphic Vbmalization withMATLAB2ndo凼廿蚰Pre嘶ceHall20022 Li Q，动ao WDatai雌and A卵IicaliofMATLABBe嵇h喾NatimmlDe6目lIndmm柚Pr，20013啦ZY，Liu R
25、z，h略z YMasm-,PmfieiemintheMATLABB咧i雌Unvemity of删cs And Aslae-namiesP嘣。19974 BurnsRSAdvancedControl功lgi嘶BeUerw,mh-Heaema啦200I5 DanielABcardHon8 Qin，ChemcalBiophysics：QuantitativeAnalysis ofCellularSystemsCambridge：Q蛐UniversityPI此0086 Huang H JPIm峨Chemical Engineer DiBi协l Simulation万方数据1342 计算机与应用化学 2
26、01 1，28(10)Technology-Applications of MATLAB in ChemicalEngineer Beijiag：Chemical lindustry Press，20047 Xu G YApplications of MATLAB in ChemistryXian：XianJiaotong University Press，20058 Millican D W and McGown L BFluorescence lifetime resolutionof spectra in the f搠uency domain using multiwayanalysis
27、Analytical Chemistry,1990，62：2242-22479 Schueler PA，Ives J T and DeLaCroix F，et a1Physical structure，Optical resonance，and surface enhanced Raman scattering ofsilver-island films on suspended polymer latexparticlesAnalytical Chemistry，1993，65：3 I 77-3 1 8610 Min Gui，Rutan S C and Armel AgbodjanKinet
28、ic detection ofoverlapped amino acids in thin-layer chromatography with a directtrilinear decomposition methodAnalytical Chemistry，1995，67：3293-3296“Bcata Walczak，Bus Van Den Bogaert，and Desire LucMassartApplication of wavelet packet transform in patternrecognition of near-IR damAnalytical Chemistry
29、，1 996，68：1742174712 Shen Hu，Le Zhang，Sergey Krylov，and Norman J DovichiCellcycle-dependent protein fingerprint from a single cancer cell：image cytometry coupled with singlecell capillary sievingelectrophoresisAnalytical Chemistry，2003，75：3495350113 Lyshevski S EEngineering and Scientific Computatio
30、ns UsingMATLABJohn WileySons，Inc，2004：9913214 Wilson H B，Tureotte L H and Halpern DAdvanced Mathenmticsand Mechanics Applications Using MatlabCRC Press，2002：15-6515 Zhang J H，Su Y Z，Song J H，Guo H YTwo dimensional graphsof the wave functions in Origin 60Computers and AppliedChemistry，2003，20(4)：467-
31、47916 Zhang J H，Su Y ZSong J H，Cmo H YThree dimensional graphsof the wave functions in 0一gin 60Computers and AppliedChemistry2003，20(5)：627-63017 Liu P，Zhang D SDrawing atomic orbital and electronic cloudusing MATLABComputers and Applied Chemistry，2003，20(4)：533-53618 Jing S S Elementary Quantum Che
32、mistry and itsApplicationsShanghai：Shanghai Science and TechnologyPublishing House，1980中文参考文献2 李强，赵伟MATLAB数据处理与应用北京：国防工业出版社。20013 张志涌，刘瑞桢，杨祖樱掌握和精通MATLAB北京：北京航空航天大学出版社，19976 黄华江实用化工计算模拟MATLAB在化学工程中的应用北京：化学工业出版社，20047 许国根MATLAB在化学中的应用西安：西安交通大学出版社，2005金松寿量子化学基础及其应用上海：上海科学技术出版社，1980Developing the Mfile
33、fibrary for visualization of spherical harmonic function basedon M【ATLABGuo Shiheng，Xu Qingling，Wang Dongyao，Su Yuzhi and Zhang Jianhua。(School ofChemistry and Chemical Engineering，Guangzhou University,Guangzhon，510006，Guangdong，China)Abstract：A accurate。fast，simple and direct m-file library for vis
34、ualization of atomic orbital spherical harmonic function Ylnl(0，p)such ass,p，d，f and sp3，d2sp3 hybrid orbimls are designedThe matrixT,P】from the discrete data of variables似币【0，2霄】)and0(0【O，兀】)in sphericalcoordinate system are used to calculate the R based on the formula expressions ofatomic orbital
35、sphedcal harmonic function Yl。(e，币)and spj，d2sp3 hybrid orbitals111e matrixX，Y,Z】for plotting are obtained from the discrete data TP,R using【X，距】=sph2cart(T,P,R)commandThusthe visualization of the atomic orbital spherical harmonic function YIm(0，p)and the hybrid orbitals are realized using the mesh(
36、X，YZ)commandProcedures of the visualization for the spherical harmonic functions ale discussedFinally 山e MATLAB readable formulaexpressions of atomic orbital spherical harmonic function Yh(0，币)and the sp，,fsphybrid orblitals based on队TLAB array elementalgorithm,the concrete m-files and graphs are gi
37、ven hereKeywords：MAnABsphericaI harmonic function,m-file(Received：201 1-03-24；Revised：2011-07-02)万方数据基于MATLAB的球谐函数M函数库作者： 郭仕恒， 许庆陵， 王东耀， 苏育志， 张建华， Guo Shiheng， Xu Qingling， Wang Dongyao， SuYuzhi， Zhang Jianhua作者单位： 广州大学化学化工学院,广东,广州,510006刊名： 计算机与应用化学英文刊名： Computers and Applied Chemistry年，卷(期)： 2011,28(10)本文链接：http:/
