【数值计算方法(黄明游)】常微分方程初值问题的数值积分法:欧拉方法(向前Euler)【理论到程序】

本文主要是介绍【数值计算方法(黄明游)】常微分方程初值问题的数值积分法:欧拉方法(向前Euler)【理论到程序】,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

文章目录

  • 一、数值积分法
    • 1. 一般步骤
    • 2. 数值方法
  • 二、欧拉方法(Euler Method)
    • 1. 向前欧拉法(前向欧拉法)
      • a. 基本理论
      • b. 典例解析
      • c. 算法实现

  常微分方程初值问题的数值积分法是一种通过数值方法求解给定初始条件下的常微分方程(Ordinary Differential Equations, ODEs)的问题。

一、数值积分法

1. 一般步骤

  • 确定微分方程:

    • 给定微分方程组 y ′ ( x ) = f ( x , y ( x ) ) y'(x) = f(x, y(x)) y(x)=f(x,y(x))
  • 确定初始条件:

    • 初值问题包含一个初始条件 y ( a ) = y 0 y(a) = y_0 y(a)=y0,其中 a a a 是定义域的起始点, y 0 y_0 y0 是初始值。
  • 选择数值方法:

    • 选择适当的数值方法来近似解(需要考虑精度、稳定性和计算效率),常见的数值方法包括欧拉方法、改进的欧拉方法、Runge-Kutta 方法等。
  • 离散化定义域:

    • 将定义域 [ a , b ] [a, b] [a,b] 分割为若干小步,即选择合适的步长 h h h。通常,较小的步长能够提高数值解的精度,但也增加计算成本。
  • 数值迭代:

    • 使用选定的数值方法进行迭代计算:根据选择的方法,计算下一个点的函数值,并更新解。
  • 判断停止条件:

    • 判断是否达到满足指定精度的近似解:可以使用某种误差估计方法,例如控制局部截断误差或全局误差。
  • 输出结果:

    • 最终得到在给定定义域上满足初值问题的近似解。

2. 数值方法

  1. 欧拉方法(Euler Method):

    • 基本思想:根据微分方程的定义,使用离散步长逼近导数,进而逼近下一个点的函数值。
    • 公式: y n + 1 = y n + h f ( t n , y n ) y_{n+1} = y_n + h f(t_n, y_n) yn+1=yn+hf(tn,yn)
      其中, y n y_n yn是第 n n n 步的函数值, h h h是步长, f ( t n , y n ) f(t_n, y_n) f(tn,yn) 是在点 ( t n , y n ) (t_n, y_n) (tn,yn) 处的导数。
  2. 改进的欧拉方法(Improved Euler Method 或梯形法 Trapezoidal Rule):

    • 基本思想:使用两次近似来提高精度,首先使用欧拉方法计算中间点,然后用该点的导数估计值来计算下一个点。
    • 公式: y n + 1 = y n + h 2 [ f ( t n , y n ) + f ( t n + 1 , y n + h f ( t n , y n ) ) ] y_{n+1} = y_n + \frac{h}{2} [f(t_n, y_n) + f(t_{n+1}, y_n + hf(t_n, y_n))] yn+1=yn+2h[f(tn,yn)+f(tn+1,yn+hf(tn,yn))]
  3. Runge-Kutta 方法:

    • 基本思想:通过多个阶段的计算来提高精度。其中最常见的是四阶 Runge-Kutta 方法。
    • 公式:
      k 1 = h f ( t n , y n ) k_1 = hf(t_n, y_n) k1=hf(tn,yn) k 2 = h f ( t n + h 2 , y n + k 1 2 ) k_2 = hf(t_n + \frac{h}{2}, y_n + \frac{k_1}{2}) k2=hf(tn+2h,yn+2k1) k 3 = h f ( t n + h 2 , y n + k 2 2 ) k_3 = hf(t_n + \frac{h}{2}, y_n + \frac{k_2}{2}) k3=hf(tn+2h,yn+2k2) k 4 = h f ( t n + h , y n + k 3 ) k_4 = hf(t_n + h, y_n + k_3) k4=hf(tn+h,yn+k3) y n + 1 = y n + 1 6 ( k 1 + 2 k 2 + 2 k 3 + k 4 ) y_{n+1} = y_n + \frac{1}{6}(k_1 + 2k_2 + 2k_3 + k_4) yn+1=yn+61(k1+2k2+2k3+k4)

  这些方法中,步长 h h h 是一个关键参数,它决定了离散化的程度,选择合适的步长对于数值解的准确性和稳定性非常重要。

二、欧拉方法(Euler Method)

1. 向前欧拉法(前向欧拉法)

