jan/fdfd_tools
2
2
Fork 1
Code Issues Pull Requests Releases Wiki Activity

Fixup poynting operators for new approach

Browse Source
This commit is contained in:
Jan Petykiewicz 2019-09-27 20:43:32 -07:00
parent a1a7aa5511
commit 487bdd61ec
1 changed files with 15 additions and 23 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
meanas/fdfd/operators.py
View File

@ -453,22 +453,18 @@ def poynting_e_cross(e: vfield_t, dxes: dx_lists_t) -> sparse.spmatrix:
""" """
shape = [len(dx) for dx in dxes[0]] shape = [len(dx) for dx in dxes[0]]
bx, by, bz = [rotation(i, shape, -1) for i in range(3)] fx, fy, fz = [rotation(i, shape, 1) for i in range(3)]
dxag = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[0], indexing='ij')] dxag = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[0], indexing='ij')]
dxbg = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[1], indexing='ij')] dxbg = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[1], indexing='ij')]
dbgx, dbgy, dbgz = [sparse.diags(dx) for dx in dxbg] Ex, Ey, Ez = [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)]
n = numpy.prod(shape) block_diags = [[ None, fx @ -Ez, fx @ Ey],
[ fy @ Ez, None, fy @ -Ex],
P = sparse.bmat( [ fz @ -Ey, fz @ Ex, None]]
[[ None, -bx @ Ez @ dbgy, bx @ Ey @ dbgz], block_matrix = sparse.bmat([[sparse.diags(x) if x is not None else None for x in row]
[ by @ Ez @ dbgx, None, -by @ Ex @ dbgz], for row in block_diags])
[-bz @ Ey @ dbgx, bz @ Ex @ dbgy, None]]) P = block_matrix @ sparse.diags(numpy.concatenate(dxag))
#TODO
P2 = sparse.block_diag((bx, by, bz)) @ cross([Ex, Ey, Ez]) @ sparse.diags(numpy.concatenate(dxbg))
print(sparse.linalg.norm((P-P2)), sparse.linalg.norm(P), sparse.linalg.norm(P2))
return P return P
@ -482,21 +478,17 @@ def poynting_h_cross(h: vfield_t, dxes: dx_lists_t) -> sparse.spmatrix:
""" """
shape = [len(dx) for dx in dxes[0]] shape = [len(dx) for dx in dxes[0]]
fx, fy, fz = [avgf(i, shape) for i in range(3)] #TODO fx, fy, fz = [rotation(i, shape, 1) for i in range(3)]
bx, by, bz = [avgb(i, shape) for i in range(3)]
dxag = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[0], indexing='ij')]
dxbg = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[1], indexing='ij')] dxbg = [dx.ravel(order='C') for dx in numpy.meshgrid(*dxes[1], indexing='ij')]
dagx, dagy, dagz = [sparse.diags(dx.ravel(order='C'))
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)]
n = numpy.prod(shape) P = (sparse.bmat(
[[ None, -Hz @ fx, Hy @ fx],
P = sparse.bmat( [ Hz @ fy, None, -Hx @ fy],
[[ None, Hz @ bx @ dagy, Hy @ bx @ dagz], [-Hy @ fz, Hx @ fz, None]])
[ Hz @ by @ dagx, None, -Hx @ by @ dagz], @ sparse.diags(numpy.concatenate(dxag)))
[-Hy @ bz @ dagx, Hx @ bz @ dagy, None]])
return P return P