Smallest and largest 100-bit square with maximum Hamming weight

Question

Introduction

In the binary representation of a \$k-\$bit \$(k>1)\$ square number \$n^2, \$ not all bits can be \$1\$. The maximum number of \$1\$-bits (Hamming weight, \$B\$ for short) is given in A357658, but only up to \$k=95\$ bits.

An example of a \$99\$-bit square is \$B(796095014224053^2)=87\$. The records in A357304 show that \$B(m^2)<88\$ for a \$100\$-bit square \$m^2\$, and that the computation is feasible.

Challenge

Your task is to write a fastest-code program which determines

• \$B_{max}(n^2)\$
• \$n_{min}\$
• \$n_{max}\$

in the range of \$​​100\$-bit squares \$n^2\$ with \$796131459065722\le n\le 1125899906842623\$

with the smallest number \$n_{min}\$ and the largest number \$n_{max}\$ such that \$B(n_{min}^2)=B(n_{max}^2)=B_{max}\$.

A brute force calculation visiting all squares \$n^2\$ in the given range is possible, but takes a very long time, even if the multiplication is implemented efficiently with the instruction set for the available hardware. See this question on Stack Overflow.

Is it possible to avoid actually performing all the multiplications? The program that determines the 3 numbers \$B_{max}, n_{min}, n_{max}\$ in the shortest run time wins.

Test

To test the order of magnitude of a program's speed, it can be applied to the interval \$1085300000000000\le n \le 1085312000000000\$. The interval contains \$2\$ numbers with \$82\$ \$1\$-bits. The length of this interval is approximately \$1/27500\$ of the total range to be searched. You can then decide if you want to try the full range.

The results will be used in a planned extension of the OEIS sequences A357658, A357659, and A357660, giving credit to the authors of the fastest programs.

Answer 1

Rust w/ Rayon -- 0.9s runtime, 6.7s cputime for whole range

(1/99) a=0 1 (best: 0)
(2/49) b=2 4 (best: 50)
(3/27) a=1 49 (best: 72)
(4/23) b=3 188 (best: 76)
(5/21) a=2 1176 (best: 78)
(6/17) b=4 4332 (best: 82)
(7/15) a=3 18424 (best: 84)
(8/15) b=5 65220 (best: 84)
(9/15) a=4 211876 (best: 84)
(10/15) b=6 721396 (best: 84)
(11/15) a=5 1906884 (best: 84)
(12/15) b=7 6249964 (best: 84)
(13/14) a=6 13983816 (best: 85)
(14/14) b=8 44156732 (best: 85)

B_max: 85
n_min: 933544295655093
n_max: 1125891114428899

Found values: 1121216655081291 1043134657624555 1120665056942829 1116514939453621 933544295655093 1125074802425067 1124798782678823 1125891114428899 1090147221387163

100 85 933544295655093 1125891114428899

________________________________________________________
Executed in  909.41 millis    fish           external
   usr time    6.70 secs      0.00 micros    6.70 secs
   sys time    0.02 secs    516.00 micros    0.02 secs

Output of time max_square 100

Tested on my laptop with 4 cores (8 virtual) @ ~4.5 GHz

use std::{str::FromStr, collections::HashSet};
use rayon::prelude::*;
use std::sync::{Arc, Mutex};
use std::sync::atomic::{AtomicU32, Ordering};
use rayon::iter;

fn solve(n: u32, verbose: bool) {
    if n < 6 {
        match n {
            1 => {println!("1 1 1 1")},
            2 => {println!("None")},
            3 => {println!("3 1 2 2")},
            4 => {println!("4 2 3 3")},
            5 => {println!("5 3 5 5")},
            _ => {println!("Error");}
        }
        return;
    }
    let mut best = Solution::new(n);
    let mut a = 0;
    let mut b = 2;
    while best.get_best() + a + b <= n {
        let cost_a = calculate_cost_a(n, a);
        let cost_b = calculate_cost_b(n, b);
        if cost_a < cost_b {
            if verbose {eprintln!("({}/{}) a={} {} (best: {})", a+b-1, n - best.get_best() - 1, a, cost_a, best.get_best());}
            step_a(n, a, &mut best);
            a += 1;
        } else {
            if verbose {eprintln!("({}/{}) b={} {} (best: {})", a+b-1, n - best.get_best() - 1, b, cost_b, best.get_best());}
            step_b(n, b, &mut best);
            b += 1;
        }
    }
    let solutions = best.solutions.lock().unwrap();
    if verbose {
        eprintln!("");
        eprintln!("B_max: {}", best.get_best());
        eprintln!("n_min: {}", solutions.iter().min().expect("No solutions"));
        eprintln!("n_max: {}", solutions.iter().max().expect("No solutions"));
        eprintln!("");
        eprint!("Found values: ");
        for z in solutions.iter() {
            eprint!("{} ", z);
        }
        eprintln!("");
        eprintln!("");
    }
    println!("{} {} {} {}", n, best.get_best(), solutions.iter().min().expect("No solutions"), solutions.iter().max().expect("No solutions"));
}

