python3 <<'PY' import math,random a=0.25; b=1.0/12 def BC(y1,y2): lo=0.0; hi=a; gr=(math.sqrt(5)-1)/2 def f(y3): beta=a+b-y3; B=math.sqrt(beta)*(b+math.sqrt(max(0,a*(a-y3)))) gam=y1+y2+y3; C=math.sqrt(gam)*(math.sqrt(max(0,b*y1))+math.sqrt(max(0,a*y2))+math.sqrt(max(0,a*y3))) return B+C for _ in range(20): c=hi-gr*(hi-lo); d=lo+gr*(hi-lo) if f(c)mu2+1e-12 or mu1<1e-12: return -1 c=nu1-a*q2/mu1 if c<-1e-10 or c>q+1e-10 or c>nu1+1e-10 or nu1-c>q2+1e-10: return -1 c2=nu1-c; x2=mu1; z2=mu2-q return math.sqrt(max(0,b*q*c))+math.sqrt(max(0,a*q*(q-c)))+math.sqrt(max(0,a*x2*c2))+math.sqrt(max(0,b*x2*(q2-c2)))+math.sqrt(max(0,b*z2*c2))+math.sqrt(max(0,a*z2*(q2-c2))) def p_col1(y1,y2,q): al=2*a+2*b-y1-y2; mu1=a+b-y1; mu2=b+a-y2; nu1=a+b if q<=0 or q>=nu1: return -1 q2=al-q; nu3=al-nu1 den=nu3 r=(a*q2-mu1*(nu1-q))/den if r<-1e-10 or r>q+1e-10: return -1 c2=nu1-q; x2=mu1-r; z2=mu2-(q-r) return math.sqrt(max(0,a*q*r))+math.sqrt(max(0,b*q*(q-r)))+math.sqrt(max(0,a*x2*c2))+math.sqrt(max(0,b*x2*nu3))+math.sqrt(max(0,b*z2*c2))+math.sqrt(max(0,a*z2*nu3)) def p_col3(y1,y2,q): al=2*a+2*b-y1-y2; mu1=a+b-y1; mu2=b+a-y2; nu1=a+b; nu3=al-nu1 if q<=0 or q>nu3+1e-12 or nu1<1e-12: return -1 q2=al-q r=mu1-a*q2/nu1 if r<-1e-10 or r>q+1e-10: return -1 x2=mu1-r; z2=mu2-(q-r); c3=nu3-q return math.sqrt(max(0,b*q*r))+math.sqrt(max(0,a*q*(q-r)))+math.sqrt(max(0,a*x2*nu1))+math.sqrt(max(0,b*x2*c3))+math.sqrt(max(0,b*z2*nu1))+math.sqrt(max(0,a*z2*c3)) for name,fn,qf in [("row2",p_row2,"mu2"),("col1",p_col1,"nu1"),("col3",p_col3,"nu3")]: best=-1; arg=None for _ in range(25000): y1=random.random()*b; y2=random.random()*a al=2*a+2*b-y1-y2; mu2=b+a-y2; nu1=a+b; nu3=al-nu1 mx={"mu2":mu2,"nu1":nu1,"nu3":nu3}[qf] q=random.uniform(1e-9,max(1e-9,mx)) A=fn(y1,y2,q) if A<0: continue tot=A+BC(y1,y2) if tot>best: best=tot; arg=(y1,y2,q,A,tot) print(name,best,arg) PY