Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

I'm trying to simulate a lightning in a programming language and I started to read how it can be done, and I found that it could be done by using the Laplacian Growing Model. My experience with mathematics is limited.

I think I understood the algorithm described here in section 4.1. My problem is that I don't know how to solve the Laplacian equation for the lightning model seen at (b) Lightning configuration.

I tried to catch something up from here but that configuration is a rectangle and the lightning configuration seems like a triangle or maybe a circular arc.

I know my boundary conditions would be the potential at the origin and the potential at the ground but I really don't know how to translate this in mathematical language and how to resolve the laplace for it.

Can you help me with a solution (along with an explanation) or a tutorial "for dummies"?

share|improve this question
    
Do you have a preferred language? Would matlab be okay? –  Omnomnomnom Jul 28 '13 at 20:20
    
@Omnomnomnom I want to implement it in Java, but I need the theoretical solution for the problem, but I think I can understand some Matlab too. But my question is theoretical, language-adnostic. –  DaJackal Jul 29 '13 at 10:51
add comment

1 Answer 1

up vote 5 down vote accepted

Here's a description of how to solve the laplacian equation for an arbitrary $2$-D grid with boundary conditions using the so-called "iteration method":


Initialization:

Create a matrix $M$ representing each point on the grid, with one entry per "pixel". So, for example, to work with a $7\times 7$ grid, we would set $$ M = \begin{bmatrix} 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 \end{bmatrix} $$ We then set permanent boundary conditions that won't be changed as we improve our approximation. To get something kind of circular, we could set $$ M = \begin{bmatrix} \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1\\ \mathbf1 & \mathbf1 & 0 & 0 & 0 & \mathbf1 & \mathbf1\\ \mathbf1 & 0 & 0 & 0 & 0 & 0 & \mathbf1\\ \mathbf1 & 0 & 0 & \mathbf0 & 0 & 0 & \mathbf1\\ \mathbf1 & 0 & 0 & 0 & 0 & 0 & \mathbf1\\ \mathbf1 & \mathbf1 & 0 & 0 & 0 & \mathbf1 & \mathbf1\\ \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 \end{bmatrix} $$ The bold entries remain constant, while the rest will update as we improve our approximation.

We then normally put in some non-zero initial guess for each point to be computed, just a very rough approximation. For instance, in this case we could put $$ M = \begin{bmatrix} \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1\\ \mathbf1 & \mathbf1 & 0.66 & 0.66 & 0.66 & \mathbf1 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & 0.33 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & \mathbf0 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & 0.33 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & \mathbf1 & 0.66 & 0.66 & 0.66 & \mathbf1 & \mathbf1\\ \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 \end{bmatrix} $$

We are now ready to begin the iterative process


Iteration:

What we use now (and will justify later, maybe) is the fact that on a $2$D plane where $\phi$ is a function satisfying $$ \nabla^2\phi=\frac{\partial^2\phi}{\partial x^2}+\frac{\partial^2\phi}{\partial y^2}=0 $$ the following is approximately true: Suppose that $\phi(x_0,y_0)=\phi_0$ Let $\Delta x$ and $\Delta y$ be "sufficiently small" and define: $$ \phi_{x^+}=\phi(x_0+\Delta x,y_0)\\ \phi_{x^-}=\phi(x_0-\Delta x,y_0)\\ \phi_{y^+}=\phi(x_0,y_0+\Delta y)\\ \phi_{y^-}=\phi(x_0,y_0-\Delta y) $$ Then we have $$ \phi_0 \approx \frac14(\phi_{x^+}+\phi_{x^-}+\phi_{y^+}+\phi_{y^-}) $$ In other words, we could calculate the potential at a point by averaging the potential at the four neighboring points. We do exactly that, using the improved values as we go.

For our matrix $M$, let's begin with the first changeable point in the upper-left corner.

$$ \begin{align} \phi_{2,3} &=\frac14(\phi_{2,4}+\phi_{2,1}+\phi_{1,3}+\phi_{3,3}) \\&= 0.25\cdot(0.66 + 1 + 1 + 0.33) \\&= 0.7475 \end{align} $$

Making the new matrix

$$ M = \begin{bmatrix} \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1\\ \mathbf1 & \mathbf1 & 0.7475 & 0.66 & 0.66 & \mathbf1 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & 0.33 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & \mathbf0 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & 0.33 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & \mathbf1 & 0.66 & 0.66 & 0.66 & \mathbf1 & \mathbf1\\ \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 \end{bmatrix} $$ In fact, if we know that we will be using a symmetrical configuration, we could make things more efficient by updating all points that will have the same outcome. That is, we could update $M$ to get

$$ M = \begin{bmatrix} \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1\\ \mathbf1 & \mathbf1 & 0.7475 & 0.66 & 0.7475 & \mathbf1 & \mathbf1\\ \mathbf1 & 0.7475 & 0.33 & 0.33 & 0.33 & 0.7475 & \mathbf1\\ \mathbf1 & 0.66 & 0.33 & \mathbf0 & 0.33 & 0.66 & \mathbf1\\ \mathbf1 & 0.7475 & 0.33 & 0.33 & 0.33 & 0.7475 & \mathbf1\\ \mathbf1 & \mathbf1 & 0.7475 & 0.66 & 0.7475 & \mathbf1 & \mathbf1\\ \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 \end{bmatrix} $$

We would then proceed to the next point $$ \phi_{2,4}=0.25\cdot(0.7475+0.7475+1+0.33)= 0.5194 $$ And update the matrix accordingly

$$ M = \begin{bmatrix} \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1\\ \mathbf1 & \mathbf1 & 0.7475 & 0.5194 & 0.7475 & \mathbf1 & \mathbf1\\ \mathbf1 & 0.7475 & 0.33 & 0.33 & 0.33 & 0.7475 & \mathbf1\\ \mathbf1 & 0.5194 & 0.33 & \mathbf0 & 0.33 & 0.5194 & \mathbf1\\ \mathbf1 & 0.7475 & 0.33 & 0.33 & 0.33 & 0.7475 & \mathbf1\\ \mathbf1 & \mathbf1 & 0.7475 & 0.5194 & 0.7475 & \mathbf1 & \mathbf1\\ \mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 &\mathbf1 \end{bmatrix} $$

Proceed in such a fashion until each changeable value is updated once (in this case, that will take $4$ steps if you take advantage of symmetry, $20$ if you make no such simplifications), then repeat. Do this enough, and eventually the updates won't significantly change the value. When that happens, you've reached the final approximation and can stop the process. You now have an approximate solution to the Laplacian equation, given the boundary conditions of your matrix.

A similar method to this one is the relaxation method, which is faster for systems with relatively few pixels and slightly more accurate.


EDIT:

I found this derivation of the iteration method (apparently, also known as the finite difference method) if you want a mathematical justification for the process.

share|improve this answer
    
Let me know if the method is clear as explained, and let me know if you want the mathematical justification –  Omnomnomnom Jul 29 '13 at 17:32
    
Thank you very much @Omnomnomnom. At the first glance, I think I udnerstood your solution, I can't wait to implement it. –  DaJackal Jul 30 '13 at 7:00
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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