karthik/Calculations
1
0
Fork 0
Code Issues Pull Requests repo.project_board Releases Wiki Activity

Compare commits

base: karthik:c2e6420243e31a3385f4b9a311aa7ca60f0eae9c
Branches Tags
karthik:master
...
compare: karthik:de683c7d6b23129e7efd4ea90667d2d2f22a19fd
Branches Tags
karthik:master

5 Commits
c2e6420243 ... de683c7d6b

Author SHA1 Message Date
Karthik Chandrashekara
de683c7d6b Minor edits to reflect changes like the new Push beam flag. 2021-07-11 21:17:49 +02:00
Karthik Chandrashekara
1788066bce Modified to reflect change in acceleration due to scattering function. 2021-07-11 21:16:53 +02:00
Karthik Chandrashekara
3a6732a4a6 Corrected scattering rate formula and removed detuning as an argument to be passed. 2021-07-11 21:16:07 +02:00
Karthik Chandrashekara
a79389c608 Cosmetic changes and included the new Push beam flag. 2021-07-11 21:14:31 +02:00
Karthik Chandrashekara
cbef467b6f Added a flag to turn on or off the Push beam. 2021-07-11 21:13:01 +02:00
5 changed files with 37 additions and 36 deletions
Show Stats Expand all files Collapse all files

4
MOT Capture Process Simulation/@MOTSimulator/MOTSimulator.m
View File

