Version 5.2.0, Revision 23132, Hash 1612051

MODIFIED *  configure include/version/version.m4 interface/INIT_activate.F interface/INIT_load.F modules/SET_defaults.F modules/mod_real_time.F real_time_drivers/RT_driver.F real_time_initialize/.objects real_time_propagation/.objects real_time_propagation/Runge_Kutta_driver.F real_time_propagation/Runge_Kutta_initialize.F

NEW *  real_time_initialize/Runge_Kutta_step_and_order_estimation.F real_time_propagation/RT_G_symmetrization.F

DELETED *  real_time_initialize/Runge_Kutta_optimization.F real_time_propagation/RT_Glesser_evolve.F

This branch is relative to issue #126

Additions:
- RK optimization procedure completed. Now the code can propose a RK order + dT based on the integration of a suitable test function.
  Notes to be added on Overleaf.
- Added the G_lesser symmetrization

Patch sent by:  Andrea Marini <[email protected]>

configure

@@ -1,6 +1,6 @@
 #! /bin/sh
 # Guess values for system-dependent variables and create Makefiles.
-# Generated by GNU Autoconf 2.71 for Yambo 5.2.0 r.23132 h.d2d5779d5.
+# Generated by GNU Autoconf 2.71 for Yambo 5.2.0 r.23133 h.161205154.
 #
 # Report bugs to <[email protected]>.
 #
@@ -610,8 +610,8 @@ MAKEFLAGS=
 # Identity of this package.
 PACKAGE_NAME='Yambo'
 PACKAGE_TARNAME='yambo'
-PACKAGE_VERSION='5.2.0 r.23132 h.d2d5779d5'
-PACKAGE_STRING='Yambo 5.2.0 r.23132 h.d2d5779d5'
+PACKAGE_VERSION='5.2.0 r.23133 h.161205154'
+PACKAGE_STRING='Yambo 5.2.0 r.23133 h.161205154'
 PACKAGE_BUGREPORT='[email protected]'
 PACKAGE_URL=''
 
@@ -1600,7 +1600,7 @@ if test "$ac_init_help" = "long"; then
   # Omit some internal or obsolete options to make the list less imposing.
   # This message is too long to be a string in the A/UX 3.1 sh.
   cat <<_ACEOF
-\`configure' configures Yambo 5.2.0 r.23132 h.d2d5779d5 to adapt to many kinds of systems.
+\`configure' configures Yambo 5.2.0 r.23133 h.161205154 to adapt to many kinds of systems.
 
 Usage: $0 [OPTION]... [VAR=VALUE]...
 
@@ -1666,7 +1666,7 @@ fi
 
 if test -n "$ac_init_help"; then
   case $ac_init_help in
-     short | recursive ) echo "Configuration of Yambo 5.2.0 r.23132 h.d2d5779d5:";;
+     short | recursive ) echo "Configuration of Yambo 5.2.0 r.23133 h.161205154:";;
    esac
   cat <<\_ACEOF
 
@@ -1876,7 +1876,7 @@ fi
 test -n "$ac_init_help" && exit $ac_status
 if $ac_init_version; then
   cat <<\_ACEOF
-Yambo configure 5.2.0 r.23132 h.d2d5779d5
+Yambo configure 5.2.0 r.23133 h.161205154
 generated by GNU Autoconf 2.71
 
 Copyright (C) 2021 Free Software Foundation, Inc.
@@ -2505,7 +2505,7 @@ cat >config.log <<_ACEOF
 This file contains any messages produced by compilers while
 running configure, to aid debugging if configure makes a mistake.
 
-It was created by Yambo $as_me 5.2.0 r.23132 h.d2d5779d5, which was
+It was created by Yambo $as_me 5.2.0 r.23133 h.161205154, which was
 generated by GNU Autoconf 2.71.  Invocation command line was
 
   $ $0$ac_configure_args_raw
@@ -3263,8 +3263,8 @@ ac_compiler_gnu=$ac_cv_c_compiler_gnu
 SVERSION="5"
 SSUBVERSION="2"
 SPATCHLEVEL="0"
-SREVISION="23132"
-SHASH="d2d5779d5"
+SREVISION="23133"
+SHASH="161205154"
 
 
 
@@ -16578,7 +16578,7 @@ cat >>$CONFIG_STATUS <<\_ACEOF || ac_write_fail=1
 # report actual input values of CONFIG_FILES etc. instead of their
 # values after options handling.
 ac_log="
-This file was extended by Yambo $as_me 5.2.0 r.23132 h.d2d5779d5, which was
+This file was extended by Yambo $as_me 5.2.0 r.23133 h.161205154, which was
 generated by GNU Autoconf 2.71.  Invocation command line was
 
   CONFIG_FILES    = $CONFIG_FILES
@@ -16642,7 +16642,7 @@ ac_cs_config_escaped=`printf "%s\n" "$ac_cs_config" | sed "s/^ //; s/'/'\\\\\\\\
 cat >>$CONFIG_STATUS <<_ACEOF || ac_write_fail=1
 ac_cs_config='$ac_cs_config_escaped'
 ac_cs_version="\\
