1

I have wrote this code to solve an equation , I know the behavior of this function has very rapid oscillations, when I RUN it gives bogus values for some "m[x]" and some "t"'s, with this error:

C:\Users\dani\anaconda3\lib\site-packages\scipy\integrate\odepack.py:247: ODEintWarning: Excess work done on this call (perhaps wrong Dfun type). Run with full_output = 1 to get quantitative information. warnings.warn(warning_msg, ODEintWarning)

I don't know what is the problem. how can I get correct results? or at least as accurate as possible? or maybe I should rewrite the code in another form? thank you.

import scipy as sio
import numpy as np
import mpmath as mp
import scipy.integrate as spi
import matplotlib.pyplot as plt
import time
from scipy.integrate import quad
initial_value=np.logspace(24,27,100)
t=np.logspace(-20,6,100)
m=np.logspace(0,6,100)
start_time=time.perf_counter()
phi_m={}
phi_m_prime={}
phi=[]
phi_prime=[]
j=0
i=np.pi*2.435*initial_value[0]
while i<(np.pi*(2.435*10**(27))):
    i=np.pi*2.435*initial_value[j]
    phi=[]
    phi_prime=[]
    for x in range (len(m)):
        def dzdt(z,T):
            return [z[1], -3*1.4441*(10**(-6))*m[x]*np.sqrt(0.69)*(mp.coth(1.5*np.sqrt(0.69)*(10**(-6))*1.4441*m[x]*T))*z[1] - z[0]]
        z0 = [i,0]
        ts = t/m[x]
        zs = spi.odeint(dzdt, z0, ts)
        phi.append(zs[99,0])
        phi_prime.append(zs[99,1])
    phi_m[j]=phi
    phi_m_prime[j]=phi_prime
    j+=1
end_time=time.perf_counter()
print(end_time-start_time,"seconds")
danial
  • 21
  • 3
  • Could you please condense down your problem to a minimal example with one equation with a minimal number of terms and coefficients? (In general the structured approach you used is better to avoid errors and compare with the corresponding abstract formulas, but the question is not about that.) Also, your question is slightly off-topic as it is not about the code used but about alternative numerical methods. A better place would be scicomp.SE. – Lutz Lehmann Jun 27 '22 at 10:11
  • @LutzLehmann actually I cant change the coefficients because the answer is very sensitive to them! I edited the code , would you please take another look? – danial Jun 27 '22 at 16:34
  • 1
    Note that your problem may arise from using very large and small values. See [this Q&A](https://stackoverflow.com/q/68500704/2127008). – Wrzlprmft Jun 27 '22 at 16:45
  • I meant for you to select some x and combine the coefficients. In the end the equation seems to look like `z''+a*t*z'+z=0`. Depending on the size of the coefficient this can degenerate rather quickly, the solution contracting to zero, but the increasing stiffness driving the step size to zero equally fast. This is also the error you get, the internal step size falling so low that the time update does not advance the floating-point time. – Lutz Lehmann Jun 27 '22 at 17:12

0 Answers0