「オイラープロジェクト71～80」の編集履歴（バックアップ）一覧に戻る
オイラープロジェクト71～80 - (2012/11/07 (水) 08:04:40) のソース
*問い71
nとdを正の整数として, 分数 n/d を考えよう. n<dかつHCF(n,d)=1のとき, 真既約分数と呼ぶ.
d ≤ 8について既約分数を大きさ順に並べると, 以下を得る:
1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
3/7のすぐ左の分数は2/5である.
d ≤ 1,000,000について真既約分数を大きさ順に並べたとき, 3/7のすぐ左の分数の分子を求めよ.
ある3/7に一番近い分数aの分母が決まった時aは3/7に一番近くなる分子が決まります。
これが決まった後aの分母を１増加した時,aの値は必ず減少します。
減少した分だけ3/7から遠ざかるので分子をインクリメントして補てんします。
この時3/7を超えない最大値の分数を求めます。
これを100万まで繰り返せば答えが出ます。
#include<stdio.h>
#include<queue>
#include<set>
int GCD(int a, int b)
{
while( 1 )
{
a = a % b;
if( a == 0 )return b;
b = b % a;
if( b == 0 )return a;
}
}
int main(){
double limit=3.0/7.0;
double d=1,u,sa=limit-3.0/8.0;
int ansA,ansB;
for(int i=9;i<=1000000;i++){
while(1){
u=(d+1.0)/i;
if(u>=limit)break;
d+=1.0;
}
u=d/i;
if(sa>limit-u){
sa=limit-u;
ansA=d;
ansB=i;
//printf("%d/%d\n",ansA,ansB);
}
}
int c=GCD(ansA,ansB);
printf("\nlast(%d %d)",ansA/c,ansB/c);
}
*問い72
nとdを正の整数として, 分数 n/d を考えよう. n<dかつHCF(n,d)=1のとき, 真既約分数と呼ぶ.
d ≤ 8について既約分数を大きさ順に並べると, 以下を得る:
1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
この集合は21個の要素をもつことが分かる.
d ≤ 1,000,000について, 真既約分数の集合は何個の要素を持つか?
オイラーのファイ関数は分母が決まった時の既約分数の数だと気がついたら一瞬で解けた問題。
しかしまあ世間ではこのレベルの問題も高校レベルなのだろう。
もっと難しい問題が解けるようになりたいものである。
#include<stdio.h>
#include<iostream>
#include<string.h>
const int up=1000000;
int data[up+1],count[up+1],dell[up+1];
__int64 calc(){
int k,add;
__int64 ans=0;
memset(data,0,sizeof(data));
memset(count,0,sizeof(count));
memset(dell,0,sizeof(dell));
for(int i=2;i<=up;i++){
//printf("\n%d ",i);
if(count[i]==0){
k=1;
for(int p=i*2;p<=up;p+=i){
data[p]+=k;
//printf("(%d %d %d)",p,k,data[p]);
k++;
count[p]++;
dell[p]++;
}
}
}
for(int i=2;i<=up;i++){
ans+=i-data[i]-1;
if(count[i]>1){
k=count[i]-1;
add=1;
}else if(count[i]==0&&dell[i]!=0){
k=-1;
add=-1;
}else{
continue;
}
for(int p=i*2;p<=up;p+=i){
data[p]-=k;
k+=add;
count[p]-=add;
}
}
return ans;
}
int main(){
std::cout<<calc();
}
*問い73
nとdを正の整数として, 分数 n/d を考えよう. n<dかつHCF(n,d)=1のとき, 真既約分数と呼ぶ.
d ≤ 8について既約分数を大きさ順に並べると, 以下を得る:
1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, 5/7, 3/4, 4/5, 5/6, 6/7, 7/8
1/3と1/2の間には3つの分数が存在することが分かる.
では, d ≤ 12,000 について真既約分数をソートした集合では, 1/3 と 1/2 の間に何個の分数があるか?
1/２～1/３の範囲に収まる分数だけ作って既約かどうか調べて足せばＯＫ。
ちょっと遅いコードだけどまあぎりぎり合格かな。
ソートという言葉が気になる。
何かもっと賢い生成方法とかあるのだろうか？
#include<stdio.h>
const int up=12000;
int GCD(int a, int b)
{
while( 1 )
{
a = a % b;
if( a == 0 )
return b;
b = b % a;
if( b == 0 )
return a;
}
}
int main(){
int c,ans=0;
for(int i=4;i<=up;i++){
for(int j=i/3+1;j<(i+1)/2;j++){
c=GCD(i,j);
if(c==1){
ans++;
}
}
}
printf("%d\n",ans);
}
*問い74
http://projecteuler.net/problem=74
各桁の階乗をとりそれを合計する操作を繰り返すと数字が循環する。
１００万以下の６０個目まで循環しない数字は何個あるか？
メモ化で少しだけ高速化。
後は愚直に計算するだけです。
#include<stdio.h>
#include<set>
#include<iostream>
__int64 perm[10];
std::set<__int64> memo;
const int up=1000000;
int lens[up]={0};
int saiki(__int64 num){
int re;
if(memo.find(num)!=memo.end()){
re=0;
}else if(num<up&&lens[(int)num]!=0){
return lens[(int)num];
}else{
memo.insert(num);
__int64 next=0,m=num;
while(num!=0){
next+=perm[(int)num%10];
num/=10;
}
//std::cout<<next<<" ";
re=saiki(next)+1;
if(m<up)lens[(int)m]=re;
}
return re;
}
int main(){
perm[0]=1;
for(int i=1;i<10;i++){
perm[i]=perm[i-1]*i;
}
int ans=0;
for(int i=1;i<up;i++){
memo.clear();
ans+=(saiki(i)==60);
}
printf("%d\n",ans);
}
*問い75
http://odz.sakura.ne.jp/projecteuler/index.php?cmd=read&page=Problem%2075
ある長さの鉄線を折り曲げたときに1通りの直角三角形を作る最短の長さは12 cmである. 他にも沢山の例が挙げられる.
12 cm: (3,4,5)
24 cm: (6,8,10)
30 cm: (5,12,13)
36 cm: (9,12,15)
40 cm: (8,15,17)
48 cm: (12,16,20)
それとは対照的に, ある長さの鉄線 (例えば20cm) は整数長さで折ることができない。 また2つ以上の折り曲げ方があるものもある. 2つ以上ある例としては、120cmの長さの鉄線を用いた場合で, 3通りの折り曲げ方がある.
120 cm: (30,40,50), (20,48,52), (24,45,51)
Lを鉄線の長さとする. 直角三角形を作るときに1通りの折り曲げ方しか存在しないような L ≤ 1,500,000 の総数を答えよ.
解法
最初aとｂの全探索でいって当然惨敗しました。
次に思いついた方法は以下のようなものでこれでまあまあの速度が出ました。
a^2+b^2=c^2を式変形すると
a^2=(c+b)(c-b)です。
素因数分解の一意性より
a^2と(c+b)(c-b)の素因数は同じとなります。
aの素因数分解をもとめればa^2の素因数は簡単にもとまります。
後はa^2の素因数分解の結果を(c+b)と(c-b)に分配することを再帰関数で全探索すれば比較的少ない計算量で答えがもとまります。
r1=(c+b),r2=(c-b)とすれば(r1+r2)/2でｃが(r1-r2)/2でｂが取り出せます。
コードはまずaを1から増加させながら試します。
a^2の素因数分解をf1で求め、その結果を使って(c-b)に分配する素因数をsaiki関数で試します。
(c-b)に分配しなかった数は(c+b)に分配されるのでc-bだけで十分です。
後はsaiki関数の計算結果から&u(){}r1,r2からbとｃの長さを求め条件を満たすなら長さmemo[a+b+c]をカウントしていくだけです。
#include<stdio.h>
#include<map>
#include<iostream>
const int up=1500000;
//const int up=100;
int memo[up+2]={0};
void saiki(int p,__int64 num,__int64 a,std::map<int,int>& nums,bool addOK){
if(num>a)return;
__int64 r1,r2,b,c;
r1=(a*a)/num;
r2=num;
c=(r1+r2)/2;
b=(r1-r2)/2;
if(a+b+c<=up&&b>=a&&addOK==true&&(r1+r2)%2==0&&(r1-r2)%2==0){
//std::cout<<a<<" "<<b<<" "<<c<<"\n";
memo[(int)(a+b+c)]++;
}
std::map<int,int>::iterator it=nums.upper_bound(p);
if(it==nums.end())return;
__int64 next=num;
for(int i=0;i<=(*it).second;i++){
saiki((*it).first,next,a,nums,i>0);
next*=(*it).first;
}
}
void f1(int a){
//素因数分解を行う
std::map<int,int> nums;
//printf("%d ",a);
int t=a;
for(int i=2;i*i<=a;){
if(a%i==0){
if(nums.find(i)==nums.end())nums[i]=0;
nums[i]+=2;
a/=i;
//printf(" %d",i);
}else{
i+=1+(i&1);
}
}
if(a>1){
//printf(" %d",a);
if(nums.find(a)==nums.end())nums[a]=0;
nums[a]+=2;
}
//printf("\n");
saiki(1,1,t,nums,true);
}
int main(){
for(int a=1;a<=up/3;a++){
f1(a);
if(a%10000==0)printf("%d ",a);
}
int ans=0;
for(int i=0;i<=up;i++){
//printf("%d %d\n",i,memo[i]);
ans+=(memo[i]==1);
}
printf("%d",ans);
}
*問い76
100の分割数を求める問題。
メモ化を使ってみたが無意味な処理だったかもしれない。
最悪のコード実行速度となっている。
今からきちんと公式の意味と式を勉強してコードを書きなおそうと思う。
#include<stdio.h>
#include<string.h>
int ans=0;
int memo[101][101];
int saiki(int n,int m){
if(memo[n][m]!=0){
ans+=memo[n][m];
return memo[n][m];
}else if(n==0){
ans++;
return 1;
}else{
int sum=0;
for(int i=m;i<=n;i++){
memo[n-i][i]=saiki(n-i,i);
sum+=memo[n][i];
}
return sum;
}
}
int main(){
memset(memo,0,sizeof(memo));
ans=0;
saiki(100,1);
printf("%d",ans-1);
}
公式を使ってみたら驚くほど速くなった。
公式すごいな。
#include<stdio.h>
#include<string.h>
int memo[101][101];
int P(int k,int n){
if(k>n)return 0;
if(k==n)return 1;
if(memo[k][n]!=-1)return memo[k][n];
memo[k][n]=P(k+1,n)+P(k,n-k);
return memo[k][n];
}
int main(){
memset(memo,-1,sizeof(memo));
printf("%d\n",P(1,100)-1);
}
*問い77
10は素数の和として5通りに表すことができる:
7 + 3
5 + 5
5 + 3 + 2
3 + 3 + 2 + 2
2 + 2 + 2 + 2 + 2
素数の和としての表し方が5000通り以上になる最初の数を求めよ。
値が小さいので全探索でもいいのですが一応メモ化で気持ち高速化しています。
値が大きくなったらメモ化の意味が出てきます。
#include<stdio.h>
#include<vector>
#include<string.h>
const int up=1000000;
std::vector<int> sosuu;
bool so[up];
int memo[10000][100];
void setSo(){
int i2;
memset(so,true,sizeof(so));
sosuu.push_back(2);
for(int i=3;i<=up;i+=2){
if(so[i]==false)continue;
sosuu.push_back(i);
i2=i*2;
for(int j=i*3;j<up;j+=i2){
so[j]=false;
}
}
}
int saiki(int n,int p){
if(n==0)return 1;
if(n==1)return 0;
if(n<10000&&p<100&&memo[n][p]!=0)return memo[n][p];
int re=0;
for(int i=p;i<sosuu.size()&&sosuu[i]<=n;i++){
re+=saiki(n-sosuu[i],i);
}
memo[n][p]=re;
return re;
}
int main(){
setSo();
memset(memo,0,sizeof(memo));
memo[2][0]=1;
for(int i=3;i<100;i++){
int t=saiki(i,0);
printf("%d %d\n",i,t);
if(t>=5000)break;
}
}
*問い78
自然数nの分割数を求める関数P(n)を考えた時p(n)が１００万で割れる最小のｎを求めよという問題。
問い76の公式を使っても遅いのでＷｉｋiの分割数の解説にあったもう一段早い処理を実装。
問い７６で解いた公式はまあなんとなくメモ化がいい具合に聞く原理は分かるのだが、こちらになると一体どういう原理なのか想像もつかない。
多分この辺でも初等整数論なんだろうけど俺には整数論は難しいすぎるな。
#include<stdio.h>
#include<vector>
#include<iostream>
const __int64 mask=1000000;
int main(){
std::vector<__int64> memo;
int s[]={1,1,2,3,5,7,11,15};
for(int i=0;i<8;i++)memo.push_back(s[i]);
int k=8;
__int64 sum;
int sa,down,t,p,base;
while(1){
sum=0;
sa=1;
t=1;
down=1;
base=1;
std::cout<<k<<"\n";
while(1){
if(k<base)break;
sum+=(t*memo[k-base])%mask;
if(k<base+sa)break;
sum+=(t*memo[k-base-sa])%mask;
down+=3;
base+=down;
sa++;
t*=-1;
//std::cout<<"("<<base<<" "<<sa<<")";
}
//std::cout<<sum<<"\n";
sum=sum%mask;
if(sum==0)break;
memo.push_back(sum);
k++;
}
std::cout<<k;
}
*問い80
1～100までの平方根のうち無理数になる数の集合を考える。
１００桁目までの各桁の和の合計はどうなるか？
このサイトを参考に問題を解く。
http://www004.upp.so-net.ne.jp/s_honma/root.htm
サイト通りだと桁数が64bit整数で溢れそうだったので苦手なＬｉｓｐで挑戦。
Ｌｉｓｐに不慣れなために構文エラーなどで恐ろしく時間がかかった。
オイラープロジェクトは学びのためのサイトだろうから調べ物しながら解くことで数学とアルゴリズムに関する知識を学ぶのが本筋だろうから調べながら勉強するのが正しい気がする。
Ｌｉｓｐで集計する方法がよくわからなかったのででてきた数字列をＥｘｃelで集計して、最後に余分な無理数でない数を引いてアセプト。
(defun mysqrt (n m sum deep)
(cond ((= deep 0)
(let ((r (floor (sqrt n))))
(mysqrt (* (- n (* r r)) 100) (* (+ r r) 10) r (+ deep 1))))
((= deep 100) sum )
(t
(dolist (x '(9 8 7 6 5 4 3 2 1 0))
(let ((r (* (+ m x) x)))
(if (<= r n)
(return (mysqrt (* (- n r) 100) (* (+ m x x) 10) (+ sum x) (+ deep 1)))))))))
mysqrt
(dotimes (x 100)
(print (mysqrt (+ x 1) 0 0 0)))
