一种解决药店卖药问题的新思路

Hellsegamosken

已于 2022-11-01 19:54:47 修改

阅读量2.9k

点赞数 21

分类专栏：实践记录文章标签：贪心算法数据结构算法搜索

于 2022-01-05 09:28:24 首次发布

本文链接：https://blog.csdn.net/DT_Kang/article/details/122316269

版权

实践记录专栏收录该内容

14 篇文章 2 订阅

订阅专栏

期末在做数字钟，在 csdn 上参考了前人的经验。遂想到把数据结构大作业的实验报告也发上来，兴许可以供学弟学妹参考。

解决这道药店卖药问题最普遍也是最容易想到的做法是贪心求解，然而本题多限制之间相互作用，关系复杂，因此简单的贪心难以找到十分优秀的解。本文提供了一种不同的思路，从新的角度对这个问题进行思考，设计了新的算法，并且详述了算法主要思路产生过程及其合理性。实际运行结果优秀，随机数据下（也许大概率）能找到全局最优解。

报告注重思路的分析和证明，实现部分稍简略，具体可见代码。

文章目录

1 实验背景及要求

你是⼀家药店的老板，这个⽉你从供货商⼿⾥收到了⼀批共50个药品，其中每个药品有独立的进价和保质期，其中进价的区间为[20,30]元，保质期的剩余天数为[1,15]天.你每天可以陈列出10个药品在商店中售卖,每天会有三个顾客各过来购买⼀个药品。药品的定价只能为商品进价加上{-1.5,-1, -0.5, 0, 2 ,4 ,6}元，不得随意定价。每位顾客购买的逻辑是，买其中最便宜的药品，如果说最便宜的药品价格⼀致则购买保质期⻓的药品。三位顾客会依次购买药品。

药品如果没有陈列在商店中，⽽是存放在仓库时，会收取管理费，其中保质期低于5天的每天收取1元管理费，其余的每天收取0.5元。每天的陈列药品可以少于10个

你的⽬标是，10天之后，你的利润（售出商品售价总和-售出商品进价总和-⽀付仓库的管理费⽤-10天内过期/丢弃商品的进价）最⼤。

2 任务分析

2.1 初步认识

首先，肯定不能通过直接枚举的方法求解，因为如果把每一种药品摆放方式及定价策略看做一个解，那么药品选择方案数为 $\prod_{i=0}^{9}{50-3i\choose 10}$ ，定价方案数为 $6^{100}$ ，因此解空间大小为 $\prod_{i=0}^{9}{50-3i\choose 10}\cdot 6^{100}$ ，巨大无比。一种解决思路是贪心，想出一些看上去靠谱的策略，然后根据这些策略做唯一选择。这样做的好处是效率较高，代码容易实现；坏处是难以考虑所有情况，不能保证结果足够优秀。

我的思路是：首先观察问题性质，根据其性质缩小解空间，后对解空间进行搜索。

不管怎么说，要解决这个问题，就要对其有一个宏观的认识。注意到存在以下几个平凡的结论：

每天都会卖出三种药。
我们倾向于先卖即将到期的那些药，因为他们会产生额外的管理费。
有时候会刻意摆不满货架，来引导顾客选择我们想卖给他们的药。
我们倾向于先卖便宜的药，因为如果要先卖贵的药，要么就通过降低售价使之被选走，要么就要减少摆上货架的药品数量增加管理费。而先卖便宜的就没有这么多顾虑。
如果一种药最后会过期或者被扔掉，那么一定是在第一天的开始就扔掉。

基于以上观察，能想到几种贪心策略，比如最简单的把所有药品按照过期时间为第一关键字、价格为第二关键字升序排序，直接按照这个顺序选择药品，并把定价都设为最高。这样确实能够得到一组不错的解，但我们的追求远不止于此。

2.2 进一步思考

经过分析，我们可以接着得出以下性质：

性质2.2.1：最优方案一定是扔掉的药数量最少的某一种方案。

