← Back to branch summary

~siesta-sip/siesta/sip-solver

« back to all changes in this revision

Viewing changes to Src/write_subs.F

  • Committer: Alberto Garcia
  • Date: 2015-06-22 09:32:02 UTC
  • Revision ID: albertog@icmab.es-20150622093202-enhjztqfku2e976q
Revert previous change affecting FA and FC files. FAC/FCC files.

* At revno 464, a change was introduced to output the constrained
  forces to the FA and FC files.  This had unintended consequences for
  some analysis tools.

  The change has been reverted, and the FC/FA files are now created as
  before. In addition, if constraints are used, new FCC/FAC files
  will be created with the constrained forces.

  The vibra utility will default to use the FCC files if the
  GeometryConstraints block exists, and in any case the
  force-constants file to use can be specified with the Vibra.FC fdf
  option.

  To support these changes, the write_forces subroutine has now 
  the geometry step as an additional argument.

  The writing of forces has been moved to its appropriate place in the
  write_subs modules.

* Some more cosmetic changes in the Vibra package.

modified:
  Docs/siesta.tex
  Src/Makefile
  Src/iofa.f
  Src/ofc.f
  Src/siesta_analysis.F
  Src/siesta_forces.F
  Src/state_analysis.F
  Src/write_subs.F
  Util/Vibra/Docs/CHANGES
  Util/Vibra/Docs/vibra.tex
  Util/Vibra/Examples/README
  Util/Vibra/Examples/si54.bands
  Util/Vibra/Examples/si54.fdf
  Util/Vibra/Src/Makefile
  Util/Vibra/Src/fcbuild.f
  Util/Vibra/Src/recoor.f
  Util/Vibra/Src/vibra.f

Show diffs side-by-side

added added

removed removed

Lines of Context:
Src/write_subs.F
65
65
      use units
66
66
      use m_energies, only: FreeE
67
67
      use m_stress, only: stress, kin_stress, mstress, cstress, tstress
68
 
      use m_forces,   only: fa, cfa   ! This does not really belong here  
69
68
 
70
69
      implicit none
71
70
 
187
186
 
188
187
          endif                  !varcel
189
188
 
190
 
        ! Write Force Constant matrix if FC calculation ...
191
 
        case(6)
192
 
          call ofc(cfa,dx,na_u)
193
 
        case(7)
194
 
!          call phonon_write_forces(fa,na_u,ucell,istep)
195
 
           if (IOnode) write(*,*) 'phonon suport deactivated'
196
189
        end select !idyn
197
190
 
198
191
      else !final
579
572
     &          siesta_write_energies, siesta_write_positions
580
573
      CONTAINS
581
574
      
582
 
      subroutine siesta_write_forces()
 
575
      subroutine siesta_write_forces(istep)
 
576
      use parallel, only: IOnode
583
577
      USE siesta_options
584
578
      use siesta_geom
585
579
      use siesta_cml
587
581
      use m_forces
588
582
      use alloc,      only : re_alloc, de_alloc
589
583
      implicit none
 
584
      integer, intent(in) :: istep
590
585
      integer :: ia, ix
591
586
      real(dp) :: cfmax, fmax, fres
592
587
      real(dp) :: ftot(3), cftot(3)
626
621
        if (writef) then
627
622
          write(6,'(i6,3f12.6)')(ia,(fa(ix,ia)*Ang/eV,ix=1,3),ia=1,na_u)
628
623
        endif
 
624
 
629
625
        ! Always write forces to .FA file
630
 
        call iofa( na_u, cfa )
 
626
        call iofa( na_u, fa , .false.)
 
627
        if ( any(fa/=cfa) ) then
 
628
           ! if any one constraint is not the same
 
629
           ! as the real forces, the FAC file 
 
630
           ! will also be created
 
631
           call iofa( na_u, cfa , .true. )
 
632
        end if
631
633
 
632
634
        write(6,'(40("-"),/,a6,3f12.6)') 'Tot',(ftot(ix)*Ang/eV,ix=1,3)
633
635
        write(6,'(40("-"),/,a6, f12.6)') 'Max',fmax*Ang/eV
635
637
     .       '    sqrt( Sum f_i^2 / 3N )'
636
638
        write(6,'(40("-"),/,a6, f12.6,a)') 'Max',cfmax*Ang/eV, 
637
639
     .       '    constrained'
 
640
 
 
641
        ! Write Force Constant matrix if FC calculation ...
 
642
        select case (idyn)
 
643
        case(6)
 
644
           ! If the istep is the first step, then it
 
645
           ! must be the first 
 
646
           call ofc(fa,dx,na_u,.false.,istep==inicoor)
 
647
           if ( any(fa/=cfa) ) then
 
648
              call ofc(cfa,dx,na_u,.true.,istep==inicoor)
 
649
           end if
 
650
        case(7)
 
651
!          call phonon_write_forces(fa,na_u,ucell,istep)
 
652
           if (IOnode) write(*,*) 'phonon support deactivated'
 
653
        end select
 
654
 
638
655
      else !not final
639
656
! In finalization, only print forces if sufficiently large.
640
657
        if (fmax .gt. ftol) then