科学的python数组语法

Scientific python array syntax

所以我在 python 中编写了这段代码。我不打算解释它,因为它是一个我似乎看不到的简单语法修复,所以解释它的用途是没有用的。

我遇到的问题是,对于给定的 d,例如 15,我得到 "cuentas" 正确的值和 "e" 正确的值。

我想做的是遍历一组 d 并获取每个 cuentas 和每个 e 的值,以便绘制 d 与 e。

我的问题是我似乎不知道如何在 python 中创建矩阵。

在 matlab 中我曾经写过两个不同的循环,比如 theese:

for i=1:1:N
    for j=1:9

      a[i,j]= and so on

a[i,j] 将是一个 N 行 9 列的矩阵,我可以访问和操作它。

在我下面的代码中,我会故意在我想要遍历距离的地方放置#条注释

import numpy as np
import matplotlib.pyplot as plt
N=100000
cos=np.zeros(N)
phi=np.zeros(N)
teta=np.zeros(N)
a=np.zeros(N)


xrn=np.zeros(N)
yrn=np.zeros(N)
zrn=np.zeros(N)

x=np.zeros(N)
y=np.zeros(N)
z=np.zeros(N)

lim1=14.7
lim2=3.35
lim3=-lim1
lim4=-lim2

#d=np.array([15,20,25,30,35,40,45,50,55])
d=15

    #for j in range(9):
for i in range(N):
    cos[i]=np.random.uniform(-1,1)
    teta[i]=np.random.uniform(-np.pi,np.pi)
    phi[i]=np.random.uniform(0,2*np.pi)

# a[i]=d[j]/cos[i]*np.cos(phi[i])
a[i]=d/cos[i]*np.cos(phi[i])


xrn[i]=np.random.uniform(-1,1)*lim1
yrn[i]=np.random.uniform(-1,1)*lim2

x[i]=a[i]*np.sin(teta[i])*np.cos(phi[i])+xrn[i]
y[i]=a[i]*np.sin(teta[i])*np.sin(phi[i])+yrn[i]

#cuentas[j]=0
cuentas=0

#for j in range(9):
for i in range(N):
    if a[i]>0 and x[i] < lim1 and x[i]>lim3 and y[i] < lim2 and y[i]>lim4:
        cuentas=cuentas+1
#e[j]=cuenta[j]/N
e=cuentas/N

非常感谢阅读的人!

您可以按以下方式使用 numpy 在 python 中创建矩阵:

n=5
k=4
a=np.zeros([n,k])
for i in range(n):
    for j in range(k):
        a[i,j]=i+j
print(a)

结果是

[[ 0.  1.  2.  3.]
 [ 1.  2.  3.  4.]
 [ 2.  3.  4.  5.]
 [ 3.  4.  5.  6.]
 [ 4.  5.  6.  7.]]

简短版本:

这与 Python

中的 MATLAB 代码完全相同
a = np.zeros([N, 9])
for i in range(N):
    for j in range(9):
        a[i,j]= and so on

唯一的主要区别是您需要预先定义数组,如果您希望代码具有合理的性能,您也确实应该在 MATLAB 中这样做。

但是,如果您事先不知道大小,可以在 Python 中使用列表,然后在最后转换为 numpy 数组。这将比你的大型数组的 MATLAB 示例快得多,因为列表和 matrices/arrays 的内部处理方式:

