euler_387_harshad.zig

euler_387_harshad.zig

// A Harshad or Niven number is a number that is divisible by the sum
// of its digits. 201 is a Harshad number because it is divisible by 3
// (the sum of its digits.)
//
// Right truncatable Harshad number: Recursively truncating the last
// digit, always results in a Harshad number. When we truncate the
// last digit from 201, we get 20, which is a Harshad number. When we
// truncate the last digit from 20, we get 2, which is also a Harshad
// number.
//
// Also: 201 / 3 = 67 which is prime.
//
// *Strong Harshad number*: A Harshad number that, when divided by the
//  sum of its digits, results in a prime.
//
// *Strong Right truncatable harshad prime*: A strong Harshad number
// that is also right truncatable. Such a number is 2011 which is
// prime and when we truncate the last digit from it we get 201.
//
// You are given that the sum of the strong, right truncatable Harshad
// primes less than 10000 is 90619.
//
// Find the sum of the strong, right truncatable Harshad primes less
// than 10^14.
const std = @import("std");
const assert = @import("std").debug.assert;
const warn = @import("std").debug.warn;

fn log(args: ...) void {
    const fmt: []const u8 = ("{} " ** args.len) ++ "\n";
    warn(fmt, args);
}

fn sum_of_digits(num: u64) u64 {
    var cur: u64 = num;
    var sum: u64 = 0;
    while (cur > 0) {
        sum += cur % 10;
        cur /= 10;
    }
    return sum;
}

test "sum_of_digits" {
    assert(sum_of_digits(201) == 3);
    assert(sum_of_digits(2011) == 4);
    assert(sum_of_digits(2010) == 3);
    assert(sum_of_digits(9999) == 9 * 4);
}

fn is_h_num(num: u64) bool {
    if (num == 0) {
        return false;
    } else if (num < 10) {
        return true;
    } else {
        var sod: u64 = sum_of_digits(num);
        return num % sod == 0;
    }
}

test "is_h_num" {
    assert(is_h_num(201) == true);
    assert(is_h_num(20) == true);
    assert(is_h_num(2) == true);
    assert(is_h_num(22) == false);
}

fn is_rth_num(num: u64) bool {
    if (num == 0) {
        return false;
    }
    var cur: u64 = num;
    while (cur > 10) {
        cur /= 10;
        if (!is_h_num(cur)) {
            return false;
        }
    }
    return true;
}

test "is_rth_num" {
    assert(is_rth_num(2011) == true);
}

fn hash_u64(x: u64) u32 {
    return @truncate(u32, x);
}

fn eql_u64(a: u64, b: u64) bool {
    return a == b;
}

const HashSet = std.HashMap(u64, bool, hash_u64, eql_u64);
const IntList = std.ArrayList(u64);

const HarshadSolver = struct {
    _alloc1: *std.heap.DirectAllocator,
    _alloc2: *std.heap.ArenaAllocator,

    prime_list: IntList,
    prime_set: HashSet,
    primes: []u64,

    pub fn init(limit: u32) !HarshadSolver {
        var direct = std.heap.DirectAllocator.init();
        var arena = std.heap.ArenaAllocator.init(&direct.allocator);
        var allocator = &arena.allocator;

        var prime_list = IntList.init(allocator);
        var prime_set = HashSet.init(allocator);

        try prime_seive(&prime_list, &prime_set, limit);

        return HarshadSolver{
            ._alloc1 = &direct,
            ._alloc2 = &arena,
            .prime_list = prime_list,
            .prime_set = prime_set,
            .primes = prime_list.toSlice(),
        };
    }

    pub fn deinit(self: HarshadSolver) void {
        log("\n\n...", "asd1", "\n");
        self.prime_list.deinit();
        log("\n\n...", "asd2", "\n");
        self.prime_set.deinit();
        log("\n\n...", "asd3", "\n");
        self._alloc2.deinit();
        log("\n\n...", "asd4", "\n");
        self._alloc1.deinit();
        log("\n\n...", "asd5", "\n");
    }

    pub fn is_srth_num(self: HarshadSolver, num: u64) bool {
        var sod: u64 = sum_of_digits(num);
        var next = num / sod;
        return self.prime_set.contains(next) and is_rth_num(num);
    }
};

