diff --git a/alib2algo/src/graph/embedding/HopcroftTarjan.cpp b/alib2algo/src/graph/embedding/HopcroftTarjan.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..ab61a0d0f1f032cc22e967c2d4c6aad1ae457fe4
--- /dev/null
+++ b/alib2algo/src/graph/embedding/HopcroftTarjan.cpp
@@ -0,0 +1,675 @@
+#include "HopcroftTarjan.h"
+
+#include <list>
+#include <stack>
+#include <algorithm>
+
+#include <exception/AlibException.h>
+
+namespace graph
+{
+
+namespace embedding
+{
+
+// Helpers
+
+using node = Node;
+using edge = DirectedEdge;
+
+template <typename T>
+using edge_array = std::unordered_map<edge, T>;
+
+template <typename T>
+using node_array = std::unordered_map<node, T>;
+
+template<typename T>
+inline void conc(T &container, T &l)
+{
+	container.splice(container.end(), l);
+}
+
+class GraphWrapper {
+public:
+	explicit GraphWrapper(DirectedGraph *g) : m_g(g) { }
+
+	std::vector<DirectedEdge> adj_edges(const node &n) const
+	{
+		std::vector<DirectedEdge> out;
+		const auto &edges = m_g->neighborEdges(n);
+		std::copy(edges.begin(), edges.end(), std::back_inserter(out));
+		sort(out);
+		return out;
+	}
+
+	std::vector<DirectedEdge> all_edges() const
+	{
+		std::vector<DirectedEdge> out;
+		const auto &edges = m_g->getEdges();
+		std::copy(edges.begin(), edges.end(), std::back_inserter(out));
+		sort(out);
+		return out;
+	}
+
+	edge first_adj_edge(const node &n) const
+	{
+		return adj_edges(n).front();
+	}
+
+	std::set<Node> all_nodes() const
+	{
+		return m_g->getNodes();
+	}
+
+	node first_node() const
+	{
+		return *m_g->getNodes().begin();
+	}
+
+	int number_of_nodes() const
+	{
+		return m_g->getNodes().size();
+	}
+
+	int number_of_edges() const
+	{
+		return m_g->getEdges().size();
+	}
+
+	void sort_edges(const edge_array<int> &sorted)
+	{
+		m_sortedEdges = sorted;
+	}
+
+	void del_edge(const edge &e)
+	{
+		m_g->removeEdge(e);
+	}
+
+private:
+	// FIXME: This kills performance!
+	// The correct way to implement this is to use adjacency lists
+	// in graph wrapper and not a DirectedGraph.
+	void sort(std::vector<DirectedEdge> &edges) const
+	{
+		if (m_sortedEdges.empty())
+			return;
+
+		std::sort(edges.begin(), edges.end(), [this](const DirectedEdge &a, const DirectedEdge &b) {
+			return m_sortedEdges.at(a) < m_sortedEdges.at(b);
+		});
+	}
+
+	DirectedGraph *m_g;
+	edge_array<int> m_sortedEdges;
+};
+
+using graph = GraphWrapper;
+
+bool createReversal(const graph &g, edge_array<edge> &reversal)
+{
+	std::vector<DirectedEdge> edges = g.all_edges();
+
+	while (!edges.empty()) {
+		DirectedEdge edge = edges.back();
+		edges.pop_back();
+		bool found = false;
+		for (size_t i = 0; i < edges.size(); ++i) {
+			const DirectedEdge &e = edges[i];
+			if (edge.getFromNode() == e.getToNode() && edge.getToNode() == e.getFromNode()) {
+				reversal.insert({edge, e});
+				reversal.insert({e, edge});
+				edges.erase(edges.begin() + i);
+				found = true;
+				break;
+			}
+		}
+		if (!found)
+			return false;
+	}
+
+	return true;
+}
+
+//  Mehlhorn, Kurt and Mutzel, Petra and Näher, Stefan (1994)
+//  An Implementation of the Hopcroft and Tarjan Planarity Test and Embedding Algorithm.
+
+const int LEFT = 1;
+const int RIGHT = 2;
+
+class Block
+{
+private:
+	std::list<int> Latt, Ratt; // list of attachments
+	std::list<edge> Lseg, Rseg; // list of segments represented by their defining edges
+
+public:
+	Block(const edge &e, std::list<int> &A)
+	{
+		Lseg.push_back(e);
+		conc(Latt, A);
+		// the other two lists are empty
+	}
+
+	void flip()
+	{
+		std::list<int> ha;
+		std::list<edge> he;
+		// we first interchange |Latt| and |Ratt| and then |Lseg| and |Rseg|
+		conc(ha, Ratt);
+		conc(Ratt, Latt);
+		conc(Latt, ha);
+		conc(he, Rseg);
+		conc(Rseg, Lseg);
+		conc(Lseg, he);
+	}
+
+	int head_of_Latt()
+	{
+		return Latt.front();
+	}
+
+	bool empty_Latt()
+	{
+		return Latt.empty();
+	}
+
+	int head_of_Ratt()
+	{
+		return Ratt.front();
+	}
+
+	bool empty_Ratt()
+	{
+		return Ratt.empty();
+	}
+
+	bool left_interlace(std::stack<Block*> &S)
+	{
+		// check for interlacing with the left side of the topmost block of |S|
+		if (Latt.empty())
+			throw exception::AlibException("HopcroftTarjan: Latt is never empty");
+
+		if (!S.empty() && !((S.top())->empty_Latt()) && Latt.back() < (S.top())->head_of_Latt())
+			return true;
+		else
+			return false;
+	}
+
+	bool right_interlace(std::stack<Block*> &S)
+	{
+		// check for interlacing with the right side of the topmost block of |S|
+		if (Latt.empty())
+			throw exception::AlibException("HopcroftTarjan: Latt is never empty");
+
+		if (!S.empty() && !((S.top())->empty_Ratt()) && Latt.back() < (S.top())->head_of_Ratt())
+			return true;
+		else
+			return false;
+	}
+
+	void combine(Block *Bprime)
+	{
+		// add block Bprime to the rear of |this| block
+		conc(Latt, Bprime->Latt);
+		conc(Ratt, Bprime->Ratt);
+		conc(Lseg, Bprime->Lseg);
+		conc(Rseg, Bprime->Rseg);
+		delete Bprime;
+	}
+
+	bool clean(int dfsnum_w, edge_array<int> &alpha)
+	{
+		// remove all attachments to |w|; there may be several
+		while (!Latt.empty() && Latt.front() == dfsnum_w)
+			Latt.pop_front();
+		while (!Ratt.empty() && Ratt.front() == dfsnum_w)
+			Ratt.pop_front();
+
+		if (!Latt.empty() || !Ratt.empty())
+			return false;
+
+		// |Latt| and |Ratt| are empty; we record the placement of the subsegments in |alpha|.
+		for (const edge &e : Lseg)
+			alpha[e] = LEFT;
+
+		for (const edge &e : Rseg)
+			alpha[e] = RIGHT;
+
+		return true;
+	}
+
+	void add_to_Att(std::list<int> &Att, int dfsnum_w0, edge_array<int> &alpha)
+	{
+		// add the block to the rear of |Att|. Flip if necessary
+		if (!Ratt.empty() && head_of_Ratt() > dfsnum_w0)
+			flip();
+
+		conc(Att, Latt);
+		conc(Att, Ratt);
+
+		// This needs some explanation. Note that |Ratt| is either empty
+		// or $\{w0\}$. Also if |Ratt| is non-empty then all subsequent sets are contained
+		// in $\{w0\}$. So we indeed compute an ordered set of attachments.
+
+		for (const edge &e : Lseg)
+			alpha[e] = LEFT;
+
+		for (const edge &e : Rseg)
+			alpha[e] = RIGHT;
+	}
+};
+
+static void dfs_in_reorder(const graph &G, std::list<edge> &Del, node v, int &dfs_count,
+	node_array<bool> &reached,
+	node_array<int> &dfsnum,
+	node_array<int> &lowpt1, node_array<int> &lowpt2,
+	node_array<node> &parent);
+
+static void reorder(graph &G, node_array<int> &dfsnum, node_array<node> &parent)
+{
+	node_array<bool> reached;
+	int dfs_count = 0;
+	std::list<edge> Del;
+	node_array<int> lowpt1, lowpt2;
+
+	dfs_in_reorder(G, Del, G.first_node(), dfs_count, reached, dfsnum, lowpt1, lowpt2, parent);
+
+	// remove forward and reversals of tree edges
+
+	for (const edge &e : Del)
+		G.del_edge(e);
+
+	// we now reorder adjacency lists as described in \cite[page 101]{Me:book}
+
+	edge_array<int> cost;
+	for (const edge &e : G.all_edges()) {
+		node v = e.getFromNode();
+		node w = e.getToNode();
+		cost[e] = ((dfsnum[w] < dfsnum[v]) ?
+		2 * dfsnum[w] :
+		((lowpt2[w] >= dfsnum[v]) ?
+		2 * lowpt1[w] : 2 * lowpt1[w] + 1));
+	}
+
+	G.sort_edges(cost);
+}
+
+static void dfs_in_reorder(const graph &G, std::list<edge> &Del, node v, int &dfs_count,
+	node_array<bool> &reached,
+	node_array<int> &dfsnum,
+	node_array<int> &lowpt1, node_array<int> &lowpt2,
+	node_array<node> &parent)
+{
+	dfsnum[v] = dfs_count++;
+	lowpt1[v] = lowpt2[v] = dfsnum[v];
+	reached[v] = true;
+
+	for (const edge &e : G.adj_edges(v)) {
+		node w = e.getToNode();
+		if (!reached[w]) { // e is a tree edge
+			parent.erase(w);
+			parent.insert({w, v});
+			dfs_in_reorder(G, Del, w, dfs_count, reached, dfsnum, lowpt1, lowpt2, parent);
+			lowpt1[v] = std::min(lowpt1[v], lowpt1[w]);
+		} // end tree edge
+		else {
+			lowpt1[v] = std::min(lowpt1[v], dfsnum[w]); // no effect for forward edges
+			if ((dfsnum[w] >= dfsnum[v]) || w == parent.at(v)) // forward edge or reversal of tree edge
+				Del.push_back(e);
+		} // end non-tree edge
+	} // end forall
+
+	// we know |lowpt1[v]| at this point and now make a second pass over all
+	// adjacent edges of |v| to compute |lowpt2|
+
+	for (const edge &e : G.adj_edges(v)) {
+		node w = e.getToNode();
+		if (parent.at(w) == v) { // tree edge
+			if (lowpt1[w] != lowpt1[v])
+				lowpt2[v] = std::min(lowpt2[v], lowpt1[w]);
+			lowpt2[v] = std::min(lowpt2[v], lowpt2[w]);
+		} // end tree edge
+		else if (lowpt1[v] != dfsnum[w]) // all other edges
+			lowpt2[v] = std::min(lowpt2[v], dfsnum[w]);
+	}
+}
+
+static bool strongly_planar(edge e0, const graph &G, std::list<int> &Att,
+	edge_array<int> &alpha,
+	node_array<int> &dfsnum,
+	node_array<node> &parent)
+{
+	node x = e0.getFromNode();
+
+	node y = e0.getToNode();
+
+	edge e = G.first_adj_edge(y);
+
+	node wk = y;
+
+	while (dfsnum[e.getToNode()] > dfsnum[wk]) { // e is a tree edge
+		wk = e.getToNode();
+		e = G.first_adj_edge(wk);
+	}
+
+	node w0 = e.getToNode();
+
+	node w = wk;
+
+	std::stack<Block*> S;
+
+	while (w != x) {
+		int count = 0;
+		for (const edge &e : G.adj_edges(w)) {
+			count++;
+			if (count != 1) { // no action for first edge
+				std::list<int> A;
+
+				if (dfsnum[w] < dfsnum[e.getToNode()]) { // tree edge
+					if (!strongly_planar(e, G, A, alpha, dfsnum, parent)) {
+						while (!S.empty()) {
+							delete S.top();
+							S.pop();
+						}
+						return false;
+					}
+				}
+				else {
+					A.push_back(dfsnum[e.getToNode()]); // a back edge
+				}
+
+				Block *B = new Block(e, A);
+
+				while (true) {
+					if (B->left_interlace(S))
+						(S.top())->flip();
+
+					if (B->left_interlace(S)) {
+						delete B;
+						while (!S.empty()) {
+							delete S.top();
+							S.pop();
+						}
+						return false;
+					};
+
+					if (B->right_interlace(S)) {
+						B->combine(S.top());
+						S.pop();
+					}
+					else
+						break;
+				} // end while
+
+				S.push(B);
+			} // end if
+		} // end forall
+
+		while (!S.empty() && (S.top())->clean(dfsnum[parent.at(w)], alpha)) {
+			delete S.top();
+			S.pop();
+		}
+
+		w = parent.at(w);
+	} // end while
+
+	Att.clear();
+	while (!S.empty()) {
+		Block *B = S.top();
+		S.pop();
+
+		if (!(B->empty_Latt()) && !(B->empty_Ratt()) &&
+		(B->head_of_Latt() > dfsnum[w0]) && (B->head_of_Ratt() > dfsnum[w0])) {
+			delete B;
+			while (!S.empty()) {
+				delete S.top();
+				S.pop();
+			}
+			return false;
+		}
+
+		B->add_to_Att(Att, dfsnum[w0], alpha);
+		delete B;
+	} // end while
+
+	// Let's not forget (as the book does) that $w0$ is an attachment of $S(e0)$
+	// except if $w0 = x$.
+
+	if (w0 != x)
+		Att.push_back(dfsnum[w0]);
+
+	return true;
+}
+
+static void embedding(edge e0, int t, const graph &G,
+	edge_array<int> &alpha,
+	node_array<int> &dfsnum,
+	std::list<edge> &T, std::list<edge> &A, int &cur_nr,
+	edge_array<int> &sort_num, node_array<edge> &tree_edge_into,
+	node_array<node> &parent, edge_array<edge> &reversal)
+{
+	node x = e0.getFromNode();
+
+	node y = e0.getToNode();
+
+	tree_edge_into.erase(y);
+	tree_edge_into.insert({y, e0});
+
+	edge e = G.first_adj_edge(y);
+
+	node wk = y;
+
+	while (dfsnum[e.getToNode()] > dfsnum[wk]) { // e is a tree edge
+		wk = e.getToNode();
+		tree_edge_into.erase(wk);
+		tree_edge_into.insert({wk, e});
+		e = G.first_adj_edge(wk);
+	}
+
+	node w0 = e.getToNode();
+	edge back_edge_into_w0 = e;
+
+	node w = wk;
+
+	std::list<edge> Al, Ar;
+	std::list<edge> Tprime, Aprime;
+
+	T.clear();
+	T.push_back(e); // |e = (wk,w0)| at this point, line (2)
+
+	while (w != x) {
+		int count = 0;
+		for (const edge &e : G.adj_edges(w)) {
+			count++;
+			if (count != 1) { // no action for first edge
+				if (dfsnum[w] < dfsnum[e.getToNode()]) { // tree edge
+					int tprime = ((t == alpha[e]) ? LEFT : RIGHT);
+					embedding(e, tprime, G, alpha, dfsnum, Tprime, Aprime, cur_nr, sort_num, tree_edge_into, parent, reversal);
+				}
+				else { // back edge
+					Tprime.push_back(e); // $e$
+					Aprime.push_back(reversal.at(e)); // reversal of $e$
+				}
+
+				if (t == alpha[e]) {
+					conc(Tprime, T);
+					conc(T, Tprime); // $T = Tprime\ conc\ T$
+					conc(Al, Aprime); // $Al = Al\ conc\ Aprime$
+				}
+				else {
+					conc(T, Tprime); // $ T\ = T\ conc\ Tprime $
+					conc(Aprime, Ar);
+					conc(Ar, Aprime); // $ Ar\ = Aprime\ conc\ Ar$
+				}
+			} // end if
+		} // end forall
+
+		T.push_back(reversal.at(tree_edge_into.at(w))); // $(w_{j-1},w_j)$
+
+		for (const edge &e : T)
+			sort_num[e] = cur_nr++;
+
+		// |w|'s adjacency list is now computed; we clear |T| and prepare for the
+		// next iteration by moving all darts incident to |parent[w]| from |Al| and
+		// |Ar| to |T|. These darts are at the rear end of |Al| and at the front end
+		// of |Ar|
+
+		T.clear();
+
+		// |parent[w]| is in |G|, |Al.tail| in |H|
+		while (!Al.empty() && Al.back().getFromNode() == parent.at(w)) {
+			T.push_front(Al.back()); // Pop removes from the rear
+			Al.pop_back();
+		}
+
+		T.push_back(tree_edge_into.at(w)); // push would be equivalent
+
+		while (!Ar.empty() && Ar.front().getFromNode() == parent.at(w)) {
+			T.push_back(Ar.front()); // pop removes from the front
+			Ar.pop_front();
+		}
+
+		w = parent.at(w);
+	} // end while
+
+	A.clear();
+
+	conc(A, Ar);
+	A.push_back(reversal.at(back_edge_into_w0));
+	conc(A, Al);
+}
+
+// |Gin| is a directed graph. |planar| decides whether |Gin| is planar.
+// If it is and |embed == true| then it also computes a
+// combinatorial embedding of |Gin| by suitably reordering
+// its adjacency lists.
+// |Gin| must be bidirected in that case.
+bool PLANAR(graph &G, bool embed, HopcroftTarjan::Result &res)
+{
+	int n = G.number_of_nodes();
+
+	if (n <= 3)
+		return true;
+
+	if (G.number_of_edges() > 6 * n - 12)
+		return false;
+
+	{
+		node_array<int> nr;
+		edge_array<int> cost;
+		int cur_nr = 0;
+		int n = G.number_of_nodes();
+
+		for (const node &v : G.all_nodes())
+			nr[v] = cur_nr++;
+
+		for (const edge &e : G.all_edges())
+			cost[e] = ((nr[e.getFromNode()] < nr[e.getToNode()]) ?
+			n * nr[e.getFromNode()] + nr[e.getToNode()] :
+			n * nr[e.getToNode()] + nr[e.getFromNode()]);
+
+		G.sort_edges(cost);
+
+		std::vector<edge> L = G.all_edges();
+		while (!L.empty()) {
+			edge e = L.front();
+			L.erase(L.begin());
+			// check whether the first edge on |L| is equal to the reversal of |e|. If so,
+			// delete it from |L|, if not, add the reversal to |G|
+			if (!L.empty() && (e.getFromNode() == L.front().getToNode()) && (L.front().getFromNode() == e.getToNode()))
+				L.erase(L.begin());
+			else
+				throw exception::AlibException("HopcroftTarjan: Graph not bidirected");
+		}
+	}
+
+	// An undirected planar graph has at most $3n-6$ edges; a directed graph may have twice as many
+
+	edge_array<edge> reversal;
+	if (!createReversal(G, reversal))
+		throw exception::AlibException("HopcroftTarjan: Graph not bidirected");
+
+	node_array<int> dfsnum;
+	node_array<node> parent;
+	for (const node &n : G.all_nodes())
+		parent.insert({n, node("HopcroftTarjan_invalid")}); // no default constructor -,-
+
+	reorder(G, dfsnum, parent);
+
+	edge_array<int> alpha;
+	for (const edge &e : G.all_edges())
+		alpha[e] = 0;
+
+	{
+		std::list<int> Att;
+
+		alpha[G.first_adj_edge(G.first_node())] = LEFT;
+
+		if (!strongly_planar(G.first_adj_edge(G.first_node()), G, Att, alpha, dfsnum, parent))
+			return false;
+	}
+
+	if (embed) {
+		std::list<edge> T, A; // lists of edges of |H|
+
+		int cur_nr = 0;
+		edge_array<int> sort_num;
+		node_array<edge> tree_edge_into;
+
+		embedding(G.first_adj_edge(G.first_node()), LEFT, G, alpha, dfsnum, T, A,
+		cur_nr, sort_num, tree_edge_into, parent, reversal);
+
+		// |T| contains all edges incident to the first node except the cycle edge into it.
+		// That edge comprises |A|
+
+		conc(T, A);
+
+		for (const edge &e : T)
+			sort_num[e] = cur_nr++;
+
+		std::vector<const DirectedEdge*> edgs(sort_num.size());
+		for (const auto &p : sort_num)
+			edgs[p.second] = &p.first;
+
+		for (const DirectedEdge *e : edgs)
+			res[e->getFromNode()].push_back(e->getToNode());
+	}
+
+	return true;
+}
+
+HopcroftTarjan::Result HopcroftTarjan::hopcrofttarjan(const Graph &graph)
+{
+	return getInstance().dispatch(graph.getData());
+}
+
+HopcroftTarjan::Result HopcroftTarjan::hopcrofttarjan(const DirectedGraph &graph)
+{
+	DirectedGraph copy = graph;
+	GraphWrapper graphw(&copy);
+	Result res;
+	if (!PLANAR(graphw, true, res))
+		throw exception::AlibException("HopcroftTarjan: Graph not planar");
+	return res;
+}
+
+HopcroftTarjan::Result HopcroftTarjan::hopcrofttarjan(const UndirectedGraph &graph)
+{
+	DirectedGraph copy;
+	for (const UndirectedEdge &e : graph.getEdges()) {
+		DirectedEdge e1(e.getFirstNode(), e.getSecondNode(), e.getName());
+		DirectedEdge e2(e.getSecondNode(), e.getFirstNode(), e.getName());
+		copy.addEdge(e1);
+		copy.addEdge(e2);
+	}
+	GraphWrapper graphw(&copy);
+	Result res;
+	if (!PLANAR(graphw, true, res))
+		throw exception::AlibException("HopcroftTarjan: Graph not planar");
+	return res;
+}
+
+} // namespace embedding
+
+} // namespace graph
diff --git a/alib2algo/src/graph/embedding/HopcroftTarjan.h b/alib2algo/src/graph/embedding/HopcroftTarjan.h
new file mode 100644
index 0000000000000000000000000000000000000000..797fdd6eed0efad3c731ba88a347c20e6b6cc8a7
--- /dev/null
+++ b/alib2algo/src/graph/embedding/HopcroftTarjan.h
@@ -0,0 +1,37 @@
+#ifndef GRAPH_HOPCROFT_TARJAN_H_
+#define GRAPH_HOPCROFT_TARJAN_H_
+
+#include <common/multipleDispatch.hpp>
+
+#include <graph/Graph.h>
+#include <graph/directed/DirectedGraph.h>
+#include <graph/undirected/UndirectedGraph.h>
+
+namespace graph
+{
+
+namespace embedding
+{
+
+// Computes combinatorial embedding of bidirected biconnected simple planar graph
+class HopcroftTarjan : public std::SingleDispatch<std::unordered_map<Node, std::vector<Node>>, graph::GraphBase>
+{
+public:
+	typedef std::unordered_map<Node, std::vector<Node>> Result;
+
+	static Result hopcrofttarjan(const Graph &graph);
+
+	static Result hopcrofttarjan(const DirectedGraph &graph);
+	static Result hopcrofttarjan(const UndirectedGraph &graph);
+
+	static HopcroftTarjan &getInstance() {
+		static HopcroftTarjan res;
+		return res;
+	}
+};
+
+} // namespace embedding
+
+} // namespace graph
+
+#endif // GRAPH_HOPCROFT_TARJAN_H_