# Q：将一个数分解为随机非零数

I have some pairs of numbers, let's say

``````[X1,Y1], [X2,Y2] ..... [Xn,Yn]    Y > X
``````

I have a another number, say Z. Now what I want is to split Z into random values Nm such that it always meets following condition:

``````N1 + N2 + N3 ...... Nn = Z
Nm > 0
Nm <=  Ym - Xm
``````

``````[X1,Y1], [X2,Y2] ..... [Xn,Yn]    Y > X
``````

``````N1 + N2 + N3 ...... Nn = Z
Nm > 0
Nm <=  Ym - Xm
``````
``````IF the series y_1 - x_1 + … + y_m - x_m < z OR x_i ≥ y_i for any i
FAIL
ELSE DO
residual := z
FOR i := 1 to m-1
n_i := (y_i - x_i) × random (0,1]
residual := residual - n_i
IF residual ≤ 0
BREAK
n_m := residual
WHILE n_m ≤ 0 OR n_m > y_m-x_m
``````

Every n_i before n_m falls into the range because it is the upper bound times a random number between 0 and 1. The last value falls into the range because, if it did not, we would have looped and tried again. The numbers add up to z because n_m is the difference of z and the sum of the rest of the series.

``````IF the series y_1 - x_1 + … + y_m - x_m < z OR x_i ≥ y_i for any i
FAIL
ELSE DO
residual := z
FOR i := 1 to m-1
n_i := (y_i - x_i) × random (0,1]
residual := residual - n_i
IF residual ≤ 0
BREAK
n_m := residual
WHILE n_m ≤ 0 OR n_m > y_m-x_m
``````

You have M intervals [ Xi, Yi ] and you want a random number in each interval so that the sum of all numbers is Z.

The sum of such numbers must needs be between SUM(Xi) and SUM(Yi). So we calculate both:

``````Zmin = 0;
Zmax = 0;
for (i = 0; i < m; i++) {
Zmin += X[i];
Zmax += Y[i];
Z[i] = X[i];
}
``````

Now we check that Z is between Zmin and Zmax. If it is not, there's no way of obtaining a Z vector that meets the constraints.

If Z is between Zmin and Zmax, we can imagine to start with a Z vector which is a copy of X, and sprinkle the difference between Z and Zmin at random among all Zi's.

``````remainder = Z - Zmin;
if ((remainder < 0) || (remainder + Zmin > Zmax)) {
// Can't do.
}

// There are several ways of equidistributing the remainder.
// This is not one of them, but it's simple.
// If Z is very near Zmax, performances will suffer.
while (remainder > 0) {
// Choose
i = rand()*m;
// Can we increase Z[i]?
if (Z[i] < Y[i]) {
Z[i]++;
remainder--;
}
}
// DONE!
``````

A more complicated way that would run in O(n) would be to equidistribute the remainder R among the m slots in proportion to each slot's width. On average, a slot (Y[i]-X[i]) wide should receive an extra of R*(Y[i]-X[i])/(Zmax - Zmin). We could use Bresenham's algorithm to do this as exactly as possible, but if we want to have a random distribution, we run the risk of under-allocating the first K slots, and having such a large remainder that even if we maxed out all the remaining slots we still wouldn't exhaust it.

To avoid this, we need to keep a running total of how much we can defer to the following intervals, starting with Zmax - Zmin. So at iteration i, if the next intervals from i+1 to m can dispatch up to D, and remainder is R, whatever we add to Z[i] must leave (remainder - Z[i] + X[i] <= D), which means that Z[i] - X[i] >= remainder - D; i.e., our random number that we add to Z[i] (currently still equal to X[i]) must start from MAX(remainder - D, 0). Ideally, this will always be 0. So we can attribute a random number from MAX(remainder-D, 0) to Y[i]-X[i]. Let this interval be of width Q. We want to assign on average Q*remainder/(D+Q) to this interval, and that means Z[i] = Y[i]-X[i]-Q+rand()*(Q*remainder/(D+Q)).

``````D = Zmax - Zmin;

for (i = 0; i < m; i++) {

V = MAX(remainder-D, 0);

Q = (Y[i] - X[i]) - V; // This is positive, or we would have failed earlier

Z[i] = Y[i]-X[i]-Q+rand()*(Q*remainder/(D+Q));

remainder -= (Z[i] - X[i]);
}
``````

