Enumeration in C#

Enumerating a collection

In C#, a collection class needs to implement IEnumerable to be used in foreach. This interface provides a function GetEnumerator that returns an IEnumerator instance (called an enumerator, which may or may not be the collection itself) supporting operations such as Reset, Current, and MoveNext. An enumerator serves as a proxy to sequentially access the elements in a collection. However, the behavior of a enumerator can be problematic if the underlying collection is modified while the enumeration is in progress. There are several approaches to deal with this problem. One is to take a snapshot of the collection for enumeration. Another approach is to locate the elements on demand, so that the ith enumerated element is always the ith element, if any, in the current collection. To guarantee this, it may take linear time to visit each element and thus quadratic time to enumerate a collection using foreach. The third approach is to throw an exception when the enumerator detects that the contents of the collection changed during enumeration. The last approach is most frequently adopted in C#'s standard libraries for performance concerns. Even so, the first approach can be realized in the call site by creating a snapshot manually using ICollection.CopyTo.
In Java, Iterator.remove is the only safe way to modify a collection during iteration; the behavior of an iterator is unspecified if the underlying collection is modified in any other way. Besides, if you happen to have a ListIterator then you can use ListIterator.add to insert elements in addition to removing them during iteration.

Enumerators are not Enumerable

While it is enumerators that are used for enumeration, foreach has to work with an IEnumerable instance from which an enumerator is obtained. Namely, one can not directly use an enumerator in a foreach loop. In Java, this is not a problem because an anonymous interface instance can be easily created on the fly when it is needed. For example, one common scenario to instantiate an interface anonymously in Java is to create a task for a thread:

new Thread (new Runnable() { public void run() { ... } }).start();

In the snippet, an anonymous Runnable class is created, used and then discarded at runtime. Anonymous classes come in handy when Java programmers want to create classes to serve as workers, listeners, adapters, etc. By contrast, C# programmers are expected to use anonymous functions (delegates or lambdas) for these purposes. While C# does supports anonymous classes, anonymous classes in C# look more like structs than classes, as they can only have public read-only properties and they are not allowed to extend classes (other than object) or implement interfaces.
Although one cannot use an enumerator with foreach, it is almost as easy to use it in a while loop:

while (enum.MoveNext()) { Process(enum.Current); }

However, in situations where an enumerable instance is explicitly required (e.g., by the constructor of a collection), it will take some tricks to convert an enumerator to an enumerable instance at runtime.

Logical relations as a proof technique — Part I

This is the part I of my lecture note based on the lecture series "Logical Relations" given by Amal Ahmed at Oregon Programming Languages Summer School, 2015. See here for Part II.

Overview. These lectures will introduce the proof method of logical relations. Logical relations have been used to prove a wide range of properties, including normalization, type safety, equivalence of programs, non-interference (in security typed languages), and compiler correctness. We will start by looking at how to prove normalization and type safety for the simply-typed lambda calculus (STLC). We will then add recursive types to STLC and develop a step-indexed logical relation, using step-indexing to ensure that the logical relation remains well-founded. We will cover polymorphic and existential types and develop a logical relation to reason about the associated notions of parametricity and representation independence. The lectures will also discuss how to establish that a logical relation corresponds to contextual equivalence.

Free theorems. In parametric polymorphism, types cannot be checked, inspected or instantiated by programs. Hence, programs with polymorphic types don't change behaviors according to their type parameters. It follows that all functions of type $\forall a.\, a \rightarrow a$ are identity functions, since nothing but the argument can be returned by the function. Similarly, all functions of type $\forall a.\, a \rightarrow Int$ are constant functions. In contrast, there doesn't exist a function of type $\forall a.\, Int \rightarrow a$, since the return value can be obtained from nowhere. Results of this kind are called free theorems, and many of them can be proved using logical relations.

Logical predicates. Unary logical relations are also called logical predicates. They are used to assert properties of one single program, such as termination. In contrast, binary logical relations are used to assert relations between two programs, such as equivalence. Note that we only put closed programs in a logical relation. By and large, to prove that a well-typed and closed term $e$ has some special property, we shall set up a predicate $P_\tau$ such that $P_\tau(e)$ holds if and only if all of the following conditions hold: i) $\vdash e: \tau$, ii) $e$ has the property of interests, and iii) the property of interests is preserved by elimination forms (or introduction forms) of type $\tau$.

Strong normalization. In history, logical relations first appeared as a technique to prove strong normalization, which means that a well-typed term always reduces to the same value regardless of the strategies of evaluation. Since typing rules in STLC are deterministic, proving strong normalization for STLC amounts to showing that every well-typed program terminates, that is, $\vdash e:\tau \implies e \Downarrow$. Observe that this statement cannot be proved by applying structural induction on typing derivations directly. The problem of induction occurs at the application rule: we cannot conclude $e_1\ e_2\Downarrow$ from $e_1 \Downarrow$ and $e_2 \Downarrow$. Nevertheless, this problem can be resolved by using a stronger induction hypothesis, which leads to the following definition of logical predicates.
We shall set up a predicate $SN_\tau(\cdot)$ such that $SN_\tau(e)$ holds only if term $e$ has type $\tau$ and is strongly normalizing. We consider two types for STLC: boolean and function types. For boolean, it suffices to define that $SN_{bool}(e)$ holds iff $\vdash e:bool \wedge e \Downarrow$, since boolean doesn't has elimination forms. For functions types, we define that $SN_{\tau\rightarrow\tau'}(e)$ holds iff

