%% Tilting of the dipoles
% Atom Number      = 1250 ppum
% System size      = [sf * unitcell_x, sf * unitcell_x]

ppum                                      = 1250; % Atom Number Density in per micrometers

%% v_z = 500, theta = 0: a_s = 75.00

a                                     = 1.8058;
scalingfactor                         = 5;
lx                                    = scalingfactor*a;
ly                                    = scalingfactor*sqrt(3)*a;

% Initialize OptionsStruct
OptionsStruct                         = struct;

% Assign values to OptionsStruct
OptionsStruct.NumberOfAtoms           = ppum * (lx*ly);     
OptionsStruct.DipolarPolarAngle       = deg2rad(0);
OptionsStruct.DipolarAzimuthAngle     = 0;
OptionsStruct.ScatteringLength        = 75.00;           

OptionsStruct.TrapFrequencies         = [0, 0, 500];
OptionsStruct.TrapPotentialType       = 'None';

OptionsStruct.NumberOfGridPoints      = [128, 128];
OptionsStruct.Dimensions              = [lx, ly];      
OptionsStruct.TimeStepSize            = 1E-3;            % in s
OptionsStruct.MinimumTimeStepSize     = 1E-5;            % in s
OptionsStruct.TimeCutOff              = 2E6;             % in s
OptionsStruct.EnergyTolerance         = 5E-10;
OptionsStruct.ResidualTolerance       = 1E-05;
OptionsStruct.NoiseScaleFactor        = 0.05;
OptionsStruct.IncludeDDICutOff        = false;

OptionsStruct.MaxIterations           = 10;      
OptionsStruct.VariationalWidth        = 1.15;       
OptionsStruct.WidthLowerBound         = 0.01;
OptionsStruct.WidthUpperBound         = 12;
OptionsStruct.WidthCutoff             = 1e-2;

OptionsStruct.PlotLive                = false;
OptionsStruct.JobNumber               = 0;
OptionsStruct.RunOnGPU                = true;
OptionsStruct.SaveData                = true;
OptionsStruct.SaveDirectory           = './Results/Data_TiltingOfDipoles/AdjustedSystemSize/Hz500';
options                               = Helper.convertstruct2cell(OptionsStruct);
clear OptionsStruct

solver                                = VariationalSolver2D.DipolarGas(options{:});
pot                                   = VariationalSolver2D.Potentials(options{:});
solver.Potential                      = pot.trap();

%-% Run Solver %-%
[Params, Transf, psi, V, VDk]         = solver.run();

%% v_z = 500, theta = 5: a_s = 75.00

a                                     = 1.795;
scalingfactor                         = 5;
lx                                    = scalingfactor*a;
ly                                    = scalingfactor*sqrt(3)*a;

% Initialize OptionsStruct
OptionsStruct                         = struct;

% Assign values to OptionsStruct
OptionsStruct.NumberOfAtoms           = ppum * (lx*ly);     
OptionsStruct.DipolarPolarAngle       = deg2rad(5);
OptionsStruct.DipolarAzimuthAngle     = 0;
OptionsStruct.ScatteringLength        = 75.00;           

OptionsStruct.TrapFrequencies         = [0, 0, 500];
OptionsStruct.TrapPotentialType       = 'None';

OptionsStruct.NumberOfGridPoints      = [128, 128];
OptionsStruct.Dimensions              = [lx, ly];      
OptionsStruct.TimeStepSize            = 1E-3;            % in s
OptionsStruct.MinimumTimeStepSize     = 1E-5;            % in s
OptionsStruct.TimeCutOff              = 2E6;             % in s
OptionsStruct.EnergyTolerance         = 5E-10;
OptionsStruct.ResidualTolerance       = 1E-05;
OptionsStruct.NoiseScaleFactor        = 0.05;
OptionsStruct.IncludeDDICutOff        = false;

