40m QIL Cryo_Lab CTN SUS_Lab TCS_Lab OMC_Lab CRIME_Lab FEA ENG_Labs OptContFac Mariner WBEEShop
 COMSOL elog, Page 2 of 3 Not logged in ID Date Author Type Category Subject
39   Thu Jun 27 17:12:33 2013 Arnaldo RodriguezOpticsGeneralDefocusing in P-Controlled Manual Simulation

To test the effectiveness of the PID loop, I've only given a nonzero value to the K_p constant and varied it over a series of runs to verify the response by the controller.

The following is a plot of the baseline compared to the P-active simulations for different proportionality gain constants. I've yet to determine whether all the trends in behavior are physical.

Each different case has an evident droop (steady-state error), and seem to be destabilizing as K_p increases. However, most cases reach steady-state quicker than the baseline. It should be noted that the ring heater is not correcting any effects from the beam, as it is turned off. The ring heater is correcting itself; or rather, its own initial heat pulse of 5 W * delta_t (600s) = 3000 J before the controller becomes active.

40   Fri Jun 28 11:54:10 2013 Emory BrownOpticsGeneralMaintaining a constant volume while modifying the radii of the test mass faces

When changing the dimensions of the mirror, we want to keep a constant volume. To do so, we must modify the height of the test mass. The area of the test mass if R is the radius of the front mirror face and r the radius of the back mirror face is (1/3)*pi*h*(R^2+R*r+r^2)=area, so using the values currently in our model, R=r=0.17 m and h=0.2 m, area=0.0181584055 m^3=pi*0.2*0.17^2 m^3, so to ensure that the volume of the test mass remains constant, h=0.01734/(R^2+R*r+r^2). Note that this does give 0.2 m when R=r=0.17 m.

41   Fri Jun 28 14:03:07 2013 Deep ChatterjeeGeneralGeneralApproach for TE noise for beam in cavity

In this post I mention of the approach that would be taken to solve for the TE noise spectrum for a beam inside the substrate cavity.

Completing the same would be my plan for the next week.

In the COMSOL simulation, a gaussian pressure will be applied on both the flat surface of the cylinder. The strain would be calculated
using the Structural mechanics module and the temperature field would be calculated using the Heat Transfer module. The approach
is similar to that followed in Liu and Thorne. The Work dissipated would be calculated using the gradient of  T. The spectral density is
then calculated using the FDT.

In the attachment, I have written briefly about the material I have studied this week along with my thoughts on the problem of the TE noise
calculation.

Attachment 1: LL_derivations.pdf    42   Mon Jul 1 13:01:58 2013 Arnaldo RodriguezOpticsGeneralBrief Parameter Study on P-active Control Loop

The following is a plot of the defocus as a function of time with only the proportional gain being nonzero. Only three cases were run before the program crashed. The dotted line indicates the setpoint.

One can clearly see the loop becomes unstable and oscillates out of control at a proportional gain of 240, and can see the characteristic steady-state offset of the P-controller from the set point develop in the other cases.

43   Mon Jul 1 17:31:40 2013 Arnaldo RodriguezOpticsGeneralParameter Study on Sampling Time for P-active Loop

The following is a plot of the defocus as a function of time for a fixed proportionality gain of 150, and different listed sampling times. The discrepancy between the behaviors is enough to warrant a proper mathematical explanation; the max. overshoot, for example, is inversely proportional to the sampling rate.

44   Tue Jul 2 13:32:30 2013 Emory BrownOpticsGeneralUsing Matlab to run our COMSOL simulations

We want to be able to test a number of values for the ratio of the radii of the frustum which will be much easier once we can get Matlab to automate running our COMSOL model.  I have been able to output the COMSOL model to a .m file and was able to figure out a way to export the Umax value we want to Matlab (see the attached .m file).

The model ran in Matlab after debugging it for a short while and returned the same result as the COMSOL model had by itself.  I setup a script to loop over some values of RatioRadii and feed them to the model.  The first value it fed to the model worked and returned 1.2099e-10 J for RatioRadii=0.1, but for RatioRadii=0.2 the model had a

meshing error.  I am going to have to see if there is some way we can modify the mesh to fix this error.  When I briefly looked into it today, I found that for several RatioRadii which I tried, the object could be meshed for certain finenesses, but not others.  I am going to spend tomorrow working on the first progress report, but after that is done I will

work on fixing this meshing error.

Attachment 1: model_7_1_13.m
```function [Umax] = model(ratioRadii)
% ratioRadii must be a scalar
% model_7_1_13.m
%
% Model exported on Jul 1 2013, 17:08 by COMSOL 4.3.2.189.

import com.comsol.model.*
import com.comsol.model.util.*

model = ModelUtil.create('Model');
```
... 220 more lines ...
Attachment 2: looper.m
```function [RatioRValues,UmaxValues] = looper
% Iterates over a range of RatioRValues

RatioRValues=[0.1:0.1:10]
UmaxValues=[]
for r=RatioRValues,
UmaxValues = [UmaxValues model_7_1_13(r)];
UmaxValues
end
```
... 4 more lines ...
45   Wed Jul 3 00:04:46 2013 KojiOpticsGeneralUsing Matlab to run our COMSOL simulations

Please use zipped files instead of uploading extremely long text files on the elog.
Or use wiki or svn instead of elog for uploading them.

And appropriate word wrapping so that we can read the entries. Please.

46   Thu Jul 4 15:23:47 2013 Deep ChatterjeeOpticsGeneralAttempt to compute TE noise following Liu and Thorne

Liu and Thorne, in their paper have analytically solved for the TE noise spectra considering the test mass to be a half infinite space.

An attempt was made to get similar results using a COMSOL model. In this post, the rough profile of the PSD plot is given and is compared
to Liu and Thorne's result.

Steps followed by Liu and Thorne:

> Oscillating pressure, scaled according to beam profile, is applied one of the surfaces.

> The stress balance equation is solved borrowing algorithm from Landau Lifshitz.

> The stress and the temperature perturbation is related using an equation given in Chapter 1 of LL

> The W_diss is calculated from the gradient of the temperature.

> PSD is calculated using FDT due to Levin

> The resulting plot looks like the following - ***Only the frequency dependence is plotted here leaving  out the prefactors which will
appear just as an additive constant in the figure. The aim was to see how the plot of the simulation looks like. Steps followed in the COMSOL model:

> 2D axisymmetric model was used since the applied was applied radially symmetrically.

> The dimensions of the test mass(cylinder) was defined larger as compared to the beam spot radius.

> The physics added was Solid Stress and Heat Transfer in Solids.

> Thermal expansion was enabled in the Solid Stress physics

> Pressure work was enabled in Heat Transfer physics which takes into consideration heat due to stresses i.e. adds this heat source to the heat eqn

> Heat insulation wsa applied to all surfaces.

> All surfaces apart from the one where pressure was applied was held fixed to mimic infinity

> A custom mesh was put.

> The model was exported as '.m'

> It was called for frequency range 10^(0 : 0.25 : 4)

> Following Liu and Thorne, the temperature gradient was integrated over the volume using mphint2() function in Matlab

> The time average was taken over the second half of the cycle to avoid transient which are assumed to die off after a couple of initial cycles.

> All the prefactors were ignored while plotting

> The plot shown below is integral{grad_T}^2} / omega^2   which is proportional to PSD

> We just wanted to see if at all we got anything looking similar to the analytic result of Liu and Thorne (shown above) The above plot is when the radius is 100 times the beam spot and height of the cylinder is 200 times the beam spot
Other values of radius and height were also tried. The one below shows similar plot as above but for radius 10 times
the beam spot and height 20 times the same While the one below is for radius = 1000 * beam spot and height = 2000 * beam spot Here, it is seen that the straight line with a negative slope is seen to some extent.