$\vdash e':\tau\rightarrow\tau'\quad\wedge\quad e' \Downarrow\quad \wedge\quad (\forall e.\, SN_{\tau}(e) \implies SN_{\tau'}(e'\ e))$.

(1)

We shall see that this definition provides the strengthened hypothesis we need to deal with the application rule. The proof of $\vdash e:\tau \implies e \Downarrow$ now breaks into two parts: A) $\vdash e:\tau \implies SN_{\tau}(e)$, and B) $SN_{\tau}(e) \implies e \Downarrow$. Part B follows immediately from the definition of $SN_\tau(e)$. It remains to prove part A, and we rely on structural induction to do this. But again, we cannot use the statement $\vdash e:\tau \implies SN_{\tau}(e)$ we are going to prove as an induction hypothesis. It only asserts closed terms, but to show that it holds for terms of function type, say, for $(\lambda x:\tau.\, e'): \tau \rightarrow \tau'$, we need the hypothesis to hold for open term $e'$ so as to perform induction.
Instead, we shall use the following stronger claim as an induction hypothesis. Given a mapping $\gamma$ from variables to values, we write $\gamma \models \Gamma$ when $dom(\gamma)=dom(\Gamma)$ and $\forall x\in dom(\Gamma).\, SN_{\Gamma(x)}(\gamma(x))$ holds. That is, $\gamma$ maps variables in $dom(\Gamma)$ to values that strongly normalize with types compatible with $\Gamma$. We claim that for any term $e$ and substitution $\gamma$,

$\Gamma \vdash e:\tau \wedge \gamma \models \Gamma  \implies SN_{\tau}(\gamma(e))$.

(2)

It is clear that part A follows from this claim by setting $\Gamma$ to be empty. Now we shall prove (2) by inducting on the derivation rules case by case.
Case $\quad\rhd\ \ \vdash true: bool\quad$ and $\quad\rhd\ \ \vdash false: bool$.
　Note that $\gamma(x)=x$ for $x\in \{\,true,\,false\,\}$, which is of course strongly normalizing.
Case $\quad\Gamma(x)=\tau\ \ \rhd\ \ \Gamma \vdash x:\tau$.
　By premises $\Gamma(x)=\tau$ and $\gamma\models\Gamma$, $\gamma(x)$ has type $\tau$. Thus $SN_{\tau}(\gamma(x))$ holds by $\gamma\models\Gamma$.
Case $\quad\vdash e': \tau\rightarrow\tau'$, $\vdash e: \tau\ \ \rhd\ \ \vdash e'\ e: \tau'$.
　Note that $\gamma(e'\ e)=\gamma(e')\ \gamma(e)$. By premise $\gamma\models\Gamma$, both $SN_{\tau\rightarrow\tau'}(\gamma(e'))$ and $SN_{\tau}(\gamma(e))$ holds by induction. Thus $SN_{\tau'}(\gamma(e')\ \gamma(e))$, which equals to $SN_{\tau'}(\gamma(e'\ e))$, holds by (1).
Case $\quad\Gamma, x:\tau \vdash e':\tau'\ \ \rhd\ \ \Gamma\vdash\lambda x:\tau.\, e': \tau\rightarrow\tau'$.
　 This is the case we got stuck with before introducing logical predicates, so it is not surprising that the proof for this case is far more complicated than others. Before we proceed, we need the following two lemmas:
Lemma 1. (substitution lemma) If $\Gamma \vdash e:\tau$ and $\gamma \models \Gamma$, then $\Gamma \vdash \gamma(e):\tau.$
Lemma 2. (SN is preserved by forward/backward reductions) Suppose that $\vdash e:\tau$ and $e \hookrightarrow e'$. Then $SN_{\tau}(e)$ holds if and only if $SN_{\tau}(e')$ holds.
To prove (2), i.e., $\Gamma \vdash \lambda x:\tau.\, e' \wedge \gamma\models\Gamma$ implies $SN_{\tau\rightarrow\tau'}(\gamma(\lambda x:\tau.\, e'))$, it suffices to prove the following three statements:
1) $\gamma(\lambda x:\tau.\, e')$ has type $\tau\rightarrow\tau'$. This follows by premise $\gamma\models \Gamma$ and Lemma 1;
2) $\gamma(\lambda x:\tau.\, e')\Downarrow$. This is obvious, since $\gamma(\lambda x:\tau.\, e')=\lambda x:\tau.\,\gamma(e')$ is a lambda value;
3) $\forall e.\,SN_\tau(e)\implies SN_{\tau'}(\gamma(\lambda x:\tau.\, e')\ e)$. By induction, we know that $\Gamma, x:\tau \vdash e':\tau'$ and $\gamma'\models \Gamma, x:\tau$ implies $SN_{\tau'}(\gamma'(e'))$. Suppose that $SN_\tau(e)$ holds for some $e$. Then $e\Downarrow v$ for some value $v$, and $SN_{\tau}(v)$ holds by Lemma 2. Given premises $\Gamma \vdash \lambda x:\tau.\, e'$ and $\gamma\models\Gamma$, define $\gamma'=\gamma[x \mapsto v]$. Note that $\gamma'\models\Gamma,x:\tau$. Thus $SN_{\tau'}(\gamma'(e'))$ holds by induction. However, note that
$\gamma'(e')=\gamma(e')[x \mapsto v]=(\lambda x:\tau.\,\gamma(e'))\ v=\gamma(\lambda x:\tau.\,e')\ v$.
It then follows from $SN_{\tau'}(\gamma'(e'))$ and Lemma 2 that $SN_{\tau'}(\gamma(\lambda x:\tau.\,e')\ e)$ holds. This concludes our proof of strong normalization for STLC.

(End of lecture note Part I. See here for Part II.)

Tips for writing Spark programs

Writing serializable closures

When you use a closure, Spark will serialize it and send it around the cluster. This means that any captured variables must be serializable. But sometimes it can't be. Consider this example:

class Foo {
  val factor = 3.14159
  val log = new Log(...) // not serializable
  def multiply(rdd: RDD[Int]) = {
    rdd.map(x => x * factor).reduce(...)
  }
}

The JVM must pass the closure to map to capture variable factor. However, since the closure contains an unserializable variable log, a NotSerializableException will be thrown at runtime. A work around here is to assign factor to a local variable:

class Foo {
  val factor = 3.14159
  val log = new Log(...) // not serializable
  def multiply(rdd: RDD[Int]) = {
    // Only factor2 will be serialized.
    val factor2 = factor 
    rdd.map(x => x * factor2).reduce(...)
  }
}

Serialization is a general issue for distributed programs written on the JVM. A future version of Scala may introduce a "serialization-safe" mechanism for defining closures for this purpose.

Preventing the driver program from OOM

1. Operations on RDDs includes transformations and actions. When the dataset is distributed over the workers, the outputs of transformations are still distributed. In contrast, actions will transfer the outputs to driver program, which can lead to OOM for large outputs. A rule of thumb here is to avoid using actions that may return unbounded outputs, such as collect, countByKey, collectAsMap, etc. If you need to probe a large RDD, try to use actions such as sample, count and take to obtain bounded results. Also, you can write outputs directly from workers to the HDFS using saveAsTextFile and then load them later.
2. Sometimes, you get an OOM Error not because your RDDs don't fit in memory, but because the working set of one of your tasks was too large. Spark’s shuffle operations (sortByKey, groupByKey, reduceByKey, join, etc) build a hash table within each task to perform the grouping, which can often be large. The simplest fix here is to increase the number of tasks, so that each task’s input set is smaller. Spark reuses one executor across many tasks with low launching costs, so you can safely increase the number of tasks to more than the number of cores in your clusters.

Avoiding imbalant join operations

1. Suppose that you want to inner join two RDDs in Spark, one is very small with $n$ partitions and the other is very large with $m\gg n$ partitions. This operation will shuffle the data in the large RDD into $n$ partitions, carried out by at most $n$ executors regardless of the number of work nodes. A better way to do inner join in this case is to first broadcasts the small RDD to the nodes on which the large RDD is distributed, and then performs the join operation at each executor. In this way, you not only avoid the need of shuffling, but also maintain the parallelism of the larger RDD (with $m$ output partitions instead of $n$). The benefits are obtained at the cost of replicating and broadcasting the smaller RDD around the cluster.
2. If you want to left join a very small RDD with a very large RDD, it may be helpful to first filter the the large RDD for only those entries that sharing keys with the small RDD. This filtering step reduces the data needed to be shuffled over the network and thus may speed up the overall execution.

Note that this kind of balancing tricks are not needed for reduce-like operations such as groupByKey and reduceByKey, where the largest parent RDD’s number of partitions will be used.

References and resources

1. The official guidelines for turning Spark
2. How-to: Tune Your Apache Spark Jobs at Cloudera
3. CME 323: Distributed Algorithms and Optimization at Stanford Univ.

Designing algorithms with Spark aggregate – Part I

In this post, I shall give a simple framework to design algorithms in terms of Spark's aggregate function. This framework is just a proof of concept and should not be applied blindly in practice. However, it does induce efficient algorithms in some special cases, as will be shown in the examples.

In general, our goal is to compute a property of any given instance from a domain. Let $G$ denote the domain of all interested instances, and $C$ the domain of the desired property. For example, if $G$ denotes graphs and the property to compute is MST, then $C$ will be the set of trees. Furthermore, we assume that there exists an associative and commutative "composition" operator $\oplus : G \rightarrow G \rightarrow G$, such that an instance of interests can be decomposed into (and composed back from) sub-instances in $G$. We shall denote the desired property of instance $g$ as $f(g)$. For example, if $g$ is a graph, then $f(g)$ can be some graph property such as diameter, MST and vertex cover. Note that the the mapping $f : G \rightarrow C$ can be non-deterministic if the property is not unique w.r.t. the given instance.

To use the aggregate function, we specify zero, seg and comb as follows:
$\mathrm{zero} \in (G, C)$ is a pair comprised of the zero values from $G$ and $C$. If there are no such zero values, we can set $\mathrm{zero} = (g, f(g))$ for an arbitrary sub-instance of $g$.
$\mathrm{seq} : (G, C) \rightarrow G \rightarrow (G, C)$ takes a partially aggregated result $(g_1, f(g_1))$ and a sub-instance $g_2$. It outputs $(g_3, f(g_1 \oplus g_2))$ such that $f(g_3 \oplus g) = f(g_1 \oplus g_2 \oplus g)$ for all $g \in G$.
$\mathrm{comb} : (G, C) \rightarrow (G, C) \rightarrow (G, C)$ takes two partially aggregated results $(g_1, f(g_1))$ and $(g_2, f(g_2))$. It outputs $(g_3, f(g_1 \oplus g_2))$ such that $f(g_3 \oplus g) = f(g_1 \oplus g_2 \oplus g)$ for all $g \in G$.
The rest of the steps are standard: we first decompose the input instance into a sequence of sub-instances. The sequence is then converted to an RDD through parallelize. Finally, we invoke RDD.aggregate(zero)(seq,comb) to compute the desired property of the instance.

Remarks

1. While functions seq and comb look very similar in spec, they can use different methods to compute $f(g_3)$. Due to the partial results provided, seq is suitable to use incremental methods and comb can use (the merge part of) divide-and-conquer methods to compute the desired property. See the convex hull example below.
2. The assumption that operator $\oplus$ is associative and commutative over $G$ can be relaxed to be so only with respect to the given property $f$. For instance, suppose that $G$ is a set of finite strings and $\oplus$ is the string concatenation operator. While $\oplus$ is not commutative per se, it is commutative w.r.t. the string length property. That is, $\mathrm{strlen}(g_1 \oplus g_2) = \mathrm{strlen}(g_2 \oplus g_1)$ for all $g_1, g_2 \in G$. The framework is thus applicable to compute the length of string. See the majority algorithm below for a more non-trivial example.
3. The framework is practical only if $|g_3|$ can be significantly less than $|g_1 \oplus g_2|$, for otherwise we will have to pass almost the whole instance from node to node. In the case that $f(g_1 \oplus g_2)$ can be computed from $f(g_1)$ and $f(g_2)$, there is no need to output $g_3$, i.e., $|g_3|=0$. This special case allows us to operate on simpler domains (e.g., the type of seq becomes $C \rightarrow G \rightarrow C$) and to obtain more efficient aggregate algorithms due to less communication overheads.

Examples

The above naive framework can yield efficient algorithms if $f(g_1 \oplus g_2)$ can be computed solely based on $f(g_1)$ and $f(g_2)$. Below we give three examples that satisfy this condition. The first is a variant of the reverse-delete algorithm that finds an MST of a graph. The second is an outline of a convex hull algorithm. The last is Boyer and Moore's algorithm that finds the majority element in a list.

Finding MST. Let $G$ be the set of all weighted undirected graphs, and $C$ be the set of trees. Given a graph $g$ from $G$, we assign a unique index to each of its edges. Let $w(e)$ and $idx(e)$ denote the weight and index of edge $e$, respectively. Given two edges $e_1$ and $e_2$ of $g$, we say $e_1 < e_2$ iff either $w(e_1) < w(e_2)$, or $w(e_1) = w(e_2)$ and $idx(e_1)$ < $idx(e_2)$. Now we define an operator msf over $G$, such that msf computes a minimum spanning forest of $g_1 \cup g_2$ from $g_1$ and $g_2$ as follows:

def msf (g1, g2) = {
  var g = g1 ∪ g2
  while (g contains a cycle c)
    remove the largest edge from c
  return g
}

Let RDD(g) denote the RDD obtained from decomposing and parallelizing graph $g$. Then RDD(g).aggregate(Nil)(msf,msf) computes a minimum spanning forest of $g$, which would be a MST if $g$ is connected. Note that the result is deterministic if the assignment of indices is deterministic.

Finding convex hull. Here $G$ is the collection of finite sets in Euclidean plane and $C$ is the set of convex polygons. The convex hull of a set of points is the smallest convex polygon containing all the points. If the points are given as a sequence, we can use aggregate to find their convex hull by computing smaller convex hulls incrementally in seq and then merging them in comb to obtain the final convex hull. Both of seq and comb can be implemented to take linear time (see the links for pseudo codes).

Finding majority. Let $C$ be a set of comparable elements. Fix some element $\bot \in C$ and let $G = (C - \{\bot\})^*$ be a set of finite sequences. The majority problem is to determine in a given sequence from $G$ whether there is an element with more occurrences than all the others, and if so, to determine this element. Click here to see a linear-time algorithm that solves the problem using aggregate. Note that the aggregate invocation only finds a candidate for majority; we need a further check to assert that the candidate is indeed the majority element.

// Input:  a list of elements and an element 
// Output: a boolean value indicating whether elem is the majority of list
def check_majority (list: G, elem: C): Boolean = {
  def zero = (elem, 0)
  def op ((a: C, c: Int), b: Int) = {
    if(a != b)
      return (a, c - 1)
    else
      return (a, c + 1)
  }
  val res = list.fold(zero)(op)
  return res._2 > 0
}
// Input:  a list of comparable elements
// Output: the majority element of the list, or ⊥ if there is no such element
def find_majority (list: G): C = {
  def zero = (⊥, 0)  
  def seq ((a: C, c: Int), b: Int) = { 
    if(c == 0) 
      return (b, 1)
    if(a == b)
      return (a, c + 1)
    else
      return (a, c - 1)
  }  
  def comb ((a1: C, c1: Int), (a2: C, c2: Int)) = {
    if(a1 == a2) 
      return (a1, c1 + c2)
    if(c1 > c2)
      return (a1, c1 - c2)
    else
      return (a2, c2 - c1)
  }
  val res = list.aggregate(zero)(seq, comb)

  if(check_majority(list, res._1))
    return res._1
  else
    return ⊥
}

Basic dependent type theory ─ Part I

Computational trinitarianism. Proof Theory, Category Theory, and Type Theory, to a significant degree, can be regarded as fields that see the same subject from different points of view. They can be related via the follows:
TT ⇔ CT: category semantics; WTT corresponds to morphisms and categories.
PT ⇔ CT: category logic; an algebraic way of thinking about logic and proofs.
PT ⇔ TT: thinking propositions as types and WTT as proofs.
When you are learning a new idea in one field, it is helpful to transpose it to the other two fields and seek its correspondence there. See here and here for more details.
Judgements. We start with three kinds of categorical judgements:
$A\ type$　　　─ 　$A$ is a type.
$a \in A$　　 　─ 　$a$ is a value of type $A.$
$a_1=a_2: A$　─ 　$a_1$ and $a_2$ are "definitionally equivalent" values of type $A.$
A judgement in form of $x_1:A_1,...,x_n:A_n \models a(x_1,...,x_n):A$ is called a hypothetical judgement. One may regard such a judgement as a mapping from $A_1,...,A_n$ to $A.$ $a(x_1,...,x_n)$ is called an open term that might involve free variables $x_1,...,x_n.$ In practice, we often consider the open term parametrized at $x_1,...,x_n.$ A closed term of type $A$ can be determined categorically after the values of $x_1,...,x_n$ are plugged in.
Definitional equality is 1) an equivalence relation, 2) a congruence (see typing rules), and 3) respected by maps. That is, $$\frac{\Gamma, x:A \models b:B, \Gamma \models a_1\equiv a_2: A}{\Gamma \models [a_1/x]b \equiv [a_2/x]b : B}.$$The last feature is called the principle of functionality property. This property and its related concepts will play a central role later and in Homotopy type theory.
Definitional equality has its limits. For example, one cannot prove definitional equivalence between $0+1$ and $1+0$ even given all the arithmetic rules about natural numbers. Definitional equality corresponds to calculating by simplification. There however are equations (such as $0+1=1+0$) that you wish to hold cannot be reached that way. These equations need not just calculation but proof, and thus cannot hold definitionally.
Typing rules. We have the following kinds of rules to reason about types. Rules 1-5 come from Proof Theory and Rule 6 comes from Category Theory.

1. Formation rules: how to construct a type.
2. Introduction rules: how to create elements of a type.
3. Elimination rules: how to use elements of a type.
4. Congruence rules: how to interact with definitional equivalence.
5. Computation rules: how to do calculation with types
6. Uniquicity/Unisersality rules: characterize precisely types

Negative types. E.g., Unit, product, and function. Elimination-oriented. Uniquicity of product: $\Gamma \vdash \langle fst(c), snd(c) \rangle \equiv c : A\times B.$ Computation of function is substitution. Elimination of function is application. Extension of function can be written as $$\frac{\Gamma, x:A \vdash a_1(x)\equiv a_2(x):B}{\Gamma \vdash a_1\equiv a_2: A\rightarrow B}$$This means that functions are characterized by their input-output behaviors.
Positive types. E.g., Void and sum. Introduction-oriented. The Uniquicity property of positive types is problematic, and cannot be shown as for the negative types as directly. Elimination rule for the Void type:$$\frac{\Gamma\vdash a:\emptyset}{\Gamma\vdash abort(a): A}$$This means if you have a value of the Void type, you can abort the program with any type. Note that one can never build a term of type $A$ in this way since there is no closed term of type $\emptyset.$ One application of this rule is that if you are in a conditional branch, then one side of the conditional may be contradictory, and you can build an $abort$ element of type $A$ to make sure both elements for the condition have the same type.
Negative types vs. positive types. Approximately, the two types differ most noticeably in how they behave under elimination. Reasoning about a negative type has noting to do with its structure. For example, a product is characterized by its behaviors on $\texttt{fst}$ and $\texttt{snd}$ without caring about what it actually is. A product can be realized by an ordered pair or just a black box that computes the $\texttt{fst}$ and $\texttt{snd}$ components on demand. Also, there is no induction principle for negative types. For example, all you know about a function is that it is a lambda term. Thus you cannot induct on its structure via the typing rules. In contrast, the structure of a positive type is enforced in the theory. The elimination form of positive types work by case analysis─the actual element in a sum is extracted by $\texttt{case}$ in the process of inference. Thus it cannot be anything else after the inference.
A correspondence btw STT and IPT. $Unit$ is $true.$ $Void$ is $false.$ Sum is disjunction. Product is conjunction. Function is implication. $A\ type$ is $A\ prop.$ $\exists e: A$ is $A\ true.$ $\exists b.\, a:A \vdash b:B$ is $A\le B.$ Thus $A\le B$ means that there is a function transforming a proof for $A$ to a proof for $B.$ Note that $\le$ is a pre-order: that is, it is reflexive and transitive. We use $A\supset B$ to denote the largest element $x$ satisfying $A\wedge x \le B$, i.e., $C \wedge A \le B$ implies $C \le A \supset B.$ If an ordered structure satisfies 1) $\forall A: A\le \top$; 2) $\forall A: \bot \le A$; 3) $A\wedge B$ (meet) is a glb for $A,\,B$; 4) $A\vee B$ (join) is a sub for $A,\,B$; 5) $\wedge$ is distributive over $\vee$; 6) $A \supset B$ exists for every $A,B$, then it is called a Heyting pre-algebra. (It is called a pre-algebra because $\le$ is not anti-reflexive.)
Negations. $\neg\ A$ is defined as $A \supset \bot.$ Thus we can get a contradiction from a value of $A.$ The fact $A\le \neg\neg\ A$ follows directly by the properties of $\supset.$ Its converse $\neg\neg\ A \le A$ however doesn't hold in general (it takes quite some efforts to check this statement). Similarly, while $A \vee \neg\ A \le \top$ is trivial, $\top \le A \vee \neg\ A$ cannot be proved in general without knowing something about $A.$
By definition, $\neg A$ is the largest element inconsistent with $A.$ That is, $A\wedge\neg A \le \bot$ and $B\wedge\neg A \le \bot \implies B\le \neg A.$ The compliment $\overline{A}$ of $A$ is defined as the smallest element $B$ such that $\top\le A\vee\overline{A}.$ Hence $\top\le A\vee B$ implies $\overline{A}\le B.$ Note that i) $\neg A$ always exists while $\overline{A}$ does not; ii) it is reasonable to add the largest/smallest requirement to the definitions, for otherwise we can trivially choose $\neg\ A=\bot$ and $\overline{A}=\top$; iii) it is easy to check $\neg\ A \le \overline{A}$ (things in $\neg\ A$ are definitely not in $A$). The converse $\overline{A} \le \neg\ A$ however does not hold in general (the negation does not fill up the whole space).
A Boolean algebra is also called a complemented Heyting algebra in that $\overline{A}$ exists for every $A$, with i) $\top \le A\vee \neg\ A$ (which implies $A\vee \neg\ A=A\vee \overline{A}$), and ii) $\overline{\overline{A}} \le A$ (note that $A \le \overline{\overline{A}}$ follows by definition). Heyting pre-algebras generalize Boolean pre-algebras in the sense that a Heyting pre-algebra satisfying $1 \le A \vee \neg A$ (law of excluded middle), or equivalently $\neg\neg A \le A$ (law of double negation), is a Boolean pre-algebra.
Heyting algebras. An exponential $y^x$ is defined as $\overline{x}\vee y.$ A Heyting algebra can also be defined as a lattice with exponentials: $H$ is called a Heyting algebra if i) it has a preorder; ii) it has finite meets and joins, and iii) satisfies $y^x$ is the largest element $z$ such that $x \wedge z \le y.$
Exercise. Show that every Heyting algebra is distributive.
We give hints about proving $x \wedge (y \vee z) = (x \wedge y) \vee (x \wedge z).$ First, note that $z \wedge x \le y$ iff $z \le y^x.$ Hence, to prove direction $x \wedge (y \vee z) \le (x \wedge y) \vee (x \wedge z)$, it suffices to show $y \vee z \le ((x \wedge y) \vee (x \wedge z))^x$, which amounts to show that $u \le ((x \wedge y) \vee (x \wedge z))^x$ for $u=y, z.$ However, $u \le ((x \wedge y) \vee (x \wedge z))^x$ iff $u \wedge y \le (x \wedge y) \vee (x \wedge z).$ For $u=y, z$ this follows by definition. As to the other direction, note that $x \le y$ iff $\forall z.\,(y\le z$ $\implies x \le z).$ Hence, to prove $(x \wedge y) \vee (x \wedge z) \le x \wedge (y \vee z)$, it suffices to show that $(x \wedge y) \vee (x \wedge z) \le u$ implies $x \wedge (y \vee z) \le u$, which is equivalent to $y \vee z \le u^x.$ Note that $(x \wedge y) \vee (x \wedge z) \le u$ implies both $x \wedge y \le u$ and $x \wedge z \le u.$ Thus we have $y \le u^x$ and $z \le u^x$, that is, $y \vee z \le u^x.$
As an alternative, you can also use the completeness theorem for IPT. Completeness assures that anything that holds in every Heyting algebra must be provable in IPT, and vice versa. Since preorder is interpreted to entailment in IPT, you can transfer e.g., $x \wedge (y \vee z) \le (x \wedge y) \vee (x \wedge z)$ to $X \wedge (Y \vee Z)\ true \vdash (X \wedge Y) \vee (X \wedge Z)\ true$, and then reason about the logical connectives to do the proof. You can further add proof terms to obtain a FP problem $\Gamma, x: A \times (B + C) \vdash M: (A \times B) + (A \times C).$ Resolving this problem amounts to writing a program $M$ of type $(A \times B) + (A \times C)$, and it is not hard to see that $M = \texttt{case}(z.\,\texttt{inl}(\langle \texttt{fst}(x), z\rangle);\, z.\,\texttt{inr}(\langle \texttt{fst}(x), z\rangle))\ \texttt{snd}(x).$ It is sort of a cyclelogical matter for people to transfer problems around the trinitarian when they are doing proofs.

