Created
October 22, 2025 19:46
-
-
Save Hermann-SW/57bb787c3a971fbaa57d96a196f1f534 to your computer and use it in GitHub Desktop.
work version for 4D in 3D display of n=p*q as sum of 4 squares
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| // https://hermann-sw.github.io/lattice_sphere_cmp/ | |
| // https://github.com/Hermann-sw/lattice_sphere_cmp | |
| // https://math.stackexchange.com/questions/4917740/which-faces-does-sphere-lattice-polyhedron-operatornamehullp-in-mathbbz | |
| // | |
| const jscad = require('@jscad/modeling') | |
| const { star } = jscad.primitives | |
| const { hull } = jscad.hulls | |
| const { geom3, poly3 } = jscad.geometries | |
| const { colorize } = jscad.colors | |
| const { cuboid, sphere, cylinder, circle } = jscad.primitives | |
| const { rotate, translate, translateX, scale:scale3d } = jscad.transforms | |
| const { angle, add, length, subtract, scale, dot } = jscad.maths.vec3 | |
| const { vec3, plane } = jscad.maths | |
| const { fromPoints } = jscad.maths.plane | |
| const { vectorText } = jscad.text | |
| const { union } = jscad.booleans | |
| const { extrudeLinear } = jscad.extrusions | |
| const { measureBoundingBox } = jscad.measurements | |
| const { hullChain } = jscad.hulls | |
| const det = ([[a,b,c],[d,e,f],[g,h,i]]) => a*(e*i-f*h)-b*(d*i-f*g)+c*(d*h-e*g) | |
| const rep = (A,b,i) => (C=[A[0].slice(),A[1].slice(),A[2].slice()], | |
| C[0][i]=b[0],C[1][i]=b[1],C[2][i]=b[2], C) | |
| const solve = (A,b) => | |
| (d=det(A), [det(rep(A,b,0))/d, det(rep(A,b,1))/d, det(rep(A,b,2))/d]) | |
| const gcd = (a,b) => b ? gcd(b, a % b) : Math.abs(a) | |
| const frac = (i,d) => (g=gcd(i,d),(d==g)?i/g:[i/g,d/g]) | |
| const solver = (A,b) => (d=det(A), | |
| [frac(det(rep(A,b,0)),d), frac(det(rep(A,b,1)),d), frac(det(rep(A,b,2)),d)]) | |
| function txt(mesg, w) { | |
| const lineRadius = w / 8 // 2 | |
| const lineCorner = circle({ radius: lineRadius }) | |
| const lineSegmentPointArrays = vectorText({ x: 0, y: 0, height: 0.125, input: mesg }) | |
| const lineSegments = [] | |
| lineSegmentPointArrays.forEach(function(segmentPoints) { | |
| const corners = segmentPoints.map((point) => translate(point, lineCorner)) | |
| lineSegments.push(hullChain(corners)) | |
| }) | |
| const message2D = union(lineSegments) | |
| ret = extrudeLinear({ height: w }, message2D) | |
| return ret | |
| } | |
| let map = new Map() | |
| function txt2(mesg, w) { | |
| let xofs = 0 | |
| let ret = [] | |
| for(i in mesg){ | |
| if(map.get(mesg[i])===undefined) { | |
| map.set(mesg[i], txt(mesg[i], w)) | |
| } | |
| ret.push(translateX(xofs, map.get(mesg[i]))) | |
| xofs += measureBoundingBox(map.get(mesg[i]))[1][0] | |
| } | |
| return ret | |
| } | |
| function vtxt(p1, num) { | |
| str = num.toString() | |
| la1 = Math.atan2(p1[1], p1[0]) | |
| ph1 = Math.PI/2-Math.acos(p1[2]/vec3.length(p1)) | |
| lab = txt2(str, 0.05) | |
| return translate([0, 0, 0.2], | |
| rotate([0, 0, la1], | |
| rotate([0, -ph1, 0], | |
| translate([vec3.length(p1)+0.1, -measureBoundingBox(lab[lab.length-1])[1][0]/2], | |
| rotate([Math.PI/2, 0, Math.PI/2], | |
| colorize([1, 1, 1], | |
| lab | |
| ) | |
| ) | |
| ) | |
| ) | |
| ) | |
| ) | |
| } | |
| segments = 10 | |
| const oneCylinder = cylinder({radius: 1, height: 1, segments: segments}) | |
| let reusedsphere = undefined | |
| function fastvertex(c) { return(translate(c, reusedsphere)) } | |
| function edge(_v, _w, r = 0.05, segs = segments) { | |
| d = [0, 0, 0] | |
| w = [0, 0, 0] | |
| subtract(d, _w, _v) | |
| add(w, _v, _w) | |
| scale(w, w, 0.5) | |
| return translate(w, | |
| rotate([0, Math.acos(d[2]/length(d)), Math.atan2(d[1], d[0])], | |
| scale3d([r,r,length(d)],oneCylinder) | |
| ) | |
| ) | |
| } | |
| let nsqrt = 0 | |
| const isPointOnPlane = (pla, pt, tol=1e-10) => | |
| Math.abs(vec3.dot(pt, pla) - pla[3]) < tol | |
| const arePointsOnPlane = (pla, pts, tol=1e-10) => | |
| pts.every((p) => isPointOnPlane(pla, p, tol)) | |
| // issue #1347 workaround: works for all points in same plane as well | |
| function fromPointsConvex(pts) { | |
| pla = plane.fromPoints(plane.create(), ...pts) | |
| if (!arePointsOnPlane(pla, pts)) { | |
| return geom3.fromPointsConvex(pts) | |
| } | |
| np = vec3.subtract(vec3.create(), | |
| pts[0], | |
| vec3.scale(vec3.create(), pla, 1e-3) | |
| ) | |
| g = geom3.fromPointsConvex([np, ...pts]) | |
| p3 = g.polygons.reduce((a,p) => | |
| (p!==undefined&&!p.vertices.some((v)=>vec3.equals(v,np)))?a=p:a,{}) | |
| g.polygons = [p3, poly3.invert(p3)] | |
| return g | |
| } | |
| const assert = (b) => {if(!b){throw("assert")}} | |
| var x=0 | |
| var sc = Math.sqrt(65-x*x) | |
| var er = 0.04 | |
| function cart2pol(p, f=sc) { | |
| return [Math.atan2(p[1],p[0]), Math.acos(p[2]/f)] | |
| } | |
| function edge2(_p1, _p2) { | |
| const v = _p1 // map3D(coords[_p1][0], coords[_p1][1]) | |
| const w = _p2 // map3D(coords[_p2][0], coords[_p2][1]) | |
| const p1 = cart2pol(v) | |
| const p2 = cart2pol(w) | |
| // al/la/ph: alpha/lambda/phi | lxy/sxy: delta lambda_xy/sigma_xy | |
| // https://en.wikipedia.org/wiki/Great-circle_navigation#Course | |
| const la1 = p1[0] | |
| const la2 = p2[0] | |
| const l12 = la2 - la1 | |
| const ph1 = Math.PI/2 - p1[1] | |
| const ph2 = Math.PI/2 - p2[1] | |
| const al1 = Math.atan2(Math.cos(ph2)*Math.sin(l12), Math.cos(ph1)*Math.sin(ph2)-Math.sin(ph1)*Math.cos(ph2)*Math.cos(l12)) | |
| // delta sigma_12 | |
| // https://en.wikipedia.org/wiki/Great-circle_distance#Formulae | |
| const s12 = Math.acos(Math.sin(ph1)*Math.sin(ph2)+Math.cos(ph1)*Math.cos(ph2)*Math.cos(l12)) | |
| return rotate([0, -ph1, la1], | |
| rotate([Math.PI/2-al1, 0, 0], | |
| makeArc2(sc, s12, 64) | |
| ) | |
| ) | |
| } | |
| // cool "makeArc2()" replaces "extrudeRotate()", from [email protected] | |
| // | |
| let cachedCylinder = rotate([-Math.PI/2,0,0],cylinder({radius:er, height:1, center:[0,0,0.5]})) | |
| function makeArc2(radius, angle, segments=64) { | |
| let correction = 0 | |
| if(angle < 0){ | |
| correction = angle | |
| angle *= -1 | |
| } | |
| // match how jscad calculates segments | |
| let stepA = Math.PI/(segments)*2 | |
| let steps = Math.ceil(angle / stepA) | |
| stepA = angle / steps | |
| let offset = 0, next=0 | |
| let out = [] | |
| while(offset < angle){ | |
| next += stepA | |
| if(next > angle) next=angle | |
| let len = (next-offset) * radius | |
| let x = Math.cos(offset) * radius | |
| let y = Math.sin(offset) * radius | |
| let part = colorize([0,0,1],translate([x,y,0],rotate([0,0,(next-(next-offset)/2)], | |
| jscad.transforms.scale([1,len,1], cachedCylinder) | |
| ))) | |
| out.push( correction ? rotate([0,0,correction], part):part) | |
| offset = next | |
| } | |
| if(!out.length) { throw("makeArc2(",radius,",",angle,",",segments,") problem") } | |
| return out | |
| } | |
| function main(params) { | |
| n = params.p*params.q | |
| nsqrt = Math.sqrt(n) | |
| reusedsphere = sphere({radius:0.014*nsqrt, segments: segments}) | |
| m = Math.floor(Math.sqrt(n)) | |
| out = [] | |
| edges = [] | |
| vertices = [] | |
| outside = [] | |
| for(k=0;k<=m;++k) | |
| { | |
| for(x=-m;x<=m;++x) | |
| for(y=-m;y<=m;++y) | |
| for(z=-m;z<=m;++z) | |
| if (x*x+y*y+z*z+k*k==n) | |
| out.push([x,y,z,k]) | |
| h = geom3.fromPointsConvex(out) | |
| out=[] | |
| sc=Math.sqrt(n-k*k) | |
| h.polygons.forEach((p)=>{ | |
| vs = p.vertices | |
| for(i=0;i<vs.length;i++){ | |
| edges.push(edge2(vs[i], vs[(i+1)%vs.length])) | |
| vertices.push(fastvertex(vs[i])) | |
| } | |
| if(params.text===true) | |
| vs.forEach((v)=>outside.push(vtxt(v,JSON.stringify(v)))) | |
| }) | |
| } | |
| return [edges, vertices, outside] | |
| } | |
| function getParameterDefinitions() { | |
| return [ | |
| { name: 'p', type: 'choice', values: [5, 13, 17, 29, 37, 41, 53, 61, 73, 89, 97], initial: 5, caption: 'p' }, | |
| { name: 'q', type: 'choice', values: [5, 13, 17, 29, 37, 41, 53, 61, 73, 89, 97], initial: 13, caption: 'q' }, | |
| { name: 'text', type: 'checkbox', checked: false, caption: 'text:' }, | |
| ] | |
| } | |
| module.exports = { main, getParameterDefinitions } |
Author
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Trying to map 4D into 3D — to be improved ...
https://jscad.app/#https://gist.github.com/Hermann-SW/57bb787c3a971fbaa57d96a196f1f534/raw/53804a9df05a7a7cdea3f2fb1b2ecc6f96205149/new.jscad