Merge branch 'wstat_restart' into 'develop'

Allow user to specify how often wstat writes its restart info

See merge request west-devel/West!32

Doc/manual.rst

@@ -159,6 +159,23 @@ wstat_control
      - Number of PDEP eigenpotentials that can be read from file.
 
 
+.. data:: n_steps_write_restart
+
+.. list-table::
+   :widths: 10 90
+   :stub-columns: 0
+
+   * - **Type**
+     - int
+   * - **Default**
+     - 1
+   * - **Description**
+     - Available options are:
+
+       - If ( n_steps_write_restart >  0 ) A checkpoint is written every n_steps_write_restart iterations in the PDEP loop.
+       - If ( n_steps_write_restart <= 0 ) A checkpoint is NEVER written in the PDEP loop. Restart will not be possible.
+
+
 .. data:: trev_pdep
 
 .. list-table::

Pytools/west_fetch_input.py

@@ -36,6 +36,7 @@ default["wstat_control"]["n_pdep_times"] = 4
 default["wstat_control"]["n_pdep_maxiter"] = 100
 default["wstat_control"]["n_dfpt_maxiter"] = 250
 default["wstat_control"]["n_pdep_read_from_file"] = 0
+default["wstat_control"]["n_steps_write_restart"] = 1
 default["wstat_control"]["trev_pdep"] = 1.e-3
 default["wstat_control"]["trev_pdep_rel"] = 1.e-1
 default["wstat_control"]["tr2_dfpt"] = 1.e-12

Tools/fetch_input.f90

@@ -50,6 +50,7 @@ SUBROUTINE add_intput_parameters_to_json_file( num_drivers, driver, json )
         CALL json%add('input.wstat_control.n_pdep_maxiter',n_pdep_maxiter)
         CALL json%add('input.wstat_control.n_dfpt_maxiter',n_dfpt_maxiter)
         CALL json%add('input.wstat_control.n_pdep_read_from_file',n_pdep_read_from_file)
+        CALL json%add('input.wstat_control.n_steps_write_restart',n_steps_write_restart)
         CALL json%add('input.wstat_control.trev_pdep',trev_pdep)
         CALL json%add('input.wstat_control.trev_pdep_rel',trev_pdep_rel)
         CALL json%add('input.wstat_control.tr2_dfpt',tr2_dfpt)
@@ -213,6 +214,7 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
         IERR = return_dict%get(n_pdep_maxiter, "n_pdep_maxiter", DUMMY_DEFAULT)
         IERR = return_dict%get(n_dfpt_maxiter, "n_dfpt_maxiter", DUMMY_DEFAULT)
         IERR = return_dict%get(n_pdep_read_from_file, "n_pdep_read_from_file", DUMMY_DEFAULT)
+        IERR = return_dict%get(n_steps_write_restart, "n_steps_write_restart", DUMMY_DEFAULT)
         IERR = return_dict%getitem(trev_pdep, "trev_pdep")
         IERR = return_dict%getitem(trev_pdep_rel, "trev_pdep_rel")
         IERR = return_dict%getitem(tr2_dfpt, "tr2_dfpt")
@@ -368,6 +370,7 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
      CALL mp_bcast(n_pdep_maxiter,root,world_comm)
      CALL mp_bcast(n_dfpt_maxiter,root,world_comm)
      CALL mp_bcast(n_pdep_read_from_file,root,world_comm)
+     CALL mp_bcast(n_steps_write_restart,root,world_comm)
      CALL mp_bcast(trev_pdep,root,world_comm)
      CALL mp_bcast(trev_pdep_rel,root,world_comm)
      CALL mp_bcast(tr2_dfpt,root,world_comm)
@@ -399,6 +402,7 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
      IF( n_pdep_maxiter == DUMMY_DEFAULT ) CALL errore('fetch_input','Err: cannot read n_pdep_maxiter')
      IF( n_dfpt_maxiter == DUMMY_DEFAULT ) CALL errore('fetch_input','Err: cannot read n_dfpt_maxiter')
      IF( n_pdep_read_from_file == DUMMY_DEFAULT ) CALL errore('fetch_input','Err: cannot read n_pdep_read_from_file')
