数据处理

遍历每一个文件夹

把处理后的数据放在原文件夹

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
import os
# 创建一个空的DataFrame，用于存储所有txt文件的数据
final_df = pd.DataFrame()
# 指定包含txt文件的文件夹路径
txt_folder = 'data/参考数据集/'

# 遍历文件夹中的每个txt文件

for subdir in range(1,11):
    subpath = os.path.join(txt_folder,str(subdir))
    for filename in os.listdir(subpath):
        if filename.endswith('.txt'):
            txt_file_path = os.path.join(subpath,filename)
            data_np = np.loadtxt(txt_file_path)
            data_np = data_np.T
            data_df = pd.DataFrame(data_np)
            data_df.to_excel(txt_file_path.replace('.txt', '.xlsx'), index=False)
print("数据已导入到Excel文件")

把数据按照原逻辑放在一个新文件夹

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
import os
import numpy as np
import pandas as pd

# 指定基础文件夹路径
base_folder = 'data/参考数据集/'

# 创建一个新的文件夹用于存放处理后的Excel文件
clean_folder = 'data/clean/参考数据集/'
if not os.path.exists(clean_folder):
    os.makedirs(clean_folder)

# 遍历基础文件夹中的子文件夹
for subdir in os.listdir(base_folder):
    subdir_path = os.path.join(base_folder, subdir)
    if os.path.isdir(subdir_path):
        # 创建一个新的子文件夹用于存放处理后的Excel文件
        new_subfolder = os.path.join(clean_folder, subdir)
        if not os.path.exists(new_subfolder):
            os.makedirs(new_subfolder)

        # 遍历当前子文件夹中的txt文件
        for filename in os.listdir(subdir_path):
            if filename.endswith('.txt'):
                txt_file_path = os.path.join(subdir_path, filename)
                data_np = np.loadtxt(txt_file_path)
                data_np = data_np.T
                data_df = pd.DataFrame(data_np)

                # 生成新的Excel文件路径，存放在新的子文件夹中
                new_excel_path = os.path.join(new_subfolder, filename.replace('.txt', 'A_new.xlsx'))

                # 将DataFrame保存为Excel文件
                data_df.to_excel(new_excel_path, sheet_name='A', index=False)

print("数据已导入到Excel文件，并按照原文件夹逻辑存放在新的文件夹中")

pandas专区

numpy专区

random.random()和random.randn()的区别

numpy.random.randn 是 NumPy 库中用于生成服从标准正态分布（均值为0，标准差为1）的随机数的函数。它返回一个具有指定形状的数组，数组中的元素是独立且服从标准正态分布的随机数。

1
2
3
4
import numpy as np

# 生成形状为 (3, 3) 的标准正态分布随机数组
random_array = np.random.randn(3, 3)

相比之下，random.random 是 Python 内置的 random 模块中用于生成0到1之间均匀分布的随机数的函数。它返回一个浮点数，该浮点数位于半开区间 [0.0, 1.0)内。

1
2
3
import random
# 生成0到1之间的均匀分布的随机数
random_number = random.random()

画图

配色

单线

颜色 线 填充
橙 #F0CC53 #FFC66B

两线

颜色
橙 #F0CC53 #F0CC53

三线

颜色
橙 #F0CC53 #00C9EF #00C9EF

模板

基本配置

1
2
3
4
5
6
7
8
# 基础设置
# 设置全局字体为Times New Roman
plt.rcParams['font.family'] = 'Times New Roman'
# 设置字号为12pt
plt.rcParams['font.size'] = 12
# 清晰
%config InlineBackend.figure_format = 'retina'
%matplotlib inline

堆叠图模板

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
import matplotlib.ticker as ticker
colors=['#F0CC53','#FF63A9']
population = population[::-1]
fig ,ax= plt.subplots(figsize=(16,4))

ax.stackplot(population['年份'],population['乡村人口(万人)'],population['城镇人口(万人)'],labels=['Rural Population','Urban Population'],colors=colors,alpha=0.6)

ax.plot(population['年份'],population['年末总人口(万人)'],'o',ls='-',color=colors[1],linewidth=3)

ax.plot(population['年份'],population['乡村人口(万人)'],'o',ls='-',color=colors[0],linewidth=3)

