具有2个连续三次样条曲线操作的问题

2024-05-29 02:09:01 发布

您现在位置:Python中文网/ 问答频道 /正文

我在2个数据列(总共301个文件)上执行第一个三次样条曲线。一切正常,数值范围等于预期值

问题来自三次样条插值的第二个操作,该操作是在python源代码中使用os.system执行一个程序之后完成的(该程序生成其他301个文件,我将其加载到kh, pk

事实上,我在第二次三次样条曲线运算中得到以下错误,以下是错误:

Traceback (most recent call last):
File "Exe_launcher.py", line 323, in <module>
P_m = CubicSpline(np.log10(kh[ppp]),np.log10(pk[ppp]))
File "/Users/fab/Library/Python/2.7/lib/python/site-packages/scipy   /interpolate/_cubic.py", line 527, in __init__
axis = axis % y.ndim
ZeroDivisionError: integer division or modulo by zero

我找不到零的除法是从哪里来的。我在下面给出产生问题的代码片段

我试图找到打印2个数据列的值及其大小的错误,这似乎是正确的,但事实并非如此,因为我得到了上述错误

zrange = range(300)
while paramo < len(fid)+1:

   ## First cubic spline interpolation (paramo = 0)
  print('paramo_final1 = ', paramo)

  if paramo == 0:
   # Test_matterpower
   FPA = fold_path_fid
   aaa = len(zrange)
   while aaa >= 0:
       if aaa == len(zrange):
       kh, pk = np.loadtxt("test_matterpower"+str(aaa)+".dat", usecols=(0,1,), unpack=True)
       elif aaa > 0:
       kh1, pk1 = np.loadtxt("test_matterpower"+str(aaa)+".dat", usecols=(0,1,), unpack=True)
       kh = np.vstack((kh,kh1))
       pk = np.vstack((pk,pk1))
       else:
       kh1, pk1 = np.loadtxt("test_matterpower"+str(len(zrange)+1)+".dat", usecols=(0,1,), unpack=True)
       kh = np.vstack((kh,kh1))
       pk = np.vstack((pk,pk1))
       aaa = aaa-1

   ooo,ppp = 0,0
   integrale = np.zeros((len(pk)))
   print('1) len(kh) = ', len(kh))
   print('len(pk) = ', len(pk))
   print('len(integrale) = ', len(integrale))
   print('kh = ', kh)
   print('kh[0] = ', kh[0])
   print('np.log10(kh[ppp] = ', np.log10(kh[ppp]))
   print('np.log10(pk[ppp] = ', np.log10(pk[ppp]))
   print('pk = ', pk)
   while ppp < len(integrale):
       P_m = CubicSpline(np.log10(kh[ppp]),np.log10(pk[ppp]))
       while ooo < len(kh[0])-1:
       integrale[ppp] = integrale[ppp] + 1./(2*np.pi**2) * (kh[ppp][ooo+1]-kh[ppp][ooo])/6. * ( (kh[ppp][ooo]**2*(3*(np.sin(8*kh[ppp][ooo]) - 8*kh[ppp][ooo]*np.cos(8*kh[ppp][ooo]))/(8*kh[ppp][ooo])**3)**2*pk[ppp][ooo]) + (kh[ppp][ooo+1]**2*(3*(np.sin(8*kh[ppp][ooo+1]) - 8*kh[ppp][ooo+1]*np.cos(8*kh[ppp][ooo+1]))/(8*kh[ppp][ooo+1])**3)**2*pk[ppp][ooo+1]) + 4.*( ((kh[ppp][ooo+1]+kh[ppp][ooo])/2.)**2*(3*(np.sin(8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2.) - 8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2*np.cos(8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2.))/(8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2.)**3)**2*10**P_m(np.log10((kh[ppp][ooo]+kh[ppp][ooo+1])/2))) )
       ooo=ooo+1
       ooo=0
       ppp=ppp+1

   print('1) kh[0] = ', kh[0])
   kh = np.delete(kh, len(kh)-1, axis=0)
   pk = np.delete(pk, len(pk)-1, axis=0)
   print('2) kh[0] = ', kh[0])

  # for paramo > 0 
  else: 

   ## Second cubic spline interpolation (paramo > 0)

     os.system('./camb Exe_launcher_modified.ini')   

     aaa = len(zrange)
     while aaa >= 0:
      if aaa == len(zrange):
          kh, pk = np.loadtxt("test_matterpower"+str(aaa)+".dat", usecols=(0,1,), unpack=True)
      elif aaa > 0:
          kh1, pk1 = np.loadtxt("test_matterpower"+str(aaa)+".dat", usecols=(0,1,), unpack=True)
          kh = np.vstack((kh,kh1))
          pk = np.vstack((pk,pk1))
      else:
          kh1, pk1 = np.loadtxt("test_matterpower"+str(len(zrange)+1)+".dat", usecols=(0,1,), unpack=True)
          kh = np.vstack((kh,kh1))
          pk = np.vstack((pk,pk1))
      aaa = aaa-1

      ooo,ppp = 0,0
      integrale = np.zeros((len(pk)))
      print('2) len(kh) = ', len(kh))
      print('len(pk) = ', len(pk))
      print('len(integrale) = ', len(integrale))
      print('kh = ', kh)
      print('kh[0] = ', kh[0])
      print('np.log10(kh[ppp] = ', np.log10(kh[ppp]))
      print('np.log10(pk[ppp] = ', np.log10(pk[ppp]))
      print('pk = ', pk)
      while ppp < len(integrale):
          P_m = CubicSpline(np.log10(kh[ppp]),np.log10(pk[ppp]))
          while ooo < len(kh[0])-1:
          integrale[ppp] = integrale[ppp] + 1./(2*np.pi**2) * (kh[ppp][ooo+1]-kh[ppp][ooo])/6. * ( (kh[ppp][ooo]**2*(3*(np.sin(8*kh[ppp][ooo]) - 8*kh[ppp][ooo]*np.cos(8*kh[ppp][ooo]))/(8*kh[ppp][ooo])**3)**2*pk[ppp][ooo]) + (kh[ppp][ooo+1]**2*(3*(np.sin(8*kh[ppp][ooo+1]) - 8*kh[ppp][ooo+1]*np.cos(8*kh[ppp][ooo+1]))/(8*kh[ppp][ooo+1])**3)**2*pk[ppp][ooo+1]) + 4.*( ((kh[ppp][ooo+1]+kh[ppp][ooo])/2.)**2*(3*(np.sin(8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2.) - 8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2*np.cos(8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2.))/(8*(kh[ppp][ooo+1]+kh[ppp][ooo])/2.)**3)**2*10**P_m(np.log10((kh[ppp][ooo]+kh[ppp][ooo+1])/2))) )
          ooo=ooo+1
          ooo=0
          ppp=ppp+1

     print('3) kh[0] = ', kh[0])
     kh = np.delete(kh, len(kh)-1, axis=0)
     pk = np.delete(pk, len(pk)-1, axis=0)

     print('4) kh[0] = ', kh[0])

在执行时,实现了第一次三次插值,但没有实现第二次插值:

Traceback (most recent call last):
File "Exe_launcher.py", line 323, in <module>
P_m = CubicSpline(np.log10(kh[ppp]),np.log10(pk[ppp]))
File "/Users/fab/Library/Python/2.7/lib/python/site-packages/scipy   /interpolate/_cubic.py", line 527, in __init__
axis = axis % y.ndim
ZeroDivisionError: integer division or modulo by zero

更新1::问题已修复,抱歉,缩进错误

问候


Tags: lennppppprintpkaxisaaakh

热门问题