TABLE OF CONTENTS
ABINIT/m_intagm_img [ Modules ]
NAME
m_intagm_img
FUNCTION
This module provides a generic interface that allows to initialize some of the geometry variables in the case of "images". Set up: acell, scalecart, rprim, angdeg, xred, xcart, vel These variables can be defined for a set of images of the cell. They also can be be defined along a path (in the configuration space). The path must be defined with its first and last points, but also with intermediate points.
COPYRIGHT
Copyright (C) 2012-2020 ABINIT group (XG, MT) 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
iimage=index of the current image jdtset=number of the dataset looked for lenstr=actual length of the input string nimage=number of images size1,size2, ...: size of array to be read (dp_data) string=character string containing 'tags' and data. token=character string for tagging the data to be read in input string typevarphys= variable type (for dimensionality purposes)
SIDE EFFECTS
dp_data(size1,size2,...)=data to be read (double precision) tread_ok=flag to be set to 1 if the data have been found in input string
NOTES
The routine is a generic interface calling subroutine according to the number of arguments of the variable to be read
PARENTS
ingeo
CHILDREN
intagm
SOURCE
48 #if defined HAVE_CONFIG_H 49 #include "config.h" 50 #endif 51 52 #include "abi_common.h" 53 54 MODULE m_intagm_img 55 56 use defs_basis 57 use m_abicore 58 use m_errors 59 60 use m_parser, only : intagm 61 62 implicit none 63 64 private
m_intagm_img/ingeo_img_1D [ Functions ]
[ Top ] [ m_intagm_img ] [ Functions ]
NAME
intagm_img_1D
FUNCTION
Read input file variables according to images path definition (1D array)
INPUTS
PARENTS
CHILDREN
intagm
SOURCE
95 subroutine intagm_img_1D(dp_data,iimage,jdtset,lenstr,nimage,size1,string,token,tread_ok,typevarphys) 96 97 !Arguments ------------------------------------ 98 !scalars 99 integer,intent(in) :: iimage,jdtset,lenstr,nimage,size1 100 integer,intent(inout) :: tread_ok 101 real(dp),intent(inout) :: dp_data(size1) 102 character(len=*),intent(in) :: typevarphys 103 character(len=*),intent(in) :: token 104 character(len=*),intent(in) :: string 105 !arrays 106 107 !Local variables------------------------------- 108 !scalars 109 integer :: iimage_after,iimage_before,marr,tread_after,tread_before,tread_current 110 real(dp) :: alpha 111 character(len=10) :: stringimage 112 character(len=3*len(token)+10) :: token_img 113 !arrays 114 integer, allocatable :: intarr(:) 115 real(dp),allocatable :: dprarr(:),dp_data_after(:),dp_data_before(:) 116 117 ! ************************************************************************* 118 119 !Nothing to do in case of a single image 120 if (nimage<=1) return 121 122 marr=size1 123 ABI_ALLOCATE(intarr,(marr)) 124 ABI_ALLOCATE(dprarr,(marr)) 125 126 !First, try to read data for current image 127 tread_current=0 128 write(stringimage,'(i10)') iimage 129 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 130 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 131 & token_img,tread_current,typevarphys) 132 if (tread_current==1)then 133 dp_data(1:size1)=dprarr(1:size1) 134 tread_ok=1 135 end if 136 if (tread_current==0.and.iimage==nimage) then 137 ! If the image is the last one, try to read data for last image (_lastimg) 138 token_img=trim(token)//'_lastimg' 139 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 140 & token_img,tread_current,typevarphys) 141 if (tread_current==1)then 142 dp_data(1:size1)=dprarr(1:size1) 143 tread_ok=1 144 end if 145 end if 146 147 if (tread_current==0) then 148 149 ! The current image is not directly defined in the input string 150 ABI_ALLOCATE(dp_data_before,(size1)) 151 ABI_ALLOCATE(dp_data_after,(size1)) 152 153 ! Find the nearest previous defined image 154 tread_before=0;iimage_before=iimage 155 do while (iimage_before>1.and.tread_before/=1) 156 iimage_before=iimage_before-1 157 write(stringimage,'(i10)') iimage_before 158 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 159 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 160 & token_img,tread_before,typevarphys) 161 if (tread_before==1) dp_data_before(1:size1)=dprarr(1:size1) 162 end do 163 if (tread_before==0) then 164 iimage_before=1 165 dp_data_before(1:size1)=dp_data(1:size1) 166 end if 167 168 ! Find the nearest following defined image 169 tread_after=0;iimage_after=iimage 170 do while (iimage_after<nimage.and.tread_after/=1) 171 iimage_after=iimage_after+1 172 write(stringimage,'(i10)') iimage_after 173 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 174 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 175 & token_img,tread_after,typevarphys) 176 if (tread_after==1) dp_data_after(1:size1)=dprarr(1:size1) 177 if (tread_after==0.and.iimage_after==nimage) then 178 token_img=trim(token)//'_lastimg' 179 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 180 & token_img,tread_after,typevarphys) 181 if (tread_after==1) dp_data_after(1:size1)=dprarr(1:size1) 182 end if 183 end do 184 if (tread_after==0) then 185 iimage_after=nimage 186 dp_data_after(1:size1)=dp_data(1:size1) 187 end if 188 189 ! Interpolate image data 190 if (tread_before==1.or.tread_after==1) then 191 alpha=real(iimage-iimage_before,dp)/real(iimage_after-iimage_before,dp) 192 dp_data(1:size1)=dp_data_before(1:size1) & 193 & +alpha*(dp_data_after(1:size1)-dp_data_before(1:size1)) 194 tread_ok=1 195 end if 196 197 ABI_DEALLOCATE(dp_data_before) 198 ABI_DEALLOCATE(dp_data_after) 199 200 end if 201 202 ABI_DEALLOCATE(intarr) 203 ABI_DEALLOCATE(dprarr) 204 205 end subroutine intagm_img_1D
m_intagm_img/ingeo_img_2D [ Functions ]
[ Top ] [ m_intagm_img ] [ Functions ]
NAME
intagm_img_2D
FUNCTION
Read input file variables according to images path definition (2D array)
INPUTS
PARENTS
CHILDREN
intagm
SOURCE
228 subroutine intagm_img_2D(dp_data,iimage,jdtset,lenstr,nimage,size1,size2,string,token,tread_ok,typevarphys) 229 230 !Arguments ------------------------------------ 231 !scalars 232 integer,intent(in) :: iimage,jdtset,lenstr,nimage,size1,size2 233 integer,intent(inout) :: tread_ok 234 real(dp),intent(inout) :: dp_data(size1,size2) 235 character(len=*),intent(in) :: typevarphys 236 character(len=*),intent(in) :: token 237 character(len=*),intent(in) :: string 238 !arrays 239 240 !Local variables------------------------------- 241 !scalars 242 integer :: iimage_after,iimage_before,marr,tread_after,tread_before,tread_current 243 real(dp) :: alpha 244 character(len=10) :: stringimage 245 character(len=3*len(token)+10) :: token_img 246 !arrays 247 integer, allocatable :: intarr(:) 248 real(dp),allocatable :: dprarr(:),dp_data_after(:,:),dp_data_before(:,:) 249 250 ! ************************************************************************* 251 252 !Nothing to do in case of a single image 253 if (nimage<=1) return 254 255 marr=size1*size2 256 ABI_ALLOCATE(intarr,(marr)) 257 ABI_ALLOCATE(dprarr,(marr)) 258 259 !First, try to read data for current image 260 tread_current=0 261 write(stringimage,'(i10)') iimage 262 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 263 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 264 & token_img,tread_current,typevarphys) 265 if (tread_current==1)then 266 dp_data(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 267 tread_ok=1 268 end if 269 if (tread_current==0.and.iimage==nimage) then 270 ! In the image is the last one, try to read data for last image (_lastimg) 271 token_img=trim(token)//'_lastimg' 272 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 273 & token_img,tread_current,typevarphys) 274 if (tread_current==1)then 275 dp_data(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 276 tread_ok=1 277 end if 278 end if 279 280 if (tread_current==0) then 281 282 ! The current image is not directly defined in the input string 283 ABI_ALLOCATE(dp_data_before,(size1,size2)) 284 ABI_ALLOCATE(dp_data_after,(size1,size2)) 285 286 ! Find the nearest previous defined image 287 tread_before=0;iimage_before=iimage 288 do while (iimage_before>1.and.tread_before/=1) 289 iimage_before=iimage_before-1 290 write(stringimage,'(i10)') iimage_before 291 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 292 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 293 & token_img,tread_before,typevarphys) 294 if (tread_before==1) & 295 & dp_data_before(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 296 end do 297 if (tread_before==0) then 298 iimage_before=1 299 dp_data_before(1:size1,1:size2)=dp_data(1:size1,1:size2) 300 end if 301 302 ! Find the nearest following defined image 303 tread_after=0;iimage_after=iimage 304 do while (iimage_after<nimage.and.tread_after/=1) 305 iimage_after=iimage_after+1 306 write(stringimage,'(i10)') iimage_after 307 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 308 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 309 & token_img,tread_after,typevarphys) 310 if (tread_after==1) & 311 & dp_data_after(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 312 if (tread_after==0.and.iimage_after==nimage) then 313 token_img=trim(token)//'_lastimg' 314 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 315 & token_img,tread_after,typevarphys) 316 if (tread_after==1) & 317 & dp_data_after(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 318 end if 319 end do 320 if (tread_after==0) then 321 iimage_after=nimage 322 dp_data_after(1:size1,1:size2)=dp_data(1:size1,1:size2) 323 end if 324 325 ! Interpolate image data 326 if (tread_before==1.or.tread_after==1) then 327 alpha=real(iimage-iimage_before,dp)/real(iimage_after-iimage_before,dp) 328 dp_data(1:size1,1:size2)=dp_data_before(1:size1,1:size2) & 329 & +alpha*(dp_data_after(1:size1,1:size2)-dp_data_before(1:size1,1:size2)) 330 tread_ok=1 331 end if 332 333 ABI_DEALLOCATE(dp_data_before) 334 ABI_DEALLOCATE(dp_data_after) 335 336 end if 337 338 ABI_DEALLOCATE(intarr) 339 ABI_DEALLOCATE(dprarr) 340 341 end subroutine intagm_img_2D