Archive for the ‘Mathematics’ Category

Digraphs in Haskell

April 3, 2012

A directed graph (or digraph) is an ordered pair D = (V, A), where V is a non-empty finite set of elements called vertices, and A is a finite set of elements called arcs, where an arc is an ordered pair of distinct vertices [1]. We call V and A the vertex set and the arc set, respectively. Note that A \subseteq V \times V. For an arc (v_i, v_j), we call v_i the tail, and v_j the head, and say that the arc (v_i, v_j) leaves vertex v_i and enters vertex v_j.

For instance, consider the digraph depicted below

This digraph has 10 vertices and 15 arcs. How would we build this digraph in Haskell? A vertex could be represented by an integer, and an arc could be represented by an ordered pair of integers. However, in Mathematics we work with sets, whereas in Haskell we work with lists. Hence, we could represent the vertex set by a list of vertices, and the arc set by a list of edges.

Here is a simplistic implementation:

-- type synonyms
type Vertex  = Int
type Arc     = (Vertex,Vertex)
type Digraph = ([Vertex],[Arc])

-- create list of vertices
vs :: [Vertex]
vs = [1..10]

-- create list of arcs
as :: [Arc]
as = [(1,4),(1,5),(1,6),(2,4),(3,4),
      (4,5),(5,6),(7,6),(8,1),(8,2),
      (8,7),(9,8),(9,10),(10,2),(10,3)]

-- create digraph
digraph :: Digraph
digraph = (vs,as)

We load this script in GHCi and perform some simple testing:

*Main> -- list of vertices
*Main> vs
[1,2,3,4,5,6,7,8,9,10]
*Main> length vs
10
*Main> -- list of arcs
*Main> as
[(1,4),(1,5),(1,6),(2,4),(3,4),(4,5),(5,6),(7,6),
(8,1),(8,2),(8,7),(9,8),(9,10),(10,2),(10,3)]
*Main> length as
15
*Main> -- digraph
*Main> digraph
([1,2,3,4,5,6,7,8,9,10],[(1,4),(1,5),(1,6),
(2,4),(3,4),(4,5),(5,6),(7,6),(8,1),(8,2),
(8,7),(9,8),(9,10),(10,2),(10,3)])

We have thus created a digraph. What can we do with it? Is there anything of interest that we could compute?

__________

Neighborhoods

Given a digraph D = (V, A), we call the following sets

N_D^{+} (v) = \{ u \in V \mid (v, u) \in A\}

N_D^{-} (v) = \{ u \in V \mid (u, v) \in A\}

the out-neighborhood of vertex v and the in-neighborhood of vertex v, respectively [1]. The vertices in N_D^{+} (v) and N_D^{-} (v) are called the out-neighbors and the in-neighbors of vertex v, respectively.

In Haskell, we compute the out- and in-neighborhoods as follows:

-- compute out-neighborbood of a vertex
outNeighbors :: Digraph -> Vertex -> [Vertex]
outNeighbors (vs,as) v = map snd (filter (p v) as)
                         where p v (t,h) | t == v = True
                                         | t /= v = False

-- compute in-neighborbood of a vertex
inNeighbors :: Digraph -> Vertex -> [Vertex]
inNeighbors (vs,as) v = map fst (filter (p v) as)
                        where p v (t,h) | h == v = True
                                        | h /= v = False

Let us do some testing in GHCi:

*Main> -- list of all out-neighborhoods
*Main> map (outNeighbors digraph) vs
[[4,5,6],[4],[4],[5],[6],[],[6],[1,2,7],[8,10],[2,3]]
*Main> -- list of all in-neighborhoods
*Main> map (inNeighbors digraph) vs
[[8],[8,10],[10],[1,2,3],[1,4],[1,5,7],[8],[9],[],[9]]

 __________

Degrees

Given a digraph D = (V, A), we call the following cardinalities

d_D^{+} (v) = |N_D^{+} (v)|

d_D^{-} (v) = |N_D^{-} (v)|

the out-degree and the in-degree of vertex v, respectively [1]. We compute the out- and in-degrees in Haskell as follows:

-- compute out-degree of a vertex
outDegree :: Digraph -> Vertex -> Int
outDegree (vs,as) v = length (outNeighbors (vs,as) v)

-- compute in-degree of a vertex
inDegree :: Digraph -> Vertex -> Int
inDegree (vs,as) v = length (inNeighbors (vs,as) v)