A logical interpretation of the function subtyping rule

We discussed the subtyping rule for function type in this post. In short, a function type is contravariant in the argument types and covariant in the result type. Several days ago, I happened to notice that the Monotonicity Theorem for implication (see p.70 in [1]) looks very similar to the aforementioned subtyping rule in form. Given propositions $p,p',q,q'$ such that $p'\Rightarrow p$ and $q'\Rightarrow q$, the theorem asserts $$p'\vee q' \Rightarrow p \vee q.$$ Besides, it is immediate to see that $p \vee q \equiv \neg q \Rightarrow p$ and $q' \Rightarrow q \equiv \neg q \Rightarrow \neg q'$, assuming the law of double negation. Now, by writing $p$, $p'$, $\neg q$ and $\neg q'$ as $P$, $P'$, $Q$ and $Q'$, respectively, we obtain \begin{equation} (Q' \Rightarrow P') \Rightarrow (Q \Rightarrow P) \label{subtyping} \end{equation} given that $P'\Rightarrow P$ and $Q\Rightarrow Q'$.

The connection between the above fact and the subtyping rule of function types can be established by treating the implication connective "$\Rightarrow$" as both a function-type constructor and a subtype-of relation. Indeed, if Apple is a subtype of Fruit then being an Apple implies being a Fruit. It is thus reasonable to use "Apple $\Rightarrow$ Fruit" to express the subtyping relation between Apple and Fruit. In view of this, the logical entailment in \eqref{subtyping} leads us to the following statement for function types: If $P'$ is a subtype of $P$ (i.e. $P' \Rightarrow P$) and $Q$ is a subtype of $Q'$ (i.e. $Q \Rightarrow Q'$), then function type $Q' \Rightarrow P'$ is a subtype of function type $Q \Rightarrow P$ (i.e. $(Q' \Rightarrow P') \Rightarrow (Q \Rightarrow P)$). In other words, a function type is contravariant in the argument types and covariant in the result type.