证明：扔掉一种药至少亏损20元成本价。记 $f_n$ 为扔 $n$ 种药的最大获利，我们采用反证法，假设有 $f_k<f_{k+1}$ ，那么 $f_{k+1}$ 对应的卖药方案中肯定存在一天可以卖最后被扔掉的药但是没卖。现在我们在这一天的货架上只摆出想要卖的三种药，其他位置为空，这样当天多赚了至少二十块成本，亏了至多7*1元的管理费，同时因为不扔药最多额外产生7.5元管理费。因此总盈利 $\geq f_{k+1}+20-7\times1+5\times1+5\times0.5>f_{k+1}$ ，与 $f_{k+1}$ 为最大获利矛盾。因此， $f_k\geq f_{k+1}$ 恒成立，证毕。

性质2.2.2：如果一种药 $v_i, t_i)$ ，另一种药 $v_j, t_j)$ ，满足 $v_i<v_j$ 并且 $t_i<t_j$ ，那么把 $i$ 放在 $j$ 前面肯定不会比把 $j$ 放在 $i$ 前面更劣。

证明：如果过期时间相同，肯定是先卖便宜的；如果价格相同，肯定先卖过期早的。

根据性质2.2.1，我们要尽量把将要过期的药都卖掉。首先把药按照过期时间升序排序，记 $s_{l, r}$ 表示过期时间在第 $l$ 天到第 $r$ 天之间的药品集合。注意到，如果 $s_{1, x}|< 3x$ ，最优策略肯定是将 $s_{1, x}$ 中的药品全部卖出；否则，过期时间在 $[1, x]$ 的药肯定不能都售出，并且从第 1 天到第 $x$ 天卖的一定都是 $s_{1, x}$ 中的药，此时第 1 天到第 $x$ 天卖的药和别的时间卖的药是相对独立的，分开计算并相加即可。

进一步地，如果 $|s_{1, x}|\geq 3x$ 并且已经处理完了 $[1, x]$ 天，那么接着把第 $x + 1$ 天作为第一天重复上述过程即可。

推论2.2.1：可以把排序后的药品划分为若干段，每段的选择相互独立。

过期时间随机生成的性质使得这个“分裂”的过程大概率会进行，从而缩小问题规模。

2.3 压缩解空间——优雅的枚举

再次从一个宏观的视角思考，该如何确定我们的策略？总的来说，为了缩小解空间，我们要对关键信息，即对结果影响大的变量进行决策。做一个不恰当的类比，就是要找到线性空间的一组基。一般来说有以下两种思路：

先确定每一种药品的定价，再确定每天摆什么药。
先确定每天摆什么药，再给要定价。

而对于每天摆什么药，还有两种思路：

直接确定摆哪些药，根据定价即可得知卖出了哪三种药。
先确定卖出哪三种药，再决定摆哪些药。

容易发现，先确定每一种药品定价是不怎么优秀的，因为最后总共只会卖出 30 瓶药，别的药的价格无关紧要（摆上货架不卖的药肯定越贵越好）；而且一种药在不同时间定价可以不同，难以预先确定。

而对于摆药，没卖出去的摆上货架的药其实没有那么重要。假设每天卖哪三种药及其定价已经确定了，那么只要再摆一些售价高于它们的药即可（当然优先摆五天内过期的）。而对于定价，只需要枚举卖出的三种药的定价，即可算出当天收益，取最优即可。

总结一下整体思路：首先确定每天卖哪三种药，然后逐天给这些药定价，顺便确定摆那些药能剩下最多的管理费。这样做的优势是，我们只需要枚举每天卖哪三种药和这三种药的售价，其他的量是可以直接计算出来的。实际上我们只需要枚举卖出的三种药里的最高售价即可。

再进一步，根据性质1.2.2，只需要枚举卖的是哪一天过期的药。粗略估算解空间一个很松的上界为 ${15\choose 3}^{10}$ ，远远小于原本的大小，而且其中仍有大量的违反性质2.2.1不优的解，可以被直接丢弃。最为关键的是这种问题规模的压缩是“无损”的，因为扔掉的都是那些不优的解，而最优解，或者说比较优的那些解一定还在剩下的解空间中。

