Added residual tolerance.

This commit is contained in:
Karthik 2024-06-20 12:16:42 +02:00
parent cafe7eeb86
commit fd5fcf6c8d
7 changed files with 65 additions and 22 deletions

View File

@ -6,10 +6,20 @@ function visualizeGSWavefunction(run_index)
Data = load(sprintf(horzcat(folder_path, '/Run_%03i/psi_gs.mat'),run_index),'psi','Params','Transf','Observ'); Data = load(sprintf(horzcat(folder_path, '/Run_%03i/psi_gs.mat'),run_index),'psi','Params','Transf','Observ');
psi = Data.psi;
Params = Data.Params; Params = Data.Params;
Transf = Data.Transf; Transf = Data.Transf;
Observ = Data.Observ; Observ = Data.Observ;
if isgpuarray(Data.psi)
psi = gather(Data.psi);
else
psi = Data.psi;
end
if isgpuarray(Data.Observ.residual)
Observ.residual = gather(Data.Observ.residual);
else
Observ.residual = Data.Observ.residual;
end
format long format long
x = Transf.x*Params.l0*1e6; x = Transf.x*Params.l0*1e6;
@ -18,8 +28,12 @@ function visualizeGSWavefunction(run_index)
dx = x(2)-x(1); dy = y(2)-y(1); dz = z(2)-z(1); dx = x(2)-x(1); dy = y(2)-y(1); dz = z(2)-z(1);
%Plotting %Plotting
figure('Position', [100, 100, 1600, 900]);
clf
subplot(2,3,1) % Subplot 1
% subplot(2,3,1)
subplot('Position', [0.05, 0.55, 0.28, 0.4])
n = abs(psi).^2; n = abs(psi).^2;
nxz = squeeze(trapz(n*dy,2)); nxz = squeeze(trapz(n*dy,2));
nyz = squeeze(trapz(n*dx,1)); nyz = squeeze(trapz(n*dx,1));
@ -27,27 +41,46 @@ function visualizeGSWavefunction(run_index)
plotxz = pcolor(x,z,nxz'); plotxz = pcolor(x,z,nxz');
set(plotxz, 'EdgeColor', 'none'); set(plotxz, 'EdgeColor', 'none');
xlabel('$x$ [$\mu$m]'); ylabel('$z$ [$\mu$m]'); colorbar
xlabel('$x$ [$\mu$m]', 'FontSize', 14); ylabel('$z$ [$\mu$m]', 'FontSize', 14);
title('$|\Psi_{zx}|^2$', 'FontSize', 14);
subplot(2,3,2) % Subplot 2
% subplot(2,3,2)
subplot('Position', [0.36, 0.55, 0.28, 0.4])
plotyz = pcolor(y,z,nyz'); plotyz = pcolor(y,z,nyz');
set(plotyz, 'EdgeColor', 'none'); set(plotyz, 'EdgeColor', 'none');
xlabel('$y$ [$\mu$m]'); ylabel('$z$ [$\mu$m]'); colorbar
xlabel('$y$ [$\mu$m]', 'FontSize', 14); ylabel('$z$ [$\mu$m]', 'FontSize', 14);
title('$|\Psi_{zy}|^2$', 'FontSize', 14);
subplot(2,3,3) % Subplot 3
% subplot(2,3,3)
subplot('Position', [0.67, 0.55, 0.28, 0.4]);
plotxy = pcolor(x,y,nxy'); plotxy = pcolor(x,y,nxy');
set(plotxy, 'EdgeColor', 'none'); set(plotxy, 'EdgeColor', 'none');
xlabel('$x$ [$\mu$m]'); ylabel('$y$ [$\mu$m]'); colorbar
xlabel('$x$ [$\mu$m]', 'FontSize', 14); ylabel('$y$ [$\mu$m]', 'FontSize', 14);
title('$ |\Psi_{xy}|^2$', 'FontSize', 14);
subplot(2,3,4) % Subplot 4
% subplot(2,3,4)
subplot('Position', [0.05, 0.05, 0.26, 0.4]);
plot(-log10(Observ.residual),'-b') plot(-log10(Observ.residual),'-b')
ylabel('$-\mathrm{log}_{10}(r)$'); xlabel('steps'); ylabel('$-\mathrm{log}_{10}(r)$', 'FontSize', 14); xlabel('Time steps', 'FontSize', 14);
title('Residual', 'FontSize', 14);
subplot(2,3,5) % Subplot 5
% subplot(2, 3, 5);
subplot('Position', [0.36, 0.05, 0.26, 0.4]);
plot(Observ.EVec,'-b') plot(Observ.EVec,'-b')
ylabel('$E$'); xlabel('steps'); ylabel('$E$', 'FontSize', 14); xlabel('Time steps', 'FontSize', 14);
title('Total Energy', 'FontSize', 14);
subplot(2,3,6) % Subplot 6
% subplot(2, 3, 6);
subplot('Position', [0.67, 0.05, 0.26, 0.4]);
plot(Observ.mucVec,'-b') plot(Observ.mucVec,'-b')
ylabel('$\mu$'); xlabel('steps'); ylabel('$\mu$', 'FontSize', 14); xlabel('Time steps', 'FontSize', 14);
title('Chemical Potential', 'FontSize', 14);
end end

