24小时热门版块排行榜    

CyRhmU.jpeg
查看: 1046  |  回复: 7

xzz1007

金虫 (小有名气)

[求助] 想请教一个有关二面角的问题!已有1人参与

我现在使用tinker程序中的xyzint程序,将笛卡尔坐标(xyz coor.)转化为内坐标(internal coor.)的形式;
然后再使用其中的intxyz程序将内坐标转化回笛卡尔坐标形式,但是却得不到原来的构型!

关键的问题出现在转化二面角上,比如:二面角ABCD

程序中是先求出angle(ABC)所在平面的法向量;然后再求出angle(BCD)所在平面的法向量,最后求出两个平面的法向量,最终求出这两个平面的夹角!

我的问题是:
1. 当这两个平面重合的情况下,是应该没有二面角,还是0度或者180度?
        因为这种情况下,程序中是处理成,二面角直接等于angle(ABD),这个平面角了,这样的处理不太理解!

2. 如果ABC三点或者BCD三点在同一条直线上,这样的二面角又该怎么处理?

3. 同样,如果两个平面重合的话,它们各自对应的法向量应该是什么关系呢?

我是学化学的,二面角以及向量这块儿,很是欠缺,求大侠帮助!
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

阳光总在风雨后,再坚持坚持吧.....
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

feixiaolin

荣誉版主 (文坛精英)

优秀版主

http://wenku.baidu.com/link?url= ... RZZOriH9pAYI98cTaCu

另外,解反三角函数需要分象限讨论。
2楼2014-09-26 15:45:37
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

polypro

木虫 (正式写手)

【答案】应助回帖

感谢参与,应助指数 +1
以前做个这个问题,有现成的fortran代码(子程)对笛卡尔和极坐标转换,需要告诉一声
泉涸,鱼相与处于陆,相呴以湿,相濡以沫,不如相忘于江湖。
3楼2014-09-28 08:28:25
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xzz1007

金虫 (小有名气)

sweety: 鼓励站内交流 2014-09-28 14:21:12
引用回帖:
3楼: Originally posted by polypro at 2014-09-28 08:28:25
以前做个这个问题,有现成的fortran代码(子程)对笛卡尔和极坐标转换,需要告诉一声

可以发给我一份吗?我想对比一下!
我总觉得tinker中的程序在转化的时候有点问题!
留个邮箱吧!
阳光总在风雨后,再坚持坚持吧.....
4楼2014-09-28 11:16:12
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xzz1007

金虫 (小有名气)

引用回帖:
4楼: Originally posted by xzz1007 at 2014-09-28 11:16:12
可以发给我一份吗?我想对比一下!
我总觉得tinker中的程序在转化的时候有点问题!
留个邮箱吧!
jane_xu@dicp.ac.cn...

你说的是笛卡尔和极坐标的转化啊?我要的是笛卡尔和内坐标的转换呀!
阳光总在风雨后,再坚持坚持吧.....
5楼2014-09-28 11:18:00
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xzz1007

金虫 (小有名气)

引用回帖:
2楼: Originally posted by feixiaolin at 2014-09-26 15:45:37
http://wenku.baidu.com/link?url=Ab7fWcvAZL6E8BYkBNnwTYZqonJWjIwrlo1ETEV9brhtg0s3BpDG-XSZRa5mAmPOCkGtGJQ-ywX8rfWi0om1IuRgRZZOriH9pAYI98cTaCu

另外,解反三角函数需要分象限讨论。

非常感谢你给我的链接!可是,我只是想知道遇到特殊角的时候应该怎么处理?:
阳光总在风雨后,再坚持坚持吧.....
6楼2014-09-28 11:21:19
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

polypro

木虫 (正式写手)

引用回帖:
4楼: Originally posted by xzz1007 at 2014-09-28 11:16:12
可以发给我一份吗?我想对比一下!
我总觉得tinker中的程序在转化的时候有点问题!
留个邮箱吧!...

极坐标转换
!this sub to transfer Cartesian coordinates (x,y,z) with Polar coordinates (rou,theta,phi) (radius)
!control number 1 for (x,y,z)=>(rou,theta,phi) and -1 is reverse
!x=rou*sin(phi)*cos(theta)
!x=rou*sin(phi)*sin(theta)
!z=rou*cos(phi)
        subroutine trans(a,b,c,iflag)
        parameter(pi=3.1415926,h=0.01745329252)  !pi/180
        real rou,theta,phi
        real a,b,c
        integer iflag !1 from Cartesan to Polar, -1 is reverse (degree)
