PDF
Stress–Optical Effects with Generalized Plane Strain
Introduction
The assumptions made for plane strain in the previous analysis of the waveguide structure (see the Application Library model Stress–Optical Effects in a Photonic Waveguide, the model name is stress_optical) do not hold in a situation where the silicon-silica laminate is free to expand in the z direction. Instead, it is necessary to use a generalized plane strain model that allows for free expansion in the z direction. The boundary conditions in the xy-plane already allow the structure to expand freely in all directions in the plane. When the different materials in a laminate expand with different expansion coefficients, the laminate bends. In this model, the silica-silicon laminate bends in both the x and z directions.
Note: This model requires the Wave Optics Module and the Structural Mechanics Module.
Model Definition
Generalized Plane Strain
One possible extension of the plane strain formulation is to assume that the strains are independent of the out-of -plane coordinate z, that is
Under the small strain assumption, the above equations have the following 3D solution:
where e0, e1, and e2 are constants. Thus, at the cross section z = 0, one has
and
The above conditions differ from the plane strain state only by the fact that the normal out-of -plane strain component can vary linearly throughout the cross section. This approximation is expected to be good when the structure is free to expand in the out-of-plane direction, and possible bending curvature is small with respect to the extents of the structure in the xy-plane. It corresponds to a small rotation that is representative of each cross section of the structure along the z-axis.
In COMSOL Multiphysics, the Generalized Plane Strain option is available under Structural Mechanics interface in 2D. Coefficients e0, e1, and e2 in the expression for the εz strain are modeled as extra degrees of freedom that are constant throughout the model (global variables). See the Structural Mechanics Module User’s Guide for more details.
Results and Discussion
Figure 1 shows the von Mises stress distribution together with the deformed shape of the waveguide.
Figure 1: von Mises equivalent stress computed with (upper) and without (lower) the generalized plane strain assumption.
A 3D representation of the generalized plane strain solution is shown in this figure:
For symmetry reasons, the strain components εx and εz should be equal. The plot in Figure 2 visualizes the area where the relative difference between εx and εz is within 5%. The model is most accurate in the regions close to the core, far from the boundaries on the far left and right.
Figure 2: In the red-colored region, the relative difference between the x and z strain components is within 5%.
Figure 3 compares the effective mode indices for the first four propagating modes using the generalized plane strain equations with those obtained from the analysis in the previous model. As the plot shows, there is a systematic shift in the propagation constants when the strain in the z direction is taken into account.
Figure 3: Effective mode indices assuming plane strain and generalized plane strain.
Figure 4 shows the stress-induced birefringence along the symmetry line within the waveguide. Nonzero out-of-plane strain leads to an increase in the birefringence effect.
Figure 4: Birefringence along the vertical symmetry line within the waveguide for plane strain (lower curve) and generalized plane strain (upper curve).
Figure 5 gives the details of the eigenmode with the lowest effective mode index. It presents the visualization of the power flow, also called the optical intensity or the Poynting vector, in the out-of-plane direction.
Figure 5: Eigenmode with lowest mode index, computed with the stress–optical effect under generalized plane strain assumption.
Application Library path: Wave_Optics_Module/Waveguides/stress_optical_generalized
Modeling Instructions
From the File menu, choose New.
New
In the New window, click  Model Wizard.
Model Wizard
1
In the Model Wizard window, click  2D.
2
In the Select Physics tree, select Structural Mechanics > Solid Mechanics (solid).
3
Click Add.
4
In the Select Physics tree, select Optics > Wave Optics > Electromagnetic Waves, Frequency Domain (ewfd).
5
Click Add.
6
Click  Study.
7
In the Select Study tree, select Preset Studies for Some Physics Interfaces > Stationary.
8
Click  Done.
Global Definitions
Parameters 1
1
In the Model Builder window, under Global Definitions click Parameters 1.
2
In the Settings window for Parameters, locate the Parameters section.
3
Click  Load from File.
4
Browse to the model’s Application Libraries folder and double-click the file stress_optical_parameters.txt.
Geometry 1
1
In the Model Builder window, under Component 1 (comp1) click Geometry 1.
2
In the Settings window for Geometry, locate the Units section.
3
From the Length unit list, choose µm.
Rectangle 1 (r1)
1
In the Geometry toolbar, click  Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type 320.
4
In the Height text field, type 83.
5
Locate the Position section. In the x text field, type -160.
6
In the y text field, type -100.
Rectangle 2 (r2)
1
In the Geometry toolbar, click  Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type 320.
4
In the Height text field, type 14.
5
Locate the Position section. In the x text field, type -160.
6
In the y text field, type -17.
Rectangle 3 (r3)
1
In the Geometry toolbar, click  Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type 320.
4
In the Height text field, type 16.
5
Locate the Position section. In the x text field, type -160.
6
In the y text field, type -3.
Rectangle 4 (r4)
1
In the Geometry toolbar, click  Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type 6.
4
In the Height text field, type 6.
5
Locate the Position section. In the x text field, type -3.
6
In the y text field, type -3.
Rectangle 5 (r5)
1
In the Geometry toolbar, click  Rectangle.
2
In the Settings window for Rectangle, locate the Size and Shape section.
3
In the Width text field, type 20.
4
In the Height text field, type 20.
5
Locate the Position section. In the x text field, type -10.
6
In the y text field, type -10.
7
Click  Build Selected.
The last rectangular region encloses the optical computational domain. It can be enlarged if needed for validating the results. The region should be chosen large enough so that the computed propagation constants do not change significantly if the region is enlarged.
Solid Mechanics (solid)
1
In the Model Builder window, under Component 1 (comp1) click Solid Mechanics (solid).
2
In the Settings window for Solid Mechanics, locate the Thickness section.
3
In the d text field, type 2*d.
Linear Elastic Material 1
In the Model Builder window, under Component 1 (comp1) > Solid Mechanics (solid) click Linear Elastic Material 1.
Thermal Expansion 1
1
In the Physics toolbar, click  Attributes and choose Thermal Expansion.
2
In the Settings window for Thermal Expansion, locate the Model Input section.
3
From the T list, choose User defined. In the associated text field, type T1.
4
Click  Go to Source for Volume reference temperature.
Global Definitions
Default Model Inputs
1
In the Model Builder window, under Global Definitions click Default Model Inputs.
2
In the Settings window for Default Model Inputs, locate the Browse Model Inputs section.
3
Find the Expression for remaining selection subsection. In the Volume reference temperature text field, type T0.
Definitions
Variables 1
1
In the Definitions toolbar, click  Local Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
N
nCore
 
