ｃｃ０２ - プログラム倉庫

#include <stdio.h>
#include <math.h>
int main(void)
{
double th[101],c[101],y[101];
double tl,d1,d2,c1,y1,l1,w1,tr;
double tr1,tr2,tr3;
double wel,maxwel,maxtl,maxtr,v1,vs,maxv;
int s,t,n;
int n1,n2,q,mx1,mx2,nx1,nx2;
double cp,lp,yp,h,u1;
double u[101][2][2],w[101][2][2],v[99][2][2][20];
double gotoc[101][2][2][20],gotoy[101][2][2][20],gotoq[101][2][2][20];
double endc[2][2],endy[2][2],endq[2][2],endv[2][2];
int m1,m2,qx;
int pp;
for (s=1;s<=100;s++){
th[s] = 0.02*s;
}
maxwel=-999;
for (n=10;n<=40;n++){
tl=0.01*n;
tr1=0.05;
tr2=0.03;
d1=0;
for (s=1;s<=100;s++){
w1=(1-tl)*th[s];
l1=(w1-tr1)/(2*w1);
if (l1<0) l1=0;
y1=th[s]*l1;
c1=w1*l1+tr1;
d1=d1+y1-c1;
}
t=0;
while(t<100){
d2=0;
for (s=1;s<=100;s++){
w1=(1-tl)*th[s];
l1=(w1-tr2)/(2*w1);
if (l1<0) l1=0;
y1=th[s]*l1;
c1=w1*l1+tr2;
d2=d2+y1-c1;
}
tr3=tr2-d2*(tr2-tr1)/(d2-d1);
tr1=tr2;
tr2=tr3;
d1=d2;
if (d2*d2<0.01) t=1000;
t=t+1;
}
tr=tr2;
wel=0;
for (s=1;s<=100;s++){
w1=(1-tl)*th[s];
l1=(w1-tr)/(2*w1);
if (l1<0) l1=0;
c1=w1*l1+tr;
wel=wel+log(c1)+log(1-l1);
}
if (wel>maxwel) maxtl=tl;
if (wel>maxwel) maxtr=tr;
if (wel>maxwel) maxwel=wel;
}
tl=maxtl;
tr=maxtr;
for (s=1;s<=100;s++){
w1=(1-tl)*th[s];
l1=(w1-tr)/(2*w1);
if (l1<0) l1=0;
y[s]=th[s]*l1;
c[s]=w1*l1+tr;
}
h=0.001;
for (s=1;s<=100;s++){
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
cp=c[s]+n1*h;
yp=y[s]+n2*h;
lp=yp/th[s];
pp=0;
if (lp<0) pp=1;
if (lp>0.99) pp=1;
if (cp<0.01) pp=1;
if (pp==1) lp=0.5;
if (pp==1) cp=0.5;
u[s][n1+1][n2+1]=log(cp)+log(1-lp);
if (pp==1) u[s][n1+1][n2+1]=-999;
}
}
}
for (s=1;s<=99;s++){
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
cp=c[s]+n1*h;
yp=y[s]+n2*h;
lp=yp/th[s+1];
pp=0;
if (lp<0) pp=1;
if (lp>0.99) pp=1;
if (cp<0.01) pp=1;
if (pp==1) lp=0.5;
if (pp==1) cp=0.5;
w[s][n1+1][n2+1]=log(cp)+log(1-lp);
if (pp==1) w[s][n1+1][n2+1]=-999;
}
}
}
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (q=-10;q<=10;q++){
v[1][n1+1][n2+1][q+10]=-999;
}
}
}
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
q=n2-n1;
v[1][n1+1][n2+1][q+10]=u[1][n1+1][n2+1];
}
}
for (s=2;s<=99;s++){
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (q=-10;q<=10;q++){
u1=u[s][n1+1][n2+1];
vs=-999;
qx=q+n1-n2;
pp=0;
if (qx<-10) pp=1;
if (qx>10) pp=1;
if (pp==1) qx=0;
for (m1=-1;m1<=1;m1++){
for (m2=-1;m2<=1;m2++){
v1=u1+v[s-1][m1+1][m2+1][qx+10];
if (w[s-1][m1][m2]>u1)v1=-999;
if (v1>vs) mx1=m1;
if (v1>vs) mx2=m2;
if (v1>vs) vs=v1;
}
}
if (pp==1) vs=-999;
v[s][n1+1][n2+1][q+10]=vs;
gotoc[s][n1+1][n2+1][q+10]=mx1;
gotoy[s][n1+1][n2+1][q+10]=mx2;
gotoq[s][n1+1][n2+1][q+10]=qx;
}
}
}
}
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
u1=u[100][n1+1][n2+1];
vs=-999;
qx=n1-n2;
for (m1=-1;m1<=1;m1++){
for (m2=-1;m2<=1;m2++){
v1=u1+v[99][m1+1][m2+1][qx+10];
if (w[99][m1][m2]>u1)v1=-999;
if (v1>vs) mx1=m1;
if (v1>vs) mx2=m2;
if (v1>vs) vs=v1;
}
}
endv[n1+1][n2+1]=vs;
endc[n1+1][n2+1]=mx1;
endy[n1+1][n2+1]=mx2;
endq[n1+1][n2+1]=qx;
}
}
vs=-999;
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
if (endv[n1+1][n2+1]>vs) nx1=n1;
if (endv[n1+1][n2+1]>vs) nx2=n2;
if (endv[n1+1][n2+1]>vs) vs=endv[n1+1][n2+1];
}
}
printf("%f",vs);
return 0;
}