View File

@ -23,6 +23,7 @@ OptionsStruct.SimulationMode = 'ImaginaryTimeEvolution'; % 'ImaginaryTi
OptionsStruct.TimeStepSize = 50E-6; % in s OptionsStruct.TimeStepSize = 50E-6; % in s
OptionsStruct.NumberOfTimeSteps = 200; % in s OptionsStruct.NumberOfTimeSteps = 200; % in s
OptionsStruct.EnergyTolerance = 5E-10; OptionsStruct.EnergyTolerance = 5E-10;
OptionsStruct.ResidualTolerance = 1E-05;
OptionsStruct.JobNumber = 1; OptionsStruct.JobNumber = 1;
OptionsStruct.SaveData = true; OptionsStruct.SaveData = true;
@ -43,3 +44,5 @@ Plotter.visualizeSpace(Transf)
Plotter.visualizeTrapPotential(V,Params,Transf) Plotter.visualizeTrapPotential(V,Params,Transf)
%% - Plot initial wavefunction %% - Plot initial wavefunction
Plotter.visualizeWavefunction(psi,Params,Transf) Plotter.visualizeWavefunction(psi,Params,Transf)
%% - Plot GS wavefunction
Plotter.visualizeGSWavefunction(Params.njob)

View File

@ -16,13 +16,14 @@ OptionsStruct.TrapDepth = 5;
OptionsStruct.BoxSize = 15; OptionsStruct.BoxSize = 15;
OptionsStruct.TrapPotentialType = 'Harmonic'; OptionsStruct.TrapPotentialType = 'Harmonic';
OptionsStruct.NumberOfGridPoints = [128, 256, 128]; OptionsStruct.NumberOfGridPoints = [256, 512, 256];
OptionsStruct.Dimensions = [50, 120, 150]; OptionsStruct.Dimensions = [50, 120, 150];
OptionsStruct.CutoffType = 'Cylindrical'; OptionsStruct.CutoffType = 'Cylindrical';
OptionsStruct.SimulationMode = 'ImaginaryTimeEvolution'; % 'ImaginaryTimeEvolution' | 'RealTimeEvolution' OptionsStruct.SimulationMode = 'ImaginaryTimeEvolution'; % 'ImaginaryTimeEvolution' | 'RealTimeEvolution'
OptionsStruct.TimeStepSize = 500E-6; % in s OptionsStruct.TimeStepSize = 500E-6; % in s
OptionsStruct.NumberOfTimeSteps = 2E6; % in s OptionsStruct.NumberOfTimeSteps = 2E6; % in s
OptionsStruct.EnergyTolerance = 5E-10; OptionsStruct.EnergyTolerance = 5E-10;
OptionsStruct.ResidualTolerance = 1E-05;
OptionsStruct.JobNumber = 1; OptionsStruct.JobNumber = 1;
OptionsStruct.RunOnGPU = true; OptionsStruct.RunOnGPU = true;

View File

