Short answer: separation of variables still works. The snag isn’t “no self-adjointness,” it’s a separation-constant mix-up.
What’s going wrong
You forced the
same eigenvalue in zzz and ρ\rhoρ (you wrote “eigenvalues must be consistent for both axial and radial directions”). That’s not how the cylindrical Laplacian splits. In cylindrical coordinates (axisymmetric or with azimuthal index mmm), the separated equations are
Z′′(z)+kn2Z(z)⏟axial S–L = 0,R′′(ρ)+1ρR′(ρ)+(λmp−m2ρ2)R(ρ)⏟radial S–L with weight ρ = 0,\underbrace{Z''(z)+k_n^2 Z(z)}_{\text{axial S–L}}\;=\;0,\qquad\underbrace{R''(\rho)+\frac1\rho R'(\rho)+\Big(\lambda_{mp}-\frac{m^2}{\rho^2}\Big)R(\rho)}_{\text{radial S–L with weight }\rho}\;=\;0,axial S–LZ′′(z)+kn2Z(z)=0,radial S–L with weight ρR′′(ρ)+ρ1R′(ρ)+(λmp−ρ2m2)R(ρ)=0,
with
two independent separation constants: knk_nkn from the axial BCs, and λmp=αmp2\lambda_{mp}=\alpha_{mp}^2λmp=αmp2 from the radial BCs. The full Laplacian eigenvalue is the
sum:
Λn,mp=kn2+αmp2.\Lambda_{n,mp}=k_n^2+\alpha_{mp}^2.Λn,mp=kn2+αmp2.
You do
not set kn=αmpk_n=\alpha_{mp}kn=αmp. The product basis
Φn,mp(ρ,ϕ,z)=Jm(αmpρ) eimϕ {cosknz, sinknz}\Phi_{n,mp}(\rho,\phi,z)=J_m(\alpha_{mp}\rho)\,e^{im\phi}\,\{ \cos k_n z,\ \sin k_n z\}Φn,mp(ρ,ϕ,z)=Jm(αmpρ)eimϕ{cosknz, sinknz}
is orthogonal in (ρ,z,ϕ)(\rho,z,\phi)(ρ,z,ϕ) with inner product ∫02π ∫0W ∫0RΦ Ψ ρ dρ dz dϕ\int_0^{2\pi}\!\int_0^W\!\int_0^R \Phi\,\Psi\,\rho\,d\rho\,dz\,d\phi∫02π∫0W∫0RΦΨρdρdzdϕ.
Self-adjointness & orthogonality (radial)
The radial operator is Sturm–Liouville on (0,R)(0,R)(0,R) with weight ρ\rhoρ. Self-adjointness holds if the boundary form
[ρ (fg′−f′g)]0R=0\big[\rho\,(f g' - f' g)\big]_{0}^{R}=0[ρ(fg′−f′g)]0R=0
vanishes for your BCs. That’s true for standard Dirichlet/Neumann/Robin at ρ=R\rho=Rρ=R with
constant coefficients; and at ρ=0\rho=0ρ=0 regularity removes the singular solution YmY_mYm (choose JmJ_mJm only). Hence the familiar orthogonality:
∫0Rρ Jm(αmpρ) Jm(αmqρ) dρ ∝ δpq.\int_0^R \rho\,J_m(\alpha_{mp}\rho)\,J_m(\alpha_{mq}\rho)\,d\rho\;\propto\; \delta_{pq}.∫0RρJm(αmpρ)Jm(αmqρ)dρ∝δpq.
If your domain is an annulus a≤ρ≤Ra\le\rho\le Ra≤ρ≤R, you keep a JmJ_mJm+YmY_mYm combination and enforce Robin/Dirichlet/Neumann at both aaa and RRR; the resulting transcendental equation defines αmp\alpha_{mp}αmp. It’s still a self-adjoint S–L pair.
When separation actually breaks
Separation fails only if your
boundary condition couples ρ\rhoρ and zzz (e.g., a Robin BC like ∂ρΦ+h ∂zΦ=0\partial_\rho \Phi + h\,\partial_z\Phi=0∂ρΦ+h∂zΦ=0 on the sidewall). Then the axial and radial problems won’t decouple. In that case:
- Use a Fourier(-series) in zzz first: Φ(ρ,z)=∑nΦ^n(ρ) eiknz\Phi(\rho,z)=\sum_n \widehat\Phi_n(\rho)\,e^{ik_n z}Φ(ρ,z)=∑nΦn(ρ)eiknz. Each mode sees a radial ODE with a parameter knk_nkn inside the Robin coefficients. Solve each mode (Bessel with parameter), then invert the Fourier series/transform. Orthogonality in ρ\rhoρ is no longer across different knk_nkn, but that’s fine—the orthogonality sits in zzz (i.e., δnn′\delta_{nn'}δnn′).
- Or build a Green’s function: G(ρ,z;ρ′,z′)=∑n,mpΦn,mp(ρ,z) Φn,mp(ρ′,z′)Λn,mpG(\rho,z;\rho',z')=\sum_{n,mp}\frac{\Phi_{n,mp}(\rho,z)\,\Phi_{n,mp}(\rho',z')}{\Lambda_{n,mp}}G(ρ,z;ρ′,z′)=∑n,mpΛn,mpΦn,mp(ρ,z)Φn,mp(ρ′,z′) if BCs decouple; with coupled BCs, use Fourier in zzz and a mode-dependent radial Green’s function Gn(ρ,ρ′)G_n(\rho,\rho')Gn(ρ,ρ′).
Practical recipe (what to do)
- Axial problem: impose your axial BCs to get knk_nkn
• Periodic of span WWW: kn=2πn/Wk_n=2\pi n/Wkn=2πn/W
• Dirichlet/Neumann at z=0,Wz=0,Wz=0,W: sines/cosines with kn=nπ/Wk_n=n\pi/Wkn=nπ/W.
- Radial problem: impose constant-coefficient BC at ρ=R\rho=Rρ=R (and regularity at ρ=0\rho=0ρ=0) to get αmp\alpha_{mp}αmp from
a R(ρ)+b ρR′(ρ)=0at ρ=R ⇒ a Jm(αR)+b αR Jm′(αR)=0.a\,R(\rho)+b\,\rho R'(\rho)=0\quad\text{at }\rho=R\;\Rightarrow\; a\,J_m(\alpha R)+b\,\alpha R\,J_m'(\alpha R)=0.aR(ρ)+bρR′(ρ)=0at ρ=R⇒aJm(αR)+bαRJm′(αR)=0.
- Expand your data in the tensor basis; coefficients are
cn,mp=∫f(ρ,ϕ,z) Φn,mp(ρ,ϕ,z) ρ dρ dϕ dz∫Φn,mp2 ρ dρ dϕ dz,c_{n,mp}=\frac{\displaystyle \int f(\rho,\phi,z)\,\Phi_{n,mp}(\rho,\phi,z)\,\rho\,d\rho\,d\phi\,dz}{\displaystyle \int \Phi_{n,mp}^2\,\rho\,d\rho\,d\phi\,dz},cn,mp=∫Φn,mp2ρdρdϕdz∫f(ρ,ϕ,z)Φn,mp(ρ,ϕ,z)ρdρdϕdz,
i.e., constants (no variable dependence), as required by Poisson/Laplace linearity.
- If axial–radial coupling exists in BCs, switch to Fourier in zzz → solve a parametrized Bessel ODE for each kkk → invert.
TL;DR
- You don’t need (and shouldn’t try) to match radial and axial eigenvalues. You get a double index and Λ=kn2+αmp2\Lambda=k_n^2+\alpha_{mp}^2Λ=kn2+αmp2.
- The radial piece is self-adjoint (with weight ρ\rhoρ) if you impose regularity at ρ=0\rho=0ρ=0 and constant-coefficient Robin/Dirichlet/Neumann at ρ=R\rho=Rρ=R. No singular YmY_mYm term on the axis.
- If your wall BC couples zzz and ρ\rhoρ, do a Fourier in zzz first (or use a mode-dependent Green’s function).
If you share the exact sidewall BC (the one behind your webp), I can write the transcendental equation for αmp\alpha_{mp}αmp explicitly and hand you the expansion formula in one go.