fn calculate_cost_a(n: u32, a: u32) -> i64 {
    let a = a as i64;
    let n = n as i64;
    binomial((n-1)/2, a)
}

fn step_a(n: u32, a: u32, best: &Solution) {
    let h = n / 2;
    iter::split(
        Combinator::with_prefix_postfix(n, a, 1, 1, (1 << h) - 1, h),
        split_combinator).flatten_iter()
        .for_each(|z|{
            best.update(isqrt(z));
        }
    );
}

fn calculate_cost_b(n: u32, b: u32) -> i64 {
    let b = b as i64;
    let n = n as i64;
    (1..b/2 + 1).map(|i| (1<<i) * binomial(n/2 - 1 - 2*i, b - 2*i)).sum::<i64>() * 2
}

fn step_b(n: u32, b: u32, best: &Solution) {
    let m = (1 << ((n + 1) / 2)) - 1;
    for k in 0..b/2 {
        let w = ((n + 1) / 2) - 3 - 2*k;
        iter::split(Combinator::new(w, b - 2 - 2*k, 1, 3), split_combinator).flatten_iter()
            .for_each(|a| {
                let x = modsqrt(w + 3, a) << k;
                let mut z = 1 << ((n + 1) / 2);
                let step_size = 1 << (((n + 1) / 2) - 1 - k);
                while z > 0 {
                    z -= step_size;
                    best.update((z + x)&m);
                    best.update((z - x)&m);
                }
            }
        );
    }
}

fn modsqrt(n: u32, a: i128) -> i128 {
    let mut x = 1;
    let m = (1<<n) - 1;
    for _ in 3..n {
        x = (((x*x - a) >> 1) + x) & m;
    }
    x
}

fn isqrt(n: i128) -> i128 {
    let mut x = (n as f64).sqrt() as i128;
    if x*x <= n {
        x += 1;
        while x*x <= n {
            x += 1;
        }
        x-1
    } else {
        x -= 1;
        while x*x > n {
            x -= 1;
        }
        x
    }
}

fn binomial(n: i64, k: i64) -> i64 {
    let k = std::cmp::min(k, n-k);
    if k < 0 {return 0;}
    if k == 0 {return 1;}
    let mut c: i64 = 1;
    for i in 0..(k as i64) {
        c = (c * (n-i)) / (i + 1);
    }
    c
}

struct Combinator {
    curr: i128,
    mask: i128,
    offset: u32,
    min: i128,
    max: i128,
    k: u32,
}

impl Combinator {
    fn new(width: u32, num_zeros: u32, mask: i128, offset: u32) -> Combinator {
        Combinator { curr: (1 << width) - (1 << num_zeros), mask: mask, offset: offset, min: 0, max: binomial(width as i64, num_zeros as i64) as i128, k: width - num_zeros }
    }

    fn with_prefix_postfix(full_width: u32, num_zeros: u32, prefix: i128, prefix_width: u32, postfix: i128, postfix_width: u32) -> Combinator {
        let width = full_width - prefix_width - postfix_width;
        Combinator::new(width, num_zeros, (prefix << (width + postfix_width)) + postfix, postfix_width)
    }
}

// calculates the next combination, in reverse lexicographic order
fn next_combination(n: i128) -> i128 {
    if n&1 == 0 {return n - (((n^(n-1))+1)>>2);}
    let a = (n^(n+1))>>1;
    let na = n-a;
    if na == 0 {return na;}
    let b = ((na ^ (na-1))>>2) - a;
    if b == 0 {return na - 1;}
    let b = b >> b.trailing_zeros();
    na - 1 - b
}

