如果代数操作失败,可以对约束进行数值求解,例如在每个时间步运行fsolve:import sys

from numpy import linspace

from scipy.integrate import odeint

from scipy.optimize import fsolve

y0 = [0, 5]

time = linspace(0., 10., 1000)

F_lon = 10.

mass = 1000.

def F_r(a, v):

return (((1 - a) / 3) ** 2 + (2 * (1 + a) / 3) ** 2) * v

def constraint(a, v):

return (F_lon - F_r(a, v)) / mass - a

def integral(y, _):

v = y[1]

a, _, ier, mesg = fsolve(constraint, 0, args=[v, ], full_output=True)

if ier != 1:

print "I coudn't solve the algebraic constraint, error:\n\n", mesg

sys.stdout.flush()

return [v, a]

dydt = odeint(integral, y0, time)

显然这会减慢你的时间整合。请始终检查fsolve是否找到了一个好的解决方案,并刷新输出,这样您就可以实现它,并停止模拟。在

关于如何在前一个时间步“缓存”变量的值,可以利用以下事实:默认参数只在函数定义中计算

^{pr2}$

请注意,为了使技巧生效,cache参数必须是可变的,这就是我使用list的原因。如果您不熟悉默认参数的工作方式,请参阅this链接。在

请注意,这两个代码不会产生相同的结果,您应该非常小心地使用上一个时间步骤中的值,以确保数值稳定性和精度。第二种方法显然要快得多。在

python解隐式方程_求解隐式ODE(微分代数方程DAE)相关推荐

  1. python解常微分方程龙格库_求解二阶常微分方程的RungeKutta四阶方法

    我试着做一个简谐振子的例子,它将用龙格-库塔四阶法求解.要求解的二阶常微分方程(ODE)和初始条件为: y''+y=0 y(0)=0和y'(0)=1/pi 范围在0到1之间,共有100步.我用u作为辅 ...

  2. python解常微分方程龙格库_求解常微分方程组初值问题的龙格库塔法分析及其C代码...

    求解常微分方程组初值问题的 龙格库塔法分析及其 C 代码 1 .概 述 由高等数学的知识可知,一些特殊类型的常微分方程(组)能够求出给定初 始值的解析解, 而在科学与工程问题中遇到的常微分方程 (组) ...

  3. 2021-01-13 Matlab求解微分代数方程 (DAE)

    Matlab求解微分代数方程 (DAE) 什么是微分代数方程? 微分代数方程是一类微分方程,其中一个或多个因变量导数未出现在方程中.方程中出现的未包含其导数的变量称为代数变量,代数变量的存在意味着不能 ...

  4. matlab求解微分代数方程组,微分代数方程(DAE)的Matlab 解法.PDF

    微分代数方程(DAE)的Matlab 解法 微分代数方程(DAE)的Matlab解法 所谓微分代数方程,是指在微分方程中,某些变量满足某些代数方程的约束.假 设微分方程的更一般形式可以写成 前面所介绍 ...

  5. python解奥数题_【新手】小白用python解小学奥数第五集

    前言 因为工作的关系,长期教小学数学奥数,最近接触到python,学到了while:for 语句这里,于是就有了这个想法:用python解奥数题 1. 订立制作目标目标: 打印题目 输入逻辑 输出答案 ...

  6. 抛物线交点式公式_抛物线交点式

    以下是范文网www zhuodaoren com 分享的6 已知抛物线与轴两交点在轴同侧,它们的距离的平方等于,则的值为( ),希望能帮助到大家!(一)6 已知抛物线与轴两交点在轴同侧,它们的距离的平 ...

  7. 交替隐式方法matlab求解,交替方向隐式CFD解法器的GPU并行计算及其优化

    0引言GPU作为一类典型的众核计算加速部件,已被广泛应用于构建高性能计算系统[1].充分发挥GPU性能,需要针对问题模型和计算方法特点,设计细粒度并行算法.对计算流体力学(ComputationalF ...

  8. python解一元二次方程复数_一元二次方程求解(包括复数各种情况)

    /*无聊时候写的,我想每个入门的编程都会拿着个练手把,多的不说,我们来谈谈正经的--此功能实现一元二次方程求解,复数情况,输入字符处理判断是否为数字*/ import os import math i ...

  9. python解三元一次方程组_关于实现求解三元一次方程组的小程序的一点小研究

    引入 在今年的暑假,作者利用空闲时间阅读了gilbert strange的introduction to linear algebra.作为享誉国际学界多年的MIT教材,该书以高斯消元法求解二三元方程 ...

最新文章

  1. AWS计算云上下载代码命令 git clone url
  2. 基于WinSvr2012共享文件夹的Hyper-V实时迁移之三实时迁移的实现及验证
  3. redis单机版安装
  4. php怎么看回调的异步通知的数据_php app支付宝回调(异步通知)详解
  5. wxWidgets:wxStatusBar类用法
  6. 【字符集UTF8】处理Toad显示乱码及Windows XP下无法插入“某些汉字”问题
  7. 《统计学习方法》—— 朴素贝叶斯方法、详细推导及其python3实现(二)
  8. 电力-101规约说明书2
  9. 电脑qq音乐显示无法代理服务器,电脑QQ音乐软件无法登录如何解决
  10. 四级单词pdf_英语单词里的字母到底有没有含义?
  11. c语言中x的n次幂程序,请教c语言关于X的N次幂问题?
  12. 自我介绍(思维导图)
  13. 用计算机管理从新分区,电脑如何分区硬盘分区_电脑怎么重新分区教程-win7之家...
  14. java8 Predicate
  15. 如何在win7下装ubuntu(硬盘版安装)
  16. Django项目之美多商城遇到的问题记录与解决(一)
  17. 写小论文之引言写什么?
  18. Node.js微服务 2 :基于Seneca和PM2构建Node.js微服务
  19. android 一天多少毫秒,Android中的日期和时间以毫秒为单位
  20. 关于 百度飞浆paddleOCR编译32位版本 的解决方案

热门文章

  1. 数据中心 CLOS 架构
  2. 联想笔记本thinkpad按F2不能直接重命名
  3. 完美国际最新服务器,完美国际2:新服“永恒”明日开启,新老朋友踏实安家...
  4. 树莓派开发文档笔记—基础配置
  5. 南加州大狼屋水上乐园部署RFID腕带系统
  6. 利用ffmpeg加水印
  7. 上海PHICOMM斐讯短短五年时间市场为何能扩大成这样?
  8. 厉害了!通过 Java 技术手段,获取女朋友定位地址...
  9. Bugku-加密-Crack it
  10. 实现类似路由器那样通过WEB页面配置设备