ax.set_ylabel('population ($1 \\times 10^4$)',weight='bold')
ax.legend(loc='upper left')
ax.set_xlabel('Year',weight='bold')
ax.grid()
ax.set_ylim(0,200000)
ax.xaxis.set_major_locator(ticker.MaxNLocator(integer=True))
ax.xaxis.set_major_locator(ticker.MultipleLocator(1))

线性图

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
# 线形图样板图
# 第0维为线颜色,第1维为填充
colors = ['#F0CC53','#FFC66B']
fig,ax=plt.subplots(figsize=(15,3))

ax.plot(population['年份'],population['年末总人口(万人)'],'o',ls='-',color=colors[0],linewidth=3,label = 'population')

ax.legend(loc='upper left')
ax.fill_between(population['年份'],0,population['年末总人口(万人)'],color=colors[1],alpha=0.3)
ax.set_xlabel('Year',weight='bold')
ax.set_ylabel('population ($1 \\times 10^4$)',weight='bold')
ax.grid()
ax.set_ylim(120000,150000)
ax.xaxis.set_major_locator(ticker.MaxNLocator(integer=True))
ax.xaxis.set_major_locator(ticker.MultipleLocator(1))

散点图

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
import matplotlib.dates as mdates
import matplotlib.ticker as ticker
# 基础设置
# 设置全局字体为Times New Roman
plt.rcParams['font.family'] = 'Times New Roman'
# 设置字号为12pt
plt.rcParams['font.size'] = 18
# 清晰
%config InlineBackend.figure_format = 'retina'
%matplotlib inline

# 画Bland-Altman图
# 设置画布大小
fig , ax = plt.subplots(figsize=(16,8))
colors=['#86AA2D','#E45C64','#B98600','#5B9BD5']
ax.plot(AQI['Date'], AQI['AQI'], 'o', color=colors[0], markersize=4,alpha=0.6,label='AQI')

# 设置y轴坐标名称为:AQI分布的英文名称
ax.set_ylabel('AQI',weight='bold')

ax.set_xlabel('Date',weight='bold')

ax.legend(loc='upper left')
# 均值
mean = 50.5368
# 可信下界
lowb= 10.5847
# 可信上界
upb= 90.4888
# 检验结果：
ans = 0.9545
# 画出三条平行于x轴的直线
ax.axhline(y=mean, color=colors[1], linestyle='-', label='Mean',linewidth=2)
ax.axhline(y=lowb, color=colors[2], linestyle='--',linewidth=2)
ax.axhline(y=upb, color=colors[3], linestyle='--',linewidth=2)
ax.xaxis.set_major_formatter(mdates.DateFormatter('%Y-%m'))

# 设置x轴的日期间隔，这里设置为每3个月
ax.xaxis.set_major_locator(mdates.MonthLocator(interval=6))

# ax.xaxis.set_major_locator(ticker.MaxNLocator(integer=True))
# ax.xaxis.set_major_locator(ticker.MultipleLocator(1))
# 把整个图的背景填充为蓝色
ax.set_facecolor('#E0E0E0')
ax.legend()

热力图

首先给出两个函数,可以辅助我们画图

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
def heatmap(data, row_labels, col_labels, ax=None,
            cbar_kw=None, cbarlabel="", **kwargs):
    """
    Create a heatmap from a numpy array and two lists of labels.

    Parameters
    ----------
    data
        A 2D numpy array of shape (M, N).
    row_labels
        A list or array of length M with the labels for the rows.
    col_labels
        A list or array of length N with the labels for the columns.
    ax
        A `matplotlib.axes.Axes` instance to which the heatmap is plotted.  If
        not provided, use current axes or create a new one.  Optional.
    cbar_kw
        A dictionary with arguments to `matplotlib.Figure.colorbar`.  Optional.
    cbarlabel
        The label for the colorbar.  Optional.
    **kwargs
        All other arguments are forwarded to `imshow`.
    """

    if ax is None:
        ax = plt.gca()

    if cbar_kw is None:
        cbar_kw = {}

    # Plot the heatmap
    im = ax.imshow(data, **kwargs)

    # Create colorbar
    cbar = ax.figure.colorbar(im, ax=ax, **cbar_kw)
    cbar.ax.set_ylabel(cbarlabel, rotation=-90, va="bottom")

    # Show all ticks and label them with the respective list entries.
    ax.set_xticks(np.arange(data.shape[1]), labels=col_labels)
    ax.set_yticks(np.arange(data.shape[0]), labels=row_labels)

    # Let the horizontal axes labeling appear on top.
    ax.tick_params(top=True, bottom=False,
                   labeltop=True, labelbottom=False)

    # Rotate the tick labels and set their alignment.
    plt.setp(ax.get_xticklabels(), rotation=-30, ha="right",
             rotation_mode="anchor")

    # Turn spines off and create white grid.
    ax.spines[:].set_visible(False)

    ax.set_xticks(np.arange(data.shape[1]+1)-.5, minor=True)
    ax.set_yticks(np.arange(data.shape[0]+1)-.5, minor=True)
    ax.grid(which="minor", color="w", linestyle='-', linewidth=3)
    ax.tick_params(which="minor", bottom=False, left=False)

    return im, cbar


