SHARE
TWEET

Untitled

a guest Jun 15th, 2019 63 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. << NDSolve`FEM`
  2.  
  3. r = 0.8;
  4. outerCirclePoints =
  5.   With[{r = 2.},
  6.    Table[{r Cos[[Theta]], r Sin[[Theta]]}, {[Theta],
  7.      Range[0, 2 [Pi], 0.05 [Pi]] // Most}]];
  8. innerCirclePoints =
  9.   With[{r = r},
  10.    Table[{r Cos[[Theta]], r Sin[[Theta]]}, {[Theta],
  11.      Range[0, 2 [Pi], 0.08 [Pi]] // Most}]];
  12.  
  13. bmesh = ToBoundaryMesh[
  14.    "Coordinates" -> Join[outerCirclePoints, innerCirclePoints],
  15.    "BoundaryElements" -> {LineElement[
  16.       Riffle[Range[Length@outerCirclePoints],
  17.         RotateLeft[Range[Length@outerCirclePoints], 1]] //
  18.        Partition[#, 2] &],
  19.      LineElement[
  20.       Riffle[Range[Length@outerCirclePoints + 1,
  21.          Length@Join[outerCirclePoints, innerCirclePoints]],
  22.         RotateLeft[
  23.          Range[Length@outerCirclePoints + 1,
  24.           Length@Join[outerCirclePoints, innerCirclePoints]], 1]] //
  25.        Partition[#, 2] &]}];
  26.  
  27. mesh = ToElementMesh[bmesh];
  28. {bmesh["Wireframe"], mesh["Wireframe"]}
  29.  
  30. glass = 1.45; air = 1.; k0 = (2 [Pi])/1.55;
  31. n[x_, y_] := If[x^2 + y^2 <= r^2, glass, air]
  32.  
  33. helm = !(
  34. *SubsuperscriptBox[([Del]), ({x, y}), (2)](u[x, y])) +
  35.    n[x, y]^2*k0^2*u[x, y];
  36. boundary = DirichletCondition[u[x, y] == 0., True];
  37.  
  38. (*region=ImplicitRegion[x^2+y^2[LessEqual]2.^2,{x,y}];*)
  39.  
  40. {vals, funs} =
  41.   NDEigensystem[{helm, boundary}, u[x, y], {x, y} [Element] mesh, 1,
  42.    Method -> {"Eigensystem" -> {"FEAST",
  43.        "Interval" -> {k0^2, glass^2 k0^2}}}];
  44. vals
  45.  
  46.  Table[Plot3D[funs[[i]], {x, y} [Element] mesh, PlotRange -> All,
  47.   PlotLabel -> vals[[i]]], {i, Length[vals]}]
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy. OK, I Understand
 
Top