In [1]:
import numpy as np
import scipy
from scipy.special import genlaguerre
from scipy.special import laguerre
In [2]:
def TSum(r):
    output = 0
    for k in range(1,100):
        output = output + ((2*r**2/r0**2)**k)*((-1)**k)/(k*np.math.factorial(k))
    return output

TSumNorm = 0

T0 = 300 #This shouldn't matter
r0 = 0.002
rmax = 5*r0
points = 10000
deltar = rmax/points

for i in range(points):
    TSumNorm = TSumNorm + 4*np.exp(-2*(rmax*i/points)**2/r0**2)*TSum(rmax*i/points)*rmax*i/points*deltar/r0**2
    
def L0(p,r):
    return np.sqrt(2/(np.pi*r0**2))*laguerre(p)(2*r**2/r0**2)*np.exp(-r**2/r0**2)    


def Aint(m):
    output = 0
    for i in range(points):
        output = output+L0(m,rmax*i/points)*L0(0,rmax*i/points)*(TSum(rmax*i/points)-0*TSumNorm)*rmax*i/points*deltar
    return output

"""
IntSum = 0
for i in range(1,10):
    IntSum = IntSum + abs(Aint(i)) 
    print(f"Aint({i}) = {abs(Aint(i))}")
print(f"IntSum = {IntSum}")
"""

IntegralNumber = 15
IntArray = np.zeros(IntegralNumber)

IntSquareSum = 0
for i in range(1,IntegralNumber):
    Int = Aint(i)
    IntArray[i-1] = Int
    IntSquareSum = IntSquareSum + (abs(Int))**2
    print(f"Aint({i}) = {abs(Aint(i))}")
print(f"IntSumSquare = {IntSquareSum}")
Aint(1) = 0.07957747155400367
Aint(2) = 0.01989436770093922
Aint(3) = 0.006631458683989795
Aint(4) = 0.0024867674697319074
Aint(5) = 0.0009949450829378457
Aint(6) = 0.00041304850518690816
Aint(7) = 0.0001847574619399611
Aint(8) = 4.848929272410496e-05
Aint(9) = 0.00013276907831016163
Aint(10) = 0.0002550577768827879
Aint(11) = 0.0006126515969918277
Aint(12) = 0.0010745552684240912
Aint(13) = 0.0014572310661290224
Aint(14) = 0.0013212553755210699
IntSumSquare = 0.006785199044440701
In [54]:
#GQuEST Silicon
betaSi = 1E-4 #Jan Meyer et al 2022 Class. Quantum Grav. 39 135001
aSi = 0.002
Power = 10E3
bulkabsoprtionSi = 200E-6
lambdal = 1550E-9
kappaSi = 380 #Lee says 739
CoatingAbsorptionGQuEST = 2*1E-6

print(f"The GQuEST bulk absorption is {aSi*bulkabsoprtionSi}")
print(f"The GQuEST total absorption is {aSi*bulkabsoprtionSi+CoatingAbsorptionGQuEST}")

#00 to HOM
FractionLostSi = ((betaSi*Power*(aSi*bulkabsoprtionSi+CoatingAbsorptionGQuEST)/(lambdal*kappaSi))**2)*(4*np.pi**2*IntSquareSum)

PowerLostSi = Power*FractionLostSi/4

print(f"The total GQuEST Silicon power loss with 00 to HOM math is {PowerLostSi}W, or {1E6*FractionLostSi}ppm")

#Lee paramaters and equations

kappaSi = 739

FractionLostSiGiven = (betaSi*(CoatingAbsorptionGQuEST)*Power/(lambdal*kappaSi))**2 + (betaSi*(1.3*aSi*bulkabsoprtionSi)*Power/(lambdal*kappaSi))**2

PowerLostSiGiven = Power*FractionLostSiGiven/4

print(f"According to Lee's numbers and math, the total GQuEST Silicon power loss is {PowerLostSiGiven}W, or {1E6*FractionLostSiGiven}ppm")
The GQuEST bulk absorption is 4.0000000000000003e-07
The GQuEST total absorption is 2.4e-06
The total GQuEST Silicon power loss with 00 to HOM math is 0.01111727438828679W, or 4.446909755314715ppm
According to Lee's numbers and math, the total GQuEST Silicon power loss is 0.008136849662447926W, or 3.2547398649791703ppm
In [21]:
##00 Overlap
output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(betaSi*Power*(aSi*bulkabsoprtionSi+CoatingAbsorptionGQuEST)/(lambdal*kappaSi))*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