impl Iterator for Combinator {
    type Item = i128;

    fn next(&mut self) -> Option<Self::Item> {
        if self.min == self.max {return None;}
        let curr = self.curr;
        self.max -= 1;
        self.curr = next_combination(curr);
        Some(self.mask | (curr << self.offset))
    }
}

// combinatorial number system
fn ntoc(n: i128, k: i128) -> i128 {
    if k == 0 {return 0;}
    if n == 0 {return (1<<k) - 1;}
    let mut x = k;
    let mut xck = 1; // x choose k
    let mut x1ck = k+1; // x+1 choose k
    while x1ck <= n {
        xck = x1ck;
        x += 1;
        x1ck = (xck * (x + 1)) / (x + 1 - k);
    }
    ntoc(n - xck, k-1) + (1 << x)
}

fn split_combinator(c: Combinator) -> (Combinator, Option<Combinator>) {
    if c.max - c.min <= 100000 {return (c, None);}
    let midpoint = c.min + (c.max - c.min) / 2;
    (Combinator {curr: ntoc(midpoint, c.k as i128), min: c.min, max: midpoint, mask: c.mask, offset: c.offset, k: c.k}, Some(Combinator {curr: c.curr, min: midpoint, max: c.max, mask: c.mask, offset: c.offset, k: c.k}))
}

struct Solution {
    best: AtomicU32,
    minz: i128,
    maxz: i128,
    solutions: Arc<Mutex<HashSet<i128>>>,
}

impl Solution {
    fn update(&self, z: i128) {
        if z < self.minz {return;}
        if z > self.maxz {return;}
        let z2 = z*z;
        let z2c = z2.count_ones();
        let prev = self.best.fetch_max(z2c, Ordering::Relaxed);
        if z2c < prev {return;}
        let mut solutions = self.solutions.lock().unwrap();
        if z2c > prev {
            solutions.clear();
        }
        solutions.insert(z);
    }
    fn new(width: u32) -> Solution {
        Solution {best: AtomicU32::new(0), minz: isqrt((1 << (width - 1)) - 1) + 1, maxz: isqrt(1 << width), solutions: Arc::new(Mutex::new(HashSet::new()))}
    }
    fn get_best(&self) -> u32 {
        self.best.load(Ordering::Relaxed)
    }
}

fn main() {
    let args: Vec<String> = std::env::args().collect();
    let i = match args.len() {
        2 => (1, true),
        3 => {
            if args[1] == "-q" {
                (2, false)
            } else if args[2] == "-q" {
                (1, false)
            } else {
                println!("Usage: max_square [-q] [number of bits]");
                return;
            }
        },
        _ => {
            println!("Usage: max_square [-q] [number of bits]");
            return;
        },
    };
    if let Ok(n) = u32::from_str(&args[i.0]) {
        solve(n, i.1);
    } else {
        println!("Error reading args");
    }
}

Explanation

We'll let \$n\$ be the number of bits (100) in the squares we're searching.

This algorithm has two types of steps that corresponding to fixing the most significant / least significant half of the square.

While there are possible solutions remaining, we chose between the two step types based on which would require fewer operations.

Step type \$a\$:

Range over all integers \$x\$ with \$\left\lceil\frac{n}{2}\right\rceil-1\$ bits of which \$a\$ are zero.

For each \$x\$ there is at most one natural number \$z\$ for which

$$ \underbrace{1x0\ldots0}_{n\text{ bits}} \le z^2 \le \underbrace{1x1\ldots1}_{n\text{ bits}} $$

After updating our solutions with those \$z^2\$, assuming we already checked the previous values of \$a\$, we now know that all unchecked squares have at least \$a+1\$ zeros in the first half.

Step type \$b\$:

Range over integers \$k\$ and \$x\$ where \$2+2k\le b\$ and \$x\$ has \$\left\lceil\frac{n}{2}\right\rceil-3-2k\$ bits of which \$b-2-2k\$ are zero.

For each \$x\$, \$k\$ there are \$4\cdot2^k\$ values of \$z\$ (mod \$2^{\left\lceil\frac{n}{2}\right\rceil}\$) satisfying

$$ z^2\equiv \underbrace{x001[00]\times k}_{\left\lceil\frac{n}{2}\right\rceil\text{ bits}}\mod2^{\left\lceil\frac{n}{2}\right\rceil} $$