-Yambo config.status 5.2.0 r.23132 h.d2d5779d5
+Yambo config.status 5.2.0 r.23133 h.161205154
 configured by $0, generated by GNU Autoconf 2.71,
   with options \\"\$ac_cs_config\\"
 

include/version/version.m4

@@ -1,9 +1,9 @@
-AC_INIT(Yambo, 5.2.0 r.23132 h.d2d5779d5, [email protected])
+AC_INIT(Yambo, 5.2.0 r.23133 h.161205154, [email protected])
 SVERSION="5"
 SSUBVERSION="2"
 SPATCHLEVEL="0"
-SREVISION="23132"
-SHASH="d2d5779d5"
+SREVISION="23133"
+SHASH="161205154"
 AC_SUBST(SVERSION)
 AC_SUBST(SSUBVERSION)
 AC_SUBST(SPATCHLEVEL)

src/interface/INIT_activate.F

@@ -532,7 +532,7 @@ subroutine INIT_activate()
    !
    call initactivate(1,'RTBands Integrator SlowDynamics LinearRegime GrKind TwoAlph RTehEny ScattTresh')
    call initactivate(1,'RADLifeTime RADmagnific PhLifeTime DephCVonly DephEThresh DephTRange')
-   call initactivate(1,'RTstep NEsteps NETime DipoleEtresh SPINprojected')
+   call initactivate(1,'RTstep RKquestTresh NEsteps NETime DipoleEtresh SPINprojected')
    call initactivate(1,'IOtime IOCachetime')
    !
    ! n_ext_fields is by default 0.  It is set to the input number of fields when a command line option is provided (below in init.F)

src/interface/INIT_load.F

@@ -68,12 +68,12 @@ subroutine INIT_load(defs,en,q,k,X,Xw,Dip)
  use stderr,        ONLY:intc
  use RT_occupations,ONLY:RT_RTA_temp,RT_RTA_chem
  use RT_lifetimes,  ONLY:RT_life_extrap_times_INPUT
- use real_time,     ONLY:RT_step, Integrator_name, RAD_LifeTime, Phase_LifeTime, &
+ use real_time,     ONLY:RT_step,Integrator_name, RAD_LifeTime, Phase_LifeTime, &
 &                        NE_tot_time,RT_bands,RT_eh_en,RT_scatt_tresh,Gr_kind, Integrator_slow_approx,     &
 &                        NE_MEM_treshold,NE_time_step_update_jump_INPUT,         &
 &                        two_alpha,DbGd_EE_percent,NE_time_step_update_jump,     &
 &                        NE_initial_time_step_update,NE_step_update_treshold,    &
-&                        RT_MAX_step,RAD_magnification,  &
+&                        RT_MAX_step,RAD_magnification,RK_quest_treshold,  &
 &                        RT_deph_deg_thresh,RT_step_manual_prefactor,NE_time_step_update_manual
  use fields,        ONLY:Efield,n_ext_fields_max,n_ext_fields
  use plasma,        ONLY:PLASMA_redux_percent
