~albertog/siesta/efield-2.5-jjunquera

« back to all changes in this revision

Viewing changes to Src/cell_broyden_optim.F

  • Committer: Alberto Garcia
  • Date: 2007-11-22 21:34:42 UTC
  • mfrom: (unknown (missing))
  • Revision ID: Arch-1:siesta@uam.es--2006%siesta-devel--reference--2.5--patch-2
Optimization enhancements (fire algorithm, step control).

* Step-size control has been add to the Broyden optimizers (Standard, ZM,
  and Cell-only). There is no theoretical justification for this procedure,
  so the Hessian is not likely to improve very much... it is probably best
  to combine a size-limited CG series with another small-Jinv Broyden series.

* New FORCE_STRESS for direct output of the energy, forces, and stress
  (in internal Siesta units).

* Re-implemented the FIRE optimization algorithm, both for geometry
  and (experimentally) for SCF.

  The advantages for the SCF mixing are:

  * Trivial algorithm (also in parallel)
  * Initial dt can be just set to alpha. 
  * Adaptive changes to alpha.
  * There should be no need for kicks.

  (Set DM.FIRE.Mixing to T to enable)
  Could also set DM.Fire.Nmin  to 2, 3, 4  instead of the default 5
  (Actually, this number could be akin to the DM.Number{Pulay,Broyden}...)

* Added .make files for opterons and matgas cluster at ICMAB in parallel.

* Compilation fixes(uncovered by a new version of g95)

  * External statements referring to obsolete names
  * Continuation lines.
  * Statement order.

(replayed patch-2 of outlier branch ref-2.4)


Show diffs side-by-side

added added

removed removed

Lines of Context:
34
34
 
35
35
      CONTAINS
36
36
      subroutine cell_broyden_optimizer( na, xa, cell, stress,
37
 
     $           tp, strtol, varcel, relaxd)
 
37
     $                tp, strtol, varcel, relaxd)
38
38
c ***************************************************************************
39
39
c Broyden geometry optimization (Cell Only)
40
40
c
79
79
 
80
80
c Internal variables and arrays
81
81
 
82
 
      real(dp)       :: volume, new_volume, trace
 
82
      real(dp)            :: new_volume, trace, volume
83
83
 
84
84
      logical           found
85
85
      integer           iu, ia, i, j, n, indi
87
87
      real(dp) ::  celli(3,3), sxx, syy, szz, sxy, sxz, syz
88
88
      real(dp) ::  stress_dif(3,3)
89
89
 
90
 
      real(dp), dimension(:), allocatable       :: gxa, gfa
 
90
      real(dp), dimension(:), allocatable       :: gxa, gfa, max_step
91
91
      real(dp), dimension(:), allocatable       :: rold, rnew, rdiff
92
92
 
93
93
 
107
107
      logical, save           :: initialization_done = .false.
108
108
      integer, save  :: numel
109
109
 
110
 
      real(dp), save :: jinv0
 
110
      real(dp), save :: jinv0, max_step_strain
111
111
      integer, save  :: maxit
112
112
      logical, save  :: cycle_on_maxit, variable_weight
113
113
      logical, save  :: broyden_debug
119
119
c ---------------------------------------------------------------------------
120
120
 
121
121
      volume = volcel(cell)
 
122
 
122
123
      if (.not. initialization_done) then
123
124
 
124
125
          call fdf_global_get(maxit,"MD.Broyden.History.Steps",5)
130
131
     $         "MD.Broyden.Debug",.false.)
131
132
          call fdf_global_get(jinv0,
132
133
     $         "MD.Broyden.Initial.Inverse.Jacobian",1.0_dp)
 
134
          call fdf_global_get(Max_step_strain,
 
135
     $         "MD.MaxCGDispl",0.02_dp,"Bohr")
133
136
 
134
137
         if (ionode) then
135
138
          write(6,*) 
237
240
         allocate(gfa(numel), stat=mem_stat)
238
241
         call memory('A','D',numel,'Broyden_optim',
239
242
     $        stat=mem_stat,id="gfa")
 
243
 
240
244
         allocate(gxa(numel), stat=mem_stat)
241
245
         call memory('A','D',numel,'Broyden_optim',
242
246
     $        stat=mem_stat,id="gxa")
 
247
 
 
248
         allocate(max_step(numel), stat=mem_stat)
 
249
         call memory('A','D',numel,'Broyden_optim',
 
250
     $        stat=mem_stat,id="max_step")
 
251
 
243
252
         allocate (rnew(numel), stat=mem_stat)
244
253
         call memory('A','D',numel,'Broyden_optim',
245
254
     $        stat=mem_stat,id="rnew")
278
287
              indi = indi + 1
279
288
              gfa(indi) = -stress_dif(i,j)*volume/precon
280
289
              gxa(indi) = strain(i,j) * precon
 
290
              max_step(indi) = Max_step_strain
281
291
           enddo
282
292
        enddo
283
293
 
296
306
 
297
307
        if (.not. broyden_is_setup(br)) then
298
308
           call broyden_reset(br,numel,maxit,cycle_on_maxit,
299
 
     $                       jinv0,0.01_dp)
 
309
     $                       jinv0,0.01_dp,max_step)
300
310
        endif
301
311
 
302
312
        weight = 1.0_dp
339
349
        deallocate (rnew, stat=mem_stat)
340
350
        call memory('D','D',numel,'Broyden_optim',
341
351
     $       stat=mem_stat,id="rnew")
 
352
 
342
353
        deallocate (gxa, stat=mem_stat)
343
354
        call memory('D','D',numel,'Broyden_optim',
344
355
     $       stat=mem_stat,id="gxa")
 
356
 
 
357
        deallocate (max_step, stat=mem_stat)
 
358
        call memory('D','D',numel,'Broyden_optim',
 
359
     $       stat=mem_stat,id="max_step")
 
360
 
345
361
        deallocate (gfa, stat=mem_stat)
346
362
        call memory('D','D',numel,'Broyden_optim',
347
363
     $       stat=mem_stat,id="gfa")