TABLE OF CONTENTS
ABINIT/m_out_spg_anal [ Modules ]
NAME
m_out_spg_anal
FUNCTION
COPYRIGHT
Copyright (C) 1998-2024 ABINIT group (DCA, XG, GMR) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_out_spg_anal 23 24 use defs_basis 25 use m_results_out 26 use m_dtset 27 use m_abicore 28 use m_errors 29 use m_xomp 30 use m_xmpi 31 #if defined HAVE_NETCDF 32 use netcdf 33 #endif 34 use m_outvar_a_h 35 use m_outvar_i_n 36 use m_outvar_o_z 37 38 use m_parser, only : ab_dimensions 39 use m_nctk, only : create_nc_file 40 use m_symfind, only : symfind_expert, symanal, symlatt 41 use m_geometry, only : metric, mkrdim 42 use m_spgdata, only : prtspgroup 43 44 implicit none 45 46 private
ABINIT/out_spg_anal [ Functions ]
NAME
out_spg_anal
FUNCTION
Perform final spacegroup analysis of the results of ABINIT, for each dataset. Compare with the initial one, and perform analysis, with adequate warning if there was a change. Possibly echo spacegroup for all dtsets and possibly all images
INPUTS
dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables iout=unit number for echoed output - the echo is done to std_out anyhow. ndtset=number of datasets ndtset_alloc=number of datasets, corrected for allocation of at least one data set. Use for most dimensioned arrays. echo_spgroup = not relevant anymore, as at present set to 1 in the calling routine. (0 => write ; 1 => echo of spacegroup for all dtsets and possibly all images) results_out(0:ndtset_alloc)=<type results_out_type>contains the results needed for outvars, including evolving variables
OUTPUT
Only writing
NOTES
Note that this routine is called only by the processor me==0 . In consequence, no use of message and wrtout routine.
SOURCE
84 subroutine out_spg_anal(dtsets,echo_spgroup,iout,ndtset,ndtset_alloc,results_out) 85 86 !Arguments ------------------------------------ 87 !scalars 88 integer,intent(in) :: echo_spgroup,iout 89 integer,intent(in) :: ndtset,ndtset_alloc 90 !arrays 91 type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc) 92 type(results_out_type),intent(in) :: results_out(0:ndtset_alloc) 93 94 !Local variables------------------------------- 95 !scalars 96 integer, save :: counter0=1, counter1=1 97 integer :: idtset,iimage,invar_z,jdtset,msym,mu,natom,nimage,nptsym,nsym 98 integer :: ptgroupma,spgroup,symmetry_changed 99 real(dp) :: tolsym,ucvol 100 character(len=500) :: msg 101 !arrays 102 integer :: bravais(11) 103 integer, allocatable :: ptsymrel(:,:,:),symafm(:),symrel(:,:,:) 104 real(dp) :: acell(3),genafm(3),gmet(3,3),gprimd(3,3),rmet(3,3),rprim(3,3),rprimd(3,3) 105 real(dp), allocatable :: tnons(:,:) 106 107 ! ************************************************************************* 108 109 !An upper bound on the number of symmetry operations is 48 (maximum of point symmetries) times 2 (for the spin flip) 110 !times the number of atoms as the latter gives the maximum number of non-symmorphic translations. 111 !Actually, a better bound might be obtained from the minimum of the numbers of atoms of the same type, 112 !but this refinement is not needed here. 113 msym=96*dtsets(1)%natom 114 if(ndtset_alloc>1)then 115 do idtset=2,ndtset_alloc 116 msym=max(96*dtsets(idtset)%natom,msym) 117 end do 118 end if 119 ABI_MALLOC(ptsymrel,(3,3,msym)) 120 ABI_MALLOC(symafm,(msym)) 121 ABI_MALLOC(symrel,(3,3,msym)) 122 ABI_MALLOC(tnons,(3,msym)) 123 124 do idtset=1,ndtset_alloc 125 126 tolsym=dtsets(idtset)%tolsym 127 natom=dtsets(idtset)%natom 128 nimage=results_out(idtset)%nimage 129 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 130 131 do iimage=1,nimage 132 133 acell=results_out(idtset)%acell(:,iimage) 134 rprim=results_out(idtset)%rprim(:,:,iimage) 135 call mkrdim(acell,rprim,rprimd) 136 call metric(gmet,gprimd,dev_null,rmet,rprimd,ucvol) 137 138 !From rprimd and tolsym, compute bravais, nptsym and ptsymrel (with maximum size msym). 139 call symlatt(bravais,dev_null,msym,nptsym,ptsymrel,rprimd,tolsym) 140 141 !DEBUG 142 ! write(std_out,*)' out_spg_data : before symfind_expert, return, msym= ',msym 143 ! write(std_out,*)' out_spg_data : before symfind_expert, continue ' 144 ! return 145 !ENDDEBUG 146 147 invar_z=0 ; if(dtsets(idtset)%jellslab/=0 .or. dtsets(idtset)%nzchempot/=0)invar_z=2 148 149 call symfind_expert(gprimd,msym,natom,nptsym,dtsets(idtset)%nspden,nsym,& 150 dtsets(idtset)%pawspnorb,dtsets(idtset)%prtvol,ptsymrel,dtsets(idtset)%spinat,symafm,symrel,& 151 tnons,tolsym,dtsets(idtset)%typat,dtsets(idtset)%usepaw,results_out(idtset)%xred(1:3,1:natom,iimage),& 152 chrgat=dtsets(idtset)%chrgat,nucdipmom=dtsets(idtset)%nucdipmom,& 153 invardir_red=dtsets(idtset)%field_red,invar_z=invar_z) 154 155 !DEBUG 156 ! write(std_out,*)' out_spg_data : before symfind_expert, return ' 157 ! write(std_out,*)' out_spg_data : after symfind_expert, return ' 158 ! return 159 !ENDDEBUG 160 161 !Set chkprim to 0, to allow detecting increase of multiplicity 162 call symanal(bravais,0,genafm,msym,nsym,ptgroupma,rprimd,spgroup,symafm,symrel,tnons,tolsym) 163 164 symmetry_changed=0 165 if( nsym/=dtsets(idtset)%nsym .or. spgroup/=dtsets(idtset)%spgroup .or. & 166 & ptgroupma/=dtsets(idtset)%ptgroupma)then 167 symmetry_changed=1 168 endif 169 170 !DEBUG 171 ! write(std_out,*)' m_out_spg_anal : determined symmetry_changed =', symmetry_changed 172 ! write(std_out,*)' m_out_spg_anal : nsym, dtsets(idtset)%nsym=',nsym, dtsets(idtset)%nsym 173 ! write(std_out,*)' m_out_spg_anal : spgroup,dtsets(idtset)%spgroup=',spgroup,dtsets(idtset)%spgroup 174 ! write(std_out,*)' m_out_spg_anal : ptgroupma,dtsets(idtset)%ptgroupma=',ptgroupma,dtsets(idtset)%ptgroupma 175 !ENDDEBUG 176 177 if(symmetry_changed==1)then 178 if(counter0==1)then 179 write(msg,'(8a)')ch10,' The spacegroup number, the magnetic point group, and/or the number of symmetries',ch10,& 180 & ' have changed between the initial recognition based on the input file',ch10,& 181 & ' and a postprocessing based on the final acell, rprim, and xred.',ch10,& 182 & ' More details in the log file.' 183 call wrtout(iout,msg,'COLL') 184 counter0=0 185 endif 186 if(echo_spgroup==1 .and. counter1==1)then 187 write(msg,'(10a)')ch10,' The spacegroup number, the magnetic point group, and/or the number of symmetries',ch10,& 188 & ' have changed between the initial recognition based on the input file',ch10,& 189 & ' and a postprocessing based on the final acell, rprim, and xred.',ch10,& 190 & ' These modifications are detailed below.',ch10,& 191 & ' The updated tnons, symrel or symrel have NOT been reported in the final echo of variables after computation.' 192 call wrtout(std_out,msg,'COLL') 193 write(msg,'(5a)')' Such change of spacegroup, or magnetic point group might happen in several cases.',ch10,& 194 & ' (1) If spgroup (+ptgroupma) defined in the input file, but the actual groups are supergroups of these; ',ch10,& 195 & ' (2) If symrel, tnons (+symafm) defined in the input file, while the system is more symmetric; ' 196 call wrtout(std_out,msg,'COLL') 197 write(msg,'(5a)')& 198 & ' (3) If the geometry has been optimized and the final structure is more symmetric than the initial one;',ch10,& 199 & ' (4) In case of GW of BSE calculation with inversion symmetry, as nsym has been reduced in such',ch10,& 200 & ' dataset, excluding the improper symmetry operations (with determinant=-1), but not in the postprocessing.' 201 call wrtout(std_out,msg,'COLL') 202 write(msg,'(5a)')' In some case, the recognition of symmetries strongly depends on the value of tolsym.',ch10,& 203 ' You might investigate its effect by restarting abinit based on the final acell, rprim and xred,',ch10,& 204 & ' and different values for tolsym.' 205 counter1=0 206 endif 207 endif 208 209 ! Echo the spacegroup (and ptgroupma) if requested, and give more information if the symmetry changed. 210 if(echo_spgroup==1)then 211 212 if(symmetry_changed==0)then 213 if(nimage==1)then 214 call prtspgroup(bravais,genafm,std_out,jdtset,ptgroupma,spgroup) 215 else 216 call prtspgroup(bravais,genafm,std_out,jdtset,ptgroupma,spgroup,iimage=iimage) 217 endif 218 else 219 write(msg,'(2a,3i8)')ch10,' Initial data. jdtset, iimage, nsym=',jdtset,iimage,dtsets(idtset)%nsym 220 call wrtout(std_out,msg,'COLL') 221 if(nimage==1)then 222 call prtspgroup(dtsets(idtset)%bravais,dtsets(idtset)%genafm,std_out,jdtset,& 223 & dtsets(idtset)%ptgroupma,dtsets(idtset)%spgroup) 224 else 225 call prtspgroup(dtsets(idtset)%bravais,dtsets(idtset)%genafm,std_out,jdtset,& 226 & dtsets(idtset)%ptgroupma,dtsets(idtset)%spgroup,iimage=iimage) 227 endif 228 write(msg,'(a,3i8)')' Final data. jdtset, iimage, nsym=',jdtset,iimage,nsym 229 call wrtout(std_out,msg,'COLL') 230 if(nimage==1)then 231 call prtspgroup(bravais,genafm,std_out,jdtset,ptgroupma,spgroup) 232 else 233 call prtspgroup(bravais,genafm,std_out,jdtset,ptgroupma,spgroup,iimage=iimage) 234 endif 235 endif 236 237 end if ! echo_spgroup==1 238 239 enddo ! iimage 240 241 enddo ! idtset 242 243 !########################################################### 244 !## Deallocations and cleaning 245 246 ABI_FREE(ptsymrel) 247 ABI_FREE(symafm) 248 ABI_FREE(symrel) 249 ABI_FREE(tnons) 250 251 if(echo_spgroup==1)then 252 write(msg,'(a,80a)')ch10,('=',mu=1,80) 253 call wrtout(std_out,msg,'COLL') 254 endif 255 256 !************************************************************************** 257 258 end subroutine out_spg_anal