@ -78,6 +78,7 @@ classdef MOTSimulator < handle & matlab.mixin.Copyable
%Flags %Flags
SpontaneousEmission; SpontaneousEmission;
Sideband; Sideband;
PushBeam;
ZeemanSlowerBeam; ZeemanSlowerBeam;
Gravity; Gravity;
BackgroundCollision; BackgroundCollision;
@ -120,6 +121,8 @@ classdef MOTSimulator < handle & matlab.mixin.Copyable
@islogical); @islogical);
addParameter(p, 'Sideband', false,... addParameter(p, 'Sideband', false,...
@islogical); @islogical);
addParameter(p, 'PushBeam', false,...
@islogical);
addParameter(p, 'ZeemanSlowerBeam', false,... addParameter(p, 'ZeemanSlowerBeam', false,...
@islogical); @islogical);
addParameter(p, 'Gravity', false,... addParameter(p, 'Gravity', false,...
@ -142,6 +145,7 @@ classdef MOTSimulator < handle & matlab.mixin.Copyable
s.SpontaneousEmission = p.Results.SpontaneousEmission; s.SpontaneousEmission = p.Results.SpontaneousEmission;
s.Sideband = p.Results.Sideband; s.Sideband = p.Results.Sideband;
s.PushBeam = p.Results.PushBeam;
s.ZeemanSlowerBeam = p.Results.ZeemanSlowerBeam; s.ZeemanSlowerBeam = p.Results.ZeemanSlowerBeam;
s.Gravity = p.Results.Gravity; s.Gravity = p.Results.Gravity;
s.BackgroundCollision = p.Results.BackgroundCollision; s.BackgroundCollision = p.Results.BackgroundCollision;

2
MOT Capture Process Simulation/@MOTSimulator/accelerationDueToPushBeam.m
View File

@ -14,7 +14,7 @@ function ret = accelerationDueToPushBeam(this, PositionVector, VelocityVector)
a_push = (Helper.PhysicsConstants.PlanckConstantReduced * this.PushBeamWaveNumber * WaveVectorEndPoint(1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.PushBeamLinewidth * 0.5) .* (s_push/(1+s_push)); a_push = (Helper.PhysicsConstants.PlanckConstantReduced * this.PushBeamWaveNumber * WaveVectorEndPoint(1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.PushBeamLinewidth * 0.5) .* (s_push/(1+s_push));
if this.SpontaneousEmission if this.SpontaneousEmission
a_scatter = this.accelerationDueToSpontaneousEmissionProcess(s_push, s_push, Detuning, this.PushBeamLinewidth, this.PushBeamWaveNumber); a_scatter = this.accelerationDueToSpontaneousEmissionProcess(s_push, s_push, this.PushBeamLinewidth, this.PushBeamWaveNumber);
else else
a_scatter = [0,0,0]; a_scatter = [0,0,0];
end end

4
MOT Capture Process Simulation/@MOTSimulator/accelerationDueToSpontaneousEmissionProcess.m
View File

@ -1,8 +1,8 @@
function ret = accelerationDueToSpontaneousEmissionProcess(this, SaturationParameter, TotalSaturationParameter, Detuning, Linewidth, WaveNumber) function ret = accelerationDueToSpontaneousEmissionProcess(this, SaturationParameter, TotalSaturationParameter, Linewidth, WaveNumber)
Vector = [2*rand(1)-1,2*rand(1)-1,2*rand(1)-1]; Vector = [2*rand(1)-1,2*rand(1)-1,2*rand(1)-1];
Vector = Vector./norm(Vector); Vector = Vector./norm(Vector);
ScatteringRate = (0.5 * SaturationParameter * Linewidth) / (1 + TotalSaturationParameter + (4 * (Detuning / Linewidth)^2)); ScatteringRate = (0.5 * SaturationParameter * Linewidth) / (1 + TotalSaturationParameter);
NumberOfScatteringEvents = floor(ScatteringRate * this.TimeStep); NumberOfScatteringEvents = floor(ScatteringRate * this.TimeStep);
if NumberOfScatteringEvents > 0 if NumberOfScatteringEvents > 0

52
MOT Capture Process Simulation/@MOTSimulator/calculateTotalAcceleration.m
View File

@ -10,12 +10,11 @@ function ret = calculateTotalAcceleration(this, PositionVector, VelocityVector)
Origin = [0,0,0]; Origin = [0,0,0];
% Calculate the Saturation Intensity at the specified point along its Gaussian Profile % Calculate the Saturation Intensity at the specified point along its Gaussian Profile
CoolingBeamLocalSaturationIntensity = [this.calculateLocalSaturationIntensity(1, PositionVector, Origin, WaveVectorEndPoint(1,:), this.CoolingBeamRadius, this.CoolingBeamWaist), ... CoolingBeamLocalSaturationIntensity = [this.calculateLocalSaturationIntensity(0.25 * this.CoolingBeamSaturationParameter, PositionVector, Origin, WaveVectorEndPoint(1,:), this.CoolingBeamRadius, this.CoolingBeamWaist), ...
this.calculateLocalSaturationIntensity(1, PositionVector, Origin, WaveVectorEndPoint(2,:), this.CoolingBeamRadius, this.CoolingBeamWaist)]; this.calculateLocalSaturationIntensity(0.25 * this.CoolingBeamSaturationParameter, PositionVector, Origin, WaveVectorEndPoint(2,:), this.CoolingBeamRadius, this.CoolingBeamWaist)];
SidebandLocalSaturationIntensity = [this.calculateLocalSaturationIntensity(1, PositionVector, Origin, WaveVectorEndPoint(1,:), this.SidebandBeamRadius, this.SidebandBeamWaist), ...
this.calculateLocalSaturationIntensity(1, PositionVector, Origin, WaveVectorEndPoint(2,:), this.SidebandBeamRadius, this.SidebandBeamWaist)];
SidebandLocalSaturationIntensity = [this.calculateLocalSaturationIntensity(0.25 * this.SidebandSaturationParameter, PositionVector, Origin, WaveVectorEndPoint(1,:), this.SidebandBeamRadius, this.SidebandBeamWaist), ...
this.calculateLocalSaturationIntensity(0.25 * this.SidebandSaturationParameter, PositionVector, Origin, WaveVectorEndPoint(2,:), this.SidebandBeamRadius, this.SidebandBeamWaist)];
TotalAcceleration = zeros(1,3); TotalAcceleration = zeros(1,3);
@ -26,29 +25,29 @@ function ret = calculateTotalAcceleration(this, PositionVector, VelocityVector)
LocalMagneticField = this.magneticFieldForMOT(PositionVector); LocalMagneticField = this.magneticFieldForMOT(PositionVector);
B = sign(LocalMagneticField(1:3) * WaveVectorEndPoint(i,1:3)') * LocalMagneticField(4); B = sign(dot(LocalMagneticField(1:3), WaveVectorEndPoint(i,:))) * LocalMagneticField(4);
ZeemanShift = this.LandegFactor * this.MagneticSubLevel * (Helper.PhysicsConstants.BohrMagneton / Helper.PhysicsConstants.PlanckConstantReduced) * B; ZeemanShift = this.LandegFactor * this.MagneticSubLevel * (Helper.PhysicsConstants.BohrMagneton / Helper.PhysicsConstants.PlanckConstantReduced) * B;
DopplerShift = (VelocityVector * WaveVectorEndPoint(i,1:3)') * this.CoolingBeamWaveNumber; DopplerShift = dot(WaveVectorEndPoint(i,:), VelocityVector) * this.CoolingBeamWaveNumber;
Delta_Cooling(i*2-1) = this.CoolingBeamDetuning + DopplerShift + ZeemanShift * Sigma(i); Delta_Cooling(i*2-1) = this.CoolingBeamDetuning + DopplerShift + (ZeemanShift * Sigma(i));
Delta_Cooling(i*2) = this.CoolingBeamDetuning - DopplerShift - ZeemanShift * Sigma(i); Delta_Cooling(i*2) = this.CoolingBeamDetuning - DopplerShift - (ZeemanShift * Sigma(i));
if this.Sideband if this.Sideband
Delta_Sideband(i*2-1) = this.SidebandDetuning + DopplerShift + ZeemanShift * Sigma(i); Delta_Sideband(i*2-1) = this.SidebandDetuning + DopplerShift + (ZeemanShift * Sigma(i));
Delta_Sideband(i*2) = this.SidebandDetuning - DopplerShift - ZeemanShift * Sigma(i); Delta_Sideband(i*2) = this.SidebandDetuning - DopplerShift - (ZeemanShift * Sigma(i));
end end
end end
SaturationParameter = [0,0,0,0,0,0,0,0]; SaturationParameter = [0,0,0,0,0,0,0,0];
for i = 1:2 for i = 1:2
SaturationParameter(2*i-1) = (0.25 * this.CoolingBeamSaturationParameter * CoolingBeamLocalSaturationIntensity(1)) /(1 + 4 * (Delta_Cooling(2*i-1)/this.CoolingBeamLinewidth)^2); SaturationParameter(2*i-1) = CoolingBeamLocalSaturationIntensity(1) /(1 + 4 * (Delta_Cooling(2*i-1)/this.CoolingBeamLinewidth)^2);
SaturationParameter(2*i) = (0.25 * this.CoolingBeamSaturationParameter * CoolingBeamLocalSaturationIntensity(1)) /(1 + 4 * (Delta_Cooling(2*i) /this.CoolingBeamLinewidth)^2); SaturationParameter(2*i) = CoolingBeamLocalSaturationIntensity(2) /(1 + 4 * (Delta_Cooling(2*i) /this.CoolingBeamLinewidth)^2);
if this.Sideband if this.Sideband
SaturationParameter(2*i-1+4) = (0.25 * this.SidebandSaturationParameter * SidebandLocalSaturationIntensity(1)) /(1 + 4 * (Delta_Sideband(2*i-1)/this.CoolingBeamLinewidth)^2); SaturationParameter(2*i-1+4) = SidebandLocalSaturationIntensity(1) /(1 + 4 * (Delta_Sideband(2*i-1)/this.CoolingBeamLinewidth)^2);
SaturationParameter(2*i+4) = (0.25 * this.SidebandSaturationParameter * SidebandLocalSaturationIntensity(2)) /(1 + 4 * (Delta_Sideband(2*i)/this.CoolingBeamLinewidth)^2); SaturationParameter(2*i+4) = SidebandLocalSaturationIntensity(2) /(1 + 4 * (Delta_Sideband(2*i)/this.CoolingBeamLinewidth)^2);
end end
end end
@ -56,28 +55,25 @@ function ret = calculateTotalAcceleration(this, PositionVector, VelocityVector)
for i = 1:2 for i = 1:2
a_1 = (Helper.PhysicsConstants.PlanckConstantReduced * this.CoolingBeamWaveNumber * WaveVectorEndPoint(i,1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.CoolingBeamLinewidth * 0.5) .* ... a_sat = (Helper.PhysicsConstants.PlanckConstantReduced * this.CoolingBeamWaveNumber * WaveVectorEndPoint(i,1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.CoolingBeamLinewidth * 0.5);
(SaturationParameter(2*i-1)/(1 + TotalSaturationParameter)); a_1 = a_sat .* (SaturationParameter(2*i-1)/(1 + TotalSaturationParameter));
a_2 = (Helper.PhysicsConstants.PlanckConstantReduced * this.CoolingBeamWaveNumber * WaveVectorEndPoint(i,1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.CoolingBeamLinewidth * 0.5) .* ... a_2 = a_sat .* (SaturationParameter(2*i) /(1 + TotalSaturationParameter));
(SaturationParameter(2*i) / (1 + TotalSaturationParameter));
if this.SpontaneousEmission if this.SpontaneousEmission
a_scattering = this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i-1), TotalSaturationParameter, Delta_Cooling(2*i-1), this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber) + ... a_scattering = this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i-1), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber) + ...
this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i), TotalSaturationParameter, Delta_Cooling(2*i), this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber); this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber);
else else
a_scattering = [0,0,0]; a_scattering = [0,0,0];
end end
if this.Sideband if this.Sideband
a_1 = a_1 + (Helper.PhysicsConstants.PlanckConstantReduced * this.CoolingBeamWaveNumber * WaveVectorEndPoint(i,1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.CoolingBeamLinewidth * 0.5) .* ... a_1 = a_1 + a_sat .* (SaturationParameter(2*i-1+4)/(1 + TotalSaturationParameter));
(SaturationParameter(2*i-1+4)/(1 + TotalSaturationParameter)); a_2 = a_2 + a_sat .* (SaturationParameter(2*i+4) /(1 + TotalSaturationParameter));
a_2 = a_2 + (Helper.PhysicsConstants.PlanckConstantReduced * this.CoolingBeamWaveNumber * WaveVectorEndPoint(i,1:3)/Helper.PhysicsConstants.Dy164Mass).*(this.CoolingBeamLinewidth * 0.5) .* ...
(SaturationParameter(2*i+4)/(1 + TotalSaturationParameter));
if this.SpontaneousEmission if this.SpontaneousEmission
a_scattering = a_scattering + ... a_scattering = a_scattering + ...
this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i-1+4), TotalSaturationParameter, Delta_Cooling(2*i-1), this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber) + ... this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i-1+4), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber) + ...
this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i+4), TotalSaturationParameter, Delta_Cooling(2*i), this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber); this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i+4), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber);
else else
a_scattering = [0,0,0]; a_scattering = [0,0,0];
end end
@ -86,7 +82,7 @@ function ret = calculateTotalAcceleration(this, PositionVector, VelocityVector)
TotalAcceleration = TotalAcceleration + (a_2 - a_1) + a_scattering; TotalAcceleration = TotalAcceleration + (a_2 - a_1) + a_scattering;
end end
if this.PushBeamRadius ~= 0 if this.PushBeam
TotalAcceleration = TotalAcceleration + this.accelerationDueToPushBeam(PositionVector, VelocityVector); TotalAcceleration = TotalAcceleration + this.accelerationDueToPushBeam(PositionVector, VelocityVector);
end end