Reference

1. Gries, David, and Fred B. Schneider. "A logical approach to discrete math." Springer, 1993.

Functors and Monads in Scala ─ Part I

A type constructor is a type that you can apply type arguments to construct a new type: given a type A and a type constructor M, you get a new type M[A]. A functor is a type constructor, say M, together with a function

    lift: (A -> B) -> M[A] -> M[B]

that lifts a function over the original types into a function over the new types. Moreover, we require that the lifting function should preserve composition and identities. That is, if function h is formed by composing functions g and f, then lift(h) should be equivalent to the composition of lift(g) and lift(f). Formally, this means

    h(a) = g(f(a))  implies  lift h a = lift g (lift f a).

Also, if h is an identity function on variables of type A, then lift(h) should be an identity function on variables of type M[A]:

    h(a) = a  implies  lift h a = a.

A monad is a functor associated with two special operations, unitizing and binding:

    unit: A -> M[A]
    bind: (A -> M[B]) -> M[A] -> M[B]

A unitizing operation takes an object of the original type and maps it to an object of the new type. A binding operation is useful in composing two monadic functions, i.e., functions with signature $*\rightarrow M[*]$. Using binding, the composition of two monadic functions f: A -> M[B] and g: B -> M[C] can be defined as

    a => bind g f(a),

