From a79389c6089e39468fa44864bb0a75683c028a6d Mon Sep 17 00:00:00 2001 From: Karthik Chandrashekara Date: Sun, 11 Jul 2021 21:14:31 +0200 Subject: [PATCH] Cosmetic changes and included the new Push beam flag. --- .../calculateTotalAcceleration.m | 56 +++++++++---------- 1 file changed, 26 insertions(+), 30 deletions(-) diff --git a/MOT Capture Process Simulation/@MOTSimulator/calculateTotalAcceleration.m b/MOT Capture Process Simulation/@MOTSimulator/calculateTotalAcceleration.m index 8042411..fc620c5 100644 --- a/MOT Capture Process Simulation/@MOTSimulator/calculateTotalAcceleration.m +++ b/MOT Capture Process Simulation/@MOTSimulator/calculateTotalAcceleration.m @@ -10,13 +10,12 @@ 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(1, PositionVector, Origin, WaveVectorEndPoint(1,:), this.CoolingBeamRadius, this.CoolingBeamWaist), ... - this.calculateLocalSaturationIntensity(1, PositionVector, Origin, WaveVectorEndPoint(2,:), this.CoolingBeamRadius, this.CoolingBeamWaist)]; + 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)]; - 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); Delta_Cooling = [0,0,0,0]; @@ -26,58 +25,55 @@ function ret = calculateTotalAcceleration(this, PositionVector, VelocityVector) 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; - 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) = 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) = (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); + 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); 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+4) = (0.25 * this.SidebandSaturationParameter * SidebandLocalSaturationIntensity(2)) /(1 + 4 * (Delta_Sideband(2*i)/this.CoolingBeamLinewidth)^2); + 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); end end - - TotalSaturationParameter = sum(SaturationParameter); + TotalSaturationParameter = sum(SaturationParameter); + for i = 1:2 - 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)); + 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)); if this.SpontaneousEmission - 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); + a_scattering = this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i-1), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber) + ... + this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber); else a_scattering = [0,0,0]; end if this.Sideband - 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)); + 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)); if this.SpontaneousEmission 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+4), TotalSaturationParameter, Delta_Cooling(2*i), this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber); + this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i-1+4), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber) + ... + this.accelerationDueToSpontaneousEmissionProcess(SaturationParameter(2*i+4), TotalSaturationParameter, this.CoolingBeamLinewidth, this.CoolingBeamWaveNumber); else a_scattering = [0,0,0]; end @@ -86,7 +82,7 @@ function ret = calculateTotalAcceleration(this, PositionVector, VelocityVector) TotalAcceleration = TotalAcceleration + (a_2 - a_1) + a_scattering; end - if this.PushBeamRadius ~= 0 + if this.PushBeam TotalAcceleration = TotalAcceleration + this.accelerationDueToPushBeam(PositionVector, VelocityVector); end