TABLE OF CONTENTS
ABINIT/m_dtset [ Modules ]
NAME
m_dtset
FUNCTION
COPYRIGHT
Copyright (C) 1992-2017 ABINIT group (XG, MG) 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 .
SOURCE
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 MODULE m_dtset 24 25 use defs_basis 26 use m_profiling_abi 27 use m_copy 28 use m_errors 29 30 use defs_abitypes, only : dataset_type 31 32 implicit none 33 34 private 35 36 public :: dtset_chkneu 37 public :: dtset_copy 38 public :: dtset_free 39 40 CONTAINS !==============================================================================
m_dtset/dtset_chkneu [ Functions ]
[ Top ] [ m_dtset ] [ Functions ]
NAME
dtset_chkneu
FUNCTION
Check neutrality of system based on band occupancies and valence charges of pseudo-atoms. Eventually initialize occ if occopt==1 or 3...8 Also return nelect, the number of valence electron per unit cell
INPUTS
charge=number of electrons missing (+) or added (-) to system (usually 0) dtset <type(dataset_type)>=all input variables in this dataset | iscf= if>0, SCF calculation ; if<=0, non SCF calculation (wtk might | not be defined) | natom=number of atoms in unit cell | nband(nkpt*nsppol)=number of bands at each k point | nkpt=number of k points | nspinor=number of spinorial components of the wavefunctions | nsppol=1 for unpolarized, 2 for spin-polarized | ntypat=number of pseudopotentials | positron=0 if electron GS calculation | 1 if positron GS calculation | 2 if electron GS calcultaion in presence of the positron | typat(natom)=atom type (integer) for each atom | wtk(nkpt)=k point weights (defined if iscf>0 or iscf==-3) | ziontypat(ntypat)=ionic charge of each pseudoatom occopt=option for occupancies
OUTPUT
Writes warning and/or aborts if error condition exists dtset <type(dataset_type)>=all input variables in this dataset | nelect=number of valence electrons per unit cell | (from counting valence electrons in psps, and taking into | account the input variable "charge")
SIDE EFFECTS
Input/Output : dtset <type(dataset_type)>=all input variables in this dataset | occ_orig(dtset%nband(1)*nkpt*nsppol)=occupation numbers for each band and k point | must be input for occopt==0 or 2, | will be an output for occopt==1 or 3 ... 8
PARENTS
invars2
CHILDREN
SOURCE
93 subroutine dtset_chkneu(charge,dtset,occopt) 94 95 96 !This section has been created automatically by the script Abilint (TD). 97 !Do not modify the following lines by hand. 98 #undef ABI_FUNC 99 #define ABI_FUNC 'dtset_chkneu' 100 use interfaces_14_hidewrite 101 !End of the abilint section 102 103 implicit none 104 105 !Arguments ------------------------------------ 106 !scalars 107 integer,intent(in) :: occopt 108 real(dp),intent(in) :: charge 109 type(dataset_type),intent(inout) :: dtset 110 111 !Local variables------------------------------- 112 !scalars 113 integer :: bantot,iatom,iband,ii,ikpt,isppol,nocc 114 real(dp) :: maxocc,nelect_occ,nelect_spin,occlast,sign_spin,zval 115 character(len=500) :: message 116 !arrays 117 real(dp),allocatable :: tmpocc(:) 118 119 ! ************************************************************************* 120 121 !(1) count nominal valence electrons according to ziontypat 122 zval=zero 123 do iatom=1,dtset%natom 124 zval=zval+dtset%ziontypat(dtset%typat(iatom)) 125 end do 126 if (dtset%positron/=1) then 127 dtset%nelect=zval-charge 128 else 129 dtset%nelect=one 130 end if 131 132 ! write(std_out,*)ch10 133 ! write(std_out,*)' chkneu : enter, dtset%nelect=',dtset%nelect 134 ! write(std_out,*)' occopt,dtset%nsppol,dtset%nspden=',occopt,dtset%nsppol,dtset%nspden 135 136 !(2) Optionally initialize occ with semiconductor occupancies 137 !(even for a metal : at this stage, the eigenenergies are unknown) 138 !Note that nband(1)=nband(2) in this section, as occopt=2 is avoided. 139 if(occopt==1 .or. (occopt>=3 .and. occopt<=8) )then 140 ! Here, initialize a real(dp) variable giving the 141 ! maximum occupation number per band 142 maxocc=2.0_dp/real(dtset%nsppol*dtset%nspinor,dp) 143 144 ! Determine the number of bands fully or partially occupied 145 nocc=(dtset%nelect-1.0d-8)/maxocc + 1 146 ! Occupation number of the highest level 147 occlast=dtset%nelect-maxocc*(nocc-1) 148 149 !write(std_out,*)' maxocc,nocc,occlast=',maxocc,nocc,occlast 150 151 152 ! The number of allowed bands must be sufficiently large 153 if( nocc<=dtset%nband(1)*dtset%nsppol .or. dtset%iscf==-2) then 154 155 if(dtset%iscf==-2 .and. nocc>dtset%nband(1)*dtset%nsppol)nocc=dtset%nband(1)*dtset%nsppol 156 157 ! First treat the case where the spin magnetization is not imposed, is zero with nspden==1, or has sufficient flexibility 158 ! for a target not to be matched at the initialisation, but later 159 if(abs(dtset%spinmagntarget+99.99_dp)<tol8 .or. (dtset%nspden==4) .or. & 160 & (abs(dtset%spinmagntarget)<tol8.and.dtset%nspden==1))then 161 162 ! Use a temporary array for defining occupation numbers 163 ABI_ALLOCATE(tmpocc,(dtset%nband(1)*dtset%nsppol)) 164 ! First do it for fully occupied bands 165 if (1<nocc) tmpocc(1:nocc-1)=maxocc 166 ! Then, do it for highest occupied band 167 if (1<=nocc) tmpocc(nocc)=occlast 168 ! Finally do it for eventual unoccupied bands 169 if ( nocc<dtset%nband(1)*dtset%nsppol ) tmpocc(nocc+1:dtset%nband(1)*dtset%nsppol)=0.0_dp 170 171 ! Now copy the tmpocc array in the occ array, taking into account the spin 172 if(dtset%nsppol==1)then 173 do ikpt=1,dtset%nkpt 174 dtset%occ_orig(1+(ikpt-1)*dtset%nband(1):ikpt*dtset%nband(1))=tmpocc(:) 175 end do 176 else 177 do ikpt=1,dtset%nkpt 178 do iband=1,dtset%nband(1) 179 do isppol=1,dtset%nsppol 180 dtset%occ_orig(iband+dtset%nband(1)*(ikpt-1+dtset%nkpt*(isppol-1))) = & 181 & tmpocc(isppol+dtset%nsppol*(iband-1)) 182 end do 183 end do 184 end do 185 end if 186 ABI_DEALLOCATE(tmpocc) 187 188 ! Second, treat the case in which one imposes the spin magnetization (only possible for nspden==2) 189 ! Also treat antiferromagnetic case (nsppol==1, nspden==2), although spinmagntarget must be zero 190 else if(abs(dtset%spinmagntarget+99.99_dp)>1.0d-10)then 191 do isppol=1,dtset%nsppol 192 sign_spin=real(3-2*isppol,dp) 193 nelect_spin=half*(dtset%nelect*maxocc+sign_spin*dtset%spinmagntarget) 194 195 !write(std_out,*)' isppol,sign_spin,nelect_spin=',isppol,sign_spin,nelect_spin 196 ! Determines the last state, and its occupation 197 if(abs(nint(nelect_spin)-nelect_spin)<tol10)then 198 nocc=nint(nelect_spin/maxocc) 199 occlast=maxocc 200 else 201 nocc=ceiling(nelect_spin/maxocc) 202 occlast=nelect_spin-(real(nocc,dp)-one)*maxocc 203 end if 204 !write(std_out,*)' dtset%nband(1),maxocc,occlast=',dtset%nband(1),maxocc,occlast 205 if(dtset%nband(1)*nint(maxocc)<nocc)then 206 write(message, '(a,i4,a, a,2i6,a, a,es16.6,a, a,es16.6,6a)' )& 207 & 'Initialization of occ, with nspden=',dtset%nspden,ch10,& 208 & 'number of bands=',dtset%nband(1:2),ch10,& 209 & 'number of electrons=',dtset%nelect,ch10,& 210 & 'and spinmagntarget=',dtset%spinmagntarget,ch10,& 211 & 'This combination is not possible, because of a lack of bands.',ch10,& 212 & 'Action : modify input file ... ',ch10,& 213 & '(you should likely increase nband, but also check nspden, nspinor, nsppol, and spinmagntarget)' 214 MSG_ERROR(message) 215 end if 216 do ikpt=1,dtset%nkpt 217 ! Fill all bands, except the upper one 218 if(dtset%nband(1)>1)then 219 do iband=1,nocc-1 220 dtset%occ_orig(iband+dtset%nband(1)*(ikpt-1+dtset%nkpt*(isppol-1)))=maxocc 221 end do 222 end if 223 ! Fill the upper occupied band 224 dtset%occ_orig(nocc+dtset%nband(1)*(ikpt-1+dtset%nkpt*(isppol-1)))=occlast 225 end do 226 !write(std_out,*)' dtset%occ_orig(1)=',dtset%occ_orig(1) 227 end do 228 229 else 230 write(message, '(a,i4,a,a,es16.6,6a)' )& 231 & 'Initialization of occ, with nspden=',dtset%nspden,ch10,& 232 & 'and spinmagntarget=',dtset%spinmagntarget,ch10,& 233 & 'This combination is not possible.',ch10,& 234 & 'Action : modify input file ... ',ch10,& 235 & '(check nspden, nspinor, nsppol and spinmagntarget)' 236 MSG_ERROR(message) 237 end if 238 239 ! Now print the values 240 if(dtset%nsppol==1)then 241 242 write(message, '(a,i4,a,a)' ) & 243 & ' chkneu : initialized the occupation numbers for occopt= ',occopt,', spin-unpolarized or antiferromagnetic case : ' 244 call wrtout(std_out,message,'COLL') 245 do ii=0,(dtset%nband(1)-1)/12 246 write(message,'(12f6.2)') dtset%occ_orig( 1+ii*12 : min(12+ii*12,dtset%nband(1)) ) 247 call wrtout(std_out,message,'COLL') 248 end do 249 250 else 251 252 write(message, '(a,i4,a,a)' ) & 253 & ' chkneu : initialized the occupation numbers for occopt= ',occopt,& 254 & ch10,' spin up values : ' 255 call wrtout(std_out,message,'COLL') 256 do ii=0,(dtset%nband(1)-1)/12 257 write(message,'(12f6.2)') dtset%occ_orig( 1+ii*12 : min(12+ii*12,dtset%nband(1)) ) 258 call wrtout(std_out,message,'COLL') 259 end do 260 write(message, '(a)' ) ' spin down values : ' 261 call wrtout(std_out,message,'COLL') 262 do ii=0,(dtset%nband(1)-1)/12 263 write(message,'(12f6.2)') & 264 & dtset%occ_orig( 1+ii*12+dtset%nkpt*dtset%nband(1) : min(12+ii*12,dtset%nband(1))+dtset%nkpt*dtset%nband(1) ) 265 call wrtout(std_out,message,'COLL') 266 end do 267 268 end if 269 270 ! Here, treat the case when the number of allowed bands is not large enough 271 else 272 write(message, '(a,i4,a,a,a,a,a,a,a,a)' )& 273 & 'Initialization of occ, with occopt=',occopt,ch10,& 274 & 'There are not enough bands to get charge balance right',ch10,& 275 & 'Action: modify input file ... ',ch10,& 276 & '(check the pseudopotential charges, the variable charge,',ch10,& 277 & 'and the declared number of bands, nband)' 278 MSG_ERROR(message) 279 end if 280 end if 281 282 !The remaining of the routine is for SCF runs and special options 283 if(dtset%iscf>0 .or. dtset%iscf==-1 .or. dtset%iscf==-3)then 284 285 ! (3) count electrons in bands (note : in case occ has just been 286 ! initialized, point (3) and (4) is a trivial test 287 nelect_occ=0.0_dp 288 bantot=0 289 do isppol=1,dtset%nsppol 290 do ikpt=1,dtset%nkpt 291 do iband=1,dtset%nband(ikpt+(isppol-1)*dtset%nkpt) 292 bantot=bantot+1 293 nelect_occ=nelect_occ+dtset%wtk(ikpt)*dtset%occ_orig(bantot) 294 end do 295 end do 296 end do 297 298 ! (4) if dtset%iscf/=-3, dtset%nelect must equal nelect_occ 299 ! if discrepancy exceeds tol11, give warning; tol8, stop with error 300 301 if (abs(nelect_occ-dtset%nelect)>tol11 .and. dtset%iscf/=-3) then 302 303 ! There is a discrepancy 304 write(message, & 305 & '(a,a,e16.8,a,e16.8,a,a,a,e22.14,a,a,a,a,a,a,a)' ) ch10,& 306 & ' chkneu: nelect_occ=',nelect_occ,', zval=',zval,',',ch10,& 307 & ' and input value of charge=',charge,',',ch10,& 308 & ' nelec_occ is computed from occ and wtk',ch10,& 309 & ' zval is nominal charge of all nuclei, computed from zion (read in psp),',ch10,& 310 & ' charge is an input variable (usually 0).' 311 call wrtout(std_out,message,'COLL') 312 313 if (abs(nelect_occ-dtset%nelect)>tol8) then 314 ! The discrepancy is severe 315 write(message,'(a,a,e9.2,a,a)')ch10,& 316 & 'These must obey zval-nelect_occ=charge to better than ',tol8,ch10,& 317 & ' This is not the case. ' 318 else 319 ! The discrepancy is not so severe 320 write(message, '(a,a,e9.2)' )ch10,'These should obey zval-nelect_occ=charge to better than ',tol11 321 end if 322 MSG_WARNING(message) 323 324 write(message, '(a,a,a,a,a,a)' ) & 325 & 'Action : check input file for occ,wtk, and charge.',ch10,& 326 & 'Note that wtk is NOT automatically normalized when occopt=2,',ch10,& 327 & 'but IS automatically normalized otherwise.',ch10 328 call wrtout(std_out,message,'COLL') 329 330 ! If the discrepancy is severe, stop 331 if (abs(nelect_occ-dtset%nelect)>tol8)then 332 MSG_ERROR(message) 333 end if 334 335 end if 336 end if ! End the condition dtset%iscf>0 or -1 or -3 . 337 338 end subroutine dtset_chkneu
m_dtset/dtset_copy [ Functions ]
[ Top ] [ m_dtset ] [ Functions ]
NAME
dtset_copy
FUNCTION
Copy all values of dataset dtin to dataset dtout. allocatables of dtout are allocated if required. Use dtset_free() to free a dataset after use.
INPUTS
dtin <type(dataset_type)>=all input variables in this dataset
OUTPUT
dtout <type(dataset_type)>
PARENTS
calc_vhxc_me,chkinp,dfpt_looppert,driver,gwls_hamiltonian,hybrid_corr m_io_kss,m_kxc,xchybrid_ncpp_cc
CHILDREN
SOURCE
365 subroutine dtset_copy(dtout, dtin) 366 367 368 !This section has been created automatically by the script Abilint (TD). 369 !Do not modify the following lines by hand. 370 #undef ABI_FUNC 371 #define ABI_FUNC 'dtset_copy' 372 !End of the abilint section 373 374 implicit none 375 376 !Arguments ------------------------------------ 377 !scalars 378 type(dataset_type),intent(in) :: dtin 379 type(dataset_type),intent(out) :: dtout 380 381 ! ************************************************************************* 382 383 DBG_ENTER("COLL") 384 385 !@dataset_type 386 387 !BEGIN VARIABLES FOR @Bethe-Salpeter 388 dtout%bs_algorithm = dtin%bs_algorithm 389 dtout%bs_haydock_niter = dtin%bs_haydock_niter 390 dtout%bs_exchange_term = dtin%bs_exchange_term 391 dtout%bs_coulomb_term = dtin%bs_coulomb_term 392 dtout%bs_calctype = dtin%bs_calctype 393 dtout%bs_coupling = dtin%bs_coupling 394 395 dtout%bs_haydock_tol = dtin%bs_haydock_tol 396 dtout%bs_hayd_term = dtin%bs_hayd_term 397 398 dtout%bs_interp_m3_width = dtin%bs_interp_m3_width 399 dtout%bs_interp_method = dtin%bs_interp_method 400 dtout%bs_interp_mode = dtin%bs_interp_mode 401 dtout%bs_interp_prep = dtin%bs_interp_prep 402 dtout%bs_interp_rl_nb = dtin%bs_interp_rl_nb 403 404 dtout%bs_interp_kmult(:) = dtin%bs_interp_kmult(:) 405 dtout%bs_eh_cutoff(:) = dtin%bs_eh_cutoff(:) 406 dtout%bs_freq_mesh(:) = dtin%bs_freq_mesh(:) 407 !END VARIABLES FOR @Bethe-Salpeter. 408 409 !Copy integers from dtin to dtout 410 dtout%iomode = dtin%iomode 411 dtout%accuracy = dtin%accuracy 412 dtout%adpimd = dtin%adpimd 413 dtout%autoparal = dtin%autoparal 414 dtout%awtr = dtin%awtr 415 dtout%bandpp = dtin%bandpp 416 dtout%bdeigrf = dtin%bdeigrf 417 dtout%berryopt = dtin%berryopt 418 dtout%berrysav = dtin%berrysav 419 dtout%berrystep = dtin%berrystep 420 dtout%brvltt = dtin%brvltt 421 dtout%bs_nstates = dtin%bs_nstates 422 dtout%builtintest = dtin%builtintest 423 dtout%cd_customnimfrqs = dtin%cd_customnimfrqs 424 dtout%cd_frqim_method = dtin%cd_frqim_method 425 dtout%cd_full_grid = dtin%cd_full_grid 426 dtout%chkdilatmx = dtin%chkdilatmx 427 dtout%chkexit = dtin%chkexit 428 dtout%chkprim = dtin%chkprim 429 dtout%chksymbreak = dtin%chksymbreak 430 dtout%cineb_start = dtin%cineb_start 431 dtout%cgtyphf = dtin%cgtyphf 432 dtout%delayperm = dtin%delayperm 433 dtout%diismemory = dtin%diismemory 434 dtout%dmatpuopt = dtin%dmatpuopt 435 dtout%dmatudiag = dtin%dmatudiag 436 dtout%dmft_dc = dtin%dmft_dc 437 dtout%dmft_entropy = dtin%dmft_entropy 438 dtout%dmft_iter = dtin%dmft_iter 439 dtout%dmft_nlambda = dtin%dmft_nlambda 440 dtout%dmft_mxsf = dtin%dmft_mxsf 441 dtout%dmft_nwlo = dtin%dmft_nwlo 442 dtout%dmft_nwli = dtin%dmft_nwli 443 dtout%dmft_read_occnd = dtin%dmft_read_occnd 444 dtout%dmft_rslf = dtin%dmft_rslf 445 dtout%dmft_solv = dtin%dmft_solv 446 dtout%dmft_t2g = dtin%dmft_t2g 447 dtout%dmft_tolfreq = dtin%dmft_tolfreq 448 dtout%dmft_tollc = dtin%dmft_tollc 449 dtout%dmftbandi = dtin%dmftbandi 450 dtout%dmftbandf = dtin%dmftbandf 451 dtout%dmftcheck = dtin%dmftcheck 452 dtout%dmftctqmc_basis = dtin%dmftctqmc_basis 453 dtout%dmftctqmc_check = dtin%dmftctqmc_check 454 dtout%dmftctqmc_correl = dtin%dmftctqmc_correl 455 dtout%dmftctqmc_gmove = dtin%dmftctqmc_gmove 456 dtout%dmftctqmc_grnns = dtin%dmftctqmc_grnns 457 dtout%dmftctqmc_meas = dtin%dmftctqmc_meas 458 dtout%dmftctqmc_mrka = dtin%dmftctqmc_mrka 459 dtout%dmftctqmc_mov = dtin%dmftctqmc_mov 460 dtout%dmftctqmc_order = dtin%dmftctqmc_order 461 dtout%dmftctqmc_triqs_nleg = dtin%dmftctqmc_triqs_nleg 462 dtout%dmftqmc_n = dtin%dmftqmc_n 463 dtout%dmftqmc_l = dtin%dmftqmc_l 464 dtout%dmftqmc_seed = dtin%dmftqmc_seed 465 dtout%dmftqmc_therm = dtin%dmftqmc_therm 466 dtout%d3e_pert1_elfd = dtin%d3e_pert1_elfd 467 dtout%d3e_pert1_phon = dtin%d3e_pert1_phon 468 dtout%d3e_pert2_elfd = dtin%d3e_pert2_elfd 469 dtout%d3e_pert2_phon = dtin%d3e_pert2_phon 470 dtout%d3e_pert3_elfd = dtin%d3e_pert3_elfd 471 dtout%d3e_pert3_phon = dtin%d3e_pert3_phon 472 dtout%efmas = dtin%efmas 473 dtout%efmas_calc_dirs = dtin%efmas_calc_dirs 474 dtout%efmas_deg = dtin%efmas_deg 475 dtout%efmas_dim = dtin%efmas_dim 476 dtout%efmas_n_dirs = dtin%efmas_n_dirs 477 dtout%efmas_ntheta = dtin%efmas_ntheta 478 dtout%enunit = dtin%enunit 479 480 ! begin eph variables 481 dtout%asr = dtin%asr 482 dtout%dipdip = dtin%dipdip 483 dtout%chneut = dtin%chneut 484 485 dtout%eph_mustar = dtin%eph_mustar 486 dtout%eph_intmeth = dtin%eph_intmeth 487 dtout%eph_extrael = dtin%eph_extrael 488 dtout%eph_fermie = dtin%eph_fermie 489 dtout%eph_fsmear = dtin%eph_fsmear 490 dtout%eph_fsewin = dtin%eph_fsewin 491 dtout%eph_ngqpt_fine = dtin%eph_ngqpt_fine 492 dtout%eph_task = dtin%eph_task 493 dtout%eph_transport = dtin%eph_transport 494 495 dtout%ph_wstep = dtin%ph_wstep 496 dtout%ph_intmeth = dtin%ph_intmeth 497 dtout%symdynmat = dtin%symdynmat 498 dtout%ph_nqshift = dtin%ph_nqshift 499 if (allocated(dtin%ph_qshift)) call alloc_copy(dtin%ph_qshift, dtout%ph_qshift) 500 dtout%ph_smear = dtin%ph_smear 501 dtout%ddb_ngqpt = dtin%ddb_ngqpt 502 dtout%ddb_shiftq = dtin%ddb_shiftq 503 504 dtout%ph_ndivsm = dtin%ph_ndivsm 505 dtout%ph_nqpath = dtin%ph_nqpath 506 dtout%ph_ngqpt = dtin%ph_ngqpt 507 if (allocated(dtin%ph_qpath)) call alloc_copy(dtin%ph_qpath, dtout%ph_qpath) 508 ! end eph variables 509 510 dtout%exchn2n3d = dtin%exchn2n3d 511 dtout%extrapwf = dtin%extrapwf 512 dtout%pawfatbnd = dtin%pawfatbnd 513 dtout%fermie_nest = dtin%fermie_nest 514 dtout%fftgw = dtin%fftgw 515 dtout%fockoptmix = dtin%fockoptmix 516 dtout%freqim_alpha = dtin%freqim_alpha 517 dtout%freqremin = dtin%freqremin 518 dtout%freqremax = dtin%freqremax 519 dtout%freqspmin = dtin%freqspmin 520 dtout%freqspmax = dtin%freqspmax 521 dtout%frzfermi = dtin%frzfermi 522 dtout%ga_algor = dtin%ga_algor 523 dtout%ga_fitness = dtin%ga_fitness 524 dtout%ga_n_rules = dtin%ga_n_rules 525 dtout%getbseig = dtin%getbseig 526 dtout%getbsreso = dtin%getbsreso 527 dtout%getbscoup = dtin%getbscoup 528 dtout%getcell = dtin%getcell 529 dtout%getddb = dtin%getddb 530 dtout%getddk = dtin%getddk 531 dtout%getdelfd = dtin%getdelfd 532 dtout%getdkdk = dtin%getdkdk 533 dtout%getdkde = dtin%getdkde 534 dtout%getden = dtin%getden 535 dtout%getgam_eig2nkq = dtin%getgam_eig2nkq 536 dtout%gethaydock = dtin%gethaydock 537 dtout%getocc = dtin%getocc 538 dtout%getpawden = dtin%getpawden 539 dtout%getqps = dtin%getqps 540 dtout%getscr = dtin%getscr 541 dtout%getsuscep = dtin%getsuscep 542 dtout%getvel = dtin%getvel 543 dtout%getwfk = dtin%getwfk 544 dtout%getwfkfine = dtin%getwfkfine 545 dtout%getwfq = dtin%getwfq 546 dtout%getxcart = dtin%getxcart 547 dtout%getxred = dtin%getxred 548 dtout%get1den = dtin%get1den 549 dtout%get1wf = dtin%get1wf 550 dtout%goprecon = dtin%goprecon 551 dtout%gpu_linalg_limit = dtin%gpu_linalg_limit 552 dtout%gwcalctyp = dtin%gwcalctyp 553 dtout%gwcomp = dtin%gwcomp 554 dtout%gwencomp = dtin%gwencomp 555 dtout%gwmem = dtin%gwmem 556 dtout%gwpara = dtin%gwpara 557 dtout%gwgamma = dtin%gwgamma 558 dtout%gwrpacorr = dtin%gwrpacorr 559 dtout%gwfockmix = dtin%gwfockmix 560 dtout%gw_customnfreqsp = dtin%gw_customnfreqsp 561 dtout%gw_nqlwl = dtin%gw_nqlwl 562 dtout%gw_nstep = dtin%gw_nstep 563 dtout%gw_frqim_inzgrid = dtin%gw_frqim_inzgrid 564 dtout%gw_frqre_inzgrid = dtin%gw_frqre_inzgrid 565 dtout%gw_frqre_tangrid = dtin%gw_frqre_tangrid 566 dtout%gw_invalid_freq = dtin%gw_invalid_freq 567 dtout%gw_qprange = dtin%gw_qprange 568 dtout%gw_sctype = dtin%gw_sctype 569 dtout%gw_sigxcore = dtin%gw_sigxcore 570 dtout%gw_toldfeig = dtin%gw_toldfeig 571 dtout%gwls_stern_kmax= dtin%gwls_stern_kmax 572 dtout%gwls_npt_gauss_quad = dtin%gwls_npt_gauss_quad 573 dtout%gwls_diel_model= dtin%gwls_diel_model 574 dtout%gwls_print_debug = dtin%gwls_print_debug 575 dtout%gwls_nseeds = dtin%gwls_nseeds 576 dtout%gwls_n_proj_freq = dtin%gwls_n_proj_freq 577 dtout%gwls_kmax_complement = dtin%gwls_kmax_complement 578 dtout%gwls_kmax_poles = dtin%gwls_kmax_poles 579 dtout%gwls_kmax_analytic = dtin%gwls_kmax_analytic 580 dtout%gwls_kmax_numeric = dtin%gwls_kmax_numeric 581 dtout%gwls_band_index = dtin%gwls_band_index 582 dtout%gwls_exchange = dtin%gwls_exchange 583 dtout%gwls_correlation = dtin%gwls_correlation 584 dtout%gwls_first_seed = dtin%gwls_first_seed 585 dtout%gwls_recycle = dtin%gwls_recycle 586 dtout%iboxcut = dtin%iboxcut 587 dtout%icoulomb = dtin%icoulomb 588 dtout%icutcoul = dtin%icutcoul 589 dtout%ieig2rf = dtin%ieig2rf 590 dtout%imgmov = dtin%imgmov 591 dtout%inclvkb = dtin%inclvkb 592 dtout%intxc = dtin%intxc 593 dtout%ionmov = dtin%ionmov 594 dtout%densfor_pred = dtin%densfor_pred 595 dtout%iprcel = dtin%iprcel 596 dtout%iprcfc = dtin%iprcfc 597 dtout%irandom = dtin%irandom 598 dtout%irdbseig = dtin%irdbseig 599 dtout%irdbsreso = dtin%irdbsreso 600 dtout%irdbscoup = dtin%irdbscoup 601 dtout%irdddb = dtin%irdddb 602 dtout%irdddk = dtin%irdddk 603 dtout%irdden = dtin%irdden 604 dtout%irdhaydock = dtin%irdhaydock 605 dtout%irdpawden = dtin%irdpawden 606 dtout%irdqps = dtin%irdqps 607 dtout%irdscr = dtin%irdscr 608 dtout%irdsuscep = dtin%irdsuscep 609 dtout%irdvdw = dtin%irdvdw 610 dtout%irdwfk = dtin%irdwfk 611 dtout%irdwfkfine = dtin%irdwfkfine 612 dtout%irdwfq = dtin%irdwfq 613 dtout%ird1den = dtin%ird1den 614 dtout%ird1wf = dtin%ird1wf 615 dtout%iscf = dtin%iscf 616 dtout%isecur = dtin%isecur 617 dtout%istatimg = dtin%istatimg 618 dtout%istatr = dtin%istatr 619 dtout%istatshft = dtin%istatshft 620 dtout%ixc = dtin%ixc 621 dtout%ixcpositron = dtin%ixcpositron 622 dtout%jdtset = dtin%jdtset 623 dtout%jellslab = dtin%jellslab 624 dtout%kptopt = dtin%kptopt 625 dtout%kssform = dtin%kssform 626 dtout%localrdwf = dtin%localrdwf 627 #if defined HAVE_LOTF 628 dtout%lotf_classic = dtin%lotf_classic 629 dtout%lotf_nitex = dtin%lotf_nitex 630 dtout%lotf_nneigx = dtin%lotf_nneigx 631 dtout%lotf_version = dtin%lotf_version 632 #endif 633 dtout%magconon = dtin%magconon 634 dtout%maxnsym = dtin%maxnsym 635 dtout%max_ncpus = dtin%max_ncpus 636 dtout%mband = dtin%mband 637 dtout%mdf_epsinf = dtin%mdf_epsinf 638 dtout%mep_solver = dtin%mep_solver 639 dtout%mem_test = dtin%mem_test 640 dtout%mffmem = dtin%mffmem 641 dtout%mgfft = dtin%mgfft 642 dtout%mgfftdg = dtin%mgfftdg 643 dtout%mkmem = dtin%mkmem 644 dtout%mkqmem = dtin%mkqmem 645 dtout%mk1mem = dtin%mk1mem 646 dtout%mpw = dtin%mpw 647 dtout%mqgrid = dtin%mqgrid 648 dtout%mqgriddg = dtin%mqgriddg 649 dtout%natom = dtin%natom 650 dtout%natrd = dtin%natrd 651 dtout%natsph = dtin%natsph 652 dtout%natsph_extra = dtin%natsph_extra 653 dtout%natpawu = dtin%natpawu 654 dtout%natvshift = dtin%natvshift 655 dtout%nbdblock = dtin%nbdblock 656 dtout%nbdbuf = dtin%nbdbuf 657 dtout%nbandhf = dtin%nbandhf 658 dtout%nberry = dtin%nberry 659 dtout%nc_xccc_gspace = dtin%nc_xccc_gspace 660 dtout%nbandkss = dtin%nbandkss 661 dtout%nconeq = dtin%nconeq 662 dtout%nctime = dtin%nctime 663 dtout%ndtset = dtin%ndtset 664 dtout%ndynimage = dtin%ndynimage 665 dtout%neb_algo = dtin%neb_algo 666 dtout%nfft = dtin%nfft 667 dtout%nfftdg = dtin%nfftdg 668 dtout%nfreqim = dtin%nfreqim 669 dtout%nfreqre = dtin%nfreqre 670 dtout%nfreqsp = dtin%nfreqsp 671 dtout%nimage = dtin%nimage 672 dtout%nkpt = dtin%nkpt 673 dtout%nkpthf = dtin%nkpthf 674 dtout%nkptgw = dtin%nkptgw 675 dtout%nline = dtin%nline 676 dtout%nnsclo = dtin%nnsclo 677 dtout%nnsclohf = dtin%nnsclohf 678 dtout%nomegasf = dtin%nomegasf 679 dtout%nomegasi = dtin%nomegasi 680 dtout%nomegasrd = dtin%nomegasrd 681 dtout%npband = dtin%npband 682 dtout%npfft = dtin%npfft 683 dtout%nphf = dtin%nphf 684 dtout%npimage = dtin%npimage 685 dtout%npkpt = dtin%npkpt 686 dtout%nppert = dtin%nppert 687 dtout%npspinor = dtin%npspinor 688 dtout%npsp = dtin%npsp 689 dtout%npspalch = dtin%npspalch 690 dtout%npulayit = dtin%npulayit 691 dtout%npvel = dtin%npvel 692 dtout%npweps = dtin%npweps 693 dtout%npwkss = dtin%npwkss 694 dtout%npwsigx = dtin%npwsigx 695 dtout%npwwfn = dtin%npwwfn 696 dtout%np_slk = dtin%np_slk 697 dtout%nqpt = dtin%nqpt 698 dtout%nqptdm = dtin%nqptdm 699 dtout%nscforder = dtin%nscforder 700 dtout%nshiftk = dtin%nshiftk 701 dtout%nshiftk_orig = dtin%nshiftk_orig 702 dtout%nspden = dtin%nspden 703 dtout%nspinor = dtin%nspinor 704 dtout%nsppol = dtin%nsppol 705 dtout%nstep = dtin%nstep 706 dtout%nsym = dtin%nsym 707 dtout%ntime = dtin%ntime 708 dtout%ntimimage = dtin%ntimimage 709 dtout%ntypalch = dtin%ntypalch 710 dtout%ntypat = dtin%ntypat 711 dtout%ntyppure = dtin%ntyppure 712 dtout%nwfshist = dtin%nwfshist 713 dtout%nzchempot = dtin%nzchempot 714 dtout%occopt = dtin%occopt 715 dtout%optcell = dtin%optcell 716 dtout%optdriver = dtin%optdriver 717 dtout%optforces = dtin%optforces 718 dtout%optnlxccc = dtin%optnlxccc 719 dtout%optstress = dtin%optstress 720 dtout%ortalg = dtin%ortalg 721 dtout%paral_atom = dtin%paral_atom 722 dtout%paral_kgb = dtin%paral_kgb 723 dtout%paral_rf = dtin%paral_rf 724 dtout%pawcpxocc = dtin%pawcpxocc 725 dtout%pawcross = dtin%pawcross 726 dtout%pawlcutd = dtin%pawlcutd 727 dtout%pawlmix = dtin%pawlmix 728 dtout%pawmixdg = dtin%pawmixdg 729 dtout%pawnhatxc = dtin%pawnhatxc 730 dtout%pawnphi = dtin%pawnphi 731 dtout%pawntheta = dtin%pawntheta 732 dtout%pawnzlm = dtin%pawnzlm 733 dtout%pawoptmix = dtin%pawoptmix 734 dtout%pawoptosc = dtin%pawoptosc 735 dtout%pawprtdos = dtin%pawprtdos 736 dtout%pawprtvol = dtin%pawprtvol 737 dtout%pawprtwf = dtin%pawprtwf 738 dtout%pawprt_k = dtin%pawprt_k 739 dtout%pawprt_b = dtin%pawprt_b 740 dtout%pawspnorb = dtin%pawspnorb 741 dtout%pawstgylm = dtin%pawstgylm 742 dtout%pawsushat = dtin%pawsushat 743 dtout%pawusecp = dtin%pawusecp 744 dtout%pawujat = dtin%pawujat 745 dtout%macro_uj = dtin%macro_uj 746 dtout%pawujrad = dtin%pawujrad 747 dtout%pawujv = dtin%pawujv 748 dtout%pawxcdev = dtin%pawxcdev 749 dtout%pimd_constraint = dtin%pimd_constraint 750 dtout%pitransform = dtin%pitransform 751 dtout%plowan_compute = dtin%plowan_compute 752 dtout%plowan_bandi = dtin%plowan_bandi 753 dtout%plowan_bandf = dtin%plowan_bandf 754 dtout%plowan_natom = dtin%plowan_natom 755 dtout%plowan_nt = dtin%plowan_nt 756 dtout%plowan_realspace = dtin%plowan_realspace 757 dtout%posdoppler = dtin%posdoppler 758 dtout%positron = dtin%positron 759 dtout%posnstep = dtin%posnstep 760 dtout%ppmodel = dtin%ppmodel 761 dtout%prepanl = dtin%prepanl 762 dtout%prepgkk = dtin%prepgkk 763 dtout%prtbbb = dtin%prtbbb 764 dtout%prtbltztrp = dtin%prtbltztrp 765 dtout%prtcif = dtin%prtcif 766 dtout%prtden = dtin%prtden 767 dtout%prtdensph = dtin%prtdensph 768 dtout%prtdipole = dtin%prtdipole 769 dtout%prtdos = dtin%prtdos 770 dtout%prtdosm = dtin%prtdosm 771 dtout%prtebands = dtin%prtebands ! TODO prteig could be replaced by prtebands... 772 dtout%prtefg = dtin%prtefg 773 dtout%prteig = dtin%prteig 774 dtout%prtelf = dtin%prtelf 775 dtout%prtfc = dtin%prtfc 776 dtout%prtfsurf = dtin%prtfsurf 777 dtout%prtgsr = dtin%prtgsr 778 dtout%prtgden = dtin%prtgden 779 dtout%prtgeo = dtin%prtgeo 780 dtout%prtgkk = dtin%prtgkk 781 dtout%prtkden = dtin%prtkden 782 dtout%prtkpt = dtin%prtkpt 783 dtout%prtlden = dtin%prtlden 784 dtout%prtnabla = dtin%prtnabla 785 dtout%prtnest = dtin%prtnest 786 dtout%prtphbands = dtin%prtphbands 787 dtout%prtphdos = dtin%prtphdos 788 dtout%prtphsurf = dtin%prtphsurf 789 dtout%prtposcar = dtin%prtposcar 790 dtout%prtpot = dtin%prtpot 791 dtout%prtpsps = dtin%prtpsps 792 dtout%prtspcur = dtin%prtspcur 793 dtout%prtsuscep = dtin%prtsuscep 794 dtout%prtstm = dtin%prtstm 795 dtout%prtvclmb = dtin%prtvclmb 796 dtout%prtvdw = dtin%prtvdw 797 dtout%prtvha = dtin%prtvha 798 dtout%prtvhxc = dtin%prtvhxc 799 dtout%prtvol = dtin%prtvol 800 dtout%prtvolimg = dtin%prtvolimg 801 dtout%prtvpsp = dtin%prtvpsp 802 dtout%prtvxc = dtin%prtvxc 803 dtout%prtwant = dtin%prtwant 804 dtout%prtwf = dtin%prtwf 805 dtout%prtwf_full = dtin%prtwf_full 806 dtout%prtxml = dtin%prtxml 807 dtout%prt1dm = dtin%prt1dm 808 dtout%ptgroupma = dtin%ptgroupma 809 dtout%qptopt = dtin%qptopt 810 dtout%random_atpos = dtin%random_atpos 811 dtout%recgratio = dtin%recgratio 812 dtout%recnpath = dtin%recnpath 813 dtout%recnrec = dtin%recnrec 814 dtout%recptrott = dtin%recptrott 815 dtout%rectesteg = dtin%rectesteg 816 dtout%rcut = dtin%rcut 817 dtout%restartxf = dtin%restartxf 818 dtout%rfasr = dtin%rfasr 819 dtout%rfddk = dtin%rfddk 820 dtout%rfelfd = dtin%rfelfd 821 dtout%rfmagn = dtin%rfmagn 822 dtout%rfmeth = dtin%rfmeth 823 dtout%rfphon = dtin%rfphon 824 dtout%rfstrs = dtin%rfstrs 825 dtout%rfuser = dtin%rfuser 826 dtout%rf2_dkdk = dtin%rf2_dkdk 827 dtout%rf2_dkde = dtin%rf2_dkde 828 dtout%rhoqpmix = dtin%rhoqpmix 829 dtout%signperm = dtin%signperm 830 dtout%slabwsrad = dtin%slabwsrad 831 dtout%slabzbeg = dtin%slabzbeg 832 dtout%slabzend = dtin%slabzend 833 dtout%smdelta = dtin%smdelta 834 dtout%spgaxor = dtin%spgaxor 835 dtout%spgorig = dtin%spgorig 836 dtout%spgroup = dtin%spgroup 837 dtout%spmeth = dtin%spmeth 838 dtout%string_algo = dtin%string_algo 839 dtout%symchi = dtin%symchi 840 dtout%symmorphi = dtin%symmorphi 841 dtout%symsigma = dtin%symsigma 842 dtout%td_mexcit = dtin%td_mexcit 843 dtout%tfkinfunc = dtin%tfkinfunc 844 dtout%timopt = dtin%timopt 845 dtout%use_gemm_nonlop = dtin%use_gemm_nonlop 846 dtout%use_gpu_cuda = dtin%use_gpu_cuda 847 dtout%use_slk = dtin%use_slk 848 dtout%usedmatpu = dtin%usedmatpu 849 dtout%usedmft = dtin%usedmft 850 dtout%useexexch = dtin%useexexch 851 dtout%usefock = dtin%usefock 852 dtout%usekden = dtin%usekden 853 dtout%use_nonscf_gkk = dtin%use_nonscf_gkk 854 dtout%usepaw = dtin%usepaw 855 dtout%usepawu = dtin%usepawu 856 dtout%usepotzero = dtin%usepotzero 857 dtout%userec = dtin%userec 858 dtout%useria = dtin%useria 859 dtout%userib = dtin%userib 860 dtout%useric = dtin%useric 861 dtout%userid = dtin%userid 862 dtout%userie = dtin%userie 863 dtout%usewvl = dtin%usewvl 864 dtout%usexcnhat_orig = dtin%usexcnhat_orig 865 dtout%useylm = dtin%useylm 866 dtout%vacnum = dtin%vacnum 867 dtout%vdw_df_acutmin = dtin%vdw_df_acutmin 868 dtout%vdw_df_aratio = dtin%vdw_df_aratio 869 dtout%vdw_df_damax = dtin%vdw_df_damax 870 dtout%vdw_df_damin = dtin%vdw_df_damin 871 dtout%vdw_df_dcut = dtin%vdw_df_dcut 872 dtout%vdw_df_dratio = dtin%vdw_df_dratio 873 dtout%vdw_df_dsoft = dtin%vdw_df_dsoft 874 dtout%vdw_df_gcut = dtin%vdw_df_gcut 875 dtout%vdw_df_ndpts = dtin%vdw_df_ndpts 876 dtout%vdw_df_ngpts = dtin%vdw_df_ngpts 877 dtout%vdw_df_nqpts = dtin%vdw_df_nqpts 878 dtout%vdw_df_nrpts = dtin%vdw_df_nrpts 879 dtout%vdw_df_nsmooth = dtin%vdw_df_nsmooth 880 dtout%vdw_df_phisoft = dtin%vdw_df_phisoft 881 dtout%vdw_df_qcut = dtin%vdw_df_qcut 882 dtout%vdw_df_qratio = dtin%vdw_df_qratio 883 dtout%vdw_df_rcut = dtin%vdw_df_rcut 884 dtout%vdw_df_rsoft = dtin%vdw_df_rsoft 885 dtout%vdw_df_tolerance = dtin%vdw_df_tolerance 886 dtout%vdw_df_threshold = dtin%vdw_df_threshold 887 dtout%vdw_df_tweaks = dtin%vdw_df_tweaks 888 dtout%vdw_df_zab = dtin%vdw_df_zab 889 dtout%vdw_nfrag = dtin%vdw_nfrag 890 dtout%vdw_xc = dtin%vdw_xc 891 dtout%wfoptalg = dtin%wfoptalg 892 dtout%wvl_bigdft_comp = dtin%wvl_bigdft_comp 893 dtout%w90iniprj = dtin%w90iniprj 894 dtout%w90prtunk = dtin%w90prtunk 895 dtout%xclevel = dtin%xclevel 896 dtout%xc_denpos = dtin%xc_denpos 897 898 !Copy allocated integer arrays from dtin to dtout 899 dtout%bdberry(:) = dtin%bdberry(:) 900 dtout%cd_subset_freq(:) = dtin%cd_subset_freq(:) 901 dtout%d3e_pert1_atpol(:) = dtin%d3e_pert1_atpol(:) 902 dtout%d3e_pert1_dir(:) = dtin%d3e_pert1_dir(:) 903 dtout%d3e_pert2_atpol(:) = dtin%d3e_pert2_atpol(:) 904 dtout%d3e_pert2_dir(:) = dtin%d3e_pert2_dir(:) 905 dtout%d3e_pert3_atpol(:) = dtin%d3e_pert3_atpol(:) 906 dtout%d3e_pert3_dir(:) = dtin%d3e_pert3_dir(:) 907 dtout%ga_rules(:) = dtin%ga_rules(:) 908 dtout%gpu_devices(:) = dtin%gpu_devices(:) 909 dtout%jfielddir(:) = dtin%jfielddir(:) 910 dtout%kptrlatt(:,:) = dtin%kptrlatt(:,:) 911 dtout%kptrlatt_orig = dtin%kptrlatt_orig 912 dtout%qptrlatt(:,:) = dtin%qptrlatt(:,:) 913 dtout%ngfft(:) = dtin%ngfft(:) 914 dtout%ngfftdg(:) = dtin%ngfftdg(:) 915 dtout%nloalg(:) = dtin%nloalg(:) 916 dtout%ngkpt(:) = dtin%ngkpt(:) 917 dtout%qprtrb(:) = dtin%qprtrb(:) 918 dtout%rfatpol(:) = dtin%rfatpol(:) 919 dtout%rfdir(:) = dtin%rfdir(:) 920 dtout%rf2_pert1_dir(:) = dtin%rf2_pert1_dir(:) 921 dtout%rf2_pert2_dir(:) = dtin%rf2_pert2_dir(:) 922 dtout%supercell(:) = dtin%supercell(:) 923 dtout%ucrpa_bands(:) = dtin%ucrpa_bands(:) 924 dtout%vdw_supercell(:) = dtin%vdw_supercell(:) 925 dtout%vdw_typfrag(:) = dtin%vdw_typfrag(:) 926 dtout%wvl_ngauss(:) = dtin%wvl_ngauss(:) 927 928 !Copy reals from dtin to dtout 929 dtout%adpimd_gamma = dtin%adpimd_gamma 930 dtout%boxcutmin = dtin%boxcutmin 931 dtout%bxctmindg = dtin%bxctmindg 932 dtout%cd_halfway_freq = dtin%cd_halfway_freq 933 dtout%cd_max_freq = dtin%cd_max_freq 934 dtout%charge = dtin%charge 935 dtout%cpus = dtin%cpus 936 dtout%ddamp = dtin%ddamp 937 dtout%diecut = dtin%diecut 938 dtout%diegap = dtin%diegap 939 dtout%dielam = dtin%dielam 940 dtout%dielng = dtin%dielng 941 dtout%diemac = dtin%diemac 942 dtout%diemix = dtin%diemix 943 dtout%diemixmag = dtin%diemixmag 944 dtout%dilatmx = dtin%dilatmx 945 dtout%dosdeltae = dtin%dosdeltae 946 dtout%dtion = dtin%dtion 947 dtout%ecut = dtin%ecut 948 dtout%ecuteps = dtin%ecuteps 949 dtout%ecutsigx = dtin%ecutsigx 950 dtout%ecutsm = dtin%ecutsm 951 dtout%ecutwfn = dtin%ecutwfn 952 dtout%effmass_free = dtin%effmass_free 953 dtout%efmas_deg_tol = dtin%efmas_deg_tol 954 dtout%elph2_imagden = dtin%elph2_imagden 955 dtout%eshift = dtin%eshift 956 dtout%esmear = dtin%esmear 957 dtout%exchmix = dtin%exchmix 958 dtout%fband = dtin%fband 959 dtout%focktoldfe = dtin%focktoldfe 960 dtout%friction = dtin%friction 961 dtout%fxcartfactor = dtin%fxcartfactor 962 dtout%ga_opt_percent = dtin%ga_opt_percent 963 dtout%gwls_model_parameter = dtin%gwls_model_parameter 964 dtout%kptnrm = dtin%kptnrm 965 dtout%kptrlen = dtin%kptrlen 966 dtout%maxestep = dtin%maxestep 967 dtout%bmass = dtin%bmass 968 dtout%magcon_lambda = dtin%magcon_lambda 969 dtout%mdwall = dtin%mdwall 970 dtout%mep_mxstep = dtin%mep_mxstep 971 dtout%nelect = dtin%nelect 972 dtout%nnos = dtin%nnos 973 dtout%noseinert = dtin%noseinert 974 dtout%omegasimax = dtin%omegasimax 975 dtout%omegasrdmax = dtin%omegasrdmax 976 dtout%pawecutdg = dtin%pawecutdg 977 dtout%pawovlp = dtin%pawovlp 978 dtout%posocc = dtin%posocc 979 dtout%postoldfe = dtin%postoldfe 980 dtout%postoldff = dtin%postoldff 981 dtout%ppmfrq = dtin%ppmfrq 982 dtout%pw_unbal_thresh = dtin%pw_unbal_thresh 983 dtout%ratsph_extra = dtin%ratsph_extra 984 dtout%recrcut = dtin%recrcut 985 dtout%recefermi = dtin%recefermi 986 dtout%rectolden = dtin%rectolden 987 dtout%dfpt_sciss = dtin%dfpt_sciss 988 dtout%mbpt_sciss = dtin%mbpt_sciss 989 dtout%spinmagntarget = dtin%spinmagntarget 990 dtout%spbroad = dtin%spbroad 991 dtout%spnorbscl = dtin%spnorbscl 992 dtout%stmbias = dtin%stmbias 993 dtout%strfact = dtin%strfact 994 dtout%strprecon = dtin%strprecon 995 dtout%tfw_toldfe = dtin%tfw_toldfe 996 dtout%tl_radius = dtin%tl_radius 997 dtout%tl_nprccg = dtin%tl_nprccg 998 dtout%td_maxene = dtin%td_maxene 999 dtout%toldfe = dtin%toldfe 1000 dtout%tolmxde = dtin%tolmxde 1001 dtout%toldff = dtin%toldff 1002 dtout%tolimg = dtin%tolimg 1003 dtout%tolmxf = dtin%tolmxf 1004 dtout%tolrde = dtin%tolrde 1005 dtout%tolrff = dtin%tolrff 1006 dtout%tolsym = dtin%tolsym 1007 dtout%tolvrs = dtin%tolvrs 1008 dtout%tolwfr = dtin%tolwfr 1009 dtout%tphysel = dtin%tphysel 1010 dtout%tsmear = dtin%tsmear 1011 dtout%ucrpa = dtin%ucrpa 1012 dtout%userra = dtin%userra 1013 dtout%userrb = dtin%userrb 1014 dtout%userrc = dtin%userrc 1015 dtout%userrd = dtin%userrd 1016 dtout%userre = dtin%userre 1017 dtout%vacwidth = dtin%vacwidth 1018 dtout%vdw_tol = dtin%vdw_tol 1019 dtout%vdw_tol_3bt = dtin%vdw_tol_3bt 1020 dtout%vis = dtin%vis 1021 dtout%wfk_task = dtin%wfk_task 1022 dtout%wtq = dtin%wtq 1023 dtout%wvl_hgrid = dtin%wvl_hgrid 1024 dtout%wvl_crmult = dtin%wvl_crmult 1025 dtout%wvl_frmult = dtin%wvl_frmult 1026 dtout%wvl_nprccg = dtin%wvl_nprccg 1027 dtout%xc_tb09_c = dtin%xc_tb09_c 1028 dtout%zcut = dtin%zcut 1029 1030 !Copy allocated real arrays from dtin to dtout 1031 dtout%boxcenter(:) = dtin%boxcenter(:) 1032 dtout%bfield(:) = dtin%bfield(:) 1033 dtout%dfield(:) = dtin%dfield(:) 1034 dtout%efield(:) = dtin%efield(:) 1035 dtout%genafm(:) = dtin%genafm(:) 1036 dtout%goprecprm(:) = dtin%goprecprm(:) 1037 dtout%mdtemp(:) = dtin%mdtemp(:) 1038 dtout%neb_spring(:) = dtin%neb_spring(:) 1039 dtout%polcen(:) = dtin%polcen(:) 1040 dtout%qptn(:) = dtin%qptn(:) 1041 dtout%pvelmax(:) = dtin%pvelmax(:) 1042 dtout%red_efield(:) = dtin%red_efield(:) 1043 dtout%red_dfield(:) = dtin%red_dfield(:) 1044 dtout%red_efieldbar(:) = dtin%red_efieldbar(:) 1045 dtout%shiftk_orig = dtin%shiftk_orig 1046 dtout%strtarget(:) = dtin%strtarget(:) 1047 dtout%ucrpa_window(:) = dtin%ucrpa_window(:) 1048 dtout%vcutgeo(:) = dtin%vcutgeo(:) 1049 dtout%vprtrb(:) = dtin%vprtrb(:) 1050 dtout%zeemanfield(:) = dtin%zeemanfield(:) 1051 1052 !Use alloc_copy to allocate and copy the allocatable arrays 1053 1054 !integer allocatables 1055 call alloc_copy( dtin%algalch, dtout%algalch) 1056 1057 call alloc_copy( dtin%bdgw, dtout%bdgw) 1058 1059 call alloc_copy( dtin%bs_loband, dtout%bs_loband) 1060 1061 call alloc_copy( dtin%dynimage, dtout%dynimage) 1062 1063 call alloc_copy( dtin%efmas_bands, dtout%efmas_bands) 1064 1065 call alloc_copy( dtin%iatfix, dtout%iatfix) 1066 1067 call alloc_copy( dtin%iatsph, dtout%iatsph) 1068 1069 call alloc_copy( dtin%istwfk, dtout%istwfk) 1070 1071 call alloc_copy( dtin%kberry, dtout%kberry) 1072 1073 call alloc_copy( dtin%lexexch, dtout%lexexch) 1074 1075 call alloc_copy( dtin%lpawu, dtout%lpawu) 1076 1077 call alloc_copy( dtin%nband, dtout%nband) 1078 1079 call alloc_copy( dtin%plowan_iatom, dtout%plowan_iatom) 1080 1081 call alloc_copy( dtin%plowan_it, dtout%plowan_it) 1082 1083 call alloc_copy( dtin%plowan_nbl, dtout%plowan_nbl) 1084 1085 call alloc_copy( dtin%plowan_lcalc, dtout%plowan_lcalc) 1086 1087 call alloc_copy( dtin%plowan_projcalc, dtout%plowan_projcalc) 1088 1089 call alloc_copy( dtin%prtatlist, dtout%prtatlist) 1090 1091 call alloc_copy( dtin%so_psp, dtout%so_psp) 1092 1093 call alloc_copy( dtin%symafm, dtout%symafm) 1094 1095 call alloc_copy( dtin%symrel, dtout%symrel) 1096 1097 call alloc_copy( dtin%typat, dtout%typat) 1098 1099 !Allocate and copy real allocatable 1100 call alloc_copy( dtin%acell_orig, dtout%acell_orig) 1101 1102 call alloc_copy( dtin%amu_orig, dtout%amu_orig) 1103 1104 call alloc_copy( dtin%atvshift, dtout%atvshift) 1105 1106 call alloc_copy( dtin%cd_imfrqs, dtout%cd_imfrqs) 1107 1108 call alloc_copy( dtin%chempot, dtout%chempot) 1109 1110 call alloc_copy( dtin%corecs, dtout%corecs) 1111 1112 call alloc_copy( dtin%densty, dtout%densty) 1113 1114 call alloc_copy( dtin%dmatpawu, dtout%dmatpawu) 1115 1116 call alloc_copy( dtin%efmas_dirs, dtout%efmas_dirs) 1117 1118 call alloc_copy( dtin%f4of2_sla, dtout%f4of2_sla) 1119 1120 call alloc_copy( dtin%f6of2_sla, dtout%f6of2_sla) 1121 1122 call alloc_copy( dtin%gw_qlwl, dtout%gw_qlwl) 1123 1124 call alloc_copy( dtin%gw_freqsp, dtout%gw_freqsp) 1125 1126 call alloc_copy( dtin%gwls_list_proj_freq, dtout%gwls_list_proj_freq) 1127 1128 call alloc_copy( dtin%jpawu, dtout%jpawu) 1129 1130 call alloc_copy( dtin%kpt, dtout%kpt) 1131 1132 call alloc_copy( dtin%kptgw, dtout%kptgw) 1133 1134 call alloc_copy( dtin%kptns, dtout%kptns) 1135 1136 call alloc_copy( dtin%mixalch_orig, dtout%mixalch_orig) 1137 1138 call alloc_copy( dtin%nucdipmom, dtout%nucdipmom) 1139 1140 call alloc_copy( dtin%occ_orig, dtout%occ_orig) 1141 1142 call alloc_copy( dtin%pimass, dtout%pimass) 1143 1144 call alloc_copy( dtin%ptcharge, dtout%ptcharge) 1145 1146 call alloc_copy( dtin%qmass, dtout%qmass) 1147 1148 call alloc_copy( dtin%qptdm, dtout%qptdm) 1149 1150 call alloc_copy( dtin%quadmom, dtout%quadmom) 1151 1152 call alloc_copy( dtin%ratsph, dtout%ratsph) 1153 1154 call alloc_copy( dtin%rprim_orig, dtout%rprim_orig) 1155 1156 call alloc_copy( dtin%rprimd_orig, dtout%rprimd_orig) 1157 1158 call alloc_copy( dtin%shiftk, dtout%shiftk) 1159 1160 call alloc_copy( dtin%spinat, dtout%spinat) 1161 1162 call alloc_copy( dtin%tnons, dtout%tnons) 1163 1164 call alloc_copy( dtin%upawu, dtout%upawu) 1165 1166 call alloc_copy( dtin%vel_orig, dtout%vel_orig) 1167 1168 call alloc_copy( dtin%vel_cell_orig, dtout%vel_cell_orig) 1169 1170 call alloc_copy( dtin%wtatcon, dtout%wtatcon) 1171 1172 call alloc_copy( dtin%wtk, dtout%wtk) 1173 1174 call alloc_copy( dtin%xred_orig, dtout%xred_orig) 1175 1176 call alloc_copy( dtin%xredsph_extra, dtout%xredsph_extra) 1177 1178 call alloc_copy( dtin%ziontypat, dtout%ziontypat) 1179 1180 call alloc_copy( dtin%znucl, dtout%znucl) 1181 1182 DBG_EXIT("COLL") 1183 1184 dtout%ndivsm = dtin%ndivsm 1185 dtout%nkpath = dtin%nkpath 1186 dtout%einterp = dtin%einterp 1187 call alloc_copy(dtin%kptbounds, dtout%kptbounds) 1188 1189 end subroutine dtset_copy
m_dtset/dtset_free [ Functions ]
[ Top ] [ m_dtset ] [ Functions ]
NAME
dtset_free
FUNCTION
Free a dataset after use.
SIDE EFFECTS
dtset <type(dataset_type)>=free all allocated allocatable.
PARENTS
calc_vhxc_me,chkinp,dfpt_looppert,driver,gwls_hamiltonian,hybrid_corr m_ab7_invars_f90,m_io_kss,m_kxc,mover_effpot,xchybrid_ncpp_cc
CHILDREN
SOURCE
1212 subroutine dtset_free(dtset) 1213 1214 1215 !This section has been created automatically by the script Abilint (TD). 1216 !Do not modify the following lines by hand. 1217 #undef ABI_FUNC 1218 #define ABI_FUNC 'dtset_free' 1219 !End of the abilint section 1220 1221 implicit none 1222 1223 !Arguments ------------------------------------ 1224 !scalars 1225 type(dataset_type),intent(inout) :: dtset 1226 1227 ! ************************************************************************* 1228 1229 !please, use the same order as the one used in the declaration of the type (see defs_abitypes). 1230 1231 !@dataset_type 1232 !integer allocatable 1233 if (allocated(dtset%algalch)) then 1234 ABI_DEALLOCATE(dtset%algalch) 1235 end if 1236 if (allocated(dtset%bdgw)) then 1237 ABI_DEALLOCATE(dtset%bdgw) 1238 end if 1239 if (allocated(dtset%bs_loband)) then 1240 ABI_DEALLOCATE(dtset%bs_loband) 1241 end if 1242 1243 if (allocated(dtset%dynimage)) then 1244 ABI_DEALLOCATE(dtset%dynimage) 1245 end if 1246 if (allocated(dtset%efmas_bands)) then 1247 ABI_DEALLOCATE(dtset%efmas_bands) 1248 end if 1249 if (allocated(dtset%iatfix)) then 1250 ABI_DEALLOCATE(dtset%iatfix) 1251 end if 1252 if (allocated(dtset%iatsph)) then 1253 ABI_DEALLOCATE(dtset%iatsph) 1254 end if 1255 if (allocated(dtset%istwfk)) then 1256 ABI_DEALLOCATE(dtset%istwfk) 1257 end if 1258 if (allocated(dtset%kberry)) then 1259 ABI_DEALLOCATE(dtset%kberry) 1260 end if 1261 if (allocated(dtset%lexexch)) then 1262 ABI_DEALLOCATE(dtset%lexexch) 1263 end if 1264 if (allocated(dtset%lpawu)) then 1265 ABI_DEALLOCATE(dtset%lpawu) 1266 end if 1267 if (allocated(dtset%nband)) then 1268 ABI_DEALLOCATE(dtset%nband) 1269 end if 1270 if (allocated(dtset%ph_qpath)) then 1271 ABI_DEALLOCATE(dtset%ph_qpath) 1272 end if 1273 if (allocated(dtset%ph_qshift)) then 1274 ABI_DEALLOCATE(dtset%ph_qshift) 1275 end if 1276 if (allocated(dtset%plowan_iatom)) then 1277 ABI_DEALLOCATE(dtset%plowan_iatom) 1278 end if 1279 if (allocated(dtset%plowan_it)) then 1280 ABI_DEALLOCATE(dtset%plowan_it) 1281 end if 1282 if (allocated(dtset%plowan_lcalc)) then 1283 ABI_DEALLOCATE(dtset%plowan_lcalc) 1284 end if 1285 if (allocated(dtset%plowan_nbl)) then 1286 ABI_DEALLOCATE(dtset%plowan_nbl) 1287 end if 1288 if (allocated(dtset%plowan_projcalc)) then 1289 ABI_DEALLOCATE(dtset%plowan_projcalc) 1290 end if 1291 if (allocated(dtset%prtatlist)) then 1292 ABI_DEALLOCATE(dtset%prtatlist) 1293 end if 1294 if (allocated(dtset%so_psp)) then 1295 ABI_DEALLOCATE(dtset%so_psp) 1296 end if 1297 if (allocated(dtset%symafm)) then 1298 ABI_DEALLOCATE(dtset%symafm) 1299 end if 1300 if (allocated(dtset%symrel)) then 1301 ABI_DEALLOCATE(dtset%symrel) 1302 end if 1303 if (allocated(dtset%typat)) then 1304 ABI_DEALLOCATE(dtset%typat) 1305 end if 1306 1307 !real allocatable 1308 if (allocated(dtset%acell_orig)) then 1309 ABI_DEALLOCATE(dtset%acell_orig) 1310 end if 1311 if (allocated(dtset%amu_orig)) then 1312 ABI_DEALLOCATE(dtset%amu_orig) 1313 end if 1314 if (allocated(dtset%atvshift)) then 1315 ABI_DEALLOCATE(dtset%atvshift) 1316 end if 1317 if (allocated(dtset%cd_imfrqs)) then 1318 ABI_DEALLOCATE(dtset%cd_imfrqs) 1319 end if 1320 if (allocated(dtset%chempot)) then 1321 ABI_DEALLOCATE(dtset%chempot) 1322 end if 1323 if (allocated(dtset%corecs)) then 1324 ABI_DEALLOCATE(dtset%corecs) 1325 end if 1326 if (allocated(dtset%densty)) then 1327 ABI_DEALLOCATE(dtset%densty) 1328 end if 1329 if (allocated(dtset%dmatpawu)) then 1330 ABI_DEALLOCATE(dtset%dmatpawu) 1331 end if 1332 if (allocated(dtset%efmas_dirs)) then 1333 ABI_DEALLOCATE(dtset%efmas_dirs) 1334 end if 1335 if (allocated(dtset%gw_qlwl)) then 1336 ABI_DEALLOCATE(dtset%gw_qlwl) 1337 end if 1338 if (allocated(dtset%gw_freqsp)) then 1339 ABI_DEALLOCATE(dtset%gw_freqsp) 1340 end if 1341 if (allocated(dtset%gwls_list_proj_freq)) then 1342 ABI_DEALLOCATE(dtset%gwls_list_proj_freq) 1343 end if 1344 if (allocated(dtset%f4of2_sla)) then 1345 ABI_DEALLOCATE(dtset%f4of2_sla) 1346 end if 1347 if (allocated(dtset%f6of2_sla)) then 1348 ABI_DEALLOCATE(dtset%f6of2_sla) 1349 end if 1350 if (allocated(dtset%jpawu)) then 1351 ABI_DEALLOCATE(dtset%jpawu) 1352 end if 1353 if (allocated(dtset%kpt)) then 1354 ABI_DEALLOCATE(dtset%kpt) 1355 end if 1356 if (allocated(dtset%kptbounds)) then 1357 ABI_DEALLOCATE(dtset%kptbounds) 1358 end if 1359 if (allocated(dtset%kptgw)) then 1360 ABI_DEALLOCATE(dtset%kptgw) 1361 end if 1362 if (allocated(dtset%kptns)) then 1363 ABI_DEALLOCATE(dtset%kptns) 1364 end if 1365 if (allocated(dtset%mixalch_orig)) then 1366 ABI_DEALLOCATE(dtset%mixalch_orig) 1367 end if 1368 if (allocated(dtset%nucdipmom)) then 1369 ABI_DEALLOCATE(dtset%nucdipmom) 1370 end if 1371 if (allocated(dtset%occ_orig)) then 1372 ABI_DEALLOCATE(dtset%occ_orig) 1373 end if 1374 if (allocated(dtset%pimass)) then 1375 ABI_DEALLOCATE(dtset%pimass) 1376 end if 1377 if (allocated(dtset%ptcharge)) then 1378 ABI_DEALLOCATE(dtset%ptcharge) 1379 end if 1380 if (allocated(dtset%qmass)) then 1381 ABI_DEALLOCATE(dtset%qmass) 1382 end if 1383 if (allocated(dtset%qptdm)) then 1384 ABI_DEALLOCATE(dtset%qptdm) 1385 end if 1386 if (allocated(dtset%quadmom)) then 1387 ABI_DEALLOCATE(dtset%quadmom) 1388 end if 1389 if (allocated(dtset%ratsph)) then 1390 ABI_DEALLOCATE(dtset%ratsph) 1391 end if 1392 if (allocated(dtset%rprim_orig)) then 1393 ABI_DEALLOCATE(dtset%rprim_orig) 1394 end if 1395 if (allocated(dtset%rprimd_orig)) then 1396 ABI_DEALLOCATE(dtset%rprimd_orig) 1397 end if 1398 if (allocated(dtset%shiftk)) then 1399 ABI_DEALLOCATE(dtset%shiftk) 1400 end if 1401 if (allocated(dtset%spinat)) then 1402 ABI_DEALLOCATE(dtset%spinat) 1403 end if 1404 if (allocated(dtset%tnons)) then 1405 ABI_DEALLOCATE(dtset%tnons) 1406 end if 1407 if (allocated(dtset%upawu)) then 1408 ABI_DEALLOCATE(dtset%upawu) 1409 end if 1410 if (allocated(dtset%vel_orig)) then 1411 ABI_DEALLOCATE(dtset%vel_orig) 1412 end if 1413 if (allocated(dtset%vel_cell_orig)) then 1414 ABI_DEALLOCATE(dtset%vel_cell_orig) 1415 end if 1416 if (allocated(dtset%wtatcon)) then 1417 ABI_DEALLOCATE(dtset%wtatcon) 1418 end if 1419 if (allocated(dtset%wtk)) then 1420 ABI_DEALLOCATE(dtset%wtk) 1421 end if 1422 if (allocated(dtset%xred_orig)) then 1423 ABI_DEALLOCATE(dtset%xred_orig) 1424 end if 1425 if (allocated(dtset%xredsph_extra)) then 1426 ABI_DEALLOCATE(dtset%xredsph_extra) 1427 end if 1428 if (allocated(dtset%ziontypat)) then 1429 ABI_DEALLOCATE(dtset%ziontypat) 1430 end if 1431 if (allocated(dtset%znucl)) then 1432 ABI_DEALLOCATE(dtset%znucl) 1433 end if 1434 1435 end subroutine dtset_free