Fundamental Algorithms

I'm starting a new series (will see how it 'll go:)) about some very basic algorithms implemented in C++. The first, order statistics . I use classic divide and conquer strategy, with random pivot (it's not the best, but reasonable choice, I think it gives nearly amortised linear time). As a bonus, there is quicksort algorithm - uses the same partition function. Code (partition and order statistic: quick_select here):

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35

int partition(int a[], int p, int r) { int t = a[r]; int i = p - 1; for (int j = p; j < r;j++){ if (a[j] <= t) { i += 1; std::swap(a[i], a[j]); } } std::swap(a[i + 1], a[r]); return i + 1; } int quick_select(int a[], int s, int n, int k) { int r = partition_random(a, s, n); if (r - s == k - 1) {return a[r];} else { if (k - 1 < r - s) { return quick_select(a, s, r - 1, k); } else { return quick_select(a, r + 1 , n, k - r + s - 1); } } } int partition_random(int a[], int p, int r) { std::random_device rd; std::mt19937 rng(rd()); std::uniform_int_distribution<int> uni(p, r); auto rn = uni(rng); std::swap(a[r], a[rn]); return partition(a, p, r); }

As seen, quick_select use partition to split initial array, according to some random index; and then recurs in certain part, depends if k - the #nth smallest is top or bottom part of the array.
Partition algorithm is commonly known, explanation can be found on the web; it returns number of elements less than the last element of the array, and mutates array in place (moving the all elements less then pivot, to the left). It could be use to compute, median:

1 2 3 4 5 6 7 8

float median(int a [], int n) { if (n % 2 == 1) return quick_select(a, 0, n - 1, n / 2 + 1); int l = n / 2; int r = n / 2 + 1; return (1.0 * (quick_select(a, 0, n - 1, l) + quick_select(a, 0, n - 1, r))) / 2; }

In C++ Standard Library, order statistics is named: nth_element. Code on github .

Thank you.

Kleisli Category by Example

My New Year resolution: studying category theory, I've started softly, reading this. From the book, there is an exercise, which is the subject of this post. I don't want to go into the details, anyone can read (strongly recommended), but shortly: 
A category is a structure consists of objects (for as they are primitives in the theory) and morphisms (not really are, but we can visualize them as functions) between them. And, for every object a there must be the identity morphism going from a to a; for every pair of morphisms going from a to b and from b to c, there must be the morphism, going from a to c which is a composition of the two morphisms and is associative. 

Thinking about programming languages, we can have category types and functions, where objects are types, functions are morphisms between types. Kleisli Category is such type and is constructed as follows: Let's say we have a category of types and functions, like Integer, Double, and morphisms, for example, int => int and so on. Let's say, that we want to make every such function (ex.: square root) total, it can be done by making the return value type of the function something like a pair: Double and String value: "Success" or "Fail", even better Boolean: True, False. So now, we have functions from int to pair. If we can make a category from these new objects and functions we have Kliesli Category. 

The language is C++, morphisms in category: T => optional<T, bool> - no surprise here:). Optional type class:

template<class A> class optional {
 bool _isValid;
 A _value;
 
 public:
 optional(): _isValid(false) {}
 optional(A e) : _isValid(true), _value(e) {}
 optional(A e, bool v) : _isValid(v), _value(e) {}
 A value() const { return _value; }
 bool validation() const {return _isValid;}
};

Pair with value and Boolean indicator, constructors, methods to get values. 
We deal with square root and reciprocal:

optional<double> safe_root(double x) {
 if (x >= 0) return optional<double>{std::sqrt(x)};
 else return optional<double>{};
}

optional<double> safe_reciprocal(double x) {
 if (! x == 0) return optional<double> {1 / x};
 else return optional<double>{};
}

As definition states, the job is to construct identity and composition. Identity seems easy:

optional<double> identity(double x) {
 return optional<double>{x};
}

If original identity would go from Double to Double, then a new must go from Double to optional and change the state of isValid to True. There is no logical alternative: True is neutral to conjuction (why conjuction in a moment) and if we don't change a value(the very heart of identity operation), then there still is a value, so indicator must be true.   

Composition, first we need a compose function (takes f and g and returns g○f - "g after f"), to do that, here is a bit of new C++ lambda syntax (thanks to Eric Niebler):

auto const compose = [](auto m1, auto m2) {
 return [m1, m2](auto x) {
 auto p1 = m1(x);
 auto p2 = m2(p1.value());
 return optional<double>(p2.value(), 
                    p1.validation() && p2.validation());
 };
};

Inside the lambda, p1 and p2 are responsible for the composition of two incoming functions. The value after the second evaluation and the logical conjunction of the two incoming optionals isValid fields is passed to the returned function. Is this really a composition? Well, must be: functions are composable, so check, the second isValid parts are composed by logical and operator which is, I hope, clear. To get a composition with a proper type, both functions validation fields must evaluate to true, if one or both are false, then the composition fails (in the term of a value type evaluation), setting validation to false. And it's also, associative - composition is associative.
All works as expected: 

auto const safe_root_reciprocal = compose(safe_reciprocal, safe_root);
auto const safe_rootIdR = compose(safe_root, identity);
auto const safe_rootIdL = compose(identity, safe_root);
std::cout << safe_root_reciprocal(0).validation() << "\n"; // -> 0
std::cout << safe_root_reciprocal(3).validation() << "\n"; // -> 1
std::cout << safe_root_reciprocal(3).value() << "\n"; // -> 0.57735
std::cout << (safe_rootIdR(2).value()) << "\n"; // -> 1.41421
std::cout << (safe_rootIdL(2).value()) << "\n"; // -> 1.41421

Left and right identities gave expected results, also composition safe_root_reciprocal does the job. Looks like it's done, one more time, Happy New Year!