當兩個物體之間存在較大的電勢差時會出現放電現象，比如生活中常見的閃電現象，閃電形成的條件就是雲層積累了大量負電荷之后與地面之間形成了強大的電勢差。目前關於閃電建模的方法比較少，下面介紹一種利用電介擊穿模型來模擬閃電的方法，電介擊穿模型可以模擬自然界許多現象，該方法通過迭代求解Laplace方程得到放電過程的中間狀態。
初始電位結構如下圖所示，首先在2維柵格正中心的單元放置一個負電荷Ф = 0（灰色），然后在其周圍放置一圈正電荷Ф = 1（黑色），而其他柵格單元可以通過求解Laplace方程得到：
Laplace方程求解完之后，我們構建一個列表記錄負電荷（Ф = 0）周圍的柵格單元，並隨機選取其中一個單元作為下一個被擊穿的單元，被選中的柵格單元設置Ф = 0，同時把它作為下一步迭代計算時的邊界條件。
相鄰柵格單元被擊穿的概率與其計算得到的電勢相關，其概率如下：
式中i代表第i個與負電荷相鄰的柵格單元，n代表與負電荷相鄰的柵格單元個數，而η是一個用戶設定的參數，通過實驗表明η值可以控制電弧的分叉密度，當η = 1時，分叉密集，當η逐漸增大時，分叉密度慢慢減小。
圖：不同結構的初始電位，其中灰色點代表負電荷Ф = 0，黑色點代表正電荷Ф = 1
圖：初始電位為（b）結構時的模擬結果
左：η = 1，中：η = 2，右：η = 3
% Laplacian Growth Model clc clear all close all w = 256; % width h = 256; % height eta = 3; % 指定正負電荷 ni = floor(w/2)*h + (1:2)'; % negative charge pi = ((1:w)' - 1)*h + h; % positive charge % 構建laplace矩陣 S1 = bsxfun(@plus, ((1:(w-1)) - 1)*h, (1:h)'); D1 = bsxfun(@plus, ((2:w) - 1)*h, (1:h)'); S2 = bsxfun(@plus, ((1:w) - 1)*h, (1:(h-1))'); D2 = bsxfun(@plus, ((1:w) - 1)*h, (2:h)'); S = [S1(:); S2(:)]; D = [D1(:); D2(:)]; A = sparse([S D], [D S], 1); n = size(A,1); % Lp = speye(h*w) - bsxfun(@rdivide, A, sum(A,2)); Lp = spdiags(sum(A,2), 0, n, n) - A; rhs = zeros(n, 1); % 循環計算下一步狀態 [X0, Y0] = ind2sub([h,w], (1:h*w)'); im = zeros(h,w); iter = 1; i = 1; while ~any(ismember(ni, pi)) phi = solve_equation(Lp, rhs, [ni;pi], [zeros(size(ni));ones(size(pi))]); [~, adj] = find(A(ni,:)); adj = unique(adj); adj = adj(~ismember(adj, ni)); k = randsample(adj, 1, true, phi(adj).^eta/sum(phi(adj).^eta)); ni = [ni; k]; if mod(iter,10) == 1 [X, Y] = ind2sub([h,w], ni); [~, dis] = knnsearch([X,Y], [X0,Y0]); im = reshape(1-dis/3, [h,w]); imshow(im, 'Border','tight'); drawnow; Frame(i) = getframe(gcf); i = i + 1; end fprintf('iteration: %d ...\n', iter); iter = iter + 1; end
本文為原創，轉載請注明出處：http://www.cnblogs.com/shushen。
參考文獻：
[1] Theodore Kim and Ming C. Lin. 2007. Fast Animation of Lightning Using an Adaptive Mesh. IEEE Transactions on Visualization and Computer Graphics 13, 2 (March 2007), 390-402.