Branch data Line data Source code
1 : : !!****m* ABINIT/m_wkk
2 : : !! NAME
3 : : !! m_wkk
4 : : !!
5 : : !! FUNCTION
6 : : !! Computation of the matrix elements of the screened interaction W_kk'.
7 : : !!
8 : : !! COPYRIGHT
9 : : !! Copyright (C) 2008-2026 ABINIT group (MG)
10 : : !! This file is distributed under the terms of the
11 : : !! GNU General Public License, see ~abinit/COPYING
12 : : !! or http://www.gnu.org/copyleft/gpl.txt.
13 : : !!
14 : : !! SOURCE
15 : :
16 : : #if defined HAVE_CONFIG_H
17 : : #include "config.h"
18 : : #endif
19 : :
20 : : #include "abi_common.h"
21 : :
22 : : module m_wkk
23 : :
24 : : !use, intrinsic :: iso_c_binding
25 : : use defs_basis
26 : : use m_abicore
27 : : use m_xmpi
28 : : use m_mpinfo
29 : : use m_errors
30 : : use m_fft
31 : : use m_wfd
32 : : use m_hdr
33 : : use netcdf
34 : : use m_nctk
35 : : use m_dtset
36 : : use m_dtfil
37 : : use m_fstab
38 : : use m_ephtk
39 : :
40 : : use defs_abitypes, only : mpi_type
41 : : use defs_datatypes, only : pseudopotential_type
42 : : use m_gwdefs, only : GW_Q0_DEFAULT, cone_gw, czero_gw
43 : : use m_hide_blas, only : xgemm, xdotc
44 : : use m_special_funcs, only : gaussian
45 : : use m_time, only : cwtime, cwtime_report, timab, sec2str
46 : : use m_fstrings, only : itoa, ftoa, sjoin, ktoa, ltoa, strcat
47 : : use m_numeric_tools, only : arth, linspace ! print_arr
48 : : use m_io_tools, only : iomode_from_fname
49 : : use m_fftcore, only : ngfft_seq, sphereboundary, get_kg, print_ngfft
50 : : use m_crystal, only : crystal_t
51 : : use m_bz_mesh, only : kmesh_t, findqg0
52 : : use m_htetra, only : htetra_t
53 : : use m_gsphere, only : gsphere_t
54 : : use m_pawtab, only : pawtab_type
55 : : use m_io_screening, only : hscr_t, get_hscr_qmesh_gsph
56 : : use m_vcoul, only : vcoul_t
57 : : !use m_occ, only : get_fact_spin_tol_empty
58 : : use m_ebands, only : ebands_t, edos_t
59 : : use m_pstat, only : pstat_proc
60 : : use m_sigtk, only : sigtk_multiply_by_vc_sqrt
61 : : use m_screening, only : epsm1_t
62 : : use m_ddk, only : ddkop_t
63 : :
64 : : implicit none
65 : :
66 : : private
67 : : !!***
68 : :
69 : : public :: wkk_run ! Main entry point to compute Wkk' matrix elements.
70 : :
71 : : !----------------------------------------------------------------------
72 : :
73 : : contains !=====================================================
74 : : !!***
75 : :
76 : : !----------------------------------------------------------------------
77 : :
78 : : !!****f* m_wkk/wkk_run
79 : : !! NAME
80 : : !! wkk_run
81 : : !!
82 : : !! FUNCTION
83 : : !! Computation of the matrix elements of the screened interaction W_kk'
84 : : !! between two Cooper pairs.
85 : : !!
86 : : !! INPUTS
87 : : !! wfk0_path=String with the path to the GS unperturbed WFK file.
88 : : !! dtfil<datafiles_type>=Variables related to files.
89 : : !! ngfft(18),ngfftf(18)=Coarse and Fine FFT meshes.
90 : : !! dtset<dataset_type>=All input variables for this dataset.
91 : : !! cryst: Crystalline structure
92 : : !! ebands<ebands_t>=The GS KS band structure (energies, occupancies, k-weights...)
93 : : !! wfk_hdr=Header of the WFK file.
94 : : !! pawtab(ntypat*usepaw)<pawtab_type>=Paw tabulated starting data.
95 : : !! psps<pseudopotential_type>=Variables related to pseudopotentials.
96 : : !! comm=MPI communicator.
97 : : !!
98 : : !! OUTPUT
99 : : !! Results are written to file in netcdf format.
100 : : !!
101 : : !! SOURCE
102 : :
103 : 0 : subroutine wkk_run(wfk0_path, dtfil, ngfft, ngfftf, dtset, cryst, ebands, wfk_hdr, &
104 : 0 : pawtab, psps, mpi_enreg, comm)
105 : :
106 : : !Arguments ------------------------------------
107 : : !scalars
108 : : character(len=*),intent(in) :: wfk0_path
109 : : type(datafiles_type),intent(in) :: dtfil
110 : : type(dataset_type),intent(inout) :: dtset
111 : : type(crystal_t),intent(in) :: cryst
112 : : type(ebands_t),target,intent(in) :: ebands
113 : : type(hdr_type),intent(in) :: wfk_hdr
114 : : type(pseudopotential_type),intent(in) :: psps
115 : : type(mpi_type),intent(inout) :: mpi_enreg
116 : : integer,intent(in) :: comm
117 : : !arrays
118 : : integer,intent(in) :: ngfft(18), ngfftf(18)
119 : : type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw)
120 : :
121 : : !Local variables ------------------------------
122 : : !scalars
123 : : integer,parameter :: LOG_MODQ = 1, LOG_MODK = 4, istwfk1 = 1, master = 0, ndat1 = 1
124 : : integer :: id_required, approx_type, ikxc, option_test, nkxc, ig
125 : : integer :: nkibz, my_rank, nsppol, iq_ibz, iq_bz, isym_qq, itim_qq
126 : : integer :: nspinor,nprocs, ii, spin, npw_x, npw_c ! max_npw_xc, min_npw_xc, my_is, cplex, ib,
127 : : integer :: bstart_k, bstop_k, nband_k, ncols, mpw, ierr !, ncerr
128 : : integer :: bstart_kp, bstop_kp, nband_kp, ik_bz, bmin, bmax, max_nb
129 : : integer :: ikp_ibz, isym_kp, trev_kp, npw_kp, istwf_kp, npw_kp_ibz, istwf_kp_ibz
130 : : integer :: ik_ibz, isym_k, trev_k, npw_k, istwf_k, npw_k_ibz, istwf_k_ibz
131 : : integer :: m_k, im_k, n_kp, in_kp, root_ncid, wkk_mode, ne
132 : : integer :: nfft,nfftf,mgfft,mgfftf,nkpg_kp,nkpg_k,edos_intmeth, nqlwl, scr_iomode ! cnt,
133 : : integer :: ikp_bz, my_ikp, my_nkp !, my_iq,
134 : : real(dp) :: cpu_all, wall_all, gflops_all, cpu_kp, wall_kp, gflops_kp ! cpu, wall, gflops,
135 : : real(dp) :: edos_step, edos_broad, e_mk, e_nkp, e_min, e_max, e_step, smear_mk, smear_nkp, faq
136 : : logical :: isirr_k, isirr_kp, qq_is_gamma, remove_exchange, print_time_kp
137 : 0 : type(wfd_t) :: wfd
138 [ # # # # : 0 : type(kmesh_t) :: qmesh, kmesh
# # # # ]
139 : 0 : type(gsphere_t),target :: gsph_x, gsph_c
140 [ # # ]: 0 : type(hscr_t),target :: hscr
141 [ # # # # : 0 : type(vcoul_t) :: vcp
# # # # #
# ]
142 : 0 : type(edos_t) :: edos
143 [ # # ]: 0 : type(epsm1_t) :: epsm1
144 [ # # # # : 0 : type(ddkop_t) :: ddkop
# # ]
145 : 0 : type(htetra_t) :: tetra
146 : : character(len=fnlen) :: screen_filepath
147 : : character(len=5000) :: msg !, qkp_string
148 : : !arrays
149 : : integer :: g0_k(3), g0_kp(3), g0_qq(3), units(2), work_ngfft(18), gmax(3), mapl_k(6), mapl_kp(6), mg0(3)
150 : 0 : integer,allocatable :: kg_kp(:,:), kg_k(:,:), gbound_kp(:,:), gbound_k(:,:), gbound_c(:,:), gbound_x(:,:)
151 : 0 : integer,allocatable :: nband(:,:), wfd_istwfk(:)
152 : 0 : integer, contiguous, pointer :: kg_c(:,:), kg_x(:,:)
153 : : real(dp) :: kk_ibz(3), kk_bz(3), kp_ibz(3), kp_bz(3), qq_bz(3), kk_diff(3) ! qq_ibz(3)
154 : : complex(gwp) :: ctmp_gwpc
155 : : character(len=fnlen) :: path
156 : : !arrays
157 : 0 : real(dp) :: n0(ebands%nsppol)
158 : 0 : real(dp),allocatable :: qlwl(:,:), kpg_kp(:,:), kpg_k(:,:), ug_kp(:,:,:), ug_k(:,:) !, cg_work(:,:)
159 : 0 : real(dp),allocatable :: work(:,:,:,:), e_mesh(:), e_args(:), wgt_mk(:), wgt_nkp(:,:) ! vcart_ibz(:,:,:,:),
160 : 0 : complex(gwp),allocatable :: cwork_ur(:), rhotwg_mn_x(:,:,:), rhotwg_mn_c(:,:,:), w_rhotwg_mn_c(:,:,:)
161 : 0 : complex(gwp),allocatable :: vc_sqrt_gx(:), ur_nkp(:,:), ur_mk(:,:), kxcg(:,:), mu_mn(:,:)
162 : 0 : complex(dp),allocatable :: w_ee(:,:) ! w_kkp(:,:,:,:)
163 : 0 : logical,allocatable :: bks_mask(:,:,:), keep_ur(:,:,:)
164 : 0 : type(fstab_t),target,allocatable :: fstab(:)
165 : : !************************************************************************
166 : :
167 [ # # ]: 0 : ABI_CHECK(psps%usepaw == 0, "PAW not implemented")
168 [ # # ]: 0 : ABI_CHECK(dtset%nsppol == 1, "nsppol 2 implemented!")
169 [ # # ]: 0 : ABI_CHECK(dtset%nspinor == 1, "nspinor 2 implemented!")
170 : :
171 : : ABI_UNUSED(mpi_enreg%nproc)
172 : :
173 [ # # ]: 0 : my_rank = xmpi_comm_rank(comm); nprocs = xmpi_comm_size(comm); units = [std_out, ab_out]
174 : 0 : call cwtime(cpu_all, wall_all, gflops_all, "start")
175 : :
176 : 0 : nsppol = ebands%nsppol; nspinor = ebands%nspinor
177 : :
178 : : ! Compute electron DOS.
179 [ # # ]: 0 : edos_intmeth = 2; if (dtset%prtdos /= 0) edos_intmeth = dtset%prtdos
180 : : edos_step = dtset%dosdeltae; edos_broad = dtset%tsmear
181 : 0 : edos_step = 0.01 * eV_Ha; edos_broad = 0.3 * eV_Ha
182 : 0 : edos = ebands%get_edos(cryst, edos_intmeth, edos_step, edos_broad, comm)
183 : :
184 : : ! Get DOS per spin channel
185 [ # # ]: 0 : n0(:) = edos%gef(1:edos%nsppol)
186 [ # # ]: 0 : if (my_rank == master) then
187 : 0 : call edos%print(units)
188 [ # # ]: 0 : path = strcat(dtfil%filnam_ds(4), "_EDOS")
189 : 0 : call edos%write(path)
190 : : end if
191 : :
192 : : ! Find Fermi surface k-points.
193 : : ! TODO: support kptopt, change setup of k-points if tetra: fist tetra weights then k-points on the Fermi surface!
194 [ # # # # : 0 : ABI_MALLOC(fstab, (nsppol))
# # # # #
# # # #
# ]
195 : 0 : call fstab_init(fstab, ebands, cryst, dtset, tetra, comm)
196 : 0 : call tetra%free()
197 : :
198 : 0 : bmin = huge(1); bmax = -1
199 [ # # ]: 0 : do spin=1,nsppol
200 : 0 : bmin = min(bmin, fstab(spin)%bmin)
201 : 0 : bmax = max(bmax, fstab(spin)%bmax)
202 : : end do
203 : 0 : max_nb = bmax - bmin + 1
204 : :
205 : : ! Build linear mesh for W(e,e').
206 : 0 : spin = 1
207 : : associate (fs => fstab(spin))
208 [ # # # # : 0 : e_min = minval(ebands%eig(fs%bmin,:,spin)) - 0.1_dp * eV_Ha
# # # # ]
209 [ # # # # : 0 : e_max = maxval(ebands%eig(fs%bmax,:,spin)) + 0.1_dp * eV_Ha
# # # # ]
210 : 0 : e_step = 0.002 * eV_Ha
211 : 0 : ne = 1 + (e_max - e_min) / e_step
212 [ # # # # : 0 : ABI_CALLOC(w_ee, (ne, ne))
# # # # #
# # # # #
# # # # #
# ]
213 [ # # # # : 0 : ABI_MALLOC(e_mesh, (ne))
# # # # #
# # # ]
214 [ # # # # : 0 : ABI_MALLOC(e_args, (ne))
# # # # ]
215 [ # # # # : 0 : ABI_MALLOC(wgt_mk, (ne))
# # # # ]
216 [ # # # # : 0 : ABI_MALLOC(wgt_nkp, (ne, max_nb))
# # # # #
# # # # #
# # ]
217 [ # # # # : 0 : e_mesh = linspace(e_min, e_max, ne)
# # ]
218 : : end associate
219 : :
220 : : ! ================
221 : : ! HANDLE SCREENING
222 : : ! ================
223 : 0 : screen_filepath = dtfil%fnameabi_scr
224 [ # # ]: 0 : ABI_CHECK(dtfil%fnameabi_scr /= ABI_NOFILE, "SCR file must be specified in input!")
225 : :
226 : : ! Read g-sphere for correlation and qmesh from SCR file.
227 : 0 : call get_hscr_qmesh_gsph(screen_filepath, dtset, cryst, hscr, qmesh, gsph_c, qlwl, comm)
228 : 0 : call hscr%print(units, dtset%prtvol, header="Header of the SCR file")
229 : :
230 : 0 : nqlwl = size(qlwl, dim=2)
231 [ # # ]: 0 : if (nqlwl == 0) then
232 : 0 : nqlwl=1
233 [ # # # # ]: 0 : ABI_MALLOC(qlwl,(3,nqlwl))
234 [ # # ]: 0 : qlwl(:,nqlwl)= GW_Q0_DEFAULT
235 : : write(msg,'(3a,i0,a,3f9.6)')&
236 : 0 : "The Header of the screening file does not contain the list of q-point for the optical limit ",ch10,&
237 : 0 : "Using nqlwl= ",nqlwl," and qlwl = ",qlwl(:,1)
238 : 0 : ABI_COMMENT(msg)
239 : : end if
240 : :
241 : : ! Init g-sphere for the exchange part from ecutsigx.
242 : 0 : call gsph_c%extend(cryst, dtset%ecutsigx, gsph_x)
243 : :
244 : : ! Initialize Coulomb term on the IBZ of the qmesh. Use largest G-sphere.
245 : 0 : call kmesh%init(cryst, wfk_hdr%nkpt, wfk_hdr%kptns, dtset%kptopt)
246 : 0 : faq = one / (cryst%ucvol*kmesh%nbz)
247 : :
248 : : ! TODO:
249 : : ! Here we sort the pp_mesh by stars so that we can split the pp wavevectors in blocks and therefore
250 : : ! reduce the number of wavevectors in the IBZ that must be stored in memory.
251 : : !call qmesh%pack_by_stars()
252 : :
253 : : ! The sum over k' is done with the k-mesh of the WFK file (normalized with faq) and q = k - k' is
254 : : ! taken from the qmesh of the SCR file. Hence the treatment of the Coulomb singularity in vcoul_init
255 : : ! (based on the k1-k2 mesh built from kmesh if qmesh%nbz /= kmesh%nbz) is consistent with this sum.
256 : : ! This is not the case in GWPT, where the sum runs over the q-mesh of the SCR file (see gwpt_run).
257 : 0 : npw_x = gsph_x%ng; npw_c = gsph_c%ng !; min_npw_xc = min(npw_x, npw_c); max_npw_xc = max(npw_x, npw_c)
258 [ # # ]: 0 : if (gsph_x%ng >= gsph_c%ng) then
259 : : call vcp%init(gsph_x, cryst, qmesh, kmesh, dtset%rcut, dtset%gw_icutcoul, dtset%vcutgeo, gsph_x%ng, &
260 : 0 : nqlwl, qlwl, comm)
261 : : else
262 : : call vcp%init(gsph_c, cryst, qmesh, kmesh, dtset%rcut, dtset%gw_icutcoul, dtset%vcutgeo, gsph_c%ng, &
263 : 0 : nqlwl, qlwl, comm)
264 : : end if
265 : :
266 [ # # ]: 0 : ABI_FREE(qlwl)
267 : :
268 [ # # ]: 0 : if (my_rank == master) then
269 : 0 : call kmesh%print(units, header="K-mesh for wavefunctions", prtvol=dtset%prtvol)
270 : 0 : call gsph_x%print(units, dtset%prtvol, header="G-sphere for exchange")
271 : 0 : call gsph_c%print(units, dtset%prtvol, header="G-sphere for correlation")
272 : 0 : call vcp%print(units, prtvol=dtset%prtvol)
273 : : end if
274 : 0 : call kmesh%free()
275 : :
276 [ # # ]: 0 : ABI_CHECK_IGE(npw_x, 1, "npw_x <= 1")
277 [ # # ]: 0 : ABI_CHECK_IGE(npw_c, 1, "npw_c <= 1")
278 : :
279 : : ! Note that in this case, the sphere is always Gamma-centered i.e. it does not depend on the qq wavevector
280 : 0 : kg_c => gsph_c%gvec(:, 1:npw_c)
281 : 0 : kg_x => gsph_x%gvec(:, 1:npw_x)
282 : :
283 : : ! Initialize the wave function descriptor.
284 : : ! Only wavefunctions on the FS are stored in wfd.
285 : : ! Need all k-points on the FS because of k+q, spin is not distributed for the time being.
286 : : ! It would be possible to reduce the memory allocated per MPI-rank via OpenMP.
287 : 0 : nkibz = wfk_hdr%nkpt
288 [ # # # # : 0 : ABI_MALLOC(nband, (nkibz, nsppol))
# # # # #
# # # # #
# # ]
289 [ # # # # : 0 : ABI_MALLOC(bks_mask, (dtset%mband, nkibz, nsppol))
# # # # #
# # # # #
# # # # #
# ]
290 [ # # # # : 0 : ABI_MALLOC(keep_ur, (dtset%mband, nkibz, nsppol))
# # # # #
# # # # #
# # ]
291 [ # # # # : 0 : nband = dtset%mband; bks_mask = .False.; keep_ur = .False.
# # # # #
# # # # #
# # ]
292 : :
293 [ # # ]: 0 : do spin=1,ebands%nsppol
294 : 0 : associate (fs => fstab(spin))
295 [ # # ]: 0 : do ik_bz=1,fs%nkfs
296 : 0 : ik_ibz = fs%indkk_fs(1, ik_bz)
297 : 0 : bstart_k = fs%bstart_cnt_ibz(1, ik_ibz); nband_k = fs%bstart_cnt_ibz(2, ik_ibz)
298 [ # # ]: 0 : bks_mask(bstart_k:bstart_k+nband_k-1, ik_ibz, spin) = .True.
299 : : end do
300 : : end associate
301 : : end do
302 : :
303 : : ! Impose istwfk=1 for all k points. This is also done in respfn (see inkpts)
304 : : ! wfd_read_wfk will handle a possible conversion if WFK contains istwfk /= 1.
305 [ # # # # : 0 : ABI_MALLOC(wfd_istwfk, (nkibz))
# # # # #
# # # ]
306 [ # # ]: 0 : wfd_istwfk = 1
307 : :
308 : : call wfd%init(cryst, pawtab, psps, keep_ur, dtset%mband, nband, nkibz, nsppol, bks_mask,&
309 : : dtset%nspden, nspinor, dtset%ecut, dtset%ecutsm, dtset%dilatmx, wfd_istwfk, ebands%kptns, ngfft,&
310 : 0 : dtset%nloalg, dtset%prtvol, dtset%pawprtvol, comm)
311 : :
312 : 0 : call pstat_proc%print(_PSTAT_ARGS_)
313 : :
314 [ # # ]: 0 : ABI_FREE(nband)
315 [ # # ]: 0 : ABI_FREE(bks_mask)
316 [ # # ]: 0 : ABI_FREE(keep_ur)
317 [ # # ]: 0 : ABI_FREE(wfd_istwfk)
318 : :
319 : : ! Read wavefunctions.
320 : 0 : call wfd%read_wfk(wfk0_path, iomode_from_fname(wfk0_path))
321 : :
322 : : ! TODO
323 : : ! mpw is the maximum number of plane-waves over k and k+q where k and k+q are in the BZ.
324 : : ! We also need the max components of the G-spheres (k, k+q) in order to allocate the workspace array work
325 : : ! that will be used to symmetrize the wavefunctions in G-space.
326 : : ! This is the maximum number of PWs for all possible k+q treated.
327 : :
328 : 0 : call ephtk_get_mpw_gmax(nkibz, wfk_hdr%kptns, dtset%ecut, cryst%gmet, mpw, gmax, comm)
329 : : ! FIXME
330 : : !gmax = 2 * gmax
331 : :
332 : 0 : call ddkop%init(dtset, cryst, pawtab, psps, wfd%mpi_enreg, mpw, wfd%ngfft)
333 : :
334 : : #if 0
335 : : call cwtime(cpu, wall, gflops, "start", msg=" Computing v_nk matrix elements for all states on the FS...")
336 : :
337 : : ABI_CALLOC(vcart_ibz, (3, bmin:bmax, nkibz, nsppol))
338 : : ABI_MALLOC(cg_work, (2, mpw * nspinor))
339 : :
340 : : cnt = 0
341 : : do spin=1,nsppol
342 : : associate (fs => fstab(spin))
343 : : do ik_ibz=1,ebands%nkpt
344 : : kk_ibz = ebands%kptns(:, ik_ibz)
345 : : npw_k = wfd%npwarr(ik_ibz); istwf_k = wfd%istwfk(ik_ibz)
346 : : ! NB: The two checks below are global --> all procs will cycle.
347 : : if (all(bks_mask(:, ik_ibz, spin) .eqv. .False.)) cycle
348 : : if (npw_k == 1) cycle
349 : : cnt = cnt + 1; if (mod(cnt, nproc) /= my_rank) cycle ! MPI parallelism.
350 : :
351 : : call ddkop%setup_spin_kpoint(dtset, cryst, psps, spin, kk_ibz, istwf_k, npw_k, wfd%kdata(ik_ibz)%kg_k)
352 : :
353 : : do band_k=fs%bmin,fs%bmax
354 : : if (.not. bks_mask(band_k, ik_ibz, spin)) cycle
355 : : call wfd%copy_cg(band_k, ik_ibz, spin, cg_work)
356 : : eig0nk = ebands%eig(band_k, ik_ibz, spin)
357 : : vk = ddkop%get_vdiag(eig0nk, istwf_k, npw_k, wfd%nspinor, cg_work, cwaveprj0)
358 : : vcart_ibz(:, band_k, ik_ibz, spin) = vk
359 : : end do
360 : : end do
361 : : end associate
362 : : end do ! spin
363 : :
364 : : call xmpi_sum(vcart_ibz, comm, ierr)
365 : : ABI_FREE(cg_work)
366 : : ABI_FREE(vcart_ibz)
367 : : call cwtime_report(" Velocities", cpu, wall, gflops)
368 : : #endif
369 : 0 : call ddkop%free()
370 : :
371 : 0 : call ngfft_seq(work_ngfft, gmax)
372 [ # # # # : 0 : ABI_MALLOC(work, (2, work_ngfft(4), work_ngfft(5), work_ngfft(6)))
# # # # #
# # # # #
# # # # #
# ]
373 : :
374 : : ! FFT meshes from input file, not necessarily equal to the ones found in the external files.
375 : : ! NB: ur arrays are always allocated with nfft and not with product(ngfft(4:6)).
376 [ # # ]: 0 : nfftf = product(ngfftf(1:3)); mgfftf = maxval(ngfftf(1:3))
377 [ # # # # ]: 0 : nfft = product(ngfft(1:3)) ; mgfft = maxval(ngfft(1:3))
378 : :
379 [ # # ]: 0 : call print_ngfft([std_out], ngfft, header="FFT mesh")
380 : :
381 : : ! Set the FFT mesh
382 : 0 : call wfd%change_ngfft(cryst, psps, ngfft)
383 : 0 : call wfd%print(units, header="Wavefunctions for GWPT calculation.")
384 : :
385 [ # # # # : 0 : ABI_MALLOC(gbound_kp, (2*mgfft+8, 2))
# # # # #
# # # # #
# # ]
386 [ # # # # : 0 : ABI_MALLOC(gbound_k, (2*mgfft+8, 2))
# # # # ]
387 [ # # # # : 0 : ABI_MALLOC(gbound_c, (2*mgfft+8, 2))
# # # # ]
388 [ # # # # : 0 : ABI_MALLOC(gbound_x, (2*mgfft+8, 2))
# # # # ]
389 : :
390 : : ! ====================================
391 : : ! This is the g-sphere for W_{gg'}(qq)
392 : : ! ====================================
393 : : ! Init g-sphere for the exchange part from ecutsigx.
394 : 0 : call sphereboundary(gbound_c, istwfk1, kg_c, mgfft, npw_c)
395 : 0 : call sphereboundary(gbound_x, istwfk1, kg_x, mgfft, npw_x)
396 : :
397 [ # # # # : 0 : ABI_MALLOC(vc_sqrt_gx, (npw_x))
# # # # #
# # # ]
398 : :
399 : : ! Allocate g-vectors centered on k, k'
400 [ # # # # : 0 : ABI_MALLOC(kg_k, (3, mpw))
# # # # #
# # # ]
401 [ # # # # : 0 : ABI_MALLOC(kg_kp, (3, mpw))
# # # # ]
402 : :
403 : 0 : spin = 1
404 : : associate (fs => fstab(spin))
405 [ # # # # : 0 : ABI_MALLOC_OR_DIE(ur_nkp, (nfft*nspinor, fs%bmax-fs%bmin+1), ierr)
# # # # #
# # # # #
# # # # ]
406 [ # # # # : 0 : ABI_MALLOC_OR_DIE(ur_mk, (nfft*nspinor, fs%bmax-fs%bmin+1), ierr)
# # # # #
# # # # #
# # # # #
# ]
407 [ # # # # : 0 : ABI_MALLOC(cwork_ur, (nfft*nspinor))
# # # # #
# # # ]
408 : :
409 : : ! Read symmetrized em1 from file
410 : 0 : id_required = 4; ikxc = 0; approx_type = 0; option_test = 0; nkxc = 0
411 [ # # # # : 0 : ABI_MALLOC(kxcg, (nfftf, nkxc))
# # # # ]
412 : 0 : scr_iomode = iomode_from_fname(screen_filepath)
413 : :
414 : 0 : call epsm1%from_file(screen_filepath, qmesh%nibz, npw_c, comm)
415 : : call epsm1%mkdump(vcp, npw_c, kg_c, nkxc, kxcg, id_required, approx_type, &
416 : 0 : ikxc, option_test, dtfil%fnameabo_scr, scr_iomode, nfftf, ngfft, comm)
417 : :
418 : 0 : call epsm1%malloc_epsm1_qbz(npw_c, epsm1%nomega)
419 [ # # ]: 0 : ABI_FREE(kxcg)
420 : :
421 : 0 : call pstat_proc%print(_PSTAT_ARGS_)
422 : :
423 : 0 : wkk_mode = 1
424 : : ! Master creates the netcdf file used to store the results of the calculation.
425 [ # # ]: 0 : if (my_rank == master) then
426 [ # # ]: 0 : NCF_CHECK(nctk_open_create(root_ncid, strcat(dtfil%filnam_ds(4), "_WKK.nc") , xmpi_comm_self))
427 [ # # ]: 0 : NCF_CHECK(cryst%ncwrite(root_ncid))
428 [ # # ]: 0 : NCF_CHECK(ebands%ncwrite(root_ncid))
429 [ # # ]: 0 : NCF_CHECK(edos%ncwrite(root_ncid))
430 : :
431 : : !ncerr = nctk_def_dims(ncid, [ &
432 : : ! nctkdim_t("smat_bsize1", smat_bsize1), nctkdim_t("smat_bsize2", smat_bsize2) &
433 : : ! ], defmode=.True.)
434 : : !NCF_CHECK(ncerr)
435 : : !ncerr = nctk_def_iscalars(ncid, [character(len=nctk_slen) :: &
436 : : ! "gwr_completed", "scf_iteration" &
437 : : !])
438 : : !NCF_CHECK(ncerr)
439 : : !ncerr = nctk_def_dpscalars(ncid, [character(len=nctk_slen) :: &
440 : : ! "wr_step", "gwr_boxcutmin", "cosft_duality_error", "regterm" &
441 : : !])
442 : : !NCF_CHECK(ncerr)
443 : :
444 : : !ncerr = nctk_def_arrays(ncid, [ &
445 : : ! nctkarr_t("gwr_task", "char", "character_string_length"), &
446 : : !])
447 : : !NCF_CHECK(ncerr)
448 : : !NCF_CHECK(nctk_set_datamode(ncid))
449 : :
450 : : !ncerr = nctk_write_iscalars(ncid, [character(len=nctk_slen) :: &
451 : : ! "gwr_completed", "sig_diago", "b1gw", "b2gw", "symsigma", "symchi", "scf_iteration"], &
452 : : ! [0, merge(1, 0, gwr%sig_diago), gwr%b1gw, gwr%b2gw, gwr%dtset%symsigma, dtset%symchi, gwr%scf_iteration])
453 : : !NCF_CHECK(ncerr)
454 : :
455 : : !ncerr = nctk_write_dpscalars(ncid, [character(len=nctk_slen) :: &
456 : : ! "wr_step", "ecuteps", "ecut", "ecutwfn", "ecutsigx", "gwr_boxcutmin", &
457 : : ! "cosft_duality_error", "regterm"], &
458 : : ! [gwr%wr_step, dtset%ecuteps, dtset%ecut, dtset%ecutwfn, dtset%ecutsigx, dtset%gwr_boxcutmin, &
459 : : ! gwr%ft_max_error(1), gwr%ft_max_error(2), gwr%ft_max_error(3), gwr%cosft_duality_error, regterm &
460 : : ! ])
461 : : !NCF_CHECK(ncerr)
462 : : !NCF_CHECK(nf90_put_var(ncid, vid("gwr_task"), trim(dtset%gwr_task)))
463 : :
464 [ # # ]: 0 : NCF_CHECK(nf90_close(root_ncid))
465 : : end if
466 : 0 : call xmpi_barrier(comm)
467 : :
468 : : ! Open WKK.nc file and go to data mode.
469 [ # # ]: 0 : NCF_CHECK(nctk_open_modify(root_ncid, strcat(dtfil%filnam_ds(4), "_WKK.nc") , comm))
470 [ # # ]: 0 : NCF_CHECK(nctk_set_datamode(root_ncid))
471 : :
472 : 0 : mg0 = [1, 1, 1]
473 : 0 : remove_exchange = .True.
474 : :
475 [ # # # # : 0 : ABI_CALLOC(mu_mn, (bmin:bmax, bmin:bmax))
# # # # #
# # # # #
# # # # #
# ]
476 : :
477 : : ! Loop over k'-points in the energy window.
478 [ # # # # ]: 0 : do ikp_bz=1,fs%nkfs
479 : 0 : my_ikp = ikp_bz; my_nkp = fs%nkfs
480 : :
481 [ # # # # : 0 : print_time_kp = my_rank == 0 .and. (ikp_bz <= LOG_MODK .or. mod(ikp_bz, LOG_MODK) == 0)
# # ]
482 : 0 : if (print_time_kp) call cwtime(cpu_kp, wall_kp, gflops_kp, "start")
483 : :
484 : : ! The k-point and the symmetries relating the BZ k-point to the IBZ.
485 [ # # ]: 0 : kp_bz = fs%kpts(:, ikp_bz)
486 : 0 : ikp_ibz = fs%indkk_fs(1, ikp_bz) ; isym_kp = fs%indkk_fs(2, ikp_bz)
487 [ # # ]: 0 : trev_kp = fs%indkk_fs(6, ikp_bz); g0_kp = fs%indkk_fs(3:5,ikp_bz)
488 [ # # # # : 0 : isirr_kp = (isym_kp == 1 .and. trev_kp == 0 .and. all(g0_kp == 0))
# # # # ]
489 [ # # ]: 0 : mapl_kp = fs%indkk_fs(:, ikp_bz)
490 : :
491 [ # # ]: 0 : kp_ibz = ebands%kptns(:,ikp_ibz)
492 : 0 : istwf_kp_ibz = wfd%istwfk(ikp_ibz); npw_kp_ibz = wfd%npwarr(ikp_ibz)
493 : :
494 : : ! Get npw_kp, kg_kp for this kp_bz.
495 : : call wfd%get_gvec_gbound(cryst%gmet, dtset%ecut, kp_bz, ikp_ibz, isirr_kp, dtset%nloalg, & ! in
496 [ # # ]: 0 : istwf_kp, npw_kp, kg_kp, nkpg_kp, kpg_kp, gbound_kp) ! out
497 : :
498 : : ! Rotate from IBZ to BZ and compute ur_nkp for all n_kp bands.
499 : 0 : bstart_kp = fs%bstart_cnt_ibz(1, ikp_ibz); nband_kp = fs%bstart_cnt_ibz(2, ikp_ibz); bstop_kp = bstart_kp + nband_kp - 1
500 [ # # # # : 0 : ABI_CALLOC(ug_kp, (2, npw_kp*nspinor, nband_kp))
# # # # #
# # # # #
# # # # #
# # # ]
501 : : call wfd%rotate_cg(bstart_kp, nband_kp, spin, kp_ibz, npw_kp, kg_kp, istwf_kp, &
502 : 0 : cryst, mapl_kp, gbound_kp, work_ngfft, work, ug_kp, urs_kbz=ur_nkp)
503 : 0 : do ii=1, nband_kp
504 : : !print *, "isirr_kp:", isirr_kp
505 : : !print *, "g:", sum(ug_kp(1,:,ii)**2 + ug_kp(2,:,ii)**2)
506 : : !print *, "r:", sum(abs(ur_nkp(:,ii)) ** 2) / nfft ! * cryst%ucvol /
507 : : end do
508 : :
509 : : ! Loop over k-points in the energy window.
510 [ # # ]: 0 : do ik_bz=1,fs%nkfs
511 : :
512 [ # # ]: 0 : mapl_k = fs%indkk_fs(:, ik_bz)
513 [ # # ]: 0 : ik_ibz = mapl_k(1); isym_k = mapl_k(2); trev_k = mapl_k(6); g0_k = mapl_k(3:5)
514 [ # # # # : 0 : isirr_k = (isym_k == 1 .and. trev_k == 0 .and. all(g0_k == 0))
# # # # ]
515 [ # # ]: 0 : kk_ibz = ebands%kptns(:, ik_ibz)
516 : 0 : istwf_k_ibz = wfd%istwfk(ik_ibz); npw_k_ibz = wfd%npwarr(ik_ibz)
517 : :
518 [ # # ]: 0 : kk_bz = fs%kpts(:, ik_bz)
519 [ # # ]: 0 : kk_ibz = ebands%kptns(:,ik_ibz)
520 : :
521 : : ! Get npw_k, kg_k for this k.
522 : : call wfd%get_gvec_gbound(cryst%gmet, dtset%ecut, kk_bz, ik_ibz, isirr_k, dtset%nloalg, & ! in
523 [ # # ]: 0 : istwf_k, npw_k, kg_k, nkpg_k, kpg_k, gbound_k) ! out
524 [ # # ]: 0 : ABI_FREE(kpg_k)
525 : :
526 : 0 : bstart_k = fs%bstart_cnt_ibz(1, ik_ibz); nband_k = fs%bstart_cnt_ibz(2, ik_ibz); bstop_k = bstart_k + nband_k - 1
527 [ # # # # : 0 : ABI_MALLOC(ug_k, (2, npw_k*nspinor))
# # # # #
# # # ]
528 [ # # # # : 0 : ABI_MALLOC(rhotwg_mn_x, (npw_x*nspinor, nband_k, nband_kp))
# # # # #
# # # # #
# # # # #
# ]
529 : :
530 [ # # ]: 0 : do m_k=bstart_k, bstop_k
531 : : ! Rotate from IBZ to BZ and compute ur_mk
532 : 0 : im_k = m_k - bstart_k + 1
533 : : call wfd%rotate_cg(m_k, ndat1, spin, kk_ibz, npw_k, kg_k, istwf_k, &
534 : 0 : cryst, mapl_k, gbound_k, work_ngfft, work, ug_k, urs_kbz=ur_mk(:,im_k))
535 : :
536 [ # # ]: 0 : do n_kp=bstart_kp, bstop_kp
537 : 0 : in_kp = n_kp - bstart_kp + 1
538 [ # # # # ]: 0 : cwork_ur = conjg(ur_nkp(:,in_kp)) * ur_mk(:,im_k)
539 : 0 : call fft_ur(npw_x, nfft, nspinor, ndat1, mgfft, ngfft, istwfk1, kg_x, gbound_x, cwork_ur, rhotwg_mn_x(:, im_k, in_kp))
540 : 0 : call sigtk_multiply_by_vc_sqrt("N", npw_x, nspinor, ndat1, vc_sqrt_gx, rhotwg_mn_x(:, im_k, in_kp))
541 : : end do
542 : : end do
543 : :
544 [ # # ]: 0 : ABI_FREE(ug_k)
545 : :
546 : : ! Identify q and G0 where q + G0 = k_GW - k_i
547 [ # # ]: 0 : kk_diff = kk_bz - kp_bz
548 : :
549 : 0 : call findqg0(iq_bz, g0_qq, kk_diff, qmesh%nbz, qmesh%bz, mG0)
550 : :
551 : : ! Find the corresponding irred qq-point in qmesh.
552 : 0 : call qmesh%get_bz_item(iq_bz, qq_bz, iq_ibz, isym_qq, itim_qq)
553 : 0 : qq_is_gamma = sum(qq_bz**2) < tol14
554 : :
555 : : ! Find k + q in the extended zone and extract symmetry info.
556 : : ! Be careful here because there are two umklapp vectors to be considered as:
557 : : !
558 : : ! k + q = k_bz + g0_bz = IS(k_ibz) + g0_ibz + g0_bz
559 : : !
560 : : ! Get Fourier components of the Coulomb interaction in the BZ
561 : : ! In 3D systems, neglecting umklapp: vc(Sq,sG) = vc(q,G) = 4pi/|q+G|**2
562 : : ! The same relation holds for 0-D systems, but not in 1-D or 2D systems. It depends on S.
563 : : ! NOTE: vc_sqrt_gx is dimensioned with npw_x --> use rottb table from gsph_x.
564 [ # # ]: 0 : do ig=1,npw_x
565 : 0 : vc_sqrt_gx(gsph_x%rottb(ig, itim_qq, isym_qq)) = vcp%vc_sqrt(ig, iq_ibz)
566 : : end do
567 : :
568 : 0 : call epsm1%rotate_iqbz(iq_bz, hscr%nomega, npw_c, gsph_c, qmesh, remove_exchange)
569 : :
570 : : ! Compute \sum_{gg'} (M_g^{kn,k'm})^* W_{gg'}(q) M_g'^{kn,k'm} with k' = k + q.
571 : : ! TODO: nspinor
572 [ # # # # : 0 : ABI_MALLOC(rhotwg_mn_c, (npw_c*nspinor, nband_k, nband_kp))
# # # # #
# # # # #
# # # # #
# ]
573 [ # # # # : 0 : rhotwg_mn_c(:,:,:) = rhotwg_mn_x(1:npw_c,:,:)
# # ]
574 : :
575 : 0 : ncols = nband_k * nband_kp
576 [ # # # # : 0 : ABI_MALLOC(w_rhotwg_mn_c, (npw_c*nspinor, nband_k, nband_kp))
# # # # #
# # # # #
# # # # #
# ]
577 : :
578 : : call xgemm("N", "N", npw_c, ncols, npw_c, cone_gw, epsm1%epsm1_qbz(:,:,1), npw_c, rhotwg_mn_c(:,:,1), npw_c*nspinor, &
579 [ # # # # : 0 : czero_gw, w_rhotwg_mn_c(:,:,1), npw_c)
# # # # #
# ]
580 : :
581 : : !do in_kp=1, nband_k
582 : : ! do im_k=1, nband_kq
583 : : ! w_rhotwg_mn_c(:,in_kp,im_k) = matmul(epsm1%epsm1_qbz(:,:,1), rhotwg_mn_c(:,im_k,in_kp))
584 : : ! end do
585 : : !end do
586 : : !print *, " rhotwg_mn_c", maxval(abs(rhotwg_mn_c)), maxloc(abs(rhotwg_mn_c))
587 : : !print *, "w_rhotwg_mn_c:", maxval(abs(w_rhotwg_mn_c))
588 : :
589 : : ! Precompute delta(e - e_nkp - ef).
590 [ # # ]: 0 : do n_kp=bstart_kp, bstop_kp
591 : 0 : in_kp = n_kp - bstart_kp + 1
592 : 0 : e_nkp = ebands%eig(n_kp, ik_ibz, spin)
593 [ # # # # ]: 0 : e_args = e_mesh - e_nkp
594 : 0 : smear_nkp = 0.1_dp * eV_Ha
595 [ # # ]: 0 : wgt_nkp(:, im_k) = gaussian(e_args, smear_nkp)
596 : : end do
597 : :
598 [ # # ]: 0 : do n_kp=bstart_kp, bstop_kp
599 : 0 : in_kp = n_kp - bstart_kp + 1
600 : 0 : e_nkp = ebands%eig(n_kp, ikp_ibz, spin)
601 : : !fs%tetra_wtk(n_kp - fs%bmin + 1, ikp_ibz)
602 [ # # ]: 0 : do m_k=bstart_k, bstop_k
603 : 0 : im_k = m_k - bstart_k + 1
604 : 0 : e_mk = ebands%eig(m_k, ik_ibz, spin)
605 : : !fs%tetra_wtk(m_k - fs%bmin + 1, ik_ibz)
606 : 0 : smear_mk = 0.1_dp * eV_Ha
607 : 0 : ctmp_gwpc = xdotc(npw_c*nspinor, rhotwg_mn_c(:,im_k,in_kp), 1, w_rhotwg_mn_c(:,im_k,in_kp), 1)
608 : : if (remove_exchange) then
609 : 0 : ctmp_gwpc = ctmp_gwpc + xdotc(npw_x*nspinor, rhotwg_mn_x(:,im_k,in_kp), 1, rhotwg_mn_x(:,im_k,in_kp), 1)
610 : : end if
611 : 0 : ctmp_gwpc = faq * ctmp_gwpc
612 : : !print *, "ctmp_gwpc:", ctmp_gwpc
613 : : !mu_mn(m_k, n_kp) = mu_mn(m_k, n_kp) + &
614 : : ! ctmp_gwpc * gaussian(e_mk - ebands%fermie, smear_mk) * gaussian(e_mp, smear_mk)
615 : : ! Computes: A := A + alpha * x * y^T
616 : : !w_ee =
617 : : !call dger(ne, ne, alpha, x, incx, y, incy, w_ee, lda)
618 : : end do
619 : : end do
620 : :
621 [ # # ]: 0 : ABI_FREE(rhotwg_mn_x)
622 [ # # ]: 0 : ABI_FREE(rhotwg_mn_c)
623 [ # # ]: 0 : ABI_FREE(w_rhotwg_mn_c)
624 : : end do ! ik_bz
625 : :
626 [ # # ]: 0 : if (print_time_kp) then
627 : 0 : call inds2str(1, "My kp-point", my_ikp, my_nkp, fs%nkfs, msg)
628 [ # # ]: 0 : call cwtime_report(msg, cpu_kp, wall_kp, gflops_kp); if (ikp_bz == LOG_MODK) call wrtout(std_out, "...", do_flush=.True.)
629 : : end if
630 : :
631 [ # # ]: 0 : ABI_FREE(ug_kp)
632 [ # # ]: 0 : ABI_FREE(kpg_kp)
633 : : end do ! ikp_bz
634 : : end associate
635 : :
636 [ # # ]: 0 : ABI_FREE(ur_nkp)
637 [ # # ]: 0 : ABI_FREE(ur_mk)
638 : : !end do ! my_is
639 : :
640 : 0 : call cwtime_report(" wkk calculation", cpu_all, wall_all, gflops_all, end_str=ch10)
641 : :
642 : 0 : call xmpi_sum(w_ee, comm, ierr)
643 : 0 : call xmpi_sum(mu_mn, comm, ierr)
644 : :
645 : : if (my_rank == master) then
646 : : !do m_k=bmin, bmax
647 : : ! do n_kp=bmin, bmax
648 : : ! call print_arr(units, mu_mn, max_r=max_nb, max_c=max_nb)
649 : : ! end do
650 : : !end do
651 : : !call wrtout(units, sjoin("Re mu: ", ftoa(sum(real(mu_mn)))))
652 : : !call wrtout(units, sjoin("Im mu: ", ftoa(sum(aimag(mu_mn)))))
653 : : end if
654 : :
655 : : ! Set vqk_completed to 1 so that we can easily check if restarted is needed.
656 : : !if (wkk_mode /= 0) then
657 : : !if (my_rank == master) then
658 : : ! NCF_CHECK(nf90_put_var(root_ncid, root_vid("vqk_completed"), 1))
659 : : !end if
660 [ # # ]: 0 : NCF_CHECK(nf90_close(root_ncid))
661 : 0 : call xmpi_barrier(comm)
662 : : !end if
663 : :
664 : : ! Output some of the results to ab_out for testing purposes
665 : : !call gstore%print_for_abitests(dtset)
666 : :
667 : : ! Free memory
668 [ # # ]: 0 : ABI_FREE(kg_kp)
669 [ # # ]: 0 : ABI_FREE(kg_k)
670 [ # # ]: 0 : ABI_FREE(cwork_ur)
671 [ # # ]: 0 : ABI_FREE(work)
672 [ # # ]: 0 : ABI_FREE(gbound_kp)
673 [ # # ]: 0 : ABI_FREE(gbound_k)
674 [ # # ]: 0 : ABI_FREE(gbound_c)
675 [ # # ]: 0 : ABI_FREE(gbound_x)
676 [ # # ]: 0 : ABI_FREE(vc_sqrt_gx)
677 [ # # ]: 0 : ABI_FREE(w_ee)
678 [ # # ]: 0 : ABI_FREE(e_mesh)
679 [ # # ]: 0 : ABI_FREE(e_args)
680 [ # # ]: 0 : ABI_FREE(wgt_mk)
681 [ # # ]: 0 : ABI_FREE(wgt_nkp)
682 [ # # ]: 0 : ABI_FREE(mu_mn)
683 : :
684 : 0 : call wfd%free(); call vcp%free(); call qmesh%free(); call gsph_x%free(); call gsph_c%free(); call hscr%free(); call edos%free()
685 : 0 : call epsm1%free()
686 : :
687 [ # # ]: 0 : do spin=1,ebands%nsppol
688 : 0 : call fstab(spin)%free()
689 : : end do
690 [ # # # # : 0 : ABI_FREE(fstab)
# # # # #
# # # # #
# # # # #
# # # #
# ]
691 : :
692 : 0 : call xmpi_barrier(comm) ! This to make sure that the parallel output of WKK is completed
693 [ # # # # : 0 : call cwtime_report(" wkk_run: MPI barrier before returning.", cpu_all, wall_all, gflops_all, end_str=ch10, comm=comm)
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # ]
694 : :
695 : : contains
696 : :
697 : 0 : subroutine inds2str(ilevel, prefix, my_ik, my_nk, nk_tot, out_str)
698 : : character(len=*),intent(in) :: prefix
699 : : integer,intent(in) :: ilevel, my_ik, my_nk, nk_tot
700 : : character(len=*),intent(out) :: out_str
701 : :
702 [ # # # # ]: 0 : out_str = sjoin(prefix, itoa(my_ik), "/", itoa(my_nk), "[", itoa(nk_tot), "]")
703 [ # # # # : 0 : out_str = repeat(' ', 4 * ilevel) // trim(out_str)
# # ]
704 : 0 : end subroutine inds2str
705 : :
706 : : integer function root_vid(var_name)
707 : : character(len=*),intent(in) :: var_name
708 : : root_vid = nctk_idname(root_ncid, var_name)
709 : : end function root_vid
710 : :
711 : : !integer function spin_vid(var_name)
712 : : ! character(len=*),intent(in) :: var_name
713 : : ! spin_vid = nctk_idname(spin_ncid, var_name)
714 : : !end function spin_vid
715 : :
716 : : end subroutine wkk_run
717 : : !!***
718 : :
719 : : end module m_wkk
720 : : !!***
|
