You have not specified the region where the problem is defined and the nature of the boundary conditions (for instance, are these conditions constant or dependent on the coordinates of the points of the boundary?). Also, you have not mentioned anything about the nature of the coefficients (for instance, are they easy to calculate?).

Thank you Professor Farid for your insight in the question.

1. My problem is defined on a plane/disc (with possible out-of-plane deformations)

2. Both U and V vanish at the edges of the plane/disc

3. The coefficients are real, continuous and differentiable over the whole plane; no discontinuities. The functions are not simple though, but some how complicated.

Please, I can send you the full equations if you don't mind.

Thank you

You are welcome Musah. You could always add the file of the full question to this page for all to see and possibly respond to.

The out-of-plane issues to which you refer are not reflected in the set of your differential equations, which explicitly depend only on x and y. This issue is therefore vague to me. Let us therefore restrict ourselves to the pure two-dimensional problem defined on a disc on the flat x-y plane. Further, since the problem involves second-order derivatives, values of U and V on the boundary are not sufficient; the normal derivatives of these on the boundary must also be specified (the derivatives tangent to the boundaries are by definition vanishing).

If the space of interest were a rectangular disc, given the zero boundary conditions for U and V, these functions would be expressible in terms of products of two sine functions, of the form sin(k_x x) sin(k_y y), where k_x and k_y both belong to discrete sets; by choosing the origin of the coordinate system at one of the corners of the rectangle (say, the lower-left corner), k_x and k_y would satisfy k_x L_x = n π and k_y L_y = n π, where n = 1,2,3, .... and L_x and L_y denote the lengths of the arms of the rectangle in the x and y directions (n=0 is excluded because sin(0) = 0, and negative n are excluded by the fact that sin(-x) = -sin(x) is not linearly independent of sin(x)). Denoting the latter k_x and k_y by k_x^n and k_y^n, one has the complete basis set {sin(k_x^m x) sin(k_y^n y) | m, n = 1,2,3, ...} in terms of which to expand the functions U(x,y) and V(x,y). Doing so, one readily obtains a matrix eigenvalue problem whose solution yields both the sought-after eigenvalues and the coefficients of the expansions of the eigenfunctions in terms of the above-mentioned set. To construct the above-mentioned matrix to be diagonalized, you will have to be able to Fourier transform the functions a(x,y), b(x,y), etc. Either this can be done analytically, or it will have to be done numerically, in which case one could preferably use fast-Fourier transformation (FFT). The matrix to be diagonalized is in practice always truncated by imposing upper bounds to the values of k_x^m and k_y^n. For sufficiently large upper bounds, the lowest-lying eigenvalues and the corresponding eigenfunctions will be accurately determined.

The above-approach is very similar to what one does (but often in one spatial dimension) in solving the time-independent Schrödinger equation for an electron inside a box with infinite walls. The discretization of the problem is natural, because k_x and k_y are discretized / quantized (see for instance Section 2.6 in Introduction to Quantum Mechanics, by David J. Griffiths (Prentice Hall, 1995), and take the limit of V_0 --> infinity).

For circular disc, you discretize the problem as above, through employing the Fourier-Bessel series expansion, instead of the above double-sine expansion. For this, you may consult the following as a starting point:

http://en.wikipedia.org/wiki/Fourier%E2%80%93Bessel_series

Essentially, what one does is using the cylindrical polar coordinate system. Then one uses separation of variables. The boundary conditions correspond to a fixed r, and in the φ direction one has periodicity, which enables one to expand the functions of φ in terms of discrete Fourier series (note that (r,φ) is defined on a rectangular region: [0,r_0] x [0,2π], where r_0 denotes the radius of the circular disc). The solution of the problem is in a way similar to solving the Schrödinger equation for the Hydrogen atom, but in two dimensions and subject to a somewhat different boundary condition in the radial direction.

Dear professor Farid, thank you so much for the thorough discussion of my problem.

Yes, it is true that I should have specified both Driclet and Neuman boundary conditions. Sorry for mentioning the the out-of-plane d.o.f, it is actually captured through a condition like z^2=x^2+y^2.

