text:
The Prisoner's Dilemma

D4. Four Operations Revisited

One of the great things about mathematics is that it’s always possible to look again at familiar things and discover something new. Sometimes this can happen when you combine well-known things in a new way. For example, you’re probably familiar with the four basic operations of addition, subtraction, multiplication, and division. And you may have heard about various number sequences, such as the powers of two: 1, 2, 4, 8, 16, and so on. In that sequence, you obtain each term by adding the previous term to itself, or in symbols $a_n = a_{n-1} + a_{n-1}$. Another similar sequence is the sum recurrence, also known as the Virahanka or Fibonacci sequence, in which you obtain each term by adding the previous term to the one before it, or $b_n = b_{n-1} + b_{n-2}$. The sum recurrence begins with 1, 2, 3, 5, 8, and so on.

This dilemma challenges you to describe the long-term behavior of four different sequences. Unlike sequences $a$ and $b$ which just add one of the last two terms unchanged to the previous one, each of these sequences uses one of the four basic operations on the last two terms to determine what to add to the previous term. But first, what is meant by “long-term behavior?” The gold standard is a closed-form formula for the sequence — an expression for the $n$th term that just involves $n,$ not previous terms. For example, the closed form for the powers of two is $a_n = 2^n$. If there doesn’t seem to be a closed form, maybe you can find a formula that works when $n$ is bigger than some number, but might be off for the first few values. Or sometimes you can find a formula that is never exactly right, but gets closer and closer as $n$ gets large. An example of this possibility is the formula $b_n \approx \phi^{n+1}/\sqrt5$ for the sum recurrence, where $\phi = (1+\sqrt5)/2$ is the golden ratio. All else failing, it could be one formula that’s always less than the sequence and another that is always greater — and hopefully you can get those two “bounding” formulas as close to each other as possible.

So all that said, can you describe the long term-behavior of each of the four following sequences? Each sequence starts with 1, then 2 (like $a$ and $b$ do), and thereafter adds the described quantity to the previous term to obtain the next term. If the described quantity is ever not a whole number, it should be rounded up to the next whole number.

1. One less than the sum of the previous two terms. So to obtain the third term of this sequence, you take 1+2-1 = 2, and add that to 2, to get 4.
2. One more than the difference of the previous two terms. This time, the amount you add to get the third term is 2-1+1 = 2, so the third term is again 4. (But don’t worry, this sequence will veer off from the first sequence soon!)
3. Half the product of the previous two terms. So in this case our next increment is 2·1/2 = 1, and the third term will be 3.
4. Twice the quotient of the previous two terms. Now the next increment is (2/1)·2 = 4, and the third term will be 6.

Which of these four sequences ultimately grows the slowest? You will see that by using each operation as part of an infinite cycle, we can get new, different behaviors. Happy sequence hunting!

Update: (2022 Sep 15) Submissions have only been received related to sequences 1, 2, and 4 above. You can see the solution to the first one below, and the solutions to 2 and 4 will be posted as soon as the Fall 2022 PMP Newsletter is out. Therefore, submissions will still be accepted for the remaining sequence, number 3, through the deadline for the other new 2022 Fall problems.

This problem originally appeared in the Prisoner’s Dilemma in the 2022 Winter issue of the PMP Newsletter. Solutions will be accepted through 2022 Dec 31.

Show solution?
Solution.

As submitter Jesse Waite noted, a good first step in attacking one of these sequences is to write out (in symbols) the rule defining it, known as its recurrence relation. For example, if we denote the first series connected with addition by $c_0, c_1, c_2,\ldots,$ its recurrence relation is $c_{n+2} = c_{n+1} + c_n – 1 + c_{n+1} = 2c_{n+1} + c_n – 1.$

If the formula doesn’t ring a bell, a next good thing to do is write out the first several terms: $c_1 = 1,$ $c_2 = 2,$ $c_3 = 4,$ $c_4=9,$ $c_5=21,$ $c_6=50$, $c_7=120,$ $c_8 = 289,$ and so on.

If the sequence still isn’t looking familiar, a good third step to try is to look at its “first difference sequence,” which just consists of the numbers you get by subtracting from each entry the one before it (you have to start from the second entry to do this). In this case, the sequence of first differences is $1, 2, 5, 12, 29, 70, 169, 408, \ldots$

And now, if you have some experience with sequences, you could recognize these numbers as the Pell Sequence $p_n$, a close cousin of the familiar sum recurrence (also known as the Virahāṅka-Fibonacci sequence), in which $p_{n+2} = 2p_{n+1} + p_n$ (instead of having $b_{n+2} = b_{n+1} +b_n$ as in the sum recurrence). Another way you could figure this out is by taking these first differences symbolically. If we let $c’_n = c_{n+1} – c_n,$ then $\begin{split}c’_{n+2} & = c_{n+3} – c_{n+2} = 2c_{n+2} + c_{n+1} – (2c_{n+1} +c_n -1) \\ & = 2(c_{n+2}-c_{n+1}) + (c_{n+1} – c_n) = 2c’_{n+1} + c’_n,\end{split}$ exactly the recurrence for the Pell sequence.

So if the first differences of our sequence $c_n$ are the Pell sequence, then our sequence must be the sum of the first $n-1$ terms of the Pell sequence (well, except possibly for an initial constant, which turns out to be one). And now you can look up or use standard linear recurrence methods to find out that the formula for $p_n$ is $p_n = \frac{(1+\sqrt2)^n-(1-\sqrt2)^n}{2\sqrt2}.$

We are only interested in the long term behavior, so note that since the absolute value of $1-\sqrt2$ is less than one, $(1-\sqrt2)^n$ goes to zero as $n$ grows without bound. Therefore we can ignore the second term in the numerator, and letting the constant $k = 1/(2\sqrt2),$ we have that $p_n\approx k(1+\sqrt2)^n.$ Therefore, $c_n\approx k\sum_{i=0}^{n-1} (1+\sqrt2)^i.$ But this sum is geometric, which has a convenient formula. So $c_n\approx k\frac{(1+\sqrt 2)^n-1}{(1+\sqrt2)-1} = \frac{(1+\sqrt2)^n-1}4.$

This expression gives a very good picture of the long-term behavior of the first sequence $c.$ (With a little more work along these same lines — don’t ignore the small term in the numerator and be careful with all of the constants — one can actually obtain an exact closed-form formula for $c_n.$)


text:
#### Definition

Angles can be measured in units of degrees or radians. A complete revolution is defined as $360^\circ$ which is equal to $2\pi$ radians

$360^\circ = 2\pi \text{ radians.}$

From this, we can derive that

$1^{\circ} = \frac{\pi}{180} \text{ radians.}$ $1 \text{ radian} = \frac{180}{\pi}^{\circ}$

Sometimes a superscript c is used to denote radians instead of degrees, though it is conventional to assume that radians are used unless otherwise specified.

#### Worked Examples

###### Example 1

Convert $83^{\circ}$ to radians.

###### Solution

Recall that $1^{\circ} = \dfrac{\pi}{180} \text{ radians}$. So multiply $\dfrac{\pi}{180}$ by $83$:

\begin{align} 83^{\circ} &= 83 \times \frac{\pi}{180}\\ &\approx 1.449 \text{ radians (to 3 d.p.)} \end{align}

###### Example 2

Convert $3 \text{ radians}$ into degrees.

###### Solution

Using the definition, if $1 \text{ radian} = \dfrac{180}{\pi} \text{ degrees}$, then multiply by $3$ to find the angle in degrees.

$3 \times \frac{180}{\pi} \approx 172^{\circ} \text{ (to 3 sig.fig.)}$

#### Video Example

Prof. Robin Johnson shows how to convert $37^{\circ}$ to radians, and $1$ radian to degrees.

#### Workbook

This workbook produced by HELM is a good revision aid, containing key points for revision and many worked examples.

#### More Support

You can get one-to-one support from Maths-Aid.

format: plain


text:
# Constructions of factorization systems in categories

@article{Bousfield1977ConstructionsOF,
title={Constructions of factorization systems in categories},
author={Aldridge Knight Bousfield},
journal={Journal of Pure and Applied Algebra},
year={1977},
volume={9},
pages={207-220}
}
• A. K. Bousfield
• Published 1977
• Mathematics
• Journal of Pure and Applied Algebra
In [2] we constructed homological localizations of spaces, groups, and 17"modules; here we generalize those constructions to give "factorization systems" and "homotopy factorization systems" for maps in categories. In Section 2 we recall the definition and basic properties of factorization systems, and in Section 3 we give our first existence theorem (3.1)for such systems. It can be viewed as a generalization of Deleanu's existence theorem [5] for localizations, and is best possible although it… Expand
155 Citations
Topological and Simplicial Models of Identity Types
• Mathematics, Computer Science
• TOCL
• 2012
By exhibiting suitable path object structures on Top and SSet, this paper endow those categories with the structure of a homotopy-theoretic model and, in this way, obtain the desired topological and simplicial models of identity types. Expand
ON TOWERS APPROXIMATING HOMOLOGICAL LOCALIZATIONS
• Mathematics
• 1997
Our object of study is the natural tower which, for any given map f : A → B and each space X, starts with the localization of X with respect to f and converges to X itself. These towers can be usedExpand
The Frobenius condition, right properness, and uniform fibrations
• Mathematics
• 2017
We develop further the theory of weak factorization systems and algebraic weak factorization systems. In particular, we give a method for constructing (algebraic) weak factorization systems whoseExpand
Localization and genus in group theory and homotopy theory
To support the development in (1) we combine the works of Adams [1], Casacuberta-Peschke-Pfenniger [7] and chapter I of the present exposition. Thus we arrive at the following pivotal point. TheExpand
On Localization and Stabilization for Factorization Systems
• Mathematics, Computer Science
• Appl. Categorical Struct.
• 1997
It is found necessary and sufficient conditions for (ε′, M*) to be another factorization system, and it is shown that a number of interesting factorization systems arise in this way. Expand
Localization with respect to a class of maps I — Equivariant localization of diagrams of spaces
Homotopical localizations with respect to a set of maps are known to exist in cofibrantly generated model categories (satisfying additional assumptions) [4, 13, 24, 35]. In this paper we expand theExpand
Adams Memorial Symposium on Algebraic Topology: On orthogonal pairs in categories and localization
• Mathematics
• 1992
X with f ∈ M, there is a unique morphism h : B → X with hf = g. The orthogonal subcategory problem [13] asks whether D is reflective in C, i.e., under which conditions the inclusion functor D → CExpand
Cellular Spaces 0. Introduction and Main Results
The aim of the present work is to study, in a systematic way, the possibility of constructing new spaces out of a given collection of persumably better understood spaces. Recently there were severalExpand
A naive diagram-chasing approach to formalisation of tame topology
• Mathematics
• 2018
We rewrite classical topological definitions using the category-theoretic notation of arrows and are led to concise reformulations in terms of simplicial categories and orthogonality of morphisms,Expand
The orthogonal subcategory problem in homotopy theory
• Mathematics
• 2005
It is known that, in a locally presentable category, localization exists with respect to every set of morphisms, while the statement that localization with respect to every (possibly proper) class ofExpand

#### References

SHOWING 1-10 OF 16 REFERENCES
The localization of spectra with respect to homology
IN [8] WE studied localizations of spaces with respect to homology, and we now develop the analogous stable theory. Let Ho” denote the stable homotopy category of CW-spectra. We show that eachExpand
On Bicategories of Relations and Pullback Spans
In this paper, we study relations in general categories. Our approach requires that these must have finite products and factorization systems. Klein [5] has obtained a condition for composition ofExpand
Monomorphisms, Epimorphisms, and Pull-Backs
As the applications of category theory increase, we find ourselves wanting to imitate in general categories much that was at first done only in abelian categories. In particular it becomes necessaryExpand
Homotopical Algebra
• in: Lecture Notes in Mathematics, Vol. 43
• 1967
Bousfield , Homoiogical localization towers for groups and wmodules
• Memoirs Amer . Math . Soc . , to appear
• 1975
Deleanu , Existence of the Adams completion for objects of cocomplete categories
• J . Pure Appl . Algebra
• 1975
Homotopy Limits
• Completions and Localizations, in: Lecture Notes in Mathematics, Vol. 304
• 1972


