# パッケージがインストール済なら次の行は不要
# import Pkg; Pkg.add("Plots")
using Plots
import Random  # 乱数を扱うためのモジュールをインポート
n = 500  # 試行回数
Random.seed!(1)  # 乱数のシード値を固定
# Pi=(xi[i], eta[i])
# U(0,1)に従う乱数をn個生成し，xiに代入
xi = rand(n)
# U(0,1)に従う乱数をn個生成し，etaに代入
eta = rand(n)
function ND(n)  # 自然数nに対してNn(D)を返す関数
    s = 0
    for i in 1:n
        if xi[i]^2 + eta[i]^2 <= 1
            s += 1
        end
    end
    return s
end
Z = zeros(n)  # Nn(D)/nを格納する配列
for i in 1:n
    Z[i] = ND(i) / i
end
n_list = collect(1:n)
plt_a = plot(n_list, Z, label = "Nn(D)/n", xlab = "n", ylab = "Nn(D)/n")
plt_a = hline!([pi / 4], label = "μ", linestyle = :dash)  # 理論値
# Pのプロット
plt_b = scatter(xi, eta, label = "P", xlab = "x", ylab = "y")
x_list = collect(0:0.01:1)
circle(x) = sqrt(1 - x^2)  # x^2 + y^2 = 1のグラフ
plt_b = plot!(x_list, circle.(x_list), label = "", linestyle = :dash)
# 図をまとめて表示
plts = plot(plt_a, plt_b, layout=(1,2), size=(800,400))