#print(output)

FractionLostSi = 1-np.abs(output)**2

PowerLostSi = Power*FractionLostSi/4

print(f"The total GQuEST Silicon power loss with 00 Overlap math is {PowerLostSi}W, or {1E6*FractionLostSi}ppm")
The total GQuEST Silicon power loss with 00 Overlap math is 0.003354190447890648W, or 1.3416761791562593ppm
In [58]:
#GQuEST Fused Slica (2mm thickness for now)
betafs = 8.5E-6
afs = 0.002
Power = 10E3
bulkabsoprtionfs = 100E-6
lambdal = 1550E-9
kappafs = 1.38
CoatingAbsorptionGQuEST = 2*1E-6

print(f"The GQuEST fused silica bulk absorption is {afs*bulkabsoprtionfs}")
print(f"The GQuEST fused silica total absorption is {afs*bulkabsoprtionfs+CoatingAbsorptionGQuEST}")

#00 to HOM
FractionLostfs = ((betafs*Power*(afs*bulkabsoprtionfs+CoatingAbsorptionGQuEST)/(lambdal*kappafs))**2)*(4*np.pi**2*IntSquareSum)

PowerLostfs = Power*FractionLostfs/4

print(f"The total GQuEST fused silica power loss with 00 to HOM math is {PowerLostfs}W, or {1E6*FractionLostfs}ppm")

#Lee paramaters and equations

FractionLostfsGiven = (betafs*(CoatingAbsorptionGQuEST+1.3*afs*bulkabsoprtionfs)*Power/(lambdal*kappafs))**2
FractionLostfsGiven = (betafs*(CoatingAbsorptionGQuEST)*Power/(lambdal*kappafs))**2 + (betafs*(1.3*afs*bulkabsoprtionfs)*Power/(lambdal*kappafs))**2

PowerLostfsGiven = Power*FractionLostfsGiven/4

print(f"According to Lee's math, the total GQuEST fused silica power loss is {PowerLostfsGiven}W, or {1E6*FractionLostfsGiven}ppm")
The GQuEST fused silica bulk absorption is 2.0000000000000002e-07
The GQuEST fused silica total absorption is 2.2e-06
The total GQuEST fused silica power loss with 00 to HOM math is 5.117620469494779W, or 2047.0481877979114ppm
According to Lee's math, the total GQuEST fused silica power loss is 16.05811373671924W, or 6423.245494687695ppm
In [4]:
betafs = 8.5E-6
aHolometer = 3*25.4/1000
PowerHolometer = 2E3
bulkabsoprtionHolometer = 30E-6
lambda1064 = 1064E-9
kappafs = 1.38
CoatingAbsorptionHolometer = 2*1E-6

print(f"The Holometer bulk absorption is {aHolometer*bulkabsoprtionHolometer}")
print(f"The Holometer total absorption is {aHolometer*bulkabsoprtionHolometer+CoatingAbsorptionHolometer}")

#00 to HOM
FractionLostHolometer = ((betafs*PowerHolometer*(aHolometer*bulkabsoprtionHolometer+CoatingAbsorptionHolometer)/(lambda1064*kappafs))**2)*(4*np.pi**2*IntSquareSum)

PowerLostHolometer = PowerHolometer*FractionLostHolometer/4

print(f"The total Holometer power loss with 00 to HOM math is {PowerLostHolometer}W, or {1E6*FractionLostHolometer}ppm")

#00 to HOM No Coating Loss
FractionLostHolometer = ((betafs*PowerHolometer*(aHolometer*bulkabsoprtionHolometer)/(lambda1064*kappafs))**2)*(4*np.pi**2*IntSquareSum)

PowerLostHolometer = PowerHolometer*FractionLostHolometer/4

print(f"The total Holometer power loss without coating with 00 to HOM math is {PowerLostHolometer}W, or {1E6*FractionLostHolometer}ppm")

#Lee paramaters and equations

