1
$\begingroup$

I'm reading a paper in which the authors solve the following equation:

$\frac{d^{2}}{dz^{2}}\hat{p}$($\bf{q}$$,z)$-$q^{2}\hat{p}$($\bf{q}$$,z)$-$\frac{iq_{y}}{(2\pi)^{2}}\delta(z-z_{2})$=0

here $\bf{q}$=$(q_{x},q_{y})$ and $q^{2}$=$q_{x}^{2}$+$q_{y}^{2}$

$\hat{p}$($\bf{q}$$,z)$ is the fourier transform of the real function $p$($\bf{s}$$,z)$:

$\hat{p}$($\bf{q}$$,z)$= $\int$$p$$(\bf{s}$$,z)$$e^{-i\bf{q}\cdot\bf{s}}$$d\bf{s}$

$z_{2}$ is a parameter in the differential equatin , $i$ is the imaginary unit.

the autors report the following solution:

$N_{1} cosh(qz)$+$N_{2} sinh(qz)$-$\frac{iq_{y}}{8\pi^{2}q}$$e^{|z-z_{2}|}$

$N_{1}$ and $N_{2}$ are coefficients that depend on the boundary conditions.

when i'm solving that equation with mathematica i'm getting different results:

I'm considering $p$ as being just a function of $z$ and treating $\bf{q}$ as a parameter

 DSolve[p''[z] - (q^2) p[z] - (I*qy/(2 Pi)^2) DiracDelta[z - z2] == 0, 
  p[z], z]

and i get:

p[z] = E^(q z) C[1] + E^(-q z) C[2] + (
 i E^(q z - q z2) qy HeavisideTheta[z - z2])/(8 Pi^2 q) - (
 i E^(-q z + q z2) qy HeavisideTheta[z - z2])/(8 Pi^2 q)

which is different from the one the authors are reporting in the paper.

Please can someone explain me how to corectly solve that equation?

thanks in advance

$\endgroup$
10
  • $\begingroup$ One thing you should do if you haven't already is compare the two answers for a specific choice of boundary conditions. It may be that the expression in Heaviside step functions is equivalent to theirs (which is plausible, given that absolute values also change their behavior discontinuously at 0) $\endgroup$ Jul 30, 2014 at 14:11
  • $\begingroup$ I have not done it yet, because i don't have the boundary conditions directly on $\hat{p}$ but i need to substitute the solution for $\hat{p}$ in another expression. I Suspect that the origin in the difference between the two solutions lies in the fact that $\hat{p}$ is a complex function while mathematica solves the differential equation as if $\hat{p}$ was a real function. $\endgroup$
    – SSC Napoli
    Jul 30, 2014 at 14:18
  • $\begingroup$ eh, maybe. but even if you don't know which boundary conditions you 'really' want, you can still pick some more or less arbitrarily just to test the correspondence $\endgroup$ Jul 30, 2014 at 14:21
  • $\begingroup$ anyway i considered the case in which all the constants are zero, the solution given by mathematica doesn't match in this case. this is because the part with heavyside adds a contribution to the solution only when $z > z_{2}$ while the abslute value is active even if $z < z_{2}$ $\endgroup$
    – SSC Napoli
    Jul 30, 2014 at 14:29
  • $\begingroup$ thank you very much, i'm really having bad times with this problem... $\endgroup$
    – SSC Napoli
    Jul 30, 2014 at 14:47

1 Answer 1

0
$\begingroup$

First, note that the first two terms of their solution are equivalent to your first two terms under the identification $(C_1,C_2)=\frac{1}{2}(N_1+N_2,N_1-N_2)$. So the only potential problem is in the second halves of the two results; even there, we can recognize a common factor of $-i q_y/8\pi q^2$, so the only issue is the exponentials themselves.

To see the issue more clearly, let's explicitly distinguish the $z<z_2$ and $z>z_2$ cases. We have $$e^{-q|z-z_2|}=\begin{cases} e^{q(z-z_2)}, & z<z_2 \\ e^{q(z_2-z)}, & z>z_2 \\ \end{cases}$$ This can be written in terms of Heaviside step functions as $\Theta(z-z_2) e^{q(z_2-z)}+\Theta(z_2-z) e^{q(z-z_2)}$, so that the chice of which sign is in the exponent flips at $z=z_2$. This is nearly the same as the Mathematica answer, except that the step functions in the Mathematica answer are identical. Consequently the Mathematica answer has the second half of the answer being zero for $z < z_2$ and equal to the wrong value for $z > z_2$. So I'd trust the paper's result more than Mathematica...

If you really want Mathematica to do it right, you probably have to introduce the Dirac delta via boundary conditions rather than in the function itself. That is, you should integrate your ODE from $z=z_2-\epsilon$ to $z=z_2+\epsilon$. This will imply a 'jump' condition for the first derivative, which can be used as an internal boundary condition for the problem.

$\endgroup$
3
  • $\begingroup$ thanks very much for the answer, i already recognized that the difference lied only in the second part of the solution, i will try what you suggest and let you know the result. $\endgroup$
    – SSC Napoli
    Jul 30, 2014 at 15:44
  • $\begingroup$ hi, i tried to solve the equation in the following way: DSolve[{p''[z] - p[z]*q^2 == 0, p[z2] == (I*qy/(4*Pi^2))}, p[z], z] however i get this solution: p[z] = 1/ 4 E^(-q z) ((i E^(q z2) qy)/Pi^2 + 4 E^(2 q z) C[1] - 4 E^(2 q z2) C[1]) which is clearly not correct, $\endgroup$
    – SSC Napoli
    Jul 31, 2014 at 8:20
  • $\begingroup$ I think the following is happening: i'm solving a differential equation in the indipentent variable z and treating $\bf{q}$ as a parameter, therefore when mathematica solves the differential equation when it crosses $z_{2}$ applies the dirac delta effects. However $\hat{p}(\bf{q},z)$ the Dirac delta is multiplied by $i q_{y}$ meaning that the Dirac delta is not applied exactly at $z_{2}$ but rather at $z_{2} + i*q_{y}$. Maybe because $\hat{p}(\bf{q},z)$ is a complex function i have to split it into the real part and the imaginary part? have i been clear? $\endgroup$
    – SSC Napoli
    Jul 31, 2014 at 8:34

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .