• bzoj1502 simpson求面积


    题目:http://www.lydsy.com/JudgeOnline/problem.php?id=1502

    题解:

    simpson积分求面积,s = (f(a)+f(b)+4*f(c))/6*Δx ,c=(a+b)/2。

    题中的树投影下来是一些圆和相邻圆的公切线组成的一个封闭图形,并且上下对称,所以可以只求上半部分。

    simpson求面积时,若f(x)代价很大,要尽量减少其的重复调用。其次尽量优化f(x)函数的计算。

    写完后还要自己出极限随机数据,将eps调到能接受的最大值。

      1 #include <cstdio>
      2 #include <cmath>
      3 #include <iostream>
      4 #define maxn 510
      5 #define eps 1e-6
      6 using namespace std;
      7 
      8 int sg( double x ) {
      9     return (x>-eps)-(x<eps);
     10 }
     11 struct Vector {
     12     double x, y;
     13     Vector(){}
     14     Vector( double x, double y ) : x(x), y(y) {}
     15     Vector operator+( Vector b ) { return Vector(x+b.x,y+b.y); }
     16     Vector operator-( Vector b ) { return Vector(x-b.x,y-b.y); }
     17     Vector operator*( double b ) { return Vector(x*b,y*b); }
     18     double len() {
     19         return sqrt( x*x+y*y );
     20     }
     21 };
     22 typedef Vector Point;
     23 
     24 struct Line {
     25     Point A, B;
     26     double k, b;
     27     double lf, rg;
     28     Line(){}
     29     Line( Point A, Point B ):A(A),B(B) {
     30         if( A.x>B.x ) swap(A,B);
     31         lf = A.x;
     32         rg = B.x;
     33         k = (B.y-A.y)/(B.x-A.x);
     34         b = A.y-A.x*k;
     35     }
     36     inline bool in( double x ) {
     37         return sg(x-lf)>=0 && sg(rg-x)>=0;
     38     }
     39     inline double f( double x ) { return k*x+b; }
     40     void out() {
     41         printf( "(%lf,%lf) (%lf,%lf) 
    ", A.x, A.y, B.x, B.y );
     42     }
     43 };
     44 struct Circle {
     45     Point C;
     46     double r;
     47     double lf, rg;
     48     Circle(){}
     49     void make() {
     50         lf = C.x-r;
     51         rg = C.x+r;
     52     }
     53     inline bool in( double x ) {
     54         return sg(x-lf)>=0 && sg(rg-x)>=0;
     55     }
     56     inline double f( double x ) {
     57         return sqrt( r*r-(C.x-x)*(C.x-x) );
     58     }
     59     void out() {
     60         printf( "(%lf,%lf) r=%lf
    ", C.x, C.y, r );
     61     }
     62 };
     63 
     64 int n; 
     65 double alpha;
     66 Line lines[maxn];
     67 Circle cirs[maxn];
     68 
     69 Line cir_tang( Circle & a, Circle & b ) {
     70     double ang = acos( (a.r-b.r)/(b.C.x-a.C.x) );
     71     Point A = a.C + Vector(cos(ang),sin(ang))*a.r;
     72     Point B = b.C + Vector(cos(ang),sin(ang))*b.r;
     73     return Line(A,B);
     74 }
     75 
     76 double F( double x ) {
     77     double y = -1.0;
     78     for( int i=0; i<=n; i++ ) {
     79         if( cirs[i].lf-eps<=x && x<=cirs[i].rg+eps ) 
     80             y = max( y, cirs[i].f(x) );
     81         if( i&& lines[i].lf-eps<=x && x<=lines[i].rg+eps ) 
     82             y = max( y, lines[i].f(x) );
     83     }
     84     return y;
     85 }
     86 inline double simpson( double a, double b, double fa, double fb, double fc ) {
     87     return (fa+4*fc+fb)*(b-a)/6;
     88 }
     89 double adapt( double a, double b, double c, double fa, double fb, double fc ) {
     90     double d = (a+c)/2, e = (c+b)/2;
     91     double fd = F(d), fe = F(e);
     92     double sa = simpson(a,c,fa,fc,fd);
     93     double sb = simpson(c,b,fc,fb,fe);
     94     double ss = simpson(a,b,fa,fb,fc);
     95     if( fabs(sa+sb-ss)<15*eps ) return sa+sb;
     96     return adapt(a,c,d,fa,fc,fd)+adapt(c,b,e,fc,fb,fe);
     97 }
     98 
     99 int main() {
    100     scanf( "%d%lf", &n, &alpha );
    101     double k = 1/tan(alpha);
    102     double hsum = 0.0, h;
    103     double lf=1e200, rg=-1e200;
    104     for( int i=0; i<=n; i++ ) {
    105         scanf( "%lf", &h );
    106         hsum += h;
    107         cirs[i].C.x = hsum*k;
    108         cirs[i].C.y = 0.0;
    109     }
    110     for( int i=0; i<n; i++ )  {
    111         scanf( "%lf", &cirs[i].r );
    112         cirs[i].make();
    113         lf = min( lf, cirs[i].lf );
    114         rg = max( rg, cirs[i].rg );
    115     }
    116     cirs[n].r = 0.0;
    117     rg = max( rg, cirs[n].C.x );
    118     for( int i=1; i<=n; i++ ) 
    119         lines[i] = cir_tang( cirs[i-1], cirs[i] );
    120     printf( "%.2lf
    ", adapt(lf,rg,(lf+rg)/2,F(lf),F(rg),F((lf+rg)/2) )*2 );
    121 }
    View Code
  • 相关阅读:
    pm2 配置
    添加项目到远程服务器(git)
    psql 命令行使用
    SQL
    iOS AFNetworking 打印从服务器返回的错误提示信息
    iOS 获取网络图片的大小
    iOS 10 常见配置的问题
    LGLTagsView
    xcode8 关闭控制台打印不用信息
    LGLProgressHUD
  • 原文地址:https://www.cnblogs.com/idy002/p/4275443.html
Copyright © 2020-2023  润新知