How I can solve the Helmholtz equation in 3D domain? According to the following post, I encounter an error when I run the program. Any help would be graet!


This section should be considered as a continuous program to avoid any possible errors

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[, 4, MaxCellMeasure -> .05, AccuracyGoal -> 2];
  • $\begingroup$ "ImplicitRegion[surf, {{x, -1, 1}, {y, -1, 1}, {z, -1, 1}}]" incorrectly defined. "surf" should be a condition. $\endgroup$ – Daniel Huber Sep 16 at 20:25
  • $\begingroup$ Thanks. How I can edit them.? Best $\endgroup$ – iman karimipour Sep 17 at 9:55
  • 1
    $\begingroup$ "surf" is a vector function of u and v. However, Implicite Region needs a condition, e.g.ImplicitRegion[x^2 + y^2 <= 1, {x, y}] $\endgroup$ – Daniel Huber Sep 17 at 10:57
  • $\begingroup$ For tours structure how I can define a condition such as you mentioned. The Coordinate in this structure is defenid as follows: en.wikipedia.org/wiki/Toroidal_coordinates $\endgroup$ – iman karimipour Sep 17 at 18:35
  • $\begingroup$ Can you please specify explicitly the error that you get? Thank you! $\endgroup$ – CA Trevillian Sep 17 at 19:28

I assume "tours structure" means "torus surface". For this I would use a parameterized torus in cart. coord. Toroidal coordinates make this unnecessary complicated. E.g.

r1 = 1; r2 = 2;
reg = ParametricRegion[{(r2 + r1 Cos[phi]) Cos[
      tau], (r2 + r1 Cos[phi]) Sin[tau], 
    r1 Sin[phi]}, {{phi, -Pi, Pi}, {tau, -Pi, Pi}}];

enter image description here

| improve this answer | |
  • $\begingroup$ Very nice! Does this solve the OP’s issue? Is it possible to show this? Thank you :) $\endgroup$ – CA Trevillian Sep 17 at 19:29
  • $\begingroup$ i Provided the updated code in the main question. Errors are: ToElementMesh::femtemnm: A mesh could not be generated. and .... $\endgroup$ – iman karimipour Sep 17 at 20:40
  • 1
    $\begingroup$ There are syntax errors: "Module[{u, x, y, z, t, pde, dirichletCondition,..." has no closing bracket. "{eigenValues, evIF, mesh}];" has no opening bracket, where does it belong? $\endgroup$ – Daniel Huber Sep 18 at 8:22
  • $\begingroup$ Please, check it again. I uploaded it again in main question. $\endgroup$ – iman karimipour Sep 18 at 18:07
  • $\begingroup$ No answer found. I do know What is the reason? This code has a solution when for example in helmholzSolve3D I replace [RR] with Cuboid[]!! $\endgroup$ – iman karimipour Sep 18 at 18:31

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.