fn is_srth_num(primes: *HashSet, num: u64) bool {
    var sod: u64 = sum_of_digits(num);
    var next = num / sod;
    return primes.contains(next) and is_rth_num(num);
}

test "is_srth_num" {
    const solver = try HarshadSolver.init(10000);

    assert(solver.is_srth_num(201) == true);

    var buf: [1024 * 1024 * 4]u8 = undefined;
    var fixed_buf = std.heap.FixedBufferAllocator.init(buf[0..]);
    var allocator = &fixed_buf.allocator;

    var prime_list = IntList.init(allocator);
    defer prime_list.deinit();

    var prime_set = HashSet.init(allocator);
    defer prime_set.deinit();

    try prime_seive(&prime_list, &prime_set, 10000);

    var primes: []u64 = prime_list.toSlice();

    assert(is_srth_num(&prime_set, 201) == true);
    log("\n\n...", "asdf", "\n");
    solver.deinit();
    log("\n\n...", "asdf", "\n");
}

// fn is_prime(primes: []u64, maybe_p: u64) bool:
//     var max_p: u64 = std.math.sqrt(maybe_p)
//     for p in primes:
//         if p > max_p:
//             return true
//         if maybe_p % p == 0:
//             return false
//     return true

fn is_prime(primes: []u64, maybe_p: u64) bool {
    var max_p: u64 = std.math.sqrt(maybe_p);
    for (primes) |p| {
        if (p > max_p) {
            return true;
        }
        if (maybe_p % p == 0) {
            return false;
        }
    }
    return true;
}

test "is_prime" {
    var primes = []u64{ 2, 3, 5, 7, 11, 13, 17, 19, 23 };
    assert(is_prime(primes[0..], 29));
    assert(!is_prime(primes[0..], 4));
    assert(!is_prime(primes[0..], 13 * 3));
}

// fn prime_seive(prime_list: *IntList, prime_set: *HashSet, limit: u64) !void:
//     assert(prime_list.count() == 0)
//     try prime_list.append(2)

//     var maybe_p := 3
//     while maybe_p < limit:
//         var primes = prime_list.toSlice()
//         if is_prime(primes, maybe_p):
//             try prime_list.append(maybe_p)
//             _ = try prime_set.put(maybe_p, true)
//         maybe_p += 2

fn prime_seive(prime_list: *IntList, prime_set: *HashSet, limit: u64) !void {
    assert(prime_list.count() == 0);
    try prime_list.append(2);

    var maybe_p: u64 = 3;
    while (maybe_p < limit) {
        if (is_prime(prime_list.toSlice(), maybe_p)) {
            try prime_list.append(maybe_p);
            _ = try prime_set.put(maybe_p, true);
        }
        maybe_p += 2;
    }
}

fn solve(limit: u64) !u64 {
    var direct = std.heap.DirectAllocator.init();
    defer direct.deinit();

    var arena = std.heap.ArenaAllocator.init(&direct.allocator);
    defer arena.deinit();

    const allocator = &arena.allocator;

    var prime_list = IntList.init(allocator);
    defer prime_list.deinit();

    var prime_set = HashSet.init(allocator);
    defer prime_set.deinit();

    try prime_seive(&prime_list, &prime_set, limit);
    var primes: []u64 = prime_list.toSlice();

    var total: u64 = 0;
    for (primes) |p| {
        if (is_srth_num(&prime_set, p)) {
            // log(total, "+=", p);
            total += p;
        }
    }
    return total;
}

pub fn main() void {
    log(solve(10000));
}