text:
506
Q:

# Find the missing number in the series?4, 18, ?, 100, 180, 294, 448

 A) 48 B) 50 C) 58 D) 60

Explanation:

Take a Look at this series :

$23$ - $22$ =  8 - 4 = 4

$33$ - $32$ = 27 - 9 = 18

$43$ - $42$ = 64 - 16 = 48

$53$ - $52$ = 125 - 25 = 100

$63$ - $62$ = 216 - 36 = 180

$73$ - $72$ = 343 - 49 = 294

$83$ - $82$ = 512 - 64 = 448

Q:

3, 6, 17.5, 64.75, ?

 A) 304.587 B) 295.875 C) 115.246 D) 101

Explanation:

3

3 x 1.5 + 1.5 = 6

6 x 2.5 + 2.5 = 17.5

17.5 x 3.5 + 3.5 = 64.75

64.75 x 4.5 + 4.5 = 295.875

7 463
Q:

Find the next number in the given sequence?

1, 2, 6, 21, 88, ?

 A) 445 B) 364 C) 296 D) 194

Explanation:

The given series is 1, 2, 6, 21, 88, ?

1

1 x 1 + 1 = 2

2 x 2 + 2 = 6

6 x 3 + 3 = 21

21 x 4 + 4 = 88

88 x 5 + 5 = 445

Hence, the next number in the given number sequence is 445.

5 1539
Q:

Find the next number in the series given 18, 90, 450, 2250, ?

 A) 4520 B) 6520 C) 7250 D) 11250

Explanation:

The given number series is 18, 90, 450, 2250, ?

18

18 x 5 = 90

90 x 5 = 450

450 x 5 = 2250

2250 x 5 = 11250

Hence, the next number in the given series is 11250.

7 1379
Q:

In following question, a number series is given with one term missing. Choose the correct alternative for these sequence.

24, 60, 120, 210, ?

 A) 336 B) 330 C) 300 D) 370

Explanation:

3 x 3 x 3 = 27 - 3 = 24

4 x 4 x 4 = 64 - 4 = 60

5 x 5 x 5 = 125 - 5 = 120

6 x 6 x 6 = 216 - 6 = 210

7 x 7 x 7 = 343 - 7 = 336

7 474
Q:

What would be the next number in the sequence 52, 56, 48, 64, 32, ?

 A) 96 B) 48 C) 54 D) 81

Explanation:

The given number series is 52, 56, 48, 64, 32, …

It follows a pattern that,

52

52 + 4 = 56

56 - 8 = 48

48 + 16 = 64

64 - 32 = 32

32 + 64 = 96

Therefore, the required next number in the given number series is 96.

4 402
Q:

Find the next number in the given series?

118, 166, 190, 202, 208, ?

 A) 211 B) 222 C) 233 D) 244

Explanation:

The given series is 118, 166, 190, 202, 208, ?

118

118 + 48 = 166

166 + 24 = 190

190 + 12 = 202

202 + 6 = 208

208 + 3 = 211

2 1275
Q:

2, 16, 3, 81, 4, ?

 A) 99 B) 125 C) 256 D) 396

Explanation:

The given number series is 2, 16, 3, 81, 4, ?

2

2^4 = 16

3

3^4 = 81

4

4^4 = 256

Hence, the next number is 256.

2 496
Q:

Find the next number in the given number sequence

200, 200, 100, 300, 75, ?

 A) 150 B) 175 C) 250 D) 375

Explanation:

The given number series is 200, 200, 100, 300, 75, ?

It follows a pattern that,

200

200 x 1 = 200

200 / 2 = 100

100 x 3 = 300

300 / 4 = 75

75 x 5 = 375

Hence, the next number in the given number sequence is 375.


text:
# Alignment of colorbar with colorbar to name

I am trying to make a group of plots with groupplots sharing the same colorbar and the same y axis. An example:

\documentclass{standalone}
\usepackage{pgfplots}
\usepgfplotslibrary{groupplots}
\begin{document}
\begin{tikzpicture}
\begin{groupplot}[%
point meta min=0, point meta max=25,
group style={group size=2 by 1,
yticklabels at=edge left, ylabels at=edge left,
},
colorbar,
colorbar to name=commonbar,
]
\nextgroupplot
\nextgroupplot
\end{groupplot}
\end{tikzpicture}
\ref{commonbar}
\end{document}

Here, the colorbar to name option takes care of using the same colorbar for all plots (there is a warning from pdflatex on "multiply defined labels" -- to be expected, I suppose, and no problem as long as the range of the colorbar is specified). The colorbar is then "placed" by the \ref{commonbar} command.

