Thorup.cs

void Main()
{
	 int[][] adj = new int[][] {
            new int[] {0, 10, 15, 10, 20},
            new int[] {10, 0, 35, 20, 25},
            new int[] {15, 35, 0, 30, 10},
			new int[] {50, 60, 30, 0, 10},
			new int[] {20, 25, 10, 10, 0}
		};

	// Generate Thorup's shortest path resolution for the graph
	int[][] parent, ancestor;
	int[][][] bucket;
	int[] middle;
	Thorups.ThorupsResolve(adj, out parent, out ancestor, out bucket, out middle);
	Console.WriteLine(parent);
	// Find the shortest path between two vertices using the resolution
	int shortestDistance = (int)Thorups.FindShortestPath(adj, parent, ancestor, bucket, middle, 0, 4);

	// Print the result
	Console.WriteLine("Shortest distance between vertices 0 and 4: " + shortestDistance);
}
public static class Thorups
{
	public static double FindShortestPath(int[][] adj, int[][] parent, int[][] ancestor, int[][][] bucket, int[] middle, int u, int v)
{
    int n = adj.Length;

		// Step 1: compute the distance from u to v
		double dist = double.MaxValue;
		int ancestorUV = -1;

		for (int w = u; w != -1; w = parent[middle[w]][w])
		{
			if (bucket[0][w].Length > 0)
			{
				int minBucket = bucket[0][w][0];
				if (adj[w][minBucket] + GetEdgeWeight(adj, w, minBucket) < dist)
				{
					dist = adj[w][minBucket] + GetEdgeWeight(adj, w, minBucket);
					ancestorUV = w;
				}
			}
		}

		for (int w = v; w != -1; w = parent[middle[w]][w])
		{
			if (bucket[0][w].Length > 0)
			{
				int minBucket = bucket[0][w][0];
				if (adj[w][minBucket] + GetEdgeWeight(adj, w, minBucket) + GetDistance(ancestor, w, ancestorUV) < dist)
				{
					dist = adj[w][minBucket] + GetEdgeWeight(adj, w, minBucket) + GetDistance(ancestor, w, ancestorUV);
					ancestorUV = w;
				}
			}
		}

		// Step 2: construct the path from u to ancestorUV
		List<int> path = new List<int>();

		for (int w = u; w != ancestorUV; w = parent[middle[w]][w])
		{
			path.Add(w);
		}

		path.Add(ancestorUV);

		// Step 3: construct the path from ancestorUV to v
		List<int> reversePath = new List<int>();

		for (int w = v; w != ancestorUV; w = parent[middle[w]][w])
		{
			reversePath.Add(w);
		}

		reversePath.Reverse();

		path.AddRange(reversePath);

		// Step 4: return the total distance of the path
		return dist;
	}
	public static double GetDistance(int[][] adj, int u, int v)
	{
		return adj[u][v];
	}

	public static void ThorupsResolve(int[][] adj, out int[][] parent, out int[][] ancestor, out int[][][] bucket, out int[] middle)
	{
		int n = adj.Length;

		// Initialize arrays
		parent = new int[n][];
		ancestor = new int[n][];
		bucket = new int[n][][];
		for (int i = 0; i < n; i++)
		{
			parent[i] = new int[2];
			ancestor[i] = new int[2];
			bucket[i] = new int[2][];
			for (int j = 0; j < 2; j++)
			{
				bucket[i][j] = new int[n];
				for (int k = 0; k < n; k++)
				{
					bucket[i][j][k] = -1;
				}
			}
		}
		middle = new int[n];

		// Run the Thorup's algorithm
		ThorupsAlgorithm(adj, parent, ancestor, bucket, middle);
	}
	public static void ThorupsAlgorithm(int[][] adj, int[][] parent, int[][] ancestor, int[][][] bucket, int[] middle)
	{
		int n = adj.Length;

		// Step 1: Initialize parent and ancestor arrays
		for (int i = 0; i < n; i++)
		{
			parent[i][0] = -1;
			parent[i][1] = -1;
			ancestor[i][0] = -1;
			ancestor[i][1] = -1;
		}

		// Step 2: Compute the middle node for each vertex
		ComputeMiddle(adj, middle);

		// Step 3: Compute the shortest path from each vertex to its middle node
		ComputeShortestPaths(adj, middle, parent, bucket);

		// Step 4: Compute the ancestor and bucket arrays
		ComputeAncestors(adj, parent, ancestor, bucket);

		// Step 5: Sort the buckets
		SortBuckets(ancestor, bucket);
	}
	public static void ComputeAncestors(int[][] adj, int[][] parent, int[][] ancestor, int[][][] bucket)
	{
		int n = adj.Length;
		for (int i = 0; i < n; i++)
		{
			for (int j = 0; j < bucket.Length; j++)
			{
				if (bucket[j][i].Length > 0)
				{
					int w = bucket[j][i][0];
					ancestor[w][i] = w;
					for (int k = 1; k < bucket[j][i].Length; k++)
					{
						int v = bucket[j][i][k];
						ancestor[v][i] = LCA(ancestor[v][i], w, parent, adj);
						w = v;
					}
				}
			}
		}
	}
	
