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.

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

Implemented Quadratic residue function #736

Merged
merged 26 commits into from Mar 23, 2016

Conversation

Projects
None yet
5 participants
@CodeMaxx
Copy link
Contributor

CodeMaxx commented Jan 2, 2016

  • Quadratic Residue function
  • is_quad_residue function
  • is_nthroot_mod_prime_power function
  • is_nthroot_mod1() function
  • is_nth_power_residue() function

@isuruf @Sumith1896 Could you please review this?

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Jan 2, 2016

There are compile errors. You can test locally by running make and then test using ctest -V.
Also, you can check here https://travis-ci.org/CodeMaxx/symengine

@CodeMaxx CodeMaxx force-pushed the CodeMaxx:quadratic_residue branch from c9c3697 to ca73632 Jan 2, 2016

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Jan 2, 2016

@isuruf Compiler errors fixed and it passes the one test I set. I'll add more tests. By that time you can review the rest of the code.

@CodeMaxx CodeMaxx force-pushed the CodeMaxx:quadratic_residue branch from ce36e32 to 7991983 Jan 2, 2016

@CodeMaxx CodeMaxx force-pushed the CodeMaxx:quadratic_residue branch from 7991983 to 20a9dc9 Jan 2, 2016

CodeMaxx added some commits Jan 2, 2016

@CodeMaxx CodeMaxx force-pushed the CodeMaxx:quadratic_residue branch from bed19fe to 13ca5ac Jan 2, 2016

@certik

This comment has been minimized.

Copy link
Contributor

certik commented Jan 3, 2016

You got a test failure on Appveyor:

test_quadratic_residue(): ntheory
-------------------------------------------------------------------------------
C:\projects\symengine\symengine\tests\ntheory\test_ntheory.cpp:640
...............................................................................

C:\projects\symengine\symengine\tests\ntheory\test_ntheory.cpp:666: FAILED:
  REQUIRE( quadratic_residue(*a100) == i100 )
with expansion:
  { 0, 1, 4, 9, 15, 16, 20, 24, 25, 28, 29, 35, 36, 40, 44, 49, 56, 60, 61, 63,
  64, 68, 69, 75, 76, 81, 83, 89, 96, 99 }
  ==
  { 0, 1, 4, 9, 16, 21, 24, 25, 29, 36, 41, 44, 49, 56, 61, 64, 69, 76, 81, 84,
  89, 96 }

===============================================================================
test cases:  25 |  24 passed | 1 failed
assertions: 230 | 229 passed | 1 failed


      Start 14: test_diophantine

Can you look into that? Is that a bug in your code?

