I have a linear program and I need to check if it has a solution.
Some coefficients entering the program are much smaller than others. E.g., some coefficients are of the order $10^{-10}$, while the others are of the order of $10^{-2}$.
When coefficients are very small, many optimization algorithms do not work properly and, to address the issue, some of them arbitrarily drop coefficients with small values.
Therefore, I need to rescale the linear program.
However, if I rescale by multiplying everything, e.g., by $10^7$, the small coefficients become acceptable, but the others become too large.
What can be a correct way of rescaling?
This is the linear program written in a compact way: $$ \begin{aligned} &(1)\quad A \text{ }x = b\\ &(2)\quad C \text{ }x\leq c\\ \end{aligned} $$
$A$ and $C$ contain the coefficients with magnitude of different order as mentioned above. $A$ is $m\times n$ with $n>m$. $C$ is $j\times n$ with $n>j$. $A$ has rank $m$. $C$ has rank $h<j$.
If I had only $A$, I think I could use the full reduced form of $A$. But, here, I also have $C$.
This is the linear program in details (if needed):
The unknown of the linear program is an array $\lambda$ of size $\overbrace{(K+1)\times (K+1) \times ...\times (K+1)}^{D \text{ times}}\times L$. I denote its $(k_1,k_2,...,k_D,y)-th$ element by $\lambda(k_1,k_2,...,k_D,y)$.
The coefficients of the linear program are three arrays, $\alpha$, $\beta$, and $\delta$.
$\alpha$ has size $\overbrace{(K+1)\times (K+1) \times ...\times (K+1)}^{D \text{ times}}\times \overbrace{L \times L}^{2 \text{ times}}$. I denote its $(k_1,k_2,...,k_D,y, y')$-th element by $\alpha(k_1,k_2,...,k_D,y, y')$.
$\beta$ has size $\overbrace{(K+1)\times (K+1) \times ...\times (K+1)}^{D \text{ times}}$. I denote its $(k_1,k_2,...,k_D)$-th element by $\beta(k_1,k_2,...,k_D)$.
$\delta$ has size $L\times 1$. I denote its $y$-th element by $\delta(y)$.
This is the linear program: $$ \begin{aligned} &(1) \quad \sum_{k_1=0}^K \sum_{k_2=0}^K... \sum_{k_D=0}^K \lambda(k_1,..., k_D,y) \alpha(k_1,...,k_D,y, y')\geq 0 \quad \forall y,y' \text{ with } y\neq y'\\ &(2) \quad \lambda(k_1,..., k_D,y) \geq 0 \quad \forall k_1,..., k_D, y\\ & (3) \sum_{y=1}^L \lambda(k_1,..., k_D,y) =1\quad \forall k_1,..., k_D\\ & (4) \sum_{k_1=0}^K \sum_{k_2=0}^K... \sum_{k_D=0}^K \lambda(k_1,..., k_D,y) \beta(k_1,..., k_D)=\delta(y) \quad \forall y \end{aligned} $$
The arrays that some very small numbers when $K$ is large are $\alpha$ and $\beta$.
It is unlikely that linearly rescaling will work, as rescaling is approximately orthogonal to the underlying double precision.
This is because IEEE764 doubles store floating point numbers as $(\pm 1) x^y$ where $x$ has up to $53$ bits of precision and $y$ has up to $11$ bits. The significand $x$ carries only about 15 base-10 digits of precision. Scaling each coefficient approximately only effects the exponents $y$, but the precision issues you're seeing are due to insufficient precision in the significands.
If the coefficients of the matrices $A$ and $C$ vary by $8$ orders of magnitude, then double precision is incapable of recognizing cancellation among the larger coefficients. For example, if the $a_i$ are each of size $10^8$ and the $b_i$ are each roughly constant, then $a_1 a_2 + b_1 \approx a_1 a_2$ in double precision arithmetic, as $b_1$ is too small to affect the precision of the significant in $a_1 a_2$.
It might be possible that using higher precision (e.g. 128-bit quad precision, or arbitrary precision) would work. It might also be necessary to try to modify your linear program to be better conditioned --- but this can be hard.