+     IF( n_steps_write_restart == DUMMY_DEFAULT ) CALL errore('fetch_input','Err: cannot read n_steps_write_restart')
      IF(gamma_only) THEN
         IF (SIZE(qlist)/=1) CALL errore('fetch_input','Err: SIZE(qlist)/=1.',1)
      ELSE 
@@ -520,13 +524,14 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
         !
         CALL io_push_title('I/O Summary : wstat_control')
         !
-        numsp=30
+        numsp = 30
         CALL io_push_value('wstat_calculation',wstat_calculation,numsp)
         CALL io_push_value('n_pdep_eigen',n_pdep_eigen,numsp)
         CALL io_push_value('n_pdep_times',n_pdep_times,numsp)
         CALL io_push_value('n_pdep_maxiter',n_pdep_maxiter,numsp)
         CALL io_push_value('n_dfpt_maxiter',n_dfpt_maxiter,numsp)
         CALL io_push_value('n_pdep_read_from_file',n_pdep_read_from_file,numsp)
+        CALL io_push_value('n_steps_write_restart',n_steps_write_restart,numsp)
         CALL io_push_es0('trev_pdep',trev_pdep,numsp)
         CALL io_push_es0('trev_pdep_rel',trev_pdep_rel,numsp)
         CALL io_push_es0('tr2_dfpt',tr2_dfpt,numsp)
@@ -547,7 +552,7 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
         !
         CALL io_push_title('I/O Summary : wfreq_control')
         !
-        numsp=40
+        numsp = 40
         CALL io_push_value('wfreq_calculation',wfreq_calculation,numsp)
         CALL io_push_value('n_pdep_eigen_to_use',n_pdep_eigen_to_use,numsp)
         CALL io_push_value('qp_bandrange(1)',qp_bandrange(1),numsp)
@@ -578,7 +583,7 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
         !
         CALL io_push_title('I/O Summary : westpp_control')
         !
-        numsp=40
+        numsp = 40
         CALL io_push_value('westpp_calculation',westpp_calculation,numsp)
         CALL io_push_value('westpp_range(1)',westpp_range(1),numsp)
         CALL io_push_value('westpp_range(2)',westpp_range(2),numsp)
@@ -602,7 +607,7 @@ SUBROUTINE fetch_input_yml( num_drivers, driver, verbose, debug )
         !
         CALL io_push_title('I/O Summary : server_control')
         !
-        numsp=40
+        numsp = 40
         CALL io_push_value('document',document,numsp)
         !
         CALL io_push_bar()

Wstat/davidson_diago.f90

 !
