TABLE OF CONTENTS


ABINIT/m_spin_hist [ Modules ]

[ Top ] [ Modules ]

NAME

 m_spin_hist

FUNCTION

 This module contains definition the type spin_hist_t
 and its related routines

 Datatypes:

 * spin_hist_t: history record of spin orientations and amplitudes

 Subroutines:

 * spin_hist_t_init
 * spin_hist_t_free
 * spin_hist_t
 * spin_hist_t_get_S
 * spin_hist_t_findIndex
 * spin_hist_t_set_vars
 * spin_hist_t_set_params

COPYRIGHT

 Copyright (C) 2001-2017 ABINIT group (hexu)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

SOURCE

35 ! TODO hexu:
36 ! sync ihist_latt when with lattice dynamics
37 ! add average , variance, etc (should they be here?)
38 ! structural information and some parameters are no longer 
39 ! used here. They should be removed form this file.
40 
41 #if defined HAVE_CONFIG_H
42 #include "config.h"
43 #endif
44 
45 #include "abi_common.h"
46 module m_spin_hist
47   use defs_basis
48   use m_abicore
49   use m_errors
50   use m_xmpi
51   implicit none
52 
53   private

m_spin_hist/spin_hist_t [ Types ]

[ Top ] [ m_spin_hist ] [ Types ]

NAME

 spin_hist_t

FUNCTION

 This type has several vectors, and index scalars to store
 a proper history of previous evaluations of forces and
 stresses,velocities,positions and energies

 It contains:
 * mxhist                  : Maximum size of history
 * ihist                   : index of history
 natom : number of atoms
 nmatom: number of magnetic atoms
 * acell(3)         : Acell (acell , rprimd, xred: only initial value kept if there is!!  no lattice dynamics. Other wise for each step, the corresponding lattice step number is kept)
 * rprimd(3,3)      : Rprimd
 * xred(3,natom)    : Xred
 * index_spin     : the index of atom in spin model, -1 if it is not in the spin model 
 * heff(3,nmatom,mxhist)   : effective magnetic field (cartesian)
 * snorm(nmatom, mxhist) : magnetitude of spin.
 * S(3,nmatom,mxhist)   : spin orientation of atoms (cartesian)
 * dSdt(3, nmatom, mxhist) : dS/dt (cartesian)
 * etot(mxhist)            : Electronic total Energy
 * entropy(mxhist)         : Entropy
 * itime(mxhist)           : index of spin dynamics step.
 * time(mxhist)            : Time (or iteration number for GO)

 * has_latt (whether lattice dynamics is also present)
 * ihist_latt(mxhist): the corresponding lattice step. 0 if none.

SOURCE

 90   type, public :: spin_hist_t
 91      ! scalars
 92      ! Index of the last element on all records
 93      integer :: ihist = 0
 94      integer :: ihist_prev = -1
 95      ! Maximun size of the historical records
 96      integer :: mxhist = 0
 97 
 98      integer :: nmatom
 99      ! whether lattice dynamics is also present
100      integer, allocatable :: ihist_latt(:)
101      logical :: has_latt
102 
103      ! arrays
104      !  placeholders for structure-related parameters. They are not used currently. 
105      integer :: natom
106      real(dp) :: acell(3)
107      real(dp) :: rprimd(3,3)
108      real(dp), allocatable :: xred(:, :)
109      integer :: ntypat
110      integer, allocatable :: typat(:)
111      real(dp), allocatable :: znucl(:)
112      integer, allocatable :: spin_index(:)
113 
114      ! spin
115      !heff(3, nmatom, mxhist)
116      real(dp), allocatable :: heff(:, :, :)
117      !snorm(nmatom, mxhist)
118      real(dp), allocatable :: snorm(:, :)
119 
120      !S(3, nmatom, mxhist)
121      real(dp), allocatable :: S(:, :, :)
122      !dSdt(3, nmatom, mxhist)
123      ! TODO hexu: is it useful?
124      real(dp), allocatable :: dSdt(:, :, :)
125 
126      ! label of spins, sublattice
127      integer(dp), allocatable :: label(:)
128      ! sc_label(3, nmatom) label of cell. (R)
129      integer(dp), allocatable :: sc_label(:,:)
130      ! etot(mxhist)
131      real(dp), allocatable :: etot(:)
132      real(dp), allocatable :: entropy(:)
133      real(dp), allocatable :: time(:)
134      integer, allocatable :: itime(:)
135      ! id of netcdf spin hist file
136      integer :: ncid
137      ! spin_nctime: interval of step for writing to netcdf hist file.
138      integer :: spin_nctime
139      real(dp) :: spin_temperature
140   end type spin_hist_t

