Programming Sandboxコミュの[お題]階乗n!の桁を求めるプログラム

/*数学的すぎる方法：スターリングの公式使用*/
#include <stdio.h>
#include <math.h>

#define LOG10e 0.43429448190325
#define TWICE_PI 6.2831853071796

int getfigure(int n);

int main()
{
int n;
for (n = 1; n <= 5; n++) {
printf("%d!->%d figures\n", n, getfigure(n));
}

printf("10!->%d figures\n", getfigure(10));
printf("100!->%d figures\n", getfigure(100));
printf("1000!->%d figures\n", getfigure(1000));

return 0;
}

int getfigure(int n)
{
double d = n * (log10((double)n) - LOG10e)
/*n>10ぐらいから、この項は抜いて計算しても変わらない*/
+ 0.5 * log10(TWICE_PI * n);

return (int)(d) + 1;
}

//intしか使わない方法についてはまた考えます＾＾；

おぉ、タイムリーですね。
是非、多倍長でお願いします!

タワーは本当にものすごいことになります。

4318285! ならまだなんとかなりますが、
1000000000000! とか、5638029384213847! とか、大変です。

#include <stdio.h>

int f(int n){
int tmp, digits = 1;
int mul = 1;
while(n > 0){
tmp = n;
while(tmp % 10 == 0){
tmp /= 10;
digits++;
}
mul *= tmp;
n--;
}
while(mul / 10 != 0){
mul /= 10;
digits++;
}
return digits;
}

int main(){
int n = 10;
int i, mul = 1;
for(i = 1; i <= n; i++)
mul *= i;
printf("%d!の桁数は%d\n", n, f(n));
printf("%d!の数値は%d\n", n, mul);

return 0;
}

ひどく簡単なものでオーバーフローとかまったく考えられていません。

これだとせいぜい10の階乗が限界みたいです・・・。
100以上はどうすればいいのか考えてみます(*'-')

なるほど。
10の倍数のところを先に処理して、桁を減らすっていう方法ですね。
その路線でいくなら、2の倍数と5の倍数の時もカウントしておいて、
後で処理するってことができそうです。

しかし、10を削っても桁が爆発的に増えて、mulがオーバーフローしてしまいますね。
いかにmulを小さくするかが課題でしょうか。