FractionLostHolometerGiven = (betafs*(CoatingAbsorptionHolometer)*PowerHolometer/(lambda1064*kappafs))**2 + (betafs*(1.3*aHolometer*bulkabsoprtionHolometer)*PowerHolometer/(lambda1064*kappafs))**2

PowerLostHolometerGiven = PowerHolometer*FractionLostHolometerGiven/4

print(f"According to Lee's math, the total Holometer power loss is {PowerLostHolometerGiven}W, or {1E6*FractionLostHolometerGiven}ppm")

r0 = 0.005
rmax = 5*r0
points = 20000
deltar = rmax/points

output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(betafs*PowerHolometer*aHolometer*bulkabsoprtionHolometer/(2*lambda1064*kappafs))*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

#print(output)

FractionLostHolometerHolger = 1-np.abs(output)**2

PowerLostHolometerHolger = PowerHolometer*FractionLostHolometerHolger/4

print(f"According to Holger's 00 to 00 math assuming a single BS Trip, the total Holometer power loss is {PowerLostHolometerHolger}W, or {1E6*FractionLostHolometerHolger}ppm")

output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(betafs*PowerHolometer*aHolometer*bulkabsoprtionHolometer/(lambda1064*kappafs))*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

#print(output)

FractionLostHolometerHolger2trip = 1-np.abs(output)**2

PowerLostHolometerHolger2trip = PowerHolometer*FractionLostHolometerHolger2trip/4

print(f"According to Holger's 00 to 00 math assuming two BS Trips, the total Holometer power loss is {PowerLostHolometerHolger2trip}W, or {1E6*FractionLostHolometerHolger2trip}ppm")





PowerLostholometerMeasured = 0.10*0.8

LambdaHolometerMeasured = 4*PowerLostholometerMeasured/PowerHolometer

print(f"Experimentally, the total Holometer power loss is {PowerLostholometerMeasured}W, or {1E6*LambdaHolometerMeasured}ppm")
The Holometer bulk absorption is 2.286e-06
The Holometer total absorption is 4.285999999999999e-06
The total Holometer power loss with 00 to HOM math is 0.3298018449059192W, or 659.6036898118383ppm
The total Holometer power loss without coating with 00 to HOM math is 0.09382124667590064W, or 187.64249335180128ppm
According to Lee's math, the total Holometer power loss is 0.8600170675031884W, or 1720.0341350063768ppm
According to Holger's 00 to 00 math assuming a single BS Trip, the total Holometer power loss is 0.023456649015096698W, or 46.913298030193396ppm
According to Holger's 00 to 00 math assuming two BS Trips, the total Holometer power loss is 0.09375741982714514W, or 187.51483965429028ppm
Experimentally, the total Holometer power loss is 0.08000000000000002W, or 160.00000000000003ppm
In [77]:
 
Out[77]:
4.042904621685007
In [201]:
#GEO600
def TSum(r):
    output = 0
    for k in range(1,100):
        output = output + ((2*r**2/r0**2)**k)*((-1)**k)/(k*np.math.factorial(k))
    return output

betafs = 8.5E-6
aGEO = 9/100
PowerGEO = 2.7E3
bulkabsoprtionGEO = 50E-6
lambda1064 = 1064E-9
kappafs = 1.38
CoatingAbsorptionGEO = 0

#OPD, nm:
(betafs*PowerGEO*(aGEO*bulkabsoprtionGEO+CoatingAbsorptionGEO)/(2*np.pi*kappafs))*(-TSum(1*r0))*1E9

#00 to HOM
FractionLostGEO = ((betafs*PowerGEO*(aGEO*bulkabsoprtionGEO+CoatingAbsorptionGEO)/(lambda1064*kappafs))**2)*(4*np.pi**2*IntSquareSum)

PowerLostGEO = PowerGEO*FractionLostGEO/4

print(f"The total GEO600 power loss with 00 to HOM math is {PowerLostGEO}W, or {1E6*FractionLostGEO}ppm")

#Lee paramaters and equations

FractionLostGEOLee = (betafs*(CoatingAbsorptionGEO)*PowerGEO/(lambda1064*kappafs))**2 + (betafs*(1.3*aGEO*bulkabsoprtionGEO)*Power/(lambda1064*kappafs))**2

PowerLostGEOLee = Power*FractionLostGEOGiven/4

print(f"According to Lee's math, the total GEO power loss is {PowerLostGEOGiven}W, or {1E6*FractionLostGEOGiven}ppm")

r0 = 0.01
rmax = 5*r0
points = 20000
deltar = rmax/points

output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(betafs*PowerGEO*aGEO*bulkabsoprtionGEO/(2*lambda1064*kappafs))*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

#print(output)

FractionLostGEOHolger = 1-np.abs(output)**2

PowerLostGEOHolger = Power*FractionLostGEOHolger/4

print(f"According to Holger's 00 to 00 math assuming a single BS Trip, the total GEO power loss is {PowerLostGEOHolger}W, or {1E6*FractionLostGEOHolger}ppm")

output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(betafs*PowerGEO*aGEO*bulkabsoprtionGEO/(lambda1064*kappafs))*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

#print(output)

FractionLostGEOHolger2trip = 1-np.abs(output)**2

PowerLostGEOHolger2trip = Power*FractionLostGEOHolger2trip/4

print(f"According to Holger's 00 to 00 math assuming two BS Trips, the total GEO power loss is {PowerLostGEOHolger2trip}W, or {1E6*FractionLostGEOHolger2trip}ppm")

PowerLostGEOGiven = 4E-4

LambdaGEOGiven = 4*PowerLostGEOGiven/PowerGEO

print(f"According to Holger's 00 to 00 math assuming two BS Trips, the total GEO power loss is {PowerLostGEOGiven}W, or {1E6*LambdaGEOGiven}ppm")
The total GEO600 power loss with 00 to HOM math is 0.8936086754931182W, or 1323.864704434249ppm
According to Lee's math, the total GEO power loss is 0.0004W, or 8360.566068450002ppm
According to Holger's 00 to 00 math assuming a single BS Trip, the total GEO power loss is 0.22343238615863958W, or 331.0109424572438ppm
According to Holger's 00 to 00 math assuming two BS Trips, the total GEO power loss is 0.893195858790341W, or 1323.2531241338386ppm
According to Holger's 00 to 00 math assuming two BS Trips, the total GEO power loss is 0.0004W, or 0.5925925925925926ppm
In [179]:
#GEO600 00 Overlap integrals
#Try to recreate their numbers
#Real version needs coating loss and factor of 2 for two trips through BS

betafs = 8.5E-6
aGEO = 9/100
PowerGEO = 2.7E3
bulkabsoprtionGEO = 50E-6
lambda1064 = 1064E-9
kappafs = 1.38
CoatingAbsorptionHolometer = 2*1E-6

r0 = 0.01
rmax = 5*r0
points = 20000
deltar = rmax/points

output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(betafs*PowerGEO*aGEO*bulkabsoprtionGEO/(2*lambda1064*kappafs))*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

print(output)
print(1-np.abs(output)**2)
(0.9998344808301465+8.089602313618409e-07j)
0.00033101094245724383
In [187]:
OPD = betafs*Power*aGEO*bulkabsoprtionGEO/(4*np.pi*kappafs)
print(OPD)

output = 0
for i in range(points):
    output = output+2*np.pi*L0(0,rmax*i/points)*L0(0,rmax*i/points)*np.exp(0+1j*(2*np.pi*OPD/lambda1064)*(TSum(rmax*i/points)-TSumNorm))*rmax*i/points*deltar

print(output)

FractionalPowerGEO = 1-np.abs(output)**2

print(FractionalPowerGEO*1E6)
5.955335778194019e-09
(0.9998344808301465+8.089602313598662e-07j)
331.0109424572438
In [184]:
print((1-np.abs(output)**2)/(4E-4/Power))
2234.323861586396
In [182]:
4E-4/Power
Out[182]:
1.4814814814814815e-07
In [204]:
(4*np.pi**2*IntSquareSum)
Out[204]:
0.2676052473213298
In [205]:
IntSquareSum
Out[205]:
0.006778520101874418
In [210]:
Power = 1
OPD = betaSi*Power*aSi*bulkabsoprtionSi/(4*np.pi*kappaSi)
print(OPD)
4.307305631715706e-15
In [220]:
Power = 2000
OPD = betafs*Power*aHolometer*bulkabsoprtionfs/(4*np.pi*kappafs)*(-TSum(r0))
print(OPD)
9.85476884704998e-09
In [ ]:
 
In [ ]: