TABLE OF CONTENTS


ABINIT/m_psps [ Modules ]

[ Top ] [ Modules ]

NAME

  m_psps

FUNCTION

  This module provides method to allocate/free/initialize the
  pseudopotential_type object.

COPYRIGHT

  Copyright (C) 2014-2018 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 .

PARENTS

CHILDREN

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 module m_psps
29 
30  use defs_basis
31  use m_abicore
32  use m_errors
33  use m_xmpi
34  use m_nctk
35  use m_copy
36 #ifdef HAVE_NETCDF
37  use netcdf
38 #endif
39 
40  use m_fstrings,      only : itoa, sjoin, yesno
41  use m_io_tools,      only : open_file
42  use m_symtk,         only : matr3inv
43  use defs_datatypes,  only : pspheader_type, pseudopotential_type, pseudopotential_gth_type, nctab_t
44  use defs_abitypes,   only : dataset_type
45  use m_paw_numeric,   only : paw_spline
46  use m_pawrad,        only : pawrad_type, pawrad_init, pawrad_free, simp_gen
47  use m_pawpsp,        only : pawpsp_cg
48  use m_parser,        only : chkint_eq
49  use m_memeval,       only : getdim_nloc, setmqgrid
50 
51  implicit none
52 
53  private
54 
55  ! Helper functions
56  public :: test_xml_xmlpaw_upf     ! Test if a pseudo potential file is in XML, XML-PAW or in UPF format.
57 
58  public :: psps_init_global        ! Allocate and init all part of psps structure that are independent of a given dataset.
59  public :: psps_init_from_dtset    ! Allocate and init all part of psps structure that are dependent of a given dataset.
60  public :: psps_free               ! Deallocate all memory of psps structure.
61  public :: psps_copy               ! Copy the psps structure.
62  public :: psps_print              ! Print info on the pseudopotentials.
63  public :: psps_ncwrite            ! Write psps data in netcdf format.
64 
65  public :: nctab_init              ! Create the object.
66  public :: nctab_free              ! Free memory.
67  public :: nctab_copy              ! Copy the object.
68  public :: nctab_eval_tvalespl     ! Evaluate spline-fit of the atomic pseudo valence charge in reciprocal space.
69  public :: nctab_eval_tcorespl     ! Evalute spline-fit of the model core charge in reciprocal space.
70  public :: nctab_mixalch           ! Mix the pseudopotential tables. Used for alchemical mixing.

m_psps/nctab_copy [ Functions ]

[ Top ] [ m_psps ] [ Functions ]

NAME

  nctab_copy

FUNCTION

PARENTS

      m_psps

CHILDREN

      nctab_free,nctab_init

SOURCE

1498 subroutine nctab_copy(nctabin, nctabout)
1499 
1500 
1501 !This section has been created automatically by the script Abilint (TD).
1502 !Do not modify the following lines by hand.
1503 #undef ABI_FUNC
1504 #define ABI_FUNC 'nctab_copy'
1505 !End of the abilint section
1506 
1507  implicit none
1508 
1509 !Arguments ------------------------------------
1510 !scalars
1511  type(nctab_t),intent(in) :: nctabin
1512  type(nctab_t),intent(out) :: nctabout
1513 
1514 ! *************************************************************************
1515 
1516  nctabout%mqgrid_vl  = nctabin%mqgrid_vl
1517  nctabout%has_tvale  = nctabin%has_tvale
1518  nctabout%has_tcore  = nctabin%has_tcore
1519  nctabout%dncdq0     = nctabin%dncdq0
1520  nctabout%d2ncdq0    = nctabin%d2ncdq0
1521  nctabout%dnvdq0     = nctabin%dnvdq0
1522 
1523  if (allocated(nctabin%tvalespl)) then
1524    call alloc_copy( nctabin%tvalespl, nctabout%tvalespl)
1525  end if
1526  if (allocated(nctabin%tcorespl)) then
1527    call alloc_copy( nctabin%tcorespl, nctabout%tcorespl)
1528  end if
1529 
1530 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

PARENTS

      pspatm

CHILDREN

      nctab_free,nctab_init

SOURCE

1639 subroutine nctab_eval_tcorespl(nctab, n1xccc, xcccrc, xccc1d, mqgrid_vl, qgrid_vl)
1640 
1641 
1642 !This section has been created automatically by the script Abilint (TD).
1643 !Do not modify the following lines by hand.
1644 #undef ABI_FUNC
1645 #define ABI_FUNC 'nctab_eval_tcorespl'
1646 !End of the abilint section
1647 
1648  implicit none
1649 
1650 !Arguments ------------------------------------
1651 !scalars
1652  integer,intent(in) :: n1xccc,mqgrid_vl
1653  real(dp),intent(in) :: xcccrc
1654  type(nctab_t),intent(inout) :: nctab
1655 !arrays
1656  real(dp),intent(in) :: xccc1d(n1xccc,6),qgrid_vl(mqgrid_vl)
1657 
1658 !Local variables-------------------------------
1659 !scalars
1660  real(dp) :: amesh,yp1,ypn
1661  type(pawrad_type) :: core_mesh
1662 
1663 ! *************************************************************************
1664 
1665  ABI_CHECK(mqgrid_vl == nctab%mqgrid_vl, "wrong mqgrid_vl")
1666 
1667  if (.not. allocated(nctab%tcorespl)) then
1668    ABI_CALLOC(nctab%tcorespl, (mqgrid_vl, 2))
1669  else
1670    ABI_CHECK(size(nctab%tcorespl, dim=1) == mqgrid_vl, "wrong mqgrid_vl")
1671  end if
1672 
1673  ! Skip loop if this atom has no core charge
1674  if (abs(xcccrc) < tol16) then
1675    nctab%has_tcore = .False.
1676    return
1677  end if
1678 
1679  nctab%has_tcore = .True.
1680  ! XCCC is given on a linear mesh.
1681  amesh = xcccrc / dble(n1xccc-1)
1682  call pawrad_init(core_mesh, mesh_size=n1xccc, mesh_type=1, rstep=amesh)
1683 
1684  ! Compute 4\pi\int[(\frac{\sin(2\pi q r)}{2\pi q r})(r^2 n(r))dr].
1685  !write(std_out,*)"xccc1d: amesh, min, max, minloc ",amesh,maxval(xccc1d(:,1)),minval(xccc1d(:,1)),minloc(xccc1d(:,1))
1686  call pawpsp_cg(nctab%dncdq0, nctab%d2ncdq0, mqgrid_vl, qgrid_vl, nctab%tcorespl(:,1), &
1687 &               core_mesh, xccc1d(:,1), yp1, ypn)
1688 
1689  ! Compute second derivative of tcorespl(q)
1690  call paw_spline(qgrid_vl, nctab%tcorespl(:,1), mqgrid_vl, yp1, ypn, nctab%tcorespl(:,2))
1691  call pawrad_free(core_mesh)
1692 
1693 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

PARENTS

      psp8in,psp9in

CHILDREN

      nctab_free,nctab_init

SOURCE