How can I align the colorbar properly with the figures? I have tried placing the \ref command in several places (the placement above gives the best result I've seen so far). You may say I'm looking for the equivalent of yticklabels at=edge left for colorbars.

-

The warning about the multiply defined labels comes from the fact that you supplied the colorbar to name key to the groupplot environment, which in turn applies it to each of the plots in the group, so it is executed twice in this case. If you just supply colormap to the last \nextgroupplot[<options>], the message will disappear and the colormap will be placed as desired, without using colorbar to name and \ref. If you do want to use colormap to name, Percusse's answer is a good way to go.

\documentclass[border=4mm]{standalone}
\usepackage{pgfplots}
\usepgfplotslibrary{groupplots}
\begin{document}
\begin{tikzpicture}
\begin{groupplot}[%
point meta min=0, point meta max=25,
group style={
group size=2 by 1,
yticklabels at=edge left,
ylabels at=edge left
}
]
\nextgroupplot
\nextgroupplot[colorbar]
\end{groupplot}
\end{tikzpicture}
\end{document}

-
Thanks, exactly the type of solution I was looking for! Does this mean that all (or perhaps a large, meaningful subset) of the options that can be applied to the axis environments also can be applied to \nextgroupplot? –  josteinb Apr 17 '12 at 7:44
@JosteinBøFløystad: Yes, all the options that can be supplied to an axis can also be applied to nextgroupplot. The groupplots environment is basically an elaborate wrapper for several full-fledged axis environments. –  Jake Apr 17 '12 at 7:50

This just works and there is manual adjustment involved. Probably, Christian or Jake would give you an automated solution. I put the commonbar into a node and placed it next to the rightmost groupplot.

\documentclass{article}
\usepackage{pgfplots}
\usepgfplotslibrary{groupplots}
\begin{document}
\begin{tikzpicture}
\begin{groupplot}[%
point meta min=0, point meta max=25,
group style={group size=2 by 1,
yticklabels at=edge left, ylabels at=edge left,
},
colorbar,
colorbar to name={commonbar},
]
\nextgroupplot

For some reason the labels of the colorbar moves too close to the colormap. I don't know why exactly this is happening.But I'll look into it. I used inner sep=0 and that is inherited by colormap. Removing that option fixes that problem.


text:
# Tag Info

Using $\mathbb E$ for the expected value symbol, $E_R$ and $V_R$ for the mean and the avriance of returns $R$, for a utility function of the form $$U(R) = \ln(1+R)$$ the second-order Taylor expansion gives (ignoring the remainder) $$U(R) \approx \ln(1+E_R) + \frac{1}{1+E_R} (R- E_R) -\frac{1}{2(1+E_R)^2}(R-E_R)^2.$$ Then \mathbb E[U(R)] \approx E_R + 0 - \...


text:
# WAN

Wide Area Network. A network that spans a great distance, usually over a public network, such as telephone. WANs are typically 2 or more connected LANs. The Internet is the largest.


text:
# Is there an example of a linear map between Banach Spaces that is not closed?

Is there a simple example of a linear map $$T: X \to Y$$ between Banach spaces that is not closed, i.e. if $$x_n \to 0$$ in $$X$$ and $$T(x_n) \to y$$ in $$Y$$, then $$y=0$$? I know the Closed Graph Theorem says that if $$T: X \to Y$$ is a linear map between Banach spaces, then \begin{align*} T \text{ is bounded } \iff T \text{ is closed } \iff \text{gr}(T) \text { is a closed subspace of } X \oplus Y. \end{align*} So I would have to find an unbounded linear map between Banach spaces, correct? What exactly would that be?

We can construct an unbounded linear functional easily. Pick your favorite basis $$e_n$$ with $$\lVert e_n \rVert = 1$$ of $$X$$ (which can be done with the axiom of choice) and define $$f$$ as $$f(x)=\sum_n n x_n,\,\,\,\,\text{where }\,\,\,\,x=\sum_n x_n e_n.$$ Note that this is well-defined, as $$x_n$$ is non-zero for finitely many $$n$$ for every $$x$$. Then pick $$y_n = \frac{e_n}{\sqrt{n}},$$ which clearly goes to $$0$$, but $$f(y_n)=\sqrt{n} \to +\infty$$.
• What exactly would $X$ and $Y$ be? Are there any concrete Banach spaces that this works for? Mar 1, 2021 at 21:37
• @MATH-LORD Pick your favorite Banach space $X$ over the field $\mathbb{F}$ and $Y=\mathbb{F}$. Note that you can't construct this function "exactly", but you can prove that it exists. Mar 1, 2021 at 21:39


text:
# Substitution for impulse in recurrence equation

1. Apr 13, 2008

### yoran

Hi,

I was wondering what substitution to try when finding a particular solution for a recurrence equation with a linear combination of impulses on the right side of the equation. I think an example will clarify this.
Given the recurrence equation
$$y[k+2]-4y[k]=\delta [k] + 2\delta [k+1]$$
with $$\delta [k]$$ being the discrete impulse, i.e $$\delta [k] = 1$$ for $$k=0$$ and $$\delta [k] = 0$$ for $$k\neq 0$$.
First you must find a solution for the homogeneous equation. I have no problems with that. But I don't know what substitution to use for finding the particular solution? Can someone help me with that?

Thanks,

Yoran


text:
#### Morphology and physical properties of polyacrylonitrile fibers spun from a spinning dope with or without zinc chloride = 염화아연을 포함 혹은 미포함한 방사용액에서 방사된 폴리아크릴로나이트릴 섬유의 형태 및 물성

Cited 0 time in Cited 0 time in
• Hit : 168
All commercial PAN homopolymers are generally manufactured by the radical polymerization of acrylonitrile(AN) in solvent with an initiator. PAN has the strong intermolecular dipolar bonds between the nitrile groups on adjacent polymer chains, and thus it has a much higher melting temperature (317$^\circ$C) than the general purpose polymers such as polyethylene(PE), polypropylene(PP), and polyvinylchloride(PVC). The solvent of PAN must have the functional groups that disrupt the strong dipolar bonds. Dimethylformamide(DMF), dimethysulfoxide(DMSO), and dimethyl acetamide(DMAc) can be used as a solvent. PAN homopolymer was prepared in the conventional was by the radical polymerization of AN in DMF solution with azobisisobutyronirile(AIBN). The molecular weight and polydispersity index were about $4\times10^5$ and 2.16, respectively. PAN fibers are generally manufactured by a wet spinning, a dry-jet wet spinning, or a dry spinning. The most successful techique among them is the dry-jet wet spinning because PAN fibers produced by this techique have better mechanical properties. $ZnCl_2$ was added to a spinning dope in order to improve the mechanical properties of PAN fibers and the addition effect was analyzed. When the amount of $ZnCl_2$ was 0.01, 0.02, 0.05, and 0.10 mole $ZnCl_2$ per mole AN units, the residual amount of $Zn^{+2}$ in PAN fiber was $6.5\times10^{-4}$, $6.8\times10^{-4}$, $7.6\times10^{-3}$, and $9.3\times10^{-3}$ wt\%, respectively. As the residual amount of $Zn^{+2}$ increased, the maximum degradation temperature(DTmax) dropped to a lower one because the residual $Zn^{+2}$ accelerated the thermal degradation of PAN fibers. The complexations between the CN groups and $Zn^{+2}$ were investigated using a Raman and Infraed Spectroscopy. Any shifts in the wavenumber $\nu_{CN}$ for CN stretching were not observed for all experiments. It was understood that the complex of CN group and $Zn^{+2}$ could be neglected because the amount of $ZnCl_2$ added to ...
Chang, In-JaeLee, Wha-Seop정인재이화설
Description
한국과학기술원 : 화학공학과,
Publisher
한국과학기술원
Issue Date
1993
Identifier
68183/325007 / 000845619
Language
eng
Description

학위논문(박사) - 한국과학기술원 : 화학공학과, 1993.8, [ xvi, 153 p. ]

URI
http://hdl.handle.net/10203/28658


text:
Men are worse drivers, but only slightly¶

This is a Jupyter Notebook for analyzing the traffic collision data available for Los Angeles County through TIMS. The analysis uses the Pandas package in Python:

In [3]:
import pandas

The data can be downloaded from TIMS with a user account. Only ~5000 collisions may be downloaded at a time. This allows for about one month worth of data for Los Angeles County.

Here, the data for the parties involved in traffic collisions fron January 2012 in Los Angeles County is used.

In [4]:
#Load the data, show data size

data.shape

Out[4]:
(8893, 34)

There were 8893 parties in recorded accidents in Los Angeles County in Jan 2012. This database records all parties involved in a collision: perpetrator, victim, affected bystanders etc.

In [5]:
# What does the data look like? limit to first five rows

data.iloc[:5,:]

Out[5]:
CASEID PARNUM PTYPE ATFAULT PSEX PAGE PSOBER PDRUG PDIRECT PSAFETY1 ... VEHYEAR VEHMAKE VEHTYPE CHPTYPE1 CHPTYPE2 PRACE INATTN TRUCKCMB TRAILRCMB YEAR_
0 5980867 1 1 Y M 57 A - S L ... 1997 HONDA A 01 00 B NaN - - 2012
1 5980827 2 1 Y M 18 A - S L ... 1996 FORD A 01 00 W NaN - - 2012
2 5980827 1 1 N F 56 A - S M ... 2001 KIA A 07 00 B NaN - - 2012
3 5980823 2 1 N F 31 A - N M ... 2005 NISSAN A 01 00 O NaN - - 2012
4 5980823 1 2 Y M 50 G G W - ... 0 - N 60 00 B NaN - - 2012

5 rows × 34 columns

Querying a Pandas DataFrame¶

Pandas DataFrames can be queried in a way which is conceptually (though not syntactically) similar to SQL databases.

First, let's count up the total number of incidents (unique CASEID), and compare that to the number of parties.

In [5]:
cases = data.CASEID.nunique()
cases

Out[5]:
4146

So there were 4146 cases (collisions) and 8893 parties. That means an average of 8893/4146 = 2.14 parties per case.

Next, we will group by sex and fault, then count the entries with each unique combination of levels on these columns.

In [8]:
#Count parties by sex & fault

faultgender = data.groupby(['ATFAULT','PSEX']).count().iloc[:,[0]].rename(columns = {'CASEID':'counts'})
faultgender

Out[8]:
counts
ATFAULT PSEX
N - 481
F 1938
M 2577
Y - 297
F 1369
M 2231
In [9]:
fg = faultgender.reset_index()
fg = fg[fg.PSEX != '-']
fg

Out[9]:
ATFAULT PSEX counts
1 N F 1938
2 N M 2577
4 Y F 1369
5 Y M 2231
In [10]:
#relative probabilites?

float(2231)/(1369+2231)/0.5

Out[10]:
1.2394444444444443

Correcting for confounders¶

This would seem to indicate that men are ~24% more likely to cause accidents, based on the assumption that men and women are equal portions of the population. However, due to the confounder of men being more likely to drive (due to driving/delivery occupations etc.) this is probably not quite true.

Luckily, the data contains enough information to correct for this, at least partly. The database includes information about both the perpetrators and victims of crashes, which should correct for the hypothesis that there are more men driving. For each gender, we can calculate the probability of being at fault. Then we can calculate the relative risk of being at fault, conditioning on being male.

$\large RR = \frac{P(Fault | Male)}{P(Fault | Female)}$

In [11]:
# Get values from dataframes
NF = float(fg.counts[(fg.PSEX == 'F')&(fg.ATFAULT == 'N')].values[0])
NM = float(fg.counts[(fg.PSEX == 'M')&(fg.ATFAULT == 'N')].values[0])
YF = float(fg.counts[(fg.PSEX == 'F')&(fg.ATFAULT == 'Y')].values[0])
YM = float(fg.counts[(fg.PSEX == 'M')&(fg.ATFAULT == 'Y')].values[0])

#probability of being at fault for male parties
PM = YM/(YM+NM)
#probability of being at fault for female parties
PF = YF/(YF+NF)
#Relative risk of being at fault, given being a party in a traffic accident and being male
RR = PM/PF
[PM,PF,RR]

Out[11]:
[0.46401830282861894, 0.4139703658905352, 1.1208973903975479]

Conclusion...¶

So men are ~12% more likely to be at fault in a collision, at least during January 2012 in Los Angeles County. We can also plot data in this notebook as follows.

In [12]:
%matplotlib inline

import matplotlib.pyplot as plt
import numpy as np

N = 2
Men = (YM, NM)

ind = np.arange(N)  # the x locations for the groups
width = 0.35       # the width of the bars

fig, ax = plt.subplots()
rects1 = ax.bar(ind, Men, width, color='r')

Women = (YF, NF)
rects2 = ax.bar(ind + width, Women, width, color='y')

# add some text for labels, title and axes ticks
ax.set_ylabel('Count')
ax.set_title('Count by Fault & Gender')
ax.set_xticks(ind + width)
ax.set_xticklabels(('Fault','Not at Fault'))

ax.legend((rects1[0], rects2[0]), ('Men', 'Women'),loc=(0.25,0.75))

def autolabel(rects):
# attach some text labels
for rect in rects:
height = rect.get_height()
ax.text(rect.get_x() + rect.get_width()/2., 1.05*height,
'%d' % int(height),
ha='center', va='bottom')

autolabel(rects1)
autolabel(rects2)

plt.show()


text:
AAS Meeting #193 - Austin, Texas, January 1999
Session 47. Supernovae
Display, Thursday, January 7, 1999, 9:20am-6:30pm, Exhibit Hall 1

## [47.04] An Analysis of the Peculiar Type IIn Supernova 1995N

M.D. Baird (Grove City College), P. M. Garnavich, E.M. Schlegel, P. M. Challis, R. P. Kirshner (SAO)

SN 1995N is a peculiar type IIn supernova. Spectroscopic and photometric data for this analysis were gathered between May 10, 1995 (two days after discovery) and July 18, 1998. A total of twenty two photometric images and eight spectra were obtained at the FLWO and MMTO. The photometric data show a broad maximum at R=17.0 occurred in late October, 1995, followed by a very slow decline at a rate of 2.39 {\rm millimag-day}-1 for R and 1.37 {\rm millimag-day}-1 for V. The R decay rate corresponds to a half life of 315 days, which is much longer than that of 56Co. The spectra show broad hydrogen (1500 km/s FWHM) and oxygen (10000 km/s FWZI) emission features along with many unresolved emission lines. Some of the more interesting narrow lines identified correspond to high ionization states for iron such as Fe~VII and Fe~X which indicate temperatures as high as 106 degrees K. These high ionization states, the X-ray detection by Lewin et al. (1996, IAUC 6445) and the slow photometric decay suggest that SN~1995N is powered by a shock propagating through a dense circumstellar environment. From the earliest observations the energy output appears dominated by the interaction and not by radioactivity, implying that the progenitor exploded well before the discovery of SN~1995N. The situation may be similar to SN~1987A, where the rise in emission from a circumstellar interaction is only now beginning and is expected to peak some 15 years after the supernova explosion.


text:
## Chirp Input Signals

A swept-frequency cosine input signal, or chirp signal, excites your system at a range of frequencies, such that the input frequency changes instantaneously.

You can use a chirp input signal for frequency-response estimation at the command line or in Model Linearizer. The estimation algorithm injects the chirp signal at the input point you specify for estimation, and measures the response at the output point. Chirp signals are useful when your system is nearly linear in the simulation range. They are also useful when you want to obtain a response quickly for a lot of frequency points. The frequency-response model that results when you use a chirp input contains only frequencies that fall within the range of the chirp.

### Create Chirp Signals Using Model Linearizer

In the Model Linearizer, to use a chirp input signal for estimation, in the Estimation tab, select Input Signal > Chirp. You can specify the frequency range and other properties of the chirp in one of two ways:

• Enter the values manually in the Create chirp input dialog box.

• Initialize the frequencies based on the dynamics of a linear model, such as a linearization of your system.

To create a chirp signal based on a linear model:

1. Obtain a linearized model, `linsys1`.

For example, see Linearize Simulink Model at Model Operating Point, which shows how to linearize a model.

2. In the Model Linearizer, in the Estimation tab, select Input Signal > Chirp.

The Create chirp input dialog box opens.

3. In the System list, select `linsys1`. Click .

The software automatically selects frequency points based on the dynamics of `linsys1`. The software also automatically determines other parameters of the chirp signal, including:

• Frequency range — Range of frequencies for the chirp, which the software chooses based on the frequencies at which the linear system has interesting dynamics.

• Amplitude — Amplitude of applied perturbation.

• Sample time — Sample time of signal. To avoid aliasing, the software chooses the sample time such that the Nyquist frequency of the signal is five times the upper end of the frequency range, $\frac{2\pi }{5*\mathrm{max}\left(FreqRange\right)}$.

• Number of samples

• Initial phase

• Sweep method

• Sweep shape

4. Click to create the chirp input signal. A new input signal `in_chirp1` appears in the Linear Analysis Workspace.

You can now select this signal in the Input Signal drop-down list for estimation.

The mapping between the parameters of the Create chirp input dialog box in the Model Linearizer and the properties of `frest.Chirp` is as follows:

Create chirp input dialog box`frest.Chirp` property
Frequency range > From First element associated with the `'FreqRange'` option
Frequency range > ToSecond element associated with the `'FreqRange'` option
Amplitude`'Amplitude'`
Sample time (sec)`'Ts'`
Number of samples`'NumSamples'`
Initial phase (deg)`'InitialPhase'`
Sweep method`'SweepMethod'`
Sweep shape`'Shape'`

### Create Chirp Signals Using MATLAB Code

To create a chirp signal for estimation at the command line with `frestimate`, use `frest.Chirp`. See that page for examples and more information about chirp signal properties.


text:
# Green's function based unparameterised multi-dimensional kernel density and likelihood ratio estimator

This paper introduces a probability density estimator based on Green's function identities. A density model is constructed under the sole assumption that the probability density is differentiable. The method is implemented as a binary likelihood estimator for classification purposes, so issues such as mis-modeling and overtraining are also discussed. The identity behind the density estimator can be interpreted as a real-valued, non-scalar kernel method which is able to reconstruct differentiable density functions.

## Authors

• 4 publications
• 2 publications
• 1 publication
• ### Strong Asymptotic Properties of Kernel Smooth Density and Hazard Function Estimation for Right Censoring NA Data

The paper considers kernel estimation of the density function together w...
01/17/2019 ∙ by Jianhua Shi, et al. ∙ 0

• ### RNADE: The real-valued neural autoregressive density-estimator

We introduce RNADE, a new model for joint density estimation of real-val...
06/02/2013 ∙ by Benigno Uria, et al. ∙ 0

• ### A simple efficient density estimator that enables fast systematic search

This paper introduces a simple and efficient density estimator that enab...
07/03/2017 ∙ by Jonathan R. Wells, et al. ∙ 0

• ### Shannon entropy estimation for linear processes

In this paper, we estimate the Shannon entropy S(f) = -[ log (f(x))] of ...
09/08/2020 ∙ by Timothy Fortune, et al. ∙ 0

• ### Condition Number Analysis of Kernel-based Density Ratio Estimation

The ratio of two probability densities can be used for solving various m...
12/15/2009 ∙ by Takafumi Kanamori, et al. ∙ 0

• ### Trimmed Density Ratio Estimation

Density ratio estimation is a vital tool in both machine learning and st...
03/09/2017 ∙ by Song Liu, et al. ∙ 0

• ### Modeling a Hidden Dynamical System Using Energy Minimization and Kernel Density Estimates

In this paper we develop a kernel density estimation (KDE) approach to m...
04/08/2019 ∙ by Trevor K. Karn, et al. ∙ 0

##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Theory of density estimation

Generally, a density estimator is an identity operator for the probability density function. If it is a linear operator, then it is an approximation of the Dirac-delta function:

 gest(x)=∫δapprox(x−x′)g(x′)dx.

Fixed kernel methods usually estimate the density with a bias, since the kernel does not converge to the delta function. The same can be seen on histograms with fixed bin sizes, and to make them unbiased certain adaptivity is necessary. Different kernels usually maintain different features of the , such as differentiability or easy recognizability of peaks [1]

. Non-adaptive methods are the estimation from the moments via Fourier integrals

[2] or the estimations that target the integral probability density. A problem with the methods using Fourier integrals, is that they may require computationally expensive volume integrals, and the result is not necessarily a density function until sufficient number of moments are known.

To avoid the complexity of finding the global minimum of the error function, most of the overtraining problems and heavy volume integrals, we tried to define a successive adaptive kernel method to estimate the probability densities, using a kernel that is suppressed at large distances. The starting point is, that functionals like linear operators are preferred, as they can be successfully approximated with a sample of a distribution, because they act like expectation values of random variables.

 limN→∞N∑i=1Lp=∫L(x,x′)p(x′)dμ(x′),

where is a linear operator,

is a parameter vector and

is a probability measure. To have such an estimator for the probability density itself, an identity operator is needed:

 μ(x)=I(μ)=P(∫L(x,x′)p(x′)dμ(x′)),

where is a necessary post-processing to convert the linear operator’s output to a scalar. Identity operators can be constructed using Green functions of operators. A particularly interesting one is the Green function of the Laplace operator in dimensions111Although the Green’s function is different for , the same procedure can be repeated and it gives the same form of functions at the end, making method valid in 2 dimensions..

 ΔG(x,x′)=δ(x−x′) (1)
 G(x,x′)=−1Sn1|x−x′|n−2.

Here is a constant, the surface of a unit sphere in dimensions. For twice differentiable functions that disappear at the boundaries, , the identity operator can be constructed:

 c(x)= −1Sn∫RnΔc(x′)|x−x′|n−2dnx′ = 1Sn∫Rn∂μ′1|x−x′|n−2⋅∂μ′c(x′)dnx′.

This can not be used directly for density estimations, since it is linear in instead of , but an identity can be constructed for the vector field of the first derivative:

 ∂μc(x)=1Sn∫Rn∂μ∂μ′1|x−x′|n−2⋅∂μ′c(x′)dnx′.

The amplitude of is still a scalar function; hence for those probability measures, , where a exists222See the subsection targeted at the existence of such a . with , there must be at least one vector field consisting of unit vectors, for which the following identity holds:

 g(x)ϕμ(x)=1Sn∫Rn∂μ∂μ′1|x−x′|n−2⋅ϕμ′(x′)g(x′)dnx′. (2)

Since the integral operator was derived from the Laplace operator, in three dimensions it is similar to the dipole interaction operator in electromagnetism, hence might be called a dipole field. What is important is that this field can be found without the knowledge of the true ; only the possibility of calculating the expectation values on the left-hand side of eq. (2) is required. In this case it is possible to calculate the left-hand side for any configuration of , until the expectation value has the same direction as at every :

 ϕμ(x)∥Eμ(x)=1Sn∫Rn∂μ∂μ′1|x−x′|n−2⋅ϕμ′(x′)g(x′)dnx′. (3)

A way to show that is to do an integration by parts first:

 Eμ(x)=1Sn∫Rn∂μ∂μ′1|x−x′|n−2⋅ϕμ′(x′)g(x′)dnx′=−1Sn∫Rn∂μ1|x−x′|n−2⋅∂μ′(ϕμ′(x′)g(x′))dnx′,

where the integral on the boundaries is zero, since the derivatives of the probability measures are also suppressed at infinity. It is possible to do an external differentiation now, calculating the divergence of , and using eq. (1):

 ∂μEμ(x)=∫Rn−1Sn∂μ∂μ1|x−x′|n−2ΔG(x,x′)=δ(x−x′)⋅∂μ′(ϕμ′(x′)g(x′))⋅dnx′=∂μ(ϕμ(x)g(x)).

For an arbitrary , and can only differ by a divergence-free vector field. In order to show , a certain measure is needed for the difference. With the square of the difference it becomes:

 ∫Rn(g(x)ϕμ(x)−Eμ(x))2dnx≥0. (4)

Expanding the square, it will have three terms:

 ∫Rndnxg2% invariant−2gϕμEμdipole energy+E2field energy.

The field energy can be expressed in terms of the dipole energy. Explicitly writing it out:

 ∫RndnxE2=1S2n∫R2ndnxdny[gxϕxμgyϕyν]∫Rndnz∂xμ∂zη1|x−z|n−2∂yν∂zη1|y−z|n−2Iint.

Performing integration by parts on the term, appears. Due to eq. (1) this can be substituted with a delta function. The remaining part has the derivative of the delta function, which can be easily evaluated:

 Iint.=−Sn∫Rndnz∂xμδ(x−z)∂yν1|y−z|n−2=Sn∂yν∂xμ1|x−y|n−2.

As a result the field energy equals the dipole energy

 ∫RndnxE2=∫RndnxgϕμEμ ,

so the eq. (4) can be modified accordingly:

 ∫Rn(g(x)ϕμ(x)−Eμ(x))2dnx=∫Rndnxg2−gϕμEμ≥0.

The minimum is not necessarily zero, but it can be reached by finding the energy minimum of the dipole system.

### 1.1 Existence of the optimum

It was assumed for eq. (2), that an identity operator exists for a given probability density with a special configuration:

 ∂μc(x)=g(x)ϕμ(x).

In other words, such a has to exist which has a gradient with specified absolute value . This only exists if can be made curl free, and all possible loop integrals should result in zero:

 ∮Lg(x)ϕμ(x)dlμ=0.

As an -dimensional cube has faces, from which share the same corner, linearly independent infinitesimal loops can be imagined. So the curl in -dimensions may have that many independent variables, which has to be zeroed out by a vector field with only independent variables. That seems to be non-trivial, but for differentiable , it turns out it can be done easily. The curl in question is:

 [curl g(x)ϕμ(x)]νη=∂(ν)(gϕμ)e(ν)μ−∂(η)(gϕμ)e(η)μ?=0, (5)

where denotes the unit basis vector in the direction, and the bracketed indices are not summed up. In the basis , the curl is simplified:

 [curl f(x)ϕμ(x)]νη=∂(ν)(gϕ(ν))−∂(η)(gϕ(η))?=0.

As only certain derivatives of appear, one possibility to anull the curl is to set every subcomponent to zero:

 ∂(ν)(gϕ(ν))=0
 ∂(ν)ϕ(ν)ϕ(ν)=−∂(ν)gg.

In this differential equation only the relative change of appears, so setting at every is still possible. As a consequence, the vector field exists for any differentiable , making curl-free, and a function exists with the property , making eq. (2) an identity.

## 2 Density estimation from a finite sample

As long as a suitable field is found that satisfies eq. (3), it can be used to estimate any differentiable probability density function from a sample :

 g(x)=limN→∞∣∣ ∣∣1NSnN∑i=1∂μ∂ν1|x−xi|n−2⋅ϕν(xi)∣∣ ∣∣. (6)

However, since the Green’s function has divergences, a finite sample may produce unnecessarily large fluctuations. Fortunately, volumes around most of the divergences can be neglected, as their integral is zero. A homogeneous , homotropic sphere surface around a point gives no contribution:

 ∫Sn∂μ∂μ′1|x−x′|n−2gϕμ′dnx′ =gϕμ′∂μ∫Sn∂μ′1|x−x′|n−2dnx′ =gϕν∂μ∫Sn−1|x−x′|n−1const.^rμ′dSn∫Sn^r′μdSn=0dr′=0.

As a result, a small sphere of radius around the point of interest can be excluded from the summation in eq. (2), since the error will go as for any differentiable and . The only exception where the error is not is for functions, so these points have to be treated independently.

In order to have control on the fluctuations, one has to choose the radius of the exclusion volume in such a way that a thick shell around this sphere contains enough sample points to estimate the integral belonging to that area. In this case the points in this layer will have a dominant contribution to the sum between , since points at larger distances will have decreasing values. To have the contribution in the shell to be of the same order of magnitude, has to be around . This approximation is derived from the behaviour of the function, which becomes nearly flat in , with a value around . From this it follows that if one needs number of points in this stable layer, approximately points have to be discriminated, which eventually determines . In practice one may choose a different for the iterative search of the parameters and for the evaluation of the density estimation.

### 2.1 Finding the parameters

With a suitable , the set of that can be used for density estimations can be found with iterative search. The stopping criterium can be formulated with discretising eq. (3):

 ϕμi∥Eμ(xi)=1NSn∑j∈{j:|xi−xj|>Ridiscr}∣∣{k:|xi−xk|

As this requirement is similar to looking for the ground state of the dipole system, a way to find it is the gradient descent on the energy function. With this direction is:

 dϕiη=(1−ϕiηϕiμ)Eiμ, (8)

which is a rotation with the field in the directions perpendicular to . A gradient descent method typically contains a line search part, but it is not necessary here. It is enough to set a limit for the angular change of the , an order or magnitude smaller than the maximal possible, , seems to be a good choice, which ensures that the field won’t change significantly during the iteration, keeping it a valid gradient. Further speed increase can be achieved, while taking care not to overshoot the angle between and , since this is what minimises the energy. However, there might be several local energy minima, but the identity in eq. (2) is only true if is curl free. Following from the inequality in eq. (4) it means that the non-global minima solutions are not satisfying the identity and they are not curl free. These solutions will hold a bias in the density estimations, which can be relatively big in the low density regions. This is because the non-optimal minima have a remnant polarisation, for example the destructive interferences do not fully cancel at zero density regions.

As can be seen on fig. 2

, the algorithm actually performs as a density estimator. A sample with a size of two thousand points from a two dimensional Gaussian distribution was fed into the parameter search, and the density estimation was evaluated in the described manner, at the position of sample points. The figure shows the average of the estimated densities for a given radius, where the theoretical value is constant, because of the rotational symmetry of the chosen Gaussian distribution. The spread was calculated from the variance of the estimations on the same radius. The uncertainty can be directly compared with the simplest, kernel-less nearest

-neighbour method [1] in fig. 2. It can be seen, that the uncertainty for the discussed density estimator in fig. 2 is approximately half of that from the nearest -neighbour method with flat kernel in fig. 2. The resolution were chosen to be the same for the two figures, given that the exclusion radii, defined in Section 2, were the same as the radii that was used for volume calculation in the nearest -neighbour method.

## 3 Classification and overtraining

Although binary classification can be done with binary regression models, their results typically depend only on the fraction of the signal and background density , which is optimal due to the Neyman-Pearson Lemma [3]

. This is also true for an optimally trained neural network with a quadratic loss function, where the response is usually the following

[1][4]:

 rNN=s(x)s(x)+b(x).

In case the network is not optimally trained, the error of the response can be related to the uncertainty of the fraction. The same fraction can be calculated from a density estimator, giving a likelihood estimation. The uncertainty in the latter case is much better under control, since the calculation of the signal and background densities are made independently. For the method described in this article, a density estimation for a single point comes from the superposition of the kernels of the most dominant sample points, hence by construction the likelihood estimation can not come from a single point, as happens during overtraining in regressions.

Figure 4 and fig. 4 show the response distribution for a network trained on a signal consisting of twelve non-overlapping Gaussians and a flat background [5]. Unlike in many regression models, here there is no need for an independent test sample for an early stopping criterium. Figure 4 is evaluated on the full training sample, while fig. 4 was made on an independent sample. The shape of the two responses are essentially same, the slight difference is at the part where either the signal or the background distribution density is very low, since these areas have the largest uncertainty in density estimators.

## 4 Conclusions

The article presents a novel class of density estimators, based on Green’s function identities. The investigated one is derived from the Green’s function of the Laplace operator which provides meaningful results for once differentiable probability density functions. It is possible to construct density estimators from other linear operators, but they are expected to require either higher order differentiability or other functional constraints on the investigated probability density. As it was shown, such a density estimator can be used for binary classification. In such case, due to the handling of the fluctuations, overtraining is not so prominent as it is in binary regression models. Therefore there is less need for independent test samples.

## Acknowledgements

We would like to thank Eckhard von Törne for giving encouragement, for pointing out the important topics needed to be touched and for the sample he prepared and gave us to use for testing the classification power of the algorithm.


text:
• A new pixel detector for the CMS experiment was built in order to cope with the instantaneous luminosities anticipated for the Phase~I Upgrade of the LHC. The new CMS pixel detector provides four-hit tracking with a reduced material budget as well as new cooling and powering schemes. A new front-end readout chip mitigates buffering and bandwidth limitations, and allows operation at low comparator thresholds. In this paper, comprehensive test beam studies are presented, which have been conducted to verify the design and to quantify the performance of the new detector assemblies in terms of tracking efficiency and spatial resolution. Under optimal conditions, the tracking efficiency is $99.95\pm0.05\,\%$, while the intrinsic spatial resolutions are $4.80\pm0.25\,\mu \mathrm{m}$ and $7.99\pm0.21\,\mu \mathrm{m}$ along the $100\,\mu \mathrm{m}$ and $150\,\mu \mathrm{m}$ pixel pitch, respectively. The findings are compared to a detailed Monte Carlo simulation of the pixel detector and good agreement is found.
• ### Adaptive orthonormal systems for matrix-valued functions(1511.07834)

Nov. 24, 2015 math.CV
In this paper we consider functions in the Hardy space $\mathbf{H}_2^{p\times q}$ defined in the unit disc of matrix-valued. We show that it is possible, as in the scalar case, to decompose those functions as linear combinations of suitably modified matrix-valued Blaschke product, in an adaptive way. The procedure is based on a generalization to the matrix-valued case of the maximum selection principle which involves not only selections of suitable points in the unit disc but also suitable orthogonal projections. We show that the maximum selection principle gives rise to a convergent algorithm. Finally, we discuss the case of real-valued signals.
• ### The $H^\infty$ functional calculus based on the $S$-spectrum for quaternionic operators and for $n$-tuples of noncommuting operators(1511.07629)

Nov. 24, 2015 math.FA
In this paper we extend the $H^\infty$ functional calculus to quaternionic operators and to $n$-tuples of noncommuting operators using the theory of slice hyperholomorphic functions and the associated functional calculus, called $S$-functional calculus. The $S$-functional calculus has two versions one for quaternionic-valued functions and one for Clifford algebra-valued functions and can be considered the Riesz-Dunford functional calculus based on slice hyperholomorphicity because it shares with it the most important properties. The $S$-functional calculus is based on the notion of $S$-spectrum which, in the case of quaternionic normal operators on a Hilbert space, is also the notion of spectrum that appears in the quaternionic spectral theorem. The main purpose of this paper is to construct the $H^\infty$ functional calculus based on the notion of $S$-spectrum for both quaternionic operators and for $n$-tuples of noncommuting operators. We remark that the $H^\infty$ functional calculus for $(n+1)$-tuples of operators applies, in particular, to the Dirac operator.
• ### The mathematics of superoscillations(1511.01938)

Nov. 5, 2015 quant-ph, math-ph, math.MP
In the past 50 years, quantum physicists have discovered, and experimentally demonstrated, a phenomenon which they termed superoscillations. Aharonov and his collaborators showed that superoscillations naturally arise when dealing with weak values, a notion that provides a fundamentally different way to regard measurements in quantum physics. From a mathematical point of view, superoscillating functions are a superposition of small Fourier components with a bounded Fourier spectrum, which result, when appropriately summed, in a shift that can be arbitrarily large, and well outside the spectrum. Purpose of this work is twofold: on one hand we provide a self-contained survey of the existing literature, in order to offer a systematic mathematical approach to superoscillations; on the other hand, we obtain some new and unexpected results, by showing that superoscillating sequences can be seen of as solutions to a large class of convolution equations and can therefore be treated within the theory of Analytically Uniform spaces. In particular, we will also discuss the persistence of the superoscillatory behavior when superoscillating sequences are taken as initial values of the Schr\"odinger equation and other equations.
• The degradation of signal in silicon sensors is studied under conditions expected at the CERN High-Luminosity LHC. 200 $\mu$m thick n-type silicon sensors are irradiated with protons of different energies to fluences of up to $3 \cdot 10^{15}$ neq/cm$^2$. Pulsed red laser light with a wavelength of 672 nm is used to generate electron-hole pairs in the sensors. The induced signals are used to determine the charge collection efficiencies separately for electrons and holes drifting through the sensor. The effective trapping rates are extracted by comparing the results to simulation. The electric field is simulated using Synopsys device simulation assuming two effective defects. The generation and drift of charge carriers are simulated in an independent simulation based on PixelAV. The effective trapping rates are determined from the measured charge collection efficiencies and the simulated and measured time-resolved current pulses are compared. The effective trapping rates determined for both electrons and holes are about 50% smaller than those obtained using standard extrapolations of studies at low fluences and suggests an improved tracker performance over initial expectations.
• ### Quantum harmonic oscillator with superoscillating initial datum(1411.4112)

Nov. 15, 2014 quant-ph, math-ph, math.MP
In this paper we study the evolution of superoscillating initial data for the quantum driven harmonic oscillator. Our main result shows that superoscillations are amplified by the harmonic potential and that the analytic solution develops a singularity in finite time. We also show that for a large class of solutions of the Schr\"odinger equation, superoscillating behavior at any given time implies superoscillating behavior at any other time.
• ### An extension of Herglotz's theorem to the quaternions(1403.0079)

July 23, 2014 math.FA
A classical theorem of Herglotz states that a function $n\mapsto r(n)$ from $\mathbb Z$ into $\mathbb C^{s\times s}$ is positive definite if and only there exists a $\mathbb C^{s\times s}$-valued positive measure $d\mu$ on $[0,2\pi]$ such that $r(n)=\int_0^{2\pi}e^{int}d\mu(t)$for $n\in \mathbb Z$. We prove a quaternionic analogue of this result when the function is allowed to have a number of negative squares. A key tool in the argument is the theory of slice hyperholomorphic functions, and the representation of such functions which have a positive real part in the unit ball of the quaternions. We study in great detail the case of positive definite functions.
• ### The quantum pigeonhole principle and the nature of quantum correlations(1407.3194)

July 11, 2014 quant-ph
The pigeonhole principle: "If you put three pigeons in two pigeonholes at least two of the pigeons end up in the same hole" is an obvious yet fundamental principle of Nature as it captures the very essence of counting. Here however we show that in quantum mechanics this is not true! We find instances when three quantum particles are put in two boxes, yet no two particles are in the same box. Furthermore, we show that the above "quantum pigeonhole principle" is only one of a host of related quantum effects, and points to a very interesting structure of quantum mechanics that was hitherto unnoticed. Our results shed new light on the very notions of separability and correlations in quantum mechanics and on the nature of interactions. It also presents a new role for entanglement, complementary to the usual one. Finally, interferometric experiments that illustrate our effects are proposed.
• ### Boundary interpolation for slice hyperholomorphic Schur functions(1404.3352)

May 25, 2014 math.FA, math.CV
A boundary Nevanlinna-Pick interpolation problem is posed and solved in the quaternionic setting. Given nonnegative real numbers $\kappa_1, \ldots, \kappa_N$, quaternions $p_1, \ldots, p_N$ all of modulus $1$, so that the $2$-spheres determined by each point do not intersect and $p_u \neq 1$ for $u = 1,\ldots, N$, and quaternions $s_1, \ldots, s_N$, we wish to find a slice hyperholomorphic Schur function $s$ so that $$\lim_{\substack{r\rightarrow 1\\ r\in(0,1)}} s(r p_u) = s_u\quad {\rm for} \quad u=1,\ldots, N,$$ and $$\lim_{\substack{r\rightarrow 1\\ r\in(0,1)}}\frac{1-s(rp_u)\overline{s_u}}{1-r}\le\kappa_u,\quad {\rm for} \quad u=1,\ldots, N.$$ Our arguments relies on the theory of slice hyperholomorphic functions and reproducing kernel Hilbert spaces.
• ### Bloch, Besov and Dirichlet spaces of slice hyperholomorphic functions(1403.7187)

March 27, 2014 math.CV
In this paper we begin the study of some important Banach spaces of slice hyperholomorphic functions, namely the Bloch, Besov and weighted Bergman spaces, and we also consider the Dirichlet space, which is a Hilbert space. The importance of these spaces is well known, and thus their study in the framework of slice hyperholomorphic functions is relevant, especially in view of the fact that this class of functions has recently found several applications in operator theory and in Schur analysis. We also discuss the property of invariance of this function spaces with respect to Moebius maps by using a suitable notion of composition.
• ### The spectral theorem for unitary operators based on the $S$-spectrum(1403.0175)

March 2, 2014 math.FA, math.SP
The quaternionic spectral theorem has already been considered in the literature, see e.g. [22], [31], [32], however, except for the finite dimensional case in which the notion of spectrum is associated to an eigenvalue problem, see [21], it is not specified which notion of spectrum underlies the theorem. In this paper we prove the quaternionic spectral theorem for unitary operators using the $S$-spectrum. In the case of quaternionic matrices, the $S$-spectrum coincides with the right-spectrum and so our result recovers the well known theorem for matrices. The notion of $S$-spectrum is relatively new, see [17], and has been used for quaternionic linear operators, as well as for $n$-tuples of not necessarily commuting operators, to define and study a noncommutative versions of the Riesz-Dunford functional calculus. The main tools to prove the spectral theorem for unitary operators are the quaternionic version of Herglotz's theorem, which relies on the new notion of $q$-positive measure, and quaternionic spectral measures, which are related to the quaternionic Riesz projectors defined by means of the $S$-resolvent operator and the $S$-spectrum. The results in this paper restore the analogy with the complex case in which the classical notion of spectrum appears in the Riesz-Dunford functional calculus as well as in the spectral theorem.
• ### Planck pre-launch status: calibration of the Low Frequency Instrument flight model radiometers(1005.2541)

May 14, 2010 astro-ph.IM
The Low Frequency Instrument (LFI) on-board the ESA Planck satellite carries eleven radiometer subsystems, called Radiometer Chain Assemblies (RCAs), each composed of a pair of pseudo-correlation receivers. We describe the on-ground calibration campaign performed to qualify the flight model RCAs and to measure their pre-launch performances. Each RCA was calibrated in a dedicated flight-like cryogenic environment with the radiometer front-end cooled to 20K and the back-end at 300K, and with an external input load cooled to 4K. A matched load simulating a blackbody at different temperatures was placed in front of the sky horn to derive basic radiometer properties such as noise temperature, gain, and noise performance, e.g. 1/f noise. The spectral response of each detector was measured as was their susceptibility to thermal variation. All eleven LFI RCAs were calibrated. Instrumental parameters measured in these tests, such as noise temperature, bandwidth, radiometer isolation, and linearity, provide essential inputs to the Planck-LFI data analysis.
• ### Dynamic validation of the Planck/LFI thermal model(1001.4646)

Jan. 26, 2010 astro-ph.CO, astro-ph.IM
The Low Frequency Instrument (LFI) is an array of cryogenically cooled radiometers on board the Planck satellite, designed to measure the temperature and polarization anisotropies of the cosmic microwave backgrond (CMB) at 30, 44 and 70 GHz. The thermal requirements of the LFI, and in particular the stringent limits to acceptable thermal fluctuations in the 20 K focal plane, are a critical element to achieve the instrument scientific performance. Thermal tests were carried out as part of the on-ground calibration campaign at various stages of instrument integration. In this paper we describe the results and analysis of the tests on the LFI flight model (FM) performed at Thales Laboratories in Milan (Italy) during 2006, with the purpose of experimentally sampling the thermal transfer functions and consequently validating the numerical thermal model describing the dynamic response of the LFI focal plane. This model has been used extensively to assess the ability of LFI to achieve its scientific goals: its validation is therefore extremely important in the context of the Planck mission. Our analysis shows that the measured thermal properties of the instrument show a thermal damping level better than predicted, therefore further reducing the expected systematic effect induced in the LFI maps. We then propose an explanation of the increased damping in terms of non-ideal thermal contacts.
• ### Planck pre-launch status: Low Frequency Instrument calibration and expected scientific performance(1001.4562)

Jan. 25, 2010 astro-ph.IM
We give the calibration and scientific performance parameters of the Planck Low Frequency Instrument (LFI) measured during the ground cryogenic test campaign. These parameters characterise the instrument response and constitute our best pre-launch knowledge of the LFI scientific performance. The LFI shows excellent $1/f$ stability and rejection of instrumental systematic effects; measured noise performance shows that LFI is the most sensitive instrument of its kind. The set of measured calibration parameters will be updated during flight operations through the end of the mission.
• ### Planck pre-launch status: Design and description of the Low Frequency Instrument(1001.3321)

Jan. 19, 2010 astro-ph.CO, astro-ph.IM
In this paper we present the Low Frequency Instrument (LFI), designed and developed as part of the Planck space mission, the ESA program dedicated to precision imaging of the cosmic microwave background (CMB). Planck-LFI will observe the full sky in intensity and polarisation in three frequency bands centred at 30, 44 and 70 GHz, while higher frequencies (100-850 GHz) will be covered by the HFI instrument. The LFI is an array of microwave radiometers based on state-of-the-art Indium Phosphide cryogenic HEMT amplifiers implemented in a differential system using blackbody loads as reference signals. The front-end is cooled to 20K for optimal sensitivity and the reference loads are cooled to 4K to minimise low frequency noise. We provide an overview of the LFI, discuss the leading scientific requirements and describe the design solutions adopted for the various hardware subsystems. The main drivers of the radiometric, optical and thermal design are discussed, including the stringent requirements on sensitivity, stability, and rejection of systematic effects. Further details on the key instrument units and the results of ground calibration are provided in a set of companion papers.
• ### Non commutative functional calculus: unbounded operators(0708.3592)

Jan. 25, 2008 math.SP
In a recent work, \cite{cgss}, we developed a functional calculus for bounded operators defined on quaternionic Banach spaces. In this paper we show how the results from \cite{cgss} can be extended to the unbounded case, and we highlight the crucial differences between the two cases. In particular, we deduce a new eigenvalue equation, suitable for the construction of a functional calculus for operators whose spectrum is not necessarily real.
• ### A new functional calculus for non-commuting operators(0708.3594)

Jan. 25, 2008 math.SP
In this paper we use the notion of slice monogenic functions \cite{slicecss} to define a new functional calculus for an $n$-tuple $T$ of not necessarily commuting operators. This calculus is different from the one discussed in \cite{jefferies} and it allows the explicit construction of the eigenvalue equation for the $n$-tuple $T$ based on a new notion of spectrum for $T$. Our functional calculus is consistent with the Riesz-Dunford calculus in the case of a single operator.
• ### Slice monogenic functions(0708.3595)

Jan. 25, 2008 math.CV
In this paper we offer a definition of monogenicity for functions defined on $\rr^{n+1}$ with values in the Clifford algebra $\rr_n$ following an idea inspired by the recent papers \cite{gs}, \cite{advances}. This new class of monogenic functions contains the polynomials (and, more in general, power series) with coefficients in the Clifford algebra $\rr_n$. We will prove a Cauchy integral formula as well as some of its consequences. Finally, we deal with the zeroes of some polynomials and power series.
• ### Non commutative functional calculus: bounded operators(0708.3591)

Aug. 27, 2007 math.SP
In this paper we develop a functional calculus for bounded operators defined on quaternionic Banach spaces. This calculus is based on the notion of slice-regularity, see \cite{gs}, and the key tools are a new resolvent operator and a new eigenvalue problem.


text:
# How to calculate “Kronecker Product” of a vector and a matrix with numpy

I need to compute the following forumula:

Fromula in TeX:

$\sum_n^N \sum_m^N a_n * a_m * C_{nm}$

Peudocode:

a = array of length N
C = NxN matrix
retval = 0
for n in range(N):
for m in range(N):
retval += a[n] * a[m] * C[n][m]

If a were a NxN matrix constructed as in the product above one could simply use np.kron for the Kronecker Matrix multiplication and then use np.sum to get the desired result. However I don't know a faster numpy way of constructing a matrix A as in the formula above.

Any ideas?

We could directly use those iterators for np.einsum string-notation -
retval = np.einsum('n,m,nm->',a,a,C)

Alternatively, with np.dot -
retval = a.dot(C).dot(a)

• Wow thanks for the quick answer. Do you know which of the two is faster? Probably np.dot right? – jaaq Jun 23 at 16:59
• Times are: Naive for loops: 6m35s dot method: 4m10s einsum: 4m4s – jaaq Jun 23 at 17:28


text:
# Counterexample of the day

For a little while we thought that there was an extension to Earnshaw’s theorem. But there isn’t: it is quite feasible to have fully confining rf traps with nonzero field at the center.


text:
Numbers as decimals, fractions, percentages

Presenting a number as a fraction, percent or decimal number is quite a simple procedure.
You do not have to do complicated calculations or use difficult patterns.

If you want to see how to do it, enter any number in the field below.
We'll show you step by step solution.

Enter number to get the solution

 Write as a decimal as a fraction as a percent


text:
# Swans

Algebra Level pending

Out of a group of swans, $$\frac { 7 }{ 2 }$$ of the the square root of the number of swans are playing on the shore of a tank. The 2 remaining ones are playing,with amorous fight, in the water. What is the total no. of swans?

×


text:
# Assignment 4 [Coding]: Conformal Parameterization (Due 4/28)

For the coding portion of your assignment on conformal parameterization, you will implement the Spectral Conformal Parameterization (SCP) algorithm as described in the course notes.Please implement the following routines in

1. core/geometry.js:
1. complexLaplaceMatrix
2. projects/parameterization/spectral-conformal-parameterization.js:
1. buildConformalEnergy
2. flatten
3. utils/solvers.js:
1. solveInversePowerMethod
2. residual

Notes

• I added a few tests for utils/solvers.js. They’re available in the git repository. If you’re familiar with git, you can pull the changes. If you’d prefer not to mess around with git, you can also just download the files and put them in tests/utils.
• Warning: I fixed a problem with the residual test. In the notes, the residual is defined as a vector $Ax – \lambda x$. But your function residual should return a float. You should return $\frac{\|Ax – \lambda x\|_2}{\|x\|_2}$. Furthermore, you should compute $\lambda$ as $\frac{x^* A x}{x^* x}$ where $x^*$ is the conjugate transpose of $x$.
• The complex version of the cotan-Laplace matrix can be built in exactly the same manner as its real counterpart. The only difference now is that the cotan values of our matrix will be complex numbers with a zero imaginary component. This time, we will work with meshes with boundary, so your Laplace matrix will have to handle boundaries properly (you just have to make sure your cotan function returns 0 for halfedges which are in the boundary).
• The buildConformalEnergy function builds a $|V| \times |V|$ complex matrix corresponding to the conformal energy $E_C(z) = E_D(z) – \mathcal A(Z)$ where $E_D$ is the Dirichlet energy (given by our complex cotan-Laplace matrix) and $\mathcal A$ is the total signed area of the region $z(M)$ in the complex plane (Please refer to Section 7.4 of the notes for more details). You may find it easiest to iterate over the halfedges of the mesh boundaries to construct the area matrix (Recall that the Mesh object has a boundaries variable which stores all the boundary loops.
• The flatten function returns a dictionary mapping each vertex to a vector (linear-algebra/vector.js) of planar coordinates by finding the eigenvector corresponding to the smallest eigenvalue of the conformal energy matrix. You can compute this eigenvector by using solveInversePowerMethod (described below).
• Your solveInversePowerMethod function should implement Algorithm 1 in Section 7.5 of the course notes with one modification – after computing $Ay_i = y_{i-1}$, center $y_i$ around the origin by subtracting its mean. You don’t have to worry about the $B$ matrix or generalized eigenvalue problem. Important: Terminate your iterations when your residual is smaller that 10^-10.
• The parameterization project directory also contains a basic implementation of the Boundary First Flattening (BFF) algorithm. Feel free to play around with it in the viewer and compare the results to your SCP implementation.

Submission Instructions

Please rename your geometry.js, spectral-conformal-parameterization.js and solvers.js files to geometry.txt, spectral-conformal-parameterization.txt and solvers.txt (respectively) and put them in a single zip file called solution.zip. These files and your solution to the written exercises should be submitted together in a single email to Geometry.Collective@gmail.com with the subject line DDG20A4.

## 17 thoughts on “Assignment 4 [Coding]: Conformal Parameterization (Due 4/28)”

1. yixinh says:

I’m confused by the calculation of $\lambda$ . Should the equation, $\frac{x^\dagger A x}{x^\dagger x}$ be $\frac{x^\dagger A x}{\|x^\dagger x\|}$? If that’s not a typo, how do you do division by a dense matrix?

1. Keenan says:

The quantity $x^T x$ is a scalar (or if you want to be really pedantic: it’s a 1×1 matrix), given by
$x^T x = \sum_{i=1}^n x_i^2.$
Likewise, $x^T A x$ is a scalar: it’s evaluation of a quadratic form, equal to
$x^T A x = \sum_{i=1}^n \sum_{j=1}^n A_{ij} x_i x_j$
(though for a sparse matrix $A$, you do not need to sum over $O(n^2)$ terms and can instead just do a sparse matrix-vector multiply $Ax$, followed by a dense vector-vector inner product $x^T(Ax)$).

Note that the symbol $^T$ denotes the matrix transpose, whereas $^\dagger$ typically denotes the pseudoinverse (a very different operation!).

2. hawaiii says:

Github shows the two files in test/util were last updated 13 months ago. And git pull told me “Already up to date.” Am I missing something?

1. Alex says:

Sounds like you’re up to date. If you have tests for the solver then you’re fine.

2. yixinh says:

It’s because the write up is copied from last year’s iteration of the class haha

3. yixinh says:

Possibly a silly question: my flatten timeout’s if I run the test from Chrome but works fine and passes if I run it in Firefox. Is that okay? Or do I need to make it faster somehow?

My solveInversePowerMethod() stops after 7 iterations, if that helps.

1. Keenan says:

If it’s *really* slow in Chrome (like, you sit there waiting for a long time) then there could be something wrong. If it’s just a little slower than in Firefox (like, you don’t really notice when you’re just sitting there watching it), then you’re probably fine.

4. Wei Dong says:

Question irrelevent to this assignment: up to now I have only received the feedback from Assignment 0… is the grading still ongoing or have I missed something? Thanks!

5. bahn says:

Question regarding some details:

1. In this homework, the Laplace-Beltrami operator is defined as -*d*d (as in Exercise 7.6), which is a negative of what we have defined previously. Is it just for this homework or can we use it interchangeably?

2. Also, I am a bit confused about when we use the mass matrix. Why can we ignore the mass matrix and use just the weak Laplacian matrix in this homework? Is our solution an approximation of the exact conformal parameterization or is it just a different normalization?

1. Daniel Li says:

1. I think the one in the homework is a sign error. My solution has would require it to be consistent with what we have before. Does yours match?

2. This is more so a matter of dimension. When we integrate the Laplacian on our mesh, we are integrating over 2-cells, so the result is a 2-form. However, applying the Laplacian to a function should return a scalar function (0-form). Thus, we must apply the dual 2-form Hodge star, which amounts to dividing by area (the Mass Matrix). In the math calculations, you will see this by the inverse mass matrices we multiply. However, in the programming, the Laplace Matrix you build is actually the zero form already. This means that we basically have $LM^{-1}$ when we build the matrix. Thus, when we multiply by M, we are just left with $L$ as desired. For example, in HW 3, for buildFlowOperator, this is precisely why you constructed $M+h\Delta$ and not $M+h\Delta M$.

1. Lqq88888 says:

I don’t quite understand. For heat equations, we have $\frac{du}{dt} = \Delta u$ and we discrete it as $\frac{u_{k+1}-u_k}{h} = Lu_{k+1}$, $(I-hL)u_{k+1} = u_k$. Since $L = M^{-1}C$, we get $(M-hC)u_{k+1} = Mu_k$ I thought the Lapacian we built is -C and shouldn’ t it be a 2-form? And if we want to convert it into 0-form, we should multiply it by $M^{-1}$?

1. Daniel Li says:

Yes, that’s what I was trying to say. Sorry if I wasn’t clear about it. I was just trying to explain why sometimes equations may “appear” incorrect.

1. bahn says:

Thanks for the answer but I don’t still understand question 2 that is about the mass matrix. As you mentioned, we need to multiply the mass matrix $M$ as well as the weak Laplacian matrix $C$ to make $L = MC$ as desired for Dirichlet energy $E_D$ (which I think would require the generalized eigenvalue problem in the end), but we didn’t multiply $M$ in this HW. Is it implicitly handled somewhere along the line or did I miss something?

6. hawner says:

I’m confused about the signed area of the region z(M), what is the value of z at each vertex? I can’t find any complex values attached to the vertices, so I’m not sure what exactly I’m supposed to be using for $$\bar{z}_i z_j – \bar{z}_j z_i$$

1. Daniel Li says:

What we are doing here is to build a constant matrix A such that $$x^TAx$$ is equal to the quadratic form. If you look at exercise 7.14, it hints at why this is the correct method of discretization. The matrix only contains information for the terms of the map z that we want to compute, the boundary edges. This amounts to placing $$\pm \frac{i}{4}$$ in entry ij, if $$e_{ij}$$ is on the boundary.


text:
Astroparticle Physics - A Joint TeVPA/ IDM Conference

Jun 23 – 28, 2014
Amsterdam
Europe/Zurich timezone

Systematic uncertainties in dark matter searches due to halo asphericity

Jun 23, 2014, 6:15 PM
15m
Room 6 (Tuschinski Theatre)

Room 6

Tuschinski Theatre

Presentation Dark Matter Indirect Detection

Nicolas Bernal

Description

We study the impact of aspherical dark matter density distribution in Milky-Way like halo on direct and indirect searches. Using data from large N-body cosmological simulation Bolshoi, we perform a complete statistical analysis and quantify the systematic uncertainties that affect the determination of local dark matter density and $J$ factors for annihilating and decaying dark matter. We find that the systematic effects can be as large as 50% at 95% CL for a spherically averaged local density of 0.3 GeV/cm^3. Similarly, systematic uncertainties for $J$ factors evaluated around the galactic center can be as large as 20% and 15% at 95% CL for decaying and annihilating dark matter, for typical NFW profiles.

 Slides


text:
#### Vol. 11, No. 2, 2020

 Recent Issues Volume 11, Issue 2 Volume 11, Issue 1
 The Journal About the Journal Editorial Board Subscriptions Submission Guidelines Submission Form Policies for Authors Ethics Statement ISSN (electronic): 2693-3004 ISSN (print): 2693-2997 To Appear Other MSP Journals
The norm of the saturation of a binomial ideal, with applications to Markov bases

### David Holmes

Vol. 11 (2020), No. 2, 169–187
##### Abstract

Let $B$ be a finite set of pure binomials in the variables ${x}_{i}$, and write ${I}_{B}$ for the ideal generated by these binomials. We define a new measure of the complexity of the saturation of the ideal ${I}_{B}$ with respect to the product of the ${x}_{i}$, which we call the norm of $B$. We give a bound on the norm in terms of easily computed invariants of $B$. We discuss statistical applications, both practical and theoretical.

##### Keywords
markov basis, saturation, toric ideals
##### Mathematical Subject Classification 2010
Primary: 13P25, 14M25


text:
# Deducing Changes in Inventories Ferguson Products Inc., a manufacturer, reported $130 million in... 1 answer below » Deducing Changes in Inventories Ferguson Products Inc., a manufacturer, reported$130 million in sales and a loss of $25 million in its absorption costing income statement provided to shareholders. According to a CVP analysis prepared for management, the company’s break-even point is$120 million in sales.

Required:

Assuming that the CVP analysis is correct, is it likely that the company’s inventory level increased, decreased, or remained unchanged during the year? Explain.

As per the question the brekeven sales is at $120 million and the sales were$130 million above the breakeven sales, than also there was loss which means there has been...


text:
# Math Help - The half life of a radioactive substance.

1. ## The half life of a radioactive substance.

I understand part a) as the equation is worked so that it is in terms of k. The solution for part b) is losing me at:

$\lambda=-\frac{ln(2)}{k}*y_0e^k^(^t^1^+ ^\lambda)$

I see that where t was before there is now the symbol lambda but the part after the * symbol is lost on me right now. Why is the equation in terms of k now being multiplied by:

$y_0e^k^(^t^1^+ ^\lambda)$

Thanks for any responses...

2. ## Re: The half life of a radioactive substance.

In part b), it is supposed to be a period, not a multiplication. They are using $\lambda$ as the value of $t$ you found in part a). So, if $y_1 = y_0 e^{kt_1}$, then $y_0 e^{k(t_1+\lambda)} = y_0e^{kt_1}e^{k\lambda}$. Now, replace $y_0e^{kt_1}$ with $y_1$ since they are equal by assumption. So, $y_0e^{kt_1}e^{k\lambda} = y_1e^{k\lambda}$. Next, substitute in the value set for $\lambda$: $y_1e^{k\lambda} = y_1e^{-\ln 2} = \dfrac{1}{2}y_1$.

