C语言——一些处理脚本(21)

# include <stdio.h>
# include <math.h>
# include <string.h>
# include <stdlib.h>

# define N   100
# define NN  200
# define NTH 181
int main(void)
{
  int i, j, k, l, m, zz=0;
  double a=4;
  double n, d, L;
  double sinp, cosp, sinth, costh;
  double sinw, cosw, cosps, om;
  double costh2,sinth2;
  double x[N], y[N], z[N];
  double xm1,xm2,ym1,ym2,zm1,zm2;
  double A[N],B,aa;
  double pi, th, phi, ps;
  double t[N],xc[N],yc[N],zc[N],r[N],e[N],f[N],g[N];
  double xn[NTH][N],yn[NTH][N],zn[NTH][N],xt[NTH][NN],yt[NTH][NN],zt[NTH][NN];
  double sig,sig2,fr,kcont=1.0;
  double sigsa=4.0;
  double rx, ry, rz, r2,r4,r6,rij,e3=0;
  double E,EE,Elow,thm,eh;
  char Istring[10];
  FILE *fp, *fp1;
  pi=acos(-1.);
  
  /*printf("th = ");
  scanf("%lf",&th);
  printf("phi = ");
  scanf("%lf",&phi);*/


  for(i=1;i<180;i++)
  {
    for(j=1;j<180;j++)
    {
      cosp=cos(j*pi/180);
      sinp=sin(j*pi/180);
 
      costh=cos(i*pi/180);
      sinth=sin(i*pi/180);
  
      costh2=cos(i*pi/2/180);
      sinth2=sin(i*pi/2/180);  

      y[4]=a*sinth2-a*costh*sinth2-a*cosp*sinth*costh2;
      x[4]=((-1)*a*cosp*sinth-y[4]*costh2)/sinth2;
      z[4]=sqrt(a*a-x[4]*x[4]-a*a*costh2*costh2-2*x[4]*a*costh2-y[4]*y[4]-a*a*sinth2*sinth2+2*y[4]*a*sinth2);

      x[0]=x[4]; y[0]=-y[4]; z[0]=-z[4];
      x[1]=(-1)*a*costh2; y[1]=(-1)*a*sinth2; z[1]=0;
      x[2]=0; y[2]=0; z[2]=0;
      x[3]=x[1]; y[3]=-y[1]; z[3]=0;

      B=sqrt( ((x[4]+x[2])/2-x[3])*((x[4]+x[2])/2-x[3])
         +((y[4]+y[2])/2-y[3])*((y[4]+y[2])/2-y[3])
         +((z[4]+z[2])/2-z[3])*((z[4]+z[2])/2-z[3]) );

      for(k=5;k<N;k++) 
      {
        x[k]=(x[k-1]+x[k-3]-2*x[k-2])*fabs(x[4])/B+2*x[k-2]-x[k-4];
        y[k]=(y[k-1]+y[k-3]-2*y[k-2])*fabs(x[4])/B+2*y[k-2]-y[k-4];
        z[k]=(z[k-1]+z[k-3]-2*z[k-2])*fabs(x[4])/B+2*z[k-2]-z[k-4];
      }

       sinw=(fabs(z[4]*sinth2))/(sqrt(z[4]*z[4]+(y[4]-2*a*sinth2)*(y[4]-2*a*sinth2)));
       cosw=sqrt(1-sinw*sinw);
       cosps=(costh+sinw*sinw)/(cosw*cosw);
       ps=acos(cosps);
       om=acos(cosw)/pi*180;

       n=2*pi/(pi-ps);
       L=a*sinw*n;
       d=a*cosw/sin(pi/n);

       xm1=(-1)*a*costh2/(2*sin(pi/n)*sin(pi/n));  ym1=z[4]/2;  zm1=(-1)*y[4]/2+a*sinth2;  
       xm2=(-1)*a*costh2/(2*sin(pi/n)*sin(pi/n));  ym2=0;  zm2=0;

       aa=sqrt((xm1-xm2)*(xm1-xm2)+(ym1-ym2)*(ym1-ym2)+(zm1-zm2)*(zm1-zm2));

       for(k=0;k<N;k++)
       {
         xn[0][k]=x[k];
         yn[0][k]=y[k];
         zn[0][k]=z[k];
       }

       for(k=0;k<N;k++)
       {
         t[k]=aa*((xm1-xm2)*(xn[0][k]-xm2)+(ym1-ym2)*(yn[0][k]-ym2)+(zm1-zm2)*(zn[0][k]-zm2))/((xm1-xm2)*(xm1-xm2)+(ym1-ym2)*(ym1-ym2)+(zm1-zm2)*(zm1-zm2));
         xc[k]=xm2+(xm1-xm2)*t[k]/aa;
         yc[k]=ym2+(ym1-ym2)*t[k]/aa;
         zc[k]=zm2+(zm1-zm2)*t[k]/aa;
         r[k]=sqrt( (xn[0][k]-xc[k])*(xn[0][k]-xc[k])+(yn[0][k]-yc[k])*(yn[0][k]-yc[k])+(zn[0][k]-zc[k])*(zn[0][k]-zc[k]) );
         e[k]=(ym1-ym2)*(zn[0][k]-zc[k])/aa/r[k]-(zm1-zm2)*(yn[0][k]-yc[k])/aa/r[k];
         f[k]=(zm1-zm2)*(xn[0][k]-xc[k])/aa/r[k]-(xm1-xm2)*(zn[0][k]-zc[k])/aa/r[k];
         g[k]=(xm1-xm2)*(yn[0][k]-yc[k])/aa/r[k]-(ym1-ym2)*(xn[0][k]-xc[k])/aa/r[k];
       }
       
       for(m=0;m<NTH;m++)
       {
         for(k=0;k<N;k++)
         {
           if(zc[k]>zm1)
           {  
             xn[m][k]=(xn[0][k]-xc[k])*cos(m*pi/180)+e[k]*r[k]*sin(m*pi/180)+xc[k];
             yn[m][k]=(yn[0][k]-yc[k])*cos(m*pi/180)+f[k]*r[k]*sin(m*pi/180)+yc[k];
             zn[m][k]=(zn[0][k]-zc[k])*cos(m*pi/180)+g[k]*r[k]*sin(m*pi/180)+zc[k];
           }
           else
           {
             xn[m][k]=(xn[0][k]-xc[k])*cos(m*pi/180)+e[k]*r[k]*sin(m*pi/180)+xc[k];
             yn[m][k]=(yn[0][k]-yc[k])*cos(m*pi/180)+f[k]*r[k]*sin(m*pi/180)+yc[k];
             zn[m][k]=(zn[0][k]-zc[k])*cos(m*pi/180)+g[k]*r[k]*sin(m*pi/180)+zc[k];
           }
         }
       }

       for(m=1;m<NTH;m++)
       {
         for(k=0;k<N;k++)
         {
           xt[m][k]=x[k];
           yt[m][k]=y[k];
           zt[m][k]=z[k];
         }
         for(k=N;k<NN;k++)
         {
           xt[m][k]=xn[m][NN-1-k];
           yt[m][k]=yn[m][NN-1-k];
           zt[m][k]=zn[m][NN-1-k];
         }
       }

       Elow=99999999;
       for(m=1;m<NTH;m++)
       {
         E=0; eh=0;
         for(k=0;k<NN;k++)
         {
           for(l=k+2;l<NN;l++)
           {
             rx=xt[m][l]-xt[m][k]; 
             ry=yt[m][l]-yt[m][k]; 
             rz=zt[m][l]-zt[m][k]; 
             r2=rx*rx+ry*ry+rz*rz;
             //if((r2>100)||(r2<1e-8)) continue;
             rij=sqrt(r2);
             sig2=sigsa*sigsa; 
             r6=sig2/r2; r4=r6*r6; r6=r4*r6;
             e3=kcont*r6*(r6-2.0);
             E+=e3;
           }
           if(k==99)
           {
             l=100;
             rx=xt[m][l]-xt[m][k]; 
             ry=yt[m][l]-yt[m][k]; 
             rz=zt[m][l]-zt[m][k]; 
             r2=rx*rx+ry*ry+rz*rz;
             //if((r2>100)||(r2<1e-8)) continue;
             rij=sqrt(r2);
             sig2=sigsa*sigsa; 
             r6=sig2/r2; r4=r6*r6; r6=r4*r6;
             e3=kcont*r6*(r6-2.0);
             E+=e3;
           }
         }
         if(m==180)
         {
           EE=E;
           for(k=0;k<NN;k++)
           {
             for(l=0;l<N;l++)
             {
               if(l!=k)
               {
                 rx=xt[m][l]-xt[m][k]; 
                 ry=yt[m][l]-yt[m][k]; 
                 rz=zt[m][l]-zt[m][k]; 
                 r2=rx*rx+ry*ry+rz*rz;
                 //if((r2>100)||(r2<1e-8)) continue;
                 rij=sqrt(r2);
                 sig2=sigsa*sigsa; 
                 r6=sig2/r2; r4=r6*r6; r6=r4*r6;
                 e3=kcont*r6*(r6-2.0);
                 if(e3>0)
                 {eh=e3;}
               }
             }
           }
         }

         if(E<Elow)
         {
           Elow=E;
           thm=m;
         }
       }
       fp=fopen("ld1","a");
       fprintf(fp,"%3d  %3d  %10g  %10g  %10g  %10g   %3.0f  %10g   180  %10g\n",i,j,n,om,L,d,thm,Elow,EE);
       fclose(fp);
       if(Elow<0)
       {
         fp=fopen("ld2","a");
         fprintf(fp,"%3d  %3d  %10g  %10g  %10g  %10g   %3.0f  %10g  %10g   180  %10g  %10g\n",i,j,n,om,L,d,thm,Elow,Elow/NN,EE,EE/NN);
         fclose(fp);
         if(thm<180)
         {
           fp=fopen("ld3","a");
           fprintf(fp,"%3d  %3d  %10g  %10g  %10g  %10g   %3.0f  %10g   180  %10g\n",i,j,n,om,L,d,thm,Elow,EE);
           fclose(fp);
         }
       } 
    }
  }
  
  return 0;
}

 1     2.01117  0.00436338  0.000612646     4.00015     0       1e+08   180  1.01501e+27
 2     2.01117  0.00872642  0.00122524     4.00015     0       1e+08   180  1.01501e+27
 3     2.01117   0.0130888  0.00183775     4.00015     0       1e+08   180  1.01501e+27
 4     2.01117   0.0174502  0.00245011     4.00015     0       1e+08   180  1.01501e+27
 5     2.01116   0.0218102  0.00306228     4.00015     0       1e+08   180  1.01501e+27
 6     2.01116   0.0261686  0.00367422     4.00015     0       1e+08   180  1.01501e+27
 7     2.01115    0.030525  0.00428587     4.00015     0       1e+08   180  1.01501e+27
 8     2.01115   0.0348791  0.00489718     4.00015     0       1e+08   180  1.01501e+27
 9     2.01114   0.0392305  0.00550812     4.00015     0       1e+08   180  1.01502e+27