After updating our solutions with those \$z^2\$, assuming we already checked the previous values of \$b\$, we now know that all unchecked squares have at least \$b+1\$ zeros in the first half.

Output for n from 3 to 126

3 1 2 2
4 2 3 3
5 3 5 5
6 3 7 7
7 5 11 11
8 4 13 15
9 6 21 21
10 6 27 27
11 8 45 45
12 8 53 53
13 9 75 89
14 9 101 117
15 13 181 181
16 11 217 235
17 13 362 362
18 12 437 491
19 14 627 723
20 15 923 949
21 16 1241 1241
22 16 1619 1773
23 17 2505 2891
24 18 3915 3915
25 19 5221 5747
26 20 6475 7093
27 21 11309 11309
28 22 15595 16203
29 23 19637 19637
30 24 31595 31595
31 25 44491 44491
32 25 61029 64747
33 26 69451 86581
34 27 113447 113447
35 28 185269 185269
36 29 244661 244661
37 30 357081 357081
38 31 453677 453677
39 31 642119 738539
40 34 980853 980853
41 33 1380917 1481453
42 34 1961706 2079669
43 37 2965685 2965685
44 37 3923411 3923411
45 38 5931189 5931189
46 38 8096813 8222581
47 39 11862197 11862965
48 39 16244453 16244453
49 41 22193965 22193965
50 41 32988235 32988235
51 42 43586515 47450235
52 44 57804981 57804981
53 44 94905541 94905541
54 44 115609962 134053275
55 46 164382389 177549685
56 47 268123979 268123979
57 47 320528203 379526189
58 49 528407335 528407335
59 50 657506485 657506485
60 51 888717493 888717493
61 52 1510322221 1510322221
62 52 2100772683 2111529939
63 53 2630119323 2630119323
64 54 4227325771 4227325771
65 55 5853487285 5853487285
66 55 8580481947 8580481947
67 57 10520476455 10520476455
68 57 14878203147 17110519093
69 58 21040953067 24295633269
70 59 31469977781 31469977781
71 62 48589783221 48589783221
72 63 66537313397 66537313397
73 62 97179566261 97179566442
74 63 133074626794 133074626794
75 64 194366549037 194366549037
76 64 259365748555 274835944267
77 66 388730126517 388730126517
78 66 475016392885 540779528373
79 68 759031429451 759031429451
80 68 1064500084555 1081632748507
81 68 1503999171957 1554929409435
82 69 2055848542027 2197923974325
83 70 2576736968885 3109840793419
84 72 3646667521867 4362198168395
85 72 5384675514549 6216643955893
86 73 8445509046091 8791729601637
87 74 11620627528805 11937419436213
88 76 17435796250805 17435796250805
89 76 21537910687315 24879106612889
90 77 34895284158283 34895284158283
91 77 39275570260917 49752141391123
92 78 60940798018325 68134217006923
93 80 93089098383179 93089098383179
94 80 121861586862261 139768858984629
95 81 184920521899109 199030966469849
96 82 222524145495085 222524145495085
97 85 398064946368587 398064946368587
98 84 461438710099403 562672855105355
99 87 796095014224053 796095014224053
100 85 933544295655093 1125891114428899
101 88 1592054863872923 1592054863872923
102 88 2233021208375059 2233021208375059
103 89 3153273599102795 3153273599102795
104 90 4498915668727861 4498915668727861
105 92 6354993066190667 6354993066190667
106 91 7790298568288075 8721180941868443
107 92 12709986132381334 12711635715085643
108 93 18005548694652121 18005548694652121
109 94 25469208479569437 25469208479569437
110 94 30874441840049387 36004547995030069
111 95 47233732013276085 50801164553247509
112 96 66646890995961483 72057405025529767
113 99 101480601046523531 101480601046523531
114 98 143542302959229109 143542302959229109
115 99 202961202093047062 202961202093047062
116 101 269614819660724085 269614819660724085
117 101 403417952942237229 403417952942237229
118 101 529472308325235947 576458535517221045
119 103 808844333869390923 808844333869390923
120 103 970104967079909557 1143803445656875189
121 104 1578696363991412555 1630401033738048053
122 106 1996776422069106507 2156918025961517429
123 106 2688433943050910795 3260133329992761163
124 107 4575515161264941861 4605489388616186187
125 108 6518622390190614475 6518622390190614475
126 110 7977185712837244873 7977185712837244873