7
MOT Capture Process Simulation/test_MOTSimulator.m
View File

@ -7,7 +7,8 @@ OptionsStruct.SimulationMode = '2D';
OptionsStruct.TimeStep = 50e-06; % in s OptionsStruct.TimeStep = 50e-06; % in s
OptionsStruct.SimulationTime = 4e-03; % in s OptionsStruct.SimulationTime = 4e-03; % in s
OptionsStruct.SpontaneousEmission = true; OptionsStruct.SpontaneousEmission = true;
OptionsStruct.Sideband = false; OptionsStruct.Sideband = true;
OptionsStruct.PushBeam = true;
OptionsStruct.Gravity = true; OptionsStruct.Gravity = true;
OptionsStruct.BackgroundCollision = true; OptionsStruct.BackgroundCollision = true;
OptionsStruct.SaveData = false; OptionsStruct.SaveData = false;
@ -115,13 +116,13 @@ QuantityOfInterestArray = LoadingRateArray;
OptionsStruct = struct; OptionsStruct = struct;
OptionsStruct.RescalingFactorForParameter = 1000; OptionsStruct.RescalingFactorForParameter = 1000;
OptionsStruct.XLabelString = 'Cooling Beam Power (mW)'; OptionsStruct.XLabelString = 'Cooling Beam Power (mW)';
OptionsStruct.RescalingFactorForYQuantity = 1e-10; OptionsStruct.RescalingFactorForYQuantity = 1e-11;
OptionsStruct.ErrorsForYQuantity = true; OptionsStruct.ErrorsForYQuantity = true;
OptionsStruct.ErrorsArray = StandardErrorArray; OptionsStruct.ErrorsArray = StandardErrorArray;
OptionsStruct.CIForYQuantity = true; OptionsStruct.CIForYQuantity = true;
OptionsStruct.CIArray = ConfidenceIntervalArray; OptionsStruct.CIArray = ConfidenceIntervalArray;
OptionsStruct.RemoveOutliers = true; OptionsStruct.RemoveOutliers = true;
OptionsStruct.YLabelString = 'Loading rate (x 10^{10} atoms/s)'; OptionsStruct.YLabelString = 'Loading rate (x 10^{11} atoms/s)';
OptionsStruct.TitleString = sprintf('Magnetic Gradient = %.0f (G/cm)', Simulator.MagneticGradient * 100); OptionsStruct.TitleString = sprintf('Magnetic Gradient = %.0f (G/cm)', Simulator.MagneticGradient * 100);
options = Helper.convertstruct2cell(OptionsStruct); options = Helper.convertstruct2cell(OptionsStruct);