Table Of Content1
HIGHPRECISIONMODELINGTOWARDSTHE10−20 LEVEL
Author:M.Andres1,L.Banz1,A.Costea1,E.Hackmann2,S.Herrmann2,C.La¨mmerzahl2,L.Nesemann1,B.Rievers2,E.P.Stephan1
1:InstituteforAppliedMathematics,LeibnizUniversityHannover,Germany
2:CenterofAppliedSpaceTechnologyandMicrogravityZARM,UniversityofBremen,Germany
Abstract
The requirements for accurate numerical simulation are increasing constantly. Modern high precision physics ex-
1 perimentsnowexceedtheachievablenumericalaccuracyofstandardcommercialandscientificsimulationtools. One
1 exampleareopticalresonatorsforwhichchangesintheopticallengtharenowcommonlymeasuredto10−15precision
0
[1].Theachievablemeasurementaccuracyforresonatorsandcavitiesisdirectlyinfluencedbychangesinthedistances
2
betweenthe opticalcomponents. If deformationsin the rangeof 10−15 occur, those effectscannotbe modeledand
n
analysedanymorewithstandardmethodsbasedondoubleprecisiondatatypes. Newexperimentalapproachespoint
a
J outthattheachievableexperimentalaccuraciesmayimprovedowntothelevelof10−17inthenearfuture[2]. Forthe
developmentandimprovementofhighprecisionresonatorsandtheanalysisofexperimentaldata,newmethodshave
5
to be developedwhich enable the needed level of simulation accuracy. Thereforewe plan the developmentof new
] highprecisionalgorithmsforthesimulationandmodelingofthermo-mechanicaleffectswithanachievableaccuracy
h of10−20. Inthispaperweanalyseatestcaseandidentifytheproblemsonthewaytothisgoal.
p
-
p
m Motivation soundway.
o Todate,thebestcavitystabilizedlasershaveachieved
c Optical high-finesse resonators are widely used as a arelativefrequencystabilityontheorderoffewparts
.
s frequencyreferenceforthestabilizationoflaserse.g. in 1016 [5]. As pointed out by Numata et al. [6,
ic in optical atomic clocks, now aiming for the 10−18 7], this precision is limited by fundamental thermal
s
y levelofprecision [3],orindirecttestsofspecialand noise, mostly affecting the highly reflective mirror
h general relativity. For example, [1] and [4] tested coatings. Consequently,low-noise,single-crystalline
p
the isotropy of the speed of light by comparing the siliconmirrorsarecurrentlybeingdeveloped[2]and
[
frequencies of two lasers referenced to crossed res- an improvement of optical cavity stability down to
v1 onators witharelativeprecision of10−15 at1sinte- the10−17levelmightbecomefeasiblealreadywithin
1 grationtime. thenextyears. Thesecavities willultimately helpto
8 Typically, an electronic feedback loop is applied to buildopticalclocksofimprovedaccuracyortocarry
9
0 actively stabilize the laser frequency to the optical out fundamental tests at the 10−17 level and below.
. resonator. The laser frequency then results directly Thisdevelopment howeverneedstobeparalleled by
1
0 from the geometrical resonator dimensions as given similar advancements in high precision modeling
1 by νL = mc = 2L, where L is the resonator length tools.
1
: and m is the longitudinal mode number. The mea- In order to react to the grown achievable accuracies
v
surement performance thus depends on the length ofopticalresonatorsandtheneedforimprovedmod-
i
X stabilityoftheimplementedcavities. Since∆ν/ν = eling methods, we intend to develop a new method
r ∆L/L, high precision frequency stabilization at the for the modeling of thermo-mechanical effects in 3
a
10−15 level and below, requires a relative cavity dimensions with an achievable numerical accuracy
length stability of the same magnitude. The domi- of 10−20. Inthis paper, weoutline the difficulties to
nant effects which cause expansions or contractions overcomebyanalysing acavity testcase.
are gravity (in particular if an experiment is assem-
bled on earth and conducted in space) and thermal
Test case configuration
expansion of the components. Therefore high preci-
sion modeling of thermo-mechanical effects is a ba-
In the following we will present the cavity test case
sic requirement for a safe design and verification of
which will be used to demonstrate the problems for
optical resonators. Only if the models and simula-
modeling of thermo-mechanical effects with a nu-
tionsofferthesameprecisionwhichcanbeachieved merical accuracy of 10−20. Theanalysis isconfined
in experiments, the optimization, design, and eval-
to a single mirror within the cavity and performed
uation of optical resonators can be performed in a
usingtwodifferentapproaches.
2
Fig.1: Designofcavitytestcaseconfiguration.
Cavity test case space and cannot be deformed. As the mirrors are
bondedtothespacer, adeformation ofthemirrorsat
An optical cavity is typically axially symmetric and
the contact areas is also not allowed. However, heat
composed of two highly reflective concave lenses
conduction atthecontactareaistakenintoaccount.
whicharebondedtoaspacer. Thedimensionsofthe
ULEisdesigned tohaveaminimum thermal expan-
test cavity analyzed here are adapted from Webster
sioncoefficientatsometemperatureTULEnearroom
et al. [8], for the design see Figure 1. The mir-
temperature,whoseexactvaluehastobedetermined
ror substrates, the mirror coatings, and the spacer
experimentally for each cavity. Even if TULE is ex-
areassumedtobemadeofultra-lowexpansionglass
actly known, in practice it is not possible to achieve
(ULE) provided by Corning with the relevant mate-
a stable temperature of TULE over the whole cav-
rial parameters listed in Tab. 1. Although ULE is
ity, and a small temperature gradient is inevitable.
knownforanisotropicthermalexpansion,weassume
However, by careful thermal control a temperature
as an approximation that all materials are homoge-
stability of a few mK can be obtained if operating
nousandisotropicandthatthespacearoundthecav-
near TULE and even about 100µK for temperatures
ityisavacuumatroomtemperature. Further,wetake
around 300K [9]. Here we assume the mirrors sub-
only thermo-elastic effects into account and neglect
stratesandcoatingstohaveastabletemperaturecor-
othererrorsources likevibration ormateriallosses.
responding to a coefficient of thermal expansion of
Materialparameter value 30 × 10−9K−1. For more rigorous treatments than
thistestcasethethermalexpansioncoefficientshould
Thermalconductivity [W/mK] 1.31
be varied in time to reflect the small disturbances of
Specificheat[J/kgK] 767
Mean coefficient of thermal ex- 0±30×10−9 thetemperature.
Also, weassume thatthis testcavity hasafinesse of
pansion [1/K]
Elasticmodulus[Pa] 67.6×109 F = 106 andatransmission coefficientoftheincou-
pling mirror of T = 10−6. For an input power of
Poisson’s ratio 0.17
Massdensity[kg/m3] 2.21×103 Pin = 50µW thisresultswiththegainfactorof
F2
Tab.1: MaterialparametersofULECorningcode7972. T ≈ 105 (1)
π2
inaneffectivepowerofaboutP = 5W storedinthe
Modeled effects, boundaries, and loads cavity. ThelossesAduetoabsorptioncorresponding
π
to the chosen values for F ≈ and T can be
T+A
Regarding the setup of the test cavity, the relevant assumedtoA≈ 10−6resultinginadissipatedpower
thermal physical effects are heat conduction and ra-
ofabout 5µW permirror. Thelaser isassumed here
diation as well as the thermal expansion of the ma-
to have a diameter of 200µm and to hit the mirror
terial. Due to the assumed room temperature, heat
perfectly atitscenter.
radiation may be neglected as a first simplification.
Once an equilibrium state is reached, the tempera-
Thespacermayradiatetothewholesurroundingspace
tureisconstant intimebutmayvaryoverthecavity.
and, thus, we assume it as a heat sink with a con-
However, due to the changed temperature the cavity
stant temperature of 300K as a second simplifica-
willdeform,whatinturnchangestheboundary con-
tion. Third, we assume that the spacer is fixed in
ditions for the heat conduction. Therefore, the final
3
static state has to be determined iteratively. For this maximal cavity temperature
300.014
test case, we neglect the influence of the changes of uniform mesh
adaptive mesh
300.013
boundaryconditionsandperformonestaticanalysis,
300.012
only.
mperature [K]30300.00.1011
FE models and results maximal te300.009
Thetestcasedescribedabovewasmodeledusingthe 300.008
Finite Element (FE) approach. Two independ pro- 300.007
grams were used for this task: the academic version 300.006
0 0.5 1 1.5 2 2.5
ofthecommercialANSYSsoftwareandthesoftware DOF x 105
packageMaiProgsdeveloped attheIfAM.
Fig.2: Max.tempvs.DOFusingANSYS
ANSYS model x 10−11 axial displacement at front of cavity
2
Asafirstapproach,onehalfoftheaxiallysymmetric
0
cavity was modeled in three dimensions. However,
due to the limited node number of 256000 nodes of m] −2
m
the academic version of ANSYS, this turned out to ment [ −4
be nonsufficient for reaching a high modeling ac- displace
curacy. Therefore, the axial symmetry of the prob- −6
lem was used to reduce the model to two dimen-
−8
sions. For the computation of heat conduction and
mechanical deformation, the eight nodes PLANE77 −1−015 −10 −5 0 5 10 15
position on diameter [mm]
and PLANE183 element types were used. On the
spacer a fixed temperature and no displacement are Fig.3: NodalsolutionplotusingANSYS
prescribed and, therefore, we modeled it in ANSYS
with a uniform mesh with a large maximal element
ofthecavityatthefrontandthebackisplottedusing
size of 1mm. On the lens we started with the same
thedataofthecomputations withMaiProgs.
element size and reduced it step by step. However,
Topointoutthesensitivityofthemathematicalmodel
the best results could be obtained by using a mesh
ontheaccuracyofthesolution,computationsontwo
adapted to the problem with maximal element sizes
different mesheswereperformed attheIfAM.Inthe
of approximatly 0.002mm in the regions of interest
firstcasetheusualh-versiononaquasiuniformmesh
andatotal nodenumberofabout247000. Themax-
was computed. The geometry is modeled via two
imum temperatures for the various degrees of free-
cylinders where the ring of the outer cylinder fits
dom are shown in figure 2. After solving the ther-
withtheringofthecavitythatisbondedtothespacer.
mal problem the resulting temperature distribution
Thetetrahedralelementshavenearlythesameshape
was applied as a body load for the following struc-
and size. Since the laser is striking the cavitiy in
tural analysis. From a preliminary analysis the dis-
placementswereexpectedintherangeof10−12mm.
Therefore, as the structural analysis isalinear prob- x 10−11 axial displacement at back of cavity
3
lem,wedecidedtorescalethethermalexpansionco-
efficientforabetter numerical performance. There- 2.5
sulting displacements atfront and back ofthe cavity
2
areshowninfigures3and4. mm]
ment [1.5
MaiProgs model displace 1
At the IfAM computations were made for linear el- 0.5
ements on tetrahedral meshes with up to seven mil-
0
−15 −10 −5 0 5 10 15
lion nodes. For the computations the self developed position on diameter [mm]
FE/BE-softwareMaiProgs[10]wasused. InFigures
Fig.4: NodalsolutionplotusingANSYS
5 and 6 the axial displacements on the cross section
4
x 10−11 axial displacement at front of cavity maximal cavity temperature
2
10−2
0
−2
displacement −4 θmax
−6
10−3
extrapolated
−8 h−version
h−adaptive version
−1−015 −10 −5 0 5 10 15 103 104 105 106 107
position on diameter DOF
Fig.5: NodalsolutionplotusingMaiProgs Fig.7: Max.tempvs.elementsize
x 10−11 axial displacement at back of cavity seemsreasonable toassumethatanadaptivecompu-
3
tationusingsomekindofaposteriorierrorestimator
2.5 (residual, hierarchical, etc.) [13, 14] would exceed-
inglyimprovetheaccuracyofthesolution.
2
ment
displace1.5 Comparison of results
1 Thetwodifferentsoftwaresandapproachesdescribed
inthissection leadtoresultscomparing welltoeach
0.5
other. However,acloserinspectionofthetworesults
0 reveal significant differences: the maximal tempera-
−15 −10 −5 0 5 10 15
position on diameter turesforthetworesultsdifferby
Fig.6: NodalsolutionplotusingMaiProgs −05
∆T = 1.904041846501059 ×10 K.
Giventhisquitelargediscrepancyitisalmostsupris-
its midpoint with a radius of 0.1mm it is obvious
ingthatthemaximaldisplacementsforbothapproaches
thatthe temperature willhavealarge gradient there.
differbyonly
This is also emphasized by the axial displacement
in Figure 5 where the displacement seems to have −12
∆D = 1.447725122807002 ×10 mm.
a singularity in the midpoint of the cavity. Using
this information the IfAM modeled a second mesh These differences are caused ontheone hand bythe
consisting of three cylinders where the inner cylin- discretizationerrorsconnectedtothedifferentmeshes,
der has the radius of the laser such that the coars- and ontheother hand bythelimited accuracy ofthe
est mesh has 400 more points in the neighborhood numerical procedures involved in the solution pro-
of the midpoint than the first version. The tetrahe- cess. Whereas the first error source can be mini-
dral elements do not have the same size any longer, mized by using finer meshes, the second one is an
sincetheelementsnearthemidpointofthecavityare inert problem of floating point arithmetic. As a rule
considerably smaller than the other elements. From of thumb, one can say that for a result with n sig-
this point of view the second mesh can be regarded nificant accurate digits anaccuracy of2nsignificant
as kind of an adaptive mesh. This leads to a bet- digitsforeachadditionandmultiplication isneeded.
ter approximation of the boundary loads and hence Asthe exponent is treated seperatly, for linear prob-
to a better approximation of the solution as can be lems with known order of magnitude for displace-
seen in Figure 7. Here the red dashed line is the ex- ment or temperature the accurate digits can be con-
trapolatedmaximaltemperatureinsidethecavity. On siderablyimprovedbyasuitablerescalingofthema-
bothmeshestheIfAMcomputedfourstepsoftheh- terial properties. However, for nonlinear problems
version [11, 12] starting with a coarsest mesh with this approach will in general not work. As a con-
2400 points and 2800 points, respectively. Figure 7 sequence, the only way to minimize the numerical
showsthatthecomputationsonthesecondmeshlead errorsistousequad(with128bits)orevenarbitrary
tomuchbetterapproximations ofthetemperature. It precision arithmetic instead of the usual double pre-
cision (with 64 bits). Since there is no appropriate
5
percentage thermal energy error percentage structural energy error
100 35
uniform mesh uniform mesh
90 adaptive mesh adaptive mesh
30
80
25
70
percentage error 456000 percentage error1250
30 10
20 5
10
0
0 0 0.5 1 1.5 2 2.5
0 0.5 1 DOF 1.5 2 x 1025.5 DOF x 105
Fig.9: Structuralenergyerrornormplot
Fig.8: Thermalenergyerrornormplot
relative errors in energy norm
hardwareavailable, thesecomputations mustbeper-
formedwithasoftware emulation. Thisslowsdown 100
thecomputation processsignificantly.
Accuracy Evaluation or
err
el.
r
ANSYS
h−version(θ)
Within the ANSYS software a posteriori error esti- 10−1 h−version(u)
h−adaptive(θ)
mator can be used to evaluate the quality of the so- h−adaptive(u)
lution. As the exact solution is not known in gen-
103 104 105 106 107
eraland,inparticular, notforourtestcase,thiserror DOF
estimator compares the (discountinous) heat flux or Fig.10:Errornormplot
stresses, respectively, in each node as computed by
ANSYStoanaveraged(continuous)function,which
sure to compare different approximations. The aim
interpolates thesevalues. IneachelementEi theen-
is to derive efficient and reliable a posteriori error
ergy norm ei of the difference of the discontinuous
estimators for the test case in order to have a better
valuesandtheinterpolated function isameasurefor
prediction of thetrue error. Ofcourse, these estima-
theerrorinthiselement. Forthewholemodelaper-
tors can be used within adaptive algorithms leading
centaged relative energy norm error can be given by
toafasterconvergence [13,14].
the sum of all ei’s divided by Piei and the thermal
dissipation energy. ThisisakindofZienkewics/Zhu
error estimator [15], which converges to the true er-
ror. The percentage error for various degree of free-
Conclusion
domsareplottedinFigures8and9forboththether-
malandstructural computations.
We analysed a cavity test case for thermal induced
changes in the cavity length, which directly influ-
ences the stability of the laser frequency used to in-
terrogatethecavity. Weusedtwodifferent andinde-
pendentsoftwaresforafiniteelementanalysisofthe
thermo-elastic deformation of the cavity lens due to
MaiProgs
thedissipated powerofthelaser.
Since the exact solution of the presented test case is Despitethenumeroussimplificationsassumedforthis
not known in advance it is impossible to compute test case the results of the two independently com-
the exact error of the computed approximations in puted values forthe maximaldisplacement differ by
the energy norm. Therefore Figure in 10 extrapo- an order of magnitude of 10−12, which is consid-
latednormsforthetemperatureandthedisplacement erably lower than the corresponding achievable fre-
arecompared withthe energy norms ofthe approxi- quencystabilitiesofultrastablecavities. Theindivid-
matedsolutions. Theerrordoes notreflecttheaccu- ual analysis have two major error sources: the dis-
racy of the solution but can be considered as a mea- cretization errorofthemeshandthenumerical error
6
ofthesolution process. Thefirsterrorcan beshown J. C. Bergquist. Visible Lasers with Subhertz
toconvergetozeroforaninfinitlyfinemesh. Succes- Linewidths,Phys.Rev.Lett.82,3799-3802.
sively finer meshes can be produced with any com-
[6] K. Numata, A. Kemery, and J. Camp. Thermal-
mercialFEsoftwareandadequatecomputingpower.
Noise Limit in the Frequency Stabilization of
Therefore, this error source can be minimized with-
Lasers with Rigid Cavities, Phys. Rev. Lett. 93,
out principle problems. The second error source is
250602 (2004).
inerttofloatingpointarithmeticandcanonlybemin-
imized if data types with more bits are used. For
[7] M.Notcutt, L-S.Ma,A.D.Ludlow,S.M.Fore-
the present hardware, these data types have in gen-
man,J.Ye,andJ.L.Hall.Contributionofthermal
eral to be software emulated, which will slow down
noise to frequency stability of rigid optical cav-
the computation process considerably. Also, to the
ity via Hertz-linewidth lasers, Phys. Rev. A 73,
knowledge of the authors, a publicly available FE-
031804(R) (2006).
solver using quad or arbitrary precision arithmetics
doesnotexist. [8] S.A. Webster, M. Oxborrow, and P.Gill. Vibra-
Thetestcase presented hereisprocessed taking into tion insensitive optical cavity, Phys.Rev. A75,
account only thermo-elastic effects. For the design 011801, 2007.
of high precision cavities, this is of course not suffi-
[9] J. Alnis, A. Matveev, N. Kolachevsky, Th.
cient. However, the principle problems for the FE-
Udem, and T. W. Ha¨nsch. Subhertz linewidth
modelling will remain the same. Therefore, cavi-
tieswithlengthstabilities onthe10−17 leveloreven diode lasers by stabilization to vibrationally and
thermally compensated ultralow-expansion glass
higher can only be simulated if the problem of nu-
Fabry-Pe´rot cavities, Phys. Rev. A 77, 053809
merical errors can be solved. We plan to solve this
(2008).
problem with the development of a stand-alone FE-
solverusingarbitrary precision arithmetics, whichis
[10] M. Maischak. Webpage of
basedonthesolveralreadyimplementedinMaiProgs,
the Software package Maiprogs,
together with an interface to ANSYS or other pre-
http://www.ifam.uni-hannover.de/∼maiprogs.
pocessors.
[11] E.P. Stephan. The h-p version of the boundary
element method for solving 2- and 3-dimensional
References
problems, Comp. Meth. Appl. Mech. Eng. 133
[1] Ch.Eisele,A.Yu.Nevsky,andS.Schiller.Labo- (1996)183-208.
ratoryTestoftheIsotropyofLightPropagationat
[12] E.P. Stephan. Hp-Finite Element Methoden,
the 10−17 Level, Phys. Rev. Letters 103, 090401
Vorlesungsskript 2009, Leibniz Universita¨t Han-
(2009).
nover.
[2] F. Bru¨ckner, D. Friedrich, T. Clausnitzer, M.
[13] E.P. Stephan. Coupling of Boundary Element
Britzger,O.Burmeister,K.Danzmann,E-B.Kley,
MethodsandFiniteElementMethods,Encyclope-
A. Tu¨nnermann, and R. Schnabel. Realization of
diaofComputationalMechanics,EditedbyErwin
aMonolithicHigh-ReflectivityCavityMirrorfrom
Stein,RenA˜(cid:13)c deBorstandThomasJ.R.Hughes.
a Single Silicon Crystal, Phys. Rev. Lett. 104,
Vol.1,Chapter13: Fundamentals. 2004JohnWi-
163903(2010).
ley&Sons.
[3] T. Rosenband et al. Frequency Ratio ofAl+ and
[14] R. Araya, A.H. Poza, and E.P. Stephan.
Hg+Single-Ion Optical Clocks; Metrology atthe
A hierarchical a posteriori error estimate for
17th Decimal Place, Science 319, 1808 – 1812
an advection-diffusion-reaction problem, Math.
(2008).
Models Methods Appl. Sci. 15 (2005), no. 7,
[4] S. Herrmann, A. Senger, K. Mo¨hle, M. Nagel, 1119-1139.
E.V. Kovalchuk, and A. Peters. Rotating optical
[15] 0.C. Zienkiewics and J.Z. Zhu. A simple error
cavity experiment testing Lorentz invariance at
estimator and adaptive procedure for practical
the 10 − 17 level, Phys. Rev. D 80, pp. 105011
engineering analysis, Int. J. Numer. Meth. Eng.
(2009).
24,337-357 (1987).
[5] B. C. Young, F. C. Cruz, W. M. Itano, and
7
[16] R. W. Fox. Temperature analysis of low-
expansion Fabry-Perot cavities, Optics Express
17,15023(2009).