1559 subroutine nctab_eval_tvalespl(nctab, zion, mesh, valr, mqgrid_vl, qgrid_vl)
1560 
1561 
1562 !This section has been created automatically by the script Abilint (TD).
1563 !Do not modify the following lines by hand.
1564 #undef ABI_FUNC
1565 #define ABI_FUNC 'nctab_eval_tvalespl'
1566 !End of the abilint section
1567 
1568  implicit none
1569 
1570 !Arguments ------------------------------------
1571 !scalars
1572  integer,intent(in) :: mqgrid_vl
1573  real(dp),intent(in) :: zion
1574  type(nctab_t),intent(inout) :: nctab
1575  type(pawrad_type),intent(in) :: mesh
1576 !arrays
1577  real(dp),intent(in) :: valr(mesh%mesh_size),qgrid_vl(mqgrid_vl)
1578 
1579 !Local variables-------------------------------
1580 !scalars
1581  real(dp) :: fact,yp1,ypn,d2nvdq0
1582 
1583 ! *************************************************************************
1584 
1585  nctab%has_tvale = .True.
1586  if (.not. allocated(nctab%tvalespl)) then
1587    ABI_MALLOC(nctab%tvalespl, (mqgrid_vl, 2))
1588  else
1589    ABI_CHECK(size(nctab%tvalespl, dim=1) == mqgrid_vl, "wrong mqgrid_vl")
1590  end if
1591 
1592  call pawpsp_cg(nctab%dnvdq0,d2nvdq0,mqgrid_vl,qgrid_vl,nctab%tvalespl(:,1),mesh,valr,yp1,ypn)
1593  call simp_gen(yp1,mesh%rad**2 * valr, mesh); write(std_out,*)"valence charge integrates to: ",four_pi*yp1
1594 
1595  ! Rescale the integral to have the correct number of valence electrons.
1596  ! In some cases, indeed, the radial mesh is not large enough and some valence charge is missing
1597  ! pawpsp_cg extrapolates the integrand beyond rmax but this is not enough.
1598  ! Remember that tvalespl is used to build an initial guess for rhor hence it's very important
1599  ! to have the correct electrostatic.
1600  fact = zion / nctab%tvalespl(1,1)
1601  nctab%tvalespl(:,1) = nctab%tvalespl(:,1) * fact
1602 
1603  ! Compute second derivative of tvalespl(q)
1604  call paw_spline(qgrid_vl,nctab%tvalespl(:,1),mqgrid_vl,yp1,ypn,nctab%tvalespl(:,2))
1605 
1606 end subroutine nctab_eval_tvalespl

m_psps/nctab_free [ Functions ]

[ Top ] [ m_psps ] [ Functions ]

NAME

  nctab_free

FUNCTION

 Free memory allocated in nctab_t

PARENTS

      m_psps,pspini

CHILDREN

      nctab_free,nctab_init

SOURCE

1452 subroutine nctab_free(nctab)
1453 
1454 
1455 !This section has been created automatically by the script Abilint (TD).
1456 !Do not modify the following lines by hand.
1457 #undef ABI_FUNC
1458 #define ABI_FUNC 'nctab_free'
1459 !End of the abilint section
1460 
1461  implicit none
1462 
1463 !Arguments ------------------------------------
1464 !scalars
1465  type(nctab_t),intent(inout) :: nctab
1466 
1467 ! *************************************************************************
1468 
1469  nctab%mqgrid_vl = 0
1470 
1471  nctab%has_tvale = .False.
1472  if (allocated(nctab%tvalespl)) then
1473    ABI_FREE(nctab%tvalespl)
1474  end if
1475 
1476  nctab%has_tcore = .False.
1477  if (allocated(nctab%tcorespl)) then
1478    ABI_FREE(nctab%tcorespl)
1479  end if
1480 
1481 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.

PARENTS

      m_psps,pspini

CHILDREN

      nctab_free,nctab_init

SOURCE

1399 subroutine nctab_init(nctab, mqgrid_vl, has_tcore, has_tvale)
1400 
1401 
1402 !This section has been created automatically by the script Abilint (TD).
1403 !Do not modify the following lines by hand.
1404 #undef ABI_FUNC
1405 #define ABI_FUNC 'nctab_init'
1406 !End of the abilint section
1407 
1408  implicit none
1409 
1410 !Arguments ------------------------------------
1411 !scalars
1412  integer,intent(in) :: mqgrid_vl
1413  logical,intent(in) :: has_tcore,has_tvale
1414  type(nctab_t),intent(inout) :: nctab
1415 
1416 ! *************************************************************************
1417 
1418  nctab%mqgrid_vl = mqgrid_vl
1419 
1420  ! The array for the model core charge is always allocated and initialized with zeros.
1421  ! This approach is similar to the one used in the PAW code.
1422  ! has_tcore tells us whether the model core charge is present or not.
1423  nctab%has_tcore = has_tcore
1424  nctab%dncdq0 = zero; nctab%d2ncdq0 = zero
1425  ABI_CALLOC(nctab%tcorespl, (mqgrid_vl, 2))
1426 
1427  ! tvalespl is allocated only if available.
1428  nctab%has_tvale = has_tvale
1429  nctab%dnvdq0 = zero
1430  if (has_tvale) then
1431    ABI_MALLOC(nctab%tvalespl, (mqgrid_vl, 2))
1432  end if
1433 
1434 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

PARENTS

      pspini

CHILDREN

      nctab_free,nctab_init

SOURCE

1721 subroutine nctab_mixalch(nctabs, npspalch, ntypalch, algalch, mixalch, mixtabs)
1722 
1723 
1724 !This section has been created automatically by the script Abilint (TD).
1725 !Do not modify the following lines by hand.
1726 #undef ABI_FUNC
1727 #define ABI_FUNC 'nctab_mixalch'
1728 !End of the abilint section
1729 
1730  implicit none
1731 
1732 !Arguments ------------------------------------
1733 !scalars
1734  integer,intent(in) :: npspalch,ntypalch
1735 !arrays
1736  integer,intent(in) :: algalch(ntypalch)
1737  real(dp),intent(in) :: mixalch(npspalch, ntypalch)
1738  type(nctab_t),intent(in) :: nctabs(npspalch)
1739  type(nctab_t),target,intent(inout) :: mixtabs(ntypalch)
1740 
1741 !Local variables-------------------------------
1742 !scalars
1743  integer :: ipspalch,itypalch
1744  logical :: has_tcore, has_tvale
1745  real(dp) :: mc
1746  type(nctab_t),pointer :: mix
1747 
1748 ! *************************************************************************
1749 
1750  ABI_CHECK(all(nctabs(:)%mqgrid_vl == nctabs(1)%mqgrid_vl), "Wrong mqgrid_vl")
1751  ABI_CHECK(all(algalch == 1), "algalch /= 1 not implemented")
1752 
1753  do itypalch=1,ntypalch
1754 
1755    ! has_tcore is true is at least one pseudo has nlcc.
1756    ! has_tvale is true if *all* mixed pseudos have the PS valence charge.
1757    has_tcore = .False.; has_tvale = .True.
1758    do ipspalch=1,npspalch
1759      if (abs(mixalch(ipspalch,itypalch)) < tol6) cycle
1760      if (nctabs(ipspalch)%has_tcore) has_tcore = .True.
1761      if (.not. nctabs(ipspalch)%has_tvale) has_tvale = .False.
1762    end do
1763    !write(std_out,*)has_tvale, has_tcore
1764 
1765    call nctab_free(mixtabs(itypalch))
1766    call nctab_init(mixtabs(itypalch), nctabs(1)%mqgrid_vl, has_tcore, has_tvale)
1767    mix => mixtabs(itypalch)
1768 
1769    do ipspalch=1,npspalch
1770      mc = mixalch(ipspalch,itypalch)
1771      if (abs(mc) < tol6) cycle
1772      ! Linear combination of the quantities
1773      ! Mix core for NLCC
1774      if (has_tcore) then
1775        mix%tcorespl = mix%tcorespl + mc * nctabs(ipspalch)%tcorespl
1776        mix%dncdq0 = mix%dncdq0 + mc * nctabs(ipspalch)%dncdq0
1777        mix%d2ncdq0 = mix%d2ncdq0 + mc * nctabs(ipspalch)%d2ncdq0
1778      end if
1779      ! Mix pseudo valence charge.
1780      if (has_tvale) then
1781        mix%tvalespl = mix%tvalespl + mc * nctabs(ipspalch)%tvalespl
1782        mix%dnvdq0 = mix%dnvdq0 + mc * nctabs(ipspalch)%dnvdq0
1783      end if
1784    end do
1785 
1786  end do
1787 
1788 end subroutine nctab_mixalch

m_psps/psp2params_copy [ Functions ]

[ Top ] [ m_psps ] [ Functions ]

NAME

 psp2params_copy

FUNCTION

INPUTS

OUTPUT

PARENTS

      m_psps

CHILDREN

      nctab_free,nctab_init

SOURCE

