Styjun

I am trying to compute the eigenfunctions of an oblate spheroid (a=75 cm and b=60 cm) using Mathematica's FEM package and Chris' answer from here. Specifically, I am looking for eigenfrequencies around 433, 893, 913 and 2400 MGHz. Is there any way I could narrow my search besides getting all eigenfrequencies initially and then looking for the desired outcome which is impractical?

Here is my code for the first 4 eigenmodes:

Needs["NDSolve`FEM`"];

helmholzSolve3D[g_, numEigenToCompute_Integer, 
opts : OptionsPattern[]] := 
Module[u, x, y, z, t, pde, dirichletCondition, mesh, boundaryMesh, 
nr, state, femdata, initBCs, methodData, initCoeffs, vd, sd, 
discretePDE, discreteBCs, load, stiffness, damping, pos, nDiri, 
numEigen, res, eigenValues, eigenVectors, 
evIF,

(*Discretize the region*)

If[Head[g] === ImplicitRegion || Head[g] === ParametricRegion, 
mesh = ToElementMesh[DiscretizeRegion[g, opts], opts], 
mesh = ToElementMesh[DiscretizeGraphics[g, opts], opts]];
boundaryMesh = ToBoundaryMesh[mesh];

(*Set up the PDE and boundary condition*)

pde = D[u[t, x, y, z], t] - Laplacian[u[t, x, y, z], x, y, z] + 
u[t, x, y, z] == 0;
dirichletCondition = DirichletCondition[u[t, x, y, z] == 0, True];
(*Pre-process the equations to obtain the FiniteElementData in 
StateData*)nr = ToNumericalRegion[mesh];
state = 
NDSolve`ProcessEquations[pde, dirichletCondition, 
u[0, x, y, z] == 0, u, t, 0, 1, Element[x, y, z, nr]];
femdata = state["FiniteElementData"];
initBCs = femdata["BoundaryConditionData"];
methodData = femdata["FEMMethodData"];
initCoeffs = femdata["PDECoefficientData"];

(*Set up the solution*)vd = methodData["VariableData"];

sd = NDSolve`SolutionData["Space" -> nr, "Time" -> 0.];

(*Discretize the PDE and boundary conditions*)

discretePDE = DiscretizePDE[initCoeffs, methodData, sd];
discreteBCs = DiscretizeBoundaryConditions[initBCs, methodData, sd];

(*Extract the relevant matrices and deploy the boundary conditions*)

load = discretePDE["LoadVector"];
stiffness = discretePDE["StiffnessMatrix"];
damping = discretePDE["DampingMatrix"];
DeployBoundaryConditions[load, stiffness, damping, discreteBCs];

(*Set the number of eigenvalues ignoring the Dirichlet positions*)

pos = discreteBCs["DirichletMatrix"]["NonzeroPositions"][[All, 2]];
nDiri = Length[pos];
numEigen = numEigenToCompute + nDiri;

(*Solve the eigensystem*)

res = Eigensystem[stiffness, damping, -numEigen];
res = Reverse /@ res;
eigenValues = res[[1, nDiri + 1 ;; Abs[numEigen]]];
eigenVectors = res[[2, nDiri + 1 ;; Abs[numEigen]]];
evIF = ElementMeshInterpolation[mesh, #] & /@ eigenVectors;

(*Return the relevant information*)

eigenValues, evIF, mesh]

ev, if, mesh = 
helmholzSolve3D[Ellipsoid[0, 0, 0, 0.75, 0.6, 0.6], 4, 
MaxCellMeasure -> 0.025]

Table[
DensityPlot[
if[[i]][x, y, 0.1], x, -1, 1, y, -1, 1,
RegionFunction -> Function[x, y, x^2/0.75^2 + y^2/0.6^2 < 1],
PlotLabel -> ev[i] ,
ColorFunction -> Hue,
PlotLegends -> Automatic
],
i, 1, 4
]

Any suggestions?

You could use something like this:

vals, funs = 
 NDEigensystem[-Laplacian[u[x, y, z], x, y, z] + u[x, y, z], 
 DirichletCondition[u[x, y, z] == 0, True], u, 
 Element[x, y, z, Ellipsoid[0, 0, 0, 0.75, 0.6, 0.6]], 4, 
 Method -> "Eigensystem" -> "FEAST", "Interval" -> 425, 500]

427.961, 428.783, 430.026, 430.156,...

And here are the density plots:

Table[DensityPlot[funs[[i]][x, y, 0.1], x, -1, 1, y, -1, 1, 
 RegionFunction -> Function[x, y, x^2/0.75^2 + y^2/0.6^2 < 1], 
 PlotLabel -> vals[[i]], ColorFunction -> Hue, 
 PlotLegends -> Automatic, PlotRange -> All], i, 1, 4]

Slice density plots:

Table[SliceDensityPlot3D[funs[[i]][x, y, z], 
 Element[ x, y, z, Ellipsoid[0, 0, 0, 0.75, 0.6, 0.6]], 
 PlotRange -> All, PlotLabel -> vals[[i]], 
 PlotTheme -> "Minimal"], i, Length[vals]]

And density plots:

Table[DensityPlot3D[funs[[i]][x, y, z], 
 Element[ x, y, z, Ellipsoid[0, 0, 0, 0.75, 0.6, 0.6]], 
 PlotRange -> All, PlotLabel -> vals[[i]], 
 PlotTheme -> "Minimal"], i, Length[vals]]

You may try Eigensystem with

Method -> "FEAST", "Interval" -> a, b

to search eigenvalue pairs within an interval. See the documentation of Eigensystem, Section "Methods", Subsection "FEAST" for more details.