这是一个基于拉普拉斯局部迭代的模拟闪电方法，但程序缺少主程序和头文件，求大神补充啊。

以下是算法关键代码：

／／模型控制模块
void chosenlist：：insertone(int px，int PY，int rllx，int my)//插入到数据结构。相应把电势更新，再
次迭代，并且把备选表中的节点删除，而且，加入相邻节点
{
node tmp(px，py，mx，my，true，false)；
tmp．1evel=g->u[mx][my]；        
g->u[mx][my]=O；
lchosenlist．push_back(tmp)；
}
void chosenlist：：addCandidate(int px，int PY，int mx，int my)//添加--个候选点，不做任何判断。
{
ModelControl(mx，my)；
node tmp(px，py，mx，my，true，false)；
lcandidates．push_back(tmp)；
}
void chosenlist：：ModelControl(int mx，int my)//对模型做控制。包括选择结束时机，区分闪电的
主干和分枝部分
{
if(my<=beginny+n_step&&isfirst)
{
if(nlift_count<4)
nlift_count++；
else
{
isfirst=false；
t_end=seconds()；
double t_dur=t_end-t_start；
cout<<”计算时间是：”<<t_dur<<endl；
rout<<”粒子数量是：”<<lchosenlist．sizeO<<endl；
for(int i=O；i<g->nx；i++)
for(intj=0；j<g->ny；j++)
{
g->u[i][j]=0；
}
list<node>：：iterator pos；
list<node>：：iterator pos_miny；
int miny = 99999;
for(pos=lchosenlist．begin()；posl=lchoserdist．End()；pos++)//求主干
{
g->u[pos->mx][pos->my]=pos->level；
if(pos->my<miny)
{
Miny = pos->my；
pos_miny = pos；
}
}．
Pos=pos_miny；
pos->isMainChanel=true；
while(1)
{
for(1ist<node>：：iterator father_pos=lchosenlist．Begin()；father_pos!=
lchosenlist．end()；father_pos++)
{
if(father_pos->mx==pos->px&&father_pos->my== pos->py)
{
father_pos->isMainChanel=true；
father_pos->level=1．2；
pos=father_pos；
break；
}
}
if(pos->px==g->nx/2&&pos->py==g->ny-1&&pos->mx==g->nx/2
&&pos->my==g->ny-2)
{
break；
}
}
}
}
}
void chosenlist：：addlncNeighbor(int mx．int my)//加入邻居队列
{
int mxleft=mx一1；
int mxright=m+l；
int mydown=my-1：
int myup = my+l；
if(mxleft>0&&mxleft<g->nx)
{
if(!bcandidate[mxleft][my])
{
bcandidate[mxleft][my]=true；
addCandidate(mx，my，mxleft．，my)；//正左方的邻居节点加入候选队列
}
if(mydown>O&&mydown<g->ny)
{
if(!bcandidate[mxleft][mydown])
{
bcandidate[mxleft][mydown]=true；
addCandidate(mx，my，mxleft，mydown)；//芝E下方的邻居节点加入候选队列
}
}
if(myup>O && myup<g->ny)
{
if(!bcandidate[mxleft][myup])
{
bcandidate[mxleft][myup]=true；
addCandidate(mx，my，mxleft，myup)；//左上方的邻居节点加入候选队列
}．
}
}
if(mxright>0 && mxright<g->nx)
{
if(!bcandidate[mxright][my])
{
bcandidate[mxright][my]=true；
addCandidate(mx，my，mxright，my)；//正右方的邻居节点加入候选队列
}
if(mydown>O&&mydown<g->ny)
{
if(!bcandidate[mxright][mydown])
{
bcandidate[mxright][mydown]=true；
addCandidate(mx，my，mxright，mydown)；//右下方的邻居节点加入候选队列
}
}
if(myup>O&&myup<g->ny)
{
if(!bcandidate[mxright][myup])
{
bcandidate[mxright][myup]=true；
addCandidate(mx，my，mxright，myup)；//右上方的邻居节点加入候选队列
}．
}
}
if(mydown>0&&mydown<g->ny)
if(!bcandidate[mx][mydown])
{
bcandidate[mx][mydown]=true；
addCandidate(mx，my，mx，mydown)；//正下方的邻居节点加入候选队列
}
if(myup>0&&myup<g->ny)
if(!bcandidate[mx][myup])
{
bcandidate[mx][myup]=true；
addCandidate(mx，my，mx，myup)；//正上方的邻居节点加入候选队列
}
}
void chosenlist：：update()
{
//选择一个最大可能的闪电粒子加入到闪电中
int RANGE_MIN=0；
int RANGE_MAX=100；
srand(time(O))；
double rand1=_twister．getDoubleLR()；
list<node>：：iterator pos；
list<node>：：iterator tmppos；
double totalPotential=0：
int k=lchosenlist．size()；
for(pos=lcandidates．begin()；pos!=lcandidates．end()；pos++)//遍历所有的候选节点，求出所有的节点电势和
{
totalPotential+=pow(g->u[pos->mx][pos一>my]，relay)；
}
double potentialSeen=0；
int index=-1；
for(pos=lcandidates．begin()；pos!=lcandidates．end()；pos++)//从所有候选点里面选出下一个闪电发生的地点
{

index++：
tmppos = pos；
if(potentialSeen<rand1)
{
potentialSeen+=pow(g->u[pos->mx][pos->my]，relay)
/totalPotential；
continue；
}
else
break；
}
//
pos=tmppos;
int firstx = g->nx/2；
int firsty = g->ny；
addlncNeighbor(pos->mx，pos->my)；
insertone(pos->px，pos->py，pos一>mx，pos一>my)；∥将生成的节点插入到节点列表k=lchosenlist．size（）；
lcandidates．erase(pos)脏候选队列中删除被选择的那个元素
laplacian(n_iter，g,&lchosenlist，a，b，n_step)；∥进行拉普拉斯局部迭代
}
∥拉普拉斯算子的迭代
Real LaplaceSolver：：timeStep(const Real dt)
{
Real dx2 = g->dx*g->dx；
Real dy2 = g->dy*g->dy；
Real dnr_inv=0．5/(dx2+dy2)；
Real tmp；
Real err=0．O：
int nx=g->nx；
int ny=g->ny；
Real**u= g->u；
nstep=1；
setltlnit()；
for(int i=firstnx+nstep；i<secondnx-nstep；i+=nstep){
for(int j=firstny+nstep；j<secondny-nstep；j+=nstep){
tmp=u[i][j]；
u[i][j]=((u[i—nstep][j]+u[i+nstepl[j])*dy2+
(u[i][j-nstep]+u[i][j+nstep])*dx2)*dnr_inv；
}
}
return 1；
}
Real LaplaceSolver：：solve(const int n_iter，const Real eps)
{
Real err=timeStep()；
int count=1；
while(1){
if(n_iter&&(count>=n_iter))
{
return 1；
}
timeStep()；
++count；
}
return Real(count)；
}
void LaplaceSolver：：setltlnit()
{
Real**u=g->u；
int nx=g->nx；
intny=g->ny；
u[nx/2][nx一1]_u[nx/2][nx-2]=u[nx/2][nx-3]=u[nx/2-1][nx-1]=u[nx/2-l][nx-2]=
u[nx/2-1][nx一3]=u[nx/2-2]nx-l]=u[nx/2-2][nx-2]=u[nx/2-2][nx-3]=0；∥上面的负
极点
for(intj=O；j<nx；++j)//底部的点，为正
{
u[j][0]=1；
}
list<node>：：iterator pos；
for(pos=1->begin()；pos!=1->end()；pos++)
{
u[pos->mx][pos->my]=O；
}
}
double laplacian(int n_iter，Grid*g,list<node>*mp_list，int mx，int my，int nstep)
{
g->setBCFunc()；初始化
LaplaceSolver s=LaplaceSolver(g，1)；
if(mx—nstep-1>0)
s．firstnx=mx-nstep-1；
else
s．firstnx=0：
if(my-nstep>0)
s．firstny = my-nstep-1；
else
s．firstny=0；
if(mx+nstep+1<g->nx)
s．secondnx=mx+nstep+1；
else
s．secondnx = g->nx；
if(my + nstep+1<g->ny)
s．secondny = my+nstep+1：
else
s．sccondny = g->ny；
s．1= mp_list；
Real t_start=seconds()；//记时
s．solve(n_iter,0)；//解拉氏方程
Real t_end=seconds()；
return(double)t_end-t_start)；
}