python3 <<'PY'
import math,random
a=0.25; b=1.0/12
def best3(y1,y2,N=500000):
    alpha=2*a+2*b-y1-y2; mu1=a+b-y1; mu2=alpha-mu1; nu1=a+b
    best=-1
    for _ in range(N):
      q=[alpha*random.random() for i in range(3)]
      s=sum(q); q=[v/s*alpha for v in q]
      # random a_i sum mu1 with ai<=qi; use stick breaking scaled
      u=[random.random() for i in range(3)]; ss=sum(u); aa=[mu1*v/ss for v in u]
      if any(aa[i]>q[i]+1e-12 for i in range(3)): continue
      # choose c1 random in range, solve c2
      c1=random.random()*min(q[0],nu1)
      a3=aa[2]; q3=q[2]; a2=aa[1]; q2=q[1]; a1=aa[0]; q1=q[0]
      if abs(a2/q2-a3/q3)<1e-14: continue
      c2=(m11 if False else a - a1*c1/q1 - a3*(nu1-c1)/q3)/(a2/q2-a3/q3)
      c3=nu1-c1-c2
      if c2<-1e-11 or c2>q2+1e-11 or c3<-1e-11 or c3>q3+1e-11: continue
      cc=[c1,c2,c3]
      val=0.0
      for i in range(3):
        qi=q[i]; ai=aa[i]; ci=cc[i]
        bi=qi-ai; di=qi-ci
        val+=math.sqrt(max(0,a*ai*ci))+math.sqrt(max(0,b*ai*di))+math.sqrt(max(0,b*bi*ci))+math.sqrt(max(0,a*bi*di))
      if val>best: best=val
    return best
for y1,y2 in [(0.02,0.22),(0.0,0.2),(0.01,0.24),(0.04,0.23)]:
  print(y1,y2,best3(y1,y2,300000))
PY