3. ## Re: The half life of a radioactive substance.

Hey sepoto.

Remember that you are looking at the time t1 + lambda and not t1 (and you also have the information about the half life at t1 and t1 + lambda from the question). Just plug in the information for the lambda + t1 and everything should be ok.

4. ## Re: The half life of a radioactive substance.

So:

$t=\frac{-ln(2)}{k}$

So t might be some number 1590 in years. So if lambda is 1590 then why this equation below. So then "t_1+\lambda" is some time t plus lambda where lambda is the half life of 1590?

$y_0e^{k(t_1+\lambda)}$

5. ## Re: The half life of a radioactive substance.

I'm not exactly sure what you are trying to say: can you please ask again but in more detail (and clarity)?


text:
# Finding expansion coefficient of a 3-d Gaussian wave packet

• I

## Main Question or Discussion Point

I'm having trouble with trying to find the expansion coefficients of a superposition of a Gaussian wave packet.

First I'm decomposing a Gaussian wave packet
$$\psi(\textbf{r},0) = \frac{1}{(2\pi)^{3/4}\sigma^{3/2}}\text{exp}\left[ -\frac{(\textbf{r} - \textbf{r}_0)^2}{4\sigma^2} + i\textbf{k}_0 \cdot \textbf{r}\right]$$
into a linear superposition of eigenfunctions of the hydrogen atom, where the initial position and momentum is chosen as $\textbf{r}_0$ and $\textbf{p}_0 = \hbar \textbf{k}_0$, and also a spatial width $\sigma$ such that the initial wave packet $\psi(\textbf{r},0)$ can to a sufficient approximation be superimposed by bound-state eigenfunctions $\varphi_\textit{nlm}(\textbf{r})$,
$$\psi(\textbf{r},0) = \sum_{n=1}^{\infty}\sum_{l=0}^{n-1}\sum_{m=-l}^{l} b_{\textit{nlm}}\varphi_\textit{nlm}(\textbf{r})$$
Having solved the Schrodinger Equation of a hydrogen atom with a Hamiltonian containing the Coulomb potential, the eigenfunctions were obtained to be a product of the Radial Wave solution and Spherical Harmonics
$$\varphi_\textit{nlm}(\textbf{r}) = R_{nl}(r) Y_{lm}(\theta, \phi)$$
where
$$R_{nl}(r) = \frac{1}{a^{3/2}}\frac{2}{n^2}\sqrt{\frac{(n-l-1)!}{(n+1)!}}\left(\frac{2r}{na}\right)^l e^{-\frac{r}{na}} L^{2l+1}_{n-l-1}\left(\frac{2r}{na}\right)$$
where $L^{2l+1}_{n-l-1}\left(\frac{2r}{na}\right)$ is the Laguerre Polynomial which has the form of
$$L^{k}_{p}\left(x\right) = \sum_{k}^{p} (-1)^s \binom {p+k}{p-s} \frac{x^s}{s!}$$
and
$$Y_{lm}(\theta, \phi) = (-1)^m \sqrt{\frac{2l+1}{4\pi}\frac{(l-m)!}{(l+m)!}} P^{m}_{l} (\cos\theta) e^{im\phi}$$
with the Legendre Polynomial contained within the Associated Legendre Functions
$$P^{m}_{l} (x) = (1 - x^2)^{\frac{m}{2}} \frac{d^m}{dx^m}(P_{l}(x))$$
$$P_{l}(x) = \frac{1}{2^l l! }\frac{d^l}{dx^l}[(x^2 - 1)^l]$$