1277 subroutine psp2params_copy(gth_paramsin, gth_paramsout)
1278 
1279 
1280 !This section has been created automatically by the script Abilint (TD).
1281 !Do not modify the following lines by hand.
1282 #undef ABI_FUNC
1283 #define ABI_FUNC 'psp2params_copy'
1284 !End of the abilint section
1285 
1286  implicit none
1287 
1288 !Arguments ------------------------------------
1289 !scalars
1290  type(pseudopotential_gth_type),intent(in) :: gth_paramsin
1291  type(pseudopotential_gth_type),intent(out) :: gth_paramsout
1292 
1293 ! *********************************************************************
1294 
1295  if (allocated(gth_paramsin%psppar)) then
1296    call alloc_copy( gth_paramsin%psppar, gth_paramsout%psppar)
1297  end if
1298  if (allocated(gth_paramsin%radii_cf)) then
1299    call alloc_copy( gth_paramsin%radii_cf, gth_paramsout%radii_cf)
1300  end if
1301  if (allocated(gth_paramsin%psp_k_par)) then
1302    call alloc_copy( gth_paramsin%psp_k_par, gth_paramsout%psp_k_par)
1303  end if
1304  if (allocated(gth_paramsin%hasGeometry)) then
1305    call alloc_copy( gth_paramsin%hasGeometry, gth_paramsout%hasGeometry)
1306  end if
1307  if (allocated(gth_paramsin%set)) then
1308    call alloc_copy( gth_paramsin%set, gth_paramsout%set)
1309  end if
1310 
1311 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.

PARENTS

      m_psps

CHILDREN

      nctab_free,nctab_init

SOURCE

1336 subroutine psp2params_free(gth_params)
1337 
1338 
1339 !This section has been created automatically by the script Abilint (TD).
1340 !Do not modify the following lines by hand.
1341 #undef ABI_FUNC
1342 #define ABI_FUNC 'psp2params_free'
1343 !End of the abilint section
1344 
1345  implicit none
1346 
1347 !Arguments ------------------------------------
1348 !scalars
1349  type(pseudopotential_gth_type),intent(inout) :: gth_params
1350 
1351 ! *********************************************************************
1352 
1353 !Check arrays.
1354  if (allocated(gth_params%set))  then
1355    ABI_DEALLOCATE(gth_params%set)
1356  end if
1357  if (allocated(gth_params%hasGeometry))  then
1358    ABI_DEALLOCATE(gth_params%hasGeometry)
1359  end if
1360 
1361 !Coefficients for local part and projectors
1362  if (allocated(gth_params%psppar))  then
1363    ABI_DEALLOCATE(gth_params%psppar)
1364  end if
1365 
1366 !Coefficients for spin orbit part
1367  if (allocated(gth_params%psp_k_par))  then
1368    ABI_DEALLOCATE(gth_params%psp_k_par)
1369  end if
1370 
1371 !Different radii
1372  if (allocated(gth_params%radii_cf))  then
1373    ABI_DEALLOCATE(gth_params%radii_cf)
1374  end if
1375 
1376 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.

PARENTS

      m_psps

CHILDREN

      nctab_free,nctab_init

SOURCE

1217 subroutine psp2params_init(gth_params, npsp)
1218 
1219 
1220 !This section has been created automatically by the script Abilint (TD).
1221 !Do not modify the following lines by hand.
1222 #undef ABI_FUNC
1223 #define ABI_FUNC 'psp2params_init'
1224 !End of the abilint section
1225 
1226  implicit none
1227 
1228 !Arguments ------------------------------------
1229 !scalars
1230  integer,intent(in) :: npsp
1231  type(pseudopotential_gth_type),intent(out) :: gth_params
1232 
1233 ! *********************************************************************
1234 
1235 !Check array, no params are currently set.
1236  ABI_ALLOCATE(gth_params%set,(npsp))
1237  gth_params%set(:) = .false.
1238 
1239 !Check array, have geometric information been filled?
1240  ABI_ALLOCATE(gth_params%hasGeometry,(npsp))
1241  gth_params%hasGeometry(:) = .false.
1242 
1243 !Coefficients for local part and projectors
1244  ABI_ALLOCATE(gth_params%psppar,(0:4, 0:6, npsp))
1245  gth_params%psppar = zero
1246 
1247 !Coefficients for spin orbit part
1248  ABI_ALLOCATE(gth_params%psp_k_par,(1:4, 1:3, npsp))
1249  gth_params%psp_k_par = zero
1250 
1251 !Different radii
1252  ABI_ALLOCATE(gth_params%radii_cf,(npsp, 3))
1253 
1254 end subroutine psp2params_init

m_psps/psps_copy [ Functions ]

[ Top ] [ m_psps ] [ Functions ]

NAME

 psps_copy

FUNCTION

 Copy the psps structure.

OUTPUT

SIDE EFFECTS

PARENTS

      m_ddb_hdr

CHILDREN

      nctab_free,nctab_init

SOURCE

677 subroutine psps_copy(pspsin, pspsout)
678 
679 
680 !This section has been created automatically by the script Abilint (TD).
681 !Do not modify the following lines by hand.
682 #undef ABI_FUNC
683 #define ABI_FUNC 'psps_copy'
684 !End of the abilint section
685 
686  implicit none
687 
688 !Arguments ------------------------------------
689 !scalars
690  type(pseudopotential_type),intent(in) :: pspsin
691  type(pseudopotential_type),intent(out) :: pspsout
692 
693 !Local variables-------------------------------
694 !scalars
695  integer :: ii
696 
697 ! *************************************************************************
698 
699  ! integer
700  pspsout%dimekb         = pspsin%dimekb
701  pspsout%lmnmax         = pspsin%lmnmax
702  pspsout%lnmax          = pspsin%lnmax
703  pspsout%mproj          = pspsin%mproj
704  pspsout%mpsang         = pspsin%mpsang
705  pspsout%mpspso         = pspsin%mpspso
706  pspsout%mpssoang       = pspsin%mpssoang
707  pspsout%mqgrid_ff      = pspsin%mqgrid_ff
708  pspsout%mqgrid_vl      = pspsin%mqgrid_vl
709  pspsout%mtypalch       = pspsin%mtypalch
710  pspsout%npsp           = pspsin%npsp
711  pspsout%npspalch       = pspsin%npspalch
712  pspsout%ntypat         = pspsin%ntypat
713  pspsout%ntypalch       = pspsin%ntypalch
714  pspsout%ntyppure       = pspsin%ntyppure
715  pspsout%n1xccc         = pspsin%n1xccc
716  pspsout%optnlxccc      = pspsin%optnlxccc
717  pspsout%positron       = pspsin%positron
718  pspsout%usepaw         = pspsin%usepaw
719  pspsout%usewvl         = pspsin%usewvl
720  pspsout%useylm         = pspsin%useylm
721  pspsout%nc_xccc_gspace = pspsin%nc_xccc_gspace
722 
723  ! logical
724  pspsout%vlspl_recipSpace = pspsin%vlspl_recipSpace
725 
726  ! integer allocatable
727  if (allocated(pspsin%algalch)) then
728    call alloc_copy( pspsin%algalch, pspsout%algalch)
729  end if
730  if (allocated(pspsin%indlmn)) then
731    call alloc_copy( pspsin%indlmn, pspsout%indlmn)
732  end if
733  if (allocated(pspsin%pspdat)) then
734    call alloc_copy( pspsin%pspdat, pspsout%pspdat)
735  end if
736  if (allocated(pspsin%pspcod)) then
737    call alloc_copy( pspsin%pspcod, pspsout%pspcod)
738  end if
739  if (allocated(pspsin%pspso)) then
740    call alloc_copy( pspsin%pspso, pspsout%pspso)
741  end if
742  if (allocated(pspsin%pspxc)) then
743    call alloc_copy( pspsin%pspxc, pspsout%pspxc)
744  end if
745 
746  ! real allocatable
747  if (allocated(pspsin%ekb)) then
748    call alloc_copy( pspsin%ekb, pspsout%ekb)
749  end if
750  if (allocated(pspsin%ffspl)) then
751    call alloc_copy( pspsin%ffspl, pspsout%ffspl)
752  end if
753  if (allocated(pspsin%mixalch)) then
754    call alloc_copy( pspsin%mixalch, pspsout%mixalch)
755  end if
756  if (allocated(pspsin%qgrid_ff)) then
757    call alloc_copy( pspsin%qgrid_ff, pspsout%qgrid_ff)
758  end if
759  if (allocated(pspsin%qgrid_vl)) then
760    call alloc_copy( pspsin%qgrid_vl, pspsout%qgrid_vl)
761  end if
762  if (allocated(pspsin%vlspl)) then
763    call alloc_copy( pspsin%vlspl, pspsout%vlspl)
764  end if
765  if (allocated(pspsin%dvlspl)) then
766    call alloc_copy( pspsin%dvlspl, pspsout%dvlspl)
767  end if
768  if (allocated(pspsin%xcccrc)) then
769    call alloc_copy( pspsin%xcccrc, pspsout%xcccrc)
770  end if
771  if (allocated(pspsin%xccc1d)) then
772    call alloc_copy( pspsin%xccc1d, pspsout%xccc1d)
773  end if
774  if (allocated(pspsin%zionpsp)) then
775    call alloc_copy( pspsin%zionpsp, pspsout%zionpsp)
776  end if
777  if (allocated(pspsin%ziontypat)) then
778    call alloc_copy( pspsin%ziontypat, pspsout%ziontypat)
779  end if
780  if (allocated(pspsin%znuclpsp)) then
781    call alloc_copy( pspsin%znuclpsp, pspsout%znuclpsp)
782  end if
783  if (allocated(pspsin%znucltypat)) then
784    call alloc_copy( pspsin%znucltypat, pspsout%znucltypat)
785  end if
786 
787  ! allocate and copy character strings
788  ABI_ALLOCATE(pspsout%filpsp,(pspsout%npsp))
789  ABI_ALLOCATE(pspsout%title,(pspsout%npsp))
790  ABI_ALLOCATE(pspsout%md5_pseudos,(pspsout%npsp))
791  do ii=1,pspsout%npsp
792    pspsout%filpsp(ii) = pspsin%filpsp(ii)
793    pspsout%title(ii) = pspsin%title(ii)
794    pspsout%md5_pseudos(ii) = pspsin%md5_pseudos(ii)
795  end do
796 
797  ! allocate and copy objects
798  if (allocated(pspsin%nctab)) then
799    ABI_DATATYPE_ALLOCATE(pspsout%nctab,(pspsout%ntypat))
800    do ii=1,pspsout%ntypat
801      call nctab_copy(pspsin%nctab(ii), pspsout%nctab(ii))
802    end do
803  end if
804 
805  call psp2params_copy(pspsin%gth_params, pspsout%gth_params)
806 
807 end subroutine psps_copy

