from ortools.linear_solver import pywraplp
import matplotlib.pyplot as plt

# データ
families = [0, 1, 2, 3, 4]  # 表現の一部として家族を表すインデックス
food_required = [3, 7, 2, 5, 4]  # 各家族が必要とする食糧の量
transport_costs = [10, 20, 5, 15, 10]  # 各家族への食糧配達のコスト

# ベクトル全体で最小化するためのスカラー乗数
scalar_multiplier = 10

# 問題を作成
solver = pywraplp.Solver.CreateSolver('SCIP')

# 変数を作成
x = []
for i in families:
    x.append(solver.IntVar(0, 1, 'x[%i]' % i))  # x[i] = 家族iに食糧を配達するかどうか (0 または 1)

# 制約を定義
budget = 30  # 予算
solver.Add(sum(transport_costs[i] * x[i] for i in families) <= budget)

# 目的関数を定義
solver.Maximize(scalar_multiplier * sum(food_required[i] * x[i] for i in families) - sum(transport_costs[i] * x[i] for i in families))

# 問題を解く
solver.Solve()

# 結果を出力
print("Total food delivered = ", sum(food_required[i] * x[i].solution_value() for i in families))
delivered = [x[i].solution_value() for i in families]
plt.figure(figsize=(10,5))
plt.bar(range(len(families)), delivered)
plt.title("Food delivered to each family")
plt.xlabel("Family")
plt.ylabel("Food quantity")
plt.xticks(range(len(families)), families)
plt.show()

from ortools.linear_solver import pywraplp
import matplotlib.pyplot as plt

families = [0, 1, 2, 3, 4]  # 家族を表すインデックス
food_required = [3, 7, 2, 5, 4]  # 各家族の食糧必要量
transport_costs = [10, 20, 5, 15, 10]  # 各家族への食糧配達コスト
scalar_multiplier = 10  # 最小化するためのスカラー乗数

solver = pywraplp.Solver.CreateSolver('SCIP')

x = []
for i in families:
   x.append(solver.IntVar(0, 1, 'x[%i]' % i))  # 家族iに食糧を配達するかどうか (0 または 1)

budget = 30
solver.Add(sum(transport_costs[i] * x[i] for i in families) <= budget)

solver.Maximize(scalar_multiplier * sum(food_required[i] * x[i] for i in families) - sum(transport_costs[i] * x[i] for i in families))
solver.Solve()

print("Total food delivered = ", sum(food_required[i] * x[i].solution_value() for i in families))
delivered = [x[i].solution_value() for i in families]
plt.figure(figsize=(10,5))
plt.bar(range(len(families)), delivered)
plt.title("Food delivered to each family")
plt.xlabel("Family")
plt.ylabel("Food quantity")
plt.xticks(range(len(families)), families)
plt.show()

from scipy.spatial.distance import cdist
import matplotlib.pyplot as plt

# 緯度経度のデータを定義
# 例: Tokyo, New York, London, Sydney
cities = {
    'Tokyo': (35.6895, 139.6917),
    'New York': (40.7128, -74.0060),
    'London': (51.5074, -0.1278),
    'Sydney': (-33.8688, 151.2093)
}

# 緯度経度をnumpy配列に変換
locations = list(cities.values())

# 距離を計算（ここではユークリッド距離を使用）
distances = cdist(locations, locations)

# 結果のグラフ化
plt.imshow(distances, cmap='viridis', interpolation='nearest')
plt.colorbar(label='Distance (degrees)')
plt.xticks(ticks=range(len(cities)), labels=cities.keys())
plt.yticks(ticks=range(len(cities)), labels=cities.keys())
plt.title('Distance between Cities')
plt.show()

from scipy.spatial.distance import cdist
import matplotlib.pyplot as plt

cities = {
    'Tokyo': (35.6895, 139.6917),
    'New York': (40.7128, -74.0060),
    'London': (51.5074, -0.1278),
    'Sydney': (-33.8688, 151.2093)
}

locations = list(cities.values())

distances = cdist(locations, locations)

plt.imshow(distances, cmap='viridis', interpolation='nearest')
plt.colorbar(label='Distance (degrees)')
plt.xticks(ticks=range(len(cities)), labels=cities.keys())
plt.yticks(ticks=range(len(cities)), labels=cities.keys())
plt.title('Distance between Cities')
plt.show()

import numpy as np
from scipy.optimize import minimize
import matplotlib.pyplot as plt

# Rosenbrock関数を定義
def rosenbrock(x):
   return (1 - x[0])**2 + 100 * (x[1] - x[0]**2)**2

