Branch data Line data Source code
1 : : !!****m* ABINIT/m_tdep_dataset
2 : : !! NAME
3 : : !! m_tdep_dataset
4 : : !!
5 : : !! FUNCTION
6 : : !! This module contains the data type holding the input variables of atdep
7 : : !! (atdep_dataset_type) and the routines reading and echoing them, as well
8 : : !! as the data type describing the MPI parallelism (MPI_enreg_type) and
9 : : !! the distribution of the MD steps over the processors.
10 : : !!
11 : : !! COPYRIGHT
12 : : !! Copyright (C) 2011-2026 ABINIT group (GA,FB,JB)
13 : : !! This file is distributed under the terms of the
14 : : !! GNU General Public License, see ~abinit/COPYING
15 : : !! or http://www.gnu.org/copyleft/gpl.txt .
16 : : !!
17 : : !! SOURCE
18 : :
19 : : #if defined HAVE_CONFIG_H
20 : : #include "config.h"
21 : : #endif
22 : :
23 : : #include "abi_common.h"
24 : :
25 : : module m_tdep_dataset
26 : :
27 : : use defs_basis
28 : : use m_errors
29 : : use m_abicore
30 : : use m_xmpi
31 : : use m_abihist
32 : : use m_parser
33 : : use netcdf
34 : : use m_nctk
35 : : use m_fstrings, only : inupper,ljust,next_token
36 : : use m_io_tools, only : file_exists
37 : : use m_ddb, only : ddb_type
38 : : use m_ddb_hdr, only : ddb_hdr_type
39 : : use m_crystal, only : crystal_t
40 : :
41 : : implicit none
42 : : !!***
43 : :
44 : : !----------------------------------------------------------------------
45 : :
46 : : !!****t* m_tdep_dataset/atdep_dataset_type
47 : : !! NAME
48 : : !! atdep_dataset_type
49 : : !!
50 : : !! FUNCTION
51 : : !! Holds the input variables of atdep (read from the input file, the
52 : : !! HIST.nc file and the DDB file), as well as a few derived quantities.
53 : : !!
54 : : !! SOURCE
55 : :
56 : : type atdep_dataset_type
57 : :
58 : : integer :: natom
59 : : ! Number of atoms in the supercell
60 : :
61 : : integer :: natom_unitcell
62 : : ! Number of atoms in the unit cell
63 : :
64 : : integer :: nshell_max
65 : : ! Maximum number of shells allowed (hard-coded)
66 : :
67 : : integer :: nstep_max
68 : : ! Last MD step used
69 : :
70 : : integer :: nstep_min
71 : : ! First MD step used
72 : :
73 : : integer :: nstep_tot
74 : : ! Total number of MD steps used: (nstep_max-nstep_min)/slice+1
75 : :
76 : : integer :: my_nstep
77 : : ! Number of MD steps held by the present processor (see tdep_init_MPIdata)
78 : :
79 : : integer :: ntypat
80 : : ! Number of atom types
81 : :
82 : : integer :: use_ideal_positions
83 : : ! Compute the displacements with respect to the average positions (0)
84 : : ! or to the ideal positions (1)
85 : :
86 : : integer :: ab_out
87 : : ! Unit number of the main output file
88 : :
89 : : integer :: std_out
90 : : ! Unit number of the standard output
91 : :
92 : : integer :: bzpath
93 : : ! Definition of the q-point path for the phonon dispersion:
94 : : ! > 0: number of special q-points given by their names (special_qpt)
95 : : ! < 0: number of q-points given by their reduced coordinates (qpt)
96 : :
97 : : integer :: ndivsm
98 : : ! Number of divisions of the smallest segment of the q-point path
99 : :
100 : : integer :: order
101 : : ! Maximum IFC order (2, 3, or 4)
102 : :
103 : : integer :: slice
104 : : ! Only one MD step every slice steps is used
105 : :
106 : : integer :: enunit
107 : : ! Unit of the phonon frequencies: 0=meV, 1=cm-1, 2=mHa, 3=THz
108 : :
109 : : integer :: readifc
110 : : ! Read the IFC from file: 0=no, 1=from _ifc_in.dat, 2=from _ifc_out.dat
111 : :
112 : : integer :: together
113 : : ! Compute the IFC of all orders simultaneously (1) or successively (0)
114 : :
115 : : integer :: alloy
116 : : ! Use the virtual crystal approximation (1) to merge the atom types
117 : : ! ityp_alloy1 and ityp_alloy2
118 : :
119 : : integer :: ityp_alloy1
120 : : ! First atom type merged when alloy=1
121 : :
122 : : integer :: ityp_alloy2
123 : : ! Second atom type merged when alloy=1
124 : :
125 : : integer :: nproc(2)
126 : : ! Number of processors for the distribution over shells (1) and over MD steps (2)
127 : :
128 : : integer :: bzlength
129 : : ! Number of segments of the q-point path whose lengths are given
130 : : ! explicitly (lgth_segments)
131 : :
132 : : integer :: ngqpt1(3)
133 : : ! q-point grid used to build the IFC and DDB objects (zero if not given in input)
134 : :
135 : : integer :: qptrlatt(3,3)
136 : : ! q-mesh lattice (as kptrlatt for k-points) of the coarse q-point grid used to build
137 : : ! the IFC and DDB objects: diag(ngqpt1) if ngqpt1 is given in input, otherwise the
138 : : ! supercell given by multiplicity (whose lines are the columns of qptrlatt).
139 : :
140 : : integer :: ngqpt2(3)
141 : : ! q-point grid used to compute the phonon density of states
142 : :
143 : : double precision :: q1shft(3)
144 : : ! Shift of the coarse q-point grid (qptrlatt), as anaddb's q1shft (with nqshft=1)
145 : :
146 : : double precision :: q2shft(3)
147 : : ! Shift of the q-point grid ngqpt2, as anaddb's q2shft
148 : :
149 : : integer :: use_weights
150 : : ! Read the weights of the MD configurations from the _weights.dat file (1),
151 : : ! or use the same weight for each of them (0)
152 : :
153 : : integer :: prtphdos
154 : : ! Print the phonon density of states (1)
155 : :
156 : : integer :: prtvol ! Verbosity level (see the ABINIT input variable of the same name)
157 : : integer :: symbec ! Symmetrize born effective charges
158 : :
159 : : integer, allocatable :: typat_unitcell(:)
160 : : ! typat_unitcell(natom_unitcell)
161 : : ! Atom type in the unitcell
162 : :
163 : : integer, allocatable :: typat(:)
164 : : ! typat(natom)
165 : : ! Atom type in the supercell
166 : :
167 : : integer, allocatable :: lgth_segments(:)
168 : : ! lgth_segments(bzlength)
169 : : ! Number of q-points along each segment of the q-point path
170 : :
171 : : logical :: dipdip
172 : : ! Whether the dipole-dipole (long-range) interaction is subtracted from
173 : : ! the dynamical matrix, a la [[anaddb:dipdip]] (same physical meaning,
174 : : ! but atdep only ever supports the on/off choice: unlike anaddb, there is
175 : : ! no dipdip=-1 fast-algorithm variant, and it is auto-enabled whenever
176 : : ! usable Born effective charges/dielectric tensor are found, see
177 : : ! tdep_read_input). Explicitly setting the input variable dipdip=0
178 : : ! forces it off even when such data is available.
179 : :
180 : : logical :: netcdf
181 : : ! Whether the MD trajectory is read from the HIST.nc netcdf file
182 : :
183 : : double precision :: dosdeltae
184 : : ! Energy step of the phonon density of states (Ha)
185 : :
186 : : double precision :: rcut
187 : : ! Cutoff radius of the second-order IFC (bohr)
188 : :
189 : : double precision :: rcut3
190 : : ! Cutoff radius of the third-order IFC (bohr)
191 : :
192 : : double precision :: rcut4
193 : : ! Cutoff radius of the fourth-order IFC (bohr)
194 : :
195 : : double precision :: temperature
196 : : ! Temperature of the MD simulation (K)
197 : :
198 : : double precision :: tolread
199 : : ! Tolerance given together with readifc=1
200 : :
201 : : double precision :: rprimd_md(3,3)
202 : : ! Dimensioned primitive vectors
203 : :
204 : : double precision :: multiplicity(3,3)
205 : : ! Supercell scaling matrix (integers)
206 : :
207 : : double precision :: dielt(3,3)
208 : : ! Dielectric tensor
209 : :
210 : : double precision, allocatable :: amu(:)
211 : : ! amu(ntypat)
212 : : ! Mass of the atoms (atomic mass unit)
213 : :
214 : : double precision, allocatable :: zeff(:,:,:)
215 : : ! zeff(3,3,natom_unitcell)
216 : : ! Born effective charge tensor
217 : :
218 : : double precision, allocatable :: znucl(:)
219 : : ! znucl(ntypat)
220 : : ! Nuclear charge of each atom type
221 : :
222 : : double precision, allocatable :: qpt(:,:)
223 : : ! qpt(3,abs(bzpath))
224 : : ! Reduced coordinates of the q-points defining the path (bzpath < 0)
225 : :
226 : : double precision, allocatable :: xred_unitcell(:,:)
227 : : ! xred_unitcell(3, natom_unitcell)
228 : : ! Reduced equilibrium position of the atoms in the unitcell, anchored so
229 : : ! that atom 1 sits at the origin (0,0,0) and wrapped into [0,1).
230 : :
231 : : double precision, allocatable :: xred_unitcell_original(:,:)
232 : : ! xred_unitcell_original(3, natom_unitcell)
233 : : ! Reduced equilibrium position of the atoms in the unitcell,
234 : : ! as given in the input file, wrapped into (-0.5,0.5].
235 : :
236 : : character (len=2), allocatable :: special_qpt(:)
237 : : ! special_qpt(bzpath)
238 : : ! Names of the special q-points defining the path (bzpath > 0)
239 : :
240 : : character (len=fnlen) :: output_prefix
241 : : ! Root name of the output files
242 : :
243 : : character (len=fnlen) :: input_prefix
244 : : ! Root name of the input files (HIST.nc, DDB, xred.dat, ...)
245 : :
246 : : character (len=fnlen) :: output_file
247 : : ! Name of the main output file
248 : :
249 : : contains
250 : :
251 : : procedure :: print => tdep_write_input
252 : : procedure :: free => tdep_destroy_invar
253 : :
254 : : end type atdep_dataset_type
255 : : !!***
256 : :
257 : : !----------------------------------------------------------------------
258 : :
259 : : !!****t* m_tdep_dataset/MPI_enreg_type
260 : : !! NAME
261 : : !! MPI_enreg_type
262 : : !!
263 : : !! FUNCTION
264 : : !! Holds the information on the MPI parallelism: the 2D grid of processors
265 : : !! (shells x MD steps), the corresponding communicators, and the
266 : : !! distribution of the MD steps over the processors.
267 : : !!
268 : : !! SOURCE
269 : :
270 : : type MPI_enreg_type
271 : :
272 : : integer :: comm_shell
273 : : ! Communicator for the distribution over shells
274 : :
275 : : integer :: comm_step
276 : : ! Communicator for the distribution over MD steps
277 : :
278 : : integer :: comm_shellstep
279 : : ! Communicator over the whole 2D grid of processors (shells x steps)
280 : :
281 : : integer :: nproc
282 : : ! Total number of processors
283 : :
284 : : integer :: nproc_shell
285 : : ! Number of processors for the distribution over shells
286 : :
287 : : integer :: nproc_step
288 : : ! Number of processors for the distribution over MD steps
289 : :
290 : : integer :: master
291 : : ! Rank of the master processor
292 : :
293 : : integer :: me_shell
294 : : ! Rank of the present processor in comm_shell
295 : :
296 : : integer :: me_step
297 : : ! Rank of the present processor in comm_step
298 : :
299 : : integer, allocatable :: my_nshell(:)
300 : : ! Currently unused
301 : :
302 : : integer :: my_nstep
303 : : ! Number of MD steps held by the present processor
304 : :
305 : : logical :: iam_master
306 : : ! Whether the present processor is the master
307 : :
308 : : integer, allocatable :: nstep_acc(:)
309 : : ! nstep_acc(nproc_step+1)
310 : : ! Accumulated number of MD steps held by the previous processors
311 : :
312 : : integer, allocatable :: nstep_all(:)
313 : : ! nstep_all(nproc_step)
314 : : ! Number of MD steps held by each processor
315 : :
316 : : integer, allocatable :: shft_step(:)
317 : : ! shft_step(nproc_step)
318 : : ! Offset of the MD steps held by each processor (used in xmpi_gatherv)
319 : :
320 : : logical, allocatable :: my_shell(:)
321 : : ! Currently unused
322 : :
323 : : logical, allocatable :: my_step(:)
324 : : ! my_step(nstep_tot)
325 : : ! Whether each MD step is held by the present processor
326 : :
327 : : contains
328 : :
329 : : procedure :: init => tdep_init_MPIdata
330 : : procedure :: free => tdep_destroy_mpidata
331 : :
332 : : end type MPI_enreg_type
333 : : !!***
334 : :
335 : : public :: tdep_read_input
336 : : public :: tdep_write_input
337 : : public :: tdep_destroy_invar
338 : : public :: tdep_init_MPIdata
339 : : public :: tdep_destroy_mpidata
340 : :
341 : : contains
342 : :
343 : : !====================================================================================================
344 : :
345 : : !!****f* ABINIT/m_tdep_dataset/tdep_read_input
346 : : !! NAME
347 : : !! tdep_read_input
348 : : !!
349 : : !! FUNCTION
350 : : !! Read the input variables of atdep and store them in the Invar object.
351 : : !! The names of the input file, of the input/output prefixes and of the
352 : : !! main output file are either given by input_path (and the input
353 : : !! variables indata_prefix, outdata_prefix and output_file), or read from
354 : : !! the standard input (deprecated "files" file). The main output file is
355 : : !! then opened.
356 : : !! If the HIST.nc netcdf file is available, the MD trajectory is read
357 : : !! into Hist, together with natom, ntypat, typat, znucl and amu.
358 : : !! If the DDB file is available, the Born effective charges and the
359 : : !! dielectric tensor are read from it (and dipdip is enabled).
360 : : !! Then, the mandatory and optional input variables are read and checked,
361 : : !! the virtual crystal approximation is applied (if alloy=1), and the
362 : : !! unitcell is anchored so that atom 1 sits at the origin. Finally,
363 : : !! the total number of MD steps (nstep_tot) is computed.
364 : : !!
365 : : !! INPUTS
366 : : !! input_path = Path of the input file. If empty, the names of the files
367 : : !! are read from the standard input.
368 : : !!
369 : : !! OUTPUT
370 : : !! Hist = abihist object containing the MD trajectory
371 : : !! (only if the HIST.nc file is available).
372 : : !! Invar = Input object containing the input variables.
373 : : !!
374 : : !! SOURCE
375 : :
376 [ - + - + : 46 : subroutine tdep_read_input(input_path,Hist,Invar)
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ ]
377 : :
378 : : ! Arguments-------------------------------
379 : : character(len=*), intent(in):: input_path
380 : : type(atdep_dataset_type),intent(out) :: Invar
381 : : type(abihist), intent(out) :: Hist
382 : :
383 : : ! Local variables-------------------------
384 : : ! scalars
385 : : integer :: ncid, ncerr, me,ierr,master
386 : : integer :: nimage, mdtime, natom_id,nimage_id,time_id,xyz_id,six_id
387 : : integer :: ntypat_id,iatcell
388 : : integer :: ii,jj,shift,iatom,itypat,sum_alloy1,sum_alloy2
389 : : integer:: lenstr, marr, jdtset, tread
390 : : logical :: has_nimage
391 : : double precision :: dtion,amu_average,born_average
392 : : character(len=500) :: msg
393 : : character(len=fnlen) :: input_filename,hist_filename,ddb_filename
394 : : character(len=strlen):: string, raw_string
395 [ + + + + : 2898 : type(crystal_t):: crystal
+ + + + +
+ + + + +
+ + + + +
+ ]
396 : 46 : type(ddb_type):: ddb
397 [ + + + + : 2898 : type(ddb_hdr_type):: ddb_hdr
+ + + + +
+ + + + +
+ + + + +
+ ]
398 : : ! arrays
399 : 46 : integer, allocatable:: intarr(:)
400 : 46 : integer, allocatable :: typat_unitcell_tmp(:)
401 : : real(dp) :: zeff_average(3,3)
402 : 46 : real(dp), allocatable:: zeff_tmp(:,:,:)
403 : 92 : real(dp), allocatable :: xred_unitcell_tmp(:,:),amu_tmp(:),born_charge_tmp(:),znucl_tmp(:)
404 [ - + ]: 46 : real(dp), allocatable:: dprarr(:)
405 : :
406 : : ! *********************************************************************
407 : :
408 : : ! Define output files
409 : 46 : Invar%ab_out=ab_out
410 : 46 : Invar%std_out=std_out
411 : :
412 : : ! Define default values
413 : 46 : Invar%nshell_max=500
414 : 46 : Invar%bzpath=0
415 : 46 : Invar%ndivsm=20
416 : 46 : Invar%order=2
417 : 46 : Invar%slice=1
418 : 46 : Invar%enunit=0
419 : 46 : Invar%together=1
420 : 46 : Invar%alloy=0
421 : 46 : Invar%ityp_alloy1=0
422 : 46 : Invar%ityp_alloy2=0
423 [ + + ]: 138 : Invar%nproc(:)=1
424 : 46 : Invar%bzlength=0
425 : 46 : Invar%readifc=0
426 : 46 : Invar%tolread=1.d-8
427 : 46 : Invar%dosdeltae=0.2_dp/Ha_cmm1
428 : 46 : Invar%prtvol=0
429 : 46 : Invar%dipdip=.false.
430 : 46 : Invar%netcdf=.false.
431 : 46 : Invar%use_ideal_positions=0
432 : 46 : Invar%use_weights=0
433 : 46 : Invar%prtphdos=1
434 : 46 : Invar%symbec=0
435 : : ! By default, the coarse q-point grid is defined by the supercell (see qptrlatt)
436 [ + + ]: 184 : Invar%ngqpt1(:)=0
437 : : ! In order to have an accuracy better than 1meV
438 [ + + ]: 184 : Invar%ngqpt2(:)=32
439 : : ! Unshifted q-point grids (they contain Gamma)
440 [ + + ]: 184 : Invar%q1shft(:)=zero
441 [ + + ]: 184 : Invar%q2shft(:)=zero
442 [ + + + + ]: 598 : Invar%dielt(:,:)=zero; Invar%dielt(1,1)=one; Invar%dielt(2,2)=one; Invar%dielt(3,3)=one
443 : :
444 : 46 : master = 0
445 : 46 : me = xmpi_comm_rank(xmpi_world)
446 [ + - ]: 46 : if (me==master) then
447 : :
448 [ - + ]: 46 : if (len_trim(input_path) == 0) then
449 : :
450 : 0 : write(std_out, "(2a)")" DeprecationWarning: ",ch10
451 : 0 : write(std_out, "(a)") " The files file has been deprecated in Abinit10 and will be removed in Abinit11."
452 : 0 : write(std_out, "(2a)")" Use the syntax `atdep t01.abi` where t01.abi is an atdep input,",ch10
453 : 0 : write(std_out, "(2a)")" and use input variables output_file, indata_prefix, outdata_prefix.",ch10
454 : :
455 [ # # ]: 0 : write(std_out,'(a)',err=10) ' Give name for input file '
456 [ # # ]: 0 : read(*, '(a)',err=10) input_filename
457 [ # # ]: 0 : if ( input_filename == "" ) input_filename='input.in'
458 [ # # ]: 0 : write(std_out, '(a)',err=10) '.'//trim(input_filename)
459 : : 10 continue
460 : : ! Check if a NetCDF file is available
461 [ # # ]: 0 : write(std_out,'(a)',err=11) ' Give root name for generic input files (NetCDF or ASCII)'
462 [ # # ]: 0 : read(*, '(a)',err=11) Invar%input_prefix
463 [ # # ]: 0 : write(std_out, '(a)',err=11) '.'//trim(Invar%input_prefix)
464 : : 11 continue
465 [ # # ]: 0 : write(std_out,'(a)', err=12)' Give root name for generic output files:'
466 [ # # ]: 0 : read (*, '(a)', err=12) Invar%output_prefix
467 [ # # ]: 0 : if ( Invar%output_prefix == "" ) Invar%output_prefix = 'atdep'
468 [ # # ]: 0 : write (std_out, '(a)', err=12 ) '.'//trim(Invar%output_prefix)
469 : : 12 continue
470 [ # # ]: 0 : Invar%output_file = trim(Invar%output_prefix)//'.abo'
471 : : else
472 [ - + ]: 46 : input_filename = input_path
473 : :
474 : : ! Read input
475 : 46 : string = repeat(" ", strlen)
476 : 46 : raw_string = repeat(" ", strlen)
477 : 46 : call instrng(input_filename, lenstr, 1, strlen, string, raw_string)
478 : : ! To make case-insensitive, map characters to upper case.
479 : 46 : call inupper(string(1:lenstr))
480 : :
481 : 46 : marr = 3
482 [ - + ]: 46 : ABI_MALLOC(intarr, (marr))
483 [ - + ]: 46 : ABI_MALLOC(dprarr, (marr))
484 : 46 : jdtset = 0
485 : :
486 : : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), &
487 : 46 : "indata_prefix", tread, 'KEY', key_value=Invar%input_prefix)
488 [ - + ]: 46 : if (tread == 0) then
489 : 0 : Invar%input_prefix = ''
490 : : end if
491 : 46 : write(std_out, "(2a)")"Root name for input files: ", trim(Invar%input_prefix)
492 : :
493 : : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), &
494 : 46 : "outdata_prefix", tread, 'KEY', key_value=Invar%output_prefix)
495 [ - + ]: 46 : if (tread == 0) then
496 : 0 : Invar%output_prefix = 'atdep'
497 : : end if
498 : 46 : write(std_out, "(2a)")"Root name for output files: ", trim(Invar%output_prefix)
499 : :
500 : : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), &
501 : 46 : "output_file", tread, 'KEY', key_value=Invar%output_file)
502 [ - + ]: 46 : if (tread == 0) then
503 [ # # ]: 0 : Invar%output_file = trim(Invar%output_prefix) // '.abo'
504 : : end if
505 : 46 : write(std_out, "(2a)")"Main output file: ", trim(Invar%output_file)
506 : :
507 [ - + ]: 46 : ABI_FREE(intarr)
508 [ - + ]: 46 : ABI_FREE(dprarr)
509 : :
510 : : end if
511 : :
512 : : ! Setup other filenames
513 [ - + ]: 46 : if ( Invar%input_prefix == "" ) then
514 : 0 : hist_filename='HIST.nc'
515 : 0 : ddb_filename='DDB'
516 : : else
517 [ + - ]: 46 : hist_filename=trim(Invar%input_prefix)//'_HIST.nc'
518 [ + - ]: 46 : ddb_filename=trim(Invar%input_prefix)//'_DDB'
519 : : end if
520 : :
521 : 46 : open(unit=ab_out,file=trim(Invar%output_file))
522 : :
523 : : end if !me
524 : :
525 : 46 : call xmpi_bcast(input_filename,master,xmpi_world,ierr)
526 : 46 : call xmpi_bcast(hist_filename,master,xmpi_world,ierr)
527 : 46 : call xmpi_bcast(ddb_filename,master,xmpi_world,ierr)
528 : 46 : call xmpi_bcast(Invar%output_prefix,master,xmpi_world,ierr)
529 : 46 : call xmpi_bcast(Invar%input_prefix,master,xmpi_world,ierr)
530 : :
531 : : #if defined HAVE_NETCDF
532 : : !Open netCDF file
533 : 46 : ncerr=nf90_open(path=trim(hist_filename),mode=NF90_NOWRITE,ncid=ncid)
534 [ + + ]: 46 : if(ncerr /= NF90_NOERR) then
535 : 39 : write(std_out,'(3a)') '-'//'Could not open ',trim(hist_filename),', starting from scratch'
536 : 39 : Invar%netcdf=.false.
537 : : else
538 : 7 : write(std_out,'(3a)') '-'//'Succesfully open ',trim(hist_filename),' for reading'
539 : 7 : write(std_out,'(a)') ' Extracting information from NetCDF file...'
540 : 7 : Invar%netcdf=.true.
541 : : end if
542 : :
543 [ + + ]: 46 : if ( Invar%netcdf) then
544 : : call get_dims_hist(ncid,Invar%natom,Invar%ntypat,nimage,mdtime,&
545 : 7 : & natom_id,ntypat_id,nimage_id,time_id,xyz_id,six_id,has_nimage)
546 [ - + + - : 30 : ABI_CALLOC(Invar%amu,(Invar%ntypat))
+ - + - -
+ - + - +
+ + ]
547 [ - + + - : 901 : ABI_CALLOC(Invar%typat,(Invar%natom))
+ - + - -
+ - + - +
+ + ]
548 [ - + + - : 30 : ABI_CALLOC(Invar%znucl,(Invar%ntypat))
+ - + - -
+ - + - +
+ + ]
549 : 7 : call read_csts_hist(ncid,dtion,Invar%typat,Invar%znucl,Invar%amu)
550 : :
551 : : ! Need to close NetCDF file because it is going to be reopened by read_md_hist
552 : 7 : ncerr = nf90_close(ncid)
553 : : ! .false. -> Velocities note used
554 : : ! .true. -> acell and rprimd may change (2017_04 only NVT/isoK used but maybe
555 : : ! .false. -> read all times
556 : : ! NPT one day ?)
557 : 21 : call read_md_hist(hist_filename,Hist,.false.,.true.,.false.)
558 : : end if
559 : : #endif
560 : :
561 : : ! =========================================================================== !
562 : : ! Read input file
563 : :
564 : 46 : string = repeat(" ", strlen)
565 : 46 : raw_string = repeat(" ", strlen)
566 : 46 : call instrng(input_filename, lenstr, 1, strlen, string, raw_string)
567 : : ! To make case-insensitive, map characters to upper case.
568 : 46 : call inupper(string(1:lenstr))
569 : :
570 : : ! marr is the aximum array size. It is thus a hard-coded maximum value
571 : : ! for 3 * (number of atoms)
572 : 46 : marr = 9000
573 [ - + ]: 46 : ABI_MALLOC(intarr, (marr))
574 [ - + ]: 46 : ABI_MALLOC(dprarr, (marr))
575 : 46 : jdtset = 0
576 : :
577 : : ! Mandatory input variables
578 : : ! -------------------------
579 : :
580 : : ! The Bravais lattice is found from rprimd, multiplicity and xred_unitcell
581 : : ! (see tdep_init_crystal): brav and angle are no longer input variables.
582 : 46 : call intagm(dprarr, intarr, jdtset, marr, 2, string(1:lenstr), 'brav', tread, 'INT')
583 [ - + ]: 46 : if (tread == 1) then
584 : 0 : ABI_WARNING('The input variable "brav" is no longer used: the Bravais lattice is found from the primitive vectors.')
585 : : end if
586 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'angle', tread, 'DPR')
587 [ - + ]: 46 : if (tread == 1) then
588 : 0 : ABI_WARNING('The input variable "angle" is no longer used: it is found from the primitive vectors.')
589 : : end if
590 : :
591 : : ! natom_unitcell
592 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'natom_unitcell', tread, 'INT')
593 [ - + ]: 46 : if (tread == 0) then
594 : : write(msg,*)&
595 : 0 : 'Variable "natom_unitcell" is mandatory, but was not found in input file.'
596 : 0 : ABI_ERROR(msg)
597 : : end if
598 : 46 : Invar%natom_unitcell = intarr(1)
599 [ - + ]: 46 : if (3*Invar%natom_unitcell.gt.marr) then
600 : : write(msg,*)&
601 : 0 : 'Maximum number of atoms exceeded. Modify source code to circumvent problem.'
602 : 0 : ABI_ERROR(msg)
603 : : end if
604 : :
605 : : ! xred_unitcell
606 [ - + + - : 582 : ABI_CALLOC(Invar%xred_unitcell,(3,Invar%natom_unitcell))
+ - + - -
+ - + - +
+ + + + ]
607 : 46 : call intagm(dprarr, intarr, jdtset, marr, 3*Invar%natom_unitcell, string(1:lenstr), 'xred_unitcell', tread, 'DPR')
608 [ - + ]: 46 : if (tread == 0) then
609 : : write(msg,*)&
610 : 0 : 'Variable "xred_unitcell" is mandatory, but was not found in input file.'
611 : 0 : ABI_ERROR(msg)
612 : : end if
613 [ + + + + : 582 : Invar%xred_unitcell(:,:) = reshape(dprarr(1:3*Invar%natom_unitcell),(/3,Invar%natom_unitcell/))
+ + ]
614 [ + + ]: 184 : do ii=1,3
615 [ + + ]: 517 : do iatcell=1,Invar%natom_unitcell
616 [ + - + + ]: 471 : if ((Invar%xred_unitcell(ii,iatcell).le.(-0.5)).or.(Invar%xred_unitcell(ii,iatcell).gt.(0.5))) then
617 [ - + ]: 3 : do while (Invar%xred_unitcell(ii,iatcell).le.(-0.5))
618 : 0 : Invar%xred_unitcell(ii,iatcell)=Invar%xred_unitcell(ii,iatcell)+1.d0
619 : : end do
620 [ + + ]: 6 : do while (Invar%xred_unitcell(ii,iatcell).gt.(0.5))
621 : 6 : Invar%xred_unitcell(ii,iatcell)=Invar%xred_unitcell(ii,iatcell)-1.d0
622 : : end do
623 : : end if
624 : : end do
625 : : end do
626 : :
627 : : ! typat_unitcell
628 [ - + + - : 249 : ABI_CALLOC(Invar%typat_unitcell,(Invar%natom_unitcell))
+ - + - -
+ - + - +
+ + ]
629 : 46 : call intagm(dprarr, intarr, jdtset, marr, Invar%natom_unitcell, string(1:lenstr), 'typat_unitcell', tread, 'INT')
630 [ - + ]: 46 : if (tread == 0) then
631 : : write(msg,*)&
632 : 0 : 'Variable "typat_unitcell" is mandatory, but was not found in input file.'
633 : 0 : ABI_ERROR(msg)
634 : : end if
635 [ + + ]: 157 : Invar%typat_unitcell(:) = intarr(1:Invar%natom_unitcell)
636 : :
637 [ + + ]: 46 : if (Invar%netcdf) then
638 [ + + + + ]: 91 : Invar%rprimd_md(:,:)=TRANSPOSE(Hist%rprimd(:,:,Hist%ihist))
639 : : else
640 : : ! ntypat
641 : 39 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'ntypat', tread, 'INT')
642 [ + - ]: 39 : if (tread == 1) then
643 : 39 : Invar%ntypat = intarr(1)
644 : : else
645 : : write(msg,*)&
646 : 0 : 'The NetCDF file .nc is not used.',ch10,&
647 : 0 : 'The variable "ntypat" is thus mandatory.',ch10,&
648 : 0 : 'ACTION : Please modify your input file'
649 : 0 : ABI_ERROR(msg)
650 : : end if
651 : :
652 : : ! amu
653 [ - + + - : 166 : ABI_MALLOC(Invar%amu,(Invar%ntypat)); Invar%amu(:)=zero
+ - + - -
+ - + - +
+ + ]
654 : 39 : call intagm(dprarr, intarr, jdtset, marr, Invar%ntypat, string(1:lenstr), 'amu', tread, 'DPR')
655 [ + - ]: 39 : if (tread == 1) then
656 [ + - - + : 127 : Invar%amu = dprarr(1:Invar%ntypat)
- - + + ]
657 : : else
658 : : write(msg,*)&
659 : 0 : 'The NetCDF file .nc is not used.',ch10,&
660 : 0 : 'The variable "amu" is thus mandatory.',ch10,&
661 : 0 : 'ACTION : Please modify your input file'
662 : 0 : ABI_ERROR(msg)
663 : : end if
664 : :
665 : : ! rprimd_md
666 : 39 : call intagm(dprarr, intarr, jdtset, marr, 9, string(1:lenstr), 'rprimd', tread, 'LEN')
667 [ + - ]: 39 : if (tread == 1) then
668 [ + + + + ]: 507 : Invar%rprimd_md(:,:) = TRANSPOSE(reshape(dprarr(1:9),(/3, 3/)))
669 : : else
670 : : write(msg,*)&
671 : 0 : 'The NetCDF file .nc is not used.',ch10,&
672 : 0 : 'The variable "rprimd" is thus mandatory.',ch10,&
673 : 0 : 'ACTION : Please modify your input file'
674 : 0 : ABI_ERROR(msg)
675 : : end if
676 : :
677 : : ! natom
678 : 39 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'natom', tread, 'INT')
679 [ + - ]: 39 : if (tread == 1) then
680 : 39 : Invar%natom = intarr(1)
681 : : else
682 : : write(msg,*)&
683 : 0 : 'The NetCDF file .nc is not used.',ch10,&
684 : 0 : 'The variable "natom" is thus mandatory.',ch10,&
685 : 0 : 'ACTION : Please modify your input file'
686 : 0 : ABI_ERROR(msg)
687 : : end if
688 : :
689 : : ! typat
690 [ - + + - : 6515 : ABI_MALLOC(Invar%typat,(Invar%natom)); Invar%typat(:)=0
+ - + - -
+ - + - +
+ + ]
691 : 39 : call intagm(dprarr, intarr, jdtset, marr, Invar%natom, string(1:lenstr), 'typat', tread, 'INT')
692 [ + - ]: 39 : if (tread == 1) then
693 [ + - - + : 6476 : Invar%typat = intarr(1:Invar%natom)
- - + + ]
694 : : else
695 : : write(msg,*)&
696 : 0 : 'The NetCDF file .nc is not used.',ch10,&
697 : 0 : 'The variable "typat" is thus mandatory.',ch10,&
698 : 0 : 'ACTION : Please modify your input file'
699 : 0 : ABI_ERROR(msg)
700 : : end if
701 : :
702 : : end if
703 : :
704 : : ! multiplicity
705 : 46 : call intagm(dprarr, intarr, jdtset, marr, 9, string(1:lenstr), 'multiplicity', tread, 'DPR')
706 [ + - ]: 46 : if (tread == 1) then
707 [ + + + + ]: 598 : Invar%multiplicity(:,:) = TRANSPOSE(reshape(dprarr(1:9),(/3, 3/)))
708 : : else
709 : : write(msg,*)&
710 : 0 : 'Variable "multiplicity" is mandatory, but was not found in input file.'
711 : 0 : ABI_ERROR(msg)
712 : : end if
713 : :
714 : : ! temperature
715 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'temperature', tread, 'DPR')
716 [ + - ]: 46 : if (tread == 1) then
717 : 46 : Invar%temperature = dprarr(1)
718 : : else
719 : : write(msg,*)&
720 : 0 : 'Variable "temperature" is mandatory, but was not found in input file.'
721 : 0 : ABI_ERROR(msg)
722 : : end if
723 : :
724 : : ! nstep_max
725 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'nstep_max', tread, 'INT')
726 [ + - ]: 46 : if (tread == 1) then
727 : 46 : Invar%nstep_max = intarr(1)
728 : : else
729 : : write(msg,*)&
730 : 0 : 'Variable "nstep_max" is mandatory, but was not found in input file.'
731 : 0 : ABI_ERROR(msg)
732 : : end if
733 : :
734 : : ! nstep_min
735 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'nstep_min', tread, 'INT')
736 [ + - ]: 46 : if (tread == 1) then
737 : 46 : Invar%nstep_min = intarr(1)
738 : : else
739 : : write(msg,*)&
740 : 0 : 'Variable "nstep_min" is mandatory, but was not found in input file.'
741 : 0 : ABI_ERROR(msg)
742 : : end if
743 : :
744 : : ! rcut
745 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'rcut', tread, 'LEN')
746 [ + - ]: 46 : if (tread == 1) then
747 : 46 : Invar%rcut = dprarr(1)
748 : : else
749 : : write(msg,*)&
750 : 0 : 'Variable "rcut" is mandatory, but was not found in input file.'
751 : 0 : ABI_ERROR(msg)
752 : : end if
753 : :
754 : : ! znucl (optional)
755 [ + + ]: 46 : if (.not.Invar%netcdf) then
756 : 39 : call intagm(dprarr, intarr, jdtset, marr, Invar%ntypat, string(1:lenstr), 'znucl', tread, 'DPR')
757 [ + - ]: 39 : if (tread == 1) then
758 [ - + + - : 166 : ABI_MALLOC(Invar%znucl,(Invar%ntypat)); Invar%znucl(:)=zero
+ - + - -
+ - + - +
+ + ]
759 [ - + + + ]: 127 : Invar%znucl = dprarr(1:Invar%ntypat)
760 : : end if
761 : : end if
762 : :
763 : : ! =========================================================================== !
764 : : ! Read DDB file if available to retrieve the Born effective charges
765 : : ! and dielectric tensor. Enforce charge neutrality by equal redistribution of excess charge
766 : : ! =========================================================================== !
767 : :
768 : : ! symbec
769 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'symbec', tread, 'INT')
770 [ - + ]: 46 : if (tread == 1) then
771 : 0 : Invar%symbec = intarr(1)
772 : : end if
773 : :
774 : :
775 [ - + + - : 138 : ABI_MALLOC(Invar%zeff, (3,3,Invar%natom_unitcell))
+ - + - -
+ - + -
+ ]
776 [ + + + + : 1489 : Invar%zeff = zero
+ + ]
777 : :
778 [ + - - + ]: 46 : if (file_exists(ddb_filename).or.file_exists(nctk_ncify(ddb_filename))) then
779 : :
780 : 0 : call ddb%from_file(ddb_filename, ddb_hdr, crystal, xmpi_world, -1)
781 : : ! TODO Check that the crystal in the DDB matches the unit cell.
782 : :
783 : 0 : ii = ddb%get_dielt_zeff(crystal,1,1,0,Invar%dielt,Invar%zeff)
784 [ # # ]: 0 : if (ii/=0) Invar%dipdip = .true.
785 : :
786 : 0 : call crystal%free()
787 : 0 : call ddb_hdr%free()
788 : 0 : call ddb%free()
789 : :
790 : : !! Take the average born effective charge for each atom.
791 [ - - ]: 46 : if (Invar%symbec==1) then
792 [ # # ]: 0 : do iatom=1,Invar%natom_unitcell
793 : : born_average = zero
794 [ # # ]: 0 : do ii=1,3
795 : 0 : born_average = born_average + Invar%zeff(ii,ii,iatom)
796 : : end do
797 [ # # # # ]: 0 : Invar%zeff(:,:,iatom) = zero
798 [ # # ]: 0 : do ii=1,3
799 : 0 : Invar%zeff(ii,ii,iatom) = born_average / 3
800 : : end do
801 : : end do
802 : : end if
803 : : end if
804 : :
805 : :
806 : : ! =========================================================================== !
807 : : ! Optional input variables
808 : : ! =========================================================================== !
809 : :
810 : : ! prtvol
811 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'prtvol', tread, 'INT')
812 [ - + ]: 46 : if (tread == 1) then
813 : 0 : Invar%prtvol = intarr(1)
814 : : end if
815 : :
816 : : ! prtphdos
817 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'prtphdos', tread, 'INT')
818 [ - + ]: 46 : if (tread == 1) then
819 : 0 : Invar%prtphdos = intarr(1)
820 : : end if
821 : :
822 : : ! dosdeltae
823 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'dosdeltae', tread, 'ENE')
824 [ - + ]: 46 : if (tread == 1) then
825 : 0 : Invar%dosdeltae = dprarr(1)
826 : : end if
827 : :
828 : : ! use_ideal_positions
829 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'use_ideal_positions', tread, 'INT')
830 [ + + ]: 46 : if (tread == 1) then
831 : 42 : Invar%use_ideal_positions = intarr(1)
832 : : end if
833 : :
834 : : ! born_charge
835 : 46 : call intagm(dprarr, intarr, jdtset, marr, Invar%ntypat, string(1:lenstr), 'born_charge', tread, 'DPR')
836 [ + + ]: 46 : if (tread == 1) then
837 : 5 : Invar%dipdip=.true.
838 [ - + + - : 25 : ABI_MALLOC(born_charge_tmp,(Invar%ntypat)) ; born_charge_tmp(:) = zero
+ - + - -
+ - + +
+ ]
839 [ + + ]: 15 : born_charge_tmp(:) = dprarr(1:Invar%ntypat)
840 [ + + + + : 148 : Invar%zeff(:,:,:) = zero
+ + ]
841 [ + + ]: 16 : do iatom=1,Invar%natom_unitcell
842 : 11 : itypat = Invar%typat_unitcell(iatom)
843 [ + + ]: 49 : do ii=1,3
844 : 44 : Invar%zeff(ii,ii,iatom) = born_charge_tmp(itypat)
845 : : end do
846 : : end do
847 [ - + ]: 5 : ABI_FREE(born_charge_tmp)
848 : : end if
849 : :
850 : : ! dielec_constant
851 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'dielec_constant', tread, 'DPR')
852 [ + + ]: 46 : if (tread == 1) then
853 [ + + ]: 20 : do ii=1,3
854 : 20 : Invar%dielt(ii,ii) = dprarr(1)
855 : : end do
856 : : end if
857 : :
858 : : ! dipdip
859 : : ! Invar%dipdip is auto-enabled above whenever usable Born effective
860 : : ! charges/dielectric tensor are found (from the DDB file or born_charge).
861 : : ! Explicitly setting dipdip=0 overrides that and forces it off.
862 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'dipdip', tread, 'INT')
863 [ - + ]: 46 : if (tread == 1) then
864 [ # # ]: 0 : if (intarr(1) == 0) Invar%dipdip = .false.
865 : : end if
866 : :
867 : : ! bzpath
868 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'bzpath', tread, 'INT')
869 [ + + ]: 46 : if (tread == 1) then
870 : 10 : Invar%bzpath = intarr(1)
871 : :
872 [ - + ]: 10 : if (Invar%bzpath.lt.0) then
873 [ # # # # : 0 : ABI_MALLOC(Invar%qpt,(3,abs(Invar%bzpath))); Invar%qpt(:,:)=zero
# # # # #
# # # # #
# # ]
874 : :
875 : 0 : call intagm(dprarr, intarr, jdtset, marr, 1+3*abs(Invar%bzpath), string(1:lenstr), 'bzpath', tread, 'DPR')
876 [ # # # # : 0 : Invar%qpt(:,:) = reshape(dprarr(2:1+3*abs(Invar%bzpath)), (/3,abs(Invar%bzpath)/))
# # ]
877 : :
878 [ + - ]: 10 : else if (Invar%bzpath.gt.0) then
879 [ - + + - : 20 : ABI_MALLOC(Invar%special_qpt,(Invar%bzpath))
- + - + -
+ ]
880 : 10 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'special_qpt', tread, 'KEY', key_value=msg)
881 [ + - ]: 10 : if (tread == 1) then
882 : 10 : jj = 1
883 [ + + ]: 64 : do ii=1,Invar%bzpath
884 : 54 : Invar%special_qpt(ii) = " "
885 : 64 : ierr = next_token(msg,jj,Invar%special_qpt(ii))
886 : : !jj = 2*ii-1
887 : : !Invar%special_qpt(ii) = msg(jj:jj)
888 : : end do
889 : : else
890 : : write(msg,*)&
891 : 0 : 'Variable bzpath > 0, but variable "special_qpt" was not found in input file.'
892 : 0 : ABI_ERROR(msg)
893 : : end if
894 : : end if
895 : : end if
896 : :
897 : : ! ndivsm
898 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'ndivsm', tread, 'INT')
899 [ - + ]: 46 : if (tread == 1) then
900 : 0 : Invar%ndivsm = intarr(1)
901 : : end if
902 : :
903 : : ! order
904 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'order', tread, 'INT')
905 [ + + ]: 46 : if (tread == 1) then
906 : 16 : Invar%order = intarr(1)
907 [ + + ]: 16 : if (Invar%order.eq.3) then
908 : 7 : call intagm(dprarr, intarr, jdtset, marr, 2, string(1:lenstr), 'order', tread, 'DPR')
909 : 7 : Invar%rcut3 = dprarr(2)
910 [ - + ]: 7 : if (Invar%rcut3.gt.Invar%rcut) then
911 : 0 : ABI_ERROR('The cutoff radius of the third order cannot be greater than the second order one.')
912 : : end if
913 [ + - ]: 9 : else if (Invar%order.eq.4) then
914 : 9 : call intagm(dprarr, intarr, jdtset, marr, 3, string(1:lenstr), 'order', tread, 'DPR')
915 : 9 : Invar%rcut3 = dprarr(2)
916 : 9 : Invar%rcut4 = dprarr(3)
917 [ - + ]: 9 : if (Invar%rcut4.gt.Invar%rcut) then
918 : 0 : ABI_ERROR('The cutoff radius of the fourth order cannot be greater than the second order one.')
919 : : end if
920 : : else
921 : 0 : ABI_ERROR('Only the 3rd and 4th orders are allowed. Change your input file.')
922 : : end if
923 : : end if
924 : :
925 : : ! slice
926 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'slice', tread, 'INT')
927 [ + + ]: 46 : if (tread == 1) then
928 : 2 : Invar%slice = intarr(1)
929 : : end if
930 : :
931 : : ! enunit
932 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'enunit', tread, 'INT')
933 [ + + ]: 46 : if (tread == 1) then
934 : 38 : Invar%enunit = intarr(1)
935 : : end if
936 : :
937 : : ! nproc
938 : 46 : call intagm(dprarr, intarr, jdtset, marr, 2, string(1:lenstr), 'nproc', tread, 'INT')
939 [ - + ]: 46 : if (tread == 1) then
940 [ # # ]: 0 : Invar%nproc(1:2) = intarr(1:2)
941 : : end if
942 : :
943 : : ! readifc
944 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'readifc', tread, 'INT')
945 [ + + ]: 46 : if (tread == 1) then
946 : 3 : Invar%readifc = intarr(1)
947 [ + + ]: 3 : if (Invar%readifc.eq.1) then
948 : 1 : call intagm(dprarr, intarr, jdtset, marr, 2, string(1:lenstr), 'readifc', tread, 'DPR')
949 : 1 : Invar%tolread = dprarr(1)
950 : : end if
951 : : end if
952 : :
953 : : ! alloy
954 : 46 : call intagm(dprarr, intarr, jdtset, marr, 3, string(1:lenstr), 'alloy', tread, 'INT')
955 [ + + ]: 46 : if (tread == 1) then
956 : 1 : Invar%alloy = intarr(1)
957 : 1 : Invar%ityp_alloy1 = intarr(2)
958 : 1 : Invar%ityp_alloy2 = intarr(3)
959 : : end if
960 : :
961 : : ! together
962 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'together', tread, 'INT')
963 [ + + ]: 46 : if (tread == 1) then
964 : 2 : Invar%together = intarr(1)
965 : : end if
966 : :
967 : : ! bzlength
968 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'bzlength', tread, 'INT')
969 [ + + ]: 46 : if (tread == 1) then
970 : 2 : Invar%bzlength = intarr(1)
971 [ - + + - : 6 : ABI_MALLOC(Invar%lgth_segments,(Invar%bzlength))
+ - + - -
+ - + -
+ ]
972 : 2 : call intagm(dprarr, intarr, jdtset, marr, 1+Invar%bzlength, string(1:lenstr), 'bzlength', tread, 'DPR')
973 [ + - - + : 14 : Invar%lgth_segments = dprarr(2:1+Invar%bzlength)
- - + + ]
974 : : end if
975 : :
976 : : ! ngqpt1
977 : 46 : call intagm(dprarr, intarr, jdtset, marr, 3, string(1:lenstr), 'ngqpt1', tread, 'INT')
978 [ + + ]: 46 : if (tread == 1) then
979 [ + + ]: 52 : Invar%ngqpt1 = intarr(1:3)
980 : : end if
981 : :
982 : : ! Coarse q-point grid: from ngqpt1 if given, otherwise from the supercell
983 [ + + + + ]: 598 : Invar%qptrlatt(:,:)=0
984 [ + + + + : 145 : if (any(Invar%ngqpt1(:) > 0)) then
+ + ]
985 [ + + ]: 52 : do ii=1,3
986 : 52 : Invar%qptrlatt(ii,ii)=Invar%ngqpt1(ii)
987 : : end do
988 : : else
989 [ + + + + : 429 : if (any(abs(Invar%multiplicity(:,:)-nint(Invar%multiplicity(:,:))) > tol8)) then
+ - - + ]
990 : 0 : ABI_ERROR('The multiplicity must be an integer matrix.')
991 : : end if
992 [ + + + + : 825 : Invar%qptrlatt(:,:)=TRANSPOSE(nint(Invar%multiplicity(:,:)))
+ + + + ]
993 : : end if
994 : :
995 : : ! q1shft
996 : 46 : call intagm(dprarr, intarr, jdtset, marr, 3, string(1:lenstr), 'q1shft', tread, 'DPR')
997 [ - + ]: 46 : if (tread == 1) then
998 [ # # ]: 0 : Invar%q1shft = dprarr(1:3)
999 : : end if
1000 : :
1001 : : ! q2shft
1002 : 46 : call intagm(dprarr, intarr, jdtset, marr, 3, string(1:lenstr), 'q2shft', tread, 'DPR')
1003 [ - + ]: 46 : if (tread == 1) then
1004 [ # # ]: 0 : Invar%q2shft = dprarr(1:3)
1005 : : end if
1006 : :
1007 : : ! ngqpt2
1008 : 46 : call intagm(dprarr, intarr, jdtset, marr, 3, string(1:lenstr), 'ngqpt2', tread, 'INT')
1009 [ + + ]: 46 : if (tread == 1) then
1010 [ + + ]: 92 : Invar%ngqpt2 = intarr(1:3)
1011 : : end if
1012 : :
1013 : : ! use_weights
1014 : 46 : call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'use_weights', tread, 'INT')
1015 [ + + ]: 46 : if (tread == 1) then
1016 : 1 : Invar%use_weights = intarr(1)
1017 : : end if
1018 : :
1019 : 46 : ABI_FREE(intarr)
1020 [ - + ]: 46 : ABI_FREE(dprarr)
1021 : :
1022 : : ! =========================================================================== !
1023 : : ! Perform some checks
1024 : :
1025 [ - + ]: 46 : if ((Invar%use_ideal_positions.ne.0).and.(Invar%use_ideal_positions.ne.1)) then
1026 : 0 : ABI_ERROR('STOP: THIS VALUE IS NOT ALLOWED FOR use_ideal_positions')
1027 : : end if
1028 : :
1029 : : ! Allowed values
1030 [ - + ]: 46 : if ((Invar%together.ne.1).and.(Invar%together.ne.0)) then
1031 : 0 : ABI_ERROR('STOP: The value of input variable TOGETHER is not allowed')
1032 : : end if
1033 [ - + ]: 46 : if ((Invar%alloy.ne.1).and.(Invar%alloy.ne.0)) then
1034 : 0 : ABI_ERROR('STOP: The value of input variable ALLOY is not allowed')
1035 : : end if
1036 [ + + ]: 46 : if (Invar%alloy.ge.1) then
1037 : : if ((Invar%ityp_alloy1.lt.1).or.(Invar%ityp_alloy2.lt.1).or.&
1038 [ + - + - : 1 : & (Invar%ityp_alloy1.gt.Invar%natom_unitcell).or.(Invar%ityp_alloy2.gt.Invar%natom_unitcell)) then
+ - - + ]
1039 : 0 : ABI_ERROR('STOP: The value of input variables IALLOY are not allowed')
1040 : : end if
1041 : : end if
1042 : :
1043 : : ! Incompatible variables :
1044 [ + + - + ]: 46 : if ((Invar%readifc.eq.1).and.(Invar%together.eq.1).and.(Invar%order.gt.2)) then
1045 : 0 : ABI_ERROR('STOP: readifc=1, together=1 and order=3 or 4 are incompatible')
1046 : : end if
1047 : :
1048 : : ! =========================================================================== !
1049 : : ! GA: TODO REMOVE ALLOYING FEATURE.
1050 : : ! We agreed this logic can be handled externally.
1051 : : !
1052 : : ! Treat virtual crystal approximation (VCA) when alloy=1.
1053 : : ! Redefine all the data depending on (n)typat(_unitcell) and natom_unitcell
1054 : :
1055 : : ! The mixing of zeff will only work for the isotropic case,
1056 : : ! that is, if it was specified with input variables born_charge
1057 : :
1058 [ + + ]: 46 : if (Invar%alloy.eq.1) then
1059 : 1 : sum_alloy1=0
1060 : 1 : sum_alloy2=0
1061 [ + + ]: 129 : do iatom=1,Invar%natom
1062 [ + + ]: 128 : if (Invar%typat(iatom).eq.Invar%ityp_alloy1) then
1063 : 96 : sum_alloy1=sum_alloy1+1
1064 : : end if
1065 [ + + ]: 129 : if (Invar%typat(iatom).eq.Invar%ityp_alloy2) then
1066 : 32 : sum_alloy2=sum_alloy2+1
1067 : : end if
1068 : : end do
1069 : : amu_average =(Invar%amu(Invar%ityp_alloy1)*sum_alloy1+&
1070 : 1 : & Invar%amu(Invar%ityp_alloy2)*sum_alloy2)/(sum_alloy1+sum_alloy2)
1071 [ - + ]: 1 : if (Invar%dipdip) then
1072 : 0 : zeff_average(:,:) = zero
1073 [ # # ]: 0 : do iatom=1,Invar%natom_unitcell
1074 [ # # ]: 0 : if (Invar%typat_unitcell(iatom).eq.Invar%ityp_alloy1) then
1075 [ # # # # ]: 0 : zeff_average(:,:) = zeff_average(:,:) + Invar%zeff(:,:,iatom) * sum_alloy1
1076 : : exit
1077 : : end if
1078 : : end do
1079 [ # # ]: 0 : do iatom=1,Invar%natom_unitcell
1080 [ # # ]: 0 : if (Invar%typat_unitcell(iatom).eq.Invar%ityp_alloy2) then
1081 [ # # # # ]: 0 : zeff_average(:,:) = zeff_average(:,:) + Invar%zeff(:,:,iatom) * sum_alloy2
1082 : : exit
1083 : : end if
1084 : : end do
1085 [ # # # # ]: 0 : zeff_average(:,:) = zeff_average(:,:) / (sum_alloy1 + sum_alloy2)
1086 [ # # ]: 0 : do iatom=1,Invar%natom_unitcell
1087 [ # # ]: 0 : if (Invar%typat_unitcell(iatom).eq.min(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1088 [ # # # # ]: 0 : Invar%zeff(:,:,iatom) = zeff_average(:,:)
1089 : : end if
1090 : : end do
1091 : : end if
1092 : 1 : shift=0
1093 [ + + ]: 3 : do iatom=1,Invar%natom_unitcell
1094 [ + + ]: 3 : if (Invar%typat_unitcell(iatom).lt.max(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1095 : 1 : Invar%typat_unitcell (iatom-shift)=Invar%typat_unitcell (iatom)
1096 [ + + ]: 4 : Invar%xred_unitcell(:,iatom-shift)=Invar%xred_unitcell(:,iatom)
1097 [ - + ]: 1 : if (Invar%dipdip) then
1098 [ # # # # ]: 0 : Invar%zeff(:,:,iatom-shift) = Invar%zeff(:,:,iatom)
1099 : : end if
1100 [ + - ]: 1 : else if (Invar%typat_unitcell(iatom).eq.max(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1101 : 1 : shift=shift+1
1102 : : else if (Invar%typat_unitcell(iatom).gt.max(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1103 : 0 : Invar%typat_unitcell (iatom-shift)=Invar%typat_unitcell (iatom) - 1
1104 [ # # ]: 0 : Invar%xred_unitcell(:,iatom-shift)=Invar%xred_unitcell(:,iatom)
1105 [ # # ]: 0 : if (Invar%dipdip) then
1106 [ # # # # ]: 0 : Invar%zeff(:,:,iatom-shift) = Invar%zeff(:,:,iatom)
1107 : : end if
1108 : : end if
1109 : : end do
1110 : 1 : Invar%natom_unitcell=Invar%natom_unitcell-shift
1111 [ + + ]: 129 : do iatom=1,Invar%natom
1112 [ + + ]: 129 : if (Invar%typat(iatom).ge.max(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1113 : 32 : Invar%typat(iatom)=Invar%typat(iatom) - 1
1114 : : end if
1115 : : end do
1116 [ + + ]: 3 : do itypat=1,Invar%ntypat
1117 [ + + ]: 3 : if (itypat.eq.min(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1118 : 1 : Invar%amu (itypat)=amu_average
1119 [ - + ]: 1 : else if (itypat.gt.max(Invar%ityp_alloy1,Invar%ityp_alloy2)) then
1120 : 0 : Invar%amu (itypat-1)=Invar%amu (itypat)
1121 : : end if
1122 : : end do
1123 : 1 : Invar%ntypat=Invar%ntypat-1
1124 [ - + + - : 3 : ABI_MALLOC(typat_unitcell_tmp,( Invar%natom_unitcell))
+ - + - -
+ - + ]
1125 [ - + + - : 3 : ABI_MALLOC(xred_unitcell_tmp ,(3,Invar%natom_unitcell))
+ - + - -
+ - + ]
1126 [ - + + - : 3 : ABI_MALLOC(amu_tmp ,( Invar%ntypat))
+ - + - -
+ - + ]
1127 [ - + ]: 1 : if (Invar%dipdip) then
1128 [ # # # # : 0 : ABI_MALLOC(zeff_tmp ,(3,3,Invar%natom_unitcell))
# # # # #
# # # ]
1129 : : end if
1130 [ + + ]: 2 : typat_unitcell_tmp (:)=Invar%typat_unitcell(1:Invar%natom_unitcell)
1131 [ + + + + ]: 5 : xred_unitcell_tmp(:,:)=Invar%xred_unitcell(:,1:Invar%natom_unitcell)
1132 [ + + ]: 2 : amu_tmp (:)=Invar%amu(1:Invar%ntypat)
1133 [ - + ]: 1 : if (Invar%dipdip) then
1134 [ # # # # : 0 : zeff_tmp(:,:,:) = Invar%zeff(:,:,1:Invar%natom_unitcell)
# # ]
1135 : : end if
1136 [ + - - + : 3 : ABI_REMALLOC(Invar%typat_unitcell,( Invar%natom_unitcell))
+ - + - +
- - + - +
- + ]
1137 [ + - - + : 3 : ABI_REMALLOC(Invar%xred_unitcell ,(3,Invar%natom_unitcell))
+ - + - +
- - + - +
- + ]
1138 [ + - - + : 3 : ABI_REMALLOC(Invar%amu ,( Invar%ntypat))
+ - + - +
- - + - +
- + ]
1139 [ + - - + : 3 : ABI_REMALLOC(Invar%zeff ,(3,3,Invar%natom_unitcell))
+ - + - +
- - + - +
- + ]
1140 [ + + ]: 2 : Invar%typat_unitcell (:)=typat_unitcell_tmp (:)
1141 [ + + + + ]: 5 : Invar%xred_unitcell(:,:)=xred_unitcell_tmp(:,:)
1142 [ + + ]: 2 : Invar%amu (:)=amu_tmp (:)
1143 [ - + ]: 1 : if (Invar%dipdip) then
1144 [ # # # # : 0 : Invar%zeff(:,:,:) =zeff_tmp(:,:,:)
# # ]
1145 : : end if
1146 [ + - ]: 1 : if (allocated(Invar%znucl)) then
1147 [ - + + - : 3 : ABI_MALLOC(znucl_tmp,(Invar%ntypat))
+ - + - -
+ - + ]
1148 [ + + ]: 2 : znucl_tmp(:) = Invar%znucl(1:Invar%ntypat)
1149 [ - + + - : 3 : ABI_REMALLOC(Invar%znucl,(Invar%ntypat))
+ - + - -
+ - + ]
1150 [ + + ]: 2 : Invar%znucl(:) = znucl_tmp(:)
1151 [ - + ]: 1 : ABI_FREE(znucl_tmp)
1152 : : end if
1153 [ - + ]: 1 : ABI_FREE(typat_unitcell_tmp)
1154 [ - + ]: 1 : ABI_FREE(xred_unitcell_tmp)
1155 [ - + ]: 1 : ABI_FREE(amu_tmp)
1156 [ - + ]: 1 : if (Invar%dipdip) then
1157 [ # # ]: 0 : ABI_FREE(zeff_tmp)
1158 : : end if
1159 : : end if
1160 : :
1161 : : ! Anchor the unitcell on atom 1: shift every atom by atom 1's own reduced
1162 : : ! position and wrap into [0,1). This is not just a local convenience: the
1163 : : ! symmetry finder (tdep_make_sym) needs atom 1 at the origin to find the
1164 : : ! space group, and this same convention is reused, unmodified, as the
1165 : : ! reference unitcell (Sym%xred_zero, Crystal%xred, SC%xred_unitcell) for the
1166 : : ! rest of the code. The original (as-given-in-input) positions are kept
1167 : : ! separately in xred_unitcell_original, so that the rigid shift between the
1168 : : ! two conventions can be applied explicitly wherever needed (see
1169 : : ! tdep_sampling_init_read).
1170 [ - + + - : 138 : ABI_MALLOC(Invar%xred_unitcell_original,(3,Invar%natom_unitcell))
+ - + - -
+ - + -
+ ]
1171 [ + + + + ]: 486 : Invar%xred_unitcell_original(:,:) = Invar%xred_unitcell(:,:)
1172 [ + + ]: 156 : do iatcell=1,Invar%natom_unitcell
1173 [ + + ]: 440 : Invar%xred_unitcell(:,iatcell) = Invar%xred_unitcell(:,iatcell) - Invar%xred_unitcell_original(:,1)
1174 [ + + ]: 486 : do ii=1,3
1175 [ + + ]: 330 : if (Invar%xred_unitcell(ii,iatcell).lt.0.d0) then
1176 : 66 : jj=int(Invar%xred_unitcell(ii,iatcell)+1.0)
1177 : : else
1178 : 264 : jj=int(Invar%xred_unitcell(ii,iatcell)+0.0)
1179 : : end if
1180 : 440 : Invar%xred_unitcell(ii,iatcell) = Invar%xred_unitcell(ii,iatcell) - real(jj)
1181 : : end do
1182 : : end do
1183 : :
1184 : : ! Compute Nstep as a function of the slice
1185 : 46 : Invar%nstep_tot=int(float(Invar%nstep_max-Invar%nstep_min)/float(Invar%slice)+1)
1186 : :
1187 [ - + - + : 92 : end subroutine tdep_read_input
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - - - -
- - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - - - -
- - - - -
- - - - -
- - - - -
- - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + -
+ ]
1188 : : !!***
1189 : :
1190 : : !====================================================================================================
1191 : :
1192 : : !!****f* ABINIT/m_tdep_dataset/tdep_write_input
1193 : : !! NAME
1194 : : !! tdep_write_input
1195 : : !!
1196 : : !! FUNCTION
1197 : : !! Echo the set of input variables to main output file.
1198 : : !!
1199 : : !! INPUTS
1200 : : !! Invar = Input object containing the input variables.
1201 : : !! unt = Unit number of the file to write in.
1202 : : !!
1203 : : !! OUTPUT
1204 : : !! Only writing.
1205 : : !!
1206 : : !! SOURCE
1207 : :
1208 : 46 : subroutine tdep_write_input(Invar, unt)
1209 : :
1210 : : class(atdep_dataset_type),intent(in) :: Invar
1211 : : integer, intent(in) :: unt
1212 : :
1213 : : integer :: ii,jj,kk
1214 : :
1215 : : ! *********************************************************************
1216 : :
1217 : 46 : write(unt,*) '#############################################################################'
1218 : 46 : write(unt,*) '######################### ECHO OF INPUT FILE ################################'
1219 : 46 : write(unt,*) '#############################################################################'
1220 : 46 : write(unt,*) ' '
1221 : 46 : write(unt,'(a)') ' ======================= Files ================================================'
1222 : 46 : write(unt,'(a1,a20,1x,a)')'-', ljust('indata_prefix',20),trim(Invar%input_prefix)
1223 : 46 : write(unt,'(a1,a20,1x,a)')'-', ljust('outdata_prefix',20),trim(Invar%output_prefix)
1224 : 46 : write(unt,'(a1,a20,1x,a)')'-', ljust('output_file',20),trim(Invar%output_file)
1225 : :
1226 : 46 : write(unt,'(a)') ' ======================= Define the unitcell ================================='
1227 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('natom_unitcell',20),Invar%natom_unitcell
1228 : 46 : write(unt,'(1x,a20)') ljust('xred_unitcell',20)
1229 [ + + ]: 156 : do ii=1,Invar%natom_unitcell
1230 [ + - + + ]: 486 : write(unt,'(22x,3(f15.10,1x))') (Invar%xred_unitcell(jj,ii), jj=1,3)
1231 : : end do
1232 : : ! write(unt,'(1x,a20)') ljust('xred_unitcell_orig',20)
1233 : : ! do ii=1,Invar%natom_unitcell
1234 : : ! write(unt,'(22x,3(f15.10,1x))') (Invar%xred_unitcell_original(jj,ii), jj=1,3)
1235 : : ! end do
1236 : : ! write(unt,'(a)') &
1237 : : !& ' NOTE: The unitcell reduced coordinates (xred_unitcell) are shifted so that atom 1 sits at'
1238 : : ! write(unt,'(a)') &
1239 : : !& ' the origin (0,0,0); this convention is used throughout the rest of the code.'
1240 : : ! write(unt,'(a)') &
1241 : : !& ' xred_unitcell_original holds the positions exactly as given in the input file.'
1242 [ + + ]: 46 : if (Invar%alloy.eq.1) then
1243 : : write(unt,'(a)') &
1244 : 1 : & ' NOTE: Virtual crystal approximation applied (alloy=1); unitcell variables above'
1245 : : write(unt,'(a)') &
1246 : 1 : & ' reflect the merged (averaged) atom type.'
1247 : : end if
1248 [ + - + + ]: 156 : write(unt,'(1x,a20,20(1x,i4))') ljust('typat_unitcell',20),(Invar%typat_unitcell(jj),jj=1,Invar%natom_unitcell)
1249 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('ntypat',20),Invar%ntypat
1250 [ + - + + ]: 103 : write(unt,'(1x,a20,20(1x,f15.10))') ljust('amu',20),(Invar%amu(jj),jj=1,Invar%ntypat)
1251 [ + - ]: 46 : if (allocated(Invar%znucl)) then
1252 [ + - + + ]: 103 : write(unt,'(1x,a20,20(1x,f15.10))') ljust('znucl',20),(Invar%znucl(jj),jj=1,Invar%ntypat)
1253 : : end if
1254 : :
1255 : 46 : write(unt,'(a)') ' ======================= Define the supercell ================================'
1256 : 46 : write(unt,'(1x,a20)') ljust('rprimd',20)
1257 [ + + ]: 184 : do ii=1,3
1258 [ + - + + ]: 598 : write(unt,'(22x,3(f15.10,1x))') (Invar%rprimd_md(ii,jj),jj=1,3)
1259 : : end do
1260 : 46 : write(unt,'(1x,a20)') ljust('multiplicity',20)
1261 [ + + ]: 184 : do ii=1,3
1262 [ + - + + ]: 598 : write(unt,'(22x,3(f15.10,1x))') (Invar%multiplicity(ii,jj),jj=1,3)
1263 : : end do
1264 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('natom',20),Invar%natom
1265 : 46 : write(unt,'(1x,a20)') ljust('typat',20)
1266 [ - + + + ]: 792 : do ii=1,Invar%natom,10
1267 [ + + ]: 792 : if (ii+9.lt.Invar%natom) then
1268 [ + - + + ]: 7700 : write(unt,'(22x,10(i4,1x))') (Invar%typat(ii+jj-1),jj=1,10)
1269 : : else
1270 [ + - + + ]: 324 : write(unt,'(22x,10(i4,1x))') (Invar%typat(jj),jj=ii,Invar%natom)
1271 : : end if
1272 : : end do
1273 : :
1274 : 46 : write(unt,'(a)') ' ======================= Define computational details ========================'
1275 : 46 : write(unt,'(1x,a20,1x,i5)') ljust('nstep_max',20),Invar%nstep_max
1276 : 46 : write(unt,'(1x,a20,1x,i5)') ljust('nstep_min',20),Invar%nstep_min
1277 : 46 : write(unt,'(1x,a20,1x,i5)') ljust('nstep_tot',20),Invar%nstep_tot
1278 : 46 : write(unt,'(1x,a20,1x,f15.10)') ljust('rcut',20),Invar%rcut
1279 : 46 : write(unt,'(1x,a20,1x,f15.10)') ljust('temperature',20),Invar%temperature
1280 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('order',20),Invar%order
1281 [ + + ]: 46 : if (Invar%order.eq.3) then
1282 : 7 : write(unt,'(1x,a20,1x,f15.10)') ljust('rcut3',20),Invar%rcut3
1283 [ + + ]: 39 : else if (Invar%order.eq.4) then
1284 : 9 : write(unt,'(1x,a20,2(1x,f15.10))') ljust('rcut3/rcut4',20),Invar%rcut3,Invar%rcut4
1285 : : end if
1286 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('prtvol',20),Invar%prtvol
1287 : :
1288 : 46 : write(unt,'(a)') ' ======================= Optional input variables ============================'
1289 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('prtphdos',20),Invar%prtphdos
1290 : 46 : write(unt,'(1x,a20,1x,f15.10)') ljust('dosdeltae',20),Invar%dosdeltae
1291 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('use_ideal_positions',20),Invar%use_ideal_positions
1292 [ + + ]: 87 : write(unt,'(1x,a20,1x,i4)') ljust('dipdip',20),merge(1,0,Invar%dipdip)
1293 [ + + ]: 46 : if (Invar%dipdip) then
1294 : 5 : write(unt,'(1x,a20)') ljust('born_charge',20)
1295 [ + + ]: 16 : do kk=1,Invar%natom_unitcell
1296 : 11 : write(unt,'(24x,a,1x,i4)') 'Atom',kk
1297 [ + + ]: 49 : do ii=1,3
1298 [ + - + + ]: 143 : write(unt,'(22x,3(f15.10,1x))') (Invar%zeff(ii,jj,kk),jj=1,3)
1299 : : end do
1300 : : end do
1301 : 5 : write(unt,'(1x,a20)') ljust('dielec_constant',20)
1302 [ + + ]: 20 : do ii=1,3
1303 [ + - + + ]: 65 : write(unt,'(22x,3(f15.10,1x))') (Invar%dielt(ii,jj),jj=1,3)
1304 : : end do
1305 : 5 : write(unt,'(1x,a20,1x,i4)') ljust('symbec',20),Invar%symbec
1306 : : end if
1307 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('bzpath',20),Invar%bzpath
1308 [ - + ]: 46 : if (Invar%bzpath.lt.0) then
1309 : 0 : write(unt,'(a)') ' Q points as given in the input file:'
1310 [ # # ]: 0 : do jj=1,abs(Invar%bzpath)
1311 : 0 : write(unt,'(22x,3(f15.10,1x))') Invar%qpt(:,jj)
1312 : : end do
1313 [ + + ]: 46 : else if (Invar%bzpath.gt.0) then
1314 : 10 : write(unt,'(a,1x,a2,10("-",a2))') ' Special q-points: ',Invar%special_qpt(:)
1315 : : end if
1316 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('ndivsm',20),Invar%ndivsm
1317 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('slice',20),Invar%slice
1318 [ + + ]: 46 : if (Invar%enunit.eq.0) write(unt,'(1x,a20,1x,i4,1x,a)') ljust('enunit',20),Invar%enunit,'(Phonon frequencies in meV)'
1319 [ + + ]: 46 : if (Invar%enunit.eq.1) write(unt,'(1x,a20,1x,i4,1x,a)') ljust('enunit',20),Invar%enunit,'(Phonon frequencies in cm-1)'
1320 [ + + ]: 46 : if (Invar%enunit.eq.2) write(unt,'(1x,a20,1x,i4,1x,a)') ljust('enunit',20),Invar%enunit,'(Phonon frequencies in mHa)'
1321 [ + + ]: 46 : if (Invar%enunit.eq.3) write(unt,'(1x,a20,1x,i4,1x,a)') ljust('enunit',20),Invar%enunit,'(Phonon frequencies in THz)'
1322 : 46 : write(unt,'(1x,a20,1x,i4,1x,i4)') ljust('nproc',20),Invar%nproc(1),Invar%nproc(2)
1323 [ + + ]: 46 : if (Invar%readifc.eq.1) then
1324 : 1 : write(unt,'(1x,a20,1x,i4,1x,f15.10)') ljust('readifc',20),Invar%readifc,Invar%tolread
1325 : : else
1326 : 45 : write(unt,'(1x,a20,1x,i4)') ljust('readifc',20),Invar%readifc
1327 : : end if
1328 : 46 : write(unt,'(1x,a20,1x,3(i4,1x))') ljust('alloy',20),Invar%alloy,Invar%ityp_alloy1,Invar%ityp_alloy2
1329 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('together',20),Invar%together
1330 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('bzlength',20),Invar%bzlength
1331 [ + + ]: 46 : if (Invar%bzlength.gt.0) then
1332 : 2 : write(unt,'(a,1x,i3,10("-",i3))') ' Length of BZ : ',Invar%lgth_segments(:)
1333 : : end if
1334 [ + + + + : 145 : if (any(Invar%ngqpt1(:) > 0)) then
+ + ]
1335 : 13 : write(unt,'(1x,a20,1x,3(i4,1x))') ljust('ngqpt1',20),Invar%ngqpt1(:)
1336 : : end if
1337 : 46 : write(unt,'(1x,a20,1x,3(f6.3,1x))') ljust('q1shft',20),Invar%q1shft(:)
1338 : 46 : write(unt,'(1x,a20,1x,3(i4,1x))') ljust('ngqpt2',20),Invar%ngqpt2(:)
1339 : 46 : write(unt,'(1x,a20,1x,3(f6.3,1x))') ljust('q2shft',20),Invar%q2shft(:)
1340 : 46 : write(unt,'(1x,a20,1x,i4)') ljust('use_weights',20),Invar%use_weights
1341 : :
1342 : 46 : end subroutine tdep_write_input
1343 : : !!***
1344 : :
1345 : : !====================================================================================================
1346 : :
1347 : : !!****f* ABINIT/m_tdep_dataset/tdep_destroy_invar
1348 : : !! NAME
1349 : : !! tdep_destroy_invar
1350 : : !!
1351 : : !! FUNCTION
1352 : : !! Free the memory allocated in the Invar object.
1353 : : !!
1354 : : !! SIDE EFFECTS
1355 : : !! Invar = Input object. All its allocated arrays are deallocated.
1356 : : !!
1357 : : !! SOURCE
1358 : :
1359 : 46 : subroutine tdep_destroy_invar(Invar)
1360 : :
1361 : : class(atdep_dataset_type), intent(inout) :: Invar
1362 : :
1363 [ - + ]: 46 : ABI_FREE(Invar%amu)
1364 [ - + ]: 46 : ABI_FREE(Invar%typat)
1365 [ - + ]: 46 : ABI_FREE(Invar%xred_unitcell)
1366 [ - + ]: 46 : ABI_FREE(Invar%xred_unitcell_original)
1367 [ - + ]: 46 : ABI_FREE(Invar%typat_unitcell)
1368 [ - + ]: 46 : if (Invar%bzpath.lt.0) then
1369 [ # # ]: 0 : ABI_FREE(Invar%qpt)
1370 [ + + ]: 46 : else if (Invar%bzpath.gt.0) then
1371 [ - + ]: 10 : ABI_FREE(Invar%special_qpt)
1372 : : end if
1373 [ + + ]: 46 : if (Invar%bzlength.gt.0) then
1374 [ - + ]: 2 : ABI_FREE(Invar%lgth_segments)
1375 : : end if
1376 [ + - ]: 46 : ABI_SFREE(Invar%zeff)
1377 [ + - ]: 46 : ABI_SFREE(Invar%znucl)
1378 : :
1379 : 46 : end subroutine tdep_destroy_invar
1380 : : !!***
1381 : :
1382 : : !====================================================================================================
1383 : :
1384 : : !!****f* ABINIT/m_tdep_dataset/tdep_init_MPIdata
1385 : : !! NAME
1386 : : !! tdep_init_MPIdata
1387 : : !!
1388 : : !! FUNCTION
1389 : : !! Initialize the MPIdata object: create the 2D cartesian grid of processors
1390 : : !! (nproc(1) processors for the shells times nproc(2) processors for the
1391 : : !! MD steps) and the corresponding communicators. If the product of
1392 : : !! nproc(1) and nproc(2) does not match the total number of processors,
1393 : : !! the parallelization is performed over the MD steps only.
1394 : : !! Then, distribute the nstep_tot MD steps over the processors, and write
1395 : : !! this distribution in the output file.
1396 : : !!
1397 : : !! OUTPUT
1398 : : !! MPIdata = Info on MPI parallelism, initialized.
1399 : : !!
1400 : : !! SIDE EFFECTS
1401 : : !! Invar = Input object containing the input variables.
1402 : : !! On exit, Invar%my_nstep is set.
1403 : : !!
1404 : : !! SOURCE
1405 : :
1406 : 46 : subroutine tdep_init_MPIdata(MPIdata,Invar)
1407 : :
1408 : : class(MPI_enreg_type), intent(out) :: MPIdata
1409 : : type(atdep_dataset_type), intent(inout) :: Invar
1410 : : integer :: ii,remain,ierr,iproc,istep
1411 : 46 : integer, allocatable :: tab_step(:)
1412 : : character(len=500) :: message
1413 : :
1414 : : #if defined HAVE_MPI
1415 : : integer :: dimcart,commcart_2d,me_cart_2d
1416 : : logical :: reorder
1417 : 46 : integer,allocatable :: coords(:),sizecart(:)
1418 [ + - ]: 46 : logical,allocatable :: periode(:), keepdim(:)
1419 : : #endif
1420 : :
1421 : : ! Check the number of processors
1422 : 46 : MPIdata%nproc_shell=Invar%nproc(1)
1423 : 46 : MPIdata%nproc_step =Invar%nproc(2)
1424 : 46 : MPIdata%nproc = xmpi_comm_size(xmpi_world)
1425 [ - + ]: 46 : if (MPIdata%nproc_step*MPIdata%nproc_shell.ne.MPIdata%nproc) then
1426 : 0 : ABI_WARNING('The parallelization is performed over steps')
1427 : 0 : MPIdata%nproc_step = xmpi_comm_size(xmpi_world)
1428 : : end if
1429 : :
1430 : 46 : MPIdata%master = 0
1431 : 46 : MPIdata%iam_master =.false.
1432 : : ! Initialize the MPIdata datastructure for sequential calculation
1433 [ + - ]: 46 : if (MPIdata%nproc.eq.1) then
1434 : 46 : MPIdata%comm_shell = xmpi_comm_null
1435 : 46 : MPIdata%comm_step = xmpi_comm_null
1436 : 46 : MPIdata%comm_shellstep = xmpi_comm_null
1437 : 46 : MPIdata%me_shell = 0
1438 : 46 : MPIdata%me_step = 0
1439 : 46 : MPIdata%my_nstep = Invar%nstep_tot
1440 : 46 : MPIdata%iam_master = (MPIdata%me_step == MPIdata%master)
1441 : : end if
1442 : :
1443 : : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1444 : : !!!!!!!!!!!!!!!!!!!!!! Parallel calculation !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1445 : : !!!!!!!!!!!!!! Definition of the processor grid !!!!!!!!!!!!!!!!!!!!!!!!!!!!
1446 : : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1447 : : #if defined HAVE_MPI
1448 : : !FB if (MPIdata%nproc.eq.1) return
1449 : :
1450 : : ! Create the global cartesian 2D-communicator
1451 : 46 : dimcart=2
1452 [ - + ]: 46 : ABI_MALLOC(sizecart,(dimcart))
1453 [ - + ]: 46 : ABI_MALLOC(periode,(dimcart))
1454 : 46 : sizecart(1)=MPIdata%nproc_shell ! MPIdata%nproc_shell
1455 : 46 : sizecart(2)=MPIdata%nproc_step ! MPIdata%nproc_step
1456 [ + + ]: 138 : periode(:)=.false.;reorder=.false.
1457 : 46 : call MPI_CART_CREATE(xmpi_world,dimcart,sizecart,periode,reorder,commcart_2d,ierr)
1458 [ - + ]: 46 : ABI_FREE(periode)
1459 : 46 : ABI_FREE(sizecart)
1460 : :
1461 : : ! Find the index and coordinates of the current processor
1462 : 46 : call MPI_COMM_RANK(commcart_2d,me_cart_2d,ierr)
1463 [ - + + - : 138 : ABI_MALLOC(coords,(dimcart))
+ - + - -
+ - + ]
1464 : 46 : call MPI_CART_COORDS(commcart_2d,me_cart_2d,dimcart,coords,ierr)
1465 : 46 : MPIdata%me_shell=coords(1)
1466 : 46 : MPIdata%me_step =coords(2)
1467 : 46 : ABI_FREE(coords)
1468 [ + - + - ]: 46 : if ((MPIdata%me_shell == MPIdata%master).and.(MPIdata%me_step == MPIdata%master)) then
1469 : 46 : MPIdata%iam_master = .true.
1470 : : end if
1471 : :
1472 [ - + + - : 138 : ABI_MALLOC(keepdim,(dimcart))
+ - + - -
+ - + ]
1473 : : ! Create the communicator for shell distribution
1474 : 46 : keepdim(1)=.true.
1475 : 46 : keepdim(2)=.false.
1476 : 46 : call MPI_CART_SUB(commcart_2d,keepdim,MPIdata%comm_shell,ierr)
1477 : : ! Create the communicator for step distribution
1478 : 46 : keepdim(1)=.false.
1479 : 46 : keepdim(2)=.true.
1480 : 46 : call MPI_CART_SUB(commcart_2d,keepdim,MPIdata%comm_step,ierr)
1481 : : ! Create the communicator for shellstep distribution
1482 : 46 : keepdim(1)=.true.
1483 : 46 : keepdim(2)=.true.
1484 : 46 : call MPI_CART_SUB(commcart_2d,keepdim,MPIdata%comm_shellstep,ierr)
1485 : 46 : ABI_FREE(keepdim)
1486 : 46 : call xmpi_comm_free(commcart_2d)
1487 : :
1488 : : ! Write some data
1489 : 46 : write(message,'(a1,a20,2(1x,i4))')'-', 'Number of processors',MPIdata%nproc_shell,MPIdata%nproc_step
1490 : 46 : call wrtout(ab_out,message,'COLL')
1491 : : !FB write(message,'(a,2i5)') 'me_shell and me_step : ',MPIdata%me_shell,MPIdata%me_step
1492 : : !FB call wrtout(ab_out,message,'COLL')
1493 : : #endif
1494 : :
1495 : : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1496 : : !!!!!!!!!!!!!! Distribution over STEP processors !!!!!!!!!!!!!!!!!!!!!!!!!!!
1497 : : !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1498 : 46 : MPIdata%my_nstep =int(Invar%nstep_tot/MPIdata%nproc_step)
1499 : 46 : remain=Invar%nstep_tot-MPIdata%nproc_step*MPIdata%my_nstep
1500 [ - + ]: 46 : do ii=1,remain
1501 [ - - ]: 46 : if ((ii-1).eq.MPIdata%me_step) MPIdata%my_nstep=MPIdata%my_nstep+1
1502 : : end do
1503 [ - + + - : 184 : ABI_MALLOC(MPIdata%nstep_all,(MPIdata%nproc_step)); MPIdata%nstep_all(:)=zero
+ - + - -
+ - + - +
+ + ]
1504 : 46 : call xmpi_allgather(MPIdata%my_nstep,MPIdata%nstep_all,MPIdata%comm_step,ierr)
1505 : 46 : write(ab_out,'(a)') ' '
1506 : 46 : write(ab_out,'(a,1x,i4)') ' All quantities are computed from nstep_min=',Invar%nstep_min
1507 : 46 : write(ab_out,'(a,1x,i4)') ' to nstep_max=',Invar%nstep_max
1508 [ + + ]: 46 : if (Invar%slice.ne.1) then
1509 : 2 : write(ab_out,'(a,1x,i4)') ' by using a slice=',Invar%slice
1510 : : end if
1511 : 46 : write(ab_out,'(a,1x,i4)') ' So, the real number of time steps is nstep=',Invar%nstep_tot
1512 [ - + ]: 46 : if (MPIdata%nproc_step.gt.1) then
1513 : 0 : write(ab_out,'(a,1000(1x,i5))') '-Distribution of number of steps wrt the number of processors=',MPIdata%nstep_all(:)
1514 : : end if
1515 : :
1516 [ + - - + : 230 : ABI_MALLOC(MPIdata%nstep_acc,(MPIdata%nproc_step+1)); MPIdata%nstep_acc(:)=zero
+ - + - +
- - + - +
- + + + ]
1517 : 46 : MPIdata%nstep_acc(1)=0
1518 [ + + ]: 92 : do ii=2,MPIdata%nproc_step+1
1519 : 92 : MPIdata%nstep_acc(ii)=MPIdata%nstep_acc(ii-1)+MPIdata%nstep_all(ii-1)
1520 : : end do
1521 [ - + ]: 46 : if (MPIdata%nstep_acc(MPIdata%nproc_step+1).ne.Invar%nstep_tot) then
1522 : 0 : write(std_out,*) 'STOP : pb in nstep_acc'
1523 : 0 : stop
1524 : : end if
1525 : :
1526 [ - + + - : 1130 : ABI_MALLOC(tab_step,(Invar%nstep_tot)); tab_step(:)=zero
+ - + - -
+ - + +
+ ]
1527 [ - + + - : 1084 : ABI_MALLOC(MPIdata%my_step ,(Invar%nstep_tot)); MPIdata%my_step (:)=.false.
- + - + -
+ + + ]
1528 [ + + ]: 92 : do iproc=1,MPIdata%nproc_step
1529 [ + + ]: 1084 : do istep=1,Invar%nstep_tot
1530 [ + - + - ]: 1038 : if ((istep.gt.MPIdata%nstep_acc(iproc)).and.(istep.le.MPIdata%nstep_acc(iproc+1))) then
1531 : 992 : tab_step(istep)=iproc-1
1532 : : end if
1533 : : end do
1534 : : end do
1535 [ + + ]: 1038 : do istep=1,Invar%nstep_tot
1536 : 1038 : MPIdata%my_step(istep) = (tab_step(istep) == MPIdata%me_step)
1537 : : end do
1538 : :
1539 [ - + + - : 184 : ABI_MALLOC(MPIdata%shft_step,(MPIdata%nproc_step)); MPIdata%shft_step(:)=zero
+ - + - -
+ - + - +
+ + ]
1540 : 46 : MPIdata%shft_step(1)=0
1541 [ - + ]: 46 : do ii=2,MPIdata%nproc_step
1542 : 46 : MPIdata%shft_step(ii)=MPIdata%shft_step(ii-1)+MPIdata%nstep_all(ii-1)
1543 : : end do
1544 [ - + ]: 46 : ABI_FREE(tab_step)
1545 : :
1546 : 46 : Invar%my_nstep=MPIdata%my_nstep
1547 : :
1548 : 46 : end subroutine tdep_init_MPIdata
1549 : : !!***
1550 : :
1551 : : !====================================================================================================
1552 : :
1553 : : !!****f* ABINIT/m_tdep_dataset/tdep_destroy_mpidata
1554 : : !! NAME
1555 : : !! tdep_destroy_mpidata
1556 : : !!
1557 : : !! FUNCTION
1558 : : !! Free the memory allocated in the MPIdata object.
1559 : : !!
1560 : : !! SIDE EFFECTS
1561 : : !! MPIdata = Info on MPI parallelism. Its allocated arrays are deallocated.
1562 : : !!
1563 : : !! SOURCE
1564 : :
1565 : 46 : subroutine tdep_destroy_mpidata(MPIdata)
1566 : :
1567 : : class(MPI_enreg_type), intent(inout) :: MPIdata
1568 : :
1569 [ - + ]: 46 : ABI_FREE(MPIdata%shft_step)
1570 [ - + ]: 46 : ABI_FREE(MPIdata%nstep_acc)
1571 [ - + ]: 46 : ABI_FREE(MPIdata%nstep_all)
1572 [ - + ]: 46 : ABI_FREE(MPIdata%my_step)
1573 : :
1574 : 46 : end subroutine tdep_destroy_mpidata
1575 : : !!***
1576 : :
1577 : : !====================================================================================================
1578 : :
1579 [ - + - - : 138 : end module m_tdep_dataset
+ + - + +
- - + - +
+ - - + -
+ - + - +
# # # # #
# # # ]
1580 : : !!***
|
