如何提高计算编程技能

How to improve programming skill for computation?

本文关键字:编程 计算 何提高      更新时间:2023-10-16

我是物理专业的学生,编写模拟电磁波的代码。我在这里发布了我的一个代码,我的问题是如何改进我的代码。(不是速度)

#include <iostream>
#include <vector>
#include <math.h>
#include <fstream>
using namespace std;
#define IE 100
#define JE 100
#define KE 100
double ***alloc3d(int l, int m, int n) {
double *data = new double [l*m*n];
double ***array = new double **[l];
for (int i=0; i<l; i++) {
    array[i] = new double *[m];
    for (int j=0; j<m; j++) {
        array[i][j] = &(data[(i*m+j)*n]);
    }
}
return array;
}

int main()
{
//////////////////////declartion/////////////////////////////
int const NFREQS=100,ia=7,ja=7,ka=7;
double const pi=3.14159;
int i,j;
double ***gax=alloc3d(IE,JE,KE);
double ***gay=alloc3d(IE,JE,KE);
double ***gaz=alloc3d(IE,JE,KE);
double ***gbx=alloc3d(IE,JE,KE);
double ***gby=alloc3d(IE,JE,KE);
double ***gbz=alloc3d(IE,JE,KE);
double ***dx=alloc3d(IE,JE,KE);
double ***dy=alloc3d(IE,JE,KE);
double ***dz=alloc3d(IE,JE,KE);
double ***ex=alloc3d(IE,JE,KE);
double ***ey=alloc3d(IE,JE,KE);
double ***ez=alloc3d(IE,JE,KE);
double ***sx=alloc3d(IE,JE,KE);
double ***sy=alloc3d(IE,JE,KE);
double ***sz=alloc3d(IE,JE,KE);
double ***hx=alloc3d(IE,JE,KE);
double ***hy=alloc3d(IE,JE,KE);
double ***hz=alloc3d(IE,JE,KE);
double ***ix=alloc3d(IE,JE,KE);
double ***iy=alloc3d(IE,JE,KE);
double ***iz=alloc3d(IE,JE,KE);
double ***idyl=alloc3d(IE,JE,KE);
double ***ihyl=alloc3d(IE,JE,KE);
double*** idxl=new double**[ia];
double*** idxh=new double**[ia];
double*** ihxl=new double**[ia];
double*** ihxh=new double**[ia];
for(i=0;i<ia;i++)
{
    idxl[i]=new double *[JE];
    idxh[i]=new double *[JE];
    ihxl[i]=new double *[JE];
    ihxh[i]=new double *[JE];
    for(j=0;j<JE;j++)
    {
        idxl[i][j]=new double [KE];
        idxh[i][j]=new double [KE];
        ihxl[i][j]=new double [KE];
        ihxh[i][j]=new double [KE];
    }
}
double***idzl=new double**[IE];
double***idzh=new double**[IE];
double***ihzl=new double**[IE];
double***ihzh=new double**[IE];
for(i=0;i<IE;i++)
{
    idzl[i]=new double *[JE];
    idzh[i]=new double *[JE];
    ihzl[i]=new double *[JE];
    ihzh[i]=new double *[JE];
    for(j=0;j<JE;j++){
        idzl[i][j]=new double [ka];
        idzh[i][j]=new double [ka];
        ihzl[i][j]=new double [ka];
        ihzh[i][j]=new double [ka];
    }
}
double *gi1=new double [IE];
double *gi2=new double [IE];
double *gi3=new double [IE];
double *gj1=new double [JE];
double *gj2=new double [JE];
double *gj3=new double [JE];
double *gk1=new double [KE];
double *gk2=new double [KE];
double *gk3=new double [KE];
double *fi1=new double [IE];
double *fi2=new double [IE];
double *fi3=new double [IE];
double *fj1=new double [JE];
double *fj2=new double [JE];
double *fj3=new double [JE];
double *fk1=new double [KE];
double *fk2=new double [KE];
double *fk3=new double [KE];
int    l,n,m,k,ic,jc,kc,nsteps,numsph,npml,ib,jb,kb,ixh,jyh,kzh,ke,percent;
double xn,xxn,xnum,xd,curl_e,curl_h,freqi,freqf;
double *ez_inc=new double [JE];
double *hx_inc=new double [JE];
double t0,spread,pulse;
double *freq=new double [NFREQS];
double *arg=new double [NFREQS];
double***real_pt=new double**[NFREQS];
double***imag_pt=new double**[NFREQS];
for(m=0;m<NFREQS;m++)
{
    real_pt[m]=new double *[IE];
    imag_pt[m]=new double *[IE];
    for(i=0;i<IE;i++){
        real_pt[m][i]=new double [JE];
        imag_pt[m][i]=new double [JE];
    }
}
double***real_pt_sx=new double**[NFREQS];
double***imag_pt_sx=new double**[NFREQS];
for(m=0;m<NFREQS;m++)
{
    real_pt_sx[m]=new double *[IE];
    imag_pt_sx[m]=new double *[IE];
    for(i=0;i<IE;i++){
        real_pt_sx[m][i]=new double [JE];
        imag_pt_sx[m][i]=new double [JE];
    }
}
double***real_pt_sy=new double**[NFREQS];
double***imag_pt_sy=new double**[NFREQS];
for(m=0;m<NFREQS;m++)
{
    real_pt_sy[m]=new double *[IE];
    imag_pt_sy[m]=new double *[IE];
    for(i=0;i<IE;i++){
        real_pt_sy[m][i]=new double [JE];
        imag_pt_sy[m][i]=new double [JE];
    }
}
double***real_pt_sz=new double**[NFREQS];
double***imag_pt_sz=new double**[NFREQS];
for(m=0;m<NFREQS;m++)
{
    real_pt_sz[m]=new double *[IE];
    imag_pt_sz[m]=new double *[IE];
    for(i=0;i<IE;i++){
        real_pt_sz[m][i]=new double [JE];
        imag_pt_sz[m][i]=new double [JE];
    }
}
double *real_in=new double[NFREQS];
double *imag_in=new double[NFREQS];
double *amp_in=new double[NFREQS];
double *phase_in=new double[NFREQS];
ic=IE/2;
jc=JE/2;
kc=KE/2;
ddx=2e-9;
dt=ddx/6e8;
epsz=8.8e-12;
ib=IE-ia-1;
jb=JE-ja-1;
kb=KE-ka-1;
/////////////////////////////initilization///////////////////////
ez_low_m2=0;
ez_high_m2=0;
ez_low_m1=0;
ez_high_m1=0;
for(j=0;j<JE;j++){
    ez_inc[j]=0.0;
    hx_inc[j]=0.0;
    for(i=0;i<IE;i++){
    }}
for (k=0;k<KE;k++){
    for(j=0;j<JE;j++){
        for(i=0;i<IE;i++){
            ex[i][j][k]=0.0;
            ey[i][j][k]=0.0;
            ez[i][j][k]=0.0;
            sx[i][j][k]=0.0;
            sy[i][j][k]=0.0;
            sz[i][j][k]=0.0;
            ix[i][j][k]=0.0;
            iy[i][j][k]=0.0;
            iz[i][j][k]=0.0;
            dx[i][j][k]=0.0;
            dy[i][j][k]=0.0;
            dz[i][j][k]=0.0;
            hx[i][j][k]=0.0;
            hy[i][j][k]=0.0;
            hz[i][j][k]=0.0;
            gax[i][j][k]=1.0;
            gay[i][j][k]=1.0;
            gaz[i][j][k]=1.0;
            gbx[i][j][k]=0.0;
            gby[i][j][k]=0.0;
            gbz[i][j][k]=0.0;
            idyl[i][j][k]=0.0;
            ihyl[i][j][k]=0.0;
        }}}
for(i=0;i<ia;i++){
    for(j=0;j<JE;j++){
        for(k=0;k<KE;k++){
            idxl[i][j][k]=0.0;
            idxh[i][j][k]=0.0;
            ihxl[i][j][k]=0.0;
            ihxh[i][j][k]=0.0;
        }}}
for(i=0;i<IE;i++){
    for(j=0;j<JE;j++){
        for(k=0;k<ka;k++){
            idzl[i][j][k]=0.0;   
            idzh[i][j][k]=0.0;
            ihzl[i][j][k]=0.0;
            ihzh[i][j][k]=0.0;
        }}}
for(i=0;i<IE;i++){
    gi1[i]=0.0;
    fi1[i]=0.0;
    gi2[i]=1.0;
    fi2[i]=1.0;
    gi3[i]=1.0;
    fi3[i]=1.0;
}
for(j=0;j<JE;j++){
    gj1[j]=0.0;
    fj1[j]=0.0;
    gj2[j]=1.0;
    fj2[j]=1.0;
    gj3[j]=1.0;
    fj3[j]=1.0;
}
for(k=0;k<KE;k++){
    gk1[k]=0.0;
    fk1[k]=0.0;
    gk2[k]=1.0;
    fk2[k]=1.0;
    gk3[k]=1.0;
    fk3[k]=1.0;
}
for(n=0;n<NFREQS;n++){
    real_in[n]=0.0;
    imag_in[n]=0.0;
    amp_in[n]=0.0;
    phase_in[n]=0.0;
    for(j=0;j<JE;j++){
        for(i=0;i<IE;i++){
            real_pt[n][i][j]=0.0;
            imag_pt[n][i][j]=0.0;
            real_pt_sx[n][i][j]=0.0;
            imag_pt_sx[n][i][j]=0.0;
            real_pt_sy[n][i][j]=0.0;
            imag_pt_sy[n][i][j]=0.0;
            real_pt_sz[n][i][j]=0.0;
            imag_pt_sz[n][i][j]=0.0;
        }}}
lambdai=200e-9;
lambdaf=700e-9;
freqi=2e8/lambdai;
freqf=7e8/lambdaf;
for(i=0;i<NFREQS;i++){
    freq[i]=freqi-i*(freqi-freqf)/NFREQS;
}
for(n=0;n<NFREQS;n++){
    arg[n]=2*pi*freq[n]*dt;
}
//////////////////////////antena//////////////////////////
t0=20.0;
spread =6.0;
T=0.0;
cout<<"enter time -->";
cin>>nsteps;
ke=0;
////////////*******************FDTD*******************///////////////////////
for (n=1;n<=nsteps;n++){
    T=T+1;
    percent=T/nsteps*100;
    if (percent>=ke+1) {
        cout<<percent<<"%"<<endl;
    }
    ke=percent;

    /////////////////////source////////////////////
    pulse=exp(-.5*(pow((t0-T)/spread,2.0)));
    ez_inc[3]=pulse;
    //////////////////ABC//////////////////////////////////
    ez_inc[0]=ez_low_m2;
    ez_low_m2=ez_low_m1;
    ez_low_m1=ez_inc[1];
    ez_inc[JE-1]=ez_high_m2;
    ez_high_m2=ez_high_m1;
    ez_high_m1=ez_inc[JE-2];
    //////////////////////////update dx///////////////////////////////
    for(i=1;i<ia;i++){
        for(j=1;j<JE;j++){
            for(k=1;k<KE;k++){
                curl_h=(hz[i][j][k]-hz[i][j-1][k]-hy[i][j][k]+hy[i][j][k-1]);
                idxl[i][j][k]=idxl[i][j][k]+curl_h;
                dx[i][j][k]=x[i][j][k]+gj2[j]*gk2[k]*.5*(curl_h+gi1[i]*idxl[i][j][k]);
            }}}
    for(i=ia;i<=ib;i++){
        for(j=1;j<JE;j++){
            for(k=1;k<KE;k++){
                curl_h=(hz[i][j][k]-hz[i][j-1][k]-hy[i][j][k]+hy[i][j][k-1]);
                dx[i][j][k]=gj3[j]*gk3[k]*dx[i][j][k]+gj2[j]*gk2[k]*.5*curl_h;
            }}}
    for(i=ib+1;i<IE;i++){
        ixh=IE-i-1;
        for(j=1;j<JE;j++){
            for(k=1;k<KE;k++){
                curl_h=(hz[i][j][k]-hz[i][j-1][k]-hy[i][j][k]+hy[i][j][k-1]);
                idxh[ixh][j][k]=idxh[ixh][j][k]+curl_h;
                dx[i][j][k]=dx[i][j][k]+gj2[j]*gk2[k]*.5*(curl_h+gi1[i]*idxh[ixh][j][k]);
            }}}
    /////////////////////////update dy/////////////////////////////////
    for(i=1;i<IE;i++){
        for(j=1;j<JE;j++){
            for(k=1;k<KE;k++){
                curl_h=(hx[i][j][k]-hx[i][j][k-1]-hz[i][j][k]+hz[i-1][j][k]);
                idyl[i][j][k]=idyl[i][j][k]+curl_h;
                dy[i][j][k]=y[i][j][k]+gi2[i]*gk2[k]*.5*(curl_h+gj1[j]*idyl[i][j][k]);
            }}}

    ////////////////////////update dz///////////////////////////////////
    for(i=1;i<IE;i++){
        for(j=1;j<JE;j++){
            for(k=0;k<ka;k++){
                curl_h=(hy[i][j][k]-hy[i-1][j][k]-hx[i][j][k]+hx[i][j-1][k]);
                idzl[i][j][k]=idzl[i][j][k]+curl_h;
                dz[i][j][k]=gi3[i]*gj3[j]*dz[i][j][k]+gi2[i]*gj2[j]*.5*(curl_h+gk1[k]*idzl[i][j][k]);
            }}}
    for(i=1;i<IE;i++){
        for(j=1;j<JE;j++){
            for(k=ka;k<=kb;k++){
                curl_h=(hy[i][j][k]-hy[i-1][j][k]-hx[i][j][k]+hx[i][j-1][k]);
                dz[i][j][k]=gi3[i]*gj3[j];
            }}}
    for(i=1;i<IE;i++){
        for(j=1;j<JE;j++){
            for(k=kb+1;k<KE;k++){
                kzh=KE-k-1;
                curl_h=(hy[i][j][k]-hy[i-1][j][k]-hx[i][j][k]+hx[i][j-1][k]);
                idzh[i][j][kzh]=idzh[i][j][kzh]+curl_h;
                dz[i][j][k]=gi3[i]*gj3[j]*dz[i][j][k]+
            }}}
    ///////////////////////////////<<<<<<<<<incident dz>>>>>>>>>>>>>////////////////////////////
    ////////////////////////////source///////////////////////////
     pulse=exp(-.5*(pow((t0-T)/spread,2.0)));
     dz[ic][jc][kc]=pulse;

    ////////////////////////////////{{{{{{{{{{{update ex,ey,ez}}}}}}}}}}}}///////////////////////
    for(i=1;i<IE-1;i++){
        for(j=1;j<JE-1;j++){
            for(k=1;k<KE-1;k++){
                ex[i][j][k]=gax[i][j][k]*(dx[i][j][k]-ix[i][j][k]);
                ix[i][j][k]=ix[i][j][k]+gbx[i][j][k]*ex[i][j][k];
                ey[i][j][k]=gay[i][j][k]*(dy[i][j][k]-iy[i][j][k]);
                iy[i][j][k]=iy[i][j][k]+gby[i][j][k]*ey[i][j][k];
                ez[i][j][k]=gaz[i][j][k]*(dz[i][j][k]-iz[i][j][k]);
                iz[i][j][k]=iz[i][j][k]+gbz[i][j][k]*ez[i][j][k];
            }}}
    //////////////////////////////update hx///////////////////////////////////
    for(i=0;i<ia;i++){
        for(j=0;j<JE-1;j++){
            for(k=0;k<KE-1;k++){
                curl_e=(ey[i][j][k+1]-ey[i][j][k]-ez[i][j+1][k]+ez[i][j][k]);
                ihxl[i][j][k]=ihxl[i][j][k]+curl_e;
                hx[i][j][k]=fj3[j]*fk3[k]*hx[i][j][k]+fj2[j]*fk2[k]*.5*(curl_e+fi1[i]*ihxl[i][j][k]);
            }}}
    for(i=ia;i<=ib;i++){
        for(j=0;j<JE-1;j++){
            for(k=0;k<KE-1;k++){
                curl_e=(ey[i][j][k+1]-ey[i][j][k]-ez[i][j+1][k]+ez[i][j][k]);
                hx[i][j][k]=.5*curl_e;
            }}}
    for(i=ib+1;i<IE;i++){
        ixh=IE-i-1;
        for(j=0;j<JE-1;j++){
            for(k=0;k<KE-1;k++){
                curl_e=(ey[i][j][k+1]-ey[i][j][k]-ez[i][j+1][k]+ez[i][j][k]);
                ihxh[ixh][j][k]=ihxh[ixh][j][k]+curl_e;
     fj3[j]*fk3[k]*hx[i][j][k]+fj2[j]*fk2[k]*.5*(curl_e+fi1[i]*ihxh[ixh][j][k]);
            }}}
    /////////////////////////////incident hx////////////////////////////////////
    /*for(i=ia;i<=ib;i++){
        for(k=ka;k<=kb;k++){
            hx[i][ja-1][k]=hx[i][ja-1][k]+.5*ez_inc[ja];
            hx[i][jb][k]=hx[i][jb][k]-.5*ez_inc[jb];
        }
    }*/
    //////////////////////////////update hy/////////////////////////////////////
    for(i=0;i<IE-1;i++){
        for(j=0;j<JE;j++){
            for(k=0;k<KE-1;k++){
                curl_e=(ez[i+1][j][k]-ez[i][j][k]-ex[i][j][k+1]+ex[i][j][k]);
                ihyl[i][j][k]=ihyl[i][j][k]+curl_e;
                hy[i][j][k]=fi3[i]*fk3[k]*hy[i][j][k]+fi2[i]*fk2[k]*.5*(curl_e+fj1[j]*ihyl[i][j][k]);
            }}}
    /*for(i=0;i<IE-1;i++){
        for(j=ja;j<=jb;j++){
            for(k=0;k<KE-1;k++){
                curl_e=(ez[i+1][j][k]-ez[i][j][k]-ex[i][j][k+1]+ex[i][j][k]);
                hy[i][j][k]=fi3[i]*fk3[k]*hy[i][j][k]+fi2[i]*fk2[k]*.5*curl_e;
            }}}
    for(i=0;i<IE-1;i++){
        for(j=jb+1;j<JE;j++){
            jyh=JE-j-1;
            for(k=0;k<KE-1;k++){
                curl_e=(ez[i+1][j][k]-ez[i][j][k]-ex[i][j][k+1]+ex[i][j][k]);
                ihyh[i][jyh][k]=ihyh[i][jyh][k]+curl_e;
                hy[i][j][k]=fi3[i]*fk3[k]*hy[i][j][k]+fi2[i]*fk2[k]*.5*(curl_e+fj1[j]*ihyh[i][jyh][k]);
            }}}*/
    /////////////////////////////incident hy////////////////////////////////////
    /*for(j=ja;j<=jb;j++){
        for(k=ka;k<=kb;k++){
            hy[ia-1][j][k]=hy[ia-1][j][k]-.5*ez_inc[j];
            hy[ib][j][k]  =hy[ib][j][k]  +.5*ez_inc[j];
        }
    }*/
    /////////////////////////////update hz////////////////////////////////////////
    for(i=0;i<IE-1;i++){
        for(j=0;j<JE-1;j++){
            for(k=0;k<ka;k++){
                curl_e=(ex[i][j+1][k]-ex[i][j][k]-ey[i+1][j][k]+ey[i][j][k]);
                ihzl[i][j][k]=ihzl[i][j][k]+curl_e;
                hz[i][j][k]=fi3[i]*fj3[j]*hz[i][j][k]+fi2[i]*fj2[j]*.5*         (curl_e+fk1[k]*ihzl[i][j][k]);
            }}}
    for(i=0;i<IE-1;i++){
        for(j=0;j<JE-1;j++){
            for(k=ka;k<=kb;k++){
                curl_e=(ex[i][j+1][k]-ex[i][j][k]-ey[i+1][j][k]+ey[i][j][k]);
                hz[i][j][k]=fi3[i]*fj3[j]*hz[i][j][k]+fi2[i]*fj2[j]*.5*curl_e;
            }}}
    for(i=0;i<IE-1;i++){
        for(j=0;j<JE-1;j++){
            for(k=kb+1;k<KE;k++){
                kzh=KE-k-1;
                curl_e=(ex[i][j+1][k]-ex[i][j][k]-ey[i+1][j][k]+ey[i][j][k]);
                ihzh[i][j][kzh]=ihzh[i][j][kzh]+curl_e;
                hz[i][j][k]=fi3[i]*fj3[j]*hz[i][j][k]+fi2[i]*fj2[j]*.5*(curl_e+fk1[k]*ihzh[i][j][kzh]);
            }}} 
}//********************************END OF FDTD*******************************////////////


}

