zoukankan      html  css  js  c++  java
  • 求解方程( 迭代法 牛顿迭代法 二分法)

    迭代的意思是反反复复地执行某一步骤、程序或者事件,在数学中用的比较常见。

    以上代码转自https://blog.csdn.net/pengwill97/article/details/77200372

    https://blog.csdn.net/akatsuki__itachi/article/details/80719686

    首先,迭代法解方程的实质是按照下列步骤构造一个序列x0,x1,…,xn,来逐步逼近方程f(x)=0的解:

    1)选取适当的初值x0;

    2)确定迭代格式,即建立迭代关系,需要将方程f(x)=0改写为x=φ(x)的等价形式;

    3)   构造序列x0,x1,……,xn,即先求得x1=φ(x0),再求x2=φ(x1),……如此反复迭代,就得到一个数列x0, x1,……,xn,若这个数列收敛,即存在极值,且函数 φ(x)连续,则很容易得到这个极限值,x*就是方程f(x)=0的根。

    #define eps 1e-8
    int main()
    {
        x0=初始近似根;
        do{
            x1=x0;
            x0=g(x1); //按特定的方程计算新的近似根
        }while(fabs(x0-x1)>Epsilon);
        printf("方程的近似根是%f
    ",x0);
    }

     注意:如果方程无解,算法求出的近似根序列就不会收敛,那么迭代过程就会变成死循环。因此,在使用迭代算法前应先考察方程是否有解,并在算法中对迭代次数给予限制。

    #include<stdio.h>
    #include<iostream>
    #include<math.h>
    #define eps 1e-8
    using namespace std;
    void NewtonIteration();//牛顿迭代法求解方程
    void DichotomySolving();//二分法求解方程
    int main ()
    {
        NewtonIteration();
        DichotomySolving();
    }
    //牛顿迭代法:xn=x(n-1)-f(x(n-1))/f'(x(n-1))
    void NewtonIteration()
    {
        int a = 3, b = 2, c = 1, d = -6;//系数
        float x1 = 1, x0;
        float f0, f1;//f0是原方程 f1是方程的一阶导
        do
        {
            x0 = x1;
            f0 = ((a * x0 + b) * x0 + c) * x0 + d;//即a*x^3+b*x^2+c*x+d
            f1 = (3 * a * x0 + 2 * b) * x0 + c;//一阶导数
            x1 = x0 - f0/f1;
        }while(fabs(x1 - x0) >= eps);
        printf("(%d)x^3+(%d)x^2+(%d)x+(%d) = 0
    ", a, b, c, d);
        cout << "方程的解为:" << x1 << endl;
    }
     
    void DichotomySolving()//二分法求解方程
    {
        float x, x1 = 0, x2 = 2;
        float f1, f2, f;
        cout << "input x1, x2 (f(x1)*f(x2)<0)" << endl;
        //假设原方程为 1/2x^3+2x^2-8
        f1 = pow(x1, 3) / 2 + pow(x1, 2) * 2 - 8;
        f2 = pow(x2, 3) / 2 + pow(x2, 2) * 2 - 8;
        if(f1 * f2 > 0)
        {
            cout << "error" << endl;
            return;
        }
        do
        {
            x = (x1 + x2) / 2;
            f = pow(x, 3) / 2 + 2 * pow(x, 2) - 8;//f为飞f1
            if(f == 0)
                break;
            if(f1 * f > 0)//已知解一定在f1与f2之间 因此以f1为切入点逐步逼近解
            {//f1 * f > 0证明解在右半部分中 因此要

    高斯消元法

    #include<stdio.h>
    #include<algorithm>
    #include<iostream>
    #include<string.h>
    #include<math.h>
    using namespace std;
    const int MAXN=50;
    int a[MAXN][MAXN];//增广矩阵
    int x[MAXN];//解集
    bool free_x[MAXN];//标记是否是不确定的变元
    int gcd(int a,int b)
    {
        if(b == 0) return a;
        else return gcd(b,a%b);
    }
    inline int lcm(int a,int b)
    {
        return a/gcd(a,b)*b;//先除后乘防溢出
    }
    // 高斯消元法解方程组(Gauss-Jordan elimination).(-2表示有浮点数解,但无整数解,
    //-1表示无解,0表示唯一解,大于0表示无穷解,并返回自由变元的个数)
    //有equ个方程,var个变元。增广矩阵行数为equ,分别为0到equ-1,列数为var+1,分别为0到var.
    int Gauss(int equ,int var)
    {
        int i,j,k;
        int max_r;// 当前这列绝对值最大的行.
        int col;//当前处理的列
        int ta,tb;
        int LCM;
        int temp;
        int free_x_num;
        int free_index;
    
        for(int i=0; i<=var; i++)
        {
            x[i]=0;
            free_x[i]=true;
        }
        //转换为阶梯阵.
        col=0; // 当前处理的列
        for(k = 0; k < equ && col < var; k++,col++) // 枚举当前处理的行.
        {
            // 找到该col列元素绝对值最大的那行与第k行交换.(为了在除法时减小误差)
            max_r=k;
            for(i=k+1; i<equ; i++)
            {
                if(abs(a[i][col])>abs(a[max_r][col])) max_r=i;
            }
            if(max_r!=k) // 与第k行交换.
            {
                for(j=k; j<var+1; j++) swap(a[k][j],a[max_r][j]);
            }
            if(a[k][col]==0) // 说明该col列第k行以下全是0了,则处理当前行的下一列.
            {
                k--;
                continue;
            }
            for(i=k+1; i<equ; i++) // 枚举要删去的行.
            {
                if(a[i][col]!=0)
                {
                    LCM = lcm(abs(a[i][col]),abs(a[k][col]));
                    ta = LCM/abs(a[i][col]);
                    tb = LCM/abs(a[k][col]);
                    if(a[i][col]*a[k][col]<0)tb=-tb;//异号的情况是相加
                    for(j=col; j<var+1; j++)
                    {
                        a[i][j] = a[i][j]*ta-a[k][j]*tb;
                    }
                }
            }
        }
        // 1. 无解的情况: 化简的增广阵中存在(0, 0, ..., a)这样的行(a != 0).
        for (i = k; i < equ; i++)  // 对于无穷解来说,如果要判断哪些是自由变元,那么初等行变换中的交换就会影响,则要记录交换.
        {
            if (a[i][col] != 0) return -1;
        }
        // 2. 无穷解的情况: 在var * (var + 1)的增广阵中出现(0, 0, ..., 0)这样的行,即说明没有形成严格的上三角阵.
        // 且出现的行数即为自由变元的个数.
        if (k < var)
        {
            return var - k; // 自由变元有var - k个.
        }
        // 3. 唯一解的情况: 在var * (var + 1)的增广阵中形成严格的上三角阵.
        // 计算出Xn-1, Xn-2 ... X0.
        for (i = var - 1; i >= 0; i--)
        {
            temp = a[i][var];
            for (j = i + 1; j < var; j++)
            {
                if (a[i][j] != 0) temp -= a[i][j] * x[j];
            }
            if (temp % a[i][i] != 0) return -2; // 说明有浮点数解,但无整数解.
            x[i] = temp / a[i][i];
        }
        return 0;
    }
    int main(void)
    {
        //    freopen("in.txt", "r", stdin);
        //    freopen("out.txt","w",stdout);
        int i, j;
        int equ,var;
        while (scanf("%d %d", &equ, &var) != EOF)
        {
            memset(a, 0, sizeof(a));
            for (i = 0; i < equ; i++)
            {
                
  • 相关阅读:
    “耐撕”团队 2016.3.25 站立会议
    “耐撕”团队 2016.03.24 站立会议
    “耐撕”团队 2016.3.22 站立会议
    windows环境下的git安装及使用
    词频统计(三)
    第二周作业
    Unity之GUI控件
    Lua的协同程序(coroutine)
    Lua与C++的交互
    Lua的元方法__newindex元方法
  • 原文地址:https://www.cnblogs.com/zcy19990813/p/9702734.html
Copyright © 2011-2022 走看看