Computes for specified edd and/or density.

This commit is contained in:
Karthik 2025-01-10 15:51:44 +01:00
parent 7be9f6ddd8
commit 6093c700ab

View File

@ -30,29 +30,27 @@ DyMagneticMoment = 9.93*BohrMagneton;
%% 2-D DDI Potential in k-space, with Gaussian ansatz width determined by constrained minimization %% 2-D DDI Potential in k-space, with Gaussian ansatz width determined by constrained minimization
w0 = 2*pi*61.6316; % Angular frequency unit [s^-1]
l0 = sqrt(PlanckConstantReduced/(Dy164Mass*w0)); % Defining a harmonic oscillator length
wz = 2 * pi * 300; % Trap frequency in the tight confinement direction wz = 2 * pi * 300; % Trap frequency in the tight confinement direction
lz = sqrt(PlanckConstantReduced/(Dy164Mass * wz)); % Defining a harmonic oscillator length lz = sqrt(PlanckConstantReduced/(Dy164Mass * wz)); % Defining a harmonic oscillator length
% Number of grid points in each direction % Number of grid points in each direction
Params.Nx = 128; Params.Nx = 128;
Params.Ny = 128; Params.Ny = 128;
% Dimensions (in units of l0)
w0 = 2*pi*61.6316; % Angular frequency unit [s^-1]
l0 = sqrt(PlanckConstantReduced/(Dy164Mass*w0)); % Defining a harmonic oscillator length
Params.Lx = 150*l0; Params.Lx = 150*l0;
Params.Ly = 150*l0; Params.Ly = 150*l0;
[Transf] = setupSpace(Params); [Transf] = setupSpace(Params);
nadd2s = 0.05:0.001:0.25; nadd2s = 0.110;
as_to_add = 0.74:0.001:0.79; as_to_add = 0.782;
Params.alpha = 10;
Params.phi = 0;
add = VacuumPermeability*DyMagneticMoment^2*Dy164Mass/(12*pi*PlanckConstantReduced^2); % Dipole length add = VacuumPermeability*DyMagneticMoment^2*Dy164Mass/(12*pi*PlanckConstantReduced^2); % Dipole length
gdd = VacuumPermeability*DyMagneticMoment^2/3; gdd = VacuumPermeability*DyMagneticMoment^2/3;
%%
var_widths = zeros(length(as_to_add), length(nadd2s));
x0 = 5; x0 = 5;
Aineq = []; Aineq = [];
Bineq = []; Bineq = [];
@ -63,42 +61,16 @@ ub = [10];
nonlcon = []; nonlcon = [];
fminconopts = optimoptions(@fmincon,'Display','off', 'StepTolerance', 1.0000e-11, 'MaxIterations',1500); fminconopts = optimoptions(@fmincon,'Display','off', 'StepTolerance', 1.0000e-11, 'MaxIterations',1500);
for idx = 1:length(nadd2s) AtomNumberDensity = nadd2s / add^2; % Areal density of atoms
for jdx = 1:length(as_to_add) as = as_to_add * add; % Scattering length
AtomNumberDensity = nadd2s(idx) / add^2; % Areal density of atoms eps_dd = add/as; % Relative interaction strength
as = (as_to_add(jdx) * add); % Scattering length gs = 4 * pi * PlanckConstantReduced^2/Dy164Mass * as; % Contact interaction strength
gs = 4 * pi * PlanckConstantReduced^2/Dy164Mass * as; % Contact interaction strength TotalEnergyPerParticle = @(x) computeTotalEnergyPerParticle(x, as, AtomNumberDensity, wz, lz, gs, add, gdd, PlanckConstantReduced);
TotalEnergyPerParticle = @(x) computeTotalEnergyPerParticle(x, as, AtomNumberDensity, wz, lz, gs, add, gdd, PlanckConstantReduced); sigma = fmincon(TotalEnergyPerParticle, x0, Aineq, Bineq, Aeq, Beq, lb, ub, nonlcon, fminconopts);
sigma = fmincon(TotalEnergyPerParticle, x0, Aineq, Bineq, Aeq, Beq, lb, ub, nonlcon, fminconopts);
var_widths(jdx, idx) = sigma;
end
end
figure(10)
clf
set(gcf,'Position',[50 50 950 750])
imagesc(nadd2s, as_to_add, var_widths); % Specify x and y data for axes
set(gca, 'YDir', 'normal'); % Correct the y-axis direction
colorbar; % Add a colorbar
xlabel('$na_{dd}^2$','fontsize',16,'interpreter','latex');
ylabel('$a_s/a_{dd}$','fontsize',16,'interpreter','latex');
%% Chosen values of interaction, density and tilt %% Chosen values of interaction, density and tilt
nadd = 0.110; MeanWidth = sigma * lz;
asadd = 0.782;
Params.alpha = 0;
Params.phi = 0;
[~, nadd2sidx] = min(abs(nadd2s - nadd));
[~, asaddidx] = min(abs(as_to_add - asadd));
AtomNumberDensity = nadd2s(nadd2sidx) / add^2; % Areal density of atoms
as = (as_to_add(asaddidx) * add); % Scattering length
eps_dd = add/as; % Relative interaction strength
gs = 4 * pi * PlanckConstantReduced^2/Dy164Mass * as; % Contact interaction strength
gdd = VacuumPermeability*DyMagneticMoment^2/3;
MeanWidth = var_widths(asaddidx, nadd2sidx)*lz;
% == 2-D DDI Potential in k-space == % % == 2-D DDI Potential in k-space == %
VDk = compute2DPotentialInMomentumSpace(Transf, Params, MeanWidth); VDk = compute2DPotentialInMomentumSpace(Transf, Params, MeanWidth);
@ -146,22 +118,10 @@ title(['$\theta = ',num2str(Params.alpha), '; \phi = ', num2str(Params.phi),'$']
%% %%
nadd = 0.110;
asadd = 0.782;
% Define values for alpha and phi % Define values for alpha and phi
alpha_values = 0:5:90; % Range of alpha values (you can modify this) alpha_values = 0:5:90; % Range of alpha values (you can modify this)
phi_values = 0:2:90; % Range of phi values (you can modify this) phi_values = 0:2:90; % Range of phi values (you can modify this)
[~, nadd2sidx] = min(abs(nadd2s - nadd));
[~, asaddidx] = min(abs(as_to_add - asadd));
AtomNumberDensity = nadd2s(nadd2sidx) / add^2; % Areal density of atoms
as = (as_to_add(asaddidx) * add); % Scattering length
eps_dd = add/as; % Relative interaction strength
gs = 4 * pi * PlanckConstantReduced^2/Dy164Mass * as; % Contact interaction strength
gdd = VacuumPermeability*DyMagneticMoment^2/3;
MeanWidth = var_widths(asaddidx, nadd2sidx)*lz;
% Set up VideoWriter object % Set up VideoWriter object
v = VideoWriter('potential_movie', 'MPEG-4'); % Create a video object v = VideoWriter('potential_movie', 'MPEG-4'); % Create a video object
v.FrameRate = 5; % Frame rate of the video v.FrameRate = 5; % Frame rate of the video