事实上，如果过期时间在 [1, 15] 内均匀分布，进一步用性质2.2.1和推论2.2.1削减后的解空间大小小于 1e8 的概率约为 70%（实测），也就是说直接枚举就可以轻松找到全局最优解。可惜的是，过期时间并不是完全随机的，因此需要另寻他法。

2.4 搜索解空间——模拟退火

如上面分析，对于解空间的缩小已经做到我能力范围内的极限了，仍然不可能通过穷举解决，因此不得不采用随机算法对其进行搜索。我采用了模拟退火的方法。

具体来说，将一张3行10列且每列无序的数表看做一个解，先按照每天卖出最早过期的三种药品的方式贪心得到初始解。每次调整时，有两种方法：1. 随机某一天的某一种药，和另一天的一种药交换；2. 随机某一天的某一种药，和某一种没被选择的药交换。注意，一定满足调整后的序列满足性质2.2.1，即不能产生新的过期的药，否则会有大量无效搜索。做如上调整 1~5 次，生成新的解。

3 任务实现

3.1 模拟部分

3.1.1 关键变量

函数/变量名称	类型	功能
read()	function	解析命令行参数并读取相应文件
simulate()	function	模拟卖药过程
a	node*	药品结构体，存放价格(v)和过期时间(t)
sold	int*	是否被卖出/扔掉
del	vector	每天扔掉哪些药

3.1.2 实现流程

写在 sim.cpp 里。逐天处理，先扣除每个药的管理费，再挨个上架药品，上架的同时把扣掉的管理费加回来。排序得到卖出的三件药品，计算收益后将 sold 设为 1。最后把当天应该扔掉的药品扔掉，扣除相应进价。

3.2 制定策略

3.2.1 关键变量

函数/变量名称	类型	功能
SA()	function	模拟退火流程
calculate()	function	计算某种策略的盈利
init()	function	贪心生成模拟退火初始解
adjust()	function	从某个解随机调整生成新的解
R	Random	随机数类
state	SAnode	模拟退火中用到的策略类（每天卖哪三种）
ksj1, ksj2	vector	计算盈利时分别存放过期时间<=5和>5的药

3.2.2 实现流程

简述模拟退火流程：

生成初始解 state，设定初始温度 T
随机调整解，计算新解的盈利
如果新解更优，就接受这个解；否则按照 $P=e^{\frac{\Delta profit}{T}}$ 的概率接受这个解
降温，若达到最低温度退出，否则回到步骤2

关于退火参数的设置本人没有什么经验，完全凭感觉乱设的。但我猜测最优解距离初始贪心解不会太远，因此初始温度设置不高。此外，在最后我还对所得最优解的附近进行爬山，企图找到一组更优的解。

这个程序不需要用到上面的 sim.cpp，而是在 calculate() 中确定每天药品摆放和定价的时候顺便求出收益，因此此部分写在 s1.cpp 中。

4 实验数据

数据	data1	data2	data3	data4	data5
教辅	24	27.5	22.5	27.5	23.5
本程序	36.5	35	30.5	39.5	31

5 缺点和不足

在扔掉那些必须要扔的药的时候，我是选择尽量便宜的几种药丢弃，这样做在绝大多数情况下都是很好的。但在构造数据下，可能会出现扔掉较贵的药反而获利更大的情况。这也是我的程序中的唯一纰漏，也是其将最优解排除在解空间外的唯一情况。综合考虑随机数据下这种情况微乎其微的概率以及完善它所需要的代码量，没有选择去修补它。

6 总结与感悟

