��վ�ܷ�������

数论——高斯消元

主要是用于解多元线性方程

快速入门我就不想多讲

普通版高斯消元

大概思路是每次一行中某一列的最大系数,

用解方程的思路消掉这一列其余行的系数,每一往后找一列,且保持上一次消过的一行不再被消掉,

直到最后形成一个集中与右上角的直角三角行,然后从最后一列只有一个系数的方程不断往上代就能求出所有的答案

洛谷模板

先上代码

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
int n,m,flag;
const int N=501;
double mp[N][N],ans[N];
double eps=1e-7;
void gauss()
{
    for(int i=1;i<=n;i++)//列数
    {
        int p=i;
        for(int j=i+1;j<=n;j++)//找最大系数
          if(fabs(mp[p][i])<fabs(mp[j][i]))p=j;
        if(fabs(mp[p][i])<eps){flag=1;return;}
        if(p!=i)swap(mp[i],mp[p]);
        for(int j=i+1;j<=n;j++)//逐步消元
        {
            double bi=mp[j][i]/mp[i][i];
            for(int k=i;k<=n+1;k++)
                mp[j][k]-=mp[i][k]*bi;

        }
    }
    for(int i=n;i>=1;i--)//从下往上依次代数
    {
        ans[i]=mp[i][n+1]/mp[i][i];
        for(int j=1;j<i;j++)
        {
            mp[j][n+1]-=ans[i]*mp[j][i];
            mp[j][i]=0;//可省

        }
    }
}
int main()
{
   scanf("%d",&n);
   for(int i=1;i<=n;i++)
     for(int j=1;j<=n+1;j++)
       scanf("%lf",&mp[i][j]);
    gauss();
    if(flag){ printf("No Solution");return 0;}
   for(int i=1;i<=n;i++)
    printf("%.2lf\n",ans[i]);
    return 0;
}

简化版约旦消元

在高斯消元的基础上简化的求法,每一次把所有除自己行的系数全部消掉,最后每一行就得到的是一元一次方程,直接把系数一除过去就行了

我悟了,差别就是jordan把从第一行开始所有非自己行的全部消除一遍)

高斯就是不消除前面的只删后面)

#include <bits/stdc++.h>
#define eps 1e-5
using namespace std;
const int N=140;
int n,m,flag;
double mp[N][N];
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while(ch>'9'||ch<'0'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
    return x*f;
}
void guass()
{
    for(int i=1;i<=n;i++)//列数
    {
         int p=i;
         for(int j=i+1;j<=n;j++)
            if(fabs(mp[j][i])>fabs(mp[p][i]))p=j;
         if(fabs(mp[p][i])<eps){flag=1;return ;}
         if(p!=i)swap(mp[i],mp[p]);
         for(int j=1;j<=n;j++)
            if(j!=i)
         {
             double ai=mp[j][i]/mp[i][i];
             for(int k=1;k<=n+1;k++)
                mp[j][k]-=ai*mp[i][k];
         }
    }
}
main()
{
    n=read();
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n+1;j++)
          scanf("%lf",&mp[i][j]);
    guass();
     if(flag){ printf("No Solution");return 0;}
     for(int i=1;i<=n;i++)
        printf("%.2lf\n",mp[i][n+1]/mp[i][i]);
    return 0;
}



例题

1.poj1222 经典关灯问题

大概是建立30个方程求解


#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
int n=29,m,flag;
int _map[500][500],ans[500];
double eps=1e-7;
void gauss()
{

   for(int i=0;i<=n;i++)//列数
    {
        int p=i;
        for(int j=i+1;j<=n;j++)//找最大系数
          if(abs(_map[p][i])<abs(_map[j][i]))p=j;
        if(p!=i)swap(_map[i],_map[p]);
        for(int j=0;j<=n;j++)//逐步消元
          if(i!=j&&_map[j][i])
          {
             int bi=_map[j][i]/_map[i][i];
              for(int k=0;k<=n+1;k++)
                _map[j][k]=( _map[j][k]+2-bi*_map[i][k])%2;//注意负数

          }
         // ans[i]=_map[i][n+1]/_map[i][i];
    }
}

int main()
{
   int Case;
   int kn,km,kx,ky;
    scanf("%d",&Case);
   for(int k=1;k<=Case;k++)
   {
       memset(_map,0,sizeof(_map));
       for(int i=0;i<30;i++)
       {
           scanf("%d",&_map[i][30]);
           ans[i]=0;
       }
       for(int i=0;i<30;i++)
       {
           kn=i/6;km=i%6;
           for(int j=0;j<30;j++)
           {
               kx=j/6;ky=j%6;
               if(abs(kx-kn)+abs(ky-km)<=1)
                 _map[i][j]=1;
               else
                 _map[i][j]=0;

           }
       }
       gauss();
       printf("PUZZLE #%d\n",k);
       for(int i=0;i<30;i++)
       {
           printf("%d",ans[i]=_map[i][n+1]/_map[i][i]);
           if((i+1)%6==0)
            printf("\n");
           else printf(" ");
       }


   }

    return 0;
}

2.[JSOI2008]球形空间产生器

题目描述

