scipy 对于给定的微分方程，如何求解库兹涅茨曲线增长模型？

flvlnr44 于 2022-11-10 发布在其他

关注(0)|答案(1)|浏览(186)

内容：

使用麦迪逊项目的真实的人均GDP数据集，我用最小二乘法推导出了以下方程：0.012406 + 0.005132ln（克）- 0.006304ln（克）*2
当我试图预测不同经济群体到2050年的人均GDP时，我参考了本文的方法“Tilman et al. 10.1073/pnas.1116437108”，以了解他们是如何将其作为微分方程来求解的：微分G/dt = G（-0.6284 +0.157lnG- 0.0093ln（G））^2
我用同样的方法转换了最小二乘法：-0.012406 克+克 0.005132 数学对数（克）-克 0.006304 数学对数（克）**2
我尝试在python中求解常微分方程的几个初始值，以获得库兹涅茨曲线和2050年的估计值。我使用了下面的代码，但我无法解决相同的问题。

Python代码：

import numpy as np
import matplotlib.pyplot as plt
import scipy as sp
from scipy.integrate import odeint
from scipy.integrate import solve_ivp 
import math

def solveit(y0):
    def gdp(g, t):
        y = g
        dgdt = [-0.012406*g  + g*0.005132*math.log(g) - g*0.006304*math.log(g)**2]
        return dgdt

# initial conditions

    #y0 = [785.60] 
    t = np.linspace(0, 60000, 1000)

# call integrator

    sol = odeint(gdp, y0, t)
    m = sol[:]
    plt.plot(t,m)
    plt.show()

ys= [[785.60],[1860],[7800]]

fig = plt.figure()
for y_ in ys:
    solveit(y_)

plt.legend(loc='best')
plt.grid()
plt.show()

错误消息：

RuntimeError: The array return by func must be one-dimensional, but got ndim=2.

关于同样的方向将是有帮助的。

scipy

来源：https://stackoverflow.com/questions/72776455/how-to-solve-kuznets-curve-growth-model-for-the-given-differential-equation

1条答案

按热度按时间

quhf5bfb1#

当odeint调用用户提供的func参数时，它将为状态向量y传入一个数组。即使系统是一个单标量方程，y也将作为长度为1的一维数组传入。在您的例子中，这意味着g将是长度为1的一维数组。您在gdp中有以下代码：

dgdt = [-0.012406*g  + g*0.005132*math.log(g) - g*0.006304*math.log(g)**2]
    return dgdt

问题是dgdt表达式中的多余括号。g已经是一维数组，因此-0.012406*g也是。结果是math.log(g)实际上返回一个标量，但是标量加上长度为1的一维NumPy会产生预期的效果，所以在这种情况下它不会引起问题。所以完整表达式-0.012406*g + g*0.005132*math.log(g) - g*0.006304*math.log(g)**2的结果也是一个一维NumPy数组（长度为1）。当您将结果括在括号中时，您增加了数据结构的“深度”。实际上，您创建了一个二维数组（具有平凡的形状（1，1））。例如，在下面的示例中，w是长度为1的一维数组：

In [23]: w = np.array([1.0])

In [24]: np.shape(w)
Out[24]: (1,)

In [25]: np.shape([w])
Out[25]: (1, 1)

请注意，[w]的形状为(1, 1)，即它是一个 2-d 数组。odeint期望用户的函数返回一个一维数组，但它从您的函数中获得了一个二维数组，并引发了错误。
解决方法很简单：删除这些括号：

def gdp(g, t):
        dgdt = -0.012406*g  + g*0.005132*math.log(g) - g*0.006304*math.log(g)**2
        return dgdt

(Note我删除了行y = g，这没有任何意义。）

赞(0）回复(0）举报 2022-11-10

我来回答

scipy 对于给定的微分方程，如何求解库兹涅茨曲线增长模型？

内容：

Python代码：

错误消息：

1条答案

相关问题

热门标签

最新问答