This is a quick implementation in PHP which appears to work.

``````<?php
\$X = [ 1, 4, 11, 3, 5, 17, 22, 35, 120, 0 ];
\$Y = [ 8, 9, 33, 9, 9, 28, 24, 36, 215, 3 ];

\$m = count(\$X);

\$ZZ= 0;
// Construct a Z that will work.
for (\$i = 0; \$i < \$m; \$i++) {
// \$Y[\$i] = rand(\$X[\$i], \$X[\$i] + 20);
\$ZZ += rand(\$X[\$i], \$Y[\$i]);
}

\$Zmin = 0;
\$Zmax = 0;

for (\$i = 0; \$i < \$m; \$i++) {
\$Zmin += \$X[\$i];
\$Zmax += \$Y[\$i];
}

\$R = \$ZZ - \$Zmin;
if ((\$R < 0) || (\$R + \$Zmin > \$Zmax)) {
die("Can't do.\n");
}

\$D = \$Zmax - \$Zmin;

for (\$i = 0; \$i < \$m; \$i++) {
\$A = max(\$R-\$D, 0); // Cannot distribute less than this.
\$B = min(\$Y[\$i]-\$X[\$i], \$R); // Nor more than this.

\$Q = (\$B - \$A); // This is positive, or we would have failed earlier
assert('\$Q > 0');

\$Z[\$i] = \$X[\$i] + rand(\$A, \$B); // floor(\$A + (\$B-\$A)*rand());

\$R -= (\$Z[\$i] - \$X[\$i]);
assert('\$R >= 0');
assert('\$X[\$i] <= \$Z[\$i]');
assert('\$Y[\$i] >= \$Z[\$i]');
}
``````

``````Zmin = 0;
Zmax = 0;
for (i = 0; i < m; i++) {
Zmin += X[i];
Zmax += Y[i];
Z[i] = X[i];
}
``````

``````remainder = Z - Zmin;
if ((remainder < 0) || (remainder + Zmin > Zmax)) {
// Can't do.
}

// There are several ways of equidistributing the remainder.
// This is not one of them, but it's simple.
// If Z is very near Zmax, performances will suffer.
while (remainder > 0) {
// Choose
i = rand()*m;
// Can we increase Z[i]?
if (Z[i] < Y[i]) {
Z[i]++;
remainder--;
}
}
// DONE!
``````

``````D = Zmax - Zmin;

for (i = 0; i < m; i++) {

V = MAX(remainder-D, 0);

Q = (Y[i] - X[i]) - V; // This is positive, or we would have failed earlier

Z[i] = Y[i]-X[i]-Q+rand()*(Q*remainder/(D+Q));

remainder -= (Z[i] - X[i]);
}
``````

``````<?php
\$X = [ 1, 4, 11, 3, 5, 17, 22, 35, 120, 0 ];
\$Y = [ 8, 9, 33, 9, 9, 28, 24, 36, 215, 3 ];

\$m = count(\$X);

\$ZZ= 0;
// Construct a Z that will work.
for (\$i = 0; \$i < \$m; \$i++) {
// \$Y[\$i] = rand(\$X[\$i], \$X[\$i] + 20);
\$ZZ += rand(\$X[\$i], \$Y[\$i]);
}

\$Zmin = 0;
\$Zmax = 0;

for (\$i = 0; \$i < \$m; \$i++) {
\$Zmin += \$X[\$i];
\$Zmax += \$Y[\$i];
}

\$R = \$ZZ - \$Zmin;
if ((\$R < 0) || (\$R + \$Zmin > \$Zmax)) {
die("Can't do.\n");
}

\$D = \$Zmax - \$Zmin;

for (\$i = 0; \$i < \$m; \$i++) {
\$A = max(\$R-\$D, 0); // Cannot distribute less than this.
\$B = min(\$Y[\$i]-\$X[\$i], \$R); // Nor more than this.

\$Q = (\$B - \$A); // This is positive, or we would have failed earlier
assert('\$Q > 0');

\$Z[\$i] = \$X[\$i] + rand(\$A, \$B); // floor(\$A + (\$B-\$A)*rand());

\$R -= (\$Z[\$i] - \$X[\$i]);
assert('\$R >= 0');
assert('\$X[\$i] <= \$Z[\$i]');
assert('\$Y[\$i] >= \$Z[\$i]');
}
``````
algorithm  math