However, the prefactors, when put in, will show the actual PSD.

Any comments or suggestions related to the model construction in COMSOL or otherwise is welcome.

Attachment 1: Liu&Thorne_type.bmp
Attachment 5: 2D_almost_infinte_TE.mph
47   Fri Jul 5 14:06:24 2013 Deep ChatterjeeOpticsGeneralAttempt to compute TE noise following Liu and Thorne

 Quote: Liu and Thorne, in their paper have analytically solved for the TE noise spectra considering the test mass to be a half infinite space. An attempt was made to get similar results using a COMSOL model. In this post, the rough profile of the PSD plot is given and is compared to Liu and Thorne's result. Steps followed by Liu and Thorne: > Oscillating pressure, scaled according to beam profile, is applied one of the surfaces. > The stress balance equation is solved borrowing algorithm from Landau Lifshitz. > The stress and the temperature perturbation is related using an equation given in Chapter 1 of LL > The W_diss is calculated from the gradient of the temperature. > PSD is calculated using FDT due to Levin > The resulting plot looks like the following - ***Only the frequency dependence is plotted here leaving  out the prefactors which will appear just as an additive constant in the figure. The aim was to see how the plot of the simulation looks like. Steps followed in the COMSOL model: > 2D axisymmetric model was used since the applied was applied radially symmetrically. > The dimensions of the test mass(cylinder) was defined larger as compared to the beam spot radius. > The physics added was Solid Stress and Heat Transfer in Solids. > Thermal expansion was enabled in the Solid Stress physics > Pressure work was enabled in Heat Transfer physics which takes into consideration heat due to stresses i.e. adds this heat source to the heat eqn > Heat insulation wsa applied to all surfaces. > All surfaces apart from the one where pressure was applied was held fixed to mimic infinity > A custom mesh was put. > The model was exported as '.m' > It was called for frequency range 10^(0 : 0.25 : 4) > Following Liu and Thorne, the temperature gradient was integrated over the volume using mphint2() function in Matlab > The time average was taken over the second half of the cycle to avoid transient which are assumed to die off after a couple of initial cycles. > All the prefactors were ignored while plotting > The plot shown below is integral{grad_T}^2} / omega^2   which is proportional to PSD > We just wanted to see if at all we got anything looking similar to the analytic result of Liu and Thorne (shown above) The above plot is when the radius is 100 times the beam spot and height of the cylinder is 200 times the beam spot Other values of radius and height were also tried. The one below shows similar plot as above but for radius 10 times the beam spot and height 20 times the same While the one below is for radius = 1000 * beam spot and height = 2000 * beam spot Here, it is seen that the straight line with a negative slope is seen to some extent. However, the prefactors, when put in, will show the actual PSD. Any comments or suggestions related to the model construction in COMSOL or otherwise is welcome.

Here is a plot showing the profiles of Liu and Thorne and the COMSOL simulation. Since the constants were not taken into account, one expects the plot to be shifted by a constant equal to the logarithm of the prefactor that has been ignored. The green 'o' is just a 1/omega^2 plot in the loglog scale which is how Liu and Thorne's result go (this is not linearized PSD, but is doesn't matter as long as we want to see the profiles). The other plot in blue is the COMSOL simulation.

The plot below is a difference between them. One expects it to be a constant but there are a few jitters. Once again, it is emphasized that the numerical values do not give any information related to the actual difference which will only be calculable
once the correct physical constants are put in the formula. I will post another entry once that is done Attachment 2: L&T_COMSOL_difference.bmp
48   Mon Jul 8 15:47:23 2013 Emory BrownOpticsGeneralMeshing error when varying geometry

I attempted to run the Matlab code built last week, but received an error that the mesh failed to generate on the inner domain/failed to respect edge element on geometry edge.  The error occurred on the second computation with ratioR=0.2, but for ratioR=0.1 the simulation completed successfully.  I reproduced the error in COMSOL and noted that for meshes with size at least as course as "course" or extremely fine, a mesh could be constructed, but for sizes in between an error was returned.  This region of sizes which throw errors varies for different values of ratioR.  By specifying an edge selection of the central line running down the frustum and another along the top edge of the inner cylinder and meshing points on those edges first the rest of the meshing is forced to use those nodes as a seed.

I also noted our COMSOL model had an error: ratioR was applied to the wrong face.  This was corrected as was the .m file.

I added an eigenfrequency study and ran it for several values of ratioR, obtaining the following:

ratioR Lowest Drumhead mode (Hz) Lowest Butterfly mode (Hz)

0.1 3624.8 3802.4

0.3 4375.9 3785.4

0.5 4366.5 3650.0

.75 3709.0 3041.85

1 3000 2186.3

1.1 2736.2 1893.0

1.3 2247.6 1426.4

The mirror's front face is 0.17 [m].  Its Gaussian beam size is 0.0156 [m].  Its height is given by 0.01734[m^3]/(R^2*(1+ratioR+ratioR^2)) such that when the two faces have the same radii it is 0.2 [m] long.

49   Mon Jul 8 18:39:31 2013 Matt A.OpticsGeneralMeshing error when varying geometry

 Quote: I attempted to run the Matlab code built last week, but received an error that the mesh failed to generate on the inner domain/failed to respect edge element on geometry edge.  The error occurred on the second computation with ratioR=0.2, but for ratioR=0.1 the simulation completed successfully.  I reproduced the error in COMSOL and noted that for meshes with size at least as course as "course" or extremely fine, a mesh could be constructed, but for sizes in between an error was returned.  This region of sizes which throw errors varies for different values of ratioR.  By specifying an edge selection of the central line running down the frustum and another along the top edge of the inner cylinder and meshing points on those edges first the rest of the meshing is forced to use those nodes as a seed. I also noted our COMSOL model had an error: ratioR was applied to the wrong face.  This was corrected as was the .m file. I added an eigenfrequency study and ran it for several values of ratioR, obtaining the following: ratioR Lowest Drumhead mode (Hz) Lowest Butterfly mode (Hz) 0.1 3624.8 3802.4 0.3 4375.9 3785.4 0.5 4366.5 3650.0 .75 3709.0 3041.85 1 3000 2186.3 1.1 2736.2 1893.0 1.3 2247.6 1426.4   The mirror's front face is 0.17 [m].  Its Gaussian beam size is 0.0156 [m].  Its height is given by 0.01734[m^3]/(R^2*(1+ratioR+ratioR^2)) such that when the two faces have the same radii it is 0.2 [m] long.

Great Emory,

Can you put these in a plot so it's easier to read?

50   Mon Jul 8 19:04:35 2013 Matt A.OpticsGeneralMeshing error when varying geometry

 Quote: I attempted to run the Matlab code built last week, but received an error that the mesh failed to generate on the inner domain/failed to respect edge element on geometry edge.  The error occurred on the second computation with ratioR=0.2, but for ratioR=0.1 the simulation completed successfully.  I reproduced the error in COMSOL and noted that for meshes with size at least as course as "course" or extremely fine, a mesh could be constructed, but for sizes in between an error was returned.  This region of sizes which throw errors varies for different values of ratioR.  By specifying an edge selection of the central line running down the frustum and another along the top edge of the inner cylinder and meshing points on those edges first the rest of the meshing is forced to use those nodes as a seed. I also noted our COMSOL model had an error: ratioR was applied to the wrong face.  This was corrected as was the .m file. I added an eigenfrequency study and ran it for several values of ratioR, obtaining the following: ratioR Lowest Drumhead mode (Hz) Lowest Butterfly mode (Hz) 0.1 3624.8 3802.4 0.3 4375.9 3785.4 0.5 4366.5 3650.0 .75 3709.0 3041.85 1 3000 2186.3 1.1 2736.2 1893.0 1.3 2247.6 1426.4   The mirror's front face is 0.17 [m].  Its Gaussian beam size is 0.0156 [m].  Its height is given by 0.01734[m^3]/(R^2*(1+ratioR+ratioR^2)) such that when the two faces have the same radii it is 0.2 [m] long.

