1.4 二分法、实数精度
在有些情况下,问题的所有数据对象为一个有序区间。二分法将这个区间等分成两个子区间,根据计算的要求决定下一步计算是在左子区间进行还是在右子区间进行;然后根据计算的要求等分所在区间,直至找到解为止。显然,对一个规模为O(n)的问题,如果采用盲目枚举的办法,则效率为O(n);若采用二分法,则计算效率可提高至O(log2(n))。
许多算法都采用了二分法,例如二分法查找、减半递推技术、快速排序、合并排序、最优二叉树、线段树等。其中比较浅显的算法是二分法查找和减半递推技术,使用这两种方法解简单计算题,可以显著提高计算时效。
假设数据是按升序排序的,二分法查找的基本思想是对于待查找值x,从序列的中间位置开始比较:
1)若当前中间位置值等于x,则查找成功。
2)若x小于当前中间位置值,则在数列的左子区间(数列的前半段)中查找。
3)若x大于当前中间位置值,则在右子区间(数列的后半段)中继续查找。
以此类推,直至找到x在序列中的位置(查找成功)或子区间不存在(查找失败)为止。若查找失败,则当前子区间右指针所指的元素是序列中大于x的最小数。
【1.4.1 Pie】
我的生日快到了,按传统,我用馅饼招待朋友。我用的不是一块馅饼,而是有很多块馅饼,口味、大小各异。有F位朋友要来参加我的生日聚会,每个朋友会得到一块馅饼,而不是几小块馅饼,因为这样的话,看起来很乱。也就是说,这一块馅饼是一整块的馅饼。
我的朋友们都很烦人,如果他们中的一个得到了比其他人更大的一块馅饼,他们就会开始抱怨。因此,他们所有人都应该得到同样大小(但不一定是同样形状)的馅饼,即使这会导致一些馅饼被浪费,但这比破坏了聚会要好。当然,我也要留一块馅饼给自己,而且那块馅饼也应该是同样大小的。
我们能得到的最大尺寸的馅饼是多少?所有的馅饼都是圆柱形的,它们都有相同的高度1,但是馅饼的半径可以是不同的。
输入
输入的第一行给出一个正整数:测试用例的数量。然后,对于每个测试用例:
1)在一行中给出两个整数N和F(1≤N,F≤10 000),分别为馅饼的数目和朋友的数目。
2)在一行中给出N个整数为ri(1≤ri≤10 000),即馅饼的半径。
输出
对于每个测试用例,在一行中输出最大可能体积的V,使我和我的朋友们都可以得到一个大小为V的馅饼。答案是一个浮点数,误差绝对值不超过10-3。
试题来源:ACM Northwestern Europe 2006
在线测试:POJ 3122,UVA 3635
试题解析
每个朋友都会得到一整块馅饼,也就是说,每个朋友得到的那块馅饼必须是从同一个馅饼上得到的,而馅饼被分后,剩下的部分就会被浪费。
有F位朋友和我得到一整块馅饼,所以,馅饼要分成F+1块。采用二分法求解馅饼的最大尺寸。
初始时,区间下界low=0,上界high=maxsize(所有馅饼里最大馅饼的体积),对当前区间的中间值mid,计算如果按照mid的尺寸分馅饼,能分给多少人?即对每个馅饼(设其体积为size),按照mid,计算能够分给的人数size/mid,并舍弃小数。如果每个朋友可以分一块,则low=mid,否则high=mid,继续二分。
本题的关键在于实数精度控制。本题要求误差绝对值不超过10-3。为此,设定实数的最小精度限制值esp=10-6,π=3.14159265359,否则,即使算法正确,也会得到Wrong Answer。
参考程序
#include<iostream> #include<iomanip> using namespace std; const double pi=3.14159265359; //最短的π长度,再短就会得到Wrong Answer const double esp=1e-6; //根据题目要求的精度,为了实数二分法设定的最小 //精度限制值 int main(void) { int test; cin>>test; while(test--) { int n,f; //馅饼数为n,朋友数为f cin>>n>>f; double* v=new double[n+1]; //每个馅饼的大小 f++; //加上自己的总人数 double maxsize=0.0; for(int i=1;i<=n;i++) { cin>>v[i]; v[i]*=v[i]; //半径平方,计算馅饼体积时先不乘π,以提高精 //度和减少时间 if(maxsize<v[i]) maxsize=v[i]; } double low=0.0; //下界,每人都分不到馅饼 double high=maxsize; //上界,每人都得到整个馅饼,而且是所有馅饼中最 //大的 double mid; while(high-low>esp) //实数精度控制作为循环结束条件 { mid=(low+high)/2; //计算按照mid的尺寸分馅饼,能分给多少人 int count_f=0; //根据mid尺寸能分给的人数 for(int i=1;i<=n;i++) //枚举每个馅饼 count_f+=(int)(v[i]/mid); //第i个馅饼按照mid的尺寸去切,最多能分的人 //数(向下取整) if(count_f < f) //当用mid尺寸分,可以分的人数小于朋友数 high=mid; //说明mid偏大,上界优化 else low=mid; //否则mid偏小,下界优化(注意'='一定要放在 //下界优化,否则精度会出错) } cout<<fixed<<setprecision(4)<<mid*pi<<endl; //之前的计算都只是利用半径平方去计 //算,最后的结果要记得乘π delete v; } return 0; }
二分法不仅可用于数据查找,也可用于函数计算。例如有变量x1、x2、x3,已知函数值是x1,x2、x3是函数的自变量,即存在函数关系式x1=f(x2,x3)。如何在函数值x1和自变量值x2确定的情况下,求满足函数式的x3值?下面介绍一种算法——减半递推技术。
所谓“减半”是指将问题的规模(例如x3的取值范围)减半,而问题的性质(例如x1=f(x2,x3))不变。假设原问题的规模为n,则可以采用与问题有关的特定方法将原问题化为c个(c是与规模无关而只与问题有关的常数)规模减半的问题,然后通过研究规模为的问题(显然与原问题性质相同,只是规模不同而已)来解决问题。问题的规模减少了,会给解决问题带来不少方便。但是,在规模减半过程中,势必增加某些辅助工作,在分析工作量时必须予以考虑。所谓“递推”是指重复上述“减半”过程。因为规模为的问题又可以转化为c个规模为的相同性质的问题。以此类推,直至问题的规模减少到最小,能够很方便地解决为止。
【1.4.2 Humidex】
湿热指数(humidex)是加拿大气象学家用来表示温度和湿度的综合影响的度量衡。它不同于美国所采用的露点(dew point),露点表示酷热指数而不表示相对的湿度(摘自维基百科)。
当温度为30℃(86℉)且露点为15℃(59℉)时,湿热指数是34(注意,湿热指数是一个没有度量单位的数,这个值表示大约的摄氏温度值)。如果温度保持在30℃并且露点上升到25℃(77℉),湿热指数就上升到42.3。
在温度和相对湿度相同的情况下,湿热指数往往比美国酷热指数高。
目前确定的湿热指数的公式是1979年由加拿大大气环境服务局的J.M.Masterton和F.A.Richardson给出的。
根据加拿大气象局的观点,湿热指数达到40会使人感到“非常不舒服”;湿热指数达到45以上就是“危险”状态;当湿热指数达到54,人就会马上中暑。
在加拿大,湿热指数的最高纪录是,1953年6月20日在加拿大安大略省温莎地区出现的,达到52.1。(温莎地区的居民并不知道,因为在当时尚未发明湿热指数。)1995年7月14日湿热指数在温莎和多伦多两地达到50。
湿热指数的公式如下:
湿热指数=温度+h
h=(0.5555)×(e-10.0)
e=6.11×exp[5417.7530×((1/273.16)-(1/(露点+273.16)))]
其中,exp(x)以2.718281828为底,x为指数。
由于湿热指数只是一个数字,电台播音员会像宣布气温一样宣布它,例如“那里气温47度……[暂停]……湿热指数”。有时天气预报会给出温度和露点,或者温度和湿热指数,但很少同时报告这三个测量值。请编写一个程序,给出其中两个测量值,计算第3个值。
本题设定,对于所有的输入,温度、露点和湿热指数的范围在-100~100°C之间。
输入
输入包含许多行,除最后一行之外,每行由空格分开的4个项组成:第一个字母,第一个数字,第二个字母,第二个数字。每个字母说明后面跟着的数字的含义:T表示温度,D表示露点,H表示湿热指数。输入结束行只有一个字母E。
输出
除最后一行之外,每行输入产生一行输出。输出的形式如下。
T数字 D数字 H数字
其中的3个数字给出温度、露点和湿热指数。每个数字都是十进制数,精确到小数点后一位,所有温度以摄氏度表示。
试题来源:Waterloo Local Contest,2007.7.14
在线测试:POJ 3299
试题解析
由题目给出的湿热指数公式可以看出,h与露点呈正比关系。显然,已知露点和温度(或湿热指数),则可先由露点推出h值,再由“湿热指数=温度+h”这一公式得出湿热指数;或者由公式“温度=湿热指数-h”得出温度。
若已知的两个测试量是温度和湿热指数,则我们采取减半递推技术计算露点值。
最初假设露点值为0,然后进入循环:露点的增量值从100开始,每次循环减半。若根据公式得出的湿热指数大于预报的湿热指数,则露点值减少一个增量(即h↘,使公式得出的湿热指数向下逼近预报值);否则露点值增加一个增量(即h↗,使公式得出的湿热指数向上逼近预报值)。这个循环过程直至增量值小于等于0.0001为止,此时得出应预报的露点值。
参考程序
#include <stdio.h> //预编译命令 #include <math.h> #include <assert.h> char a,b; //定义两个测试标志字符 double A,B,temp,hum,dew; double dohum(double tt,double dd){ //根据温度tt和露点dd计算湿热指数 double e = 6.11 * exp (5417.7530 * ((1/273.16) - (1/(dd+273.16)))); double h = (0.5555)*(e - 10.0); return tt + h; //返回湿热指数 } double dotemp(){ //根据露点dew和湿热指数hum计算温度 double e = 6.11 * exp (5417.7530 * ((1/273.16) - (1/(dew+273.16)))); double h = (0.5555)*(e - 10.0); return hum - h; //返回温度 } double dodew(){ //根据温度temp和湿热指数hum计算露点 double x = 0; //露点值及其增量初始化 double delta=100; //循环:增量值从100开始,每次循环减少一半,若根据当前温度temp和露点x得出的湿热指数大于hum,//则露点x减少一个增量值,否则露点x增加一个增量值,这个循环过程直至增量值小于等于0.0001为止 for (delta=100;delta>.00001;delta *=.5) { if (dohum(temp,x)>hum) x -= delta; else x += delta; } return x; //返回露点x } int main() //主函数 { //循环:依次输入每次天气预报的两个测试量,直至结束标志'E' while (4 == scanf(" %c %lf %c %lf",&a,&A,&b,&B) && a != 'E'){ temp = hum = dew = -99999; //温度、湿热指数和露点值初始化 if (a == 'T') temp = A; //预报的第1个测量值是温度 if (a == 'H') hum = A; //预报的第1个测量值是湿热指数 if (a == 'D') dew = A; //预报的第1个测量值是露点 if (b == 'T') temp = B; //预报的第2个测量值是温度 if (b == 'H') hum = B; //预报的第2个测量值是湿热指数 if (b == 'D') dew = B; //预报的第2个测量值是露点 if (hum== -99999) hum=dohum(temp,dew); //若缺失湿热指数,则根据温度和露点计算 if (dew == -99999) dew=dodew(); //若缺失露点,则根据温度和湿热指数计算 if (temp == -99999) temp = dotemp(); //若缺失温度,则根据湿热指数和露点计算 printf("T %0.1lf D %0.1lf H %0.1lf\n",temp,dew,hum); //输出温度、露点和湿热 //指数 } assert(a == 'E'); //遇到结束标志'E'时使用断言 }
在实数运算中,有时需要判断实数x和实数y是否相等,如果编程者把判断条件简单设成y-x是否等于0,就有可能产生精度误差。避免精度误差的办法是设一个精度常量delta。若y-x的实数值与0之间的区间长度小于delta,则认定x和y相等,这样就可将误差控制在delta范围内,如图1.4-1所示。显然,判断实数x和实数y是否相等的条件应设成|y-x|≤delta。
图 1.4-1
【1.4.3 Hangover】的解答涉及离线计算、二分法和实数精度的控制。
【1.4.3 Hangover】
你能使一叠在桌子上的卡片向桌子外伸出多远?如果是一张卡片,这张卡片能向桌子外伸出卡片的一半长度(卡片以直角伸出桌子)。如果有两张卡片,就让上面一张卡片向外伸出下面那张卡片的一半长度,而下面的那张卡片向桌子外伸出卡片的三分之一长度,所以两张卡片向桌子外延伸的总长度是1/2+1/3=5/6卡片长度。依次类推,n张卡片向桌子外延伸的总长度是1/2+1/3+1/4+…+1/(n+1)卡片长度:最上面的卡片向外延伸1/2,第二张卡片向外延伸1/3,第三张卡片向外延伸1/4,……,最下面一张卡片向桌子外延伸1/(n+1),如图1.4-2所示。
图 1.4-2
输入
输入由一个或多个测试用例组成,最后一行用0.00表示输入结束,每个测试用例一行,是一个3位正浮点数c,最小值为0.01,最大值为5.20。
输出
对每个测试数据c,输出要伸出卡片长度c最少要用的卡片的数目,输出形式见样例输出。
试题来源:ACM Mid-Central USA 2001
在线测试:POJ 1003,UVA 2294
试题解析
由于数据范围很小,因此先离线计算向桌子外延伸的卡片长度不超过5.20所需的最少卡片数。设total为卡片数,len[i]为前i张卡片向桌子外延伸的长度,即len[i]=len[i-1]+,i≥1,len[0]=0。显然len[]为递增序列。
注意:由于len的元素和被查找的要伸出卡片长度x为实数,因此要严格控制精度误差。设精度delta=1e-8。zero(x)为实数x是正负数和0的标志。
初始时len[0]=0,通过结构为for(total=1;zero(len[total-1]-5.20)<0;total++)len[total]=len[total-1]+1.0/double(total+1)的循环,递推计算len序列。
在计算出len数组后,先输入第1个测试用例x,并进入结构为while(zero(x))的循环,每一次循环,使用二分法在len表中查找伸出卡片长度x最少要用的卡片数,并输入下一个测试用例x。这个循环过程直至输入测试数据x=0.00为止。
二分查找的过程如下。
初始时区间[l,r]=[1,total],区间的中间指针。若zero(len[mid]-x)<0,则所需的卡片数在右区间(l=mid);否则所需的卡片数在左区间(r=mid)。继续二分区间[l,r],直至l+1≥r为止。此时得出的r即最少要用的卡片数。
参考程序:
#include <iostream> //预编译命令 using namespace std; //使用 C++标准程序库中的所有标识符 const int maxn = 300; //len数组容量 const double delta = 1e-8; //设定精度 int zero(double x) //在精度delta的范围内,若x是小于0的负实数, //则返回-1;若x是大于0的正实数,则返回1;若 //x为0,则返回0 { if (x < -delta) return -1; return x > delta; } int main(void) //主函数 { //主函数开始 double len[maxn]; //定义len数组和数组长度 int total; len[0] = 0.0; //直接计算出截止长度不超过5.20所需的最少卡片数 for (total = 1; zero(len[total - 1] - 5.20) < 0; total++) len[total] = len[total - 1] + 1.0 / double(total + 1); double x; cin >> x; //输入第1个测试用例x while (zero(x)) { // 用二分法在len表中查找不小于x的最少卡片数 int l, r; l = 0; //设定查找区间的左右指针 r = total; while (l + 1 < r) { //循环条件是查找区间存在 int mid = (l + r) / 2; //计算查找区间的中间指针 if (zero(len[mid] - x) < 0) //若中间元素值小于x,则在右区间查找;否则在左区 //间查找 l = mid; else r = mid; } cout << r << "card(s)" << endl; //输出至少伸出x长度最少要用的卡片数 cin >> x; //读下一个测试用例 } return 0; }