从结果看，运行多次程序结果相同，这令我非常吃惊，因为我的程序里有随机化成分，在我预想中最终结果应该有微小浮动，至少不会总是一样。对此有两种可能的解释：一种解释是这个解就是全局最优解，自然不会有更好的解出现；另一种解释是这是局部最优解而全局最优解极为刁钻，甚至不能近其身。我个人倾向于第一种。因为每多一种药过期在我的解法里相当于多了一条限制，会帮助我缩小解空间。而一开始下发的data1中没有药过期，理论上这种情况下我的算法找到最优解的概率最小，但它仅用很短的时间就找到了一个解，并且通过手算可以证明它就是最优解。

当然，这也无法完全验证或者证伪，因为目前还没有别人跑出更优的解。究其原因，可能是2.3中给出的解空间上界确实非常非常松，实际数据下解空间大小远没有这么大，从而比较容易找到最优解；也有可能是最优解距离初始贪心解的距离不太远。

回过头来看，整个算法最为关键的部分是前面对于压缩解空间的分析，后面的模拟退火反而没那么重要，只是顺势而为罢了，或者只能算一种妥协。事实上，我把模拟退火的部分去掉后，直接从初始贪心解周围随机若干解的效果也不错，甚至也超过了教辅的收益。将随机换成直接爬山效果更佳，距离模拟退火退出来的结果只有几元钱的差距。

最终的 sim.cpp 和 s1.cpp 加上实验报告总共花了三个晚上的时间，比预期短，是因为我吸取了银行的教训，没有过于自信直接上手而花费很多时间在不可行的方案上。这次我在前期花了大量的时间准备思考，写了大量测试数据的代码（虽然这部分最后没派上用场），全部考虑清楚后再着手写代码，避免了很多无用功，也让我的思路更加清楚，代码一气呵成，没花多少时间调试。

/**
 * TODO:
 * 
 * 
 * */

# include <bits/stdc++.h>
# include <getopt.h>
# define pb push_back
# define ll long long
using namespace std;
const double deltav[] = {-1.5, -1, -0.5, 0, 2, 4, 6};  // 定价表
struct node {
    double v;
    int t;
} a[100];

string fpath = "testdata\\data";  // 路径，测试用

// 策略结构体，包含strategy和delete
struct Koutput {
    int s[11][11][2], d[51][2], n;
    Koutput() { n = 0; }
    Koutput& operator=(const Koutput &b) {
        memcpy(s, b.s, sizeof(s));
        memcpy(d, b.d, sizeof(d));
        n = b.n;
        return *this;
    }

    void print() {
        ofstream fout;
        fout.open((fpath + "\\strategy.my").c_str());
        for (int i = 1; i <= 10; i++) {
            for (int j = 0; j < 10; j++) {
                fout << s[i][j][0] << ',' << s[i][j][1] << ' ';
            }
            fout << '\n';
        }
        fout.close();
        fout.open((fpath + "\\delete.my").c_str());
        for (int i = 0; i < n; i++) {
            fout << d[i][0] << ' ' << d[i][1] <<'\n';
        }
        fout.close();
    }
} beststr;

// 退火用到的3*10的数表结构体
struct SAnode{
    int A[11][3], B[30], n;
    SAnode(): n(0) {}
} initstate, bestSAstr;

int n = 50;  // 药品总数
int s_a[100];  // 过期时间升的药品编号
int initDel[100], initDelcnt;  // 一开始扔掉哪些药及其数量
double bestans = -10000;  // 最大获利
double initCost = 0;  // 扔药的损失
vector<int> ksj1, ksj2;  // 计算获利时维护的集合，ksj1是即将过期的药，ksj2是其余药

bool cmp_greater(int x, int y) { return a[x].v > a[y].v; }
bool cmp_id(int x, int y) {
    if (a[x].t ^ a[y].t) return a[x].t < a[y].t;
    return a[x].v > a[y].v;
}

// 生成[0,1]浮点数
double frand() { return 1.0 * rand() / RAND_MAX; }

