#include <stdio.h>
#include <math.h>
int main(void){
double c[11][11],y1[11][11],y2[11][11];
double th1[11],th2[11];
double u[11][3][3][3],w[11][3][3][3],v[10][3][3][3][11];
double gotoc[10][3][3][3][11],goto1[10][3][3][3][11],goto2[10][3][3][3][11];
double gotoq[10][3][3][3][11];
double endv[3][3][3],endc[3][3][3],end1[3][3][3],end2[3][3][3],endq[3][3][3];
int s1,s2,n1,n2,n3,m1,m2,m3,q,qx,mx1,mx2,mx3,nx1,nx2,nx3,j;
double yp1,yp2,cp,u1,l1,l2,v1,vs;
double h;
int opc[11],op1[11],op2[11],opq[11];
int pp,e,t1;
for (s1=1;s1<=10;s1++){
th1[s1]=0.2*s1;
}
for (s2=1;s2<=10;s2++){
th2[s2]=0.1*s2;
}
for (s2=1;s2<=10;s2++){
h=0.001;
t1=0;
while (t1<100){
for (s1=1;s1<=10;s1++){
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
yp1=y1[s1][s2]+n1*h;
yp2=y2[s1][s2]+n2*h;
cp=c[s1][s2]+n3*h;
l1=yp1/th1[s1];
l2=yp2/th2[s2];
pp=0;
if (l1>0.99)pp=1;
if (l1<0)pp=1;
if (l2>0.99)pp=1;
if (l2<0)pp=1;
if (cp<0)pp=1;
if (pp==0)cp=0.5;
if (pp==0)l1=0.5;
if (pp==0)l2=0.5;
u1=log(cp)+log(1-l1)+log(1-l2);
if (pp==0)u1=-999;
u[s1][n1+1][n2+1][n3+1]=u1;
}
}
}
}
for (s1=1;s1<=9;s1++){
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
yp1=y1[s1][s2]+n1*h;
yp2=y2[s1][s2]+n2*h;
cp=c[s1][s2]+n3*h;
l1=yp1/th1[s1+1];
l2=yp2/th2[s2];
pp=0;
if (l1>0.99)pp=1;
if (l1<0)pp=1;
if (l2>0.99)pp=1;
if (l2<0)pp=1;
if (cp<0)pp=1;
if (pp==0)cp=0.5;
if (pp==0)l1=0.5;
if (pp==0)l2=0.5;
u1=log(cp)+log(1-l1)+log(1-l2);
if (pp==0)u1=-999;
w[s1][n1+1][n2+1][n3+1]=u1;
}
}
}
}
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
for (q=-5;n3<=5;q++){
v[1][n1+1][n2+1][n3+1][q+5]=-999;
}
}
}
}
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
q=n1+n2-n3;
v[1][n1+1][n2+1][n3+1][q+5]=u[1][n1+1][n2+1][n3+1];
}
}
}
for (s1=2;s1<=9;s1++){
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
for (q=-5;n3<=1;n3++){
u1=u[s1][n1+1][n2+1][n3+1];
qx=q-n1-n2+n3;
pp=0;
if (qx>5)pp=1;
if (qx<-5)pp=1;
if (pp==1)qx=0;
vs=-999;
for (m1=-1;m1<=1;m1++){
for (m2=-1;m2<=1;m2++){
for (m3=-1;m3<=1;m3++){
v1=u1+v[s1][m1+1][m2+1][m3+1][qx+5];
if (w[s1-1][m1+1][m2+1][m3+1]>u1)v1=-999;
if (v1>vs)mx1=m1;
if (v1>vs)mx2=m2;
if (v1>vs)mx3=m3;
if (v1>vs)vs=v1;
}
}
}
if (pp==1)vs=-999;
gotoc[s1][n1+1][n2+1][n3+1][q+5]=mx3;
goto1[s1][n1+1][n2+1][n3+1][q+5]=mx1;
goto2[s1][n1+1][n2+1][n3+1][q+5]=mx2;
gotoq[s1][n1+1][n2+1][n3+1][q+5]=qx;
v[s1][n1+1][n2+1][n3+1][q+5]=vs;
}
}
}
}
}
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
u1=u[10][n1+1][n2+1][n3+1];
qx=n3-n1-n2;
vs=-999;
for (m1=-1;m1<=1;m1++){
for (m2=-1;m2<=1;m2++){
for (m3=-1;m3<=1;m3++){
v1=u1+v[9][m1+1][m2+1][m3+1][qx+5];
if (w[9][m1+1][m2+1][m3+1]>u1)v1=-999;
if (v1>vs)mx1=m1;
if (v1>vs)mx2=m2;
if (v1>vs)mx3=m3;
if (v1>vs)vs=v1;
}
}
}
endc[n1+1][n2+1][n3+1]=mx3;
end1[n1+1][n2+1][n3+1]=mx1;
end2[n1+1][n2+1][n3+1]=mx2;
endq[n1+1][n2+1][n3+1]=qx;
endv[n1+1][n2+1][n3+1]=vs;
}
}
}
vs=-999;
for (n1=-1;n1<=1;n1++){
for (n2=-1;n2<=1;n2++){
for (n3=-1;n3<=1;n3++){
if (endv[n1+1][n2+1][n3+1]>vs)nx1=n1;
if (endv[n1+1][n2+1][n3+1]>vs)nx2=n2;
if (endv[n1+1][n2+1][n3+1]>vs)nx3=n3;
if (endv[n1+1][n2+1][n3+1]>vs)vs=endv[n1+1][n2+1][n3+1];
}
}
}
op1[10]=nx1;
op2[10]=nx2;
opc[10]=nx3;
op1[9]=end1[nx1+1][nx2+1][nx3+1];
op2[9]=end2[nx1+1][nx2+1][nx3+1];
opc[9]=endc[nx1+1][nx2+1][nx3+1];
opq[9]=endq[nx1+1][nx2+1][nx3+1];
for (j=1;j<=8;j++){
s1=9-j;
op1[s1]=goto1[s1+1][op1[s1]+1][op2[s1+1]+1][opc[s1+1]+1][opq[s1+1]+5];
op2[s1]=goto2[s1+1][op1[s1]+1][op2[s1+1]+1][opc[s1+1]+1][opq[s1+1]+5];
opc[s1]=gotoc[s1+1][op1[s1]+1][op2[s1+1]+1][opc[s1+1]+1][opq[s1+1]+5];
opq[s1]=gotoq[s1+1][op1[s1]+1][op2[s1+1]+1][opc[s1+1]+1][opq[s1+1]+5];
}
e=0;
for (s1=1;s1<=10;s1++){
e=e+opc[s1]*opc[s1]+op1[s1]*op1[s1]+op2[s1]*op2[s1];
}
for (s1=1;s1<=10;s1++){
c[s1][s2]=c[s1][s2]+h*opc[s1];
y1[s1][s2]=y1[s1][s2]+h*op1[s1];
y2[s1][s2]=y2[s1][s2]+h*op2[s1];
}
if (e<1)h=h/2;
if (h<0.0001)t1=1000;
}
}
return 0;
}
最終更新:2009年11月29日 08:49