Fastest algorithm to take the product of all subsets

23

8

Given n numbers in an array (you can't assume they are integers), I would like to compute the product of all subsets of size n-1.

You can do this by multiplying all the numbers together and then dividing by each one in turn, as long as none of the numbers is zero. However, how quickly can you do this with no division?

If you don't allow division, what is the minimum number of arithmetic operations (e.g. multiplication and addition) needed to compute the product of all subsets of size n-1?

Clearly you can do it in (n-1)*n multiplications.

To clarify, the output is n different products and the only operations apart from reading and writing to memory allowed are multiplication, addition and subtraction.

Example

If the input has three numbers 2,3,5, then the output is three numbers 15 = 3*5, 10 = 2*5 and 6 = 2*3.

Winning criterion

Answers should give an exact formula for the number of arithmetic operations their code will use in terms of n. To make life simple, I will just plug n = 1000 into your formula to judge its score. The lower the better.

If it's too hard to produce an exact formula for your code, you can just run it for n = 1000 and count the arithmetic operations in the code. An exact formula would be best however.

You should add your score for n=1000 to your answer for easy comparison.

Arthur

Posted 2017-06-18T17:45:19.093

Reputation: 473

Welcome to PPCG! Nice first question! – Erik the Outgolfer – 2017-06-18T17:47:37.403

@EriktheOutgolfer Thank you! – Arthur – 2017-06-18T17:49:31.333

Do the products have to be listed in any particular order? – xnor – 2017-06-18T17:56:40.940

@xnor No. Any order will do. – Arthur – 2017-06-18T17:58:54.987

4May we count multiplying by 1 as free? Otherwise I'd define a custom-multiplication function that does this. – xnor – 2017-06-18T18:22:20.393

For simplicity, all multiplications should be counted in the score, even by 1. – Arthur – 2017-06-18T18:29:27.053

3Would it be against the rules to do a whole bunch of multiplications in parallel by concatenating numbers together with sufficiently many spacer 0 digits? – xnor – 2017-06-18T18:37:50.107

@xnor Yes it would! – Arthur – 2017-06-18T18:39:15.347

Are operations like >= or max also banned? You might want to make this an atomic code golf so people don't try to use language features to get around your restrictions.

– xnor – 2017-06-18T18:48:04.340

@xnor I basically don't want to allow any operations that are going to substitute for the arithmetic operations via some trickery. I really want to know how you can do it using a small number of multiplications, subtractions and additions. – Arthur – 2017-06-18T18:50:50.307

1Do the operations such as + on indices count? If this is the case, does array indexing count as well? (since it is after all syntactic sugar for addition and dereferencing). – nore – 2017-06-18T18:56:44.087

I removed the restricted-code tag since my understanding is that atomic-code-golf supersedes it. – xnor – 2017-06-18T18:58:09.230

@nore You can look up/access/read the input elements for free. All the arithmetic operations in your code should be counted however. – Arthur – 2017-06-18T19:00:06.993

@Arthur does this mean a loop has a cost of the number of operations for the counter? – nore – 2017-06-18T19:03:27.987

Are architecture-specific tricks allowed in our answers? – Govind Parmar – 2017-06-18T19:56:13.033

2@nore OK I give in :) Just count arithmetic operations that involve the input in some way. – Arthur – 2017-06-18T20:12:11.813

1Clearly you can do it in (n-1)*n multiplications You mean (n-2)*n, right? – Luis Mendo – 2017-06-18T21:16:04.723

Does a SIMD operation count as a single operation? – Veedrac – 2017-06-19T02:13:19.317

@Veedrac Not for this challenge, sorry. – Arthur – 2017-06-19T08:25:45.000

Maybe I'm being thick, but what does it mean to multiply two or more sets? No definition of set that I know has a multiply operation. Perhaps you could give a worked example for what the result would be for a small input array? – Toby Speight – 2017-06-19T17:05:45.213

@TobySpeight If the input has three number 1,2,3, then the output is three numbers 6 = 2*3,3 = 1*3 and 2 = 1*2. Does that help? – Arthur – 2017-06-19T17:15:41.933

1Oh, so you're asking for a set of outputs - no wonder I was confused. The way it's worded in the question suggests you are looking for (2,3) * (1,3) * (1,2) - perhaps it could be edited to include your helpful example? – Toby Speight – 2017-06-19T17:21:31.257

