博客
关于我
强烈建议你试试无所不能的chatGPT,快点击我
三分法——求解凸性函数的极值问题
阅读量:4703 次
发布时间:2019-06-10

本文共 6431 字,大约阅读时间需要 21 分钟。

       二分法作为分治中最常见的方法,适用于单调函数,逼近求解某点的值。但当函数是凸性函数时,二分法就无法适用,这时三分法就可以“大显身手”~~

       如图,类似二分的定义Left和Right,mid = (Left + Right) / 2,midmid = (mid + Right) / 2; 如果mid靠近极值点,则Right = midmid;否则(即midmid靠近极值点),则Left = mid;
程序模版如下:

double Calc(Type a){    /* 根据题目的意思计算 */}void Solve(void){    double Left, Right;    double mid, midmid;    double mid_value, midmid_value;    Left = MIN; Right = MAX;    while (Left + EPS < Right)    {        mid = (Left + Right) / 2;        midmid = (mid + Right) / 2;        mid_area = Calc(mid);        midmid_area = Calc(midmid);        // 假设求解最大极值.        if (mid_area >= midmid_area) Right = midmid;        else Left = mid;    }}

现根据几道的OJ题目来分析三分法的具体实现。

描述

In this problem, you're to calculate the distance between a point P(xp, yp, zp) and a segment (x1, y1, z1) ? (x2, y2, z2), in a 3D space, i.e. the minimal distance from P to any point Q(xq, yq, zq) on the segment (a segment is part of a line).
输入
The first line contains a single integer T (1 ≤ T ≤ 1000), the number of test cases. Each test case is a single line containing 9 integers xp, yp, zp, x1, y1, z1, x2, y2, z2. These integers are all in [-1000,1000].
输出
For each test case, print the case number and the minimal distance, to two decimal places.
样例输入
3
0 0 0 0 1 0 1 1 0
1 0 0 1 0 1 1 1 0
-1 -1 -1 0 1 0 -1 0 -1
样例输出
Case 1: 1.00
Case 2: 0.71
Case 3: 1.00题意为在一条线段上找到一点,与给定的P点距离最小。很明显的凸性函数,用三分法来解决。

Calc函数即为求某点到P点的距离。

#include
#include
#include
#include
#include
#define eps 1e-8using namespace std;struct point { double x,y,z; }; point l,r,p; point MID(point a,point b){ point k; k.x=(a.x+b.x)*0.5; k.y=(a.y+b.y)*0.5; k.z=(a.z+b.z)*0.5; return k;}double dist(point a,point b){ return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)+(a.z-b.z)*(a.z-b.z));}int sgn(double a){ return (a>eps)-(a<-eps);}point work(){ point mid,midmid; while(sgn(dist(l,r))>0) { mid=MID(l,r); midmid=MID(mid,r); if(dist(mid,p)

如图,人左右走动,求影子L的最长长度。
根据图,很容易发现当灯,人的头部和墙角成一条直线时(假设此时人站在A点),此时的长度是影子全在地上的最长长度。当人再向右走时,影子开始投影到墙上,当人贴着墙,影子长度即为人的高度。所以当人从A点走到墙,函数是先递增再递减,为凸性函数,所以我们可以用三分法来求解。
下面只给出Calc函数,其他直接套模版即可。
double Calc(double x)
{
    return (h * D - H * x) / (D - x) + x;
}

//-------common header---------------#include 
#include
#include
#include
#include
typedef unsigned long u32;using namespace std;int main(){ int caseNum = 0; scanf("%d",&caseNum); while(caseNum--) { double H,h,D; scanf("%lf %lf %lf",&H,&h,&D); double result = 0.0; if(H-h>=D) { result = h; } else if(H*H<=D*(H-h)) { result = h*D/H; } else { double a = -sqrt((H-h)*D)+H; result = D*(h-a)/(H-a)+a; } printf("%.3f/n",(float)result); } return 1;}

汽车拐弯问题,给定X, Y, l, d判断是否能够拐弯。首先当X或者Y小于d,那么一定不能。
其次我们发现随着角度θ的增大,最大高度h先增长后减小,即为凸性函数,可以用三分法来求解。
这里的Calc函数需要比较繁琐的推倒公式:
s = l * cos(θ) + w * sin(θ) - x;
h = s * tan(θ) + w * cos(θ);
其中s为汽车最右边的点离拐角的水平距离, h为里拐点最高的距离, θ范围从0到90。

#include 
#include
#include
using namespace std;double pi = acos(-1.0);double x,y,l,w,s,h;double cal(double a){ s = l*cos(a)+w*sin(a)-x; h = s*tan(a)+w*cos(a); return h;}int main(){ double left,right,mid,midmid; while(scanf("%lf%lf%lf%lf",&x,&y,&l,&w)!=EOF) { left = 0.0; right = pi/2; while(fabs(right-left)>1e-8) { mid = (left+right)/2; midmid = (mid+right)/2; if(cal(mid)>=cal(midmid))right = midmid; else left = mid; } if(cal(mid)<=y)printf("yes\n"); else printf("no\n"); } return 0;}

题意为给定n(n <= 30)个点,求出饱含这些点的面积最小的正方形。
有两种解法,一种为逼近法,就是每次m分角度,求出最符合的角度,再继续m分,如此进行times次,即可求出较为精确的解。(m 大概取10, times取30即可)
第二种解法即为三分法,首先旋转的角度只要在0到180度即可,超过180度跟前面的相同的。坐标轴旋转后,坐标变换为:
X’ = x * cosa - y * sina;
y’ = y * cosa + x * sina;
至于这题的函数是否是凸性的,为什么是凸性的,我也无法给出准确的证明,希望哪位路过的大牛指点一下~~

#include 
#include
using namespace std;#define MAX 33#define eps 0.00000005#define max(a,b) (a>b?a:b)struct { double x, y; } p[MAX];double cal1 ( int n, double d ){ double dis1 = 0.0, temp; for ( int i = 1; i < n; i++ ) { for ( int j = i + 1; j <= n; j++ ) { temp = fabs( (p[i].y-p[j].y) * sin(d) + (p[i].x-p[j].x) * cos(d)); dis1 = max ( dis1, temp ); } } return dis1;}double cal2 ( int n, double d ){ double dis2 = 0.0, temp; for ( int i = 1; i < n; i++ ) { for ( int j = i + 1; j <= n; j++ ) { temp = fabs( (p[i].y-p[j].y) * cos(d) - (p[i].x-p[j].x) * sin(d)); dis2 = max ( dis2, temp ); } } return dis2;} int main(){ int cs, n; scanf("%d",&cs); while ( cs-- ) { scanf("%d",&n); for ( int i = 1; i <= n; i++ ) scanf("%lf%lf", &p[i].x, &p[i].y); double ans = 1000, low = 0, high = asin(1.0); while ( high - low > eps ) { double mid1 = low + ( high - low ) / 3; double mid2 = high - ( high - low ) / 3; double len1 = max ( cal1 ( n, mid1 ), cal2 ( n, mid1 )); double len2 = max ( cal1 ( n, mid2 ), cal2 ( n, mid2 )); if ( len1 < len2 ) { ans = len1; high = mid2; } else { ans = len2; low = mid1; } } printf("%.2lf\n",ans*ans); } return 0;}

典型的三分法,先三分第一条线段,找到一个点,然后根据这个点再三分第二条线段即可,想出三分的思路基本就可以过了。

对于求解一些实际问题,当公式难以推导出来时,二分、三分法可以较为精确地求解出一些临界值,且效率也是令人满意的。

#include
#include
using namespace std;#define eps 1e-6struct POINT{ double x,y;}a,b,c,d;double p,q,r;double dis(POINT A,POINT B){ return sqrt((A.x-B.x)*(A.x-B.x)+(A.y-B.y)*(A.y-B.y)+eps);}double ftime(POINT t,double distd)//time between t and D{ POINT temp; temp.x=d.x+(c.x-d.x)*((q*distd)/dis(c,d)); temp.y=d.y+(c.y-d.y)*((q*distd)/dis(c,d)); return dis(t,temp)/r+distd;}double time(double posa)//posa->time between t and A{ POINT t; t.x=a.x+(b.x-a.x)*((posa*p)/dis(a,b)); t.y=a.y+(b.y-a.y)*((posa*p)/dis(a,b)); //test distance on CD double low=0,high=dis(c,d)/q,mid,m; while(high-low>eps) { m=(high+low)/2; mid=(high+m)/2; if(ftime(t,m)<=ftime(t,mid)) high=mid; else low=m; } return ftime(t,low)+posa;}int main(){ freopen("1.txt","w",stdout); int T; cin>>T; while(T--) { cin>>a.x>>a.y>>b.x>>b.y>>c.x>>c.y>>d.x>>d.y; cin>>p>>q>>r; double low=0,high=dis(a,b)/p,mid,m; while(high-low>eps) { m=(low+high)/2; mid=(low+m)/2; if(time(mid)>=time(m)) low=mid; else high=m; } printf("%.2lf/n",time(low)); }}

转载于:https://www.cnblogs.com/bo-jwolf/archive/2013/04/14/3033103.html

你可能感兴趣的文章
Android-事件分发(ViewGroup)
查看>>
Djngo 请求的生命周期
查看>>
[kuangbin带你飞]专题1-23题目清单总结
查看>>
EF6.0批量插入
查看>>
Knockout v3.4.0 中文版教程-1-入门和安装
查看>>
JavaScript 基础一
查看>>
数据库范式1NF 2NF 3NF BCNF
查看>>
ComboBox控件隐藏fieldLabel不能隐藏问题解决
查看>>
Hibernate自动生成DO手写DAO的注意事项
查看>>
IDEA中打包Spark项目提示Error:(16, 48) java: -source 1.5 中不支持 lambda 表达式
查看>>
poj 3013 Big Christmas Tree 最短路 dijkstra算法
查看>>
Java 多线程,Thread,Runnable,Callable
查看>>
HTTP长连接、短连接究竟是什么?
查看>>
这是二零一四年十点整的广州
查看>>
linux学习笔记(十四)
查看>>
Groovy 与 DSL
查看>>
统一异常处理@ExceptionHandler
查看>>
线程安全的简单理解
查看>>
[Silverlight] Visual Studio2010不能安装Silverlight4_Tools,提示语言不一致
查看>>
Java基础(十)数据结构
查看>>