0
$\begingroup$

I want to solve the following differential equation from a paper with the boundary condition: enter image description here

The paper used the Crank–Nicolson method for solving it. I think I understand the method after googling it, but most websites discussing it use the heat equation as an example. Here we can replace the usual t variable with Xi, and the usual x as Rho. Now unlike the heat equation, this equation has an additional first derivative of Xi term (the first term gives both first and second derivative of Xi). Does it matter what finite difference scheme I use to approximate this term?

Another thing is that the equation has nonlinear terms involving |A|^2 and |A|^4. How can I incorporate this to the matrix equation I am going to solve? (PS: A is complex)

I am pretty new to numerical methods so forgive me if I have neglected anything significant.

Improve this question
$\endgroup$
1
  • $\begingroup$ There are similar questions with Crank Nicolson and nonlinear PDEs, e.g. this, where I have explained the concept and how it is applied. What you need first is a solver for the linear problem and then use the Newton-Kantorovich scheme on every "time step" with the Frechet Derivative of your functional. $\endgroup$
    – Bort
    Feb 12, 2018 at 14:47

1 Answer 1

Reset to default
1
$\begingroup$

You need at least one initial condition for $\xi$ and two BC for $\rho$. Where are they? Your equation looks like the heat equation in cylindrical coordinates assuming angular and plane symmetry, with nonlinear heat sources: $$ 2ik\partial_tA-\triangle A=N(A)$$.

I would discretise the spatial derivative with finite differences (which is essentialy the same as finite elements in this case) and discretise the dynamic system with what you've proposed. If you discretise in space you will obtain: $$ \partial_t\vec{A}-L(\rho) \vec{A}=\vec{N}(\vec{A})$$ Where $L(\rho)$ is the matrix operator for the Laplacian.

This operator looks like $$L(\rho)=L+RG$$ Where $L$ is the famous Laplacian matrix operator in cartesian coordinates and $G$ the gradient operator. Without any BC these operators look like: $$L(i,i)=-2/\Delta \rho^2, L(i-1,i)=L(i,i+1)=1/\Delta \rho^2 $$ $$ G(i,i+1)=-G(i,i-1)=1/\Delta\rho$$ $$ R(i,i)=1/\rho_i$$

The resulting scheme will be sth. like this: $$ [2ki\mathbb{I}-\Delta \xi/2L(\rho)]\vec{A}^{n+1}+\Delta\xi/2\vec{N}(\vec{A}^{n+1})= 2ki\vec{A}^n+\Delta\xi/2 L(\rho)\vec{A}^{n} +\Delta\xi/2\vec{N}(\vec{A}^{n})$$

Has been obtained from the general Crank-Nicolson scheme: $$\frac{W^{n+1}-W^n}{\Delta t}=\frac{1}{2}[F(W^{n+1})+F(W^n)]+\mathcal{O}(\Delta t^2)$$ for the differential equation: $$ \frac{dW}{dt}=F(t,W)$$

For each $n$ you will solve a nonlinear system of equations...

Improve this answer
$\endgroup$
10
  • $\begingroup$ Using finite differences in cylindrical coordinate systems is difficult, because you need to do something with the node at $\rho=0$ where the equation is singular. It's simpler to do this with finite elements. $\endgroup$
    – Wolfgang Bangerth
    Feb 6, 2018 at 18:33
  • $\begingroup$ Also, you really only need one boundary condition at $\rho=R$. No boundary condition is necessary at $\rho=0$ -- it is implied by the symmetry of the solution. $\endgroup$
    – Wolfgang Bangerth
    Feb 6, 2018 at 18:34
  • $\begingroup$ The simplest method (in cartesian grid) is FD whatever the equation, in my opinion. Difficult would be if used for solving problems with curved computational domain. The singularity is avoidable if Neuman BC is provided, and it does not matter if Dirichlet, so there is not any problem with that issue. You are assuming that is a solid cylinder, I am considerig a more general case with a hole. If you compute the matrix operators you will find that actually two conditions must be given for the uniqueness of the solution. Name them as you wish. $\endgroup$
    – HBR
    Feb 7, 2018 at 11:04
  • $\begingroup$ If it's a cylindrical shell, you are correct. If the cylinder is solid, mathematically you can't pose a Neumann-type boundary condition. $\endgroup$
    – Wolfgang Bangerth
    Feb 7, 2018 at 14:43
  • $\begingroup$ So the second equation A_entrance is basically the initial condition for ξ at ξ=0 (not to be confused with z, which is just a parameter that I input). The equation describes the propagation of electric-field in nonlinear medium, so A tends to zero as rho tends to infinity as the A has a Gaussian profile initially (and shouldn't deviate from Gaussian shape too much after the whole propagation). Are there anymore BC/IC I need? $\endgroup$
    – Physicist
    Feb 9, 2018 at 9:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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