hdu5017 Ellipsoid 三分套三分

最新推荐文章于 2021-04-27 21:58:58 发布

原创最新推荐文章于 2021-04-27 21:58:58 发布 · 644 阅读

0 ·

CC 4.0 BY-SA版权

math 专栏收录该内容

10 篇文章

订阅专栏

该博客介绍了如何解决一个数学问题，即给定椭球方程，利用三分法求解椭球面上任意点到原点的最短距离。博主在实现过程中遇到一些困惑，如区间选择和精度控制，并表示最终解决方案依赖于参数调整。

题意大概是给一个椭球方程，求椭球面上到原点（0,0,0）的最小距离。

三分x然后三分y,用求根公式求出z返回距离..虽说是过来，但还是感觉有很多不科学的地方...比如三分区间（我感觉可以直接算出来啊，但直接算出区间结果不对..），还有eps的大小...感觉最后都是调参数调过去的......

/*=============================================================================
#  Author:Erich
#  FileName:
=============================================================================*/
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#include <cstring>
#include <string>
#include <vector>
#include <map>
#include <queue>
#include <stack>
#define lson id<<1,l,m
#define rson id<<1|1,m+1,r
using namespace std;
typedef long long ll;
const double inf=1e9;
const double PI=acos(-1.0);
const double eps=1e-12;
int n,m;
double a,b,c,d,e,f;
inline int dcmp(double x)
{
	if (fabs(x)<eps) return 0;
	if (x<0) return -1;
	return 1;
}
double slove(double x,double y)
{
	double z;
	double aa,bb,cc;
	aa=c;
	bb=e*x+d*y;
	cc=a*x*x+b*y*y+f*x*y-1;
	double dlt=bb*bb-4*aa*cc;
	if (dcmp(dlt)==-1) return 1e10+x*x+y*y;

	else if (dcmp(dlt)==0)
	{
		z=-bb/(2.0*a);
		return (x*x+y*y+z*z);
	}
	else
	{
		z=-bb+sqrt(dlt);
		z/=(2.0*aa);
		double mi=x*x+y*y+z*z;

		z=-bb-sqrt(dlt);
		z/=(2.0*aa);
		mi=min(mi,x*x+y*y+z*z);

		return mi;
	}
}
double search(double x)
{
	double L,R;
	double tc=1-a*x*x;
	L=-sqrt(tc);
	R=sqrt(tc);
	L=-inf; R=inf;

	for (int i=0; i<200; i++)
	{
		double m1=L+(R-L)/3;
		double m2=R-(R-L)/3;
		if (slove(x,m1)>slove(x,m2)) L=m1; else R=m2;
	}
	return slove(x,(L+R)/2.0);
}
int main()
{
//	freopen("in.txt","r",stdin);
//	freopen("out.txt","w",stdout);
	while(~scanf("%lf",&a))
	{
		scanf("%lf%lf%lf%lf%lf",&b,&c,&d,&e,&f);
//		double L=-sqrt(1.0/a),R=sqrt(1.0/a);
        double L=-inf,R=inf;
		for (int i=0; i<200; i++)
		{
			double m1=L+(R-L)/3;
			double m2=R-(R-L)/3;
//			printf("%.10f %.10lf\n",L,R);
//			printf("%.10f %.10lf\n",search(m1),search(m2));
//			puts("----------");
			if (search(m1)>search(m2)) L=m1; else  R=m2;
		}
		printf("%.10lf\n",sqrt(search((L+R)/2.0)));
	}
	return 0;
}