Go4Expert

Go4Expert (http://www.go4expert.com/)
-   C (http://www.go4expert.com/forumdisplay.php?f=60)
-   -   Generating prime set which add upto particular number (http://www.go4expert.com/showthread.php?t=18865)

SaswatPadhi 7Aug2009 20:48
Generating prime set which add upto particular number
 

Introduction



On 4th August, shabbir posted this question : (at the $1 Daily contest thread)
List a set of prime numbers where the sum total of the complete set is 100.
So, I thought why not generate try generating *all* such sets of primes, instead of *any* such set ? And, I sat down to code my own generator .. :)

About the program



The program conforms to the latest ANSI C++ standards (compiled and tested with MinGW on Windows XP). On my machine (Celeron D 2.6 GHz, 512 MB RAM), the program generates 2473838 different sets of prime numbers, all of which add to 500; in 38.8 seconds (on average).

The main logic of the program is :
(1) Generate primes till some upper limit : call it PMax.
(2) Get the target number from the user : call it UserInput.
(3) Compute the approximate number of primes less than UserInput : call it LastIndex.
(4) Through a loop, "try to" break UserInput into k different prime numbers by the function RPS.
--> Vary k from 1 to sqrt(UserInput), 'coz at max UserInput can be broken into approximately sqrt(UserInput) increasing prime numbers.
( sqrt(UserInput) is a loose upper limit, but the gap is negligible within 0~500, which we are concerned with. )
--> After the RPS func return the lists, sort them if required and remove duplicate lists.

The backbone of the program -- RPS function works as follows :
(1) Accept for arguments : N, K, Index, PIndex.
--> N is the number we are trying to break.
--> K is the number of parts we want.
--> Index is the list position where we store currently.
--> PIndex is the index of the prime number from which we start breaking.
(2) If K is 1, we just need to check if N is prime or not. So, we do a binary search from PIndex to LastIndex and return check if N is prime.
--> If N is prime, then add N to the list and return 0 value.
[ Let me mention here that, the function always returns X-1 where X is the number of elements added to the list. ]
--> If N is not prime, then return -1. ( 'coz 0 values added. )
(3) If K is not equal to one, then ... (the important portion),
--> Loop through all primes starting from index PIndex till we reach N/K ('coz if we break N into K parts, at least one of them will be less than N/K) : call it P, and it's index Pi.
--> Recursively check if N-P can be broken into K-1 parts : RPS(N-1, K, Index, Pi+1).
--> If N-P could be broken, add the lists along with P to the current list and increase Index.

That's all -- a highly recursive function that breaks N into K parts recursively.

Now, we have some additional thing(s) to take care of :
(*) LastIndex -- The last index of the prime numbers array, till which we conduct the binary search. Calculating this can be omitted, and we can always check till the last index, but that will reduce speed.
To calculate LastIndex approximately, I referred the Wikipedia article on Prime-counting function ( Π(x) ). There, you can find a list that compares Π(x), x/ln x and li(x). From the list, you can see that, for x = 1000, Π(x) - (x/ln x) = 23. So, we can approximately set our LastIndex = 23 + (UserInput/log(UserInput));

The code



So, here is the code ---
Please do give me some credit, if you want to use my code. Just mention my name somewhere :)

Code: C++
/**
      Program written by SASWAT PADHI.

      This program generates a list of sets of prime numbers which add
    up to the number entered by the user.

      There are no repititions in any sets and no two sets contain all
    elements identical.The program detects equivalent sets even if the
    order of elements differ.

**/
#include <cmath>
#include <ctime>
#include <vector>
#include <fstream>
#include <iostream>
#include <algorithm>

using namespace std;

#define VectorList vector< vector<int> >

#define MINI_OUTPUT
#define MEGA_OUTPUT
//#define ASC_ORDERED
//#define FILE_OUTPUT
#define TIME_CALC

#define MaxLimit 600

VectorList MyList;
int Primes[200];
int LastIndex;
vector<int> TList;
#ifdef FILE_OUTPUT
    ofstream fout;
#endif

void PrimeGen();
int RPS(int,int,int,int);
int BS(int,int,int);
void PrintVector(vector<int>);
#ifdef FILE_OUTPUT
    void PrintFVector(vector<int>);
#endif
bool VectorCompare(vector<int>,vector<int>);

int main()
{
    PrimeGen();
    int TOTAL, Found, UserInput, Parts;
    UserInput = 0;
    TList.clear();

    #ifdef FILE_OUTPUT
        fout.open("C:\\MyList.txt");
    #endif

    while ( ( UserInput < 1 ) || ( UserInput > MaxLimit ) )
    {
        cout << "(*) Enter the number (<= " << MaxLimit << "): ";
        cin >> UserInput;
        #ifdef FILE_OUTPUT
            fout << "(*) Enter the number (<= " << MaxLimit << "): " << UserInput << endl;
        #endif
    }

    #ifdef TIME_CALC
        clock_t t1, t2;
        t1 = clock();
    #endif
    LastIndex = 23+(UserInput/log(UserInput));
    TOTAL = 0;
    for (Parts=1; Parts <= sqrt(UserInput); ++Parts)
    {
        MyList.clear();
        Found = RPS(UserInput,Parts,0,0);

        if (++Found)
        {
            #ifdef ASC_ORDERED
                for (int i=0; i<Found; ++i)
                    sort(MyList[i].begin(), MyList[i].end());
            #endif

            Found = MyList.end()-MyList.begin();

            #ifdef MINI_OUTPUT
                cout << "\n(~) Found " << Found << " ways of expressing " << UserInput << " as a sum of " << Parts << " primes.";
                #ifdef FILE_OUTPUT
                    fout << "\n(~) Found " << Found << " ways of expressing " << UserInput << " as a sum of " << Parts << " primes.";
                #endif
                #ifdef MEGA_OUTPUT
                    for (int i=0; i<Found; ++i)
                    {
                        cout << "\n--> ";
                        PrintVector(MyList[i]);
                        #ifdef FILE_OUTPUT
                            fout << "\n--> ";
                            PrintFVector(MyList[i]);
                        #endif
                    }
                    cout << endl;
                #endif
            #endif
            TOTAL += Found;
        }
    }
    cout << "\n\n----< TOTAL " << TOTAL << " UNIQUE SETS OF PRIMES FOUND, WHICH ADD TO " << UserInput << ". >----\n";
    #ifdef FILE_OUTPUT
        fout << "\n\n----< TOTAL " << TOTAL << " UNIQUE SETS OF PRIMES FOUND, WHICH ADD TO " << UserInput << ". >----\n";
    #endif
    #ifdef TIME_CALC
        t2 = clock();
        cout << "\n----< TIME OF EXECUTION = " << ((float)t2 - (float)t1)/1000.0F <<" sec. >----\n";
        #ifdef FILE_OUTPUT
            fout << "\n----< TIME OF EXECUTION = " << ((float)t2 - (float)t1)/1000.0F <<" sec. >----\n";
        #endif
    #endif
    #ifdef FILE_OUTPUT
        fout.close();
    #endif
    return 0;
}

void PrimeGen()
{
    int flag, Cntr=0;
    Primes[0] = 2;
    for (int i=3; i < MaxLimit+1; ++i)
    {
        flag=0;
        for (int j=2; j <= sqrt(i); ++j)
        {
            if ( ! (i%j) )
            {
                flag=-1;
                break;
            }
        }
        if (!flag)
            Primes[++Cntr]=i;
    }
    LastIndex = Cntr;
}

int BS(int first, int last, int key)
{
    while (first <= last)
    {
        int mid = (first + last) / 2;
        if (key > Primes[mid])
            first = mid + 1;
        else if (key < Primes[mid])
            last = mid - 1;
        else
            return 0;
    }
    return -1;
}

void PrintVector( vector<int> MyVector )
{
    int i, j = MyVector.size();
    cout << " { ";
    for (i=0; i<j-1; ++i)
        cout << MyVector[i] << ", ";
    cout << MyVector[i] << " }";
}

#ifdef FILE_OUTPUT
    void PrintFVector( vector<int> MyVector )
    {
        int i, j = MyVector.size();
        fout << " { ";
        for (i=0; i<j-1; ++i)
            fout << MyVector[i] << ", ";
        fout << MyVector[i] << " }";
    }
#endif

bool VectorCompare( vector<int> i, vector<int> j)
{
    return equal(i.begin(), i.end(), j.begin());
}

int RPS(int n, int k, int MIndex, int PStart)
{
    int R, BMI = MIndex+1;
    if (k==1)
    {
        if ( BS(PStart,LastIndex,n) )
            return -1;
        else
        {
            if( MyList.size() <= MIndex )
                MyList.push_back(TList);
            MyList[MIndex].push_back(n);
            return 0;
        }
    }
    else
    {
        for (int i=PStart; Primes[i]<= n/k; ++i)
        {
            R = RPS(n-Primes[i], k-1, MIndex,i+1);
            if ( R != -1)
            {
                for (int j=0; j <= R; ++j)
                {
                    if( MyList.size() <= MIndex )
                        MyList.push_back(TList);
                    MyList[MIndex].push_back(Primes[i]);
                    ++MIndex;
                }
            }
        }
        return MIndex-BMI;
    }
}

Using the code


Using (modifying) the code is simple. Near the top of the program, you can see some #defines :

Code: C++
#define MINI_OUTPUT
#define MEGA_OUTPUT
#define ASC_ORDERED
#define FILE_OUTPUT
#define TIME_CALC
#define MaxLimit 600

Let me explain what their purpose are :
(1) MINI_OUTPUT :: If you un-comment this #define, you can see how many sets are generated of a particular length K. So, verbosity of output increases.
(2) MEGA_OUTPUT :: The most verbose setting. If you un-comment this #define, you can see the whole list of all different sets generated by the program.
(3) ASC_ORDERED :: If you un-comment this #define, the elements inside the sets will be ordered in ascending order.
(4) FILE_OUTPUT :: If you un-comment this #define, the program writes the output to C:\MyFile.txt.
(5) TIME_CALC :: If you un-comment this #define, the program shows the time elapsed to generate the list of sets. (does not count the time you take to input values)

(6) MaxLimit :: Set this to 600 only. Higher values will eat more CPU and memory. Don't try entering 1000 as the target :p

A plot of number of sets vs UserInput




I tried to find a best fit function for the point-series and found the blue one ( e^(x^(13/30)) ) as a satisfactory best fit till UserInput <= 500.

So ... we have come to the end of this article then :)

Thanx for reading this.
I hope this was useful for you.

Good bye and,
Take Care.

SaswatPadhi 8Aug2009 10:48
Re: CODE :: Generating sets of primes which add upto a particular number.
 
A sample run of the code for UserInput = 75.

Program options :

#define MINI_OUTPUT
#define MEGA_OUTPUT
#define ASC_ORDERED
#define TIME_CALC
#define MaxLimit 600

Output :

Code: OUTPUT
(*) Enter the number (<= 600): 75

(~) Found 1 ways of expressing 75 as a sum of 2 primes.
-->  { 2, 73 }
(~) Found 14 ways of expressing 75 as a sum of 3 primes.
-->  { 3, 5, 67 }
-->  { 3, 11, 61 }
-->  { 3, 13, 59 }
-->  { 3, 19, 53 }
-->  { 3, 29, 43 }
-->  { 3, 31, 41 }
-->  { 5, 11, 59 }
-->  { 5, 17, 53 }
-->  { 5, 23, 47 }
-->  { 5, 29, 41 }
-->  { 7, 31, 37 }
-->  { 11, 17, 47 }
-->  { 11, 23, 41 }
-->  { 13, 19, 43 }
(~) Found 17 ways of expressing 75 as a sum of 4 primes.
-->  { 2, 3, 11, 59 }
-->  { 2, 3, 17, 53 }
-->  { 2, 3, 23, 47 }
-->  { 2, 3, 29, 41 }
-->  { 2, 5, 7, 61 }
-->  { 2, 5, 31, 37 }
-->  { 2, 7, 13, 53 }
-->  { 2, 7, 19, 47 }
-->  { 2, 7, 23, 43 }
-->  { 2, 7, 29, 37 }
-->  { 2, 11, 19, 43 }
-->  { 2, 13, 17, 43 }
-->  { 2, 13, 19, 41 }
-->  { 2, 13, 23, 37 }
-->  { 2, 13, 29, 31 }
-->  { 2, 17, 19, 37 }
-->  { 2, 19, 23, 31 }
(~) Found 22 ways of expressing 75 as a sum of 5 primes.
-->  { 3, 5, 7, 13, 47 }
-->  { 3, 5, 7, 17, 43 }
-->  { 3, 5, 7, 19, 41 }
-->  { 3, 5, 7, 23, 37 }
-->  { 3, 5, 7, 29, 31 }
-->  { 3, 5, 11, 13, 43 }
-->  { 3, 5, 11, 19, 37 }
-->  { 3, 5, 13, 17, 37 }
-->  { 3, 5, 13, 23, 31 }
-->  { 3, 5, 17, 19, 31 }
-->  { 3, 7, 11, 13, 41 }
-->  { 3, 7, 11, 17, 37 }
-->  { 3, 7, 11, 23, 31 }
-->  { 3, 7, 13, 23, 29 }
-->  { 3, 7, 17, 19, 29 }
-->  { 3, 11, 13, 17, 31 }
-->  { 3, 11, 13, 19, 29 }
-->  { 3, 13, 17, 19, 23 }
-->  { 5, 7, 11, 23, 29 }
-->  { 5, 7, 13, 19, 31 }
-->  { 5, 11, 13, 17, 29 }
-->  { 5, 11, 17, 19, 23 }
(~) Found 16 ways of expressing 75 as a sum of 6 primes.
-->  { 2, 3, 5, 7, 11, 47 }
-->  { 2, 3, 5, 7, 17, 41 }
-->  { 2, 3, 5, 11, 13, 41 }
-->  { 2, 3, 5, 11, 17, 37 }
-->  { 2, 3, 5, 11, 23, 31 }
-->  { 2, 3, 5, 13, 23, 29 }
-->  { 2, 3, 5, 17, 19, 29 }
-->  { 2, 3, 7, 11, 23, 29 }
-->  { 2, 3, 7, 13, 19, 31 }
-->  { 2, 3, 11, 13, 17, 29 }
-->  { 2, 3, 11, 17, 19, 23 }
-->  { 2, 5, 7, 11, 13, 37 }
-->  { 2, 5, 7, 11, 19, 31 }
-->  { 2, 5, 7, 13, 17, 31 }
-->  { 2, 5, 7, 13, 19, 29 }
-->  { 2, 7, 11, 13, 19, 23 }
(~) Found 1 ways of expressing 75 as a sum of 7 primes.
-->  { 3, 5, 7, 11, 13, 17, 19 }

----< TOTAL 71 UNIQUE SETS OF PRIMES FOUND, WHICH ADD TO 75. >----

----< TIME OF EXECUTION = 0.172 sec. >----

shabbir 2Sep2009 18:56
Re: Generating prime set which add upto particular number
 
Nomination for Article of the month - Aug 2009 Started.

shabbir 19Sep2009 11:41
Re: Generating prime set which add upto particular number
 
Vote for this article for Article of the month - August 2009

SaswatPadhi 20Sep2009 18:20
Re: Generating prime set which add upto particular number
 
I had changed (improved the performance of) the above code but had forgotten to post it here.
When I saw this threads' name in the recent threads list, I remembered to post it. :)

Here it is, the new code :
Code: CPP
#include <cmath>
#include <ctime>
#include <vector>
#include <fstream>
#include <iostream>
#include <algorithm>

using namespace std;

#define VectorList vector< vector<int> >

#define MINI_OUTPUT
#define MEGA_OUTPUT
//#define ASC_ORDERED
//#define FILE_OUTPUT
#define TIME_CALC

#define MaxLimit 600

VectorList MyList;
int Primes[200];
int LastIndex;
vector<int> TList;
#ifdef FILE_OUTPUT
    ofstream fout;
#endif

void PrimeGen();
int RPS(int,int,int,int);
int BS(int,int,int);
void PrintVector(vector<int>);
#ifdef FILE_OUTPUT
    void PrintFVector(vector<int>);
#endif
bool VectorCompare(vector<int>,vector<int>);

int main()
{
    PrimeGen();
    int TOTAL, Found, UserInput, Parts;
    UserInput = 0;
    TList.clear();

    #ifdef FILE_OUTPUT
        fout.open("C:\\MyList.txt");
    #endif

    while ( ( UserInput < 1 ) || ( UserInput > MaxLimit ) )
    {
        cout << "(*) Enter the number (<= " << MaxLimit << "): ";
        cin >> UserInput;
        #ifdef FILE_OUTPUT
            fout << "(*) Enter the number (<= " << MaxLimit << "): " << UserInput << endl;
        #endif
    }

    #ifdef TIME_CALC
        clock_t t1, t2;
        t1 = clock();
    #endif
    LastIndex = 23+(UserInput/log(UserInput));
    TOTAL = 0;
    for (Parts=1; Parts <= sqrt(UserInput); ++Parts)
    {
        MyList.clear();
        Found = RPS(UserInput,Parts,0,0);

        if (++Found)
        {
            #ifdef ASC_ORDERED
                for (int i=0; i<Found; ++i)
                    sort(MyList[i].begin(), MyList[i].end());
            #endif

            Found = MyList.end()-MyList.begin();

            #ifdef MINI_OUTPUT
                cout << "\n(~) Found " << Found << " ways of expressing " << UserInput << " as a sum of " << Parts << " primes.";
                #ifdef FILE_OUTPUT
                    fout << "\n(~) Found " << Found << " ways of expressing " << UserInput << " as a sum of " << Parts << " primes.";
                #endif
                #ifdef MEGA_OUTPUT
                    for (int i=0; i<Found; ++i)
                    {
                        cout << "\n--> ";
                        PrintVector(MyList[i]);
                        #ifdef FILE_OUTPUT
                            fout << "\n--> ";
                            PrintFVector(MyList[i]);
                        #endif
                    }
                    cout << endl;
                #endif
            #endif
            TOTAL += Found;
        }
    }
    cout << "\n\n----< TOTAL " << TOTAL << " UNIQUE SETS OF PRIMES FOUND, WHICH ADD TO " << UserInput << ". >----\n";
    #ifdef FILE_OUTPUT
        fout << "\n\n----< TOTAL " << TOTAL << " UNIQUE SETS OF PRIMES FOUND, WHICH ADD TO " << UserInput << ". >----\n";
    #endif
    #ifdef TIME_CALC
        t2 = clock();
        cout << "\n----< TIME OF EXECUTION = " << ((float)t2 - (float)t1)/1000.0F <<" sec. >----\n";
        #ifdef FILE_OUTPUT
            fout << "\n----< TIME OF EXECUTION = " << ((float)t2 - (float)t1)/1000.0F <<" sec. >----\n";
        #endif
    #endif
    #ifdef FILE_OUTPUT
        fout.close();
    #endif
    return 0;
}

void PrimeGen()
{
    int flag, Cntr=0;
    Primes[0] = 2;
    for (int i=3; i < MaxLimit+1; ++i)
    {
        flag=0;
        for (int j=2; j <= sqrt(i); ++j)
        {
            if ( ! (i%j) )
            {
                flag=-1;
                break;
            }
        }
        if (!flag)
            Primes[++Cntr]=i;
    }
    LastIndex = Cntr;
}

int BS(int first, int last, int key)
{
    while (first <= last)
    {
        int mid = (first + last) / 2;
        if (key > Primes[mid])
            first = mid + 1;
        else if (key < Primes[mid])
            last = mid - 1;
        else
            return 0;
    }
    return -1;
}

void PrintVector( vector<int> MyVector )
{
    int i, j = MyVector.size();
    cout << " { ";
    for (i=0; i<j-1; ++i)
        cout << MyVector[i] << ", ";
    cout << MyVector[i] << " }";
}

#ifdef FILE_OUTPUT
    void PrintFVector( vector<int> MyVector )
    {
        int i, j = MyVector.size();
        fout << " { ";
        for (i=0; i<j-1; ++i)
            fout << MyVector[i] << ", ";
        fout << MyVector[i] << " }";
    }
#endif

bool VectorCompare( vector<int> i, vector<int> j)
{
    return equal(i.begin(), i.end(), j.begin());
}

int RPS(int n, int k, int MIndex, int PStart)
{
    int R, BMI = MIndex+1;
    if (k==1)
    {
        if ( BS(PStart,LastIndex,n) )
            return -1;
        else
        {
            if( MyList.size() <= MIndex )
                MyList.push_back(TList);
            MyList[MIndex].push_back(n);
            return 0;
        }
    }
    else
    {
        for (int i=PStart; Primes[i]<= n/k; ++i)
        {
            R = RPS(n-Primes[i], k-1, MIndex,i+1);
            if ( R != -1)
            {
                for (int j=0; j <= R; ++j)
                {
                    if( MyList.size() <= MIndex )
                        MyList.push_back(TList);
                    MyList[MIndex].push_back(Primes[i]);
                    ++MIndex;
                }
            }
        }
        return MIndex-BMI;
    }
}

If you notice closely, there isn't much change in the code.
Actually, what I had done earlier was, I was checking all vectors (which were the prime-sets) for their possible duplicates.
But I just noticed (few days after posting the article here) that the way I generate the primes-sets, all of them must be unique. And they actually were !
So, I just removed the duplicate check :
Code: C++
Found = unique(MyList.begin(), MyList.end(), VectorCompare)-MyList.begin();
to :
Code: C++
Found = MyList.end()-MyList.begin();

And, you get a 20% performance improvement. :D

[[ I request shabbir to replace the code in the article with the one above. ]]

shabbir 2Oct2009 10:22
Re: Generating prime set which add upto particular number
 
Winner of Article of the month - August 2009 Jointly


All times are GMT +5.5. The time now is 21:23.