计算方法考试

第 1 章 Lagrange 插值方法

#include <iostream>

using namespace std;

const long MAXSIZE = 100;

int main() {

    double x[MAXSIZE], y[MAXSIZE], xx, yy, t;

    int n;

    cout << "\n请输入插值节点的个数: ";

    cin >> n;

    for (int i = 0; i <= n - 1; i++) {

        cout << "\n请输入插值节点x[" << i << "], y[" << i << "]: ";

        cin >> x[i] >> y[i];

    }

    cout << "\n请输入插值点：";

    cin >> xx;

    yy = 0;

    for (int i = 0; i <= n - 1; i++) {

        t = 1;

        for (int j = 0; j <= n - 1; j++) {

            if (j != i) {

                t *= (xx - x[j]) / (x[i] - x[j]);

            }

        }

        yy += t * y[i];

    }

    cout << "\n插值点(x,y) = (" << xx << ", " << yy << ")。" << endl;

    return 0;

}

第 2 章 牛顿插值方法

#include <iostream>

#include <vector>

using namespace std;

double Difference(vector<double>& x, vector<double>& y, int start, int end) {

    if (start == end) {

        return y[start];

    } else {

        return (Difference(x, y, start + 1, end) - Difference(x, y, start, end - 1)) / (x[end] - x[start]);

    }

}

double newton(double x, vector<double>& xx, vector<double>& yy) {

    double result = 0.0;

    for (int i = 0; i < xx.size(); i++) {

        double term = Difference(xx, yy, 0, i);

        for (int j = 0; j < i; j++) {

            term *= (x - xx[j]);

        }

        result += term;

    }

    return result;

}

int main() {

    vector<double> xx (5);

    vector<double> yy (5);

    cout<<"输入五个插值点："<<endl;

    for(int i=0;i<5;i++)

    {

        cout<<"输入插值点x["<<i<<"] y["<<i<<"]: ";

        cin>>xx[i]>>yy[i];

    }

    double x_to_interpolate ;

    cout<<"请输入要插值的点： ";

    cin>>x_to_interpolate;

    double interpolated_value = newton(x_to_interpolate, xx, yy);

    cout << "结果为 x = " << x_to_interpolate << ": " << interpolated_value << endl;

    return 0;

}

第 3 章 Newton-Cotes 方法

#include <iostream>

#include <vector>

#include<bits/stdc++.h>

#define MAXSIZE 7

using namespace std;

// 定义新的被积函数

double customFunction(double x) {

    return sin(x)/x;

}

// 添加函数参数，接收用户定义的被积函数

void input(vector<double>& f, double a, double b, long n, double (*func)(double)) {

    long i;

    double h = (b - a) / (n - 1);

    cout << "\n请输入求积节点纵坐标:";

    // f[0]=1;

    // cout<<"x[0] = 0, f[0] = 1"<<endl;

    for (i = 0; i <= n - 1; i++) {

        double x = a + i * h;

        f[i] = func(x);  // 将用户定义的被积函数应用于节点值

        cout << "\nx[" << i << "] = " << x << ", f[" << i << "] = " << f[i];

    }

}

int main(void) {

    long c[MAXSIZE][MAXSIZE/2+2] = {{2, 1}, {6, 1, 4}, {8, 1, 3}, {90, 7, 32, 12}, {288, 19, 75, 50}, {840, 41, 216, 27, 272}, {17280, 751, 3577, 1323, 2989}};

    double a, b, integral;

    long n, i;

    cout << "\n请输入积分区间边界 a, b:";

    cin >> a >> b;

    cout << "\n请输入积分节点的个数 (2~8):";

    cin >> n;

    vector<double> f(n + 1);

    // 调用新的输入函数

    input(f, a, b, n, customFunction);

    integral = 0;

    for (i = 0; i < n / 2; i++)

        integral += (f[i] + f[n - i - 1]) * c[n - 2][i + 1] / c[n - 2][0];

    if (n % 2)

        integral += f[n / 2] * c[n - 2][n / 2 + 1] / c[n - 2][0];

    integral *= b - a;

    cout << "\n积分值为= " << integral << endl;

    return 0;

}

实验 4 求非线性方程根的牛顿法

#include <iostream>

#include <cmath>

using namespace std;

double f(double x) {

    return x * exp(x) - 1;

}

double df(double x) {

    return exp(x) + x * exp(x);//导数

}

int main(void) {

    double epsilon, x0, x1, fx0, dfx0;

    int i, maxi;

    cout << "\n请输入x的精度要求：";

    cin >> epsilon;

    cout << "\n请输入迭代初值：";

    cin >> x1;

    cout << "\n请输入最大迭代次数：";

    cin >> maxi;

    for (i = 0; i < maxi; i++) {

        x0 = x1;

        fx0 = f(x0);

        dfx0 = df(x0);

        x1 = x0 - fx0 / dfx0;

        if (fabs(x1 - x0) <= epsilon)

            break;

    }

    if (i < maxi)

        cout << "\n方程 f(x)=0 的根 x=" << x1 << "。" << endl;

    else

        cout << "\n迭代次数已超过上限。" << endl;

    return 0;

}

实验 5 解线性方程组的迭代法

#include<bits/stdc++.h>

#define MAXSIZE 100

using namespace std;