Let us do some testing in GHCi:

*Main> -- out-degrees of all vertices
*Main> map (outDegree digraph) vs
[3,1,1,1,1,0,1,3,2,2]
*Main> -- in-degrees of all vertices
*Main> map (inDegree digraph) vs
[1,2,1,3,2,3,1,1,0,1]

__________

Concluding remarks

Here is the whole Haskell script:

-- type synonyms
type Vertex  = Int
type Arc     = (Vertex,Vertex)
type Digraph = ([Vertex],[Arc])

-- create list of vertices
vs :: [Vertex]
vs = [1..10]

-- create list of arcs
as :: [Arc]
as = [(1,4),(1,5),(1,6),(2,4),(3,4),
      (4,5),(5,6),(7,6),(8,1),(8,2),
      (8,7),(9,8),(9,10),(10,2),(10,3)]

-- create digraph
digraph :: Digraph
digraph = (vs,as)

-- compute out-neighborbood of a vertex
outNeighbors :: Digraph -> Vertex -> [Vertex]
outNeighbors (vs,as) v = map snd (filter (p v) as)
                         where p v (t,h) | t == v = True
                                         | t /= v = False

-- compute in-neighborbood of a vertex
inNeighbors :: Digraph -> Vertex -> [Vertex]
inNeighbors (vs,as) v = map fst (filter (p v) as)
                        where p v (t,h) | h == v = True
                                        | h /= v = False

-- compute out-degree of a vertex
outDegree :: Digraph -> Vertex -> Int
outDegree (vs,as) v = length (outNeighbors (vs,as) v)

-- compute in-degree of a vertex
inDegree :: Digraph -> Vertex -> Int
inDegree (vs,as) v = length (inNeighbors (vs,as) v)

We now know how to create digraphs in Haskell and compute simple things, like neighborhoods and degrees. What else would it be of interest to compute?

__________

References

[1] Jørgen Bang-Jensen, Gregory Z. Gutin, Digraphs: Theory, Algorithms and Applications (2nd edition), Springer, 2009.

Newton’s method using rational arithmetic

April 2, 2012

The floating-point number system puts further restraints on what can be expected. We can at best hope to keep the relative error under control, and we cannot expect to find zeros far from the origin with great absolute accuracy.

Richard W. Hamming [1]

I am somewhat disenchanted with floating-point arithmetic. We all have implemented numerical methods using floating-point numbers. But, how many of us have implemented numerical methods using rational numbers of arbitrary precision? I certainly have not. Hence, let us give it a try. In this post I will implement Newton’s method [1] (also known as the Newton-Raphson method) in Haskell using arbitrary-precision rational numbers.

Given a continuous (real) function f : \mathbb{R} \to \mathbb{R} and x_0, an initial approximation of the zero of function f, we then obtain the next approximation using the following recurrence relation