def annotate_heatmap(im, data=None, valfmt="{x:.2f}",
                     textcolors=("black", "white"),
                     threshold=None, **textkw):
    """
    A function to annotate a heatmap.

    Parameters
    ----------
    im
        The AxesImage to be labeled.
    data
        Data used to annotate.  If None, the image's data is used.  Optional.
    valfmt
        The format of the annotations inside the heatmap.  This should either
        use the string format method, e.g. "$ {x:.2f}", or be a
        `matplotlib.ticker.Formatter`.  Optional.
    textcolors
        A pair of colors.  The first is used for values below a threshold,
        the second for those above.  Optional.
    threshold
        Value in data units according to which the colors from textcolors are
        applied.  If None (the default) uses the middle of the colormap as
        separation.  Optional.
    **kwargs
        All other arguments are forwarded to each call to `text` used to create
        the text labels.
    """

    if not isinstance(data, (list, np.ndarray)):
        data = im.get_array()

    # Normalize the threshold to the images color range.
    if threshold is not None:
        threshold = im.norm(threshold)
    else:
        threshold = im.norm(data.max())/2.

    # Set default alignment to center, but allow it to be
    # overwritten by textkw.
    kw = dict(horizontalalignment="center",
              verticalalignment="center")
    kw.update(textkw)

    # Get the formatter in case a string is supplied
    if isinstance(valfmt, str):
        valfmt = matplotlib.ticker.StrMethodFormatter(valfmt)

    # Loop over the data and create a `Text` for each "pixel".
    # Change the text's color depending on the data.
    texts = []
    for i in range(data.shape[0]):
        for j in range(data.shape[1]):
            kw.update(color=textcolors[int(im.norm(data[i, j]) > threshold)])
            text = im.axes.text(j, i, valfmt(data[i, j], None), **kw)
            texts.append(text)

    return texts
1
2
3
4
5
6
7
8
9
10
# 读取一个excel文件,并利用该excel文件画出热力图
df = pd.read_excel('./heatmap.xlsx',index_col=0,header=0)

fig,ax = plt.subplots()
im,cbar = heatmap(df.values,df.index,df.columns,ax=ax)
for i in range(len(df.index)):
    for j in range(len(df.columns)):
        # 注意这里的顺序为(j,i),因为第一个参数代表的是x的坐标,也即列的位置
        text = ax.text(j,i,"%.2f"%df.values[i,j],
                       ha="center",va="center",color="w")

画地图

画世界人口地图

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
from pyecharts import options as opts
from pyecharts.charts import Map
from pyecharts.faker import Faker

c = (
    Map()
    .add("Population",[list(z) for z in zip(df['Country'],df['Population'])], "world")
    .set_series_opts(label_opts=opts.LabelOpts(is_show=False))#不显示国家名字
    .set_global_opts(
        title_opts=opts.TitleOpts(title="Population"),
        visualmap_opts=opts.VisualMapOpts(max_=1500000000),
    )
    .render("../figure/map_world_population.html")
)

用folium画地图,且栅格化

第一步先生成一个图:

1
2
3
4
5
6
7
8
9
10
m = folium.Map(location = [np.mean(latitudes),np.mean(longitudes)],
               width = '100%',
               height = '100%',
               zoom_start = 5,#缩放比例
                min_zoom = 0,
                max_zoom=18,
                tiles='OpenStreetMap',
                attr = '高德地图'
              )
m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
import folium
import pandas as pd
import webbrowser as wb
import numpy as np
from folium.plugins import HeatMap
archive = pd.read_csv('dataset_of_vctria.csv') # 导入数据集 （数据集在项目的文件夹下）