a. 基本理论

  1. 等距节点组:

    • { X n } \{X_n\} {Xn}被定义为区间 [ a , b ] [a, b] [a,b] 上的等距节点组,其中 X n = a + n h X_n = a + nh Xn=a+nh h h h 是步长, n n n 是节点索引,这样的离散化有助于数值计算。
  2. 向前差商近似微商:

    • 在节点 X n X_n Xn 处,通过向前差商 y ( X n + 1 ) − y ( X n ) h \frac{y(X_{n+1}) - y(X_n)}{h} hy(Xn+1)y(Xn) 近似替代微分方程 y ′ ( x ) = f ( x , y ( x ) ) y'(x) = f(x, y(x)) y(x)=f(x,y(x)) 中的导数项,得到 y ′ ( X n ) ≈ y ( X n + 1 ) − y ( X n ) h = f ( X n , y ( X n ) ) y'(X_n) \approx \frac{y(X_{n+1}) - y(X_n)}{h} = f(X_n, y(X_n)) y(Xn)hy(Xn+1)y(Xn)=f(Xn,y(Xn))
    • 这个近似通过将差商等于导数的思想,将微分方程转化为递推关系式。
  3. 递推公式:

    • 将上述近似公式改为等式,得到递推公式 y n + 1 = y n + h f ( X n , y n ) y_{n+1} = y_n + hf(X_n, y_n) yn+1=yn+hf(Xn,yn)
    • 这个公式是 Euler 方法的核心,通过这个公式可以逐步计算得到近似解的数值。
  4. 步骤解释:

    • n = 0 n=0 n=0 时,使用初始条件 y 0 y_0 y0 计算 y 1 y_1 y1
    • 然后,利用 y 1 y_1 y1 计算 y 2 y_2 y2,以此类推,得到 y n y_n yn,直到 n = N n=N n=N,其中 N N N 是节点数。
    • 这个过程形成了一个逐步逼近微分方程解的序列。
  5. 几何解释:

    • 在几何上,Euler 方法的求解过程可以解释为在积分曲线上通过连接相邻点的折线来逼近微分方程的解,因而被称为折线法
  6. 截断误差:

    • 通过 Taylor 展开,可以得到 Euler 方法的截断误差公式(忽略 h 2 h^2 h2 项) y ( x n + 1 ) = y ( x n ) + h f ( X n , y n ) + O ( h 2 ) y(x_{n+1}) = y(x_n) + hf(X_n, y_n) + O(h^2) y(xn+1)=y(xn)+hf(Xn,yn)+O(h2)
    • 这表明 Euler 方法的误差主要来自于 h h h 的一阶项,因此选择较小的步长可以提高方法的精度。

b. 典例解析

在这里插入图片描述

计算过程:

  1. 初始化: x 0 = 0 x_0 = 0 x0=0, y 0 = 1 y_0 = 1 y0=1.
  2. 计算 x 1 x_1 x1 y 1 y_1 y1
    x 1 = x 0 + h = 0.1 x_1 = x_0+h=0.1 x1=x0+h=0.1 y 1 = y 0 + h f ( x 0 , y 0 ) = 1 + 0.1 ⋅ ( y 0 − 2 x 0 y 0 ) = 1 + 0.1 ⋅ 1 = 1.1 y_1 = y_0 + h f(x_0, y_0) = 1 + 0.1 \cdot (y_0-\frac{2x_0}{y_0}) = 1 + 0.1 \cdot 1 = 1.1 y1=y0+hf(x0,y0)=1+0.1(y0y02x0)=1+0.11=1.1.
  3. 计算 x 2 x_2 x2 y 2 y_2 y2
    x 2 = x 1 + h = 0.2 x_2 = x_1+h=0.2 x2=x1+h=0.2 y 2 = y 1 + h f ( x 1 , y 1 ) = 1.1 + 0.1 ⋅ ( y 1 − 2 x 1 y 1 ) = 1.1 + 0.1 ⋅ ( 1.1 − 0.181819 ) = 1.191818 y_2 = y_1 + h f(x_1, y_1) = 1.1 + 0.1 \cdot (y_1-\frac{2x_1}{y_1}) = 1.1 + 0.1 \cdot (1.1-0.181819)= 1.191818 y2=y1+hf(x1,y1)=1.1+0.1(y1y12x1)=1.1+0.1(1.10.181819)=1.191818.
  4. 计算 x 3 x_3 x3 y 3 y_3 y3
    … … … … … … ……………… ………………

c. 算法实现

import numpy as np
import matplotlib.pyplot as pltdef forward_euler(f, y0, a, b, h):"""使用向前欧拉法求解一阶常微分方程初值问题Parameters:- f: 函数,表示微分方程的右侧项,形式为 f(x, y)- y0: 初始条件,表示在 x=a 处的函数值- a: 区间起点- b: 区间终点- h: 步长Returns:- x_values: 区间 [a, b] 上的离散节点- y_values: 对应节点上的函数值的近似解"""num_steps = int((b - a) / h) + 1  # 计算步数x_values = np.linspace(a, b, num_steps)  # 生成离散节点y_values = np.zeros(num_steps)  # 初始化结果数组y_values[0] = y0  # 设置初始条件# 使用向前欧拉法进行逐步迭代for i in range(1, num_steps):x = x_values[i - 1]y = y_values[i - 1]y_values[i] = y + h * f(x, y)return x_values, y_values# 示例:求解 y' = y - 2x/y,初始条件 y(0) = 1 在区间 [0, 1] 上的近似解
def example_function(x, y):return y - 2*x/ya, b = 0, 1  # 区间 [a, b]
y0 = 1  # 初始条件 y(0) = 1
h = 0.1  # 步长x_values, y_values = forward_euler(example_function, y0, a, b, h)# 绘制结果
plt.plot(x_values, y_values, label='Forward Euler')
plt.plot(x_values, np.sqrt(1+2*x_values), label='Exact Solution')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.show()