Refractive index for stress-free material
4
Locate the Geometric Entity Selection section. From the Geometric entity level list, choose Domain.
5
Select Domain 6 only.
Variables 2
1
In the Definitions toolbar, click  Local Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
N
nSiO2
 
Refractive index for stress-free material
4
Locate the Geometric Entity Selection section. From the Geometric entity level list, choose Domain.
5
Select Domains 4 and 5 only.
Variables 3
1
In the Definitions toolbar, click  Local Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
Nx
N-(B1*solid.sx+B2*(solid.sy+solid.sz))
 
Refractive index, x-component
Ny
N-(B1*solid.sy+B2*(solid.sx+solid.sz))
 
Refractive index, y-component
Nz
N-(B1*solid.sz+B2*(solid.sx+solid.sy))
 
Refractive index, z-component
4
Locate the Geometric Entity Selection section. From the Geometric entity level list, choose Domain.
5
Select Domains 4–6 only.
Solid Mechanics (solid)
Rigid Motion Suppression 1
1
In the Physics toolbar, click  Domains and choose Rigid Motion Suppression.
2
Select Domain 1 only.
Electromagnetic Waves, Frequency Domain (ewfd)
The computational domain is reduced significantly for the optical mode analysis.
1
In the Model Builder window, under Component 1 (comp1) click Electromagnetic Waves, Frequency Domain (ewfd).
2
Select Domains 4–6 only.
Wave Equation, Electric 1
1
In the Model Builder window, under Component 1 (comp1) > Electromagnetic Waves, Frequency Domain (ewfd) click Wave Equation, Electric 1.
2
In the Settings window for Wave Equation, Electric, locate the Electric Displacement Field section.
3
From the n list, choose User defined. From the list, choose Diagonal.
4
Specify the n matrix as
Nx
0
0
0
Ny
0
0
0
Nz
5
From the k list, choose User defined.
Materials
Si
1
In the Model Builder window, under Component 1 (comp1) right-click Materials and choose Blank Material.
2
In the Settings window for Material, type Si in the Label text field.
3
Select Domain 1 only.
4
Locate the Material Contents section. In the table, enter the following settings:
Property
Variable
Value
Unit
Property group
Young’s modulus
E
ESi
Pa
Young’s modulus and Poisson’s ratio
Poisson’s ratio
nu
nuSi
1
Young’s modulus and Poisson’s ratio
Density
rho
rhoSi
kg/m³
Basic
Coefficient of thermal expansion
alpha_iso ; alphaii = alpha_iso, alphaij = 0
alphaSi
1/K
Basic
SiO2
1
Right-click Materials and choose Blank Material.
2
In the Settings window for Material, type SiO2 in the Label text field.
3
Select Domains 2–6 only.
4
Locate the Material Contents section. In the table, enter the following settings:
Property
Variable
Value
Unit
Property group
Young’s modulus
E
ESiO2
Pa
Young’s modulus and Poisson’s ratio
Poisson’s ratio
nu
nuSiO2
1
Young’s modulus and Poisson’s ratio
Density
rho
rhoSiO2
kg/m³
Basic
Coefficient of thermal expansion
alpha_iso ; alphaii = alpha_iso, alphaij = 0
alphaSiO2
1/K
Basic
Mesh 1
Free Triangular 1
In the Mesh toolbar, click  Free Triangular.
Size 1
1
Right-click Free Triangular 1 and choose Size.
2
In the Settings window for Size, locate the Geometric Entity Selection section.
3
From the Geometric entity level list, choose Domain.
4
Select Domains 4–6 only.
5
Locate the Element Size section. Click the Custom button.
6
Locate the Element Size Parameters section.
7
Select the Maximum element size checkbox. In the associated text field, type 0.2.
Size
1
In the Model Builder window, under Component 1 (comp1) > Mesh 1 click Size.
2
In the Settings window for Size, locate the Element Size section.
3
From the Predefined list, choose Extra fine.
4
Click  Build All.
Study 1
Step 2: Mode Analysis
1
In the Study toolbar, click  Study Steps and choose Other > Mode Analysis.
2
In the Settings window for Mode Analysis, locate the Study Settings section.
3
From the Transform list, choose Effective mode index.
4
In the Search for modes around shift text field, type 1.46.
5
Select the Desired number of modes checkbox. In the associated text field, type 4.
6
In the Mode analysis frequency text field, type c_const/lambda0_ewfd.
These settings make the eigenmode solver search for the 4 eigenmodes with effective mode indices closest to the value 1.46. This value is an estimate of the effective mode index for the fundamental mode.
Exclude Solid Mechanics from the Mode Analysis step.
7
Locate the Physics and Variables Selection section. In the Solve for column of the table, under Component 1 (comp1), clear the checkbox for Solid Mechanics (solid).
Step 1: Stationary
In the Study toolbar, click  Compute.
Results
Stress (solid)
1
In the Settings window for 2D Plot Group, locate the Data section.
2
From the Dataset list, choose Study 1/Solution Store 1 (sol2).
3
In the Stress (solid) toolbar, click  Plot.
4
Click the  Go to Default View button in the Graphics toolbar.
5
Click the  Zoom Extents button in the Graphics toolbar.
6
In the Study toolbar, click  Create Solution Copy.
Study 1/Solution 1 - Plane Strain
1
In the Model Builder window, expand the Results > Datasets node, then click Study 1/Solution 1 - Copy 1 (sol3).
2
In the Settings window for Solution, type Study 1/Solution 1 - Plane Strain in the Label text field.
Solid Mechanics (solid)
1
In the Model Builder window, under Component 1 (comp1) click Solid Mechanics (solid).
2
In the Settings window for Solid Mechanics, locate the 2D Approximation section.
3
From the list, choose Generalized plane strain.
Definitions
Variables 4
1
In the Model Builder window, under Component 1 (comp1) right-click Definitions and choose Variables.
2
In the Settings window for Variables, locate the Variables section.
3
In the table, enter the following settings:
Name
Expression
Unit
Description
u1
u-solid.wxZ/2*Z^2
m
Displacement, x-component
v1
v-solid.wyZ/2*Z^2
m
Displacement, y-component
w1
(solid.wxZ*X+solid.wyZ*Y+solid.w0Z)*Z
m
Displacement, z-component
Study 1
In the Study toolbar, click  Compute.
Results
Stress (solid)
1
In the Model Builder window, under Results click Stress (solid).
2
In the Stress (solid) toolbar, click  Plot.
3
Click the  Go to Default View button in the Graphics toolbar.
4
Click the  Zoom Extents button in the Graphics toolbar.
Strain level
To visualize the area, where the relative difference between the x and z strain components is within 5%, follow these steps:
1
In the Model Builder window, under Results click Electric Field (ewfd).
2
In the Settings window for 2D Plot Group, type Strain level in the Label text field.
3
Locate the Data section. From the Dataset list, choose Study 1/Solution Store 1 (sol2).
Surface 1
1
In the Model Builder window, expand the Strain level node, then click Surface 1.
2
In the Settings window for Surface, locate the Expression section.
3
In the Expression text field, type abs((solid.eX-solid.eZ)/solid.eZ)<0.05.
4
Locate the Coloring and Style section. Clear the Color legend checkbox.
5
Click to expand the Quality section. From the Resolution list, choose Extra fine.
6
In the Strain level toolbar, click  Plot.
7
Click the  Zoom Extents button in the Graphics toolbar.
Cut Line 2D 1
Next, plot the stress-induced birefringence along the symmetry line within the waveguide.
1
In the Results toolbar, click  Cut Line 2D.
2
In the Settings window for Cut Line 2D, locate the Line Data section.
3
In row Point 1, set Y to -10.
4
In row Point 2, set X to 0 and Y to 10.
5
Click  Plot.
Cut Line 2D 2
1
Right-click Cut Line 2D 1 and choose Duplicate.
2
In the Settings window for Cut Line 2D, locate the Data section.
3
From the Dataset list, choose Study 1/Solution 1 - Plane Strain (sol3).
1D Plot Group 3
In the Results toolbar, click  1D Plot Group.
Line Graph 1
1
Right-click 1D Plot Group 3 and choose Line Graph.
2
In the Settings window for Line Graph, locate the Data section.
3
From the Dataset list, choose Cut Line 2D 1.
4
Locate the y-Axis Data section. In the Expression text field, type Nx-Ny.
5
Locate the x-Axis Data section. From the Parameter list, choose Expression.
6
In the Expression text field, type Y.
Line Graph 2
1
Right-click Line Graph 1 and choose Duplicate.
2
In the Settings window for Line Graph, locate the Data section.
3
From the Dataset list, choose Cut Line 2D 2.
1D Plot Group 3
1
In the Model Builder window, click 1D Plot Group 3.
2
In the Settings window for 1D Plot Group, click to expand the Title section.
3
From the Title type list, choose Manual.
4
In the Title text area, type Birefringence.
5
Locate the Plot Settings section.
6
Select the x-axis label checkbox. In the associated text field, type Y.
7
Select the y-axis label checkbox. In the associated text field, type Nx-Ny.
8
In the 1D Plot Group 3 toolbar, click  Plot.
To visualize the details of the eigenmode with the lowest effective mode index, you first set up a view that includes the optical computation domain only.
9
Click the  Show More Options button in the Model Builder toolbar.
10
In the Show More Options dialog, in the tree, select the checkbox for the node Results > Views.
11
Click OK.
View 2D 2
1
In the Model Builder window, under Results right-click Views and choose View 2D.
2
In the Settings window for View 2D, locate the View section.
3
Select the Lock axis checkbox.
Axis
1
In the Model Builder window, expand the View 2D 2 node, then click Axis.
2
In the Settings window for Axis, locate the Axis section.
3
In the x minimum text field, type -10.
4
In the x maximum text field, type 10.
5
In the y minimum text field, type -10.
6
In the y maximum text field, type 10.
2D Plot Group 4
1
In the Results toolbar, click  2D Plot Group.
2
In the Settings window for 2D Plot Group, locate the Data section.
3
From the Effective mode index (1) list, choose 1.4514.
4
Locate the Plot Settings section. From the View list, choose View 2D 2.
Surface 1
1
Right-click 2D Plot Group 4 and choose Surface.
2
In the Settings window for Surface, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1) > Electromagnetic Waves, Frequency Domain > Energy and power > Power flow, time average (spatial frame) - W/m² > ewfd.Poavz - Power flow, time average, z-component.
3
In the 2D Plot Group 4 toolbar, click  Plot. This creates a visualization of the power flow, also called optical intensity or the Poynting vector, in the z direction (out-of-plane direction).
Global Evaluation 1
To collect all computed effective mode indices in a table, follow these steps:
1
In the Results toolbar, click  Global Evaluation.
2
In the Settings window for Global Evaluation, click Replace Expression in the upper-right corner of the Expressions section. From the menu, choose Component 1 (comp1) > Electromagnetic Waves, Frequency Domain > Global > ewfd.neff - Effective mode index - 1.
3
Locate the Expressions section. In the table, enter the following settings:
Expression
Unit
Description
ewfd.neff
1
Effective mode index (Generalized plane strain)
4
Click  Evaluate.
Table 1
1
Go to the Table 1 window.
If you see too few digits in the table, click the Full Precision toolbar button.
Global Evaluation 2
1
Right-click Global Evaluation 1 and choose Duplicate.
2
In the Settings window for Global Evaluation, locate the Data section.
3
From the Dataset list, choose Study 1/Solution 1 - Plane Strain (sol3).
4
Locate the Expressions section. In the table, enter the following settings:
Expression
Unit
Description
ewfd.neff
1
Effective mode index (Plane strain)
5
Click  Evaluate.
Create a table plot to visualize the shift of the effective mode indices.
1D Plot Group 5
In the Results toolbar, click  1D Plot Group.
Table Graph 1
1
Right-click 1D Plot Group 5 and choose Table Graph.
2
In the Settings window for Table Graph, locate the Data section.
3
From the Plot columns list, choose Manual.
4
From the x-axis data list, choose Row index.
5
In the Columns list, choose Effective mode index (Generalized plane strain) (1) and Effective mode index (Plane strain) (1).
6
Locate the Coloring and Style section. Find the Line style subsection. From the Line list, choose None.
7
Find the Line markers subsection. From the Marker list, choose Circle.
8
Click to expand the Legends section. Select the Show legends checkbox.
1D Plot Group 5
1
In the Model Builder window, click 1D Plot Group 5.
2
In the Settings window for 1D Plot Group, locate the Title section.
3
From the Title type list, choose Manual.
4
In the Title text area, type Eigenvalue shift.
5
Locate the Plot Settings section. Select the x-axis label checkbox.
6
Select the y-axis label checkbox.
7
In the x-axis label text field, type Mode number.
8
In the y-axis label text field, type Effective mode index.
9
Locate the Legend section. From the Layout list, choose Outside graph axis area.
10
From the Position list, choose Bottom.
11
In the Number of rows text field, type 2.
12
In the 1D Plot Group 5 toolbar, click  Plot.
Finally, visualize the out-of-plane warping of the waveguide.
Extrusion 2D 1
1
In the Results toolbar, click  More Datasets and choose Extrusion 2D.
2
In the Settings window for Extrusion 2D, locate the Data section.
3
From the Dataset list, choose Study 1/Solution Store 1 (sol2).
4
Locate the Extrusion section. In the z maximum text field, type d.
5
In the z variable text field, type Z.
3D Plot Group 6
1
In the Results toolbar, click  3D Plot Group.
2
In the Settings window for 3D Plot Group, click to expand the Number Format section.
3
Select the Manual grid settings checkbox.
4
In the Precision text field, type 4.
Surface 1
1
Right-click 3D Plot Group 6 and choose Surface.
2
In the Settings window for Surface, locate the Expression section.
3
In the Expression text field, type v1.
Deformation 1
1
Right-click Surface 1 and choose Deformation.
2
In the Settings window for Deformation, locate the Expression section.
3
In the x-component text field, type u1.
4
In the y-component text field, type v1.
5
In the z-component text field, type w1.
6
Locate the Scale section.
7
Select the Scale factor checkbox. In the associated text field, type 1.
8
In the 3D Plot Group 6 toolbar, click  Plot.
9
Click the  Go to Default View button in the Graphics toolbar.
10
Click the  Zoom Extents button in the Graphics toolbar.