TABLE OF CONTENTS
- ABINIT/m_psps
- m_psps/nctab_copy
- m_psps/nctab_eval_tcorespl
- m_psps/nctab_eval_tvalespl
- m_psps/nctab_free
- m_psps/nctab_init
- m_psps/nctab_mixalch
- m_psps/psp2params_copy
- m_psps/psp2params_free
- m_psps/psp2params_init
- m_psps/psps_copy
- m_psps/psps_free
- m_psps/psps_init_from_dtset
- m_psps/psps_init_global
- m_psps/psps_ncwrite
- m_psps/psps_print
- m_psps/test_xml_xmlpaw_upf
ABINIT/m_psps [ Modules ]
NAME
m_psps
FUNCTION
This module provides method to allocate/free/initialize the pseudopotential_type object.
COPYRIGHT
Copyright (C) 2014-2022 ABINIT group (XG,DC,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 .
SOURCE
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 module m_psps 24 25 use defs_basis 26 use m_abicore 27 use m_errors 28 use m_xmpi 29 use m_nctk 30 use m_copy 31 use m_dtset 32 #ifdef HAVE_NETCDF 33 use netcdf 34 #endif 35 36 use m_fstrings, only : itoa, sjoin, yesno, atoi 37 use m_io_tools, only : open_file 38 use m_symtk, only : matr3inv 39 use defs_datatypes, only : pspheader_type, pseudopotential_type, pseudopotential_gth_type, nctab_t 40 use m_paw_numeric, only : paw_spline 41 use m_pawrad, only : pawrad_type, pawrad_init, pawrad_free, simp_gen 42 use m_pawpsp, only : pawpsp_cg 43 use m_parser, only : chkint_eq 44 use m_memeval, only : getdim_nloc, setmqgrid 45 46 implicit none 47 48 private 49 50 ! Helper functions 51 public :: test_xml_xmlpaw_upf ! Test if a pseudo potential file is in XML, XML-PAW or in UPF format. 52 53 public :: psps_init_global ! Allocate and init all part of psps structure that are independent of a given dataset. 54 public :: psps_init_from_dtset ! Allocate and init all part of psps structure that are dependent of a given dataset. 55 public :: psps_free ! Deallocate all memory of psps structure. 56 public :: psps_copy ! Copy the psps structure. 57 public :: psps_print ! Print info on the pseudopotentials. 58 public :: psps_ncwrite ! Write psps data in netcdf format. 59 60 public :: nctab_init ! Create the object. 61 public :: nctab_free ! Free memory. 62 public :: nctab_copy ! Copy the object. 63 public :: nctab_eval_tvalespl ! Evaluate spline-fit of the atomic pseudo valence charge in reciprocal space. 64 public :: nctab_eval_tcorespl ! Evalute spline-fit of the model core charge in reciprocal space. 65 public :: nctab_mixalch ! Mix the pseudopotential tables. Used for alchemical mixing.
m_psps/nctab_copy [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
nctab_copy
FUNCTION
Copy the object.
SOURCE
1174 subroutine nctab_copy(nctabin, nctabout) 1175 1176 !Arguments ------------------------------------ 1177 class(nctab_t),intent(in) :: nctabin 1178 class(nctab_t),intent(out) :: nctabout 1179 1180 ! ************************************************************************* 1181 1182 nctabout%mqgrid_vl = nctabin%mqgrid_vl 1183 nctabout%has_tvale = nctabin%has_tvale 1184 nctabout%has_tcore = nctabin%has_tcore 1185 nctabout%dncdq0 = nctabin%dncdq0 1186 nctabout%d2ncdq0 = nctabin%d2ncdq0 1187 nctabout%dnvdq0 = nctabin%dnvdq0 1188 1189 if (allocated(nctabin%tvalespl)) call alloc_copy(nctabin%tvalespl, nctabout%tvalespl) 1190 if (allocated(nctabin%tcorespl)) call alloc_copy(nctabin%tcorespl, nctabout%tcorespl) 1191 1192 end subroutine nctab_copy
m_psps/nctab_eval_tcorespl [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
nctab_eval_tcorespl
FUNCTION
Evalute spline-fit of the model core charge in reciprocal space.
INPUTS
xcccrc=maximum radius of the pseudo-core charge n1xccc=Number of radial points for the description of the pseudo-core charge (in the framework of the non-linear XC core correction) mqgrid_vl=Number of points in the reciprocal space grid qgrid_vl(mqgrid_vl)=The coordinates of all the points of the radial q-grid xccc1d(n1xccc,6)= The component xccc1d(n1xccc,1) is the pseudo-core charge on the radial grid. The components xccc1d(n1xccc,ideriv) give the ideriv-th derivative of the pseudo-core charge with respect to the radial distance.
SIDE EFFECTS
nctabl%tcorespl(mqgrid_vl,2) nctab%d2ncdq0 nctab%dncdq0
SOURCE
1279 subroutine nctab_eval_tcorespl(nctab, n1xccc, xcccrc, xccc1d, mqgrid_vl, qgrid_vl) 1280 1281 !Arguments ------------------------------------ 1282 !scalars 1283 class(nctab_t),intent(inout) :: nctab 1284 integer,intent(in) :: n1xccc,mqgrid_vl 1285 real(dp),intent(in) :: xcccrc 1286 !arrays 1287 real(dp),intent(in) :: xccc1d(n1xccc,6),qgrid_vl(mqgrid_vl) 1288 1289 !Local variables------------------------------- 1290 real(dp) :: amesh,yp1,ypn 1291 type(pawrad_type) :: core_mesh 1292 1293 ! ************************************************************************* 1294 1295 ABI_CHECK(mqgrid_vl == nctab%mqgrid_vl, "wrong mqgrid_vl") 1296 1297 if (.not. allocated(nctab%tcorespl)) then 1298 ABI_CALLOC(nctab%tcorespl, (mqgrid_vl, 2)) 1299 else 1300 ABI_CHECK(size(nctab%tcorespl, dim=1) == mqgrid_vl, "wrong mqgrid_vl") 1301 end if 1302 1303 ! Skip loop if this atom has no core charge 1304 if (abs(xcccrc) < tol16) then 1305 nctab%has_tcore = .False. 1306 return 1307 end if 1308 1309 nctab%has_tcore = .True. 1310 ! XCCC is given on a linear mesh. 1311 amesh = xcccrc / dble(n1xccc-1) 1312 call pawrad_init(core_mesh, mesh_size=n1xccc, mesh_type=1, rstep=amesh) 1313 1314 ! Compute 4\pi\int[(\frac{\sin(2\pi q r)}{2\pi q r})(r^2 n(r))dr]. 1315 ! write(std_out,*)"xccc1d: amesh, min, max, minloc ",amesh,maxval(xccc1d(:,1)),minval(xccc1d(:,1)),minloc(xccc1d(:,1)) 1316 call pawpsp_cg(nctab%dncdq0, nctab%d2ncdq0, mqgrid_vl, qgrid_vl, nctab%tcorespl(:,1), & 1317 core_mesh, xccc1d(:,1), yp1, ypn) 1318 1319 ! Compute second derivative of tcorespl(q) 1320 call paw_spline(qgrid_vl, nctab%tcorespl(:,1), mqgrid_vl, yp1, ypn, nctab%tcorespl(:,2)) 1321 call pawrad_free(core_mesh) 1322 1323 end subroutine nctab_eval_tcorespl
m_psps/nctab_eval_tvalespl [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
nctab_eval_tvalespl
FUNCTION
Evalute spline-fit of the atomic pseudo valence charge in reciprocal space.
INPUTS
zion=nominal valence of atom as specified in psp file. Used to rescale the f(q=0) component mesh<pawrad_type>Radial mesh (r-space) used for the valence denity. valr(mesh%mesh_size)=Valence density in real space. mqgrid_vl=Number of points in the reciprocal space grid qgrid_vl(mqgrid_vl)=The coordinates of all the points of the radial q-grid
SIDE EFFECTS
nctabl%tvalspl(mqgrid_vl,2) nctab%d2ncdq0
SOURCE
1215 subroutine nctab_eval_tvalespl(nctab, zion, mesh, valr, mqgrid_vl, qgrid_vl) 1216 1217 !Arguments ------------------------------------ 1218 class(nctab_t),intent(inout) :: nctab 1219 integer,intent(in) :: mqgrid_vl 1220 real(dp),intent(in) :: zion 1221 type(pawrad_type),intent(in) :: mesh 1222 !arrays 1223 real(dp),intent(in) :: valr(mesh%mesh_size),qgrid_vl(mqgrid_vl) 1224 1225 !Local variables------------------------------- 1226 real(dp) :: fact,yp1,ypn,d2nvdq0 1227 1228 ! ************************************************************************* 1229 1230 nctab%has_tvale = .True. 1231 if (.not. allocated(nctab%tvalespl)) then 1232 ABI_MALLOC(nctab%tvalespl, (mqgrid_vl, 2)) 1233 else 1234 ABI_CHECK(size(nctab%tvalespl, dim=1) == mqgrid_vl, "wrong mqgrid_vl") 1235 end if 1236 1237 call pawpsp_cg(nctab%dnvdq0, d2nvdq0, mqgrid_vl, qgrid_vl, nctab%tvalespl(:,1), mesh, valr, yp1, ypn) 1238 call simp_gen(yp1, mesh%rad**2 * valr, mesh) 1239 write(std_out,*)" valence charge (before rescaling) integrates to: ",four_pi*yp1 1240 1241 ! Rescale the integral to have the correct number of valence electrons. 1242 ! In some cases, indeed, the radial mesh is not large enough and some valence charge is missing 1243 ! pawpsp_cg extrapolates the integrand beyond rmax but this is not enough. 1244 ! Remember that tvalespl is used to build an initial guess for rhor hence it's very important 1245 ! to have the correct electrostatic. 1246 fact = zion / nctab%tvalespl(1,1) 1247 nctab%tvalespl(:,1) = nctab%tvalespl(:,1) * fact 1248 1249 ! Compute second derivative of tvalespl(q) 1250 call paw_spline(qgrid_vl,nctab%tvalespl(:,1),mqgrid_vl,yp1,ypn,nctab%tvalespl(:,2)) 1251 1252 end subroutine nctab_eval_tvalespl
m_psps/nctab_free [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
nctab_free
FUNCTION
Free memory allocated in nctab_t
SOURCE
1147 subroutine nctab_free(nctab) 1148 1149 !Arguments ------------------------------------ 1150 class(nctab_t),intent(inout) :: nctab 1151 1152 ! ************************************************************************* 1153 1154 ABI_SFREE(nctab%tvalespl) 1155 ABI_SFREE(nctab%tcorespl) 1156 ABI_SFREE(nctab%tphi_qspl) 1157 ABI_SFREE(nctab%tphi_n) 1158 ABI_SFREE(nctab%tphi_l) 1159 ABI_SFREE(nctab%tphi_jtot) 1160 ABI_SFREE(nctab%tphi_occ) 1161 1162 end subroutine nctab_free
m_psps/nctab_init [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
nctab_init
FUNCTION
Create nctab_t.
INPUTS
mqgrid_vl=Number of q-points has_tcore=True if the pseudo has NLCC. has_tvale=True if the atomic valence density is available.
SOURCE
1110 subroutine nctab_init(nctab, mqgrid_vl, has_tcore, has_tvale) 1111 1112 !Arguments ------------------------------------ 1113 class(nctab_t),intent(inout) :: nctab 1114 integer,intent(in) :: mqgrid_vl 1115 logical,intent(in) :: has_tcore, has_tvale 1116 1117 ! ************************************************************************* 1118 1119 nctab%mqgrid_vl = mqgrid_vl 1120 1121 ! The array for the model core charge is always allocated and initialized with zeros. 1122 ! This approach is similar to the one used in the PAW code. 1123 ! has_tcore tells us whether the model core charge is present or not. 1124 nctab%has_tcore = has_tcore 1125 nctab%dncdq0 = zero; nctab%d2ncdq0 = zero 1126 ABI_CALLOC(nctab%tcorespl, (mqgrid_vl, 2)) 1127 1128 ! tvalespl is allocated only if available. 1129 nctab%has_tvale = has_tvale 1130 nctab%dnvdq0 = zero 1131 if (has_tvale) then 1132 ABI_CALLOC(nctab%tvalespl, (mqgrid_vl, 2)) 1133 end if 1134 1135 end subroutine nctab_init
m_psps/nctab_mixalch [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
nctab_mixalch
FUNCTION
Mix the pseudopotential tables. Used for alchemical mixing.
INPUTS
nctabs(npspalch)=NC tables to be mixed npspalch=Number of alchemical pseudos. ntypalch=Number of types of alchemical pseudoatoms algalch(ntypalch)=For each type of pseudo atom, the algorithm to mix the pseudopotentials mixalch(npspalch,ntypalch)=Mixing coefficients to generate alchemical pseudo atoms
OUTPUT
mixtabs(ntypalch)=NC tables describing the alchemical pseudos
SOURCE
1345 subroutine nctab_mixalch(nctabs, npspalch, ntypalch, algalch, mixalch, mixtabs) 1346 1347 !Arguments ------------------------------------ 1348 !scalars 1349 integer,intent(in) :: npspalch,ntypalch 1350 !arrays 1351 integer,intent(in) :: algalch(ntypalch) 1352 real(dp),intent(in) :: mixalch(npspalch, ntypalch) 1353 type(nctab_t),intent(in) :: nctabs(npspalch) 1354 type(nctab_t),target,intent(inout) :: mixtabs(ntypalch) 1355 1356 !Local variables------------------------------- 1357 !scalars 1358 integer :: ipspalch,itypalch 1359 logical :: has_tcore, has_tvale 1360 real(dp) :: mc 1361 type(nctab_t),pointer :: mix 1362 1363 ! ************************************************************************* 1364 1365 ABI_CHECK(all(nctabs(:)%mqgrid_vl == nctabs(1)%mqgrid_vl), "Wrong mqgrid_vl") 1366 ABI_CHECK(all(algalch == 1), "algalch /= 1 not implemented") 1367 1368 do itypalch=1,ntypalch 1369 1370 ! has_tcore is true is at least one pseudo has nlcc. 1371 ! has_tvale is true if *all* mixed pseudos have the PS valence charge. 1372 has_tcore = .False.; has_tvale = .True. 1373 do ipspalch=1,npspalch 1374 if (abs(mixalch(ipspalch,itypalch)) < tol6) cycle 1375 if (nctabs(ipspalch)%has_tcore) has_tcore = .True. 1376 if (.not. nctabs(ipspalch)%has_tvale) has_tvale = .False. 1377 end do 1378 !write(std_out,*)has_tvale, has_tcore 1379 1380 call nctab_free(mixtabs(itypalch)) 1381 call nctab_init(mixtabs(itypalch), nctabs(1)%mqgrid_vl, has_tcore, has_tvale) 1382 mix => mixtabs(itypalch) 1383 1384 do ipspalch=1,npspalch 1385 mc = mixalch(ipspalch,itypalch) 1386 if (abs(mc) < tol6) cycle 1387 ! Linear combination of the quantities 1388 ! Mix core for NLCC 1389 if (has_tcore) then 1390 mix%tcorespl = mix%tcorespl + mc * nctabs(ipspalch)%tcorespl 1391 mix%dncdq0 = mix%dncdq0 + mc * nctabs(ipspalch)%dncdq0 1392 mix%d2ncdq0 = mix%d2ncdq0 + mc * nctabs(ipspalch)%d2ncdq0 1393 end if 1394 ! Mix pseudo valence charge. 1395 if (has_tvale) then 1396 mix%tvalespl = mix%tvalespl + mc * nctabs(ipspalch)%tvalespl 1397 mix%dnvdq0 = mix%dnvdq0 + mc * nctabs(ipspalch)%dnvdq0 1398 end if 1399 end do 1400 1401 end do 1402 1403 end subroutine nctab_mixalch
m_psps/psp2params_copy [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psp2params_copy
FUNCTION
INPUTS
OUTPUT
SOURCE
1030 subroutine psp2params_copy(gth_paramsin, gth_paramsout) 1031 1032 !Arguments ------------------------------------ 1033 class(pseudopotential_gth_type),intent(in) :: gth_paramsin 1034 class(pseudopotential_gth_type),intent(out) :: gth_paramsout 1035 1036 ! ********************************************************************* 1037 1038 if (allocated(gth_paramsin%psppar)) then 1039 call alloc_copy( gth_paramsin%psppar, gth_paramsout%psppar) 1040 end if 1041 if (allocated(gth_paramsin%radii_cf)) then 1042 call alloc_copy( gth_paramsin%radii_cf, gth_paramsout%radii_cf) 1043 end if 1044 if (allocated(gth_paramsin%psp_k_par)) then 1045 call alloc_copy( gth_paramsin%psp_k_par, gth_paramsout%psp_k_par) 1046 end if 1047 if (allocated(gth_paramsin%hasGeometry)) then 1048 call alloc_copy( gth_paramsin%hasGeometry, gth_paramsout%hasGeometry) 1049 end if 1050 if (allocated(gth_paramsin%set)) then 1051 call alloc_copy( gth_paramsin%set, gth_paramsout%set) 1052 end if 1053 1054 end subroutine psp2params_copy
m_psps/psp2params_free [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psp2params_free
FUNCTION
Deallocate a previously allocated data structure for storage of GTH parameters.
INPUTS
SIDE EFFECTS
gth_params <type (pseudopotential_gth_type)>=the values to deallocate.
SOURCE
1073 subroutine psp2params_free(gth_params) 1074 1075 !Arguments ------------------------------------ 1076 !scalars 1077 type(pseudopotential_gth_type),intent(inout) :: gth_params 1078 1079 ! ********************************************************************* 1080 1081 ABI_SFREE(gth_params%set) 1082 ABI_SFREE(gth_params%hasGeometry) 1083 1084 ! Coefficients for local part and projectors 1085 ABI_SFREE(gth_params%psppar) 1086 1087 ! Coefficients for spin orbit part 1088 ABI_SFREE(gth_params%psp_k_par) 1089 1090 ! Different radii 1091 ABI_SFREE(gth_params%radii_cf) 1092 1093 end subroutine psp2params_free
m_psps/psp2params_init [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psp2params_init
FUNCTION
Allocate and initialise the data structure holding parameters for the GTH pseudo-potentials. MJV note: this should be renamed: psp2 suggests it relates to pspcod 2, whereas it is actually 3 the parameters would also be better off separated into C and h arrays
INPUTS
npsp=number of true pseudo used (not alchemy).
OUTPUT
gth_params <type (pseudopotential_gth_type)>=the values to allocate and initialise.
SOURCE
986 subroutine psp2params_init(gth_params, npsp) 987 988 !Arguments ------------------------------------ 989 class(pseudopotential_gth_type),intent(out) :: gth_params 990 integer,intent(in) :: npsp 991 992 ! ********************************************************************* 993 994 !Check array, no params are currently set. 995 ABI_MALLOC(gth_params%set,(npsp)) 996 gth_params%set(:) = .false. 997 998 !Check array, have geometric information been filled? 999 ABI_MALLOC(gth_params%hasGeometry,(npsp)) 1000 gth_params%hasGeometry(:) = .false. 1001 1002 !Coefficients for local part and projectors 1003 ABI_MALLOC(gth_params%psppar,(0:4, 0:6, npsp)) 1004 gth_params%psppar = zero 1005 1006 !Coefficients for spin orbit part 1007 ABI_MALLOC(gth_params%psp_k_par,(1:4, 1:3, npsp)) 1008 gth_params%psp_k_par = zero 1009 1010 !Different radii 1011 ABI_MALLOC(gth_params%radii_cf,(npsp, 3)) 1012 1013 end subroutine psp2params_init
m_psps/psps_copy [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psps_copy
FUNCTION
Copy the psps structure.
SOURCE
536 subroutine psps_copy(pspsin, pspsout) 537 538 !Arguments ------------------------------------ 539 class(pseudopotential_type),intent(in) :: pspsin 540 class(pseudopotential_type),intent(out) :: pspsout 541 542 !Local variables------------------------------- 543 integer :: ii 544 545 ! ************************************************************************* 546 547 ! integer 548 pspsout%dimekb = pspsin%dimekb 549 pspsout%lmnmax = pspsin%lmnmax 550 pspsout%lnmax = pspsin%lnmax 551 pspsout%mproj = pspsin%mproj 552 pspsout%mpsang = pspsin%mpsang 553 pspsout%mpspso = pspsin%mpspso 554 pspsout%mpssoang = pspsin%mpssoang 555 pspsout%mqgrid_ff = pspsin%mqgrid_ff 556 pspsout%mqgrid_vl = pspsin%mqgrid_vl 557 pspsout%mtypalch = pspsin%mtypalch 558 pspsout%npsp = pspsin%npsp 559 pspsout%npspalch = pspsin%npspalch 560 pspsout%ntypat = pspsin%ntypat 561 pspsout%ntypalch = pspsin%ntypalch 562 pspsout%ntyppure = pspsin%ntyppure 563 pspsout%n1xccc = pspsin%n1xccc 564 pspsout%optnlxccc = pspsin%optnlxccc 565 pspsout%positron = pspsin%positron 566 pspsout%usepaw = pspsin%usepaw 567 pspsout%usewvl = pspsin%usewvl 568 pspsout%useylm = pspsin%useylm 569 pspsout%nc_xccc_gspace = pspsin%nc_xccc_gspace 570 571 ! logical 572 pspsout%vlspl_recipSpace = pspsin%vlspl_recipSpace 573 574 ! integer allocatable 575 if (allocated(pspsin%algalch)) call alloc_copy(pspsin%algalch, pspsout%algalch) 576 if (allocated(pspsin%indlmn)) call alloc_copy(pspsin%indlmn, pspsout%indlmn) 577 if (allocated(pspsin%pspdat)) call alloc_copy(pspsin%pspdat, pspsout%pspdat) 578 if (allocated(pspsin%pspcod)) call alloc_copy(pspsin%pspcod, pspsout%pspcod) 579 if (allocated(pspsin%pspso)) call alloc_copy(pspsin%pspso, pspsout%pspso) 580 if (allocated(pspsin%pspxc)) call alloc_copy(pspsin%pspxc, pspsout%pspxc) 581 582 ! real allocatable 583 if (allocated(pspsin%ekb)) call alloc_copy( pspsin%ekb, pspsout%ekb) 584 if (allocated(pspsin%ffspl)) call alloc_copy( pspsin%ffspl, pspsout%ffspl) 585 if (allocated(pspsin%mixalch)) call alloc_copy(pspsin%mixalch, pspsout%mixalch) 586 if (allocated(pspsin%qgrid_ff)) call alloc_copy(pspsin%qgrid_ff, pspsout%qgrid_ff) 587 if (allocated(pspsin%qgrid_vl)) call alloc_copy(pspsin%qgrid_vl, pspsout%qgrid_vl) 588 if (allocated(pspsin%vlspl)) call alloc_copy(pspsin%vlspl, pspsout%vlspl) 589 if (allocated(pspsin%dvlspl)) call alloc_copy(pspsin%dvlspl, pspsout%dvlspl) 590 if (allocated(pspsin%xcccrc)) call alloc_copy(pspsin%xcccrc, pspsout%xcccrc) 591 if (allocated(pspsin%xccc1d)) call alloc_copy(pspsin%xccc1d, pspsout%xccc1d) 592 if (allocated(pspsin%zionpsp)) call alloc_copy(pspsin%zionpsp, pspsout%zionpsp) 593 if (allocated(pspsin%ziontypat)) call alloc_copy(pspsin%ziontypat, pspsout%ziontypat) 594 if (allocated(pspsin%znuclpsp)) call alloc_copy(pspsin%znuclpsp, pspsout%znuclpsp) 595 if (allocated(pspsin%znucltypat)) call alloc_copy(pspsin%znucltypat, pspsout%znucltypat) 596 597 ! allocate and copy character strings 598 ABI_MALLOC(pspsout%filpsp,(pspsout%npsp)) 599 ABI_MALLOC(pspsout%title,(pspsout%npsp)) 600 ABI_MALLOC(pspsout%md5_pseudos,(pspsout%npsp)) 601 do ii=1,pspsout%npsp 602 pspsout%filpsp(ii) = pspsin%filpsp(ii) 603 pspsout%title(ii) = pspsin%title(ii) 604 pspsout%md5_pseudos(ii) = pspsin%md5_pseudos(ii) 605 end do 606 607 ! allocate and copy objects 608 if (allocated(pspsin%nctab)) then 609 ABI_MALLOC(pspsout%nctab,(pspsout%ntypat)) 610 do ii=1,pspsout%ntypat 611 call nctab_copy(pspsin%nctab(ii), pspsout%nctab(ii)) 612 end do 613 end if 614 615 call psp2params_copy(pspsin%gth_params, pspsout%gth_params) 616 617 end subroutine psps_copy
m_psps/psps_free [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psps_free
FUNCTION
Deallocate all memory of psps structure.
SOURCE
478 subroutine psps_free(psps) 479 480 !Arguments ------------------------------------ 481 type(pseudopotential_type),intent(inout) :: psps 482 483 !Local variables------------------------------- 484 integer :: ii 485 486 ! ************************************************************************* 487 488 !Allocation of some arrays independent of the dataset 489 ABI_SFREE(psps%filpsp) 490 ABI_SFREE(psps%pspcod) 491 ABI_SFREE(psps%pspdat) 492 ABI_SFREE(psps%pspso) 493 ABI_SFREE(psps%pspxc) 494 ABI_SFREE(psps%title) 495 ABI_SFREE(psps%zionpsp) 496 ABI_SFREE(psps%znuclpsp) 497 ABI_SFREE(psps%algalch) 498 ABI_SFREE(psps%mixalch) 499 ABI_SFREE(psps%ekb) 500 ABI_SFREE(psps%indlmn) 501 ABI_SFREE(psps%ffspl) 502 ABI_SFREE(psps%qgrid_ff) 503 ABI_SFREE(psps%qgrid_vl) 504 ABI_SFREE(psps%vlspl) 505 ABI_SFREE(psps%dvlspl) 506 ABI_SFREE(psps%xccc1d) 507 ABI_SFREE(psps%xcccrc) 508 ABI_SFREE(psps%ziontypat) 509 ABI_SFREE(psps%znucltypat) 510 ABI_SFREE(psps%md5_pseudos) 511 512 ! Free types. 513 call psp2params_free(psps%gth_params) 514 515 if (allocated(psps%nctab)) then 516 do ii=1,size(psps%nctab) 517 call nctab_free(psps%nctab(ii)) 518 end do 519 ABI_FREE(psps%nctab) 520 end if 521 522 end subroutine psps_free
m_psps/psps_init_from_dtset [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psps_init_from_dtset
FUNCTION
Allocate and initialise all part of psps structure that are dependent of a given dataset.
INPUTS
dtset=<type dataset_type>a given dataset pspheads(npsp)=<type pspheader_type>all the important information from the pseudopotential file header, as well as the psp file name
SIDE EFFECTS
psps=<type pseudopotential_type>the pseudopotentials description
SOURCE
251 subroutine psps_init_from_dtset(psps, dtset, idtset, pspheads) 252 253 !Arguments ------------------------------------ 254 !scalars 255 class(pseudopotential_type),intent(inout) :: psps 256 integer,intent(in) :: idtset 257 type(dataset_type),intent(in) :: dtset 258 !arrays 259 type(pspheader_type),intent(in) :: pspheads(psps%npsp) 260 261 !Local variables------------------------------- 262 !scalars 263 integer,save :: dimekb_old=-1,lmnmax_old=-1,lnmax_old=-1,mqgridff_old=0 264 integer,save :: mqgridvl_old=0,ntypat_old=-1,usepaw_old=-1 265 integer :: ipsp,lmnmax,lmnmaxso,lnmax,lnmaxso,newmqgrid,newmqgriddg,nptsgvec 266 integer :: changed,ii,itypat 267 real(dp) :: gprimd_orig(3,3) 268 269 ! ************************************************************************* 270 271 psps%optnlxccc = dtset%optnlxccc 272 !Determine the number of points needed in reciprocal space to represent the 273 !pseudopotentials (either set by hand from input variable or set automatically by abinit) 274 nptsgvec = 200 !This has to be chosen one and for all or else ?? 275 newmqgrid = dtset%mqgrid 276 newmqgriddg = dtset%mqgriddg 277 278 !JB:Which image to use ? I guess 1 always works 279 call matr3inv(dtset%rprimd_orig(:,:,1),gprimd_orig) 280 if ( dtset%usewvl == 0) then 281 call setmqgrid(newmqgrid,newmqgriddg,dtset%ecut*dtset%dilatmx**2,& 282 & dtset%pawecutdg*dtset%dilatmx**2,gprimd_orig,nptsgvec,psps%usepaw) 283 else 284 call setmqgrid(newmqgrid,newmqgriddg,one,one,gprimd_orig,nptsgvec,psps%usepaw) 285 end if 286 psps%mqgrid_ff = newmqgrid 287 if (psps%usepaw == 1) then 288 psps%mqgrid_vl = newmqgriddg 289 else 290 psps%mqgrid_vl = newmqgrid 291 end if 292 293 !Determine the maximum number of projectors, for the set of pseudo atom 294 call getdim_nloc(lmnmax,lmnmaxso,lnmax,lnmaxso,dtset%mixalch_orig,dtset%nimage,psps%npsp,dtset%npspalch,& 295 & dtset%ntypat,dtset%ntypalch,pspheads) 296 297 psps%npspalch = dtset%npspalch 298 psps%ntypat = dtset%ntypat 299 psps%ntypalch = dtset%ntypalch 300 psps%ntyppure = dtset%ntyppure 301 302 !Set the flag for reciprocal space or real space calculations 303 psps%vlspl_recipSpace = (dtset%icoulomb /= 1) 304 psps%positron = dtset%positron 305 psps%useylm = dtset%useylm 306 psps%usewvl = dtset%usewvl 307 308 ! Define treatment of the model core density for NC pseudos. 309 psps%nc_xccc_gspace = dtset%nc_xccc_gspace 310 311 if (idtset > 1) then 312 ABI_SFREE(psps%algalch) 313 ABI_SFREE(psps%mixalch) 314 end if 315 316 ABI_MALLOC(psps%algalch,(psps%ntypalch)) 317 ABI_MALLOC(psps%mixalch,(psps%npspalch,psps%ntypalch)) 318 psps%algalch(1:psps%ntypalch)=dtset%algalch(1:psps%ntypalch) 319 !This value will be overwritten elsewhere in case there are different images ... 320 psps%mixalch(1:psps%npspalch,1:psps%ntypalch)=dtset%mixalch_orig(1:psps%npspalch,1:psps%ntypalch,1) 321 322 !Set mpspso and psps%pspso 323 !Warning: mpspso might be different for each dataset. 324 ! mpspso not relevant in case of PAW. 325 psps%mpspso=1 326 do ipsp=1,dtset%npsp 327 if(dtset%nspinor==1)then 328 psps%pspso(ipsp)=0 329 ! This is needed to treate SOC perturbatively in sigma. 330 !if (dtset%optdriver == RUNL_SIGMA .and. dtset%so_psp(ipsp) /= 0) then 331 ! ABI_WARNING("Setting pspso to 2 although nspinor == 1") 332 ! psps%pspso(ipsp) = 2 333 !end if 334 335 ! Ideally the following line should not exist, but at present, the space has to be booked 336 if(pspheads(ipsp)%pspso/=0)psps%mpspso=2 337 else if (psps%usepaw==0) then 338 if(dtset%so_psp(ipsp)/=1)then 339 psps%pspso(ipsp)=dtset%so_psp(ipsp) 340 else 341 psps%pspso(ipsp)=pspheads(ipsp)%pspso 342 end if 343 if(psps%pspso(ipsp)/=0)psps%mpspso=2 344 if(pspheads(ipsp)%pspso/=0)psps%mpspso=2 345 else 346 psps%pspso(ipsp)=1+dtset%pawspnorb 347 end if 348 end do 349 350 !Set mpssoang, lmnmax, lnmax 351 if(psps%mpspso==1)then 352 psps%mpssoang=psps%mpsang 353 psps%lmnmax =lmnmax 354 psps%lnmax =lnmax 355 else 356 psps%mpssoang=2*psps%mpsang-1 357 psps%lmnmax=lmnmaxso 358 psps%lnmax=lnmaxso 359 end if 360 361 !T. Rangel: for wvl + paw do not change psps%lmnmax 362 if (psps%useylm==0 .and. psps%usepaw/=1 ) then 363 psps%lmnmax=psps%lnmax 364 end if 365 366 !Set dimekb 367 if (psps%usepaw==0) then 368 psps%dimekb=psps%lnmax 369 else 370 psps%dimekb=psps%lmnmax*(psps%lmnmax+1)/2 371 end if 372 373 !The following arrays are often not deallocated before the end of the dtset loop 374 !and might keep their content from one dataset to the other, if the conditions are fulfilled 375 changed = 0 376 377 if(dimekb_old/=psps%dimekb .or. ntypat_old/=dtset%ntypat .or. usepaw_old/=psps%usepaw) then 378 changed = changed + 1 379 if(idtset/=1) then 380 ABI_SFREE(psps%ekb) 381 end if 382 ABI_MALLOC(psps%ekb,(psps%dimekb,dtset%ntypat*(1-psps%usepaw))) 383 dimekb_old=psps%dimekb 384 end if 385 386 if(lmnmax_old/=psps%lmnmax .or. ntypat_old/=dtset%ntypat)then 387 changed = changed + 1 388 if(idtset/=1) then 389 ABI_SFREE(psps%indlmn) 390 end if 391 ABI_MALLOC(psps%indlmn,(6,psps%lmnmax,dtset%ntypat)) 392 lmnmax_old=psps%lmnmax 393 end if 394 395 if(mqgridff_old/=psps%mqgrid_ff .or. lnmax_old/=psps%lnmax .or. ntypat_old/=dtset%ntypat)then 396 changed = changed + 1 397 if(idtset/=1) then 398 ABI_SFREE(psps%ffspl) 399 ABI_SFREE(psps%qgrid_ff) 400 end if 401 ABI_MALLOC(psps%ffspl,(psps%mqgrid_ff,2,psps%lnmax,dtset%ntypat)) 402 ABI_MALLOC(psps%qgrid_ff,(psps%mqgrid_ff)) 403 mqgridff_old=psps%mqgrid_ff 404 lnmax_old=psps%lnmax 405 end if 406 407 if(mqgridvl_old/=psps%mqgrid_vl .or. ntypat_old/=dtset%ntypat)then 408 changed = changed + 1 409 if(idtset/=1) then 410 ABI_SFREE(psps%qgrid_vl) 411 ABI_SFREE(psps%vlspl) 412 if (allocated(psps%nctab)) then 413 do ii=1,size(psps%nctab) 414 call nctab_free(psps%nctab(ii)) 415 end do 416 ABI_FREE(psps%nctab) 417 end if 418 end if 419 if (idtset/=1 .and. .not.psps%vlspl_recipSpace) then 420 ABI_SFREE(psps%dvlspl) 421 end if 422 423 ABI_MALLOC(psps%qgrid_vl,(psps%mqgrid_vl)) 424 ABI_MALLOC(psps%vlspl,(psps%mqgrid_vl,2,dtset%ntypat)) 425 426 if (psps%usepaw == 0) then 427 ! If you change usepaw in the input, you will get what you deserve! 428 ABI_MALLOC(psps%nctab, (dtset%ntypat)) 429 do itypat=1,dtset%ntypat 430 call nctab_init(psps%nctab(itypat), psps%mqgrid_vl, .False., .False.) 431 end do 432 end if 433 434 if (.not.psps%vlspl_recipSpace) then 435 ABI_MALLOC(psps%dvlspl,(psps%mqgrid_vl,2,dtset%ntypat)) 436 end if 437 mqgridvl_old=psps%mqgrid_vl 438 end if 439 440 if(ntypat_old/=dtset%ntypat.or. usepaw_old/=psps%usepaw)then 441 changed = changed + 1 442 if(idtset/=1) then 443 ABI_SFREE(psps%xccc1d) 444 end if 445 ABI_MALLOC(psps%xccc1d,(psps%n1xccc*(1-psps%usepaw),6,dtset%ntypat)) 446 usepaw_old=psps%usepaw 447 end if 448 449 if(ntypat_old/=dtset%ntypat)then 450 changed = changed + 1 451 if(idtset/=1) then 452 ABI_SFREE(psps%xcccrc) 453 ABI_SFREE(psps%ziontypat) 454 ABI_SFREE(psps%znucltypat) 455 end if 456 ABI_MALLOC(psps%xcccrc,(dtset%ntypat)) 457 ABI_MALLOC(psps%znucltypat,(dtset%ntypat)) 458 ABI_MALLOC(psps%ziontypat,(dtset%ntypat)) 459 ntypat_old=dtset%ntypat 460 end if 461 462 psps%ziontypat(:)=dtset%ziontypat(:) 463 464 end subroutine psps_init_from_dtset
m_psps/psps_init_global [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psps_init_global
FUNCTION
Allocate and initialise all part of psps structure that are independent of a given dataset.
INPUTS
npsp=the number of read pseudo files. pspheads(npsp)=<type pspheader_type>all the important information from the pseudopotential file header, as well as the psp file name
OUTPUT
SIDE EFFECTS
psps=<type pseudopotential_type>the pseudopotentials description
SOURCE
174 subroutine psps_init_global(psps, mtypalch, npsp, pspheads) 175 176 !Arguments ------------------------------------ 177 !scalars 178 class(pseudopotential_type),intent(inout) :: psps 179 integer,intent(in) :: mtypalch,npsp 180 !arrays 181 type(pspheader_type),intent(in) :: pspheads(npsp) 182 183 !Local variables------------------------------- 184 !scalars 185 integer :: ii, mpsang, n1xccc 186 187 ! ************************************************************************* 188 189 !Allocation of some arrays independent of the dataset 190 ABI_MALLOC(psps%filpsp,(npsp)) 191 ABI_MALLOC(psps%pspcod,(npsp)) 192 ABI_MALLOC(psps%pspdat,(npsp)) 193 ABI_MALLOC(psps%pspso,(npsp)) 194 ABI_MALLOC(psps%pspxc,(npsp)) 195 ABI_MALLOC(psps%title,(npsp)) 196 ABI_MALLOC(psps%zionpsp,(npsp)) 197 ABI_MALLOC(psps%znuclpsp,(npsp)) 198 call psp2params_init(psps%gth_params, npsp) 199 200 psps%filpsp(1:npsp)=pspheads(1:npsp)%filpsp 201 psps%pspcod(1:npsp)=pspheads(1:npsp)%pspcod 202 psps%pspdat(1:npsp)=pspheads(1:npsp)%pspdat 203 psps%pspso(1:npsp)=pspheads(1:npsp)%pspso 204 psps%pspxc(1:npsp)=pspheads(1:npsp)%pspxc 205 psps%title(1:npsp)=pspheads(1:npsp)%title 206 psps%zionpsp(1:npsp)=pspheads(1:npsp)%zionpsp 207 psps%znuclpsp(1:npsp)=pspheads(1:npsp)%znuclpsp 208 209 ! Transfer md5 checksum 210 ABI_MALLOC(psps%md5_pseudos, (npsp)) 211 psps%md5_pseudos = pspheads(1:npsp)%md5_checksum 212 !Set values independant from dtset 213 psps%npsp = npsp 214 !Note that mpsang is the max of 1+lmax, with minimal value 1 (even for local psps, at present) 215 mpsang=1 216 n1xccc=pspheads(1)%xccc 217 do ii=1,psps%npsp 218 mpsang=max(pspheads(ii)%lmax+1,mpsang) 219 n1xccc=max(pspheads(ii)%xccc,n1xccc) 220 end do 221 psps%mpsang = mpsang 222 psps%n1xccc = n1xccc 223 ! Determine here whether the calculation is PAW 224 ! If paw, all pspcod necessarily are 7 or 17 (see iofn2) 225 psps%usepaw =0 226 if (pspheads(1)%pspcod==7.or.pspheads(1)%pspcod==17) psps%usepaw=1 227 psps%mtypalch = mtypalch 228 229 end subroutine psps_init_global
m_psps/psps_ncwrite [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psps_ncwrite
FUNCTION
Writes on file the most important arrays defined in the pseudopotential_type for futher post-processing. This function should be called by master node only.
INPUTS
path=File name.
SOURCE
841 subroutine psps_ncwrite(psps, path) 842 843 !Arguments ------------------------------------ 844 class(pseudopotential_type),intent(in) :: psps 845 character(len=*),intent(in) :: path 846 847 !Local variables------------------------------- 848 integer :: ipsp,itypat,ncid,ncerr 849 850 ! ************************************************************************* 851 852 #ifdef HAVE_NETCDF 853 NCF_CHECK(nctk_open_create(ncid, path, xmpi_comm_self)) 854 855 ! Define dimensions 856 ncerr = nctk_def_dims(ncid, [ & 857 nctkdim_t("fnlen", fnlen + 1), & 858 nctkdim_t("md5_slen", md5_slen + 1), & 859 nctkdim_t("ntypat", psps%ntypat), & 860 nctkdim_t("npsp", psps%npsp), & 861 nctkdim_t("lnmax", psps%lnmax), & 862 nctkdim_t("lmnmax", psps%lnmax), & 863 nctkdim_t("dimekb", psps%dimekb), & 864 nctkdim_t("mqgrid_vl", psps%mqgrid_vl), & 865 nctkdim_t("mqgrid_ff", psps%mqgrid_ff) & 866 ]) 867 NCF_CHECK(ncerr) 868 869 if (psps%n1xccc /= 0) then ! 0 means unlimited! 870 NCF_CHECK(nctk_def_dims(ncid, nctkdim_t("n1xccc", psps%n1xccc))) 871 end if 872 873 ! Define variables 874 ncerr = nctk_def_iscalars(ncid, [character(len=nctk_slen) :: "usepaw", "useylm"]) 875 NCF_CHECK(ncerr) 876 877 ncerr = nctk_def_arrays(ncid, [& 878 nctkarr_t("ziontypat", "dp", "ntypat"), & 879 nctkarr_t("znucltypat", "dp", "ntypat"), & 880 nctkarr_t("qgrid_vl", "dp", "mqgrid_vl"), & 881 nctkarr_t("qgrid_ff", "dp", "mqgrid_ff"), & 882 nctkarr_t("vlspl", "dp", "mqgrid_vl, two, ntypat"), & 883 nctkarr_t("xcccrc", "dp", "ntypat"), & 884 nctkarr_t("indlmn", "int", "six, lmnmax, ntypat"), & 885 nctkarr_t("ffspl", "dp", "mqgrid_ff, two, lnmax, ntypat"), & 886 nctkarr_t("filpsp", "char", "fnlen, npsp"), & 887 nctkarr_t("md5_pseudos", "char", "md5_slen, npsp") & 888 ]) 889 NCF_CHECK(ncerr) 890 891 if (psps%usepaw == 0) then 892 NCF_CHECK(nctk_def_arrays(ncid, nctkarr_t("ekb", "dp", "dimekb, ntypat"))) 893 if (psps%n1xccc /= 0) then 894 NCF_CHECK(nctk_def_arrays(ncid, nctkarr_t("xccc1d", "dp", "n1xccc, six, ntypat"))) 895 end if 896 ncerr = nctk_def_arrays(ncid, [& 897 nctkarr_t("nc_tvalespl", "dp", "mqgrid_vl, two, ntypat"), & 898 nctkarr_t("nc_tcorespl", "dp", "mqgrid_vl, two, ntypat") & 899 ]) 900 NCF_CHECK(ncerr) 901 end if 902 903 ! Write data 904 NCF_CHECK(nf90_enddef(ncid)) 905 NCF_CHECK(nf90_put_var(ncid, vid("usepaw"), psps%usepaw)) 906 NCF_CHECK(nf90_put_var(ncid, vid("useylm"), psps%useylm)) 907 NCF_CHECK(nf90_put_var(ncid, vid("ziontypat"), psps%ziontypat)) 908 NCF_CHECK(nf90_put_var(ncid, vid("znucltypat"), psps%znucltypat)) 909 do ipsp=1,psps%npsp 910 NCF_CHECK(nf90_put_var(ncid, vid("filpsp"), trim(psps%filpsp(ipsp)), start=[1, ipsp])) 911 NCF_CHECK(nf90_put_var(ncid, vid("md5_pseudos"), trim(psps%md5_pseudos(ipsp)), start=[1, ipsp])) 912 end do 913 NCF_CHECK(nf90_put_var(ncid, vid("qgrid_vl"), psps%qgrid_vl)) 914 NCF_CHECK(nf90_put_var(ncid, vid("qgrid_ff"), psps%qgrid_ff)) 915 NCF_CHECK(nf90_put_var(ncid, vid("indlmn"), psps%indlmn)) 916 917 ! Local part in q-space and second derivative 918 if (allocated(psps%vlspl)) then 919 NCF_CHECK(nf90_put_var(ncid, vid("vlspl"), psps%vlspl)) 920 end if 921 922 ! Form factors for each type of atom 923 ! for each type and each (l,n) channel, ffnl(q) and second derivative 924 if (allocated(psps%ffspl)) then 925 NCF_CHECK(nf90_put_var(ncid, vid("ffspl"), psps%ffspl)) 926 end if 927 928 ! Pseudo-core charge for each type of atom, on the real-space radial 929 NCF_CHECK(nf90_put_var(ncid, vid("xcccrc"), psps%xcccrc)) 930 if (psps%usepaw == 0 .and. allocated(psps%xccc1d)) then 931 NCF_CHECK(nf90_put_var(ncid, vid("xccc1d"), psps%xccc1d)) 932 end if 933 934 ! NC-only: add tcore_spl and tvalespl in q-space 935 if (psps%usepaw == 0) then 936 NCF_CHECK(nf90_put_var(ncid, vid("ekb"), psps%ekb)) 937 do itypat=1,psps%ntypat 938 if (psps%nctab(itypat)%has_tvale) then 939 ncerr = nf90_put_var(ncid, vid("nc_tvalespl"), psps%nctab(itypat)%tvalespl, start=[1,1,itypat]) 940 NCF_CHECK(ncerr) 941 end if 942 if (psps%nctab(itypat)%has_tcore) then 943 ncerr = nf90_put_var(ncid, vid("nc_tcorespl"), psps%nctab(itypat)%tcorespl, start=[1,1,itypat]) 944 NCF_CHECK(ncerr) 945 end if 946 end do 947 end if 948 949 NCF_CHECK(nf90_close(ncid)) 950 951 #else 952 ABI_WARNING("netcdf support not activated. psps file cannot be created!") 953 #endif 954 955 contains 956 integer function vid(vname) 957 character(len=*),intent(in) :: vname 958 vid = nctk_idname(ncid, vname) 959 end function vid 960 961 end subroutine psps_ncwrite
m_psps/psps_print [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
psps_print
FUNCTION
Print the content of a pseudopotential_type derived type
INPUTS
psps=<type pseudopotential_type>=Info on the pseudopotentials. unit(optional)=unit number for output prtvol(optional)=verbosity level mode_paral(optional): either "COLL" or "PERS"
OUTPUT
Only writing
SOURCE
640 subroutine psps_print(psps, unit, prtvol, mode_paral) 641 642 !Arguments ------------------------------------ 643 !scalars 644 class(pseudopotential_type),intent(in) :: psps 645 integer,intent(in),optional :: prtvol,unit 646 character(len=4),intent(in),optional :: mode_paral 647 648 !Local variables------------------------------- 649 !scalars 650 integer :: ierr,ips,ipsp_alch,ityp_alch,itypat,unt,my_prtvol 651 character(len=4) :: mode 652 character(len=500) :: msg 653 !arrays 654 integer :: cond_values(4) 655 character(len=9) :: cond_string(4) 656 657 ! ************************************************************************* 658 659 ! Provide defaults 660 my_prtvol=0; if (present(prtvol)) my_prtvol=prtvol 661 unt=std_out; if (present(unit)) unt=unit 662 mode='COLL'; if (present(mode_paral)) mode=mode_paral 663 ierr=0; cond_string(1:4)=' '; cond_values(:)=0 664 665 ! General info including spin-orbit 666 call wrtout(unt,' ==== Info on pseudopotentials ==== ', mode) 667 668 SELECT CASE (psps%usepaw) 669 CASE (0) 670 call wrtout(unt,' Norm-conserving pseudopotentials ', mode) 671 !call wrtout(unt, sjoin(' Max number of Kleinman-Bylander energies ', itoa(psps%dimekb)), mode) 672 !do itypat=1,psps%ntypat 673 ! write(msg,'(a,i4,a,f9.4)')' Type ',itypat,' K-B energies ',(psps%ekb(ikbe,itypat),ikbe=1,psps%dimekb) 674 !end do 675 CASE (1) 676 write(msg,'(a)') 677 call wrtout(unt,' PAW calculation', mode) 678 !call wrtout(unt,sjoin(' Max number of D_ij coefficients ', itoa(psps%dimekb)), mode) 679 CASE DEFAULT 680 call chkint_eq(0,0,cond_string,cond_values,ierr,'usepaw',psps%usepaw,2,[0,1],unt) 681 END SELECT 682 683 !SELECT CASE (psps%positron) 684 !CASE (0) 685 ! call wrtout(unt, ' Standard Electron Calculation ', mode) 686 !CASE (1,2) 687 ! write(msg,'(a,i0)')' Positron Calculation with positron .. ',psps%positron 688 ! call wrtout(unt,msg,mode) 689 !CASE DEFAULT 690 ! call chkint_eq(0,0,cond_string,cond_values,ierr,'positron',psps%positron,3,[0,1,2],unt) 691 !END SELECT 692 693 write(msg,'(a,i4,2a,i4)')& 694 ' Number of pseudopotentials .. ',psps%npsp,ch10,& 695 ' Number of types of atoms .. ',psps%ntypat 696 call wrtout(unt,msg,mode) 697 698 if (psps%usepaw==0) then 699 SELECT CASE (psps%mpspso) 700 CASE (1) 701 call wrtout(unt,' Scalar calculation (no spin-orbit term) ',mode) 702 CASE (2) 703 write(msg,'(3a,i3)')& 704 ' Calculation with spin-orbit coupling ',ch10,& 705 ' Max number of channels (spin-orbit included) ',psps%mpssoang 706 call wrtout(unt,msg,mode) 707 do itypat=1,psps%ntypat 708 if (psps%pspso(itypat) /= 1) then 709 write(msg,'(a,i4,a,i2,a)')& 710 ' - Atom type ',itypat,' has spin-orbit characteristics (pspso= ',psps%pspso(itypat),")" 711 call wrtout(unt,msg,mode) 712 end if 713 end do 714 CASE DEFAULT 715 call chkint_eq(0,0,cond_string,cond_values,ierr,'mpspso',psps%mpspso,2,[1,2],unt) 716 END SELECT 717 else 718 SELECT CASE (maxval(psps%pspso)) 719 CASE (0,1) 720 msg=' Scalar calculation (no spin-orbit term) ' 721 CASE (2) 722 msg=' Calculation with spin-orbit coupling ' 723 END SELECT 724 call wrtout(unt,msg,mode) 725 end if 726 727 ! Info on nonlocal part 728 SELECT CASE (psps%useylm) 729 CASE (0) 730 msg = ' Nonlocal part applied using Legendre polynomials ' 731 CASE (1) 732 msg = ' Nonlocal part applied using real spherical harmonics ' 733 CASE DEFAULT 734 call chkint_eq(0,0,cond_string,cond_values,ierr,'psps%useylm',psps%useylm,2,(/0,1/),unt) 735 END SELECT 736 call wrtout(unt,msg,mode) 737 738 write(msg,'(a,i3)')' Max number of non-local projectors over l and type ',psps%mproj 739 call wrtout(unt,msg,mode) 740 741 write(msg,'(a,i3,2a,i3,2a,i3)')& 742 ' Highest angular momentum +1 ....... ',psps%mpsang,ch10,& 743 ' Max number of (l,n) components .. ',psps%lnmax, ch10,& 744 ' Max number of (l,m,n) components .. ',psps%lmnmax 745 call wrtout(unt,msg,mode) 746 747 !FIXME for paw n1xccc==1 748 ! Non-linear Core correction 749 if (psps%n1xccc/=0) then 750 write(msg,'(3a,2(a,i4,a),2a)')ch10,& 751 ' Pseudo-Core Charge Info: ',ch10,& 752 ' Number of radial points for pseudo-core charge .. ',psps%n1xccc,ch10,& 753 ' XC core-correction treatment (optnlxccc) ........ ',psps%optnlxccc,ch10,& 754 ' Radius for pseudo-core charge for each type ..... ',ch10 755 call wrtout(unt,msg,mode) 756 do itypat=1,psps%ntypat 757 write(msg,'(a,i4,a,f7.4)')' - Atom type ',itypat,' has pseudo-core radius .. ',psps%xcccrc(itypat) 758 call wrtout(unt,msg,mode) 759 end do 760 end if 761 762 ! Alchemical mixing 763 if (psps%mtypalch/=0) then 764 write(msg,'(3a,3(a,i4,a))')ch10,& 765 ' Calculation with alchemical mixing:',ch10,& 766 ' Number of pure pseudoatoms .... ',psps%ntyppure,ch10,& 767 ' Number of pseudos for mixing .. ',psps%npspalch,ch10,& 768 ' Alchemical pseudoatoms ........ ',psps%ntypalch,ch10 769 call wrtout(unt,msg,mode) 770 do ipsp_alch=1,psps%npspalch 771 do ityp_alch=1,psps%ntypalch 772 write(std_out,*)' mixalch ',psps%mixalch(ipsp_alch,ityp_alch) 773 end do 774 end do 775 do ityp_alch=1,psps%ntypalch 776 write(msg,'(a,i4,a,i4)')' For alchemical atom no. ',ityp_alch,' algalch is .. ',psps%algalch(ityp_alch) 777 call wrtout(unt,msg,mode) 778 end do 779 end if 780 781 ! Info in Q-grid for spline of form factors 782 write(msg,'(3a,a,i6,a,a,i6)')ch10,& 783 ' Info on the Q-grid used for form factors in spline form: ',ch10,& 784 ' Number of q-points for radial functions ffspl .. ',psps%mqgrid_ff,ch10,& 785 ' Number of q-points for vlspl ................... ',psps%mqgrid_vl 786 call wrtout(unt,msg,mode) 787 788 if (psps%vlspl_recipSpace) then 789 call wrtout(unt,' vloc is computed in Reciprocal Space ',mode) 790 else 791 call wrtout(unt,' vloc is computed in Real Space ',mode) 792 end if 793 if (psps%usepaw == 0) then 794 if (psps%nc_xccc_gspace == 0) call wrtout(unt,' model core charge treated in real-space', mode) 795 if (psps%nc_xccc_gspace == 1) call wrtout(unt,' model core charge treated in G-space', mode) 796 end if 797 798 !TODO additional stuff that might be printed 799 call wrtout(unt, "", mode) 800 do itypat=1,psps%ntypat 801 write(msg,'(a,i0,a,i0)')' XC functional for type ',itypat,' is ',psps%pspxc(itypat) 802 call wrtout(unt,msg,mode) 803 !write(std_out,*)psps%ziontypat(itypat),psps%znucltypat(itypat) 804 if (psps%usepaw == 0) then 805 call wrtout(unt, sjoin(" Pseudo valence available: ", yesno(psps%nctab(itypat)%has_tvale)), mode) 806 end if 807 end do 808 809 !integer, pointer :: pspxc(:) 810 ! pspxc(ntypat) 811 ! For each type of psp, the XC functional that was used to generate it, as given by the psp file 812 if (my_prtvol>=3) then 813 do ips=1,psps%npsp 814 write(std_out,*)' Pseudo number ',ips,' read from ',trim(psps%filpsp(ips)) 815 write(std_out,*)' Format or code ',psps%pspcod(ips) 816 write(std_out,*)' Generation date ',psps%pspdat(ips) 817 write(std_out,*)' Content of first line: ', trim(psps%title(ips)) 818 end do 819 end if 820 821 call wrtout(unt, "", mode) 822 823 end subroutine psps_print
m_psps/test_xml_xmlpaw_upf [ Functions ]
[ Top ] [ m_psps ] [ Functions ]
NAME
test_xml_xmlpaw_upf
FUNCTION
Test if a pseudo potential file is in XML, XML-PAW or in UPF format.
INPUTS
path=Pseudopotential file
OUTPUT
usexml=1 if XML file xmlpaw=1 if PAW file in XML format useupf=1 or 2 if UPF file.
SOURCE
87 subroutine test_xml_xmlpaw_upf(path, usexml, xmlpaw, useupf) 88 89 !Arguments ------------------------------------ 90 !scalars 91 character(len=*),intent(in) :: path 92 integer,intent(out) :: usexml, xmlpaw, useupf 93 94 !Local variables------------------------------- 95 !scalars 96 integer :: temp_unit, ii 97 character(len=500) :: msg,errmsg 98 character(len=70) :: testxml 99 100 ! ************************************************************************* 101 102 ! Check if the file pseudopotential file is written in XML 103 usexml = 0; xmlpaw = 0; useupf = 0 104 105 if (open_file(path,msg,newunit=temp_unit,form='formatted',status='old') /= 0) then 106 ABI_ERROR(msg) 107 end if 108 rewind (unit=temp_unit,err=10,iomsg=errmsg) 109 110 read(temp_unit, "(a)",err=10,iomsg=errmsg) testxml 111 if(testxml(1:5)=='<?xml')then 112 usexml = 1 113 read(temp_unit,*,err=10,iomsg=errmsg) testxml 114 if(testxml(1:4)=='<paw') xmlpaw = 1 115 else 116 usexml = 0 117 if (testxml(1:4) == '<UPF') then 118 ! Make sure this is not UPF version >= 2 119 ! "<UPF version="2.0.1"> 120 ii = index(testxml, '"') 121 if (ii /= 0) then 122 useupf = atoi(testxml(ii+1:ii+1)) 123 !if (useupf >= 2) then 124 ! ABI_ERROR(sjoin("UPF version >= 2 is not supported by Abinit. Use psp8 or psml format.", ch10, "Pseudo:", path)) 125 !end if 126 else 127 ABI_ERROR(sjoin("Cannot find version attributed in UPF file:", path)) 128 end if 129 130 end if 131 end if 132 133 ! Check if pseudopotential file is a Q-espresso UPF1 file 134 if (useupf == 0) then 135 rewind (unit=temp_unit,err=10,iomsg=errmsg) 136 read(temp_unit,*,err=10,iomsg=errmsg) testxml ! just a string, no relation to xml. 137 if(testxml(1:9)=='<PP_INFO>')then 138 useupf = 1 139 else 140 useupf = 0 141 end if 142 end if 143 144 close(unit=temp_unit,err=10,iomsg=errmsg) 145 146 return 147 148 ! Handle IO error 149 10 continue 150 ABI_ERROR(errmsg) 151 152 end subroutine test_xml_xmlpaw_upf