When we talked today, I thought that these frequencies were a little low, so I ran a really simple COMSOL model of a cylinder with radius 0.17 m and height 0.2 m (attached).

I solved for the first 10 modes (the first 6 being trivial), and the frequencies I got are listed below:

Mode #          Freq [Hz]          Type

7                      5950                 Butterfly

8                      5950                 Butterfly

10                   8226                 Shear

These are more in the frequency range that I expected.

**What material are you using? I used fused silica.

**I used 'finer' meshing. Perhaps your more complicated meshing is screwing with the solutions, although I'd expect it to make it stiffer, not softer.

**Are you applying force during the eigenvalue calculation? That might make it softer.

**Perhaps some of your constraints are moving the modes around. My model has zero constraints: it's just a cylinder floating in space.

Attachment 1: SimpleCylinder.mph
51   Tue Jul 9 14:45:04 2013 Emory BrownOpticsGeneralMeshing error when varying geometry

I tried removing the constraints and applied force and was able to get close to the same frequencies for the butterfly modes as the simple model.  By also removing the interior cylinder from our model I was able to restore the butterfly modes as well.  I can see why removing the constraints and forces would make a difference, but I don't understand how the interior region which was only used to aid in making the meshing in the central region more dense would have such a strong effect. With that removed as well, we get the following results using a finer mesh (plot attached):

ratioR  Lowest Drumhead mode (Hz) Lowest Butterfly mode (Hz) Any lower modes (Hz)

0.1 8777.8 7253.3 3059.7

0.3 5636.9 7264.7 4202.3

0.5 6068.5 7308.0 5980.9

0.75 7341.5 7183.0 no

1 8110.2 5951.3 no

1.1 7302.8 5308.8 no

1.3 5744.3 4087.9 no

This data is more in line with what we would expect, and indicates that the optimal region is probably between ratioR=0.6 and ratioR=0.8.  I will run the looping script we designed to map out how Umax changes when varying ratioR, particularly in this interesting region.

This line was used to generate the plot from data.

figure(1); plot(ratio,drum); hold on; plot(ratio,butterfly,'-rs'); plot(ratio(1:3),others,'g'); hleg1=legend('Drumhead','Butterfly','Other'); title('Eigenfrequencies of mode types for different ratioR values'); xlabel('ratioR'); ylabel('lowest eigenfrequency'); shg;

Quote:

 Quote: I attempted to run the Matlab code built last week, but received an error that the mesh failed to generate on the inner domain/failed to respect edge element on geometry edge.  The error occurred on the second computation with ratioR=0.2, but for ratioR=0.1 the simulation completed successfully.  I reproduced the error in COMSOL and noted that for meshes with size at least as course as "course" or extremely fine, a mesh could be constructed, but for sizes in between an error was returned.  This region of sizes which throw errors varies for different values of ratioR.  By specifying an edge selection of the central line running down the frustum and another along the top edge of the inner cylinder and meshing points on those edges first the rest of the meshing is forced to use those nodes as a seed. I also noted our COMSOL model had an error: ratioR was applied to the wrong face.  This was corrected as was the .m file. I added an eigenfrequency study and ran it for several values of ratioR, obtaining the following: ratioR Lowest Drumhead mode (Hz) Lowest Butterfly mode (Hz) 0.1 3624.8 3802.4 0.3 4375.9 3785.4 0.5 4366.5 3650.0 .75 3709.0 3041.85 1 3000 2186.3 1.1 2736.2 1893.0 1.3 2247.6 1426.4   The mirror's front face is 0.17 [m].  Its Gaussian beam size is 0.0156 [m].  Its height is given by 0.01734[m^3]/(R^2*(1+ratioR+ratioR^2)) such that when the two faces have the same radii it is 0.2 [m] long.

When we talked today, I thought that these frequencies were a little low, so I ran a really simple COMSOL model of a cylinder with radius 0.17 m and height 0.2 m (attached).

I solved for the first 10 modes (the first 6 being trivial), and the frequencies I got are listed below:

Mode #          Freq [Hz]          Type

7                      5950                 Butterfly

8                      5950                 Butterfly

10                   8226                 Shear

These are more in the frequency range that I expected.

**What material are you using? I used fused silica.

**I used 'finer' meshing. Perhaps your more complicated meshing is screwing with the solutions, although I'd expect it to make it stiffer, not softer.

**Are you applying force during the eigenvalue calculation? That might make it softer.

**Perhaps some of your constraints are moving the modes around. My model has zero constraints: it's just a cylinder floating in space.

Attachment 1: Eigenfrequency.png 52   Wed Jul 10 14:15:20 2013 Emory BrownOpticsGeneralComparison of Umax and the frequency of the lowest eigenmode against ratioR

I was able to run the looper program designed last week over the interval ratioR=[0.5, 0.9] with a step size of 0.01 and obtained the plot below with a minimal value of  Umax=1.4723e-10 J obtained for ratioR=0.71, which is about 3.2% less than the value of Umax=1.5202e-10 J when ratioR=1.  This is not very significant, but the shifting upward of the frequencies of the lowest eigenmodes as described in the previous post is.  Additionally, these curves both have their optimal points in the same region, with a ratioR value of about 0.7.

I also ran a higher resolution test to see how the frequency of the lowest real eigenmode responded to changes in ratioR over this region.  The plot generated is also attached with a maximal (most favorable) value of 7210.0 Hz at ratioR=0.74.

Data from both of these runs is included in the attached tarball stored as txt files delineated by tabs.

The next few things to try are substituting silicon for fused silica as the material in the model and seeing how the curves change and seeing the response from reversing the orientation of the optimal test mass (ratioR~0.72).  I expect that the orientation change will result in much higher noise and a lower minimal real eigenfrequency.

Attachment 1: UmaxValuesVsRatioR.png Attachment 2: FrequencyLowestEigenmodes.png Attachment 3: data.tar.gz
53   Wed Jul 10 17:52:18 2013 KojiGeneralGeneralHow to force time steps for a time dependent simulation

A common trap in COMSOL:

When you run a time series simulation, COMSOL tries to be "clever" by skipping some time steps
and returns you an interpolated result in order to accelerate the calculation.
But for certain types of applications, such as time series response analysis by applying periodic force/heat,
this "clever" adaptive algorithm messes up our calculation.

You can apply forced time steps by giving a specific option:
- In Model Builder, open your study. Then open your Solver
- Click "Time-Dependent Solver" in Model Builder
- You get bunch of options in Time-Dependent Solver pane.
Open "Time Stepping".
- Change "Method" from "Generalized alpha" to "BDF".
- Change "Steps taken by solver" to "Strict"

Then, you should be OK.

This can be done by

`model.sol('sol1').feature('t1').set('tstepsbdf', 'strict');`

in COMSOL with MATLAB. (Of course `sol1` and `t1` should be changed appropriately.)

54   Thu Jul 11 19:16:12 2013 Deep ChatterjeeGeneralGeneralTE noise problem - Matlab script

The present aim is the calculation of the TE noise like Liu and Thorne. Although, the results of Liu and Thorne are for the case of an infinite mirror, in our case, we try
to model a test mass with large dimensions in an attempt to get closer to the result due to Liu and Thorne

Following my previous post - The presents results given by the COMSOL simulations gave the same profile as Liu and Thorne's result which has a frequency dependence
of 1/omega^2 but was displaced from the analytic result by a constant factor which was dependent on the applied pressure on the face of the test mass which should not be
the case.

Following the steps followed Liu and Thorne, we had constructed the test mass as a cylinder, large compared to the beam spot size. An oscillating pressure was applied
to one of the faces and the from the temperature gradient generated in the process due to strains, one can calculate the work dissipated. The process involves integrating
the gradient of temperature over the entire geometry and taking the time avrage. However, this did not give the correct results, so it was decide to extract the Fourier coefficient
of the signal and perform the integration on it, as done in the case of TR niose by Koji Arai.

I mention the steps it was done -

> The data extracted was stored in a 3D array in matlab. using mphinterp - each dimension for r, z, t

> The number of radial slices and z slices is defined by the user previously

> For each r,z value there was a time signal.

> The fourier coefficient corresponding to the pressure oscillation frequency was extracted from the last three cycles
(In theory only one cycle should suffice)

> The above step was done so that random data( which is small in magnitude) generated by COMSOL is avoided.

> The fourier coefficients are stored in a 2D array corresponding to r and z values.

> The integration was performed using trapz twice once on each dimension to get the total volume integral

> The rest of the calculation is same as the previous script - plugging in the prefactors in the formulae and plotting

**The primary problem with this script is that the data extraction takes significant amount of time - For a mirror with 200 times the radius of the beam spot and 200 radial slices, it takes close to 6 minutes
to evaluate work dissipated for each frequency. The script solves for 16 frequency values. The results for small number of radial slices does not follow the straight line profile. For larger number of slices
the program takes a longer time. The results have not been checked yet.

Attachment 1: TE_model.m
```% TE calculation using 2D axisymmetric model. This script calls the model
% file, extracts grad_T and does the calculations for the different
% frequency
% All parameter are defined in this script and passed to the COMSOL model
% via arguments
% Parameters, Quantitites in [m]
%Define the radius of the mirror as 100 times that of the beam radius
param.mirror.height = 200 * param.beam.radius;
```
... 94 more lines ...
Attachment 2: almost_infinite_2D_axisym_TE_4.m
```function out = almost_infinite_2D_axisym_TE_4(param)
% Putting all the parameters in the variables
H = param.mirror.height;
dt = param.COMSOL.simulation_time_step;
fmod = param.COMSOL.fmod;
t_end = param.COMSOL.simulation_end_time;
p0 = param.constant.F0/(pi*r0^2);% The amplitude of Osc. pressure
T_amb = param.constant.T_amb;
```
... 121 more lines ...
55   Fri Jul 12 00:09:17 2013 Emory BrownOpticsGeneralAnalysis of substrate Brownian noise in a silicon test mass with changing ratioR value

I performed the same analysis we previously did on a test mass made of silicon (using the built-in silicon within COMSOL) and obtained the attached plots.  The plots show more of a difference between the two solutions optimal points than for fused silica, but they are still quite close.  All data is included in the attached tarball.  The values for ratioR=1 or a cylindrical test mass place the frequency of the lowest real eigenmode at 8491 Hz and the Umax value of  6.5484*10^-11 J.

Attachment 1: SiliconEigenfrequencyPlot.png Attachment 2: UmaxPlotSilicon.png Attachment 3: Silicon.tar.gz
56   Fri Jul 12 00:31:11 2013 Emory BrownOpticsGeneralReversed face mirror design

We had noticed that when we ran an eigenfrequency analysis, there was often a node (see attached screenshot) on the face opposite the applied force.  This indicated that the noise might be lowered by reversing the face the force was applied to, but it also seems like this would be equivalent to increasing the ratioR value and making the primary face smaller while retaining the mass of the test mass.  It seems worthwhile to see how the Umax values compare for a beam incident on either face of a test mass in the optimal region with Umax=0.71.  For the initial configuration, Umax=1.47228*10^-10 J and changing the side of the mirror the force is applied to we get Umax=1.47148*10^-10 J, which is about the same value.

Attachment 1: node.png 57   Fri Jul 12 19:43:32 2013 Deep ChatterjeeGeneralGeneralTE noise problem - Matlab script

 Quote: The present aim is the calculation of the TE noise like Liu and Thorne. Although, the results of Liu and Thorne are for the case of an infinite mirror, in our case, we try to model a test mass with large dimensions in an attempt to get closer to the result due to Liu and Thorne Following my previous post - The presents results given by the COMSOL simulations gave the same profile as Liu and Thorne's result which has a frequency dependence of 1/omega^2 but was displaced from the analytic result by a constant factor which was dependent on the applied pressure on the face of the test mass which should not be the case. Following the steps followed Liu and Thorne, we had constructed the test mass as a cylinder, large compared to the beam spot size. An oscillating pressure was applied to one of the faces and the from the temperature gradient generated in the process due to strains, one can calculate the work dissipated. The process involves integrating the gradient of temperature over the entire geometry and taking the time avrage. However, this did not give the correct results, so it was decide to extract the Fourier coefficient of the signal and perform the integration on it, as done in the case of TR niose by Koji Arai. I mention the steps it was done - > The data extracted was stored in a 3D array in matlab. using mphinterp - each dimension for r, z, t > The number of radial slices and z slices is defined by the user previously > For each r,z value there was a time signal. > The fourier coefficient corresponding to the pressure oscillation frequency was extracted from the last three cycles (In theory only one cycle should suffice) > The above step was done so that random data( which is small in magnitude) generated by COMSOL is avoided. > The fourier coefficients are stored in a 2D array corresponding to r and z values. > The integration was performed using trapz twice once on each dimension to get the total volume integral > The rest of the calculation is same as the previous script - plugging in the prefactors in the formulae and plotting   **The primary problem with this script is that the data extraction takes significant amount of time - For a mirror with 200 times the radius of the beam spot and 200 radial slices, it takes close to 6 minutes to evaluate work dissipated for each frequency. The script solves for 16 frequency values. The results for small number of radial slices does not follow the straight line profile. For larger number of slices the program takes a longer time. The results have not been checked yet.

I am replying to report of the modification made to the script from the last one. The problem of too many radial and z slices have been avoided by using a slicing based on a gaussian
distribution both along the r and z directions - i.e. the slicing is heavy in the region around the surface is applied and much thinner towards the edges of the cylinder. The slicing has
been separately handled by a function named giveSlices() which does the slicing and returns the values as two arrays corresponding to r and z.
The extent to which the fineness of slicing is controlled in the geometry is by controlling the parameter SD in the code which is the standard deviation of the gaussian according to which
the slicing is controlled as being fine throughout or fine in the centre and thin with increasing r.

This method reduces the number of slices by almost 2 orders for significant large values of the dimensions of the test mass while it is expected to give similar results since slicing the
test mass finely at the edge is not required as the gradient of temperature almost falls to zero.

It is suspected that the long time of simulation is attributed to calling mphinterp() command a large number of times. In my next modification, it will be tried to use the command just once
followed by the proper data extraction from the output given by the command.

Attachment 1: codes_TE_calc.zip
58   Sun Jul 14 15:13:46 2013 Deep ChatterjeeGeneralGeneralTE noise problem - Matlab script