What I ultimately would like to find is the closed form expansion coefficients $b_{nlm}$ of the eigenfunction superposition, now from my knowledge of finding the expansion coefficients of a single-indexed superposition of a wavefunction, the same method can be applied on this superposition which involves 3 different indices, such that
$$\begin{equation*} \begin{split} \langle\varphi_\textit{nlm}(\textbf{r})|\psi(\textbf{r},0)\rangle &= \langle\varphi_\textit{nlm}(\textbf{r})|\sum_{n=1}^{\infty}\sum_{l=0}^{n-1}\sum_{m=-l}^{l} b_{\textit{nlm}}\varphi_\textit{nlm}(\textbf{r})\rangle\\ &= b_{\textit{nlm}}\sum_{n=1}^{\infty}\sum_{l=0}^{n-1}\sum_{m=-l}^{l}\langle\varphi_\textit{n'l'm'}(\textbf{r})| \varphi_\textit{nlm}(\textbf{r})\rangle\\ &= b_{\textit{nlm}}\sum_{n=1}^{\infty}\sum_{l=0}^{n-1}\sum_{m=-l}^{l}\delta_{n'n}\delta_{l'l}\delta_{m'm}\\ &= b_{nlm} \end{split} \end{equation*}$$

My attempt:
To avoid doing the scalar product directly and having to compute the whole integral, I was thinking of using the decomposition of stationary harmonic plane wave into partial waves, as such
$$e^{\textbf{k}\cdot\textbf{r}} = e^{ikz} = e^{ikr\cos\theta} = \sum^{l=0}_{\infty}(2l+1)i^l j_l(kr) P_l(\cos\theta)$$
where $i$ is the imaginary number, $\textbf{k}$ and $\textbf{r}$ are the wave and position vectors, $j_l$ are the spherical Bessel functions, $P_l$ are Legendre Polynomials. But according to Wikipedia (https://en.wikipedia.org/wiki/Plane_wave_expansion), it can also be rewritten as
$$e^{\textbf{k}\cdot\textbf{r}} = 4\pi \sum_{l=0}^{\infty}\sum_{m=-l}^{l} i^l j_l(kr)Y^m_l(\hat{\textbf{k}})Y^{m*}_l(\hat{\textbf{r}})$$

In theory, this should now mean that I can rewrite the scalar product in terms of these sum of spherical waves, so that

$$\begin{equation*} \begin{split} b_{nlm} &= \langle\varphi_\textit{nlm}(\textbf{r})|\psi(\textbf{r},0)\rangle\\ &= \langle R_{nl}(r) Y_{lm}(\theta, \phi)|\frac{1}{(2\pi)^{3/4}\sigma^{3/2}}\text{exp}\left[ -\frac{(\textbf{r} - \textbf{r}_0)^2}{4\sigma^2}\right]\cdot 4\pi\sum_{l=0}^{\infty}\sum_{m=-l}^{l} i^l j_l(kr)Y^m_l(\hat{\textbf{k}})Y^{m*}_l(\hat{\textbf{r}})\rangle\\ \end{split} \end{equation*}$$
and so this is where I'm at. I'm not sure how to go from here. I hope I've clarified my question well.

## Answers and Replies

Related Quantum Physics News on Phys.org
hilbert2
Science Advisor
Gold Member
Is this related to "Kepler wavepackets", where an electron wavepacket orbits a nucleus like in classical central-force motion?

Is this related to "Kepler wavepackets", where an electron wavepacket orbits a nucleus like in classical central-force motion?
Yes how did you know?Ultimately I'm trying to investigate the correspondence principle of the hydrogen atom, so in the limit of high quantum number how it'd would collapse from a quantum mechanical regime to a classical one.

hilbert2
Science Advisor
Gold Member
Someone I worked with many years ago made his master's thesis about that subject.

There seems to be some related material that can be found with a Google search, but I'm not sure if it's useful to you. I suggest you represent the Gaussian wavepacket in spherical polar coordinates, as the hydrogenic eigenstates are, and then try to calculate the expansion coefficients by integration with Wolfram Mathematica. If that application can't solve the integrals exactly, it's very unlikely that they could be found with pen and paper, either.

Someone I worked with many years ago made his master's thesis about that subject.

There seems to be some related material that can be found with a Google search, but I'm not sure if it's useful to you. I suggest you represent the Gaussian wavepacket in spherical polar coordinates, as the hydrogenic eigenstates are, and then try to calculate the expansion coefficients by integration with Wolfram Mathematica. If that application can't solve the integrals exactly, it's very unlikely that they could be found with pen and paper, either.
My supervisor did his master's thesis on it too but he did it on a Quantum Harmonic Oscillator. Most of the stuff can be found from the book "The Picture Book of Quantum Mechanics" but they didn't show how to find the expansion coefficients. This project is waaay above my knowledge of QM which is why I'm always stuck at every step of progress that I'm trying to make lol.

hilbert2
Science Advisor
Gold Member
Based on what I can find about it problem, it doesn't seem to be possible to find a closed form solution for the overlap integral of a gaussian wavepacket and the hydrogenic orbitals. It may be possible for a wavepacket in a 3d harmonic oscillator potential $V(x,y,z)=\frac{1}{2}k(x^2 + y^2 + z^2 )$.

Based on what I can find about it problem, it doesn't seem to be possible to find a closed form solution for the overlap integral of a gaussian wavepacket and the hydrogenic orbitals. It may be possible for a wavepacket in a 3d harmonic oscillator potential $V(x,y,z)=\frac{1}{2}k(x^2 + y^2 + z^2 )$.
It seems that the Legendre functions are related to the Spherical Harmonics by $$P_l(\cos\theta) = \sqrt{\frac{4\pi}{2l+1}} Y_{l,0}(\theta, \phi)$$ when m = 0, so the decomposition of stationary harmonic plane wave into partial waves can now be expressed as $$e^{\textbf{k}\cdot\textbf{r}} = e^{ikz} = e^{ikr\cos\theta} = \sum^{\infty}_{l=0}(2l+1)i^l j_l(kr) P_l(\cos\theta) = \sqrt{\frac{4\pi}{2l+1}}\sum^{\infty}_{l=0}(2l+1)i^l j_l(kr) Y_{l,0}(\theta, \phi)$$
which I can then incorporate into my Gaussian wave packet, but how can I go about changing the $(\textbf{r} - \textbf{r}_0)^2$ into spherical polar coordinates? I'm only trying this approach because my supervisor said to "expand wavepackets into plane waves, then expand plane waves into spherical harmonics", and I think hopefully the spherical harmonics would cancel out nicely due to them being orthonormal when I do the scalar product.

Which also sounds right because partial waves decomposition works well with a central potential, which is exactly the case with this problem because the potential in the Hamiltonian is central in the Hamiltonian of the SE that I had to solve earlier.

hilbert2
Science Advisor
Gold Member
but how can I go about changing the $(\textbf{r} - \textbf{r}_0)^2$ into spherical polar coordinates?
That's not a problem. If the position vectors are $\mathbf{r}=(x,y,z)$ and $\mathbf{r}_0=(x_0 ,y_0 ,z_0 )$, then

$(\mathbf{r}-\mathbf{r}_0 )^2 = (x-x_0 )^2 + (y-y_0 )^2 + (z-z_0 )^2$,

and then you just use

$x = r\sin\theta \cos\phi$
$y = r\sin\theta \sin\phi$
$z = r\cos\theta$.

Only the $x,y,z$ are changed to these polar coordinates, the $x_0 ,y_0 ,z_0$ stay as they are, because they're just numbers and not variables.

That's not a problem. If the position vectors are $\mathbf{r}=(x,y,z)$ and $\mathbf{r}_0=(x_0 ,y_0 ,z_0 )$, then

$(\mathbf{r}-\mathbf{r}_0 )^2 = (x-x_0 )^2 + (y-y_0 )^2 + (z-z_0 )^2$,

and then you just use

$x = r\sin\theta \cos\phi$
$y = r\sin\theta \sin\phi$
$z = r\cos\theta$.

Only the $x,y,z$ are changed to these polar coordinates, the $x_0 ,y_0 ,z_0$ stay as they are, because they're just numbers and not variables.
Just in case that you're still interested in helping me, all of the equations above were followed closely from Chapter 10.1, 10.9, and 10.10 of the book "The Picture Book of Quantum Mechanics 4th edition", so it might be useful if you have a look. Meanwhile I'm trying all the approach I could lol

hilbert2
Science Advisor
Gold Member
I found the 3rd edition of that book as a pdf document, and it looks quite interesting. The section about Kepler wavepackets doesn't seem to go very far in explaining how to actually do the calculations.

If you have Wolfram Mathematica available, you can try to calculate the integral of the wavepacket multiplied by a hydrogenic wavefunction $\psi_{n,l,m_l}(r,\phi ,\theta )$ over all 3d space, representing all position variables in spherical polar coordinates. It's likely that it's not possible, and the integrations have to be done by some numerical method.

Just used Mathematica for the first time and had to figure out how to do integration. So doing the whole integral of the scalar product gave me this, not sure if it's right or what to make of it. Any opinions would be appreciated!

#### Attachments

• 19.2 KB Views: 91
That's not a problem. If the position vectors are $\mathbf{r}=(x,y,z)$ and $\mathbf{r}_0=(x_0 ,y_0 ,z_0 )$, then

$(\mathbf{r}-\mathbf{r}_0 )^2 = (x-x_0 )^2 + (y-y_0 )^2 + (z-z_0 )^2$,

and then you just use

$x = r\sin\theta \cos\phi$
$y = r\sin\theta \sin\phi$
$z = r\cos\theta$.

Only the $x,y,z$ are changed to these polar coordinates, the $x_0 ,y_0 ,z_0$ stay as they are, because they're just numbers and not variables.
So following your advice of converting my wavepacket wavefunction into spherical polar coordinates, I get
$$\begin{equation*} \begin{split} \psi(\textbf{r},0) &= \frac{1}{(2\pi)^{3/4}\sigma^{3/2}}\text{exp}\left[ -\frac{(\textbf{r} - \textbf{r}_0)^2}{4\sigma^2} + i\textbf{k}_0 \cdot \textbf{r}\right]\\ &= \frac{1}{(2\pi)^{3/4}\sigma^{3/2}}\text{exp}\left[ -\frac{1}{4\sigma^2}[(x - x_0)^2 + (y - y_0)^2 + (z - z_0)^2]\right] \text{exp}\left[i\textbf{k}_0 \cdot \textbf{r}\right]\\ &= \frac{1}{(2\pi)^{3/4}\sigma^{3/2}}\text{exp}\left[ -\frac{1}{4\sigma^2}[(r\sin\theta\cos\phi - x_0)^2 + (r\sin\theta\sin\phi - y_0)^2 + (r\cos\theta - z_0)^2]\right] \sqrt{\frac{4\pi}{2l+1}}\sum^{\infty}_{l=0}(2l+1)i^l j_l(kr) Y_{l,0}(\theta, \phi)\\ &= \eta(r, \theta, \phi) \sqrt{4\pi}\sum^{\infty}_{l=0}\sqrt{(2l+1)}i^l j_l(kr) Y_{l,0}(\theta, \phi) \end{split} \end{equation*}$$
So then I thought I'd use this expression to do the inner product instead, so I got to

$$\begin{equation*} \begin{split} b_{nlm} &= \langle\varphi_\textit{nlm}(\textbf{r})|\psi(\textbf{r},0)\rangle\\ &= \langle R_{nl}(r) Y_{lm}(\theta, \phi)|\eta(r, \theta, \phi) \sqrt{4\pi}\sum^{\infty}_{l=0}\sqrt{(2l+1)}i^l j_l(kr) Y_{l,0}(\theta, \phi)\rangle\\ \end{split} \end{equation*}$$
Then I realised here it's starting to look like that I should be able to use the orthonormalization relations of Spherical Harmonics, where
$$\int_{\theta = 0}^{\pi}\int_{\phi = 0}^{2\pi} Y^m_l {Y^{m'}_{ l'}}^{*} \sin\theta d\phi d\theta = \delta_{ll'}\delta_{mm'}$$
since $m$ for one of the Spherical Harmonics is already 0, we can get an expression for using that relation for $l=l'$ and $m=m'=0$, but I can't really do that with $\eta(r, \theta, \phi)$ in place, since I've converted it into using spherical polar coordinates. The way I see it could only work if I can somehow change $\frac{1}{(2\pi)^{3/4}\sigma^{3/2}}\text{exp}\left[ -\frac{(\textbf{r} - \textbf{r}_0)^2}{4\sigma^2}\right]$ to a scalar function that only depends on $r$.