If this was a simple code golf, I would submit this Perl 6 one {.combinations(.end).map: *.reduce: &prefix:<[*]>} – Brad Gilbert b2gills – 2017-06-19T18:26:12.273

Answers

25

Python, 3(n-2) operations, score = 2994

l = list(map(float, input().split()))
n = len(l)

left = [0] * len(l)
right = [0] * len(l)
left[0] = l[0]
right[-1] = l[-1]
for i in range(1,len(l)-1):
  left[i] = l[i] * left[i - 1]
  right[-i-1] = l[-i-1] * right[-i]

result = [0] * len(l)
result[-1] = left[-2]
result[0] = right[1]
for i in range(1, len(l) - 1):
  result[i] = left[i - 1] * right[i+1]

print(result)

The arrays left and right contain the cumulated products of the array from the left and from the right, respectively.

EDIT: Proof that 3(n-2) is the optimal number of operations needed for n >= 2, if we only use multiplication.

We will do that by induction; by the above algorithm, we just have to prove that for n >= 2, 3(n-2) is a lower bound on the number of multiplications needed.

For n = 2, we need at least 0 = 3(2-2) multiplications, so the result is trivial.

Let n > 2, and suppose for n - 1 elements, we need at least 3(n-3) multiplications. Consider a solution for n elements with k multiplications. Now, we remove the last of these elements as if it was 1, and simplify all multiplications directly by it. We also remove the multiplication leading to the product of all the other elements, since that one is not needed as it can never be used as an intermediate value to get the product of n-2 of the other elements, since division is not allowed. This leaves us with l multiplications, and a solution for n - 1 elements.

By induction hypothesis, we have l >= 3(n-3).

Now, let's have a look at how many multiplications were removed. One of them was the one leading to the product of all elements except the last. Moreover, the last element was used directly in at least two multiplications: if it was used in only one, then it was used when multiplying by an intermediate result consisting in some product of the other elements; let's say, without loss of generality, this this intermediate result included the first element in the product. Then, there is no way to get the product of all the elements but the first, since every product that contains the last element is either the last element, or contains the first element.

We thus have k >= l+3 >= 3(n-2), proving the claimed theorem.

nore

Posted 2017-06-18T17:45:19.093

Reputation: 436

8

This turns out very clean in Haskell: f l = zipWith (*) (scanl (*) 1 l) (scanr (*) 1 $ tail l).

– xnor – 2017-06-18T22:25:30.670

Comments are not for extended discussion; this conversation has been moved to chat.

– Dennis – 2017-06-21T03:19:12.437

12

Haskell, score 2994

group :: Num a => [a] -> [[a]]
group (a:b:t) = [a,b] : group t
group [a] = [[a]]
group [] = []

(%) :: (Num a, Eq a) => a -> a -> a
a % 1 = a
1 % b = b
a % b = a * b

prod_one_or_two :: (Num a, Eq a) => [a] -> a
prod_one_or_two [a, b] = a % b
prod_one_or_two [x] = x

insert_new_value :: (Num a, Eq a) => ([a], a) -> [a]
insert_new_value ([a, b], c) = [c % b, c % a]
insert_new_value ([x], c) = [c]

products_but_one :: (Num a, Eq a) => [a] -> [a]
products_but_one [a] = [1]
products_but_one l = 
    do combination <- combinations ; insert_new_value combination
    where 
        pairs = group l
        subresults = products_but_one $ map prod_one_or_two pairs
        combinations = zip pairs subresults

Try it online!

Say we're given the list [a,b,c,d,e,f,g,h]. We first group it into pairs [[a,b],[c,d],[e,f],[g,h]]. Then, we recurse on the half-size list pairs of their products to get subresults

[a*b, c*d, e*f, g*h] -> [(c*d)*(e*f)*(g*h), (a*b)*(e*f)*(g*h), (a*b)*(c*d)*(g*h), (a*b)*(c*d)*(e*f)]

If we take the first element (c*d)*(e*f)*(g*h), and multiply it by b and a respectively, we get the product of all but a and all but b. Doing this for every pair and recursive result with that pair missing, we get out final result. The odd-length case is specially handled by having the odd element passed unpaired to the recursive step, and the product of the remaining elements returned is the product without it.

The number of multiplications t(n) is n/2 for the pairing product, t(n/2) for the recursive call, and another n for the products with individual elements. This gives t(n) = 1.5 * n + t(n/2) for odd n. Using a more precise counts for odd n and ignoring multiplying with 1 for the base case gives score 2997 for n=1000.