Quote:

 Quote: The present aim is the calculation of the TE noise like Liu and Thorne. Although, the results of Liu and Thorne are for the case of an infinite mirror, in our case, we try to model a test mass with large dimensions in an attempt to get closer to the result due to Liu and Thorne Following my previous post - The presents results given by the COMSOL simulations gave the same profile as Liu and Thorne's result which has a frequency dependence of 1/omega^2 but was displaced from the analytic result by a constant factor which was dependent on the applied pressure on the face of the test mass which should not be the case. Following the steps followed Liu and Thorne, we had constructed the test mass as a cylinder, large compared to the beam spot size. An oscillating pressure was applied to one of the faces and the from the temperature gradient generated in the process due to strains, one can calculate the work dissipated. The process involves integrating the gradient of temperature over the entire geometry and taking the time avrage. However, this did not give the correct results, so it was decide to extract the Fourier coefficient of the signal and perform the integration on it, as done in the case of TR niose by Koji Arai. I mention the steps it was done - > The data extracted was stored in a 3D array in matlab. using mphinterp - each dimension for r, z, t > The number of radial slices and z slices is defined by the user previously > For each r,z value there was a time signal. > The fourier coefficient corresponding to the pressure oscillation frequency was extracted from the last three cycles (In theory only one cycle should suffice) > The above step was done so that random data( which is small in magnitude) generated by COMSOL is avoided. > The fourier coefficients are stored in a 2D array corresponding to r and z values. > The integration was performed using trapz twice once on each dimension to get the total volume integral > The rest of the calculation is same as the previous script - plugging in the prefactors in the formulae and plotting   **The primary problem with this script is that the data extraction takes significant amount of time - For a mirror with 200 times the radius of the beam spot and 200 radial slices, it takes close to 6 minutes to evaluate work dissipated for each frequency. The script solves for 16 frequency values. The results for small number of radial slices does not follow the straight line profile. For larger number of slices the program takes a longer time. The results have not been checked yet.

I am replying to report of the modification made to the script from the last one. The problem of too many radial and z slices have been avoided by using a slicing based on a gaussian
distribution both along the r and z directions - i.e. the slicing is heavy in the region around the surface is applied and much thinner towards the edges of the cylinder. The slicing has
been separately handled by a function named giveSlices() which does the slicing and returns the values as two arrays corresponding to r and z.
The extent to which the fineness of slicing is controlled in the geometry is by controlling the parameter SD in the code which is the standard deviation of the gaussian according to which
the slicing is controlled as being fine throughout or fine in the centre and thin with increasing r.

This method reduces the number of slices by almost 2 orders for significant large values of the dimensions of the test mass while it is expected to give similar results since slicing the
test mass finely at the edge is not required as the gradient of temperature almost falls to zero.

It is suspected that the long time of simulation is attributed to calling mphinterp() command a large number of times. In my next modification, it will be tried to use the command just once
followed by the proper data extraction from the output given by the command.

Following the previous reply regarding the use of mphinterp(), the appropriate changes in code were made. The overhead time due to the call to COMSOL was responsible for the longer
simulation time. Right now, for relatively large no. of total slices (about 27000) just the extraction of data happens within a few seconds. However, for the larger number of slices, the codes
will probably still take some time to complete. The codes are attached.

Attachment 1: codes_TE_calc_(2).zip
59   Mon Jul 15 13:14:30 2013 Emory BrownOpticsGeneralUmax curves for varying height test masses

I ran sweeps over the ratioR range from 0.5 to 0.9 for test masses with modified height, and thus mass, which were otherwise equivalent to our previous test mass.  The attached plot was generated and indicates that the optimal ratioR value for the test mass was shifted significantly, for the half mass test mass into lower frequency regions, and for the 10x mass test mass into higher frequency ranges.  I took the data collected from the half mass run and fit a polynomial of degree 2 to it in Matlab using the command p=polyfit(x,y,2) and obtained p=1.089*10^-10*x^2 - 1.024*10^-10*x +1.774*10^-10 which has a minimum value for ratioR=0.47.  This fits with the plot generated as it appears to be near its minumum value, but I should run the simulation more around this value to see if it fits the expected behaviour.  I also did this fitting for the 10x mass curve, but based on the greater variance in the curve over this region of ratioR values I don't expect the predicted minimum to be as accurate, so if we want to find it we will have to search that region of ratioR values. For the 10x mass test mass, p=4.459*10^-10*x^2 - 8.859*10^-10*x + 6.289*10^-10 which has a minimum of about 1.

I will run simulations in these regions overnight, and edit this post with the plots/results.

Attachment 1: MassComparison.png 60   Mon Jul 15 16:05:29 2013 Deep ChatterjeeGeneralGeneralTE noise problem - Matlab script

Quote:

Quote:

 Quote: The present aim is the calculation of the TE noise like Liu and Thorne. Although, the results of Liu and Thorne are for the case of an infinite mirror, in our case, we try to model a test mass with large dimensions in an attempt to get closer to the result due to Liu and Thorne Following my previous post - The presents results given by the COMSOL simulations gave the same profile as Liu and Thorne's result which has a frequency dependence of 1/omega^2 but was displaced from the analytic result by a constant factor which was dependent on the applied pressure on the face of the test mass which should not be the case. Following the steps followed Liu and Thorne, we had constructed the test mass as a cylinder, large compared to the beam spot size. An oscillating pressure was applied to one of the faces and the from the temperature gradient generated in the process due to strains, one can calculate the work dissipated. The process involves integrating the gradient of temperature over the entire geometry and taking the time avrage. However, this did not give the correct results, so it was decide to extract the Fourier coefficient of the signal and perform the integration on it, as done in the case of TR niose by Koji Arai. I mention the steps it was done - > The data extracted was stored in a 3D array in matlab. using mphinterp - each dimension for r, z, t > The number of radial slices and z slices is defined by the user previously > For each r,z value there was a time signal. > The fourier coefficient corresponding to the pressure oscillation frequency was extracted from the last three cycles (In theory only one cycle should suffice) > The above step was done so that random data( which is small in magnitude) generated by COMSOL is avoided. > The fourier coefficients are stored in a 2D array corresponding to r and z values. > The integration was performed using trapz twice once on each dimension to get the total volume integral > The rest of the calculation is same as the previous script - plugging in the prefactors in the formulae and plotting   **The primary problem with this script is that the data extraction takes significant amount of time - For a mirror with 200 times the radius of the beam spot and 200 radial slices, it takes close to 6 minutes to evaluate work dissipated for each frequency. The script solves for 16 frequency values. The results for small number of radial slices does not follow the straight line profile. For larger number of slices the program takes a longer time. The results have not been checked yet.

I am replying to report of the modification made to the script from the last one. The problem of too many radial and z slices have been avoided by using a slicing based on a gaussian
distribution both along the r and z directions - i.e. the slicing is heavy in the region around the surface is applied and much thinner towards the edges of the cylinder. The slicing has
been separately handled by a function named giveSlices() which does the slicing and returns the values as two arrays corresponding to r and z.
The extent to which the fineness of slicing is controlled in the geometry is by controlling the parameter SD in the code which is the standard deviation of the gaussian according to which
the slicing is controlled as being fine throughout or fine in the centre and thin with increasing r.

This method reduces the number of slices by almost 2 orders for significant large values of the dimensions of the test mass while it is expected to give similar results since slicing the
test mass finely at the edge is not required as the gradient of temperature almost falls to zero.

It is suspected that the long time of simulation is attributed to calling mphinterp() command a large number of times. In my next modification, it will be tried to use the command just once
followed by the proper data extraction from the output given by the command.

Following the previous reply regarding the use of mphinterp(), the appropriate changes in code were made. The overhead time due to the call to COMSOL was responsible for the longer
simulation time. Right now, for relatively large no. of total slices (about 27000) just the extraction of data happens within a few seconds. However, for the larger number of slices, the codes
will probably still take some time to complete. The codes are attached.

The dependence of the Power Spectral Density on the parameters F0, which is the pressure amplitude, has been somewhat corrected for in the codes attached. The mistake lay in the fact that
extraction of the Fourier coefficients was being done for the ht.gradTmag, the Fourier coefficients of which are not the correct ones to be used in this case since ht.gradTmag is not the correct time signal.
In the present case the Fourier coefficients are extracted from the r and z components, squared and added to carry out the calculation. However, the results are still off from the L&T results by order of 5.
The codes and the plots are attached. The red values are the L&T result while the blue sketch is the simulation

