Calculations/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/setupParameters.m

103 lines
4.0 KiB
Matlab

function [Params] = setupParameters(this)
CONSTANTS = Helper.PhysicsConstants;
hbar = CONSTANTS.PlanckConstantReduced; % [J.s]
kbol = CONSTANTS.BoltzmannConstant; % [J/K]
mu0 = CONSTANTS.VacuumPermeability; % [N/A^2]
muB = CONSTANTS.BohrMagneton; % [J/T]
a0 = CONSTANTS.BohrRadius; % [m]
m0 = CONSTANTS.AtomicMassUnit; % [kg]
w0 = 2*pi*61.6582; % Angular frequency unit [s^-1]
mu0factor = 0.3049584233607396; % =(m0/me)*pi*alpha^2 -- me=mass of electron, alpha=fine struct. const.
% mu0=mu0factor *hbar^2*a0/(m0*muB^2)
% Number of points in each direction
Params.Nx = this.NumberOfGridPoints(1);
Params.Ny = this.NumberOfGridPoints(2);
% Dimensions (in units of l0)
Params.Lx = this.Dimensions(1);
Params.Ly = this.Dimensions(2);
% Mass, length scale
Params.m = CONSTANTS.Dy164Mass;
l0 = sqrt(hbar/(Params.m*w0)); % Defining a harmonic oscillator length
% Atom numbers
Params.N = this.NumberOfAtoms;
% Dipole angle
Params.theta = this.DipolarPolarAngle; % pi/2 dipoles along x, theta=0 dipoles along z
Params.phi = this.DipolarAzimuthAngle;
% Dipole lengths (units of muB)
Params.mu = CONSTANTS.DyMagneticMoment;
% Scattering lengths
Params.as = this.ScatteringLength*a0;
% Trapping frequencies
Params.wx = 2*pi*this.TrapFrequencies(1);
Params.wy = 2*pi*this.TrapFrequencies(2);
Params.wz = 2*pi*this.TrapFrequencies(3);
% Tolerances
Params.Etol = this.EnergyTolerance;
Params.rtol = this.ResidualTolerance;
Params.sim_time_cut_off = this.TimeCutOff; % sometimes the imaginary time gets a little stuck
% even though the solution is good, this just stops it going on forever
Params.mindt = this.MinimumTimeStepSize; % Minimum size for a time step using adaptive dt
Params.njob = this.JobNumber;
% ================ Variational method parameters ================ %
% FMinCon Settings
Params.SelfConIter = this.MaxIterations; % Max number of iterations to perform self-consistent calculation
Params.ell = this.VariationalWidth; % initial [ell], ell is the "width" - psi ~ e^(z^2/ell^2)
% Window of optimization
Params.ell_lower = this.WidthLowerBound;
Params.ell_upper = this.WidthUpperBound;
% Relative cutoffs
Params.ellcutoff = this.WidthCutoff;
% ================ Parameters defined by those above ================ %
% Contact interaction strength (units of l0/m)
Params.gs = 4*pi*Params.as/l0;
% Dipole lengths
Params.add = mu0*Params.mu^2*Params.m/(12*pi*hbar^2);
% DDI strength
Params.gdd = 12*pi*Params.add/l0; %sometimes the 12 is a 4 --> depends on how Vdk (DDI) is defined
% Trap gamma
Params.gx = (Params.wx/w0)^2;
Params.gy = (Params.wy/w0)^2;
Params.gz = (Params.wz/w0)^2;
% == Calculate LHY correction to account for quantum fluctuations == %
eps_dd = Params.add/Params.as;
if eps_dd == 0
Q5 = 1;
elseif eps_dd == 1
Q5 = 3*sqrt(3)/2;
else
yeps = (1-eps_dd)/(3*eps_dd);
Q5 = (3*eps_dd)^(5/2)*( (8+26*yeps+33*yeps^2)*sqrt(1+yeps) + 15*yeps^3*log((1+sqrt(1+yeps))/sqrt(yeps)) )/48;
Q5 = real(Q5);
end
Params.gammaQF = 128/3*sqrt(pi*(Params.as/l0)^5)*Q5;
% Loading the rest into Params
Params.hbar = hbar;
Params.kbol = kbol;
Params.mu0 = mu0;
Params.muB = muB;
Params.a0 = a0;
Params.w0 = w0;
Params.l0 = l0;
end