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
 The observables are also calculated. 

 Datatypes:

 * spin_hist_t: history record of spin orientations and amplitudes

 Subroutines:

 * init
 * free
 * spin_hist_t
 * get_S
 * findIndex
 * set_vars
 * set_params

COPYRIGHT

 Copyright (C) 2001-2022 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

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

m_spin_hist/finalize [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 finalize

FUNCTION

 free memory for spin_hist_t

INPUTS

OUTPUT

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

SOURCE

336   subroutine finalize(self)
337 
338     class(spin_hist_t) , intent(inout) :: self 
339 
340     if (allocated(self%xred)) then
341        ABI_FREE(self%xred)
342     end if
343     if (allocated(self%typat)) then
344        ABI_FREE(self%typat)
345     end if
346     if (allocated(self%znucl)) then
347        ABI_FREE(self%znucl)
348     end if
349     if (allocated(self%spin_index)) then
350        ABI_FREE(self%spin_index)
351     end if
352     if (allocated(self%heff)) then
353        ABI_FREE(self%heff)
354     end if
355     if (allocated(self%snorm)) then
356        ABI_FREE(self%snorm)
357     end if
358     if (allocated(self%S)) then
359        ABI_FREE(self%S)
360     end if
361     if (allocated(self%dSdt)) then
362        ABI_FREE(self%dSdt)
363     end if
364     if (allocated(self%etot)) then
365        ABI_FREE(self%etot)
366     end if
367     if (allocated(self%entropy)) then
368        ABI_FREE(self%entropy)
369     end if
370     if (allocated(self%time)) then
371        ABI_FREE(self%time)
372     end if
373     if (allocated(self%itime)) then
374        ABI_FREE(self%itime)
375     end if
376     if (allocated(self%ihist_latt)) then
377        ABI_FREE(self%ihist_latt)
378     end if
379 
380   end subroutine finalize

m_spin_hist/get_S [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 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, nspin)=spin orientations at step

SOURCE

400   function get_S(self, step) result(S)
401     class(spin_hist_t), intent(inout) :: self
402     integer, intent(in), optional:: step
403     real(dp) :: S(3, self%nspin)
404     integer :: i, j
405     if (.not. present(step)) then
406        j=0
407     else
408        j=step
409     end if
410     i=self%findIndex(step=j)
411     S(:,:)=self%S(:,:,i)
412   end function get_S

m_spin_hist/inc1 [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 inc1

FUNCTION

 time counter increase

INPUTS

OUTPUT

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

SOURCE

429   subroutine inc1(self)
430 
431     class(spin_hist_t), intent(inout) :: self
432     if(self%ihist_prev ==0 ) then
433         self%itime(self%ihist)=1
434     else
435         self%itime(self%ihist)=self%itime(self%ihist_prev)+1
436     endif
437     self%ihist_prev=self%ihist
438     self%ihist=self%findIndex(1)
439   end subroutine inc1

m_spin_hist/initialize [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 initialize

FUNCTION

 initialize spin hist

INPUTS

 nspin = 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

SOURCE

182   subroutine initialize(self, nspin, mxhist, has_latt)
183 
184     implicit none
185     class(spin_hist_t), intent(inout) :: self 
186     integer, intent(in) :: nspin, mxhist
187     logical, intent(in) :: has_latt
188     !integer, optional,  intent(in) :: calc_traj_obs, calc_thermo_obs, calc_correlation_obs
189 
190     self%nspin=nspin
191     self%ntypat=0
192     self%ihist=1
193     self%ihist_prev=0
194     self%mxhist=mxhist
195     self%natoms=0
196     self%has_latt=has_latt
197 
198     ABI_MALLOC(self%heff, (3, nspin, mxhist))
199     ABI_MALLOC(self%snorm, (nspin, mxhist))
200     ABI_MALLOC(self%S, (3, nspin, mxhist))
201     ABI_MALLOC(self%dSdt, (3, nspin, mxhist))
202 
203     ABI_MALLOC(self%etot, (mxhist))
204     ABI_MALLOC(self%entropy, (mxhist))
205     ABI_MALLOC(self%time, (mxhist))
206     ABI_MALLOC(self%itime, (mxhist))
207 
208     ABI_MALLOC(self%ihist_latt, (mxhist))
209 
210     ! TODO: add observable allocation here.
211 
212     self%etot(1) =zero
213     self%entropy(1) =zero
214     self%time(1) =zero
215 
216     !self%acell(:)=zero
217     !self%rprimd(:, :)=zero
218     !self%xred(:,:) =zero
219     self%heff(:,:,:)=zero
220     self%S(:,:,:)=zero
221     self%dSdt(:,:,:)=zero
222     self%snorm(:,:)=zero
223   end subroutine initialize

m_spin_hist/set_atomic_structure [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 set_atomic_structure

FUNCTION

 set atomic structure 

INPUTS

 acell(3) = acell
 rprimd(3, 3) = 
 xred(3, natoms) = positions in reduced coordinates
 spin_index(3, natoms) = 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

SOURCE

272   subroutine set_atomic_structure(self, acell, rprimd, xred, spin_index, ntypat,  typat, znucl)
273 
274     class(spin_hist_t), intent(inout) :: self
275     real(dp), intent(in) :: acell(3), rprimd(3,3), xred(:,:), znucl(:)
276     integer, intent(in):: spin_index(:), ntypat, typat(:)
277     integer :: natoms
278     natoms=size(typat)
279     ABI_MALLOC(self%xred, (3, natoms))
280     ABI_MALLOC(self%spin_index, (natoms))
281     ABI_MALLOC(self%typat,(ntypat))
282     ABI_MALLOC(self%znucl, (ntypat))
283 
284     self%acell(:)=acell(:)
285     self%rprimd(:,:)=rprimd(:,:)
286     self%xred(:,:)=xred(:,:)
287     self%spin_index(:)=spin_index(:)
288     self%ntypat=ntypat
289     self%typat(:)=typat(:)
290     self%znucl(:)=znucl(:)
291   end subroutine set_atomic_structure

m_spin_hist/set_params [ Functions ]

[ Top ] [ m_spin_hist ] [ Functions ]

NAME

 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

SOURCE

312   subroutine set_params(self, spin_nctime, spin_temperature)
313 
314     class(spin_hist_t), intent(inout) :: self
315     integer, intent(in) :: spin_nctime
316     real(dp), intent(in) :: spin_temperature
317     self%spin_nctime= spin_nctime
318     self%spin_temperature=spin_temperature
319   end subroutine set_params

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
 natoms : number of atoms
 nspin: 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,natoms)    : Xred
 * index_spin     : the index of atom in spin model, -1 if it is not in the spin model 
 * heff(3,nspin,mxhist)   : effective magnetic field (cartesian)
 * snorm(nspin, mxhist) : magnetitude of spin.
 * S(3,nspin,mxhist)   : spin orientation of atoms (cartesian)
 * dSdt(3, nspin, 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

 89   type, public :: spin_hist_t
 90      ! scalars
 91      ! Index of the last element on all records
 92      integer :: ihist = 0
 93      integer :: ihist_prev = -1
 94      ! Maximun size of the historical records
 95      integer :: mxhist = 0
 96 
 97      integer :: nspin, nspin_prim
 98      ! whether lattice dynamics is also present
 99      integer, allocatable :: ihist_latt(:)
100      logical :: has_latt
101 
102      ! arrays
103      !  placeholders for structure-related parameters. They are not used currently. 
104      integer :: natoms
105      real(dp) :: acell(3)
106      real(dp) :: rprimd(3,3)
107      real(dp), allocatable :: xred(:, :)
108      integer :: ntypat
109      integer, allocatable :: typat(:)
110      real(dp), allocatable :: znucl(:)
111      integer, allocatable :: spin_index(:)
112 
113      ! spin
114      !heff(3, nspin, mxhist)
115      real(dp), allocatable :: heff(:, :, :)
116      !snorm(nspin, mxhist)
117      real(dp), allocatable :: snorm(:, :)
118 
119      !S(3, nspin, mxhist)
120      real(dp), allocatable :: S(:, :, :)
121      !dSdt(3, nspin, mxhist)
122      ! TODO hexu: is it useful?
123      real(dp), allocatable :: dSdt(:, :, :)
124 
125      ! etot(mxhist)
126      real(dp), allocatable :: etot(:)
127      real(dp), allocatable :: entropy(:)
128      real(dp), allocatable :: time(:)
129      integer, allocatable :: itime(:)
130 
131      ! spin_nctime: interval of step for writing to netcdf hist file.
132      integer :: spin_nctime
133      real(dp) :: spin_temperature
134 
135      ! observables
136      integer:: calc_thermo_obs, calc_traj_obs, calc_correlation_obs
137 
138      real(dp), allocatable :: ms_sub(:,:)   ! staggered M. 
139      real(dp), allocatable :: Cv(:) ! specfic heat
140      real(dp), allocatable :: binderU4_sub(:,:), binderU4(:)
141      real(dp), allocatable :: chi_sub(:, :), chi(:) ! magnetic susceptibility
142      real(dp), allocatable :: rcorr(:,:)
143      real(dp), allocatable :: sp_corr_func(:,:,:)
144    contains
145      procedure :: initialize
146      procedure :: finalize
147      procedure :: reset
148      procedure :: set_vars
149      procedure :: get_S => get_S
150      procedure :: findIndex => findIndex
151      procedure :: set_params => set_params
152      procedure :: inc1
153   end type spin_hist_t