add infinite volume residual gas damping

gwinc/ifo/aLIGO/__init__.py

@@ -30,9 +30,12 @@ class ExcessGas(nb.Budget):
     )
 
     noises = [
-        ExcessGasH2,
-        ExcessGasN2,
-        ExcessGasH2O,
+        ExcessGasScatteringH2,
+        ExcessGasScatteringN2,
+        ExcessGasScatteringH2O,
+        ExcessGasDampingH2,
+        ExcessGasDampingN2,
+        ExcessGasDampingH2O,
     ]
 
 

gwinc/ifo/aLIGO/ifo.yaml

@@ -37,11 +37,13 @@ Infrastructure:
       pressure: 4.0e-7              # Pa
       mass: 3.35e-27                # kg; Mass of H_2 (ref. 10)
       polarizability: 7.8e-31       # m^3
+      TestMassPressure: 4e-7
 
     N2:
       pressure: 1.33e-8
       mass: 4.65e-26
       polarizability: 1.79e-30
+      TestMassPressure: 4.65e-8
 
     H2O:
       pressure: 1.33e-8

gwinc/ifo/noises.py

@@ -764,7 +764,7 @@ class ExcessGas(nb.Noise):
 
     def calc(self):
         cavity = arm_cavity(self.ifo)
-        n = noise.residualgas.residual_gas_cavity(
+        n = noise.residualgas.residual_gas_scattering(
             self.freq, self.ifo, cavity, self.ifo.Infrastructure.ResidualGas)
         # FIXME HACK: it's unclear if a phase noise in the arms like
         # the excess gas noise should get the same dhdL strain
@@ -779,55 +779,109 @@ class ExcessGas(nb.Noise):
         return n * 2 / sinc_sqr
 
 
-class ExcessGasH2(nb.Noise):
-    """Excess gas for H2
+class ExcessGasScatteringH2(nb.Noise):
+    """Excess gas scattering for H2
 
     """
     style = dict(
-        label='H2',
+        label='H2 scattering',
         color='xkcd:red orange'
     )
 
     def calc(self):
         cavity = arm_cavity(self.ifo)
         species = self.ifo.Infrastructure.ResidualGas.H2
-        n = noise.residualgas.residual_gas_cavity(
+        n = noise.residualgas.residual_gas_scattering(
             self.freq, self.ifo, cavity, species)
         dhdl_sqr, sinc_sqr = dhdl(self.freq, self.ifo.Infrastructure.Length)
         return n * 2 / sinc_sqr
 
 
-class ExcessGasN2(nb.Noise):
-    """Excess gas for N2
+class ExcessGasScatteringN2(nb.Noise):
+    """Excess gas scattering for N2
 
     """
     style = dict(
-        label='N2',
+        label='N2 scattering',
         color='xkcd:emerald'
     )
 
     def calc(self):
         cavity = arm_cavity(self.ifo)
         species = self.ifo.Infrastructure.ResidualGas.N2
-        n = noise.residualgas.residual_gas_cavity(
+        n = noise.residualgas.residual_gas_scattering(
             self.freq, self.ifo, cavity, species)
         dhdl_sqr, sinc_sqr = dhdl(self.freq, self.ifo.Infrastructure.Length)
         return n * 2 / sinc_sqr
 
 
-class ExcessGasH2O(nb.Noise):
-    """Excess gas for H2O
+class ExcessGasScatteringH2O(nb.Noise):
+    """Excess gas scattering for H2O
 
     """
     style = dict(
-        label='H2O',
+        label='H2O scattering',
         color='xkcd:water blue'
     )
 
     def calc(self):
         cavity = arm_cavity(self.ifo)
         species = self.ifo.Infrastructure.ResidualGas.H2O
-        n = noise.residualgas.residual_gas_cavity(
+        n = noise.residualgas.residual_gas_scattering(
             self.freq, self.ifo, cavity, species)
         dhdl_sqr, sinc_sqr = dhdl(self.freq, self.ifo.Infrastructure.Length)
         return n * 2 / sinc_sqr
+
+
+class ExcessGasDampingH2(nb.Noise):
+    """Residual gas damping for H2
+
+    """
+    style = dict(
+        label='H2 damping',
+        color='xkcd:red orange',
+        linestyle='--',
+    )
+
+    @nb.precomp(sustf=precomp_suspension)
+    def calc(self, sustf):
+        species = self.ifo.Infrastructure.ResidualGas.H2
+        n = noise.residualgas.residual_gas_damping(
+            self.freq, self.ifo, species, sustf)
+        return 4 * n
+
+
+class ExcessGasDampingN2(nb.Noise):
+    """Excess gas damping for N2
+
+    """
+    style = dict(
+        label='N2 damping',
+        color='xkcd:emerald',
+        linestyle='--',
+    )
+
+    @nb.precomp(sustf=precomp_suspension)
+    def calc(self, sustf):
+        species = self.ifo.Infrastructure.ResidualGas.N2
+        n = noise.residualgas.residual_gas_damping(
+            self.freq, self.ifo, species, sustf)
+        return 4 * n
+
+
+class ExcessGasDampingH2O(nb.Noise):
+    """Excess gas damping for H2O
+
+    """
+    style = dict(
+        label='H2O damping',
+        color='xkcd:water blue',
+        linestyle='--',
+    )
+
+    @nb.precomp(sustf=precomp_suspension)
+    def calc(self, sustf):
+        species = self.ifo.Infrastructure.ResidualGas.H2O
+        n = noise.residualgas.residual_gas_damping(
+            self.freq, self.ifo, species, sustf)
+        return 4 * n

gwinc/noise/residualgas.py

@@ -7,18 +7,18 @@ from numpy import sqrt, log, pi
 from .. import const
 
 
-def residual_gas_cavity(f, ifo, cavity, species):
-    """Residual gas noise strain spectrum in an evacuated cavity.
+def residual_gas_scattering(f, ifo, cavity, species):
+    """Residual gas noise strain spectrum due to scattering.
 
     Noise caused by the passage of residual gas molecules through the
-    laser beams in a cavity.
+    laser beams in a cavity due to scattering.
 
     :f: frequency array in Hz
     :ifo: gwinc IFO structure
     :cavity: arm cavity structure
     :species: molecular species structure
 
-    :returns: arm displacement noise power spectrum at :f:
+    :returns: arm strain noise power spectrum at :f:
 
     The method used here is presented by Rainer Weiss, Micheal
     E. Zucker, and Stanley E. Whitcomb in their paper Optical
@@ -55,3 +55,38 @@ def residual_gas_cavity(f, ifo, cavity, species):
     zint[zint < 0] = 0
 
     return zint
+
+
+def residual_gas_damping(f, ifo, species, sustf):
+    """Noise due to residual gas damping for one test mass
+
+    :f: frequency array in Hz
+    :ifo: gwinc IFO structure
+    :species: molecular species structure
+    :sustf: suspension transfer function structure
+
+    :returns: displacement noise
+    """
+    kT = ifo.Infrastructure.Temp * const.kB
+    mass = species.mass
+    radius = ifo.Materials.MassRadius
+    thickness = ifo.Materials.MassThickness
+    thermal_vel = sqrt(kT/mass)  # thermal velocity
+
+    # pressure near the test mass
+    # possibly different from the pressure in the arms due to outgassing near
+    # the test mass
+    if 'TestMassPressure' in species:
+        pressure = species.TestMassPressure
+    else:
+        pressure = species.pressure
+
+    # infinite volume viscous damping coefficient for a cylinder
+    # table 1 of https://doi.org/10.1016/j.physleta.2010.06.041
+    beta_inf = pi * radius**2 * pressure/thermal_vel
+    beta_inf *= (1 + thickness/(2*radius) + pi/4)
+
+    # convert force to displacement noise using the suspension susceptibility
+    noise = 4 * kT * beta_inf * abs(sustf.tst_suscept)**2
+
+    return noise