Branch data Line data Source code
1 : : !!****m* ABINIT/m_tdep_sampling
2 : : !! NAME
3 : : !! m_tdep_sampling
4 : : !!
5 : : !! FUNCTION
6 : : !! This module contains the TDEP Sampling data type
7 : : !! which holds the set of configurations from which the IFC will be fit.
8 : : !!
9 : : !! COPYRIGHT
10 : : !! Copyright (C) 2011-2026 ABINIT group (GA,FB,JB)
11 : : !! This file is distributed under the terms of the
12 : : !! GNU General Public License, see ~abinit/COPYING
13 : : !! or http://www.gnu.org/copyleft/gpl.txt .
14 : : !!
15 : : !! SOURCE
16 : :
17 : : #if defined HAVE_CONFIG_H
18 : : #include "config.h"
19 : : #endif
20 : :
21 : : #include "abi_common.h"
22 : :
23 : : module m_tdep_sampling
24 : :
25 : : use defs_basis
26 : : use m_errors
27 : : use m_abicore
28 : : use m_xmpi
29 : : use m_io_tools
30 : : use m_abihist, only : abihist
31 : : use m_tdep_dataset, only : atdep_dataset_type, MPI_enreg_type
32 : : use m_tdep_latt, only : Lattice_type, tdep_make_inbox
33 : : use m_tdep_sym, only : Symmetries_type, tdep_SearchS_1at
34 : :
35 : : implicit none
36 : :
37 : : type tdep_Sampling_type
38 : :
39 : : integer :: natom
40 : : ! Number of atoms in the supercell
41 : :
42 : : integer :: natom_unitcell
43 : : ! Number of atoms in the unitcell
44 : :
45 : : integer :: my_nstep
46 : : ! Number of MD steps held locally
47 : :
48 : : integer :: nstep_tot
49 : : ! Total of MD steps
50 : :
51 : : integer, allocatable :: typat_unitcell(:)
52 : : ! typat_unitcell(natom_unitcell)
53 : : ! Atom type in the unitcell.
54 : :
55 : : integer, allocatable :: typat(:)
56 : : ! typat(natom)
57 : : ! Atom type in the supercell.
58 : :
59 : : double precision, allocatable :: xred_unitcell(:,:)
60 : : ! xred_unitcell(3, natom_unitcell)
61 : : ! Reduced equilibrium position of the atoms in the unitcell.
62 : :
63 : : double precision, allocatable :: xred_ideal(:,:)
64 : : ! xred_ideal(3, natom)
65 : : ! Reduced equilibrium position of the atoms in the supercell.
66 : :
67 : : double precision, allocatable :: xred(:,:,:)
68 : : ! xred(3, natom, my_nstep)
69 : : ! Reduced position of the atom at each step.
70 : :
71 : : double precision, allocatable :: fcart(:,:,:)
72 : : ! fcart(3, natom, my_nstep)
73 : : ! Cartesian forces at each step.
74 : :
75 : : double precision, allocatable :: etot(:)
76 : : ! etot(my_nstep)
77 : : ! Total energy at each step.
78 : :
79 : : double precision, allocatable :: weights(:)
80 : : ! weights(my_nstep)
81 : : ! The weight of each configuration for the fitting.
82 : : ! By default, these will be 1 / nstep_tot.
83 : :
84 : : double precision, allocatable :: Rlatt_cart(:,:,:)
85 : : ! Rlatt_cart(3, natom_unitcell, natom)
86 : : ! Cartesian coordinate of the lattice vectors of the unitcell
87 : : ! within the supercell, for each atom.
88 : :
89 : : double precision, allocatable :: Rlatt_scaled(:,:,:)
90 : : ! Rlatt_scaled(3, natom_unitcell, natom)
91 : : ! Rlatt_cart divided by acell_unitcell.
92 : :
93 : : double precision, allocatable :: ucart(:,:,:)
94 : : ! ucart(3, natom, my_nstep)
95 : : ! Cartesian displacements of the atoms with respect to their equilibrium
96 : : ! positions at each time step.
97 : :
98 : : double precision, allocatable :: distance(:,:,:)
99 : : ! distance(natom, natom, 4)
100 : : ! Distance between the ideal positions of the atoms in the supercell,
101 : : ! (norm, and cartesian components).
102 : :
103 : : double precision, allocatable :: Forces(:)
104 : : ! Forces(3*natom*my_nstep)
105 : : ! The cartesian forces for all configurations, as a flat array.
106 : :
107 : : end type tdep_Sampling_type
108 : :
109 : : public :: tdep_sampling_init_read
110 : : public :: tdep_sampling_free
111 : : public :: tdep_sampling_shift_xred
112 : : public :: tdep_sampling_rotate
113 : : public :: tdep_MatchIdeal2Average
114 : : public :: tdep_write_xred_average
115 : :
116 : : contains
117 : :
118 : : !=====================================================================================================
119 : :
120 : 44 : subroutine tdep_sampling_init_read(MD,Invar,MPIdata,Hist)
121 : :
122 : : type(tdep_Sampling_type), intent(inout) :: MD
123 : : type(atdep_dataset_type), intent(in) :: Invar
124 : : type(MPI_enreg_type), intent(in) :: MPIdata
125 : : type(abihist), intent(in) :: Hist
126 : :
127 : : integer :: this_istep,istep,iatom,jstep
128 : : double precision :: tmp1,tmp2,tmp3
129 : :
130 : 44 : MD%natom = Invar%natom
131 : 44 : MD%natom_unitcell = Invar%natom_unitcell
132 : 44 : MD%my_nstep = Invar%my_nstep
133 : 44 : MD%nstep_tot = Invar%nstep_tot
134 : :
135 [ - + + - : 238 : ABI_CALLOC(MD%typat_unitcell, (MD%natom_unitcell))
+ - + - -
+ - + - +
+ + ]
136 [ - + + - : 7014 : ABI_CALLOC(MD%typat, (MD%natom))
+ - + - -
+ - + - +
+ + ]
137 [ - + + - : 556 : ABI_CALLOC(MD%xred_unitcell, (3,MD%natom_unitcell))
+ - + - -
+ - + - +
+ + + + ]
138 [ - + + - : 27660 : ABI_CALLOC(MD%xred_ideal, (3,MD%natom))
+ - + - -
+ - + - +
+ + + + ]
139 [ - + + - : 557736 : ABI_CALLOC(MD%xred, (3,MD%natom,MD%my_nstep))
- + + - +
- + - - +
- + - + +
+ + + +
+ ]
140 [ - + + - : 557692 : ABI_CALLOC(MD%fcart, (3,MD%natom,MD%my_nstep))
- + + - -
+ - + - +
+ + + + +
+ ]
141 [ - + + - : 1084 : ABI_CALLOC(MD%etot, (MD%my_nstep))
+ - + - -
+ - + - +
+ + ]
142 [ - + + - : 1040 : ABI_CALLOC(MD%weights, (MD%my_nstep))
- + - + -
+ + + ]
143 [ - + + - : 5917156 : ABI_CALLOC(MD%distance, (MD%natom,MD%natom,4))
- + + - +
- + - + -
- + - + -
+ + + + +
+ + ]
144 [ - + + - : 66602 : ABI_CALLOC(MD%Rlatt_scaled, (3,MD%natom_unitcell,MD%natom))
- + + - +
- + - - +
- + - + +
+ + + +
+ ]
145 [ - + + - : 66558 : ABI_CALLOC(MD%Rlatt_cart, (3,MD%natom_unitcell,MD%natom))
- + + - -
+ - + - +
+ + + + +
+ ]
146 [ - + + - : 557692 : ABI_CALLOC(MD%ucart, (3,MD%natom,MD%my_nstep))
- + + - -
+ - + - +
+ + + + +
+ ]
147 [ + - - + : 417588 : ABI_CALLOC(MD%Forces, (3*MD%natom*MD%my_nstep))
+ - + - +
- - + - +
- + + + ]
148 : :
149 [ + + ]: 150 : MD%typat_unitcell(:) = Invar%typat_unitcell(:)
150 [ + + + + ]: 468 : MD%xred_unitcell(:,:) = Invar%xred_unitcell(:,:)
151 [ + + ]: 6926 : MD%typat(:) = Invar%typat(:)
152 : :
153 : : ! Read xred.dat, fcart.dat and etot.dat ASCII files or extract them from the HIST.nc netcdf file.
154 : 44 : write(Invar%stdout,'(a)') ' '
155 : 44 : this_istep=0
156 : 44 : jstep=0
157 [ + + ]: 44 : if (Invar%use_weights.eq.1) then
158 : 1 : open(unit=30,file=trim(Invar%input_prefix)//'_weights.dat')
159 [ + - ]: 43 : else if (Invar%use_weights.eq.0) then
160 [ + + ]: 975 : MD%weights=1.0d0/real(MD%nstep_tot)
161 : : endif
162 [ + + ]: 44 : if (Invar%netcdf) then
163 [ + + ]: 217 : do istep=Invar%nstep_min,Invar%nstep_max
164 [ + + ]: 217 : if (mod(istep-Invar%nstep_min,Invar%slice).ne.0) then
165 : : cycle
166 : : else
167 : 201 : jstep=jstep+1
168 [ - + ]: 201 : if (.not.MPIdata%my_step(jstep)) cycle
169 : 201 : this_istep=this_istep+1
170 [ + + + + ]: 62441 : MD%xred(:,:,this_istep) =Hist%xred (:,:,istep)
171 [ + + + + ]: 62441 : MD%fcart(:,:,this_istep)=Hist%fcart(:,:,istep)
172 : 211 : MD%etot(this_istep) =Hist%etot (istep)
173 : : end if
174 : : end do !istep
175 : 6 : write(Invar%stdout,'(a)') ' The positions, forces and energies are extracted from the NetCDF file: HIST.nc'
176 : : else
177 : 38 : open(unit=60,file=trim(Invar%input_prefix)//'_fcart.dat')
178 : 38 : open(unit=50,file=trim(Invar%input_prefix)//'_xred.dat')
179 : 38 : open(unit=40,file=trim(Invar%input_prefix)//'_etot.dat')
180 [ - + ]: 38 : do istep=1,Invar%nstep_min-1
181 [ # # ]: 0 : if (Invar%use_weights.eq.1) then
182 : 0 : read(30,*) tmp1
183 : : endif
184 : 0 : read(40,*) tmp1
185 [ - - ]: 38 : do iatom=1,MD%natom
186 : 0 : read(50,*) tmp1,tmp2,tmp3
187 : 0 : read(60,*) tmp1,tmp2,tmp3
188 : : end do
189 : : end do
190 [ + + ]: 799 : do istep=Invar%nstep_min,Invar%nstep_max
191 [ + + ]: 799 : if (mod(istep-Invar%nstep_min,Invar%slice).ne.0) then
192 [ - + ]: 10 : if (Invar%use_weights.eq.1) then
193 : 0 : read(30,*) tmp1
194 : : endif
195 : 10 : read(40,*) tmp1
196 [ + + ]: 1290 : do iatom=1,Invar%natom
197 : 1280 : read(50,*) tmp1,tmp2,tmp3
198 : 1290 : read(60,*) tmp1,tmp2,tmp3
199 : : end do
200 : : else
201 : 751 : jstep=jstep+1
202 [ - + ]: 751 : if (.not.MPIdata%my_step(jstep)) then
203 [ # # ]: 0 : if (Invar%use_weights.eq.1) then
204 : 0 : read(30,*) tmp1
205 : : endif
206 : 0 : read(40,*) tmp1
207 [ # # ]: 0 : do iatom=1,Invar%natom
208 : 0 : read(50,*) tmp1,tmp2,tmp3
209 : 0 : read(60,*) tmp1,tmp2,tmp3
210 : : end do
211 : : else
212 : 751 : this_istep=this_istep+1
213 [ + + ]: 751 : if (Invar%use_weights.eq.1) then
214 : 20 : read(30,*) MD%weights(this_istep)
215 : : endif
216 : 751 : read(40,*) MD%etot(this_istep)
217 [ + + ]: 124343 : do iatom=1,MD%natom
218 : 123592 : read(50,*) MD%xred (1,iatom,this_istep),MD%xred (2,iatom,this_istep),MD%xred (3,iatom,this_istep)
219 : 124343 : read(60,*) MD%fcart(1,iatom,this_istep),MD%fcart(2,iatom,this_istep),MD%fcart(3,iatom,this_istep)
220 : : end do
221 : : end if !my_step
222 : : end if !slice
223 : : end do !istep
224 : 38 : close(40)
225 : 38 : close(50)
226 : 38 : close(60)
227 : 38 : write(Invar%stdout,'(2a)') ' The positions, forces and energies are extracted from the ASCII files:',&
228 : 76 : & ' xred.dat, fcart.dat & etot.dat'
229 : : end if !netcdf
230 [ + + ]: 44 : if (Invar%use_weights.eq.1) then
231 : 1 : close(30)
232 : : end if
233 : :
234 : 44 : end subroutine tdep_sampling_init_read
235 : :
236 : : !=====================================================================================================
237 : :
238 : 44 : subroutine tdep_sampling_free(MD)
239 : :
240 : : type(tdep_Sampling_type), intent(inout) :: MD
241 : :
242 [ - + ]: 44 : ABI_FREE(MD%typat_unitcell)
243 [ - + ]: 44 : ABI_FREE(MD%typat)
244 [ - + ]: 44 : ABI_FREE(MD%xred_unitcell)
245 [ - + ]: 44 : ABI_FREE(MD%xred_ideal)
246 [ - + ]: 44 : ABI_FREE(MD%xred)
247 [ - + ]: 44 : ABI_FREE(MD%fcart)
248 [ - + ]: 44 : ABI_FREE(MD%etot)
249 [ - + ]: 44 : ABI_FREE(MD%weights)
250 [ - + ]: 44 : ABI_FREE(MD%distance)
251 [ - + ]: 44 : ABI_FREE(MD%Rlatt_scaled)
252 [ - + ]: 44 : ABI_FREE(MD%Rlatt_cart)
253 [ - + ]: 44 : ABI_FREE(MD%ucart)
254 [ - + ]: 44 : ABI_FREE(MD%Forces)
255 : :
256 : 44 : end subroutine tdep_sampling_free
257 : :
258 : : !=====================================================================================================
259 : :
260 : : ! Shift xred to keep atoms in the same unit cell at each step.
261 : 44 : subroutine tdep_sampling_shift_xred(MD,MPIdata)
262 : :
263 : : type(tdep_Sampling_type), intent(inout) :: MD
264 : : type(MPI_enreg_type), intent(in) :: MPIdata
265 : : integer :: natom,ii,iatom,istep,ierr
266 : : integer :: shift,shift_max,shift_best
267 : : double precision :: xi, dist, best_dist
268 : 44 : double precision, allocatable :: x0(:,:)
269 : :
270 : 44 : natom = MD%natom
271 [ - + + - : 132 : ABI_MALLOC(x0,(3,natom))
+ - + - -
+ - + ]
272 : :
273 : : ! Communicate xred at the first step
274 [ + + + + ]: 27572 : x0(:,:) = zero
275 [ + - ]: 44 : if (MPIdata%my_step(1)) then
276 [ + + + + ]: 27572 : x0(:,:) = MD%xred(:,:,1)
277 : : end if
278 : 44 : call xmpi_sum(x0,MPIdata%comm_step,ierr)
279 : :
280 : : ! Shift xred from all steps in the same unitcell as the first step
281 : 44 : shift_max = 1
282 [ + + ]: 996 : do istep=1, MD%my_nstep
283 [ + + ]: 140148 : do iatom=1,natom
284 [ + + ]: 557560 : do ii=1,3
285 : 417456 : best_dist = abs(MD%xred(ii,iatom,istep) - x0(ii,iatom))
286 : 417456 : shift_best = 0
287 [ + + ]: 1669824 : do shift=-shift_max,shift_max
288 : 1252368 : xi = MD%xred(ii,iatom,istep) + shift
289 : 1252368 : dist = abs(xi - x0(ii,iatom))
290 [ + + ]: 1669824 : if (dist < best_dist) then
291 : 366 : best_dist = dist
292 : 366 : shift_best = shift
293 : : end if
294 : : end do
295 : 556608 : MD%xred(ii,iatom,istep) = MD%xred(ii,iatom,istep) + shift_best
296 : : end do
297 : : end do
298 : : end do
299 : :
300 [ - + ]: 44 : ABI_FREE(x0)
301 : :
302 : 44 : end subroutine tdep_sampling_shift_xred
303 : :
304 : : !=====================================================================================================
305 : :
306 : 44 : subroutine tdep_sampling_rotate(MD,rotation_cart)
307 : :
308 : : type(tdep_Sampling_type), intent(inout) :: MD
309 : : integer :: iatom,istep
310 : : double precision :: rotation_cart(3,3)
311 : :
312 : : ! Apply rotation to fcart
313 [ + + ]: 996 : do istep=1,MD%my_nstep
314 [ + + ]: 140148 : do iatom=1,MD%natom
315 [ + + + + : 2783992 : MD%fcart(:,iatom,istep) = MATMUL(rotation_cart, MD%fcart(:,iatom,istep))
+ + + + +
- + - ]
316 : : end do
317 : : end do
318 : :
319 : 44 : end subroutine tdep_sampling_rotate
320 : :
321 : : !=====================================================================================================
322 : :
323 : : !!***
324 : : !!****f* ABINIT/m_tdep_sampling/tdep_MatchIdeal2Average
325 : : !! NAME
326 : : !! tdep_MatchIdeal2Average
327 : : !!
328 : : !! FUNCTION
329 : : !! Find the mapping between the atoms in the ideal (equilibrium) supercell,
330 : : !! and the atoms of the input moledular dynamics using their average positions.
331 : : !! Then compute the atom displacements with respect to the equilibrium positions
332 : : !! at each time step of the MD.
333 : : !!
334 : : !! INPUTS
335 : : !! MD = TDEP Sampling object containing the input positions and forces.
336 : : !! Invar = Input object containing the input variables.
337 : : !! Lattice = Lattice object describing the ideal structure.
338 : : !! Sym = Symetries object describing all the symmetry operations of the crystal.
339 : : !! MPIdata = Info on MPI parallelism.
340 : : !!
341 : : !! OUTPUT
342 : : !!
343 : : !! SIDE EFFECTS
344 : : !! The following quantities in MD are computed:
345 : : !!
346 : : !! distance = Distance between the ideal positions of the atoms in the supercell,
347 : : !! (norm, and cartesian components).
348 : : !! Forces = Cartesian forces on the atoms at each time steps, as a flat array.
349 : : !! ucart = Cartesian displacements of the atoms with respect to their equilibrium
350 : : !! positions at each time step.
351 : : !! Rlatt_cart = Cartesian coordinate of the lattice vectors of the unitcell
352 : : !! within the supercell, for each atom.
353 : : !! This array seems to have an extra dimension, for algorithmic simplicity.
354 : : !! Rlatt_scaled = Rlatt_cart divided by acell_unitcell.
355 : : !! These are used when reading an IFC file, to compare with the R vectors
356 : : !! that are stored in the file.
357 : : !!
358 : : !! Some of the reduced positions of the atoms MD%xred are shifted by a supercell
359 : : !! lattice vector in order to re-center the crystal.
360 : : !!
361 : : !! NOTES
362 : : !!
363 : : !! SOURCE
364 : :
365 : 44 : subroutine tdep_MatchIdeal2Average(MD,Invar,Lattice,Sym,MPIdata)
366 : :
367 : : type(tdep_Sampling_type),intent(inout) :: MD
368 : : type(atdep_dataset_type),intent(inout) :: Invar
369 : : type(Lattice_type),intent(in) :: Lattice
370 : : type(Symmetries_type),intent(inout) :: Sym
371 : : type(MPI_enreg_type),intent(in) :: MPIdata
372 : :
373 : : integer :: ii,jj,kk,max_ijk,iatcell,jatcell,iatom,jatom,eatom,fatom,istep
374 : : integer :: iatom_ref,ierr
375 : : integer :: ndir_match,natom_match
376 : : double precision :: tmp(3),tmp1(3),tmp2(3),Rlatt(3),xred_tmp(3),rprimd_md_tmp(3,3),distance_tmp(3)
377 : 44 : double precision, allocatable :: dist_unitcell(:,:,:),xcart_average(:,:)
378 : 44 : double precision, allocatable :: fcart_tmp(:,:,:),ucart_tmp(:,:,:)
379 : 44 : double precision, allocatable :: xred_average(:,:)
380 : 44 : double precision, allocatable :: xred_center(:,:)
381 : 44 : double precision, allocatable :: Rlatt_red (:,:,:)
382 : 44 : double precision, allocatable :: xred_ideal(:,:)
383 : : ! double precision, allocatable :: distance_average(:,:,:)
384 : 44 : integer, allocatable :: FromIdeal2Average(:)
385 : 44 : double precision, allocatable :: xcart(:,:,:)
386 : 44 : double precision, allocatable :: xcart_ideal(:,:)
387 : : logical :: ok,must_shift,discard_R
388 : : character(len=500) :: msg
389 : :
390 : 44 : write(Invar%stdout,*)' '
391 : 44 : write(Invar%stdout,*) '#############################################################################'
392 : 44 : write(Invar%stdout,*) '###### Find the matching between ideal and average positions ###############'
393 : 44 : write(Invar%stdout,*) '#############################################################################'
394 : :
395 : : !==========================================================================================
396 : : !======== 1/ Determine ideal positions and distances ======================================
397 : : !==========================================================================================
398 : 44 : write(Invar%stdout,*)' Determine ideal positions and distances...'
399 : : ! Define the bigbox with ideal positions
400 [ - + + - : 66602 : ABI_CALLOC(Rlatt_red ,(3,MD%natom_unitcell,MD%natom))
- + + - +
- + - - +
- + + + +
+ + + ]
401 [ - + + - : 27660 : ABI_CALLOC(xred_ideal,(3,MD%natom))
+ - + - -
+ - + + +
+ + ]
402 : 44 : max_ijk=20
403 : 44 : iatom=1
404 [ + + ]: 1848 : do ii=-max_ijk,max_ijk
405 [ + + ]: 75812 : do jj=-max_ijk,max_ijk
406 [ + + ]: 3108292 : do kk=-max_ijk,max_ijk
407 : :
408 : 3032524 : Rlatt(1)=real(ii-1)
409 : 3032524 : Rlatt(2)=real(jj-1)
410 : 3032524 : Rlatt(3)=real(kk-1)
411 : :
412 : 3032524 : discard_R = .false.
413 [ + + ]: 10412114 : do iatcell=1,MD%natom_unitcell
414 : :
415 [ + + ]: 7305626 : if (discard_R) cycle
416 : :
417 : : ! Compute the reduced positions
418 [ + + ]: 12141444 : tmp(:) = Rlatt(:) + MD%xred_unitcell(:,iatcell)
419 : 3035361 : call DGEMV('T',3,3,1.d0,Lattice%multiplicitym1(:,:),3,tmp(:),1,0.d0,xred_tmp(:),1)
420 : :
421 : : ! If the first atom of the pattern is in the [0;1[ range then keep all the
422 : : ! atoms of the pattern (even if the others are outside the box). Else,
423 : : ! none are taken.
424 [ + + ]: 3035361 : if (iatcell==1) then
425 [ - + + - : 27292716 : if (minval(xred_tmp(:)).lt.0.d0.or.maxval(xred_tmp(:)).ge.(1.d0-tol12)) then
+ + - + +
- + + + +
+ + ]
426 : : discard_R = .true.
427 : : cycle
428 : : end if
429 : : end if
430 : :
431 : : !GA: Why natom+1 ?
432 [ - + ]: 6882 : if (iatom.gt.(MD%natom+1)) then
433 : 0 : ABI_ERROR('The number of atoms found in the bigbox exceeds natom' )
434 : : end if
435 : :
436 [ + + ]: 27528 : xred_ideal(:,iatom) = xred_tmp(:)
437 : 6882 : call DGEMV('T',3,3,1.d0,Lattice%multiplicitym1(:,:),3,Rlatt(:),1,0.d0,Rlatt_red(:,1,iatom),1)
438 : 7309671 : iatom = iatom + 1
439 : : end do
440 : : end do
441 : : end do
442 : : end do
443 : :
444 [ - + ]: 44 : if (iatom.lt.MD%natom+1) then
445 : 0 : ABI_ERROR('The number of atoms found in the big box is smaller than natom')
446 : : end if
447 : :
448 : : ! Compute the distances between ideal positions in the SUPERcell
449 [ + + ]: 6926 : do eatom=1,MD%natom
450 [ + + ]: 1479234 : do fatom=1,MD%natom
451 [ + + ]: 5889232 : tmp(:)=xred_ideal(:,fatom)-xred_ideal(:,eatom)
452 : 1472308 : call tdep_make_inbox(tmp,1,1d-4)
453 [ + + + + ]: 19140004 : rprimd_md_tmp(:,:)=Lattice%rprimd_md(:,:)
454 [ + + ]: 5889232 : distance_tmp(:)=MD%distance(eatom,fatom,2:4)
455 : 1472308 : call DGEMV('T',3,3,1.d0,rprimd_md_tmp,3,tmp,1,0.d0,distance_tmp,1)
456 [ + + ]: 5889232 : MD%distance(eatom,fatom,2:4)=distance_tmp(:)
457 [ + + ]: 5889232 : do ii=1,3
458 : : ! Remove the rounding errors before writing (for non regression testing purposes)
459 [ + + ]: 4416924 : if (abs(MD%distance(eatom,fatom,ii+1)).lt.tol8) MD%distance(eatom,fatom,ii+1)=zero
460 : 5889232 : MD%distance(eatom,fatom,1)=MD%distance(eatom,fatom,1)+(MD%distance(eatom,fatom,ii+1))**2
461 : : end do
462 : 1472308 : MD%distance(eatom,fatom,1)=MD%distance(eatom,fatom,1)**0.5
463 : 1479190 : MD%distance(eatom,fatom,1)=tol12 * dint(MD%distance(eatom,fatom,1) / tol12)
464 : : end do
465 : : end do
466 : :
467 : : ! Compute the distances between ideal positions in the UNITcell
468 [ - + + - : 2212 : ABI_MALLOC(dist_unitcell,(MD%natom_unitcell,MD%natom_unitcell,3)); dist_unitcell(:,:,:)=zero
- + + - +
- + - + -
- + - + +
+ + + +
+ ]
469 [ + + ]: 150 : do iatcell=1,MD%natom_unitcell
470 [ + + ]: 664 : do jatcell=1,MD%natom_unitcell
471 [ + + ]: 2056 : tmp(:) = xred_ideal(:,jatcell)-xred_ideal(:,iatcell)
472 : 514 : call tdep_make_inbox(tmp,1,tol8)
473 [ + + ]: 2162 : dist_unitcell(iatcell,jatcell,:) = tmp(:)
474 : : end do
475 : : end do
476 : :
477 : : !==========================================================================================
478 : : !======== 2/ Find the matching between the ideal and average ==============================
479 : : !======== (from the MD simulations) positions. ==========================================
480 : : !==========================================================================================
481 : : ! NOTE: - xred_center is used to find the matching with the ideal positions
482 : : ! - xred_average is used to compute the displacements (from MD trajectories)
483 : : ! The difference between those two is that xred_center will be shifted to bring
484 : : ! one of the average positions at the origin, for an easier comparison with
485 : : ! xred_ideal. Some shifts by a supercell lattice vector will be computed
486 : : ! from the difference between xred_center and xred_ideal, and those shifts
487 : : ! will be applied to xred_average and xred at all steps.
488 : :
489 : 44 : write(Invar%stdout,*)' Compute average positions...'
490 [ - + + - : 27660 : ABI_CALLOC(xred_average,(3,MD%natom))
+ - + - -
+ - + + +
+ + ]
491 [ - + + - : 27616 : ABI_CALLOC(xred_center,(3,MD%natom))
- + - + +
+ + + ]
492 : : ! Average positions from MD (on nstep steps)
493 [ + + ]: 996 : do istep=1,MD%my_nstep
494 [ + + ]: 140148 : do iatom=1,MD%natom
495 [ + + ]: 557560 : xred_average(:,iatom)=xred_average(:,iatom)+MD%xred(:,iatom,istep)
496 : : end do
497 : : end do
498 : 44 : call xmpi_sum(xred_average,MPIdata%comm_step,ierr)
499 [ + + + + ]: 27572 : xred_average(:,:) = xred_average(:,:) / real(MD%nstep_tot)
500 : :
501 : : ! Search the basis of atoms in the supercell
502 : : ! in order to find iatom_ref
503 : 44 : write(Invar%stdout,*)' Search the unitcell basis of atoms in the MD trajectory...'
504 : 44 : ok=.false.
505 [ + + + + ]: 27572 : xred_center(:,:) = xred_average(:,:)
506 : 44 : iatcell=1
507 [ + - ]: 51 : do iatom=1,MD%natom
508 [ + - ]: 51 : if (MD%typat(iatom).ne.MD%typat_unitcell(iatcell)) cycle
509 : 51 : natom_match = 0
510 [ + + ]: 7797 : do jatom=1,MD%natom
511 : :
512 [ + + ]: 30984 : tmp(:)=xred_center(:,jatom)-xred_center(:,iatom)
513 : 7746 : call tdep_make_inbox(tmp,1,Invar%tolinbox)
514 : :
515 [ + + ]: 24087 : do jatcell=1,MD%natom_unitcell
516 [ + + ]: 16403 : if (MD%typat(jatom).ne.MD%typat_unitcell(jatcell)) cycle
517 : 14263 : ndir_match = 0
518 [ + + ]: 57052 : do ii=1,3
519 [ + + ]: 57052 : if (abs(tmp(ii)-dist_unitcell(iatcell,jatcell,ii)).le.Invar%tolmotif) then
520 : 5723 : ndir_match=ndir_match+1
521 : : end if
522 : : end do
523 [ + + ]: 21896 : if (ndir_match==3) then
524 : 113 : natom_match = natom_match + 1
525 : 113 : exit
526 : : end if
527 : : end do
528 : : end do
529 [ + + ]: 51 : if (natom_match.eq.MD%natom_unitcell) then
530 : : iatom_ref = iatom
531 : : ok=.true.
532 : : exit
533 [ - + ]: 7 : else if (natom_match.gt.MD%natom_unitcell) then
534 : 0 : write(msg,'(5a)') 'Too many atoms match the unit cell.',ch10,&
535 : 0 : 'Perhaps the value of tolmotif is too large,',ch10,&
536 : 0 : 'or the value of tolinbox is too small.'
537 : 0 : ABI_ERROR(msg)
538 : : endif
539 : : end do
540 : : if (.not.ok) then
541 : 0 : call tdep_write_xred_average(Invar,MPIdata,Lattice,xred_ideal,xred_center)
542 : 0 : write(msg,'(3a)') 'The basis of atoms written in input.in file does not appear in the MD trajectory.',ch10,&
543 : 0 : 'Perhaps, you can adjust the tolerance (tolmotif).'
544 : 0 : ABI_ERROR(msg)
545 : : end if
546 [ - + ]: 44 : ABI_FREE(dist_unitcell)
547 : :
548 : 44 : write(Invar%stdout,*)' Compare ideal and average positions using PBC...'
549 : : ! Modification of xred and Rlatt tabs
550 : : ! for averaged quantities: xred_center, xred_average, xred
551 : : ! 1/ The "iatom_ref" atom is put in (0.0;0.0;0.0)
552 [ + + ]: 176 : tmp(:) = xred_center(:,iatom_ref)
553 [ + + ]: 6926 : do jatom=1,MD%natom
554 [ + + ]: 27572 : xred_center(:,jatom) = xred_center(:,jatom) - tmp(:)
555 : : end do
556 : : ! 2/ All the atoms are put in the range [-0.5;0.5[ (use of PBC)
557 [ + + ]: 6926 : do jatom=1,MD%natom
558 [ + + ]: 27528 : tmp(:)=xred_center(:,jatom)
559 : 6882 : call tdep_make_inbox(tmp,1,Invar%tolinbox,xred_center(:,jatom))
560 : 6882 : call tdep_make_inbox(tmp,1,Invar%tolinbox,xred_average(:,jatom))
561 [ + + ]: 146078 : do istep=1,MD%my_nstep
562 : 146034 : call tdep_make_inbox(tmp,1,Invar%tolinbox,MD%xred(:,jatom,istep))
563 : : end do
564 : : end do
565 : : ! Modification of xred and Rlatt tabs
566 : : ! for ideal quantities: Rlatt_red et xred_ideal
567 : : ! 1/ The atom 1 is put in (0.0;0.0;0.0)
568 [ + + ]: 176 : tmp1(:)=xred_ideal(:,1)
569 [ + + ]: 176 : tmp2(:)=Rlatt_red(:,1,1)
570 [ + + ]: 6926 : do jatom=1,MD%natom
571 [ + + ]: 27528 : xred_ideal(:,jatom)= xred_ideal(:,jatom) -tmp1(:)
572 [ + + ]: 27572 : Rlatt_red (:,1,jatom)=Rlatt_red (:,1,jatom)-tmp2(:)
573 : : end do
574 : : ! 2/ All the atoms are put in the range [-0.5;0.5[ (use of PBC)
575 [ + + ]: 6926 : do jatom=1,MD%natom
576 [ + + ]: 27528 : tmp(:)=xred_ideal(:,jatom)
577 : 6882 : call tdep_make_inbox(tmp,1,tol8,xred_ideal(:,jatom))
578 : 6926 : call tdep_make_inbox(tmp,1,tol8,Rlatt_red(:,1,jatom))
579 : : !FB call tdep_make_inbox(Rlatt_red(:,1,jatom),1,tol8)
580 : : end do
581 : :
582 : : ! When the multiplicity equals 1 along one direction, there is some trouble
583 : : ! To clean!!!!!!!
584 [ + + ]: 176 : do ii=1,3
585 : : if ((Invar%multiplicity(ii,ii).eq.1).and.(Invar%multiplicity(ii,mod(ii ,3)+1).eq.0)&
586 [ + + + + : 176 : & .and.(Invar%multiplicity(ii,mod(ii+1,3)+1).eq.0)) then
+ - ]
587 [ + + ]: 65 : Rlatt_red(ii,1,:)=0.d0
588 : 1 : write(Invar%stdout,*) 'WARNING: multiplicity=1 for ii=',ii
589 : : end if
590 : : end do
591 : :
592 : : ! Define Rlatt for all the atoms in the basis (Rlatt_red varies as a function of iatcell)
593 [ + + ]: 44 : if (MD%natom_unitcell.gt.1) then
594 [ + + ]: 93 : do iatcell=2,MD%natom_unitcell
595 [ + + + + ]: 32109 : Rlatt_red(:,iatcell,:)=Rlatt_red(:,1,:)
596 : : end do
597 : : end if
598 [ + + ]: 6926 : do iatom=1,MD%natom
599 [ + + ]: 21812 : do iatcell=1,MD%natom_unitcell
600 [ + + ]: 59544 : tmp(:)=xred_ideal(:,iatom)-xred_ideal(:,iatcell)
601 : 21768 : call tdep_make_inbox(tmp,1,tol8,Rlatt_red(:,iatcell,iatom))
602 : : end do
603 : : end do
604 [ + + ]: 44 : if (Invar%debug) then
605 [ + + ]: 12 : do iatcell=1,MD%natom_unitcell
606 : 10 : write(Invar%stdout,*) 'For iatcell=',iatcell
607 [ + + ]: 1036 : do jatom=1,MD%natom
608 : 1034 : write(Invar%stdout,'(a,i4,a,3(f16.10,1x))') 'For jatom=',jatom,', Rlatt=',Rlatt_red(1:3,iatcell,jatom)
609 : : end do
610 : : end do
611 : : end if
612 : :
613 : : ! Matching between Ideal and Average positions: xred_ideal and xred_center
614 : : ! Then, write them in the xred_average.xyz file.
615 : 44 : write(Invar%stdout,*)' Write the xred_average.xyz file with ideal and average positions...'
616 [ - + + - : 7014 : ABI_CALLOC(FromIdeal2Average,(MD%natom))
+ - + - -
+ - + +
+ ]
617 [ + + ]: 6926 : do iatom=1,MD%natom
618 : 6882 : ok =.false.
619 [ + - ]: 739595 : do jatom=1,MD%natom
620 [ + + ]: 739595 : if (MD%typat(iatom).ne.MD%typat_unitcell(mod(jatom-1,MD%natom_unitcell)+1)) cycle
621 : 624459 : must_shift=.false.
622 : 624459 : ndir_match=0
623 [ + + ]: 2497836 : do ii=1,3
624 [ + + ]: 2497836 : if (abs(xred_center(ii,iatom)-xred_ideal(ii,jatom)).le.Invar%tolmatch) then
625 : 246347 : ndir_match=ndir_match+1
626 : : else if ((abs(xred_center(ii,iatom)-xred_ideal(ii,jatom)-1.d0).le.Invar%tolmatch) &
627 [ + + + + ]: 1627030 : & .or.(abs(xred_center(ii,iatom)-xred_ideal(ii,jatom)+1.d0).le.Invar%tolmatch)) then
628 : 1114 : ndir_match=ndir_match+1
629 : 1114 : must_shift=.true.
630 : : endif
631 : : end do
632 [ + + ]: 624459 : if (ndir_match==3.and..not.must_shift) then
633 : 6882 : FromIdeal2Average(jatom)=iatom
634 : : ok=.true.
635 : : exit
636 [ - + ]: 617577 : else if (ndir_match==3.and.must_shift) then
637 [ # # ]: 0 : do ii=1,3
638 [ # # ]: 0 : if (abs(xred_center(ii,iatom)-xred_ideal(ii,jatom)-1.d0).le.Invar%tolmatch) then
639 : 0 : xred_center(ii,iatom)=xred_center(ii,iatom)-1d0
640 : 0 : xred_average(ii,iatom)=xred_average(ii,iatom)-1d0
641 [ # # ]: 0 : do istep=1,MD%my_nstep
642 : 0 : MD%xred(ii,iatom,istep)=MD%xred(ii,iatom,istep)-1d0
643 : : end do
644 : 0 : FromIdeal2Average(jatom)=iatom
645 [ # # ]: 0 : else if (abs(xred_center(ii,iatom)-xred_ideal(ii,jatom)+1.d0).le.Invar%tolmatch) then
646 : 0 : xred_center(ii,iatom)=xred_center(ii,iatom)+1d0
647 : 0 : xred_average(ii,iatom)=xred_average(ii,iatom)+1d0
648 [ # # ]: 0 : do istep=1,MD%my_nstep
649 : 0 : MD%xred(ii,iatom,istep)=MD%xred(ii,iatom,istep)+1d0
650 : : end do
651 : 0 : FromIdeal2Average(jatom)=iatom
652 : : end if
653 : : end do
654 : : ok=.true.
655 : : exit
656 : : end if
657 : : end do
658 [ - - ]: 44 : if (.not.ok) then
659 : 0 : write(Invar%stdlog,*) 'Problem to find the average position for iatom=',iatom
660 : 0 : write(Invar%stdlog,*) ' Reasons:'
661 : 0 : write(Invar%stdlog,*) ' 1/ One atom jump to another equilibrium position'
662 : 0 : write(Invar%stdlog,*) ' 2/ The system is no more solid'
663 : 0 : write(Invar%stdlog,*) ' 3/ Perhaps, you can adjust the tolerance (tolmatch)'
664 : 0 : write(Invar%stdlog,*) ' xred_center=',(xred_center(ii,iatom),ii=1,3)
665 [ # # ]: 0 : do eatom=1,MD%natom
666 : 0 : write(Invar%stdlog,'(a,1x,3(f10.6,1x))') 'I',xred_ideal (:,eatom)
667 : 0 : write(Invar%stdlog,'(a,1x,3(f10.6,1x))') 'C',xred_center(:,eatom)
668 : : end do
669 : 0 : ABI_ERROR('Problem to find the average position')
670 : : end if
671 : : end do
672 : :
673 : : ! WARNING: VERY IMPORTANT: The positions are displayed/sorted
674 : : ! (and used in the following) according to ideal positions xred_ideal.
675 : 44 : call tdep_write_xred_average(Invar,MPIdata,Lattice,xred_ideal,xred_center,FromIdeal2Average)
676 [ - + ]: 44 : ABI_FREE(xred_center)
677 : :
678 : : !====================================================================================
679 : : !====================== END OF REDUCED COORDINATES ==================================
680 : : !====================================================================================
681 : : ! a/ Get cartesian coordinates from reduced ones
682 : : ! b/ Compute ucart and fcart tabs
683 : : ! c/ The atoms are sorted according the IDEAL arrangement
684 : : ! The correspondance function is contained in: FromIdeal2Average
685 : : ! WARNING : Consequently the arrangement of the xcart* tabs is not modified.
686 : 44 : write(Invar%stdout,*)' Compute cartesian coordinates and forces...'
687 [ - + + - : 557736 : ABI_MALLOC(xcart ,(3,MD%natom,MD%my_nstep)); xcart(:,:,:)=0.d0
- + + - +
- + - - +
- + + + +
+ + + ]
688 [ - + + - : 27660 : ABI_MALLOC(xcart_ideal ,(3,MD%natom)) ; xcart_ideal(:,:)=0.d0
+ - + - -
+ - + + +
+ + ]
689 [ - + + - : 27616 : ABI_MALLOC(xcart_average,(3,MD%natom)) ; xcart_average(:,:)=0.d0
- + - + +
+ + + ]
690 [ - + + - : 557692 : ABI_MALLOC(ucart_tmp ,(3,MD%natom,MD%my_nstep)); ucart_tmp(:,:,:)=0.d0
- + + - -
+ - + + +
+ + + + ]
691 [ + + ]: 6926 : do iatom=1,MD%natom
692 : 6882 : call DGEMV('T',3,3,1.d0,Lattice%rprimd_md(:,:),3,xred_ideal (:,iatom),1,0.d0,xcart_ideal (:,iatom),1)
693 : 6882 : call DGEMV('T',3,3,1.d0,Lattice%rprimd_md(:,:),3,xred_average(:,iatom),1,0.d0,xcart_average(:,iatom),1)
694 [ + + ]: 21812 : do iatcell=1,MD%natom_unitcell
695 : 21768 : call DGEMV('T',3,3,1.d0,Lattice%rprimd_md(:,:),3,Rlatt_red(:,iatcell,iatom),1,0.d0,MD%Rlatt_cart(:,iatcell,iatom),1)
696 : : end do
697 : : end do
698 [ + + ]: 996 : do istep=1,MD%my_nstep
699 [ + + ]: 140148 : do iatom=1,MD%natom
700 : 139152 : jatom = FromIdeal2Average(iatom)
701 : : call DGEMV('T',3,3,1.d0,Lattice%rprimd_md(:,:),3,MD%xred(:,jatom,istep),&
702 : 139152 : & 1,0.d0,xcart(:,jatom,istep),1)
703 [ + + ]: 140104 : if (Invar%use_ideal_positions.eq.0) then
704 [ + + ]: 115648 : ucart_tmp(:,iatom,istep) = xcart(:,jatom,istep) - xcart_average(:,jatom)
705 : : else
706 [ + + ]: 440960 : ucart_tmp(:,iatom,istep) = xcart(:,jatom,istep) - xcart_ideal(:,iatom)
707 : : end if
708 : : end do
709 : : end do
710 [ - + ]: 44 : ABI_FREE(xred_average)
711 [ - + ]: 44 : ABI_FREE(xcart)
712 [ - + ]: 44 : ABI_FREE(xcart_ideal)
713 [ - + ]: 44 : ABI_FREE(xcart_average)
714 : :
715 : : ! Rearrangement of the fcart tabs in column --> MD%Forces
716 [ - + + - : 557736 : ABI_CALLOC(fcart_tmp,(3,MD%natom,MD%my_nstep))
- + + - +
- + - - +
- + + + +
+ + + ]
717 [ + + ]: 996 : do istep=1,MD%my_nstep
718 [ + + ]: 140148 : do iatom=1,MD%natom
719 [ + + ]: 557560 : fcart_tmp(:,iatom,istep) = MD%fcart(:,FromIdeal2Average(iatom),istep)
720 : : end do
721 : : end do
722 [ + + ]: 996 : do istep=1,MD%my_nstep
723 [ + + ]: 140148 : do jatom=1,MD%natom
724 [ + + ]: 557560 : do ii=1,3
725 : 417456 : jj = ii + 3*(jatom-1) + 3*MD%natom*(istep-1)
726 : 417456 : MD%Forces(jj) = fcart_tmp(ii,jatom,istep)
727 : 556608 : MD%ucart(ii,jatom,istep) = ucart_tmp(ii,jatom,istep)
728 : : enddo
729 : : enddo
730 : : enddo
731 [ - + ]: 44 : ABI_FREE(FromIdeal2Average)
732 [ - + ]: 44 : ABI_FREE(ucart_tmp)
733 [ - + ]: 44 : ABI_FREE(fcart_tmp)
734 : :
735 : : ! Define Rlatt_scaled, fulfilling the definition of mkphdos (ABINIT routine)
736 [ + + ]: 176 : do ii=1,3
737 [ + + ]: 572 : rprimd_md_tmp(ii,:) = Lattice%rprimd_md(ii,:) / Lattice%acell_unitcell(ii)
738 : : end do
739 [ + + ]: 6926 : do iatom=1,MD%natom
740 [ + + ]: 21812 : do iatcell=1,MD%natom_unitcell
741 : 21768 : call DGEMV('T',3,3,1.d0,rprimd_md_tmp,3,Rlatt_red(:,iatcell,iatom),1,0.d0,MD%Rlatt_scaled(:,iatcell,iatom),1)
742 : : end do
743 : : end do
744 : :
745 : : ! Find the symetry operation between 2 atoms
746 : 44 : call tdep_SearchS_1at(Invar,MPIdata,Sym,xred_ideal)
747 [ + + + + ]: 27572 : MD%xred_ideal(:,:)=xred_ideal(:,:)
748 [ - + ]: 44 : ABI_FREE(xred_ideal)
749 [ - + ]: 44 : ABI_FREE(Rlatt_red)
750 : :
751 : 44 : end subroutine tdep_MatchIdeal2Average
752 : :
753 : : !====================================================================================================
754 : :
755 : 44 : subroutine tdep_write_xred_average(Invar,MPIdata,Lattice,&
756 : 44 : xred_ideal,xred_center,&
757 : 44 : FromIdeal2Average)
758 : : type(atdep_dataset_type), intent(in) :: Invar
759 : : type(MPI_enreg_type),intent(in) :: MPIdata
760 : : type(Lattice_type),intent(in) :: Lattice
761 : : double precision,intent(in) :: xred_ideal(3,Invar%natom)
762 : : double precision,intent(in) :: xred_center(3,Invar%natom)
763 : : integer,intent(in),optional :: FromIdeal2Average(Invar%natom)
764 : :
765 : : integer :: unt
766 : : !integer :: natom,natom_unitcell
767 : : integer :: iatom,jatom,ii,jj
768 : : !logical :: with_xcart
769 : 44 : integer,allocatable :: ideal2average(:)
770 : : double precision :: rprimd(3,3)
771 : : double precision :: xred_C(3),xred_I(3),xcart_C(3),xcart_I(3)
772 : :
773 [ + - ]: 44 : if (MPIdata%iam_master) then
774 : :
775 [ + + + + ]: 572 : rprimd(:,:) = Lattice%rprimd_md(:,:)
776 : :
777 [ - + + - : 132 : ABI_MALLOC(ideal2average,(Invar%natom))
+ - + - -
+ - + ]
778 [ + + ]: 6926 : ideal2average(:)=0
779 [ + - ]: 44 : if (present(FromIdeal2Average)) then
780 [ + + ]: 6926 : ideal2average(:) = FromIdeal2Average(:)
781 : : else
782 [ # # ]: 0 : do iatom=1,Invar%natom
783 : 0 : ideal2average(iatom) = iatom
784 : : end do
785 : : end if
786 : :
787 : 44 : unt=31
788 : 44 : open(unit=unt,file=trim(Invar%output_prefix)//'_xred_average.xyz')
789 : 44 : write(unt,'(a,i4)') '# natom = ',Invar%natom
790 : 44 : write(unt,'(a,i4)') '# natom_unitcell = ',Invar%natom_unitcell
791 [ + - + + : 572 : write(unt,'(a,9(f4.1,1x))') '# multiplicity = ',((Lattice%multiplicity(ii,jj),jj=1,3),ii=1,3 )
+ - + + ]
792 : 44 : write(unt,'(a)') '#'
793 : :
794 : 44 : write(unt,'(a1,1x,a8,2x,a6,2x,2(a5,30x))') '#', 'position', 'iatom', 'xred ', 'xcart'
795 : 44 : write(unt,'(a)')''
796 : :
797 : 44 : xred_I = zero
798 : 44 : xred_C = zero
799 [ + + ]: 6926 : do iatom=1,Invar%natom
800 : 6882 : jatom = ideal2average(iatom)
801 [ + + ]: 27528 : xred_I = xred_ideal (:,iatom)
802 [ + + ]: 27528 : xred_C = xred_center(:,jatom)
803 : :
804 : 6882 : xcart_I(:)=zero
805 : 6882 : xcart_C(:)=zero
806 : 6882 : call DGEMV('T',3,3,1.d0,rprimd(:,:),3,xred_I,1,0.d0,xcart_I,1)
807 : 6882 : call DGEMV('T',3,3,1.d0,rprimd(:,:),3,xred_C,1,0.d0,xcart_C,1)
808 : :
809 : 6882 : write(unt,'(2x,a6,4x,i6,2x,3(f10.6,1x),2x,3(f10.6,1x))')'Ideal ',iatom,xred_I,xcart_I
810 : 6882 : write(unt,'(2x,a6,4x,i6,2x,3(f10.6,1x),2x,3(f10.6,1x))')'Center',jatom,xred_C,xcart_C
811 : 6926 : write(unt,'(a)')''
812 : :
813 : : end do
814 : :
815 : 44 : close(unt)
816 [ - + ]: 44 : ABI_FREE(ideal2average)
817 : : end if
818 : :
819 [ + - ]: 44 : end subroutine tdep_write_xred_average
820 : :
821 : : !====================================================================================================
822 : :
823 [ # # # # : 139152 : end module m_tdep_sampling
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # ]
824 : : !!***
|