@@ -411,15 +411,16 @@ subroutine INIT_load(defs,en,q,k,X,Xw,Dip)
  call it(defs,'TwoAlpha',   '[RT] C_nk ~ alpha*Gamma_nk^2 ',two_alpha,Verb_level=V_real_time)
  call it(defs,'GrKind',     '[RT] G-ret kind: Lorentzian (QP)/ Hyperbolic QP_secant (HS)',Gr_kind,Verb_level=V_real_time)
  call it(defs,'RADLifeTime','[RT] Radiative life-time (if negative Yambo sets it equal to Phase_LifeTime in NL)', &
-&         RAD_LifeTime,unit=Time_unit(1))
+&                            RAD_LifeTime,unit=Time_unit(1),Verb_level=V_real_time)
  call it(defs,'RADmagnific','[RT] Radiative life-time magnification',RAD_magnification,Verb_level=V_real_time)
- call it(defs,'PhLifeTime', '[RT] Constant Dephasing Time',Phase_LifeTime,unit=Time_unit(1))
+ call it(defs,'PhLifeTime', '[RT] Constant Dephasing Time',Phase_LifeTime,unit=Time_unit(1),Verb_level=V_real_time)
  call it(defs,'DephEThresh', '[RT] Threshold on the energy difference between two states to dephase them',&
-         &   RT_deph_deg_thresh,unit=E_unit,Verb_level=V_real_time)
+&                            RT_deph_deg_thresh,unit=E_unit,Verb_level=V_real_time)
  !
  ! Dynamics
  !
- call it(defs,'RTstep',      '[RT] Real Time step length',RT_step,unit=Time_unit(3))
+ call it(defs,'RTstep',      '[RT] Real Time step: >0 absolute, <0 used as seed for the Runge-Kutta optimization procedure.',&
+&                             RT_step,unit=Time_unit(3))
  call it(defs,'NETime',      '[RT] Simulation Time',NE_tot_time,unit=Time_unit(2))
  call it(defs,'dTupdateTimeSet','[RT] Time for manual deltaT update',NE_time_step_update_manual,&
 &                            unit=Time_unit(1),verb_level=V_real_time)
@@ -428,11 +429,12 @@ subroutine INIT_load(defs,en,q,k,X,Xw,Dip)
  call it(defs,'dTupdateJump','[RT] Time betweem two deltaT updates',  NE_time_step_update_jump_INPUT,&
 &                            unit=Time_unit(1),verb_level=V_real_time)
  call it(defs,'dTupdateTresh','[RT][o/o] Treshold of deltaT updates',  NE_step_update_treshold ,verb_level=V_real_time)
+ call it(defs,'RKquestTresh','[RT][o/oo] Treshold used to identify the best dT/integrator',RK_quest_treshold,verb_level=V_real_time)
  call it(defs,'dT_MAX',      '[RT] Maximum value for the time-dependent dT',  RT_MAX_step ,verb_level=V_real_time,unit=Time_unit(1))
  call it(defs,'dT_SET',      '[RT] Prefactor for manual dT update',  RT_step_manual_prefactor,verb_level=V_real_time)
- call it(defs,'Integrator',  '[RT] Time-evolution (TE) integrator: "EULER", "RK2/RK4", "HEUN"',Integrator_name)
+ call it(defs,'Integrator',  '[RT] Time-evolution RK integrator: "EULER/EMP/RK2/RK3/H3/RK42/RK41/best"',Integrator_name)
  call it(defs,'SlowDynamics','[RT] Slow-TE approximations: "EXPn" (n<=6), "EXP" (n=3), "INV+ACC", "EXP+ACC", "INV+DIAG"',&
-&        Integrator_slow_approx,verb_level=V_real_time)
+&                             Integrator_slow_approx,verb_level=V_real_time)
  call it('f',defs,'LinearRegime','[RT] Adopt specific tricks in the EULER case assuming the linear-regime',verb_level=V_real_time)
  !
  IO_times=(/CARR_RT_IO_t%INTERVAL_time_INPUT,Gless_RESTART_RT_IO_t%INTERVAL_time_INPUT,OUTPUT_RT_IO_t%INTERVAL_time_INPUT/)