Attachment 1: codes_TE_calc.zip
61   Tue Jul 23 18:26:04 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses.

The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other
geometries.

The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands
like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB.

Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time
are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne).

The errors expected are
> Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors

> The numerical errors by COMSOL are therefore not filtered off.

> The plot differs from the analytic solution for larger frequencies over 3000 Hz.

> It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In
fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent
in the log scale plots

The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Attachment 2: TE_calc_Jul23.zip
Attachment 3: Jul_23.eps 67   Tue Jul 23 20:53:37 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 66   Tue Jul 23 20:53:37 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 65   Tue Jul 23 20:53:37 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 64   Tue Jul 23 20:53:37 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 63   Tue Jul 23 20:53:37 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 62   Tue Jul 23 20:53:37 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 68   Tue Jul 23 20:53:45 2013 Deep ChatterjeeOpticsGeneralComparison between Liu and Thorne Results and COMSOL results for TE noise

 Quote: In this post I report of the results of TE noise simulated by COMSOL for the TE noise of Infinite test masses. The aim was to follow the procedure by Liu and Thorne in their analytic calculations so that the same model could be used for the other geometries. The simulation is done in a different way than the TR simulations. It was observed that the output given by COMSOL by the use of commands like mphinterp() or taking an export resulted in certain discrepancies between the results computed in COMSOL and that read by MATLAB. Thus, the volume integration of the temperature gradient is performed in COMSOL itself and the results of the integration for each time are sent to files. Matlab read these values and time averages them to get the result as in the paper (Sec. 2 of Liu and Thorne). The errors expected are > Fourier analysis is not done at all. This would have involved exporting data which, as mentioned before is giving errors > The numerical errors by COMSOL are therefore not filtered off. > The plot differs from the analytic solution for larger frequencies over 3000 Hz. > It is to be noted from the paper by Liu and Thorne that the TE noise for the finite and infinite case are not very different. In fact the correction factor goes as O(1). Thus, differences between finite and inifinte cases are unlikely to be prominent in the log scale plots    The codes are put as a zip file. Corrections made to the codes will be uploaded as a reply. Here is another plot with a mesh size slightly finer than the default Extra fine mesh in COMSOL.

One may notice that the value for the final frequency i.e. 10000Hz is different from the previous plot.
It maybe that the error for the higher frequencies is a result of the FEA. However, it may also be that
the appropriate boundary conditions required for an infinite model break down at high frequencies.

Attachment 1: Jul_23_bettermesh2.eps 69   Wed Jul 24 21:08:24 2013 Deep ChatterjeeOpticsGeneralTR results for different dimensions

In this post I simulate the procedure of calculating the TR noise for finite cavities as proposed by Heinert and check for a
match.

The technique of performing all necessary calculations in COMSOL and exporting the results was applied to the TR codes.
It was seen that the codes gives similar output as the technique of extraction of Fourier coefficients in place of time averaging
as has been done in the codes of Koji Arai. One can see the output as the present code runs to be similar to the previous ones
found in the SVN.

However, the results in the present case were off by a constant factor close to 100. This maybe due to some 'm' - 'cm' or similar difference between
analytic calculations and COMSOL values of parameters. Although, it has not been found yet, the correction is hopeful to be
found soon.

The codes give results similar to the analytic result for other values of the mirror radius and beam radii (apart from the constant
factor I have mentioned above). One may have a look at the trend of the graphs between analytic and simulated values in the plots
attached. These plots are for the case when the mirror radius = 25m while the beam radius = 9 cm i.e. the original radii were 0.25m
and 9cm respectively i.e. the ratio has been changed by a order of 2.

As mentioned before the reason for the constant factor difference will be looked into.

Attachment 1: Jul24.pdf Attachment 2: Jul24_acalc.pdf 70   Thu Jul 25 11:36:17 2013 ranaOpticsGeneralTR results for different dimensions

PDF please instead of EPS or BMP or JFIF or TARGA or GIF or ascii art.

71   Thu Jul 25 13:21:46 2013 Deep ChatterjeeOpticsGeneralTR results for different dimensions

 Quote: In this post I simulate the procedure of calculating the TR noise for finite cavities as proposed by Heinert and check for a match. The technique of performing all necessary calculations in COMSOL and exporting the results was applied to the TR codes. It was seen that the codes gives similar output as the technique of extraction of Fourier coefficients in place of time averaging as has been done in the codes of Koji Arai. One can see the output as the present code runs to be similar to the previous ones found in the SVN. However, the results in the present case were off by a constant factor close to 100. This maybe due to some 'm' - 'cm' or similar difference between analytic calculations and COMSOL values of parameters. Although, it has not been found yet, the correction is hopeful to be found soon. The codes give results similar to the analytic result for other values of the mirror radius and beam radii (apart from the constant factor I have mentioned above). One may have a look at the trend of the graphs between analytic and simulated values in the plots attached. These plots are for the case when the mirror radius = 25m while the beam radius = 9 cm i.e. the original radii were 0.25m and 9cm respectively i.e. the ratio has been changed by a order of 2. As mentioned before the reason for the constant factor difference will be looked into.

The discrepancy related to the difference between the analytic and COMSOL results has been partially addressed. Attached is another
plot showing the comparison. The ratio this time between the COMSOL results and the analytic results is between 0.7 - 0.8. This difference
will be looked into. It is, however, observed that the difference is not a constant factor - it has to do with the model file.

Attachment 1: Jul_25.pdf 72   Thu Jul 25 15:54:58 2013 Deep ChatterjeeOpticsGeneralTR results for different dimensions

Quote:

 Quote: In this post I simulate the procedure of calculating the TR noise for finite cavities as proposed by Heinert and check for a match. The technique of performing all necessary calculations in COMSOL and exporting the results was applied to the TR codes. It was seen that the codes gives similar output as the technique of extraction of Fourier coefficients in place of time averaging as has been done in the codes of Koji Arai. One can see the output as the present code runs to be similar to the previous ones found in the SVN. However, the results in the present case were off by a constant factor close to 100. This maybe due to some 'm' - 'cm' or similar difference between analytic calculations and COMSOL values of parameters. Although, it has not been found yet, the correction is hopeful to be found soon. The codes give results similar to the analytic result for other values of the mirror radius and beam radii (apart from the constant factor I have mentioned above). One may have a look at the trend of the graphs between analytic and simulated values in the plots attached. These plots are for the case when the mirror radius = 25m while the beam radius = 9 cm i.e. the original radii were 0.25m and 9cm respectively i.e. the ratio has been changed by a order of 2. As mentioned before the reason for the constant factor difference will be looked into.

The discrepancy related to the difference between the analytic and COMSOL results has been partially addressed. Attached is another
plot showing the comparison. The ratio this time between the COMSOL results and the analytic results is between 0.7 - 0.8. This difference
will be looked into. It is, however, observed that the difference is not a constant factor - it has to do with the model file.

