research(R3): floor ownership on an overlapping helical passage
Read-only study for L38: a synthetic FI Style 2 helix (3 m radius, 3.31 m pitch) shows Euclidean closest-centerline ownership is wrong when bore radius > pitch/2; a floor-near (+-0.125 m, within bore radius) nearest-floor-point owner gives 2.9 cm p-p on helix and style-like paths. The L38 0.94 m ripple period equals the 20-segment chord. L38's continuation applied it: floor ripple 0; wall ripple open. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_016dhxM43dnWvUu2B5Ph1iG1
This commit is contained in:
Fr0zkaandClaude Opus 5.5 committed
1 parent
7f2d84d2c7
commit
945225d942
3 files changed
+309
No files matched your search
@@ -0,0 +1,212 @@
|
||||
import math
|
||||
import bisect
|
||||
|
||||
DS = 0.125
|
||||
PITCH = 3.31 # L27 measured FI Style 2 main helix pitch (m/turn)
|
||||
HELIX_R = 3.0 # measured mean centerline radius (m)
|
||||
BORE_R = 4.5 # conservative authored near-mouth bore radius (m)
|
||||
LIFT = 0.15
|
||||
TURNS = 12.95
|
||||
SEGMENTS_PER_TURN = 20
|
||||
GRID = 0.125
|
||||
|
||||
|
||||
def snap(v):
|
||||
return round(v / GRID) * GRID
|
||||
|
||||
|
||||
def catmull(p0, p1, p2, p3, t):
|
||||
t2, t3 = t*t, t*t*t
|
||||
return 0.5*((2*p1)+(-p0+p2)*t+(2*p0-5*p1+4*p2-p3)*t2+(-p0+3*p1-3*p2+p3)*t3)
|
||||
|
||||
|
||||
def helix_run(radius=HELIX_R, pitch=PITCH, turns=TURNS, bore=BORE_R):
|
||||
n = math.ceil(turns * SEGMENTS_PER_TURN)
|
||||
theta_end = turns * 2*math.pi
|
||||
k = pitch/(2*math.pi)
|
||||
# At control i, z is the centerline height; its floor is z - bore + LIFT.
|
||||
pts = []
|
||||
for i in range(n+1):
|
||||
th = theta_end*i/n
|
||||
pts.append((radius*math.cos(th), radius*math.sin(th), k*th))
|
||||
length_per_rad = math.hypot(radius, k)
|
||||
count = math.ceil(theta_end*length_per_rad/DS)
|
||||
samples = []
|
||||
for j in range(count+1):
|
||||
th = theta_end*j/count
|
||||
true_floor = k*th-bore+LIFT
|
||||
q = (snap(radius*math.cos(th)), snap(radius*math.sin(th)), snap(true_floor))
|
||||
samples.append((th, true_floor, q))
|
||||
return pts, samples, bore, (lambda th: th)
|
||||
|
||||
|
||||
def other_run(kind, bore=1.5):
|
||||
# Monotone-grade non-helix controls: straight, bounded worm, and scalloped cascade.
|
||||
length = 80.0
|
||||
n = math.ceil(length / 0.75)
|
||||
pts = []
|
||||
xy = []
|
||||
for i in range(n+1):
|
||||
s = length*i/n
|
||||
if kind == 'straight':
|
||||
x, y = s, 0.0
|
||||
elif kind == 'worm':
|
||||
x, y = s, 2.0*math.sin(s/5.0)+0.35*math.sin(s/1.7)
|
||||
elif kind == 'cascading':
|
||||
th = s/5.0
|
||||
r = 10.0-1.3*abs(math.sin(4.0*th))
|
||||
x, y = r*math.cos(th), r*math.sin(th)
|
||||
else:
|
||||
raise ValueError(kind)
|
||||
xy.append((x,y))
|
||||
arc=[0.0]
|
||||
for i in range(1,len(xy)):
|
||||
arc.append(arc[-1]+math.hypot(xy[i][0]-xy[i-1][0],xy[i][1]-xy[i-1][1]))
|
||||
for i,(x,y) in enumerate(xy):
|
||||
pts.append((x,y,-0.19*arc[i]))
|
||||
samples=[]
|
||||
total=arc[-1]
|
||||
count=math.ceil(total/DS)
|
||||
for j in range(count+1):
|
||||
s=total*j/count
|
||||
i=0
|
||||
while i+1<len(arc) and arc[i+1]<s: i+=1
|
||||
if i+1>=len(arc): i=len(arc)-2
|
||||
t=(s-arc[i])/(arc[i+1]-arc[i])
|
||||
x=pts[i][0]*(1-t)+pts[i+1][0]*t
|
||||
y=pts[i][1]*(1-t)+pts[i+1][1]*t
|
||||
z=pts[i][2]*(1-t)+pts[i+1][2]*t
|
||||
f=z-bore+LIFT
|
||||
samples.append((s,f,(snap(x),snap(y),snap(f))))
|
||||
return pts,samples,bore,None
|
||||
|
||||
|
||||
def project_candidates(pts, q, bore):
|
||||
qx,qy,qz=q
|
||||
out=[]
|
||||
for i in range(len(pts)-1):
|
||||
a,b=pts[i],pts[i+1]
|
||||
dx,dy=b[0]-a[0],b[1]-a[1]
|
||||
h2=dx*dx+dy*dy
|
||||
if h2<1e-12:
|
||||
dz=b[2]-a[2]
|
||||
den=dx*dx+dy*dy+dz*dz
|
||||
t=max(0.0,min(1.0,((qx-a[0])*dx+(qy-a[1])*dy+(qz-a[2])*dz)/den))
|
||||
else:
|
||||
t=max(0.0,min(1.0,((qx-a[0])*dx+(qy-a[1])*dy)/h2))
|
||||
cx=a[0]+t*dx
|
||||
cy=a[1]+t*dy
|
||||
# Match ShapeWorldInternal.h: smooth only the lateral coordinates with Catmull-Rom.
|
||||
if i>2 and i+12<len(pts)-1:
|
||||
p0,p1,p2,p3=pts[i-1],a,b,pts[i+2]
|
||||
cx=catmull(p0[0],p1[0],p2[0],p3[0],t)
|
||||
cy=catmull(p0[1],p1[1],p2[1],p3[1],t)
|
||||
cz=a[2]+t*(b[2]-a[2])
|
||||
floor=cz-bore+LIFT
|
||||
hx,hy=qx-cx,qy-cy
|
||||
h2=hx*hx+hy*hy
|
||||
fd2=h2+(qz-floor)*(qz-floor)
|
||||
cd2=h2+(qz-cz)*(qz-cz)
|
||||
out.append((floor,fd2,cd2,h2,i,t,cx,cy))
|
||||
return out
|
||||
|
||||
|
||||
def metrics(label, pts, samples, bore, helix_info=None):
|
||||
names=['nearest centerline','nearest floor point','below+within-radius (strict)',
|
||||
'below+within-radius (+0.125m)','topmost floor below',
|
||||
'near-floor band (±0.125m)',
|
||||
'height-bin local segment (±1)','angle+height']
|
||||
if helix_info:
|
||||
names.insert(1,'analytic closest helix centerline')
|
||||
errors={n:[] for n in names}
|
||||
owners={n:[] for n in names}
|
||||
for loc,true_floor,q in samples:
|
||||
c=project_candidates(pts,q,bore)
|
||||
# Stable tie: lower segment index after score.
|
||||
a=min(c,key=lambda x:(x[2],x[4]))
|
||||
b=min(c,key=lambda x:(x[1],x[4]))
|
||||
filt=[x for x in c if x[0]<=q[2]+1e-12 and math.sqrt(x[3])<=bore]
|
||||
fc=min(filt,key=lambda x:(x[1],x[4])) if filt else b
|
||||
filt_tol=[x for x in c if x[0]<=q[2]+0.125+1e-12 and math.sqrt(x[3])<=bore]
|
||||
ft=min(filt_tol,key=lambda x:(x[1],x[4])) if filt_tol else b
|
||||
below=[x for x in c if x[0]<=q[2]+1e-12 and math.sqrt(x[3])<=bore]
|
||||
top=max(below,key=lambda x:(x[0],-x[4])) if below else b
|
||||
near=[x for x in c if abs(x[0]-q[2])<=0.125+1e-12 and math.sqrt(x[3])<=bore]
|
||||
nf=min(near,key=lambda x:(x[1],x[4])) if near else b
|
||||
floors=[p[2]-bore+LIFT for p in pts]
|
||||
if floors[0] < floors[-1]:
|
||||
i=max(0,min(len(pts)-2,bisect.bisect_right(floors,q[2])-1))
|
||||
else:
|
||||
neg=[-v for v in floors]
|
||||
i=max(0,min(len(pts)-2,bisect.bisect_right(neg,-q[2])-1))
|
||||
local=[x for x in c if abs(x[4]-i)<=1]
|
||||
hb=min(local,key=lambda x:(x[1],x[4])) if local else b
|
||||
if helix_info:
|
||||
radius,pitch,theta_end=helix_info
|
||||
phi=math.atan2(q[1],q[0])
|
||||
# Select the angle-aligned turn whose analytic floor is closest to voxel height.
|
||||
theta_z=(q[2]+bore-LIFT)/(pitch/(2*math.pi))
|
||||
k0=round((theta_z-phi)/(2*math.pi))
|
||||
th=max(0.0,min(theta_end,phi+2*math.pi*k0))
|
||||
hfloor=(pitch/(2*math.pi))*th-bore+LIFT
|
||||
howner=k0
|
||||
else:
|
||||
hfloor=b[0]
|
||||
howner=b[4]
|
||||
selected=[(names[0],a[0],a[4])]
|
||||
if helix_info:
|
||||
radius,pitch,theta_end=helix_info
|
||||
qx,qy,qz=q
|
||||
def dist2(th):
|
||||
cx=radius*math.cos(th); cy=radius*math.sin(th); cz=pitch*th/(2*math.pi)
|
||||
return (qx-cx)**2+(qy-cy)**2+(qz-cz)**2
|
||||
step=math.pi/16.0
|
||||
count=math.ceil(theta_end/step)
|
||||
vals=[dist2(min(theta_end,i*step)) for i in range(count+1)]
|
||||
best=(vals[0],0.0)
|
||||
if (vals[-1],theta_end)<best: best=(vals[-1],theta_end)
|
||||
gr=(math.sqrt(5.0)-1.0)/2.0
|
||||
for i in range(1,count):
|
||||
if vals[i]>vals[i-1] or vals[i]>vals[i+1]: continue
|
||||
lo=max(0.0,(i-1)*step); hi=min(theta_end,(i+1)*step)
|
||||
x1=hi-gr*(hi-lo); x2=lo+gr*(hi-lo)
|
||||
f1,f2=dist2(x1),dist2(x2)
|
||||
for _ in range(30):
|
||||
if f1<=f2:
|
||||
hi,x2,f2=x2,x1,f1; x1=hi-gr*(hi-lo); f1=dist2(x1)
|
||||
else:
|
||||
lo,x1,f1=x1,x2,f2; x2=lo+gr*(hi-lo); f2=dist2(x2)
|
||||
th=(lo+hi)*0.5; value=dist2(th)
|
||||
if value<best[0]: best=(value,th)
|
||||
exact_floor=pitch*best[1]/(2*math.pi)-bore+LIFT
|
||||
selected.append((names[1],exact_floor,best[1]))
|
||||
selected += [(names[2],b[0],b[4]),(names[3],fc[0],fc[4]),
|
||||
(names[4],ft[0],ft[4]),(names[5],top[0],top[4]),
|
||||
(names[6],nf[0],nf[4]),(names[7],hb[0],hb[4]),
|
||||
(names[8],hfloor,howner)]
|
||||
else:
|
||||
selected += [(names[1],b[0],b[4]),(names[2],fc[0],fc[4]),
|
||||
(names[3],ft[0],ft[4]),(names[4],top[0],top[4]),
|
||||
(names[5],nf[0],nf[4]),(names[6],hb[0],hb[4]),
|
||||
(names[7],hfloor,howner)]
|
||||
for name,f,owner in selected:
|
||||
errors[name].append(f-true_floor)
|
||||
owners[name].append(owner)
|
||||
print(f'[{label}] probes={len(samples)} grid={GRID:.3f} m bore_radius={bore:.2f} m')
|
||||
for name in names:
|
||||
es=errors[name]
|
||||
pp=max(es)-min(es)
|
||||
rms=math.sqrt(sum(e*e for e in es)/len(es))
|
||||
mae=sum(abs(e) for e in es)/len(es)
|
||||
flips=sum(owners[name][i]!=owners[name][i-1] for i in range(1,len(owners[name])))
|
||||
print(f' {name:39s} ripple_pp={pp:8.4f} max_abs={max(map(abs,es)):8.4f} rms={rms:8.4f} owner_steps={flips}')
|
||||
|
||||
|
||||
pts,samples,bore,_=helix_run()
|
||||
metrics('FI-style2 helix',pts,samples,bore,(HELIX_R,PITCH,TURNS*2*math.pi))
|
||||
# Mid-run bore from the published route is 3 m; this checks the mouth-radius stress case above.
|
||||
pts3,samples3,bore3,_=helix_run(bore=3.0)
|
||||
metrics('FI-style2 mid-bore control',pts3,samples3,bore3,(HELIX_R,PITCH,TURNS*2*math.pi))
|
||||
for style in ('straight','worm','cascading'):
|
||||
p,s,b,_=other_run(style)
|
||||
metrics(style,p,s,b)
|
||||
Reference in new issue
Block a user