径向分布函数求解配位数？

Aristotler

各位老师，同学，想向大家求教一个问题：我想通过径向分布函数求解配位数，目前看到两种计算方法：第一种是利用图1的方法，那么如何求解这个积分函数呢.r^2需要积分吗？（数学渣渣）以上公式单位如何转换？


第二种是如下：

哪种是对的呢？第一种是怎么积分呢？感谢大家

wbn

...你说的第二种方法用的也是第一种的公式...

我以前写过一个f90 code通过gmx的rdf输出算coordination num，你可以拿去改改。不保证对啊

1. ! the program to calculate coordination number from g(r)
   
2.       program gr2coln
   
3.       implicit none
   
4.       type rdf
   
5.       real*8 r,gr
   
6.       endtype
   
7.       character rdffile*50,head
   
8.       real*8 peakhead,peaktail,total,boxlength,nip
   
9.       integer i
   
10.       type(rdf) g_r(2000)
    
11.       call getarg(1,rdffile)
    
12.       write(*,*) "input where the peak begins in nm"
    
13.       read(*,*) peakhead
    
14.       write(*,*) "input where the peak ends in nm"
    
15.       read(*,*) peaktail
    
16.       write(*,*) "the box length in nm"
    
17.       read(*,*) boxlength
    
18.       write(*,*) "number of ion pairs in box"
    
19.       read(*,*) nip
    
20.       open(10,file=rdffile,status='old')
    
21.       read(10,*) head
    
22.       do while(head=='@'.or.head=='#')
    
23.           read(10,*) head
    
24.       enddo
    
25.       i=1
    
26.       total=0.0
    
27.       read(10,*) g_r(i)%r,g_r(i)%gr
    
28.       do while(g_r(i)%r<peakhead)
    
29.           i=i+1
    
30.           read(10,*) g_r(i)%r,g_r(i)%gr
    
31.       enddo
    
32.       do while(g_r(i)%r<peaktail)
    
33.           i=i+1
    
34.           read(10,*) g_r(i)%r,g_r(i)%gr
    
35.           total=total+g_r(i)%gr*4*3.14159*g_r(i)%r*g_r(i)%r*(0.002)
    
36.       enddo
    
37.       total=total*nip/(boxlength*boxlength*boxlength)
    
38.       write(*,*) "The coordination number is : ",total
    
39.       close(10)
    
40.       end
    

复制代码

sobereva

稍微学点编程，用梯形积分很容易就能做
如果不编程，利用origin的曲线积分功能也完全可以

Aristotler

sobereva 发表于 2018-5-15 06:48
稍微学点编程，用梯形积分很容易就能做
如果不编程，利用origin的曲线积分功能也完全可以

谢谢老师，想向您确认下，第二种方法对吗？我目前只会第二种的

Aristotler

wbn 发表于 2018-5-15 04:26
...你说的第二种方法用的也是第一种的公式...

我以前写过一个f90 code通过gmx的rdf输出算coordination n ...

非常感谢你，想问你第二种方法对吗？

fhh2626

第二种方法也是用第一种方法的公式积分。。。

wbn

Aristotler 发表于 2018-5-15 07:53
非常感谢你，想问你第二种方法对吗？

根本就没有两种方法，只有这一种积分公式，你要在积分里把r^2去掉了那肯定不对啊

PS：我的code里积分直接用的累加，你需要比较精确的结果的话需要梯形积分

Aristotler

sobereva 发表于 2018-5-15 06:48
稍微学点编程，用梯形积分很容易就能做
如果不编程，利用origin的曲线积分功能也完全可以

老师，今天又弄了一天怎么求解配位数。还是有些问题，想求问您：我只要将g(r)用origin积分,取计算第一个峰谷的g(r)积分值就是第一壳层的配位数吗？但是配位数与径向分布函数的关系式里面又包括对r*r的积分，那么以上做法对吗？

sobereva

Aristotler 发表于 2018-5-16 22:53
老师，今天又弄了一天怎么求解配位数。还是有些问题，想求问您：我只要将g(r)用origin积分,取计算第一个 ...

你在origin里增加一列，把g(r)乘上r^2转化成被积函数之后再用origin积分，然后再乘上前面的系数即可

Aristotler

sobereva 发表于 2018-5-17 09:23
你在origin里增加一列，把g(r)乘上r^2转化成被积函数之后再用origin积分，然后再乘上前面的系数即可

非常感谢老师，非常

gpp201013

你好，我目前也在计算径向分布函数，我想知道具体算法，求前辈不吝赐教，拜托啦

Jack

楼主能分享下含有上面公式的文献么

Aristotler

gpp201013 发表于 2019-7-4 21:00
你好，我目前也在计算径向分布函数，我想知道具体算法，求前辈不吝赐教，拜托啦

我做的就是vasp第一性原理的动力学，先用MS的AMORPHOUS CELL 模块建模，之后用vasp进行计算，之后用VMD导出给g（r）之后就是用sob老师上面提到的方法，计算配位数。希望对你有帮助

Aristotler

Jack 发表于 2019-7-5 07:03
楼主能分享下含有上面公式的文献么

文献是在小*虫上截的

highlight

Aristotler 发表于 2019-7-5 08:52
我做的就是vasp第一性原理的动力学，先用MS的AMORPHOUS CELL 模块建模，之后用vasp进行计算，之后用VMD导 ...

有个问题啊，你既然用的是vmd，为啥要自己算积分呢？
http://www.ks.uiuc.edu/Research/vmd/plugins/gofrgui/

The number integrals are computed directly and thus provide accurate coordination numbers.