This is my code in python, the dimension of sx should be of 100X4 and sy 100X1 by the multiplication (sx)(B)(sy).
import numpy as np
B= [[-6.08066634428988e-10, -8.61023850910464e-11, 5.48222828615260e-12, -9.49229025004441e-14],
[-3.38148313553674e-11, 6.47759097087283e-12, 1.14900158474371e-13, -5.70078947874486e-15],
[-2.55893304237669e-13, -1.40941560399352e-13, 5.76510238931847e-15, -5.52980385181738e-17],
[3.39795122177475e-15, 7.95704191204353e-16, -5.31260642039813e-17, 7.83532802015832e-19]]
[X, Y] = np.meshgrid(np.arange(0, 3, 0.01*3),np.arange(0, 15, 0.01*(15)))
sx=[]
sy=[]
F=[]
for i in range(len(X)):
for j in range(len(X)):
for k in range(len(B)):
sx[i,k].append(X[i,j]**k)
for l in range(len(B)):
sy[l].append((Y[i,j]**l))
F[i,j] = sx*B*sy
The error:
sx[i,k].append(X[i,j]**k) TypeError: list indices must be integers or slices, not tuple
MATLAB code copied from comment (guess as to formatting)
[x,y]=meshgrid(0:0.01*3:3,0:0.01*15:15);
for i=1:size(x)
for j=1:size(x)
for k=0:size(B) -1
sx(1,k 1)=(x(i,j)^k);
end
for k=0:size(B) -1
sy(k 1,1)=(y(i,j)^k);
end
G(i,j)=sx*B*sy;
end
end
CodePudding user response：
If sx or X is a 2D list then indices must be [i][j]. If you're trying to append to two indices i and j then it should be separate calls to append.
CodePudding user response：
In an Octave session:
B =
-6.0807e-10 -8.6102e-11 5.4822e-12 -9.4923e-14
-3.3815e-11 6.4776e-12 1.1490e-13 -5.7008e-15
-2.5589e-13 -1.4094e-13 5.7651e-15 -5.5298e-17
3.3980e-15 7.9570e-16 -5.3126e-17 7.8353e-19
>>
>> [x,y]=meshgrid(0:0.01*3:3,0:0.01*15:15);
>> for i=1:size(x)
for j=1:size(x)
for k=0:size(B) -1
sx(1,k 1)=(x(i,j)^k);
end
for k=0:size(B) -1
sy(k 1,1)=(y(i,j)^k);
end
G(i,j)=sx*B*sy;
end
end
produces
x, y, G (101 x 101)
>> sx (1,4)
sx =
1 3 9 27
>> sy (4,1)
sy =
1
15
225
3375
So the G element is (1,4) * (4,4) * (4,1) => (1,1)
Looks like I should be able to make a
In [100]: B= [[-6.08066634428988e-10, -8.61023850910464e-11, 5.48222828615260e-12, -9.49229025004441e-14],
...: [-3.38148313553674e-11, 6.47759097087283e-12, 1.14900158474371e-13, -5.70078947874486e-15],
...: [-2.55893304237669e-13, -1.40941560399352e-13, 5.76510238931847e-15, -5.52980385181738e-17],
...: [3.39795122177475e-15, 7.95704191204353e-16, -5.31260642039813e-17, 7.83532802015832e-19]]
...:
In [101]: B = np.array(B)
In [106]: [X, Y] = np.meshgrid(np.linspace(0, 3, 101),np.linspace(0, 15, 101),indexing='ij')
In [107]: X.shape
Out[107]: (101, 101)
In [108]: k = np.arange(0,4)
In [109]: k
Out[109]: array([0, 1, 2, 3])
In [110]: SX = X[:,:,None]**k # (101,101,4)
In [111]: SY = Y[:,:,None]**k
In [114]: G = np.einsum('ijk,kl,ijl->ij',SX,B,SY)
In [115]: G.shape
Out[115]: (101, 101)
Allowing for the "F" order of MATLAB (ie. transpose), looks like these results match:
>> G(1,1)
ans = -0.00000000060807
In [118]: G[0,0]
Out[118]: -6.08066634428988e-10
>> G(50,23)
ans = -0.00000000097117
In [119]: G[22,49]
Out[119]: -9.71172989297259e-10
With broadcasting I don't to make the meshgrid arrays
In [121]: x, y = np.linspace(0,3,101), np.linspace(0,15,101)
In [124]: sx = x[:,None]**k
In [125]: sy = y[:,None]**k
In [126]: sx.shape
Out[126]: (101, 4)
In [129]: g = sx@[email protected]
In [130]: g.shape
Out[130]: (101, 101)
In [131]: np.allclose(G,g)
Out[131]: True
Here I'm doing a matrix product of
(101,4) (4,4) (4,100) => (101,101)