Really odd maths bug

Hello,

I've just come across a strange bug in my program. Firstly, here's the code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
 /*

Outputs the prime factorisation of a number.

*/

#include <iostream>
#include <cmath>

using namespace std;

int checkprime(int k);
int smallestprimedivisor(int d);

int main(){

int n;

for( ; ; ){

cout << "Which number would you like to decompose into primes?:";
cin >> n;

while(!checkprime(n)){ 	cout << smallestprimedivisor(n) << 'x';
						n = n/smallestprimedivisor(n);	}

if(!(n==1)){cout << n;}

cout << '\n';
}

return 0;

}

int checkprime(int k){

//Checks if a number is prime.  Outputs 1 if it is, and 0 if it isn't.
				
				int i;
				if( k==2){ return 1;}
				if(!(k%2)){return 0;}

				for(i=3; i<=pow(k, 0.5); i+=2){if(!(k%i)){	return 0;}}
				return 1;

}

int smallestprimedivisor(int d){

//Finds the smallest prime that divides a number

				int h;
				if(checkprime(d)){ return d;}
				if(!(d%2)){ return 2; }
				
				for(int i=3; i<=d; i+=2){	if(checkprime(i)){ 	h = d%i;
																if(!h) {return i;}}}

}

				
					


This outputs the prime factorisation of a number. I'm on school computers right now, and I can't compile or copy and paste from the terminal, but if I enter 987 I get:

3 x 7 x 47

Which is right from here:

http://www.mathwarehouse.com/answered-questions/prime-factorizations/what-is-the-prime-factorization-of-987-solved.php

I tried it with some other numbers too, and they work perfectly.

However, this is where the problems start.

If you enter 169 (which is 13^2), it should output 13 x 13.

However, it spits out 169 again.

The same problem occurs with 841 (which is 29^2). It ouputs 841 instead of 29 x 29.

Can anyone see why this is? I have checked my checkprime function by getting it to output every prime number less than 500, and cross checking it with a website. The problem isn't there. I think the problem is with my while loop, but I can't see why.

I'll have a go myself when i'm home, but i'm working quite late tonight.

Any help would be appreciated.
> i<=pow(k, 0.5)

Floating point computations and comparisons are inexact.

1
2
3
4
5
6
7
8
9
10
11
/*int*/ bool checkprime( int k ) {

    if( k < 2 ) return false ;
    if( k == 2 ) return true ;
    if( k%2 == 0 ) return false ;

    const int ul = std::ceil( std::sqrt(k) ) + 1.1 ; // 
    for( int i = 3 ; i < ul ; i += 2 ) if( k%i == 0 ) return false ;

    return true ;
}
Wow, that looks crazy.

I'm afraid i'm a beginner, so I haven't come across a problem like this before.

Can you explain this bit:

const int ul = std::ceil( std::sqrt(k) ) + 1.1 ; //

My understanding is that ul is an integer whose "ceiling" is the square root of k. I don't really see why you're adding 1.1 though.

By "ceiling" I mean something like ceiling(1.1)=2 or ceiling (5.6)=6. It's like an opposite of an integer part.

For example, if k was 169, then the ceiling(sqrt(169))= 13. Then you add 1.1 to get 14.1. The integer loop then kicks in checking numbers from 3 to 14.1.

How is this different to what I had written? I checked my checkprime function when I wrote it, and was able to find all the primes up to 500. I cross checked these with an internet base, so what went wrong here?

EDIT: Thanks for reminding me of the "sqrt" term, and also that the function should have been a boolean. I should have noticed those!

EDIT: I've messed around with it a bit and why have you put in std:ceil instead of just ceil? It still works if you get rid of the std. Is it standard practice to use std?

Also, why does ul have to be a const int? It works without a const in there.

I should probably post the amended program:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
/*

Outputs the prime factorisation of a number.

*/

#include <iostream>
#include <cmath>

using namespace std;

bool prime(long long int k);
long long int smallestprimedivisor(long long int d);

int main(){

long long int n;

for( ; ; ){

cout << "Which number would you like to decompose into primes?:";
cin >> n;

while(!prime(n)){ 	cout << smallestprimedivisor(n) << 'x';
						n = n/smallestprimedivisor(n);	}

if(!(n==1)){cout << n;}

cout << '\n';

}

return 0;

}

bool prime( long long int k ) {

    if( k < 2 ) return 0;
    if( k == 2 ) return 1;
    if( !(k%2) ) return 0;

    int ub = ceil( sqrt(k) ) ; 
    for( int i = 3 ; i <= ub ; i += 2 ) if( !(k%i) ) return 0 ;

    return 1 ;
}