int main() {

    double A[MAXSIZE][MAXSIZE], x[MAXSIZE], b[MAXSIZE];

    int n;

    double e;

    cout << " 请输入原方程的阶数：";

    cin >> n;

    cout << "请输入原方程的增广矩阵：";

    for (int i = 0; i < n; ++i) {

        for (int j = 0; j < n ; ++j) {

            cin >> A[i][j];

        }

        cin >> b[i];

    }

    cout << "请输入初始迭代向量值";

    for (int i = 0; i < n; ++i) {

        cin >> x[i];

    }

    cout << "请输入误差上限";

    cin >> e;

    int count = 0;

    while (true) {

        cout << "第" << ++count << "次迭代：";

        int flag = 0;

        for (int i = 0; i < n; ++i) {

            double tmp = x[i];

            x[i] = 0;

            for (int j = 0; j < n; ++j) {

                if (i != j) {

                    x[i] += - A[i][j] * x[j];

                }

            }

            x[i] = (x[i] + b[i]) / A[i][i];

            if (fabs(x[i] - tmp) < e)

                ++flag;

               

            cout << x[i]<<" "<<fabs(x[i] - tmp)<<" ";;

        }

        cout << endl;

        if (flag == n) break;

    }

    for (int i = 0; i < n; ++i) {

        cout << x[i] << " ";

    }

}

实验 6 线性方程组的高斯消元法

#include <iostream>

#include <iomanip>

using namespace std;

void inputMatrix(double matrix[50][51], int n) {

    cout << "请输入增广矩阵（每行元素以空格分隔）：" << endl;

    for (int i = 0; i < n; ++i) {

        for (int j = 0; j <= n; ++j) {

            cin >> matrix[i][j];

        }

    }

}

void gaussElimination(double matrix[50][51], int n) {

    for (int i = 0; i < n; ++i) {

        double divisor = matrix[i][i];

        for (int j = i; j <= n; ++j) {

            matrix[i][j] /= divisor;

        }

       

        for (int k = 0; k < n; ++k) {

            if (k != i) {

                double factor = matrix[k][i];

                for (int j = i; j <= n; ++j) {

                    matrix[k][j] -= factor * matrix[i][j];

                }

            }

        }

    }

}

void outputSolution(double matrix[50][51], int n) {

    for (int i = 0; i < n; ++i) {

        cout << "x" << i + 1 << " = " << fixed << setprecision(6) << matrix[i][n] << endl;

    }

}

int main() {

    int n;

    double matrix[50][51];

    cout << "请输入线性方程组的阶数：";

    cin >> n;

    inputMatrix(matrix, n);

    gaussElimination(matrix, n);

    cout << "高斯消元法的解为：" << endl;

    outputSolution(matrix, n);

    return 0;

}

实验 7 线性方程组的矩阵分解法

#include <iostream>

#include <iomanip>

#define MAX_N 20

using namespace std;

int main() {

    int n;

    static double a[MAX_N][MAX_N], b[MAX_N], x[MAX_N], y[MAX_N];

    static double l[MAX_N][MAX_N], u[MAX_N][MAX_N];

    cout << "\n请输入 n 的值（方程组的阶数）: ";

    cin >> n;

    cout << "现在输入矩阵 a(i, j)，其中 i, j = 0, ..., " << n - 1 << ":\n";

    for (int i = 0; i < n; i++)

        for (int j = 0; j < n; j++)

            cin >> a[i][j];

    cout << "现在输入矩阵 b(i)，其中 i = 0, ..., " << n - 1 << ":\n";

    for (int i = 0; i < n; i++)

        cin >> b[i];

    for (int i = 0; i < n; i++)

        l[i][i] = 1;

    for (int k = 0; k < n; k++) {

        for (int j = k; j < n; j++) {

            u[k][j] = a[k][j];

            for (int i = 0; i <= k - 1; i++)

                u[k][j] -= (l[k][i] * u[i][j]);

        }

        for (int i = k + 1; i < n; i++) {

            l[i][k] = a[i][k];

            for (int j = 0; j <= k - 1; j++)

                l[i][k] -= (l[i][j] * u[j][k]);

            l[i][k] /= u[k][k];

        }

    }

    for (int i = 0; i < n; i++) {

        y[i] = b[i];

        for (int j = 0; j <= i - 1; j++)

            y[i] -= (l[i][j] * y[j]);

    }

    for (int i = n - 1; i >= 0; i--) {

        x[i] = y[i];

        for (int j = i + 1; j < n; j++)

            x[i] -= (u[i][j] * x[j]);

        x[i] /= u[i][i];

    }

    cout << "求解结果...x_i = \n";

    for (int i = 0; i < n; i++)

        cout << fixed << setprecision(6) << x[i] << "\n";

    return 0;

}

实验 8 常微分方程求解算法

#include <iostream>

#define N 10

using namespace std;

float f1(float x, float y) {

    return 3 * x - 2 * y * y - 12;

}

void modEuler(float (*f)(float, float), float x0, float y0, float xn, int n) {

    float h = (xn - x0) / n;

    float x = x0, y = y0;

    cout << "x[0] = " << x << "\ty[0] = " << y << endl;

    for (int i = 1; i <= n; i++) {

        float k1 = h * f(x, y);

        float k2 = h * f(x + h, y + k1);

       

        y = y + 0.5 * (k1 + k2);

        x = x0 + i * h;

        cout << "x[" << i << "] = " << x << "\ty[" << i << "] = " << y << endl;

    }

}

int main() {

    float xn = 0.0, x0 = -1.0, y0 = 2.0;

    modEuler(f1, x0, y0, xn, N);

    return 0;

}