RでProject Euler (Problem 31~40)

[Problem 21~30](https://qiita.com/Rigel_/items/215f1e27d7c59384f5ee)
[Problem 41~50](https://qiita.com/Rigel_/items/484196dbf1dfa2e391c6)

# Problem 31 「硬貨の和」

>イギリスでは硬貨はポンド£とペンスpがあり，一般的に流通している硬貨は以下の8種類である.
1p, 2p, 5p, 10p, 20p, 50p, £1 (100p) and £2 (200p).
以下の方法で£2を作ることが可能である．
1×£1 + 1×50p + 2×20p + 1×5p + 1×2p + 3×1p
これらの硬貨を使って£2を作る方法は何通りあるか?

1pと2pで$m$p支払う方法は、2pを0枚, 1枚, ..., $\lfloor m/2 \rfloor$枚があるので、合計$\lfloor m/2 \rfloor+1$通りある。
よって、5p以上を考えて、あとは再帰を使って求める。

```R
count = function(m, coins){
  if(length(coins)==0) return(m%/%2+1)
  p = 0
  for(i in 0:(m/coins[1]))
    p = p + count(m-i*coins[1], coins[-1])
  return(p)
}
ans = count(200, c(200, 100, 50, 20, 10, 5))
print(ans)
```

# Problem 32 「パンデジタル積」

>すべての桁に$1$から$n$が一度だけ使われている数を$n$桁の数がパンデジタル(pandigital)であるということにしよう: 例えば$5$桁の数$15234$は$1$から$5$のパンデジタルである.
$7254$は面白い性質を持っている.$39 × 186 = 7254$と書け, 掛けられる数, 掛ける数, 積が$1$から$9$のパンデジタルとなる.
掛けられる数/掛ける数/積が$1$から$9$のパンデジタルとなるような積の総和を求めよ.
HINT: いくつかの積は, 1通り以上の掛けられる数/掛ける数/積の組み合わせを持つが1回だけ数え上げよ.

掛けられる数(左の数) < 掛ける数(右の数)とする。
掛けられる数, 掛ける数, 積が9桁になるのは、
・掛けられる数:1桁, 掛ける数:4桁, 積:4桁
・掛けられる数:2桁, 掛ける数:3桁, 積:4桁
のどちらかしかない。
それぞれ考えられる全パターン作り、条件を満たす数を探す。

```R
data = rbind(expand.grid(i=1:9, j=1234:9876), expand.grid(i=12:98, j=123:987))
data$ij = data$i*data$j
data = data[(data$ij<10000) & !grepl(0, data$i) & !grepl(0, data$j) & !grepl(0, data$ij), ]
pas = as.numeric(apply(data, 1, paste, collapse=""))
sum(unique(data$ij[Reduce("&", lapply(1:9, function(x) grepl(x, pas)))]))

```

# Problem 33 「桁消去分数」

>$49/98$は面白い分数である.「分子と分母からそれぞれ$9$を取り除くと,$49/98 = 4/8$となり, 簡単な形にすることができる」と経験の浅い数学者が誤って思い込んでしまうかもしれないからである. (方法は正しくないが，$49/98$の場合にはたまたま正しい約分になってしまう．)
我々は$30/50 = 3/5$のようなタイプは自明な例だとする.
このような分数のうち, $1$より小さく分子・分母がともに2桁の数になるような自明でないものは, 4個ある.
その4個の分数の積が約分された形で与えられたとき, 分母の値を答えよ.

この正しくない計算を桁消去と呼ぶことにする。
桁消去する前の数を$(10a+b)/(10c+d)$と表す($a, b, c, d$は$1$から$9$の自然数)。
このとき、$a, b, c, d$の組み合わせは$9^4$通りあり、全通り作成し$10a+b<10c+d$かつ以下のどちらか条件を満たす数を探す。
・$b=c$かつ$(10a+b)/(10c+d)=a/d$
・$a=d$かつ$(10a+b)/(10c+d)=b/c$

```R
data = expand.grid(a=1:9, b=1:9, c=1:9, d=1:9)
data = data[10*data$a+data$b<10*data$c+data$d, ]
data = data[(data$b==data$c & (10*data$a+data$b)/(10*data$c+data$d)==data$a/data$d) | (data$a==data$d & (10*data$a+data$b)/(10*data$c+data$d)==data$b/data$c), ]
gcd = function(a, b){
  if(b==0) return(a)
  return(gcd(b, a%%b))
}
ans = prod(10*data$c + data$d)/gcd(prod(10*data$a + data$b), prod(10*data$c + data$d))
print(ans)
```

# Problem 34 「桁の階乗」

>$145$は面白い数である.$1! + 4! + 5! = 1 + 24 + 120 = 145$となる.
各桁の数の階乗の和が自分自身と一致するような数の和を求めよ.
注:$1! = 1$と$2! = 2$は総和に含めてはならない.

$9! = 362880$なので、各桁の数の階乗の和が自分自身と一致するような数は7桁以下。
仮に8桁だとすると、各桁すべて9だとしても$8\cdot9! = 2903040$で7桁までしかいかないため、ありえない。
各桁の数の階乗の和を全パターン作成し、自分自身の番号と一致する数を探す。
以下のRの1行目では例えば$1$は、$0!+0!+0!+0!+0!+0!+1!$と計算されている。
そのため2行目で$0!$の部分を引いている。

```R
d = Reduce("+", lapply(1:7, function(x) rep(rep(gamma(1:10), each=10^(x-1)), times=10^(7-x))))
d = d - c(0, unlist(lapply(1:7, function(x) rep(7-x, 9*10^(x-1)))))
ans = d[d==0:(length(d)-1)]
ans = sum(ans[-(1:2)])
print(ans)
```

# Problem 35 「巡回素数」

>$197$は巡回素数と呼ばれる. 桁を回転させたときに得られる数$197, 971, 719$が全て素数だからである.
100未満には巡回素数が13個ある:$2, 3, 5, 7, 11, 13, 17, 31, 37, 71, 73, 79,$および$97$である.
100万未満の巡回素数はいくつあるか?

2桁以上の巡回素数には、$0, 2, 4, 5, 6, 8$が含まれない。
仮に含まれていると、1の位がそのような数になるように巡回したとき、素数ではなくなる。
従って以下のステップで巡回素数を探す。

ステップ1. エラトステネスの篩で100万未満の素数を探す。
~~ステップ2. その中から、(1桁)or(0, 2, 4, 5, 6, 8が含まれていない)素数に絞る。~~
ステップ3. 最後に巡回させても素数のものを探す。
(ステップ2は行わないほうが計算が早かった。)

```R
getprimes = function(limit){
  p = rep(TRUE, limit)
  p[1] = FALSE
  for(i in c(2, seq(3, sqrt(limit), 2)))
    if(p[i])
      p[seq(i*2, limit, i)] = FALSE
  return(as.numeric(which(p)))
}
limit = 1000000
p = getprimes(limit)
for(i in 1:10){
  pp = (p%%10)*10^floor(log10(p)) + p%/%10
  p = pp[pp%in%p]
}
ans = length(p)
print(ans)
```

# Problem 36 「二種類の基数による回文数」

>$585 = 1001001001_2$(2進) は10進でも2進でも回文数である.
100万未満で10進でも2進でも回文数になるような数の総和を求めよ.
(注: 先頭に$0$を含めて回文にすることは許されない.)

100万未満の数の中から、10進表記で(上位1桁=下位1桁)and(上位2桁=下位2桁)and(上位3桁=下位3桁)を満たす数を探す。
その後に、2進表記で(上位1桁=下位1桁)and(上位2桁=下位2桁)and...and(上位10桁=下位10桁)を満たす数を探す。
($2^{20}>10^6$なので、10桁目までで十分)

ある自然数$n$を$p$進数表記したときに、
・上位$d$桁の数は$\lfloor n/p^{\lfloor \log_{p}{n} \rfloor+1-d} \rfloor \bmod p$
・下位$d$桁の数は$\lfloor n/p^{d-1} \rfloor \bmod p$
で求められる。

```R
v = 1:1000000
for(d in 1:3)
  v = v[(v%/%(10^(d-1))%%10)==(v%/%(10^(floor(log10(v))+1-d))%%10)]
for(d in 1:10)
  v = v[(v%/%(2^(d-1))%%2)==(v%/%(2^(floor(log2(v))+1-d))%%2)]
ans = sum(v)
print(ans)
```

# Problem 37 「切り詰め可能素数」

>$3797$は面白い性質を持っている. まずそれ自身が素数であり, 左から右に桁を除いたときに全て素数になっている($3797, 797, 97, 7$). 同様に右から左に桁を除いたときも全て素数である($3797, 379, 37, 3$).
右から切り詰めても左から切り詰めても素数になるような素数は11個しかない. 総和を求めよ.
注: $2, 3, 5, 7$を切り詰め可能な素数とは考えない.

1桁の素数$2, 3, 5, 7$の右側に数字を加えていき、新たな素数を作っていくことを考える。
$v = (2, 3, 5, 7)$とし、以下の操作を$v$がからになるまで繰り返す。
"$v$のそれぞれの右側に$1, 3, 7, 9$を加えた数のうち素数である数を新たな$v$とする"
この操作を繰り返していく過程で、$v$を保存しておく。
こうして得られた数は、右から左に桁を除いたとき全て素数である。
そして次は、これらの数の中から左から右に桁を除いたときに全て素数であるものを探す。


```R
getprimes = function(limit){
  p = rep(TRUE, limit)
  p[1] = FALSE
  for(i in c(2, seq(3, sqrt(limit), 2)))
    if(p[i])
      p[seq(i*2, limit, i)] = FALSE
  return(as.numeric(which(p)))
}
x = c()
v = c(2, 3, 5, 7)
while(length(v)>0){
  v = c(outer(v*10, c(1, 3, 7, 9), "+"))
  v = v[!Reduce("|", lapply(seq(3, sqrt(max(v))), function(x) (v%%x==0)))]
  x = c(x, v)
}
limit = max(x)%%(10^floor(log10(max(x))))
p = c(0, getprimes(limit))
for(i in 1:floor(log10(max(x))))
  x = x[x%%(10^pmax(floor(log10(x)-i+1), 0)) %in% p]
ans = sum(x)
print(ans)
```

# Problem 38 「パンデジタル倍数」

>$192$に$1, 2, 3$を掛けてみよう.
$192 × 1 = 192$
$192 × 2 = 384$
$192 × 3 = 576$
積を連結することで$1$から$9$の パンデジタル数$192384576$が得られる.$192384576$を$192$と$(1,2,3)$の連結積と呼ぶ.
同じようにして,$9$を$1,2,3,4,5$と掛け連結することでパンデジタル数$918273645$が得られる. これは$9$と$(1,2,3,4,5)$との連結積である.
整数と$(1,2,...,n) (n > 1)$との連結積として得られる9桁のパンデジタル数の中で最大のものはいくつか?

掛けられる数(左の数)は4桁以下を調べれば良い。
行列を使い、全パターン作成し、条件で絞っていく。

```R
m = (1:9999)%o%(1:9)
m = m[apply(t(apply(floor(log10(m))+1, 1, cumsum))==9, 1, any), ]
class(m) = "character"
ans = substr(apply(m, 1, paste, collapse=""), 1, 9)
ans = as.numeric(ans)
ans = ans[!grepl(0, ans)]
ans = max(ans[Reduce("&", lapply(1:9, function(x) grepl(x, ans)))])
print(ans)
```

# Problem 39 「整数の直角三角形」

>辺の長さが$\\{a,b,c\\}$と整数の3つ組である直角三角形を考え, その周囲の長さを$p$とする.$p = 120$のときには3つの解が存在する:
${20,48,52}, {24,45,51}, {30,40,50}$
$p ≤ 1000$のとき解の数が最大になる$p$はいくつか?

Problem 9と同じ考え方で、解の数を探す。
$a+b+c=p$ $(a<b<c, p\geq 12)$と$a^2+b^2=c^2$を満たす自然数$\\{a,b,c\\}$の組の数は、
$$\frac{p(p-2a)}{2(p-a)}\ \ (3 \leq a \leq \lfloor p/3 \rfloor-1)$$が整数となるものの数と一致する。


```R
limit = 1000
p = seq(12, limit, 2)
count = function(p){
  a = 3:((p/3)-1)
  return(sum((p*(p-2*a))%%(2*(p-a))==0))
}
ans = p[which.max(unlist(lapply(p,count)))]
print(ans)
```

# Problem 40 「チャンパーノウン定数」

>正の整数を順に連結して得られる以下の10進の無理数を考える:
$0.123456789101112131415161718192021...$
小数第$12$位は$1$である.
$d_{n}$で小数第$n$位の数を表す.$d_{1} × d_{10} × d_{100} × d_{1000} × d_{10000} × d_{100000} × d_{1000000}$を求めよ.

$1$桁の数は$1\sim9$まであり、合計$9$桁。
$2$桁の数は$10\sim99$まであり、合計$180$桁。
一般に$m$桁の数は$10^{m-1}\sim10^m-1$まであり、合計$m(10^m-10^{m-1})$桁。
あとは、$n$から何桁の数字の何番目の数の上位何桁の数かがわかる。
あまりで上位何桁かを求める際上位$0$桁となった場合は、1つ前の数字の下位$1$桁にずらす。

```R
g = function(n){
  m = 0:(floor(log10(n))+1)
  s = cumsum(m*(10^m-10^(m-1)))
  k = max(which(s<n))
  q = 10^(k-1) + (n-s[k])%/%(k)
  r = (n-s[k])%%k
  if(r==0){
    r = r + k
    q = q - 1
  }
  return(as.numeric(substr(as.character(q), r, r)))
}
ans = prod(unlist(lapply(10^(0:6), g)))
print(ans)
```