TABLE OF CONTENTS
ABINIT/symanal [ Functions ]
NAME
symanal
FUNCTION
Find the space group, Bravais lattice, including Shubnikov characteristics from the list of symmetries (including magnetic characteristics), and lattice parameters Warning: the recognition of the space group might not yet work for the Shubnikov group of type IV
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (XG, 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
chkprim= if 1 then stop if the cell is not primitive msym=default maximal number of symmetries nsym=actual number of symmetries rprimd(3,3)=dimensional primitive translations for real space (bohr) symafm(1:msym)=(anti)ferromagnetic part of symmetry operations symrel(3,3,1:msym)=symmetry operations in real space in terms of primitive translations tnons(3,1:msym)=nonsymmorphic translations for symmetry operations tolsym=tolerance for the symmetry operations
OUTPUT
bravais(11)=characteristics of Bravais lattice (see symlatt.F90) genafm(3)=magnetic translation generator (in case of Shubnikov group type IV) ptgroupma = magnetic point group number spgroup=symmetry space group
SIDE EFFECTS
NOTES
PARENTS
ingeo,m_ab7_symmetry,m_tdep_sym,m_use_ga
CHILDREN
chkprimit,getptgroupma,symbrav,symlatt,symptgroup,symspgr,wrtout
SOURCE
50 #if defined HAVE_CONFIG_H 51 #include "config.h" 52 #endif 53 54 #include "abi_common.h" 55 56 subroutine symanal(bravais,chkprim,genafm,msym,nsym,ptgroupma,rprimd,spgroup,symafm,symrel,tnons,tolsym) 57 58 use defs_basis 59 use m_profiling_abi 60 61 !This section has been created automatically by the script Abilint (TD). 62 !Do not modify the following lines by hand. 63 #undef ABI_FUNC 64 #define ABI_FUNC 'symanal' 65 use interfaces_14_hidewrite 66 use interfaces_41_geometry, except_this_one => symanal 67 !End of the abilint section 68 69 implicit none 70 71 !Arguments ------------------------------------ 72 !scalars 73 integer,intent(in) :: chkprim,msym,nsym 74 integer,intent(out) :: ptgroupma,spgroup 75 real(dp),intent(in) :: tolsym 76 !arrays 77 integer,intent(out) :: bravais(11) 78 integer,intent(in) :: symafm(msym),symrel(3,3,msym) 79 real(dp),intent(in) :: rprimd(3,3) 80 real(dp),intent(inout) :: tnons(3,msym) 81 real(dp),intent(out) :: genafm(3) 82 83 !Local variables------------------------------- 84 !scalars 85 integer :: iholohedry_nomagn,isym,isym_nomagn,multi 86 integer :: nptsym,nsym_nomagn,shubnikov 87 character(len=5) :: ptgroup,ptgroupha 88 character(len=500) :: message 89 !arrays 90 integer :: identity(3,3) 91 integer,allocatable :: ptsymrel(:,:,:),symrel_nomagn(:,:,:) 92 real(dp),allocatable :: tnons_nomagn(:,:) 93 94 ! ************************************************************************* 95 96 !DEBUG 97 !write(std_out,*)' symanal : enter ' 98 !call flush(6) 99 !stop 100 !ENDDEBUG 101 102 !This routine finds the Bravais characteristics, without actually 103 !looking at the symmetry operations. 104 ABI_ALLOCATE(ptsymrel,(3,3,msym)) 105 call symlatt(bravais,msym,nptsym,ptsymrel,rprimd,tolsym) 106 ABI_DEALLOCATE(ptsymrel) 107 108 !Check whether the cell is primitive or not. 109 call chkprimit(chkprim,multi,nsym,symafm,symrel) 110 111 spgroup=0 ; ptgroupma=0 ; genafm(:)=zero 112 113 if(multi>1)then ! Modify bravais if the cell is not primitive ; no determination of the space group 114 bravais(1)=-bravais(1) 115 else 116 117 ! The cell is primitive, so that the space group can be 118 ! determined. Need to distinguish Fedorov and Shubnikov groups. 119 ! Do not distinguish Shubnikov types I and II. 120 ! Also identify genafm, in case of Shubnikov type IV 121 identity(:,:)=reshape((/1,0,0,0,1,0,0,0,1/),(/3,3/)) 122 shubnikov=1 123 do isym=1,nsym 124 if(symafm(isym)==-1)then 125 shubnikov=3 126 if(sum(abs(symrel(:,:,isym)-identity(:,:)))==0)then 127 shubnikov=4 128 genafm(:)=tnons(:,isym) 129 ! DEBUG 130 ! write(std_out,*)' isym=',isym 131 ! write(std_out,*)' symrel(:,:,isym)',symrel(:,:,isym) 132 ! write(std_out,*)' tnons(:,isym)',tnons(:,isym) 133 ! write(std_out,*)' symafm(isym)',symafm(isym) 134 ! ENDDEBUG 135 exit 136 end if 137 end if 138 end do 139 140 if(shubnikov/=1)then 141 if(shubnikov==3)write(message, '(a)' )' Shubnikov space group type III' 142 if(shubnikov==4)write(message, '(a)' )' Shubnikov space group type IV' 143 call wrtout(std_out,message,'COLL') 144 end if 145 146 if(shubnikov==1 .or. shubnikov==3)then 147 ! Find the correct Bravais characteristics and point group 148 ! Should also be used for Shubnikov groups of type IV ... 149 call symbrav(bravais,msym,nsym,ptgroup,rprimd,symrel,tolsym) 150 ! Find the space group 151 call symspgr(bravais,nsym,spgroup,symrel,tnons,tolsym) 152 end if 153 154 if(shubnikov/=1)then 155 156 ! Determine nonmagnetic symmetry operations 157 nsym_nomagn=nsym/2 158 ABI_ALLOCATE(symrel_nomagn,(3,3,nsym_nomagn)) 159 ABI_ALLOCATE(tnons_nomagn,(3,nsym_nomagn)) 160 isym_nomagn=0 161 do isym=1,nsym 162 if(symafm(isym)==1)then 163 isym_nomagn=isym_nomagn+1 164 symrel_nomagn(:,:,isym_nomagn)=symrel(:,:,isym) 165 tnons_nomagn(:,isym_nomagn)=tnons(:,isym) 166 end if 167 end do 168 169 if(shubnikov==3)then 170 171 ! DEBUG 172 ! write(std_out,*)' symanal : will enter symbrav with halved symmetry set' 173 ! write(std_out,*)' Describe the different symmetry operations (index,symrel,tnons,symafm)' 174 ! do isym=1,nsym_nomagn 175 ! write(std_out,'(i3,2x,9i3,3es12.2,i3)')isym,symrel_nomagn(:,:,isym),tnons_nomagn(:,isym) 176 ! end do 177 ! ENDDEBUG 178 179 ! Find the point group of the halved symmetry set 180 call symptgroup(iholohedry_nomagn,nsym_nomagn,ptgroupha,symrel_nomagn) 181 182 ! Deduce the magnetic point group (ptgroupma) from ptgroup and ptgroupha 183 call getptgroupma(ptgroup,ptgroupha,ptgroupma) 184 185 else if(shubnikov==4)then 186 187 ! Find the Fedorov space group of the halved symmetry set 188 call symspgr(bravais,nsym_nomagn,spgroup,symrel_nomagn,tnons_nomagn,tolsym) 189 190 ! The magnetic translation generator genafm has already been determined 191 192 ! DEBUG 193 ! write(std_out,*)' genafm =',genafm 194 ! write(std_out,*)' spgroup=',spgroup 195 ! ENDDEBUG 196 197 end if 198 199 ABI_DEALLOCATE(symrel_nomagn) 200 ABI_DEALLOCATE(tnons_nomagn) 201 202 end if ! Shubnikov groups 203 204 end if 205 206 end subroutine symanal