TABLE OF CONTENTS
ABINIT/m_ingeo_img [ Modules ]
NAME
m_ingeo_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, xangst, 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-2018 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
49 #if defined HAVE_CONFIG_H 50 #include "config.h" 51 #endif 52 53 #include "abi_common.h" 54 55 MODULE m_ingeo_img 56 57 use defs_basis 58 use m_profiling_abi 59 use m_errors 60 61 use m_parser, only : intagm 62 63 implicit none 64 65 private
m_ingeo_img/ingeo_img_1D [ Functions ]
[ Top ] [ m_ingeo_img ] [ Functions ]
NAME
ingeo_img_1D
FUNCTION
Read input file variables according to images path definition (1D array)
INPUTS
PARENTS
CHILDREN
intagm
SOURCE
96 subroutine ingeo_img_1D(dp_data,iimage,jdtset,lenstr,nimage,size1,string,token,tread_ok,typevarphys) 97 98 99 !This section has been created automatically by the script Abilint (TD). 100 !Do not modify the following lines by hand. 101 #undef ABI_FUNC 102 #define ABI_FUNC 'ingeo_img_1D' 103 !End of the abilint section 104 105 implicit none 106 107 !Arguments ------------------------------------ 108 !scalars 109 integer,intent(in) :: iimage,jdtset,lenstr,nimage,size1 110 integer,intent(inout) :: tread_ok 111 real(dp),intent(inout) :: dp_data(size1) 112 character(len=*),intent(in) :: typevarphys 113 character(len=*),intent(in) :: token 114 character(len=*),intent(in) :: string 115 !arrays 116 117 !Local variables------------------------------- 118 !scalars 119 integer :: iimage_after,iimage_before,marr,tread_after,tread_before,tread_current 120 real(dp) :: alpha 121 character(len=10) :: stringimage 122 character(len=3*len(token)+10) :: token_img 123 !arrays 124 integer, allocatable :: intarr(:) 125 real(dp),allocatable :: dprarr(:),dp_data_after(:),dp_data_before(:) 126 127 ! ************************************************************************* 128 129 !Nothing to do in case of a single image 130 if (nimage<=1) return 131 132 marr=size1 133 ABI_ALLOCATE(intarr,(marr)) 134 ABI_ALLOCATE(dprarr,(marr)) 135 136 !First, try to read data for current image 137 tread_current=0 138 write(stringimage,'(i10)') iimage 139 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 140 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 141 & token_img,tread_current,typevarphys) 142 if (tread_current==1)then 143 dp_data(1:size1)=dprarr(1:size1) 144 tread_ok=1 145 end if 146 if (tread_current==0.and.iimage==nimage) then 147 ! If the image is the last one, try to read data for last image (_lastimg) 148 token_img=trim(token)//'_lastimg' 149 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 150 & token_img,tread_current,typevarphys) 151 if (tread_current==1)then 152 dp_data(1:size1)=dprarr(1:size1) 153 tread_ok=1 154 end if 155 end if 156 157 if (tread_current==0) then 158 159 ! The current image is not directly defined in the input string 160 ABI_ALLOCATE(dp_data_before,(size1)) 161 ABI_ALLOCATE(dp_data_after,(size1)) 162 163 ! Find the nearest previous defined image 164 tread_before=0;iimage_before=iimage 165 do while (iimage_before>1.and.tread_before/=1) 166 iimage_before=iimage_before-1 167 write(stringimage,'(i10)') iimage_before 168 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 169 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 170 & token_img,tread_before,typevarphys) 171 if (tread_before==1) dp_data_before(1:size1)=dprarr(1:size1) 172 end do 173 if (tread_before==0) then 174 iimage_before=1 175 dp_data_before(1:size1)=dp_data(1:size1) 176 end if 177 178 ! Find the nearest following defined image 179 tread_after=0;iimage_after=iimage 180 do while (iimage_after<nimage.and.tread_after/=1) 181 iimage_after=iimage_after+1 182 write(stringimage,'(i10)') iimage_after 183 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 184 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 185 & token_img,tread_after,typevarphys) 186 if (tread_after==1) dp_data_after(1:size1)=dprarr(1:size1) 187 if (tread_after==0.and.iimage_after==nimage) then 188 token_img=trim(token)//'_lastimg' 189 call intagm(dprarr,intarr,jdtset,marr,size1,string(1:lenstr),& 190 & token_img,tread_after,typevarphys) 191 if (tread_after==1) dp_data_after(1:size1)=dprarr(1:size1) 192 end if 193 end do 194 if (tread_after==0) then 195 iimage_after=nimage 196 dp_data_after(1:size1)=dp_data(1:size1) 197 end if 198 199 ! Interpolate image data 200 if (tread_before==1.or.tread_after==1) then 201 alpha=real(iimage-iimage_before,dp)/real(iimage_after-iimage_before,dp) 202 dp_data(1:size1)=dp_data_before(1:size1) & 203 & +alpha*(dp_data_after(1:size1)-dp_data_before(1:size1)) 204 tread_ok=1 205 end if 206 207 ABI_DEALLOCATE(dp_data_before) 208 ABI_DEALLOCATE(dp_data_after) 209 210 end if 211 212 ABI_DEALLOCATE(intarr) 213 ABI_DEALLOCATE(dprarr) 214 215 end subroutine ingeo_img_1D
m_ingeo_img/ingeo_img_2D [ Functions ]
[ Top ] [ m_ingeo_img ] [ Functions ]
NAME
ingeo_img_2D
FUNCTION
Read input file variables according to images path definition (2D array)
INPUTS
PARENTS
CHILDREN
intagm
SOURCE
238 subroutine ingeo_img_2D(dp_data,iimage,jdtset,lenstr,nimage,size1,size2,string,token,tread_ok,typevarphys) 239 240 241 !This section has been created automatically by the script Abilint (TD). 242 !Do not modify the following lines by hand. 243 #undef ABI_FUNC 244 #define ABI_FUNC 'ingeo_img_2D' 245 !End of the abilint section 246 247 implicit none 248 249 !Arguments ------------------------------------ 250 !scalars 251 integer,intent(in) :: iimage,jdtset,lenstr,nimage,size1,size2 252 integer,intent(inout) :: tread_ok 253 real(dp),intent(inout) :: dp_data(size1,size2) 254 character(len=*),intent(in) :: typevarphys 255 character(len=*),intent(in) :: token 256 character(len=*),intent(in) :: string 257 !arrays 258 259 !Local variables------------------------------- 260 !scalars 261 integer :: iimage_after,iimage_before,marr,tread_after,tread_before,tread_current 262 real(dp) :: alpha 263 character(len=10) :: stringimage 264 character(len=3*len(token)+10) :: token_img 265 !arrays 266 integer, allocatable :: intarr(:) 267 real(dp),allocatable :: dprarr(:),dp_data_after(:,:),dp_data_before(:,:) 268 269 ! ************************************************************************* 270 271 !Nothing to do in case of a single image 272 if (nimage<=1) return 273 274 marr=size1*size2 275 ABI_ALLOCATE(intarr,(marr)) 276 ABI_ALLOCATE(dprarr,(marr)) 277 278 !First, try to read data for current image 279 tread_current=0 280 write(stringimage,'(i10)') iimage 281 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 282 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 283 & token_img,tread_current,typevarphys) 284 if (tread_current==1)then 285 dp_data(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 286 tread_ok=1 287 end if 288 if (tread_current==0.and.iimage==nimage) then 289 ! In the image is the last one, try to read data for last image (_lastimg) 290 token_img=trim(token)//'_lastimg' 291 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 292 & token_img,tread_current,typevarphys) 293 if (tread_current==1)then 294 dp_data(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 295 tread_ok=1 296 end if 297 end if 298 299 if (tread_current==0) then 300 301 ! The current image is not directly defined in the input string 302 ABI_ALLOCATE(dp_data_before,(size1,size2)) 303 ABI_ALLOCATE(dp_data_after,(size1,size2)) 304 305 ! Find the nearest previous defined image 306 tread_before=0;iimage_before=iimage 307 do while (iimage_before>1.and.tread_before/=1) 308 iimage_before=iimage_before-1 309 write(stringimage,'(i10)') iimage_before 310 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 311 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 312 & token_img,tread_before,typevarphys) 313 if (tread_before==1) & 314 & dp_data_before(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 315 end do 316 if (tread_before==0) then 317 iimage_before=1 318 dp_data_before(1:size1,1:size2)=dp_data(1:size1,1:size2) 319 end if 320 321 ! Find the nearest following defined image 322 tread_after=0;iimage_after=iimage 323 do while (iimage_after<nimage.and.tread_after/=1) 324 iimage_after=iimage_after+1 325 write(stringimage,'(i10)') iimage_after 326 token_img=trim(token)//'_'//trim(adjustl(stringimage))//'img' 327 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 328 & token_img,tread_after,typevarphys) 329 if (tread_after==1) & 330 & dp_data_after(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 331 if (tread_after==0.and.iimage_after==nimage) then 332 token_img=trim(token)//'_lastimg' 333 call intagm(dprarr,intarr,jdtset,marr,size1*size2,string(1:lenstr),& 334 & token_img,tread_after,typevarphys) 335 if (tread_after==1) & 336 & dp_data_after(1:size1,1:size2)=reshape( dprarr(1:size1*size2),(/size1,size2/) ) 337 end if 338 end do 339 if (tread_after==0) then 340 iimage_after=nimage 341 dp_data_after(1:size1,1:size2)=dp_data(1:size1,1:size2) 342 end if 343 344 ! Interpolate image data 345 if (tread_before==1.or.tread_after==1) then 346 alpha=real(iimage-iimage_before,dp)/real(iimage_after-iimage_before,dp) 347 dp_data(1:size1,1:size2)=dp_data_before(1:size1,1:size2) & 348 & +alpha*(dp_data_after(1:size1,1:size2)-dp_data_before(1:size1,1:size2)) 349 tread_ok=1 350 end if 351 352 ABI_DEALLOCATE(dp_data_before) 353 ABI_DEALLOCATE(dp_data_after) 354 355 end if 356 357 ABI_DEALLOCATE(intarr) 358 ABI_DEALLOCATE(dprarr) 359 360 end subroutine ingeo_img_2D