________________________________________________________
Executed in  498.22 secs    fish           external
   usr time   65.63 mins   15.29 millis   65.63 mins
   sys time    0.03 mins   19.41 millis    0.03 mins

Answer 2

Going further

49 41 22193965 22193965
50 41 32988235 32988235
51 42 43586515 47450235
52 44 57804981 57804981
53 44 94905541 94905541
54 44 115609962 134053275
55 46 164382389 177549685
56 47 268123979 268123979
57 47 320528203 379526189
58 49 528407335 528407335
59 50 657506485 657506485
60 51 888717493 888717493
61 52 1510322221 1510322221
62 52 2100772683 2111529939
63 53 2630119323 2630119323
64 54 4227325771 4227325771
65 55 5853487285 5853487285
66 55 8580481947 8580481947
67 57 10520476455 10520476455
68 57 14878203147 17110519093
69 58 21040953067 24295633269
70 59 31469977781 31469977781
71 62 48589783221 48589783221
72 63 66537313397 66537313397
73 62 97179566261 97179566442
74 63 133074626794 133074626794
75 64 194366549037 194366549037
76 64 259365748555 274835944267
77 66 388730126517 388730126517
78 66 475016392885 540779528373
79 68 759031429451 759031429451
80 68 1064500084555 1081632748507
81 68 1503999171957 1554929409435
82 69 2055848542027 2197923974325
83 70 2576736968885 3109840793419
84 72 3646667521867 4362198168395
85 72 5384675514549 6216643955893
86 73 8445509046091 8791729601637
87 74 11620627528805 11937419436213
88 76 17435796250805 17435796250805
89 76 21537910687315 24879106612889
90 77 34895284158283 34895284158283
91 77 39275570260917 49752141391123
92 78 60940798018325 68134217006923
93 80 93089098383179 93089098383179
94 80 121861586862261 139768858984629
95 81 184920521899109 199030966469849
96 82 222524145495085 222524145495085
97 85 398064946368587 398064946368587
98 84 461438710099403 562672855105355
99 87 796095014224053 796095014224053
100 85 933544295655093 1125891114428899
101 88 1592054863872923 1592054863872923
102 88 2233021208375059 2233021208375059
103 89 3153273599102795 3153273599102795
104 90 4498915668727861 4498915668727861
105 92 6354993066190667 6354993066190667
106 91 7790298568288075 8721180941868443
107 92 12709986132381334 12711635715085643
108 93 18005548694652121 18005548694652121
109 94 25469208479569437 25469208479569437
110 94 30874441840049387 36004547995030069
111 95 47233732013276085 50801164553247509
112 96 66646890995961483 72057405025529767
113 99 101480601046523531 101480601046523531
114 98 143542302959229109 143542302959229109
115 99 202961202093047062 202961202093047062
116 101 269614819660724085 269614819660724085
117 101 403417952942237229 403417952942237229
118 101 529472308325235947 576458535517221045
119 103 808844333869390923 808844333869390923
120 103 970104967079909557 1143803445656875189
121 104 1578696363991412555 1630401033738048053
122 106 1996776422069106507 2156918025961517429
123 106 2688433943050910795 3260133329992761163
124 107 4575515161264941861 4605489388616186187
125 108 6518622390190614475 6518622390190614475
126 110 7977185712837244873 7977185712837244873