m_psps/psps_free [ Functions ]

[ Top ] [ m_psps ] [ Functions ]

NAME

 psps_free

FUNCTION

 Deallocate all memory of psps structure.

OUTPUT

SIDE EFFECTS

 psps=<type pseudopotential_type>the pseudopotentials description

PARENTS

      driver,m_ddb_hdr

CHILDREN

      nctab_free,nctab_init

SOURCE

554 subroutine psps_free(psps)
555 
556 
557 !This section has been created automatically by the script Abilint (TD).
558 !Do not modify the following lines by hand.
559 #undef ABI_FUNC
560 #define ABI_FUNC 'psps_free'
561 !End of the abilint section
562 
563  implicit none
564 
565 !Arguments ------------------------------------
566 !scalars
567  type(pseudopotential_type),intent(inout) :: psps
568 
569 !Local variables-------------------------------
570 !scalars
571  integer :: ii
572 
573 ! *************************************************************************
574 
575 !Allocation of some arrays independent of the dataset
576  if (allocated(psps%filpsp))  then
577    ABI_DEALLOCATE(psps%filpsp)
578  end if
579  if (allocated(psps%pspcod))  then
580    ABI_DEALLOCATE(psps%pspcod)
581  end if
582  if (allocated(psps%pspdat))  then
583    ABI_DEALLOCATE(psps%pspdat)
584  end if
585  if (allocated(psps%pspso))  then
586    ABI_DEALLOCATE(psps%pspso)
587  end if
588  if (allocated(psps%pspxc))  then
589    ABI_DEALLOCATE(psps%pspxc)
590  end if
591  if (allocated(psps%title))  then
592    ABI_DEALLOCATE(psps%title)
593  end if
594  if (allocated(psps%zionpsp))  then
595    ABI_DEALLOCATE(psps%zionpsp)
596  end if
597  if (allocated(psps%znuclpsp))  then
598    ABI_DEALLOCATE(psps%znuclpsp)
599  end if
600  if (allocated(psps%algalch))  then
601    ABI_DEALLOCATE(psps%algalch)
602  end if
603  if (allocated(psps%mixalch))  then
604    ABI_DEALLOCATE(psps%mixalch)
605  end if
606  if (allocated(psps%ekb))  then
607    ABI_DEALLOCATE(psps%ekb)
608  end if
609  if (allocated(psps%indlmn))  then
610    ABI_DEALLOCATE(psps%indlmn)
611  end if
612  if (allocated(psps%ffspl))  then
613    ABI_DEALLOCATE(psps%ffspl)
614  end if
615  if (allocated(psps%qgrid_ff))  then
616    ABI_DEALLOCATE(psps%qgrid_ff)
617  end if
618  if (allocated(psps%qgrid_vl))  then
619    ABI_DEALLOCATE(psps%qgrid_vl)
620  end if
621  if (allocated(psps%vlspl))  then
622    ABI_DEALLOCATE(psps%vlspl)
623  end if
624  if (allocated(psps%dvlspl))  then
625    ABI_DEALLOCATE(psps%dvlspl)
626  end if
627  if (allocated(psps%xccc1d))  then
628    ABI_DEALLOCATE(psps%xccc1d)
629  end if
630  if (allocated(psps%xcccrc))  then
631    ABI_DEALLOCATE(psps%xcccrc)
632  end if
633  if (allocated(psps%ziontypat))  then
634    ABI_DEALLOCATE(psps%ziontypat)
635  end if
636  if (allocated(psps%znucltypat))  then
637    ABI_DEALLOCATE(psps%znucltypat)
638  end if
639  if (allocated(psps%md5_pseudos))  then
640    ABI_DEALLOCATE(psps%md5_pseudos)
641  end if
642 
643  ! Free types.
644  call psp2params_free(psps%gth_params)
645 
646  if (allocated(psps%nctab)) then
647    do ii=1,size(psps%nctab)
648      call nctab_free(psps%nctab(ii))
649    end do
650    ABI_DT_FREE(psps%nctab)
651  end if
652 
653 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

OUTPUT

SIDE EFFECTS

 psps=<type pseudopotential_type>the pseudopotentials description

PARENTS

      driver

CHILDREN

      nctab_free,nctab_init

SOURCE

