How to randomly generate two integer matrices $A$ and $B$, so that entries of 3 matrices $A$, $B$, and $AB$ are within certain range? I ran into this question when writing a program. I need to generate two matrices, and calculate their product.
However, I must ensure all entries are within 8-bit signed integer range, i.e. $[-128, 128)$. Is there a way to algorithmically achieve this?
Furthermore, what if I need every intermediate result during calculation is also within such range?
Besides, I think the algorithm needs not to be deterministic, but with high probability is enough. For example, is there a way to randomly generate the initial $A$ and $B$, so that the entries of product $AB$ are highly likely to within $[-128, 127)$? If so, I can run the program several times to be lucky.

Update
The dimensions of $A$ and $B$ are inputs of the program, which are typically 1K ~ 5K.
 A: Let $D$ be the dimension of the matrices. For each row of $A$, simply choose $127$ of them at random to equal $1$ and set the rest to $0$. Do the same for each column of $B$. You are then guaranteed that the entries of $AB$ are nonnegative and at most $127$.
For a more random approach, choose each entry of $A$ and $B$ independently, taking the value $-1$ with probability $p$, $1$ with probability $p$, and $0$ with probability $1-2p$. The variance of each entry of $AB$ is $4p^2D$, which will be quite small if $p$ is sufficiently small as a function of $D$. Thus with high probability, the entries of $AB$ will all be less than $128$ in absolute value.
A: With Python...
import numpy as np

low = -128
high = 128
size = (5,5)
dt = np.dtype('i1')
A = np.random.randint(low, high,size, dtype= dt)
B = np.random.randint(low, high, size, dtype=dt)

C = np.dot(A,B)

The entries can't be outside your range because the representation doesn't exist.
array([[-124,  -95,   93,  -61,   40],
       [ -54, -117,   -7,   35,   68],
       [  79,   68,  -82,   83,  -68],
       [  91,  -14,   28,  102,   92],
       [ -40,   30, -122,  -30,  -40]], dtype=int8)

If you change the dtype you'd see that it will go outside the range..
dt = np.dtype('i2')


array([[ -7161,   3352,  13377, -12262,  -2379],
       [ -7873,   8113,  -5097,  -7848,   3062],
       [ -8811,   9435,   9431, -11306,   1616],
       [ 12141, -24141,  24527,   5043,   -622],
       [  4488,  -6618,  -5862,  15442,   1308]], dtype=int16)

For a large matrix you can generate a sparse matrix. The library is scipy sparse. Then determine the level of sparsity so you can fit it in memory.
