我是C和MATLAB的用户。当我开始学习Python时(一周前),我注意到我没有充分利用MATLAB的潜力,特别是数组运算。我经常使用for循环,可能是因为我学会了用C语言编程。
在上一篇技巧文章中,我学习了使用cumsum和其他高效的数组操作,例如:
alpha = [1e-4,1e-3,1e-4,1e-1,1e-2,1e-3,1e-6,1e-3];
zeta = alpha / (dz*dz)
nz = 101
l=[0.3,0.1,0.2,0.1,0.1,0.1,0.2];
wz = cumsum(l*(nz-1));
nl = lenght(l); 有没有可能在Python (Numpy)或MATLAB中简化以下代码?
A = zeros(nz,nz);
i=1;
for j = 2:wz(i)-1
A(j,j-1) = zeta(1,1);
A(j,j) = -2*zeta(1,1);
A(j,j+1) = zeta(1,1); % layer 1 nodes
end
%cicle to n-layers
for i=2:nl
for j=wz(i-1):wz(i-1)
A(j,j-1) = zeta(1,i-1);
A(j,j) = -zeta(1,i-1)-zeta(1,i);
A(j,j+1) = zeta(1,i);
end
for j=wz(i-1)+1:wz(i)
A(j,j-1) = zeta(1,i);
A(j,j) = -2*zeta(1,i);
A(j,j+1) = zeta(1,i);
end
end
end发布于 2012-01-18 15:47:12
在有机会在我的机器和你的机器上并排运行后,我修改了下面的代码。还有几个问题(A在最后的循环中会变得更大吗?dz是什么?)你在运行这个之前遇到的问题是,我忘记了idx_matrix必须符合逻辑。
dz=0.1;
alpha = [1e-4,1e-3,1e-4,1e-1,1e-2,1e-3,1e-6,1e-3];
zeta = alpha / (dz*dz);
nz = 101;
l=[0.3,0.1,0.2,0.1,0.1,0.1,0.2];
wz = cumsum(l*(nz-1));
nl = length(l);
A = zeros(nz);
i=1;
%replaces 1st loop
j_start = 2;
j_end = wz(i)-1;
idx_matrix = false(size(A));
idx_matrix(j_start:j_end,j_start:j_end) = eye(j_end-j_start+1);
A(idx_matrix) = -2*zeta(1,1);
idx_matrix(idx_matrix) = false;
idx_matrix(j_start:j_end,j_start-1:j_end-1) = eye(j_end-j_start+1);
A(idx_matrix) = zeta(1,1);
idx_matrix(idx_matrix) = false;
idx_matrix(j_start:j_end,j_start+1:j_end+1) = eye(j_end-j_start+1);
A(idx_matrix) = zeta(1,1);
%cicle to n-layers
for i=2:nl
%replaces 3rd loop
j_start = wz(i-1);
A(j_start,j_start) = -zeta(1,i-1)-zeta(1,i);
A(j_start,j_start-1) = zeta(1,i-1);
A(j_start,j_start+1) = zeta(1,i);
%replaces 4th loop
j_start = wz(i-1)+1;
j_end = min(wz(i),size(A,2)-1);
idx_matrix = false(size(A));
idx_matrix(j_start:j_end,j_start:j_end) = eye(j_end-j_start+1);
A(idx_matrix) = -2*zeta(1,i);
idx_matrix(idx_matrix) = false;
idx_matrix(j_start:j_end,j_start-1:j_end-1) = eye(j_end-j_start+1);
A(idx_matrix) = zeta(1,i);
idx_matrix(idx_matrix) = false;
idx_matrix(j_start:j_end,j_start+1:j_end+1) = eye(j_end-j_start+1);
A(idx_matrix) = zeta(1,i);
end发布于 2012-01-16 22:08:28
要简化循环,可以使用函数spdiags。
http://www.mathworks.fr/help/techdoc/ref/spdiags.html
例如,你的第一个循环可以写成:
A=full(spdiags(repmat([zeta(1,1),-2*zeta(1,1),zeta(1,1)],wz(i),1),[-1 0 1],wz(i),wz(i)))https://stackoverflow.com/questions/8878413
复制相似问题