280 subroutine psps_init_from_dtset(dtset, idtset, psps, pspheads)
281 
282 
283 !This section has been created automatically by the script Abilint (TD).
284 !Do not modify the following lines by hand.
285 #undef ABI_FUNC
286 #define ABI_FUNC 'psps_init_from_dtset'
287 !End of the abilint section
288 
289  implicit none
290 
291 !Arguments ------------------------------------
292 !scalars
293  integer,intent(in) :: idtset
294  type(dataset_type),intent(in) :: dtset
295  type(pseudopotential_type),intent(inout) :: psps
296 !arrays
297  type(pspheader_type),intent(in) :: pspheads(psps%npsp)
298 
299 !Local variables-------------------------------
300 !scalars
301  integer,save :: dimekb_old=-1,lmnmax_old=-1,lnmax_old=-1,mqgridff_old=0
302  integer,save :: mqgridvl_old=0,ntypat_old=-1,usepaw_old=-1
303  integer :: ipsp,lmnmax,lmnmaxso,lnmax,lnmaxso,newmqgrid,newmqgriddg,nptsgvec
304  integer :: changed,ii,itypat
305  real(dp) :: gprimd_orig(3,3)
306 
307 ! *************************************************************************
308 
309  psps%optnlxccc   = dtset%optnlxccc
310 !Determine the number of points needed in reciprocal space to represent the
311 !pseudopotentials (either set by hand from input variable or set automatically by abinit)
312  nptsgvec         = 200 !This has to be chosen one and for all or else ??
313  newmqgrid        = dtset%mqgrid
314  newmqgriddg      = dtset%mqgriddg
315  !JB:Which image to use ? I guess 1 always works
316  call matr3inv(dtset%rprimd_orig(:,:,1),gprimd_orig)
317  if ( dtset%usewvl == 0) then
318    call setmqgrid(newmqgrid,newmqgriddg,dtset%ecut*dtset%dilatmx**2,&
319 &       dtset%pawecutdg*dtset%dilatmx**2,gprimd_orig,nptsgvec,psps%usepaw)
320  else
321    call setmqgrid(newmqgrid,newmqgriddg,one,one,gprimd_orig,nptsgvec,psps%usepaw)
322  end if
323  psps%mqgrid_ff   = newmqgrid
324  if (psps%usepaw == 1) then
325    psps%mqgrid_vl = newmqgriddg
326  else
327    psps%mqgrid_vl = newmqgrid
328  end if
329 
330 !Determine the maximum number of projectors, for the set of pseudo atom
331  call getdim_nloc(lmnmax,lmnmaxso,lnmax,lnmaxso,dtset%mixalch_orig,dtset%nimage,psps%npsp,dtset%npspalch,&
332 & dtset%ntypat,dtset%ntypalch,pspheads)
333 
334  psps%npspalch = dtset%npspalch
335  psps%ntypat   = dtset%ntypat
336  psps%ntypalch = dtset%ntypalch
337  psps%ntyppure = dtset%ntyppure
338 
339 !Set the flag for reciprocal space or real space calculations
340  psps%vlspl_recipSpace = (dtset%icoulomb /= 1)
341 !changed by RShaltaf
342  psps%positron = dtset%positron
343  psps%useylm   = dtset%useylm
344 
345 !Added by T. Rangel for WVL+PAW
346  psps%usewvl   = dtset%usewvl
347 
348 ! Define treatment of the model core density for NC pseudos.
349  psps%nc_xccc_gspace = dtset%nc_xccc_gspace
350 
351  if (idtset > 1) then
352    if (allocated(psps%algalch))  then
353      ABI_DEALLOCATE(psps%algalch)
354    end if
355    if (allocated(psps%mixalch))  then
356      ABI_DEALLOCATE(psps%mixalch)
357    end if
358  end if
359  ABI_ALLOCATE(psps%algalch,(psps%ntypalch))
360  ABI_ALLOCATE(psps%mixalch,(psps%npspalch,psps%ntypalch))
361  psps%algalch(1:psps%ntypalch)=dtset%algalch(1:psps%ntypalch)
362 !This value will be overwritten elsewhere in case there are different images ...
363  psps%mixalch(1:psps%npspalch,1:psps%ntypalch)=dtset%mixalch_orig(1:psps%npspalch,1:psps%ntypalch,1)
364 
365 !Set mpspso and psps%pspso
366 !Warning: mpspso might be different for each dataset.
367 !         mpspso not relevant in case of PAW.
368  psps%mpspso=1
369  do ipsp=1,dtset%npsp
370    if(dtset%nspinor==1)then
371      psps%pspso(ipsp)=0
372      ! This is needed to treate SOC perturbatively in sigma.
373      !if (dtset%optdriver == RUNL_SIGMA .and. dtset%so_psp(ipsp) /= 0) then
374      !  MSG_WARNING("Setting pspso to 2 although nspinor == 1")
375      !  psps%pspso(ipsp) = 2
376      !end if
377 !    Ideally the following line should not exist,
378 !      but at present, the space has to be booked
379      if(pspheads(ipsp)%pspso/=0)psps%mpspso=2
380    else if (psps%usepaw==0) then
381      if(dtset%so_psp(ipsp)/=1)then
382        psps%pspso(ipsp)=dtset%so_psp(ipsp)
383      else
384        psps%pspso(ipsp)=pspheads(ipsp)%pspso
385      end if
386      if(psps%pspso(ipsp)/=0)psps%mpspso=2
387      if(pspheads(ipsp)%pspso/=0)psps%mpspso=2
388    else
389      psps%pspso(ipsp)=1+dtset%pawspnorb
390    end if
391  end do
392 
393 !Set mpssoang, lmnmax, lnmax
394  if(psps%mpspso==1)then
395    psps%mpssoang=psps%mpsang
396    psps%lmnmax  =lmnmax
397    psps%lnmax   =lnmax
398  else
399    psps%mpssoang=2*psps%mpsang-1
400    psps%lmnmax=lmnmaxso
401    psps%lnmax=lnmaxso
402  end if
403 
404 !T. Rangel: for wvl + paw do not change psps%lmnmax
405  if (psps%useylm==0 .and. psps%usepaw/=1 ) then
406    psps%lmnmax=psps%lnmax
407  end if
408 
409 !Set dimekb
410  if (psps%usepaw==0) then
411    psps%dimekb=psps%lnmax
412  else
413    psps%dimekb=psps%lmnmax*(psps%lmnmax+1)/2
414  end if
415 
416 !The following arrays are often not deallocated before the end of the dtset loop
417 !and might keep their content from one dataset to the other, if the conditions are fulfilled
418  changed = 0
419 
420  if(dimekb_old/=psps%dimekb .or. ntypat_old/=dtset%ntypat .or. usepaw_old/=psps%usepaw) then
421    changed = changed + 1
422    if(idtset/=1) then
423      if (allocated(psps%ekb))  then
424        ABI_DEALLOCATE(psps%ekb)
425      end if
426    end if
427    ABI_ALLOCATE(psps%ekb,(psps%dimekb,dtset%ntypat*(1-psps%usepaw)))
428    dimekb_old=psps%dimekb
429  end if
430 
431  if(lmnmax_old/=psps%lmnmax .or. ntypat_old/=dtset%ntypat)then
432    changed = changed + 1
433    if(idtset/=1) then
434      if (allocated(psps%indlmn))  then
435        ABI_DEALLOCATE(psps%indlmn)
436      end if
437    end if
438    ABI_ALLOCATE(psps%indlmn,(6,psps%lmnmax,dtset%ntypat))
439    lmnmax_old=psps%lmnmax
440  end if
441 
442  if(mqgridff_old/=psps%mqgrid_ff .or. lnmax_old/=psps%lnmax .or. ntypat_old/=dtset%ntypat)then
443    changed = changed + 1
444    if(idtset/=1) then
445      if (allocated(psps%ffspl))  then
446        ABI_DEALLOCATE(psps%ffspl)
447      end if
448      if (allocated(psps%qgrid_ff))  then
449        ABI_DEALLOCATE(psps%qgrid_ff)
450      end if
451    end if
452    ABI_ALLOCATE(psps%ffspl,(psps%mqgrid_ff,2,psps%lnmax,dtset%ntypat))
453    ABI_ALLOCATE(psps%qgrid_ff,(psps%mqgrid_ff))
454    mqgridff_old=psps%mqgrid_ff
455    lnmax_old=psps%lnmax
456  end if
457 
458  if(mqgridvl_old/=psps%mqgrid_vl .or. ntypat_old/=dtset%ntypat)then
459    changed = changed + 1
460    if(idtset/=1) then
461      if (allocated(psps%qgrid_vl))  then
462        ABI_DEALLOCATE(psps%qgrid_vl)
463      end if
464      if (allocated(psps%vlspl))  then
465        ABI_DEALLOCATE(psps%vlspl)
466      end if
467      if (allocated(psps%nctab)) then
468        do ii=1,size(psps%nctab)
469          call nctab_free(psps%nctab(ii))
470        end do
471        ABI_DT_FREE(psps%nctab)
472      end if
473    end if
474    if (idtset/=1 .and. .not.psps%vlspl_recipSpace) then
475      if (allocated(psps%dvlspl))  then
476        ABI_DEALLOCATE(psps%dvlspl)
477      end if
478    end if
479 
480    ABI_ALLOCATE(psps%qgrid_vl,(psps%mqgrid_vl))
481    ABI_ALLOCATE(psps%vlspl,(psps%mqgrid_vl,2,dtset%ntypat))
482 
483    if (psps%usepaw == 0) then
484      ! If you change usepaw in the input, you will get what you deserve!
485      ABI_DT_MALLOC(psps%nctab, (dtset%ntypat))
486      do itypat=1,dtset%ntypat
487        call nctab_init(psps%nctab(itypat), psps%mqgrid_vl, .False., .False.)
488      end do
489    end if
490 
491    if (.not.psps%vlspl_recipSpace) then
492      ABI_ALLOCATE(psps%dvlspl,(psps%mqgrid_vl,2,dtset%ntypat))
493    end if
494    mqgridvl_old=psps%mqgrid_vl
495  end if
496 
497  if(ntypat_old/=dtset%ntypat.or. usepaw_old/=psps%usepaw)then
498    changed = changed + 1
499    if(idtset/=1) then
500      if (allocated(psps%xccc1d))  then
501        ABI_DEALLOCATE(psps%xccc1d)
502      end if
503    end if
504    ABI_ALLOCATE(psps%xccc1d,(psps%n1xccc*(1-psps%usepaw),6,dtset%ntypat))
505    usepaw_old=psps%usepaw
506  end if
507 
508  if(ntypat_old/=dtset%ntypat)then
509    changed = changed + 1
510    if(idtset/=1) then
511      if (allocated(psps%xcccrc))  then
512        ABI_DEALLOCATE(psps%xcccrc)
513      end if
514      if (allocated(psps%ziontypat))  then
515        ABI_DEALLOCATE(psps%ziontypat)
516      end if
517      if (allocated(psps%znucltypat))  then
518        ABI_DEALLOCATE(psps%znucltypat)
519      end if
520    end if
521    ABI_ALLOCATE(psps%xcccrc,(dtset%ntypat))
522    ABI_ALLOCATE(psps%znucltypat,(dtset%ntypat))
523    ABI_ALLOCATE(psps%ziontypat,(dtset%ntypat))
524    ntypat_old=dtset%ntypat
525  end if
526 
527  psps%ziontypat(:)=dtset%ziontypat(:)
528 
529 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

