Line data Source code
1 : !
2 : ! Copyright (C) 2020 Quantum ESPRESSO group
3 : ! This file is distributed under the terms of the
4 : ! GNU General Public License. See the file `License'
5 : ! in the root directory of the present distribution,
6 : ! or http://www.gnu.org/copyleft/gpl.txt .
7 : !
8 : !--------------------------------------------------------
9 :
10 :
11 : #if defined HAVE_CONFIG_H
12 : #include "config.h"
13 : #endif
14 :
15 : #include "abi_common.h"
16 :
17 : MODULE read_upf_new_module
18 : !-----------------------------------------------------
19 : !! this module contains the simplified code for reading
20 : !! pseudopotential files in either UPF v.2 or xml
21 : !
22 : USE xmltools
23 : !USE upf_kinds, ONLY: dp
24 : use defs_basis, only : dp, std_out
25 : USE pseudo_types, ONLY: pseudo_upf, pseudo_config
26 :
27 : IMPLICIT NONE
28 :
29 : !
30 : LOGICAL :: v2
31 : !! true if UPF v.2 version, false if new UPF with xml schema
32 : INTEGER :: iun
33 : !! unit for reading data
34 : !
35 : PUBLIC
36 : !
37 : CONTAINS
38 : !
39 : !------------------------------------------------+
40 16 : SUBROUTINE read_upf_new (filename, upf, ierr) !
41 : !---------------------------------------------+
42 : !! Reads pseudopotential in UPF format (either v.2 or upf_schema).
43 : !! Derived-type variable *upf* store in output the data read from file.
44 : !! File *filename* is opened and closed inside the routine
45 : !
46 : !IMPLICIT NONE
47 : CHARACTER(len=*), INTENT(IN) :: filename
48 : !! i/o filename
49 : TYPE(pseudo_upf),INTENT(OUT) :: upf
50 : !! the derived type storing the pseudo data
51 : INTEGER, INTENT(OUT) :: ierr
52 : !! ierr=0 : xml schema, ierr=-2: UPF v.2
53 : !! ierr=-81: error reading PP file
54 : !
55 16 : iun = xml_open_file ( filename )
56 16 : IF ( iun == -1 ) CALL upf_error('read_upf', 'cannot open file',1)
57 16 : call xmlr_opentag ( 'qe_pp:pseudo', IERR = ierr )
58 16 : if ( ierr == 0 ) then
59 0 : v2 =.false.
60 16 : else if ( ierr == 1 ) then
61 16 : rewind (iun)
62 16 : call xmlr_opentag ( 'UPF', IERR = ierr )
63 16 : if ( ierr == 0 ) then
64 16 : v2 =.true.
65 16 : ierr = -2
66 16 : CALL get_attr ( 'version', upf%nv )
67 : end if
68 : end if
69 16 : if ( ierr /= 0 .and. ierr /= -2 ) then
70 0 : call xml_closefile( )
71 0 : ierr = -81
72 0 : return
73 : end if
74 : !
75 : ! The header sections differ a lot between UPF v.2 and UPF with schema
76 : !
77 16 : IF ( v2 ) THEN
78 16 : CALL read_pp_header_v2 ( upf )
79 : ELSE
80 0 : CALL read_pp_header_schema ( upf )
81 : END IF
82 : ! compatibility
83 16 : upf%is_gth = .false.
84 16 : upf%is_multiproj = .true.
85 16 : upf%with_metagga_info = .false.
86 : !
87 : ! From here on the format of v2 and schema do not differ much:
88 : ! the most frequent difference is capitalization of tags
89 : ! (see function capitalize_if_v2)
90 : !
91 16 : CALL read_pp_mesh ( upf )
92 : !
93 48 : allocate ( upf%rho_atc(upf%mesh) )
94 16 : IF(upf%nlcc) then
95 : CALL xmlr_readtag( capitalize_if_v2('pp_nlcc'), &
96 16 : upf%rho_atc(:) )
97 : else
98 0 : upf%rho_atc(:) = 0.0_dp
99 : end if
100 16 : IF( .NOT. upf%tcoulombp) then
101 48 : allocate ( upf%vloc(upf%mesh) )
102 : CALL xmlr_readtag( capitalize_if_v2('pp_local'), &
103 16 : upf%vloc(:), ierr )
104 : !
105 : ! existing PP files may have pp_nlcc first, pp_local later,
106 : ! but also the other way round - check that everything was right
107 : !
108 16 : if ( ierr /= 0 ) then
109 0 : ierr = -81
110 0 : return
111 : end if
112 : end if
113 : !
114 16 : CALL read_pp_semilocal ( upf )
115 : !
116 16 : CALL read_pp_nonlocal ( upf )
117 : !
118 16 : CALL read_pp_pswfc ( upf )
119 : !
120 16 : CALL read_pp_full_wfc ( upf )
121 : !
122 48 : allocate( upf%rho_at(1:upf%mesh) )
123 : CALL xmlr_readtag( capitalize_if_v2('pp_rhoatom'), &
124 16 : upf%rho_at(1:upf%mesh) )
125 : !
126 48 : allocate( upf%tau_at(1:upf%mesh) )
127 : CALL xmlr_readtag( capitalize_if_v2('pp_tauatom'), &
128 16 : upf%tau_at(1:upf%mesh) )
129 : !
130 48 : allocate( upf%tau_mod(1:upf%mesh) )
131 : CALL xmlr_readtag( capitalize_if_v2('pp_taumod'), &
132 16 : upf%tau_mod(1:upf%mesh) )
133 : !
134 16 : CALL read_pp_spinorb ( upf )
135 : !
136 16 : CALL read_pp_paw ( upf )
137 : !
138 16 : CALL read_pp_gipaw ( upf )
139 : !
140 : ! close initial tag, qe_pp:pseudo or UPF
141 : !
142 16 : CALL xmlr_closetag ( )
143 : !
144 16 : CALL xml_closefile ( )
145 : !
146 16 : END SUBROUTINE read_upf_new
147 : !
148 192 : FUNCTION capitalize_if_v2 ( strin ) RESULT ( strout )
149 : !
150 : ! returns a capitalized string for UPF v.2, the same string otherwise
151 : ! (UPF v.2 uses capitalized tags, UPF with schema use lowercase)
152 : !
153 : USE upf_utils, ONLY: capital
154 : !IMPLICIT NONE
155 : CHARACTER(LEN=*) :: strin
156 : !
157 : INTEGER :: n
158 : CHARACTER(LEN=:), ALLOCATABLE :: strout
159 : !
160 192 : IF ( v2 ) THEN
161 192 : strout = ''
162 1712 : DO n = 1,LEN_TRIM(strin)
163 1712 : strout = strout // capital(strin(n:n))
164 : END DO
165 : ELSE
166 0 : strout = TRIM(strin)
167 : END IF
168 : !
169 192 : END FUNCTION capitalize_if_v2
170 : !--------------------------------------------------------
171 0 : SUBROUTINE read_pp_header_schema ( upf )
172 : !--------------------------------------------------------
173 : !
174 : !IMPLICIT NONE
175 : TYPE(pseudo_upf), INTENT(INOUT) :: upf ! the pseudo data
176 : !
177 0 : CALL xmlr_opentag( capitalize_if_v2('pp_header') )
178 : !
179 0 : CALL xmlr_readtag( 'element', upf%psd )
180 0 : CALL xmlr_readtag( 'z_valence', upf%zp )
181 0 : CALL xmlr_readtag( 'type', upf%typ )
182 0 : CALL xmlr_readtag( 'functional', upf%dft )
183 0 : call remove_non_ascii(upf%dft)
184 :
185 0 : CALL xmlr_readtag( 'relativistic', upf%rel )
186 0 : CALL xmlr_readtag( 'is_ultrasoft', upf%tvanp )
187 0 : CALL xmlr_readtag( 'is_paw', upf%tpawp )
188 0 : CALL xmlr_readtag( 'is_coulomb', upf%tcoulombp )
189 0 : CALL xmlr_readtag( 'with_metagga_info', upf%with_metagga_info )
190 0 : CALL xmlr_readtag( 'has_so', upf%has_so )
191 0 : CALL xmlr_readtag( 'has_wfc', upf%has_wfc )
192 0 : CALL xmlr_readtag( 'has_gipaw', upf%has_gipaw )
193 0 : CALL xmlr_readtag( 'paw_as_gipaw', upf%paw_as_gipaw)
194 0 : CALL xmlr_readtag( 'core_correction', upf%nlcc)
195 0 : CALL xmlr_readtag( 'total_psenergy', upf%etotps )
196 0 : CALL xmlr_readtag( 'wfc_cutoff', upf%ecutwfc )
197 0 : CALL xmlr_readtag( 'rho_cutoff', upf%ecutrho )
198 0 : CALL xmlr_readtag( 'l_max', upf%lmax )
199 0 : CALL xmlr_readtag( 'l_max_rho', upf%lmax_rho )
200 0 : CALL xmlr_readtag( 'l_local', upf%lloc )
201 0 : CALL xmlr_readtag( 'mesh_size', upf%mesh )
202 0 : CALL xmlr_readtag( 'number_of_wfc', upf%nwfc )
203 0 : CALL xmlr_readtag( 'number_of_proj', upf%nbeta )
204 : !
205 0 : CALL xmlr_closetag( )
206 : !
207 0 : END SUBROUTINE read_pp_header_schema
208 : !
209 : !--------------------------------------------------------
210 16 : SUBROUTINE read_pp_header_v2 ( upf )
211 : !--------------------------------------------------------
212 : !
213 : !IMPLICIT NONE
214 : TYPE(pseudo_upf), INTENT(INOUT) :: upf ! the pseudo data
215 : !
216 : CHARACTER(LEN=1) :: dummy
217 : !
218 16 : CALL xmlr_readtag ( capitalize_if_v2('pp_header'), dummy )
219 16 : CALL get_attr ('generated', upf%generated)
220 16 : CALL get_attr ('author', upf%author)
221 16 : CALL get_attr ('date', upf%date)
222 16 : CALL get_attr ('comment', upf%comment)
223 16 : CALL get_attr ('element', upf%psd)
224 16 : CALL get_attr ('pseudo_type', upf%typ)
225 16 : CALL get_attr ('relativistic', upf%rel)
226 16 : CALL get_attr ('is_ultrasoft', upf%tvanp)
227 16 : CALL get_attr ('is_paw', upf%tpawp)
228 16 : CALL get_attr ('is_coulomb', upf%tcoulombp)
229 16 : CALL get_attr ('with_metagga_info', upf%with_metagga_info)
230 16 : CALL get_attr ('has_so', upf%has_so)
231 16 : CALL get_attr ('has_wfc', upf%has_wfc)
232 16 : CALL get_attr ('has_gipaw', upf%has_gipaw)
233 16 : CALL get_attr ('paw_as_gipaw', upf%paw_as_gipaw)
234 16 : CALL get_attr ('core_correction', upf%nlcc)
235 16 : CALL get_attr ('functional', upf%dft)
236 16 : call remove_non_ascii(upf%dft)
237 16 : CALL get_attr ('z_valence', upf%zp)
238 16 : CALL get_attr ('total_psenergy', upf%etotps)
239 16 : CALL get_attr ('wfc_cutoff', upf%ecutwfc)
240 16 : CALL get_attr ('rho_cutoff', upf%ecutrho)
241 16 : CALL get_attr ('l_max', upf%lmax)
242 16 : CALL get_attr ('l_max_rho', upf%lmax_rho)
243 16 : CALL get_attr ('l_local', upf%lloc)
244 16 : CALL get_attr ('mesh_size', upf%mesh)
245 16 : CALL get_attr ('number_of_wfc', upf%nwfc)
246 16 : CALL get_attr ('number_of_proj', upf%nbeta )
247 : !
248 16 : END SUBROUTINE read_pp_header_v2
249 : !
250 : !--------------------------------------------------------
251 16 : SUBROUTINE read_pp_mesh ( upf )
252 : !--------------------------------------------------------
253 : !
254 : !IMPLICIT NONE
255 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
256 : integer :: mesh
257 : !
258 16 : CALL xmlr_opentag( capitalize_if_v2('pp_mesh') )
259 16 : CALL get_attr ( 'mesh', mesh )
260 16 : if ( mesh == 0 ) THEN
261 : #if defined (__debug)
262 : call upf_error('read_pp_mesh',&
263 : 'mesh size missing, using the one in header',-1)
264 : #else
265 : continue
266 : #endif
267 0 : else if ( mesh /= upf%mesh ) THEN
268 : call upf_error('read_pp_mesh',&
269 0 : 'mismatch in mesh size, discarding the one in header',-1)
270 0 : upf%mesh = mesh
271 : end if
272 16 : CALL get_attr ( 'dx' , upf%dx )
273 16 : CALL get_attr ( 'xmin', upf%xmin )
274 16 : CALL get_attr ( 'rmax', upf%rmax )
275 16 : CALL get_attr ( 'zmesh', upf%zmesh )
276 48 : allocate ( upf%r(1:upf%mesh) )
277 16 : CALL xmlr_readtag( capitalize_if_v2('pp_r'), upf%r(1:upf%mesh) )
278 48 : allocate ( upf%rab(1:upf%mesh) )
279 16 : CALL xmlr_readtag( capitalize_if_v2('pp_rab'), upf%rab(1:upf%mesh) )
280 : !
281 16 : CALL xmlr_closetag( ) ! end pp_mesh
282 : !
283 16 : END SUBROUTINE read_pp_mesh
284 : !
285 : !--------------------------------------------------------
286 16 : SUBROUTINE read_pp_semilocal ( upf )
287 : !--------------------------------------------------------
288 : !
289 : IMPLICIT NONE
290 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
291 : !
292 : INTEGER :: nb, ind, l, j, ierr
293 : CHARACTER(LEN=8) :: tag
294 16 : real(dp), allocatable :: vnl(:)
295 : !
296 16 : IF ( upf%typ == "SL" ) THEN
297 : !
298 0 : IF ( upf%has_so ) then
299 0 : ALLOCATE(upf%vnl(upf%mesh,0:upf%lmax,2))
300 : else
301 0 : ALLOCATE(upf%vnl(upf%mesh,0:upf%lmax,1))
302 : end if
303 0 : allocate ( vnl(1:upf%mesh) )
304 0 : CALL xmlr_opentag( capitalize_if_v2('pp_semilocal') )
305 : !
306 0 : tag = 'vnl'
307 0 : DO nb = 1,upf%nbeta
308 0 : IF ( v2 ) THEN
309 : ! NOTA BENE: v2 format follows available PP files, written
310 : ! using original write_upf_v2; not FoX-based write_upf_v2
311 0 : IF ( nb - 1 == upf%lloc ) CYCLE
312 0 : tag = 'PP_VNL.'//i2c(nb-1)
313 : END IF
314 0 : CALL xmlr_readtag( tag, vnl, ierr )
315 0 : if ( ierr /= 0 ) &
316 0 : call upf_error('read_pp_semilocal','error reading SL PPs',1)
317 0 : CALL get_attr ( 'l', l)
318 0 : ind = 1
319 0 : IF ( upf%has_so ) then
320 0 : CALL get_attr ( 'j', j)
321 0 : IF ( l > 0 .AND. ABS(j-l-0.5_dp) < 0.001_dp ) ind = 2
322 : ! FIXME: what about spin-orbit case for v.2 upf?
323 0 : if ( v2 ) &
324 0 : call upf_error('read_pp_semilocal','check spin-orbit',1)
325 : END IF
326 0 : upf%vnl(:,l,ind) = vnl(:)
327 : END DO
328 0 : deallocate ( vnl )
329 : !
330 0 : CALL xmlr_closetag( ) ! end pp_semilocal
331 : !
332 : END IF
333 : !
334 16 : END SUBROUTINE read_pp_semilocal
335 : !
336 : !--------------------------------------------------------
337 16 : SUBROUTINE read_pp_nonlocal ( upf )
338 : !--------------------------------------------------------
339 : !
340 : IMPLICIT NONE
341 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
342 : !
343 : LOGICAL :: isnull
344 : INTEGER :: nb, l, l_, ln, lm, mb, nmb !ind,
345 : CHARACTER(LEN=15) :: tag
346 16 : REAL(dp), ALLOCATABLE :: aux(:)
347 : !
348 16 : nb = upf%nbeta
349 16 : IF ( nb == 0 ) nb = 1
350 64 : ALLOCATE (upf%beta(upf%mesh,nb) )
351 0 : ALLOCATE (upf%els_beta(nb), &
352 0 : upf%lll(nb), &
353 0 : upf%kbeta(nb), &
354 0 : upf%rcut(nb), &
355 0 : upf%rcutus(nb), &
356 0 : upf%dion(nb,nb), &
357 208 : upf%qqq(nb,nb) )
358 : !
359 28 : IF (upf%has_so) ALLOCATE( upf%jjj(upf%nbeta))
360 : !
361 16 : IF ( upf%nbeta == 0 ) THEN
362 0 : upf%nqf = 0
363 0 : upf%nqlc= 0
364 0 : upf%kkbeta = 0
365 0 : upf%qqq_eps=-1.0_dp
366 : RETURN
367 : END IF
368 : !
369 16 : CALL xmlr_opentag( capitalize_if_v2('pp_nonlocal') )
370 : !
371 133 : DO nb = 1,upf%nbeta
372 : !
373 117 : IF ( v2 ) THEN
374 117 : tag = 'PP_BETA.'//i2c(nb)
375 : ELSE
376 0 : tag = 'pp_beta'
377 : END IF
378 117 : CALL xmlr_readtag( tag, upf%beta(1:upf%mesh,nb) )
379 117 : CALL get_attr('index', mb)
380 : ! not-so-strict test: index is absent or incorrect in some UPF v.2 files
381 117 : IF ( .NOT. v2 .AND. nb /= mb ) &
382 0 : CALL upf_error('read_pp_nonlocal','mismatch',nb)
383 117 : CALL get_attr('label', upf%els_beta(nb))
384 117 : CALL get_attr('angular_momentum', upf%lll(nb))
385 117 : IF ( .NOT. v2 .AND. upf%has_so ) &
386 0 : CALL get_attr('tot_ang_mom', upf%jjj(nb))
387 117 : CALL get_attr('cutoff_radius_index', upf%kbeta(nb))
388 117 : CALL get_attr('cutoff_radius', upf%rcut(nb))
389 250 : CALL get_attr('ultrasoft_cutoff_radius', upf%rcutus(nb))
390 : !
391 : END DO
392 : !
393 : ! pp_dij (D_lm matrix)
394 : !
395 16 : CALL xmlr_readtag ( capitalize_if_v2 ('pp_dij'), upf%dion )
396 : !
397 : ! pp_augmentation
398 : !
399 16 : IF (upf%tvanp .or. upf%tpawp) THEN
400 0 : CALL xmlr_opentag( capitalize_if_v2('pp_augmentation') )
401 : !
402 0 : IF ( v2 ) THEN
403 0 : CALL get_attr ( 'q_with_l', upf%q_with_l )
404 0 : CALL get_attr ( 'nqf', upf%nqf )
405 0 : CALL get_attr ( 'nqlc', upf%nqlc )
406 0 : IF (upf%tpawp) THEN
407 0 : CALL get_attr ( 'shape', upf%paw%augshape )
408 0 : CALL get_attr ( 'cutoff_r', upf%paw%raug )
409 0 : CALL get_attr ( 'cutoff_r_index', upf%paw%iraug )
410 0 : CALL get_attr ( 'augmentation_epsilon', upf%qqq_eps )
411 0 : CALL get_attr ( 'l_max_aug', upf%paw%lmax_aug )
412 : ENDIF
413 : ELSE
414 0 : CALL xmlr_readtag( 'q_with_l', upf%q_with_l )
415 0 : CALL xmlr_readtag( 'nqf', upf%nqf )
416 0 : CALL xmlr_readtag( 'nqlc', upf%nqlc )
417 0 : IF (upf%tpawp) THEN
418 0 : CALL xmlr_readtag( 'shape', upf%paw%augshape )
419 0 : CALL xmlr_readtag( 'cutoff_r', upf%paw%raug )
420 0 : CALL xmlr_readtag( 'cutoff_r_index', upf%paw%iraug )
421 0 : CALL xmlr_readtag( 'augmentation_epsilon', upf%qqq_eps )
422 0 : CALL xmlr_readtag( 'l_max_aug', upf%paw%lmax_aug )
423 : ENDIF
424 : ENDIF
425 : !
426 0 : CALL xmlr_readtag( capitalize_if_v2('pp_q'), upf%qqq )
427 : !
428 0 : IF ( upf%tpawp ) THEN
429 0 : ALLOCATE ( upf%paw%augmom(1:upf%nbeta,1:upf%nbeta,0:2*upf%lmax) )
430 0 : CALL xmlr_readtag( capitalize_if_v2('pp_multipoles'), upf%paw%augmom )
431 : ENDIF
432 : !
433 : ! read polinomial coefficients for Q_ij expansion at small radius
434 : !
435 0 : IF ( upf%nqlc == 0 ) upf%nqlc = 2*upf%lmax+1
436 0 : ALLOCATE( upf%rinner( upf%nqlc ) )
437 0 : IF ( v2 .AND. upf%nqf > 0) THEN
438 0 : ALLOCATE ( upf%qfcoef(upf%nqf, upf%nqlc, upf%nbeta, upf%nbeta) )
439 0 : CALL xmlr_opentag('PP_QFCOEF')
440 0 : READ(iun,*) upf%qfcoef
441 0 : CALL xmlr_closetag ()
442 0 : CALL xmlr_readtag('PP_RINNER',upf%rinner)
443 0 : ELSE IF ( upf%nqf == 0 ) THEN
444 0 : ALLOCATE( upf%qfcoef(1,1,1,1) )
445 0 : upf%qfcoef =0.0_dp
446 : ENDIF
447 : !
448 : ! Read augmentation charge Q_ij
449 : !
450 0 : IF( upf%q_with_l ) THEN
451 0 : ALLOCATE( upf%qfuncl(upf%mesh,upf%nbeta*(upf%nbeta+1)/2,0:2*upf%lmax) )
452 0 : upf%qfuncl(:,:,:) = 0.0_dp
453 : ! NOTE: it would be wiser to dimension qfuncl as (:,:,0:upf%lmax)
454 : ! and store the q_l(r) with index l=L/2 (see loop_on_l below)
455 : ! This would save some storage and avoid "holes" in the array
456 : ! that may be a source of trouble if not initialized to zero
457 : ELSE
458 0 : ALLOCATE ( upf%qfunc(upf%mesh,upf%nbeta*(upf%nbeta+1)/2) )
459 0 : upf%qfunc (:,:) = 0.0_dp
460 : END IF
461 0 : ALLOCATE ( aux(upf%mesh) )
462 0 : loop_on_nb: DO nb = 1,upf%nbeta
463 0 : ln = upf%lll(nb)
464 0 : loop_on_mb: DO mb = nb,upf%nbeta
465 0 : lm = upf%lll(mb)
466 0 : IF( upf%q_with_l ) THEN
467 0 : loop_on_l: DO l = abs(ln-lm),ln+lm,2 ! only even terms
468 0 : isnull = .FALSE.
469 0 : IF( upf%tpawp ) isnull = (abs(upf%paw%augmom(nb,mb,l)) < upf%qqq_eps)
470 0 : IF(isnull) CYCLE loop_on_l
471 0 : IF ( v2 ) THEN
472 0 : tag = 'PP_QIJL.'//i2c(nb)//'.'//i2c(mb)//'.'//i2c(l)
473 : ELSE
474 0 : tag = 'pp_qijl'
475 : END IF
476 0 : CALL xmlr_readtag( tag, aux )
477 0 : CALL get_attr ('composite_index', nmb)
478 0 : IF ( nmb /= mb*(mb-1)/2 + nb ) &
479 0 : CALL upf_error ('read_pp_nonlocal','mismatch',1)
480 0 : CALL get_attr ('angular_momentum', l_)
481 0 : IF ( l /= l_ ) CALL upf_error ('read_pp_nonlocal','mismatch',2)
482 0 : upf%qfuncl(:,nmb,l) = aux(:)
483 0 : IF (upf%tpawp) upf%qfuncl(upf%paw%iraug+1:,nmb,l) = 0._DP
484 : ENDDO loop_on_l
485 : ELSE
486 0 : isnull = .FALSE.
487 0 : IF ( upf%tpawp ) isnull = ( abs(upf%qqq(nb,mb)) < upf%qqq_eps )
488 0 : IF (isnull) CYCLE loop_on_mb
489 0 : IF ( v2 ) THEN
490 0 : tag = 'PP_QIJ.'//i2c(nb)//'.'//i2c(mb)
491 : ELSE
492 0 : tag = 'pp_qij'
493 : END IF
494 0 : CALL xmlr_readtag( tag, aux )
495 0 : CALL get_attr ('composite_index', nmb)
496 0 : IF ( nmb /= mb*(mb-1)/2 + nb ) &
497 0 : CALL upf_error ('read_pp_nonlocal','mismatch',3)
498 0 : upf%qfunc(:,nmb) = aux(:)
499 : !
500 : ENDIF
501 : ENDDO loop_on_mb
502 : ENDDO loop_on_nb
503 : !
504 0 : DEALLOCATE (aux)
505 0 : CALL xmlr_closetag( ) ! end pp_augmentation
506 : !
507 : END IF
508 16 : CALL xmlr_closetag( ) ! end pp_nonlocal
509 : !
510 : ! Maximum radius of beta projector: outer radius to integrate
511 133 : upf%kkbeta = MAXVAL(upf%kbeta(1:upf%nbeta))
512 : ! For PAW, augmentation charge may extend a bit further:
513 16 : IF(upf%tpawp) upf%kkbeta = MAX(upf%kkbeta, upf%paw%iraug)
514 : !
515 16 : END SUBROUTINE read_pp_nonlocal
516 : !
517 : !--------------------------------------------------------
518 16 : SUBROUTINE read_pp_pswfc ( upf )
519 : !--------------------------------------------------------
520 : !
521 : !IMPLICIT NONE
522 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
523 : !
524 : INTEGER :: nw, ind !, l
525 : CHARACTER(LEN=8) :: tag
526 : !
527 64 : allocate ( upf%chi(1:upf%mesh,upf%nwfc) )
528 0 : allocate ( upf%els(upf%nwfc), &
529 0 : upf%oc(upf%nwfc), &
530 0 : upf%lchi(upf%nwfc), &
531 0 : upf%nchi(upf%nwfc), &
532 0 : upf%rcut_chi(upf%nwfc), &
533 0 : upf%rcutus_chi(upf%nwfc), &
534 160 : upf%epseu(upf%nwfc) )
535 16 : IF ( upf%has_so ) THEN
536 12 : allocate ( upf%nn(upf%nwfc) )
537 12 : allocate ( upf%jchi(upf%nwfc) )
538 : END IF
539 : !
540 16 : CALL xmlr_opentag( capitalize_if_v2('pp_pswfc') )
541 70 : DO nw=1,upf%nwfc
542 54 : IF ( v2 ) THEN
543 54 : tag = 'PP_CHI.'//i2c(nw)
544 : ELSE
545 0 : tag = 'pp_chi'
546 : END IF
547 54 : CALL xmlr_readtag( tag, upf%chi(1:upf%mesh,nw) )
548 54 : call get_attr('index', ind)
549 : ! not-so-strict test: index is absent or incorrect in some UPF v.2 files
550 54 : if ( .NOT. v2 .AND. ind /= nw ) &
551 0 : call upf_error('read_pp_pswfc','mismatch reading PSWFC', nw)
552 54 : call get_attr( 'label', upf%els(nw) )
553 54 : call get_attr( 'l', upf%lchi(nw) )
554 54 : IF ( .not. v2 .and. upf%has_so ) THEN
555 0 : call get_attr( 'nn', upf%nn(nw) )
556 0 : call get_attr( 'jchi', upf%jchi(nw) )
557 : END IF
558 54 : call get_attr( 'occupation', upf%oc(nw) )
559 54 : call get_attr( 'n', upf%nchi(nw) )
560 54 : call get_attr( 'pseudo_energy', upf%epseu(nw) )
561 54 : call get_attr( 'cutoff_radius', upf%rcut_chi(nw) )
562 124 : call get_attr( 'ultrasoft_cutoff_radius', upf%rcutus_chi(nw) )
563 : END DO
564 16 : CALL xmlr_closetag( ) ! end pp_pswfc
565 : !
566 16 : END SUBROUTINE read_pp_pswfc
567 : !
568 : !--------------------------------------------------------
569 16 : SUBROUTINE read_pp_full_wfc ( upf )
570 : !--------------------------------------------------------
571 : !
572 : !IMPLICIT NONE
573 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
574 : !
575 : INTEGER :: nb, mb
576 : CHARACTER(LEN=15) :: tag
577 : !
578 16 : IF ( upf%has_wfc ) THEN
579 : !
580 0 : ALLOCATE (upf%aewfc(1:upf%mesh,upf%nbeta) )
581 0 : CALL xmlr_opentag( capitalize_if_v2('pp_full_wfc') )
582 : !
583 0 : DO nb = 1, upf%nbeta
584 0 : IF ( v2 ) THEN
585 0 : tag = 'PP_AEWFC.'//i2c(nb)
586 : ELSE
587 0 : tag = 'pp_aewfc'
588 : END IF
589 0 : CALL xmlr_readtag( tag, upf%aewfc(1:upf%mesh,nb) )
590 0 : CALL get_attr ('index',mb)
591 : ! not-so-strict test (and two more below):
592 : ! index may be absent or incorrect in some UPF v.2 files
593 0 : IF ( .NOT. v2 .AND. nb /= mb ) CALL upf_error('read_pp_full_wfc','mismatch',1)
594 : END DO
595 : !
596 0 : IF ( upf%has_so .AND. upf%tpawp ) THEN
597 0 : ALLOCATE (upf%paw%aewfc_rel(1:upf%mesh,upf%nbeta) )
598 0 : DO nb = 1, upf%nbeta
599 0 : IF ( v2 ) THEN
600 0 : tag = 'PP_AEWFC_REL.'//i2c(nb)
601 : ELSE
602 0 : tag = 'pp_aewfc_rel'
603 : END IF
604 0 : CALL xmlr_readtag(tag, upf%paw%aewfc_rel(1:upf%mesh,nb) )
605 0 : CALL get_attr ('index',mb)
606 0 : IF ( .NOT. v2 .AND. nb /= mb ) CALL upf_error('read_pp_full_wfc','mismatch',2)
607 : END DO
608 : END IF
609 : !
610 0 : ALLOCATE (upf%pswfc(1:upf%mesh,upf%nbeta) )
611 0 : DO nb = 1, upf%nbeta
612 0 : IF ( v2 ) THEN
613 0 : tag = 'PP_PSWFC.'//i2c(nb)
614 : ELSE
615 0 : tag = 'pp_pswfc'
616 : END IF
617 0 : CALL xmlr_readtag(tag, upf%pswfc(1:upf%mesh,nb) )
618 0 : CALL get_attr ('index',mb)
619 0 : IF ( .NOT. v2 .AND. nb /= mb ) CALL upf_error('read_pp_full_wfc','mismatch',3)
620 : END DO
621 : !
622 0 : CALL xmlr_closetag( )
623 : !
624 : END IF
625 : !
626 16 : END SUBROUTINE read_pp_full_wfc
627 : !
628 : !--------------------------------------------------------
629 16 : SUBROUTINE read_pp_spinorb ( upf )
630 : !--------------------------------------------------------
631 : !
632 : !IMPLICIT NONE
633 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
634 : INTEGER :: nw, nb, ierr
635 : CHARACTER(LEN=1) :: dummy
636 : !
637 16 : IF ( .NOT. v2 .OR. .NOT. upf%has_so ) RETURN
638 : !
639 6 : CALL xmlr_opentag( 'PP_SPIN_ORB' )
640 36 : DO nw = 1,upf%nwfc
641 30 : CALL xmlr_readtag( 'PP_RELWFC.'//i2c(nw), dummy )
642 30 : CALL get_attr( 'index' , nb )
643 : ! not-so-strict test: index absent or incorrect in some UPF v.2 files
644 30 : IF ( .NOT. v2 .AND. nb /= nw ) CALL upf_error('read_pp_spinorb','mismatch',1)
645 30 : CALL get_attr( 'nn', upf%nn(nw) )
646 66 : CALL get_attr( 'jchi', upf%jchi(nw) )
647 : !
648 : ! the following data is already known and was not read in old versions
649 : ! of UPF-reading code. upf%oc is actually missing in some UPF files:
650 : ! reading it here may spoil the value read earlier and break DFT+U
651 : !
652 : ! CALL get_attr( 'lchi', upf%lchi(nw) )
653 : ! CALL get_attr( 'els', upf%els(nw) )
654 : ! CALL get_attr( 'oc', upf%oc(nw) )
655 : ENDDO
656 : !
657 66 : DO nb = 1,upf%nbeta
658 60 : CALL xmlr_readtag( 'PP_RELBETA.'//i2c(nb), dummy, ierr )
659 : !
660 : ! existing PP files may have pp_relbeta first, pp_relwfc later,
661 : ! but also the other way round - check that everything was right
662 : !
663 60 : if ( ierr > 0 ) then
664 : ierr = -81
665 : return
666 : end if
667 60 : CALL get_attr( 'index' , nw )
668 60 : IF ( .NOT.v2 .AND. nb /= nw ) CALL upf_error('read_pp_spinorb','mismatch',2)
669 60 : CALL get_attr( 'lll', upf%lll(nb) )
670 126 : CALL get_attr( 'jjj', upf%jjj(nb) )
671 : ENDDO
672 6 : CALL xmlr_closetag () ! end pp_spin_orb
673 : !
674 : END SUBROUTINE read_pp_spinorb
675 : !
676 : !--------------------------------------------------------
677 16 : SUBROUTINE read_pp_paw ( upf )
678 : !--------------------------------------------------------
679 : !
680 : !IMPLICIT NONE
681 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
682 : INTEGER :: nb, mb
683 : !
684 16 : IF ( .NOT. upf%tpawp ) RETURN
685 : !
686 0 : CALL xmlr_opentag( capitalize_if_v2('pp_paw') )
687 0 : CALL get_attr ('paw_data_format', upf%paw_data_format)
688 0 : CALL get_attr ('core_energy', upf%paw%core_energy)
689 : ! Full occupation (not only > 0 ones)
690 0 : ALLOCATE (upf%paw%oc(upf%nbeta) )
691 0 : ALLOCATE (upf%paw%ae_rho_atc(upf%mesh) )
692 0 : ALLOCATE (upf%paw%ae_vloc(upf%mesh) )
693 : CALL xmlr_readtag( capitalize_if_v2('pp_occupations'), &
694 0 : upf%paw%oc(1:upf%nbeta) )
695 : ! All-electron core charge
696 : CALL xmlr_readtag( capitalize_if_v2('pp_ae_nlcc'), &
697 0 : upf%paw%ae_rho_atc(1:upf%mesh) )
698 : ! All-electron local potential
699 : CALL xmlr_readtag( capitalize_if_v2('pp_ae_vloc'), &
700 0 : upf%paw%ae_vloc(1:upf%mesh) )
701 0 : CALL xmlr_closetag () ! end pp_paw
702 : !
703 0 : ALLOCATE(upf%paw%pfunc(upf%mesh, upf%nbeta,upf%nbeta) )
704 0 : upf%paw%pfunc(:,:,:) = 0._dp
705 0 : IF (upf%has_so) THEN
706 0 : ALLOCATE(upf%paw%pfunc_rel(upf%mesh, upf%nbeta,upf%nbeta) )
707 0 : upf%paw%pfunc_rel(:,:,:) = 0._dp
708 : ENDIF
709 0 : DO nb=1,upf%nbeta
710 0 : DO mb=1,nb
711 : upf%paw%pfunc (1:upf%mesh, nb, mb) = &
712 0 : upf%aewfc(1:upf%mesh, nb) * upf%aewfc(1:upf%mesh, mb)
713 0 : IF (upf%has_so) THEN
714 : upf%paw%pfunc_rel (1:upf%paw%iraug, nb, mb) = &
715 : upf%paw%aewfc_rel(1:upf%paw%iraug, nb) * &
716 0 : upf%paw%aewfc_rel(1:upf%paw%iraug, mb)
717 : !
718 : ! The small component is added to pfunc. pfunc_rel is useful only
719 : ! to add a small magnetic contribution
720 : !
721 : upf%paw%pfunc (1:upf%paw%iraug, nb, mb) = &
722 : upf%paw%pfunc (1:upf%paw%iraug, nb, mb) + &
723 0 : upf%paw%pfunc_rel (1:upf%paw%iraug, nb, mb)
724 : ENDIF
725 0 : upf%paw%pfunc(upf%paw%iraug+1:,nb,mb) = 0._dp
726 : !
727 0 : upf%paw%pfunc (1:upf%mesh, mb, nb) = upf%paw%pfunc (1:upf%mesh, nb, mb)
728 0 : IF (upf%has_so) upf%paw%pfunc_rel (1:upf%mesh, mb, nb) = &
729 0 : upf%paw%pfunc_rel (1:upf%mesh, nb, mb)
730 : ENDDO
731 : ENDDO
732 : !
733 : ! Pseudo wavefunctions (not only the ones for oc > 0)
734 : ! All-electron wavefunctions
735 0 : ALLOCATE(upf%paw%ptfunc(upf%mesh, upf%nbeta,upf%nbeta) )
736 0 : upf%paw%ptfunc(:,:,:) = 0._dp
737 0 : DO nb=1,upf%nbeta
738 0 : DO mb=1,upf%nbeta
739 : upf%paw%ptfunc (1:upf%mesh, nb, mb) = &
740 0 : upf%pswfc(1:upf%mesh, nb) * upf%pswfc(1:upf%mesh, mb)
741 0 : upf%paw%ptfunc(upf%paw%iraug+1:,nb,mb) = 0._dp
742 : !
743 0 : upf%paw%ptfunc (1:upf%mesh, mb, nb) = upf%paw%ptfunc (1:upf%mesh, nb, mb)
744 : ENDDO
745 : ENDDO
746 : !
747 : END SUBROUTINE read_pp_paw
748 : !--------------------------------------------------------
749 16 : SUBROUTINE read_pp_gipaw ( upf )
750 : !--------------------------------------------------------
751 : !
752 : !IMPLICIT NONE
753 : TYPE(pseudo_upf),INTENT(INOUT) :: upf ! the pseudo data
754 : !
755 : INTEGER :: nb, mb, ierr
756 : CHARACTER(LEN=24) :: tag
757 : !
758 16 : IF (.NOT. upf%has_gipaw) RETURN
759 : !
760 0 : CALL xmlr_opentag( capitalize_if_v2('pp_gipaw') )
761 0 : CALL get_attr ('gipaw_data_format', upf%gipaw_data_format )
762 0 : IF ( v2 ) THEN
763 0 : CALL xmlr_opentag( 'PP_GIPAW_CORE_ORBITALS', IERR=ierr )
764 : ! ierr= 0 for case <PP_GIPAW_CORE_ORBITALS>...</PP_GIPAW_CORE_ORBITALS>
765 : ! ierr=-1 for case <PP_GIPAW_CORE_ORBITALS ... />
766 0 : CALL get_attr ('number_of_core_orbitals', upf%gipaw_ncore_orbitals)
767 : ELSE
768 0 : CALL xmlr_readtag ('number_of_core_orbitals', upf%gipaw_ncore_orbitals)
769 0 : IF ( .NOT. upf%paw_as_gipaw) &
770 0 : CALL xmlr_readtag( 'number_of_valence_orbitals', upf%gipaw_wfs_nchannels)
771 : END IF
772 0 : ALLOCATE ( upf%gipaw_core_orbital(upf%mesh,upf%gipaw_ncore_orbitals) )
773 0 : ALLOCATE ( upf%gipaw_core_orbital_n(upf%gipaw_ncore_orbitals) )
774 0 : ALLOCATE ( upf%gipaw_core_orbital_el(upf%gipaw_ncore_orbitals) )
775 0 : ALLOCATE ( upf%gipaw_core_orbital_l(upf%gipaw_ncore_orbitals) )
776 0 : DO nb = 1,upf%gipaw_ncore_orbitals
777 0 : IF ( v2 ) THEN
778 0 : tag = "PP_GIPAW_CORE_ORBITAL."//i2c(nb)
779 : ELSE
780 0 : tag = 'pp_gipaw_core_orbital'
781 : END IF
782 0 : CALL xmlr_readtag( tag, upf%gipaw_core_orbital(1:upf%mesh,nb) )
783 0 : CALL get_attr ('index', mb)
784 0 : IF ( nb /= mb ) CALL upf_error('read_pp_gipaw','mismatch',1)
785 0 : CALL get_attr ('label', upf%gipaw_core_orbital_el(nb) )
786 0 : CALL get_attr ('n', upf%gipaw_core_orbital_n(nb) )
787 0 : CALL get_attr ('l', upf%gipaw_core_orbital_l(nb) )
788 : END DO
789 : ! close only for case <PP_GIPAW_CORE_ORBITALS> ... </PP_GIPAW_CORE_ORBITALS>
790 0 : IF ( v2 .AND. ierr == 0 ) CALL xmlr_closetag ( )
791 : !
792 0 : IF ( upf%paw_as_gipaw) THEN
793 : !
794 : ! PAW as GIPAW case: all-electron and pseudo-orbitals not read here
795 : !
796 0 : upf%gipaw_wfs_nchannels = upf%nbeta
797 0 : ALLOCATE ( upf%gipaw_wfs_el(upf%gipaw_wfs_nchannels) )
798 0 : ALLOCATE ( upf%gipaw_wfs_ll(upf%gipaw_wfs_nchannels) )
799 0 : ALLOCATE ( upf%gipaw_wfs_rcut(upf%gipaw_wfs_nchannels) )
800 0 : ALLOCATE ( upf%gipaw_wfs_rcutus(upf%gipaw_wfs_nchannels) )
801 0 : ALLOCATE ( upf%gipaw_wfs_ae(upf%mesh,upf%gipaw_wfs_nchannels) )
802 0 : ALLOCATE ( upf%gipaw_wfs_ps(upf%mesh,upf%gipaw_wfs_nchannels) )
803 0 : DO nb = 1,upf%gipaw_wfs_nchannels
804 0 : upf%gipaw_wfs_el(nb) = upf%els_beta(nb)
805 0 : upf%gipaw_wfs_ll(nb) = upf%lll(nb)
806 0 : upf%gipaw_wfs_ae(:,nb) = upf%aewfc(:,nb)
807 : ENDDO
808 0 : DO nb = 1,upf%gipaw_wfs_nchannels
809 0 : upf%gipaw_wfs_ps(:,nb) = upf%pswfc(:,nb)
810 : ENDDO
811 0 : ALLOCATE ( upf%gipaw_vlocal_ae(upf%mesh) )
812 0 : ALLOCATE ( upf%gipaw_vlocal_ps(upf%mesh) )
813 0 : upf%gipaw_vlocal_ae(:)= upf%paw%ae_vloc(:)
814 0 : upf%gipaw_vlocal_ps(:)= upf%vloc(:)
815 0 : DO nb = 1,upf%gipaw_wfs_nchannels
816 0 : upf%gipaw_wfs_rcut(nb)=upf%rcut(nb)
817 0 : upf%gipaw_wfs_rcutus(nb)=upf%rcutus(nb)
818 : ENDDO
819 : !
820 : ELSE
821 : !
822 : ! Read valence all-electron and pseudo orbitals
823 : !
824 0 : IF ( v2 ) THEN
825 0 : CALL xmlr_opentag( 'PP_GIPAW_ORBITALS' )
826 : CALL get_attr( 'number_of_valence_orbitals', &
827 0 : upf%gipaw_wfs_nchannels )
828 : END IF
829 0 : ALLOCATE ( upf%gipaw_wfs_el(upf%gipaw_wfs_nchannels) )
830 0 : ALLOCATE ( upf%gipaw_wfs_ll(upf%gipaw_wfs_nchannels) )
831 0 : ALLOCATE ( upf%gipaw_wfs_rcut(upf%gipaw_wfs_nchannels) )
832 0 : ALLOCATE ( upf%gipaw_wfs_rcutus(upf%gipaw_wfs_nchannels) )
833 0 : ALLOCATE ( upf%gipaw_wfs_ae(upf%mesh,upf%gipaw_wfs_nchannels) )
834 0 : ALLOCATE ( upf%gipaw_wfs_ps(upf%mesh,upf%gipaw_wfs_nchannels) )
835 0 : DO nb = 1,upf%gipaw_wfs_nchannels
836 0 : IF ( v2 ) THEN
837 0 : tag = "PP_GIPAW_ORBITAL."//i2c(nb)
838 : ELSE
839 0 : tag = 'pp_gipaw_orbital'
840 : END IF
841 0 : CALL xmlr_opentag( tag )
842 0 : CALL get_attr ('index', mb)
843 0 : IF ( nb /= mb ) CALL upf_error('read_pp_gipaw','mismatch',2)
844 0 : CALL get_attr ('label', upf%gipaw_wfs_el(nb) )
845 0 : CALL get_attr ('l', upf%gipaw_wfs_ll(nb) )
846 0 : CALL get_attr ('cutoff_radius', upf%gipaw_wfs_rcut(nb) )
847 0 : CALL get_attr ('ultrasoft_cutoff_radius', upf%gipaw_wfs_rcutus(nb) )
848 : CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_wfs_ae'), &
849 0 : upf%gipaw_wfs_ae(1:upf%mesh,nb) )
850 : CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_wfs_ps'),&
851 0 : upf%gipaw_wfs_ps(1:upf%mesh,nb) )
852 0 : CALL xmlr_closetag ()
853 : END DO
854 0 : IF ( v2 ) CALL xmlr_closetag( )
855 : !
856 : ! Read all-electron and pseudo local potentials
857 : !
858 0 : ALLOCATE ( upf%gipaw_vlocal_ae(upf%mesh) )
859 0 : ALLOCATE ( upf%gipaw_vlocal_ps(upf%mesh) )
860 0 : CALL xmlr_opentag( capitalize_if_v2('pp_gipaw_vlocal') )
861 : CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_vlocal_ae'), &
862 0 : upf%gipaw_vlocal_ae(1:upf%mesh) )
863 : CALL xmlr_readtag( capitalize_if_v2('pp_gipaw_vlocal_ps'), &
864 0 : upf%gipaw_vlocal_ps(1:upf%mesh) )
865 0 : CALL xmlr_closetag ()
866 : END IF
867 0 : CALL xmlr_closetag () ! end pp_gipaw
868 : !
869 : END SUBROUTINE read_pp_gipaw
870 : !
871 :
872 0 : SUBROUTINE upf_error( calling_routine, message, ierr )
873 : use m_abicore
874 : use m_errors
875 :
876 : CHARACTER(LEN=*), INTENT(IN) :: calling_routine, message
877 : ! the name of the calling calling_routine
878 : ! the output message
879 : INTEGER, INTENT(IN) :: ierr
880 :
881 0 : ABI_ERROR(message)
882 : if (.False.) write(std_out,*)trim(calling_routine), ierr
883 :
884 0 : END SUBROUTINE upf_error
885 :
886 : !
887 16 : subroutine remove_non_ascii(in_str)
888 :
889 : character(len=*),intent(inout) :: in_str
890 :
891 16 : character(len=len(in_str)) :: tmp_str
892 : integer :: ii, cnt
893 :
894 : ! in UTF8 character value c2 a0 (194 160) is defined as NO-BREAK SPACE.
895 : ! According to ISO/IEC 8859 this is a space that does not allow a line break to be inserted.
896 16 : cnt = 0
897 16 : tmp_str = ""
898 416 : do ii=1,len(in_str)
899 400 : if ((iachar(in_str(ii:ii)) >= 127)) cycle
900 400 : cnt = cnt + 1 !; write(std_out, *) "copying`", in_str(ii:ii), "` with iachar:", iachar(in_str(ii:ii))
901 416 : tmp_str(cnt:cnt) = in_str(ii:ii)
902 : end do
903 :
904 416 : do ii=1,len(in_str)
905 416 : in_str(ii:ii) = " "
906 : end do
907 16 : in_str = tmp_str
908 :
909 16 : end subroutine remove_non_ascii
910 :
911 1520 : END MODULE read_upf_new_module
|