xnor

Posted 2017-06-18T17:45:19.093

Reputation: 115 687

This is very nice. – Arthur – 2017-06-18T20:14:05.180

I think the reason why the score is 2995 and not 2994 as in my answer is that it computes the product of all numbers as well in the non-power of two case, which is later truncated. Maybe a careful handling in products_but_one' could avoid that by returning something of the correct length. – nore – 2017-06-18T20:19:59.740

@nore I found I had an extra multiplication in my count because I forgot the base case has a 1 that is free to multiply. I think the padding 1 didn't affect things, but I cleaned up my algorithm to not use them. – xnor – 2017-06-18T20:38:08.313

Does this code assume the input is integers? – None – 2017-06-19T03:52:55.143

@Lembik It does, but only in the optional type annotations. I'll change them all to float. – xnor – 2017-06-19T03:57:54.457

Surely you mean Num a => a, or at least Double. – Ørjan Johansen – 2017-06-19T04:02:26.413

@ØrjanJohansen Is this right? – xnor – 2017-06-19T04:11:28.987

Looks good, I forgot about the need for Eq. Although group doesn't need either constraint. – Ørjan Johansen – 2017-06-19T04:54:43.873

@ØrjanJohansen I was surprised Num doesn't imply Eq. Do you know why that is? – xnor – 2017-06-19T05:27:56.443

@xnor It used to, but was removed in GHC 7.4 because it rules out (or makes awkward) useful new instances. See libraries list thread.

– Ørjan Johansen – 2017-06-19T06:04:19.913

9

Haskell, score 9974

partition :: [Float] -> ([Float], [Float])
partition = foldr (\a (l1,l2) -> (l2, a:l1)) ([],[])

(%) :: Float -> Float -> Float
a % 1 = a
1 % b = b
a % b = a*b

merge :: (Float, [Float]) -> (Float, [Float]) -> (Float, [Float])
merge (p1,r1) (p2, r2) = (p1%p2, map(%p1)r2 ++ map(%p2)r1)

missing_products' :: [Float] -> (Float, [Float])
missing_products' [a] = (a,[1])
missing_products' l = merge res1 res2
    where
        (l1, l2) = partition l
        res1 = missing_products' l1
        res2 = missing_products' l2

missing_products :: [Float] -> [Float]
missing_products = snd . missing_products'

Try it online!

A divide-and-conquer strategy, very reminiscent of merge sort. Doesn't do any indexing.

The function partition splits the list into as-equal-as-possible halves by putting alternating elements on opposite sides of the partition. We recursively merge the results (p,r) for each of the halves, with r the list of products-with-one-missing, and p the overall product.

For the output for the full list, the missing element must be in one of the halves. The product with that element missing is a one-missing-product for the half it's in, multiplied by the full product for the other half. So, we multiply each product-with-one-missing by the full product of the other half and make a list of the results, as map(*p1)r2 ++ map(*p2)r1). This takes n multiplications, where n is the length. We also need to make a new full-product p1*p2 for future use, costing 1 more multiplication.

This gives the general recursion for for the number of operations t(n) with n even: t(n) = n + 1 + 2 * t(n/2). The odd one is similar, but one of the sublists is 1 larger. Doing out the recursion, we get n*(log_2(n) + 1) multiplications, though the odd/even distinction affects that exact value. The values up to t(3) are improved by not multiplying by 1 by defining a variant (%) of (*) that shortcuts the _*1 or 1*_ cases.

This gives 9975 multiplications for n=1000. I believe Haskell's laziness means the unused overall product in the outer layer is not computed for 9974; if I'm mistaken, I could omit it explicitly.

xnor

Posted 2017-06-18T17:45:19.093

Reputation: 115 687

You beat me by the timestamp one minute earlier. – nore – 2017-06-18T18:15:09.053

If it's hard to work out the formula exactly, feel free to just run it for n = 1000 and count the arithmetic operations in the code. – Arthur – 2017-06-18T18:18:07.673

Since our code is basically the same, I don't understand how you got to 9974 and not 9975 multiplications for n = 1000 (in the case of computing the overall product in the outer layer). Did you include a 1 in the input you used to test it? – nore – 2017-06-18T18:39:03.507