!print*,"a b cini",a,b,c,iflag
    if(iflag==1)then !Cartesan to Polar
          rou=sqrt(a*a+b*b+c*c)
          if(rou.eq.0)then
          a=0
          b=0
          c=0
          goto 80
          endif
          phi=acos(c/rou)
          if(abs(rou-c).le.1e-3.or.abs(b).le.1e-3)then
             theta=0  !point on Z aixs or Y axis, thets==0 fixed!
          else
             if(abs(a).le.1e-3)then
                theta=pi/2.0
             else
                sinphi=sqrt(1-c*c/rou/rou)
                theta=acos(a/rou/sinphi)
             endif
          endif
          if(b.lt.0)theta=-theta+pi   !all is in degree
          a=rou
          b=theta/h
          c=phi/h !in radius
    else  !polar to Cartesan
          rou=a
          if(rou.eq.0)then
          a=0
          b=0
          c=0          
          goto 80
          endif
          theta=b*h
          phi=c*h
          c=rou*cos(phi) !z
          b=rou*sin(phi)*sin(theta) !y
          a=rou*sin(phi)*cos(theta) !x         
        endif
80   continue
!print*,"a b c",a,b,c,iflag
    END subroutine trans
泉涸,鱼相与处于陆,相呴以湿,相濡以沫,不如相忘于江湖。
7楼2014-09-28 15:07:01
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

polypro

木虫 (正式写手)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
xzz1007: 金币+20, ★★★很有帮助 2014-09-28 16:20:28
另一个是计算二面角的,不需要极坐标转换(慢)
你说二个平面重合时,二面角等于0;也可以设为pi
!!this sub to calculate torsion angle
!!calculate torsion angle  REF:wiki dihedral angle
!!      ¦Ó (a, b, c) = arg (¡Ýa N7 c + (a N7 b)(b N7 c), a N7 (b NW c)) .
!!¦È = arg(x, y)
!!cos ¦È = x/ (x2 + y 2 ) and sin ¦È = y/ (x2 + y 2 ).
  function torsion_angle(tia1,tia2,tia3,tia4)   !4个atom series
    implicit none
    real torsion_angle
    real ax,ay,az,bx,by,bz,cx,cy,cz,r1,r2,r3
    real tab,tac,tbc,t1,ctbcx,ctbcy,ctbcz,ctbax,ctbay,ctbaz,t2
    integer tia1,tia2,tia3,tia4
    ax=xx(tia2)-xx(tia1)
    ay=yy(tia2)-yy(tia1)
    az=zz(tia2)-zz(tia1)
    bx=xx(tia3)-xx(tia2)
    by=yy(tia3)-yy(tia2)
    bz=zz(tia3)-zz(tia2)
    cx=xx(tia4)-xx(tia3)
    cy=yy(tia4)-yy(tia3)
    cz=zz(tia4)-zz(tia3)
    r1=sqrt(ax*ax+ay*ay+az*az)
    r2=sqrt(bx*bx+by*by+bz*bz)
    r3=sqrt(cx*cx+cy*cy+cz*cz)
    if(r1*r2*r3<0.001)then
       torsion_angle=0.0
       goto 33
    endif
    ax=ax/r1
    ay=ay/r1
    az=az/r1
    bx=bx/r2
    by=by/r2
    bz=bz/r2
    cx=cx/r3
    cy=cy/r3
    cz=cz/r3
    tab=ax*bx+ay*by+az*bz  !dot product
    tac=ax*cx+ay*cy+az*cz
    tbc=bx*cx+by*cy+bz*cz
    t1=tab*tbc-tac
    ctbcx=by*cz-bz*cy      !cross product  bxc
    ctbcy=bz*cx-bx*cz
    ctbcz=bx*cy-by*cx
    ctbax=by*az-bz*ay      !cross product bxa
    ctbay=bz*ax-bx*az
    ctbaz=bx*ay-by*ax
    t2=ax*ctbcx+ay*ctbcy+az*ctbcz  !last term
       torsion_angle=acos(t1/sqrt(t1*t1+t2*t2))   !!0~pi
!!judge the r2xr3 angle with r1, if acute (-180,0),else(-180,0)
    if(t2<0)torsion_angle=-1.0*torsion_angle    !-pi~0
33       return
  end function torsion_angle
泉涸,鱼相与处于陆,相呴以湿,相濡以沫,不如相忘于江湖。
8楼2014-09-28 15:10:25
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 xzz1007 的主题更新
信息提示
请填处理意见