which has type A -> M[C]. Composition of monadic functions is supposed to be associative and has unit as a left and right identity. For this purpose, definitions of bind and unit must satisfy three "monad laws". That is, the following three equations

1.    bind (a => bind g f(a)) la = bind g (bind (a => f(a)) la)
2.    bind f (unit a) = f(a)
3.    bind unit la = la

should hold for any values a with type A, la with type M[A], functions f with type A -> M[B] and g with type B -> M[C].

A monad can be thought of as representing computations just like a type represents values. These computations are operations waiting to be carried out. Presumably, the programmer has some interest in not carrying them out directly. For example, the computations to be done may not be fully known at compile time. One can find similar ideas of representing computations as objects in the Command design pattern of GoF. Regarding monads as computations, the purpose of unit is to coerce a value into a computation; the purpose of bind is to construct a computation that evaluates another computation and yields a value. Informally, unit gets us into a monad and bind gets us around the monad. To get out of the monad, one usually uses an operation with type M[A] -> A, though the precise type depends on the purpose of the monad at hand.

In Scala, a monad is conventionally treated like a collection associated with list operations: the unit operation corresponds to a constructor of singleton lists; lift and bind are named as map and flatMap, respectively, in the sense that lift maps a list to another list and bind maps a list of lists to another list of lists and then flattens the result to a one-layered list. Map and flatMap operations kind of subsume each other given the flatten operation, as

    flatMap g la = flatten (map g la)
    map (a => f(a)) la = flatMap (a => unit f(a)) la

Hence, the minimal set of operations a monad must provide can be {unit, flatMap} or {unit, map, flatten}. Scala doesn't define a base trait for Monad and leaves the decision of interface to the programmer. On the other hand, if a monad implements methods map, flapMap and filter, then it can be composed with other monads via for-comprehension to streamline monadic computations.

Examples

Optional values. Option (aka. Maybe) monads are used extensively in many functional programming languages: see here for a more general treatment of option monads, and here for a detailed explanation of their usages in Scala. Below is a minimalistic implementation of Option type, and its companion subtypes Some and None:

sealed trait Option[T] {
    def map[U](f: T => U):Option[U]
    def flatMap[U](f: T => Option[U]): Option[U]
}
case class Some[T](val t:T) extends Option[T] { 
    override def map[U](f: T => U): Option[U] = new Some[U](f(t))
    override def flatMap[U](f: T => Option[U]): Option[U] = f(t)
}
case class None[T] extends Option[T] {
    override def map[U](f: T => U): Option[U] = new None[U]
    override def flatMap[U](f: T => Option[U]): Option[U] = new None[U]
}

The canonical advantage of the option monad is that it allows composition of partial functions, i.e., functions that return a null pointer when they fail to give correct results. Instead, a partial function returns None (representing an "empty" Option) and this result may be composed with other functions. In other words, it allows transformation of partial functions into total ones. For example, suppose g1, g2 and g3 are functions that return either a null pointer or an addable value, and we want to print g3(g2(x) + g1(x)) given that the result is meaningful. Without Option, one has to use nested if-then null checks or a try-catch block to handle possible null inputs/outputs, for example,

try{ 
  val y = g3(g2(x) + g1(x))
  if(y != null) println(y)
}catch {
  case e : NullPointerException => ; 
}

If we make g1, g2 and g3 return an option monad via a wrapper, we can achieve the same purpose using for-comprehension, so that nothing would be printed if any of g1, g2 or g3 returns null:

def monadify[T](f: T => T) =
  (x: T) => { 
    val y = if(x!=null) f(x) else x
    if(y != null) Some(y) else None
  }
for {
  y <- monadify(g1)(x)
  z <- monadify(g2)(x)
  r <- monadify(g3)(y + z)
} yield println(r)

Whether the functional ("idiomatic") way of handling null pointers is preferable to the traditional one is still under debate, especially among sophisticated users of imperative programming languages. Since these two approaches differ very little in performance and Scala supports both of them, the choice is really only a matter of taste than truth IMHO.

An important thing to notice about this example is that the monad is used here as a control structure. We have got a sequence of operations chained together in a monad via for-comprehension. When we evaluate the monadic sequence, the implementation of the binding operator that combines monads actually does control flow management. This trick of using the chaining operator to insert control flow into monadic sequences is an incredibly useful technique and is used in Scala by a variety of monads beside Option, such as Try, Either, Future and Promise.

State transformers. A stateful computation can be modeled as a function that takes inputs with an initial state and produces an output paired with the new state. Due to the presence of states, a function may produce different results even though it is given the the same inputs. Let Env be a type representing the possible states of the environment. An Effect is a mapping from Env to Env that describes changes of states, i.e., the side-effects caused by a computation, and the monad type M adjoins an extra side-effect to any type:

    Effect = Env -> Env
    M[A] = (A, Effect)

The unit function pairs up a value with the identity on side-effects; the map function propagates a state leaving it unchanged:

    unit a = (a, e => e)
    (a, e) map f = (f(a), e)

