scnc
/
quantum-espresso
mirror of https://gitlab.com/QEF/q-e.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Bug fix: yet another problem with D_2h. 

git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@5358 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
dalcorso 2009-01-30 09:14:51 +00:00
parent 082fbb1bd5
commit 7fd3fa7ea7
2 changed files with 60 additions and 47 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

66
PW/divide_class.f90
View File

@ -48,8 +48,8 @@ REAL(DP) :: smat(3,3,nrot), cmat(3,3), ax(3), ars
INTEGER :: done(48), irot, jrot, krot, iclass, i
INTEGER :: tipo_sym, ipol, axis, axis1, axis2, ts
REAL(DP), PARAMETER :: eps = 1.d-7
REAL(DP) :: angle_rot, angle_rot_s
LOGICAL :: compare_mat, is_axis, first, first1, done_ax(6)
REAL(DP) :: angle_rot, angle_rot_s, ax_save(3,2:4)
LOGICAL :: compare_mat, is_axis, is_parallel, first, first1, done_ax(6)
!
! Divide the group in classes.
!
@ -434,22 +434,11 @@ ELSEIF (code_group==20) THEN
which_irr(iclass)=4
done_ax(3)=.FALSE.
END IF
ax_save(:,which_irr(iclass))=ax(:)
ELSEIF (ts==2) THEN
which_irr(iclass)=5
ELSEIF (ts==5) THEN
CALL mirror_axis(smat(1,1,elem(1,iclass)),ax)
IF (is_axis(ax,3)) THEN
which_irr(iclass)=6
done_ax(4)=.FALSE.
ELSE IF (is_axis(ax,2)) THEN
which_irr(iclass)=7
done_ax(5)=.FALSE.
ELSE IF (is_axis(ax,1)) THEN
which_irr(iclass)=8
done_ax(6)=.FALSE.
END IF
END IF
END DO
ENDIF
ENDDO
!
! Otherwise choose the first free axis
!
@ -465,20 +454,27 @@ ELSEIF (code_group==20) THEN
END IF
END DO
100 CONTINUE
ELSEIF (ts==5) THEN
DO i=4,6
IF (done_ax(i)) THEN
which_irr(iclass)=i+2
done_ax(i)=.FALSE.
GOTO 120
END IF
END DO
120 CONTINUE
CALL versor(smat(1,1,elem(1,iclass)),ax)
ax_save(:,which_irr(iclass))=ax(:)
ENDIF
ENDIF
ENDDO
!
! Finally it orders the mirror planes. The perpendicular to the plane
! must be parallel to one of the C_2 axis.
!
!
DO iclass=2,nclass
ts=tipo_sym(smat(1,1,elem(1,iclass)))
IF (ts==5) THEN
CALL mirror_axis(smat(1,1,elem(1,iclass)),ax)
DO i=2,4
IF (is_parallel(ax,ax_save(1,i))) which_irr(iclass)=i+4
ENDDO
END IF
IF (which_irr(iclass)==0) CALL errore('divide_class',&
'something wrong D_2h',1)
ENDDO
END DO
ELSEIF (code_group==21) THEN
!
! D_3h
@ -2633,3 +2629,21 @@ is_complex= complex_aux(code)
RETURN
END FUNCTION is_complex
FUNCTION is_parallel(a,b)
!
! This function returns true if a(3) and b(3) are parallel vectors
!
USE kinds, ONLY : DP
IMPLICIT none
LOGICAL :: is_parallel
REAL(DP) :: a(3), b(3)
REAL(DP) :: cross
cross=(a(2)*b(3)-a(3)*b(2))**2+(a(3)*b(1)-a(1)*b(3))**2+(a(1)*b(2)-a(2)*b(1))**2
is_parallel=(ABS(cross)< 1.d-6)
RETURN
END FUNCTION is_parallel

41
PW/divide_class_so.f90
View File

@ -55,8 +55,8 @@ COMPLEX(DP) :: d_spin(2,2,48), d_spine(2,2,96), c_spin(2,2)
INTEGER :: done(96), irot, jrot, krot, iclass, i
INTEGER :: tipo_sym, set_e, ipol, axis, axis1, axis2, ts
REAL(DP), PARAMETER :: eps = 1.d-7
REAL(DP) :: angle_rot, angle_rot_s, ars
LOGICAL :: compare_mat_so, is_axis, first, first1
REAL(DP) :: angle_rot, angle_rot_s, ars, ax_save(3,3:5)
LOGICAL :: compare_mat_so, is_axis, first, first1, is_parallel
LOGICAL :: done_ax(6)
!
! Divide the group in classes.
@ -514,21 +514,13 @@ ELSEIF (code_group==20) THEN
which_irr(iclass)=5
done_ax(3)=.FALSE.
END IF
ax_save(:,which_irr(iclass))=ax(:)
ELSEIF (ts==2) THEN
IF (has_e(1,iclass)==-1) THEN
which_irr(iclass)=7
ELSE
which_irr(iclass)=6
END IF
ELSEIF (ts==5) THEN
CALL mirror_axis(smat(1,1,elem(1,iclass)),ax)
IF (is_axis(ax,3)) THEN
which_irr(iclass)=8
ELSE IF (is_axis(ax,2)) THEN
which_irr(iclass)=9
ELSE IF (is_axis(ax,1)) THEN
which_irr(iclass)=10
END IF
END IF
END DO
!
@ -546,17 +538,24 @@ ELSEIF (code_group==20) THEN
END IF
END DO
100 CONTINUE
ELSEIF (ts==5) THEN
DO i=4,6
IF (done_ax(i)) THEN
which_irr(iclass)=i+4
done_ax(i)=.FALSE.
GOTO 110
END IF
CALL versor(smat(1,1,elem(1,iclass)),ax)
ax_save(:,which_irr(iclass))=ax(:)
END IF
END IF
END DO
!
! Finally consider the mirror planes
!
DO iclass=2,nclass
IF (which_irr(iclass)==0) THEN
ts=tipo_sym(smat(1,1,elem(1,iclass)))
IF (ts==5) THEN
CALL mirror_axis(smat(1,1,elem(1,iclass)),ax)
DO i=3,5
IF (is_parallel(ax,ax_save(:,i))) which_irr(iclass)=i+5
END DO
110 CONTINUE
ENDIF
ENDIF
END IF
END IF
IF (which_irr(iclass)==0) CALL errore('divide_class_so',&
'something wrong D_2h',1)
ENDDO