~roberto-robles/siesta/trunk-RRR

« back to all changes in this revision

Viewing changes to Src/m_energies.F90

  • Committer: Nick Papior
  • Date: 2016-08-18 05:56:45 UTC
  • Revision ID: nickpapior@gmail.com-20160818055645-judotnqdl75x5y9b
Run all tests, solved merge problems and updated convergence criteria

- Re-run all tests.

  This showed a couple of interesting things.

  1. There was a mistake in the compute_dm code after
  the PEXSI merge when spin-orbit coupling was introduced
  2. I have removed Hprev as it essentially is the same as
  Hold. Either way it shouldn't produce a huge difference
  in the tracking of the dEbs, etc. (after all they are
  not physically used other than for convergence criteria)
  3. The change to SCF.Mix Hamiltonian resulted in a huge
  number of changes in the output. This is because the first
  step prints out the energies at INIT. However, the Hamiltonian
  is different because it is initialized after the compute_dm step.

- Changed the logic in convergence criteria.
  Now the convergence criterias are additive and may be fully controlled.
  However, at least one convergence criteria must be used.

  Now the default convergence criteria is both the Hamiltonian and the density
  matrix.

  This is updated in the manual and the compatibility note.

- Initially I thought the above differences in the energies was
  due to inconsistencies after r538. Hence I have created some
  simple routines in the m_energies.f90 code which updates a selected
  few of the energies.
  I think this should be adopted in the future to ensure that all
  calls to update energies are consistent.
  This will make changes to energy calculations less error-prone.

- Implemented the spin-type in the following routines:
   compute_dm
   final_H_f_stress
   state_init

- Changed the m_compute_max_diff to an interface code with appropriate
  size calculations. There was no reason for explicitly using the
  sparse pattern.

- Fixed a bug in the molecularmechanics code (introduced by Nick r542)

- Added a cyclediffs.sh script which loops on OUT.diffs files and
  it lets one easily cycle the diffs, simply do:

  cd Tests
  make check
  ./cyclediffs.sh

  and answer all the questions. Basically it makes deletes OUT.diffs
  which you have agreed isn't really a change.

- The tests may now be runned via:

   make MPI="mpirun -np 4"
   
  which then uses the default SIESTA location.

  Currently the script checks whether mpirun/mpiexec
  is in SIESTA variable, and if so, does not use MPI variable.
  This makes it easier to decide on the number of cores without
  writing the full path.

Show diffs side-by-side

added added

removed removed

Lines of Context:
44
44
  real(dp):: Ebs        ! Band-structure energy, Tr(DM*H), calculated in compute_dm
45
45
  real(dp):: Eso        ! Spin-orbit energy
46
46
  real(dp):: Eldau      
47
 
  real(dp):: DEldau      
 
47
  real(dp):: DEldau
 
48
 
 
49
contains
 
50
 
 
51
  !> Initialize ALL energies to 0.
 
52
  subroutine init_Energies()
 
53
 
 
54
    DEharr = 0._dp
 
55
    DEna = 0._dp
 
56
    DUext = 0._dp
 
57
    DUscf = 0._dp
 
58
    Dxc = 0._dp
 
59
    Ecorrec = 0._dp
 
60
    ef = 0._dp
 
61
    Eharrs = 0._dp
 
62
    Eions = 0._dp
 
63
    Ekin = 0._dp
 
64
    Ekinion = 0._dp
 
65
    Emad = 0._dp
 
66
    Ena = 0._dp
 
67
    Enaatm = 0._dp
 
68
    Enascf = 0._dp
 
69
    Enl = 0._dp
 
70
    Emeta = 0._dp
 
71
    Entropy = 0._dp
 
72
    Etot = 0._dp
 
73
    Exc = 0._dp
 
74
    E0 = 0._dp
 
75
    Emm = 0._dp
 
76
    FreeE = 0._dp
 
77
    FreeEharris = 0._dp
 
78
    Uatm = 0._dp
 
79
    Uscf = 0._dp
 
80
    Ebs = 0._dp
 
81
    Eso = 0._dp
 
82
    Eldau = 0._dp      
 
83
    DEldau = 0._dp
 
84
 
 
85
  end subroutine init_Energies
 
86
 
 
87
  !> To ease the computation of specific deferred
 
88
  !> quantites we allow the computation of these quantites
 
89
  !> In ONE place
 
90
 
 
91
  subroutine update_DEna()
 
92
 
 
93
    DEna = Enascf - Enaatm
 
94
    
 
95
  end subroutine update_DEna
 
96
 
 
97
  subroutine update_E0()
 
98
 
 
99
    E0 = Ena + Ekin + Enl + Eso - Eions
 
100
    
 
101
  end subroutine update_E0
 
102
  
 
103
  subroutine update_Etot()
 
104
    
 
105
    ! DUext (external electric field) -- should it be in or out?
 
106
    Etot = Ena + Ekin + Enl + Eso - Eions + &
 
107
         DEna + DUscf + DUext + Exc + &
 
108
         Ecorrec + Emad + Emm + Emeta + Eldau
 
109
    
 
110
  end subroutine update_Etot
 
111
 
 
112
  !> @param kBT the temperature in energy
 
113
  subroutine update_FreeE( kBT )
 
114
    real(dp), intent(in) :: kBT
 
115
 
 
116
    FreeE = Etot - kBT * Entropy
 
117
 
 
118
  end subroutine update_FreeE
 
119
 
 
120
  !> @param kBT the temperature in energy
 
121
  subroutine update_FreeEHarris( kBT )
 
122
    real(dp), intent(in) :: kBT
 
123
 
 
124
    FreeEHarris = Eharrs - kBT * Entropy
 
125
 
 
126
  end subroutine update_FreeEHarris
48
127
 
49
128
end module m_energies
50
129