-! Copyright (C) 2015-2017 M. Govoni 
+! Copyright (C) 2015-2017 M. Govoni
 ! This file is distributed under the terms of the
 ! GNU General Public License. See the file `License'
 ! in the root directory of the present distribution,
@@ -7,7 +7,7 @@
 !
 ! This file is part of WEST.
 !
-! Contributors to this file: 
+! Contributors to this file:
 ! Marco Govoni
 !
 #define ZERO ( 0.D0, 0.D0 )
@@ -17,17 +17,17 @@
 SUBROUTINE davidson_diago ( )
   !----------------------------------------------------------------------------
   !
-  USE control_flags,   ONLY : gamma_only 
+  USE control_flags,   ONLY : gamma_only
   !
   IMPLICIT NONE
   !
-  IF( gamma_only ) THEN 
-    CALL davidson_diago_gamma( ) 
+  IF( gamma_only ) THEN
+    CALL davidson_diago_gamma( )
   ELSE
     CALL davidson_diago_k( )
-  ENDIF 
+  ENDIF
   !
-END SUBROUTINE 
+END SUBROUTINE
 !
 !----------------------------------------------------------------------------
 SUBROUTINE davidson_diago_gamma ( )
@@ -51,7 +51,7 @@ SUBROUTINE davidson_diago_gamma ( )
   USE mp,                   ONLY : mp_sum
   USE gvect,                ONLY : gstart
   USE wstat_tools,          ONLY : diagox,serial_diagox,build_hr,symm_hr_distr,redistribute_vr_distr,&
-                                   & update_with_vr_distr,refresh_with_vr_distr 
+                                   & update_with_vr_distr,refresh_with_vr_distr
   USE types_coulomb,        ONLY : pot3D
   !
   IMPLICIT NONE
@@ -98,8 +98,8 @@ SUBROUTINE davidson_diago_gamma ( )
   !
   pert = idistribute()
   CALL pert%init(nvecx,'i','nvecx',.TRUE.)
-  !CALL index_distr_init(pert,nvecx,'i','nvecx',.TRUE.) 
-  CALL wstat_memory_report() ! Before allocating I report the memory required. 
+  !CALL index_distr_init(pert,nvecx,'i','nvecx',.TRUE.)
+  CALL wstat_memory_report() ! Before allocating I report the memory required.
   !
   ! ... MEMORY ALLOCATION
   !
@@ -141,7 +141,7 @@ SUBROUTINE davidson_diago_gamma ( )
   dvg  = 0._DP
   hr_distr(:,:) = 0._DP
   vr_distr(:,:) = 0._DP
-  notcnv  = nvec 
+  notcnv  = nvec
   dav_iter = -2
   !
   CALL pot3D%init(fftdriver,.FALSE.,'default')
@@ -183,7 +183,7 @@ SUBROUTINE davidson_diago_gamma ( )
      ! Apply operator with DFPT
      !
      mloc = 0
-     mstart = 1 
+     mstart = 1
      DO il1 = 1, pert%nloc
         ig1 = pert%l2g(il1)
         IF( ig1 < 1 .OR. ig1 > nvec ) CYCLE
@@ -192,17 +192,17 @@ SUBROUTINE davidson_diago_gamma ( )
      ENDDO
      !
      pccg_res_tr2 = -1._DP
-     CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, 1 ) 
+     CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, 1 )
      dav_iter = -1
-     CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr)
+     IF(n_steps_write_restart == 1) CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr )
      !
   ENDIF
   !
-  IF( dav_iter == -2 ) CALL errore( 'chidiago','Cannot find the 1st starting loop',1) 
+  IF( dav_iter == -2 ) CALL errore( 'chidiago','Cannot find the 1st starting loop',1)
   !
   IF( dav_iter == -1 ) THEN
      !
-     ! < EXTRA STEP > 
+     ! < EXTRA STEP >
      !
      dvg = dng
      CALL do_mgs( dvg, 1, nvec)
@@ -218,7 +218,7 @@ SUBROUTINE davidson_diago_gamma ( )
      ! Apply operator with DFPT
      !
      mloc = 0
-     mstart = 1 
+     mstart = 1
      DO il1 = 1, pert%nloc
         ig1 = pert%l2g(il1)
         IF( ig1 < 1 .OR. ig1 > nvec ) CYCLE
@@ -228,12 +228,12 @@ SUBROUTINE davidson_diago_gamma ( )
      !
      pccg_res_tr2 = MIN(0.01_DP,1000000._DP*tr2_dfpt)
      CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, 1 )
-     ! 
+     !
      ! </ EXTRA STEP >
      !
      ! hr = <dvg|dng>
      !
-     CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, 1, nvec ) 
+     CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, 1, nvec )
      !
      ! ... diagonalize the reduced hamiltonian
      !
@@ -246,7 +246,7 @@ SUBROUTINE davidson_diago_gamma ( )
      !
      CALL output_ev_and_time(nvec,ev,conv,time_spent,sternop_ncalls,pccg_res_tr2,dfpt_dim,diago_dim,dav_iter,notcnv)
      dav_iter = 0
-     CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr)
+     IF(n_steps_write_restart == 1) CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr )
      !
   ENDIF
   !
@@ -279,19 +279,19 @@ SUBROUTINE davidson_diago_gamma ( )
         !
         IF ( .NOT. conv(n) ) THEN
            !
-           ! ... this root not yet converged ... 
+           ! ... this root not yet converged ...
            !
            np = np + 1
            !
            ! ... reorder eigenvectors so that coefficients for unconverged
-           ! ... roots come first. This allows to use quick matrix-matrix 
+           ! ... roots come first. This allows to use quick matrix-matrix
            ! ... multiplications to set a new basis vector (see below)
            !
            !IF ( np /= n ) vr(:,np) = vr(:,n)
            ishift(nbase+np) = n
            !
            ew(nbase+np) = ev(n)
-           !   
+           !
         END IF
         !
      END DO
@@ -310,7 +310,7 @@ SUBROUTINE davidson_diago_gamma ( )
      !
      ! determine image that actually compute dng first
      !
-     mloc = 0 
+     mloc = 0
      mstart = 1
      DO il1 = 1, pert%nloc
         ig1 = pert%l2g(il1)
@@ -322,14 +322,14 @@ SUBROUTINE davidson_diago_gamma ( )
      ! Apply operator with DFPT
      !
      pccg_res_tr2 = tr2_dfpt
-     CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, 1) 
+     CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, 1)
      !
      ! ... update the reduced hamiltonian
      !
      !
      ! hr = <dvg|dng>
      !
-     CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, nbase+1, nbase+notcnv ) 
+     CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, nbase+1, nbase+notcnv )
      !
      nbase = nbase + notcnv
      !
@@ -378,7 +378,7 @@ SUBROUTINE davidson_diago_gamma ( )
            !CALL output_a_report(-1)
            !
            WRITE(iter_label,'(i8)') kter
-           CALL io_push_title("Convergence achieved !!! in "//TRIM(iter_label)//" steps") 
+           CALL io_push_title("Convergence achieved !!! in "//TRIM(iter_label)//" steps")
            l_is_wstat_converged = .TRUE.
            !
            EXIT iterate
@@ -412,7 +412,7 @@ SUBROUTINE davidson_diago_gamma ( )
         !
         DO il1 = 1, pert%nloc
            ig1 = pert%l2g(il1)
-           IF( ig1 > nbase ) CYCLE 
+           IF( ig1 > nbase ) CYCLE
            hr_distr(ig1,il1) = ev(ig1)
            vr_distr(ig1,il1) = 1._DP
         ENDDO
@@ -421,7 +421,8 @@ SUBROUTINE davidson_diago_gamma ( )
         !
      END IF
      !
-     CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr)
+     IF(n_steps_write_restart > 0 .AND. MOD(dav_iter,n_steps_write_restart) == 0) &
+        CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr )
      !CALL output_a_report(dav_iter)
      !
   END DO iterate
@@ -435,7 +436,7 @@ SUBROUTINE davidson_diago_gamma ( )
   !
   !
   DEALLOCATE( dng )
-  DEALLOCATE( dvg )  
+  DEALLOCATE( dvg )
   !
   CALL stop_clock( 'chidiago' )
   !
@@ -467,7 +468,7 @@ SUBROUTINE davidson_diago_k ( )
   USE gvect,                ONLY : gstart, g, ngm
   USE gvecw,                ONLY : gcutw
   USE wstat_tools,          ONLY : diagox,serial_diagox,build_hr,symm_hr_distr,redistribute_vr_distr,&
-                                   & update_with_vr_distr,refresh_with_vr_distr 
+                                   & update_with_vr_distr,refresh_with_vr_distr
   USE types_bz_grid,        ONLY : q_grid
   USE types_coulomb,        ONLY : pot3D
   !
@@ -518,8 +519,8 @@ SUBROUTINE davidson_diago_k ( )
   !
   pert=idistribute()
   CALL pert%init(nvecx,'i','nvecx',.TRUE.)
-  !CALL index_distr_init(pert,nvecx,'i','nvecx',.TRUE.) 
-  CALL wstat_memory_report() ! Before allocating I report the memory required. 
+  !CALL index_distr_init(pert,nvecx,'i','nvecx',.TRUE.)
+  CALL wstat_memory_report() ! Before allocating I report the memory required.
   !
   ! ... MEMORY ALLOCATION
   !
@@ -568,7 +569,7 @@ SUBROUTINE davidson_diago_k ( )
      dvg  = 0._DP
      hr_distr(:,:) = 0._DP
      vr_distr(:,:) = 0._DP
-     notcnv  = nvec 
+     notcnv  = nvec
      dav_iter = -2
      !
      ! set local number of G vectors for perturbation at q
@@ -602,7 +603,7 @@ SUBROUTINE davidson_diago_k ( )
      !
      ! KIND OF CALCULATION
      !
-     IF( wstat_calculation(1:1) == "R" .OR. wstat_calculation(2:2) == "R" ) THEN 
+     IF( wstat_calculation(1:1) == "R" .OR. wstat_calculation(2:2) == "R" ) THEN
         !
         IF ( .NOT. l_restart_q_done ) THEN
            !
@@ -626,9 +627,9 @@ SUBROUTINE davidson_diago_k ( )
         !
         ! ... Eventually read from file
         !
-        IF (iq==1) THEN 
+        IF (iq==1) THEN
            l_print_pdep_read = .TRUE.
-        ELSE 
+        ELSE
            l_print_pdep_read = .FALSE.
         ENDIF
         IF(n_pdep_read_from_file>0) CALL pdep_db_read( n_pdep_read_from_file, iq, l_print_pdep_read)
@@ -652,7 +653,7 @@ SUBROUTINE davidson_diago_k ( )
         ! Apply operator with DFPT
         !
         mloc = 0
-        mstart = 1 
+        mstart = 1
         DO il1 = 1, pert%nloc
            ig1 = pert%l2g(il1)
            IF( ig1 < 1 .OR. ig1 > nvec ) CYCLE
@@ -662,17 +663,17 @@ SUBROUTINE davidson_diago_k ( )
         !
         pccg_res_tr2 = -1._DP
         !
-        CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, iq ) 
+        CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, iq )
         dav_iter = -1
-        CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr, iq)
+        IF(n_steps_write_restart == 1) CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr, iq )
         !
      ENDIF
      !
-     IF( dav_iter == -2 ) CALL errore( 'chidiago','Cannot find the 1st starting loop',1) 
+     IF( dav_iter == -2 ) CALL errore( 'chidiago','Cannot find the 1st starting loop',1)
      !
      IF( dav_iter == -1 ) THEN
         !
-        ! < EXTRA STEP > 
+        ! < EXTRA STEP >
         !
         dvg = dng
         CALL do_mgs( dvg, 1, nvec)
@@ -688,7 +689,7 @@ SUBROUTINE davidson_diago_k ( )
         ! Apply operator with DFPT
         !
         mloc = 0
-        mstart = 1 
+        mstart = 1
         DO il1 = 1, pert%nloc
            ig1 = pert%l2g(il1)
            IF( ig1 < 1 .OR. ig1 > nvec ) CYCLE
@@ -698,13 +699,13 @@ SUBROUTINE davidson_diago_k ( )
         !
         pccg_res_tr2 = MIN(0.01_DP,1000000._DP*tr2_dfpt)
         !
-        CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, iq ) 
-        ! 
+        CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, iq )
+        !
         ! </ EXTRA STEP >
         !
         ! hr = <dvg|dng>
         !
-        CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, 1, nvec ) 
+        CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, 1, nvec )
         !
         ! ... diagonalize the reduced hamiltonian
         !
@@ -717,7 +718,7 @@ SUBROUTINE davidson_diago_k ( )
         !
         CALL output_ev_and_time_q( nvec,ev,conv,time_spent,sternop_ncalls,pccg_res_tr2,dfpt_dim,diago_dim,dav_iter,notcnv,iq )
         dav_iter = 0
-        CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr, iq )
+        IF(n_steps_write_restart == 1) CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr, iq )
         !
      ENDIF
      !
@@ -750,19 +751,19 @@ SUBROUTINE davidson_diago_k ( )
            !
            IF ( .NOT. conv(n) ) THEN
               !
-              ! ... this root not yet converged ... 
+              ! ... this root not yet converged ...
               !
               np = np + 1
               !
               ! ... reorder eigenvectors so that coefficients for unconverged
-              ! ... roots come first. This allows to use quick matrix-matrix 
+              ! ... roots come first. This allows to use quick matrix-matrix
               ! ... multiplications to set a new basis vector (see below)
               !
               !IF ( np /= n ) vr(:,np) = vr(:,n)
               ishift(nbase+np) = n
               !
               ew(nbase+np) = ev(n)
-              !   
+              !
            END IF
            !
         END DO
@@ -781,7 +782,7 @@ SUBROUTINE davidson_diago_k ( )
         !
         ! determine image that actually compute dng first
         !
-        mloc = 0 
+        mloc = 0
         mstart = 1
         DO il1 = 1, pert%nloc
            ig1 = pert%l2g(il1)
@@ -794,14 +795,14 @@ SUBROUTINE davidson_diago_k ( )
         !
         pccg_res_tr2 = tr2_dfpt
         !
-        CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, iq ) 
+        CALL apply_operator ( mloc, dvg(1,mstart), dng(1,mstart), pccg_res_tr2, iq )
         !
         ! ... update the reduced hamiltonian
         !
         !
         ! hr = <dvg|dng>
         !
-        CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, nbase+1, nbase+notcnv ) 
+        CALL build_hr( dvg, dng, mstart, mstart+mloc-1, hr_distr, nbase+1, nbase+notcnv )
         !
         nbase = nbase + notcnv
         !
@@ -850,7 +851,7 @@ SUBROUTINE davidson_diago_k ( )
               !CALL output_a_report(-1)
               !
               WRITE(iter_label,'(i8)') kter
-              CALL io_push_title("Convergence achieved !!! in "//TRIM(iter_label)//" steps") 
+              CALL io_push_title("Convergence achieved !!! in "//TRIM(iter_label)//" steps")
               l_is_wstat_converged = .TRUE.
               !
               EXIT iterate
@@ -893,13 +894,14 @@ SUBROUTINE davidson_diago_k ( )
            !
         END IF
         !
-        CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr, iq)
-        !!CALL output_a_report(dav_iter)
+        IF(n_steps_write_restart > 0 .AND. MOD(dav_iter,n_steps_write_restart) == 0) &
+           CALL wstat_restart_write( dav_iter, notcnv, nbase, ew, hr_distr, vr_distr, iq )
+        !CALL output_a_report(dav_iter)
         !
      END DO iterate
-  !
-  CALL stop_clock( 'chidiago' )
-  !
+     !
+     CALL stop_clock( 'chidiago' )
+     !
   ENDDO QPOINTS_LOOP ! iq
   !
   DEALLOCATE( conv )
@@ -920,19 +922,19 @@ END SUBROUTINE
 ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 !
-SUBROUTINE do_mgs(amat,m_global_start,m_global_end) 
+SUBROUTINE do_mgs(amat,m_global_start,m_global_end)
   !
   ! MGS of the vectors beloging to the interval [ m_global_start, m_global_end ]
   !    also with respect to the vectors belonging to the interval [ 1, m_global_start -1 ]
   !
   USE kinds,                  ONLY : DP
   USE io_global,              ONLY : stdout
-  USE mp_global,              ONLY : intra_bgrp_comm,inter_image_comm,my_image_id,nimage,world_comm 
+  USE mp_global,              ONLY : intra_bgrp_comm,inter_image_comm,my_image_id,nimage,world_comm
   USE gvect,                  ONLY : gstart
   USE mp,                     ONLY : mp_sum,mp_barrier,mp_bcast
   USE westcom,                ONLY : npwq,npwqx
   USE control_flags,          ONLY : gamma_only
-  USE io_push,                ONLY : io_push_title 
+  USE io_push,                ONLY : io_push_title
   USE distribution_center,    ONLY : pert
   !
   IMPLICIT NONE
@@ -942,7 +944,7 @@ SUBROUTINE do_mgs(amat,m_global_start,m_global_end)
   INTEGER, INTENT(IN) :: m_global_start,m_global_end
   COMPLEX(DP) :: amat(npwqx,pert%nlocx)
   !
-  ! Workspace 
+  ! Workspace
   !
   INTEGER :: ig, ip, j, ncol
   INTEGER :: k_global, k_local, j_local, k_id
@@ -970,11 +972,11 @@ SUBROUTINE do_mgs(amat,m_global_start,m_global_end)
   !
   ! 2) Localize m_global_start
   !
-  m_local_start = 1 
+  m_local_start = 1
   DO ip = 1, pert%nloc
      ig = pert%l2g(ip)
      IF( ig < m_global_start ) CYCLE
-     m_local_start = ip 
+     m_local_start = ip
      EXIT
   ENDDO
   !
@@ -1046,8 +1048,8 @@ SUBROUTINE do_mgs(amat,m_global_start,m_global_end)
            DO ip = j_local,m_local_end !pert%nloc
               braket(ip) = 2._DP * DDOT(2*npwq,vec(1),1,amat(1,ip),1)
            ENDDO
-           !IF (gstart==2) FORALL( ip=j_local:pert%nloc ) braket(ip) = braket(ip) - REAL(vec(1),KIND=DP)*REAL(amat(1,ip),KIND=DP) 
-           IF (gstart==2) FORALL( ip=j_local:m_local_end ) braket(ip) = braket(ip) - REAL(vec(1),KIND=DP)*REAL(amat(1,ip),KIND=DP) 
+           !IF (gstart==2) FORALL( ip=j_local:pert%nloc ) braket(ip) = braket(ip) - REAL(vec(1),KIND=DP)*REAL(amat(1,ip),KIND=DP)
+           IF (gstart==2) FORALL( ip=j_local:m_local_end ) braket(ip) = braket(ip) - REAL(vec(1),KIND=DP)*REAL(amat(1,ip),KIND=DP)
            !CALL mp_sum(braket(j_local:pert%nloc),intra_bgrp_comm)
            CALL mp_sum(braket(j_local:m_local_end),intra_bgrp_comm)
            !FORALL(ip=j_local:pert%nloc) zbraket(ip) = CMPLX( braket(ip), 0._DP, KIND=DP)
@@ -1080,7 +1082,7 @@ END SUBROUTINE
 ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 !
-SUBROUTINE do_randomize ( amat, mglobalstart, mglobalend  ) 
+SUBROUTINE do_randomize ( amat, mglobalstart, mglobalend  )
   !
   ! Randomize in dvg the vectors belonging to [ mglobalstart, mglobalend ]
   !
@@ -1154,7 +1156,7 @@ END SUBROUTINE
 ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 ! %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 !
-SUBROUTINE do_randomize_q (amat, mglobalstart, mglobalend, iq) 
+SUBROUTINE do_randomize_q (amat, mglobalstart, mglobalend, iq)
   !
   ! Randomize in dvg the vectors belonging to [ mglobalstart, mglobalend ]
   !
@@ -1214,7 +1216,7 @@ SUBROUTINE do_randomize_q (amat, mglobalstart, mglobalend, iq)
         IF ( qgnorm2 < eps8 ) CYCLE
         rr = random_num_debug(1,ig_l2g(igq_q(ig,iq)))
         arg = tpi * random_num_debug(2,ig_l2g(igq_q(ig,iq)))
-        amat(ig,il1) = CMPLX( rr*COS( arg ), rr*SIN( arg ), KIND=DP) / ( qgnorm2 + 1.0_DP ) 
+        amat(ig,il1) = CMPLX( rr*COS( arg ), rr*SIN( arg ), KIND=DP) / ( qgnorm2 + 1.0_DP )
      ENDDO
 !$OMP ENDDO
 !$OMP END PARALLEL
@@ -1254,20 +1256,20 @@ END SUBROUTINE
 !   converged = 0
 !   DO ip=1,n_pdep_eigen
 !      !out_tab(ip,1) = REAL(ip,DP)
-!      !ipert(ip) = ip 
+!      !ipert(ip) = ip
 !      !out_tab(ip,2) = ev(ip)
 !      !out_tab(ip,3) = 0._DP
 !      !IF(conv(ip)) out_tab(ip,3) = 1._DP
 !      IF(conv(ip)) converged(ip) = 1
 !   ENDDO
-!   IF(iteration>=0) THEN 
+!   IF(iteration>=0) THEN
 !      WRITE(pref,"('itr_',i5.5)") iteration
 !   ELSE
 !      pref='converged'
 !   ENDIF
 !   !CALL serial_table_output(mpime==root,4000,'wstat.'//TRIM(ADJUSTL(pref)),out_tab,n_pdep_eigen,3,(/'   iprt','eigenv.','  conv.'/))
 !   !
-!   IF( mpime == root ) THEN 
+!   IF( mpime == root ) THEN