text:
/-
Copyright (c) 2019 Johannes Hölzl. All rights reserved.
Released under Apache 2.0 license as described in the file LICENSE.
Authors: Johannes Hölzl, Patrick Massot, Casper Putz, Anne Baanen
-/
import linear_algebra.matrix.determinant
import tactic.fin_cases

/-!
# Block matrices and their determinant

This file defines a predicate `matrix.block_triangular_matrix` saying a matrix
is block triangular, and proves the value of the determinant for various
matrices built out of blocks.

## Main definitions

 * `matrix.block_triangular_matrix` expresses that a `o` by `o` matrix is block triangular,
   if the rows and columns are ordered according to some order `b : o → ℕ`

## Main results
  * `det_of_block_triangular_matrix`: the determinant of a block triangular matrix
    is equal to the product of the determinants of all the blocks
  * `det_of_upper_triangular` and `det_of_lower_triangular`: the determinant of
    a triangular matrix is the product of the entries along the diagonal

## Tags

matrix, diagonal, det, block triangular

-/

open_locale big_operators

universes v

variables {m n : Type*} [decidable_eq n] [fintype n] [decidable_eq m] [fintype m]
variables {R : Type v} [comm_ring R]

namespace matrix

lemma det_to_block (M : matrix m m R) (p : m → Prop) [decidable_pred p] :
  M.det = (matrix.from_blocks (to_block M p p) (to_block M p (λ j, ¬p j))
    (to_block M (λ j, ¬p j) p) (to_block M (λ j, ¬p j) (λ j, ¬p j))).det :=
begin
  rw ← matrix.det_reindex_self (equiv.sum_compl p).symm M,
  rw [det_apply', det_apply'],
  congr, ext σ, congr, ext,
  generalize hy : σ x = y,
  cases x; cases y;
  simp only [matrix.reindex_apply, to_block_apply, equiv.symm_symm,
    equiv.sum_compl_apply_inr, equiv.sum_compl_apply_inl,
    from_blocks_apply₁₁, from_blocks_apply₁₂, from_blocks_apply₂₁, from_blocks_apply₂₂,
    matrix.minor_apply],
end

lemma det_to_square_block (M : matrix m m R) {n : nat} (b : m → fin n) (k : fin n) :
  (to_square_block M b k).det = (to_square_block_prop M (λ i, b i = k)).det :=
by simp

