#include <cassert>
#include <cstddef>
#include <iostream>
#include <limits>
#include <new>

// Previously, you wrote a function that calculates
// the sum of the kth powers of the positive divisors
// of an integer n.
//
// For example, the positive divisors of 12 are
//
//     1, 2, 3, 4, 6, 12.
//
// Thus,
//
//     sigma(12, 0) = 6
//
//     sigma(12, 1) = 1 + 2 + 3 + 4 + 6 + 12
//                  = 28
//
//     sigma(12, 2) = 1 + 4 + 9 + 16 + 36 + 144
//                  = 210
//
// Suppose we want to calculate these sums for
// every integer from 1 through n.
//
// We could call our previous function n times,
// but this would repeat a considerable amount
// of work.
//
// Instead, create an array a_sigma with n + 1
// entries, initially all zero.
//
// For each divisor d from 1 through n:
//
//   1. Calculate d raised to the power k.
//
//   2. Add this value to every array entry
//      whose index is a multiple of d.
//
// For example, when d == 3, add 3^k to
//
//     a_sigma[3], a_sigma[6], a_sigma[9], ...
//
// At the end, a_sigma[m] contains the sum of
// the kth powers of all positive divisors of m.
//
// We will use a_sigma[0] = 0 by convention.
// The usual divisor-sum function is defined
// only for positive integers.
//
// The function must dynamically allocate an
// array, populate it, and return its address.
//
// The calling function is responsible for
// releasing that memory using delete[].
//
// Assume that all calculated powers and sums
// fit in an unsigned long.

// Function declarations
unsigned long *sigma_array(
    std::size_t   const n,
    unsigned long const k = 1
);

unsigned long pow(
    unsigned long const base,
    unsigned long       exponent
);

int main();

// Function definitions
int main() {
    std::size_t const N{ 40 };

    unsigned long *a_sigma{
        sigma_array( N, 0 )
    };

    for ( std::size_t m{ 0 }; m <= N; ++m ) {
        std::cout << "sigma(" << m << ", 0) = "
                  << a_sigma[m] << std::endl;
    }

    std::cout << std::endl;

    // Release the first array before allocating
    // another one. The pointer is immediately
    // assigned the address of the new array.

    delete[] a_sigma;

    a_sigma = sigma_array( N );

    for ( std::size_t m{ 0 }; m <= N; ++m ) {
        std::cout << "sigma(" << m << ", 1) = "
                  << a_sigma[m] << std::endl;
    }

    delete[] a_sigma;
    a_sigma = nullptr;

    return 0;
}

// Allocate and populate an array containing
// sigma(m, k) for every m from 0 through n.
//
// The caller owns the returned array and
// must eventually release it using delete[].

unsigned long *sigma_array(
    std::size_t   const n,
    unsigned long const k
) {
    // Enter your implementation here.

    return nullptr;
}

// Calculate base raised to exponent using
// exponentiation by squaring.
//
// Assume that the result fits in an unsigned long.

unsigned long pow(
    unsigned long const base,
    unsigned long       exponent
) {
    unsigned long product{ 1 };
    unsigned long current_power{ base };

    while ( exponent > 0 ) {
        if ( (exponent & 1UL) != 0 ) {
            product *= current_power;
        }

        exponent >>= 1;

        // No further squaring is needed after
        // the last bit has been processed.
        if ( exponent > 0 ) {
            current_power *= current_power;
        }
    }

    return product;
}


















































// Hint:
//
// Start by allocating and initializing an array:
//
//     unsigned long *a_values{
//         new unsigned long[n + 1]{}
//     };
//
// The empty braces initialize every entry to zero.
//
// Then consider each possible divisor d.
//
// Calculate
//
//     unsigned long const contribution{ pow(d, k) };
//
// once, and add it to all entries whose indices
// are multiples of d:
//
//     d, 2*d, 3*d, ...
//
// Rather than multiplying a counter by d,
// start at d and repeatedly add d.
//
// Be careful not to increment the index beyond
// the largest value representable by std::size_t.


















































// Here is one solution:

/*
unsigned long *sigma_array(
    std::size_t   const n,
    unsigned long const k
) {
    assert(
        n < std::numeric_limits<std::size_t>::max()
    );

    unsigned long *a_values{
        new unsigned long[n + 1]{}
    };

    for ( std::size_t divisor{ 1 };
          divisor <= n;
          ++divisor ) {
        unsigned long const contribution{
            pow(
                static_cast<unsigned long>( divisor ),
                k
            )
        };

        for ( std::size_t multiple{ divisor };
              multiple <= n; ) {
            a_values[multiple] += contribution;

            // Stop before an addition that could
            // overflow or exceed n.

            if ( divisor > n - multiple ) {
                break;
            }

            multiple += divisor;
        }
    }

    return a_values;
}
*/