正如你所看到的,我的代码中充满了一些作用于一些矩阵的for循环,我怎么能忽略这个for循环呢?例如,我能用函数而不是为循环写函数吗。我真的很感激你的建议。

下载一个现有的库,如GLM、Eigen、Boost.UBLAS,它为您提供矩阵函数。

另外,使用std::vector

请首先使用单词表示变量,而不是单字母或双字母。我的直觉是,如果我遇到这样的代码,我会小心地将事物重命名为原来的样子。或者如果不可能的话,把它扔掉重写。

脑海中立刻浮现出三件事。

首先是更长、更具描述性的变量名和/或常量。

lmn这样的东西并不能告诉我关于这个变量的任何信息。

第二,更多的评论。每个做了重要事情的块都应该有一个注释,在高层说明它做了什么(而不是如何,除了非常普遍的情况外,我可以阅读代码来找到它)。

第三,我们并不缺乏函数原型,所以您不必将几乎每一条语句都放入main():-)中

如果一个函数占用了大约50行以上(例如),那么它可能太长了,应该调用更多的子函数)。

您的代码有多种想法。

我脑海中浮现的第一个问题是,你为什么不使用matlab进行这样的数学模拟?

第二,我是一名电气工程专业的学生,因此我对电磁波有相当的了解。我真的想,如果我看到模拟波浪的代码,我至少会识别出它的巨大轮廓。但我没有在你的代码中看到任何告诉我的东西。改进:

  • 更改变量名,有时使用长的变量名可能会很困难,但如果您想稍后重新阅读代码并进行更改,则完全值得
  • double ***gax=alloc3d(IE,JE,KE);你真的不应该使用三分球
  • 为了摆脱一些循环,编写为您执行循环的函数,给这些函数一个不同的名称,告诉它到底做什么。这将清理大部分代码