0

我正在使用 JiTCDDE 在修改后的 Oregantor 模型上解决 DDE。我遇到的问题在分岔点附近,它将返回负值。虽然我知道这些是数学上有效的解决方案,但 Oregantor 代表了一个化学系统。因此,否定答案对于真实系统来说是不现实的。有没有办法设置代码以在 <= 0 时返回变量的最小值。以下是我到目前为止的代码的主要部分。

def P1(k):
        return(
            ((H*y(k))/(k01+H*y(k)+kl*H*H*A))*phi
        )

def C(i,j):

    return(
        M1 * ( y(j,t-tau1)-y(i) ) 
      + M2 * ( y(j,t-tau2)-y(i) ) 
    )




MO4 = [
 k1*A*y(1)-k2*y(0)*y(1)+ k3*A*y(0)-2.0*k4*y(0)*y(0)-(y(0)-xsur)*kf,           #HBrO2
 -k1*A*y(1)-k2*y(0)*y(1)+f1*k5*y(2)-(y(1)-ysur)*kf+P1(3)+C(2,6),              #Bromide
 2*k3*A*y(0)-k5*y(2)+P1(3)+C(2,6),                                            #Cataylst
 k1*A*y(1)+2*k2*y(0)*y(1)+k4*y(0)*y(0)-k6*y(3)-(y(3)-vsur)*kf-P1(3)-C(2,6),   #BrMa
 k1*A*y(5)-k2*y(4)*y(5)+ k3*A*y(4)-2.0*k4*y(4)*y(4)-(y(4)-xsur)*kf,           #HBrO2
 -k1*A*y(5)-k2*y(4)*y(5)+f2*k5*y(6)-(y(5)-ysur)*kf+P1(7)+C(6,2),              #Bromide
 2*k3*A*y(4)-k5*y(6)+P1(7)+C(6,2),                                            #Cataylst
 k1*A*y(5)+2*k2*y(4)*y(5)+k4*y(4)*y(4)-k6*y(7)-(y(7)-vsur)*kf-P1(7)-C(6,2),   #BrMa
]

I = jitcdde(MO4)
I.set_integration_parameters(rtol=1e-7,atol=1e-7)
I.constant_past ([0,1.0e-6,0,0,1.0e-6,1.0e-6,1.0e-6,1.0e-6], time=0.0)
I.step_on_discontinuities(max_step=.00001)


data=[]
for time in times:
    data.append( I.integrate(time))
np.savetxt('peaks_%d.dat'%(i), data,)

data1=np.loadtxt('peaks_%d.dat'%(i),dtype = float,delimiter=' ',skiprows=200,usecols=(2,6)).T #,skiprows=80
plt.plot(data1[0],'r')
plt.plot(data1[1],'-.b')
plt.title( 'Catalyst ' )
plt.xlabel('time(sec)')
plt.ylabel('Amplitude')
plt.show()
print('DONE')
4

1 回答 1

0

在 JiTCDDE 中没有直接的方法来限制解的符号。(来源:我是作者;我会知道的。)

你需要做的是找出为什么你的解决方案变得消极。到目前为止,我可以想到三个可能的原因(如果不知道您的确切设置,我无法告诉您更多信息):

  • 数值噪声将解从数学上正确的正解踢到数学上正确的负解(但解之间的转换在数学上不正确)。在这种情况下,您应该能够通过减小绝对容差 ( atol) 或可能的最大步长 ( max_step) 来避免这种情况。另一种方法是在对数域中工作(另请参见类似问题的答案),即,您通过构造使零无法达到(此外,您更强烈地权衡接近零的错误)。

  • 您的模型允许在不应该发生的情况下发生这种情况。在这种情况下,如果您确实需要避免负值,请避免动态变量的负导数已经非常接近于零,方法是在它周围包装一个完全执行此操作的函数 - 通过使用 sigmoids 来实现此逻辑(另请参阅此我对另一个问题的回答)。抵制使用阶跃函数的诱惑,因为集成商不喜欢这样(进一步阅读)。

  • 您的初始条件是积极的,但不合理。在这种情况下,使用不同的初始条件。如果没有直接的方法来判断哪些初始条件是合理的,那么选择随机条件可能是合适的,直到找到一个好的条件。

最后,我建议您考虑一下如果解决方案变得消极,您希望发生什么,以及这是否真的会更好。

于 2018-11-12T15:19:16.510 回答