# Rosenbrock関数の勾配を定義
def rosenbrock_der(x):
   return np.array([-2 * (1 - x[0]) - 400 * x[0] * (x[1] - x[0]**2), 200 * (x[1] - x[0]**2)])

# Rosenbrock関数のHessianを定義
def rosenbrock_hess(x):
   return np.array([[2 + 1200 * x[0]**2 - 400 * x[1], -400 * x[0]],
                   [-400 * x[0], 200]])

# 初期値を設定
x0 = np.array([0.5, 0])

# 最適化問題を解く
res = minimize(rosenbrock, x0, method='trust-constr', jac=rosenbrock_der, hess=rosenbrock_hess)

# 結果を表示
print(res.x)

# 結果をグラフ化
x = np.linspace(-1, 2, 100)
y = np.linspace(-1, 3, 100)
X, Y = np.meshgrid(x, y)
Z = rosenbrock([X, Y])

plt.contour(X, Y, Z, levels=np.logspace(-2, 2, 20))
plt.xlabel('x')
plt.ylabel('y')
plt.title('Rosenbrock function')
plt.show()

import numpy as np
from scipy.optimize import minimize
import matplotlib.pyplot as plt

def rosenbrock(x):
   return (1 - x[0])**2 + 100 * (x[1] - x[0]**2)**2

def rosenbrock_der(x):
   return np.array([-2 * (1 - x[0]) - 400 * x[0] * (x[1] - x[0]**2), 200 * (x[1] - x[0]**2)])

def rosenbrock_hess(x):
   return np.array([[2 + 1200 * x[0]**2 - 400 * x[1], -400 * x[0]],
                   [-400 * x[0], 200]])

x0 = np.array([0.5, 0])
res = minimize(rosenbrock, x0, method='trust-constr', jac=rosenbrock_der, hess=rosenbrock_hess)

print(res.x)

x = np.linspace(-1, 2, 100)
y = np.linspace(-1, 3, 100)
X, Y = np.meshgrid(x, y)
Z = rosenbrock([X, Y])

plt.contour(X, Y, Z, levels=np.logspace(-2, 2, 20))
plt.xlabel('x')
plt.ylabel('y')
plt.title('Rosenbrock function')
plt.show()

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np

# 定数定義
G = 6.67430e-11  # 万有引力定数
m1 = 5.972e24  # 地球の質量(kg)
m2 = 7.35e22  # 月の質量(kg)
r = 384400000  # 地球から月までの距離(m)
v1 = 0  # 地球の初期速度(m/s)
v2 = 1022  # 月の初期速度(m/s)
timestep = 1000  # タイムステップ
num_steps = 10000  # ステップ数

# 初期条件
x1, y1, z1 = 0, 0, 0  # 地球の初期位置
x2, y2, z2 = r, 0, 0  # 月の初期位置

# 速度成分
vx1, vy1, vz1 = v1, 0, 0
vx2, vy2, vz2 = 0, v2, 0

# 位置、速度の更新
x1_values, y1_values, z1_values = [], [], []
x2_values, y2_values, z2_values = [], [], []

for _ in range(num_steps):
    # 地球と月の距離
    dx = x2 - x1
    dy = y2 - y1
    dz = z2 - z1
    r = np.sqrt(dx ** 2 + dy ** 2 + dz ** 2)

    # 万有引力の計算
    F = G * m1 * m2 / r ** 2

    # 重力の影響による加速度
    ax1 = F * dx / r / m1
    ay1 = F * dy / r / m1
    az1 = F * dz / r / m1
    ax2 = -F * dx / r / m2
    ay2 = -F * dy / r / m2
    az2 = -F * dz / r / m2

    # 速度の更新
    vx1 += ax1 * timestep
    vy1 += ay1 * timestep
    vz1 += az1 * timestep
    vx2 += ax2 * timestep
    vy2 += ay2 * timestep
    vz2 += az2 * timestep

    # 位置の更新
    x1 += vx1 * timestep
    y1 += vy1 * timestep
    z1 += vz1 * timestep
    x2 += vx2 * timestep
    y2 += vy2 * timestep
    z2 += vz2 * timestep

    # 位置を保存
    x1_values.append(x1)
    y1_values.append(y1)
    z1_values.append(z1)
    x2_values.append(x2)
    y2_values.append(y2)
    z2_values.append(z2)