m_spin_hist/spin_hist_t_free [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 spin_hist_t_free

FUNCTION

 free memory for spin_hist_t

INPUTS

OUTPUT

 hist <type(spin_hist_t)()> = spin hist type

PARENTS

      m_spin_hist

CHILDREN

SOURCE

345   subroutine spin_hist_t_free(hist)
346 
347 
348 !This section has been created automatically by the script Abilint (TD).
349 !Do not modify the following lines by hand.
350 #undef ABI_FUNC
351 #define ABI_FUNC 'spin_hist_t_free'
352 !End of the abilint section
353 
354     class(spin_hist_t) , intent(inout) :: hist
355 
356     if (allocated(hist%xred)) then
357        ABI_DEALLOCATE(hist%xred)
358     end if
359     if (allocated(hist%typat)) then
360        ABI_DEALLOCATE(hist%typat)
361     end if
362     if (allocated(hist%znucl)) then
363        ABI_DEALLOCATE(hist%znucl)
364     end if
365     if (allocated(hist%label)) then
366        ABI_DEALLOCATE(hist%label)
367     end if
368     if (allocated(hist%spin_index)) then
369        ABI_DEALLOCATE(hist%spin_index)
370     end if
371     if (allocated(hist%heff)) then
372        ABI_DEALLOCATE(hist%heff)
373     end if
374     if (allocated(hist%snorm)) then
375        ABI_DEALLOCATE(hist%snorm)
376     end if
377     if (allocated(hist%S)) then
378        ABI_DEALLOCATE(hist%S)
379     end if
380     if (allocated(hist%dSdt)) then
381        ABI_DEALLOCATE(hist%dSdt)
382     end if
383     if (allocated(hist%etot)) then
384        ABI_DEALLOCATE(hist%etot)
385     end if
386     if (allocated(hist%entropy)) then
387        ABI_DEALLOCATE(hist%entropy)
388     end if
389     if (allocated(hist%time)) then
390        ABI_DEALLOCATE(hist%time)
391     end if
392     if (allocated(hist%itime)) then
393        ABI_DEALLOCATE(hist%itime)
394     end if
395     if (allocated(hist%ihist_latt)) then
396        ABI_DEALLOCATE(hist%ihist_latt)
397     end if
398   end subroutine spin_hist_t_free

m_spin_hist/spin_hist_t_get_S [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 spin_hist_t_get_S

FUNCTION

 get the S for step. step=0 is current. step=-1 is last...

INPUTS

 hist <type(spin_hist_t)()> = spin hist type
 step = index of step. current step is 0. last step is -1. 

OUTPUT

 S(3, nmatom)=spin orientations at step

PARENTS

      m_spin_hist

CHILDREN

SOURCE

423   function spin_hist_t_get_S(hist, step) result(S)
424 
425 
426 !This section has been created automatically by the script Abilint (TD).
427 !Do not modify the following lines by hand.
428 #undef ABI_FUNC
429 #define ABI_FUNC 'spin_hist_t_get_S'
430 !End of the abilint section
431 
432     class(spin_hist_t), intent(inout) :: hist
433     integer, intent(in), optional:: step
434     real(dp) :: S(3, hist%nmatom)
435     integer :: i, j
436     if (.not. present(step)) then
437        j=0
438     else
439        j=step
440     end if
441     i=spin_hist_t_findIndex(hist,step=j)
442     S(:,:)=hist%S(:,:,i)
443   end function spin_hist_t_get_S

m_spin_hist/spin_hist_t_inc [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 spin_hist_t_get_inc

FUNCTION

 time counter increase

INPUTS

OUTPUT

   hist <type(spin_hist_t)()> = spin hist type

PARENTS

      m_spin_hist

CHILDREN

SOURCE

465   subroutine spin_hist_t_inc(hist)
466 
467 
468 !This section has been created automatically by the script Abilint (TD).
469 !Do not modify the following lines by hand.
470 #undef ABI_FUNC
471 #define ABI_FUNC 'spin_hist_t_inc'
472 !End of the abilint section
473 
474     class(spin_hist_t), intent(inout) :: hist
475     if(hist%ihist_prev ==0 ) then
476         hist%itime(hist%ihist)=1
477     else
478         hist%itime(hist%ihist)=hist%itime(hist%ihist_prev)+1
479     endif
480     hist%ihist_prev=hist%ihist
481     hist%ihist=spin_hist_t_findIndex(hist, 1)
482   end subroutine spin_hist_t_inc

m_spin_hist/spin_hist_t_init [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 spin_hist_t_init

FUNCTION

 initialize spin hist

INPUTS

 nmatom = number of magnetic atoms 
 mxhist = maximum number of hist steps
 has_latt = whether spin dynamics in with lattice dynamics

OUTPUT

 hist <type(spin_hist_t)()> = spin hist type

PARENTS

      m_spin_hist

CHILDREN

SOURCE

180   subroutine spin_hist_t_init(hist, nmatom, mxhist, has_latt)
181 
182 
183 !This section has been created automatically by the script Abilint (TD).
184 !Do not modify the following lines by hand.
185 #undef ABI_FUNC
186 #define ABI_FUNC 'spin_hist_t_init'
187 !End of the abilint section
188 
189     implicit none
190     class(spin_hist_t), intent(inout) :: hist
191     integer, intent(in) :: nmatom, mxhist
192     logical, intent(in) :: has_latt
193     hist%nmatom=nmatom
194     hist%ntypat=0
195     hist%ihist=1
196     hist%ihist_prev=0
197     hist%mxhist=mxhist
198     hist%natom=0
199     hist%has_latt=has_latt
200 
201     !print *, "initialize HIST spin"
202     ABI_ALLOCATE(hist%heff, (3, nmatom, mxhist))
203     ABI_ALLOCATE(hist%snorm, (nmatom, mxhist))
204     ABI_ALLOCATE(hist%S, (3, nmatom, mxhist))
205     ABI_ALLOCATE(hist%dSdt, (3, nmatom, mxhist))
206     ABI_ALLOCATE(hist%label, (nmatom))
207 
208     ABI_ALLOCATE(hist%etot, (mxhist))
209     ABI_ALLOCATE(hist%entropy, (mxhist))
210     ABI_ALLOCATE(hist%time, (mxhist))
211     ABI_ALLOCATE(hist%itime, (mxhist))
212 
213     ABI_ALLOCATE(hist%ihist_latt, (mxhist))
214 
215     hist%etot(1) =zero
216     hist%entropy(1) =zero
217     hist%time(1) =zero
218 
219     !hist%acell(:)=zero
220     !hist%rprimd(:, :)=zero
221     !hist%xred(:,:) =zero
222     hist%heff(:,:,1)=zero
223     hist%S(:,:,1)=zero
224     hist%dSdt(:,:,1)=zero
225     hist%snorm(:,1)=zero
226     !print *, "Initialization spin hist finished"
227   end subroutine spin_hist_t_init

m_spin_hist/spin_hist_t_set_atomic_structure [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 spin_hist_t_set_atomic_structure

FUNCTION

 set atomic structure 

INPUTS

 acell(3) = acell
 rprimd(3, 3) = 
 xred(3, natom) = positions in reduced coordinates
 spin_index(3, natom) = index of atom in spin hamiltonian
 ntypat = number of types of atoms
 typat(ntypat)=types of atoms
 znucl=z of atoms

OUTPUT

 hist <type(spin_hist_t)()> = spin hist type

PARENTS

      m_spin_hist

CHILDREN

SOURCE

256   subroutine spin_hist_t_set_atomic_structure(hist, acell, rprimd, xred, spin_index, ntypat,  typat, znucl)
257 
258 
259 !This section has been created automatically by the script Abilint (TD).
260 !Do not modify the following lines by hand.
261 #undef ABI_FUNC
262 #define ABI_FUNC 'spin_hist_t_set_atomic_structure'
263 !End of the abilint section
264 
265     class(spin_hist_t), intent(inout) :: hist
266     real(dp), intent(in) :: acell(3), rprimd(3,3), xred(:,:), znucl(:)
267     integer, intent(in):: spin_index(:), ntypat, typat(:)
268     integer :: natom
269     natom=size(typat)
270     ABI_ALLOCATE(hist%xred, (3, natom))
271     ABI_ALLOCATE(hist%spin_index, (natom))
272     ABI_ALLOCATE(hist%typat,(ntypat))
273     ABI_ALLOCATE(hist%znucl, (ntypat))
274 
275     hist%acell(:)=acell(:)
276     hist%rprimd(:,:)=rprimd(:,:)
277     hist%xred(:,:)=xred(:,:)
278     hist%spin_index(:)=spin_index(:)
279     hist%ntypat=ntypat
280     hist%typat(:)=typat(:)
281     hist%znucl(:)=znucl(:)
282   end subroutine spin_hist_t_set_atomic_structure

m_spin_hist/spin_hist_t_set_params [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 spin_hist_t_set_params

FUNCTION

 set parameters for spin_hist_t

INPUTS

 spin_nctime=number of step between two write to netcdf hist file
 spin_temperate= temperature of spin

OUTPUT

 hist <type(spin_hist_t)()> = spin hist type

PARENTS

      m_spin_hist

CHILDREN

SOURCE

309   subroutine spin_hist_t_set_params(hist, spin_nctime, spin_temperature)
310 
311 
312 !This section has been created automatically by the script Abilint (TD).
313 !Do not modify the following lines by hand.
314 #undef ABI_FUNC
315 #define ABI_FUNC 'spin_hist_t_set_params'
316 !End of the abilint section
317 
318     class(spin_hist_t), intent(inout) :: hist
319     integer, intent(in) :: spin_nctime
320     real(dp), intent(in) :: spin_temperature
321     hist%spin_nctime= spin_nctime
322     hist%spin_temperature=spin_temperature
323   end subroutine spin_hist_t_set_params