In fact, I have no doubt your answer will make an impact on my ongoing research. More expecially the last part on cylindrical coordinates. Unfortunately seperation of variables will be difficult to achieve (analytically) in my case. for instance, the coefficients contain such factors as cos[r sin( φ)]

I have read Griff's Quantum book before but with ur advice I need to look at it again.

Please, may I get your lecture notes on this topic.

Thank you.

Dear Musah, you are welcome. I have no lecture notes of my own on the subject matter. Regarding the Fourier-Bessel expansion, you could consult the book Mathematical Methods for Physicists by G.B. Arfken & H.J. Weber. Regarding separation of variables, it is not necessary for solving your particular problem. The double-sine basis I presented earlier, has the advantage that it automatically satisfies the boundary conditions, is countable, is orthogonal, and the matrix representation of your eigenvalue problem can be determined with the aid the FFT. Complete separability would be convenient and welcome, in particular if you aimed at solving your problem analytically, or if your problem were so computationally demanding that without complete separability, if would be impractical to solve. As I see it, your problem is not computationally very demanding. Very similar arguments apply to the use of the Fourier-Bessel expansion (for dealing with a circular disc).

Dear Professor Farid, I have put my equations in a pdf file which is attached here.

I am sorry for the delay in putting the equations here, I had to travel from my country to ICTP-Italy for a program. Could you let me make further clarification if the need arises.

I thank you in advance.

Dear Rabiu, thank you for the details. In order to make progress, I recommend you to use bra-ket notation. That will make the calculations more transparent.

First some observations. Eqs (2) and (3) constitute an eigenvalue problem. After having solved this problem, Ψ1 and Ψ2 are to be used to evaluate the second term on the left-hand side of Eq. (1), after which one has a partial differential equation to solve for Φ. Further, since you have used cylindrical polar coordinates, some of the "boudary conditions" you present lack meaning. For instance, Ψ(0,φ) = ψ(φ) is meaningless, since φ at r=0 is undefined. Also, for instance, Ψ(r,π) = Ψ(r,-π) is redundant, since the partial derivatives with respect to φ in your differential equations necessite continuity in the φ direction. Etc. Importantly, the Bessel-Fourier series for each of Ψ1(r,φ) and Ψ2(r,φ) must involve a double sum, one for the Bessel functions decribing the r-dependence of the latter functions, and one for describing their φ dependence; the latter series is just a Fourier series (compare with the Hydrogen problem and the way one arrives at the spherical harmonics); it is not for nothing that the expansion under discussion is called Bessel-Fourier expansion.

After correcting all the above-mentioned problems, you will have to substitute the expression for Ψ1(r,φ) and Ψ2(r,φ) in Eqs (2) and (3), multiply on the left with the complex conjugate of the basis functions and integrate with respect to r and φ. By the orthogonality of the Bessel functions (as presented in Eq. (11.49) of the book by Arfken and Web) and the exponential functions corresponding to the underlying Fourier series, in this way you will obtain two coupled matrix eigenvalue problems. Let now the Bessel-Fourier expansions of Ψ1(r,φ) and Ψ2(r,φ) involve N terms each (equal to the totality of the indeces used in the series involving the Bessel functions and the series involving the exponential functions of the Fourier series). Let further the column N-vector u1 contain the expansion coefficients of Ψ1(r,φ) and the column N-vector u2 the expansion coefficients of Ψ2(r,φ). You will obtain two coupled matrix eigenvalue problems of the form:

A u1 + B u2 = ε u1

C u1 + A u2 = ε u2

where A, B and C are N x N matrices. Now, if you define the 2N-vector w in such a way that its first N components coincide with those of u1 and the next N components with those of u2, the above equations can be cast into the following standard matrix eigenvalue problem:

D w = ε w,

where D is a 2N x 2N matrix, consisting of 4 blocks: its (1,1) and (2,2) blocks coincide with A, its (1,2) block with B and its (2,1) block with C.