# 3D軌道のプロット
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot(x2_values, y2_values, z2_values, label='Moon')
ax.plot([0], [0], [0], marker='o', markersize=8, color='blue', label='Earth')
ax.set_xlabel('X coordinate (m)')
ax.set_ylabel('Y coordinate (m)')
ax.set_zlabel('Z coordinate (m)')
ax.set_title('Moon Orbit around Earth')
ax.legend()
plt.show()

import plotly.graph_objects as go
import numpy as np

# ヘリックスの方程式
t = np.linspace(0, 10, 50)
x, y, z = np.cos(t), np.sin(t), t

fig = go.Figure(data=[go.Scatter3d(x=x, y=y, z=z, mode='markers')])
fig.show()

x = np.outer(np.linspace(-2, 2, 30), np.ones(30))
y = x.copy().T
z = np.cos(x ** 2 + y ** 2)

fig = go.Figure(data=[go.Surface(x=x, y=y, z=z)])
fig.show()

import plotly.graph_objects as go
import numpy as np

# ヘリックスの方程式
t = np.linspace(0, 10, 50)
x, y, z = np.cos(t), np.sin(t), t

fig = go.Figure(data=[go.Scatter3d(x=x, y=y, z=z, mode='markers')])
fig.show()

x = np.outer(np.linspace(-2, 2, 30), np.ones(30))
y = x.copy().T
z = np.cos(x ** 2 + y ** 2)

fig = go.Figure(data=[go.Surface(x=x, y=y, z=z)])
fig.show()

import numpy as np
import matplotlib.pyplot as plt

def kernel(x, y):
 return np.exp(-0.5*(x - y)**2 / 2)

def f(x):
 return 2 * x + 5 * np.cos(x)

np.random.seed(seed=4)
n = 13
m = 100
X = np.linspace(0, 10, m)
y_truth = f(X)
x = 0.5 + 8.5 * np.random.rand(n)
y = f(x)

K = np.array([[kernel(x[i], x[j]) for j in range(n)] for i in range(n)])

c = np.linalg.inv(K) @ y
y_pred = []
for i in range(m):
 sum = 0
 for j in range(n):
   sum += c[j] * kernel(x[j], X[i])
 y_pred.append(sum)
y_pred = np.array(y_pred)

fig, ax = plt.subplots()
ax.scatter(x, y, label="data")
ax.plot(X, y_pred, color="red", label="predict")
ax.plot(X, y_truth, label="ground truth")
ax.legend()
plt.show()

from constraint import Problem

def find_numbers():
    problem = Problem()

    # 変数の定義
    problem.addVariable('x', range(1, 11))
    problem.addVariable('y', range(1, 11))
    problem.addVariable('z', range(1, 11))

    # 制約の設定
    def constraint_function(x, y, z):
        if x + y + z == 20:
            return True

    problem.addConstraint(constraint_function, ['x', 'y', 'z'])

    solutions = problem.getSolutions()
    return solutions

solutions = find_numbers()
if solutions:
    print("整数の組み合わせ:")
    for solution in solutions:
        print(f"x: {solution['x']}, y: {solution['y']}, z: {solution['z']}")
else:
    print("条件を満たす組み合わせは見つかりませんでした。")

import cvxpy as cp
import numpy as np
import matplotlib.pyplot as plt

# 変数を定義
x = cp.Variable()
y = cp.Variable()

# 目的関数と制約条件を定義
objective = cp.Minimize(3*x + 4*y)
constraints = [x + 2*y >= 10, 3*x - y <= 12, x >= 0, y >= 0]

# 最適化問題の定義
problem = cp.Problem(objective, constraints)

# 問題を解く
problem.solve()

# 最適な解の表示
print("Optimal value:", problem.value)
print("x =", x.value)
print("y =", y.value)

# グラフ描画
x_values = np.linspace(0, 10, 100)
y1_values = (10 - x_values) / 2
y2_values = 3 * x_values - 12

plt.figure(figsize=(6, 6))
plt.plot(x_values, y1_values, label='x + 2y >= 10')
plt.plot(x_values, y2_values, label='3x - y <= 12')
plt.fill_between(x_values, np.maximum(0, y1_values), np.minimum(10, (12 + y2_values) / 3), color='gray', alpha=0.3)
plt.plot(x.value, y.value, 'ro')  # 最適解のプロット
plt.xlabel('x')
plt.ylabel('y')
plt.title('Optimization with CVXPY')
plt.legend()
plt.grid(True)
plt.show()

import cvxpy as cp
import numpy as np
import matplotlib.pyplot as plt

x = cp.Variable()
y = cp.Variable()

objective = cp.Minimize(3*x + 4*y)
constraints = [x + 2*y >= 10, 3*x - y <= 12, x >= 0, y >= 0]

