作业

3.2

点击查看代码

def X(n):  
    # 差分方程的解  
    return 2 * (-1)**(n + 1)  
  
 
n_values = [0, 1, 2, 3, 4, 5]  
for n in n_values:  
    print(f"X({n}) = {X(n)}")
 
print("学号：3127")

3.3

点击查看代码

import numpy as np
 
from scipy.sparse.linalg import eigs
 
import pylab as plt
 
  
 
w = np.array([[0, 1, 0, 1, 1, 1],
              [0, 0, 0, 1, 1, 1],
              [1, 1, 0, 1, 0, 0],
              [0, 0, 0, 0, 1, 1],
              [0, 0, 1, 0, 0, 1],
              [0, 0, 1, 0, 0, 0]])
 
r = np.sum(w,axis=1,keepdims=True)
 
n = w.shape[0]
 
d = 0.85
 
P = (1-d)/n+d*w/r #利用矩阵广播
 
w,v = eigs(P.T,1) #求最大特征值及对应的特征向量
 
v = v/sum(v)
 
v = v.real
 
print("最大特征值为：",w.real)
 
print("归一化特征向量为：\n",np.round(v,4))
 
plt.bar(range(1,n+1),v.flatten(),width=0.6)
 
plt.show()
 
 
print("学号：3127")

4.3

点击查看代码

import matplotlib.pyplot as plt
import numpy as np
import cvxpy as cp
 
x=cp.Variable(6,pos=True)
obj=cp.Minimize(x[5])
a1=np.array([0.025, 0.015, 0.055, 0.026])
a2=np.array([0.05, 0.27, 0.19, 0.185, 0.185])
a3=np.array([1, 1.01, 1.02, 1.045, 1.065])
k=0.05; kk=[]; qq=[]
while k<0.27:
    con=[cp.multiply(a1,x[1:5])-x[5]<=0,a2@x[:-1]>=k, a3@x[:-1]==1]
    prob=cp.Problem(obj,con)
    prob.solve(solver='GLPK_MI')
    kk.append(k); qq.append(prob.value)
    k=k+0.005
 
plt.rc('text',usetex=False); plt.rc('font',size=16); plt.rc('font',family='SimHei')
plt.plot(kk,qq,'k')
plt.plot(kk,qq,'b.')
plt.xlabel("收益 k"); plt.ylabel("风险 Q",rotation=0)
plt.show()
 
 
print("学号：3127")

4.4

点击查看代码