// 解析命令行参数并读入
void read(int argc, char* argv[]) {
    if (argc == 1) {
        int tmp;
        scanf("%d", &tmp);
        fpath = fpath + to_string(tmp);
        FILE* fp = fopen((fpath + "\\medicine.txt").c_str(), "r");
        for (int i = 0; i < n; i++) {
            fscanf(fp, "%lf%d", &a[i].v, &a[i].t);
            s_a[i] = i;
        }
        fclose(fp);
        return;
    }

    int opt;
    const char* optstring = "m:s:d:";

    while ((opt = getopt(argc, argv, optstring)) != -1) {
        FILE* fp = fopen(optarg, "r");
        if (opt == 'm') {
            for (int i = 0; i < n; i++) {
                fscanf(fp, "%lf%d", &a[i].v, &a[i].t);
                s_a[i] = i;
            }
        }
        fclose(fp);
    }
}

// 计算退火解的盈利
double calculate(SAnode s) {
    Koutput str;
    double profit = 0;
    ksj1.clear(), ksj2.clear();
    for (int i = 0; i < s.n; i++) {
        int u = s.B[i];
        if (a[u].t <= 10) str.d[str.n][0] = 0, str.d[str.n][1] = u, str.n++, profit -= a[u].v;
        else ksj2.pb(u);
    }
    // 初始化
    for (int i = 1; i <= 10; i++) {
        for (int j = 0; j < 3; j++) {
            if (a[s.A[i][j]].t <= 5) ksj1.pb(s.A[i][j]);  //=
            else ksj2.pb(s.A[i][j]);
        }
    }
    vector<int>::iterator it;
    // 开始遍历每一天
    for (int i = 1; i <= 10; i++) {
        sort(s.A[i], s.A[i] + 3, cmp_greater);
        double maxv = a[s.A[i][0]].v;

        // 把临近过期的药转移到 ksj1
        for (it = ksj2.begin(); it != ksj2.end(); ) {
            int u = *it;
            if (a[u].t - i + 1 <= 5) ksj1.pb(u), it = ksj2.erase(it);  //=
            else it++;
        }

        // 上架卖的三种药
        for (int j = 0; j < 3; j++) {
            int u = s.A[i][j];
            str.s[i][j][0] = u;
            if (a[u].t - i + 1 <= 5) ksj1.erase(find(ksj1.begin(), ksj1.end(), u));  //=
            else ksj2.erase(find(ksj2.begin(), ksj2.end(), u));
        }
        sort(ksj1.begin(), ksj1.end(), cmp_greater);
        sort(ksj2.begin(), ksj2.end(), cmp_greater);

        profit -= ksj1.size() + ksj2.size() * 0.5;  // 先减去所有药的管理费，这样可以把上架药品的管理费看做收益

        int p = 0, q = 0;
        int k0 = 6, k1 = 6, k2 = 6, u = s.A[i][0], v = s.A[i][1], w = s.A[i][2];
        int maxk0, maxk1, maxk2, maxp, maxq;
        double Max = -10000;

        // 枚举卖出药品的定价不大于 maxv+k，计算此时当天收益（卖药利润+上架不卖的7件药品的管理费）
        for (int k = 6; k >= -1; k--) {
            while (p < ksj1.size() && a[ksj1[p]].v + 6 > maxv + k) p++;
            while (q < ksj2.size() && a[ksj2[q]].v + 6 > maxv + k) q++;
            while (a[u].v + deltav[k0] > maxv + k) k0--;
            while (a[v].v + deltav[k1] > maxv + k) k1--;
            while (a[w].v + deltav[k2] > maxv + k) k2--;

            double tmp = deltav[k0] + deltav[k1] + deltav[k2] + min(p, 7) + 0.5 * min(q, 7 - min(p, 7));
            if (tmp > Max) {
                Max = tmp;
                maxp = p, maxq = q;
                maxk0 = k0, maxk1 = k1, maxk2 = k2;
            }
        }
        int t = 0;
        str.s[i][0][1] = maxk0;
        str.s[i][1][1] = maxk1;
        str.s[i][2][1] = maxk2;

        while (t < 7 && maxp) str.s[i][3 + t][0] = ksj1[--maxp], str.s[i][3 + t][1] = 6, t++;
        while (t < 7 && maxq) str.s[i][3 + t][0] = ksj2[--maxq], str.s[i][3 + t][1] = 6, t++;
        while (t < 7) str.s[i][3 + t][0] = -1, str.s[i][3 + t][1] = 6, t++;
        profit += Max;
    }
    
    if (profit > bestans) {
        bestans = profit;
        beststr = str;
        bestSAstr = s;
    }

    return profit;
}