If you want to adjoin two side-effects e and e', you can use bind to compose them and get a single effect e' ○ e, where ○ denotes function composition:

    (a, e) bind (a => (f(a), e')) = (f(a), e' ○ e)    

In Scala, the state monad described above may be written as

case class M[+A](a: A, e: Env => Env) {
    private def bind[B](f: A => M[B]): M[B] = {
      val fa = f(a); M(fa.a, env => fa.e(e(env)))
    }
    override def toString = showM(this)
    def map[B](f: A => B): M[B] = bind(x => unitM(f(x)))
    def flatMap[B](f: A => M[B]): M[B] = bind(f)
}
def unitM[A](a: A) = M(a, env => env)

The following example shows a state monad that increments a counter each time an arithmetic operation is invoked. (Note that one can actually declare the counter as a field member, avoiding the use of monads. This may explain the viewpoint that state monads are arguably less important in impure FPs like Scala than in pure FPs like Haskell.) The add and div operations demonstrate how to handle a monad wrapped in another monad using nested for-comprehensions.

type Env = Int
// Ad hoc operations for this particular monad
def showM[A](m: M[A]) = "Value: " + m.a + "; Counter: " + m.e(0)
val tickS = M(0, e => e + 1)

// Following code demonstrates the use of monads
type Value = Double
def add(a: M[Option[Value]], b: M[Option[Value]]): M[Option[Value]] = 
  for {
    _  <- tickS
    aa <- a 
    bb <- b
    c = for {
      a <- aa if aa.nonEmpty
      b <- bb if bb.nonEmpty
    } yield a + b     // yield None if aa is None or bb is None
  } yield c

def div(a: M[Option[Value]], b: M[Option[Value]]): M[Option[Value]] = 
  for {
    _  <- tickS
    aa <- a
    bb <- b
    c = for {
      a <- aa if aa.nonEmpty
      b <- bb if bb.nonEmpty
      b <- Some(b) if b != 0
    } yield a / b     // yield None if b==0 or aa is None or bb is None
  } yield c

val const0 = unitM(Some(0.0))
val const1 = unitM(Some(1.0))
val const2 = unitM(Some(2.0))
val const3 = unitM(Some(3.0))
val expr = div(add(const3, add(const1, const2)), const2)  
println(expr)                            // Value: Some(3.0); Counter: 2
println(add(div(expr, const0), const1))  // Value: None; Counter: 5

The syntax sugar of for-comprehensions abstracts the details and help the user manipulate monads without getting involved in the underlying machinery of binding and lifting. The fact that the abstraction works as expected is a consequence of the monad laws, and it is up to the designer of the monad to assure that the implementation does keep the laws.

We just stop here as this post is getting too long. More examples, techniques and discussions will be introduced in the Part II of this tutorial.

References and resources

1. Monads are Elephants, Part 1 ~ Part 4, is a gentle and pragmatic tutorial that explains exactly what a Scala programmer needs to know to exploit monads. Scala Monads: Declutter Your Code With Monadic Design is one of the videos on YouTube that teaches you how to write modular and extensible programs with monadically structured code.
2. Many excellent papers and slides about monads could be found in the homepage of Philip Wadler, one of the most enthusiastic advocates of monadic design in the FP community. Among others, The Essence of Functional Programming and Monads for Functional Programming are two must-reads for anyone who is struggling to turn monadic reasoning into his intuitions.
3. For a general and rigid mathematical treatment of monads, see e.g., Arrows, Structures, and Functors: The Categorical Imperative, where Section 10.2 describes monads as generalized monoids in a way that basically makes everything that we can write down or model using abstract syntax a monad. Mathematically mature computer science readers will find everything they need to know about the subject in the celebrated book Categories for the Working Mathematician by Mac Lane, the co-founder of category theory.
4. There are also books available online that introduce programmers with moderate math background into the subject, including Category Theory for Computing Science by M. Barr and C. Wells, and Computational Category Theory by D.E. Rydeheard and R.M. Burstall. Also check this thread and references therein about the relations between category theory, type theory, and functional programming.

Equivalence Checking in Lambda Calculus

Definitional equivalence checking is the problem of determining whether two λ-terms in λ-calculus are equivalent in definition. (By "definitional", we mean that equivalence is defined directly by rules, rather than by an appeal to an operational semantics.) Equivalence of λ-terms in untyped λ-calculus is undecidable, see Scott's Theorem in [3]. This result is not surprising, since untyped λ-calculus is Turing complete. On the other hand, it is decidable for simply typed λ-calculus (aka. simple type theory, STT) with a single base type, see Sec. 6 of [4]. Equivalence checking would become undecidable once we extend STT with dependent types such as lists, or recursive data types such as natural numbers, see Chap 4 and 5 of [5].

It is interesting to check the weaker relation that whether two λ-terms are equivalent in operational semantics, i.e., they reduce to the same value given the same context. For example, $\lambda x.x:\;Nat\rightarrow Nat$ and $\lambda x.pred(suc\ x): Nat\rightarrow Nat$ are not definitionally equivalent but they are operationally equivalent. For untyped λ-calculus this problem is clearly undecidable, but for STT I am not sure. I haven't found enough materials to say anything on this subject by far.

References

1. Notes on Simply Typed Lambda Calculus (another lecture note)
2. Simply Typed Lambda Calculus on Wiki
3. Practical Foundations for Programming Languages (Preview of the 1st ed.)
4. Advanced Topics in Types and Programming Languages (available online)
5. http://www.cs.cmu.edu/~rwh/courses/logic/www/handouts/logic.pdf

Scala: Subtyping

When a class inherits from another, the first class is said to be a nominal subtype of the other one. It is a nominal subtype because each type has a name, and the names are explicitly declared to have a subtyping relationship. On the other hand, a structural subtype relationship holds if two types have the same members. While the compiler only verifies this property at the level of type checking, a subtype should correctly implement the contracts of its supertypes so that Liskov's Substitution Principle applies, which essentially says that beside supporting the same operations, the operations in the subtype should require no more and provide no less than the corresponding operations in the supertype.

The Scala compiler allows a type’s subtypes to be used as a substitute wherever that type is required. For classes and traits that take no type parameters, the subtype relationship simply mirrors the subclass relationship. For classes and traits that take type parameters, however, variance comes into play. For example, class List is covariant in its type parameter, and thus List[Cat] is a subtype of List[Animal] if Cat is a subtype of Animal. These subtype relationships exist even though the class of each of these types is List. By contrast, class Array is nonvariant in its type parameter, so Array[Cat] is not a subtype of Array[Animal]. See here for a more detailed explanation of Scala's variation rules.

Subtyping rules of function types

Suppose that Cat and Dog are subtypes of type Animal. Consider function types A and B defined as follows:
type A := Animal => Dog
type B := Cat => Animal
What should be the subtyping relation between A and B? According to Liskov's Substitution Principle, $T$ is a subtype of $T'$ iff whatever we can do with $T'$, we can do it with $T$. In this example, we can pass a cat to a function of type B and get some animal in return. We can do the same thing to A: if we pass a cat to a function of type A, it returns a dog, which is an animal. Hence, A satisfies the same contract with B (in fact its contract says more than that, as it guarantees that the returned animal is a dog), and thus A is a subtype of B.

In general, if we denote the "subtype of" relation by $<$, then we can express the subtyping rule as $(A \Rightarrow B) < (C \Rightarrow D)$ iff $C < A$ and $B < D$. To see this, note that arguments are something that’s required, whereas return values are something that’s provided. To satisfy Liskov's Substitution Principle, a function subtype must require less and provide more than its supertype. Since a subtype contains more information than its supertype, a reasonable function type should be contravariant in the argument types and covariant in the result type.

Behavioral subtyping

Let A and B be two types. We say A is a behavioral subtype of B iff the set of properties an object of type A satisfies is a superset of the properties satisfied by an object of type B. In terms of Liskov's principle, this means that a property provable by objects of type B should also be provable by objects of type A.

Behavioral subtyping is a stronger notion than classical subtyping rules of function types, which rely merely on type signatures. It is a semantic rather than syntactic relation because it intends to guarantee semantic interoperability of types in a hierarchy. Behavioral subtyping is thus trivially undecidable: if the property to prove is "a function always terminates for arguments of type A", then it is impossible in general for a program (e.g. a compiler) to verify that property holds for any subtype of A.

Behavioral subtyping of function types imposing requirements on signatures just like the classical subtyping rules. In addition, there are a number of behavioral conditions a function subtype must meet:

• Preconditions cannot be strengthened in a subtype.
• Postconditions cannot be weakened in a subtype.
• Invariants of the supertype must be preserved in a subtype.
• History constraint: objects are regarded as being modifiable only through their methods.

A violation of history constraint can be illustrated by defining a mutable point as a subtype of an immutable point. On the other hand, fields added to the subtype may be safely modified because they are not observable through the supertype methods. Thus, one may derive a circle with fixed center but mutable radius from immutable point without violating subtyping rules.

Loop Invariants & Loop Variants

Loop invariant

In Hoare logic, the partial correctness of a while loop is governed by the following rule of inference:$$\frac{\{C\land I\}\;S\;\{I\}} {\{I\}\;\mathbf{while}\; C \; \mathbf{do}\;S\;\{\lnot C\land I\}}$$ The rule above is a deductive step where the premise states that the loop's body does not change $I$ from true to false given the condition $C$, and the conclusion that if the loop terminates then it leads from a state satisfying $I$ to a state satisfying $\lnot C\land I$. The formula $I$ in this rule is known as the loop invariant. In words, a loop invariant

• holds before the loop begins
• establishes the desired result (i.e. the postcondition) given that $\lnot C$ holds
• holds before and after the execution of each iteration
• describes what has been done so far and what remains to be done at each iteration

Finding a useful invariant is often a tricky task requiring sophisticated heuristics. As an illustration, consider the following program that computes the factorial: $$\{x=n \land y=1\}\quad\mathbf{while}\;x\neq 0\;\mathbf{do}\;y:=y\cdot x;\ x:=x-1\quad\{y=n!\}$$ In each iteration, $y$ holds the result so far, $x!$ is what remains to be computed, and $n!$ is the desired results. Also note that $y$ is increasing while $x$ is decreasing. Hence, we can put the decrement of $x$ and the increment of $y$ together to obtain an invariant $x!\cdot y=n!$, as desired.
Now consider a similar program that also computes $n!$: $$\{x=1 \land y=1\}\quad\mathbf{while}\; x<n\; \mathbf{do}\;x:=x+1;\ y:=y\cdot x\quad\{y=n!\}$$This time, both $x$ and $y$ increase as the loop iterates. One can immediately find an invariant $y=x!$ according to rules of thumb. This invariant, however, turns out to be too weak to prove the desired result: we need $y=n!$ after the loop terminates, while the iteration rules only gives $x\ge n \land y=x!$. To fix this, we can strengthen the invariant by conjuncting another invariant $x\le n$ with it. The new invariant $x\le n \land y=x!$ yields $x\ge n\land x\le n \land y=x!$ on termination, which is precisely what we need. (A similar problem would occur in the first example if the guard was $x\ge 1$ instead of $x\neq 0$. Again, this can be solved by strengthening the invariant with $x\ge 0$.)
We refer the reader to this note for more interesting examples of loop invariants.

Loop variant

In Hoare logic, the termination of a while loop can be ensured by the condition $[V=z]\;S\;[V<z]$, where $<$ is a well-founded relation, $V$ is called the loop variant, and $z$ is a unbound symbol that is universally quantified. Combining it with the aforementioned rule of partial correctness, we obtain a rule that expresses the total correctness of the program:$$\frac{[I \land C \land V=z ]\;S\;[I \land V < z]}{[I]\;\mathbf{while}\;C\; \mathbf{do}\;S\;[I\land\lnot C]}.$$ The existence of a variant implies that a while loop terminates. It may seem surprising that the converse is true, as long as we assume the Axiom of Choice: every while loop that terminates (given its invariant) has a variant. To prove this, we need the following theorem:
Theorem. Assuming the Axiom of Choice, a partially ordered set does not possess any infinite descending chains if and only if the corresponding strict order is well-founded.
Assume that the while loop terminates with invariant $I$ and the total correctness assertion $[C\land I]\;S\;[I]$ holds. Consider a "successor" relation $>$ on the state space $\Sigma$ induced by the execution of statement $S$, such that $s>s'$ iff (i) $s$ satisfies $C\land I$, and (ii) $s'$ is the state reached from state $s$ after executing statement $S$. Note that $s\neq s'$, for otherwise the loop would not terminate. Next, we define the "descendant" relation, denoted by $\ge$, as the reflexive and transitive closure of the successor relation. That is, $s\ge s'$ if either $s=s'$, or there exist states $s_1, ..., s_n$ such that $s>s_1>...>s_n>s'$. Note that $\ge$ is antisymmetric, since $s\ge s'$ and $s'\ge s$ implies $s=s'$ under the termination assumption. Therefore, $(\Sigma,\ge)$ is a partially ordered set. Observe that each state has only finitely many distinct descendants, and thus there is no infinite descending chain. Assuming the Axiom of Choice, it follows from above theorem that the successor relation we defined for the loop is well-founded on the state space, since it is strict (irreflexive) and contained in the $\ge$ relation. Thus the identity function on this state space is a variant for the while loop, as we have shown that the state must strictly decrease (to its successor) each time the body $S$ is executed given the invariant $I$ and the condition $C$.

A fixed-point characterization

In predicate transformers semantics, invariant and variant are built by mimicking the Kleene fixed-point theorem. Below, this construction is sketched in set theory. Let $\Sigma$ denote the state space as before. We first define a family $\{A_k\}_{k\ge1}$ of subsets of $\Sigma$ by induction over natural number $k$. Informally, $A_k$ represents the set of initial states that makes the postcondition $P$ satisfied after less than $k$ iterations of the loop: \begin{array}{rcl} A_0 & = & \emptyset \\ A_{k+1} & = & \left\{\ y \in \Sigma\ | \ ((C \Rightarrow wp(S, x \in A_k)) \wedge (\neg C \Rightarrow P))[x \leftarrow y]\ \right\} \\ \end{array}We can then define invariant $I$ as the predicate $\exists k. x \in A_k$ and variant $V$ as the identity function on the state space with a well-founded partial order $<$ such that for $y,z\in\Sigma$, $y<z$ iff $\exists i.\forall j.(y \in A_i \wedge (z \in A_j \Rightarrow i < j))$.

While the theory looks elegant, such an abstract construction cannot be handled efficiently by theorem provers in practice. Hence, loop invariants and variants are often provided by human users, or are inferred by some abstract interpretation procedure.

Deriving the Y Combinator in Javascript

A fixpoint combinator is a higher-order function that computes a fixed point of other functions. Fix-point computation is essential in the theory of programming languages, because it allows the implementation of recursion in the form of rewrite rules, without explicit support from the language's runtime engine. Fixpoint combinators are also used to show that typed lambda calculus is Turing complete, which is an important result in the theory of computation that provides a theoretical foundation for functional programming.

The Y combinator is one of the many fixpoint combinators in lambda calculus. These fixpoint combinators implement recursion without requiring a function to call itself. The requirement, though, is that the language supports anonymous and call-by-name functions. If you want to use this in a strict call-by-value language, you will need in addition the related Z combinator; otherwise the Y combinator derived here will lead to an infinite call stack.

Three derivations of fixpoint combinators are introduced in this post:

• Deriving a fixpoint a via the Fixpoint Theorem
• A step-by-step derivation of the Y combinator
• A step-by-step derivation of the U combinator

Derivation based on the Fixpoint Theorem

Let $D$ be a chain-complete partially ordered set (abbr. ccpo) and $FIX$ be a fixpoint operator for Scott continuous function on $D$. Hence, if $g:D\rightarrow D$ is Scott continuous then $FIX(g)$ is the least upper bound of chain $\{g^n(\perp):n\ge0\}$. The fact that operator $FIX$ is well-defined follows from the Kleene Fixpoint Theorem, see here or [2] for more details. It turns out that the property $g(FIX(g))=FIX(g)$ of a fixpoint suffices to characterize the computation of $FIX$:

function FIX(g) = { return function(x){ g(FIX(g))(x); } }

This function effectively finds a fixpoint of functional $g$ if there exists one. Note the return value is eta-abstracted so as to avoid infinite recursion due to strict evaluation. By exploiting a construct called the U combinator, we can eliminate explicit recursion from the definition of $FIX$ (see here for a step-by-step example of recursion elimination using the U combinator) and obtain

function FIX(g) { 
    return (function(h) { return h(h,g); })
    (function (h,g) { 
        return function(x) { return g(h(h,g))(x); }
    }) 
}

Now fix a function $f$ on non-negative integers. Let $f_m$ be a partial function such that $f_m(n)=f(n)$ for $0\le n \le m$ and $f_m(n)=\perp$ otherwise. Suppose that $g$ is a functional such that $g(f_m)=f_{m+1}$ for all $m\ge0$. Since $g$ is Scott continuous on chain $\{f_m:m\ge0\}$, it follows from the Fixpoint Theorem that $FIX(g)$ exists and $FIX(g)=\bigsqcup \{f_m\}=f$. For example, if $f$ is the factorial function and $g$ is

function g(f) { return function(n) { return n<=1 ? 1 : n*f(n-1); } }

then $FIX(g)(n)=f(n)=n!$ for $n=0,1,2,...$.

An ad hoc derivation of the Y combinator

In the following, I shall start from an ordinary recursive factorial function, apply transformations over the initial function, eliminate the recursive structure step by step, and derive the Y combinator in the end. A recursive factorial function is

function fact (n) {
    if(n<2) return 1; else return fact(n-1)*n;
}

Now consider a function, denoted by makeFact, that takes a factorial function as parameter and generates another factorial function that "unwind" the computation of factorial for one step:

function makeFact (fact) {
    return function(n) {
        if(n<2) return 1; else return fact(n-1)*n;
    }
}

Here we have makeFact(fact)(n) = fact(n). It seems that we cannot use makeFact to generate factorial without already having a factorial function at hand. However, it turns out that makeFact does not really need a factorial function. Given that makeFact takes a factorial function and returns a factorial function, the following function will take makeFact as input and generates a factorial function:

function createFact (makeFact) {
    var fact = function(n) { return makeFact(fact)(n); }
    return fact;
}

Instead of using a factorial function, createFact uses a function that generates a factorial function to generate another factorial function, ie., createFact(makeFact)(n) = fact(n). Thus it avoids the need of a predefined factorial function. However, the definition of variable fact (in the function body of makeFact) is recursive, since it needs to call itself to compute the return value. Such a reference can be avoided if we define another function that makes the reference for it:

function createFact (makeFact) {
    var genFact = function() {
        var fact = function(n) { return genFact()(n); };
        return makeFact(fact);
    };
    return genFact();
}

Now, the definition of fact is not recursive, but we have a new function genFact that refers to itself. This reference can be made via a parameter, rather than calling the function directly:

function createFact (makeFact) {
    var genFact = function(genFactRef) {
        var fact = function(n) { return genFactRef(genFactRef)(n); };
        return makeFact(fact);
    };
    return genFact(genFact);
}

Now, all the functions (ie. genFact and fact) declared in the function body of createFact are non-recursive. We proceed to show that the declarations can be further eliminated and the result is a version of the famous Y combinator. First, we avoid declaring fact by making it an anonymous function:

function createFact (makeFact) {
    var genFact = function(genFactRef) {
        return makeFact(function(n) { return genFactRef(genFactRef)(n) });
    };
    return genFact(genFact);
}

Also, we use an anonymous helper function function(f){ return f(f) } to return genFact(genFact), avoiding the need of declaring genFact:

function createFact (makeFact) {
    return function(f){ return f(f); }(function(genFactRef){
        return makeFact(function(n){ return genFactRef(genFactRef)(n) });
    });
}

Now we reach a function which has no declarations, only parameters. By renaming the parameters to general names, we call this function the Y combinator (for single-parameter functions):

function Y (rec) {
    return function(f){ return f(f) }(function(f) {
        return rec(function(x) { return f(f)(x) });
    });
} // In untyped lambda calculus: Y = λf.(λx.x x) (λx.f (x x))

Take a close look at the definition of Y. It uses only three kinds of expression: anonymous functions, variable reference and function application. It also makes no explicit reference to an outside variable or to itself. The Y combinator allows a concise transformation from a recursive function to a non-recursive function. For example, we can use the Y combinator function to define non-recursive factorial function and a Fibonacci function. To compute factorial,

Y(function(f) { return function(n) { return n<2 ? 1 : n*f(n-1); }})

To compute Fibonacci number,

Y(function(f) { return function(n) { return n<=2 ? 1 : f(n-1)*f(n-2); }})

Finally, if you need a Y combinator for functions that takes more than one argument, you can use Javascript's built-in argument object to obtain all the arguments indefinitely:

function Y(rec) {
    return function(f) { return f(f) }(function(f) {
        return rec(function() { 
            return f(f).apply(null, arguments) 
        });
    });
}

An ad hoc derivation of the U combinator

The Y combinator it not the only way to eliminate recursive calls from a program. The U combinator also does the job, and it is much simpler than the Y combinator. The derivation of the U combinator is also more intuitive and straightforward, at least for those who are comfortable with currying. The ease of derivation, however, comes at the cost that the input function of the U combinator has to make self-application explicit, which may make the user program (which is the factorial function here) look more artificial.

First note that a reference to function itself can be provided to the function as a parameter, avoiding the self reference:

function fact0(f,n) { 
    return n<2 ? 1 : n*f(f,n-1);
}

Hence we have fact0(fact0, n) = fact(n). If we define a function U as

function U(f,n) { return f(f,n); }

then U(fact0, n) yields fact(n). Here we reach a function with two inputs, but we are looking for a combinator that should take only one input. Can we work around this problem? Well a smart guy by the name of Haskell Curry has a neat trick: if you have good higher-order functions then you only need functions of one variable. The proof is that you can get from functions of two (or more in the general case) variables to one variable with a purely mechanical text transformation called currying. For example, by currying f(x,y) we can obtain a function g such that g(x)(y) = f(x,y):

var f = function(x,y) { ... } // original
var g = function(x) { return function(y) { return f(x,y) } } // curried

Now we just apply this transformation in the right places to get our final version of the U combinator.

function fact0(f,n) { // fact0(fact,n) = fact(n)
    return n<2 ? 1 : n*f(n-1);
}
function fact0(f) { // curried, now fact0(fact0)(n) = fact(n)
    return function(n) { return n<2 ? 1 : n*f(f)(n-1); } 
}
function U(f,n) { // U(fact0,n) = fact(n)
    return f(f)(n);
}
function U(f) { // curried, reaching the famous U combinator λf.(f f)
    return function(n) { return f(f)(n); }
}
function U(f) { // allowing an arbitrary number of parameters
    return function() {
        return f(f).apply(null, arguments); 
    }
}

One can check that U(fact0)(n)=fact(n), where U(fact0) can be expanded to

U(function(f) { return function(n) { return n<2 ? 1 : n*f(f)(n-1); }})

References

1. Many faces of the fixed-point combinator: This web page also describes in detail the fix-point combinators for mutual recursion, called poly-variadic fixpoint combinators. It derives the simplest combinator I've seen so far.
2. Semantics with Applications: An Appetizer, Section 5.2.
3. Fixed-point combinators in JavaScript: Memoizing recursive functions
4. Also check Y combinators on Rosetta Code, and many other fixed-point combinators on Wiki.

Tossing coins

Simulate a biased coin with a fair coin

The following algorithm simulates a biased coin with head probability 0<p<1, using a random number generator toss_coin that returns 0 and 1 with equal probability.

// Pre-condition: p is a rational number and 0 < p < 1
while (true) {
    var p = p * 2,
    var p2 = p - floor(p);
    var msb = p - p2;
    var random_bit = toss_coin(); // 0 or 1
    if (random_bit < msb) {
        print "head";
        break;
    }
    if (random_bit > msb) {
        print "tail";
        break;
    }
    p = p2;
}
// Post-condition: "head" is outputted with probability p
//                 "tail" is outputted with probability 1-p

To see why it works, consider a random number 0 < r < 1 generated by enumerating its binary representation with a fair coin. Observe that given 0 < p < 1, the probability that 0 < r < p is exactly p. The algorithm then enumerates p and r, digit by digit, and determines the return value according to the following fact: r < p if and only if r_i < p_i for some i and r_j = p_j for all j < i, where x_k denotes the kth digit in the binary representation of x.

The following code fragment is adapted from the book Abstraction, Refinement and Proof for Probabilistic Systems, A McIver, C Morgan, 2005. Succinct as it is, its behavior is not as straightforward to understand as the above one and we skip the analysis of its correctness.

// Pre-condition: p is a rational number and 0 < p < 1
var x = p;
while (toss_coin() == 1) {
    x = 2 * x;
    if(x > 1){
      x = x - 1;
    }
}
// Post-condition: x > 0.5 with probability p
print x > 0.5 ? "head" : "tail";

Simulate a fair coin with a biased coin

while (true) {
    var a = toss_biased_coin();
    var b = toss_biased_coin();
    if(a == b){
      print a == 1 ? "head" : "tail";
      break;
    }
}

The expected running time is $\frac{2}{1-2p(1-p)}$.

Checking whether a coin is fair

First, note that it is impossible to rule out arbitrarily small deviations from fairness such as might be expected to affect only one flip in a lifetime of flipping. Also, it is always possible for a biased coin to happen to turn up exactly 10 heads in 20 flips. As such, any fairness test must only establish a certain degree of confidence in a certain degree of fairness (a certain maximum bias). In more rigorous terminology, this problem can be viewed as a special case of determining the parameters of a Bernoulli process, given only a limited sample of Bernoulli trials.
Suppose that there are $h$ heads and $t$ tails in an experiment of $n=h+t$ tosses. The best estimator for the actual head probability $p$ is $$q=\frac{h}{h+t}.$$ One can associate this estimator with a margin of error $\epsilon$ at a particular confidence level, which is given by the Z-value of a standard normal distribution in the following sense: $$P\{|p-q|\le z\sigma\}=P\{|Z|\le z\},$$ where $\sigma$ is the standard deviation of error of the estimate and $Z$ is a standard-normal distributed random variable. For $n$ Bernoulli trials, the standard deviation is $$\sigma=\sqrt{\frac{p(1-p)}{n}}\le\frac{1}{2\sqrt n}.$$ Hence, suppose that $P\{|Z|\le z\}=c$ then we have $$P\{|p-q|\le \frac{z}{2\sqrt n}\}\ge c.$$ For example, by looking up the table we see that $z=2.5759$ for $c=99\%$. Suppose now $h=54$ and $t=46$, then the actual head probability $p$ falls into $0.54\pm 0.12879$ at 99% level of confidence.