long long int smallestprimedivisor(long long int d){

//Finds the smallest prime that divides a number

				int h;
				if(prime(d)){ return d;}
				if(!(d%2)){ return 2; }
				
				for(int i=3; i<=d; i+=2){	if(prime(i)){ 	h = d%i;
																if(!h) {return i;}}}

}

				
					

Last edited on
On my computer, 169 does print out 13x13.

I did notice if you type in 0, or a letter it goes in to a endless loop.
> Can you explain this bit:
> const int ul = std::ceil( std::sqrt(k) ) + 1.1 ; //

On a particular implementation:
std::sqrt(169) may yield 12.9999999999999999
and std::ceil( std::sqrt(169) ) may yield 13.0000000000000001
and 13.0000000000000001 + 1.0 may yield 13.9999999999999999


> I don't really see why you're adding 1.1 though.

I favoured the idiomatic loop with i < ul over a less idiomatic one with i <= maxv


> Also, why does ul have to be a const int?

It does not have to be a const int. The const came out of force of habit.

15. Use const proactively.

const is your friend: Immutable values are easier to understand, track, and reason about, so prefer constants over variables wherever it is sensible and make const your default choice when you define a value: It's safe, it's checked at compile time, and it's integrated with C++'s type system. ...

Constants simplify code because you only have to look at where the constant is defined to know its value everywhere. Consider this code:

1
2
3
4
5
void Fun( vector<int>& v)
{  //...
    const size_t len = v.size();
    //... 30 more lines ... 
}


When seeing len's definition above, you gain instant confidence about len's semantics throughout its scope (assuming the code doesn't cast away const, which it should not do): It's a snapshot of v's length at a specific point. Just by looking up one line of code, you know len's semantics over its whole scope. Without the const, len might be later modified, either directly or through an alias. Best of all, the compiler will help you ensure that this truth remains true. ...

- Alexandrescu and Sutter in 'C++ Coding Standards: 101 Rules, Guidelines, and Best Practices'



> It still works if you get rid of the std. Is it standard practice to use std?

Yes; in using a name in the standard C++ library.
Cubbi:
C++ source code is much easier to read with std::. When someone reads std::copy in your code, they know you're referring to the C++ standard algorithm and not boost::copy or some company::team::project::copy() function you've pulled in from some obscure header file. Granted, it doesn't really promise that: the only certain way is to write ::std::copy, but almost nobody goes that far. http://www.cplusplus.com/forum/general/72248/#msg385442


Many programmers avoid dumping using namespace std; at global namespace scope.
http://www.cplusplus.com/forum/general/72248/#msg385659


You may want to consider using a prime number sieve:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
#include <iostream>
#include <cmath>
#include <vector>

std::vector<std::size_t> prime_numbers_till( std::size_t ubound )
{
    // generate a prime number sieve up to ubound
    std::vector<bool> sieve( ubound, true ) ;
    for( std::size_t i = 2 ; i < ubound ; ++i )
        if( sieve[i] )
            for( auto j = i*i ; j < ubound ; j += i ) sieve[j] = false ;

    std::vector<std::size_t> primes( 1, 2 ) ;
    for( std::size_t i = 3 ; i < ubound ; i += 2 ) if( sieve[i] ) primes.push_back(i) ;
    return primes ;
}

std::vector<long long> prime_factors_of( long long n )
{
    if( n < 0 ) return prime_factors_of(-n) ;

    std::vector<long long> factors ;
    for( int p : prime_numbers_till( std::sqrt(n) + 2 ) ) while( n%p == 0 )
    {
         factors.push_back(p) ;
         n /= p ;
         if( n == 1 ) break ; // no more factors
    }

    if( n > 1 ) factors.push_back(n) ; // n is a prime number

    return factors ;
}

int main()
{
    const long long n = 2LL * 19*19*19*19*19 * 41*41 * 1151 * 2239 ;
    for( int v : prime_factors_of(n) ) std::cout << v << ' ' ;
    std::cout << '\n' ;
}

http://coliru.stacked-crooked.com/a/4a98c6ae846bf358
It would be more efficient do prime() in terms of smallestprimedivisor. This works because the first factor that your current prime() code finds is itself prime. It has to be because otherwise you would have found one of its factors earlier in the loop. Something like this (I'm away from my compiler so it might have some errors):

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
long long smallestprimedivisor(long long int k ) {

    if( k < 2 ) return 0;
    if( k == 2 ) return 2;
    if( !(k%2) ) return 2;

    long long ub = ceil( sqrt(k) ) ; 
    for( long long i = 3 ; i <= ub ; i += 2 ) if( !(k%i) ) return i;

    return k;
}

bool prime(long long int d){
    return smallestprimedivisor(d) == d;
}
Topic archived. No new replies allowed.