Calculation of TripleLite2 model with violin mode stuff (for T1200418-v3)

Note: this notebook needs to be in a calculation directory of case {“mark.barton”,”20120120hstsMC2damp”} of the TripleLite2 model to be run.

Setup

Switches to enable loading of previously saved results instead of recalculating from scratch

useprecomputed = True; (* set to True to use saved results from precomputed subdirectory *)
If[useprecomputed,
	exceptdamping = True, (* False by default, True to recalculate just damping-dependent stuff*)
	exceptdamping = True (* DON"T CHANGE *)
];

loadcasefromuser["ASUS3L2ModelCaseDefn.m"];

overridesorig = overrides;

modelcase

modelcasecomment

Reset[All]

Calculate[constval];

Parameters of interest

Wire length

l3/.constval

Young’s modulus

Y3/.constval

Wire radius

r3/.constval

Mass (nominal value from model)

m3/.constval

Mass of bare IMCC-03 (L1, MC2) from nebula

MC2bare =  2.898;

Mass of primary prism

prismmass1 = 3980*(0.197*0.0254)*(0.223*0.0254)*(0.79*0.0254)/2

Mass of secondary prism

prismmass2 = rhosteel*(0.079*0.0254)*(0.091*0.0254)*(0.787*0.0254)/2/.constval

Mass of MC2 with prisms

MC2mass =MC2bare+2*prismmass1+2*prismmass2

Flexure length (nominal)

flex3/.constval

Twice flexure length as a fraction of total (cf. beta=2a/L)

2*flex3/l3/.constval

Half front-back wire separation

sl/.constval

Half left-right wire separation

n5/.constval

COM offset due to 0.5° horizontal wedge

COMy = 0.17/1000;

Wire tension  (nominal value from model)

tension = m3*g/nw3/.constval

Roll torque due to wedge

rolltorque = m3*COMy*g/.constval

rolltorque*n5/.constval

Extra tension to offset roll torque

deltatensionroll = rolltorque/n5/4/.constval

Extra tension as a percentage

100*deltatensionroll/tension

Pitch torque due to  0.1 mm prism offset (gluing tolerance)

pitchtorque = 0.0001*m3*g/.constval

Extra tension to offset pitch torque

deltatensionpitch = pitchtorque/sl/.constval

0.528*sl/pitchtorque/.constval

0.528*sl/.constval

1.201*sl/pitchtorque/.constval

100*deltatensionpitch/tension

Wire second moment of area

M31/.constval (* a.k.a. I *)

Fundamental violin mode for wires above mass 3 (the optic)

f3n=n*Sqrt[(m3*g/c3/nw3)/(rhosteel*A3)]/(2*l3*(1-2*flex3/l3-n^2*Pi^2*flex3^2/2/l3^2))

f31 = f3n/.n->1/.constval

f32 = f3n/.n->2/.constval

f33 = f3n/.n->3/.constval

Giles’ formula for dissipation dilution, D

D1 = (2/l3)*Sqrt[Y2*M31/tension]*(1+(1/2/l3)*(Pi*n)^2*Sqrt[Y2*M31/tension])/.constval/.n->1

D2 = (2/l3)*Sqrt[Y2*M31/tension]*(1+(1/2/l3)*(Pi*n)^2*Sqrt[Y2*M31/tension])/.constval/.n->2

D3 = (2/l3)*Sqrt[Y2*M31/tension]*(1+(1/2/l3)*(Pi*n)^2*Sqrt[Y2*M31/tension])/.constval/.n->3

Structural component of wire phi (before dissipation dilution)

phisteel/.constval

Expected Q at fundamental violin mode neglecting thermoelastic damping

Qstruct = 1/(D1*phisteel)/.constval

Total phi(f) including thermoelastic (note damping tag is “fibreatype” but model really does have wire parameters)

LogLinearPlot[damping[imag,fibreatype][f]/.constval,{f,1,10000},PlotRange->{0,Automatic},AxesLabel->{"f (Hz)","ϕ"}]

Total phi at fundamental violin mode (before dissipation dilution)

damping[imag,fibreatype][f31]/.constval

Total phi at fundamental violin mode (after dissipation dilution)

D1 * damping[imag,fibreatype][f31]/.constval

Q at fundamental violin mode (after dissipation dilution)

Q1 = 1/(D1*damping[imag,fibreatype][f31])/.constval

Q2 = 1/(D2*damping[imag,fibreatype][f32])/.constval

Q3 = 1/(D3*damping[imag,fibreatype][f33])/.constval

Keiko’s data

Raw data from https://alog.ligo-la.caltech.edu/aLOG/index.php?callRep=5097

keikodata20121029 = {
{{631.56,247105},{1263.3,163943},{1895.19,141496}},
{{644.977,251198},{1290.09,172007},{1935.45,150543}},
{{660.305,249889},{1320.62,187364},{1981.46,177172}},
{{671.023,229827},{1342.16,175060},{2013.56,169019}}
};

Refitted data from https://alog.ligo-la.caltech.edu/aLOG/index.php?callRep=5280

keikodata20121114 = {
{{631.56,234598},{1263.3,160933},{1895.19,195142}},
{{644.977,238817},{1290.09,165250},{1935.45,178511}},
{{660.305,239867},{1320.62,167053},{1981.46,174990}},
{{671.023,229072},{1342.16,164167},{2013.56,177825}}
};

Analysis of tensions

The frequency expressed as a function of tension

f1s = keikodata20121029[[All,1,1]];

foftension = ( Sqrt[(m3*g/nw3)/(rhosteel*A3)]/(2*(l3-2*flexoften)))/.flexoften->(Sqrt[((M31 nw3 Y3)/(g m3))])/.g*m3/nw3->ten/.nw3/(g*m3)->1/ten/.constval

The tensions required to account for the observed frequencies (N)

tensions = Table[ten/.NSolve[f1s[[i]]==foftension,ten][[1]],{i,4}]

Deviations from the mean (N)

tensions-Mean[tensions]

Percentage deviations from the mean

100*(tensions-Mean[tensions])/Mean[tensions]

Net tension

netten = Plus@@tensions

Net load force

load = MC2mass*g/.constval

Percent tension-load discrepancy

100*(netten/load-1)

Plot of tensions

ListPlot[tensions,PlotLabel->"Histogram of tensions",AxesLabel->{"wire number","tension (N)"},AxesOrigin->{0.5,6.5}]

Functions giving the three linear combinations representing tension imbalances

pitch[{{FL_,FR_},{BL_,BR_}}] := BL+BR-FL-FR

roll[{{FL_,FR_},{BL_,BR_}}] := FL+BL-FR-BR

pringle[{{FL_,FR_},{BL_,BR_}}] := FL+BR-FR-BL

Most likely mapping of tension values to the four corners of the optic, in {{FL,FR},{BL,BR}} order

corners  = {{tensions[[1]],tensions[[3]]},{tensions[[4]],tensions[[2]]}};

The imbalances for the likely mapping

{pitch[corners], roll[corners], pringle[corners]}

The torque corresponding to the pitch imbalance

pitch[corners]*sl/.constval

The roll torque corresponding to the roll imbalance

roll[corners]*n5/.constval

Damping comparison using mostly stock parameters but correcting for tension differences

Doftension[n_] = (2/l3)*Sqrt[Y3*M31/ten]*(1+(1/2/l3)*(Pi*n)^2*Sqrt[Y3*M31/ten])/.constval

Doftension[1]/.ten->tensions[[1]]

The dissipation dilution factors allowing for different tensions, for the four wires and three harmonics

Ds = Table[(2/l3)*Sqrt[Y3*M31/tensions[[i]]]*(1+(1/2/l3)*(Pi*n)^2*Sqrt[Y3*M31/tensions[[i]]])/.constval,{i,4},{n,3}]

The thermoelastic tau expressed as a function of tension (actually tension-independent)

tauoftension = taufibre/.constval

taus = Table[tauoftension/.ten->tensions[[i]],{i,4}]

The thermoelastic delta parameter expressed as a function of tension

deltaoftension=deltafibre/.overrides/.g*m3/nw3->ten/.constval

The delta values for the four wires, allowing for the tension differences

deltas = Table[deltaoftension/.ten->tensions[[i]],{i,4}]

The phi function from the standard damping model in terms of the three constants to be fitted

phifn[f_,phi_,tau_,delta_] = damping[imag,fibreatype][f]/.overrides/.taufibre->tau/.deltafibre->delta/.phisteel->phi

Dissipation diluted phi as a function of tension

fitfn[f_,n_,ten_,phi_,tau_,delta_]=Doftension[n]*phifn[f,phi,tau,delta]

Analysis of “improved” Q’s

keikodata = keikodata20121114;

Equations expressing the requirement that the fitted curves pass through the data points

Do[
wireeqns[w] = Table[fitfn[keikodata[[w,n,1]],n,tensions[[w]],phi,tau,delta]==1/keikodata[[w,n,2]],{n,3}],
{w,1,4}
]

Convert to expressions, LHS-RHS

Do[
wireexp[w] = wireeqns[w]/.x_==y_:>x-y,
{w,1,4}
]

Solve using FindRoot[] (which turns out to do a better job than NSolve[] and Minimize[]

wireroot[1] = FindRoot[wireexp[1],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];
wireroot[2] = FindRoot[wireexp[2],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];
wireroot[3] = FindRoot[wireexp[3],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];
wireroot[4] = FindRoot[wireexp[4],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];

Quick summary table

TableForm[Table[{phi,tau,delta}/.wireroot[w],{w,1,4}],TableHeadings->{keikodata[[All,1,1]],{phi,"tau (s)",delta}}]

A function to give the expected phi from the stock model, except with corrections for observed tension.

netphitencorronly[w_,n_]:= Ds[[w,n]]*damping[imag,fibreatype][keikodata[[w,n,1]]]/.overrides/.taufibre->taus[[w]]/.deltafibre->deltas[[w]]/.constval

Summary plot

Show[
LogLinearPlot[fitfn[f,f/f31,m3*g/4,phisteel,taufibre,deltafibre]/.constval,{f,20,2000},PlotStyle->{Black,Dotted}],
ListLogLinearPlot[{{f31,1/Q1},{f32,1/Q2},{f33,1/Q3}},PlotStyle->{Black,PointSize[Medium]}],

ListLogLinearPlot[keikodata[[1]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Red,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[1,n,1]],netphitencorronly[1,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Red}],
LogLinearPlot[fitfn[f,f/keikodata[[1,1,1]],tensions[[1]],phi,tau,delta]/.wireroot[1],{f,20,2000},PlotStyle->{Red,Dotted}],

ListLogLinearPlot[keikodata[[2]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Blue,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[2,n,1]],netphitencorronly[2,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Blue}],
LogLinearPlot[fitfn[f,f/keikodata[[2,1,1]],tensions[[2]],phi,tau,delta]/.wireroot[2],{f,20,2000},PlotStyle->{Blue,Dotted}],

ListLogLinearPlot[keikodata[[3]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Green,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[3,n,1]],netphitencorronly[3,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Green}],
LogLinearPlot[fitfn[f,f/keikodata[[3,1,1]],tensions[[3]],phi,tau,delta]/.wireroot[3],{f,20,2000},PlotStyle->{Green,Dotted}],

ListLogLinearPlot[keikodata[[4]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Brown,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[4,n,1]],netphitencorronly[4,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Brown}],
LogLinearPlot[fitfn[f,f/keikodata[[4,1,1]],tensions[[4]],phi,tau,delta]/.wireroot[4],{f,20,2000},PlotStyle->{Brown,Dotted}],

AxesOrigin->{Automatic,0},
PlotRange->All,
AxesLabel->{"f (Hz)","1/Q"}
]

Averages and standard deviations

fmean=Mean[{keikodata[[1,1,1]],keikodata[[2,1,1]],keikodata[[3,1,1]],keikodata[[4,1,1]]}]

fstdev=StandardDeviation[{keikodata[[1,1,1]],keikodata[[2,1,1]],keikodata[[3,1,1]],keikodata[[4,1,1]]}]

phimean=Mean[{phi/.wireroot[1],phi/.wireroot[2],phi/.wireroot[3],phi/.wireroot[4]}]

phistdev=StandardDeviation[{phi/.wireroot[1],phi/.wireroot[2],phi/.wireroot[3],phi/.wireroot[4]}]

deltamean=Mean[{delta/.wireroot[1],delta/.wireroot[2],delta/.wireroot[3],delta/.wireroot[4]}]

deltastdev=StandardDeviation[{delta/.wireroot[1],delta/.wireroot[2],delta/.wireroot[3],delta/.wireroot[4]}]

taumean=Mean[{tau/.wireroot[1],tau/.wireroot[2],tau/.wireroot[3],tau/.wireroot[4]}]

taustdev=StandardDeviation[{tau/.wireroot[1],tau/.wireroot[2],tau/.wireroot[3],tau/.wireroot[4]}]

ffmean=Mean[{1/(2*Pi*tau)/.wireroot[1],1/(2*Pi*tau)/.wireroot[2],1/(2*Pi*tau)/.wireroot[3],1/(2*Pi*tau)/.wireroot[4]}]

ffstdev=StandardDeviation[{1/(2*Pi*tau)/.wireroot[1],1/(2*Pi*tau)/.wireroot[2],1/(2*Pi*tau)/.wireroot[3],1/(2*Pi*tau)/.wireroot[4]}]

Summary table for publication

fitted = {
{f31,phisteel,deltafibre,taufibre,1/taufibre/2/Pi}/.constval,
{keikodata[[1,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[1],
{keikodata[[2,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[2],
{keikodata[[3,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[3],
{keikodata[[4,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[4],
{fmean,phimean,deltamean,taumean, ffmean},
{fstdev,phistdev,deltastdev,taustdev,ffstdev}
};

TableForm[fitted,TableHeadings->{{"stock","#1","#2","#3","#4","mean","stdev"},{"f1 (Hz)","phi","delta","tau (s)","1/(2*pi*tau) (Hz)"}}]

Average fitted delta as percentage

100*deltamean/deltafibre/.constval

Analysis of older Q’s

keikodata = keikodata20121029;

Equations expressing the requirement that the fitted curves pass through the data points

Do[
wireeqns[w] = Table[fitfn[keikodata[[w,n,1]],n,tensions[[w]],phi,tau,delta]==1/keikodata[[w,n,2]],{n,3}],
{w,1,4}
]

Convert to expressions, LHS-RHS

Do[
wireexp[w] = wireeqns[w]/.x_==y_:>x-y,
{w,1,4}
]

Solve using FindRoot[] (which turns out to do a better job than NSolve[] and Minimize[]

wireroot[1] = FindRoot[wireexp[1],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];
wireroot[2] = FindRoot[wireexp[2],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];
wireroot[3] = FindRoot[wireexp[3],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];
wireroot[4] = FindRoot[wireexp[4],{phi, phisteel/.constval},{delta,deltafibre/.constval},{tau,taufibre/.constval}];

Quick summary table

TableForm[Table[{phi,tau,delta}/.wireroot[w],{w,1,4}],TableHeadings->{keikodata[[All,1,1]],{phi,"tau (s)",delta}}]

A function to give the expected phi from the stock model, except with corrections for observed tension.

netphitencorronly[w_,n_]:= Ds[[w,n]]*damping[imag,fibreatype][keikodata[[w,n,1]]]/.overrides/.taufibre->taus[[w]]/.deltafibre->deltas[[w]]/.constval

Summary plot

Show[
LogLinearPlot[fitfn[f,f/f31,m3*g/4,phisteel,taufibre,deltafibre]/.constval,{f,20,2000},PlotStyle->{Black,Dotted}],
ListLogLinearPlot[{{f31,1/Q1},{f32,1/Q2},{f33,1/Q3}},PlotStyle->{Black,PointSize[Medium]}],

ListLogLinearPlot[keikodata[[1]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Red,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[1,n,1]],netphitencorronly[1,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Red}],
LogLinearPlot[fitfn[f,f/keikodata[[1,1,1]],tensions[[1]],phi,tau,delta]/.wireroot[1],{f,20,2000},PlotStyle->{Red,Dotted}],

ListLogLinearPlot[keikodata[[2]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Blue,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[2,n,1]],netphitencorronly[2,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Blue}],
LogLinearPlot[fitfn[f,f/keikodata[[2,1,1]],tensions[[2]],phi,tau,delta]/.wireroot[2],{f,20,2000},PlotStyle->{Blue,Dotted}],

ListLogLinearPlot[keikodata[[3]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Green,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[3,n,1]],netphitencorronly[3,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Green}],
LogLinearPlot[fitfn[f,f/keikodata[[3,1,1]],tensions[[3]],phi,tau,delta]/.wireroot[3],{f,20,2000},PlotStyle->{Green,Dotted}],

ListLogLinearPlot[keikodata[[4]]/.{f_,Q_}:>{f,1/Q},PlotStyle->{Brown,PointSize[Medium]}],
ListLogLinearPlot[Table[{keikodata[[4,n,1]],netphitencorronly[4,n]},{n,3}],AxesOrigin->{20,0},PlotStyle->{Brown}],
LogLinearPlot[fitfn[f,f/keikodata[[4,1,1]],tensions[[4]],phi,tau,delta]/.wireroot[4],{f,20,2000},PlotStyle->{Brown,Dotted}],

AxesOrigin->{Automatic,0},
PlotRange->All,
AxesLabel->{"f (Hz)","1/Q"}
]

Averages and standard deviations

fmean=Mean[{keikodata[[1,1,1]],keikodata[[2,1,1]],keikodata[[3,1,1]],keikodata[[4,1,1]]}]

fstdev=StandardDeviation[{keikodata[[1,1,1]],keikodata[[2,1,1]],keikodata[[3,1,1]],keikodata[[4,1,1]]}]

phimean=Mean[{phi/.wireroot[1],phi/.wireroot[2],phi/.wireroot[3],phi/.wireroot[4]}]

phistdev=StandardDeviation[{phi/.wireroot[1],phi/.wireroot[2],phi/.wireroot[3],phi/.wireroot[4]}]

deltamean=Mean[{delta/.wireroot[1],delta/.wireroot[2],delta/.wireroot[3],delta/.wireroot[4]}]

deltastdev=StandardDeviation[{delta/.wireroot[1],delta/.wireroot[2],delta/.wireroot[3],delta/.wireroot[4]}]

taumean=Mean[{tau/.wireroot[1],tau/.wireroot[2],tau/.wireroot[3],tau/.wireroot[4]}]

taustdev=StandardDeviation[{tau/.wireroot[1],tau/.wireroot[2],tau/.wireroot[3],tau/.wireroot[4]}]

ffmean=Mean[{1/(2*Pi*tau)/.wireroot[1],1/(2*Pi*tau)/.wireroot[2],1/(2*Pi*tau)/.wireroot[3],1/(2*Pi*tau)/.wireroot[4]}]

ffstdev=StandardDeviation[{1/(2*Pi*tau)/.wireroot[1],1/(2*Pi*tau)/.wireroot[2],1/(2*Pi*tau)/.wireroot[3],1/(2*Pi*tau)/.wireroot[4]}]

Summary table for publication

fitted = {
{f31,phisteel,deltafibre,taufibre,1/taufibre/2/Pi}/.constval,
{keikodata[[1,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[1],
{keikodata[[2,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[2],
{keikodata[[3,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[3],
{keikodata[[4,1,1]],phi,delta,tau,1/tau/2/Pi}/.wireroot[4],
{fmean,phimean,deltamean,taumean, ffmean},
{fstdev,phistdev,deltastdev,taustdev,ffstdev}
};

TableForm[fitted,TableHeadings->{{"stock","#1","#2","#3","#4","mean","stdev"},{"f1 (Hz)","phi","delta","tau (s)","1/(2*pi*tau) (Hz)"}}]

Average fitted delta as percentage

100*deltamean/deltafibre/.constval

Plots with damping

tablefs = Table[N[10^logf],{logf,-1,3,1/2}]

Calculate/plot thermal noise with stock damping

useprecomputed = True; (* set to True to use saved results from precomputed subdirectory *)
If[useprecomputed,
	exceptdamping = True, (* False by default, True to recalculate just damping-dependent stuff*)
	exceptdamping = True (* DON"T CHANGE *)
];

overrides = overridesorig;

Reset[All]; Calculate[Stage2]

stocktn= Table[Sqrt[noise2[eom2,opticxoutput,tablefs[[i]]]],{i,Length[tablefs]}]

DoWithStatus["Plotting stage 2 x thermal noise",
stocktnplot = plotTN[eom2,opticxoutput,0.1,100,1,PlotStyle->{Black}]]

Calculate/plot thermal noise with fitted damping parameters

overrides = Override[
overridesorig,
{
phisteellf ->  phi/.wireroot[1],
deltafibrelf ->  delta/.wireroot[1],
taufibrelf ->  tau/.wireroot[1],
phisteellb ->  phi/.wireroot[3],
deltafibrelb ->  delta/.wireroot[3],
taufibrelb ->  tau/.wireroot[3],
phisteelrf ->  phi/.wireroot[4],
deltafibrerf ->  delta/.wireroot[4],
taufibrerf ->  tau/.wireroot[4],
phisteelrb ->  phi/.wireroot[2],
deltafibrerb ->  delta/.wireroot[2],
taufibrerb ->  tau/.wireroot[2]
}];

useprecomputed = True; (* set to True to use saved results from precomputed subdirectory *)
If[useprecomputed,
	exceptdamping = True, (* False by default, True to recalculate just damping-dependent stuff*)
	exceptdamping = True (* DON"T CHANGE *)
];

Reset[All]; Calculate[Stage2]

fittedtn = Table[Sqrt[noise2[eom2,opticxoutput,tablefs[[i]]]],{i,Length[tablefs]}]

fittedtnplot = plotTN[eom2,opticxoutput,0.1,100,1,PlotStyle->{Red}]

Comparison

Show[stocktnplot,fittedtnplot,PlotLabel->"MC2 thermal noise with stock damping (black) and fitted damping (red)",AxesLabel->{"f (Hz)","x (m/rtHz)"}]

TableForm[Transpose[{tablefs,stocktn,fittedtn,fittedtn/stocktn}],TableHeadings->{None,{"f (Hz)","x (m/rtHz), stock","x (m/rtHz), fitted","fitted/stock"}}]