OptionsStruct.MaxIterations           = 10;      
OptionsStruct.VariationalWidth        = 1.15;       
OptionsStruct.WidthLowerBound         = 0.01;
OptionsStruct.WidthUpperBound         = 12;
OptionsStruct.WidthCutoff             = 1e-2;

OptionsStruct.PlotLive                = false;
OptionsStruct.JobNumber               = 1;
OptionsStruct.RunOnGPU                = true;
OptionsStruct.SaveData                = true;
OptionsStruct.SaveDirectory           = './Results/Data_TiltingOfDipoles/AdjustedSystemSize/Hz500';
options                               = Helper.convertstruct2cell(OptionsStruct);
clear OptionsStruct

solver                                = VariationalSolver2D.DipolarGas(options{:});
pot                                   = VariationalSolver2D.Potentials(options{:});
solver.Potential                      = pot.trap();

%-% Run Solver %-%
[Params, Transf, psi, V, VDk]         = solver.run();

%% v_z = 500, theta = 7.5: a_s = 75.00

a                                     = 2.055;
scalingfactor                         = 5;
lx                                    = scalingfactor*a;
ly                                    = scalingfactor*sqrt(3)*a;

% Initialize OptionsStruct
OptionsStruct                         = struct;

% Assign values to OptionsStruct
OptionsStruct.NumberOfAtoms           = ppum * (lx*ly);     
OptionsStruct.DipolarPolarAngle       = deg2rad(7.5);
OptionsStruct.DipolarAzimuthAngle     = 0;
OptionsStruct.ScatteringLength        = 75.00;           

OptionsStruct.TrapFrequencies         = [0, 0, 500];
OptionsStruct.TrapPotentialType       = 'None';

OptionsStruct.NumberOfGridPoints      = [128, 128];
OptionsStruct.Dimensions              = [lx, ly];      
OptionsStruct.TimeStepSize            = 1E-3;            % in s
OptionsStruct.MinimumTimeStepSize     = 1E-5;            % in s
OptionsStruct.TimeCutOff              = 2E6;             % in s
OptionsStruct.EnergyTolerance         = 5E-10;
OptionsStruct.ResidualTolerance       = 1E-05;
OptionsStruct.NoiseScaleFactor        = 0.05;
OptionsStruct.IncludeDDICutOff        = false;

OptionsStruct.MaxIterations           = 10;      
OptionsStruct.VariationalWidth        = 1.15;       
OptionsStruct.WidthLowerBound         = 0.01;
OptionsStruct.WidthUpperBound         = 12;
OptionsStruct.WidthCutoff             = 1e-2;

OptionsStruct.PlotLive                = false;
OptionsStruct.JobNumber               = 2;
OptionsStruct.RunOnGPU                = true;
OptionsStruct.SaveData                = true;
OptionsStruct.SaveDirectory           = './Results/Data_TiltingOfDipoles/AdjustedSystemSize/Hz500';
options                               = Helper.convertstruct2cell(OptionsStruct);
clear OptionsStruct

solver                                = VariationalSolver2D.DipolarGas(options{:});
pot                                   = VariationalSolver2D.Potentials(options{:});
solver.Potential                      = pot.trap();

%-% Run Solver %-%
[Params, Transf, psi, V, VDk]         = solver.run();
%% v_z = 500, theta = 10: a_s = 75.00

a                                     = 2.055;
scalingfactor                         = 5;
lx                                    = scalingfactor*a;
ly                                    = scalingfactor*sqrt(3)*a;

% Initialize OptionsStruct
OptionsStruct                         = struct;

% Assign values to OptionsStruct
OptionsStruct.NumberOfAtoms           = ppum * (lx*ly);     
OptionsStruct.DipolarPolarAngle       = deg2rad(10);
OptionsStruct.DipolarAzimuthAngle     = 0;
OptionsStruct.ScatteringLength        = 75.00;           

OptionsStruct.TrapFrequencies         = [0, 0, 500];
OptionsStruct.TrapPotentialType       = 'None';

OptionsStruct.NumberOfGridPoints      = [128, 128];
OptionsStruct.Dimensions              = [lx, ly];      
OptionsStruct.TimeStepSize            = 1E-3;            % in s
OptionsStruct.MinimumTimeStepSize     = 1E-5;            % in s
OptionsStruct.TimeCutOff              = 2E6;             % in s
OptionsStruct.EnergyTolerance         = 5E-10;
OptionsStruct.ResidualTolerance       = 1E-05;
OptionsStruct.NoiseScaleFactor        = 0.05;
OptionsStruct.IncludeDDICutOff        = false;

OptionsStruct.MaxIterations           = 10;      
OptionsStruct.VariationalWidth        = 1.15;       
OptionsStruct.WidthLowerBound         = 0.01;
OptionsStruct.WidthUpperBound         = 12;
OptionsStruct.WidthCutoff             = 1e-2;

OptionsStruct.PlotLive                = false;
OptionsStruct.JobNumber               = 3;
OptionsStruct.RunOnGPU                = true;
OptionsStruct.SaveData                = true;
OptionsStruct.SaveDirectory           = './Results/Data_TiltingOfDipoles/AdjustedSystemSize/Hz500';
options                               = Helper.convertstruct2cell(OptionsStruct);
clear OptionsStruct

solver                                = VariationalSolver2D.DipolarGas(options{:});
pot                                   = VariationalSolver2D.Potentials(options{:});
solver.Potential                      = pot.trap();

%-% Run Solver %-%
[Params, Transf, psi, V, VDk]         = solver.run();

%% v_z = 500, theta = 15: a_s = 75.00

a                                     = 2.0875;
scalingfactor                         = 5;
lx                                    = scalingfactor*a;
ly                                    = scalingfactor*sqrt(3)*a;

% Initialize OptionsStruct
OptionsStruct                         = struct;

% Assign values to OptionsStruct
OptionsStruct.NumberOfAtoms           = ppum * (lx*ly);     
OptionsStruct.DipolarPolarAngle       = deg2rad(15);
OptionsStruct.DipolarAzimuthAngle     = 0;
OptionsStruct.ScatteringLength        = 75.00;           

OptionsStruct.TrapFrequencies         = [0, 0, 500];
OptionsStruct.TrapPotentialType       = 'None';

OptionsStruct.NumberOfGridPoints      = [128, 128];
OptionsStruct.Dimensions              = [lx, ly];      
OptionsStruct.TimeStepSize            = 1E-3;            % in s
OptionsStruct.MinimumTimeStepSize     = 1E-5;            % in s
OptionsStruct.TimeCutOff              = 2E6;             % in s
OptionsStruct.EnergyTolerance         = 5E-10;
OptionsStruct.ResidualTolerance       = 1E-05;
OptionsStruct.NoiseScaleFactor        = 0.05;
OptionsStruct.IncludeDDICutOff        = false;

OptionsStruct.MaxIterations           = 10;      
OptionsStruct.VariationalWidth        = 1.15;       
OptionsStruct.WidthLowerBound         = 0.01;
OptionsStruct.WidthUpperBound         = 12;
OptionsStruct.WidthCutoff             = 1e-2;

OptionsStruct.PlotLive                = false;
OptionsStruct.JobNumber               = 4;
OptionsStruct.RunOnGPU                = true;
OptionsStruct.SaveData                = true;
OptionsStruct.SaveDirectory           = './Results/Data_TiltingOfDipoles/AdjustedSystemSize/Hz500';
options                               = Helper.convertstruct2cell(OptionsStruct);
clear OptionsStruct

solver                                = VariationalSolver2D.DipolarGas(options{:});
pot                                   = VariationalSolver2D.Potentials(options{:});
solver.Potential                      = pot.trap();

%-% Run Solver %-%
[Params, Transf, psi, V, VDk]         = solver.run();