QInchworm.qmc_integrate

Quasi Monte Carlo integration routines.

qmc_integral(f; ...)
qmc_integral(f, init; trans_f, jacobian_f, seq, N)

Compute a quasi Monte Carlo estimate of a $d$-dimensional integral

\[F = \int_\mathscr{D} d^d \mathbf{u}\ f(\mathbf{u}).\]

The domain $\mathscr{D}$ is defined by a variable change $\mathbf{u} = \mathbf{u}(\mathbf{x}): [0,1]^d \mapsto \mathscr{D}$,

\[F = \int_{[0,1]^d} d^d\mathbf{x}\ f(\mathbf{u}(\mathbf{x})) \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right| \approx \frac{1}{N} \sum_{i=1}^N f(\mathbf{u}(\mathbf{x}_i)) J(\mathbf{u}(\mathbf{x}_i)).\]

Parameters

• f: Integrand $f(\mathbf{u})$.
• init: Initial (zero) value used in qMC summation.
• trans_f: Domain transformation function $\mathbf{u}(\mathbf{x})$.
• jacobian_f: Jacobian $J(\mathbf{u}) = \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right|$.
• seq: Quasi-random sequence generator.
• N: Number of points to be taken from the quasi-random sequence.

Returns

Estimated value of the integral.

qmc_integral_n_samples(f; ...)
qmc_integral_n_samples(
    f,
    init;
    trans_f,
    jacobian_f,
    seq,
    N_samples
)

Compute a quasi Monte Carlo estimate of a $d$-dimensional integral

\[F = \int_\mathscr{D} d^d \mathbf{u}\ f(\mathbf{u}).\]

The domain $\mathscr{D}$ is defined by a variable change $\mathbf{u} = \mathbf{u}(\mathbf{x}): [0,1]^d \mapsto \mathscr{D}$,

\[F = \int_{[0,1]^d} d^d\mathbf{x}\ f(\mathbf{u}(\mathbf{x})) \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right| \approx \frac{1}{N} \sum_{i=1}^N f(\mathbf{u}(\mathbf{x}_i)) J(\mathbf{u}(\mathbf{x}_i)).\]

Unlike qmc_integral(), this function performs qMC summation until a given number of valid (non-nothing) samples of $f(\mathbf{u})$ are taken.

Parameters

• f: Integrand $f(\mathbf{u})$.
• init: Initial (zero) value used in qMC summation.
• trans_f: Domain transformation function $\mathbf{u}(\mathbf{x})$.
• jacobian_f: Jacobian $J(\mathbf{u}) = \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right|$.
• seq: Quasi-random sequence generator.
• N_samples: Number of valid samples of $f(\mathbf{u})$ to be taken.

Returns

Estimated value of the integral.

abstract type AbstractDomainTransform{D}

Abstract domain transformation $[0, 1]^d \mapsto \mathscr{D}$.

ndims(
    s::QInchworm.scrambled_sobol.ScrambledSobolSeq{D}
) -> Int64

Dimension D of a scrambled Sobol sequence.

ndims(s::QInchworm.utility.RandomSeq{RNG, D}) -> Int64

Dimension D of a random sequence.

ndims(
    _::QInchworm.qmc_integrate.AbstractDomainTransform{D}
) -> Any

Return the number of dimensions $d$ of a domain transformation $[0, 1]^d \mapsto \mathscr{D}$.

struct ExpModelFunctionTransform{D} <: QInchworm.qmc_integrate.AbstractDomainTransform{D}

Domain transformation

\[\mathbf{x}\in[0,1]^d \mapsto \{\mathbf{u}: u_f \geq u_1 \geq u_2 \geq \ldots \geq u_d > -\infty\}\]

induced by the implicit variable change