x_{k+1} = x_k - \displaystyle\frac{f (x_k)}{f' (x_k)}

where f' is the first derivative of f. We will use arbitrary-precision rational numbers (represented as a ratio of two Integer values)

type Rational = Ratio Integer

which are defined in the Data.Ratio library.

__________

Example

Suppose we would like to find a rational approximation of \sqrt{2}. We create a real function defined by f (x) = x^2 - 2, whose zeros are  \sqrt{2} and -\sqrt{2}. Note that the first derivative of f is f' (x) = 2 x. The initial approximation is, say, x_0 = 1 / 1, which is closer to \sqrt{2} than to -\sqrt{2}. Hopefully, the sequence \{x_k\} will converge to \sqrt{2} instead of converging to -\sqrt{2}. We are performing a numerical experiment, and success is not guaranteed!

The following Haskell script implements Newton’s method:

import Data.Ratio

-- define function f
f :: Rational -> Rational
f x = x*x - 2

-- define 1st derivative of function f
f' :: Rational -> Rational
f' x = 2*x

-- define Newton map
g :: Rational -> Rational
g x = x - (f x / f' x)

-- initial approximation
x0 :: Rational
x0 = 1 % 1

-- list of rational approximations
xs :: [Rational]
xs = iterate g x0

where we used (higher-order) function iterate to create the list of rational approximations. We run the script above in GHCi and use take to compute the first few rational approximations in the list:

*Main> take 6 xs
[1 % 1,3 % 2,17 % 12,577 % 408,665857 % 470832,
886731088897 % 627013566048]

Hence, we have, for instance, the following rational approximation

x_5 = \displaystyle\frac{886731088897}{627013566048} \approx \sqrt{2}.

How good is this rational approximation? Let us see:

*Main> let x5 = xs !! 5
*Main> x5*x5 - 2
1 % 393146012008229658338304
*Main> 886731088897 / 627013566048
1.4142135623730951
*Main> sqrt 2
1.4142135623730951

Fairly good indeed! And all it took was five iterations!

__________

References

[1] Richard W. Hamming, Numerical Methods for Scientists and Engineers, Dover Publications, 1973.

Deciding the surjectivity of finite functions

March 19, 2012

Consider a function f : \mathcal{X} \to \mathcal{Y}. We say that f is a finite function if and only if both \mathcal{X} and \mathcal{Y} are finite sets. As mentioned two weeks ago, function f is surjective if and only if

\forall y \exists x \left( f (x) = y \right)

where x and y range over sets \mathcal{X} and \mathcal{Y}, respectively. We now introduce the surjectivity predicate \varphi : \mathcal{X} \times \mathcal{Y} \to \{\text{True}, \text{False}\}, defined by \varphi (x, y) =\left( f (x) = y\right) so that function f is surjective if and only if \forall y \exists x \, \varphi (x, y) evaluates to \text{True}.

Let m := |\mathcal{X}| and n := |\mathcal{Y}| denote the cardinalities of sets \mathcal{X} and \mathcal{Y}, respectively. Deciding \forall y \exists x \, \varphi (x, y) should then require a total of m n evaluations of predicate \varphi. Let us write \mathcal{X} = \{x_1, x_2, \dots, x_m\} and \mathcal{Y} = \{y_1, y_2, \dots, y_n\}, and introduce \Phi \in \{\text{True}, \text{False}\}^{m \times n}, a Boolean matrix defined as follows

\Phi = \left[\begin{array}{cccc} \varphi (x_1,y_1) & \varphi (x_1,y_2) & \ldots & \varphi (x_1,y_n)\\ \varphi (x_2,y_1) & \varphi (x_2,y_2) & \ldots & \varphi (x_2,y_n)\\ \vdots & \vdots & \ddots & \vdots\\ \varphi (x_m,y_1) & \varphi (x_m,y_2) & \ldots & \varphi (x_m,y_n)\\\end{array}\right]

Stating that f is surjective, i.e., \forall y \exists x \, \varphi (x, y) evaluates to \text{True}, is equivalent to saying that every column of matrix \Phi contains at least one entry that evaluates to \text{True}. We now introduce a new predicate, \psi : \mathcal{Y} \to \{\text{True}, \text{False}\}, defined by

\psi (y) = \exists x \, \varphi (x,y)

so that \forall y \exists x \, \varphi (x, y) can be rewritten in the form \forall y \, \psi (y), where y ranges over set \mathcal{Y}. Note that \psi (y) can be viewed as the disjunction of the m atoms \varphi (x_1, y), \varphi (x_2, y), \dots, \varphi (x_m, y), i.e.,

\psi (y) = \displaystyle\bigvee_{i =1}^m \varphi (x_i, y).

Note also that \forall y \, \psi (y) can be viewed as the conjunction of the n atoms \psi (y_1), \psi (y_2), \dots, \psi (y_n). Thus, we can decide surjectivity using the following procedure:

  1. Compute \psi (y) = \exists x \, \varphi (x,y) for all y \in \mathcal{Y}.
  2. Compute the conjunction \bigwedge_{j=1}^n \psi (y_j). If its truth value is \text{True}, then the function is surjective.

Let us now consider two simple examples and implement this decision procedure in Haskell.

__________

Example #1

Let us define \mathbb{Z}_4 := \{0, 1, 2, 3\}. Consider the finite function

\begin{array}{rl} f : \mathbb{Z}_4 &\to \mathbb{Z}_4\\ n &\mapsto 3 n \pmod{4}\end{array}

where \pmod{4} denotes modulo 4 arithmetic. Note that f (0) = 0, f (1) = 3, f (2) = 2, and f (3) = 1. Hence, we easily conclude that f is both injective and surjective (and, thus, is bijective). In Haskell:

Prelude> let f n = (3 * n) `mod` 4
Prelude> let xs = [0..3]
Prelude> let ys = [0..3]
Prelude> let psi y = or [(f x == y) | x <- xs]
Prelude> all psi ys
True

which allows us to conclude that f is surjective. Note that we used function or to perform disjunction, and function all to perform universal quantification.

__________

Example #2

Consider now the following finite function

\begin{array}{rl} g : \mathbb{Z}_4 &\to \mathbb{Z}_4\\ n &\mapsto 2 n \pmod{4}\end{array}

Note that g (0) = g (2) = 0 and g (1) = g (3) = 2. Hence, function g is neither injective nor surjective. In Haskell, we have the following:

Prelude> let g n = (2 * n) `mod` 4
Prelude> let xs = [0..3]
Prelude> let ys = [0..3]
Prelude> let psi y = or [(g x == y) | x <- xs]
Prelude> all psi ys
False

which allows us to conclude that g is not surjective. Why? The following code returns a list of lists of Booleans where each sublist is a column of matrix \Phi (where each entry is \Phi_{ij} = (g(x_i) = y_j))

Prelude> [[(g x == y) | x <- xs] | y <- ys]
[[True,False,True,False],[False,False,False,False],
[False,True,False,True],[False,False,False,False]]

Note that the 2nd and 4th sublists contain only \text{False} truth values. Hence, it is not the case that every column of matrix \Phi contains at least one entry that evaluates to \text{True}. Thus, g is not surjective.

Manin on abbreviated notation

March 17, 2012

Yuri Manin (Юрий Манин) on abbreviated notation:

If written down, most of the interesting expressions and texts in a formal language either would be physically extremely long, or else would be psychologically difficult to decipher and learn in an acceptable amount of time, or both.

They are therefore replaced by “abbreviated notation” (which can sometimes turn out to be physically longer). The expression “x x x x x x” can be briefly written “x \dots x (six times)” or “x^6.” The expression “\forall z \left( z \in x \Leftrightarrow z \in y\right)” can be briefly written “x = y.” Abbreviated notation can also be a way of denoting any expression of a definite type, not only a single such expression; (any expression 101010 \dots 10 can be briefly written “the sequence of length 2 n with ones in odd places and zeros in even places” or “the binary expansion of \frac{2}{3} (4^n - 1).”)

Ever since our tradition started with Vieta, Descartes, and Leibniz, abbreviated notation has served as an inexhaustible source of inspiration and errors. There is no sense in, or possibility of, trying to systematize its devices; they bear the indelible imprint of the fashion and spirit of the times, the artistry and pedantry of the authors. The symbols \sum, \int, \in are classical models worthy of imitation. Frege’s notation, now forgotten, for P and Q (…) shows what should be avoided.

Translated from the original Russian by Neal Koblitz.

__________

Source:

Yuri Manin, A Course in Mathematical Logic, Springer, 1977.

Deciding the injectivity of finite functions II

March 17, 2012

Last week we considered the following finite function

\begin{array}{rl} f : \mathbb{Z}_4 &\to \mathbb{Z}_4\\ n &\mapsto 3 n \pmod{4}\end{array}

where \mathbb{Z}_4 := \{0, 1, 2, 3\} and \pmod{4} denotes modulo 4 arithmetic. We know that function f is injective if and only if

\forall x_1 \forall x_2 \left( f (x_1) = f (x_2) \implies x_1 = x_2\right)

where x_1 and x_2 range over set \mathbb{Z}_4. We introduce the injectivity predicate \varphi : \mathbb{Z}_4 \times \mathbb{Z}_4 \to \{\text{True}, \text{False}\}, defined by

\varphi (x_1, x_2) =\left( f (x_1) = f (x_2) \implies x_1 = x_2\right)

so that function f is injective if and only if the universally-quantified sentence \forall x_1 \forall x_2 \, \varphi (x_1, x_2) evaluates to \text{True}.

__________

A naive decision procedure

We can determine the truth value of \forall x_1 \forall x_2 \, \varphi (x_1, x_2) by evaluating the predicate \varphi at all (x_1, x_2) \in \mathbb{Z}_4 \times \mathbb{Z}_4, which will require a total of |\mathbb{Z}_4|^2 = 4^2 = 16 evaluations. We implemented this decision procedure in Haskell last week. Here’s a rather brief GHCi session:

Prelude> let f n = (3 * n) `mod` 4
Prelude> let phi (x1,x2) = (f x1 /= f x2) || (x1 == x2)
Prelude> all phi [(x1,x2) | x1 <- [0..3], x2 <- [0..3]]
True

which allows us to conclude that f is injective. Note that we used the following equivalence

\left( f (x_1) = f (x_2) \implies x_1 = x_2\right) \equiv \left( f (x_1) \neq f (x_2) \lor x_1 = x_2\right)

to implement \varphi. Can we do better than this? In fact, we can.

__________

A better decision procedure

We introduce matrix \Phi \in \{\text{True}, \text{False}\}^{4 \times 4}, defined by

\Phi = \left[\begin{array}{cccc} \varphi (0,0) & \varphi (0,1) & \varphi (0,2) & \varphi (0,3)\\ \varphi (1,0) & \varphi (1,1) & \varphi (1,2) & \varphi (1,3)\\ \varphi (2,0) & \varphi (2,1) & \varphi (2,2) & \varphi (2,3)\\ \varphi (3,0) & \varphi (3,1) & \varphi (3,2) & \varphi (3,3)\\\end{array}\right]

Stating that f is injective, i.e., \forall x_1 \forall x_2 \, \varphi (x_1, x_2) evaluates to \text{True}, is equivalent to saying that all the 4^2 = 16 entries of matrix \Phi evaluate to \text{True}. Note, however, that the four entries on the main diagonal of \Phi will always evaluate to \text{True}, as the implication

f (x) = f (x) \implies x = x

evaluates to \text{True} for every x. Moreover, from the equivalence

\left( f (x_1) = f (x_2) \implies x_1 = x_2\right) \equiv \left( f (x_2) = f (x_1) \implies x_2 = x_1\right)

we can conclude that \varphi (x_1, x_2) = \varphi (x_2, x_1) for all x_1 and x_2, which tells us that matrix \Phi is symmetric (i.e., \Phi = \Phi^T). Therefore, we do not need to evaluate all entries of matrix \Phi, only the ones above the main diagonal, as illustrated below

\left[\begin{array}{cccc} \ast & \varphi (0,1) & \varphi (0,2) & \varphi (0,3)\\ \ast & \ast & \varphi (1,2) & \varphi (1,3)\\ \ast & \ast & \ast & \varphi (2,3)\\ \ast & \ast & \ast & \ast\\\end{array}\right]

where the symbol \ast denotes “don’t care”, i.e., entries we do not need to evaluate. Hence, instead of evaluating the predicate \varphi a total of 4^2 = 16 times, we now only need to evaluate it {4 \choose 2} = 6 times, which is 37.5 \% of the original total. We can generate all pairs (x_1, x_2) with x_2 > x_1 using the following Haskell code:

Prelude> [(x1,x2) | x1 <- [0..3], x2 <- [(x1+1)..3]]
[(0,1),(0,2),(0,3),(1,2),(1,3),(2,3)]

We thus implement a more efficient decision procedure as follows:

Prelude> let f n = (3 * n) `mod` 4
Prelude> let phi (x1,x2) = (f x1 /= f x2) || (x1 == x2)
Prelude> all phi [(x1,x2) | x1 <- [0..3], x2 <- [(x1+1)..3]]
True

which requires less than half of the number of evaluations required by the previous (naive) procedure.

__________

Conclusions

Given a finite function f : \mathcal{X} \to \mathcal{Y}, we introduce the injectivity predicate \varphi : \mathcal{X} \times \mathcal{X} \to \{\text{True}, \text{False}\}, defined by

\varphi (x_1, x_2) =\left( f (x_1) = f (x_2) \implies x_1 = x_2\right)

so that (finite) function f is injective if and only if \forall x_1 \forall x_2 \, \varphi (x_1, x_2) evaluates to \text{True}. Let n := |\mathcal{X}| denote the cardinality of set \mathcal{X}. Deciding \forall x_1 \forall x_2 \, \varphi (x_1, x_2) using the naive procedure will require a total of n^2  evaluations of the predicate \varphi. Since \varphi (x,x) = \text{True} for every x, and taking into account that \varphi (x_1,x_2) = \varphi (x_2,x_1), using the improved procedure we can decide injectivity at a cost of only {n \choose 2} = \frac{n (n-1)}{2} evaluations of the predicate \varphi. As n approaches infinity, the ratio {n \choose 2} / n^2 approaches 1/2. Hence, for “large” problems, the improved procedure will require approximately half the number of evaluations required by the naive procedure.


Follow

Get every new post delivered to your Inbox.

Join 47 other followers