10     2.01113    0.043579  0.00611863     4.00015     0       1e+08   180  1.01502e+27
11     2.01112   0.0479241  0.00672867     4.00015     0       1e+08   180  1.01502e+27
12     2.01111   0.0522655  0.00733819     4.00015     0       1e+08   180  1.01502e+27
13      2.0111    0.056603  0.00794714     4.00015     0       1e+08   180  1.01502e+27
14     2.01109   0.0609362  0.00855547     4.00015     0       1e+08   180  1.01502e+27
15     2.01108   0.0652647  0.00916314     4.00015     0       1e+08   180  1.01503e+27
16     2.01106   0.0695883   0.0097701     4.00015     0       1e+08   180  1.01503e+27
17     2.01105   0.0739065   0.0103763     4.00015     0       1e+08   180  1.01503e+27
18     2.01103   0.0782192   0.0109817     4.00014     0       1e+08   180  1.01503e+27
19     2.01102   0.0825258   0.0115863     4.00014     0       1e+08   180  1.01503e+27
20       2.011   0.0868262   0.0121899     4.00014     0       1e+08   180  1.01504e+27
21     2.01099     0.09112   0.0127926     4.00014     0       1e+08   180  1.01504e+27
22     2.01097   0.0954068   0.0133943     4.00014     0       1e+08   180  1.01504e+27
23     2.01095   0.0996864    0.013995     4.00014     0       1e+08   180  1.01505e+27
24     2.01093    0.103958   0.0145946     4.00014     0       1e+08   180  1.01505e+27
25     2.01091    0.108222   0.0151931     4.00014     0       1e+08   180  1.01505e+27
26     2.01089    0.112478   0.0157904     4.00014     0       1e+08   180  1.01506e+27
27     2.01086    0.116725   0.0163865     4.00014     0       1e+08   180  1.01506e+27
28     2.01084    0.120964   0.0169813     4.00013     0       1e+08   180  1.01506e+27
29     2.01082    0.125193   0.0175748     4.00013     0       1e+08   180  1.01507e+27
30     2.01079    0.129413   0.0181669     4.00013     0       1e+08   180  1.01507e+27
31     2.01076    0.133622   0.0187576     4.00013     0       1e+08   180  1.01507e+27
32     2.01074    0.137822   0.0193469     4.00013     0       1e+08   180  1.01508e+27
33     2.01071    0.142011   0.0199346     4.00013     0       1e+08   180  1.01508e+27
34     2.01068    0.146189   0.0205209     4.00013     0       1e+08   180  1.01509e+27
35     2.01065    0.150356   0.0211055     4.00012     0       1e+08   180  1.01509e+27
36     2.01062    0.154512   0.0216885     4.00012     0       1e+08   180  1.01509e+27
37     2.01059    0.158656   0.0222698     4.00012     0       1e+08   180  1.0151e+27
38     2.01056    0.162788   0.0228495     4.00012     0       1e+08   180  1.0151e+27
39     2.01053    0.166907   0.0234273     4.00012     0       1e+08   180  1.01511e+27
40      2.0105    0.171014   0.0240033     4.00012     0       1e+08   180  1.01511e+27
41     2.01046    0.175108   0.0245775     4.00011     0       1e+08   180  1.01512e+27
42     2.01043    0.179188   0.0251498     4.00011     0       1e+08   180  1.01512e+27
43     2.01039    0.183255   0.0257201     4.00011     0       1e+08   180  1.01513e+27
44     2.01036    0.187307   0.0262884     4.00011     0       1e+08   180  1.01513e+27
45     2.01032    0.191346   0.0268547     4.00011     0       1e+08   180  1.01514e+27
46     2.01028     0.19537   0.0274189     4.00011     0       1e+08   180  1.01514e+27
47     2.01024    0.199379    0.027981      4.0001     0       1e+08   180  1.01515e+27
48      2.0102    0.203373    0.028541      4.0001     0       1e+08   180  1.01515e+27
49     2.01016    0.207351   0.0290987      4.0001     0       1e+08   180  1.01516e+27
50     2.01012    0.211314   0.0296542      4.0001     0       1e+08   180  1.01516e+27
51     2.01008     0.21526   0.0302074      4.0001     0       1e+08   180  1.01517e+27
52     2.01004     0.21919   0.0307582     4.00009     0       1e+08   180  1.01517e+27
53     2.00999    0.223103   0.0313067     4.00009     0       1e+08   180  1.01518e+27
54     2.00995       0.227   0.0318527     4.00009     0       1e+08   180  1.01518e+27
55      2.0099    0.230879   0.0323963     4.00009     0       1e+08   180  1.01519e+27
56     2.00986     0.23474   0.0329374     4.00009     0       1e+08   180  1.01519e+27
57     2.00981    0.238584    0.033476     4.00008     0       1e+08   180  1.0152e+27
58     2.00977     0.24241   0.0340119     4.00008     0       1e+08   180  1.0152e+27
59     2.00972    0.246217   0.0345452     4.00008     0       1e+08   180  1.01521e+27
60     2.00967    0.250005   0.0350759     4.00008     0       1e+08   180  1.01521e+27
61     2.00962    0.253774   0.0356038     4.00007     0       1e+08   180  1.01522e+27
62     2.00957    0.257524    0.036129     4.00007     0       1e+08   180  1.01522e+27
63     2.00952    0.261254   0.0366514     4.00007     0       1e+08   180  1.01523e+27
64     2.00947    0.264964    0.037171     4.00007     0       1e+08   180  1.01523e+27
65     2.00942    0.268655   0.0376877     4.00006     0       1e+08   180  1.01524e+27
66     2.00936    0.272324   0.0382015     4.00006     0       1e+08   180  1.01524e+27
67     2.00931    0.275973   0.0387123     4.00006     0       1e+08   180  1.01525e+27
68     2.00925    0.279601   0.0392202     4.00006     0       1e+08   180  1.01525e+27
69      2.0092    0.283208    0.039725     4.00005     0       1e+08   180  1.01526e+27
70     2.00914    0.286793   0.0402268     4.00005     0       1e+08   180  1.01526e+27
71     2.00909    0.290356   0.0407254     4.00005     0       1e+08   180  1.01527e+27
72     2.00903    0.293898   0.0412209     4.00005     0       1e+08   180  1.01527e+27
73     2.00897    0.297416   0.0417132     4.00004     0       1e+08   180  1.01528e+27
74     2.00891    0.300912   0.0422023     4.00004     0       1e+08   180  1.01528e+27
75     2.00885    0.304386   0.0426882     4.00004     0       1e+08   180  1.01529e+27
76     2.00879    0.307836   0.0431707     4.00004     0       1e+08   180  1.01529e+27
77     2.00873    0.311262     0.04365     4.00003     0       1e+08   180  1.0153e+27
78     2.00867    0.314665   0.0441258     4.00003     0       1e+08   180  1.0153e+27
79     2.00861    0.318044   0.0445983     4.00003     0       1e+08   180  1.01531e+27
80     2.00855    0.321399   0.0450673     4.00003     0       1e+08   180  1.01531e+27
81     2.00848    0.324729   0.0455328     4.00002     0       1e+08   180  1.01532e+27
82     2.00842    0.328034   0.0459948     4.00002     0       1e+08   180  1.01532e+27
83     2.00836    0.331315   0.0464533     4.00002     0       1e+08   180  1.01533e+27
84     2.00829     0.33457   0.0469082     4.00002     0       1e+08   180  1.01533e+27
85     2.00823      0.3378   0.0473595     4.00001     0       1e+08   180  1.01534e+27
86     2.00816    0.341004   0.0478071     4.00001     0       1e+08   180  1.01534e+27
87     2.00809    0.344182    0.048251     4.00001     0       1e+08   180  1.01535e+27
88     2.00802    0.347334   0.0486912     4.00001     0       1e+08   180  1.01535e+27
89     2.00796    0.350459   0.0491277           4     0       1e+08   180  1.01535e+27
90     2.00789    0.353558   0.0495604           4     0       1e+08   180  1.01536e+27
91     2.00782     0.35663   0.0499892           4     0       1e+08   180  1.01536e+27
92     2.00775    0.359674   0.0504142     3.99999     0       1e+08   180  1.01537e+27
93     2.00768    0.362692   0.0508354     3.99999     0       1e+08   180  1.01537e+27
94     2.00761    0.365681   0.0512526     3.99999     0       1e+08   180  1.01537e+27
95     2.00753    0.368643   0.0516658     3.99999     0       1e+08   180  1.01538e+27
96     2.00746    0.371577   0.0520751     3.99998     0       1e+08   180  1.01538e+27
97     2.00739    0.374482   0.0524804     3.99998     0       1e+08   180  1.01538e+27
98     2.00732    0.377359   0.0528816     3.99998     0       1e+08   180  1.01539e+27
99     2.00724    0.380207   0.0532788     3.99998     0       1e+08   180  1.01539e+27
  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值