#提取经纬度数据
latitude = archive['latitude'].head(30000)
longitude = archive['longitude'].head(30000)
# 绘制基准地图
m = folium.Map(location=[np.mean(latitude), np.mean(longitude)],  #表示开始时地图的起始位置
               width='100%', #窗口的宽高限制（此处指打开窗口时占据显示屏的百分比
               height='100%',
               zoom_start=5,   #初始放大倍数
               min_zoom=0,    #最小缩放倍数，由于folium地图控件有良好的交互效果，因此可以进行自由的缩放操作
               max_zoom=18,
               tiles='OpenStreetMap',  #地图类型
               attr="高德地图"  #直接用这个就行了，表示的是地图来源
               )

#在地图上标识出数据点
def pointpoint(m, latitude, lontitude):
#这里创建了一个点集一类的东西用来存储新增加的点
    point = folium.map.FeatureGroup()
    for lat, lon in zip(latitude, longitude):
        point.add_child(
            folium.CircleMarker(
                (lat, lon),
                radius=2,   # 圆圈半径
                color='red', # 半径颜色
                fill=True,     # 内部填充
                fill_color='red',
                fill_opacity=1,
            )
        )
    m.add_child(point)

#生成栅格
#选定起始点 不断生成四方格 由于经度、纬度的增量与距离的转换公式不同 因此需要进行增量的调整，这里为了简便直接取整了
#栅格大小 5km*5km
#生成400个格子
for i in range(20):
    for j in range(20):
    #folium.Polygon 函数用来生成封闭的图形，我们需要提供四个顶点的信息
    # 生成方式为向纬度减小方向生成格子，把j前面的+改为-就是向纬度增大方生成格子
        folium.Polygon(
            [
                [-38.02 + 0.04 * j, 142.91 + i * 0.05],            # 第一个参数表示四个顶点，需要画什么图形就打多少个起始点
                [-38.02 + 0.04 * j, 142.91 + (i + 1) * 0.05],
                [-38.02 + 0.04 * (j + 1), 142.91 + (i + 1) * 0.05],
                [-38.02 + 0.04 * (j + 1), 142.91 + i * 0.05]
            ], color='orange', weight=1，fill=True, fill_color='blue', fill_opacity=0.2 # 后面添加的是 线的颜色，线的粗细，填充与否，填充颜色，透明度
        ).add_to(m)                                                # add_to是将这些生成的图形加入到原图层中
m.add_child(folium.LatLngPopup())       # 这个地图组件是可以通过点击地图在地图上选点的   

m.save('../figure/map.html')

同时写了一个函数来计算某个方块应该填上什么样的颜色:

1
2
3
4
5
6
7
8
9
def color_func(i,j):
    a = [-38.02 +0.04*j,142.91 +i*0.05]
    b = [-38.02 +0.04*j,142.91 +(i+1)*0.05]
    c = [-38.02 +0.04*(j+1),142.91 +(i+1)*0.05]
    d = [-38.02 +0.04*(j+1),142.91 +i*0.05]
    if (d[1] < 143.0859)&(b[1]>143.0859)&(d[0]<-37.3547)&(b[0]>-37.7605):
        return 'orange'
    else:
        return 'blue'

使用的时候,把(-38.02,143.91)替换成与画方块的循环内对应的开始坐标,这个函数是把某一列涂成橙色,其余为蓝色.
下面封装一个函数,传入开始坐标,经度方向的大小,纬度方向的大小,经度方向的个数,纬度方向的个数,地图m

1
2
3
4
5
6
7
8
9
10
11
12
13
14
# 生成栅格
# 选定起始点 不断生成四方格
def grid_func(latitude_start,longitude_start,latitude_size,longitude_size,latitude_num,longitude_num,m):
    for i in range(longitude_num):
        for j in range(latitude_num):
            folium.Polygon(
                [
                    [latitude_start +latitude_size*j,longitude_start +i*longitude_size],
                    [latitude_start +latitude_size*j,longitude_start +(i+1)*longitude_size],
                    [latitude_start +latitude_size*(j+1),longitude_start +(i+1)*longitude_size],
                    [latitude_start +latitude_size*(j+1),longitude_start +i*longitude_size],
                ],
                color = 'orange',weight=1,fill=True,fill_color='orange',fill_opacity=0.1,
            ).add_to(m)

最优化问题

gurobi求解器

gurobi可以用来求解线性和二次型最优问题
教程

基本代码框架:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
import gurobipy
# 创建模型
MODEL = gurobipy.Model()