PARENTS

      driver

CHILDREN

      nctab_free,nctab_init

SOURCE

184 subroutine psps_init_global(mtypalch, npsp, psps, pspheads)
185 
186 
187 !This section has been created automatically by the script Abilint (TD).
188 !Do not modify the following lines by hand.
189 #undef ABI_FUNC
190 #define ABI_FUNC 'psps_init_global'
191 !End of the abilint section
192 
193  implicit none
194 
195 !Arguments ------------------------------------
196 !scalars
197  integer,intent(in) :: mtypalch,npsp
198  type(pseudopotential_type),intent(inout) :: psps
199 !arrays
200  type(pspheader_type),intent(in) :: pspheads(npsp)
201 
202 !Local variables-------------------------------
203 !scalars
204  integer :: ii, mpsang, n1xccc
205 
206 ! *************************************************************************
207 
208 !Allocation of some arrays independent of the dataset
209  ABI_ALLOCATE(psps%filpsp,(npsp))
210  ABI_ALLOCATE(psps%pspcod,(npsp))
211  ABI_ALLOCATE(psps%pspdat,(npsp))
212  ABI_ALLOCATE(psps%pspso,(npsp))
213  ABI_ALLOCATE(psps%pspxc,(npsp))
214  ABI_ALLOCATE(psps%title,(npsp))
215  ABI_ALLOCATE(psps%zionpsp,(npsp))
216  ABI_ALLOCATE(psps%znuclpsp,(npsp))
217  call psp2params_init(psps%gth_params, npsp)
218 
219  psps%filpsp(1:npsp)=pspheads(1:npsp)%filpsp
220  psps%pspcod(1:npsp)=pspheads(1:npsp)%pspcod
221  psps%pspdat(1:npsp)=pspheads(1:npsp)%pspdat
222  psps%pspso(1:npsp)=pspheads(1:npsp)%pspso
223  psps%pspxc(1:npsp)=pspheads(1:npsp)%pspxc
224  psps%title(1:npsp)=pspheads(1:npsp)%title
225  psps%zionpsp(1:npsp)=pspheads(1:npsp)%zionpsp
226  psps%znuclpsp(1:npsp)=pspheads(1:npsp)%znuclpsp
227 
228  ! Transfer md5 checksum
229  ABI_ALLOCATE(psps%md5_pseudos, (npsp))
230  psps%md5_pseudos = pspheads(1:npsp)%md5_checksum
231 !Set values independant from dtset
232  psps%npsp   = npsp
233 !Note that mpsang is the max of 1+lmax, with minimal value 1 (even for local psps, at present)
234 !mpsang=max(maxval(pspheads(1:npsp)%lmax)+1,1) ! might not work with HP compiler
235 !n1xccc=maxval(pspheads(1:npsp)%xccc)
236  mpsang=1
237  n1xccc=pspheads(1)%xccc
238  do ii=1,psps%npsp
239    mpsang=max(pspheads(ii)%lmax+1,mpsang)
240    n1xccc=max(pspheads(ii)%xccc,n1xccc)
241  end do
242  psps%mpsang = mpsang
243  psps%n1xccc = n1xccc
244 ! Determine here whether the calculation is PAW
245 ! If paw, all pspcod necessarily are 7 or 17 (see iofn2)
246  psps%usepaw  =0
247  if (pspheads(1)%pspcod==7.or.pspheads(1)%pspcod==17) psps%usepaw=1
248  psps%mtypalch = mtypalch
249 
250 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.

PARENTS

      pspini

CHILDREN

      nctab_free,nctab_init

SOURCE

1052 subroutine psps_ncwrite(psps, path)
1053 
1054 
1055 !This section has been created automatically by the script Abilint (TD).
1056 !Do not modify the following lines by hand.
1057 #undef ABI_FUNC
1058 #define ABI_FUNC 'psps_ncwrite'
1059 !End of the abilint section
1060 
1061  implicit none
1062 
1063 !Arguments ------------------------------------
1064 !scalars
1065  character(len=*),intent(in) :: path
1066  type(pseudopotential_type),intent(in) :: psps
1067 
1068 !Local variables-------------------------------
1069 !scalars
1070  integer :: ipsp,itypat,ncid,ncerr
1071 
1072 ! *************************************************************************
1073 
1074 #ifdef HAVE_NETCDF
1075  NCF_CHECK(nctk_open_create(ncid, path, xmpi_comm_self))
1076 
1077  ! Define dimensions
1078  ncerr = nctk_def_dims(ncid, [&
1079      nctkdim_t("fnlen", fnlen+1), nctkdim_t("md5_slen", md5_slen+1), nctkdim_t("ntypat", psps%ntypat), &
1080      nctkdim_t("npsp", psps%npsp), nctkdim_t("lnmax", psps%lnmax), &
1081      nctkdim_t("lmnmax", psps%lnmax), nctkdim_t("dimekb", psps%dimekb), &
1082      nctkdim_t("mqgrid_vl", psps%mqgrid_vl), nctkdim_t("mqgrid_ff", psps%mqgrid_ff) &
1083  ])
1084  NCF_CHECK(ncerr)
1085 
1086  if (psps%n1xccc /= 0) then ! 0 means unlimited!
1087    NCF_CHECK(nctk_def_dims(ncid, nctkdim_t("n1xccc", psps%n1xccc)))
1088  end if
1089 
1090  ! Define variables
1091  ncerr = nctk_def_iscalars(ncid, [character(len=nctk_slen) :: "usepaw", "useylm"])
1092  NCF_CHECK(ncerr)
1093 
1094  ncerr = nctk_def_arrays(ncid, [&
1095    nctkarr_t("ziontypat", "dp", "ntypat"), &
1096    nctkarr_t("znucltypat", "dp", "ntypat"), &
1097    nctkarr_t("qgrid_vl", "dp", "mqgrid_vl"), &
1098    nctkarr_t("qgrid_ff", "dp", "mqgrid_ff"), &
1099    nctkarr_t("vlspl", "dp", "mqgrid_vl, two, ntypat"), &
1100    nctkarr_t("xcccrc", "dp", "ntypat"), &
1101    nctkarr_t("indlmn", "int", "six, lmnmax, ntypat"), &
1102    nctkarr_t("ffspl", "dp", "mqgrid_ff, two, lnmax, ntypat"), &
1103    nctkarr_t("filpsp", "char", "fnlen, npsp"), &
1104    nctkarr_t("md5_pseudos", "char", "md5_slen, npsp") &
1105  ])
1106  NCF_CHECK(ncerr)
1107 
1108  if (psps%usepaw == 0) then
1109    NCF_CHECK(nctk_def_arrays(ncid, nctkarr_t("ekb", "dp", "dimekb, ntypat")))
1110    if (psps%n1xccc /= 0) then
1111      NCF_CHECK(nctk_def_arrays(ncid, nctkarr_t("xccc1d", "dp", "n1xccc, six, ntypat")))
1112    end if
1113    ncerr = nctk_def_arrays(ncid, [&
1114      nctkarr_t("nc_tvalespl", "dp", "mqgrid_vl, two, ntypat"), &
1115      nctkarr_t("nc_tcorespl", "dp", "mqgrid_vl, two, ntypat")  &
1116    ])
1117    NCF_CHECK(ncerr)
1118  end if
1119 
1120  ! Write data
1121  NCF_CHECK(nf90_enddef(ncid))
1122  NCF_CHECK(nf90_put_var(ncid, vid("usepaw"), psps%usepaw))
1123  NCF_CHECK(nf90_put_var(ncid, vid("useylm"), psps%useylm))
1124  NCF_CHECK(nf90_put_var(ncid, vid("ziontypat"), psps%ziontypat))
1125  NCF_CHECK(nf90_put_var(ncid, vid("znucltypat"), psps%znucltypat))
1126  do ipsp=1,psps%npsp
1127    NCF_CHECK(nf90_put_var(ncid, vid("filpsp"), trim(psps%filpsp(ipsp)), start=[1, ipsp]))
1128    NCF_CHECK(nf90_put_var(ncid, vid("md5_pseudos"), trim(psps%md5_pseudos(ipsp)), start=[1, ipsp]))
1129  end do
1130  NCF_CHECK(nf90_put_var(ncid, vid("qgrid_vl"), psps%qgrid_vl))
1131  NCF_CHECK(nf90_put_var(ncid, vid("qgrid_ff"), psps%qgrid_ff))
1132  NCF_CHECK(nf90_put_var(ncid, vid("indlmn"), psps%indlmn))
1133 
1134  ! Local part in q-space and second derivative
1135  if (allocated(psps%vlspl)) then
1136    NCF_CHECK(nf90_put_var(ncid, vid("vlspl"), psps%vlspl))
1137  end if
1138 
1139  ! Form factors for each type of atom
1140  ! for each type and each (l,n) channel, ffnl(q) and second derivative
1141  if (allocated(psps%ffspl)) then
1142    NCF_CHECK(nf90_put_var(ncid, vid("ffspl"), psps%ffspl))
1143  end if
1144 
1145  ! Pseudo-core charge for each type of atom, on the real-space radial
1146  NCF_CHECK(nf90_put_var(ncid, vid("xcccrc"), psps%xcccrc))
1147  if (psps%usepaw == 0 .and. allocated(psps%xccc1d)) then
1148    NCF_CHECK(nf90_put_var(ncid, vid("xccc1d"), psps%xccc1d))
1149  end if
1150 
1151  ! NC-only: add tcore_spl and tvalespl in q-space
1152  if (psps%usepaw == 0) then
1153    NCF_CHECK(nf90_put_var(ncid, vid("ekb"), psps%ekb))
1154    do itypat=1,psps%ntypat
1155      if (psps%nctab(itypat)%has_tvale) then
1156        ncerr = nf90_put_var(ncid, vid("nc_tvalespl"), psps%nctab(itypat)%tvalespl, start=[1,1,itypat])
1157        NCF_CHECK(ncerr)
1158      end if
1159      if (psps%nctab(itypat)%has_tcore) then
1160        ncerr = nf90_put_var(ncid, vid("nc_tcorespl"), psps%nctab(itypat)%tcorespl, start=[1,1,itypat])
1161        NCF_CHECK(ncerr)
1162      end if
1163    end do
1164  end if
1165 
1166  NCF_CHECK(nf90_close(ncid))
1167 
1168 #else
1169  MSG_WARNING("netcdf support not activated. psps file cannot be created!")
1170 #endif
1171 
1172 contains
1173  integer function vid(vname)
1174 
1175 
1176 !This section has been created automatically by the script Abilint (TD).
1177 !Do not modify the following lines by hand.
1178 #undef ABI_FUNC
1179 #define ABI_FUNC 'vid'
1180 !End of the abilint section
1181 
1182    character(len=*),intent(in) :: vname
1183    vid = nctk_idname(ncid, vname)
1184  end function vid
1185 
1186 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

PARENTS

      pspini

CHILDREN

      nctab_free,nctab_init

