|
TRIQS/TRIQS 4.0.0
Researching Interacting Quantum Systems
|
#include <triqs/mesh/tail_fitter.hpp>
Fit the high- and low-frequency tail of a function \( f \) defined on a triqs::mesh::refreq or a triqs::mesh::imfreq frequency mesh.
For the moment, we assume that the function \( f \) to be fitted only depends on a single frequency argument \( z \) and is matrix-valued, i.e. \( f(z) \in \mathbb{C}^{d_1 \times d_2} \). Then we can write the high- and low-frequency expansion as
\[ f(z) = \sum_{n=0}^{q} \frac{A_n}{z^n} + \mathcal{O}(z^{-q-1}) \; , \]
where \( A_n \in \mathbb{C}^{d_1 \times d_2} \) is a matrix containing the n-th order expansion coefficients.
Since the function is defined on discrete mesh points, its values can be stored in a data array \( D \) of rank \( 3 \) with shape \( (N, d_1, d_2) \), where \( N \) is the number of points (or discrete frequencies).
Let \( \Omega = \{ z_l, \dots, z_h \} \) be the total set of \( N \) mesh points, with the lowest frequency \(z_l \) and the highest frequency \( z_h \). The two subranges of mesh points, \( \widetilde{\Omega}_l = \{ z_l, \dots, z_{l+p_r-1} \} \) and \( \widetilde{\Omega}_h = \{ z_{h-p_r+1}, \dots, z_h \} \), relevant for the low- and high-frequency tail fit are determined by \( p_r = \lfloor r N / 2 + 0.5 \rfloor \) with \( 0 < r \leq 1 \). The actual number of frequencies used for the fitting procedure is \( p = \min\{ p_r, p_\text{max} \} \), where \( p_\text{max} > 0\) is a given maximum number of points. \( r \) and \( p_\text{max} \) are set in the constructor of the tail fitter.
Once we have the \( p \) low-frequencies \( \Omega_l = \{ z_{l_0}, \dots, z_{l_{p-1}} \} \subseteq \widetilde{\Omega}_l \) and the \( p \) high-frequencies \( \Omega_h = \{ z_{h_0}, \dots, z_{h_{p-1}} \} \subseteq \widetilde{\Omega}_h \), we can write the equations for the unknown coefficients \( A_n \) as a system of linear equations
\[ \begin{bmatrix} 1 & z_{l_0} & z_{l_0}^{-2} & \cdots & z_{l_0}^{-q} \\ 1 & z_{h_0} & z_{h_0}^{-2} & \cdots & z_{h_0}^{-q} \\ 1 & z_{l_1} & z_{l_1}^{-2} & \cdots & z_{l_1}^{-q} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & z_{h_{p-1}} & z_{h_{p-1}}^{-2} & \cdots & z_{h_{p-1}}^{-q} \\ \end{bmatrix} \begin{bmatrix} [A_0]_{11} & [A_0]_{12} & [A_0]_{21} & [A_0]_{22} \\ [A_1]_{11} & [A_1]_{12} & [A_1]_{21} & [A_1]_{22} \\ [A_2]_{11} & [A_2]_{12} & [A_2]_{21} & [A_2]_{22} \\ \vdots & \vdots & \vdots & \vdots \\ [A_{q-1}]_{11} & [A_{q-1}]_{12} & [A_{q-1}]_{21} & [A_{q-1}]_{22} \\ \end{bmatrix} = \begin{bmatrix} [f(z_{l_0})]_{11} & [f(z_{l_0})]_{12} & [f(z_{l_0})]_{21} & [f(z_{l_0})]_{22} \\ [f(z_{h_0})]_{11} & [f(z_{h_0})]_{12} & [f(z_{h_0})]_{21} & [f(z_{h_0})]_{22} \\ [f(z_{l_1})]_{11} & [f(z_{l_1})]_{12} & [f(z_{l_1})]_{21} & [f(z_{l_1})]_{22} \\ \vdots & \vdots & \vdots & \vdots \\ [f(z_{h_{p-1}})]_{11} & [f(z_{h_{p-1}})]_{12} & [f(z_{h_{p-1}})]_{21} & [f(z_{h_{p-1}})]_{22} \\ \end{bmatrix} \; , \]
or in a more compact form
\[ V A = B \; . \]
Here, \( V \) is the \( 2p \times (q+1) \) Vandermonde matrix (see also triqs::mesh::vander), \( A \) is the \( (q+1) \times (d_1 \times d_2) \) matrix of unknown coefficients and \( B \) is the \( 2p \times (d_1 \times d_2) \) matrix of function values at the frequencies in \( \Omega_l \) and \( \Omega_h \).
This system of equations can then be solved with the linear least squares worker classes nda::lapack::gelss_worker or nda::lapack::gelss_worker_hermitian. Hermitian means that the coefficient matrices \( A_n \) are required to be hermitian, i.e. \( [A_n]_{ij} = [A_n]_{ji}^* \), which enforces the symmetry \([f(z)]_{ij} = [f(-z)]_{ji}^* \) in the function values.
Definition at line 198 of file tail_fitter.hpp.
Public Member Functions | |
| tail_fitter (double r, int p_max, std::optional< int > q={}) | |
| Construct a tail fitter for a given fraction \( r \) of the mesh, the maximum number of mesh points \( p_{\text{max}} \) to use in the fit and an optional expansion order \( q \). | |
| template<int P, bool enforce_hermiticity = false, typename M, int R> | |
| auto | fit (M const &m, nda::array_const_view< std::complex< double >, R > D, bool rescale, nda::array_const_view< std::complex< double >, R > C, std::optional< long > d={}) |
| Perform a linear least squares fit and return the coefficients \( A_n \) of the tail expansion together with the error of the fit. | |
| template<int P, typename M, int R> | |
| auto | fit_hermitian (M const &m, nda::array_const_view< std::complex< double >, R > D, bool rescale, nda::array_const_view< std::complex< double >, R > C, std::optional< long > d={}) |
| Perform a linear least squares fit and return the coefficients \( A_n \) of the tail expansion together with the error of the fit. | |
| template<bool enforce_hermiticity = false> | |
| auto & | get_lss () |
| Get linear least squares workers. | |
| template<typename M> | |
| auto | get_tail_fit_indices (M const &m) |
| Get a vector containing the indices of all mesh points to use in the tail fit, i.e. \( (l_0, h_0, l_1,
h_1, \dots, l_{p-1}, h_{p-1}) \). | |
| double | get_tail_fraction () const |
| Get the fraction \( r \) of the mesh to be considered for the tail fit. | |
| template<typename M> | |
| int | n_pts_in_tail (M const &m) const |
| Get the number of mesh points used in the tail fit for the given mesh. | |
| template<bool enforce_hermiticity = false, typename M> | |
| void | setup_lss (M const &m, int n_A) |
| Set up the linear least squares workers for a given mesh and a given number \( n_A \) of known coefficient arrays \( A_n \). | |
|
inline |
Construct a tail fitter for a given fraction \( r \) of the mesh, the maximum number of mesh points \( p_{\text{max}} \) to use in the fit and an optional expansion order \( q \).
| r | Fraction of the mesh to consider in the tail fit ( \( 0 < r \leq 1 \)). |
| p_max | Maximum number of points to use in the tail fit ( \( p_\text{max} > 0 \)). |
| q | Optional expansion order \( q \leq q_{\text{max}} = 9 \). If not set, it will be adjusted automatically. |
Definition at line 215 of file tail_fitter.hpp.
|
inline |
Perform a linear least squares fit and return the coefficients \( A_n \) of the tail expansion together with the error of the fit.
If the function \( f \) depends on more arguments, we first permute the indices of the given data array \( D \) such that the frequency mesh we are fitting on corresponds to the first dimension. Then we combine all other dimensions to form a matrix of size \( N \times M \), where \( N \) is the size of the frequency mesh and \( M \) is the product of the remaining dimensions. The same is done with the array \( C \) containing the known coefficient arrays \( A_n \) for \( n < n_A \).
Before the fit is performed, the expansion terms corresponding to the known coefficients (contained in \( C \)) are subtracted from the function values (contained in \( D \)), i.e.
\[ \tilde{f}(z_i) = f(z_i) - \sum_{n=0}^{n_A-1} \tilde{A}_n \left( \frac{|z_{\text{max}}|}{z_i} \right)^n \; , \]
where \( n_A \) is the number of known \( A_n \) and \( \tilde{A}_n = A_n / |z_{\text{max}}|^n \) (see setup_lss() for more information).
After the least squares procedure, the original coefficients \( A_n \) can be recovered from \( \tilde{A}_n \) by setting rescale to true. Otherwise, the scaled coefficients \( \tilde{A}_n \) are returned.
Furthermore, if enforce_hermiticity is set to true, the calculated coefficient matrices \( A_n \) are forced to be hermitian. This setting only makes sense for triqs::mesh::imfreq Matsubara meshes and requires an inner matrix dimension \( d \) to be specified, i.e. the number of rows/columns of the square matrix that \( f(z) \) returns.
| P | Position of the frequency mesh in case of a product mesh. |
| enforce_hermiticity | Enforce hermiticity in the coefficient matrices \( A_n \). |
| M | Frequency mesh type (either triqs::mesh::imfreq or triqs::mesh::refreq). |
| R | Rank of the data array and the array containing the known moments. |
| m | Frequency mesh. |
| D | Data array \( D \) containing the function values on the (product) mesh points. |
| rescale | Should we rescale the calculated coefficients \( \tilde{A}_n \) by \( |z_{\text{max}}|^q \)? |
| C | Data array \( C \) containing the known coefficient arrays. |
| d | Inner matrix dimensions \( d \) (only needed if enforce_hermiticity is true). |
Definition at line 398 of file tail_fitter.hpp.
|
inline |
Perform a linear least squares fit and return the coefficients \( A_n \) of the tail expansion together with the error of the fit.
It simply calls fit() with enforce_hermiticity set to true.
| P | Position of the frequency mesh in case of a product mesh. |
| M | Frequency mesh type (either triqs::mesh::imfreq or triqs::mesh::refreq). |
| R | Rank of the data array and the array containing the known moments. |
| m | Frequency mesh. |
| D | Data array \( D \) containing the function values on the (product) mesh points. |
| rescale | Should we rescale the calculated coefficients \( \tilde{A}_n \) by \( |z_{\text{max}}|^q \)? |
| C | Data array \( C \) containing the known coefficient arrays. |
| d | Inner matrix dimensions \( d \). |
Definition at line 498 of file tail_fitter.hpp.
|
inline |
Get linear least squares workers.
| enforce_hermiticity | Enforce hermiticity in the coefficient matrices \( A_n \). |
Definition at line 279 of file tail_fitter.hpp.
|
inline |
Get a vector containing the indices of all mesh points to use in the tail fit, i.e. \( (l_0, h_0, l_1, h_1, \dots, l_{p-1}, h_{p-1}) \).
| M | Frequency mesh type (either triqs::mesh::imfreq or triqs::mesh::refreq). |
| m | Frequency mesh. |
Definition at line 246 of file tail_fitter.hpp.
|
inline |
Get the number of mesh points used in the tail fit for the given mesh.
The number of points \( p \) depends on the size \( N \) of the underlying mesh, the fraction \( r \) of the mesh that should be considered and the maximum number of points \( p_\text{max} \) we want to use:
\[ p = \min \left\{ \left\lfloor \frac{N}{2} r + 0.5 \right\rfloor, p_\text{max} \right\} \; . \]
| M | Frequency mesh type (either triqs::mesh::imfreq or triqs::mesh::refreq). |
| m | Frequency mesh. |
Definition at line 233 of file tail_fitter.hpp.
|
inline |
Set up the linear least squares workers for a given mesh and a given number \( n_A \) of known coefficient arrays \( A_n \).
To set up the workers, we first determine the frequencies for the tail fit with get_tail_fit_indices() and then build the Vandermonde matrix with triqs::mesh::vander. For numerical reasons, we scale the frequencies by the absolute value of the maximum frequency in the mesh such that \( V_{ij} = \left(|z_{\text{max}}| / z_i \right)^{j} \). Note that this will also scale the coefficients \( \tilde{A}_n = A_n / |z_{\text{max}}|^{n} \).
The least squares workers are then initialized with the Vandermonde matrix. If the given number \( n_A \) of known coefficient arrays is \( > 0 \), then only the \( 2p \times (q - n_A + 1) \) submatrix of \( V \) is used. That means that only \( q - n_A + 1 \) coefficients \( A_n \) with \( n \geq n_A \) will be calculated with the fitting procedure.
Furthermore, if no expansion order \( q \) was specified during construction, this function tries to find the largest \( q \) such that \( n_A \leq q \leq q_{\text{max}} = 9 \) and for which the smallest singular value of the Vandermonde matrix is larger than some threshold \( r_{\text{cond}} = 10^{-8} \). If no suitable \( q \) is found, an exception is thrown.
| enforce_hermiticity | Enforce hermiticity in the coefficient matrices \( A_n \). |
| M | Frequency mesh type (either triqs::mesh::imfreq or triqs::mesh::refreq). |
| m | Frequency mesh. |
| n_A | Number of known coefficient arrays. |
Definition at line 311 of file tail_fitter.hpp.