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

Compare commits

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

No commits in common. "de683c7d6b23129e7efd4ea90667d2d2f22a19fd" and "c2e6420243e31a3385f4b9a311aa7ca60f0eae9c" have entirely different histories.
de683c7d6b ... c2e6420243

5 changed files with 36 additions and 37 deletions
Show Stats Expand all files Collapse all files

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

@ -78,7 +78,6 @@ classdef MOTSimulator < handle & matlab.mixin.Copyable
%Flags
SpontaneousEmission;
Sideband;
PushBeam;
ZeemanSlowerBeam;
Gravity;
BackgroundCollision;
@ -121,8 +120,6 @@ classdef MOTSimulator < handle & matlab.mixin.Copyable
@islogical);
addParameter(p, 'Sideband', false,...
@islogical);
addParameter(p, 'PushBeam', false,...
@islogical);
addParameter(p, 'ZeemanSlowerBeam', false,...
@islogical);
addParameter(p, 'Gravity', false,...
@ -145,7 +142,6 @@ classdef MOTSimulator < handle & matlab.mixin.Copyable
s.SpontaneousEmission = p.Results.SpontaneousEmission;
s.Sideband = p.Results.Sideband;
s.PushBeam = p.Results.PushBeam;
s.ZeemanSlowerBeam = p.Results.ZeemanSlowerBeam;
s.Gravity = p.Results.Gravity;
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));
if this.SpontaneousEmission
a_scatter = this.accelerationDueToSpontaneousEmissionProcess(s_push, s_push, this.PushBeamLinewidth, this.PushBeamWaveNumber);
a_scatter = this.accelerationDueToSpontaneousEmissionProcess(s_push, s_push, Detuning, this.PushBeamLinewidth, this.PushBeamWaveNumber);
else
a_scatter = [0,0,0];
end

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

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

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

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

7
MOT Capture Process Simulation/test_MOTSimulator.m
View File

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