SOURCE

 836 subroutine psps_print(psps,unit,prtvol,mode_paral)
 837 
 838 
 839 !This section has been created automatically by the script Abilint (TD).
 840 !Do not modify the following lines by hand.
 841 #undef ABI_FUNC
 842 #define ABI_FUNC 'psps_print'
 843 !End of the abilint section
 844 
 845  implicit none
 846 
 847 !Arguments ------------------------------------
 848 !scalars
 849  integer,intent(in),optional :: prtvol,unit
 850  character(len=4),intent(in),optional :: mode_paral
 851  type(pseudopotential_type),intent(in) :: psps
 852 
 853 !Local variables-------------------------------
 854 !scalars
 855  integer :: ierr,ips,ipsp_alch,ityp_alch,itypat,unt,my_prtvol
 856  character(len=4) :: mode
 857  character(len=500) :: msg
 858 !arrays
 859  integer :: cond_values(4)
 860  character(len=9) :: cond_string(4)
 861 
 862 ! *************************************************************************
 863 
 864  ! Provide defaults
 865  my_prtvol=0; if (present(prtvol)) my_prtvol=prtvol
 866  unt=std_out; if (present(unit)) unt=unit
 867  mode='COLL'; if (present(mode_paral)) mode=mode_paral
 868  ierr=0; cond_string(1:4)=' '; cond_values(:)=0
 869 
 870  ! General info including spin-orbit
 871  call wrtout(unt,' ==== Info on pseudopotentials ==== ', mode)
 872 
 873  SELECT CASE (psps%usepaw)
 874  CASE (0)
 875    call wrtout(unt,'  Norm-conserving pseudopotentials ', mode)
 876    !call wrtout(unt, sjoin('  Max number of Kleinman-Bylander energies ', itoa(psps%dimekb)), mode)
 877    !do itypat=1,psps%ntypat
 878    ! write(msg,'(a,i4,a,f9.4)')' Type ',itypat,' K-B energies ',(psps%ekb(ikbe,itypat),ikbe=1,psps%dimekb)
 879    !end do
 880  CASE (1)
 881    write(msg,'(a)')
 882    call wrtout(unt,'  PAW calculation', mode)
 883    !call wrtout(unt,sjoin('  Max number of D_ij coefficients ', itoa(psps%dimekb)), mode)
 884  CASE DEFAULT
 885    call chkint_eq(0,0,cond_string,cond_values,ierr,'usepaw',psps%usepaw,2,[0,1],unt)
 886  END SELECT
 887 
 888  !SELECT CASE (psps%positron)
 889  !CASE (0)
 890  !  call wrtout(unt, '  Standard Electron Calculation ', mode)
 891  !CASE (1,2)
 892  !  write(msg,'(a,i0)')'  Positron Calculation with positron .. ',psps%positron
 893  !  call wrtout(unt,msg,mode)
 894  !CASE DEFAULT
 895  !  call chkint_eq(0,0,cond_string,cond_values,ierr,'positron',psps%positron,3,[0,1,2],unt)
 896  !END SELECT
 897 
 898  write(msg,'(a,i4,2a,i4)')&
 899 &  '  Number of pseudopotentials .. ',psps%npsp,ch10,&
 900 &  '  Number of types of atoms   .. ',psps%ntypat
 901  call wrtout(unt,msg,mode)
 902 
 903  if (psps%usepaw==0) then
 904    SELECT CASE (psps%mpspso)
 905    CASE (1)
 906      call wrtout(unt,'  Scalar calculation (no spin-orbit term) ',mode)
 907    CASE (2)
 908      write(msg,'(3a,i3)')&
 909 &      '  Calculation with spin-orbit coupling ',ch10,&
 910 &      '  Max number of channels (spin-orbit included) ',psps%mpssoang
 911      call wrtout(unt,msg,mode)
 912      do itypat=1,psps%ntypat
 913        if (psps%pspso(itypat) /= 1) then
 914          write(msg,'(a,i4,a,i2,a)')&
 915 &          '  - Atom type ',itypat,' has spin-orbit characteristics (pspso= ',psps%pspso(itypat),")"
 916          call wrtout(unt,msg,mode)
 917        end if
 918      end do
 919    CASE DEFAULT
 920      call chkint_eq(0,0,cond_string,cond_values,ierr,'mpspso',psps%mpspso,2,[1,2],unt)
 921    END SELECT
 922  else
 923    SELECT CASE (maxval(psps%pspso))
 924    CASE (0,1)
 925      msg='  Scalar calculation (no spin-orbit term) '
 926    CASE (2)
 927      msg='  Calculation with spin-orbit coupling '
 928    END SELECT
 929    call wrtout(unt,msg,mode)
 930  end if
 931 
 932  ! Info on nonlocal part
 933  SELECT CASE (psps%useylm)
 934  CASE (0)
 935    msg = '  Nonlocal part applied using Legendre polynomials '
 936  CASE (1)
 937    msg = '  Nonlocal part applied using real spherical harmonics '
 938  CASE DEFAULT
 939    call chkint_eq(0,0,cond_string,cond_values,ierr,'psps%useylm',psps%useylm,2,(/0,1/),unt)
 940  END SELECT
 941  call wrtout(unt,msg,mode)
 942 
 943  write(msg,'(a,i3)')'  Max number of non-local projectors over l and type ',psps%mproj
 944  call wrtout(unt,msg,mode)
 945 
 946  write(msg,'(a,i3,2a,i3,2a,i3)')&
 947 & '  Highest angular momentum +1 ....... ',psps%mpsang,ch10,&
 948 & '  Max number of (l,n)   components .. ',psps%lnmax, ch10,&
 949 & '  Max number of (l,m,n) components .. ',psps%lmnmax
 950  call wrtout(unt,msg,mode)
 951 
 952  !FIXME for paw n1xccc==1
 953  ! Non-linear Core correction
 954  if (psps%n1xccc/=0) then
 955    write(msg,'(3a,2(a,i4,a),2a)')ch10,&
 956 &    ' Pseudo-Core Charge Info: ',ch10,&
 957 &    '   Number of radial points for pseudo-core charge .. ',psps%n1xccc,ch10,&
 958 &    '   XC core-correction treatment (optnlxccc) ........ ',psps%optnlxccc,ch10,&
 959 &    '   Radius for pseudo-core charge for each type ..... ',ch10
 960    call wrtout(unt,msg,mode)
 961    do itypat=1,psps%ntypat
 962      write(msg,'(a,i4,a,f7.4)')'  - Atom type ',itypat,' has pseudo-core radius .. ',psps%xcccrc(itypat)
 963      call wrtout(unt,msg,mode)
 964    end do
 965  end if
 966 
 967  ! Alchemical mixing
 968  if (psps%mtypalch/=0) then
 969    write(msg,'(3a,3(a,i4,a))')ch10,&
 970 &    ' Calculation with alchemical mixing:',ch10,&
 971 &    '   Number of pure pseudoatoms .... ',psps%ntyppure,ch10,&
 972 &    '   Number of pseudos for mixing .. ',psps%npspalch,ch10,&
 973 &    '   Alchemical pseudoatoms ........ ',psps%ntypalch,ch10
 974    call wrtout(unt,msg,mode)
 975    do ipsp_alch=1,psps%npspalch
 976      do ityp_alch=1,psps%ntypalch
 977        write(std_out,*)' mixalch ',psps%mixalch(ipsp_alch,ityp_alch)
 978      end do
 979    end do
 980    do ityp_alch=1,psps%ntypalch
 981      write(msg,'(a,i4,a,i4)')' For alchemical atom no. ',ityp_alch,' algalch is .. ',psps%algalch(ityp_alch)
 982      call wrtout(unt,msg,mode)
 983    end do
 984  end if
 985 
 986  ! Info in Q-grid for spline of form factors
 987  write(msg,'(3a,a,i6,a,a,i6)')ch10,&
 988 &  ' Info on the Q-grid used for form factors in spline form: ',ch10,&
 989 &  '   Number of q-points for radial functions ffspl .. ',psps%mqgrid_ff,ch10,&
 990 &  '   Number of q-points for vlspl ................... ',psps%mqgrid_vl
 991  call wrtout(unt,msg,mode)
 992 
 993  if (psps%vlspl_recipSpace) then
 994    call wrtout(unt,'   vloc is computed in Reciprocal Space ',mode)
 995  else
 996    call wrtout(unt,'   vloc is computed in Real Space ',mode)
 997  end if
 998  if (psps%usepaw == 0) then
 999    if (psps%nc_xccc_gspace == 0) call wrtout(unt,'   model core charge treated in real-space', mode)
1000    if (psps%nc_xccc_gspace == 1) call wrtout(unt,'   model core charge treated in G-space', mode)
1001  end if
1002 
1003  !TODO additional stuff that might be printed
1004  call wrtout(unt, "", mode)
1005  do itypat=1,psps%ntypat
1006    write(msg,'(a,i0,a,i0)')'  XC functional for type ',itypat,' is ',psps%pspxc(itypat)
1007    call wrtout(unt,msg,mode)
1008    !write(std_out,*)psps%ziontypat(itypat),psps%znucltypat(itypat)
1009    if (psps%usepaw == 0) then
1010      call wrtout(unt, sjoin("  Pseudo valence available: ", yesno(psps%nctab(itypat)%has_tvale)), mode)
1011    end if
1012  end do
1013 
1014  !integer, pointer :: pspxc(:)
1015  ! pspxc(ntypat)
1016  ! For each type of psp, the XC functional that was used to generate it, as given by the psp file
1017  if (my_prtvol>=3) then
1018    do ips=1,psps%npsp
1019      write(std_out,*)' Pseudo number   ',ips,' read from ',trim(psps%filpsp(ips))
1020      write(std_out,*)' Format or code  ',psps%pspcod(ips)
1021      write(std_out,*)' Generation date ',psps%pspdat(ips)
1022      write(std_out,*)' Content of first line: ', trim(psps%title(ips))
1023    end do
1024  end if
1025 
1026  call wrtout(unt, "", mode)
1027 
1028 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 if UPF file.

PARENTS

      pspatm

CHILDREN

      nctab_free,nctab_init

SOURCE

 98 subroutine test_xml_xmlpaw_upf(path, usexml, xmlpaw, useupf)
 99 
100 
101 !This section has been created automatically by the script Abilint (TD).
102 !Do not modify the following lines by hand.
103 #undef ABI_FUNC
104 #define ABI_FUNC 'test_xml_xmlpaw_upf'
105 !End of the abilint section
106 
107  implicit none
108 
109 !Arguments ------------------------------------
110 !scalars
111  character(len=*),intent(in) :: path
112  integer,intent(out) :: usexml, xmlpaw, useupf
113 
114 !Local variables-------------------------------
115 !scalars
116  integer :: temp_unit
117  character(len=500) :: msg,errmsg
118  character(len=70) :: testxml
119 
120 ! *************************************************************************
121 
122 !  Check if the file pseudopotential file is written in XML
123  usexml = 0; xmlpaw = 0; useupf = 0
124 
125  if (open_file(path,msg,newunit=temp_unit,form='formatted',status='old') /= 0) then
126    MSG_ERROR(msg)
127  end if
128  rewind (unit=temp_unit,err=10,iomsg=errmsg)
129 
130  read(temp_unit,*,err=10,iomsg=errmsg) testxml
131  if(testxml(1:5)=='<?xml')then
132    usexml = 1
133    read(temp_unit,*,err=10,iomsg=errmsg) testxml
134    if(testxml(1:4)=='<paw') xmlpaw = 1
135  else
136    usexml = 0
137  end if
138 
139  !Check if pseudopotential file is a Q-espresso UPF file
140  rewind (unit=temp_unit,err=10,iomsg=errmsg)
141  read(temp_unit,*,err=10,iomsg=errmsg) testxml ! just a string, no relation to xml.
142  if(testxml(1:9)=='<PP_INFO>')then
143    useupf = 1
144  else
145    useupf = 0
146  end if
147 
148  close(unit=temp_unit,err=10,iomsg=errmsg)
149 
150  return
151 
152  ! Handle IO error
153 10 continue
154  MSG_ERROR(errmsg)
155 
156 end subroutine test_xml_xmlpaw_upf