Lastly, I strongly recommend you to solve a simpler problem along the above lines, before tackling your actual problem. Take the time-independent Schrödinger equation for a particle on a two-dimensional plane, subject to a potential that is constant (let us say, zero) inside a disc, and infinite outside. Solve this problem numerically by using the Bessel-Fourier series expansion for the eigenfunctions. Once you have managed solving this problem, solving your actual problem will be far more easy.

Dear Professor Farid, thank you once more for your very comprehensive and detail insight. Please,

1. In your N-by-N matrix equation for the coeffients, you left out "Phi" term. Did you assume its fourier coefficients are orthogonal to the Ψ(r,φ) coefficients?

2. Sorry for my redundancies in the boundary conditions (BCs). Trully, I was initially working with explicit time-dependence in the equations. I converted them to Omega and ε and forgot to remove them in the BCs, particularly the Ψ(0,φ) = ψ(φ). The zero was meant for time. So this condition should be ignored. I also understand that I was very wrong with the BC: Ψ(r,π) = Ψ(r,-π), since this gives discrete BCs. The right one to employ should be Ψ(r,φ) = Ψ(r,φ+2π) or Ψ(r,0) = Ψ(r, 2π).

3. Eventhough I am familiar with Dirac Bra-Ket notation, I am not very sure how I will rewrite my PDEs in the notation. Or you meant writing them like;

Del^2 Ψ(r,φ) = L Ψ(r,φ)  as  < r, φ|Del^2|Ψ> = L < r, φ|Ψ>.

Please, I  don't understand your last statement!!!

Thank you

Dear Musah, you are welcome. As for your questions: 1. I did not leave out anything; while your problem is not fully separable, the basis functions (along the r and φ directions) separate the two variables. In Eqs (2) and (3) the function Φ is assumed as given, and as I wrote in my previous response, it has to be solved separately. This of course implies that the actual problem is to be solved self-consistently (through using some iterative process). Of course, the function Φ(r,φ) can be expressed in terms of the same set of basis functions as used for expressing Ψ1(r,φ) and  Ψ2(r,φ). 2. Agree. 3. Almost yes. "Del^2" is in the coordinate representation (some writers ignore this fact). You should employ the abstract p-operator, taking into account the h-bar involved in it, or employ an alternative symbol. I proposed use of the bra-ket notation by the fact that deducing matrix equations with the aid of this notation does not suffer from unnecessary integral symbols and thus results in quite compact expressions. This use is by no means mandatory. You could therefore just use the coordinate representation if you feel more comfortable with it.

Thank you for your prompt reply. At least, I can continue from here now.

I am most grateful for all your contributions.

You are most welcome Musah. I would strongly recommend you to begin with dealing with a simpler eigenvalue problem in two dimesnions, with cylindrical symmetry. Once you have that under control (both algebrically and numerically), dealing with your present problem will be much easier. In general, the simpler the problem, the more insight one gains by treating them in detail.

Thank you professor Farid.

In fact, I am looking at hydrogen problem now.

You are welcome. Hydrogen is fine. There is one thing, and that is that in hydrogen the boundary condition in the radial direction corresponds to infinity. Of course, for the bound states of hydrogen one can always choose a sufficiently large but finite value of r, say R, at which to required the radial wave function and its derivative to be vanishing (this due to the exponential decay of the bound states for r --> infinity).

Dear professor,

Now, when do we use Bessel or Lagurre Solutions? I have seen these two types of solutions for the same system.

The Bessel function in the Fourier-Bessel expansion is the Bessel function of the first kind. This has to do with the analytic property expected of the solution and the analytic property of the Bessel function of the first kind. Similarly as regards the Lagendre function that one encounters in the solution of the Schrödinger equation for Hydrogen; the choice of Pl(x) (the Legendre function of the first kind), to be contrasted with Ql(x) (the Legenrde function of the second kind), is related to the analytic property of Pl(x); the function Ql(x) is logarithmically singular at x=+/- 1. For details, consult any textbook on quantum mechanics. See also the book Handbook of Mathematical Functions, edited by Abramowitz and Stegun (Dover Publications).