龙格库塔解一阶微分方程c语言,四阶龙格库塔法解微分方程
四阶龙格库塔法解微分方程
精品资料 欢迎下载 四阶龙格库塔法解微分方程 一、四阶龙格库塔法解一阶微分方程 ①一阶微分方程: ,初始值y(0)=0,求解区间[0 10]。 MATLAB程序: %%%%%%%%%%% 四阶龙哥库塔法解一阶微分方程 %%%%%%%%%%% y =cost %%%%%%%%%%% y(0)=0, 0≤t≤10,h=0.01 %%%%%%%%%%% y=sint h=0.01; hf=10; t=0:h:hf; y=zeros(1,length(t)); y(1)=0; F=@(t,y)(cos(t)); for i=1:(length(t)-1) k1=F(t(i),y(i)); k2=F(t(i)+h/2,y(i)+k1*h/2); k3=F(t(i)+h/2,y(i)+k2*h/2); k4=F(t(i)+h,y(i)+k3*h); y(i+1)=y(i)+1/6*(k1+2*k2+2*k3+k4)*h; end subplot(211) plot(t,y, - ) xlabel( t ); ylabel( y ); title( Approximation ); span=[0,10]; p=y(1); [t1,y1]=ode45(F,span,p); subplot(212) plot(t1,y1) xlabel( t ); ylabel( y ); title( Exact ); 图1 ②一阶微分方程: ,初始值x(1)=2,求解区间[1 3]。 MATLAB程序: %%%%%%%%%%% 四阶龙哥库塔法解微分方程 %%%%%%%%%%% x (t)=(t*x-x^2)/t^2 %%%%%%%%%%% x(1)=2, 1≤t≤3, h=1/128 %%%%%%%%%%% 精确解:x(t)=t/(0.5+lnt) h=1/128; %%%%% 步长 tf=3; t=1:h:tf; x=zeros(1,length(t)); x(1)=2; %%%%% 初始值 F_tx=@(t,x)(t.*x-x.^2)./t.^2; for i=1:(length(t)-1) k_1=F_tx(t(i),x(i)); k_2=F_tx(t(i)+0.5*h,x(i)+0.5*h*k_1); k_3=F_tx((t(i)+0.5*h),(x(i)+0.5*h*k_2)); k_4=F_tx((t(i)+h),(x(i)+k_3*h)); x(i+1)=x(i)+(1/6)*(k_1+2*k_2+2*k_3+k_4)*h; end subplot(211) plot(t,x, - ); xlabel( t ); ylabel( x ); legend( Approximation ); %%%%%%%%%%%%%%%%%%%%%%%%%%%% ode45求精确解 t0=t(1);x0=x(1); xspan=[t0 tf]; [x_ode45,y_ode45]=ode45(F_tx,xspan,x0); subplot(212) plot(x_ode45,y_ode45, -- ); xlabel( t ); ylabel( x ); legend( Exact ); 图2 二、四阶龙格库塔法解二阶微分方程 ①二阶微分方程: ,初始值y(0)=0,y (0)=-1,求解区间[0 10]。 MATLAB程序: %%%%%%%%% 龙格库塔欧拉方法解二阶微分方程 %%%%%%%%% y =cost %%%%%%%%% y (0)=-1, y(0)=0 %%%%%%%%% 转换为一阶微分方程 h=0.01; ht=10; t=0:h:ht; %%%%%%%%% y =z1, z1 =y =cost y=zeros(1,length(t)); z=zeros(1,length(t)); y(1)=0; z(1)=-1; f=@(t,y,z)(z); g=@(t,y,z)(sin(t)); for i=1:(length(t)-1) K1=f(t(i),y(i),z(i)); L1=g(t(i),y(i),z(i)); K2=f(t(i)+h/2,y(i)+h/2*K1,z(i)+h/2*L1); L2=g(t(i)+h/2,y(i)+h/2*K1,z(i)+h/2*L1); K3=f(t(i)+h/2,y(i)+h/2*K2,z(i)+h/2*L2); L3=g(t(i)+h/2,y(i)+h/2*K2,z(i)+h/2*L2); K4=f(t(i)+h,y(i)+h*K3,z(i)+h*L3); L4=g(t(i)+h,y(i)+h*K3,z(i)+h*L3); y(i+1)=y(i)+h/6*(K1+2*K2+2*K3+K4); z(i+1)=z(i)+h/6*(L1+2*L2+2*L3+L4); end plot(t,y) xlabel( t ); ylabel( y ); title( y =sin(t) ); legend( y =sint ); 图3 ②二阶微分方程: ,初始值y(1)=0,y (1)=0,求解区间[0 25]。 MATLAB程序: %%%%%%%%% 龙格库塔欧拉方法解二阶微分方程 %%%%%%%%% y =7.35499*cosx %set(0, RecursionLimit ,500) h=0.01; a=0; b=25; x=a:h:b; RK_y(1)=0; %初值 RK_z(1)=0; for i=1:length(x)-1 K1=RK_z(i); L1=rightf_sys1(x(i),RK_y(i),RK_z(i)); % K1 and L1 K2=RK_z(i)+0.5*h*L1; L2=rightf_sys1(x(i)+0.5*h,RK_y(i)+0.5*h*K1,RK_z(i)+0.5*h*L1); K3=RK_z(i)+0.5*h*L2; L3=rightf_sys1(x(i)+0.5*h,RK_y(i)+0.5*h*K2,RK_z(i)+0.5*h*L2); K4=RK_z(i)+h*L3; L4=rightf_sys1(x(i)+h,RK_y(i)+h*K3,RK_z(i)+h*L3); % K4 and L4 RK_
龙格库塔解一阶微分方程c语言,四阶龙格库塔法解微分方程相关推荐
- 一阶欧拉近似matlab,MATLAB改进欧拉法与四阶龙格-库塔求解一阶常微分方程.doc
MATLAB改进欧拉法与四阶龙格-库塔求解一阶常微分方程 姓名:樊元君 学号:2012200902 日期:2012.11.06 一.实验目的 掌握改进欧拉法与四阶龙格-库塔求解一阶常微分方程的初值问题 ...
- 四阶龙格库塔法的基本思想_经典四阶龙格库塔法解一阶微分方程组讲义.doc
1.经典四阶龙格库塔法解一阶微分方程组 1.1运用四阶龙格库塔法解一阶微分方程组算法分析 , 经过循环计算由 推得 -- 每个龙格-库塔方法都是由一个合适的泰勒方法推导而来,使得其最终全局误差为,一种 ...
- 二阶偏微分方程组 龙格库塔法_1、经典四阶龙格库塔法解一阶微分方程组
1.经典四阶龙格库塔法解一阶微分方程组 陕 西 科 技 大 学 数值计算课程设计任务书 理学院信息与计算科学/应用数学专业信息08/数学08 班级 学生: 题目:典型数值算法的C++语言程序设计 课程 ...
- 四阶龙格库塔法的基本思想_请问用四阶龙格库塔法解二阶微分方程的思想是什么?...
默认y的单位是弧度 k=1000; t=0:0.001:1; Y=[]; err=1 K=[]; Ymax=[]; xishu=1.01; while err X=[0 0]; k=xishu*k; ...
- [计算机数值分析]四阶龙格-库塔经典格式解常微分方程的初值问题
龙格-库塔方法的设计思想: 四阶龙格-库塔方法的经典格式: 程序设计框图: 例:设取步长h=0.2,从x=0到x=1用四阶经典格式解决以下常微分方程的初值问题. 运行示例: 源码: #include& ...
- 龙格库塔公式法解微分方程组初值问题实例
用四阶龙格库塔公式(取h=0.1)解下列微分方程组初值问题: { y 1 ′ = 120 − 2 y 1 + 2 y 2 ( 0 ≤ x ≤ 1.0 ) , y 2 ′ = 2 y 1 − 5 y 2 ...
- Ode45以及龙格-库塔算法
Ode45及龙格-库塔算法 ODE45 ODE45函数描述: 求解非刚性微分方程 - 中阶方 ODE是Matlab专门用于解微分方程的功能函数.该求解器有变步长(variable-step)和定步长( ...
- 【算法基础学习 6】龙格库塔法 求微分方程
龙格库塔法 龙格库塔法对于一个工科学生来说估计并不陌生,因为我们做项目的时候或多或少会接触到微分方程,解微分方程的数值解是就需要用到龙格库塔方程了.我第一次接触到该方法是对四轴姿态角进行四元数解算的时 ...
- 龙格库塔法求微分方程
龙格库塔法 龙格库塔法对于一个工科学生来说估计并不陌生,因为我们做项目的时候或多或少会接触到微分方程,解微分方程的数值解是就需要用到龙格库塔方程了.我第一次接触到该方法是对四轴姿态角进行四元数解算的时 ...
- 阿当姆斯方法求解微分方程初值问题:四阶龙格库塔提供出发值(C语言)
#include<stdio.h>float f1(float x,float y) //第一问 {return y*y;//方程 }float f2(float x,float y)// ...
最新文章
- ui产品小结 - 包含小程序 前端等
- java获得项目绝对路径
- 每日一皮:努力寻找Bug的程序员
- FZU - 2218 Simple String Problem(状压dp)
- 我的世界java怎么玩起床战争_我的世界怎么玩起床战争_我的世界起床战争怎么玩_52pk单机游戏...
- LeetCode3:Longest Substring Without Repeating Characters
- 前端学习(3311):redux的state hook对象
- vasp软件全名是什么_Materials Studio软件常见问题与解答
- mysql查询不超过19_mysql45讲 19.为什么我只查一行的语句,也执行这么慢?
- FinSpy 发布 Mac 和 Linux OS 版本攻击埃及组织机构
- linux编译c文件for循环,Linux C 循环队列的实现
- Python实现对给定的列表中连续数字的寻找
- chrome 模拟点击_详解爬虫模拟登陆的三种方法
- atitit.团队建设总结o6o fix
- GP数据库初始化失败定位
- 一个新手对软件开发的理解(写自第一个项目--Linpop之后)
- 云原生数据中台:架构、方法论与实践
- 掷骰子java程序_掷骰子游戏窗体实现--Java初级小项目
- 电脑开机密码忘记了怎么办
- java教材管理系统,基于web的教材管理系统