jan/meanas
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

randomize imaginary part of starting vector

Browse Source
This commit is contained in:
jan 2018-01-09 00:00:58 -08:00
parent e02040c709
commit 0e47fdd5fb
1 changed files with 8 additions and 8 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

16
fdfd_tools/bloch.py
View File

@ -422,24 +422,24 @@ def eigsolve(num_modes: int,
prev_E = 0 prev_E = 0
d_scale = 1 d_scale = 1
prev_traceGtKG = 0 prev_traceGtKG = 0
prev_theta = 0.5 #prev_theta = 0.5
D = numpy.zeros(shape=y_shape, dtype=complex) D = numpy.zeros(shape=y_shape, dtype=complex)
y0 = None y0 = None
if y0 is None: if y0 is None:
Z = numpy.random.rand(*y_shape).astype(complex) Z = numpy.random.rand(*y_shape) + 1j * numpy.random.rand(*y_shape)
else: else:
Z = y0 Z = y0
while True: while True:
Z2 = Z.conj().T @ Z ZtZ = Z.conj().T @ Z
Z_norm = numpy.sqrt(real(trace(Z2))) / num_modes Z_norm = numpy.sqrt(real(trace(ZtZ))) / num_modes
Z /= Z_norm Z /= Z_norm
Z2 /= Z_norm * Z_norm ZtZ /= Z_norm * Z_norm
try: try:
U = numpy.linalg.inv(Z2) U = numpy.linalg.inv(ZtZ)
except numpy.linalg.LinAlgError: except numpy.linalg.LinAlgError:
Z = numpy.random.rand(*y_shape).astype(complex) Z = numpy.random.rand(*y_shape) + 1j * numpy.random.rand(*y_shape)
continue continue
trace_U = real(trace(U)) trace_U = real(trace(U))
@ -565,7 +565,7 @@ def eigsolve(num_modes: int,
Z += D * numpy.sin(theta) Z += D * numpy.sin(theta)
prev_traceGtKG = traceGtKG prev_traceGtKG = traceGtKG
prev_theta = theta #prev_theta = theta
prev_E = E prev_E = E
''' '''