python3 <<'PY' import math,random a=0.25; b=1.0/12 def G(y1,y2,q,r): alpha=2*a+2*b-y1-y2; mu1=a+b-y1; mu2=b+a-y2; nu1=a+b q2=alpha-q if y1<-1e-9 or y1>b+1e-9 or y2<-1e-9 or y2>a+1e-9 or q<=1e-12 or q2<=1e-12: return -1 lo=max(0.0,q-mu2); hi=min(q,mu1) if rhi+1e-9: return -1 den=r*q2/q-(mu1-r) if abs(den)<1e-14: return -1 c=(a*q2-(mu1-r)*nu1)/den if c<-1e-8 or c>q+1e-8 or c>nu1+1e-8 or nu1-c>q2+1e-8: return -1 r2=mu1-r; c2=nu1-c t1=math.sqrt(max(0,a*r*c))+math.sqrt(max(0,b*r*(q-c)))+math.sqrt(max(0,b*(q-r)*c))+math.sqrt(max(0,a*(q-r)*(q-c))) t2=math.sqrt(max(0,a*r2*c2))+math.sqrt(max(0,b*r2*(q2-c2)))+math.sqrt(max(0,b*(q2-r2)*c2))+math.sqrt(max(0,a*(q2-r2)*(q2-c2))) C=math.sqrt(y1+y2+a)*(math.sqrt(b*y1)+math.sqrt(a*y2)+a) return t1+t2+b**1.5+C t0=0.23045484476568612; y3=0.22112258331133355 y1c=3*t0-2/3 G0=G(y1c,t0,2*a+2*b-y1c-t0,(a+b-y1c)) print("G0",G0) eps=1e-6 bestder=-1e9; bestdir=None for _ in range(20000): a1=random.uniform(-20,20); a2=random.uniform(-20,20) d=random.random(); e=random.random() y1=y1c+a1*eps; y2=t0+a2*eps al=2*a+2*b-y1-y2; mu1=a+b-y1 q=al-eps r=mu1-d*eps g=G(y1,y2,q,r) if g<0: continue der=(g-G0)/eps if der>bestder: bestder=der; bestdir=(a1,a2,d,e,g) print("bestder",bestder,bestdir) PY