実は６で上の考え方の限界を感じて対数を使うコード書いたんですが、dc1394さんが書いたコードが既に適切だったので別のアプローチ考えてみてます(ﾟωﾟ；

ちなみにそのときのコードはやっぱり書き込むことにしますｗ

#include <stdio.h>
#include <math.h>


int f(int n){
    double var = 0.0;
    while(n > 0){
        var += log10((double) n);
        n--;
    }
    return (int) floor(var) + 1;
}

int main(){
    int n = 12;
    unsigned long long i, mul = 1;
    for(i = 1; i <= n; i++)
        mul *= i;
    printf("%d!の桁数は%d\n", n, f(n));
    printf("%lu!の数値は%lu\n", n, mul);

    return 0;
}

printfの中身が変なので一部修正です。
printf("%d!の数値は%llu\n", n, mul);

iをunsigned long longにするならば、nも同じくunsigned long longにしたいですね。
しかし、それだとlog10((double) n)が使えなくなってしまうので、さてどうしたものか。

遅くなりましたが＾＾；
多倍長演算で作ってみました。
参考にしたのはこちらのサイトです
http://www5.airnet.ne.jp/tomy/cpro/longint.htm
上記のサイトを参考にしつつ、多倍長演算クラスを作成してみました。

//=============================================================
// MyLargeInt.h
// MyLargeIntクラスの定義
//=============================================================

#include <iostream>
#include <iomanip>
#include <stdexcept>
#include <vector>
#include <cstring>

typedef unsigned int uint;
static const uint BASE = 10000;
static const int MAX_FIG = 4;

class MyLargeInt {
private:
std::vector<uint> num;

public:
MyLargeInt(int num);
MyLargeInt(const char * num);
MyLargeInt(const MyLargeInt & LInt);
~MyLargeInt(){};

uint getFigure() const;

void set(uint num);
void set(const char * num);
void operator*=(uint num);
void operator*=(const char * num);

friend MyLargeInt operator*(const MyLargeInt & left, const MyLargeInt & right);
friend std::ostream & operator<<(std::ostream & os, const MyLargeInt & LInt);
};

uint charToUInt(const char * num, int startfig, int endfig);
MyLargeInt myFactorial(int n);

// _MyLargeInt_h_

//=============================================================
// MyLargeInt.cpp
// MyLargeIntクラスの実装
//=============================================================

#include "MyLargeInt.h"

MyLargeInt::MyLargeInt(int num)
{
if (num < 0) {
throw new std::range_error("");
}

set(num);
}

MyLargeInt::MyLargeInt(const char * num)
{
set(num);
}

MyLargeInt::MyLargeInt(const MyLargeInt & LInt)
{
this->num = LInt.num;
}

void MyLargeInt::set(uint num)
{
int i = 0;

do {
this->num.resize(this->num.size() + 1);

this->num[i++] = num % BASE;
num /= BASE;
} while (num > 0);
}

void MyLargeInt::set(const char * num)
{
int fig = static_cast<int>(strlen(num));

if (fig <= MAX_FIG) {
this->num.resize(1);
this->num[0] = charToUInt(num, fig - 1, 0);
} else {
int size = (fig % MAX_FIG) ? (fig / MAX_FIG + 1) : fig / MAX_FIG;
this->num.resize(size);

for (int i = 0; i < size; i++) {
int n = fig - MAX_FIG * i;
if (i != size - 1) {
this->num[i] = charToUInt(num, n - 1, n - MAX_FIG);
} else {
this->num[i] = charToUInt(num, n - 1, 0);
}
}
}
}

uint MyLargeInt::getFigure() const
{
uint fig = (this->num.size() - 1) * MAX_FIG;

uint i = num[this->num.size() - 1];

while (i) {
fig++;
i /= 10;
}

return fig;
}

void MyLargeInt::operator*=(uint n)
{
if (n == 0) {
num.resize(1);
num[0] = 0;

return;
}

if (n >= BASE) {
MyLargeInt tmp(n);
*this = *this * tmp;
} else {
uint t = 0;

this->num.resize(this->num.size() + 1);
this->num[this->num.size() - 1] = 0;

for (uint i = 0; i < this->num.size(); i++) {
t += (this->num[i] * n);
this->num[i] = t % BASE;
t /= BASE;
}

if (!this->num[this->num.size() - 1]) {
this->num.resize(this->num.size() - 1);
}
}
}

MyLargeInt operator*(const MyLargeInt & left, const MyLargeInt & right)
{
std::vector<uint>::const_iterator citr;
std::vector<uint>::iterator itr;

MyLargeInt tmp(0);

if (!(left.num[0] * right.num[0])) {
return tmp;
}

if (left.num.size() == 1) {
MyLargeInt tmp(right);
tmp *= left.num[0];

return tmp;
}

if (right.num.size() == 1) {
MyLargeInt tmp(left);
tmp *= right.num[0];

return tmp;
}

tmp.num.resize(left.num.size() + right.num.size());

for (uint j = 0; j < right.num.size(); j++) {
uint u;

if ((u = right.num[j]) == 0) {
continue;
}

uint t = 0;
citr = left.num.begin() + 1;
itr = tmp.num.begin() + j;
uint i = left.num.size();

while (i--) {
*itr = (t += (*itr + u * *citr++)) % BASE;
t /= BASE;
itr++;
}
*itr += t;
}

if (!tmp.num[tmp.num.size() - 1]) {
tmp.num.resize(tmp.num.size() - 1);
}

return tmp;
}

std::ostream & operator<<(std::ostream & os, const MyLargeInt & LInt)
{
for (int i = LInt.num.size() - 1; i >= 0; i--) {
if (i == LInt.num.size() - 1) {
os << LInt.num[i];
continue;
}

os << std::setw(MAX_FIG) << std::setfill('0') << LInt.num[i];
}

return os;
}

uint charToUInt(const char * num, int startfig, int endfig)
{
uint t = 0;
int fi = 1;

for (int i = startfig; i >= endfig; i--) {
int n = num[i] - '0';

if (n < 0 || n > 9) {
throw new std::invalid_argument("");
}

n *= fi;
t += n;
fi *= 10;
}

return t;
}

MyLargeInt myFactorial(int n)
{
MyLargeInt ml(n);

for (int i = n - 1; i > 1; i--) {
ml *= i;
}

return ml;
}

// _MyLargeInt_cpp_

#include "MyLargeInt.h"

int main()
{
for (int n = 1; n < 10; n++) {
std::cout << n << "!->" << myFactorial(n).getFigure() << " figures" << std::endl;
}

std::cout << 10 << "!->" << myFactorial(10).getFigure() << " figures" << std::endl;
std::cout << 100 << "!->" << myFactorial(100).getFigure() << " figures" << std::endl;
std::cout << 1000 << "!->" << myFactorial(1000).getFigure() << " figures" << std::endl;

return EXIT_SUCCESS;
}

今気づいたのですが、operator*=(const char * num)は、宣言はしていますが実装していません＾＾；

バグを発見しました＾＾；
13の
>citr = left.num.begin() + 1;
は、
citr = left.num.begin();
の間違いです＾＾；

またまたバグを発見しました＾＾；
C++では例外をthrowするときにnewしなくていいんですね。
（newするとメモリリークする…？）
C#やJavaの書き方に慣れてしまっていました＾＾；

throw new std::range_error("");
→throw std::range_error("");

throw new std::invalid_argument("");
→throw std::invalid_argument("");

です。失礼しました＾＾；

ループを減らしてC++っぽいソースに直したので投稿し直します＾＾；

//=============================================================
// MyLargeInt.h
// MyLargeIntクラスの定義
//=============================================================

#include <iostream>
#include <iomanip>
#include <algorithm>
#include <functional>
#include <stdexcept>
#include <vector>
#include <cstring>

typedef unsigned int uint;

class MyLargeInt {
private:
std::vector<uint> num;
void removezero();

public:
static const uint BASE = 10000;
static const int MAX_FIG = 4;
explicit MyLargeInt(int num);
explicit MyLargeInt(const char * num);
MyLargeInt(const MyLargeInt & LInt);
~MyLargeInt(){};

uint getFigure();

void set(uint num);
void set(const char * num);
void operator*=(uint num);

friend const MyLargeInt operator*(const MyLargeInt & left, const MyLargeInt & right);
friend std::ostream & operator<<(std::ostream & os, MyLargeInt & LInt);
};

uint charToUInt(const char * num, int startfig, int endfig);
MyLargeInt myFactorial(int n);

class numset :
public std::unary_function<uint, uint> {
private:
uint _num;
public:
numset(uint num) : _num(num){};
uint operator()(uint n) {
uint tmp = _num;
_num /= MyLargeInt::BASE;

return tmp % MyLargeInt::BASE;
}
};

class charnumset:
public std::unary_function<uint, uint> {
private:
const char * _num;
int _size, _fig, cnt;
public:
charnumset(const char * num, int size, int fig) : _num(num), _size(size), _fig(fig), cnt(0){};
uint operator()(uint) {
int n = _fig - MyLargeInt::MAX_FIG * cnt;
uint result = (cnt != _size - 1) ?
charToUInt(_num, n - 1, n - MyLargeInt::MAX_FIG) :
charToUInt(_num, n - 1, 0);
cnt++;

return result;

}
};

class multiply :
public std::unary_function<uint, uint> {
private:
uint _n;
uint result;
public:
multiply(uint n) : _n(n), result(0){};
uint operator()(uint num) {
result += num * _n;
uint tmp = result;
result /= MyLargeInt::BASE;

return tmp % MyLargeInt::BASE;
}
};

struct iszero :
public std::unary_function<uint, bool> {
public:
bool operator()(uint num) const{
return num ? true : false;
}
};

// _MyLargeInt_h_

//=============================================================
// MyLargeInt.cpp
// MyLargeIntクラスの実装
//=============================================================

#include "MyLargeInt.h"

MyLargeInt::MyLargeInt(int num)
{
if (num < 0) {
throw std::range_error("");
}

set(num);
}

MyLargeInt::MyLargeInt(const char * num)
{
set(num);
}

MyLargeInt::MyLargeInt(const MyLargeInt & LInt)
{
this->num = LInt.num;
}

void MyLargeInt::set(uint num)
{
int len = 0;
uint tmp = num;
while (tmp) {
len++;
tmp /= BASE;
}

this->num.resize(len);

std::transform(this->num.begin(), this->num.end(), this->num.begin(), numset(num));
}

void MyLargeInt::set(const char * num)
{
int fig = static_cast<int>(strlen(num));

if (fig <= MAX_FIG) {
this->num.reserve(1);
this->num.push_back(charToUInt(num, fig - 1, 0));
} else {
int size = (fig % MAX_FIG) ? (fig / MAX_FIG + 1) : fig / MAX_FIG;
this->num.resize(size);

std::transform(this->num.begin(), this->num.end(), this->num.begin(), charnumset(num, size, fig));
}
}

void MyLargeInt::removezero()
{
this->num.erase(std::find_if(
this->num.rbegin(), this->num.rend(),
iszero()).base(),
this->num.end());
}

uint MyLargeInt::getFigure()
{
removezero();

uint fig = (this->num.size() - 1) * MAX_FIG;
uint i = num[this->num.size() - 1];

while (i) {
fig++;
i /= 10;
}

return fig;
}

void MyLargeInt::operator*=(uint n)
{
if (n == 0) {
num.resize(1);
num[0] = 0;

return;
}

if (n >= BASE) {
MyLargeInt tmp(n);
*this = *this * tmp;
} else {
uint t = 0;

this->num.reserve(this->num.size() + 1);
this->num.push_back(0);

std::transform(this->num.begin(), this->num.end(), this->num.begin(), multiply(n));
}
}

const MyLargeInt operator*(const MyLargeInt & left, const MyLargeInt & right)
{
std::vector<uint>::const_iterator citr;
std::vector<uint>::iterator itr;

MyLargeInt tmp(0);

if (!(left.num[0] * right.num[0])) {
return tmp;
}

if (left.num.size() == 1) {
MyLargeInt tmp(right);
tmp *= left.num[0];

return tmp;
}

if (right.num.size() == 1) {
MyLargeInt tmp(left);
tmp *= right.num[0];

return tmp;
}

tmp.num.resize(left.num.size() + right.num.size());

for (uint j = 0; j < right.num.size(); j++) {
uint u;

if ((u = right.num[j]) == 0) {
continue;
}

uint t = 0;
citr = left.num.begin() + 1;
itr = tmp.num.begin() + j;
uint i = left.num.size();

while (i--) {
*itr = (t += (*itr + u * *citr++)) % MyLargeInt::BASE;
t /= MyLargeInt::BASE;
itr++;
}
*itr += t;
}

return tmp;
}

std::ostream & operator<<(std::ostream & os, MyLargeInt & LInt)
{
LInt.removezero();

for (int i = LInt.num.size() - 1; i >= 0; i--) {
if (i == LInt.num.size() - 1) {
os << LInt.num[i];
continue;
}

os << std::setw(MyLargeInt::MAX_FIG) << std::setfill('0') << LInt.num[i];
}

return os;
}

uint charToUInt(const char * num, int startfig, int endfig)
{
uint t = 0;
int fi = 1;

for (int i = startfig; i >= endfig; i--) {
int n = num[i] - '0';

if (n < 0 || n > 9) {
throw std::invalid_argument("");
}

n *= fi;
t += n;
fi *= 10;
}

return t;
}

MyLargeInt myFactorial(int n)
{
MyLargeInt ml(n);

for (int i = n - 1; i > 1; i--) {
ml *= i;
}

return ml;
}

// _MyLargeInt_cpp_

Haskellの練習がてら書いてみました。

単純に階乗を計算→文字列に変換→桁数を数えるをしています。
今までやった言語の中で、一番学習コストがかかってそうな気がしますが、
CやC++と違ってメモリ管理を考えなくて楽です。


「ここから」
-- 数値演算用モジュールの読み込み.(数値→文字列変換用)
import Numeric

-- 計算する階乗n!のnを指定
n = 10

-- エントリーポイント
main = print $ f n

-- 階乗の計算
f = length . flip showInt "" . product . flip take ( iterate (+1) 1 )