Mode solving#
Caution
This example is under construction As Julia-Dicts are not ordered, the mesh might become incorrect when adjusted (for now, better do the meshing in python)
using PyCall
np = pyimport("numpy")
shapely = pyimport("shapely")
shapely.affinity = pyimport("shapely.affinity")
clip_by_rect = pyimport("shapely.ops").clip_by_rect
OrderedDict = pyimport("collections").OrderedDict
mesh_from_OrderedDict = pyimport("femwell.mesh").mesh_from_OrderedDict
wg_width = 0.5
wg_thickness = 0.22
slab_width = 3
slab_thickness = 0.11
core = shapely.geometry.box(-wg_width / 2, 0, wg_width / 2, wg_thickness)
slab = shapely.geometry.box(slab_width / 2, 0, slab_width / 2, slab_thickness)
env = shapely.affinity.scale(core.buffer(3, resolution = 8), xfact = 0.5)
polygons = OrderedDict(
core = core,
slab = slab,
box = clip_by_rect(env, -np.inf, -np.inf, np.inf, 0),
clad = clip_by_rect(env, -np.inf, 0, np.inf, np.inf),
)
resolutions = Dict(
"core" => Dict("resolution" => 0.03, "distance" => 0.5),
"slab" => Dict("resolution" => 0.03, "distance" => 0.5),
)
mesh = mesh_from_OrderedDict(
polygons,
resolutions,
default_resolution_max = 1,
filename = "mesh.msh",
)
Show code cell output
PyObject <meshio mesh object>
Number of points: 338
Number of cells:
vertex: 1
vertex: 1
line: 4
line: 17
line: 8
line: 17
line: 8
line: 6
line: 8
line: 3
triangle: 342
triangle: 116
triangle: 178
Cell sets: slab_points, slab, core___box, core___clad, box___clad, core, box, clad, gmsh:bounding_entities
Point data: gmsh:dim_tags
Cell data: gmsh:physical, gmsh:geometrical
Field data: slab_points, slab, core___box, core___clad, box___clad, core, box, clad
using Gridap
using Gridap.Geometry
using Gridap.Visualization
using Gridap.ReferenceFEs
using GridapGmsh
using GridapMakie, CairoMakie
using Femwell.Maxwell.Waveguide
CairoMakie.inline!(true)
Show code cell output
The following 1 direct dependency failed to precompile:
GLMakie [e9467ef8-e4e7-5192-8a1a-b1aee30e663a]
Failed to precompile GLMakie [e9467ef8-e4e7-5192-8a1a-b1aee30e663a] to "/home/runner/.julia/compiled/v1.10/GLMakie/jl_ugYon7".
┌ Warning: OpenGL/GLFW wasn't loaded correctly or couldn't be initialized.
│ This likely means, you're on a headless server without having OpenGL support setup correctly.
│ Have a look at the troubleshooting section in the readme:
│ https://github.com/MakieOrg/Makie.jl/tree/master/GLMakie#troubleshooting-opengl.
└ @ GLMakie ~/.julia/packages/GLMakie/TH3rf/src/gl_backend.jl:4
ERROR: LoadError: InitError: Exception[GLFW.GLFWError(65550, "X11: The DISPLAY environment variable is missing"), ErrorException("glfwInit failed")]
Stacktrace:
[1] __init__()
@ GLFW ~/.julia/packages/GLFW/wmoTL/src/GLFW.jl:35
[2] run_module_init(mod::Module, i::Int64)
@ Base ./loading.jl:1134
[3] register_restored_modules(sv::Core.SimpleVector, pkg::Base.PkgId, path::String)
@ Base ./loading.jl:1122
[4] _include_from_serialized(pkg::Base.PkgId, path::String, ocachepath::String, depmods::Vector{Any})
@ Base ./loading.jl:1067
[5] _require_search_from_serialized(pkg::Base.PkgId, sourcepath::String, build_id::UInt128)
@ Base ./loading.jl:1581
[6] _require(pkg::Base.PkgId, env::String)
@ Base ./loading.jl:1938
[7] __require_prelocked(uuidkey::Base.PkgId, env::String)
@ Base ./loading.jl:1812
[8] #invoke_in_world#3
@ ./essentials.jl:926 [inlined]
[9] invoke_in_world
@ ./essentials.jl:923 [inlined]
[10] _require_prelocked(uuidkey::Base.PkgId, env::String)
@ Base ./loading.jl:1803
[11] macro expansion
@ ./loading.jl:1790 [inlined]
[12] macro expansion
@ ./lock.jl:267 [inlined]
[13] __require(into::Module, mod::Symbol)
@ Base ./loading.jl:1753
[14] #invoke_in_world#3
@ ./essentials.jl:926 [inlined]
[15] invoke_in_world
@ ./essentials.jl:923 [inlined]
[16] require(into::Module, mod::Symbol)
@ Base ./loading.jl:1746
[17] top-level scope
@ ~/.julia/packages/GLMakie/TH3rf/src/gl_backend.jl:2
[18] include(mod::Module, _path::String)
@ Base ./Base.jl:495
[19] include(x::String)
@ GLMakie ~/.julia/packages/GLMakie/TH3rf/src/GLMakie.jl:1
[20] top-level scope
@ ~/.julia/packages/GLMakie/TH3rf/src/GLMakie.jl:87
[21] include
@ ./Base.jl:495 [inlined]
[22] include_package_for_output(pkg::Base.PkgId, input::String, depot_path::Vector{String}, dl_load_path::Vector{String}, load_path::Vector{String}, concrete_deps::Vector{Pair{Base.PkgId, UInt128}}, source::Nothing)
@ Base ./loading.jl:2222
[23] top-level scope
@ stdin:3
during initialization of module GLFW
in expression starting at /home/runner/.julia/packages/GLMakie/TH3rf/src/gl_backend.jl:1
in expression starting at /home/runner/.julia/packages/GLMakie/TH3rf/src/GLMakie.jl:1
in expression starting at stdin:
Stacktrace:
[1] pkgerror(msg::String)
@ Pkg.Types /opt/hostedtoolcache/julia/1.10.3/x64/share/julia/stdlib/v1.10/Pkg/src/Types.jl:70
[2] precompile(ctx::Pkg.Types.Context, pkgs::Vector{Pkg.Types.PackageSpec}; internal_call::Bool, strict::Bool, warn_loaded::Bool, already_instantiated::Bool, timing::Bool, _from_loading::Bool, kwargs::@Kwargs{io::IJulia.IJuliaStdio{Base.PipeEndpoint}})
@ Pkg.API /opt/hostedtoolcache/julia/1.10.3/x64/share/julia/stdlib/v1.10/Pkg/src/API.jl:1659
[3] precompile(pkgs::Vector{Pkg.Types.PackageSpec}; io::IJulia.IJuliaStdio{Base.PipeEndpoint}, kwargs::@Kwargs{_from_loading::Bool})
@ Pkg.API /opt/hostedtoolcache/julia/1.10.3/x64/share/julia/stdlib/v1.10/Pkg/src/API.jl:159
[4] precompile
@ /opt/hostedtoolcache/julia/1.10.3/x64/share/julia/stdlib/v1.10/Pkg/src/API.jl:147 [inlined]
[5] #precompile#114
@ /opt/hostedtoolcache/julia/1.10.3/x64/share/julia/stdlib/v1.10/Pkg/src/API.jl:146 [inlined]
[6] #invokelatest#2
@ ./essentials.jl:894 [inlined]
[7] invokelatest
@ ./essentials.jl:889 [inlined]
[8] _require(pkg::Base.PkgId, env::String)
@ Base ./loading.jl:1963
[9] __require_prelocked(uuidkey::Base.PkgId, env::String)
@ Base ./loading.jl:1812
[10] #invoke_in_world#3
@ ./essentials.jl:926 [inlined]
[11] invoke_in_world
@ ./essentials.jl:923 [inlined]
[12] _require_prelocked(uuidkey::Base.PkgId, env::String)
@ Base ./loading.jl:1803
[13] macro expansion
@ ./loading.jl:1790 [inlined]
[14] macro expansion
@ ./lock.jl:267 [inlined]
[15] __require(into::Module, mod::Symbol)
@ Base ./loading.jl:1753
[16] #invoke_in_world#3
@ ./essentials.jl:926 [inlined]
[17] invoke_in_world
@ ./essentials.jl:923 [inlined]
[18] require(into::Module, mod::Symbol)
@ Base ./loading.jl:1746
model = GmshDiscreteModel("mesh.msh")
Ω = Triangulation(model)
#fig = plot(Ω)
#fig.axis.aspect=DataAspect()
#wireframe!(Ω, color=:black, linewidth=1)
#display(fig)
labels = get_face_labeling(model)
epsilons = ["core" => 3.5^2, "slab" => 1.44^2, "box" => 1.444^2, "clad" => 1.44^2]
ε(tag) = Dict(get_tag_from_name(labels, u) => v for (u, v) in epsilons)[tag]
#dΩ = Measure(Ω, 1)
τ = CellField(get_face_tag(labels, num_cell_dims(model)), Ω)
modes = calculate_modes(model, ε ∘ τ, λ = 1.55, num = 1, order = 1)
println(n_eff(modes[1]))
# write_mode_to_vtk("mode", modes[2])
plot_mode(modes[1])
#plot_mode(modes[2])
modes
Perturbations#
Let’s add a minor perturbation and calculate the effective refractive index using perturbation theory:
epsilons_p = ["core" => -0.1im, "box" => -0.0im, "slab" => -0.0im, "clad" => -0.0im]
ε_p(tag) = Dict(get_tag_from_name(labels, u) => v for (u, v) in epsilons_p)[tag]
println(perturbed_neff(modes[1], ε_p ∘ τ))
For comparison, we also calculate directly the effective refractive index:
modes_p = calculate_modes(model, ε ∘ τ + ε_p ∘ τ, λ = 1.55, num = 2, order = 1)
println(n_eff(modes_p[1]))
plot_mode(modes_p[1])