problem = cp.Problem(objective, constraints)

problem.solve()

print("Optimal value:", problem.value)
print("x =", x.value)
print("y =", y.value)

x_values = np.linspace(0, 10, 100)
y1_values = (10 - x_values) / 2
y2_values = 3 * x_values - 12

plt.figure(figsize=(6, 6))
plt.plot(x_values, y1_values, label='x + 2y >= 10')
plt.plot(x_values, y2_values, label='3x - y <= 12')
plt.fill_between(x_values, np.maximum(0, y1_values), np.minimum(10, (12 + y2_values) / 3), color='gray', alpha=0.3)
plt.plot(x.value, y.value, 'ro')  # 最適解のプロット
plt.xlabel('x')
plt.ylabel('y')
plt.title('Optimization with CVXPY')
plt.legend()
plt.grid(True)
plt.show()

from pulp import *

# ビンの数とアイテムの数を設定
bin_count = 2
item_count = 5

# 各ビンの容量、アイテムの重さを設定
bin_capacity = [10, 10]
item_weight = [2, 5, 4, 7, 1]

# 問題設定
prob = LpProblem("BinPacking", LpMinimize)

# アイテムをビンに入れるかどうかを示すノードを作成（バイナリ変数）
x = [[LpVariable("x(%s,%s)" % (i, j), cat="Binary") for j in range(item_count)] for i in range(bin_count)]

# 目的関数: 使用するビンの数を最小化する
prob += lpSum(x[i][j] for i in range(bin_count) for j in range(item_count))

# 制約条件: 各アイテムは1つのビンにしか入らない
for j in range(item_count):
    prob += lpSum(x[i][j] for i in range(bin_count)) == 1

# 制約条件: 各ビンの容量を超えてアイテムを入れられない
for i in range(bin_count):
    prob += lpSum(item_weight[j] * x[i][j] for j in range(item_count)) <= bin_capacity[i]

# 最適化計算
prob.solve()

# 結果表示
for i in range(bin_count):
    print("*****Bin", i, "*****")
    total_weight = 0
    for j in range(item_count):
        if value(x[i][j])==1:
            print("Item", j)
            total_weight+=item_weight[j]
    print("Total weight:", total_weight)

from pulp import *

bin_count = 2
item_count = 5
bin_capacity = [10, 10]
item_weight = [2, 5, 4, 7, 1]

prob = LpProblem("BinPacking", LpMinimize)

x = [[LpVariable("x(%s,%s)" % (i, j), cat="Binary") for j in range(item_count)] for i in range(bin_count)]

prob += lpSum(x[i][j] for i in range(bin_count) for j in range(item_count))

for j in range(item_count):
    prob += lpSum(x[i][j] for i in range(bin_count)) == 1

for i in range(bin_count):
    prob += lpSum(item_weight[j] * x[i][j] for j in range(item_count)) <= bin_capacity[i]

prob.solve()

for i in range(bin_count):
    print("*****Bin", i, "*****")
    total_weight = 0
    for j in range(item_count):
        if value(x[i][j])==1:
            print("Item", j)
            total_weight+=item_weight[j]
    print("Total weight:", total_weight)

import pandas as pd
import numpy as np
from prophet import Prophet
import matplotlib.pyplot as plt

# ランダムな株価データを生成
np.random.seed(42)
dates = pd.date_range(start='2023-01-01', periods=365, freq='D')
prices = np.cumsum(np.random.randn(365)) + 100  # ランダムな価格変動

data = {
    'ds': dates,
    'y': prices
}

df = pd.DataFrame(data)

# Prophetモデルの作成と学習
model = Prophet()
model.fit(df)

# 予測用データの作成
future = model.make_future_dataframe(periods=30)  # 30日間の予測
forecast = model.predict(future)

# 予測結果をグラフ化
fig = model.plot(forecast)
plt.xlabel('Date')
plt.ylabel('Stock Price')
plt.title('Random Stock Price Forecast')
plt.show()

np.random.seed(42)
dates = pd.date_range(start='2023-01-01', periods=365, freq='D')
prices = np.cumsum(np.random.randn(365)) + 100  # ランダムな価格変動

data = {
    'ds': dates,
    'y': prices
}

df = pd.DataFrame(data)

model = Prophet()
model.fit(df)

future = model.make_future_dataframe(periods=30)  # 30日間の予測
forecast = model.predict(future)

fig = model.plot(forecast)
plt.xlabel('Date')
plt.ylabel('Stock Price')
plt.title('Random Stock Price Forecast')
plt.show()