• 洛谷P4035 [JSOI2008]球形空间产生器(高斯消元)


    洛谷题目传送门

    球啊球 @xzz_233 qaq

    高斯消元模板题,关键在于将已知条件转化为方程组。

    可以发现题目要求的未知量有(n)个,题目却给了我们(n+1)个点的坐标,这其中必有玄机。

    由高中数学知识可以知道,三点定圆(二维),四点定球(三维)······以此类推,应该是(n+1)个点才能确定一个(n)维空间下的球。

    那么隐藏的另一个关键未知量在哪里呢?

    想想圆的标准方程((x-x_0)^2+(y-y_0)^2=r^2),除了圆心坐标,半径不也对这个圆起到决定性作用么?

    接下来,额外设一个未知量——球的半径(r),开始试着对条件式进行变换。

    对于(n+1)个点,它们与球心的距离是定值(r),那么我们可以得到形式如下的(n+1)个方程(a为球面点坐标,x为球心坐标)

    [(a_1-x_1)^2+(a_2-x_2)^2+···+(a_n-x_n)^2=r^2 ]

    显然我们要把已知量和未知量分开,于是展开,移项

    [a_1^2-2a_1x_1+x_1^2+a_2^2-2a_2x_2+x_2^2+···+a_n^2-2a_nx_n+x_n^2=r^2 ]

    [2a_1x_1+2a_2x_2+···+2a_nx_n+r^2-x_1^2-x_2^2-···-x_n^2=a_1^2+a_2^2+···+a_n^2 ]

    发现(r^2-x_1^2-x_2^2-···-x_n^2)(a)无关,所以考虑换元,设(t=r^2-x_1^2-x_2^2-···-x_n^2)(实际上我们并不用求(r)

    终于,我们可以看到一个关于(x_1,x_2,···,x_n,t)((n+1))元方程组了,上高斯消元

    具体实现看代码

    #include<cmath>
    #include<cstdio>
    #define R register
    #define init for(i=1;i<=n;++i)ne[i-1]=pr[i+1]=i
    const int N=19;
    int p[N],pr[N],ne[N];
    double a[N][N];
    int main(){
        R int n,i,j,k,x;
        R double mx,d;
        scanf("%d",&n);++n;
        init;//链表初始化,为了实现交换行
        for(i=1;i<=n;++i){
            for(j=1;j<n;++j){
                scanf("%lf",&d);//处理系数
                a[i][j]=d*2;a[i][n+1]+=d*d;
            }
            a[i][n]=1;//t的系数为1
        }
        for(k=1;k<=n;++k){
            mx=0;//蒟蒻没有交换主元,而是交换行
                    //这样做防掉精度的效果可能不如交换主元
            for(i=ne[0];i;i=ne[i])
                if(mx<fabs(a[i][k]))
                    mx=fabs(a[i][k]),x=i;
            d=a[p[k]=x][k];//选择当前a最大的一行
            pr[ne[pr[x]]=ne[x]]=pr[x];
            for(j=1;j<=n+1;++j)
                a[x][j]/=d;
            for(i=ne[0];i;i=ne[i])
                for(d=a[i][j=k];j<=n+1;++j)
                    a[i][j]-=d*a[x][j];
        }//高斯消元
        init;
        for(k=n;k;--k){
            d=a[x=p[k]][n+1];
            pr[ne[pr[x]]=ne[x]]=pr[x];
            for(i=ne[0];i;i=ne[i])
                a[i][n+1]-=d*a[i][k];
        }//回代
        for(k=1;k<n;++k)
            printf("%.3f ",a[p[k]][n+1]);
        puts("");
        return 0;
    }
    
  • 相关阅读:
    使用Delphi自带的TDockTabSet组件实现停靠功能(Jeremy North)
    揭秘换肤技术(转载)
    cdecl、stdcall、fastcall函数调用约定区别(转)
    Delphi XE的RTTI增强,动态Hook某些内部事件
    sizeof和strlen解析
    由swap引出的局部变量,形参和指针的小问题
    单链表的逆置算法
    关于C++中的虚拟继承的一些总结
    关于"引用"的几点说明
    求最长公共子序列(子序列在原串中可以不连续)
  • 原文地址:https://www.cnblogs.com/flashhu/p/9193318.html
Copyright © 2020-2023  润新知