TABLE OF CONTENTS


ABINIT/m_psxml2ab [ Modules ]

[ Top ] [ Modules ]

NAME

 m_psxml2ab

FUNCTION

  From a SIESTA XML format pseudopotential file
  convert to abinit internal datastructures for pspheader.

COPYRIGHT

 Copyright (C) 2005-2018 ABINIT group (MJV).
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

 psxml  = pseudopotential data structure

OUTPUT

 psphead = psp information structure
 atmsymb = atomic symbol

PARENTS

      inpspheads,pspatm_abinit

CHILDREN

      wrtout

SOURCE

 32 #if defined HAVE_CONFIG_H
 33 #include "config.h"
 34 #endif
 35 
 36 #include "abi_common.h"
 37 
 38 module m_psxml2ab
 39 
 40  use defs_basis
 41  use defs_datatypes
 42  use m_abicore
 43  use m_errors
 44 #ifdef HAVE_PSML
 45  use m_psml
 46  use m_psml_api
 47 #endif
 48 
 49  use m_fstrings,     only : yesno
 50 
 51 implicit none
 52 
 53 private
 54 
 55 #ifdef HAVE_PSML
 56 public :: psxml2abheader
 57 !public :: psxml2abfull
 58 
 59 CONTAINS
 60 
 61 subroutine psxml2abheader(psxmlfile, psphead, atmsymb, creator, iwrite)
 62 
 63 
 64 !This section has been created automatically by the script Abilint (TD).
 65 !Do not modify the following lines by hand.
 66 #undef ABI_FUNC
 67 #define ABI_FUNC 'psxml2abheader'
 68 !End of the abilint section
 69 
 70  implicit none
 71 
 72 !Arguments ------------------------------------
 73 !scalars
 74  integer, intent(in) :: iwrite
 75  character(len=fnlen), intent(in) :: psxmlfile
 76  type(pspheader_type),intent(inout) :: psphead
 77  character(len=3), intent(out) :: atmsymb
 78  character(len=30), intent(out) :: creator
 79 
 80 !Local variables-------------------------------
 81 !scalars
 82  character(len=500) :: flavor, frmt_str, label, message, relat, xc_name
 83  character(len=1) :: g1,g2
 84  integer :: dd,dm,dy
 85  integer :: il,lmm,ii
 86  integer :: nxc
 87  integer :: ishell, nvshells, shell_l, shell_n
 88  integer :: iproj, ll, ll_previous, nprojs, nprojsr, nprojso
 89  integer,parameter :: n1xccc_default=2501
 90  logical :: has_nlcc, has_spin
 91  real(dp) :: ekb
 92  type(ps_t) :: psxml
 93 !arrays
 94  integer, allocatable :: idx_sr(:), idx_so(:)
 95  real(dp),allocatable :: zeld(:),zelu(:)
 96 
 97 ! *********************************************************************
 98 
 99  call ps_destroy(psxml)