@nore You're right, I was off by one. I wrote code to do the recursion for the number of multiplication function calls. Counting calls directly would be more reliable -- does anyone know how I'd do that in Haskell? – xnor – 2017-06-18T18:43:38.470

@xnor From my very limited Haskell experience, I know no easy way to do that except by using unsafePerformIO. Maybe someone has a better idea for that than me, however.

– nore – 2017-06-18T18:46:45.670

1@xnor you can use trace from Debug.Trace with a catch-all | trace "call!" False = undefined guard, I think. But this uses unsafePerformIO under the hood, so it's not really that much of an improvement. – Soham Chowdhury – 2017-06-19T15:11:55.827

6

Haskell, score 2994

group :: [a] -> Either [(a, a)] (a, [(a, a)])
group [] = Left []
group (a : l) = case group l of
  Left pairs -> Right (a, pairs)
  Right (b, pairs) -> Left ((a, b) : pairs)

products_but_one :: Num a => [a] -> [a]
products_but_one [_] = [1]
products_but_one [a, b] = [b, a]
products_but_one l = case group l of
  Left pairs ->
    let subresults =
          products_but_one [a * b | (a, b) <- pairs]
    in do ((a, b), c) <- zip pairs subresults; [c * b, c * a]
  Right (extra, pairs) ->
    let subresult : subresults =
          products_but_one (extra : [a * b | (a, b) <- pairs])
    in subresult : do ((a, b), c) <- zip pairs subresults; [c * b, c * a]

Try it online!

How it works

This is a cleaned up version of xnor’s algorithm that deals with the odd case in a more straightforward way (edit: it looks like xnor has cleaned it up the same way):

[a, b, c, d, e, f, g] ↦
[a, bc, de, fg] ↦
[(bc)(de)(fg), a(de)(fg), a(bc)(fg), a(bc)(de)] by recursion ↦
[(bc)(de)(fg), a(de)(fg)c, a(de)(fg)b, a(bc)(fg)e, a(bc)(fg)d, a(bc)(de)g, a(bc)(de)f]

[a, b, c, d, e, f, g, h] ↦
[ab, cd, ef, gh] ↦
[(cd)(ef)(gh), (ab)(ef)(gh), (ab)(cd)(gh), (ab)(cd)(ef)] by recursion ↦
[(cd)(ef)(gh)b, (cd)(ef)(gh)a, (ab)(ef)(gh)d, (ab)(ef)(gh)c, (ab)(cd)(gh)f, (ab)(cd)(gh)e, (ab)(cd)(ef)h, (ab)(cd)(ef)g].

Anders Kaseorg

Posted 2017-06-18T17:45:19.093

Reputation: 29 242

"Given n numbers in an array (you can't assume they are integers), " We can't assume they are integers – None – 2017-06-19T05:15:03.467

5

O(n log n) operations, score = 9974

Works with a binary tree.

Python

l = list(map(int, input().split()))
n = len(l)

p = [0] * n + l
for i in range(n - 1, 1, -1):
  p[i] = p[i + i] * p[i + i+1]

def mul(x, y):
  if y == None:
    return x
  return x * y

r = [None] * n + [[None]] * n
for i in range(n - 1, 0, -1):
  r[i] = [mul(p[i + i + 1], x) for x in r[i + i]] + [mul(p[i + i], x) for x in r[i + i + 1]]

u = r[1]
j = 1
while j <= n:
  j += j
print(u[n+n-j:] + u[:n+n-j])

This also requires list addition operations, and some arithmetic on numbers that are not the input values; not sure if that counts. The mul function is there to save n operations for the base case, to avoid wasting them by multiplying by 1. In any case, this is O(n log n) operations. The exact formula is, if only counting arithmetic operations on input numbers, with j = floor(log_2(n)): j * (2^(j + 1) - n) + (j + 1) * (2 * n - 2^(j + 1)) - 2.

Thanks to @xnor for saving one operation with the idea of not computing the outer product!

The last part is to output the products in order of the missing term.

nore

Posted 2017-06-18T17:45:19.093

Reputation: 436

If it's hard to work out the formula exactly, feel free to just run it for n = 1000 and count the arithmetic operations in the code. – Arthur – 2017-06-18T18:17:44.903

I counted 10975 operations...? – HyperNeutrino – 2017-06-18T18:29:09.810

p[i] = p[i + i] * p[i + i+1] isn't counted – HyperNeutrino – 2017-06-18T18:29:49.127

