Skip to main content
2 of 2
compute vi using formula, better computation of b, switched to lowercase variables (mathematica standard)
2012rcampion
  • 7.9k
  • 26
  • 44

We can start by eliminating $C$ from your equations, as it's trivial to solve for once you have the value of $v$. I'm also going to simplify some of your variable names:

$$ \begin{align} bt_ixt'_i &= a_i \\ bt'_i &= b_i \end{align} $$

Now we simply have:

$$ c \sum_i \left|a_i - vb_i\right| + v = \sum_i a_i $$

Defining $s_i$ as:

$$ s_i = \left\{ \begin{array}{rl} 1 & : a_i - vb_i > 0 \\ -1 & : a_i - vb_i < 0 \end{array} \right. $$

we can rewrite the absolute value as:

$$ c \sum_i s_i\left(a_i - vb_i\right) + v = \sum_i a_i \\ f(v) = c \sum_i s_i a_i - \sum_i a_i - v \left(c \sum_i s_i b_i - 1\right) = 0 $$

Now if we know the $s_i$, we can solve for $v$ explicitly:

$$ v = \frac{c \sum_i s_i a_i - \sum_i a_i}{c \sum_i s_i b_i - 1} = \frac{\sum_i \left(c s_i - 1\right) a_i}{\sum_i cs_i b_i - 1} $$

This matches your solution, with the $s_i$ negated.

We can rewrite $s_i$ to make something a little more clear:

$$ \begin{align} s_i &= \left\{ \begin{array}{rl} 1 & : vb_i < a_i \\ -1 & : vb_i > a_i \end{array} \right. \\ &= \left\{ \begin{array}{rl} 1 & : v < a_i/b_i \oplus b_i < 0 \\ -1 & : v > a_i/b_i \oplus b_i < 0 \end{array} \right. \end{align} $$

We can now see that we don't have to try every combination of the $s_i$, since they will change one by one as $v$ increases. We can now write a simple $\text{O}(n)$ algorithm to find the solution.


The algorithm takes three inputs: two arrays a and b of equal length, and a real number c. Our first step is to sort a and b by a/b:

{a, b} = Transpose[SortBy[Transpose[{a,b}], Divide @@ # &]]

Next, we calculate the breakpoints for v (where the s[[i]] change sign):

vb = {-Infinity} ~Join~ (a/b) ~Join~ {+Infinity}

Now we loop over each interval in vb. For each one we will compute s and solve for v:

m = Length[a]

vi = Table[
  With[{s = Sign[b] Join @@ ConstantArray @@@ {{-1, i}, {1, m - i}} },
   (c s - 1).a/(c s.b - 1)
  ], {i, 0, m}]

Now we select only the vi that are in the proper range:

solns = Pick[vi, Thread[Most[vb] <= vi <= Rest[vb]]]

If you have very large a and b lists, there is a trick that you can use to speed up finding the solutions. Note that the equation $f(v)=0$ is piecewise linear and piecewise continuous. When $v=-\infty$, $s_i$ is just $\text{sign}(b_i)$. This means that all of the linear terms $-v s_i b_i)$ are negative.

As $v$ increases, the signs flip one by one until they are all positive. This means that the slope is always nondecreasing; therefore $f(v)$ is concave up (unless $c$ is negative, in which case it is concave down), and there can be only two solutions. (There is a special case where one segment of the function is colinear with the x-axis, but when a and b are real numbers the probability of this is zero, i.e. it is almost certain not to happen.) You can do a modified binary search to quickly find the two by looking at the sign of the residual to figure out where in the function you are.

2012rcampion
  • 7.9k
  • 26
  • 44