Decompose the domain

Question

In scientific computing, it is often desirable to decompose some large domain into smaller "subdomains", then parallelize the computation by performing the work on each subdomain in parallel and the interaction of a subdomain with its neighboring subdomains. Since these interactions are typically data sync points, minimizing their size improves the scalability of the code to very large machines.

Consider a domain consisting of an N-dimensional Cartesian grid, where N is between 1 up to 6, inclusive. The minimum extent of this grid is at the coordinate (0,0,...0), and the maximum extent is likewise (n,n,...,n) (note that n may be very large!).

Your goal is to partition the domain into m non-overlapping hyperslabs (N-D rectangles) such that the sum of the surface area of all subdomains is minimized, and the ratio of surface areas of any two subdomains is approximately unity (within 10% of an optimal solution or as close as possible).

Inputs

Your program should take from any input source 3 positive integer numbers: n, N, and m. You may assume that m is smaller than the total volume of the domain (i.e. a valid solution can be found). As mentioned earlier, n can be quite large. For the purposes of this challenge, you may assume it to fit in your language's natural integer type, or a signed 32-bit integer (whichever is larger). Error checking of inputs is not required.

Outputs

Your program should output to any output sink a set of m hyperslabs describing the produced partition (order does not matter).

A hyperslab description must contain enough information to reconstruct the minimum and maximum extents of the hyperslab in the coordinates of the global domain, i.e. you can either output these values directly, or output the centroid+width in each dimension, etc.

None of the output hyperslabs can be empty (i.e. contain 0 volume). Note that every hyperslab edge must be at an integer increment.

The union of all outputted hyperslabs must cover the entire domain.

The output does not need to be deterministic (i.e. concurrent runs with the same inputs need not give the same outputs, so long as all outputs are valid).

Examples

Note that the example outputs for these inputs may not necessarily match what your code produces. Hyperslabs are specified by the min/max extent in each dimension.

input:
4, 1, 1
output:
[(0),(4)]

~

input:
4, 1, 2
output:
[(0),(1)], [(1),(4)]
sidenote:
typically in a "real" partitioner, this would be a bad partition since one subdomain contains more work than the other.
However, for the purposes of this challenge this is a perfectly valid output.

~

input:
4, 1, 2
output:
[(0),(2)], [(2),(4)]
sidenote:
an alternative "better" valid output for the same inputs as above.

~

input:
4, 2, 2
output:
[(0,0),(2,4)], [(2,0),(4,4)]

~

input:
1000000, 6, 7
output:
[(0, 0, 0, 0, 0, 0),(428571, 333333, 1000000, 1000000, 1000000, 1000000)],
[(428571, 0, 0, 0, 0, 0),(1000000, 500000, 500000, 1000000, 1000000, 1000000)],                                                                                              
[(0, 333333, 0, 0, 0, 0),(428571, 1000000, 500000, 1000000, 1000000, 1000000)],                                                                                              
[(0, 333333, 500000, 0, 0, 0),(428571, 1000000, 1000000, 1000000, 1000000, 1000000)],
[(428571, 500000, 0, 0, 0, 0),(1000000, 1000000, 500000, 1000000, 1000000, 1000000)],
[(428571, 0, 500000, 0, 0, 0),(1000000, 500000, 1000000, 1000000, 1000000, 1000000)],
[(428571, 500000, 500000, 0, 0, 0),(1000000, 1000000, 1000000, 1000000, 1000000, 1000000)]

Scoring

This is code golf; shortest code wins. Standard loopholes apply. You may use any built-ins so long as they weren't specifically designed to solve this problem. Your code should be able to run in a reasonable amount of time (say, 10 minutes on whatever hardware you have on hand, assuming m isn't too large)