	public static int LCA(int u, int v, int[][] parent, int[][] adj)
	{
		if (u == v)
		{
			return u;
		}
		if (adj[u].Length > adj[v].Length)
		{
			int temp = u;
			u = v;
			v = temp;
		}
		int diff = adj[v].Length - adj[u].Length;
		for (int i = 0; i < adj[u].Length; i++)
		{
			if (adj[u][i] == v)
			{
				return u;
			}
		}
		for (int i = 0; i < adj[u].Length; i++)
		{
			int w = adj[u][i];
			if (parent[w][diff] != -1 && parent[w][diff] != u)
			{
				int p = LCA(parent[w][diff], v, parent, adj);
				if (p != -1)
				{
					return p;
				}
			}
		}
		return -1;
	}
	public static void SortBuckets(int[][] ancestor, int[][][] bucket)
	{
		for (int i = 0; i < ancestor.Length; i++)
		{
			for (int j = 0; j < bucket[i].Length; j++)
			{
				if (bucket[i][j] == null)
				{
					continue;
				}
				int[] b = bucket[i][j];
				Array.Sort(b, (a, c) => ancestor[i][a] - ancestor[i][c]);
			}
		}
	}
	public static int GetEdgeWeight(int[][] adj, int u, int v)
{
		
		return adj[u][v];
	}
	public static void ComputeShortestPaths(int[][] adj, int[] middle, int[][] parent, int[][][] bucket)
	{
		int n = adj.Length;

		// Initialize parent array to -1 (i.e., no parent)
		for (int i = 0; i < n; i++)
		{
			for (int j = 0; j < n; j++)
			{
				parent[i][j] = -1;
			}
		}

		// Initialize bucket array to empty buckets
		for (int i = 0; i < n; i++)
		{
			for (int j = 0; j < n; j++)
			{
				bucket[i][j] = new int[0];
			}
		}

		// Initialize distance array to infinity
		double[][] dist = new double[n][];
		for (int i = 0; i < n; i++)
		{
			dist[i] = new double[n];
			for (int j = 0; j < n; j++)
			{
				dist[i][j] = double.PositiveInfinity;
			}
		}

		// Compute shortest paths using middle nodes
		for (int i = 0; i < n; i++)
		{
			// Initialize the distance to the middle node to 0
			dist[i][middle[i]] = 0;

			// Add vertex i to the bucket of the middle node
			bucket[i][middle[i]] = new int[] { i };
		}

		// Iterate over all distances
		for (int k = 0; k < n; k++)
		{
			for (int i = 0; i < n; i++)
			{
				for (int j = 0; j < n; j++)
				{
					if (dist[i][j] > dist[i][k] + GetEdgeWeight(adj, k, j))
					{
						// Update the distance and parent arrays
						dist[i][j] = dist[i][k] + GetEdgeWeight(adj, k, j);
						parent[i][j] = k;

						// Update the bucket of vertex i and middle node j
						int middleNode = middle[j];
						bucket[i][middleNode] = AddToBucket(bucket[i][middleNode], j);
					}
				}
			}
		}
	}
	private static int[] AddToBucket(int[] bucket, int node)
{
    int[] newBucket = new int[bucket.Length + 1];
    int i = 0;
    bool added = false;

		foreach (int item in bucket)
		{
			if (!added && node < item)
			{
				newBucket[i++] = node;
				added = true;
			}

			if (item != node)
			{
				newBucket[i++] = item;
			}
		}

		if (!added)
		{
			newBucket[newBucket.Length - 1] = node;
		}

		return newBucket;
	}

	public static void ComputeMiddle(int[][] adj, int[] middle)
	{
		int n = adj.Length;

		// Compute the sum of edge weights for each vertex
		long[] sumWeights = new long[n];
		for (int i = 0; i < n; i++)
		{
			foreach (int j in adj[i])
			{
				sumWeights[i] += GetEdgeWeight(adj, i, j);
			}
		}

		// Compute the middle node for each vertex
		for (int i = 0; i < n; i++)
		{
			int m = adj[i].Length;
			if (m == 0)
			{
				// Vertex i has no neighbors, so its middle node is itself
				middle[i] = i;
			}
			else
			{
				// Vertex i has neighbors, so its middle node is the one with the highest sum of edge weights
				int maxNeighbor = adj[i][0];
				for (int j = 1; j < m; j++)
				{
					int neighbor = adj[i][j];
					if (sumWeights[neighbor] > sumWeights[maxNeighbor])
					{
						maxNeighbor = neighbor;
					}
				}
				middle[i] = maxNeighbor;
			}
		}
	}
}
// You can define other methods, fields, classes and namespaces here