Switch to C-ordered arrays

This commit is contained in:
Jan Petykiewicz 2017-03-26 18:22:12 -07:00
parent 1e80a66b50
commit 48ddd9f512
2 changed files with 11 additions and 13 deletions

View File

@ -302,11 +302,9 @@ def rotation(axis: int, shape: List[int], shift_distance: int=1) -> sparse.spmat
n = numpy.prod(shape) n = numpy.prod(shape)
i_ind = numpy.arange(n) i_ind = numpy.arange(n)
j_ind = ijk[0] + ijk[1] * shape[0] j_ind = numpy.ravel_multi_index(ijk, shape, order='C')
if len(shape) == 3:
j_ind += ijk[2] * shape[0] * shape[1]
vij = (numpy.ones(n), (i_ind, j_ind.flatten(order='F'))) vij = (numpy.ones(n), (i_ind, j_ind.flatten(order='C')))
d = sparse.csr_matrix(vij, shape=(n, n)) d = sparse.csr_matrix(vij, shape=(n, n))
@ -350,7 +348,7 @@ def shift_with_mirror(axis: int, shape: List[int], shift_distance: int=1) -> spa
if len(shape) == 3: if len(shape) == 3:
j_ind += ijk[2] * shape[0] * shape[1] j_ind += ijk[2] * shape[0] * shape[1]
vij = (numpy.ones(n), (i_ind, j_ind.flatten(order='F'))) vij = (numpy.ones(n), (i_ind, j_ind.flatten(order='C')))
d = sparse.csr_matrix(vij, shape=(n, n)) d = sparse.csr_matrix(vij, shape=(n, n))
return d return d
@ -371,7 +369,7 @@ def deriv_forward(dx_e: List[numpy.ndarray]) -> List[sparse.spmatrix]:
def deriv(axis): def deriv(axis):
return rotation(axis, shape, 1) - sparse.eye(n) return rotation(axis, shape, 1) - sparse.eye(n)
Ds = [sparse.diags(+1 / dx.flatten(order='F')) @ deriv(a) Ds = [sparse.diags(+1 / dx.flatten(order='C')) @ deriv(a)
for a, dx in enumerate(dx_e_expanded)] for a, dx in enumerate(dx_e_expanded)]
return Ds return Ds
@ -392,7 +390,7 @@ def deriv_back(dx_h: List[numpy.ndarray]) -> List[sparse.spmatrix]:
def deriv(axis): def deriv(axis):
return rotation(axis, shape, -1) - sparse.eye(n) return rotation(axis, shape, -1) - sparse.eye(n)
Ds = [sparse.diags(-1 / dx.flatten(order='F')) @ deriv(a) Ds = [sparse.diags(-1 / dx.flatten(order='C')) @ deriv(a)
for a, dx in enumerate(dx_h_expanded)] for a, dx in enumerate(dx_h_expanded)]
return Ds return Ds
@ -463,8 +461,8 @@ def poynting_e_cross(e: vfield_t, dxes: dx_lists_t) -> sparse.spmatrix:
fx, fy, fz = [avgf(i, shape) for i in range(3)] fx, fy, fz = [avgf(i, shape) for i in range(3)]
bx, by, bz = [avgb(i, shape) for i in range(3)] bx, by, bz = [avgb(i, shape) for i in range(3)]
dxag = [dx.flatten(order='F') for dx in numpy.meshgrid(*dxes[0], indexing='ij')] dxag = [dx.flatten(order='C') for dx in numpy.meshgrid(*dxes[0], indexing='ij')]
dbgx, dbgy, dbgz = [sparse.diags(dx.flatten(order='F')) dbgx, dbgy, dbgz = [sparse.diags(dx.flatten(order='C'))
for dx in numpy.meshgrid(*dxes[1], indexing='ij')] for dx in numpy.meshgrid(*dxes[1], indexing='ij')]
Ex, Ey, Ez = [sparse.diags(ei * da) for ei, da in zip(numpy.split(e, 3), dxag)] Ex, Ey, Ez = [sparse.diags(ei * da) for ei, da in zip(numpy.split(e, 3), dxag)]
@ -492,8 +490,8 @@ def poynting_h_cross(h: vfield_t, dxes: dx_lists_t) -> sparse.spmatrix:
fx, fy, fz = [avgf(i, shape) for i in range(3)] fx, fy, fz = [avgf(i, shape) for i in range(3)]
bx, by, bz = [avgb(i, shape) for i in range(3)] bx, by, bz = [avgb(i, shape) for i in range(3)]
dxbg = [dx.flatten(order='F') for dx in numpy.meshgrid(*dxes[1], indexing='ij')] dxbg = [dx.flatten(order='C') for dx in numpy.meshgrid(*dxes[1], indexing='ij')]
dagx, dagy, dagz = [sparse.diags(dx.flatten(order='F')) dagx, dagy, dagz = [sparse.diags(dx.flatten(order='C'))
for dx in numpy.meshgrid(*dxes[0], indexing='ij')] for dx in numpy.meshgrid(*dxes[0], indexing='ij')]
Hx, Hy, Hz = [sparse.diags(hi * db) for hi, db in zip(numpy.split(h, 3), dxbg)] Hx, Hy, Hz = [sparse.diags(hi * db) for hi, db in zip(numpy.split(h, 3), dxbg)]

View File

@ -27,7 +27,7 @@ def vec(f: field_t) -> vfield_t:
""" """
if numpy.any(numpy.equal(f, None)): if numpy.any(numpy.equal(f, None)):
return None return None
return numpy.hstack(tuple((fi.flatten(order='F') for fi in f))) return numpy.hstack(tuple((fi.flatten(order='C') for fi in f)))
def unvec(v: vfield_t, shape: numpy.ndarray) -> field_t: def unvec(v: vfield_t, shape: numpy.ndarray) -> field_t:
@ -45,5 +45,5 @@ def unvec(v: vfield_t, shape: numpy.ndarray) -> field_t:
""" """
if numpy.any(numpy.equal(v, None)): if numpy.any(numpy.equal(v, None)):
return None return None
return [vi.reshape(shape, order='F') for vi in numpy.split(v, 3)] return [vi.reshape(shape, order='C') for vi in numpy.split(v, 3)]