Skip to content

fix sqrt with bits > limb_bits #17919

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Nov 8, 2023
Merged

fix sqrt with bits > limb_bits #17919

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
9a67b89
fix sqrt with bits > limb_bits
tisonkun Nov 7, 2023
214c146
add tests
tisonkun Nov 7, 2023
a3ec694
Merge branch 'master' into fix-sqrt-oversize
tisonkun Nov 7, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions lib/std/math/big/int.zig
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -69,7 +69,7 @@ pub fn calcPowLimbsBufferLen(a_bit_count: usize, y: usize) usize {
pub fn calcSqrtLimbsBufferLen(a_bit_count: usize) usize {
const a_limb_count = (a_bit_count - 1) / limb_bits + 1;
const shift = (a_bit_count + 1) / 2;
const u_s_rem_limb_count = 1 + ((shift - 1) / limb_bits + 1);
const u_s_rem_limb_count = 1 + ((shift / limb_bits) + 1);
return a_limb_count + 3 * u_s_rem_limb_count + calcDivLimbsBufferLen(a_limb_count, u_s_rem_limb_count);
}

Expand Down Expand Up @@ -1380,7 +1380,7 @@ pub const Mutable = struct {
var u = b: {
const start = buf_index;
const shift = (a.bitCountAbs() + 1) / 2;
buf_index += 1 + ((shift - 1) / limb_bits + 1);
buf_index += 1 + ((shift / limb_bits) + 1);
Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  1. 1 - The initial value of m is 1, it takes at lest one limb_bit.
  2. shift / limb_bits - For each the shift bits == limb_bits, it takes one more limb_bit.
  3. 1 - The sentinel.

Copy link
Contributor Author

@tisonkun tisonkun Nov 7, 2023 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

and this is how we calculate in the next shiftLeft:

r.normalize(a.limbs.len + (shift / limb_bits) + 1);

notice that we don't shift - 1 here.

Copy link
Contributor Author

@tisonkun tisonkun Nov 8, 2023 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually I'm thinking of that we may not need the so-called sentinel + 1 in both code place. But this is safe while may be waste one usize space.

var m = Mutable.init(limbs_buffer[start..buf_index], 1);
m.shiftLeft(m.toConst(), shift); // u must be >= ⌊√a⌋, and should be as small as possible for efficiency
break :b m;
Expand Down
19 changes: 19 additions & 0 deletions lib/std/math/big/int_test.zig
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3161,6 +3161,7 @@ test "big.int.Managed sqrt(0) = 0" {
try a.setString(10, "0");

try res.sqrt(&a);
try testing.expectEqual(@as(i32, 0), try res.to(i32));
}

test "big.int.Managed sqrt(-1) = error" {
Expand All @@ -3175,3 +3176,21 @@ test "big.int.Managed sqrt(-1) = error" {

try testing.expectError(error.SqrtOfNegativeNumber, res.sqrt(&a));
}

test "big.int.Managed sqrt(n) succeed with res.bitCountAbs() >= usize bits" {
const allocator = testing.allocator;
var a = try Managed.initSet(allocator, 1);
defer a.deinit();

var res = try Managed.initSet(allocator, 1);
defer res.deinit();

// a.bitCountAbs() = 127 so the first attempt has 64 bits >= usize bits
try a.setString(10, "136036462105870278006290938611834481486");
try res.sqrt(&a);

var expected = try Managed.initSet(allocator, 1);
defer expected.deinit();
try expected.setString(10, "11663466984815033033");
try std.testing.expectEqual(std.math.Order.eq, expected.order(res));
}