Yu-Xuan's Blog

在平时的计算过中,通常会遇到多层的循环嵌套求和,这是一个非常消耗时间的事情,自己现在除了在大型矩阵对角化的时候还是利用Fortran在写程序,平时一些涉及到矢量运算或者动量空间中的一些计算都是在利用Julia做的,那么将这个很耗时的事情变成多线程自然是可以节省很多时间,所以这里就折腾一下如何把一个求和分发到多个线程,然后将所有的结果在放到一起,实现多线程并行,可以节省很多的时间.{:.info} 简单测试1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071using ProgressMeterusing Distributed@everywhere using SharedArrays, LinearAlgebra# ----------------------------------------------------------------@everywhere function fun1(qx::Float64,qy::Float64)::Float64 kn = 15000 re::Float64 = 0.0 for i1 in 0:kn for i2 in 0:kn kx = i1*pi/kn ky = i2*pi/kn re = re + kx + ky + qy + qx end end return re/kn^2end#-----------------------------------------function fun2(kn::Int64)::Float64 re::Float64 = 0.0 for i1 in 1:kn for i2 in 1:kn qx = i1*pi/kn qy = i2*pi/kn re = re + fun1(qx,qy) end end return...

因为对格林函数很熟悉,在看文章的时候也有遇到计算体系准粒子干涉的文章,最后只明白可以通过干涉的分布样式来判断不同费米面之间的散射问题,从而来确定体系费米面上的一些性质,还有就是这些干涉的图案还是很好看的,所以一直也就想自己动手去算一下,这里我就想从一个很简单的例子来学习一些如何编程计算准粒子干涉的图样,同时也算是对文章的进一步理解.{:.info}这里主要是想重复Quasiparticle scattering interference in superconducting iron pnictides这篇文章中的主要结果,因为自己就是搞超导方向的(好吧,我承认我对超导其实理解并没有很深入),所以就利用这篇铁基超导的准粒子干涉计算来练练手.{:.warning} Fortran123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181 module pub implicit none integer N,ne,kn parameter(N = 4,ne = 100,kn = 200) real,parameter::pi = 3.14159265358 complex,parameter::im = (0.0,1.0) complex ham(N,N) complex gf0(N,N) ! bare...

虽然之前也整理了如何计算Chern数和$Z_2$拓扑不变量,但是对于最简单的Chern Insulator却没有认真的研究过,最近在做一些和反常量子霍尔相关的一些内容,就正好把这个最简单的Chern绝缘体模型的边界态以及Chern数计算的结果整理到一起.{:.info} 边界态计算Chern Insulator的哈密顿量非常简单 \begin{equation} H(\mathbf{k})=(m_0-t_x\cos(k_x)+t_y\cos(k_y))\sigma_z+\lambda_x\sin(k_x)+\lambda_y\sin(k_y) \end{equation}在这里写成$(m_0-t_x\cos(k_x)+t_y\cos(k_y))$的形式是为了破坏旋转对称性,这样的画采用cylinder geometry(一个方向开边界,一个方向周期)画能带图的时候,就会发现$x$方向和$y$方向的边界态一个出现的$k_{x/y}=0$,另外一个方向的边界态出现在$k_{y/x}=\pi$的位置上,如下图所示 计算代码如下,之前算这东西吸怪用Fortran了,所以就一直沿用原来的习惯了,不好改过来123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192! Author: YuXuanLi! ...

虽然一直在计算和学习凝聚态中的拓扑,不过感觉对Berry Phase和Chern Number的一些认识还不是很到位,这里正好自己最近开组会讨论到这个问题,就将自己整理的PPT分享出来,具体的理解差不多在PPT内容中也写了出来.{:.info} Main 至于PPT的源代码,可以在这里下载 公众号相关内容均会在公众号进行同步，若对该Blog感兴趣，欢迎关注微信公众号。{:.info} Email yxliphy@gmail.com

之前在Julia,Python,Fortran,Mathematica循环计算速度比较博客中，我简单的对集中编程语言的循环进行了对比，虽然没什么太大的使用价值，不过对我对自己写代码时候的速度考虑开了一个不错的头，所以这里就把自己利用循环加速改写了一个实例展示出来，看看到底效率如何。{:.info} 循环加速并行矩阵运算在前面的博客中，我只是对多重循环中进行了简单的求和，但是在实际的运算中可定不会那么的简单，所以我首先想测试的是如果循环里面是矩阵的运算，速度到底会不会提升的很好。12345678910111213141516171819202122232425262728293031323334from numba import jitimport timeimport numpy as np@jit # 函数闭包#@jit(nopython=True, parallel=True)def f1(): c = 0 cont = 100000 for i in range(cont): for j in range(cont): for k in range(cont): c = c + i + j + k return c# @jit#@jit(nopython=True, parallel=True)@jit(parallel=True)def f2(): ndim = 100 cont = 1000 mat1 = np.random.rand(ndim,ndim) matre = np.zeros((ndim,ndim)) for i in range(cont): for j in range(cont): for k in range(cont): matre = np.dot(mat1,mat1) return matret1 = time.time()# print(f1()) # 计算循环求和print(f2())t2 = time.time()print('Timing cost is(s): ',t2 -...

在之前的博客中Python稀疏矩阵对角化库,我虽然介绍了一个python的库,可以用来计算稀疏矩阵的一些特定要求的本征矢量和本征值,但是python通常情况下是比较慢的,幸运的是这个功能最近正好的Julia中科成功的实现了,而且用法和python中几乎相同,这里正好就拿Julia来实验一下.{:.info} 前言对于julia同样是需要安装一个库才可以Arpack12import PkgPkg.add("Arpack")安装完成之后,就可以使用了 实例(Julia)123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114using Arpack,LinearAlgebra,DelimitedFiles# ====================================function matset(xn::Int64,ky::Float64) m0 = 1.5 # Effect mass mu = 0 # Chemical potential tx = 1.0 ty = 1.0 ax = 1.0 ay = 1.0 del0 = 0 delx = 0. dely = -0. h = 0 N = xn*8 ham = zeros(ComplexF64,N,N) #------------------------------------------------------------------------------- g1 = zeros(ComplexF64,8,8) g2 = zeros(ComplexF64,8,8) g3 = zeros(ComplexF64,8,8) g4 =...