lemma det_to_square_block' (M : matrix m m R) (b : m → ℕ) (k : ℕ) :
  (to_square_block' M b k).det = (to_square_block_prop M (λ i, b i = k)).det :=
by simp

lemma two_block_triangular_det (M : matrix m m R) (p : m → Prop) [decidable_pred p]
  (h : ∀ i (h1 : ¬p i) j (h2 : p j), M i j = 0) :
  M.det = (to_square_block_prop M p).det * (to_square_block_prop M (λ i, ¬p i)).det :=
begin
  rw det_to_block M p,
  convert det_from_blocks_zero₂₁ (to_block M p p) (to_block M p (λ j, ¬p j))
    (to_block M (λ j, ¬p j) (λ j, ¬p j)),
  ext,
  exact h ↑i i.2 ↑j j.2
end

lemma equiv_block_det (M : matrix m m R) {p q : m → Prop} [decidable_pred p] [decidable_pred q]
  (e : ∀x, q x ↔ p x) : (to_square_block_prop M p).det = (to_square_block_prop M q).det :=
by convert matrix.det_reindex_self (equiv.subtype_equiv_right e) (to_square_block_prop M q)

lemma to_square_block_det'' (M : matrix m m R) {n : nat} (b : m → fin n) (k : fin n) :
  (to_square_block M b k).det = (to_square_block' M (λ i, ↑(b i)) ↑k).det :=
begin
  rw [to_square_block_def', to_square_block_def],
  apply equiv_block_det,
  intro x,
  apply (fin.ext_iff _ _).symm
end

/-- Let `b` map rows and columns of a square matrix `M` to `n` blocks. Then
  `block_triangular_matrix' M n b` says the matrix is block triangular. -/
def block_triangular_matrix' {o : Type*} (M : matrix o o R) {n : ℕ}
  (b : o → fin n) : Prop :=
∀ i j, b j < b i → M i j = 0

lemma upper_two_block_triangular' {m n : Type*}
  (A : matrix m m R) (B : matrix m n R) (D : matrix n n R) :
  block_triangular_matrix' (from_blocks A B 0 D) (sum.elim (λ i, (0 : fin 2)) (λ j, 1)) :=
begin
  intros k1 k2 hk12,
  have h0 : ∀ (k : m ⊕ n), sum.elim (λ i, (0 : fin 2)) (λ j, 1) k = 0 → ∃ i, k = sum.inl i,
  { simp },
  have h1 : ∀ (k : m ⊕ n), sum.elim (λ i, (0 : fin 2)) (λ j, 1) k = 1 → ∃ j, k = sum.inr j,
  { simp },
  set mk1 := (sum.elim (λ i, (0 : fin 2)) (λ j, 1)) k1 with hmk1,
  set mk2 := (sum.elim (λ i, (0 : fin 2)) (λ j, 1)) k2 with hmk2,
  fin_cases mk1 using h; fin_cases mk2 using h_1; rw [h, h_1] at hk12,
  { exact absurd hk12 (nat.not_lt_zero 0) },
  { exact absurd hk12 (by norm_num) },
  { rw hmk1 at h,
    obtain ⟨i, hi⟩ := h1 k1 h,
    rw hmk2 at h_1,
    obtain ⟨j, hj⟩ := h0 k2 h_1,
    rw [hi, hj], simp },
  { exact absurd hk12 (irrefl 1) }
end

/-- Let `b` map rows and columns of a square matrix `M` to blocks indexed by `ℕ`s. Then
  `block_triangular_matrix M n b` says the matrix is block triangular. -/
def block_triangular_matrix {o : Type*} (M : matrix o o R) (b : o → ℕ) : Prop :=
∀ i j, b j < b i → M i j = 0

lemma upper_two_block_triangular {m n : Type*}
  (A : matrix m m R) (B : matrix m n R) (D : matrix n n R) :
  block_triangular_matrix (from_blocks A B 0 D) (sum.elim (λ i, 0) (λ j, 1)) :=
begin
  intros k1 k2 hk12,
  have h01 : ∀ (k : m ⊕ n), sum.elim (λ i, 0) (λ j, 1) k = 0 ∨ sum.elim (λ i, 0) (λ j, 1) k = 1,
  { simp },
  have h0 : ∀ (k : m ⊕ n), sum.elim (λ i, 0) (λ j, 1) k = 0 → ∃ i, k = sum.inl i, { simp },
  have h1 : ∀ (k : m ⊕ n), sum.elim (λ i, 0) (λ j, 1) k = 1 → ∃ j, k = sum.inr j, { simp },
  cases (h01 k1) with hk1 hk1; cases (h01 k2) with hk2 hk2; rw [hk1, hk2] at hk12,
  { exact absurd hk12 (nat.not_lt_zero 0) },
  { exact absurd hk12 (nat.not_lt_zero 1) },
  { obtain ⟨i, hi⟩ := h1 k1 hk1,
    obtain ⟨j, hj⟩ := h0 k2 hk2,
    rw [hi, hj], simp },
  { exact absurd hk12 (irrefl 1) }
end

lemma det_of_block_triangular_matrix (M : matrix m m R) (b : m → ℕ)
  (h : block_triangular_matrix M b) :
  ∀ (n : ℕ) (hn : ∀ i, b i < n), M.det = ∏ k in finset.range n, (to_square_block' M b k).det :=
begin
  intros n hn,
  unfreezingI { induction n with n hi generalizing m M b },
  { rw finset.prod_range_zero,
    apply det_eq_one_of_card_eq_zero,
    apply fintype.card_eq_zero_iff.mpr,
    exact ⟨λ i, nat.not_lt_zero (b i) (hn i)⟩ },
  { rw [finset.prod_range_succ_comm],
    have h2 : (M.to_square_block_prop (λ (i : m), b i = n.succ)).det =
      (M.to_square_block' b n.succ).det,
    { dunfold to_square_block', dunfold to_square_block_prop, refl },
    rw two_block_triangular_det M (λ i, ¬(b i = n)),
    { rw mul_comm,
      apply congr (congr_arg has_mul.mul _),
      { let m' := {a // ¬b a = n },
        let b' := (λ (i : m'), b ↑i),
        have h' :
          block_triangular_matrix (M.to_square_block_prop (λ (i : m), ¬b i = n)) b',
        { intros i j, apply h ↑i ↑j },
        have hni : ∀ (i : {a // ¬b a = n}), b' i < n,
        { exact λ i, (ne.le_iff_lt i.property).mp (nat.lt_succ_iff.mp (hn ↑i)) },
        have h1 := hi (M.to_square_block_prop (λ (i : m), ¬b i = n)) b' h' hni,
        rw ←fin.prod_univ_eq_prod_range at h1 ⊢,
        convert h1,
        ext k,
        simp only [to_square_block_def', to_square_block_def],
        let he : {a // b' a = ↑k} ≃ {a // b a = ↑k},
        { have hc : ∀ (i : m), (λ a, b a = ↑k) i → (λ a, ¬b a = n) i,
          { intros i hbi, rw hbi, exact ne_of_lt (fin.is_lt k) },
          exact equiv.subtype_subtype_equiv_subtype hc },
        exact matrix.det_reindex_self he (λ (i j : {a // b' a = ↑k}), M ↑i ↑j) },
      { rw det_to_square_block' M b n,
        have hh : ∀ a, b a = n ↔ ¬(λ (i : m), ¬b i = n) a,
        { intro i, simp only [not_not] },
        exact equiv_block_det M hh }},
    { intros i hi j hj,
      apply (h i), simp only [not_not] at hi,
      rw hi,
      exact (ne.le_iff_lt hj).mp (nat.lt_succ_iff.mp (hn j)) }}
end

lemma det_of_block_triangular_matrix'' (M : matrix m m R) (b : m → ℕ)
  (h : block_triangular_matrix M b) :
  M.det = ∏ k in finset.image b finset.univ, (to_square_block' M b k).det :=
begin
  let n : ℕ := (Sup (finset.image b finset.univ : set ℕ)).succ,
  have hn : ∀ i, b i < n,
  { have hbi : ∀ i, b i ∈ finset.image b finset.univ, { simp },
    intro i,
    dsimp only [n],
    apply nat.lt_succ_iff.mpr,
    exact le_cSup (finset.bdd_above _) (hbi i) },
  rw det_of_block_triangular_matrix M b h n hn,
  refine (finset.prod_subset _ _).symm,
  { intros a ha, apply finset.mem_range.mpr,
    obtain ⟨i, ⟨hi, hbi⟩⟩ := finset.mem_image.mp ha,
    rw ←hbi,
    exact hn i },
  { intros k hk hbk,
    apply det_eq_one_of_card_eq_zero,
    apply fintype.card_eq_zero_iff.mpr,
    constructor,
    simp only [subtype.forall],
    intros a hba, apply hbk,
    apply finset.mem_image.mpr,
    use a,
    exact ⟨finset.mem_univ a, hba⟩ }
end

lemma det_of_block_triangular_matrix' (M : matrix m m R) {n : ℕ} (b : m → fin n)
  (h : block_triangular_matrix' M b) :
  M.det = ∏ (k : fin n), (to_square_block M b k).det :=
begin
  let b2 : m → ℕ := λ i, ↑(b i),
  simp_rw to_square_block_det'',
  rw fin.prod_univ_eq_prod_range (λ (k : ℕ), (M.to_square_block' b2 k).det) n,
  apply det_of_block_triangular_matrix,
  { intros i j hij, exact h i j (fin.coe_fin_lt.mp hij) },
  { intro i, exact fin.is_lt (b i) }
end

lemma det_of_upper_triangular {n : ℕ} (M : matrix (fin n) (fin n) R)
  (h : ∀ (i j : fin n), j < i → M i j = 0) :
  M.det = ∏ i : (fin n), M i i :=
begin
  convert det_of_block_triangular_matrix' M id h,
  ext i,
  have h2 : ∀ (j : {a // id a = i}), j = ⟨i, rfl⟩ :=
    λ (j : {a // id a = i}), subtype.ext j.property,
  haveI : unique {a // id a = i} := ⟨⟨⟨i, rfl⟩⟩, h2⟩,
  simp [h2 default]
end

lemma det_of_lower_triangular {n : ℕ} (M : matrix (fin n) (fin n) R)
  (h : ∀ (i j : fin n), i < j → M i j = 0) :
  M.det = ∏ i : (fin n), M i i :=
begin
  rw ← det_transpose,
  exact det_of_upper_triangular _ (λ (i j : fin n) (hji : j < i), h j i hji)
end

end matrix

subset_name: curated


TITLE: Generating a random number from a given distribution
QUESTION [2 upvotes]: I have a problem (a part of a Monte Carlo simulation) where I'm given the energy of an incoming particle, $\varepsilon$ and want to split this energy in two parts, randomly generating the fraction that goes into each part according to a given function
$P(\varepsilon_s) \propto %\frac{1}{\bar E \arctan\left(\frac{\varepsilon}{\bar E}\right)}
\frac{1}{1+\left(\frac{\varepsilon_s}{\bar E}\right)^2}$
$\bar E$ is a known shape parameter. Now, the incoming particle can have energies ranging from 0 to a known maximum, $\varepsilon \in \left[0,\varepsilon_{\text{max}}\right]$, but I do not know the energy until the collision event, when I want to split it up, and therefore I can't figure out how to normalize the distribution function properly.
I have a pseudo-random number generator at my disposal which generates uniformly distributed numbers on $[0,1)$. How do I generate random numbers that properly describe how the two collisional products share the available energy?
I've tried reading up on methods for this in available books, but they're all a little to theoretical for me to be able to apply them to my scenario. It's been way too long since I took a proper statistics course, and now it seems I need it :)
Update: As requested, a more thorough description of $P$ (which will get into detail about the physics of this experiment...) follows here:
The physics of this problem is ionization - I have an incoming electron with kinetic energy $\varepsilon$, that collides with a molecule and ionizes it. After the collision event, I have two electrons with energies $\varepsilon_p$ (primary) and $\varepsilon_s$ (secondary), with $\varepsilon_p+\varepsilon_s=\varepsilon-\varepsilon_i$, the latter being the (known) ionization energy. Now, I don't really care which electron is which after the collision, so determining the two energies is simply a problem of dividing $\varepsilon-\varepsilon_i$ into two parts.
To do this, I am using the differential cross section of the collision, as suggested by Y. Tzeng and E. E. Kunhardt1, which is given by them as
$q(\varepsilon, \varepsilon_s) = \frac{C(\varepsilon)}{\bar E \arctan\left((\varepsilon-\varepsilon_i)/\bar E\right)}\cdot \frac{1}{1+\left({\varepsilon_s/\bar E}\right)^2}$
$C(\varepsilon)$ in this epxression is the total collision cross section, which is also known. I'm not 100% sure on this, but if I've understood the physics correctly then normalizing this function over the available energies should give me a probability density. As it turn out, this funciton is normalized so that
$\int_0^{\varepsilon-\varepsilon_i} q(\varepsilon, \varepsilon_s)d\varepsilon_s = C(\varepsilon)$
so normalizing w.r.t. the entire interval is just a matter of removing the factor $C(\varepsilon)$. However, if I insert $\varepsilon=\varepsilon_{\text{max}}$, that gives me a probability distribution which allows for energies higher than the available energy, and if I don't then I don't know how to generate random numbers from the distribution correctly.

1 I don't think I'm at liberty to relay the entire paper here, unfortunately...

REPLY [3 votes]: If my calculations are correct, the normalizing constant is $c \left(
\varepsilon_{\max}, \bar{E} \right) = \bar{E} \arctan \left(
\frac{\varepsilon_{\max}}{\bar{E}} \right)$, the CDF is $F \left( \varepsilon_s
\right) = \frac{\arctan \left( \frac{\varepsilon_s}{\bar{E}} \right)}{\arctan
\left( \frac{\varepsilon_{\max}}{\bar{E}} \right)}$ and thus the inverse CDF
is
$$ \varepsilon_s = F^{- 1} \left( u \right) = \bar{E} \tan \left( \arctan \left(
   \frac{\varepsilon_{\max}}{\bar{E}} \right) u \right) $$
Thus draw a uniform number on $[0,1]$, set it to $u$ and compute the random draw $\varepsilon_s$ from it using the function $F^{-1}$.


TITLE: Fourier transforming between time space and wavelength space
QUESTION [3 upvotes]: Recently, we talked about a homework assignment, which seemed rather odd to me.
Let the Fourier transformed $ F(\xi) $ of a function $f(x)$ be defined via:
$$F(\xi) = \int_{-\infty}^{\infty}f(x) e^{-\mathrm{i} x\xi}\ \mathrm{d} x $$
Now, the assignment is similar to:
Let
$$
S(\lambda) = A \mathrm{e}^{-\frac{1}{2} \left(\frac{\lambda-\lambda_0}{\Delta\lambda}\right)^2}
$$
be the spectrum of a light pulse. Calculate the Forier transformed of $S(\lambda)$ ($\rightarrow$ time space).
It seems to imply that I should plug $S(\lambda)$ into the definition of the Fourier transformed (or the reverse transformation), so that in the above defintion:
$$x = \lambda $$
$$\xi = t$$
I noticed some problems with this approach. For example, the exponent $\mathrm{e}^{-\mathrm{i}\lambda t}$ is not dimensionless.
As a response to the question, whether this is allowed, the response was simply "Why shouldn't it"?
My question is: Is this approach valid?
EDIT
We finally tried to solve the following integral:
$$F(t) = \int_{-\infty}^{\infty}S(\lambda)\ e^{-\mathrm{i} \lambda t}\ \mathrm{d} \lambda $$
To rephrase the question: Is this the correct time space function corresponding to $S(\lambda)$? Or do you need $\omega$ to get the time space function?

REPLY [3 votes]: Your approach is vaguely valid, but you should have a big red flashing light that says "dimensional analysis" on it: $\lambda t$ has dimensions $[\lambda t]=[LT]$ and it has no business being exponentiated as is.
What should be on that exponent is something like $\omega t$, and you need to transform your spectrum from wavelength to frequency space via $\omega = 2\pi c/\lambda$, but what this really brings up is a flaw in your assignment: normally, when we talk about pulses having a gaussian spectrum, we make them gaussian in frequency space: they're gaussians in the wavevector $k$,
$$
S(k) = B \exp\left(-\frac12\left(\frac{k-k_0}{\Delta k}\right)^2\right), \tag1
$$
not over its inverse, the wavelength.
As such, if you want to plow on with the pulse you've been given, it's definitely something you can attempt (once you fix that non-dimensionless exponent), but it won't be particularly pretty.
What does normally get the time on the spotlight, on the other hand, is the $k$-space spectrum above, which does have a nice relationship with the corresponding 1D waveform
$$
F(x,t) = \int S(k) e^{ik(x-ct)} \mathrm dk, \tag2
$$
which you can integrate explicitly, and whose relationship with $S(k)$ (and the equivalent $S(\omega)$) is important to keep in mind in many contexts.

Edit: perhaps it's worth going into a little bit more depth into the relationship between the two spectra at hand: the $k$-space gaussian
$$
S_k(k) = B \exp\left(-\frac12\left(\frac{k-k_0}{\Delta k}\right)^2\right), \tag1
$$
from above, and the $\lambda$-space gaussian
$$
S_\lambda(\lambda) = A \exp\left(-\frac12\left(\frac{\lambda-\lambda_0}{\Delta \lambda}\right)^2\right), \tag3
$$
which when translated to $k$-space via $k=2\pi/\lambda$ gives
$$
S_\lambda(k) = A \frac{\sqrt{2\pi}}{k}\exp\left(-\frac{1}{2(\Delta \lambda)^2}\left(\frac{1}{k}-\frac{1}{k_0}\right)^2\right).\tag4
$$
The spectrum in $(4)$ can in principle be used as in $(2)$, and it will give you some pulse, but you're highly unlikely to be able to integrate it effectively. In principle, it is a different function to the $k$-space gaussian $(1)$, but in some parameter regimes it can look similar. To emphasize how this works, consider the following plots: here we have $S_\lambda(k)$ on the left and $S_k(k)$ on the right, for different values of $\Delta k$ - a broadband pulse on the top, with $\Delta k = \frac{1}{5}k_0$, and one with a narrower bandwidth on the bottom, with $\Delta k = \frac{1}{15}k_0$.

Mathematica graphics through Import["http://halirutan.github.io/Mathematica-SE-Tools/decode.m"]["http://i.stack.imgur.com/Ad7vo.png"]
As you can see, for broadband pulses the functions are obviously different, but as the bandwidth becomes narrower, you just don't have that much bandwidth to explore the nonlinearities of the $k\propto 1/\lambda$ transformation.
So, what is the take-home lesson here? I want to emphasize that gaussian pulses are ultimately just models, and they are useful models in as much as they are both realistic and solvable. In the real world, you are unlikely to encounter a pulse that has a truly gaussian spectrum (unless you work very hard at it) with the kind of broad bandwidth in the top row.
It helps to put in some numbers: say that you have a pulse at optical frequencies that lasts about a nanosecond, which sounds really short. In this case, the bandwidth will be in the GHz regime, but the frequency will be in the order of $10^{14}\: \mathrm{Hz}$, and the bandwidth will look something like $\delta k \sim 10^{-5}k_0$; as far as the plots above go, they would both be indistinguishable delta-like spikes at that resolution. To get spectra that are as broad as the top row you would need an octave-spanning pulse that was shorter than a few tens of femtoseconds; this is nowadays doable but requires lots of dedicated effort.
Ultimately, then, what I'm saying is this: it can be worthwhile to consider what would happen with the time-domain waveform of $S_\lambda(\lambda)$, but if you really want to do it then first you need to be answering some hard questions about what bandwidth you're actually interested in, whether you're not served equally well by a standard $k$-space gaussian, and (if the two spectra are actually different at your bandwidths of interest) just how realistic a model $S_\lambda(\lambda)$ is for the phenomenon of interest, and whether there isn't some better model that you should be calculating.


TITLE: Using interpolation in Sobolev Space
QUESTION [1 upvotes]: good afternoon. I'm studying an article and the author uses a notation that I don't know and this is making it difficult for me to understand an inequality.
Consider a function $\eta$ in $H^{1}(0,1)$, that is, $\eta \in L^{2}(0,1)$ and $\eta_{x} \in L^{2} (0.1)$. So, the author claims that by interpolation, we have
$$
\|\eta\|_{L^{2}(0,1)} \leq C\|\eta\|_{-1,2}\|\eta\|_{1,2}
$$
where $C >0$.
I don't know what the meaning of $\| \bullet\|_{-1,2}$. Does anyone know a reference that addresses this notation? And what about inequality? Is it some classical theorem?

REPLY [2 votes]: I am following Adams and Fournier's Sobolev Spaces.
In general, $W^{-k,p}(\Omega)=(W_{0}^{-k,p}(\Omega))'$ (well, this is not the precise definition, but you can think this space in this form because they are isometrically isomorphic, see Theorem 3.12 from the book that I mentioned or Theorem 11.62 of Leoni's Sobolev Spaces 2nd ed), where $'$ is used to denote the dual space. In your case, $H^{-1,2}(0,1)=(H_{0}^{1,2}(0,1))'$.

What is the norm here? and there make sense to compute $\|\eta\|_{-1,2}$ for $\eta \in L^{2}$?: for each $v \in L^2(0,1)$ define
$$L_{v}(u)=\langle u,v\rangle=\int_{\Omega}u(x)v(x)dx,$$
where $u \in H_{0}^{1,2}(0,1)$
Note that $L_{v}$ makes sense as an element of $H^{-1,2}$ because
$$|L_{v}(u)|=|\langle u,v \rangle| \underbrace{\leq}_{\text{Holder}} \|u\|_{L^{2}(0,1)} \|v\|_{L^{2}(0,1)} \leq \|u\|_{H^{1,2}(0,1)} \|v\|_{L^{2}(0,1)}<\infty.$$
Now, define the norm of $v$ in this space as the operator norm of $L_{v}$, that is,
$$ \|v\|_{H^{-1,2}(0,1)}=\sup_{{u \in H_{0}^{1,2}(0,1) \\ \|u\|_{H_{0}^{1,2}(0,1)} \leq 1}} |L_{v}(u)|.$$

How do you obtain the required inequality?: Let $u \in H_{0}^{1,2}(0,1)$, $v \in L^{2}(0,1)$ we have
$$|\langle u,v \rangle|=\|u\|_{H_{0}^{1,2}(0,1)} \Bigg| \left\langle \frac{u}{\|u\|_{H_{0}^{1,2}(0,1)}},v \right\rangle \Bigg| \leq \|u\|_{H_{0}^{1,2}(0,1)} \|v\|_{H^{-1,2}(0,1)}.$$
So, taking $u=v=\eta$, you obtain the desire inequality.

Finally, let me comment that this inequality is "well-known" (for people that works on PDE), and sometimes is called "Sobolev duality inequality" or "generalized Hölder inequality".


TITLE: Equation with a sum for the prime-counting function involving the Mobius function
QUESTION [1 upvotes]: I have come across the statement that
$$
 \sum_{n\leq x}\sum_{d\mid(n,P_z)}\mu(d) = \sum_{d\mid P_z}\mu(d)
\left[\frac{x}{d}\right],
$$
where $P_z=\prod_{p\leq z}p$ where $p$ is prime, $\mu(d)$ is the Mobius function, and $[x]$ denotes the greatest integer less than or equal to $x$. I would appreicate some help in seeing the above. (Both of the above are equal to $\pi(x,z)$, the number of (positive) integers less than or equal to $x$ which are prime to all primes less than or equal to $z$.)
I am aware of the property that $\sum_{d\mid n}\mu(d)=0$ unless $n=1$ in which case, the sum is equal to $1$. I can see how for a given $d$, there are $[x/d]$ integers less than or equal to $x$ with this divisor.
Thank you!

REPLY [2 votes]: For every $d$, count how often $\mu(d)$ occurs in the sum on the left.
If it occurs at all, there is an $n\leqslant x$ with $d\mid (n,P_z)$. Since $(n,P_z) \leqslant n \leqslant x$, it follows that $d\leqslant x$. Also, since $d$ divides $(n,P_z)$, it divides $P_z$.
Now, if $d \leqslant x$ and $d\mid P_z$, it occurs once in the inner sum for every $n \leqslant x$ such that $d \mid n$. There are exactly $\left[\frac{x}{d}\right]$ such $n$. For the divisors of $P_z$ that are greater than $x$, we have $\left[\frac{x}{d}\right] = 0$, and above we saw that $\mu(d)$ does not occur in the sum on the left for these $d$, so $\mu(d)$ occurs $\left[\frac{x}{d}\right]$ times in the sum on the left for all divisors of $P_z$, whether $\leqslant x$ or $> x$, i.e.
$$\sum_{n\leqslant x} \sum_{d\mid (n,P_z)} \mu(d) = \sum_{d\mid P_z} \mu(d)\left[\frac{x}{d}\right].$$


text:
\section{Table of notation}
\label{app:notation}
\begin{tabular} {l l}
    \hline
    $\aunc$         & uncertain parameter\\
    $F_n^{1-\alpha}$& CDF of the $\chi^2$\\
    $u$             & dual variable\\
    $\x, \y$        & decision variable vectors\\
    $\y_i$          & subset of decision variables $\y$\\
    $\z$            & observation vector in original space\\
    $f(\cdot)$      & black-box objective function\\
    $g(\cdot)$      & black-box constraint\\
    $h(\cdot)$      & warping function\\
    $K(\cdot,\cdot)$& kernel function of GP\\
    $w(\cdot)$      & constraint defining $\U$\\
    $\A$            & $\alpha$-confidence ellipsoid\\
    $\U$            & (warped) uncertainty set\\
    $\X$            & deterministic feasible set\\
    $\alpha$        & probability of constraint violation\\
    $\vec{\delta}$  & disturbances vector\\
    $\epsilon$      & estimated feasibility\\
    $\epsilon_0$    & target feasibility\\
    $\xivec$        & observation vector in latent space\\
    $\vec{\psi} = \{a_j, b_j, c_j\}$ & parameters of warping function\\
    $\muvec$        & mean of GP at $\y_i$\\
    $\sigma_{ij}^2$ & $ij$-element of covariance matrix\\
    $\Sigma$        & covariance matrix of GP at $\y_i$\\
    \hline
                    & \emph{Production planning}\\
    \hline
    $c_t$           & production cost in period $t$\\
    $\punc_t$       & uncertain price in period $t$\\
    $x_t$           & production amount in period $t$\\
    \hline
                    & \emph{Drill scheduling}\\
    \hline
    $W$             & weight on bit\\
    $\dot{N}$       & rotational speed\\
    $V$             & rate of penetration\\
    $\Dp$           & differential pressure\\
    $R$             & degradation indicator\\
    $M$             & set of maintenance depths\\
    \hline
\end{tabular}

config: arxiv