// 调整退火解
SAnode adjust(SAnode s) {
    for (int i = 0; i <= rand() % 5; i++) {
        while (1) {
            int opt = rand() % 4;
            if (opt == 0) {
                int x = rand() % 10 + 1, y = rand() % 3, z = rand() % s.n;
                if (a[s.B[z]].t >= x && a[s.A[x][y]].t > 10) {
                    swap(s.B[z], s.A[x][y]);
                    break;
                }
            } else {
                int x = rand() % 10 + 1, y = rand() % 3;
                int u = rand() % 10 + 1, v = rand() % 3;
                if (a[s.A[x][y]].t >= u && a[s.A[u][v]].t >= x) {
                    swap(s.A[x][y], s.A[u][v]);
                    break;
                }
            }
        }
    }
    return s;
}

// 贪心生成初始解
SAnode greedy() {
    int p = 0;
    SAnode ret;
    for (int i = 1; i <= 10; i++) {
        for (int j = 0; j < 3; j++) {
            while (a[s_a[p]].t < i) {
                int u = s_a[p];
                for (int ii = 1; ii <= a[u].t; ii++) {
                    for (int jj = 0; jj < 3; jj++) {
                        if (a[ret.A[ii][jj]].v < a[u].v) swap(u, ret.A[ii][jj]);
                    }
                }
                initDel[initDelcnt++] = u;
                initCost += a[u].v;
                p++;
            }
            ret.A[i][j] = s_a[p++];
        }
    }
    while (p < 50) {
        ret.B[ret.n++] = s_a[p++];
    }
    return ret;
}

// 进行模拟退火
double SA() {
    double T = 10000, alpha = 0.99, lst = -10000;
    SAnode state = initstate;
    lst = calculate(state);
    while (T > 0.000001) {
        SAnode newState = adjust(state);
        double tmp = calculate(newState), delta = tmp - lst;
        if (delta > 0 || exp(delta/T) > frand()) lst = tmp, state = newState;
        T *= alpha;
    }
    return lst;
}

int main(int argc, char* argv[]) {
    srand(time(0) + 9927);
    //cout<<"#";
    read(argc, argv);

    sort(s_a, s_a + n, cmp_id);
    initstate = greedy();
    for (int i = 1; i <= 300; i++) SA();

    // 在所得解附近爬山
    SAnode state = bestSAstr;
    double lst = bestans;
    for (int i = 1; i <= 5000; i++) {
        SAnode newState = adjust(state);
        double tmp = calculate(newState);
        if (tmp > lst) {
            lst = tmp;
            state = newState;
        }
    }
    
    for (int i = 0; i < initDelcnt; i++) beststr.d[beststr.n][0] = 0, beststr.d[beststr.n][1] = initDel[i], beststr.n++; 
    beststr.print();
    cout << bestans - initCost;
    return 0;
}
/* by Hellsegamosken */

Hellsegamosken

关注

21
点赞
踩
33

收藏

觉得还不错? 一键收藏
7
评论
一种解决药店卖药问题的新思路

解决这道药店卖药问题最普遍也是最容易想到的做法是贪心求解，然而本题多限制之间相互作用，关系复杂，因此简单的贪心难以找到十分优秀的解。本文提供了一种不同的思路，从新的角度对这个问题进行思考，设计了新的算法，并且详述了算法主要思路产生过程及其合理性。实际运行结果优秀，随机数据下（也许大概率）能找到全局最优解。报告注重思路的分析和证明，实现部分稍简略，具体可见代码。
复制链接

扫一扫