import mmf_setup;mmf_setup.nbinit()
from gpe.imports import *
This cell adds /home/docs/checkouts/readthedocs.org/user_builds/gpe/checkouts/latest/src to your path, and contains some definitions for equations and some CSS for styling the notebook. If things look a bit strange, please try the following:
- Choose "Trust Notebook" from the "File" menu.
- Re-execute this cell.
- Reload the notebook.
[I 19:59:35 numexpr.utils] NumExpr defaulting to 2 threads.
---------------------------------------------------------------------------
KeyboardInterrupt Traceback (most recent call last)
Cell In[1], line 2
1 import mmf_setup;mmf_setup.nbinit()
----> 2 from gpe.imports import *
File ~/checkouts/readthedocs.org/user_builds/gpe/checkouts/latest/src/gpe/imports.py:28
26 from mmfutils.plot import imcontourf # noqa: E402
27 from gpe.minimize import MinimizeState # noqa: E402
---> 28 from gpe.utils import evolve_to, evolve, evolves # noqa: E402
29 from gpe.plot_utils import MPLGrid # noqa: E402
30 from pytimeode.evolvers import EvolverSplit, EvolverABM # noqa: E402
File ~/checkouts/readthedocs.org/user_builds/gpe/checkouts/latest/src/gpe/utils.py:35
32 from persist.objects import Archivable
33 from persist.archive import Archive
---> 35 from pytimeode.evolvers import EvolverABM
36 from pytimeode.mixins import ArrayStateMixin
37 from pytimeode.interfaces import implementer, IStateForABMEvolvers
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/pytimeode/__init__.py:2
1 from . import interfaces
----> 2 from . import mixins
3 from . import evolvers
5 __all__ = ["interfaces", "mixins", "evolvers"]
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/pytimeode/mixins.py:39
30 from .interfaces import (
31 IState,
32 IStateApply,
(...) 36 implementer,
37 )
38 from . import interfaces
---> 39 from .utils import expr
41 __all__ = [
42 "StateMixin",
43 "StatesMixin",
(...) 47 "ArraysStateWithBraketMixin",
48 ]
51 @implementer(IState)
52 class StateMixin:
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/pytimeode/utils/expr.py:18
15 import numpy as np
17 try:
---> 18 import sympy
19 except ImportError:
20 sympy = None
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/__init__.py:77
70 from .logic import (to_cnf, to_dnf, to_nnf, And, Or, Not, Xor, Nand, Nor,
71 Implies, Equivalent, ITE, POSform, SOPform, simplify_logic, bool_map,
72 true, false, satisfiable)
74 from .assumptions import (AppliedPredicate, Predicate, AssumptionsContext,
75 assuming, Q, ask, register_handler, remove_handler, refine)
---> 77 from .polys import (Poly, PurePoly, poly_from_expr, parallel_poly_from_expr,
78 degree, total_degree, degree_list, LC, LM, LT, pdiv, prem, pquo,
79 pexquo, div, rem, quo, exquo, half_gcdex, gcdex, invert,
80 subresultants, resultant, discriminant, cofactors, gcd_list, gcd,
81 lcm_list, lcm, terms_gcd, trunc, monic, content, primitive, compose,
82 decompose, sturm, gff_list, gff, sqf_norm, sqf_part, sqf_list, sqf,
83 factor_list, factor, intervals, refine_root, count_roots, all_roots,
84 real_roots, nroots, ground_roots, nth_power_roots_poly, cancel,
85 reduced, groebner, is_zero_dimensional, GroebnerBasis, poly,
86 symmetrize, horner, interpolate, rational_interpolate, viete, together,
87 BasePolynomialError, ExactQuotientFailed, PolynomialDivisionFailed,
88 OperationNotSupported, HeuristicGCDFailed, HomomorphismFailed,
89 IsomorphismFailed, ExtraneousFactors, EvaluationFailed,
90 RefinementFailed, CoercionFailed, NotInvertible, NotReversible,
91 NotAlgebraic, DomainError, PolynomialError, UnificationFailed,
92 GeneratorsError, GeneratorsNeeded, ComputationFailed,
93 UnivariatePolynomialError, MultivariatePolynomialError,
94 PolificationFailed, OptionError, FlagError, minpoly,
95 minimal_polynomial, primitive_element, field_isomorphism,
96 to_number_field, isolate, round_two, prime_decomp, prime_valuation,
97 galois_group, itermonomials, Monomial, lex, grlex,
98 grevlex, ilex, igrlex, igrevlex, CRootOf, rootof, RootOf,
99 ComplexRootOf, RootSum, roots, Domain, FiniteField, IntegerRing,
100 RationalField, RealField, ComplexField, PythonFiniteField,
101 GMPYFiniteField, PythonIntegerRing, GMPYIntegerRing, PythonRational,
102 GMPYRationalField, AlgebraicField, PolynomialRing, FractionField,
103 ExpressionDomain, FF_python, FF_gmpy, ZZ_python, ZZ_gmpy, QQ_python,
104 QQ_gmpy, GF, FF, ZZ, QQ, ZZ_I, QQ_I, RR, CC, EX, EXRAW,
105 construct_domain, swinnerton_dyer_poly, cyclotomic_poly,
106 symmetric_poly, random_poly, interpolating_poly, jacobi_poly,
107 chebyshevt_poly, chebyshevu_poly, hermite_poly, hermite_prob_poly,
108 legendre_poly, laguerre_poly, apart, apart_list, assemble_partfrac_list,
109 Options, ring, xring, vring, sring, field, xfield, vfield, sfield)
111 from .series import (Order, O, limit, Limit, gruntz, series, approximants,
112 residue, EmptySequence, SeqPer, SeqFormula, sequence, SeqAdd, SeqMul,
113 fourier_series, fps, difference_delta, limit_seq)
115 from .functions import (factorial, factorial2, rf, ff, binomial,
116 RisingFactorial, FallingFactorial, subfactorial, carmichael,
117 fibonacci, lucas, motzkin, tribonacci, harmonic, bernoulli, bell, euler,
(...) 138 Znm, elliptic_k, elliptic_f, elliptic_e, elliptic_pi, beta, mathieus,
139 mathieuc, mathieusprime, mathieucprime, riemann_xi, betainc, betainc_regularized)
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/polys/__init__.py:79
3 __all__ = [
4 'Poly', 'PurePoly', 'poly_from_expr', 'parallel_poly_from_expr', 'degree',
5 'total_degree', 'degree_list', 'LC', 'LM', 'LT', 'pdiv', 'prem', 'pquo',
(...) 65 'field', 'xfield', 'vfield', 'sfield'
66 ]
68 from .polytools import (Poly, PurePoly, poly_from_expr,
69 parallel_poly_from_expr, degree, total_degree, degree_list, LC, LM,
70 LT, pdiv, prem, pquo, pexquo, div, rem, quo, exquo, half_gcdex, gcdex,
(...) 76 nth_power_roots_poly, cancel, reduced, groebner, is_zero_dimensional,
77 GroebnerBasis, poly)
---> 79 from .polyfuncs import (symmetrize, horner, interpolate,
80 rational_interpolate, viete)
82 from .rationaltools import together
84 from .polyerrors import (BasePolynomialError, ExactQuotientFailed,
85 PolynomialDivisionFailed, OperationNotSupported, HeuristicGCDFailed,
86 HomomorphismFailed, IsomorphismFailed, ExtraneousFactors,
(...) 91 MultivariatePolynomialError, PolificationFailed, OptionError,
92 FlagError)
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/polys/polyfuncs.py:10
8 from sympy.polys.polyoptions import allowed_flags, build_options
9 from sympy.polys.polytools import poly_from_expr, Poly
---> 10 from sympy.polys.specialpolys import (
11 symmetric_poly, interpolating_poly)
12 from sympy.polys.rings import sring
13 from sympy.utilities import numbered_symbols, take, public
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/polys/specialpolys.py:298
294 return dmp_mul(f, h, n, K), dmp_mul(g, h, n, K), h
296 # A few useful polynomials from Wang's paper ('78).
--> 298 from sympy.polys.rings import ring
300 def _f_0():
301 R, x, y, z = ring("x,y,z", ZZ)
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/polys/rings.py:31
27 from sympy.polys.polyoptions import (Domain as DomainOpt,
28 Order as OrderOpt, build_options)
29 from sympy.polys.polyutils import (expr_from_dict, _dict_reorder,
30 _parallel_dict_from_expr)
---> 31 from sympy.printing.defaults import DefaultPrinting
32 from sympy.utilities import public, subsets
33 from sympy.utilities.iterables import is_sequence
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/printing/__init__.py:11
7 from .mathml import mathml, print_mathml
9 from .python import python, print_python
---> 11 from .pycode import pycode
13 from .codeprinter import print_ccode, print_fcode
15 from .codeprinter import ccode, fcode, cxxcode, rust_code # noqa:F811
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/printing/pycode.py:11
9 from sympy.core.mod import Mod
10 from .precedence import precedence
---> 11 from .codeprinter import CodePrinter
13 _kw = {
14 'and', 'as', 'assert', 'break', 'class', 'continue', 'def', 'del', 'elif',
15 'else', 'except', 'finally', 'for', 'from', 'global', 'if', 'import', 'in',
16 'is', 'lambda', 'not', 'or', 'pass', 'raise', 'return', 'try', 'while',
17 'with', 'yield', 'None', 'False', 'nonlocal', 'True'
18 }
20 _known_functions = {
21 'Abs': 'abs',
22 'Min': 'min',
23 'Max': 'max',
24 }
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/printing/codeprinter.py:13
11 from sympy.core.sorting import default_sort_key
12 from sympy.core.symbol import Symbol
---> 13 from sympy.functions.elementary.complexes import re
14 from sympy.printing.str import StrPrinter
15 from sympy.printing.precedence import precedence, PRECEDENCE
File ~/checkouts/readthedocs.org/user_builds/gpe/conda/latest/lib/python3.14/site-packages/sympy/functions/__init__.py:28
25 from sympy.functions.elementary.integers import floor, ceiling, frac
26 from sympy.functions.elementary.piecewise import (Piecewise, piecewise_fold,
27 piecewise_exclusive)
---> 28 from sympy.functions.special.error_functions import (erf, erfc, erfi, erf2,
29 erfinv, erfcinv, erf2inv, Ei, expint, E1, li, Li, Si, Ci, Shi, Chi,
30 fresnels, fresnelc)
31 from sympy.functions.special.gamma_functions import (gamma, lowergamma,
32 uppergamma, polygamma, loggamma, digamma, trigamma, multigamma)
33 from sympy.functions.special.zeta_functions import (dirichlet_eta, zeta,
34 lerchphi, polylog, stieltjes, riemann_xi)
File <frozen importlib._bootstrap>:1371, in _find_and_load(name, import_)
1367 getattr(getattr(module, "__spec__", None), "_initializing", False)):
1368 with _ModuleLockManager(name):
1369 module = sys.modules.get(name, _NEEDS_LOADING)
1370 if module is _NEEDS_LOADING:
-> 1371 return _find_and_load_unlocked(name, import_)
1372
1373 # Optimization: only call _bootstrap._lock_unlock_module() if
1374 # module.__spec__._initializing is True.
File <frozen importlib._bootstrap>:1345, in _find_and_load_unlocked(name, import_)
1341 try:
1342 module = _load_unlocked(spec)
1343 finally:
1344 if parent_spec:
-> 1345 parent_spec._uninitialized_submodules.pop()
1346 if parent:
1347 # Set the module as an attribute on its parent.
1348 parent_module = sys.modules[parent]
File <frozen importlib._bootstrap>:953, in _load_unlocked(spec)
949 module = sys.modules.pop(spec.name)
950 sys.modules[spec.name] = module
951 _verbose_message('import {!r} # {!r}', spec.name, spec.loader)
952 finally:
--> 953 spec._initializing = False
954
955 return module
File <frozen importlib._bootstrap_external>:755, in _LoaderBasics.exec_module(self, module)
753 def exec_module(self, module):
754 """Execute the module."""
--> 755 code = self.get_code(module.__name__)
756 if code is None:
757 raise ImportError(f'cannot load module {module.__name__!r} when '
758 'get_code() returns None')
File <frozen importlib._bootstrap_external>:852, in SourceLoader.get_code(self, fullname)
848 else:
849 source_mtime = int(st['mtime'])
850 try:
851 data = self.get_data(bytecode_path)
--> 852 except OSError:
853 pass
854 else:
855 exc_details = {
File <frozen importlib._bootstrap_external>:951, in FileLoader.get_data(self, path)
947 def get_data(self, path):
948 """Return the data from path as raw bytes."""
949 if isinstance(self, (SourceLoader, SourcelessFileLoader, ExtensionFileLoader)):
950 with _io.open_code(str(path)) as file:
--> 951 return file.read()
952 else:
953 with _io.FileIO(path, 'r') as file:
954 return file.read()
KeyboardInterrupt:
SOC Supersolids#
Here we consider a two-component BEC with SOC of equal Rashba and Dresselhaus couplings.
\[\begin{split}
\op{H} = \begin{pmatrix}
K(p) - \mu + \frac{\delta}{2} & \frac{\Omega}{2} e^{2\I k_Rx}\\
\frac{\Omega}{2} e^{-2\I k_Rx} & K(p) - \mu - \frac{\delta}{2}\\
\end{pmatrix}, \qquad
\mat{R} = \begin{pmatrix}
e^{\I k_Rx}\\
& e^{-\I k_Rx}
\end{pmatrix},\\
\mat{R}^\dagger\cdot\op{H}\cdot\mat{R} = \mat{H}_R
= \begin{pmatrix}
K(p+k_R) - \mu + \frac{\delta}{2} & \frac{\Omega}{2}\\
\frac{\Omega}{2}& K(p-k_R) - \mu - \frac{\delta}{2}\\
\end{pmatrix}\\
\I\hbar \dot{\Psi} = \op{H}\Psi, \qquad
\Psi = \mat{R}\Psi_R,\\
\I\hbar \dot{\Psi}_R
= (\op{H}_R - \I\hbar\mat{R}^\dagger\dot{\mat{R}})\Psi_R,\\
\end{split}\]
(In this case, \(\dot{\mat{R}} = 0\), but we will consider a more general case below.)
Natural scales are expressed in terms \(E_R\) and \(k_R\) with the dimensionless parameters \(k\), \(d\), and \(w\):
\[
4E_R = 4\frac{\hbar^2k_R^2}{2m}, \qquad
k = \frac{p}{\hbar k_R}, \qquad
d = \frac{\delta}{4E_R}, \qquad
w = \frac{\Omega}{4E_R}.
\]
One can realize this by choosing units so that \(\hbar = m/2 = k_R = 1\). In terms of these, the single-particle dispersion bands are:
\[
\frac{E_{\pm}(p)}{2E_R} = \frac{k^2 + 1}{2} \pm \sqrt{(k - d)^2 + w^2}.
\]
In addition to these single-particle properties, one has non-linear couplings with an interaction energy of the form:
\[
E_{\mathrm{int}} = g_{aa} \frac{n_a^2}{2} + g_{bb} \frac{n_a^2}{2} + g_{ab}n_a n_b.
\]
The relevant dimensionless quantities are:
\[
e_1 = \frac{\bar{n}}{4E_R}\frac{g_{aa} + g_{bb} + 2g_{ab}}{8}, \qquad
e_2 = \frac{\bar{n}}{4E_R}\frac{g_{aa} + g_{bb} - 2g_{ab}}{8}, \qquad
e_3 = \frac{\bar{n}}{4E_R}\frac{g_{bb} - g_{aa}}{4}.
\]
Homogeneous States#
Excluding the super
%pylab inline --no-import-all
k = np.linspace(-1.1,1.1,100)
def E_m(k, d=0, w=1.0):
return (k**2+1)/2 - np.sqrt((k-d)**2+w**2)
plt.plot(k, E_m(k, d=0.01, w=0.9))
Supersolid Stripe Phase#
As discussed in Martone:2015, a supersolid stripe phase can be realized if \(\delta = 0\) and \(w\) is sufficiently small if the superfluid is miscible \(g=\sqrt{g_{aa}g_{bb}} > g_{ab}\). (The miscible \(^{87}\)Rb hyperfine states \(a=\ket{1, -1}\) and \(b=\ket{1,0}\) states have \(a_{aa} = 100.40a_B\), \(a_{bb}=100.86a_B\), and \(a_{ab}= 100.41a_B\). The relevant quantity is the geometric mean \(a = 100.63a_B\), so the miscibility criterion is satisfied.)
They use the following Anzatz:
\[\begin{split}
\Psi = \sqrt{\bar{n}}\left[
C_+
\begin{pmatrix}
\cos\theta_+\\
-\sin\theta_+
\end{pmatrix}
e^{\I k_+ x}
+
C_-
\begin{pmatrix}
\sin\theta_-\\
-\cos\theta_-
\end{pmatrix}
e^{-\I k_- x}
\right], \qquad
\abs{C_+}^2+\abs{C_-}^2 = 1.
\end{split}\]
The supersolid phase occurs when \(\abs{C_\pm} = 1/\sqrt{2}\) with \(k_{\pm} = \pm k_R\) with density
\[
\frac{n(x)}{\bar{n}} =
1 + \frac{\Omega}{2(2E_R + G_1)}\cos(2k_R x + \phi)
= 1 + \frac{w}{1 + 2e_1}\cos(2k_R x + \phi)
\]
(Note: in their notations \(e_i = G_i/4E_R\).) .
In the limits of small \(w\) and small \(d\), the conditions for a stripe phase is:
\[
\lim_{w\rightarrow 0}: \quad \abs{d + 2e_3} \leq 4e_2, \qquad
\lim_{d \rightarrow 0}: \quad w < \sqrt{\frac{1}{1+\frac{e_1}{e_2}}}
\approx 0.033.
\]
The actual region is approximately triangular with these vertices. The bound on \(w\) is for the \(a=\ket{1, -1}\) and \(b=\ket{1,0}\) states of \(^{87}\)Rb, which also sets an upper bound on the maximum relative density contrast is \(\delta n / \bar{n} < w/(1 + 2e_1) < w\).
Demler#
To compare with [Kasper:2020], we note the following correspondence:
\[\begin{split}
g \equiv g_{aa} = g_{bb}, \qquad
0 \equiv g_{ab},\qquad
0 \equiv \delta, \qquad
\frac{Q}{2} \equiv k_R, \qquad
M \equiv m, \\
2\hbar J \equiv \Omega, \qquad
x_0 = \sqrt{\frac{\hbar}{2MJ}} \equiv \frac{\hbar}{\sqrt{m\Omega}},\qquad
E_0 = \frac{\hbar^2}{2Mx_0^2} \equiv \frac{\Omega}{2},\\
\bar{Q} = x_0Q \equiv \sqrt{\frac{8 E_R}{\Omega}} = \sqrt{\frac{2}{w}}, \qquad
\bar{\mu} = \frac{\mu}{E_0} \equiv \frac{2 \mu}{\Omega}, \qquad
\bar{g} = \frac{g}{x_0E_0} \equiv \frac{2g\sqrt{m}}{\hbar\sqrt{\Omega}}.
\end{split}\]
We thus have:
\[\begin{split}
\mu = g \bar{n} = \frac{\Omega}{2}\bar{\mu}\\
e_1 = e_2 = \frac{w}{8}\bar{\mu}, \qquad
e_3 = 0, \qquad
w = \frac{2}{\bar{Q}^2}, \qquad
\end{split}\]
In Figure 2. they have \(\bar{\mu}=30\) and \(\bar{g}=0.6\). Since \(d=e_3=0\), and \(e_1=e_2\), there should be a SS phase for \(w < \sqrt{1/2} = 0.707\) corresponding to \(\bar{Q} > \sqrt[4]{8} = 1.68\). This is consistent with their homogeneous (commensurate) phase for \(\bar{Q}=0.7\) and a inhomogeneous (incommensurate) SS-phase for \(\bar{Q}=1.9\).
np.random.seed(2)
g, m, k_R, Omega, hbar, n = np.random.random(6)
E_R = hbar**2*k_R**2/2/m
w = Omega/4/E_R
mu = g*n
Q = 2*k_R
J = Omega/2/hbar
x0 = np.sqrt(hbar/2/m/J)
E0 = hbar**2/x0**2/2/m
Qbar = x0*Q
mubar = mu/E0
gbar = g/x0/E0
assert np.allclose(E0, Omega/2)
assert np.allclose(Qbar, np.sqrt(2/w))
assert np.allclose(mubar, 2*mu/Omega)
assert np.allclose(gbar, 2*g*np.sqrt(m/Omega)/hbar)
Hamner Phase Winding#
In Hamner:2013, Peter looks at phase winding from a coupled system with the following Hamiltonian:
\[\begin{split}
\op{H} = \begin{pmatrix}
K(\op{p}) + \Delta_0 + \delta x& \Omega_0\\
\Omega_0 & K(\op{p}),
\end{pmatrix}, \qquad
K(\op{p}) = \frac{\op{p}^2}{2m},
\end{split}\]
with the usual non-linear couplings and \(\Omega_0/2\pi \hbar = 7.4\)kHz. Here they couple the \(\ket{1,-1}\) and \(\ket{2,0}\) states with \(a_{aa} = 100.40a_B\), \(a_{bb}=94.57a_B\), and \(a_{ab}= 98.13a_B\). These states are weakly immiscible, so supersolidity is not relevant, but perhaps there is an analogy?
Consider now the transformation
\[\begin{split}
\mat{R}_{a} = \begin{pmatrix}
e^{-\I \delta t x /\hbar}\\
& 1
\end{pmatrix}, \qquad
\Psi = \mat{R}_a\Psi_r.
\end{split}\]
The effective Hamiltonian for \(\Psi_R\) is:
\[\begin{split}
\I\hbar\dot{\Psi}_R = \mat{H}_R\dot{\Psi}_R, \qquad
\mat{H}_R = \mat{R}_a^\dagger\cdot\op{H}\cdot\mat{R}_a - \I\hbar \mat{R}_a^\dagger\dot{\mat{R}}_a
= \begin{pmatrix}
K(\op{p}_x - \delta t) + \Delta_0 & \Omega_0 e^{\I \delta t x/\hbar}\\
\Omega_0 e^{-\I \delta t x/\hbar} & K(\op{p}_x)
\end{pmatrix}.
\end{split}\]
This looks somewhat like a SOC Hamiltonian with a time-dependent \(k_R = \delta t x/2\hbar\).
ga, gb, gab = (100.4, 94.57, 98.13)
#ga, gb, gab = (100.4, 100.86, 100.41)
print(np.sqrt(1/(1+(ga+gb + 2*gab)/(ga+gb - 2*gab))))
from gpe.bec import u
delta = u.hbar * 2*np.pi *2.3*u.kHz
wc = 0.033
Omega0 = u.hbar * 2*np.pi * u.kHz
t = 2*np.sqrt(u.m * Omega0/wc)/delta
t/u.ms
If we identify this with an SOC system, then we have a time-varying \(k_R = \delta t/2\hbar\), hence we might obtain a transition to the supersolid phase when:
\[
w = \frac{m\Omega_0}{\hbar^2 k_R^2} < 0.033,\qquad
t > \frac{2}{\delta}\sqrt{\frac{m\Omega_0}{0.033}}.
\]
We consider this system now
\[\begin{split}
\op{H} = \begin{pmatrix}
K(p-\kappa t) & \frac{\Omega}{2} e^{2\I \kappa t x}\\
\frac{\Omega}{2} e^{-2\I \kappa t x} & K(p+\kappa t)\\
\end{pmatrix}, \qquad
\mat{R} = \begin{pmatrix}
e^{\I \kappa t x}\\
& e^{-\I \kappa t x}
\end{pmatrix},\qquad
-\I\hbar\mat{R}^\dagger\dot{\mat{R}} = \begin{pmatrix}
\hbar\kappa x \\
& -\hbar\kappa x
\end{pmatrix}\\
(\op{H}_R - \I\hbar\mat{R}^\dagger\dot{\mat{R}}) =
\begin{pmatrix}
K(p) + \hbar \kappa x & \frac{\Omega}{2}\\
\frac{\Omega}{2}& K(p) - \hbar \kappa x\\
\end{pmatrix}\\
\end{split}\]
References#
Approach for making visible and stable stripes in a spin-orbit-coupled Bose-Einstein superfluid:
Martone:2015: Characterization of the phase diagram.
[Kasper:2020]: Simulating a quantum commensurate-incommensurate phase transition using two Raman coupled one dimensional condensates: Paper by Demler’s group that produces a similar Hamiltonian with a time-varying field coupling two tubes that can tunnel. A SOC-like Hamiltonian emerges (they call it a Pokrovsky-Talapov (PT) model). They do not mention the SOC case.
Rotating a Supersolid Dipolar Gas: Different context, but similar supersolid physics.
Hamner:2013: Peter’s phase-winding experiment.
Trap#
Here we use our code to construct the supersolid strip phase in a small external trap. We use this trap to break the degeneracy of the ground state to help our minimizer find the solution.
import mmf_setup;mmf_setup.nbinit()
from gpe.imports import *
import soc_supersolid;reload(soc_supersolid)
from soc_supersolid import u
e = soc_supersolid.Supersolid(
#lattice_height_V_TF=0.1*0.5, T__x=np.inf,
lattice_height_V_TF=0.0, T__x=np.inf,
#xi_micron=0.125/2,
xi_micron=0.1,
#dx=0.01 * u.micron,
#w=2.7/4.0, d=0.0,
#w=2.7/4.0, d=0./16,lattice_k_k_r=0.78,
w=0.1/4.0, d=-0.001223419561528528,
lattice_k_k_r=1.0,
#lattice_height_V_TF=0.0, T__x=2.0,
barrier_x=0*u.micron,
cells_x=20)
s = e.get_state()
s = e.get_initial_state()
clear_output(wait=True)
res = e.plot(s)
e1, e2, e3 = res['es']
w, d = res['w'], res['d']
w_c = np.sqrt(1/(1+e1/e2))
print(f"SS window: w={w:.4f} < {w_c:.4f}, d={d:.4f} in {-2*e3:.4f}+-{4*e2:.4f}")
if not(w < w_c):
print(f"w={w:.4f} > {w_c:.4f} too large for SS")
if not(abs(d + 2*e3) < 4*e2):
print(f"d = {d:.4f} outside window ({-2*e3-4*e2:.4f}, {-2*e3+4*e2:.4f}) for SS")
n = s.get_density().sum(axis=0)
(n.max() - n.min())/(n.max() + n.min()), w/(1+2*e1)
E = e.get_dispersion()
ks = np.linspace(-3,3,100)
plt.plot(ks, E(ks)[0])
E.get_k0()
s = e.get_initial_state()
clear_output()
s.plot(show_momenta=True);
history = [s]
dhistory = [s - s]
ev = EvolverABM(history[-1], dt=0.3*s.t_scale)
pe = None
with NoInterrupt() as interrupted:
while not interrupted:
for n in range(10):
ev.evolve(100)
history.append(ev.get_y())
dhistory.append(history[-1] - history[0])
pe = ev.y.plot(plot_elements=None, history=dhistory)
display(pe.fig)
clear_output(wait=True)
plt.close('all')
import mmfutils.plot.colors
cm = mmfutils.plot.colors.Colormaps.diverging
s = history[0]
n0 = s.get_density()
xs = s.xyz[0].ravel()
ts = [_s.t/_s.t_unit for _s in history]
plt.subplot(211)
ns = np.array([(_s.get_density()-n0)[0] for _s in history])
imcontourf(ts, xs, ns, cmap=cm, vmax=0.4, vmin=-0.4);plt.xlabel('t [ms]')
plt.colorbar()
plt.subplot(212)
ns = np.array([(_s.get_density()-n0)[1] for _s in history])
imcontourf(ts, xs, ns, cmap=cm, vmax=0.4, vmin=-0.4);plt.xlabel('t [ms]')
plt.colorbar()
Lattice#
str(b"dijf".decode()