HDU 5017 Ellipsoid

Ellipsoid

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/65536 K (Java/Others)
Total Submission(s): 1017    Accepted Submission(s): 361
Special Judge

Problem Description

Given a 3-dimension ellipsoid(椭球面)

your task is to find the minimal distance between the original point (0,0,0) and points on the ellipsoid. The distance between two points (x ₁,y ₁,z ₁) and (x ₂,y ₂,z ₂) is defined as 

 

Input

There are multiple test cases. Please process till EOF.

For each testcase, one line contains 6 real number a,b,c(0 < a,b,c,< 1),d,e,f (0 ≤ d,e,f < 1), as described above.  It is guaranteed that the input data forms a ellipsoid.All numbers are fit in double.

 

Output

For each test contains one line. Describes the minimal distance. Answer will be considered as correct if their absolute error is less than 10 ^-5.

 

Sample Input

  
  

   
   1 0.04 0.01 0 0 0
  
  

 

Sample Output

  
  

   
   1.0000000
  
  

 

Source

2014 ACM/ICPC Asia Regional Xi'an Online

 

学了下模拟退火，代码都是照人敲的。。

不过这题正解好像是拉氏函数

#include <cstdlib>
#include <cctype>
#include <cstring>
#include <cstdio>
#include <cmath>
#include <algorithm>
#include <vector>
#include <string>
#include <iostream>
#include <map>
#include <set>
#include <queue>
#include <stack>
#include <bitset>

using namespace std;

#define PB push_back
#define MP make_pair
#define REP(i,n) for(int i=0;i<(n);++i)
#define FOR(i,l,h) for(int i=(l);i<=(h);++i)
#define DWN(i,h,l) for(int i=(h);i>=(l);--i)
#define CLR(vis,pos) memset(vis,pos,sizeof(vis))
#define PI acos(-1.0)
#define INF 0x3f3f3f3f
#define LINF 1000000000000000000LL
#define eps 1e-8

typedef long long ll;

double r=0.99;//降温速度
int dx[8]={0,0,1,-1,1,-1,1,-1};
int dy[8]={1,-1,0,0,-1,1,1,-1};
double a,b,c,d,e,f;

double dis(double x,double y,double z){
    return sqrt(x*x+y*y+z*z);
}

double getz(double x, double y){
    double A=c,B=e*x+d*y,
        C=a*x*x+b*y*y+f*x*y-1;
    double pos=B*B-4*A*C;
    if (pos<0) return 1e60;
    double z1,z2;
    z1=(-B+sqrt(pos))/2.0/A;
    z2=(-B-sqrt(pos))/2.0/A;
    if (z1*z1<z2*z2) return z1;
    else             return z2;
}

void solve(){
    double step=1;
    double x,y,z,nx,ny,nz;
    x=0,y=0;
    while(step>eps){
        z=getz(x,y);
        REP(i,8){
            nx=x+step*dx[i];
            ny=y+step*dy[i];
            nz=getz(nx,ny);
            if(nz > 1e30) continue;
            if(dis(nx,ny,nz)<dis(x,y,z)){
                x=nx,y=ny,z=nz;
            }
        }
        step*=r;
    }
    printf("%.8lf\n",dis(x,y,z));
}

int main()
{
    //freopen("in.txt","r",stdin);
    //freopen("out.txt","w",stdout);
    while(scanf("%lf%lf%lf%lf%lf%lf",&a,&b,&c,&d,&e,&f)!=EOF){
        solve();
    }

    return 0;
}