The issue related to the difference between the analytic and simulated values has been resolved. The codes seems to give reasonable match
between the analytic and simulated case. There is, however, a difference between the formulas being used from the previous cases. Note that
the 1/2 in front of Eq.(15) of Heinert is a because the time average has already been considered. However, in the present codes, the volume
integral of grad_T is evaluated in COMSOL and exported as a function of time. It is then averaged in MATLAB. Thus the factor of 1/2 is to be
omitted in this case(see Liu and Thorne Eq.(5). The presence of this extra factor of 1/2 was giving error in the last upoaded plots. From the
relative difference plot, one can see the maximum difference between COMSOL and analytic results go upto 7% but for most of the graph
it is close to 1% which is a fair result.

Attachment 1: Jul_25_1.pdf Attachment 2: Jul_25_2.pdf Attachment 3: relative_plot.pdf 73   Mon Jul 29 22:42:57 2013 Deep ChatterjeeOpticsGeneralAvoiding transient solutions in the Computation of TE/TR noise

An error was being encountered in the computation of the TR noise lately. It was observed that while running the simulations in the case of the materials which have a lower value of the thermal diffusivity (silicon / sapphire at room temperature), the simulated result were slightly off from the analytic result. On the other hand, if the simulation was run with a material of higher diffusivity(same materials at lower temperature), the match would be better. The reason being the transient solution not dying off significantly during the period of the simulation. Since a time average was being taken of the quantity integral{grad_T ^2},  the transient contributed to the integral. To get the correct value, the fourier coefficient of the time signal of integral{grad_T ^2},  was extracted at twice the frequency of the pressure oscillation. The reason being that the signal was squared. Extracting the response at this frequency after the integration is logical since the integration is over space while the response we extract is over time.

The same procedure was also applied to the TE noise calculation. However, this time we obtained similar result as the case where this procedure was not applied but a simple time averaging was performed. The tail of the plot, at high frequency, is still seen to deviate from the analytic result of Liu and Thorne as was the case previously. A plot is attached showing the spectrum for Fused silica at 290K. The conclusion being that the transients do not affect the TE noise calculation - the plot stayed the same even after filtering them out. This is probably because unlike the TR case which has a heat source present along the cylinder axis, the TE noise calculation involves applying pressure only on the face of the cylinder, and the transient do not contribute much to the volume integral.

Attachment 1: Jul29_2.pdf 74   Wed Jul 31 15:39:11 2013 Deep ChatterjeeOpticsGeneralFirst try with paramter optimization for TE and TR noise profiles

After the simulations have been found to match to a fair extent with the analytic results by Heinert and Liu and Thorne, the attempt is check out the parameters for which the TE and TR noise are
close to each other. This was done with the analytic results. The frequency range 10 - 1000 Hz was looked into. In between this range, the quantity that was minimized was the absolite value
of the logarithm of the ratio between the TR and TE noise. The fminsearch function was used to minimize the mentioned quantity. The parameters that were changed were - conductivity, thermo
optic coefficient and coefficient of linear expansion. The reason for choosing these three were -

> TR noise is independent of coefficient of linear expansion

> TE noise is independent of thermo optic coefficient

> The power law dependence on conductivity is different for the TE and TR cases as can be seen from the analytic expressions

once the code returned the optimized parameters, these values were plugged in and the results were plotted.

**Note that the minimization was done for frequencies between 10 to 1000 Hz

Attachment 1: Jul31_param_opt.pdf 75   Thu Aug 8 17:17:19 2013 Deep ChatterjeeOpticsGeneralSomething like cancellation

For the material parameters of Sapphire at 300K, the TE and TR Noise profiles, though not very close, lie close to within an order of magnitude. Sapphire has a positive coefficient of linear expansion. We just inverted the sign of this quantity
and the ran the codes that puts heat and pressure simultaneously to the test mass. The total noise looks to be lower than the TR noise which is greater.

Mirror height 0.46[m]

If we have physical parameters which make the two Noise sources come closer to each other and then flip the sign of alpha, we may be able to see some noise reduction to a greater extent.

Attachment 1: suspected_cance_AUg8l.pdf 86   Tue Jun 24 14:35:42 2014 Sam MooreOpticsGeneralTrying to Verify the Heinert Model
Attachment 1: 06_23_14.pdf 87   Tue Jun 24 17:05:24 2014 Sam MooreOpticsGeneralTrying to Verify the Heinert Model

It does appear that the simplified model is only relevant for the simulations.  To quote Heinert: "An efficient computation is only possible for the simple model, as the advanced model would require an element of size more than 106 ."  I have run Koji's code that replicates Heinert's figure 3.  I have attached the resulting temperature distribution and noise amplitude curve.  In the noise amplitude curve, the red line is the analytical result, while the dots are from COMSOL.

The next step is to convert this code to an efficient complex time-independent solution.  As stated before, my main concern here is whether COMSOL actually solves the right equation in the stationary case.

Attachment 1: temperature.png Attachment 2: noiseAmplitude_agreement.png 88   Sat Jun 28 21:59:11 2014 Sam MooreOpticsGeneralDifficulty with the COMSOL stationary module; Test Cases

Here, I describe some test cases to see if COMSOL's solutions are agreeing with some simple analytical solutions.  Right now, I have two plots showing COMSOL's solution and my analytical solution on separate plots.  I will be plotting there difference to see if they really match up.

Attachment 1: 6_27_14.pdf  89   Sun Jun 29 15:37:18 2014 Sam MooreOpticsGeneralDifficulty with the COMSOL stationary module; Test Cases

 Quote: Here, I describe some test cases to see if COMSOL's solutions are agreeing with some simple analytical solutions.  Right now, I have two plots showing COMSOL's solution and my analytical solution on separate plots.  I will be plotting there difference to see if they really match up.

The following document shows the relative difference between these two plots.

Attachment 1: 6_29_14.pdf 90   Sun Jun 29 20:25:44 2014 KojiOpticsGeneralDifficulty with the COMSOL stationary module; Test Cases

What about this example? The result is easier to understand intuitively.

Consider a bar with the length of L.

Let's say there is no body heat applied, but the temperature of the bar at x=L is kept at T=0
and at x=0 is kept at T=T0 Exp[I w t].

The equation for the bar is ...(1)

Consider the solution with the form of T(x, t) = T(x) T0 Exp(I w t), where T(x) is the position dependent transfer function.
T(x) is a complex function.

Eq.1 is modified with T(x) as With the boundary condition of This can be analytically solved in the following form where alpha is defined by So kappa/Cp is the characteristic (angular) frequency of the system.
Here is the example plot for L=1 and alpha = 1 (red), 10 (yellow green), 100 (turquoise), 1000 (blue)  If the oscillation is slow enough, the temperature decay length is longer than the bar length and thus the temperature is linear to the position.
If the oscillation is fast, the decay length is significantly shorter  than the bar length and the temperature dependence on the position is exponential.

Now what we need is to solve this in COMSOL

93   Thu Jul 10 16:51:14 2014 Sam MooreOpticsGeneralDuan and Heinert Comparison

(See Plots in attached document)

My plan has been to replicate Duan's numerical thermoconductive (TE + TR) phase noise plot presented in his paper (section V).  I am trying to match Duan's analytical expression with Heinert's analytical expression.  This requires some rescaling of Heinert's TR displacement noise. (I also needed to divide Heinert's expression by 4 pi to match the Fourier Transform convention. )   Duan's analytical expression for the phase noise is obtained by evaluating the triple integral given in equation 13 of the Duan paper "General Treatment of Thermal Noise in Optical Fibers".

It turns out that an additional factor of 2 multiplies the phase noise because Duan's Fourier Transform only takes into account positive frequencies; there are also negative frequencies that occur in equal amplitude.
This integral was evaluated in Mathematica due to numerical noise in MATLAB's calculation.  The calculation in Mathematica was very slow, so the upper limits on the integral were truncated.  The following plots in the attached document show the resulting noise profile agreements for two different upper limits.

If the residual for the highest upper limit is considered acceptable for a match between the two plots, then I will use Heinert's plot as a reference when using the COMSOL steady-state method for Duan's numerical case (Heinert's plot runs much faster).
Attachment 1: 7_10_14.pdf  95   Mon Jul 14 19:09:14 2014 Sam MooreOpticsGeneralUsing Heinert's Solution for Duan's Parameters

I have plotted Heinert's analytical solution for TR noise using Duan's parameters.  Since TO and TE noise can be found by simply rescaling TR noise, these have been included in the plot as well.  The solid curve represents the analytical solution, while the tick marks represent COMSOL's solution.  I have used COMSOL for both a 1D axisymmetric and a 3D model.  Since Duan's cylinder has a radius of 125 microns, but a length of 1 m, the meshing was difficult for the 3D model.  I ended up shortening the length of the cylinder, converting to the actual length when finally calculating the thermal noise.

Attachment 1: oneD_duanParams_residual-eps-converted-to.pdf Attachment 2: threeD_duanParams_residual-eps-converted-to.pdf 96   Mon Jul 14 19:14:31 2014 Sam MooreOpticsGeneralDuan and Heinert Comparison

 Quote: (See Plots in attached document)   My plan has been to replicate Duan's numerical thermoconductive (TE + TR) phase noise plot presented in his paper (section V).  I am trying to match Duan's analytical expression with Heinert's analytical expression.  This requires some rescaling of Heinert's TR displacement noise. (I also needed to divide Heinert's expression by 4 pi to match the Fourier Transform convention. )   Duan's analytical expression for the phase noise is obtained by evaluating the triple integral given in equation 13 of the Duan paper "General Treatment of Thermal Noise in Optical Fibers". It turns out that an additional factor of 2 multiplies the phase noise because Duan's Fourier Transform only takes into account positive frequencies; there are also negative frequencies that occur in equal amplitude.   This integral was evaluated in Mathematica due to numerical noise in MATLAB's calculation.  The calculation in Mathematica was very slow, so the upper limits on the integral were truncated.  The following plots in the attached document show the resulting noise profile agreements for two different upper limits.   If the residual for the highest upper limit is considered acceptable for a match between the two plots, then I will use Heinert's plot as a reference when using the COMSOL steady-state method for Duan's numerical case (Heinert's plot runs much faster).

I have now plotted the Duan-Heinert comparison for the case of an infinite upper bound.  It turns out that the curves differ by a maximum of 10 percent for low frequencies.  Such a discrepancy has been attributed to lack of experimental investigation into this regime (according to Duan). For our purposes, such a discrepancy is acceptable. We will therefore use the Heinert curve for subsequent calculations due to its faster computation time.

Attachment 1: duan_heinert_comparisonInfinite-eps-converted-to.pdf 97   Thu Jul 31 20:55:38 2014 Sam MooreOpticsGeneralFinding the Right Meshing for the TIR cavity

In this document, I try to identify I good mesh by comparing the numerical solution from that mesh with my analytical model.  Since there are problems with carrying out the analytical calculation, it is still not entirely clear which mesh should be used.

Attachment 1: 7_30_14.pdf     98   Sat Aug 2 00:22:34 2014 Sam MooreOpticsGeneralFinding the Right Meshing for the TIR cavity

 Quote: In this document, I try to identify I good mesh by comparing the numerical solution from that mesh with my analytical model.  Since there are problems with carrying out the analytical calculation, it is still not entirely clear which mesh should be used.

I have refined the analytical calculation procedure, as outlined in this new document.  The procedure indicates that the discrepancy between the analytical and numerical solutions are more likely attributed to meshing inaccuracies.

Attachment 1: 8_1_14.pdf   108   Fri Jul 29 14:33:41 2016 Joy WestlandMechanicsGeneralA Tutorial in Importing SolidWorks Files and Applying a Gaussian Force in ANSYS for a LIGO Test Mass

Here is a tutorial in importing SolidWorks into ANSYS and the steps needed to apply a Gaussian Force to the LIGO test mass that’s imported.

Using SolidWorks:

1. Download the SolidWorks Zip Folders from the DCC
2. Open/export the zip folders into a different file
2. For only the test mass and ears, open the file called: D0902456 ITM OPTIC WITH EARS ASSEMBLY. This will represent the bottom mass.
3. An assembly can be imported into ANSYS but it’s easier to convert an assembly into a part before importing it into ANSYS. Also this will allow you to perform a split line.
1. To do this, resave the assembly but before pressing save, change the file to PRT extension instead of Assembly.
4. Open the PRT file of the test mass, make a new sketch on the face of the test mass.
5. Draw a circle (Or whatever shape) on the test mass. In this case a 0.1 m circle was made in the middle of the test mass.
6. From there go to “Insert” and then press “Curve” and then go to “Split Line”. Split line allows a user to project the sketch onto the part. This makes it possible to press the sketch separately from the part. This is a useful tool because in ANSYS, once split line has been used on a sketch, a user can press that shape separately on the face of the object. The importance of this tool allows a user to click on the circle separately and apply a force in that particular area.
1. In “Split Line” press the “sketch” you want to project onto the “face” of the object.

Importing SolidWorks into ANSYS:

1. Once the part has been made, choose the “static structural analysis”.
2. In the geometry: Import the test mass that was converted into a part
3. Once imported, press the geometry again and “edit”
4. In the edit module, “generate” the part and once the part has loaded notice the circle can be pressed such that it is a different part of the test mass surface.
5. Exit geometry module

Applying Gaussian Force:

1. Enter the “model” module to edit the setup of the analysis
2. Assign the material that you want for the test mass under “Geometry” à “Solid” under the “Material” and “Assignment”
3. Under the “Coordinate Systems”, make a new coordinate system
4. “Coordinate System:
1. Click on the face, in this case the circle made on the test mass for the “geometry” and “apply”
2. Under “Definition” à “Type” à Change from Cartesian to Cylindrical. This allows the Gaussian Force to be distributed correctly
5. Under “Analysis Settings” Turn on “Large Deflection”
6. Add “gravity” under “Static Structural”
1. Treat the test mass is if it were the bottom mass so make sure gravity is pointing the correct way, such that the wires on the ears would extend upward.
7. Insert “Fix” to the top surfaces of the ears
8. Insert “Pressure”
1. Press on the circle made to apply for the pressure
2. Under “Definition” change the “Magnitude” to “Function” (the arrow at the end of the Magnitude entry bar)
3. Once function has been activated: A Gaussian force of: 1/((3.141592654*0.0156^2)*2.718281828^((x/0.0156)^2)) was used in SI units. Note that ANSYS does not use symbols. Once that’s entered into the “Magnitude”, under “Function” the “Coordinate System” will appear. Change that to the coordinate system that was made in step 4 of the cylindrical system.
4. Under “Graph Controls” Make sure that the X-Axis is changed from Time to X.
5. Pick a range for the graph and the number of segments that you want to look at.
9. “Solve” the system and right click on “Solution” and evaluate the different results you need.
Attachment 1: Bottom_Test_Mass_With_Ears.PNG Attachment 2: Applying_Parameters.PNG Attachment 3: Applying_Gaussian_Force.PNG 114   Mon Jul 31 22:18:57 2017 ranaGeneralGeneralusing more than 12 cores in matlab

Since 2014, the limit of 12 workers using the matlab parallel computing toolbox has been lifted. Today, I was able to get this to work. There's a trick.

Usually, when you start up matlab and run a parallel thing like 'parfor', it just uses a default profile 'local' which limits you to 12 workers. You can try to ask for more by doing 'parpool(40)' for 40 workers, but it will tell you that NumWorkers = 12 and you're out of luck. So instead:

`myCluster = parcluster('local')`

`myCluster.NumWorkers = 40;`

`saveProfile(myCluster);`

`parpool('myCluster', 40)`

It seems that it needs the max # of workers and the requested number of workers to be 40 to use 40, otherwise you'll just get 12 (as of matlab 2016a).

Attachment 1: Screen_Shot_2017-07-31_at_10.11.35_PM.png 127   Sat Mar 17 15:27:48 2018 ranaGeneralGeneralfile size >> small

When saving your COMSOL files do these two things to make the files much smaller (good for saving in version control and sharing):

1. File -> Compact History
2. Preferences -> Files -> Optimize for File Size (not speed)
ELOG V3.1.3-