如何使用python计算常微分方程？

常用形式

创新互联公司主营平舆网站建设的网络公司,主营网站建设方案,重庆APP开发,平舆h5小程序开发搭建,平舆网站营销推广欢迎平舆等地区企业咨询

odeint(func, y0, t,args,Dfun)

一般这种形式就够用了。

下面是官方的例子，求解的是

D(D(y1))-t*y1=0

为了方便，采取D=d/dt。如果我们令初值

y1(0) = 1.0/3**(2.0/3.0)/gamma(2.0/3.0)

D(y1)(0) = -1.0/3**(1.0/3.0)/gamma(1.0/3.0)

这个微分方程的解y1=airy(t)。

令D(y1)=y0，就有这个常微分方程组。

D(y0)=t*y1

D(y1)=y0

Python求解该微分方程。

from scipy.integrate import odeint

from scipy.special import gamma, airy

y1_0 = 1.0/3**(2.0/3.0)/gamma(2.0/3.0)

y0_0 = -1.0/3**(1.0/3.0)/gamma(1.0/3.0)

y0 = [y0_0, y1_0]

def func(y, t):

... return [t*y[1],y[0]]

def gradient(y,t):

... return [[0,t],[1,0]]

x = arange(0,4.0, 0.01)

t = x

ychk = airy(x)[0]

y = odeint(func, y0, t)

y2 = odeint(func, y0, t, Dfun=gradient)

print ychk[:36:6]

[ 0.355028 0.339511 0.324068 0.308763 0.293658 0.278806]

print y[:36:6,1]

[ 0.355028 0.339511 0.324067 0.308763 0.293658 0.278806]

print y2[:36:6,1]

[ 0.355028 0.339511 0.324067 0.308763 0.293658 0.278806]

得到的解与精确值相比，误差相当小。

=======================================================================================================

args是额外的参数。

用法请参看下面的例子。这是一个洛仑兹曲线的求解，并且用matplotlib绘出空间曲线图。（来自《python科学计算》）

from scipy.integrate import odeint

import numpy as np

def lorenz(w, t, p, r, b):

# 给出位置矢量w，和三个参数p, r, b 计算出

# dx/dt, dy/dt, dz/dt 的值

x, y, z = w

# 直接与lorenz 的计算公式对应

return np.array([p*(y-x), x*(r-z)-y, x*y-b*z])

t = np.arange(0, 30, 0.01) # 创建时间点

# 调用ode 对lorenz 进行求解, 用两个不同的初始值

track1 = odeint(lorenz, (0.0, 1.00, 0.0), t, args=(10.0, 28.0, 3.0))

track2 = odeint(lorenz, (0.0, 1.01, 0.0), t, args=(10.0, 28.0, 3.0))

# 绘图

from mpl_toolkits.mplot3d import Axes3D

import matplotlib.pyplot as plt

fig = plt.figure()

ax = Axes3D(fig)

ax.plot(track1[:,0], track1[:,1], track1[:,2])

ax.plot(track2[:,0], track2[:,1], track2[:,2])

plt.show()

===========================================================================

scipy.integrate.odeint(func, y0, t, args=(), Dfun=None, col_deriv=0, full_output=0, ml=None, mu=None, rtol=None, atol=None, tcrit=None, h0=0.0, hmax=0.0, hmin=0.0, ixpr=0, mxstep=0, mxhnil=0, mxordn=12, mxords=5, printmessg=0)

计算常微分方程（组）

使用 FORTRAN库odepack中的lsoda解常微分方程。这个函数一般求解初值问题。

参数：

func : callable(y, t0, ...) 计算y在t0 处的导数。

y0 : 数组 y的初值条件（可以是矢量）

t : 数组 为求出y，这是一个时间点的序列。初值点应该是这个序列的第一个元素。

args : 元组 func的额外参数

Dfun : callable(y, t0, ...) 函数的梯度(Jacobian)。即雅可比多项式。

col_deriv : boolean. True，Dfun定义列向导数（更快），否则Dfun会定义横排导数

full_output : boolean 可选输出，如果为True 则返回一个字典，作为第二输出。

printmessg : boolean 是否打印convergence 消息。

返回： y : array, shape (len(y0), len(t))

数组，包含y值，每一个对应于时间序列中的t。初值y0 在第一排。

infodict : 字典，只有full_output == True 时，才会返回。

字典包含额为的输出信息。

键值：

‘hu’ vector of step sizes successfully used for each time step.

‘tcur’ vector with the value of t reached for each time step. (will always be at least as large as the input times).

‘tolsf’ vector of tolerance scale factors, greater than 1.0, computed when a request for too much accuracy was detected.

‘tsw’ value of t at the time of the last method switch (given for each time step)

‘nst’ cumulative number of time steps

‘nfe’ cumulative number of function evaluations for each time step

‘nje’ cumulative number of jacobian evaluations for each time step

‘nqu’ a vector of method orders for each successful step.

‘imxer’index of the component of largest magnitude in the weighted local error vector (e / ewt) on an error return, -1 otherwise.

‘lenrw’ the length of the double work array required.

‘leniw’ the length of integer work array required.

‘mused’a vector of method indicators for each successful time step: 1: adams (nonstiff), 2: bdf (stiff)

其他参数，官方网站和文档都没有明确说明。相关的资料，暂时也找不到。

Python怎么做最优化

一、概观

