fgrNRG

Usage: [[[om,] a0], Ifgr] = fgrNRG(NRG0, NRG2, Lambda, B [,C], [,opts])

   full density matrix (FDM) Numerical Renormalization Group (NRG)
   calculation for Fermi Golden Rule (FGR) rates, with extension
   to resonanant ineleastic x-ray scattering (RIXS).

FGR mode (default)

   This computes the FGR spectral data (om,a0) witin the fdm-NRG
   framework at arbitrary temperature (default T=TN^+ with N the
   length of the Wilson chain). When Fourier tansformed, it corresponds
   to the time-dependent correlation function [see Ifgr.Itdm.at]

      a(t) = < B(t)C'> = tr[ rho_0 * B(t) C' ]

   Here B and C are local operators acting at or in close vicinity
   of the impurity. They are assumed to switch back and forth
   between two disconnected sectors of Hilbert space (e.g. by the
   presence of an (otherwise unaccounted for) valence core hole.
   The time evolution C(t) = exp(iH0*t) C exp(-iH*t) therefore becomes
   `mixed' with `initial' Hamiltonian H0 and final Hamiltonian,
   iteratively diagonalized for both, NRG0 and NRG2, respetively.
   The frequency correponds to om=E_final-E_initial.

   If C is not specified or empy, it is automatically assumed
   equal to B. This is the typical situation in FGR.
   The inital state is determined by the full density matrix (FDM)
   rho_0 of NRG0 at given temperature (default T ~ T_N ~ 0^+K).
   The dynamics occurs within NRG2.

Input

   NRG*       (string pointer to) data of previous NRG runs;
              if this contains <path/file>, then the data will be
              read from the files <path/file>_##.mat files as
              generated by NRGWilsonQS.

   Lambda     NRG discretization parameter used in NRG*.
              While EScale and E0 will be read from NRG0_info,
              Lambda is needed for simple estimates of energy scales.

   B,C        local operators, e.g. acting on the impurity.  *)
              If operators B or C are specified by name, they will be
              read from workspace and stored in NRG2_info
              (this is disabled by `nostore').

   NB! Note that exchanging NRG0 <> NRG2 and B<>C corresponds
   to changing from emission to absorption spectra or vice versa.

Options / Flags

  'Eoffset',. additional energy offset which affects the omega
              binning (om-binning) only. Eoffset acts *on top* of
              E0_final-E0_initial, which ensures that by default,
              (Eoffset=0) omega=0 corresponds to minimal excitation
              energy for a transition from groundstate in NRG0
              to groundstate in NRG2 (with NRG2 considered (much)
              higher in energy in x-ray absorption).

  'T',...     temperature at which correlation function is calculated (0.)
  'NRho', k   builds density matrix from iteration k only (cf. `DM-NRG')

  'emin',...  min. energy for omega binning (TN/1000)  **)
  'emax',...  max. energy for omega binning (10.)
  'nlog',...  number of bins per decade for om-binning (256)

  'nostore'   do not store operator matrix elements of C[12] with NRG data
  'locRho'    if Rho needs to be calculated it is stored in RAM only
              and will be deleted after the program is finished

  'sigma',..  applies broadening in frequency domain (log.gauss)
              NB! either alpha or sigma can be specified but not both.
              default sigma=0.3 (with alpha ignored)
  'alpha',..  exponential decay of energies exp(-alpha dE) [0.]
              this applies alpha in time domain (F. Anders 2005)

  'partial'   return partial data from every NRG iteration
  'raw'       raw data only, i.e. no smoothing of data
  'RAW'       calculates time dependent data from exact frequency data
              i.e. without binning of energies on log. scale
  'disp'      more detaild output (0)

  'Z0',..     Z0 operator for first site (impurity) required for
              fermionic setup where Z0=(-1)^(number of particles)
  'zflags',.. on which operators in B.*C to apply fermionic signs Z0
              ex. 'zflags', [ 1 0 ... ]  would apply Z0 to the first
              but not to the second operator  (default: all 1s,
              as typically the case for absorption/emission setups.

RIXS mode

   With the flag '--rixs', this routine computes disrete, i.e., binned,
   RIXS (resonant ineleastic x-ray scattering) energy loss spetra
   for given incomming photon energy E_inc.
   This uses the Kramers-Heisenberg higher-order generalization of
   Fermi golden rule.

   Note on semantics: within RIXS, i(nitial) and f(inal) refer to the
   *same* Hamiltonian, wheras it is the intermediate Hamiltonian takes
   the role of the former H2:=H_final. With

       V_if := <i| B 1/[E_inc^* + E_i - H2] C' |f>

   the RIXS specral data would correspond, when Fourier transformed
   into time-domain, to

       R(t) = < V(t)*V'> = tr[ rho_0 * V(t) V' ]

   thus recovering the same final form as for absorption/emission
   spectra above. The RIXS mode requires the additional parameter

   'E_inc',.. incoming photon energy E_inc = complex(omega_in, gamma);
              this must already also include the broadening gamma
              related to finite life-time of the core-hole due to other
              unaccounted for relaxation processes such as Auger etc.

   The returned spectra correspond to the RIXS intensity,
   i.e., spectral intensity a0 vs. energy loss om.

Output

   om, a0     binned spectral data (according to emin,emax,nlog)
   Ifgr       info structure

Notes / comments

 *) Alternatively, if B and C have composite structure that oct on
    multiple sites/indices simultaneously, the operators B and C can
    also be specified as cell arrays of QSpace vectors.
    A given QSpace vector then, e.g., B{i}, is interpreted as the
    sequence of operators that act like a tensor product on indices
    L00,s00,s01,... where L00 and s00 are the `left' and local (`sigma)
    indices of A0, respectively, and s01, etc. the local indices
    of A01, etc. In this usage, Z0 must be already properly included
    in B and C upon input, i.e. Z0 and zflags cannot be used.

**) The bin a0(end) collects all data with omega >= om(end) = emax.
    For emin, the smallest bin a0(1) (given the logarthmic grid!)
    contains all data where omega is exactly zero (e.g. for diagonal
    contributions in Lehmann representation). The next smallest bin
    a0(2) then contains all data for 0 < omega <= om(2).
% 
Alternative usage:
fgrNRG('nrgdata0','nrgdata2','setup.mat','B','C')

    This reads all input parameters and options from given setup.mat
    file. Hence the operators B and C must be specified by their name.
    Output will be written / updated to NRG0 and NRG2,
    with results also appended to setup.mat.

See also NRGWilsonQS, fdmNRG_QS, tdmNRG

AW (C) 2008 ; RIXS 2019