Lenstra.cc

#include <iostream>
#include <vector>
#include <cmath>
#include <algorithm>
using namespace std;
struct punt{
	int x;
	int y;
	int z;
};

//Comparar punts:
bool operator ==(const punt& A, const punt& B){ 
	return (A.x == B.x) && (A.y == B.y) && (A.z == B.z);
}

//Màxim comú divisor:
int mcd(int a, int b){
	if (b == 0){
		return a;
	}
	if (a == 0){
		return b;
	}
	else{
		return mcd(b, a%b);
	}
}

//Mínim comú múltiple dels enters des de 1 a n:
long long mcmdelsprimers(int n){
	long long P = 1;
	vector<int> Y;
	vector<bool> T(n + 1, true);
	for (int i = 2; i < n + 1; ++i){
		if (T[i]){
			Y.push_back(i);
			for (int j = i*i; j < n + 1; j += i){
				T[j] = false;
			}
		}
	}
	for (int q = 0; q < Y.size(); ++q){
		P *= pow(Y[q], floor(log(n) / log(Y[q])));
	}
	return P;
}

//Resol la identitat de Bézout xa+by=mcd(a, b):
vector<int> bezout(int a, int b){
	vector<int> bez(3);
	int q, r, x, y, g;
	if (b == 0){
		bez[0] = 1, bez[1] = 0, bez[2] = 0;
		return bez;
	}
	q = floor(a / b);
	r = a % b;
	vector<int> k = bezout(b, r);
	x = k[0];
	y = k[1];
	g = k[2];
	bez[0] = y;
	bez[1] = x - q*y;
	bez[2] = mcd(a, b);
	return bez;
}

//Suma a la Llei de Grup:
punt suma(punt P, punt Q, int A, int& N, vector<int>& L){
	int num, den, m, invden; // m=num/den, invden=1/den
	punt R;
	R.x = 0, R.y = 0, R.z = 0;
	if (((P.x == Q.x) && (P.y != Q.y)) || ((P.y == 0) && (Q == P))){
		R.z = 1; //punt a l'infinit
		return R;
	}
	if (!(P == Q)){
		num = (Q.y - P.y) % N;
		den = (Q.x - P.x) % N;
	}
	if (P == Q){
		num = (3 * P.x*P.x + A) % N;
		den = (2 * P.y) % N;
	}
	if (den < 0){
		while (den < 0){
			den += N;
		}
		den = den % N;
	}
	vector<int> j = bezout(den, N);
	if (j[2] > 1){
		if (j[2] != N && j[2] != 1){
			N /= j[2];
			L.push_back(j[2]);
		}
	}
	invden = j[0];
	m = num*invden;
	R.x = (m*m - P.x - Q.x) % N;
	R.y = (m*(P.x - Q.x) - P.y) % N;
	return R;
}

//Muliplica un punt P per un escalar amb la Llei de Grup:
punt multiplica(long long k, punt P, int A, int& N, vector<int>& L){ 
	punt R;
	R.x = 0, R.y = 0, R.z = 0;
	punt Q = P;
	while (k > 0){
		if (k % 2 == 1){
			R = suma(R, Q, A, N, L);
		}
		Q = suma(Q, Q, A, N, L);
		k = floor(k / 2);
	}
	return R;
}

//Retorna un nombre aleatori que pertanyi a Z/NZ i  4A^3+27 mod N pertanyi a (Z/NZ)*
int aleat(int N){
	int A = rand() % N;
	if (mcd(((4 * A*A*A + 27) % N), N) == 1){
		return A;
	}
	else{
		aleat(N);
	}
	
}

//Factors de N:
vector<int> Lenstra(int N){
	int limit = 42;
	vector<int> factors;
	punt P, R;
	P.x = 0, P.y = 1, P.z = 0;
	int A = aleat(N);
	if (round(exp(sqrt(2 * log(N)*log(log(N))) / 2)) < 42){
		limit = round(exp(sqrt(2 * log(N)*log(log(N))) / 2));
	}
	R = multiplica(mcmdelsprimers(limit), P, A, N, factors);
	cout << mcmdelsprimers(limit) << "P = (" << R.x << ", " << R.y << ", " << R.z << ")" << endl;
	factors.push_back(N);
	factors.erase(remove(factors.begin(), factors.end(), 1), factors.end());
	sort(factors.begin(), factors.end());
	return factors;
}


int main(){
	int entradaoriginal;
	int N;
	while (cin >> entradaoriginal){
		vector<int> L;
		vector<int> U;
		N = entradaoriginal;
		while (N >= 2 && N >= 3  && (N % 2 == 0 || N % 3 == 0 || N == round(sqrt(N))*round(sqrt(N)))){
			if (N % 2 == 0){
				N /= 2;
				U.push_back(2);
			}
			else if (N % 3 == 0){
				N /= 3;
				U.push_back(3);
			}
			else if (N == round(sqrt(N))*round(sqrt(N))){
				N = sqrt(N);
				U.push_back(N);
			}
		}
		if (N>3){
            L = Lenstra(N);
		}
		L.insert(L.begin(), U.begin(), U.end());
		sort(L.begin(), L.end());
		for (int j = 0; j < L.size(); ++j){
			if (j == 0){
				cout << entradaoriginal << " = ";
			}
			if (j == L.size() - 1){
				cout << L[j];
			}
			else{
				cout << L[j] << "*";
			}
		}
		cout << endl << endl;
	}
}