Branch data Line data Source code
1 : : !!****m*ABINIT/m_lwf_primitive_potential
2 : : !! NAME
3 : : !! m_lwf_primitive_potential
4 : : !!
5 : : !! FUNCTION
6 : : !! This module contains the example lwf primitive potential. It is not the
7 : : !! one really used. The purpose of this is to show how to extend multibinit
8 : : !! to new degree of freedom
9 : : !!
10 : : !! Datatypes:
11 : : !! LWF_primitive_potential_t
12 : : !!
13 : : !! Subroutines:
14 : : !!
15 : : !! COPYRIGHT
16 : : !! Copyright (C) 2001-2026 ABINIT group (hexu)
17 : : !! This file is distributed under the terms of the
18 : : !! GNU General Public License, see ~abinit/COPYING
19 : : !! or http://www.gnu.org/copyleft/gpl.txt .
20 : : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
21 : : !!
22 : : !! SOURCE
23 : :
24 : :
25 : : #if defined HAVE_CONFIG_H
26 : : #include "config.h"
27 : : #endif
28 : : #include "abi_common.h"
29 : :
30 : : #define HAVE_NETCDF 1
31 : :
32 : : module m_lwf_primitive_potential
33 : : use, intrinsic :: iso_c_binding
34 : : use defs_basis
35 : : use m_abicore
36 : : use m_errors
37 : : use m_nctk
38 : : #if defined HAVE_NETCDF
39 : : use netcdf
40 : : #endif
41 : : use m_xmpi
42 : : use m_mathfuncs, only: eigensh
43 : : use m_multibinit_dataset, only: multibinit_dtset_type
44 : : use m_multibinit_cell, only: mbcell_t, mbsupercell_t
45 : : use m_primitive_potential, only: primitive_potential_t
46 : : use m_abstract_potential, only: abstract_potential_t
47 : : use m_dynamic_array, only: int2d_array_type
48 : : use m_supercell_maker, only: supercell_maker_t
49 : : use m_spmat_ndcoo, only: ndcoo_mat_t
50 : : use m_lwf_potential, only: lwf_potential_t
51 : : implicit none
52 : : private
53 : : !!***
54 : :
55 : : !-------------------------------------------------------------------!
56 : : ! An harmonic potential which has only the IFC (with no dipole-dipole)
57 : : !
58 : : ! IFC is written in a coefficient matrix M(R, i, j)=val
59 : : !-------------------------------------------------------------------!
60 : : type, public, extends(primitive_potential_t):: lwf_primitive_potential_t
61 : : integer:: nlwf = 0, natom = 0, nR = 0 ! number of LWF
62 : : type(ndcoo_mat_t):: coeff ! A N-dimensional COO matrix.
63 : : ! The indices are (ind_R, i, j). Note that xyz is included in i and j.
64 : : real(dp), allocatable:: lwf_masses(:)
65 : :
66 : : integer:: onebody_nterm = 0
67 : :
68 : : real(dp), allocatable:: lattice_coeffs(:,:, :) ! (natom*3, nlwf, nR)
69 : :
70 : : integer, allocatable:: onebody_i(:), onebody_order(:)
71 : : real(dp), allocatable:: onebody_val(:)
72 : :
73 : : type(ndcoo_mat_t):: coeff_twobody ! Two body interaction parameters of higher order.
74 : : ! The indices are (indR, i, j, orderi, orderj), H(R, i, j)= val A_i^(orderi) A_j^(orderj)
75 : : integer, allocatable:: Rlist(:,:) ! The list of R points (3, number of R-points)
76 : :
77 : : ! For testing the bounding.
78 : : !logical:: has_self_bound_term= .False.
79 : : !integer:: self_bound_order = 0
80 : : !real(dp):: self_bound_coeff = 0.0_dp
81 : :
82 : : !. Rlist(:, ind_R) is a R-vector.
83 : : real(dp):: ref_energy = 0.0 ! reference energy
84 : :
85 : : logical:: as_lattice_anharmonic=.False.
86 : :
87 : : contains
88 : : procedure:: initialize
89 : : procedure:: finalize
90 : : procedure:: load_from_files ! load potential from files listed in files file
91 : : procedure:: load_from_netcdf ! load potential from a netcdf file
92 : : procedure:: fill_supercell ! fill a supercell potential
93 : : procedure:: get_hamk ! generate hamiltonian for one k point.
94 : : procedure:: get_eigen ! eigen values and eigen vectors
95 : : !procedure:: add_self_bound_term ! add a self bound term
96 : : end type lwf_primitive_potential_t
97 : :
98 : :
99 : : contains
100 : :
101 : : !-------------------------------------------------------------------!
102 : : ! Initialize
103 : : ! Input:
104 : : ! nlwf: number of Wannier functions
105 : : ! primcell: the reference primitive cell.
106 : : !-------------------------------------------------------------------!
107 : 1 : subroutine initialize(self, primcell)
108 : : class(lwf_primitive_potential_t), intent(inout):: self
109 : : type(mbcell_t), target, intent(inout):: primcell
110 : : !integer, intent(in):: nlwf
111 : 1 : self%primcell=>primcell
112 : 1 : self%label="lwf_primitive_potential"
113 : 1 : self%has_spin=.False.
114 : 1 : self%has_displacement=.False.
115 : 1 : self%has_strain=.False.
116 : 1 : self%has_lwf=.True.
117 : 1 : end subroutine initialize
118 : :
119 : :
120 : : !-------------------------------------------------------------------!
121 : : ! Finalize
122 : : !-------------------------------------------------------------------!
123 : 1 : subroutine finalize(self)
124 : : class(lwf_primitive_potential_t), intent(inout):: self
125 : 1 : call self%coeff%finalize()
126 : 1 : call self%coeff_twobody%finalize()
127 [ + - ]: 1 : ABI_SFREE(self%Rlist)
128 : 1 : nullify(self%primcell)
129 : 1 : self%nlwf = 0
130 : 1 : self%label="Destroyed lwf_primitive_potential"
131 : 1 : call self%primitive_potential_t%finalize()
132 : :
133 [ + - ]: 1 : ABI_SFREE(self%lattice_coeffs)
134 [ + - ]: 1 : ABI_SFREE(self%lwf_masses)
135 : 1 : self%natom = 0
136 : 1 : self%nR = 0
137 : :
138 [ + - ]: 1 : if (self%onebody_nterm /= 0) then
139 [ + - ]: 1 : ABI_SFREE(self%onebody_i)
140 [ + - ]: 1 : ABI_SFREE(self%onebody_val)
141 [ + - ]: 1 : ABI_SFREE(self%onebody_order)
142 : : end if
143 : 1 : self%onebody_nterm = 0
144 : 1 : end subroutine finalize
145 : :
146 : :
147 : : !-------------------------------------------------------------------!
148 : : ! Load from files:
149 : : !> params: parameters
150 : : !> fnames: file names from files file
151 : : !-------------------------------------------------------------------!
152 : 1 : subroutine load_from_files(self, params, fnames)
153 : : class(lwf_primitive_potential_t), intent(inout):: self
154 : : type(multibinit_dtset_type), intent(in):: params
155 : : character(len = fnlen), intent(in):: fnames(:)
156 : 1 : call self%load_from_netcdf(fnames(1))
157 : : !call self%add_self_bound_term(params%lwf_self_bound_order, &
158 : : ! & params%lwf_self_bound_coeff)
159 : 1 : self%as_lattice_anharmonic= (params%latt_lwf_anharmonic == 1)
160 : 1 : end subroutine load_from_files
161 : :
162 : : !-------------------------------------------------------------------!
163 : : ! load potential from netcdf file
164 : : ! Note that the lattic part of the primitive cell is also loaded.
165 : : ! Input:
166 : : ! fname: filename
167 : : !-------------------------------------------------------------------!
168 : 1 : subroutine load_from_netcdf(self, fname)
169 : : class(lwf_primitive_potential_t), intent(inout):: self
170 : : character(len = fnlen), intent(in):: fname
171 : : integer:: ncid, ierr
172 : : integer:: iR, nR, nlwf, natom, twobody_nterm, onebody_nterm
173 : : real(dp):: cell(3, 3)
174 : 1 : real(dp), allocatable:: ifc(:, :, :), xcart(:,:), masses(:)
175 : 1 : real(dp), allocatable :: twobody_val(:)
176 : 1 : integer, allocatable :: zion(:), twobody_iR(:), twobody_i(:), twobody_j(:), twobody_orderi(:), twobody_orderj(:)
177 : : integer:: varid, i, j
178 : : #if defined HAVE_NETCDF
179 : 1 : ierr = nf90_open(trim(fname), NF90_NOWRITE, ncid)
180 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "Open netcdf file")
181 : :
182 : : ! read primcell info
183 : 1 : ierr = nctk_get_dim(ncid, "wann_natom", natom)
184 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "getting natom in lwf potential file")
185 : :
186 : 1 : ierr = nctk_get_dim(ncid, "nR" , nR)
187 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "getting nR in lwf potential file")
188 : 1 : self%nR = nR
189 : :
190 : 1 : ierr = nctk_get_dim(ncid, "wann_nwann", nlwf)
191 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "getting wann_nwann in lwf potential file")
192 : :
193 : : ! TODO : add lattice
194 : : !call self%primcell%set_lattice()
195 : :
196 : 1 : call self%primcell%set_lwf(nlwf)
197 : 1 : self%nlwf = nlwf
198 [ + + ]: 4 : call self%coeff%initialize([ nR, nlwf, nlwf])
199 [ + + ]: 6 : call self%coeff_twobody%initialize([nR, nlwf, nlwf, -1, -1])
200 : :
201 [ - + + - : 5 : ABI_MALLOC(ifc, (nlwf, nlwf, nR))
- + + - -
+ + - + -
+ - - + -
+ ]
202 [ - + + - : 3 : ABI_MALLOC(self%Rlist, (3, nR))
+ - + - -
+ - + -
+ ]
203 [ - + + - : 3 : ABI_MALLOC(xcart, (3, natom))
+ - + - -
+ - + ]
204 [ - + + - : 3 : ABI_MALLOC(zion, (natom))
+ - + - -
+ - + ]
205 [ - + + - : 3 : ABI_MALLOC(masses, (natom))
+ - + - -
+ - + ]
206 : :
207 [ - + + - : 3 : ABI_MALLOC(self%lwf_masses, (nlwf))
+ - + - -
+ - + -
+ ]
208 : :
209 : 1 : ierr = nf90_inq_varid(ncid, "wann_lwf_masses", varid)
210 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "lwf_masses")
211 : 1 : ierr = nf90_get_var(ncid, varid, self%lwf_masses)
212 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "lwf_masses")
213 [ + + ]: 3 : self%lwf_masses(:) = self%lwf_masses(:) * amu_emass
214 [ + - - + : 4 : self%primcell%lwf%lwf_masses = self%lwf_masses
- - - - -
- + + ]
215 : :
216 : 1 : ierr = nf90_inq_dimid(ncid, "wann_onebody_nterm", onebody_nterm)
217 [ - + ]: 1 : if (ierr /= nf90_noerr) then
218 : 0 : onebody_nterm = 0
219 : : else
220 : 1 : ierr = nctk_get_dim(ncid, "wann_onebody_nterm", onebody_nterm)
221 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "getting wann_onebody_nterm in lwf potential file")
222 : : end if
223 : 1 : self%onebody_nterm = onebody_nterm
224 : :
225 : :
226 : 1 : ierr = nf90_inq_dimid(ncid, "wann_twobody_nterm", twobody_nterm)
227 [ + - ]: 1 : if (ierr /= nf90_noerr) then
228 : 1 : twobody_nterm = 0
229 : : else
230 : 0 : ierr = nctk_get_dim(ncid, "wann_twobody_nterm", twobody_nterm)
231 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "getting wann_twobody_nterm in lwf potential file")
232 : : end if
233 : :
234 : 1 : self%natom = natom
235 : :
236 : :
237 : 1 : ierr = nf90_inq_varid(ncid, "wann_cell", varid)
238 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_cell")
239 : 1 : ierr = nf90_get_var(ncid, varid, cell)
240 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_cell")
241 [ + + + + ]: 13 : cell(:,:)=cell(:,:)/ Bohr_Ang
242 : :
243 : :
244 : 1 : ierr = nf90_inq_varid(ncid, "wann_atomic_xcart", varid)
245 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_atomic_xcart")
246 : 1 : ierr = nf90_get_var(ncid, varid, xcart)
247 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_atomic_xcart")
248 [ + + + + ]: 25 : xcart(:,:)=xcart(:,:)/ Bohr_Ang
249 : :
250 : 1 : ierr = nf90_inq_varid(ncid, "wann_atomic_numbers", varid)
251 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_atomic_numbers")
252 : 1 : ierr = nf90_get_var(ncid, varid, zion)
253 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_atomic_numbers")
254 : :
255 : 1 : ierr = nf90_inq_varid(ncid, "wann_atomic_masses", varid)
256 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_atomic_masses")
257 : 1 : ierr = nf90_get_var(ncid, varid, masses)
258 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_atomic_masses")
259 : :
260 : 1 : call self%primcell%set_lattice(natom, cell, xcart, masses, zion)
261 : :
262 : :
263 [ - + + - : 5 : ABI_MALLOC(self%lattice_coeffs, (nlwf, natom*3, nR))
+ - - + +
- - + + -
+ - + - -
+ - + -
+ ]
264 : :
265 : 1 : ierr = nf90_inq_varid(ncid, "wann_wannier_function_real", varid)
266 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_wannier_function_real")
267 : 1 : ierr = nf90_get_var(ncid, varid, self%lattice_coeffs)
268 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_wannier_function_real")
269 : :
270 : 1 : ierr = nf90_inq_varid(ncid, "wann_Rlist", varid)
271 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_Rlist")
272 : 1 : ierr = nf90_get_var(ncid, varid, self%Rlist)
273 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_Rlist")
274 : :
275 : 1 : ierr = nf90_inq_varid(ncid, "wann_HamR_real", varid)
276 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_HamR_real")
277 : 1 : ierr = nf90_get_var(ncid, varid, ifc)
278 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_HamR_real")
279 : :
280 : :
281 [ + + + + : 316 : ifc(:,:,:) = ifc(:,:,:) * eV_Ha * (Bohr_Ang*Bohr_Ang)
+ + ]
282 : :
283 [ + + ]: 46 : do iR = 1, nR
284 [ + + ]: 136 : do i = 1, nlwf
285 [ + + ]: 315 : do j = 1, nlwf
286 [ + + ]: 270 : if (abs(ifc(j, i, iR))>1e-4) then
287 : : ! NOTE: in fortran the order of index in reversed when reading netcdf array.
288 [ + + ]: 528 : call self%coeff%add_entry([iR, i, j], ifc(j, i, iR))
289 : : end if
290 : : end do
291 : : end do
292 : : end do
293 : :
294 : :
295 [ + - ]: 1 : if(self%onebody_nterm /= 0) then
296 [ - + + - : 3 : ABI_MALLOC(self%onebody_i, (self%onebody_nterm))
+ - + - -
+ - + -
+ ]
297 [ - + + - : 2 : ABI_MALLOC(self%onebody_order, (self%onebody_nterm))
- + - + -
+ ]
298 [ - + + - : 3 : ABI_MALLOC(self%onebody_val, (self%onebody_nterm))
+ - + - -
+ - + -
+ ]
299 : :
300 : 1 : ierr = nf90_inq_varid(ncid, "wann_onebody_i", varid)
301 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_onebody_i")
302 : 1 : ierr = nf90_get_var(ncid, varid, self%onebody_i)
303 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_onebody_i")
304 : :
305 [ + + ]: 7 : self%onebody_i = self%onebody_i+1
306 : :
307 : 1 : ierr = nf90_inq_varid(ncid, "wann_onebody_val", varid)
308 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_onebody_val")
309 : 1 : ierr = nf90_get_var(ncid, varid, self%onebody_val)
310 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_onebody_val")
311 : :
312 : 1 : ierr = nf90_inq_varid(ncid, "wann_onebody_order", varid)
313 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_onebody_order")
314 : 1 : ierr = nf90_get_var(ncid, varid, self%onebody_order)
315 [ - + ]: 1 : NCF_CHECK_MSG(ierr, "wann_onebody_order")
316 : :
317 [ + + ]: 7 : do i = 1, self%onebody_nterm
318 : 7 : self%onebody_val(i)=self%onebody_val(i) * eV_Ha * (Bohr_Ang**self%onebody_order(i))
319 : : end do
320 : : end if
321 : :
322 [ - + ]: 1 : if(twobody_nterm /= 0) then
323 [ # # # # : 0 : ABI_MALLOC(twobody_iR, (twobody_nterm))
# # # # #
# # # ]
324 [ # # # # : 0 : ABI_MALLOC(twobody_i, (twobody_nterm))
# # # # ]
325 [ # # # # : 0 : ABI_MALLOC(twobody_j, (twobody_nterm))
# # # # ]
326 [ # # # # : 0 : ABI_MALLOC(twobody_orderi, (twobody_nterm))
# # # # ]
327 [ # # # # : 0 : ABI_MALLOC(twobody_orderj, (twobody_nterm))
# # # # ]
328 [ # # # # : 0 : ABI_MALLOC(twobody_val, (twobody_nterm))
# # # # #
# # # ]
329 : :
330 : 0 : ierr = nf90_inq_varid(ncid, "wann_twobody_iR", varid)
331 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_iR")
332 : 0 : ierr = nf90_get_var(ncid, varid, twobody_iR)
333 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_iR")
334 : :
335 : 0 : ierr = nf90_inq_varid(ncid, "wann_twobody_i", varid)
336 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_i")
337 : 0 : ierr = nf90_get_var(ncid, varid, twobody_i)
338 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_i")
339 : :
340 : 0 : ierr = nf90_inq_varid(ncid, "wann_twobody_j", varid)
341 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_j")
342 : 0 : ierr = nf90_get_var(ncid, varid, twobody_j)
343 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_j")
344 : :
345 : 0 : ierr = nf90_inq_varid(ncid, "wann_twobody_orderi", varid)
346 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_orderi")
347 : 0 : ierr = nf90_get_var(ncid, varid, twobody_orderi)
348 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_i")
349 : :
350 : 0 : ierr = nf90_inq_varid(ncid, "wann_twobody_orderj", varid)
351 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_orderj")
352 : 0 : ierr = nf90_get_var(ncid, varid, twobody_orderj)
353 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_orderj")
354 : :
355 : 0 : ierr = nf90_inq_varid(ncid, "wann_twobody_val", varid)
356 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_val")
357 : 0 : ierr = nf90_get_var(ncid, varid, twobody_val)
358 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "wann_twobody_val")
359 : :
360 : 0 : ierr = nf90_close(ncid)
361 [ # # ]: 0 : NCF_CHECK_MSG(ierr, "Close netcdf file")
362 : :
363 [ # # ]: 0 : do i = 1, twobody_nterm
364 : : twobody_val(i)=twobody_val(i) * eV_Ha * &
365 : 0 : & (Bohr_Ang**(twobody_orderi(i)+ twobody_orderj(i)))
366 : : call self%coeff_twobody%add_entry([twobody_iR(i), twobody_i(i), &
367 [ # # ]: 0 : & twobody_j(i), twobody_orderi(i), twobody_orderj(i)], twobody_val(i))
368 : : end do
369 : : end if
370 : :
371 [ + - ]: 1 : ABI_SFREE(xcart)
372 [ + - ]: 1 : ABI_SFREE(zion)
373 [ + - ]: 1 : ABI_SFREE(masses)
374 [ + - ]: 1 : ABI_SFREE(ifc)
375 : :
376 [ - + ]: 1 : if(twobody_nterm /= 0) then
377 [ # # ]: 0 : ABI_SFREE(twobody_i)
378 [ # # ]: 0 : ABI_SFREE(twobody_j)
379 [ # # ]: 0 : ABI_SFREE(twobody_orderi)
380 [ # # ]: 0 : ABI_SFREE(twobody_orderj)
381 [ # # ]: 0 : ABI_SFREE(twobody_val)
382 : : endif
383 : : #else
384 : : NETCDF_NOTENABLED_ERROR()
385 : : #endif
386 : :
387 [ - + - + : 1 : end subroutine load_from_netcdf
- + - + -
+ - + - +
- + - + -
+ ]
388 : :
389 : : !-------------------------------------------------------------------!
390 : : !Fill supercell
391 : : ! Inputs:
392 : : ! scmaker : supercell_maker_t
393 : : ! Output:
394 : : ! scpot: a class pointer to an ABSTRACT potential.
395 : : ! Not a potential because this is inherited from
396 : : ! an abstract_primitive_potential_t, which doesn't know
397 : : ! the type of the supercell potential.
398 : : !
399 : : !-------------------------------------------------------------------!
400 : 1 : subroutine fill_supercell(self, scmaker, params, scpot, supercell)
401 : : use m_spmat_convert, only: COO_to_dense
402 : :
403 : : class(lwf_primitive_potential_t), intent(inout):: self
404 : : type(supercell_maker_t), intent(inout):: scmaker
405 : : type(multibinit_dtset_type), intent(inout):: params
406 : : class(abstract_potential_t), pointer, intent(inout):: scpot
407 : : type(mbsupercell_t), target:: supercell
408 : :
409 : : integer:: nlwf, sc_nlwf
410 : : integer:: inz, iR, R(3), i, j, icell
411 : 1 : integer, allocatable:: ilist_sc(:), jlist_sc(:), Rlist_sc(:,:)
412 : : real(dp):: val
413 : :
414 : 1 : ABI_UNUSED_A(params)
415 : :
416 : :
417 : 1 : nlwf = self%nlwf
418 : 1 : sc_nlwf = nlwf*scmaker%ncells
419 : :
420 : : ! Step 1: allocate the scpot as a corresponding supercell potential
421 [ - + ]: 1 : ABI_MALLOC_TYPE_SCALAR(lwf_potential_t, scpot)
422 : : ! Fortran does not know the functions specific to the derived class pointer.
423 : : ! Only the ones inheritated from abstract class,
424 : : ! unless select type is used:
425 : : select type(scpot)
426 : : type is (lwf_potential_t)
427 : 1 : call scpot%initialize(sc_nlwf)
428 : 1 : call scpot%set_supercell(supercell)
429 : : ! IFC is an COO_mat_t, which has the index of R1, R2, R3, i, j and the value of val
430 : : ! list of index R: coeff%ind%data(1, 1:coeff%nnz)
431 : : ! list of i: coeff%ind%data(2, 1:coeff%nnz)
432 : : ! list of j: coeff%ind%data(3, 1:coeff%nnz)
433 : : ! IFC: (ind) = val
434 [ - + ]: 1 : if(self%as_lattice_anharmonic) then
435 : 0 : call scpot%use_as_lattice_anharmonic()
436 : : else
437 : : ! harmonic terms
438 [ + + ]: 133 : do inz = 1, self%coeff%nnz
439 : : ! For each non-zero entry in the coeff matrix
440 : : ! get the R, i, and j, val
441 : 132 : iR = self%coeff%ind%data(1, inz)
442 [ + + ]: 528 : R(:) = self%Rlist(:, iR)
443 : 132 : i = self%coeff%ind%data(2, inz)
444 : 132 : j = self%coeff%ind%data(3, inz)
445 : 132 : val = self%coeff%val%data(inz)
446 : : ! translate i to i in supercell.
447 : : ! No need to allocate, it is done by trans_i . but remember to deallocate!
448 : : ! nbasis is the number in one primitive cell.
449 : : ! e.g. there are 3*natom possible i (3: x, y, z) in each primitive cell.
450 : 132 : call scmaker%trans_i(nbasis = self%nlwf, i = i, i_sc = ilist_sc )
451 : : ! translate j, Rj to supercell.
452 : 132 : call scmaker%trans_j_and_Rj(nbasis = self%nlwf, j = j, Rj = R, j_sc = jlist_sc, Rj_sc = Rlist_sc)
453 : : ! values are repeated in cells
454 [ + + ]: 67716 : do icell = 1, scmaker%ncells
455 : 67716 : call scpot%add_term(ilist_sc(icell), jlist_sc(icell), val )
456 : : end do
457 [ + - ]: 132 : ABI_SFREE(ilist_sc)
458 [ + - ]: 132 : ABI_SFREE(jlist_sc)
459 [ + - ]: 133 : ABI_SFREE(Rlist_sc)
460 : : end do
461 : : end if
462 : :
463 : : ! coefficients of atomic displacements in supercell
464 : 1 : call scpot%supercell%lwf%lwf_latt_coeffs%initialize(self%nlwf*scmaker%ncells, self%natom*3*scmaker%ncells )
465 [ + + ]: 3 : do i = 1, self%nlwf
466 : 2 : call scmaker%trans_i(nbasis = nlwf, i = i, i_sc = ilist_sc)
467 : :
468 [ + + ]: 92 : do iR = 1, self%nR
469 [ + + ]: 360 : R(:) = self%Rlist(:, iR)
470 [ + + ]: 1712 : do j = 1, self%natom*3
471 : 1620 : val = self%lattice_coeffs(i, j, iR) ! ilwf, iatom3, iR
472 [ + + ]: 1710 : if (abs(val) > 1e-4) then
473 : 1348 : call scmaker%trans_j_and_Rj(nbasis = self%natom*3, j = j, Rj = R, j_sc = jlist_sc, Rj_sc = Rlist_sc)
474 [ + + ]: 691524 : do icell = 1, scmaker%ncells
475 : : !call scpot%supercell%lwf%lwf_latt_coeffs(ilist_sc(icell))%push(jlist_sc(icell), val)
476 [ + + ]: 2071876 : call scpot%supercell%lwf%lwf_latt_coeffs%coeffs%add_entry([jlist_sc(icell), ilist_sc(icell)], val)
477 : : end do
478 [ + - ]: 1348 : ABI_SFREE(jlist_sc)
479 [ + - ]: 1348 : ABI_SFREE(Rlist_sc)
480 : : end if
481 : : end do
482 : : end do
483 [ + - ]: 3 : ABI_SFREE(ilist_sc)
484 : : end do
485 : :
486 : : ! anharmonic terms
487 : : !call scpot%set_ref_energy(self%ref_energy*scmaker%ncells)
488 : : !if (self%has_self_bound_term) then
489 : : ! call scpot%add_self_bound_term(self%self_bound_order, self%self_bound_coeff)
490 : : !end if
491 : :
492 [ + - ]: 2 : if (self%onebody_nterm /= 0) then
493 [ + + ]: 7 : do i = 1, self%onebody_nterm
494 : 6 : call scmaker%trans_i(nbasis = self%nlwf, i = self%onebody_i(i), i_sc = ilist_sc)
495 [ + + ]: 3078 : do icell = 1, scmaker%ncells
496 : 3078 : call scpot%add_onebody_term(i = ilist_sc(icell), order = self%onebody_order(i), val = self%onebody_val(i))
497 : : end do
498 [ + - ]: 7 : ABI_SFREE(ilist_sc)
499 : : end do
500 : : endif
501 : : ! Test the phonon energy
502 : : !call COO_to_dense(scpot%coeff, real_sc_evecs)
503 : : !sc_evecs(:,:) = real_sc_evecs(:,:)
504 : : !call eigensh(sc_evals, sc_evecs)
505 : :
506 : : end select
507 : :
508 : :
509 [ - + - + : 1 : end subroutine fill_supercell
- + ]
510 : :
511 : : !-------------------------------------------------------------------!
512 : : ! calculate hamiltonian at k point (or more precisely, q point)
513 : : ! H_ij(k) = sum_R H_ij(R) exp( i2pi k.dot.R)
514 : : !-------------------------------------------------------------------!
515 [ # # ]: 0 : subroutine get_hamk(self, kpoint, hamk)
516 : : class(lwf_primitive_potential_t), intent(in):: self
517 : : real(dp), intent(in):: kpoint(3)
518 : : complex(dp), intent(inout):: hamk(:,:)
519 : : integer:: inz, iR, R(3), i, j
520 : : real(dp):: val
521 [ # # # # ]: 0 : hamk(:,:) = cmplx(0.0, 0.0)
522 [ # # ]: 0 : do inz = 1, self%coeff%nnz
523 : 0 : iR = self%coeff%ind%data(1, inz)
524 [ # # ]: 0 : R(:) = self%Rlist(:, iR)
525 : 0 : i = self%coeff%ind%data(2, inz)
526 : 0 : j = self%coeff%ind%data(3, inz)
527 : 0 : val = self%coeff%val%data(inz)
528 [ # # ]: 0 : hamk(i, j) =hamk(i, j) + val*exp(cmplx(0.0, 2.0) *pi*dot_product(kpoint, R))
529 : : end do
530 : 0 : end subroutine get_hamk
531 : :
532 : : !-------------------------------------------------------------------!
533 : : !calculate eigenvalue and eigen vector
534 : : !-------------------------------------------------------------------!
535 [ # # # # ]: 0 : subroutine get_eigen(self, kpoint, evals, evecs)
536 : : class(lwf_primitive_potential_t), intent(in):: self
537 : : real(dp), intent(in):: kpoint(3)
538 : : real(dp), intent(inout):: evals(:)
539 : : complex(dp), intent(inout):: evecs(:,:)
540 : 0 : call self%get_hamk(kpoint, evecs)
541 : : ! The evecs array is reused both as the matrix and eigenvectors.
542 : 0 : call eigensh(evals, evecs)
543 : 0 : end subroutine get_eigen
544 : :
545 : : !subroutine add_self_bound_term(self, order, coeff)
546 : : ! class(lwf_primitive_potential_t), intent(inout):: self
547 : : ! integer, intent(in):: order
548 : : ! real(dp), intent(in):: coeff
549 : : ! if (order /= 0) then
550 : : ! self%has_self_bound_term=.True.
551 : : ! self%self_bound_order = order
552 : : ! self%self_bound_coeff = coeff
553 : : ! end if
554 : :
555 : : !end subroutine add_self_bound_term
556 : :
557 : :
558 [ + + - + : 4 : end module m_lwf_primitive_potential
+ + - + -
+ - + - +
- + + + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
|