Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_types_lapw
8 : USE m_judft
9 : use m_types_fleurinput
10 : use m_types_nococonv
11 : IMPLICIT NONE
12 : PRIVATE
13 : !These dimensions should be set once per call of FLEUR
14 : !They can be queried by the functions lapw%dim_nvd,...
15 : !You probably should avoid using the variables directly
16 : integer, save :: lapw_dim_nvd
17 : integer, save :: lapw_dim_nv2d
18 : integer, save :: lapw_dim_nbasfcn
19 :
20 : TYPE t_lapw
21 : INTEGER :: nv(2)
22 : INTEGER :: num_local_cols(2)
23 : INTEGER :: nv_tot
24 : INTEGER :: nmat
25 : INTEGER :: nlotot
26 : INTEGER, ALLOCATABLE:: k1(:, :)
27 : INTEGER, ALLOCATABLE:: k2(:, :)
28 : INTEGER, ALLOCATABLE:: k3(:, :)
29 : INTEGER, ALLOCATABLE:: gvec(:, :, :) !replaces k1,k2,k3
30 : INTEGER, ALLOCATABLE:: kp(:, :)
31 : REAL, ALLOCATABLE::rk(:, :)
32 : REAL, ALLOCATABLE::gk(:, :, :)
33 : REAL, ALLOCATABLE::vk(:, :, :)
34 : INTEGER, ALLOCATABLE::index_lo(:, :)
35 : INTEGER, ALLOCATABLE::kvec(:, :, :)
36 : INTEGER, ALLOCATABLE::nkvec(:, :)
37 : REAL :: bkpt(3)
38 : REAL :: qphon(3)
39 : CONTAINS
40 : procedure :: lapw_init => t_lapw_init
41 : procedure :: lapw_init_fi => t_lapw_init_fi
42 : GENERIC :: init => lapw_init, lapw_init_fi
43 : PROCEDURE, PASS :: alloc => lapw_alloc
44 : PROCEDURE, PASS :: phase_factors => lapw_phase_factors
45 : procedure, pass :: hyb_num_bas_fun => hyb_num_bas_fun
46 : PROCEDURE, NOPASS:: dim_nvd
47 : PROCEDURE, NOPASS:: dim_nv2d
48 : PROCEDURE, NOPASS:: dim_nbasfcn
49 : PROCEDURE, NOPASS:: init_dim => lapw_init_dim
50 : END TYPE t_lapw
51 : PUBLIC :: t_lapw, lapw_dim_nbasfcn, lapw_dim_nvd, lapw_dim_nv2d
52 :
53 : CONTAINS
54 160 : function hyb_num_bas_fun(lapw, fi) result(nbasfcn)
55 : implicit NONE
56 : class(t_lapw), intent(in) :: lapw
57 : type(t_fleurinput), intent(in) :: fi
58 :
59 : integer :: nbasfcn
60 160 : if (fi%noco%l_noco) then
61 0 : nbasfcn = lapw%nv(1) + lapw%nv(2) + 2*fi%atoms%nlotot
62 : else
63 160 : nbasfcn = lapw%nv(1) + fi%atoms%nlotot
64 : endif
65 160 : end function hyb_num_bas_fun
66 :
67 160 : subroutine lapw_init_dim(nvd_in, nv2d_in, nbasfcn_in)
68 : IMPLICIT NONE
69 : INTEGER, INTENT(IN) :: nvd_in, nv2d_in, nbasfcn_in
70 160 : lapw_dim_nvd = nvd_in
71 160 : lapw_dim_nv2d = nv2d_in
72 160 : lapw_dim_nbasfcn = nbasfcn_in
73 160 : end subroutine
74 :
75 795 : PURE INTEGER function dim_nvd()
76 795 : dim_nvd = lapw_dim_nvd
77 795 : end function
78 39526 : PURE INTEGER function dim_nv2d()
79 39526 : dim_nv2d = lapw_dim_nv2d
80 39526 : end function
81 694 : PURE INTEGER function dim_nbasfcn()
82 694 : dim_nbasfcn = lapw_dim_nbasfcn
83 694 : end function
84 :
85 16632 : SUBROUTINE lapw_alloc(lapw, cell, input, noco, nococonv)
86 : !
87 : !*********************************************************************
88 : ! determines dimensions of the lapw basis set with |k+G|<rkmax.
89 : ! bkpt is the k-point given in internal units
90 : !*********************************************************************
91 : USE m_boxdim
92 : USE m_types_fleurinput
93 : USE m_types_nococonv
94 : IMPLICIT NONE
95 : TYPE(t_cell), INTENT(IN) :: cell
96 : TYPE(t_input), INTENT(IN) :: input
97 : TYPE(t_noco), INTENT(IN) :: noco
98 : TYPE(t_nococonv), INTENT(IN) :: nococonv
99 : CLASS(t_lapw), INTENT(INOUT) :: lapw
100 :
101 : INTEGER j1, j2, j3, mk1, mk2, mk3, nv, addX, addY, addZ
102 : INTEGER ispin, nvh(2)
103 :
104 : REAL arltv1, arltv2, arltv3, rkm, rk2, r2, s(3), sonlyg(3)
105 : ! ..
106 : !
107 : !-------> ABBREVIATIONS
108 : !
109 :
110 : ! rkmax : cut-off for |g+k|
111 : ! arltv(i) : length of reciprical lattice vector along
112 : ! direction (i)
113 : !
114 : !---> Determine rkmax box of size mk1, mk2, mk3,
115 : ! for which |G(mk1,mk2,mk3) + (k1,k2,k3)| < rkmax
116 : !
117 16632 : CALL boxdim(cell%bmat, arltv1, arltv2, arltv3)
118 :
119 : ! (add 1+1 due to integer rounding, strange k_vector in BZ)
120 16632 : mk1 = int(input%rkmax/arltv1) + 2
121 16632 : mk2 = int(input%rkmax/arltv2) + 2
122 16632 : mk3 = int(input%rkmax/arltv3) + 2
123 :
124 16632 : rkm = input%rkmax
125 16632 : rk2 = rkm*rkm
126 : !---> obtain vectors
127 : !---> in a spin-spiral calculation different basis sets are used for
128 : !---> the two spin directions, because the cutoff radius is defined
129 : !---> by |G + k +/- qss/2| < rkmax.
130 16632 : nvh(2) = 0
131 49896 : DO ispin = 1, MERGE(2, 1, noco%l_ss)
132 16734 : addX = abs(NINT((lapw%bkpt(1) + (2*ispin - 3)/2.0*nococonv%qss(1)+lapw%qphon(1))/arltv1))
133 16734 : addY = abs(NINT((lapw%bkpt(2) + (2*ispin - 3)/2.0*nococonv%qss(2)+lapw%qPhon(2))/arltv2))
134 16734 : addZ = abs(NINT((lapw%bkpt(3) + (2*ispin - 3)/2.0*nococonv%qss(3)+lapw%qphon(3))/arltv3))
135 16734 : nv = 0
136 179032 : DO j1 = -mk1 - addX, mk1 + addX
137 1793598 : DO j2 = -mk2 - addY, mk2 + addY
138 19588354 : DO j3 = -mk3 - addZ, mk3 + addZ
139 71245960 : s = lapw%bkpt + (/j1, j2, j3/) + (2*ispin - 3)/2.0*nococonv%qss + lapw%qphon
140 17811490 : sonlyg = (/j1, j2, j3/)
141 284983840 : r2 = dot_PRODUCT(MATMUL(s, cell%bbmat), s)
142 : !r2 = dot_PRODUCT(MATMUL(sonlyg, cell%bbmat), sonlyg)
143 19426056 : IF (r2 .LE. rk2) nv = nv + 1
144 : END DO
145 : END DO
146 : END DO
147 33366 : nvh(ispin) = nv
148 : END DO
149 16632 : nv = MAX(nvh(1), nvh(2))
150 :
151 16632 : IF (ALLOCATED(lapw%rk)) THEN
152 42813 : IF (SIZE(lapw%rk) == nv) THEN
153 129 : RETURN !
154 : ELSE
155 14142 : DEALLOCATE (lapw%rk, lapw%gvec, lapw%vk, lapw%gk)
156 14142 : DEALLOCATE (lapw%k1, lapw%k2, lapw%k3)
157 : ENDIF
158 : ENDIF
159 66012 : ALLOCATE (lapw%k1(nv, input%jspins)) !should be removed
160 49509 : ALLOCATE (lapw%k2(nv, input%jspins)) !
161 49509 : ALLOCATE (lapw%k3(nv, input%jspins)) !
162 66012 : ALLOCATE (lapw%rk(nv, input%jspins))
163 66012 : ALLOCATE (lapw%gvec(3, nv, input%jspins))
164 66012 : ALLOCATE (lapw%vk(3, nv, input%jspins))
165 49509 : ALLOCATE (lapw%gk(3, nv, input%jspins))
166 :
167 18051693 : lapw%rk = 0; lapw%gvec = 0; lapw%nv = 0
168 : END SUBROUTINE lapw_alloc
169 :
170 136 : subroutine t_lapw_init_fi(lapw, fi, nococonv, nk, mpi, dfpt_q)
171 : USE m_types_mpi
172 : use m_types_fleurinput
173 : implicit none
174 : CLASS(t_lapw), INTENT(INOUT) :: lapw
175 : type(t_fleurinput), intent(in) :: fi
176 : TYPE(t_nococonv), INTENT(IN) :: nococonv
177 : INTEGER, INTENT(IN) :: nk
178 : TYPE(t_mpi), INTENT(IN), OPTIONAL:: mpi
179 : REAL, INTENT(IN), OPTIONAL :: dfpt_q(3)
180 :
181 136 : IF (PRESENT(mpi)) THEN
182 0 : IF (PRESENT(dfpt_q)) THEN
183 0 : call lapw%lapw_init(fi%input, fi%noco, nococonv, fi%kpts, fi%atoms, fi%sym, nk, fi%cell, mpi, dfpt_q)
184 : ELSE
185 0 : call lapw%lapw_init(fi%input, fi%noco, nococonv, fi%kpts, fi%atoms, fi%sym, nk, fi%cell, mpi)
186 : END IF
187 : ELSE
188 136 : IF (PRESENT(dfpt_q)) THEN
189 0 : call lapw%lapw_init(fi%input, fi%noco, nococonv, fi%kpts, fi%atoms, fi%sym, nk, fi%cell, dfpt_q=dfpt_q)
190 : ELSE
191 136 : call lapw%lapw_init(fi%input, fi%noco, nococonv, fi%kpts, fi%atoms, fi%sym, nk, fi%cell)
192 : END IF
193 : END IF
194 136 : end subroutine t_lapw_init_fi
195 :
196 16632 : SUBROUTINE t_lapw_init(lapw, input, noco, nococonv, kpts, atoms, sym, &
197 : nk, cell, mpi, dfpt_q)
198 : USE m_types_mpi
199 : USE m_sort
200 : USE m_boxdim
201 : USE m_types_fleurinput
202 : USE m_types_kpts
203 : USE m_types_nococonv
204 : IMPLICIT NONE
205 :
206 :
207 : TYPE(t_input), INTENT(IN) :: input
208 : TYPE(t_noco), INTENT(IN) :: noco
209 : TYPE(t_nococonv), INTENT(IN) :: nococonv
210 : TYPE(t_cell), INTENT(IN) :: cell
211 : TYPE(t_atoms), INTENT(IN) :: atoms
212 : TYPE(t_sym), INTENT(IN) :: sym
213 : TYPE(t_kpts), INTENT(IN) :: kpts
214 : TYPE(t_mpi), INTENT(IN), OPTIONAL:: mpi
215 : CLASS(t_lapw), INTENT(INOUT) :: lapw
216 :
217 : REAL, INTENT(IN), OPTIONAL :: dfpt_q(3)
218 : ! ..
219 : ! .. Scalar Arguments ..
220 : INTEGER, INTENT(IN) :: nk
221 : !LOGICAL, INTENT(IN) :: l_zref
222 : ! ..
223 : ! .. Array Arguments ..
224 : ! ..
225 : ! .. Local Scalars ..
226 : REAL arltv1, arltv2, arltv3, r2, rk2, rkm, r2q, gla, eps, t, r2g, r2phon
227 : INTEGER i, j, j1, j2, j3, k, l, mk1, mk2, mk3, n, ispin, gmi, m, nred, n_inner, n_bound, itt(3), addX, addY, addZ
228 : ! ..
229 : ! .. Local Arrays ..
230 : REAL :: s(3), sq(3), sg(3), qphon(3), sphon(3)
231 16632 : REAL, ALLOCATABLE :: rk(:), rkq(:), rkqq(:), rg(:)
232 16632 : INTEGER, ALLOCATABLE :: gvec(:, :), index3(:)
233 :
234 16632 : call timestart("t_lapw_init")
235 : ! ..
236 : !---> in a spin-spiral calculation different basis sets are used for
237 : !---> the two spin directions, because the cutoff radius is defined
238 : !---> by |G + k +/- qss/2| < rkmax.
239 :
240 66528 : lapw%qphon = [0.0,0.0,0.0]
241 16632 : IF (PRESENT(dfpt_q)) lapw%qphon = dfpt_q
242 :
243 16632 : IF (nk > kpts%nkpt) THEN
244 1616 : lapw%bkpt(:) = kpts%bkf(:, nk)
245 : ELSE
246 64912 : lapw%bkpt(:) = kpts%bk(:, nk)
247 : ENDIF
248 :
249 16632 : CALL lapw%alloc(cell, input, noco, nococonv)
250 :
251 49896 : ALLOCATE (gvec(3, SIZE(lapw%gvec, 2)))
252 83160 : ALLOCATE (rk(SIZE(lapw%gvec, 2)), rkq(SIZE(lapw%gvec, 2)), rkqq(SIZE(lapw%gvec, 2)))
253 33264 : ALLOCATE (rg(SIZE(lapw%gvec, 2)))
254 49896 : ALLOCATE (index3(SIZE(lapw%gvec, 2)))
255 :
256 : !---> Determine rkmax box of size mk1, mk2, mk3,
257 : ! for which |G(mk1,mk2,mk3) + (k1,k2,k3)| < rkmax
258 : ! arltv(i) length of reciprical lattice vector along direction (i)
259 : !
260 16632 : CALL boxdim(cell%bmat, arltv1, arltv2, arltv3)
261 :
262 : ! (add 1+1 due to integer rounding, strange k_vector in BZ)
263 16632 : mk1 = int(input%rkmax/arltv1) + 4
264 16632 : mk2 = int(input%rkmax/arltv2) + 4
265 16632 : mk3 = int(input%rkmax/arltv3) + 4
266 :
267 16632 : rk2 = input%rkmax*input%rkmax
268 : !---> if too many basis functions, reduce rkmax
269 18834 : spinloop: DO ispin = 1, input%jspins
270 16734 : addX = abs(NINT((lapw%bkpt(1) + (2*ispin - 3)/2.0*nococonv%qss(1)+lapw%qphon(1))/arltv1))
271 16734 : addY = abs(NINT((lapw%bkpt(2) + (2*ispin - 3)/2.0*nococonv%qss(2)+lapw%qphon(2))/arltv2))
272 16734 : addZ = abs(NINT((lapw%bkpt(3) + (2*ispin - 3)/2.0*nococonv%qss(3)+lapw%qphon(3))/arltv3))
273 : !---> obtain vectors
274 16734 : n = 0
275 245968 : DO j1 = -mk1 - addX, mk1 + addX
276 3429862 : DO j2 = -mk2 - addY, mk2 + addY
277 50907506 : DO j3 = -mk3 - addZ, mk3 + addZ
278 189977512 : s = lapw%bkpt + (/j1, j2, j3/) + (2*ispin - 3)/2.0*nococonv%qss + lapw%qphon
279 189977512 : sq = lapw%bkpt + (/j1, j2, j3/)
280 189977512 : sg = (/j1, j2, j3/)
281 759910048 : r2 = dot_PRODUCT(s, MATMUL(s, cell%bbmat))
282 759910048 : r2q = dot_PRODUCT(sq, MATMUL(sq, cell%bbmat))
283 759910048 : r2g = dot_PRODUCT(sg, MATMUL(sg, cell%bbmat))
284 50678272 : IF (r2 .LE. rk2) THEN
285 : !IF (r2g .LE. rk2) THEN
286 2043850 : n = n + 1
287 8175400 : gvec(:, n) = (/j1, j2, j3/)
288 2043850 : rk(n) = SQRT(r2)
289 2043850 : rkq(n) = SQRT(r2q)
290 2043850 : rg(n) = SQRT(r2g)
291 : END IF
292 : ENDDO
293 : ENDDO
294 : ENDDO
295 16734 : lapw%nv(ispin) = n
296 :
297 : !Sort according to k+g, first construct secondary sort key
298 2060584 : DO k = 1, lapw%nv(ispin)
299 : rkqq(k) = (mk1 + gvec(1, k)) + (mk2 + gvec(2, k))*(2*mk1 + 1) + &
300 2060584 : (mk3 + gvec(3, k))*(2*mk1 + 1)*(2*mk2 + 1)
301 : ENDDO
302 16734 : CALL sort(index3(:lapw%nv(ispin)), rkq, rkqq)
303 : !CALL sort(index3(:lapw%nv(ispin)), rg, rkqq)
304 2060584 : DO n = 1, lapw%nv(ispin)
305 8175400 : lapw%gvec(:, n, ispin) = gvec(:, index3(n))
306 2060584 : lapw%rk(n, ispin) = rk(index3(n))
307 : ENDDO
308 : !---> determine pairs of K-vectors, where K_z = K'_-z to use
309 : !---> z-reflection
310 2060584 : DO k = 1, lapw%nv(ispin)
311 8175400 : lapw%vk(:, k, ispin) = lapw%bkpt + lapw%gvec(:, k, ispin) + (ispin - 1.5)*nococonv%qss + lapw%qphon
312 14323684 : lapw%gk(:, k, ispin) = MATMUL(TRANSPOSE(cell%bmat), lapw%vk(:, k, ispin))/MAX(lapw%rk(k, ispin), 1.0e-30)
313 : ENDDO
314 :
315 18834 : IF (.NOT. noco%l_ss .AND. input%jspins == 2) THEN
316 : !Second spin is the same
317 14532 : lapw%nv(2) = lapw%nv(1)
318 6278380 : lapw%gvec(:, :, 2) = lapw%gvec(:, :, 1)
319 1580494 : lapw%rk(:, 2) = lapw%rk(:, 1)
320 6278380 : lapw%vk(:, :, 2) = lapw%vk(:, :, 1)
321 6278380 : lapw%gk(:, :, 2) = lapw%gk(:, :, 1)
322 : EXIT spinloop
323 : END IF
324 :
325 : ENDDO spinloop
326 : !should be removed later...
327 7315420 : lapw%k1 = lapw%gvec(1, :, :)
328 7315420 : lapw%k2 = lapw%gvec(2, :, :)
329 7315420 : lapw%k3 = lapw%gvec(3, :, :)
330 :
331 : !Count No of lapw distributed to this PE
332 49896 : lapw%num_local_cols = 0
333 47898 : DO ispin = 1, input%jspins
334 47898 : IF (PRESENT(mpi)) THEN
335 29698 : DO k = mpi%n_rank + 1, lapw%nv(ispin), mpi%n_size
336 2053065 : lapw%num_local_cols(ispin) = lapw%num_local_cols(ispin) + 1
337 : END DO
338 : ELSE
339 1568 : lapw%num_local_cols(ispin) = lapw%nv(ispin)
340 : END IF
341 : END DO
342 :
343 19324 : IF (ANY(atoms%nlo > 0)) CALL priv_lo_basis_setup(lapw, atoms, input, sym, noco, nococonv, cell)
344 :
345 16632 : lapw%nv_tot = lapw%nv(1)
346 16632 : lapw%nmat = lapw%nv(1) + atoms%nlotot
347 16632 : IF (noco%l_noco) lapw%nv_tot = lapw%nv_tot + lapw%nv(2)
348 16632 : IF (noco%l_noco) lapw%nmat = lapw%nv_tot + 2*atoms%nlotot
349 :
350 :
351 16632 : call timestop("t_lapw_init")
352 : CONTAINS
353 :
354 15192 : SUBROUTINE priv_lo_basis_setup(lapw, atoms, input, sym, noco, nococonv, cell)
355 : USE m_types_fleurinput
356 :
357 : IMPLICIT NONE
358 : TYPE(t_lapw), INTENT(INOUT):: lapw
359 : TYPE(t_atoms), INTENT(IN) :: atoms
360 : TYPE(t_input), INTENT(IN) :: input
361 : TYPE(t_sym), INTENT(IN) :: sym
362 : TYPE(t_cell), INTENT(IN) :: cell
363 : TYPE(t_noco), INTENT(IN) :: noco
364 : TYPE(t_nococonv), INTENT(IN) :: nococonv
365 :
366 : INTEGER:: n, na, nn, np, lo, nkvec_sv, nkvec(atoms%nlod, 2), iindex
367 15192 : IF (.NOT. ALLOCATED(lapw%kvec)) THEN
368 8730 : ALLOCATE (lapw%kvec(2*(2*atoms%llod + 1), atoms%nlod, atoms%nat))
369 15164 : ALLOCATE (lapw%nkvec(atoms%nlod, atoms%nat));lapw%nkvec=0
370 5238 : ALLOCATE (lapw%index_lo(atoms%nlod, atoms%nat))
371 : ENDIF
372 15192 : iindex = 0
373 15192 : na = 0
374 15192 : nkvec_sv = 0
375 48254 : DO n = 1, atoms%ntype
376 81724 : DO nn = 1, atoms%neq(n)
377 33470 : na = na + 1
378 33470 : if (sym%invsat(na) > 1) cycle
379 33288 : np = sym%invtab(sym%ngopr(na))
380 33288 : CALL priv_vec_for_lo(atoms, input, sym, na, n, np, noco, nococonv, lapw, cell)
381 103496 : DO lo = 1, atoms%nlo(n)
382 37146 : lapw%index_lo(lo, na) = iindex
383 70616 : iindex = iindex + lapw%nkvec(lo, na)
384 : ENDDO
385 : ENDDO
386 : ENDDO
387 15192 : END SUBROUTINE priv_lo_basis_setup
388 :
389 : END SUBROUTINE t_lapw_init
390 :
391 60895 : SUBROUTINE lapw_phase_factors(lapw, iintsp, tau, qss, cph)
392 : USE m_constants
393 : USE m_types_fleurinput
394 : IMPLICIT NONE
395 : CLASS(t_lapw), INTENT(in):: lapw
396 : INTEGER, INTENT(IN) :: iintsp
397 : REAL, INTENT(in) :: tau(3), qss(3)
398 : COMPLEX, INTENT(out) :: cph(:)
399 :
400 : INTEGER:: k
401 : REAL:: th
402 :
403 : !$OMP PARALLEL DO DEFAULT(none) &
404 : !$OMP& SHARED(lapw,iintsp,tau,qss,cph)&
405 60895 : !$OMP& PRIVATE(k,th)
406 : DO k = 1, lapw%nv(iintsp)
407 : th = DOT_PRODUCT(lapw%gvec(:, k, iintsp) + (iintsp - 1.5)*qss + lapw%bkpt + lapw%qphon, tau)
408 : cph(k) = CMPLX(COS(tpi_const*th), SIN(tpi_const*th))
409 : END DO
410 : !$OMP END PARALLEL DO
411 60895 : END SUBROUTINE lapw_phase_factors
412 :
413 : SUBROUTINE priv_vec_for_lo_old(atoms, input, sym, na, n, np, noco, nococonv, lapw, cell)
414 :
415 : USE m_constants
416 : USE m_orthoglo
417 : USE m_ylm
418 : USE m_types_fleurinput
419 :
420 : IMPLICIT NONE
421 :
422 : TYPE(t_noco), INTENT(IN) :: noco
423 : TYPE(t_nococonv), INTENT(IN):: nococonv
424 : TYPE(t_sym), INTENT(IN) :: sym
425 : TYPE(t_cell), INTENT(IN) :: cell
426 : TYPE(t_atoms), INTENT(IN) :: atoms
427 : TYPE(t_input), INTENT(IN) :: input
428 : TYPE(t_lapw), INTENT(INOUT):: lapw
429 : ! ..
430 : ! .. Scalar Arguments ..
431 : INTEGER, INTENT(IN) :: na, n, np
432 : ! ..
433 : ! .. Array Arguments ..
434 : ! ..
435 : ! .. Local Scalars ..
436 : COMPLEX term1
437 : REAL th, con1
438 : INTEGER l, lo, mind, ll1, lm, iintsp, k, nkmin, ntyp, lmp, m, nintsp, k_start
439 : LOGICAL linind, enough, l_lo1
440 : ! ..
441 : ! .. Local Arrays ..
442 : INTEGER :: nkvec(atoms%nlod, 2)
443 : REAL qssbti(3), bmrot(3, 3), v(3), vmult(3)
444 : REAL :: gkrot(3, SIZE(lapw%gk, 2), 2)
445 : REAL :: rph(SIZE(lapw%gk, 2), 2)
446 : REAL :: cph(SIZE(lapw%gk, 2), 2)
447 : COMPLEX ylm((atoms%lmaxd + 1)**2)
448 : COMPLEX cwork(-2*atoms%llod:2*atoms%llod + 1, 2*(2*atoms%llod + 1), atoms%nlod, 2)
449 : ! ..
450 : ! .. Data statements ..
451 : REAL, PARAMETER :: eps = 1.0E-8
452 : REAL, PARAMETER :: linindq = 1.0e-6
453 :
454 : con1 = fpi_const/SQRT(cell%omtil)
455 : ntyp = n
456 : nintsp = MERGE(2, 1, noco%l_ss)
457 : DO iintsp = 1, nintsp
458 : IF (iintsp .EQ. 1) THEN
459 : qssbti = -nococonv%qss/2
460 : ELSE
461 : qssbti = +nococonv%qss/2
462 : ENDIF
463 :
464 : !---> set up phase factors
465 : DO k = 1, lapw%nv(iintsp)
466 : th = tpi_const*DOT_PRODUCT((/lapw%k1(k, iintsp), lapw%k2(k, iintsp), lapw%k3(k, iintsp)/) + qssbti + lapw%qphon, atoms%taual(:, na))
467 : rph(k, iintsp) = COS(th)
468 : cph(k, iintsp) = -SIN(th)
469 : END DO
470 :
471 : IF (np .EQ. 1) THEN
472 : gkrot(:, :, iintsp) = lapw%gk(:, :, iintsp)
473 : ELSE
474 : bmrot = MATMUL(1.*sym%mrot(:, :, np), cell%bmat)
475 : DO k = 1, lapw%nv(iintsp)
476 : !--> apply the rotation that brings this atom into the
477 : !--> representative (this is the definition of ngopr(na))
478 : !--> and transform to cartesian coordinates
479 : v(:) = lapw%vk(:, k, iintsp)
480 : gkrot(:, k, iintsp) = MATMUL(v, bmrot)
481 : END DO
482 : END IF
483 : !---> end loop over interstitial spin
484 : ENDDO
485 :
486 : nkvec(:, :) = 0
487 : cwork(:, :, :, :) = CMPLX(0.0, 0.0)
488 : enough = .FALSE.
489 :
490 : IF (noco%l_ss) THEN
491 : k_start = 2 ! avoid k=1 !!! GB16
492 : ELSE
493 : k_start = 1
494 : ENDIF
495 :
496 : DO k = k_start, MIN(lapw%nv(1), lapw%nv(nintsp))
497 : ! IF (ANY(lapw%rk(k,:nintsp).LT.eps)) CYCLE
498 : IF (.NOT. enough) THEN
499 : DO iintsp = 1, nintsp
500 :
501 : !--> generate spherical harmonics
502 : vmult(:) = gkrot(:, k, iintsp)
503 : CALL ylm4(atoms%lnonsph(ntyp), vmult, ylm)
504 : l_lo1 = .false.
505 : IF ((lapw%rk(k, iintsp) .LT. eps) .AND. (.not. noco%l_ss)) THEN
506 : l_lo1 = .true.
507 : ELSE
508 : l_lo1 = .false.
509 : ENDIF
510 : ! --> here comes a part of abccoflo()
511 : IF (l_lo1) THEN
512 : DO lo = 1, atoms%nlo(ntyp)
513 : IF ((nkvec(lo, iintsp) .EQ. 0) .AND. (atoms%llo(lo, ntyp) .EQ. 0)) THEN
514 : enough = .false.
515 : nkvec(lo, iintsp) = 1
516 : lapw%kvec(nkvec(lo, iintsp), lo, na) = k
517 : term1 = con1*((atoms%rmt(ntyp)**2)/2)
518 : cwork(0, 1, lo, iintsp) = term1/sqrt(2*tpi_const)
519 : IF ((sym%invsat(na) .EQ. 1) .OR. (sym%invsat(na) .EQ. 2)) THEN
520 : cwork(1, 1, lo, iintsp) = conjg(term1)/sqrt(2*tpi_const)
521 : ENDIF
522 : ENDIF
523 : ENDDO
524 : ELSE
525 : enough = .TRUE.
526 : term1 = con1*((atoms%rmt(ntyp)**2)/2)*CMPLX(rph(k, iintsp), cph(k, iintsp))
527 : DO lo = 1, atoms%nlo(ntyp)
528 : IF (sym%invsat(na) .EQ. 0) THEN
529 : IF ((nkvec(lo, iintsp)) .LT. (2*atoms%llo(lo, ntyp) + 1)) THEN
530 : enough = .FALSE.
531 : nkvec(lo, iintsp) = nkvec(lo, iintsp) + 1
532 : l = atoms%llo(lo, ntyp)
533 : ll1 = l*(l + 1) + 1
534 : DO m = -l, l
535 : lm = ll1 + m
536 : cwork(m, nkvec(lo, iintsp), lo, iintsp) = term1*ylm(lm)
537 : END DO
538 : CALL orthoglo(input%l_real, atoms, nkvec(lo, iintsp), lo, l, linindq, .FALSE., cwork(-2*atoms%llod, 1, 1, iintsp), linind)
539 : IF (linind) THEN
540 : lapw%kvec(nkvec(lo, iintsp), lo, na) = k
541 : ELSE
542 : nkvec(lo, iintsp) = nkvec(lo, iintsp) - 1
543 : ENDIF
544 : ENDIF
545 : ELSE
546 : IF ((sym%invsat(na) .EQ. 1) .OR. (sym%invsat(na) .EQ. 2)) THEN
547 : IF (nkvec(lo, iintsp) .LT. 2*(2*atoms%llo(lo, ntyp) + 1)) THEN
548 : enough = .FALSE.
549 : nkvec(lo, iintsp) = nkvec(lo, iintsp) + 1
550 : l = atoms%llo(lo, ntyp)
551 : ll1 = l*(l + 1) + 1
552 : DO m = -l, l
553 : lm = ll1 + m
554 : mind = -l + m
555 : cwork(mind, nkvec(lo, iintsp), lo, iintsp) = term1*ylm(lm)
556 : mind = l + 1 + m
557 : lmp = ll1 - m
558 : cwork(mind, nkvec(lo, iintsp), lo, iintsp) = ((-1)**(l + m))*CONJG(term1*ylm(lmp))
559 : END DO
560 : CALL orthoglo(input%l_real, atoms, nkvec(lo, iintsp), lo, l, linindq, .TRUE., cwork(-2*atoms%llod, 1, 1, iintsp), linind)
561 : IF (linind) THEN
562 : lapw%kvec(nkvec(lo, iintsp), lo, na) = k
563 : ! write(*,*) nkvec(lo,iintsp),k,' <- '
564 : ELSE
565 : nkvec(lo, iintsp) = nkvec(lo, iintsp) - 1
566 : END IF
567 : END IF
568 : END IF
569 : END IF
570 : END DO
571 : IF ((k .EQ. lapw%nv(iintsp)) .AND. (.NOT. enough)) THEN
572 : WRITE (oUnit, FMT=*) 'vec_for_lo did not find enough linearly independent'
573 : WRITE (oUnit, FMT=*) 'clo coefficient-vectors. the linear independence'
574 : WRITE (oUnit, FMT=*) 'quality, linindq, is set: ', linindq
575 : WRITE (oUnit, FMT=*) 'this value might be to large.'
576 : WRITE (*, *) na, k, lapw%nv
577 : CALL juDFT_error("not enough lin. indep. clo-vectors", calledby="vec_for_lo")
578 : END IF
579 : ! -- > end of abccoflo-part
580 : ENDIF
581 : ENDDO
582 : ENDIF
583 :
584 : ! --> check whether we have already enough k-vecs
585 : enough = .TRUE.
586 : DO lo = 1, atoms%nlo(ntyp)
587 : IF (nkvec(lo, 1) .EQ. nkvec(lo, nintsp)) THEN ! k-vec accepted by both spin channels
588 : IF (sym%invsat(na) .EQ. 0) THEN
589 : IF (nkvec(lo, 1) .LT. (2*atoms%llo(lo, ntyp) + 1)) THEN
590 : enough = .FALSE.
591 : ENDIF
592 : ELSE
593 : IF (nkvec(lo, 1) .LT. (2*(2*atoms%llo(lo, ntyp) + 1))) THEN
594 : enough = .FALSE.
595 : ENDIF
596 : ENDIF
597 : ELSE
598 : nkmin = MIN(nkvec(lo, 1), nkvec(lo, nintsp)) ! try another k-vec
599 : nkvec(lo, 1) = nkmin; nkvec(lo, nintsp) = nkmin
600 : enough = .FALSE.
601 : ENDIF
602 : ENDDO
603 : IF (enough) THEN
604 : lapw%nkvec(:atoms%nlo(ntyp), na) = nkvec(:atoms%nlo(ntyp), 1)
605 : RETURN
606 : ENDIF
607 : ENDDO
608 :
609 : END SUBROUTINE priv_vec_for_lo_old
610 :
611 33288 : SUBROUTINE priv_vec_for_lo(atoms, input, sym, na, ntype, np, noco, nococonv, lapw, cell)
612 :
613 : USE m_constants
614 : USE m_orthoglo
615 : USE m_ylm
616 : USE m_types_fleurinput
617 :
618 : IMPLICIT NONE
619 :
620 : TYPE(t_noco), INTENT(IN) :: noco
621 : TYPE(t_nococonv), INTENT(IN):: nococonv
622 : TYPE(t_sym), INTENT(IN) :: sym
623 : TYPE(t_cell), INTENT(IN) :: cell
624 : TYPE(t_atoms), INTENT(IN) :: atoms
625 : TYPE(t_input), INTENT(IN) :: input
626 : TYPE(t_lapw), INTENT(INOUT):: lapw
627 : ! ..
628 : ! .. Scalar Arguments ..
629 : INTEGER, INTENT(IN) :: na, ntype, np
630 : ! ..
631 : ! .. Array Arguments ..
632 : ! ..
633 : ! .. Local Scalars ..
634 : COMPLEX term1, norm
635 : REAL th, con1, linindq, linindqStart, linindqEnd, numK, stepSize
636 : INTEGER l, lo, mind, ll1, lm, iintsp, k, nkmin, lmp, m, nintsp, k_start, k_end, minIndex, maxIndex, increment
637 : INTEGER nApproach, nApproachEnd
638 : LOGICAL linind, enough, l_lo1, l_norm
639 : ! ..
640 : ! .. Local Arrays ..
641 33288 : INTEGER :: nkvec(atoms%nlod, 2)
642 : REAL qssbti(3), bmrot(3, 3), v(3), vmult(3)
643 33288 : REAL :: gkrot(3, SIZE(lapw%gk, 2), 2)
644 33288 : REAL :: rph(SIZE(lapw%gk, 2), 2)
645 33288 : REAL :: cph(SIZE(lapw%gk, 2), 2)
646 33288 : COMPLEX ylm((atoms%lmaxd + 1)**2)
647 33288 : COMPLEX cwork(-2*atoms%llod:2*atoms%llod + 1, 2*(2*atoms%llod + 1), atoms%nlod, 2)
648 : ! ..
649 : ! .. Data statements ..
650 : REAL, PARAMETER :: eps = 1.0E-8
651 :
652 33288 : con1 = fpi_const/SQRT(cell%omtil)
653 33288 : nintsp = MERGE(2, 1, noco%l_ss)
654 66576 : DO iintsp = 1, nintsp
655 33288 : IF (iintsp .EQ. 1) THEN
656 133152 : qssbti = -nococonv%qss/2
657 : ELSE
658 0 : qssbti = +nococonv%qss/2
659 : ENDIF
660 :
661 : !---> set up phase factors
662 4130628 : DO k = 1, lapw%nv(iintsp)
663 16389360 : th = tpi_const*DOT_PRODUCT((/lapw%k1(k, iintsp), lapw%k2(k, iintsp), lapw%k3(k, iintsp)/) + qssbti + lapw%qphon, atoms%taual(:, na))
664 4097340 : rph(k, iintsp) = COS(th)
665 4130628 : cph(k, iintsp) = -SIN(th)
666 : END DO
667 :
668 66576 : IF (np .EQ. 1) THEN
669 16073352 : gkrot(:, :, iintsp) = lapw%gk(:, :, iintsp)
670 : ELSE
671 11648 : bmrot = MATMUL(1.*sym%mrot(:, :, np), cell%bmat)
672 87492 : DO k = 1, lapw%nv(iintsp)
673 : !--> apply the rotation that brings this atom into the
674 : !--> representative (this is the definition of ngopr(na))
675 : !--> and transform to cartesian coordinates
676 349072 : v(:) = lapw%vk(:, k, iintsp)
677 1396512 : gkrot(:, k, iintsp) = MATMUL(v, bmrot)
678 : END DO
679 : END IF
680 : !---> end loop over interstitial spin
681 : ENDDO
682 :
683 194328 : nkvec(:, :) = 0
684 4402536 : cwork(:, :, :, :) = CMPLX(0.0, 0.0)
685 33288 : enough = .FALSE.
686 :
687 : ! Typically the search for linearly independent K vectors starts from the
688 : ! top and then goes down. This seems to be more stable than the opposite
689 : ! direction. The exception is a calculation with a spin spiral. There it
690 : ! is important that the K vectors for both spins feature the same G
691 : ! vectors. When starting from the bottom this is typically simple to check
692 : ! by just comparing the indices of the two K vectors. From the top this
693 : ! would require a more elaborate comparison.
694 :
695 33288 : nApproachEnd = 8
696 :
697 33288 : IF (noco%l_ss) THEN
698 0 : nApproachEnd = 4
699 : END IF
700 :
701 70434 : DO lo = 1, atoms%nlo(ntype)
702 : enough = .FALSE.
703 : nApproach = 0
704 107580 : DO WHILE (.NOT.enough)
705 37146 : nApproach = nApproach + 1
706 111438 : nkvec(lo, :) = 0
707 3468342 : cwork(:,:,lo,:) = CMPLX(0.0, 0.0)
708 37146 : k_start = 2
709 37146 : k_end = 2
710 37146 : increment = 1
711 37146 : SELECT CASE (nApproach)
712 : CASE (1)
713 37146 : k_start = 2 ! avoid k=1 !!! GB16
714 37146 : k_end = MIN(lapw%nv(1), lapw%nv(nintsp))
715 37146 : increment = 1
716 37146 : l_norm = .FALSE.
717 37146 : linindqStart = 1.0e-6
718 37146 : linindqEnd = 1.0e-6
719 : CASE (2)
720 0 : k_start = 2 ! avoid k=1 !!! GB16
721 0 : k_end = MIN(lapw%nv(1), lapw%nv(nintsp))
722 0 : increment = 1
723 0 : l_norm = .TRUE.
724 0 : linindqStart = 1.0e-6
725 0 : linindqEnd = 1.0e-6
726 : CASE (3)
727 0 : k_start = 2 ! avoid k=1 !!! GB16
728 0 : k_end = MIN(lapw%nv(1), lapw%nv(nintsp))
729 0 : increment = 1
730 0 : l_norm = .FALSE.
731 0 : linindqStart = 2.0e-5
732 0 : linindqEnd = 1.0e-6
733 : CASE (4)
734 0 : k_start = 2 ! avoid k=1 !!! GB16
735 0 : k_end = MIN(lapw%nv(1), lapw%nv(nintsp))
736 0 : increment = 1
737 0 : l_norm = .TRUE.
738 0 : linindqStart = 2.0e-5
739 0 : linindqEnd = 1.0e-6
740 : CASE (5)
741 0 : k_start = MIN(lapw%nv(1), lapw%nv(nintsp))
742 0 : k_end = 1
743 0 : increment = -1
744 0 : l_norm = .FALSE.
745 0 : linindqStart = 1.0e-6
746 0 : linindqEnd = 1.0e-6
747 : CASE (6)
748 0 : k_start = MIN(lapw%nv(1), lapw%nv(nintsp))
749 0 : k_end = 1
750 0 : increment = -1
751 0 : l_norm = .TRUE.
752 0 : linindqStart = 1.0e-6
753 0 : linindqEnd = 1.0e-6
754 : CASE (7)
755 0 : k_start = MIN(lapw%nv(1), lapw%nv(nintsp))
756 0 : k_end = 1
757 0 : increment = -1
758 0 : l_norm = .FALSE.
759 0 : linindqStart = 2.0e-5
760 0 : linindqEnd = 1.0e-6
761 : CASE (8)
762 0 : k_start = MIN(lapw%nv(1), lapw%nv(nintsp))
763 0 : k_end = 1
764 0 : increment = -1
765 0 : l_norm = .TRUE.
766 0 : linindqStart = 2.0e-5
767 37146 : linindqEnd = 1.0e-6
768 : END SELECT
769 37146 : k = k_start - increment
770 37146 : DO WHILE (.NOT.enough)
771 98442 : enough = .FALSE.
772 98442 : k = k + increment
773 98442 : IF ((k.GT.MAX(k_start,k_end)).OR.(k.LT.MIN(k_start,k_end))) THEN
774 0 : IF (nApproach.GT.nApproachEnd) THEN
775 0 : WRITE (oUnit, FMT=*) 'vec_for_lo did not find enough linearly independent'
776 0 : WRITE (oUnit, FMT=*) 'clo coefficient-vectors. the linear independence'
777 0 : WRITE (oUnit, FMT=*) 'quality, linindq, is set: ', linindqEnd
778 0 : WRITE (oUnit, FMT=*) 'this value might be to large.'
779 0 : WRITE (*, *) 'Atom: ', na, 'type: ', ntype, 'nv: ', lapw%nv, 'lo: ', lo, 'l: ', atoms%llo(lo, ntype)
780 0 : CALL juDFT_error("not enough lin. indep. clo-vectors", calledby="vec_for_lo")
781 : END IF
782 : EXIT
783 : END IF
784 :
785 196884 : DO iintsp = 1, nintsp
786 :
787 393768 : vmult(:) = gkrot(:, k, iintsp)
788 98442 : CALL ylm4(atoms%lnonsph(ntype), vmult, ylm)
789 98442 : l_lo1 = .false.
790 98442 : IF ((lapw%rk(k, iintsp) .LT. eps) .AND. (.not. noco%l_ss)) THEN
791 : l_lo1 = .true.
792 : ELSE
793 98442 : l_lo1 = .false.
794 : END IF
795 :
796 98442 : numK = MERGE(2*atoms%llo(lo, ntype)+1,2*(2*atoms%llo(lo, ntype)+1),sym%invsat(na).EQ.0)
797 :
798 196884 : IF (l_lo1) THEN
799 0 : IF ((nkvec(lo, iintsp) .EQ. 0) .AND. (atoms%llo(lo, ntype) .EQ. 0)) THEN
800 0 : nkvec(lo, iintsp) = 1
801 0 : lapw%kvec(nkvec(lo, iintsp), lo, na) = k
802 0 : term1 = con1 * (atoms%rmt(ntype)**2.0) / 2.0
803 0 : cwork(0, 1, lo, iintsp) = term1 / sqrt(2.0*tpi_const)
804 0 : norm = DOT_PRODUCT(cwork(0:0, 1, lo, iintsp),cwork(0:0, 1, lo, iintsp))
805 0 : IF (l_norm) cwork(0:0, 1, lo, iintsp) = cwork(0:0, 1, lo, iintsp) / SQRT(norm)
806 0 : IF ((sym%invsat(na) .EQ. 1) .OR. (sym%invsat(na) .EQ. 2)) THEN
807 0 : cwork(1, 1, lo, iintsp) = conjg(term1) / sqrt(2.0*tpi_const)
808 : END IF
809 :
810 : END IF
811 : ELSE
812 98442 : term1 = con1 * ((atoms%rmt(ntype)**2.0)/2.0) * CMPLX(rph(k, iintsp), cph(k, iintsp))
813 98442 : IF (sym%invsat(na) .EQ. 0) THEN
814 97240 : IF ((nkvec(lo, iintsp)) .LT. (2*atoms%llo(lo, ntype) + 1)) THEN
815 97240 : nkvec(lo, iintsp) = nkvec(lo, iintsp) + 1
816 97240 : l = atoms%llo(lo, ntype)
817 97240 : ll1 = l*(l + 1) + 1
818 389724 : DO m = -l, l
819 292484 : lm = ll1 + m
820 389724 : cwork(m, nkvec(lo, iintsp), lo, iintsp) = term1 * ylm(lm)
821 : END DO
822 389724 : norm = DOT_PRODUCT(cwork(-l:l, 1, lo, iintsp),cwork(-l:l, 1, lo, iintsp))
823 97240 : IF (l_norm) cwork(-l:l, 1, lo, iintsp) = cwork(-l:l, 1, lo, iintsp) / SQRT(norm)
824 97240 : stepSize = (REAL(linindqStart - linindqEnd)) / numK
825 97240 : linindq = (numK - (REAL(nkvec(lo, iintsp) + 1))) * stepSize + linindqEnd
826 97240 : CALL orthoglo(input%l_real, atoms, nkvec(lo, iintsp), lo, l, linindq, .FALSE., cwork(-2*atoms%llod, 1, 1, iintsp), linind)
827 97240 : IF (linind) THEN
828 87752 : lapw%kvec(nkvec(lo, iintsp), lo, na) = k
829 : ELSE
830 9488 : nkvec(lo, iintsp) = nkvec(lo, iintsp) - 1
831 : END IF
832 : END IF
833 1202 : ELSE IF ((sym%invsat(na) .EQ. 1) .OR. (sym%invsat(na) .EQ. 2)) THEN
834 1202 : IF (nkvec(lo, iintsp) .LT. 2*(2*atoms%llo(lo, ntype) + 1)) THEN
835 1202 : nkvec(lo, iintsp) = nkvec(lo, iintsp) + 1
836 1202 : l = atoms%llo(lo, ntype)
837 1202 : ll1 = l*(l + 1) + 1
838 4264 : DO m = -l, l
839 3062 : lm = ll1 + m
840 3062 : mind = -l + m
841 3062 : cwork(mind, nkvec(lo, iintsp), lo, iintsp) = term1 * ylm(lm)
842 3062 : mind = l + 1 + m
843 3062 : lmp = ll1 - m
844 4264 : cwork(mind, nkvec(lo, iintsp), lo, iintsp) = ((-1)**(l + m))*CONJG(term1*ylm(lmp))
845 : END DO
846 1202 : minIndex = -l - l
847 1202 : maxIndex = l + l + 1
848 7326 : norm = DOT_PRODUCT(cwork(-2*l:2*l+1, 1, lo, iintsp),cwork(-2*l:2*l+1, 1, lo, iintsp))
849 1202 : IF (l_norm) cwork(-2*l:2*l+1, 1, lo, iintsp) = cwork(-2*l:2*l+1, 1, lo, iintsp) / SQRT(norm)
850 1202 : stepSize = (REAL(linindqStart - linindqEnd)) / numK
851 1202 : linindq = (numK - (REAL(nkvec(lo, iintsp) + 1))) * stepSize + linindqEnd
852 1202 : CALL orthoglo(input%l_real, atoms, nkvec(lo, iintsp), lo, l, linindq, .TRUE., cwork(-2*atoms%llod, 1, 1, iintsp), linind)
853 1202 : IF (linind) THEN
854 1172 : lapw%kvec(nkvec(lo, iintsp), lo, na) = k
855 : ELSE
856 30 : nkvec(lo, iintsp) = nkvec(lo, iintsp) - 1
857 : END IF
858 : END IF
859 :
860 : END IF
861 : END IF
862 : END DO
863 :
864 98442 : enough = .TRUE.
865 98442 : IF (nkvec(lo, 1) .EQ. nkvec(lo, nintsp)) THEN ! k-vec accepted by both spin channels
866 98442 : IF (sym%invsat(na) .EQ. 0) THEN
867 97240 : IF (nkvec(lo, 1) .LT. (2*atoms%llo(lo, ntype) + 1)) THEN
868 : enough = .FALSE.
869 : ENDIF
870 : ELSE
871 1202 : IF (nkvec(lo, 1) .LT. (2*(2*atoms%llo(lo, ntype) + 1))) THEN
872 : enough = .FALSE.
873 : ENDIF
874 : ENDIF
875 : ELSE
876 0 : nkmin = MIN(nkvec(lo, 1), nkvec(lo, nintsp)) ! try another k-vec
877 0 : nkvec(lo, 1) = nkmin
878 0 : nkvec(lo, nintsp) = nkmin
879 0 : enough = .FALSE.
880 : END IF
881 : END DO
882 : END DO
883 : END DO
884 :
885 70434 : lapw%nkvec(:atoms%nlo(ntype), na) = nkvec(:atoms%nlo(ntype), 1)
886 :
887 33288 : END SUBROUTINE priv_vec_for_lo
888 :
889 0 : END MODULE m_types_lapw
|