# 创建变量
X = MODEL.addVar(vtype=gurobipy.GRB.INTEGER,name="X")

# 更新变量环境
MODEL.update()

# 创建目标函数
MODEL.setObjective('目标函数表达式', gurobipy.GRB.MINIMIZE)

# 创建约束条件
MODEL.addConstr('约束表达式，逻辑运算')

# 执行线性规划模型
MODEL.optimize()

# 输出模型结果
print("Obj:", MODEL.objVal)
for x in X:
    print(f"{x.varName}：{round(x.X,3)}")

实战代码

结合copilot快速求解

  • 把所有变量名,取值范围,以及最优函数先列在注释中:
  • 根据copilot的提示,快速写代码

注释样例

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
# 在一家航空公司，每个月需要大量管理者进行飞机调度以安排货运订单。现在出现了一种AI可以替他们进行调度。由于AI本身的缺陷，
# 每使用1个单元的AI，就会发生er次调度错误。
# 因此，使用AI后，公司需要投入人力进行修正，称为“修正者”。如果无法及时修正，每一笔调度错误的订单要蒙受cc元的损失。
# 已知本月订单数为d ，为了使航空公司调度成本最低，请制定雇佣管理者、人工智能、修正者的数量，来进行本月的订单调度。
# 编程符号	参数大小	含义	
	 
# D_zon 121200	该公司当月的订单数
	 
# D_m 200	每个管理者最多调度 
	 
# D_a5000	每单元AI最多调度 
	 
# Err_a 50	每单元AI会有Err_a 次调度错误
	 
# Fix 24	每个修正者最多修正Fix 个订单
	   
# C_m 5000 每个管理者的花费
# C_a 1000 每个AI的花费
# C_r 8000	每个修正者花费 
	 
# F_m 10000
# F_a 50000
# F_r 10000	投入固定的启用成本
	 
# P 3000	AI税系数
# 模型符号	变量类型	  含义
# m,a,r INT	管理者、AI、修正者数量
 
# xm[i],xa[i],xr[i],0-1	是否启用r i个管理者、AI、修正者
 
# xmm,xaa,xrr,0-1	是否启用管理者、AI、修正者
   
# dm[i],da[i] INT	管理者和ai调配的订单数
# dr[i] INT	修正者修正的订单数
# err INT	错误的订单数

# 约束:
# M 代表无穷大
# gurobipy.quicksum(xm[i])=m,gurobipy.quicksum(xa[i])=a,gurobipy.quicksum(xr[i])=r
# xmm<=gurobipy.quicksum(xm[i])<=M*xmm,xaa<=gurobipy.quicksum(xa[i])<=M*xaa,xrr<=gurobipy.quicksum(xr[i])<=M*xrr
# xm[i]<=dm[i]<=M*xm[i],xa[i]<=da[i]<=M*xa[i],xr[i]<=dr[i]<=M*xr[i]
# sum(dm[i]+da[i])=D_zon
# err=a*Err_a - sum(dr[i])
# err>=0
# dm[i]<=D_m,da[i]<=D_,dr[i]<=Fix

# C_err = 10000
# 目标函数:
# min C=[m*C_m+a*C_a+r*C_r]+[xmm*F_m+xaa*F_a+xrr*F_r]+err*C_err+a*a*P

有了上面的注释,copilot基本上就能开始辅助我们求解了

下面为代码部分:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
import gurobipy
N = 10000
M = 1000
D_ = 121200
D_m = 200
D_a = 5000
Err_a = 50
Fix = 24
C_m = 5000
C_a = 1000
C_r = 8000
F_m = 10000
F_a = 50000
F_r = 10000
C_err = 10000
P = 3000

# 创建模型
MODEL = gurobipy.Model()
# 创建变量
m = MODEL.addVar(vtype=gurobipy.GRB.INTEGER, name='m')
a = MODEL.addVar(vtype=gurobipy.GRB.INTEGER, name='a')
r = MODEL.addVar(vtype=gurobipy.GRB.INTEGER, name='r')
xm = MODEL.addVars(range(1,M),vtype=gurobipy.GRB.BINARY, name='xm')
xa = MODEL.addVars(range(1,M),vtype=gurobipy.GRB.BINARY, name='xa')
xr = MODEL.addVars(range(1,M),vtype=gurobipy.GRB.BINARY, name='xr')
xmm = MODEL.addVar(vtype=gurobipy.GRB.BINARY, name='xmm')
xaa = MODEL.addVar(vtype=gurobipy.GRB.BINARY, name='xaa')
xrr = MODEL.addVar(vtype=gurobipy.GRB.BINARY, name='xrr')
dm = MODEL.addVars(range(1,M),vtype=gurobipy.GRB.INTEGER, name='dm')
da = MODEL.addVars(range(1,M),vtype=gurobipy.GRB.INTEGER, name='da')
dr = MODEL.addVars(range(1,M),vtype=gurobipy.GRB.INTEGER, name='dr')
err = MODEL.addVar(vtype=gurobipy.GRB.INTEGER, name='err')
# 更新变量环境
MODEL.update()
# 设置模型类型为二次规划
# 改为二次规划模型,作者的解决方法:
MODEL.setParam("NonConvex", 2)
MODEL.setObjective(m*C_m+a*C_a+r*C_r+xmm*F_m+xaa*F_a+xrr*F_r+err*C_err+a*a*P,sense = gurobipy.GRB.MINIMIZE)

# 创建约束条件
# m = gurobipy.quicksum(xm)
# a = gurobipy.quicksum(xa)
# r = gurobipy.quicksum(xr)
MODEL.addConstr(gurobipy.quicksum(xm) == m,name='c0')
MODEL.addConstr(gurobipy.quicksum(xa) == a,name='c1')
MODEL.addConstr(gurobipy.quicksum(xr) == r,name='c2')
MODEL.addConstr(xmm <= gurobipy.quicksum(xm),name='c3')
MODEL.addConstr(gurobipy.quicksum(xm) <= 10000*xmm,name='c4')
MODEL.addConstr(xaa <= gurobipy.quicksum(xa),name='c5')
MODEL.addConstr(gurobipy.quicksum(xa) <= 10000*xaa,name='c6')
MODEL.addConstr(xrr <= gurobipy.quicksum(xr),name='c7')
MODEL.addConstr(gurobipy.quicksum(xr) <= 10000*xrr,name='c8')
for i in range(1,M):
    MODEL.addConstr(xm[i] <= dm[i],name=f'c9_{i}')
for i in range(1,M):
    MODEL.addConstr(dm[i] <= 10000*xm[i],name=f'c10_{i}')
for i in range(1,M):
    MODEL.addConstr(xa[i] <= da[i],name=f'c11_{i}')
for i in range(1,M):
    MODEL.addConstr(da[i] <= 10000*xa[i],name=f'c12_{i}')
for i in range(1,M):
    MODEL.addConstr(xr[i] <= dr[i],name=f'c13_{i}')
for i in range(1,M):
    MODEL.addConstr(dr[i] <= 10000*xr[i],name=f'c14_{i}')

MODEL.addConstr(gurobipy.quicksum(dm)+gurobipy.quicksum(da) == D_,name='c15')
MODEL.addConstr(err == a*50 - gurobipy.quicksum(dr),name='c16')
MODEL.addConstr(err >= 0,name='c17')
for i in range(1,M):
    MODEL.addConstr(dm[i] <= D_m ,name=f'c18_{i}')
for i in range(1,M):
    MODEL.addConstr(da[i] <= D_a,name=f'c19_{i}')
for i in range(1,M):
    MODEL.addConstr(dr[i] <= Fix,name=f'c20_{i}')
# # 设置此次执行的为二次优化问题
# MODEL.setParam(gurobipy.GRB.Param.ModelSense, -1)
# 执行最优化
MODEL.optimize()
# 输出结果
print("Obj:", MODEL.objVal)
for v in MODEL.getVars():
    print(f"{v.varName}:{round(v.x,3)}")
print("m",m)
print("a",a)
print("r",r)
print("最少花费:", MODEL.objVal)

如果想求解二次型问题,需要加上这段代码:

1
MODEL.setParam("NonConvex", 2)

MODEL.setParam(“NonConvex”, 2) 这行代码用于将 Gurobi 模型的 NonConvex 参数设置为 2。这个参数控制 Gurobi 是否允许处理非凸优化问题。

具体来说，“NonConvex” 参数的取值有以下几种：

  • 0: 默认值，表示 Gurobi 将尝试解决凸优化问题。
  • 1: 表示 Gurobi 将尝试解决非凸优化问题，但它可能无法找到全局最优解。
  • 2: 表示 Gurobi 将尝试解决非凸优化问题，且不会进行全局最优性的验证，这可能导致找到局部最优解。