有一个球形空间产生器能够在 $n$维空间中产生一个坚硬的球体。现在,你被困在了这个$n$维球体中,你只知道球面上$n+1$个点的坐标,你需要以最快的速度确定这个 $n$维球体的球心坐标,以便于摧毁这个球形空间产生器。

题解

这也是一个比较水的板子题。

首先提示:

距离:设两个n为空间上的点A, B的坐标为$(a_1, a_2, \cdots , a_n), (b_1, b_2, \cdots , b_n)$,

则AB的距离定义为:
$$
dist = \sqrt{ (a_1-b_1)^2 + (a_2-b_2)^2 + \cdots + (a_n-b_n)^2 }
$$

利用这个式子得到
$$
-2 a_0x_0+a_0^2-2a_1x_1+a_1^2····-2a_nx_n+a_n^2=-2 b_0x_0+b_0^2-2b_1x_1+b_1^2····-2b_nx_n+b_n^2
$$
移相得到
$$
x_0 \times 2(b_0-a_0)+x_1 \times 2(b_1-a_1)+···+x_n \times 2(b_n-a_n)=b_0^2-a_0^2+b_1^2-a_1^2+···+b_n^2-a_n^2
$$
于是列n个方程就显而易见了。。

代码

#include <bits/stdc++.h>
#define eps 1e-5
using namespace std;
const int N=140;
int n,m,flag;
double mp[N][N],a[N][N];
inline int read()
{
    int x=0,f=1;char ch=getchar();
    while(ch>'9'||ch<'0'){if(ch=='-')f=-1;ch=getchar();}
    while(ch>='0'&&ch<='9'){x=x*10+ch-'0';ch=getchar();}
    return x*f;
}
void gauss()
{
    for(int i=1;i<=n;i++)//列数
    {
         int p=i;
         for(int j=i+1;j<=n;j++)
            if(fabs(mp[j][i])>fabs(mp[p][i]))p=j;
         if(fabs(mp[p][i])<eps){flag=1;return ;}
         if(p!=i)swap(mp[i],mp[p]);
         for(int j=1;j<=n;j++)
            if(j!=i)
         {
             double ai=mp[j][i]/mp[i][i];
             for(int k=1;k<=n+1;k++)
                mp[j][k]-=ai*mp[i][k];
         }
    }
}
int main()
{
    n=read();
    for(int i=1;i<=n+1;i++)
        for(int j=1;j<=n;j++)
        scanf("%lf",&a[i][j]);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n;j++)
        mp[i][j]=2*(a[i+1][j]-a[i][j]);
    for(int i=1;i<=n;i++)
    {
        double tmp=0;
        for(int j=1;j<=n;j++)
            tmp+=a[i+1][j]*a[i+1][j]-a[i][j]*a[i][j];
        mp[i][n+1]=tmp;
    }
    gauss();
     for(int i=1;i<=n;i++)
        printf("%.3lf ",mp[i][n+1]/mp[i][i]);
    return 0;
}


3.poj1681画家问题

也是一个模板题,多年之前做的东西就不讲解了。。。

当年卡了很久

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>

using namespace std;
const int N=1e6+110;
const int M=1010;
int n,m,flag;
int mp[M][M],ans[N];
double eps=1e-7;
bool f[N];
int gauss()
{
   for(int i=0;i<n*n;i++)//列数
    {
        int p=i;
        for(int j=i+1;j<n*n;j++)//找最大系数
          if(abs(mp[p][i])<abs(mp[j][i]))p=j;

        if(p!=i)swap(mp[i],mp[p]);
        if(mp[i][i]==0) f[i] = true;
        for(int j=0;j<n*n;j++)//逐步消元
          if(i!=j&&mp[j][i]&&mp[i][i])
          {
             int bi=mp[j][i]/mp[i][i];
              for(int k=0;k<=n*n;k++)
                mp[j][k]=( mp[j][k]+2-bi*mp[i][k])%2;

          }

    }
    for(int i = 0; i < n * n; i++)
    {
        if(f[i] && mp[i][n * n]) return 0;

    }
    return 1;
}

int main()
{
   int Case=1;
   int kn,km,kx,ky;
   for(int k=1;k<=Case;k++)
   {
       int sum=0,cnt=0,flag=0;
       memset(mp,0,sizeof(mp));
       memset(ans,0,sizeof(ans));
       scanf("%d",&n);
       for(int i=0;i<n;i++)
       {
           char s[20];
           scanf("%s",s);
           for(int j=0;j<n;j++)
           {
               if(s[j]=='w')mp[cnt][n*n]=1;
               else mp[cnt][n*n]=0;
               cnt++;
           }
       }
       for(int i=0;i<n*n;i++)
       {
           kn=i/n;km=i%n;
           for(int j=0;j<n*n;j++)
           {
               kx=j/n;ky=j%n;
               if(abs(kx-kn)+abs(ky-km)<=1)
                 mp[i][j]=1;
               else
                 mp[i][j]=0;

           }
       }
       int t=gauss();
       if(!t){printf("inf\n");return 0;}

       for(int i=0;i<n*n;i++)
       {
           if(!mp[i][i])ans[i]=0;
           else ans[i]=mp[i][n*n]/mp[i][i];
           if(ans[i]==1)sum++;
       }

       printf("%d\n",sum);

   }

    return 0;
}