127 overflows basic type X( (maybe) Odd

Try it online!

All files.tar.xz

Answer 3

c++ - 6s whole range (single threaded)

Attempt this online!

This algorithm has two stages where the first stage checks some very unlikely possibilities for optimal squares. Then stage 2 checks the rest of the possibilities, and gains a massive performance boost from the knowledge of stage one.

More specifically we are interested in where the zeroes are in the square. In stage 1 we assume all but at most max_upper_zeros are located in the lower half (bottom 50 bits) of the square. We just brute force trough all possible upper halves of the square, check if they have a square root and if so, remember that value as a potential solution.

Now in stage 2 we enumerate all possibilities, except the ones checked in stage 1. Thus in stage 2 we can assume that the upper half contains at least max_upper_zeros + 1 zeros.

This algorithm has "steps". In step b we enumerate all possible lower segments of b bits of the root. At step 50 we have enumerated all possible roots and can calculate how many ones the corresponding squares have.

As stated, this would not make the calculation any faster than brute force, however we can use the following trick: the b lowest bits of the root determine the b lowest bits of the square. Therefore we can get an upper bound to the number of ones in the square from the lowest b bits, which allows us to discard numbers at early steps of stage 2.

The stage 2 algorithm is run with a "target hamming weight" tg. This is the (minimum) number of ones that we want in the square. We start with tg=100 and decrement until we find solutions.

The stage 2 algorithm in the following code is quite optimized since I invented it first. It has manual tail-call elimination and it only considers odd numbers (odd solutions can be trivially converted to even solutions). The variable p stores the number of bits in the initial segment with a bias given by tg and max_upper_zeros.

The code outputs some stats to stderr and writes the solution to stdout.

#include <vector>
#include <algorithm>
#include <iostream>
#include <bit>
#include <cmath>

const int LIM = 50;
const int TB = LIM*2;
const unsigned long long RANGE_START = 796131459065722;
const unsigned long long RANGE_END = 1125899906842623;

std::vector<unsigned long long> solutions = {};

// Stage 2
void count(int tg, int b, int p, unsigned long long v) {
    unsigned __int128 mul = v * (unsigned __int128)v;
    while (b < LIM) {
        int bit = (mul >> b) & 1;
        unsigned long long vp = v | (1ull << b);
        p+= bit;
        b+= 1;
        // Here we discard solutions that don't have enough ones
        // in the first b bits
        if (p < b) {
            return;
        }
        count(tg, b, p, vp);
    }
    if (__builtin_expect(std::__popcount(mul) >= tg, 0)) {
        // Account for taking into account only odd numbers
        while (RANGE_START > v) {
            v<<= 1;
        }
        if (RANGE_START <= v && v <= RANGE_END) {
            std::cerr << "Found solution " << v << std::endl;
            solutions.push_back(v);
        }
    }
}

unsigned long long isqrt(unsigned __int128 x) {
    unsigned long long apprx = (unsigned long long) sqrt((double)x);
    while (apprx * (unsigned __int128) apprx > x) {
        apprx-= 1;
    }
    return apprx;
}

int main() {
    // First try to find solutions with at most max_upper_zeros zeros in the upper half
    const int max_upper_zeros = 7;

    std::vector<unsigned long long> initial_solutions = {}; 
    int best_hamming_weight = 0;

    unsigned long long upper_neg = 0; 
    while (upper_neg < (1ull << LIM)) {
        unsigned long long upper = ((1ull << LIM)-1)^upper_neg;
        unsigned __int128 m = ((unsigned __int128)upper << LIM);

        unsigned long long v = isqrt(m);
        do {
            unsigned __int128 mul = v * (unsigned __int128)v;
            unsigned long long upp = mul >> LIM;
            if (upp > upper) {
                break;
            }
            if (upp < upper) {
                continue;
            }
            int ones = std::__popcount(mul);
            if (ones >= best_hamming_weight && RANGE_START <= v && v <= RANGE_END) {
                if (ones > best_hamming_weight) {
                    initial_solutions.clear();
                    best_hamming_weight = ones;
                }
                initial_solutions.push_back(v);
            }
        } while (v++);
        if (std::__popcount(upper_neg) == max_upper_zeros) {
            upper_neg = (upper_neg | (upper_neg - 1)) + 1;
        } else {
            upper_neg+= 1;
        }
    }

    std::cerr
    << "First stage found the following solutions with hamming weight "
    << best_hamming_weight << ":" << std::endl;
    for (auto &s: initial_solutions) {
        std::cerr << s << std::endl;
    }

    // Then use this knowledge to calculate solutions more efficiently
    int tg = TB;;
    for (; tg >= best_hamming_weight; --tg) {
        count(tg, 1, 1 + TB - tg - max_upper_zeros, 1);
        std::cerr << "Calculated target " << tg << std::endl;
        if (!solutions.empty()) {
            break;
        }
    }

    if (tg + 1 == best_hamming_weight) {
        solutions.insert(solutions.end(),initial_solutions.begin(),initial_solutions.end());
    }

    unsigned long long min = *std::min_element(solutions.begin(), solutions.end());
    unsigned long long max = *std::max_element(solutions.begin(), solutions.end());
    int ones = std::__popcount(min * (unsigned __int128)min);

    

    std::cout << ones << std::endl << min << std::endl << max << std::endl;

    return 0;
}