欧拉法求解常微分方程（c++）【转载】

 摘自《c++和面向对象数值计算》，代码简洁明快，采用类进行封装实现代码，增强代码的重用性，通过继承可实现代码的重用，采用函数指针，通用性增强，在函数改变时只需要单独改变函数部分的代码，无需对所有代码进行重写，对其中代码稍加改动，源代码有缺陷，没有释放内存，会造成内存泄露，在构造函数当中，将源代码的在函数体中赋值，改为列表赋值，可以稍微提高代码的执行效率。 

#include<iostream>

#include <cmath>

#include <algorithm>

using namespace std;

   

class ode

{

      double tini;      //初始时间

      double ison;      //初始解

      double tend;      //结束时间

      double(*sfn)(double t, double x);      //源函数，此处采用函数指针的方式，更具通用性

public:

      ode(double t0, double x0, double T, double(*f)(double, double)) :      //类的构造函数

            tini(t0), ison(x0), tend(T), sfn(f){}

      double* euler(int n) const;            //n子区间欧拉方法

};

   

double f(double t, double x)      //源函数

{

      return x*(1 - exp(t)) / (1 + exp(t));

}

   

double exact(double t)      //真实解

{

      return 12 * exp(t) / pow(1 + exp(t), 2);

}

   

int main()

{

      ode exmp(0, 3, 2, f);      //x(0)=3,T=2

      double* soln = exmp.euler(100);            //100个子区间

      double norm = 0;

      double h = 2.0 / 100;      //计算步长

      for (int k = 1; k <= 100; k++)

            norm = max(norm, fabs(exact(k*h) - soln[k]));

      cout << "max norm of error by euler's method = "

            << norm << endl;

      delete[] soln;

      return 0;

}

   

double* ode::euler(int n) const            //显式欧拉方法

{

      double* x = new double[n + 1];      //近似解

      double h = (tend - tini) / n;      //步长

      x[0] = ison;            //初始解

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

            x[k + 1] = x[k] + h*sfn(tini + k*h, x[k]);

      return x;

}
 

转载于:https://www.cnblogs.com/liusuanyatong/p/11259909.html