From 4909a0ee7b8c92aeb0e00cafc5a74582b4d36b56 Mon Sep 17 00:00:00 2001
From: Christopher Wipf <wipf@ligo.mit.edu>
Date: Tue, 14 Apr 2020 15:44:50 -0700
Subject: [PATCH] Revise suspension model
Upper joint temperature comes from the stage above
Drop support for BQuad suspension type
Use the fiber end cross-section consistently in wireTELoss()
Remove the "questionable approximation" in continuumWireKh
Correct the spring constant calculation in continuumWireKv
---
gwinc/suspension.py | 22 ++--------------------
1 file changed, 2 insertions(+), 20 deletions(-)
diff --git a/gwinc/suspension.py b/gwinc/suspension.py
index 16b1d51..6434d00 100644
--- a/gwinc/suspension.py
+++ b/gwinc/suspension.py
@@ -100,7 +100,6 @@ def getJointParams(sus, n):
TempUpper = stages[n-1].Temp
else:
TempUpper = sus.Temp
- TempUpper = Temp # FIXME: reproduces the old calculation
##############################
# material parameters
@@ -114,10 +113,6 @@ def getJointParams(sus, n):
elif 'WireMaterial' in stage:
wireMatUpper = stage.WireMaterial
wireMatLower = stage.WireMaterial
- # FIXME: reproduces the old calculation
- elif n == last_stage and sus.Type == 'BQuad':
- wireMatUpper = 'Silicon'
- wireMatLower = 'Silicon'
elif n == last_stage:
wireMatUpper = 'Silica'
wireMatLower = 'Silica'
@@ -132,9 +127,6 @@ def getJointParams(sus, n):
# support blade (upper joint only)
if 'BladeMaterial' in stage:
bladeMat = stage.BladeMaterial
- # FIXME: reproduces the old calculation
- elif n == last_stage and sus.Type == 'BQuad':
- bladeMat = 'Silicon'
else:
bladeMat = 'MaragingSteel'
@@ -206,7 +198,7 @@ def wireGeometry(r, N, RibbonThickness=None, TaperedEndRadius=None, **kwargs):
return (xsect, xsectEnd, xII, mu_v, mu_h)
-def wireTELoss(w, tension, xsectEnd, xII, Temp, alpha, beta, rho, C, K, Y, xsect,
+def wireTELoss(w, tension, xsectEnd, xII, Temp, alpha, beta, rho, C, K, Y,
RibbonThickness=None, **kwargs):
"""Thermoelastic calculations for wires
@@ -222,9 +214,6 @@ def wireTELoss(w, tension, xsectEnd, xII, Temp, alpha, beta, rho, C, K, Y, xsect
# horizontal TE time constant, wires
# The constant 7.37e-2 is 1/(4*q0^2) from eq 12, C. Zener 10.1103/PhysRev.53.90 (1938)
tau = 7.37e-2 * 4 * (rho * C * xsectEnd) / (pi * K)
- # FIXME: reproduces the old calculation
- # then xsect can be dropped from the argument list
- tau = 7.37e-2 * 4 * (rho * C * xsect) / (pi * K)
# deal with ribbon geometry
if RibbonThickness is not None:
@@ -306,16 +295,12 @@ def continuumWireKh(w, N, length, tension, xsect, xII, rho, Y, phi):
# = T k (cos(k L) + k delta sin(k L))
# for w -> 0, this reduces to N_w * T * k
khnum = N * tension * k * (coskl + dk * sinkl)
- # FIXME: reproduces the old calculation
- khnum = N * tension * k * (1 + dk**2) * (coskl + dk * sinkl)
# denominator, horiz spring constant
# D after equation 8 in GG
# D = sin(k L) - 2 k delta cos(k L)
# for w -> 0, this reduces to k (L - 2 delta)
khden = sinkl - 2 * dk * coskl
- # FIXME: reproduces the old calculation
- khden = (1 - dk**2) * sinkl - 2 * dk * coskl
return khnum/khden
@@ -336,9 +321,6 @@ def continuumWireKv(w, N, length, xsect, xsectEnd, rho, Y, phi,
kv_mid = N * xsect * Y * k / tan(k * l_mid)
kv_end = N * xsectEnd * Y * k / tan(k * l_end)
kv = 1/(2/kv_end + 1/kv_mid)
- # FIXME: reproduces the old calculation
- kv_end = N * xsectEnd * Y * k / tan(k * 2*l_end)
- kv = 1/(1/kv_end + 1/kv_mid)
return kv
@@ -491,7 +473,7 @@ def suspQuad(f, sus):
# nominally this term becomes zero for steel wires as ds_w is zero
# The last term is for TE loss
phih_TE = wireTELoss(w, tension, xsectEnd, xII, Temp,
- alpha_w, beta_w, rho_w, C_w, K_w, Y_w, xsect,
+ alpha_w, beta_w, rho_w, C_w, K_w, Y_w,
**wireShape)
phih = phi_w * (1 + mu_h * ds_w) + phih_TE
--
GitLab