2024年高教社杯全国大学生数学建模C题-农作物的种植策略详解+思路+Python源码(一)

前言

简单介绍一下我自己：博主专注建模五年，参与过大大小小数十来次数学建模，理解各类模型原理以及每种模型的建模流程和各类题目分析方法。现提供免费的思路和部分源码，以后的数模比赛只要我还有时间肯定会第一时间写出免费开源思路。博主紧跟各类数模比赛，每场数模竞赛博主都会将最新的思路和代码写进此专栏以及详细思路和完全代码且完全免费。希望有需求的小伙伴不要错过笔者精心打造的文章。

数学建模的基本步骤大致如下：

1. 问题定义：明确你需要解决的问题。这包括对问题的背景进行研究，理解问题的实际意义，以及确定问题的边界条件和假设。
2. 建立数学模型：将实际问题转换为数学问题。这通常涉及到定义变量、参数、约束条件以及目标函数（如果是优化问题的话）。
3. 求解数学模型：选择合适的数学工具和方法求解模型。这可能包括解析方法、数值方法、仿真等。
4. 模型验证与分析：使用实际数据测试模型的有效性，分析模型结果，对模型进行必要的调整和优化。
5. 撰写模型报告：详细记录模型的建立过程、求解步骤、结果分析以及模型的局限性等。

我们依旧按如上步骤完成对C题的解析作答。

C题背景

某乡村地处华北山区，常年温度偏低，大多数耕地每年只能种植一季农作物。该乡村现有露天耕地 1201 亩，分散为 34 个大小不同的地块，包括平旱地、梯田、山坡地和水浇地 4 种类型。平旱地、梯田和山坡地适宜每年种植一季粮食类作物；水浇地适宜每年种植一季水稻或两季蔬菜。该乡村另有 16 个普通大棚和 4 个智慧大棚，每个大棚耕地面积为 0.6 亩。普通大棚适宜每年种植一季蔬菜和一季食用菌，智慧大棚适宜每年种植两季蔬菜。同一地块（含大棚）每季可以合种不同的作物。

详见附件 1。

根据农作物的生长规律，每种作物在同一地块（含大棚）都不能连续重茬种植，否则会减产；因含有豆类作物根菌的土壤有利于其他作物生长，从 2023 年开始要求每个地块（含大棚）的所有土地三年内至少种植一次豆类作物。同时，种植方案应考虑到方便耕种作业和田间管理，譬如：每种作物每季的种植地不能太分散，每种作物在单个地块（含大棚）种植的面积不宜太小，等等。2023年的农作物种植和相关统计数据见附件 2。

请建立数学模型，研究下列问题：

问题 1

假定各种农作物未来的预期销售量、种植成本、亩产量和销售价格相对于 2023 年保持稳定，每季种植的农作物在当季销售。如果某种作物每季的总产量超过相应的预期销售量，超过部分不能正常销售。请针对以下两种情况，分别给出该乡村 2024~2030 年农作物的最优种植方案，将结果分别填入 result1_1.xlsx 和 result1_2.xlsx 中（模板文件见附件 3）。

(1) 超过部分滞销，造成浪费；

(2) 超过部分按 2023 年销售价格的 50%降价出售。

第一步：问题定义

本问题针对华北山区低温环境下的农作物种植规划与布局。在该地区，常规露天耕地每年只能种一季粮食类作物，而水浇地可以种水稻或两季蔬菜；大棚（普通与智慧大棚）则可多季种植特定作物。如何在有限的耕地和设施条件下，合理安排不同作物的轮作、搭配以及分布，从而达到最大化经济收益、满足作物轮作要求、保证土壤肥力与可持续发展成为一个重要的问题。

梳理条件

这是一道典型的农业经济数学建模问题，涉及多目标优化、种植规划和约束条件设计。首先，让我们梳理问题的关键信息和约束条件：

耕地资源：

1. 露天耕地：1201亩，分34个地块
   • 包括平旱地、梯田、山坡地、水浇地
2. 大棚：
   • 16个普通大棚（每个0.6亩）
   • 4个智慧大棚（每个0.6亩）

种植约束：

1. 轮作限制
   • 同一地块不能连续重茬
   • 三年内每块地必须种植一次豆类作物.
2. 种植多样性要求
   • 每种作物每季种植地不能太分散
   • 单个地块种植面积不宜太小

数据解读与整理

作物信息整理

附件1中列出作物编号、名称、类型及适宜种植耕地类型、季节限制，以及若干说明性约束。下面对关键信息进行提炼：

作物分类及种植季节限制：

• 粮食类作物（含豆类）： 可在平旱地、梯田和山坡地每年种一季（春夏季或秋季，视题意），水浇地（单季为水稻，多季为蔬菜，不适宜粮食类除豆类外的多季种植），大棚（不涉及粮食类种植，除豆类蔬菜可能在智慧大棚中种植）。
  • 豆类粮食作物（黄豆、黑豆、红豆、绿豆、爬豆）编号：1-5
  • 其他粮食作物（小麦、玉米、谷子、高粱、黍子、荞麦、南瓜、红薯、莜麦、大麦）编号：6-15
  • 水稻（编号16）仅适合水浇地单季。

蔬菜类作物： 分为豆类蔬菜（豇豆、刀豆、芸豆）和其他蔬菜类（土豆、西红柿、茄子、菠菜、青椒、菜花、包菜、油麦菜、小青菜、黄瓜、生菜、辣椒、空心菜、黄心菜、芹菜、大白菜、白萝卜、红萝卜）。

• 豇豆、刀豆、芸豆（17-19号）为蔬菜（豆类）。
  • 豇豆（17）：可在水浇地第一季、普通大棚第一季、智慧大棚两季种植。
  • 刀豆（18）、芸豆（19）未在说明中列出特定季节限制（默认与一般蔬菜相同，可能在水浇地第一季、大棚第一季、智慧大棚两季可行，但需参考通用规则）。
• 其他蔬菜（20-34号）：
  • 水浇地可种第一季多种蔬菜（除大白菜、白萝卜、红萝卜外），第二季只能种大白菜（35）、白萝卜（36）、红萝卜（37）中的一种。
  • 普通大棚第一季可种多种蔬菜（不含大白菜、白萝卜、红萝卜），第二季只能种食用菌。
  • 智慧大棚两季皆可种蔬菜，但同样不包括大白菜、白萝卜、红萝卜。

食用菌类作物（38-41号：榆黄菇、香菇、白灵菇、羊肚菌）： 只能在普通大棚的第二季种植。

种植季节定义（参考注释）：

• 水浇地：
  • 第一季：3月~6月前后
  • 第二季：7月~10月前后
• 普通大棚：
  • 第一季：5月~9月前后
  • 第二季：9月~次年4月前后（适合食用菌）
• 智慧大棚：
  • 第一季：3月~7月前后
  • 第二季：8月~次年2月前后

特殊要求：

• 若在同一水浇地种植两季作物，第一季可以种植多种蔬菜（不含大白菜、白萝卜、红萝卜），第二季只能在大白菜、白萝卜、红萝卜中选一种。
• 大白菜、白萝卜、红萝卜只能在水浇地的第二季种植，不适用其他地块类型或季节。
• 普通大棚第二季只能种食用菌（第一季种蔬菜，但不包括大白菜、白萝卜、红萝卜）。
• 智慧大棚每年两季都可以种蔬菜（不包含大白菜、白萝卜、红萝卜）。
• 平旱地、梯田、山坡地每年一季仅适合种粮食类（豆类或非豆类），不可种水稻、水浇地专属作物或蔬菜。
• 水稻只适合水浇地种植一季（替代蔬菜种植方案的选择），不适合其他地块类型。

地块信息整理

(1) 平旱地、梯田、山坡地：

• 每年只能种一季粮食类作物（包括豆类粮食，但不含水稻）。

(2) 水浇地：

• 每年可以选择只种一季水稻或选择种两季蔬菜。
• 若种两季蔬菜：第一季可种各类蔬菜（除大白菜、白萝卜、红萝卜），第二季只能从大白菜、白萝卜、红萝卜中选一种。

(3) 普通大棚：

• 每年两季：第一季可种多种蔬菜（不含大白菜、白萝卜、红萝卜），第二季只能种食用菌类。

(4) 智慧大棚：

• 每年两季都可种蔬菜（不含大白菜、白萝卜、红萝卜）。

(5) 豆类轮作要求（来自上一轮描述，不在附件1中直接重复，但已知题意中提及三年内至少有一季豆类作物）

(6) 不连续重茬要求（同一地块不能连续种同一种作物）

(7) 同一季同一地块可以混合种不同作物（合种），但需满足最低面积分配要求等（题中可能会有后续约束）。

地块与作物关联

根据附件1与附件2，我们有以下信息：

地块类型与适宜作物（见附件1解析）：

• 平旱地、梯田、山坡地：一年一季，仅种粮食类作物（豆类粮食可计入）。不得种植水稻、蔬菜和菌类。
• 水浇地：可选择一年一季水稻，或两季蔬菜（第一季各类蔬菜不含大白菜(35)、白萝卜(36)、红萝卜(37)，第二季必须从大白菜、白萝卜、红萝卜中选一种）。
• 普通大棚：每年两季。第一季可种蔬菜（不包括35-37号），第二季只能种食用菌（38-41）。
• 智慧大棚：每年两季可种蔬菜（不包括35-37号）。无食用菌或水稻或粮食种植。

产量、成本、价格（见附件2）：

• 附件2给出每种作物在不同地块类型和种植季次下的亩产量（斤/亩）、单位种植成本（元/亩）、销售价格区间（元/斤）。
• 同一种作物在不同地块类型、不同季节、不同大棚类型下产量和成本可能不同。
• 由于问题假设2023年为基准年份，未来价格、产量和成本保持稳定，可从附件2的数据中直接取2023年的数据作为未来参数：
  • 若价格为区间，可取区间中值或参考价作为确定值。例如，若价格为2.50-4.00元/斤，可取中值3.25元/斤作为稳定价格。
  • 本模型中记为：

$$

Y{j,i,s},C{j,i,s},P_{j}

$$

分别表示某作物j在地块类型与季节条件下的亩产量、成本和价格。

实际实现中，需先根据地块类型和季节，将附件2中的数据匹配到对应的i和s，然后再求解。

数学模型描述(问题1)

集合与索引

• 年份集合： $$ T={2024,2025,…,2030},t∈TT = {2024,2025,\dots,2030}, \quad t \in T $$
• 季节集合（根据题意，一年内有 1 或 2 个种植季次，具体以题中定义为准。假设有两个季次，则： $$ S={{1,2}} \quad s \in S $$ 若部分地块类型一年只有一季，则在相应的决策中会有对应的约束令第二季决策变量为0。
• 地块集合（含露天地和大棚）： $$ I={{1,2,…,N}}, \quad i \in I $$ 不同地块有不同类型，如平旱地、梯田、山坡地、水浇地、普通大棚和智慧大棚。
• 作物集合： $$ J={{1,2,…,M}}, \quad j \in J $$ 作物包括粮食类（豆类粮食和其他粮食）、水稻、蔬菜类（豆类蔬菜和其他蔬菜）、食用菌类。

参数

• $A_{i}$：地块 $i$ 的面积（亩）。
• $Y_{j,i,s}$：在地块类型和季节条件下，种植作物 $j$ 在地块 $i$ 上的单位亩产量（斤/亩）。该值可根据附件2中地块类型与季节的适用数据确定。
• $ C_{j,i,s}$：在地块类型和季节条件下，作物 $j$ 的单位种植成本（元/亩）。
• $P_j$：作物 $j$ 的销售价格（元/斤）。从价格区间中取平均值或中位值作为稳定价格。
• $D_{j,s,t}$：第 $t$ 年第 ss 季作物 jj 的预期销售量上限（斤）。超出此量时，在情景(1)中无收益，情景(2)中以半价出售。

决策变量

• $x_{i,j,s,t}≥0$：在第 $t$ 年第 $s$ 季于地块 $i$ 上种植作物 $j$ 的面积（亩）。
• 辅助变量（用于处理超过或低于预期销售量的部分）：

$Uj,s,t≥0$:该季该年作物 j 的产量超过预期销售量部分（斤）$U{j,s,t} \geq 0$: $\text{该季该年作物 } j \text{ 的产量超过预期销售量部分（斤）}$

$Wj,s,t≥0$:该季该年作物 j 的产量低于预期销售量部分（斤）$W{j,s,t} \geq 0$: $\text{该季该年作物 } j \text{ 的产量低于预期销售量部分（斤）}$

在最优解中，$U{j,s,t} $和 $W{j,s,t} $不会同时为正值。

约束条件

1. 土地利用约束： 每块地在每季的种植面积总和不超过地块面积： $$ \sum{j} x{i,j,s,t} \leq A_i, \quad \forall i,s,t $$ .
2. 适宜性与季节性约束： 根据地块类型和作物季节要求（附件1所述），对不允许的作物-地块类型-季节组合施加限制：

$xi,j,s,t=0$若组合$ (i,j,s) $不可行$x_{i,j,s,t} = 0$ $\quad \text{若组合 }(i,j,s)\text{ 不可行}$

例如：平旱地、梯田、山坡地只种一季粮食类作物；水浇地两季蔬菜时第一季不种 35-37号作物，第二季只能种 35-37号之一；普通大棚第二季只能种食用菌；智慧大棚两季皆可种蔬菜（除35-37号）等。

1. 不连续重茬约束： 同一地块不可在连续季节（包括跨年）种同一种作物： $$ y{i,j,s,t} + y{i,j,s',t'} \leq 1 $$ 其中 $(s',t')$ 是$(s,t)$的上一季期。 利用二元变量 $y{i,j,s,t} $表示该地块该季是否种植作物 $j$，并通过逻辑关系 $(x{i,j,s,t}>0⇒y{i,j,s,t}=1)(x{i,j,s,t}>0 )$ 实现。
2. 三年豆类轮作约束： 每个地块在任意连续三年内至少有一季种豆类作物：

产量定义与超产/不足部分：

总产量：

超出预期量：

未达预期量：

1. 通过这两个约束，可线性刻画 $min⁡(Q{j,s,t},D{j,s,t})=D{j,s,t}−U_j$和 $max⁡(Q{j,s,t}−D{j,s,t},0)=U{j,s,t}$。

目标函数

问题1要求求解两种情景下的最优种植方案。

记：

情形(1)：超过部分滞销无收益

收益 = 正常价出售 $\min(Q{j,s,t},D{j,s,t})$ 的产量，超过部分无收益。

因此，

目标函数：

情形(2)：超过部分按半价出售

若超产量为$U_{j,s,t}$，则该部分以$0.5P_j $出售。

同时，当 $Q{j,s,t}<D{j,s,t}$ 时，有不足部分 $W{j,s,t}$，正常价出售的数量为 $D{j,s,t}-W_{j,s,t}$。

综合可得：

目标函数：

Python建模

这里将给出建模的基础代码结构和思路。后续可根据附件数据整理出相应的参数字典或数据表，再读入模型。

# 设定跨季约束需要知道上一季对应关系，如一年两季则上一季可能为(s-1,t)或跨年
# 对于季的衔接，这里假设S=[1,2], 则(s=1,t=2025)上一季为(s=2,t=2024)等
def prev_season(s, t):
    # 返回上一季的(s', t')
    if s == 1:
        return (2, t-1)
    else:
        return (1, t)
# -----------------------------------------
# 建立模型
model = pyo.ConcreteModel("Farming_Optimization")

# 定义变量
model.x = pyo.Var(I, J, S, T, domain=pyo.NonNegativeReals) # x_{i,j,s,t}
model.U = pyo.Var(J, S, T, domain=pyo.NonNegativeReals)    # U_{j,s,t}
model.W = pyo.Var(J, S, T, domain=pyo.NonNegativeReals)    # W_{j,s,t}

# 二元变量 y_{i,j,s,t} 用于表示是否种植，用大M法或逻辑约束
# 为简化，这里直接定义y为二元变量
model.y = pyo.Var(I, J, S, T, domain=pyo.Binary)

# -----------------------------------------
# 约束条件

# 1. 土地利用约束
def land_use_rule(model, i, s, t):
    return sum(model.x[i,j,s,t] for j in J) <= A[i]
model.LandUse = pyo.Constraint(I, S, T, rule=land_use_rule)

# 2. 适宜性与季节性约束
def feasibility_rule(model, i, j, s, t):
    if no_plant[(i,j,s)]:
        return model.x[i,j,s,t] == 0
    return pyo.Constraint.Skip
model.Feasibility = pyo.Constraint(I, J, S, T, rule=feasibility_rule)

# 3. 不连续重茬约束
# 需要确保若上一季同地块同作物种了，则本季不种同作物
def no_repeated_rule(model, i, j, s, t):
    s_prev, t_prev = prev_season(s, t)
    if t_prev < min(T):
        return pyo.Constraint.Skip
    return model.y[i,j,s,t] + model.y[i,j,s_prev,t_prev] <= 1
model.NoRepeated = pyo.Constraint(I, J, S, T, rule=no_repeated_rule)

# x与y的关联：当x>0时 y=1，否则y=0
# 可以用如下方式： x[i,j,s,t] <= A[i]*y[i,j,s,t]
def x_y_link_rule(model, i, j, s, t):
    return model.x[i,j,s,t] <= A[i]*model.y[i,j,s,t]
model.XYLink = pyo.Constraint(I, J, S, T, rule=x_y_link_rule)

# 4. 三年豆类轮作约束
# 对任意 i,t，连续三年内至少一季种豆类
# 以3年为一个窗口 t, t+1, t+2
def bean_rotation_rule(model, i, start_t):
    # start_t从2024到2028
    # 三年窗口 = start_t, start_t+1, start_t+2
    # bean类作物某一季y=1则满足
    years = [start_t, start_t+1, start_t+2]
    return sum(model.y[i,j,s,t_] for t_ in years for s in S for j in bean_crops) >= 1
model.BeanRotation = pyo.Constraint(I, [t for t in T if t <= 2028], rule=bean_rotation_rule)

# 5. 产量与超/不足定义
# Q_{j,s,t} = sum_i Y_{j,i,s} x_{i,j,s,t}
def Q_expr(model, j, s, t):
    return sum(Y[(j,i,s)]*model.x[i,j,s,t] for i in I)
model.Q = pyo.Expression(J, S, T, rule=Q_expr)

# 超出部分: U_{j,s,t} >= Q_{j,s,t}-D_{j,s,t}
def U_rule(model, j, s, t):
    return model.U[j,s,t] >= model.Q[j,s,t] - D[(j,s,t)]
model.UConstraint = pyo.Constraint(J, S, T, rule=U_rule)

# 不足部分: W_{j,s,t} >= D_{j,s,t}-Q_{j,s,t}
def W_rule(model, j, s, t):
    return model.W[j,s,t] >= D[(j,s,t)] - model.Q[j,s,t]
model.WConstraint = pyo.Constraint(J, S, T, rule=W_rule)


# -----------------------------------------
# 目标函数
# 对于问题1，我们有两种情景：
# 情形(1): 超出部分无收益
#   Revenue = sum_{t,s,j} P_j*(D_{j,s,t}-U_{j,s,t})
# 情形(2): 超出部分半价出售
#   Revenue = sum_{t,s,j} [P_j(D_{j,s,t}-W_{j,s,t}) + 0.5P_j U_{j,s,t}]

# 我们可以先定义一个参数 scenario 表示选择场景
scenario = 1  # 改为2即可切换情形(2)

def obj_rule(model):
    if scenario == 1:
        # 情形(1)
        revenue = sum(P[j]*(D[(j,s,t)] - model.U[j,s,t]) for j in J for s in S for t in T)
    else:
        # 情形(2)
        revenue = sum(P[j]*(D[(j,s,t)] - model.W[j,s,t]) + 0.5*P[j]*model.U[j,s,t] for j in J for s in S for t in T)
    cost = sum(C[(j,i,s)]*model.x[i,j,s,t] for i in I for j in J for s in S for t in T)
    return revenue - cost

model.OBJ = pyo.Objective(rule=obj_rule, sense=pyo.maximize)

# -----------------------------------------
# 求解模型
solver = pyo.SolverFactory('gurobi')  # 或 cplex/glpk
results = solver.solve(model, tee=True)

# 输出结果到文件
# 假设需要将 x_{i,j,s,t} 的结果输出到 result1_1.xlsx 或 result1_2.xlsx
# 使用pandas写入excel是常规手段，此处仅示意：
import pandas as pd

solution_x = []
for i in I:
    for j in J:
        for s in S:
            for t in T:
                val = pyo.value(model.x[i,j,s,t])
                if val > 1e-6: # 若有正种植面积
                    solution_x.append([i,j,s,t,val])

目前此题先到此，后续再整合完成其他问题解答，内容有点过多需要深刻理解。