diff --git a/Dipolar-Gas-Simulator/+Scripts/run_on_cluster.m b/Dipolar-Gas-Simulator/+Scripts/run_on_cluster.m index d100b79..558dc31 100644 --- a/Dipolar-Gas-Simulator/+Scripts/run_on_cluster.m +++ b/Dipolar-Gas-Simulator/+Scripts/run_on_cluster.m @@ -5,7 +5,7 @@ OptionsStruct = struct; OptionsStruct.NumberOfAtoms = 5E5; OptionsStruct.DipolarPolarAngle = 0; OptionsStruct.DipolarAzimuthAngle = 0; -OptionsStruct.ScatteringLength = 81.0; +OptionsStruct.ScatteringLength = 78.0; AspectRatio = 2.8; HorizontalTrapFrequency = 125; @@ -15,7 +15,7 @@ OptionsStruct.TrapPotentialType = 'Harmonic'; OptionsStruct.NumberOfGridPoints = [256, 256]; OptionsStruct.Dimensions = [18, 18]; -OptionsStruct.TimeStepSize = 0.001; % in s +OptionsStruct.TimeStepSize = 1E-3; % in s OptionsStruct.MinimumTimeStepSize = 1E-5; % in s OptionsStruct.TimeCutOff = 2E6; % in s OptionsStruct.EnergyTolerance = 5E-10; @@ -23,7 +23,7 @@ OptionsStruct.ResidualTolerance = 1E-04; OptionsStruct.NoiseScaleFactor = 0.05; OptionsStruct.MaxIterations = 10; -OptionsStruct.VariationalWidth = 1.8; +OptionsStruct.VariationalWidth = 2.0; OptionsStruct.WidthLowerBound = 0.01; OptionsStruct.WidthUpperBound = 12; OptionsStruct.WidthCutoff = 1e-2; @@ -43,14 +43,14 @@ solver.Potential = pot.trap(); %-% Run Solver %-% [Params, Transf, psi, V, VDk] = solver.run(); -%% AR = 4 +%% AR = 4.0 OptionsStruct = struct; OptionsStruct.NumberOfAtoms = 5E5; OptionsStruct.DipolarPolarAngle = 0; OptionsStruct.DipolarAzimuthAngle = 0; -OptionsStruct.ScatteringLength = 81.0; +OptionsStruct.ScatteringLength = 78.0; AspectRatio = 4.0; HorizontalTrapFrequency = 125; @@ -60,7 +60,7 @@ OptionsStruct.TrapPotentialType = 'Harmonic'; OptionsStruct.NumberOfGridPoints = [256, 256]; OptionsStruct.Dimensions = [18, 18]; -OptionsStruct.TimeStepSize = 0.001; % in s +OptionsStruct.TimeStepSize = 1E-3; % in s OptionsStruct.MinimumTimeStepSize = 1E-5; % in s OptionsStruct.TimeCutOff = 2E6; % in s OptionsStruct.EnergyTolerance = 5E-10; @@ -68,7 +68,7 @@ OptionsStruct.ResidualTolerance = 1E-04; OptionsStruct.NoiseScaleFactor = 0.05; OptionsStruct.MaxIterations = 10; -OptionsStruct.VariationalWidth = 1.5; +OptionsStruct.VariationalWidth = 2.0; OptionsStruct.WidthLowerBound = 0.01; OptionsStruct.WidthUpperBound = 12; OptionsStruct.WidthCutoff = 1e-2; diff --git a/Dipolar-Gas-Simulator/+Scripts/run_on_cluster_optimal_system_size.m b/Dipolar-Gas-Simulator/+Scripts/run_on_cluster_optimal_system_size.m index 580802c..37a4df6 100644 --- a/Dipolar-Gas-Simulator/+Scripts/run_on_cluster_optimal_system_size.m +++ b/Dipolar-Gas-Simulator/+Scripts/run_on_cluster_optimal_system_size.m @@ -31,6 +31,7 @@ for Lx = 2:0.4:6 OptionsStruct.EnergyTolerance = 5E-10; OptionsStruct.ResidualTolerance = 1E-05; OptionsStruct.NoiseScaleFactor = 0.05; + OptionsStruct.IncludeDDICutOff = false; OptionsStruct.MaxIterations = 10; OptionsStruct.VariationalWidth = 1.00; diff --git a/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVDkWithoutCutoff.m b/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVDkWithoutCutoff.m new file mode 100644 index 0000000..4f15d39 --- /dev/null +++ b/Dipolar-Gas-Simulator/+VariationalSolver2D/@Calculator/calculateVDkWithoutCutoff.m @@ -0,0 +1,19 @@ +function VDk = calculateVDkWithoutCutoff(~, Transf, Params, ell) +% == Calculating the DDI potential in Fourier space with appropriate cutoff == % + + % Interaction in K space + QX = Transf.KX*ell/sqrt(2); + QY = Transf.KY*ell/sqrt(2); + + Qsq = QX.^2 + QY.^2; + absQ = sqrt(Qsq); + QDsq = QX.^2*cos(Params.eta)^2 + QY.^2*sin(Params.eta)^2; + + % Bare interaction + Fpar = -1 + 3*sqrt(pi)*QDsq.*erfcx(absQ)./absQ; % Scaled complementary error function erfcx(x) = e^(x^2) * erfc(x) + Fperp = 2 - 3*sqrt(pi).*absQ.*erfcx(absQ); + Fpar(absQ == 0) = -1; + + % Full DDI + VDk = Fpar*sin(Params.theta)^2 + Fperp*cos(Params.theta)^2; +end \ No newline at end of file diff --git a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/DipolarGas.m b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/DipolarGas.m index 8790e49..9bccb66 100644 --- a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/DipolarGas.m +++ b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/DipolarGas.m @@ -28,6 +28,7 @@ classdef DipolarGas < handle & matlab.mixin.Copyable Calculator; SimulationParameters; + IncludeDDICutOff; PlotLive; JobNumber; RunOnGPU; @@ -84,13 +85,15 @@ classdef DipolarGas < handle & matlab.mixin.Copyable @(x) assert(isnumeric(x) && isscalar(x) && (x > 0))); addParameter(p, 'JobNumber', 0,... @(x) assert(isnumeric(x) && isscalar(x) && (x >= 0))); + addParameter(p, 'IncludeDDICutOff', true,... + @islogical); addParameter(p, 'PlotLive', false,... @islogical); addParameter(p, 'RunOnGPU', false,... @islogical); - addParameter(p, 'DebugMode', false,... + addParameter(p, 'DebugMode', false,... @islogical); - addParameter(p, 'SaveData', false,... + addParameter(p, 'SaveData', false,... @islogical); addParameter(p, 'SaveDirectory', './Data',... @ischar); @@ -119,6 +122,7 @@ classdef DipolarGas < handle & matlab.mixin.Copyable this.VariationalWidthTolerance = p.Results.VariationalWidthTolerance; this.VariationalEnergyTolerance = p.Results.VariationalEnergyTolerance; + this.IncludeDDICutOff = p.Results.IncludeDDICutOff; this.PlotLive = p.Results.PlotLive; this.JobNumber = p.Results.JobNumber; this.RunOnGPU = p.Results.RunOnGPU; diff --git a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m index 37dc116..6e63bee 100644 --- a/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m +++ b/Dipolar-Gas-Simulator/+VariationalSolver2D/@DipolarGas/initialize.m @@ -4,7 +4,11 @@ function [psi,V,VDk] = initialize(this,Params,VParams,Transf) assert(~anynan(V), 'Potential not defined! Specify as .Potential = .trap() + .'); % == Calculating the DDIs == % - VDk = this.Calculator.calculateVDkWithCutoff(Transf, Params, VParams.ell); + if this.IncludeDDICutOff + VDk = this.Calculator.calculateVDkWithCutoff(Transf, Params, VParams.ell); + else + VDk = this.Calculator.calculateVDkWithoutCutoff(Transf, Params, VParams.ell); + end % == Setting up the initial wavefunction == % psi = this.setupWavefunction(Params,Transf);