src/modules/SET_defaults.F

@@ -604,7 +604,7 @@ subroutine SET_defaults(INSTR,IND,OD,COM_DIR)
  NE_time_step_update_jump_INPUT  = 0._SP
  NE_initial_time_step_update    = 0._SP
  NE_step_update_treshold        = 0.1_SP
- NE_time_step_update_manual     =-1.0_SP
+ NE_time_step_update_manual     =-1.0_SP*FS2AUT
  !
  CACHE_OBS_steps      = 0
  CACHE_OBS_steps_now  = 0
@@ -624,7 +624,7 @@ subroutine SET_defaults(INSTR,IND,OD,COM_DIR)
  RAD_LifeTime = -1._SP*FS2AUT  ! Automatic set equal to the dephasing
 #endif
  Phase_LifeTime = 0._SP
- RT_deph_deg_thresh = 1.E-5_SP
+ RT_deph_deg_thresh = 1.E-5_SP/HA2EV
  RT_deph_cv_only =.FALSE.
  !
  NE_MEM_treshold=0._SP

src/modules/mod_real_time.F

@@ -224,6 +224,7 @@ module real_time
  integer, parameter     :: N_RK_integrators=7
  type(RK_integrator_t)  :: RK_integrators(N_RK_integrators)
  type(RK_integrator_t)  :: RK_integrator
+ real(SP)               :: RK_quest_treshold=1._SP
  !
  ! ... Slowly varying dynamics
  logical  :: l_slow_TE=.FALSE. ! Slow time-evolution => INV/EXPn/LINEAR_RESP approximations
@@ -239,12 +240,6 @@ module real_time
  ! ... INV related
  logical  :: l_RT_DIAG=.FALSE. ! Used in INV (RT_INV_step_diago)
  !
- !
- ! ... Coefficients in the implicit/explicit propagation procedure
- real(SP) :: a_tableau(4)
- real(SP) :: b_tableau(4)
- real(SP) :: c_tableau(4)
- !
  ! ... names (referring to RKn) and approximation (Physics based)
  character(schlen)            :: Integrator_name
  character(schlen)            :: Integrator_slow_approx='no'
@@ -368,6 +363,7 @@ subroutine RK_integrator_alloc(I,N,short,name)
      I%short=short
      I%name=name
      I%N=N
+     if (allocated(I%a)) return
      YAMBO_ALLOC(I%a,(I%N,I%N))
      YAMBO_ALLOC(I%b,(I%N))
      YAMBO_ALLOC(I%c,(I%N))

src/real_time_drivers/RT_driver.F

@@ -131,7 +131,7 @@ subroutine RT_driver(E,X,k,q,Dip)
  !
  call RT_apply(RT_bands,E,k,what="G",VERBOSE=.true.)
  !
- if (STRING_same(Integrator_name,"best")) call Runge_Kutta_optimization(E)
+ if (RT_step<0._SP.or.STRING_same(Integrator_name,"best")) call Runge_Kutta_step_and_order_estimation(E)
  !
  call section('=','Initialization')
  !#################################

src/real_time_initialize/.objects

@@ -2,7 +2,7 @@
 ELPH_objs = RT_ELPH_initialize.o 
 #endif
 #if defined _RT
-RT_head_objs = RT_initialize.o RT_start_and_restart.o Runge_Kutta_optimization.o
+RT_head_objs = RT_initialize.o RT_start_and_restart.o Runge_Kutta_step_and_order_estimation.o
 RT_foot_objs = RT_Field_Commensurable_Frequencies.o RT_Dephasing_Matrix.o RT_G_lesser_init.o RT_occupations_update.o
 #endif
 objs = $(RT_head_objs) RT_occupations_and_levels_init.o $(RT_foot_objs) $(ELPH_objs)