scipy中的optimize子包中提供了常用的最优化算法函数实现。我们可以直接调用这些函数完成我们的优化问题。optimize中函数最典型的特点就是能够从函数名称上看出是使用了什么算法。下面optimize包中函数的概览：

1.非线性最优化

fmin -- 简单Nelder-Mead算法

fmin_powell -- 改进型Powell法

fmin_bfgs -- 拟Newton法

fmin_cg -- 非线性共轭梯度法

fmin_ncg -- 线性搜索Newton共轭梯度法

leastsq -- 最小二乘

2.有约束的多元函数问题

fmin_l_bfgs_b ---使用L-BFGS-B算法

fmin_tnc ---梯度信息

fmin_cobyla ---线性逼近

fmin_slsqp ---序列最小二乘法

nnls ---解|| Ax - b ||_2 for x=0

3.全局优化

anneal ---模拟退火算法

brute --强力法

4.标量函数

fminbound

brent

golden

bracket

5.拟合

curve_fit-- 使用非线性最小二乘法拟合

6.标量函数求根

brentq ---classic Brent (1973)

brenth ---A variation on the classic Brent（1980）ridder ---Ridder是提出这个算法的人名

bisect ---二分法

newton ---牛顿法

fixed_point

7.多维函数求根

fsolve ---通用

broyden1 ---Broyden’s first Jacobian approximation.

broyden2 ---Broyden’s second Jacobian approximationnewton_krylov ---Krylov approximation for inverse Jacobiananderson ---extended Anderson mixing

excitingmixing ---tuned diagonal Jacobian approximationlinearmixing ---scalar Jacobian approximationdiagbroyden ---diagonal Broyden Jacobian approximation8.实用函数

line_search ---找到满足强Wolfe的alpha值

check_grad ---通过和前向有限差分逼近比较检查梯度函数的正确性二、实战非线性最优化

fmin完整的调用形式是：

fmin(func, x0, args=(), xtol=0.0001, ftol=0.0001, maxiter=None, maxfun=None, full_output=0, disp=1, retall=0, callback=None)不过我们最常使用的就是前两个参数。一个描述优化问题的函数以及初值。后面的那些参数我们也很容易理解。如果您能用到，请自己研究。下面研究一个最简单的问题，来感受这个函数的使用方法：f(x)=x**2-4*x+8，我们知道，这个函数的最小值是4，在x=2的时候取到。

from scipy.optimize import fmin #引入优化包def myfunc(x):

return x**2-4*x+8 #定义函数

x0 = [1.3] #猜一个初值

xopt = fmin(myfunc, x0) #求解

print xopt #打印结果

运行之后，给出的结果是：

Optimization terminated successfully.

Current function value: 4.000000

Iterations: 16

Function evaluations: 32

[ 2.00001953]

程序准确的计算得出了最小值，不过最小值点并不是严格的2，这应该是由二进制机器编码误差造成的。

除了fmin_ncg必须提供梯度信息外，其他几个函数的调用大同小异，完全类似。我们不妨做一个对比：

from scipy.optimize import fmin,fmin_powell,fmin_bfgs,fmin_cgdef myfunc(x):

return x**2-4*x+8

x0 = [1.3]

xopt1 = fmin(myfunc, x0)

print xopt1

print

xopt2 = fmin_powell(myfunc, x0)

print xopt2

print

xopt3 = fmin_bfgs(myfunc, x0)

print xopt3

print

xopt4 = fmin_cg(myfunc,x0)

print xopt4

给出的结果是：

Optimization terminated successfully.

Current function value: 4.000000

Iterations: 16

Function evaluations: 32

[ 2.00001953]

Optimization terminated successfully.

Current function value: 4.000000

Iterations: 2

Function evaluations: 53

1.99999999997

Optimization terminated successfully.

Current function value: 4.000000

Iterations: 2

Function evaluations: 12

Gradient evaluations: 4

[ 2.00000001]

Optimization terminated successfully.

Current function value: 4.000000

Iterations: 2

Function evaluations: 15

Gradient evaluations: 5

[ 2.]

我们可以根据给出的消息直观的判断算法的执行情况。每一种算法数学上的问题，请自己看书学习。个人感觉，如果不是纯研究数学的工作，没必要搞清楚那些推导以及定理云云。不过，必须了解每一种算法的优劣以及能力所及。在使用的时候，不妨多种算法都使用一下，看看效果分别如何，同时，还可以互相印证算法失效的问题。

在from scipy.optimize import fmin之后，就可以使用help(fmin)来查看fmin的帮助信息了。帮助信息中没有例子，但是给出了每一个参数的含义说明，这是调用函数时候的最有价值参考。

有源码研究癖好的，或者当你需要改进这些已经实现的算法的时候，可能需要查看optimize中的每种算法的源代码。在这里：https:/ / github. com/scipy/scipy/blob/master/scipy/optimize/optimize.py聪明的你肯定发现了，顺着这个链接往上一级、再往上一级，你会找到scipy的几乎所有源码！

python如何实现求函数的在一个连续区间的最值？

如果函数是确定的，可以用导数的方法进行计算，但是如果函数是不确定的，就需要用优化的方法来处理了，比如常用的梯度上升法，模拟退火等，希望可以帮到你。

本文标题：python函数的梯度 python 阶跃函数
URL分享：http://kswjz.com/article/dogjscj.html