@ -14,6 +14,7 @@ classdef DipolarGas < handle & matlab.mixin.Copyable
TimeStepSize; TimeStepSize;
NumberOfTimeSteps; NumberOfTimeSteps;
EnergyTolerance; EnergyTolerance;
ResidualTolerance;
MinimumTimeStepSize; MinimumTimeStepSize;
Calculator; Calculator;
@ -59,6 +60,8 @@ classdef DipolarGas < handle & matlab.mixin.Copyable
@(x) assert(isnumeric(x) && isscalar(x) && (x > 0))); @(x) assert(isnumeric(x) && isscalar(x) && (x > 0)));
addParameter(p, 'EnergyTolerance', 1e-10,... addParameter(p, 'EnergyTolerance', 1e-10,...
@(x) assert(isnumeric(x) && isscalar(x) && (x > 0))); @(x) assert(isnumeric(x) && isscalar(x) && (x > 0)));
addParameter(p, 'ResidualTolerance', 1e-10,...
@(x) assert(isnumeric(x) && isscalar(x) && (x > 0)));
addParameter(p, 'MinimumTimeStepSize', 1e-6,... addParameter(p, 'MinimumTimeStepSize', 1e-6,...
@(x) assert(isnumeric(x) && isscalar(x) && (x > 0))); @(x) assert(isnumeric(x) && isscalar(x) && (x > 0)));
addParameter(p, 'JobNumber', 1,... addParameter(p, 'JobNumber', 1,...
@ -86,6 +89,7 @@ classdef DipolarGas < handle & matlab.mixin.Copyable
this.TimeStepSize = p.Results.TimeStepSize; this.TimeStepSize = p.Results.TimeStepSize;
this.NumberOfTimeSteps = p.Results.NumberOfTimeSteps; this.NumberOfTimeSteps = p.Results.NumberOfTimeSteps;
this.EnergyTolerance = p.Results.EnergyTolerance; this.EnergyTolerance = p.Results.EnergyTolerance;
this.ResidualTolerance = p.Results.ResidualTolerance;
this.MinimumTimeStepSize = p.Results.MinimumTimeStepSize; this.MinimumTimeStepSize = p.Results.MinimumTimeStepSize;
this.JobNumber = p.Results.JobNumber; this.JobNumber = p.Results.JobNumber;

View File

@ -45,10 +45,12 @@ function [Params] = setupParameters(this)
Params.gamma_S = 7.5*10^(-3); % gamma for the stochastic GPE Params.gamma_S = 7.5*10^(-3); % gamma for the stochastic GPE
Params.muchem = 12.64*Params.wz/w0; % fixing the chemical potential for the stochastic GPE Params.muchem = 12.64*Params.wz/w0; % fixing the chemical potential for the stochastic GPE
Params.Etol = this.EnergyTolerance; % Tolerances % Tolerances
Params.Etol = this.EnergyTolerance;
Params.rtol = this.ResidualTolerance;
Params.sim_time_cut_off = this.NumberOfTimeSteps; % sometimes the imaginary time gets a little stuck Params.sim_time_cut_off = this.NumberOfTimeSteps; % sometimes the imaginary time gets a little stuck
% even though the solution is good, this just stops it going on forever % 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.mindt = this.MinimumTimeStepSize; % Minimum size for a time step using adaptive dt
Params.njob = this.JobNumber; Params.njob = this.JobNumber;

View File

@ -8,7 +8,7 @@
#SBATCH --cpus-per-task=10 #SBATCH --cpus-per-task=10
#SBATCH --mem=24G #SBATCH --mem=24G
# Estimated wallclock time for job # Estimated wallclock time for job
#SBATCH --time=15:00:00 #SBATCH --time=12:00:00
#SBATCH --job-name=simulation #SBATCH --job-name=simulation
#SBATCH --error=simulation.err #SBATCH --error=simulation.err
#SBATCH --output=simulation.out #SBATCH --output=simulation.out

View File

@ -8,7 +8,7 @@
#SBATCH --gres=gpu:4 #SBATCH --gres=gpu:4
#SBATCH --mem=24G #SBATCH --mem=24G
# Estimated wallclock time for job # Estimated wallclock time for job
#SBATCH --time=15:00:00 #SBATCH --time=12:00:00
#SBATCH --job-name=simulation #SBATCH --job-name=simulation
#SBATCH --error=simulation.err #SBATCH --error=simulation.err
#SBATCH --output=simulation.out #SBATCH --output=simulation.out