MAX_A = 15 
MAX_B = 24  
MAX_DEBUG = 5 
  
 
products = [  
    {"name": "Ⅰ", "A_hours": 1, "B_hours": 6, "debug_hours": 1, "profit": 2},  # 假设产品Ⅰ至少使用1小时设备A  
    {"name": "Ⅱ", "A_hours": 5, "B_hours": 2, "debug_hours": 1, "profit": 1}  
]  
  
 
max_profit = 0  
best_plan = {}  
  
 
for i in range(MAX_A // products[0]["A_hours"] + 1):  
    for j in range(MAX_B // products[1]["B_hours"] + 1):  
        # 计算调试时间是否足够  
        if (i + j) * max(products[0]["debug_hours"], products[1]["debug_hours"]) > MAX_DEBUG:  
            continue   
  
 
        total_A_hours = i * products[0]["A_hours"] + j * products[1]["A_hours"]  
        total_B_hours = i * products[0]["B_hours"] + j * products[1]["B_hours"]  
  
 
        if total_A_hours > MAX_A or total_B_hours > MAX_B:  
            continue  
  
 
        total_profit = i * products[0]["profit"] + j * products[1]["profit"]  
  
 
        if total_profit > max_profit:  
            max_profit = total_profit  
            best_plan = {"Ⅰ": i, "Ⅱ": j}  
  
 
print(f"最优生产计划：产品Ⅰ生产{best_plan['Ⅰ']}件，产品Ⅱ生产{best_plan['Ⅱ']}件")  
print(f"最大利润为：{max_profit}元")
 
 
print("学号：3127")

5.4

点击查看代码

import numpy as np  
from scipy.optimize import minimize  
  
def objective(x):  
    return -np.sum(np.sqrt(x) * np.arange(1, 101))  
  
def constraint1(x):  
    return x[1] - 10  
  
def constraint2(x):  
    return 20 - (x[1] + 2*x[2])  
  
def constraint3(x):  
    return 30 - (x[1] + 2*x[2] + 3*x[3])  
  
def constraint4(x):  
    return 40 - (x[1] + 2*x[2] + 3*x[3] + 4*x[4])  
  
def constraint5(x):  
    return 1000 - np.dot(x, np.arange(1, 101))  
  
constraints = [  
    {'type': 'ineq', 'fun': constraint1},  
    {'type': 'ineq', 'fun': constraint2},  
    {'type': 'ineq', 'fun': constraint3},  
    {'type': 'ineq', 'fun': constraint4},  
    {'type': 'ineq', 'fun': constraint5}  
]  
  
bounds = [(0, None)] * 100  
x0 = np.ones(100) * 0.1 
  
result = minimize(objective, x0, method='SLSQP', constraints=constraints, bounds=bounds)  
  
print('Optimal solution:', result.x)  
print('Objective function value at optimal solution:', -result.fun)
 
print("学号：3127")

5.5

点击查看代码

import numpy as np  
from scipy.optimize import minimize  
  
def objective(x):  
    return 2*x[0] + 3*x[0]**2 + 3*x[1] + x[1]**2 + x[2]  
  
def constraint1(x):  
    return 10 - (x[0] + 2*x[0]**2 + x[1] + 2*x[1]**2 + x[2])  
  
def constraint2(x):  
    return 50 - (x[0] + x[0]**2 + x[1] + x[1]**2 - x[2])  
  
def constraint3(x):  
    return 40 - (2*x[0] + x[0]**2 + 2*x[1] + x[2])  
  
 
def constraint4(x):  
    return x[0]**2 + x[2] - 2  
  
def constraint5(x):  
    return 1 - (x[0] + 2*x[1])  
  
constraints = [  
    {'type': 'ineq', 'fun': constraint1},  
    {'type': 'ineq', 'fun': constraint2},  
    {'type': 'ineq', 'fun': constraint3},  
    # {'type': 'eq', 'fun': constraint4}, 
    {'type': 'ineq', 'fun': constraint5}  
]  
  
 
bounds = [(0, None)] * 3  
x0 = np.array([0.1, 0.1, 0.1])  
  
result = minimize(objective, x0, method='SLSQP', constraints=constraints, bounds=bounds)  
  
print('Optimal solution:', result.x)  
print('Objective function value at optimal solution:', result.fun)  
 
print("学号：3127")

5.7

点击查看代码

import numpy as np  
  
demands = [40, 60, 80]  
max_production = 100  
total_demand = sum(demands)  
  
dp = np.full((4, total_demand + 1), float('inf'))  
dp[0][0] = 0  
  
prev_production = np.full((4, total_demand + 1), -1) 
  
for i in range(1, 4):  
    prev_demand = sum(demands[:i-1])  
    for j in range(total_demand + 1):  
        if j < prev_demand + demands[i-1]:  
            
            continue  
        for x in range(max(0, j - prev_demand - demands[i-1] + 1), min(max_production + 1, j - prev_demand + 1)):  
            production_cost = 50 * x + 0.2 * x**2  
            storage_cost = 4 * (j - prev_demand - x)  
            total_cost = dp[i-1][j-x] + production_cost + storage_cost  
            if total_cost < dp[i][j]:  
                dp[i][j] = total_cost  
                prev_production[i][j] = x  
   
min_cost = float('inf')  
final_state = -1  
for j in range(total_demand, total_demand + 1):  
    if dp[3][j] < min_cost:  
        min_cost = dp[3][j]  
        final_state = j  
  
production_plan = [0] * 3  
current_state = final_state  
for i in range(3, 0, -1):  
    production_plan[i-1] = prev_production[i][current_state]  
    current_state -= prev_production[i][current_state]  
 
print(f"最小总费用为: {min_cost} 元")  
print("生产计划为:")  
for i, plan in enumerate(production_plan, 1):  
    print(f"第{i}季度生产: {plan} 台")
 
print("学号：3127)

6.2

点击查看代码

edges = [  
    ("Pe", "T", 13),  
    ("Pe", "N", 68),  
    ("Pe", "M", 78),  
    ("Pe", "L", 51),  
    ("Pe", "Pa", 51),  
    ("T", "N", 68),  
    ("T", "M", 70),  
    ("T", "L", 60),  
    ("T", "Pa", 61),  
    ("N", "M", 21),  
    ("N", "L", 35),  
    ("N", "Pa",36),  
    ("M", "L", 56),  
    ("M", "Pa", 57),  
    ("L", "Pa", 21),  
]
class UnionFind:  
    def __init__(self, n):  
        self.parent = list(range(n))  
        self.rank = [0] * n  
  
    def find(self, u):  
        if self.parent[u] != u:  
            self.parent[u] = self.find(self.parent[u])  
        return self.parent[u]  
  
    def union(self, u, v):  
        root_u = self.find(u)  
        root_v = self.find(v)  
        if root_u != root_v:  
            if self.rank[root_u] > self.rank[root_v]:  
                self.parent[root_v] = root_u  
            elif self.rank[root_u] < self.rank[root_v]:  
                self.parent[root_u] = root_v  
            else:  
                self.parent[root_v] = root_u  
                self.rank[root_u] += 1  
  
def kruskal(edges):  
    edges.sort(key=lambda edge: edge[2]) 
    uf = UnionFind(len(set(city for edge in edges for city in edge[:2])))  
    mst = []  
      
    for u, v, weight in edges:  
        if uf.find(u) != uf.find(v):  
            uf.union(u, v)  
            mst.append((u, v, weight))  
      
    return mst  
 
city_map = {  
    "Pe": 0, "T": 1, "N": 2, "M": 3, "L": 4, "Pa": 5  
}  
 
edges_with_indices = [(city_map[u], city_map[v], weight) for u, v, weight in edges]  
 
mst_edges = [(list(city_map.keys())[u], list(city_map.keys())[v], weight) for u, v, weight in kruskal(edges_with_indices)]  
  
print("最小生成树的边:")  
for u, v, weight in mst_edges:  
    print(f"{u} - {v}: {weight}")
 
print("学号：3127")

6.3

点击查看代码

import numpy as np  
   
distances = np.array([  
    [0, 2, 7, np.inf, np.inf, np.inf],  
    [2, 0, 4, 6, 8, np.inf],           
    [7, 4, 0, 1, 3, np.inf],          
    [np.inf, 6, 1, 0, 1, 6],           
    [np.inf, 8, 3, 1, 0, 3],          
    [np.inf, np.inf, np.inf, 6, 3, 0]  
], dtype=float)  
 
students = np.array([50, 40, 60, 20, 70, 90])  
  
hospital_distances_sum = np.zeros(6)  
for i in range(6):   
    connected_distances = distances[i, :i+1].copy()  
    connected_distances = connected_distances[connected_distances != np.inf]  
    hospital_distances_sum[i] = np.sum(connected_distances)  
hospital_location = np.argmin(hospital_distances_sum)  
print(f"医院应该建在村庄 {chr(65 + hospital_location)} 处，使得最远村庄的人到医院看病所走的路最短。")  
  
school_total_distances = np.zeros(6)  
for i in range(6):  
   
    weighted_distances = 0  
    for j in range(6):  
        if distances[j, i] != np.inf:  
            weighted_distances += students[j] * distances[j, i]  
    school_total_distances[i] = weighted_distances  
  
school_location = np.argmin(school_total_distances)  
print(f"小学应该建在村庄 {chr(65 + school_location)} 处，使得所有学生上学走的总路程最短。")
 
print("学号：3127")

6.4

点击查看代码

import heapq  
  
def prim(graph, start):  
    num_nodes = len(graph)  
    visited = [False] * num_nodes  
    min_heap = [(0, start, -1)]   
    mst_cost = 0  
    mst_edges = []  
  
    while min_heap:  
        weight, u, parent = heapq.heappop(min_heap)  
        if visited[u]:  
            continue  
        visited[u] = True  
        mst_cost += weight  
        if parent != -1:  
            mst_edges.append((parent, u, weight))  
  
        for v in range(num_nodes):  
            if not visited[v] and graph[u][v] != 0:  
                heapq.heappush(min_heap, (graph[u][v], v, u))  
  
    return mst_cost, mst_edges  
    
graph = [  
    [0,20,0,0,15,0],  
    [20,0,20,60,25,0],  
    [0,20,0,30,18,0],  
    [0,60,30,0,35,10],  
    [0,0,0,10,15,0]  
]  
  
mst_cost, mst_edges = prim(graph, 0)  
print("Prim's MST Cost:", mst_cost)  
print("Prim's MST Edges:", mst_edges)
 
print("学号：3127")

6.5

点击查看代码

import numpy as np  
  
matches = np.array([  
    [0, 1, 0, 1, 1, 1],  # 1队  
    [0, 0, 0, 1, 1, 1],  # 2队  
    [1, 1, 0, 1, 0, 0],  # 3队  
    [0, 0, 0, 0, 1, 1],  # 4队  
    [0, 0, 1, 0, 0, 1],  # 5队  
    [0, 0, 1, 0, 0, 0]   # 6队  
], dtype=int)  
  
n = matches.shape[0]  
closure = matches.copy()  
for k in range(n):  
    for i in range(n):  
        for j in range(n):  
            closure[i, j] = closure[i, j] or (closure[i, k] and closure[k, j])  
  
strength = closure.sum(axis=1)  
   
ranking = np.argsort(-strength) 
  
for i, rank in enumerate(ranking):  
    print(f"{chr(65 + rank)}队 排名 {i + 1}")
    
    
import numpy as np  
from scipy.sparse import csr_matrix  
  
edges = [  
    (0, 1), (0, 3), (0, 4), (0, 5),  # 1队胜  
    (1, 3), (1, 4), (1, 5),          # 2队胜  
    (2, 0), (2, 1), (2, 3),          # 3队胜  
    (3, 4), (3, 5),                  # 4队胜  
    (4, 2), (4, 5),                  # 5队胜  
    (5, 2)                           # 6队胜  
]  
  
 
num_teams = 6  
  
 
row_ind = []  
col_ind = []  
data = []  
for u, v in edges:  
    row_ind.append(u)  
    col_ind.append(v)  
    data.append(1)  
adj_matrix = csr_matrix((data, (row_ind, col_ind)), shape=(num_teams, num_teams))  
  
 
adj_matrix_T = adj_matrix.T  
  
 
d = 0.85  
out_degree = np.array(adj_matrix_T.sum(axis=1)).flatten()  
out_degree[out_degree == 0] = 1  
M = adj_matrix_T.multiply(1.0 / out_degree).tocsr()  
M = M + (1 - d) / num_teams * csr_matrix(np.ones((num_teams, num_teams)))  
  
 
R = np.ones(num_teams) / num_teams  
  
 
num_iterations = 100  
for _ in range(num_iterations):  
    R = R.dot(M.toarray())  
  
 
pagerank_ranking = np.argsort(-R) 
  
 
for i, rank in enumerate(pagerank_ranking):  
    print(f"{chr(65 + rank)}队 PageRank排名 {i + 1}")
 
print("学号：3127")