TABLE OF CONTENTS
ABINIT/bldgrp [ Functions ]
FUNCTION
Yields all the symmetry operations starting from the generators. Applies all the generators onto themselves, and obtains all the other operations. Iterates until it reaches nsym.
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (RC) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
msym = default number of symmetry operations nsym = number of symmetry operations symafm(msym)=(anti)ferromagnetic part of symmetry operations symrel(3,3,msym) = 3D matrix containg symmetry operations tnons(3,msym) = 2D matrix containing translations of the symmery operations
OUTPUT
symafm(msym)=(anti)ferromagnetic part of symmetry operations symrel(3,3,msym) = 3D matrix containg symmetry operations tnons(3,msym) = 2D matrix containing translations of the symmery operations
SIDE EFFECTS
nogen = number of generators, number of operations to be applied onto themselves
PARENTS
symsgcube,symsghexa,symsgortho,symsgtetra
CHILDREN
SOURCE
42 #if defined HAVE_CONFIG_H 43 #include "config.h" 44 #endif 45 46 #include "abi_common.h" 47 48 49 subroutine bldgrp(msym,nogen,nsym,symafm,symrel,tnons) 50 51 use defs_basis 52 use m_errors 53 use m_profiling_abi 54 55 !This section has been created automatically by the script Abilint (TD). 56 !Do not modify the following lines by hand. 57 #undef ABI_FUNC 58 #define ABI_FUNC 'bldgrp' 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 integer,intent(in) :: msym,nsym 66 integer,intent(inout) :: nogen 67 !arrays 68 integer,intent(inout) :: symafm(msym),symrel(3,3,msym) 69 real(dp),intent(inout) :: tnons(3,msym) 70 71 !Local variables ------------------------------ 72 !matrintoper(3,3) & matrinttransl(3) are intermediate arrays of the new 73 ! symmetry operations obtained, in order to check their uniqueness. 74 !flagop,flagtr = flags used during the checking of the similarity between 75 ! the obtained operation and the already existent ones 76 !ii,ijk,ijkl,jjj,kk = counters in the cycles 77 !scalars 78 integer :: flagop,flagtr,ii,ijk,ijkl,jjj,kk,matrintsymafm,nogen_new 79 real(dp) :: nastyzero 80 character(len=500) :: message 81 !arrays 82 integer :: bcksymafm(2*msym),bcksymrel(3,3,2*msym),matrintoper(3,3) 83 real(dp) :: bcktnons(3,2*msym),matrinttransl(3) 84 85 ! ************************************************************************* 86 87 nastyzero=0.1 88 89 !DEBUG 90 ! write(std_out,*)' bldgrp : enter, builds the space group symmetry ' 91 ! write(std_out,*)' bldgrp : number of generators : ',nogen 92 ! write(std_out,*)' bldgrp : nsym,msym=',nsym,msym 93 !ENDDEBUG 94 95 if (nogen<1) then 96 write(message, '(a,i4,a,a,a,a,a)' )& 97 & 'The number of generators nogen is ',nogen,& 98 & 'and it should be greater than one',ch10,& 99 & 'This is not allowed. ',ch10,& 100 & 'Action: Contact ABINIT group ' 101 MSG_ERROR(message) 102 end if 103 104 !Transfer the generators to bcksymrel 105 do ii=1,nogen 106 bcksymrel(:,:,ii)=symrel(:,:,ii) 107 bcktnons(:,ii)=tnons(:,ii) 108 bcksymafm(ii)=symafm(ii) 109 end do 110 111 !DEBUG 112 write(std_out,*)' Describe the different generators (index,symrel,tnons,symafm)' 113 do ii=1,nogen 114 write(std_out,'(i3,2x,9i3,3es12.2,i3)')ii,symrel(:,:,ii),tnons(:,ii),symafm(ii) 115 end do 116 !ENDDEBUG 117 118 !Simply iterate until the group is complete 119 do ijkl=1,nsym 120 121 ! DEBUG 122 ! write(std_out,*)' bldgrp : in loop, ijkl,nogen=',ijkl,nogen 123 ! ENDDEBUG 124 125 nogen_new=nogen 126 127 do jjj=2,nogen 128 do kk=2,nogen 129 130 ! Computing block of the new symmetry operation according to: 131 ! ! $ { R1 | v1 }{ R2 | v2 } = { R1.R2 | v1+R1.v2 } $ 132 matrintoper(:,:) = matmul(bcksymrel(:,:,jjj),bcksymrel(:,:,kk)) 133 matrinttransl(:) = bcktnons(:,jjj)+matmul(bcksymrel(:,:,jjj),bcktnons(:,kk)) 134 matrintsymafm = bcksymafm(jjj)*bcksymafm(kk) 135 136 ! Rescaling translation between 0 and 1 137 do ii=1,3 138 if (matrinttransl(ii)>=0.9) then 139 do while (matrinttransl(ii)>=0.9) 140 matrinttransl(ii)=matrinttransl(ii)-1.0 141 end do 142 end if 143 if (matrinttransl(ii)<0.0) then 144 do while (matrinttransl(ii)<0.0) 145 matrinttransl(ii)=matrinttransl(ii)+1.0 146 end do 147 end if 148 if ( abs(matrinttransl(ii))<nastyzero) matrinttransl(ii)=0.0 149 if ( abs(matrinttransl(ii)-1.0)<nastyzero) matrinttransl(ii)=0.0 150 end do 151 152 ! Cheking block to validate the new symmetry operation 153 do ijk=1,nogen_new 154 155 flagop=0 ; flagtr=0 156 157 ! Check for rotation similarity 158 if(sum((matrintoper-bcksymrel(:,:,ijk))**2)==0)flagop=1 159 160 ! Check for translation similarity 161 if(maxval((matrinttransl-bcktnons(:,ijk))**2)<nastyzero**2)flagtr=1 162 163 if(flagop+flagtr==2)exit 164 165 end do 166 167 ! Add the new determined symmetry if it is unique 168 if (flagtr+flagop<2) then 169 nogen_new=nogen_new+1 170 ! DEBUG 171 ! write(std_out,*)' added one more symmetry : nogen_new=',nogen_new 172 ! ENDDEBUG 173 bcksymrel(:,:,nogen_new)=matrintoper(:,:) 174 bcktnons(:,nogen_new)=matrinttransl(:) 175 bcksymafm(nogen_new)=matrintsymafm 176 end if 177 178 end do 179 end do 180 181 nogen=nogen_new 182 183 if(nogen==nsym)exit 184 185 end do 186 187 !Transfer of the calculated symmetry to the routine output 188 if (nogen==nsym) then 189 symrel(:,:,1:nsym)=bcksymrel(:,:,1:nsym) 190 tnons(:,1:nsym)=bcktnons(:,1:nsym) 191 symafm(1:nsym)=bcksymafm(1:nsym) 192 else 193 ! Problem with the generation of the symmetry operations 194 write(message, '(a,i7,a,a,i7)' )& 195 & 'The symmetries obtained are ',nogen,ch10,& 196 & 'and they should be ',nsym 197 MSG_BUG(message) 198 end if 199 200 !DEBUG 201 !write(std_out,*)' bldgrp : exit with ',nogen,' operation symmetries' 202 !ENDDEBUG 203 204 end subroutine bldgrp