This is about n log2 n + n operations (which is O(nlogn) btw – HyperNeutrino – 2017-06-18T18:30:21.437

@HyperNeutrino the operations in p[i] = p[i + i] * p[i + i + 1] should be saved by the multiplication optimization. I might have counted one too many, however. – nore – 2017-06-18T18:34:41.463

Can't you assign a variable to i + i in the second for loop so you save on operations? – Zacharý – 2017-06-18T18:51:36.033

@ZacharyT I could, but I didn't include the operations on indices (to oppose to input numbers) in the count. I could do that of course, but I don't think this is the spirit of the challenge. – nore – 2017-06-18T18:54:38.403

Did you count the n-1 in the range calls? – Zacharý – 2017-06-18T19:06:50.123

@ZacharyT same, I didn't count this operation as it is only an operation on indices; I am going to wait for an explicit answer on this point. – nore – 2017-06-18T19:08:10.800

3

O((n-2)*n) = O(n2): Trivial Solution

This is just the trivial solution that multiplies together each of the subsets:

Python

def product(array): # Requires len(array) - 1 multiplication operations
    if not array: return 1
    result = array[0]
    for value in array[1:]:
        result *= value
    return result

def getSubsetProducts(array):
    products = []
    for index in range(len(array)): # calls product len(array) times, each time calling on an array of size len(array) - 1, which means len(array) - 2 multiplication operations called len(array) times
        products.append(product(array[:index] + array[index + 1:]))
    return products

Note that this also requires n list-addition operations; not sure if that counts. If that is not allowed, then product(array[:index] + array[index + 1:]) can be replaced to product(array[:index]) * product(array[index + 1:]), which changes the formula to O((n-1)*n).

HyperNeutrino

Posted 2017-06-18T17:45:19.093

Reputation: 26 575

You can add your own score to the answer. 998*1000 in this case. – Arthur – 2017-06-18T18:03:41.527

Doesn't need your product function O(n) operations? one for each element in the array (althought this can easily be changed to O(n-1)) – Roman Gräf – 2017-06-18T18:34:41.073

@RomanGräf True. I'll change it to O(n-1) but thanks for pointing that out. – HyperNeutrino – 2017-06-18T18:42:35.523

This has been changed to [tag:atomic-code-golf]... – Erik the Outgolfer – 2017-06-18T19:33:22.573

@EriktheOutgolfer What does that make my score now? Unless I am being blatantly stupid, don't the tag and the specs contradict each other now? – HyperNeutrino – 2017-06-18T19:44:32.183

3

Python, 7540

A tripartite merge strategy. I think I can do even better than this, with a yet large merge. It's O(n log n).

EDIT: Fixed a miscount.

count = 0
def prod(a, b):
    if a == 1: return b
    if b == 1: return a
    global count
    count += 1
    return a * b

def tri_merge(subs1, subs2, subs3):
    total1, missing1 = subs1
    total2, missing2 = subs2
    total3, missing3 = subs3

    prod12 = prod(total1, total2)
    prod13 = prod(total1, total3)
    prod23 = prod(total2, total3)

    new_missing1 = [prod(m1, prod23) for m1 in missing1]
    new_missing2 = [prod(m2, prod13) for m2 in missing2]
    new_missing3 = [prod(m3, prod12) for m3 in missing3]

    return prod(prod12, total3), new_missing1 + new_missing2 + new_missing3

def tri_partition(nums):
    split_size = len(nums) // 3
    a = nums[:split_size]
    second_split_length = split_size + (len(nums) % 3 == 2)
    b = nums[split_size:split_size + second_split_length]
    c = nums[split_size + second_split_length:]
    return a, b, c

def missing_products(nums):
    if len(nums) == 1: return nums[0], [1]
    if len(nums) == 0: return 1, []
    subs = [missing_products(part) for part in tri_partition(nums)]
    return tri_merge(*subs)

def verify(nums, res):
    actual_product = 1
    for num in nums:
        actual_product *= num
    actual_missing = [actual_product // num for num in nums]
    return actual_missing == res[1] and actual_product == res[0]

nums = range(2, int(input()) + 2)
res = missing_products(nums)

print("Verified?", verify(nums, res))
if max(res[1]) <= 10**10: print(res[1])

print(len(nums), count)

The relevant function is missing_products, which gives the overall product and all of the ones with a missing element.

isaacg

Posted 2017-06-18T17:45:19.093

Reputation: 39 268

did you count the multiplications in tri_merge? Also you can replace the 2 * split_size + ... in tri_partition with split_size + split_size + .... – Roman Gräf – 2017-06-18T19:51:54.927

@RomanGräf I restructured it as per your suggestion. – isaacg – 2017-06-19T00:03:46.163

1

C++, score: 5990, O([2NlogN]/3)

This implementation uses a binary tree look up table. My first implementation was O(NlogN), but a last minute optimization, which looks up the product of all array elements minus a pair, + 2 multiplications saved the day. I think this could still be optimized a bit further, maybe another 16%...

I've left some debugging traces in, only because it's easier to delete them than to rewrite them :)

[Edit] the actual complexity is measured at O([2NlogN]/3) for 100. It is actually a bit worse than O(NlogN) for small sets, but tends towards O([NlogN]/2) as the array grows very large O(0.57.NlogN) for a a set of 1 million elements.

#include "stdafx.h"
#include <vector>
#include <iostream>
#include <random>
#include <cstdlib>

using DataType = long double;

using DataVector = std::vector<DataType>;

struct ProductTree
{
    std::vector<DataVector> tree_;
    size_t ops_{ 0 };

    ProductTree(const DataVector& v) : ProductTree(v.begin(), v.end()) {}
    ProductTree(DataVector::const_iterator first, DataVector::const_iterator last)
    {
        Build(first, last);
    }

    void Build(DataVector::const_iterator first, DataVector::const_iterator last)
    {
        tree_.emplace_back(DataVector(first, last));

        auto size = std::distance(first, last);
        for (auto n = size; n >= 2; n >>= 1)
        {
            first = tree_.back().begin();
            last = tree_.back().end();

            DataVector v;
            v.reserve(n);
            while (first != last) // steps in pairs
            {
                auto x = *(first++);
                if (first != last)
                {
                    ++ops_;
                    x *= *(first++); // could optimize this out,small gain
                }
                v.push_back(x);
            }
            tree_.emplace_back(v);
        }
    }

    // O(NlogN) implementation... 
    DataVector Prod()
    {
        DataVector result(tree_[0].size());
        for (size_t i = 0; i < tree_[0].size(); ++i)
        {
            auto depth = tree_.size() - 1;
            auto k = i >> depth;
            result[i] = ProductAtDepth(i, depth);
        }
        return result;
    }

    DataType ProductAtDepth(size_t index, size_t depth) 
    {
        if (depth == 0)
        {
            return ((index ^ 1) < tree_[depth].size())
                ? tree_[depth][index ^ 1]
                : 1;
        }
        auto k = (index >> depth) ^ 1;

        if ((k < tree_[depth].size()))
        {
            ++ops_;
            return tree_[depth][k] * ProductAtDepth(index, depth - 1);
        }
        return ProductAtDepth(index, depth - 1);
    }    

    // O([3NlogN]/2) implementation... 
    DataVector Prod2()
    {
        DataVector result(tree_[0].size());
        for (size_t i = 0; i < tree_[0].size(); ++i)    // steps in pairs
        {
            auto depth = tree_.size() - 1;
            auto k = i >> depth;
            auto x = ProductAtDepth2(i, depth);
            if (i + 1 < tree_[0].size())
            {
                ops_ += 2;
                result[i + 1] = tree_[0][i] * x;
                result[i] = tree_[0][i + 1] * x;
                ++i;
            }
            else
            {
                result[i] = x;
            }
        }
        return result;
    }

    DataType ProductAtDepth2(size_t index, size_t depth)
    {
        if (depth == 1)
        {
            index = (index >> 1) ^ 1;
            return (index < tree_[depth].size())
                ? tree_[depth][index]
                : 1;
        }
        auto k = (index >> depth) ^ 1;

        if ((k < tree_[depth].size()))
        {
            ++ops_;
            return tree_[depth][k] * ProductAtDepth2(index, depth - 1);
        }
        return ProductAtDepth2(index, depth - 1);
    }

};


int main()
{
    //srand(time());

    DataVector data;
    for (int i = 0; i < 1000; ++i)
    {
        auto x = rand() & 0x3;          // avoiding overflow and zero vaolues for testing
        data.push_back((x) ? x : 1);
    }

    //for (int i = 0; i < 6; ++i)
    //{
    //  data.push_back(i + 1);
    //}

    //std::cout << "data:[";
    //for (auto val : data)
    //{
    //  std::cout << val << ",";
    //}
    //std::cout << "]\n";

    ProductTree pt(data);
    DataVector result = pt.Prod2();

    //std::cout << "result:[";
    //for (auto val : result)
    //{
    //  std::cout << val << ",";
    //}
    //std::cout << "]\n";
    std::cout << "N = " << data.size() << " Operations :" << pt.ops_ << '\n';

    pt.ops_ = 0;
    result = pt.Prod();

    //std::cout << "result:[";
    //for (auto val : result)
    //{
    //  std::cout << val << ",";
    //}
    //std::cout << "]\n";

    std::cout << "N = " << data.size() << " Operations :" << pt.ops_ << '\n';

    return 0;
}

I'm adding @nore's algorithm, for completeness. It's really nice, and is the fastest.

class ProductFlat
{
private:
    size_t ops_{ 0 };

    void InitTables(const DataVector& v, DataVector& left, DataVector& right)
    {
        if (v.size() < 2)
        {
            return;
        }

        left.resize(v.size() - 1);
        right.resize(v.size() - 1);

        auto l = left.begin();
        auto r = right.rbegin();
        auto ol = v.begin();
        auto or = v.rbegin();

        *l = *ol++;
        *r = *or++;
        if (ol == v.end())
        {
            return;
        }

        while (ol + 1 != v.end())
        {
            ops_ += 2;
            *l = *l++ * *ol++;
            *r = *r++ * *or++;
        }
    }

public:
    DataVector Prod(const DataVector& v)
    {
        if (v.size() < 2)
        {
            return v;
        }

        DataVector result, left, right;
        InitTables(v, left, right);

        auto l = left.begin();
        auto r = right.begin();
        result.push_back(*r++);
        while (r != right.end())
        {
            ++ops_;
            result.push_back(*l++ * *r++);
        }
        result.push_back(*l++);
        return result;
    }

    auto Ops() const
    {
        return ops_;
    }
};

Michaël Roy

Posted 2017-06-18T17:45:19.093

Reputation: 111

1

dc, score 2994

#!/usr/bin/dc -f

# How it works:
# The required products are
#
#   (b × c × d × e × ... × x × y × z)
# (a) × (c × d × e × ... × x × y × z)
# (a × b) × (d × e × ... × x × y × z)
# ...
# (a × b × c × d × e × ... × x) × (z)
# (a × b × c × d × e × ... × x × y)
#
# We calculate each parenthesised term by
# multiplying the one above (on the left) or below
# (on the right), for 2(n-2) calculations, followed
# by the n-2 non-parenthesised multiplications
# giving a total of 3(n-2) operations.

# Read input from stdin
?

# We will store input values into stack 'a' and
# accumulated product into stack 'b'.  Initialise
# stack b with the last value read.
sb

# Turnaround function at limit of recursion: print
# accumulated 'b' value (containing b..z above).
[Lbn[ ]nq]sG

# Recursive function - on the way in, we stack up
# 'a' values and multiply up the 'b' values.  On
# the way out, we multiply up the 'a' values and
# multiply each by the corresponding 'b' value.
[dSalb*Sb
z1=G
lFx
dLb*n[ ]n
La*]dsFx

# Do the last a*b multiplication
dLb*n[ ]n

# And we have one final 'a' value that doesn't have a
# corresponding 'b':
La*n

I'm assuming that the integer comparison z1= (which terminates the recursion when we reach the last value) is free. This is equivalent to the likes of foreach in other languages.

Demonstrations

for i in '2 3 5' '2 3 5 7' '0 2 3 5' '0 0 1 2 3 4'
do printf '%s => ' "$i"; ./127147.dc <<<"$i"; echo
done
2 3 5 => 15 10 6
2 3 5 7 => 105 70 42 30
0 2 3 5 => 30 0 0 0
0 0 1 2 3 4 => 0 0 0 0 0 0

A demo with large and small inputs:

./127147.dc <<<'.0000000000000000000542101086242752217003726400434970855712890625 1 18446744073709551616'
18446744073709551616 1.0000000000000000000000000000000000000000000000000000000000000000 .0000000000000000000542101086242752217003726400434970855712890625

Toby Speight

Posted 2017-06-18T17:45:19.093

Reputation: 5 058