\[x_n(v_n) = \frac{\int_0^{v_n} d\bar v_n h(\bar v_n)} {\int_0^\infty d\bar v_n h(\bar v_n)}, \quad v_n = \left\{ \begin{array}{ll} u_f - u_1, &n=1,\\ u_{n-1} - u_n, &n>1, \end{array} \right.\]

where $h(v)$ is an exponential model function parametrized by decay rate $\tau$, $h(v) = e^{-v/\tau}$.

The corresponding Jacobian is $J(\mathbf{u}) = \tau^d / e^{-(u_f - u_d) / \tau}$.

Fields

• u_f::Float64: Upper bound of the transformed domain $u_f$

• τ::Float64: Decay rate parameter of the exponential model function

ExpModelFunctionTransform(
    d::Integer,
    c::Keldysh.AbstractContour,
    t_f::Keldysh.BranchPoint,
    τ::Real
) -> QInchworm.qmc_integrate.ExpModelFunctionTransform

Make an ExpModelFunctionTransform object suitable for time contour integration over the domain

\[\{(t_1, \ldots, t_d) \in \mathcal{C}^d: t_f \succeq t_1 \succeq t_2 \succeq \ldots \succeq t_d \succeq \text{starting point of }\mathcal{C}\}\]

N.B. ExpModelFunctionTransform describes an infinite domain where integration variables $u_n$ can approach $-\infty$. Negative values of $u_n$ cannot be mapped onto time points on $\mathcal{C}$ and will be discarded by contour_integral().

Parameters

• d: Number of dimensions $d$.
• c: Time contour $\mathcal{C}$.
• t_f: Upper bound $t_f$.
• τ: Decay rate parameter $\tau$.

struct RootTransform{D} <: QInchworm.qmc_integrate.AbstractDomainTransform{D}

Hypercube-to-simplex domain transformation

\[\mathbf{x}\in[0,1]^d \mapsto \{\mathbf{u}: u_f \geq u_1 \geq u_2 \geq \ldots \geq u_d \geq u_i\}\]

induced by the Root mapping^[2]

\[\left\{ \begin{array}{ll} \tilde u_1 &= x_1^{1/d},\\ \tilde u_2 &= \tilde u_1 x_2^{1/(d-1)},\\ &\ldots\\ \tilde u_d &= \tilde u_{d-1} x_d, \end{array}\right.\]

\[\mathbf{u} = u_i + (u_f - u_i) \tilde{\mathbf{u}}.\]

The corresponding Jacobian is $J(\mathbf{u}) = (u_f - u_i)^d / d!$.

Fields

• u_i::Float64: Lower bound of the transformed domain $u_i$

• u_diff::Float64: Difference $u_f - u_i$

RootTransform(
    d::Integer,
    c::Keldysh.AbstractContour,
    t_i::Keldysh.BranchPoint,
    t_f::Keldysh.BranchPoint
) -> QInchworm.qmc_integrate.RootTransform

Make a RootTransform object suitable for time contour integration over the domain

\[\{(t_1, \ldots, t_d) \in \mathcal{C}^d: t_f \succeq t_1 \succeq t_2 \succeq \ldots \succeq t_d \succeq t_i\}\]

Parameters

• d: Number of dimensions $d$.
• c: Time contour $\mathcal{C}$.
• t_i: Lower bound $t_i$.
• t_f: Upper bound $t_f$.

struct SortTransform{D} <: QInchworm.qmc_integrate.AbstractDomainTransform{D}

Hypercube-to-simplex domain transformation

\[\mathbf{x}\in[0,1]^d \mapsto \{\mathbf{u}: u_f \geq u_1 \geq u_2 \geq \ldots \geq u_d \geq u_i\}\]

induced by the Sort mapping^[2]

\[\mathbf{u} = u_i + (u_f - u_i) \mathrm{sort}(x_1, \ldots, x_d).\]

The corresponding Jacobian is $J(\mathbf{u}) = (u_f - u_i)^d / d!$.

Fields

• u_i::Float64: Lower bound of the transformed domain $u_i$

• u_diff::Float64: Difference $u_f - u_i$

SortTransform(
    d::Integer,
    c::Keldysh.AbstractContour,
    t_i::Keldysh.BranchPoint,
    t_f::Keldysh.BranchPoint
) -> QInchworm.qmc_integrate.SortTransform

Make a SortTransform object suitable for time contour integration over the domain

\[\{(t_1, \ldots, t_d) \in \mathcal{C}^d: t_f \succeq t_1 \succeq t_2 \succeq \ldots \succeq t_d \succeq t_i\}\]

Parameters

• d: Number of dimensions $d$.
• c: Time contour $\mathcal{C}$.
• t_i: Lower bound $t_i$.
• t_f: Upper bound $t_f$.

struct DoubleSimplexRootTransform{D} <: QInchworm.qmc_integrate.AbstractDomainTransform{D}

Hypercube-to-double-simplex domain transformation

\[\mathbf{x}\in[0,1]^d \mapsto \{\mathbf{u}: u_f \geq u_1 \geq u_2 \geq \ldots \geq u_{d^>} \geq u_w \geq u_{d^>+1} \geq \ldots \geq u_d \geq u_i\}\]

induced by the Root mapping^[2] applied independently to two sets of variables $\{x_1, \ldots, x_{d^>} \}$ and $\{x_{d^>+1}, \ldots, x_d \}$ (cf. RootTransform).

The corresponding Jacobian is $J(\mathbf{u}) = \frac{(u_w - u_i)^{d^<}}{d^<!} \frac{(u_f - u_w)^{d^>}}{d^>!}$.

Here, $d^<$ and $d^>$ denote numbers of variables in the 'lesser' and 'greater' simplex respectively, and $d = d^< + d^>$. The two simplices are separated by a fixed boundary located at $u_w$.

Fields

• d_lesser::Int64: Number of variables in the 'lesser' simplex, $d^<$

• d_greater::Int64: Number of variables in the 'greater' simplex, $d^>$

• u_i::Float64: Lower bound of the transformed domain $u_i$

• u_w::Float64: Boundary $u_w$ separating the 'lesser' and the 'greater' simplices

• u_diff_wi::Float64: Difference $u_w - u_i$

• u_diff_fw::Float64: Difference $u_f - u_w$

DoubleSimplexRootTransform(
    d_before::Int64,
    d_after::Int64,
    c::Keldysh.AbstractContour,
    t_i::Keldysh.BranchPoint,
    t_w::Keldysh.BranchPoint,
    t_f::Keldysh.BranchPoint
) -> QInchworm.qmc_integrate.DoubleSimplexRootTransform

Make a DoubleSimplexRootTransform object suitable for time contour integration over the domain

\[\{(t_1, \ldots, t_d) \in \mathcal{C}^d: t_f \succeq t_1 \succeq t_2 \succeq \ldots \succeq t_{d_\text{after}} \succeq t_w \succeq t_{d_\text{after}+1} \succeq \ldots \succeq t_d \succeq t_i\},\]

where $d = d_\text{before} + d_\text{after}$.

Parameters

• d_before: Number of time variables in the 'before' region, $d_\text{before}$.
• d_after: Number of time variables in the 'after' region, $d_\text{after}$.
• c: Time contour $\mathcal{C}$.
• t_i: Lower bound $t_i$.
• t_w: Boundary point $t_w$.
• t_f: Upper bound $t_f$.

make_trans_f(
    t::QInchworm.qmc_integrate.ExpModelFunctionTransform
) -> QInchworm.qmc_integrate.var"#3#4"{QInchworm.qmc_integrate.ExpModelFunctionTransform{D}} where D

Return the function $\mathbf{u}(\mathbf{x})$ corresponding to a given ExpModelFunctionTransform object t.

make_trans_f(
    t::QInchworm.qmc_integrate.RootTransform
) -> QInchworm.qmc_integrate.var"#7#8"{QInchworm.qmc_integrate.RootTransform{D}} where D

Return the function $\mathbf{u}(\mathbf{x})$ corresponding to a given RootTransform object t.

make_trans_f(
    t::QInchworm.qmc_integrate.SortTransform
) -> QInchworm.qmc_integrate.var"#11#12"{QInchworm.qmc_integrate.SortTransform{D}} where D

Return the function $\mathbf{u}(\mathbf{x})$ corresponding to a given SortTransform object t.

make_trans_f(
    t::QInchworm.qmc_integrate.DoubleSimplexRootTransform
) -> QInchworm.qmc_integrate.var"#15#16"{QInchworm.qmc_integrate.DoubleSimplexRootTransform{D}} where D

Return the function $\mathbf{u}(\mathbf{x})$ corresponding to a given DoubleSimplexRootTransform object t.

make_jacobian_f(
    t::QInchworm.qmc_integrate.ExpModelFunctionTransform
) -> QInchworm.qmc_integrate.var"#5#6"{QInchworm.qmc_integrate.ExpModelFunctionTransform{D}} where D

Return the Jacobian $J(\mathbf{u}) = \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right|$ corresponding to a given ExpModelFunctionTransform object t.

make_jacobian_f(
    t::QInchworm.qmc_integrate.RootTransform
) -> QInchworm.qmc_integrate.var"#9#10"{Float64}

Return the Jacobian $J(\mathbf{u}) = \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right|$ corresponding to a given RootTransform object t.

make_jacobian_f(
    t::QInchworm.qmc_integrate.SortTransform
) -> QInchworm.qmc_integrate.var"#13#14"{Float64}

Return the Jacobian $J(\mathbf{u}) = \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right|$ corresponding to a given SortTransform object t.

make_jacobian_f(
    t::QInchworm.qmc_integrate.DoubleSimplexRootTransform
) -> QInchworm.qmc_integrate.var"#17#18"{Float64}

Return the Jacobian $J(\mathbf{u}) = \left|\frac{\partial\mathbf{u}}{\partial\mathbf{x}}\right|$ corresponding to a given DoubleSimplexRootTransform object t.

contour_integral(
    f,
    c::Keldysh.AbstractContour,
    dt::QInchworm.qmc_integrate.AbstractDomainTransform;
    init,
    seq,
    N
)

Compute a quasi Monte Carlo estimate of a contour integral over a $d$-dimensional domain $\mathscr{D}$.

Parameters

• f: Integrand.
• c: Time contour to integrate over.
• dt: Domain transformation $[0, 1]^d \mapsto \mathscr{D}$.
• init: Initial (zero) value of the integral.
• seq: Quasi-random sequence generator.
• N: The number of points to be taken from the quasi-random sequence.

Returns

Estimated value of the contour integral.

contour_integral_n_samples(
    f,
    c::Keldysh.AbstractContour,
    dt::QInchworm.qmc_integrate.AbstractDomainTransform;
    init,
    seq,
    N_samples
)

Compute a quasi Monte Carlo estimate of a contour integral over a $d$-dimensional domain $\mathscr{D}$.

Unlike contour_integral(), this function performs qMC summation until a given number of valid (non-nothing) samples of the integrand are taken.

Parameters

• f: Integrand.
• c: Time contour to integrate over.
• dt: Domain transformation $[0, 1]^d \mapsto \mathscr{D}$.
• init: Initial (zero) value of the integral.
• seq: Quasi-random sequence generator.
• N_samplex: Number of valid samples of the integrand to be taken.

Returns

• Estimated value of the contour integral.
• The total number of points taken from the quasi-random sequence.

Dictionary mapping branches of the Keldysh contour to their unitary direction coefficients in the complex time plane.

contour_function_return_type(f::Function) -> Any

Detect the return type of a function applied to a vector of Keldysh.BranchPoint.