100  call psml_reader(psxmlfile, psxml, debug=.true.)
101 
102  psphead%pspdat = 0
103  call ps_Provenance_Get(psxml, 1, creator=creator, date=message)
104  read (message(1:10), '(I4,A1,I2,A1,I2)', err=10) &
105 &  dy, g1, dm, g2, dd
106  psphead%pspdat = MODULO(dy,100) * 10000 + dm * 100 + dd
107 
108 10 continue
109 
110  call ps_PseudoAtomSpec_Get(psxml, &
111 &  atomic_symbol=atmsymb, atomic_label=label, &
112 &  atomic_number=psphead%znuclpsp, z_pseudo=psphead%zionpsp, &
113 &  pseudo_flavor=flavor, relativity=relat, &
114 &  spin_dft=has_spin, core_corrections=has_nlcc)
115 
116  psphead%pspcod = 9
117 
118 ! impose libxc coding for pspxc
119  psphead%pspxc = 0
120  call ps_ExchangeCorrelation_Get(psxml, n_libxc_functionals=nxc)
121  do ii=1, min(2,nxc)
122    call ps_LibxcFunctional_Get(psxml, ii, code=dd)
123    psphead%pspxc = dd + psphead%pspxc*1000
124  end do
125  psphead%pspxc = -psphead%pspxc
126 
127  call ps_ValenceConfiguration_Get(psxml, nshells=nvshells)
128 
129  if (iwrite == 1) then
130    write (message,'(a,a)') '- psxml2ab: ps_PseudoFlavor ', trim(message)
131 !  call wrtout(ab_out,  message,'COLL')
132    call wrtout(std_out,  message,'COLL')
133    write (message,'(a,I5)') '- psxml2ab: ps_NLibxcFunctionals ', nxc
134 !  call wrtout(ab_out,  message,'COLL')
135    call wrtout(std_out,  message,'COLL')
136    write (message,'(a,I5)') '- psxml2ab: ps_NValenceShells ', nvshells
137 !  call wrtout(ab_out,  message,'COLL')
138    call wrtout(std_out,  message,'COLL')
139  end if
140 
141  ABI_ALLOCATE(zeld, (nvshells))
142  ABI_ALLOCATE(zelu, (nvshells))
143  zeld = zero
144  zelu = zero
145  frmt_str="(a"
146  do ishell = 1, nvshells
147    if (has_spin) then
148      call ps_ValenceShell_Get(psxml, ishell, n=shell_n, l=shell_l, &
149 &    occ_up=zelu(ishell), occ_down=zeld(ishell))
150    else
151      call ps_ValenceShell_Get(psxml, ishell, n=shell_n, l=shell_l, &
152 &    occupation=zeld(ishell))
153    end if
154    if (iwrite == 1) then
155      write (message,'(a,I5)') '- psxml2ab: ps_ValenceShellN ', shell_n
156 !    call wrtout(ab_out,  message,'COLL')
157      call wrtout(std_out,  message,'COLL')
158      write (message,'(a,I5)') '- psxml2ab: ps_ValenceShellL ', shell_l
159 !    call wrtout(ab_out,  message,'COLL')
160      call wrtout(std_out,  message,'COLL')
161    end if
162    frmt_str = trim(frmt_str) // ", F10.3"
163  end do
164  frmt_str = trim(frmt_str) // ")"
165  if (iwrite == 1) then
166    write (message,frmt_str) '- psxml2ab: zeld ', zeld
167 !  call wrtout(ab_out,  message,'COLL')
168    call wrtout(std_out,  message,'COLL')
169    write (message,frmt_str) '- psxml2ab: zelu ', zelu
170 !  call wrtout(ab_out,  message,'COLL')
171    call wrtout(std_out,  message,'COLL')
172  end if
173 
174 !    Find the number of projectors per angular momentum shell
175  nprojs = 0
176  nprojsr = 0
177  nprojso = 0
178  psphead%nproj(:)=0
179  call ps_NonlocalProjectors_Filter(psxml, set=SET_NONREL, number=nprojs)
180  if (nprojs > 0) then
181    call ps_NonlocalProjectors_Filter(psxml, set=SET_NONREL, indexes=idx_sr)
182    do iproj = 1, nprojs
183      if (iwrite == 1) then
184        write (message,'(a,2I5)') '- psxml2ab: iproj, idx for nonrel ', iproj, idx_sr(iproj)
185 !      call wrtout(ab_out,  message,'COLL')
186        call wrtout(std_out,  message,'COLL')
187      end if
188      call ps_Projector_Get(psxml, idx_sr(iproj), l=il)
189      psphead%nproj(il) = psphead%nproj(il) + 1
190    end do
191  else
192    call ps_NonlocalProjectors_Filter(psxml, set=SET_SREL, number=nprojsr)
193    if (nprojsr > 0) then
194      call ps_NonlocalProjectors_Filter(psxml, set=SET_SREL, indexes=idx_sr)
195      do iproj = 1, nprojsr
196        if (iwrite == 1) then
197          write (message,'(a,2I5)') '- psxml2ab: iproj, idx for srel ', iproj, idx_sr(iproj)
198 !        call wrtout(ab_out,  message,'COLL')
199          call wrtout(std_out,  message,'COLL')
200        end if
201        call ps_Projector_Get(psxml, idx_sr(iproj), l=il)
202        psphead%nproj(il) = psphead%nproj(il) + 1
203      end do
204    else
205      MSG_BUG('Your psml potential should have either scalar- or non- relativistic projectors')
206    end if
207  end if
208 
209  psphead%nprojso(:)=0
210  call ps_NonlocalProjectors_Filter(psxml, set=SET_SO, number=nprojso, &
211 &  indexes=idx_so)
212  do iproj = 1, nprojso
213    if (iwrite == 1) then
214      write (message,'(a,2I5)') '- psxml2ab: iproj, idx for soc ', iproj, idx_so(iproj)
215 !    call wrtout(ab_out,  message,'COLL')
216      call wrtout(std_out,  message,'COLL')
217    end if
218    call ps_Projector_Get(psxml, idx_so(iproj), l=il)
219    psphead%nprojso(il) = psphead%nprojso(il) + 1
220  end do
221  if (iwrite == 1) then
222    write (message,'(a,5I5)') '- psxml2ab: nproj ', psphead%nproj
223 !  call wrtout(ab_out,  message,'COLL')
224    call wrtout(std_out,  message,'COLL')
225    write (message,'(a,5I5)') '- psxml2ab: nprojso ', psphead%nprojso
226 !  call wrtout(ab_out,  message,'COLL')
227    call wrtout(std_out,  message,'COLL')
228  end if
229 
230  if( has_nlcc) then
231    psphead%xccc  = n1xccc_default
232  else
233    psphead%xccc  = 0
234  end if
235 
236  psphead%pspso = 0
237  if (sum(abs(psphead%nprojso(:))) > 0) psphead%pspso = 2
238 
239 
240  psphead%lmax = 0
241  if (iwrite == 1) then
242    write (message,'(a,I5)') '- psxml2ab: ps_Number_of_Projectors not relativistic ',&
243 &        nprojs
244    call wrtout(ab_out,  message,'COLL')
245    call wrtout(std_out,  message,'COLL')
246    write (message,'(a,I5)') '- psxml2ab: ps_Number_of_Projectors scalar relativistic ', &
247 &        nprojsr
248    call wrtout(ab_out,  message,'COLL')
249    call wrtout(std_out,  message,'COLL')
250  end if
251  ll_previous=-1
252  do iproj = 1, max(nprojs, nprojsr)
253    call ps_Projector_Get(psxml, idx_sr(iproj), l=ll)
254    if (iwrite == 1) then
255      if(ll/=ll_previous)then
256        write (message,'(a,I5)') '- psxml2ab: ps_Projector_L ', ll
257        call wrtout(ab_out,  message,'COLL')
258        call wrtout(std_out,  message,'COLL')
259        ll_previous=ll
260      endif
261      call ps_Projector_Get(psxml, idx_sr(iproj), ekb=ekb)
262      write (message,'(a,E20.10)') '- psxml2ab: ps_Projector_Ekb ', ekb
263      call wrtout(ab_out,  message,'COLL')
264      call wrtout(std_out,  message,'COLL')
265    end if
266    psphead%lmax = max( psphead%lmax, ll)
267  end do
268 
269  if (iwrite == 1) then
270    write (message,'(a,I5)') '- psxml2ab: ps_Number_of_Projectors SOC ', nprojso
271    call wrtout(ab_out,  message,'COLL')
272    call wrtout(std_out,  message,'COLL')
273  end if
274  ll_previous=-1
275  do iproj = 1, nprojso
276    call ps_Projector_Get(psxml, idx_so(iproj), l=ll)
277    if (iwrite == 1) then
278      if(ll/=ll_previous)then
279        write (message,'(a,I5)') '- psxml2ab: ps_Projector_L ', ll
280        call wrtout(ab_out,  message,'COLL')
281        call wrtout(std_out,  message,'COLL')
282        ll_previous=ll
283      endif
284      call ps_Projector_Get(psxml, idx_so(iproj), ekb=ekb)
285      write (message,'(a,E20.10)') '- psxml2ab: ps_Projector_Ekb ', ekb
286      call wrtout(ab_out,  message,'COLL')
287      call wrtout(std_out,  message,'COLL')
288    end if
289    psphead%lmax = max( psphead%lmax, ll)
290  end do
291 
292  lmm     = max( nprojs, nprojsr)
293  lmm     = max( lmm, nprojso)
294 
295  write(std_out,'(a,f5.1,a,i4,a,i4)' ) '- psxml2ab:  read the values zionpsp=',&
296 & psphead%zionpsp,' , pspcod=',psphead%pspcod,' , lmax=',psphead%lmax
297 
298  if ( iwrite .eq. 1 ) then
299 
300    write(message,'(a,a)') &
301 &   '- psxml2ab: Atomic Label:                      ', &
302 &  trim(label)
303    call wrtout(ab_out,  message,'COLL')
304    call wrtout(std_out,  message,'COLL')
305 
306    write(message,'(a,f12.5)') &
307 &   '- psxml2ab: Atomic Number:                     ', &
308 &   psphead%znuclpsp
309    call wrtout(ab_out,  message,'COLL')
310    call wrtout(std_out,  message,'COLL')
311 
312    write(message,'(a,f12.5)') &
313 &   '- psxml2ab: Valence charge:                    ', &
314 &   psphead%zionpsp
315    call wrtout(ab_out,  message,'COLL')
316    call wrtout(std_out,  message,'COLL')
317 
318    write(message,'(a,a)') &
319 &   '- psxml2ab: Pseudopotential generator code :    ', &
320 &   creator
321    call wrtout(ab_out,  message,'COLL')
322    call wrtout(std_out,  message,'COLL')
323 
324    write(message,'(a,i8)') &
325 &   '- psxml2ab: Date of pseudopotential generation: ', &
326 &   psphead%pspdat
327    call wrtout(ab_out,  message,'COLL')
328    call wrtout(std_out,  message,'COLL')
329 
330    write(message,'(a,a)') &
331 &   '- psxml2ab: Pseudopotential flavor:             ', &
332 &   trim(flavor)
333    call wrtout(ab_out,  message,'COLL')
334    call wrtout(std_out,  message,'COLL')
335 
336    do ii = 1, nxc
337      call ps_LibxcFunctional_Get(psxml, ii, name=xc_name, code=dd)
338      write(message,'(a,I4,2a," (",I4,")")') &
339 &     '- psxml2ab: Exchange-correlation functional ',ii,' :    ', &
340 &     trim(xc_name), dd
341      call wrtout(ab_out,  message,'COLL')
342      call wrtout(std_out,  message,'COLL')
343    end do
344 
345    write(message,'(2a)') &
346 &   '- psxml2ab: Relativistically generated pseudopotential (not necessarily SOC!):   ', &
347 &   trim(relat)
348    call wrtout(ab_out,  message,'COLL')
349    call wrtout(std_out,  message,'COLL')
350 
351    write(message,'(2a)') &
352 &   '- psxml2ab: Spin-polarized pseudopotential:     ', &
353 &   yesno(has_spin)
354    call wrtout(ab_out,  message,'COLL')
355    call wrtout(std_out,  message,'COLL')
356 
357    select case(has_nlcc)
358      case(.true.)
359        write(message, '(a)' ) &
360 &       '- psxml2ab: XC core correction read in from XML file.'
361      case(.false.)
362        write(message, '(a)' ) &
363 &       '- psxml2ab: No core corrections.'
364        call wrtout(ab_out,message,'COLL')
365        call wrtout(std_out,  message,'COLL')
366    end select
367  end if
368 
369  ABI_DEALLOCATE(zeld)
370  ABI_DEALLOCATE(zelu)
371 
372  if (allocated(idx_sr)) then
373    ABI_FREE_NOCOUNT(idx_sr)
374  end if
375  if (allocated(idx_so)) then
376    ABI_FREE_NOCOUNT(idx_so)
377  end if
378 
379  call ps_destroy(psxml)
380 
381 end subroutine psxml2abheader

ABINIT/psml_die [ Functions ]

[ Top ] [ Functions ]

NAME

 psml_die

FUNCTION

  auxiliary die function for calling libPSML. Needed at link time
  allows calling software to decide how fatal the PSML die call actually is.

COPYRIGHT

 Copyright (C) 2005-2018 ABINIT group (MJV).
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

   str = string with error message

OUTPUT

PARENTS

CHILDREN

SOURCE

416 subroutine psml_die(str)
417 
418   use m_errors
419 
420 !This section has been created automatically by the script Abilint (TD).
421 !Do not modify the following lines by hand.
422 #undef ABI_FUNC
423 #define ABI_FUNC 'psml_die'
424 !End of the abilint section
425 
426   implicit none
427 
428   character(len=*), intent(in) :: str
429 
430   MSG_BUG(str)
431 
432 end subroutine psml_die