@@ -169,6 +169,9 @@ bool powermod(const Ptr<RCP<const Integer>> &powm, const RCP<const Integer> &a,
//! All solutions to x**s == a**r mod m where b = r / s. Return false if none exists.
void powermod_list(std::vector<RCP<const Integer>> &pows, const RCP<const Integer> &a,
const RCP<const Number> &b, const RCP<const Integer> &m);

//! Finds all Quadratic Residue of a Positive Integer

This comment has been minimized.

@certik

certik Jan 3, 2016

Contributor

Residue -> Residues?

This comment has been minimized.

@CodeMaxx

CodeMaxx Jan 3, 2016

Contributor

Corrected.

std::vector<int> residue;
for( int i = 0; i <= a.as_int()/2; i++)
{
residue.push_back( (int)pow(i,2) % a.as_int());

This comment has been minimized.

@certik

certik Jan 3, 2016

Contributor

Could the Appveyor failure be caused by using a 32bit int that overflows here? If so, you need to use an mpz_class.

This comment has been minimized.

@certik

certik Jan 3, 2016

Contributor

If this is the case, then we just got lucky. I created #737 to test this robustly.

This comment has been minimized.

@CodeMaxx

CodeMaxx Jan 3, 2016

Contributor

I guess that could be the problem since tests pass on 64 bit version. I'll look into it.

This comment has been minimized.

@CodeMaxx

CodeMaxx Jan 3, 2016

Contributor

@certik Actually the value here goes no more than 2500 so it shouldn't be an overflow.

throw std::runtime_error("quadratic_residue: Input must be > 0");
}

std::vector<int> residue;

This comment has been minimized.

@certik

certik Jan 3, 2016

Contributor

You might need to use mpz_class instead of an int here, see my other comment.

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Jan 3, 2016

To see if it works I changed all of the ints to mpz_class

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Jan 3, 2016

@certik Its working with mpz_class. But I still don't beleive it was an int overflow. The largest number in that test case was just 2500.

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Jan 3, 2016

You were using (int) pow(i, 2) instead of just i*i. Casting double to int leads to incorrect results.

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Jan 3, 2016

@isuruf But pow(i,2) always gave an integral answer that too <=2500 .

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Jan 5, 2016

@CodeMaxx, the compiler is free to inline pow(i, 2) call to i*i when optimizations are on. That is why you get integral answers in Release mode.

pow(i, 2) is not inlined in Debug mode and in 32-bit MinGW, std::pow implementation returns an inexact answer. A compiler is not required to give the exact answer in floating point functions nor can a compiler represent numbers exactly.

You are right about overflow though, for 100, it is not possible to overflow.

@Sumith1896

This comment has been minimized.

Copy link
Member

Sumith1896 commented Jan 5, 2016

@CodeMaxx Let us know if you need more help.

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Jan 5, 2016

@Sumith1896 I'm having some trouble using mpz_class instead of int. Currently I'm pushing with some parts commented out which I can't get to work.

@certik

This comment has been minimized.

Copy link

certik commented on symengine/ntheory.cpp in 26e066f Jan 5, 2016

You can't assign to a const (i.e. an input) argument. What are you trying to do?

This comment has been minimized.

Copy link
Owner

CodeMaxx replied Jan 5, 2016

According to the specifications of the program, I need to check if a % p is a quadratic residue. To keep the code clean I was trying to omit the use of a new variable to store a % p. I guess i'll have to define a new mpz_class variable then.

return false;

RCP<const Integer> x;
const RCP<const Integer> a1 = integer(a_final);

This comment has been minimized.

@Sumith1896

Sumith1896 Jan 5, 2016

Member

Take care of tabs and spaces.

This comment has been minimized.

@CodeMaxx

CodeMaxx Jan 5, 2016

Contributor

Everything is spaces. Is there an inconsistency anywhere?
Got it.

@Sumith1896

This comment has been minimized.

Copy link
Member

Sumith1896 commented Jan 5, 2016

Can you add tests for the second function?

@CodeMaxx CodeMaxx force-pushed the CodeMaxx:quadratic_residue branch from 3a6b147 to e76fc7e Jan 5, 2016

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Jan 5, 2016

I'll add tests by tomorrow.
---------------------------append
or maybe in about an hour or something.

// Returns whether Solution for x**n == a mod p**k exists or not
bool is_nthroot_mod_prime_power(const integer_class &a, const integer_class &n, const integer_class &p, const unsigned k)
{
integer_class pk;

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

Fix formatting

*/

integer_class p2 = p.as_mpz();
if (p2 < 1)

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

p can be negative right?

This comment has been minimized.

@CodeMaxx

CodeMaxx Mar 10, 2016

Contributor

No p is a positive number. Otherwise its prime factorisation and stuff will not be defined.

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

Why not? A negative number can be factored. What else is undefined?

This comment has been minimized.

@CodeMaxx

CodeMaxx Mar 10, 2016

Contributor

Its actually how the quadratic residues are defined....See http://mathworld.wolfram.com/QuadraticResidue.html.

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

In all other places including Mathematica, SageMath, it is defined

This comment has been minimized.

@CodeMaxx

CodeMaxx Mar 10, 2016

Contributor

So do I just need to exclude p = 0 then?

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

Yes. And let p = -p for p < 0

integer_class a_final = a.as_mpz();
if (a.as_mpz() >= p2 || a.as_mpz() < 0)
mp_fdiv_r(a_final, a.as_mpz(), p2);
if (a_final < 2 || p2 < 3)

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

If a_final < 2 then p2 < 3 case doesn't need to be handled right?

This comment has been minimized.

@CodeMaxx

CodeMaxx Mar 10, 2016

Contributor

Its not necessary but since no computation was required for them I decided to add this statement and return true straight away.

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

I mean, a_final < 2 is fine. What I'm saying is that p2 < 3 is redundant

This comment has been minimized.

@CodeMaxx

CodeMaxx Mar 10, 2016

Contributor

Oh right. Fixed it.

return is_nthroot_mod1(a, n, p, k);
}
} else {
integer_class _a;

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

See this from wikipedia

In general, to determine if a quadratic congruence with composite modulus is solvable use the following theorem:[37]

Let n > 1, and gcd(a,n) =1. Then x2 ≡ a (mod n) is solvable if and only if:

The Legendre symbol \left(\frac{a}{p}\right)=1 for all odd prime divisors p of n.
a ≡ 1 (mod 4) if 4|n, but 8 does not divide n; a ≡ 1 (mod 8) if 8|n.

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

Ignore that

@@ -159,6 +159,8 @@ int legendre(const Integer &a, const Integer &n);
int jacobi(const Integer &a, const Integer &n);
//! Kronecker Function
int kronecker(const Integer &a, const Integer &n);
// Returns whether Solution for x**n == a mod p**k exists or not
bool is_nthroot_mod_prime_power(const integer_class &a, const integer_class &n, const integer_class &p, const unsigned k);

This comment has been minimized.

@isuruf

isuruf Mar 10, 2016

Member

This is not needed in the header

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Mar 10, 2016

@CodeMaxx, since the code is generic for nthroot, can you add a new method for nthroot as well?

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Mar 10, 2016

@isuruf Do you mean two functions? ( nth_residue and is_nth_residue)

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Mar 10, 2016

Or just is_nth_residue.

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Mar 10, 2016

@isuruf I'm not sure if conditions like 1 and 2 hold for nth_residue or not. Do they?

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Mar 10, 2016

legendre is only for quadratic residues.

@isuruf So code is not generic for nthroot.

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Mar 14, 2016

@isuruf Can review this.

*/
{
if (mod->as_mpz() <= 0) {
return false;

This comment has been minimized.

@isuruf

isuruf Mar 22, 2016

Member

mod can be negative


bool is_nth_power_residue(const RCP<const Integer> &a, const RCP<const Integer> &n, const RCP<const Integer> &mod)
/*
Returns true if ``a`` (mod ``p``) is in the set of squares mod ``p``,

This comment has been minimized.

@isuruf

isuruf Mar 22, 2016

Member

p is mod in the function signature

return mp_legendre(a_final, p2) == 1;
}

bool is_nth_power_residue(const RCP<const Integer> &a, const RCP<const Integer> &n, const RCP<const Integer> &mod)

This comment has been minimized.

@isuruf

isuruf Mar 22, 2016

Member

This signature differs significantly from that of is_quad_residue

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Mar 22, 2016

@CodeMaxx, I'm sorry I didn't see this message for review. I left a few comments. When you are done ping me.

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Mar 22, 2016

@isuruf Ping


bool is_nth_residue(const RCP<const Integer> &a, const RCP<const Integer> &n, const RCP<const Integer> &mod)
/*
Returns true if ``a`` (mod ``mod``) is in the set of squares mod ``mod``,

This comment has been minimized.

@isuruf

isuruf Mar 22, 2016

Member

set of nth powers

return true;
}

if(_mod < 0)

This comment has been minimized.

@isuruf

isuruf Mar 22, 2016

Member

Space after if

return mp_legendre(a_final, p2) == 1;
}

bool is_nth_residue(const RCP<const Integer> &a, const RCP<const Integer> &n, const RCP<const Integer> &mod)

This comment has been minimized.

@isuruf

isuruf Mar 22, 2016

Member

You don't have to use RCP here

This comment has been minimized.

@CodeMaxx

CodeMaxx Mar 22, 2016

Contributor

@isuruf Oh right. Fixed.

@isuruf

This comment has been minimized.

Copy link
Member

isuruf commented Mar 22, 2016

+1 to merge when tests pass

@CodeMaxx

This comment has been minimized.

Copy link
Contributor

CodeMaxx commented Mar 22, 2016

@isuruf ping

@isuruf isuruf merged commit 3912c27 into symengine:master Mar 23, 2016

3 checks passed

codecov/project 74.15% (+0.14%) compared to b281b8f at 74.01%
Details
continuous-integration/appveyor/pr AppVeyor build succeeded
Details
continuous-integration/travis-ci/pr The Travis CI build passed
Details
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment