gpe.soc ======= .. py:module:: gpe.soc .. autoapi-nested-parse:: Experimental setup for Peter Engels' 2-component Rb 87 counterflow experiments with spin-orbit coupling. The results are summarized in the SOC Soliton.ipynb notebook. The design of this module is as follows: State1, State2: These are two state objects for single component and two component SOC states respectively. The single component state uses the Dispersion class to provide a modified dispersion. These maintain a reference to the relevant Experiment object (see below) and overload the get_Vext() function of the base class to use the appropriate experimental protocol. Experiment*: Each experimental setup is characterized by a subclass of Experiment which specifies the experimental parameters and protocol in terms of physical dimensions. These experimental objects are responsible for returning an appropriate State* object to implement the dynamics. Experiment objects can define a set of time-dependent parameters by defining methods with the name `*_t_()` that take an argument `t_` and returns the value at a specified time. To facilitate plotting, another method `*_info()` can be defined which returns `(unit, label)` with the associated numerical unit, and a text label. These are used for plotting. For example, for a varying Omega, one might define:: def Omega_t_(self, t_): return np.sin(t_) @property def Omega_info(self): return (self.E_R, r'$\Omega/E_R$') The underlying Experiment class designed here allows for the following manipulations. * Time dependent parameters: * Omega: SOC coupling strength * detuning: Detuning * B_gradient: counterflow * xyz0: Location of trap center * ws: Trap frequencies * Functions (which may contain time dependence): * get_Vext(state, fiducial, expt, single_band): External trapping potential. The default version of this combines the following pieces: * Vtrap: (Always) * Vt: (If not fiducial) * Magnetic field gradient: (If not single band) * SOC: (If not single band) * get_Vtrap_expt(state, xyz): Trapping potential used in the experiment. This is used to compute the fiducial_V_TF. Unless you need to modify this structure you should customize the following pieces instead. * get_Vtrap(state, xyz): Trapping potential for your simulation. * get_Vt(state): Time dependent trap. This is only used for initial state preparation but not computing the fiducial_V_TF. The need for this is the case where the initial state is specified with an x_TF without this potential, but then the system is prepared by cooling with this potential. Unless you are matching an experiment where they report x_TF this way, you probably don't need to modify this. * barrier_xyz0: Position of a barrier potential * barrier_depth: Depth of a barrier potential Initial States -------------- The initial state is defined in terms of the chemical potentials which define the initial wavefunction through the Thomas Fermi approximation. From the experimental standpoint, the initial state is typically expressed in terms of the Thomas Fermi "radius" `x_TF` at which point in space the "chemical potential" should match the external potential `V_TF = V_ext(x_TF)`. In the presence of axial degrees of freedom and the various tube approximations, the chemical potential required to reproduce this extent of the cloud can differ from `V_TF` (for example, due to a shift $\hbar\omega_\perp$) so we provide the following functions to convert between the two:: mu = State.get_mu_from_V_TF(V_TF) V_TF = State.get_V_TF_from_mu(mu) To specify the initial state then one should define `V_TF`. This is the fundamental quantity determining the initial state. The alternative - used by default for the experiments - is to specify `x_TF` and then to compute `V_TF` from the potential. This is typically determined as the extent of a cloud without any spin-orbit coupling, barrier potential, etc. The actual initial state is usually obtained by then evolving this state adiabatically into the ground state with SOC, possibly including the barrier potential. This adiabatic preparation should be thought of as roughly performed at constant total particle number (though there may be some particle losses). Thus, the external chemical potential that should be used to compute `V_TF` from `x_TF` may differ from the potentials defined for evolution. Finally, the state may actually be too small to represent the entire and may not even include `x_TF` in the range of abscissa. This leads us to define a `fiducial_V_TF` which is the value of `V_TF` that should be used to initialize the state. This may differ from the `V_TF = V_ext(x_TF)` because of these factors. The function `get_fiducial_V_TF()` is responsible for computing the appropriate `V_TF` to match the experiment and may need to construct another larger state in order to solve for the value of `V_TF` that will appropriately match the experiment. For these reasons, `Experiment.get_Vext(state, fiducial, expt)` and `State.set_initial_state(V_TF, fiducial, expt)` take two flags `fiducial` and `expt` in order to accommodate this functionality (see their docstrings for details). Attributes ---------- .. autoapisummary:: gpe.soc._EPS Classes ------- .. autoapisummary:: gpe.soc.Dispersion gpe.soc.Dispersion3 gpe.soc.SOCMixin gpe.soc.State2 gpe.soc.State2Tube gpe.soc.State2Axial gpe.soc.State1Mixin gpe.soc.State1 gpe.soc.State1Tube gpe.soc.State1Axial gpe.soc.Experiment gpe.soc.ExperimentStub gpe.soc.ExperimentBarrier gpe.soc.Bloch Module Contents --------------- .. py:data:: _EPS .. py:class:: Dispersion Tools for computing properties of the lower band dispersion. Everything is expressed in dimensionless units with $k/k_r$, $d=\delta/4E_R$, $w=\Omega/4E_R$, and $E_\pm/2E_R$. .. rubric:: Examples >>> E = Dispersion(w=1.5/4, d=0.543/4) >>> ks = np.linspace(-3, 3, 500) >>> ks_ = (ks[1:] + ks[:-1])/2.0 >>> ks__ = ks[1:-1] >>> Es = E(ks).T >>> dEs = E(ks, d=1).T >>> ddEs = E(ks, d=2).T >>> d3Es = E(ks, d=3).T >>> d4Es = E(ks, d=4).T >>> dEs_ = np.diff(Es, axis=0)/np.diff(ks)[:, None] >>> ddEs_ = np.diff(dEs, axis=0)/np.diff(ks)[:, None] >>> d3Es_ = np.diff(ddEs, axis=0)/np.diff(ks)[:, None] >>> d4Es_ = np.diff(d3Es, axis=0)/np.diff(ks)[:, None] >>> np.allclose((dEs[1:] + dEs[:-1])/2, dEs_, rtol=0.01) True >>> np.allclose((ddEs[1:] + ddEs[:-1])/2, ddEs_, rtol=0.02) True >>> np.allclose((d3Es[1:] + d3Es[:-1])/2, d3Es_, rtol=0.01) True >>> np.allclose((d4Es[1:] + d4Es[:-1])/2, d4Es_, rtol=0.1) True Here are some regression tests >>> E = Dispersion(w=1.0, d=0) >>> float(E.get_k0()) 0.0 .. py:attribute:: d .. py:attribute:: w .. py:attribute:: E0 .. py:attribute:: k0 .. py:method:: Es(k, d=0) .. py:attribute:: __call__ .. py:method:: newton(k) .. py:method:: get_k0(N=5) Return the minimum of the lower dispersion branch. .. py:class:: Dispersion3 Tools for computing properties of the lower band dispersion for a three-component system. Everything is expressed in dimensionless units with $k/k_r$, $d=\delta/4E_R$, $w=\Omega/4E_R$, $e=\epsilon/2E_R$ and $E_\pm/2E_R$. .. py:attribute:: d .. py:attribute:: w .. py:attribute:: e .. py:attribute:: E0 .. py:attribute:: k0 .. py:method:: Es(k, d=0) .. py:attribute:: __call__ .. py:method:: newton(k) .. py:method:: get_k0(N=5) Return the minimum of the lower dispersion branch. .. py:class:: SOCMixin Bases: :py:obj:`gpe.utils.GPUHelper` State for spin-orbit coupled experiments with a counterflow induced by a magnetic field gradient. This class depends on an Experiment() instance to define the time-dependent properties such as ramping on of the magnetic field gradient, Raman lasers, and ramping off everything for expansion imaging. To do this, the experiment can define the following methods, which will be used if they exist. Note: for performance reasons, these should return `NotImplemented` if they are not used. `experiment.delta_t_(t_)`: `experiment.delta` `experiment.Omega_t_(t_)`: `experiment.Omega` `experiment.B_gradient_t_(t_)`: `experiment.B_gradient` In addition, the following attributes are expected in experiments: `experiment.t_unit`: All experiment times are specified in this unit `t_ = t/t_unit`. `experiment.t__final`: For `t_ > self.experiment.t__final`, all the potentials are set to zero. `experiment.t__image`: If imaging is performed (see `image_ts_` in the Simulation class), then expansion occurs for this length of time (previously this was `t__expand`). .. py:attribute:: experiment :value: None .. py:method:: get_mu_from_V_TF(V_TF) Return the corrected chemical potential from V_TF. :param V_TF: External potential at the Thomas Fermi "radius". (The external potential is evaluated at this position and this is used to get `mu`.) :type V_TF: float .. py:method:: get_V_TF_from_mu(mu) Return V_TF from the chemical potential mu. :param mu: Physical chemical potential (i.e. what you would pass to the minimizer). :type mu: float .. py:method:: get_Vext_GPU(mu=None, fiducial=False, expt=False, single_band=None) Return the external potentials `(V_a, V_b, V_ab)` (or `V_ext` in the single-band case). This method just delegates to the experiment, and provides some simple memoization for performance. :param mu: If None, then subtract `self.mu`. This will minimize the phase rotations during time evolution, but is undesirable when computing initial states. In the latter case, set `mu=0`. :type mu: float, None .. py:method:: get_ns_TF(Vs_TF, fiducial=False, expt=False) Return the TF densities. Assumes gaa = gbb = gab and populates only the lower band. .. py:method:: set_initial_state(V_TF=None, fiducial=False, expt=False) Set the state using the Thomas Fermi (TF) approximation. This will overload the base class versions to pass the `fiducial` and `expt` arguments to the experiments. See `get_fiducial_V_TF()` for details about the meaning of these parameters. This also sets more appropriate phases for the state. :param V_TF: Value of the external potential at the Thomas-Fermi "radius". Note: this is the chemical potential without any corrections from the tube or SOC and is applied directly to the external potential to get the densities. These are then transformed into components using the spin-quasimomentum map. :type V_TF: float :param fiducial: If True, then use the fiducial initial state - i.e. do not include the barrier potential. :type fiducial: bool :param expt: If True, then initialize the state using experimental parameters. This is needed, for example, if the desired state is homogeneous (`wx=0`) but we want to initialize the state with the experimental `x_TF` parameter. the `expt=True` flag should be used when initializing a (typically) larger state with `Lx=Lx_expt` and frequencies `wx=wx_expt` to compute the `fiducial_V_TF` in the presence of a barrier potential. :type expt: bool .. py:method:: get_k0(t_=0, branch=-1) Return the momentum on the lattice that has minimum dispersion. :param branch: 1 for upper branch, -1 for lower branch (default). :type branch: -1 or 1 .. py:method:: get_ws_perp(t=None) Return the frequencies at time `t`. Required by tube2 classes. .. py:property:: ms .. py:property:: bcast Return a set of indices suitable for broadcasting masses etc. .. py:attribute:: rotating_phases :value: False .. py:property:: t_unit Get the time unit from the experiment. .. py:property:: t_final .. py:property:: xyz_trap Return the abscissa for use in the trapping potentials. .. py:property:: k_r .. py:property:: gs Get gs from experiment so it can manipulate them. .. py:property:: ws Get ws from experiment so it can manipulate them. .. py:method:: get(name, t_=None) Return the value of the parameter `name` at time `self.t`. This just delegates to `self.experiment.get()`. .. py:method:: get_density_x() Return the density along x. For periodic boxes, this is the average density in the transverse plane (which is the central density if the transverse plane is translationally invariant) and the integrated 1D line density for tube and axial codes. .. py:method:: get_dispersion(t_=None) .. py:method:: get_homogeneous_ground_state() Return the properties of the homogeneous ground state. .. py:method:: get_branch_mixture() Return the population [n_+, n_-] of the two energy branches assuming the TF approximation is valid. .. py:method:: get_branch_psi_ab(n0, k=None, branch=-1, tol=10 * _EPS) Return the wavefunction factors `(psi_a, psi_b)` for the specified band. In some sense, this is the inverse of `get_branch_mixture()` but for purely homogeneous states. These include the appropriate phase factors of `exp(1j*(k+k_r)*x)` and `exp(1j*(k-k_r)*x)`. :param k: Momentum of homogeneous state (in the single band model) If None, then get the minimium of the dispersion. :type k: float :param n0: Background density (3D) :type n0: float :param experiment: Defines `Omega`, `detuning`, `gs`, `E_R` and `k_r`. :type experiment: Experiment :param branch: 1 for upper branch, -1 for lower branch (default). :type branch: -1 or 1 .. py:class:: PlotElements Collection of elements for redrawing a plot. .. py:attribute:: fig .. py:attribute:: densities .. py:attribute:: history .. py:attribute:: momenta .. py:attribute:: master_grid .. py:method:: _repr_html_() Provide display capability in IPython. .. py:method:: plot(plot_elements=None, fig=None, data=None, grid=None, show_n=True, show_na=True, show_nb=True, show_log_n=True, show_momenta=False, show_mixtures=False, show_V=False, show_history_ab=True, show_history_a=False, show_history_b=False, show_log_history=False, combine_momenta_and_history=False, dynamic_range=100, log_dynamic_range=5, parity=False, clip_momenta_modes=1, k_range_k_r=(-2.5, 2.5), nk_amplitude_factor=0.5, mu_background=None, symmetric=True, fig_width=15, pane_height=2, space=0.1, history=None, split=False, rescale=False, plot_dim=None, subplot_spec=None) Plot the state. :param data: Instance returned by self.get_plot_data() of the data to be plotted. :type data: PlotData :param plot_elements: If defined, then update the plot here (faster, but will not rescale). :type plot_elements: PlotElements :param history: List of states defining the frames. :type history: [State] :param parity: If `True`, then plot abs(x) on the abscissa of the density plots to show parity violations. :type parity: bool :param clip_momenta_modes: Clip this many momenta modes (peaks). If most of the cloud occupies a single background momentum state, then there will be a large peak that obscures features. This many peaks will be clipped. :type clip_momenta_modes: int :param plot_dim: If provided, then try to reduce plots to this dimension. If None, then use self.dim. This is particularly useful for plot_dim = 1 which will plot the line density along the x axis, averaging or integrating over the other dimensions. :type plot_dim: None, int :param combine_momenta_and_history: If True, then combine the momenta and history plots in a single pane. :type combine_momenta_and_history: bool :param The following "boolean" flags can actually be numbers which will: :param specify the relative height of the corresponding pane in the plot: :param show_n: :param show_mixtures: :param show_momenta: :param show_history_ab: :param show_history_a: :param : :param show_history_b: .. py:method:: plot_densities(data, plot_elements, grid=None, split=False, show_n=True, show_na=True, show_nb=True, show_mixtures=False, show_log_n=False, show_V=False, parity=False, symmetric=False, dynamic_range=100, log_alpha=0.3, log_dynamic_range=5) .. py:property:: time_dependent_quantities .. py:method:: plot_time_dependent_parameters(fig=None, subplot_spec=None, share=False) Plot all the time-dependent parameters stacked on top of each other. Used as part of other plotting routines. .. py:method:: plot_mixtures(data, plot_elements, grid=None, dynamic_range=100, show_V=False) Plot the occupation of the two bands. .. py:method:: plot_history(data, plot_elements, grid=None, rescale=False, show_history_ab=True, show_history_a=False, show_history_b=False, show_log_history=False, dynamic_range=100, log_dynamic_range=5, **kw) Show an imcountorf plot of of the states. :param states: List of states to show. :type states: list :param rescale: If True, then scale each time-slice so the maximum is 1. :type rescale: bool .. py:method:: get_momenta_data(clip_momenta_modes=1, log_dynamic_range=5, mu_background=None, k_range_k_r=(-2.5, 2.5)) Return (ks, Em, Ep, E_ph, nk, log_nk). :param clip_momenta_modes: Clip this many momenta modes (peaks). If most of the cloud occupies a single background momentum state, then there will be a large peak that obscures features. This many peaks will be clipped. :type clip_momenta_modes: int :param log_dynamic_range: Return the log of the density scaled so that 0 corresponds to a 10**log_dynamic_range suppression. :type log_dynamic_range: int :param k_range_k_r: Range of momenta (in units of k_r). :type k_range_k_r: (float, float) :param mu_background: Background chemical potential used for computing the phonon dispersion. If this is None, the the maximum density will be used. :type mu_background: float, None :returns: * **ks** (*array*) -- Momenta abscissa. * **Em, Ep** (*array*) -- Dispersion bands. * **Ph_m, Ph_p** (*array*) -- Phonon dispersions * **nk** (*array*) -- Momentum distribution ranging from 0 to the `Em.max()-Em.min()`. This is intended to be added to Em for plotting. * **log_nk** (*array*) -- Log of the momentum distribution ranging from 0 to the `Em.max()-Em.min()`. This is intended to be added to Em for plotting. .. py:class:: MomentaData Bases: :py:obj:`tuple` .. py:attribute:: ks .. py:attribute:: Em .. py:attribute:: Ep .. py:attribute:: E_ph .. py:attribute:: nk .. py:attribute:: log_nk .. py:attribute:: d .. py:attribute:: w .. py:class:: PlotData Bases: :py:obj:`tuple` .. py:attribute:: frame .. py:attribute:: frames .. py:attribute:: xs .. py:attribute:: ts .. py:attribute:: t .. py:attribute:: na .. py:attribute:: nb .. py:attribute:: n .. py:attribute:: nas .. py:attribute:: nbs .. py:attribute:: ns .. py:attribute:: log_na .. py:attribute:: log_nb .. py:attribute:: log_n .. py:attribute:: mask .. py:attribute:: log_mask .. py:attribute:: n_lim .. py:attribute:: state .. py:attribute:: ks .. py:attribute:: Em .. py:attribute:: Ep .. py:attribute:: E_ph .. py:attribute:: nk .. py:attribute:: log_nk .. py:attribute:: d .. py:attribute:: w .. py:method:: get_plot_data(states=None, frame=None, t_unit=u.ms, x_unit=u.micron, n_unit=1.0 / u.micron**3, dynamic_range=100, log_dynamic_range=5, get_momenta_data=False, **get_momenta_kw) Return a PlotData namedtuple with the data needed for plotting. :param states: If a list of states is provided, then this is used to generate data for a history vs time contour plot. :type states: None, [State] :param frame: Frame when making movies. :type frame: int :param \*\*get_momenta_kw: Additional keywords are passed to get_momenta_data. .. py:method:: plot_momenta(data, plot_elements, grid=None, amplitude_factor=0.5, ax=None) Plot the dispersion and momentum distribution. Return plot_elements. :param amplitude_factor: The occupation is shown added on top of the dispersion such that the maximum occupation is this fraction of the dispersion range. I.e. If amplitude_factor=0.5, then the maximum momentum occupation will be half of the dispersion range. :type amplitude_factor: float .. py:method:: make_movie(filename, states, skip=1, fps=20, skip_frames=10, show_frames=2, dpi=None, display=True, extra_args=('-vcodec', 'libx264', '-pix_fmt', 'yuv420p'), **kw) Make a movie of the states. :param filename: Name of movie file such as "movie.mp4". :type filename: str :param states: List of states defining the frames. :type states: [State] :param skip: If this is >1, then skip this many states when making the movie. This allows one to have more states for smooth history plots, but speeds making the movie. :type skip: int :param fps: Frames-per-second for output movie :type fps: int :param show_frames: Display this many initial frames to aid debugging. Does not affect the generated movie. :type show_frames: int :param skip_frames: After the initial frames, draw only these frames. (Drawing frames greatly slows down the making of the movie so we limit the display). Does not affect the generated movie. :type skip_frames: int :param display: If `False`, then do not display anything. (Use for offline generation of movies for example.) :type display: bool .. py:class:: State2(experiment, **kw) Bases: :py:obj:`SOCMixin`, :py:obj:`gpe.bec2.State` State for spin-orbit coupled experiments with a counterflow induced by a magnetic field gradient. This class depends on an Experiment() instance to define the time-dependent properties such as ramping on of the magnetic field gradient, Raman lasers, and ramping off everything for expansion imaging. To do this, the experiment can define the following methods, which will be used if they exist. Note: for performance reasons, these should return `NotImplemented` if they are not used. `experiment.delta_t_(t_)`: `experiment.delta` `experiment.Omega_t_(t_)`: `experiment.Omega` `experiment.B_gradient_t_(t_)`: `experiment.B_gradient` In addition, the following attributes are expected in experiments: `experiment.t_unit`: All experiment times are specified in this unit `t_ = t/t_unit`. `experiment.t__final`: For `t_ > self.experiment.t__final`, all the potentials are set to zero. `experiment.t__image`: If imaging is performed (see `image_ts_` in the Simulation class), then expansion occurs for this length of time (previously this was `t__expand`). .. py:attribute:: single_band :value: False .. py:attribute:: experiment .. py:attribute:: _time_independent_Vext :value: False .. py:attribute:: _Vext :value: [None, None] .. py:class:: State2Tube(experiment, **kw) Bases: :py:obj:`SOCMixin`, :py:obj:`gpe.tube2.StateGPEdrZ` State for spin-orbit coupled experiments with a counterflow induced by a magnetic field gradient. This class depends on an Experiment() instance to define the time-dependent properties such as ramping on of the magnetic field gradient, Raman lasers, and ramping off everything for expansion imaging. To do this, the experiment can define the following methods, which will be used if they exist. Note: for performance reasons, these should return `NotImplemented` if they are not used. `experiment.delta_t_(t_)`: `experiment.delta` `experiment.Omega_t_(t_)`: `experiment.Omega` `experiment.B_gradient_t_(t_)`: `experiment.B_gradient` In addition, the following attributes are expected in experiments: `experiment.t_unit`: All experiment times are specified in this unit `t_ = t/t_unit`. `experiment.t__final`: For `t_ > self.experiment.t__final`, all the potentials are set to zero. `experiment.t__image`: If imaging is performed (see `image_ts_` in the Simulation class), then expansion occurs for this length of time (previously this was `t__expand`). .. py:attribute:: single_band :value: False .. py:attribute:: experiment .. py:attribute:: _time_independent_Vext :value: False .. py:attribute:: _Vext :value: [None, None] .. py:class:: State2Axial(experiment, **kw) Bases: :py:obj:`SOCMixin`, :py:obj:`gpe.axial.State2Axial` State for spin-orbit coupled experiments with a counterflow induced by a magnetic field gradient. This class depends on an Experiment() instance to define the time-dependent properties such as ramping on of the magnetic field gradient, Raman lasers, and ramping off everything for expansion imaging. To do this, the experiment can define the following methods, which will be used if they exist. Note: for performance reasons, these should return `NotImplemented` if they are not used. `experiment.delta_t_(t_)`: `experiment.delta` `experiment.Omega_t_(t_)`: `experiment.Omega` `experiment.B_gradient_t_(t_)`: `experiment.B_gradient` In addition, the following attributes are expected in experiments: `experiment.t_unit`: All experiment times are specified in this unit `t_ = t/t_unit`. `experiment.t__final`: For `t_ > self.experiment.t__final`, all the potentials are set to zero. `experiment.t__image`: If imaging is performed (see `image_ts_` in the Simulation class), then expansion occurs for this length of time (previously this was `t__expand`). .. py:attribute:: single_band :value: False .. py:attribute:: experiment .. py:attribute:: _time_independent_Vext :value: False .. py:attribute:: _Vext :value: [None, None] .. py:class:: State1Mixin Mixin for single-band states. .. py:attribute:: single_band :value: True .. py:method:: laplacian(y, factor, exp=False) Return the laplacian applied to the state `y` multiplied by `factor`. :param y: The laplacian will be applied to this state (not `self`). :type y: State :param factor: The result will be multiplied by this factor. :type factor: array-like :param exp: If `True` then `exp(factor*laplacian)(y)` will be computed instead. :type exp: bool .. py:property:: kx2 Return the effective `kx**2` for use in the laplacian. .. py:class:: State1(experiment, **kw) Bases: :py:obj:`State1Mixin`, :py:obj:`SOCMixin`, :py:obj:`gpe.bec.StateTwist_x` Mixin for single-band states. .. py:attribute:: experiment .. py:attribute:: _time_independent_Vext :value: False .. py:attribute:: _Vext :value: None .. py:class:: State1Tube(experiment, **kw) Bases: :py:obj:`State1Mixin`, :py:obj:`SOCMixin`, :py:obj:`gpe.tube.StateGPEdrZ` Tube state for single band model. .. py:attribute:: experiment .. py:attribute:: _time_independent_Vext :value: False .. py:attribute:: _Vext :value: [None, None] .. py:class:: State1Axial(experiment, **kw) Bases: :py:obj:`State1Mixin`, :py:obj:`SOCMixin`, :py:obj:`gpe.axial.StateAxial` Axial state for single_band models. .. py:attribute:: experiment .. py:attribute:: _time_independent_Vext :value: False .. py:attribute:: _Vext :value: None .. py:class:: Experiment(_local_dict=None, **kw) Bases: :py:obj:`gpe.utils.ExperimentBase` Represents a set of experimental parameters. 1. Responsible for translating these into a proper State object. This generally amounts to converting experimental parameters into values useful to the state class. These are stored in self.state_args for use in `self.get_state` and `self.get_initial_state`. 2. Provides get_Vext(state) with a (possibly time dependent) external potential. 3. Implements an experimental protocol. To do this, subclass and define the run_experiment method. This can be used to get the initial state, then to manipulate it. (Not implemented yet.) Note: these experiments are designed for SOC along the x-axis only. The following basis options are supported: '1D' : Periodic basis along the x-axis. `tube==False` : Assume a constant density along y and z. This should be thought of as modeling a 3D gas with translation symmetry along the transverse directions. (The coupling constants are not modified for modeling 1D systems.) The density here will be the central density. `tube==True` : Use the the tube code to model the transverse trap. The density here will be the line density obtained by integrating over the transverse directions. This will effectively modify the coupling constants. 'axial' : Periodic basis along the x-axis with radial excitations assuming axial symmetry. The transverse trapping frequency will be the average of the y and z frequencies. The radial extent and lattice spacing will be inferred from `Nx`, `dx` and `Lx` using the ratio of the trapping frequencies at time `t_=0`. '2D', '3D' : Full 2D or 3D simulation assuming no symmetry. The transverse extents and lattice spacing should be provided in `Nxyz` and `Lxyz`, otherwise it will be inferred from `Nx`, `dx` and `Lx` using the ratio of the trapping frequencies at time `t_=0`. Note: All times in the Experiment classes are specified in units of `t_unit`. In arguments this is `t_ = t/t_unit` and we use `t_` as a kwarg to ensure that this is properly used. :param cells_x: If specified (not `None`), then the length of the box in the x direction `Lxyz[0]` will be computed to be an integral number of cells to hold `cells_x` periods of a wave with frequency `k_r`. In addition, the trap will be reset to `'cos'` which will replace the HO potential with a compatible, but periodic cosine. This allows one to simulate the center of the trapped system with matching gradients and densities. If this is defined, then the original `Lx` will be stored as `Lx_expt` which is used to obtain initial states. :type cells_x: None, int :param dx: If specified (not `None`), then `Nx` will be computed such that his minimum lattice spacing is realized, guaranteeing a certain UV cutoff. :type dx: None, (float,) :param initial_imbalance: Initial state is prepared in (1,-1) then this fraction is transfered to the (1, 0) state. :type initial_imbalance: float :param Omega: This is the strength of the spin-orbit coupling. :type Omega: float :param delta: This is the value of detuning in energy units. :type delta: float :param B_gradient: This is the strength of the magnetic field gradient which induces a counterflow. :type B_gradient: float :param x0: Location of the center of the trapping potential. :type x0: float :param The following attributes are guaranteed to exist for use by the State: :param class.: .. attribute:: k_r Recoil wave-vector. :type: float .. attribute:: E_R Recoil energy `E_R = (hbar*k_r)**2/2/m`. :type: float .. py:attribute:: State :value: None .. py:attribute:: t_unit :value: 63.50779925891489 .. py:attribute:: t_name :value: 'ms' .. py:attribute:: t__final Return `t__final = t_final/t_unit`. .. py:attribute:: t__image :value: 7 .. py:attribute:: species .. py:attribute:: gs :value: None .. py:attribute:: trapping_frequencies_Hz :value: (3, 273, 277) .. py:attribute:: trapping_wavelength :value: 1.064 .. py:attribute:: x_TF :value: 234.6 .. py:attribute:: mu :value: None .. py:attribute:: detuning_kHz :value: 0.1 .. py:attribute:: detuning_E_R :value: None .. py:attribute:: recoil_frequency_Hz :value: 3683.8 .. py:attribute:: rabi_frequency_E_R :value: 2.9 .. py:attribute:: rabi_frequency_kHz :value: None .. py:attribute:: B_gradient_mG_cm :value: 10 .. py:attribute:: Nx :value: None .. py:attribute:: Lx :value: 1000.0 .. py:attribute:: dx :value: 0.061 .. py:attribute:: Nr :value: None .. py:attribute:: R :value: None .. py:attribute:: cells_x :value: None .. py:attribute:: old_cells_x :value: True .. py:attribute:: initial_imbalance :value: 0.5 .. py:attribute:: cooling_phase :value: 1.0 .. py:attribute:: tube :value: True .. py:attribute:: twist :value: 0 .. py:attribute:: x0 :value: 0.0 .. py:attribute:: basis_type :value: '1D' .. py:attribute:: gab_attenuation_t_ :value: None .. py:attribute:: single_band :value: False .. py:attribute:: active_components :value: 2 .. py:attribute:: epsilon_E_R :value: 3.8 .. py:property:: dim .. py:method:: init() Overload this to perform any initial computations. .. py:method:: get_Vext(state, fiducial=False, expt=False) Return (V_a, V_b, V_ab), the external potentials. For `t_ > self.t__final`, all the potentials are set to zero. :param state: Current state. Use this to get the time `state.t` and abscissa `state.basis.xyz`. :type state: IState :param fiducial: If `True`, then return the potential that should be used to define the initial state in terms of the Thomas Fermi radius of the cloud `x_TF`. :type fiducial: bool :param expt: If `True`, then return the proper experimental potential rather than the potential used in the simulation. :type expt: bool .. seealso:: * interface.IExperiment.get_Vext .. py:method:: get_state(expt=False, initialize=True) Quickly return an appropriate initial state. .. py:method:: get_initial_state(perturb=0.0, E_tol=1e-12, psi_tol=1e-12, disp=1, tries=20, cool_steps=100, cool_dt_t_scale=0.1, minimize=True, **kw) Return an initial state with the specified population fractions. This initial state is prepared in state[0] with the potentials as they are at time `t=0`, then the `initial_imbalance` is transferred as specified simulating an RF pulse by simply the appropriate fraction in each state. Phases are kept the same as in the state[0]. .. py:method:: get_Vtrap_expt(state, xyz, ws=None) Return the external trapping potential used in the experiment. This is used to determine the chemical potential and `fidcucial_V_TF` from a value of `x_TF`. It uses the experimental trapping frequencies stored in `ws_expt`. :param state: Current state. :type state: State :param xyz: Use these abscissa (not state.xyz) to compute the potential. This is likely being applied to a much larger state in order to compute the `fiducial_V_TF`. :type xyz: [array] :param ws: These are the trapping frequencies. This is just a convenience if your get_Vtrap needs to compute a harmonic potential. :type ws: [float] .. py:method:: get_Vtrap(state, xyz) Return the external trapping potential used in the simulation. This should return the trapping potential that is used in cases, both when determining x_TF, and when cooling into the initial state. It should not include time-dependent potential like the bucket that are used for preparing the initial state but NOT used when determining x_TF. :param state: Current state. :type state: State :param xyz: Use these abscissa (not state.xyz) to compute the potential. This is likely being applied to a much larger state in order to compute the `fiducial_V_TF`. :type xyz: [array] .. py:method:: get_Vt(state) Optional time-dependent trapping potentials. These potentials are not included in the fiducial state used to determine the initial conditions, however, if `Vt` is non-zero at time `t=0`, then this *will* be included in the initial state preparation. .. seealso:: * interface.IExperiment.get_Vext .. py:method:: _xyz_trap(state, x_periodic=False) Return the abscissa for use in the trapping potentials. This has the abscissa shifted by x0 so that this value is zero, with periodic wrapping. I.e. it assumes that the trapping potentials are periodic with center at zero. :param x_periodic: If `True`, then transform the abscissa so that any potential applied is smooth and periodic as discussed in the Utils.ipynb notebook. This assumes additionally that the potentials are even functions about zero. :type x_periodic: bool .. py:method:: _get_Lx(state) Return the actual length of the box represented by the state. This is needed because self.Lx might not be the same as self.Lx_expt and the latter is sometimes needed for fiducial states. .. py:property:: fiducial_V_TF This may be slow to calculate, so we defer calculation until we really need it. .. py:method:: get_fiducial_V_TF(t_=0.0, Nx=2**12, Lx_factor=1.1) Return the V_TF required to initialize the state. If V_TF is None or not defined, then compute the V_TF that defines the state in terms of the Thomas-Fermi radius x_TF along the x axis using a Harmonic trapping potential with frequencies `ws_expt` as follows: 1. Construct an axially symmetric state using a CylindricalBasis with enough points to model the harmonically trapped cloud in the TF approximation. The assumption is made here that the extremities of the cloud are harmonic. Inside, the potential need not be harmonic. This construction uses `self.ws_expt`. 2. Compute `V_TF` from `self.x_TF` and `self.ws_expt` assuming harmonic extremities. 3. Get the external potential using this state by calling `get_Vext(expt_state, fiducial=True, expt=True)`. This allows the function `get_Vext()` to customize what is meant by the fiducial state (i.e. including or excluding parts of the time-varying potential.) 4. Set the initial state and compute the total particle number in the TF approximation using: `expt_state.set_initial_state(V_TF=V_TF, fiducial=True, expt=True)` 5. Adjust V_TF to reproduce the total particle number in the TF approximation using: `expt_state.set_initial_state(V_TF=V_TF, fiducial=False, expt=True)` This allows one to simulate a long time cooling into the ground with a different potential, while fixing `x_TF` without that potential. .. py:method:: get_dispersion(t_=0.0, active_components=None) Return Dispersion instance. :param t_: Time at which to get parameters. :type t_: float :param active_components: Specify which dispersion to use - Dispersion or Dispersion3 :type active_components: 2, 3 .. py:method:: gs_t_(t_) .. py:method:: plotter(fig=None, plot=True) .. py:method:: run(get_data, fig=None, figsize=(20, 5), plot=True, factor=8.0, filename='generated_plots.mp4') .. py:method:: plot_dispersion(t_=0.0, klims=(-2, 2), hv=None) .. py:class:: ExperimentStub .. py:attribute:: Omega :value: 1.0 .. py:attribute:: delta :value: 1.0 .. py:attribute:: t_unit :value: 1.0 .. py:attribute:: xyz0 :value: (0.0, 0.0, 0.0) .. py:method:: get(key, t_) .. py:class:: ExperimentBarrier(_local_dict=None, **kw) Bases: :py:obj:`Experiment` Adds to the basic Experiment class a moving barrier potential (along the `x` axis only) with the following parameters: .. attribute:: barrier_x Position of barrier potential. :type: float .. attribute:: barrier_depth Strength of the barrier potential (energy). :type: float .. attribute:: barrier_width Width of the barrier potential (length) :type: float .. attribute:: barrier_type Type of barrier potential. The 'constant' barrier will just provide the cosine modulation. :type: 'gaussian', 'harmonic', 'constant' .. attribute:: barrier_k If this is non-zero, then the barrier has an underlying cosine modulation with this wavelength. :type: float .. attribute:: These use the function get_barrier_potential(state). .. py:attribute:: barrier_type :value: 'gaussian' .. py:attribute:: barrier_k :value: 0.0 .. py:method:: _barrier_potential(x_barrier_width) Return the barrier potential without the power factor. x_barrier_width : x/barrier_width .. py:method:: get_Vt(state) Optional time-dependent trapping potentials. These potentials are not included in the fiducial state used to determine the initial conditions, however, if `Vt` is non-zero at time `t=0`, then this *will* be included in the initial state preparation. .. seealso:: * interface.IExperiment.get_Vext .. py:class:: Bloch(experiment=ExperimentStub()) Classes and Glasses .. py:attribute:: hbar :value: 1.0 .. py:attribute:: m :value: 1.0 .. py:attribute:: experiment .. py:property:: k_r .. py:property:: t_unit .. py:property:: Omega .. py:property:: delta .. py:method:: pack(psi) .. py:method:: unpack(q) .. py:method:: get_su2_generators() .. py:property:: beta .. py:method:: get_beta_norm(beta=None) .. py:method:: get_Bloch_period(t_ms=None) .. py:method:: rhs(q, t) .. py:method:: get_initial_state(theta=np.pi, phi=np.pi / 2.0) returns an initial state packed for odeint. theta is the azimuthal angle, phi is in the xy-plane in the bloch sphere .. py:method:: get_psi_t(t_=20.0, Nt_=5000) .. py:method:: plot() .. py:method:: get_Bloch_vector(psi) returns the theta and phi of the Bloch vector theta is the azimuthal angle, phi is in the xy-plane .. py:method:: get_psi(a) :staticmethod: Return the two-component state from the vector `a` where the total density is the length of `a` and the state points in the direction specified with the top component real and positive.