a = []
for i in range(N):
    a.append([])
    for j in range(9):
        a[-1].append( and so on
a = np.array(a)

[-1] 表示("last element of a",append 将括号内的内容放在列表末尾。所以。a[-1].append(foo) 表示“将fooa.

的最后一个元素中的任何内容

长版:

您的 MATLAB 代码在 Python 中的工作方式大致相同,但您需要考虑一些显着差异。

首先,分配给大于现有 array/matrix 的索引在 MATLAB 中有效,但在 numpy 中无效。所以如果你有一个大小 [5, 5] array/matrix,在 MATLAB 中你可以分配给元素 [5, 6],但你不能在 numpy 中。这意味着在 MATLAB 中您可以从一个空数组开始,而在 numpy 中您必须事先设置数组大小。请注意,MATLAB 矩阵实际上不能调整大小,实际上每次循环都在创建一个新矩阵并将所有数据复制到它。这非常慢,这就是 MATLAB 警告您预先分配数组的原因。 Numpy 只是不会假装能够调整数组大小,因此您需要更明确地进行复制、预分配或使用列表(可调整大小)。

其次,同样,MATLAB 不要求您在使用前定义矩阵,而 numpy 则需要。这是因为传统上 MATLAB 具有三种数据结构(矩阵、元胞数组和结构体),每种都有自己的索引方式。因此 MATLAB 可以根据您的索引方式确定您想要创建哪种数据结构。 Python 只有一种索引方式,所以无法进行这种猜测。

第三,在 MATLAB 中使用单个数组大小和一些(但不是全部)函数创建一个二维方阵,每个维度都具有该大小,而在 numpy 中它创建一个一维数组。我不确定您的代码是否符合您的期望。坦率地说,我不知道为什么 MATLAB 会这样工作。

第四,numpy 数组可以有任意维数,0(标量)、1(向量)、2、3、4 等。另一方面,MATLAB 矩阵必须至少有两个维度。这可能会导致一些意想不到的差异,比如转置对 numpy 向量什么都不做。

至于您的 Python 代码,如果您不说出了什么问题,我就无法告诉您如何修复它。但希望我已经给了你足够的信息让你自己做。

所以我采纳了你的答案并且成功了!

如果有人想知道这是一个 Monte Carlo 模拟有多少粒子会通过两个矩形探测器。当我从探测器 1 扔出粒子时,它们通过它是微不足道的,我计算通过探测器 2 的粒子数。

更正后的代码是

N=100000
"La dirección viene dada por v=[rsin(teta)*cos(phi),rsin(teta)sin(phi),rcos(teta)]"
"Los vectores que vamos a usar debemos inicializarlos como un vector de ceros"
cos=np.zeros([10,N])
phi=np.zeros([10,N])
teta=np.zeros([10,N])
a=np.zeros([10,N])

xrn=np.zeros(N)
yrn=np.zeros(N)
zrn=np.zeros(N)

x=np.zeros([10,N])
y=np.zeros([10,N])
z=np.zeros([10,N])

lim1=14.7
lim2=3.35
lim3=-lim1
lim4=-lim2
"d son las disversas distancias a las que colocamos la fuente con respecto al detector"

d=np.array([0.00001,15,20,25,30,35,40,45,50,55])

"e es la eficiencia geométrica simulada"
e=np.zeros(10)

"Debemos definir el coseno como números aleatorios en vez de el ángulo teta, debido a que queremos"
"que se distribuyan uniformemente por toda la esfera"
for j in range(10):
    for i in range(N):
        cos[j,i]=np.random.uniform(0,1)

        phi[j,i]=np.random.uniform(0,2*np.pi)

        a[j,i]=d[j]/cos[j,i]

        xrn[i]=np.random.uniform(-1,1)*lim1
        yrn[i]=np.random.uniform(-1,1)*lim2

        x[j,i]=a[j,i]*np.sin(math.acos(cos[j,i]))*np.cos(phi[j,i])+xrn[i]
        y[j,i]=a[j,i]*np.sin(math.acos(cos[j,i]))*np.sin(phi[j,i])+yrn[i]


cuentas=np.zeros(10)
for j in range(10):
    for i in range(N):
        if a[j,i]>0 and x[j,i] < lim1 and x[j,i]>lim3 and y[j,i] < lim2 and y[j,i]>lim4:
            cuentas[j]=cuentas[j]+1

    e[j]=cuentas[j]/N

感谢大家!