在这里插入图片描述
在这里插入图片描述

这篇关于【数值计算方法(黄明游)】常微分方程初值问题的数值积分法:欧拉方法(向前Euler)【理论到程序】的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/432289

相关文章

Springboot控制反转与Bean对象的方法

《Springboot控制反转与Bean对象的方法》文章介绍了SpringBoot中的控制反转(IoC)概念,描述了IoC容器如何管理Bean的生命周期和依赖关系,它详细讲解了Bean的注册过程,包括... 目录1 控制反转1.1 什么是控制反转1.2 SpringBoot中的控制反转2 Ioc容器对Bea

C++实现回文串判断的两种高效方法

《C++实现回文串判断的两种高效方法》文章介绍了两种判断回文串的方法:解法一通过创建新字符串来处理,解法二在原字符串上直接筛选判断,两种方法都使用了双指针法,文中通过代码示例讲解的非常详细,需要的朋友... 目录一、问题描述示例二、解法一:将字母数字连接到新的 string思路代码实现代码解释复杂度分析三、

mysql8.0无备份通过idb文件恢复数据的方法、idb文件修复和tablespace id不一致处理

《mysql8.0无备份通过idb文件恢复数据的方法、idb文件修复和tablespaceid不一致处理》文章描述了公司服务器断电后数据库故障的过程,作者通过查看错误日志、重新初始化数据目录、恢复备... 周末突然接到一位一年多没联系的妹妹打来电话,“刘哥,快来救救我”,我脑海瞬间冒出妙瓦底,电信火苲马扁.

SpringBoot使用Jasypt对YML文件配置内容加密的方法(数据库密码加密)

《SpringBoot使用Jasypt对YML文件配置内容加密的方法(数据库密码加密)》本文介绍了如何在SpringBoot项目中使用Jasypt对application.yml文件中的敏感信息(如数... 目录SpringBoot使用Jasypt对YML文件配置内容进行加密(例:数据库密码加密)前言一、J

Spring Boot 中正确地在异步线程中使用 HttpServletRequest的方法

《SpringBoot中正确地在异步线程中使用HttpServletRequest的方法》文章讨论了在SpringBoot中如何在异步线程中正确使用HttpServletRequest的问题,... 目录前言一、问题的来源:为什么异步线程中无法访问 HttpServletRequest?1. 请求上下文与线

解读为什么@Autowired在属性上被警告,在setter方法上不被警告问题

《解读为什么@Autowired在属性上被警告,在setter方法上不被警告问题》在Spring开发中,@Autowired注解常用于实现依赖注入,它可以应用于类的属性、构造器或setter方法上,然... 目录1. 为什么 @Autowired 在属性上被警告?1.1 隐式依赖注入1.2 IDE 的警告:

SpringBoot快速接入OpenAI大模型的方法(JDK8)

《SpringBoot快速接入OpenAI大模型的方法(JDK8)》本文介绍了如何使用AI4J快速接入OpenAI大模型,并展示了如何实现流式与非流式的输出,以及对函数调用的使用,AI4J支持JDK8... 目录使用AI4J快速接入OpenAI大模型介绍AI4J-github快速使用创建SpringBoot

Android开发中gradle下载缓慢的问题级解决方法

《Android开发中gradle下载缓慢的问题级解决方法》本文介绍了解决Android开发中Gradle下载缓慢问题的几种方法,本文给大家介绍的非常详细,感兴趣的朋友跟随小编一起看看吧... 目录一、网络环境优化二、Gradle版本与配置优化三、其他优化措施针对android开发中Gradle下载缓慢的问

python 3.8 的anaconda下载方法

《python3.8的anaconda下载方法》本文详细介绍了如何下载和安装带有Python3.8的Anaconda发行版,包括Anaconda简介、下载步骤、安装指南以及验证安装结果,此外,还介... 目录python3.8 版本的 Anaconda 下载与安装指南一、Anaconda 简介二、下载 An

Java中将异步调用转为同步的五种实现方法

《Java中将异步调用转为同步的五种实现方法》本文介绍了将异步调用转为同步阻塞模式的五种方法:wait/notify、ReentrantLock+Condition、Future、CountDownL... 目录异步与同步的核心区别方法一:使用wait/notify + synchronized代码示例关键