// Optimized Cloth Simulation
// Original by Reinder Nijhoff 2019
// Optimized version with:
// - Mutable vector operations (no array allocations)
// - Cached calculations (normals, distances)
// - Squared distance comparisons (avoiding sqrt)
// - Pre-allocated arrays
// - Reduced object creation

const turtle = new Slowpoke();
const cameraPos = [-7,3,8];
const cameraLookAt = [-0.1,-1.3,0];
const viewProjectionMatrix = setupCamera(cameraPos, cameraLookAt);
const polygons = new Polygons();
const shape = 1; // min=0, max=1, step=1 (Sphere, Box)
const simulationSteps = 175; // min=1, max=250, step=1

// Pre-allocated vectors for reuse
const tempVec1 = [0,0,0];
const tempVec2 = [0,0,0];
const tempVec3 = [0,0,0];

// Mutable vector operations (modify first argument in place)
function scale3m(out, a, b) { 
    out[0] = a[0]*b; out[1] = a[1]*b; out[2] = a[2]*b; 
    return out;
}
function add3m(out, a, b) { 
    out[0] = a[0]+b[0]; out[1] = a[1]+b[1]; out[2] = a[2]+b[2]; 
    return out;
}
function sub3m(out, a, b) { 
    out[0] = a[0]-b[0]; out[1] = a[1]-b[1]; out[2] = a[2]-b[2]; 
    return out;
}
function len3(a) { return Math.sqrt(a[0]*a[0]+a[1]*a[1]+a[2]*a[2]); }
function len3sq(a) { return a[0]*a[0]+a[1]*a[1]+a[2]*a[2]; }
function dot3(a,b) { return a[0]*b[0]+a[1]*b[1]+a[2]*b[2]; }
function normalize3m(out, a) { 
    const l = 1/Math.sqrt(a[0]*a[0]+a[1]*a[1]+a[2]*a[2]);
    out[0] = a[0]*l; out[1] = a[1]*l; out[2] = a[2]*l;
    return out;
}
function cross3m(out, a, b) {
    out[0] = a[1]*b[2]-a[2]*b[1];
    out[1] = a[2]*b[0]-a[0]*b[2];
    out[2] = a[0]*b[1]-a[1]*b[0];
    return out;
}

// Standard verlet integration with optimizations
class point {
    constructor(p, mass) {
        this.p = [...p];
        this.p_prev = [...p];
        this.mass = mass / dot3(p, p);
        this.acc = [0, -9.8, 0];
        this.links = [];
    }
    attach(p2, stiffness) {
        this.links.push(new link(this, p2, stiffness));
    }
    update(dt) {
        // Reuse temp vectors instead of allocating
        sub3m(tempVec1, this.p, this.p_prev); // vel
        scale3m(tempVec1, tempVec1, .9); // damping
        scale3m(tempVec2, this.acc, dt*dt/2);
        add3m(tempVec2, tempVec2, tempVec1);
        add3m(tempVec2, tempVec2, this.p);
        
        // Update positions
        this.p_prev[0] = this.p[0];
        this.p_prev[1] = this.p[1];
        this.p_prev[2] = this.p[2];
        this.p[0] = tempVec2[0];
        this.p[1] = tempVec2[1];
        this.p[2] = tempVec2[2];
    }
    solveCollisions() {
        if (shape == 0) { // collision with sphere
            if (len3sq(this.p) <= 1) {
                normalize3m(this.p, this.p_prev);
                this.p_prev[0] = this.p[0];
                this.p_prev[1] = this.p[1];
                this.p_prev[2] = this.p[2];
            }
        } else { // collision with block
            if (Math.abs(this.p[0]) < 1 && Math.abs(this.p[2]) < 1 && 
                this.p[1] <= 1 && this.p[1] >= .9) {
                const p = this.p_prev;
                this.p[0] = p[0];
                this.p[1] = Math.max(.9, Math.min(1, p[1]+.01));
                this.p[2] = p[2];
                this.p_prev[0] = this.p[0];
                this.p_prev[1] = this.p[1];
                this.p_prev[2] = this.p[2];
            }
        }
    }
    solve() {
        this.solveCollisions();
        for (let i = 0; i < this.links.length; i++) {
            this.links[i].solve();
        }
    }
}

class link {
    constructor(p1, p2, stiffness) {
        this.p1 = p1;
        this.p2 = p2;
        this.stiffness = stiffness;
        // Cache the initial distance
        sub3m(tempVec1, p1.p, p2.p);
        this.distance = len3(tempVec1);
        this.massRatio = (1/p1.mass) / (1/p1.mass + 1/p2.mass);
    }
    solve() {
        sub3m(tempVec1, this.p1.p, this.p2.p);
        const d = len3(tempVec1);
        const diff = (this.distance - d)/d;
        const stiffDiff = this.stiffness * diff;
        
        // Apply corrections directly
        const s = this.massRatio * stiffDiff;
        this.p1.p[0] += tempVec1[0] * s;
        this.p1.p[1] += tempVec1[1] * s;
        this.p1.p[2] += tempVec1[2] * s;
        
        const t = (1-this.massRatio) * stiffDiff;
        this.p2.p[0] -= tempVec1[0] * t;
        this.p2.p[1] -= tempVec1[1] * t;
        this.p2.p[2] -= tempVec1[2] * t;
    }
}

class world {
    constructor() {
        this.points = [];
    }
    addPoint(p) {
        this.points.push(p);
        return p;
    }
    update(dt, relax) {
        const pts = this.points;
        const len = pts.length;
        
        // Update all points
        for (let i = 0; i < len; i++) {
            pts[i].update(dt);
        }
        
        // Constraint solving iterations
        for (let j = 0; j < relax; j++) {
            for (let i = 0; i < len; i++) {
                pts[i].solve();
            }
        }
    }
}

// Cloth with optimizations
class cloth {
    constructor(world, w, h, s, center) {
        this.width = w;
        this.height = h;
        const ps = this.points = [];
        const tr = this.triangles = [];
        const trgrid = [];
        
        const mass = 1;
        
        // Create points
        for (let x = 0; x < w; x++) {
            ps[x] = [];
            for (let y = 0; y < h; y++) {
                const p = new point(
                    [(x-(w-1)/2)*s + center[0], center[1], (y-(h-1)/2)*s + center[2]], 
                    mass
                );
                ps[x][y] = world.addPoint(p);
            }
        }
        
        const stiffness = .95;
        const stiffness2 = .1;
        
        // Create links
        for (let x = 0; x < w; x++) {
            for (let y = 0; y < h; y++) {
                const p = ps[x][y];
                if (y < h-1) p.attach(ps[x][y+1], stiffness);
                if (x < w-1) p.attach(ps[x+1][y], stiffness);
                if (x < w-1 && y < h-1) p.attach(ps[x+1][y+1], stiffness);
                if (y < h-2) p.attach(ps[x][y+2], stiffness2);
                if (x < w-2) p.attach(ps[x+2][y], stiffness2);
                if (x > 1 && y < h-2) p.attach(ps[x-2][y+2], stiffness2);
            }
        }
        
        // Create triangles with neighbor references
        for (let x = 0; x < w-1; x++) {
            trgrid[x] = [];
            for (let y = 0; y < h-1; y++) {
                const t1 = new triangle(0,x,y,ps[x][y],ps[x+1][y],ps[x+1][y+1]);
                const t2 = new triangle(1,x,y,ps[x][y],ps[x][y+1],ps[x+1][y+1]);
                trgrid[x][y*2+0] = t1; 
                trgrid[x][y*2+1] = t2;
                tr.push(t1, t2);
            }
        }
        
        const trg = (x, y) => (x >= 0 && x < w-1 && y >= 0 && y < h*2-1) ? trgrid[x][y] : false;
        for (let x = 0; x < w-1; x++) {
            for (let y = 0; y < h-1; y++) {
                trgrid[x][y*2+0].setNeighbours(trg(x, y*2-1), trg(x+1, y*2+1), trg(x, y*2+1));
                trgrid[x][y*2+1].setNeighbours(trg(x-1, y*2+0), trg(x, y*2+2), trg(x, y*2+0));
            }
        }
    }
}

// Optimized triangle rendering
class triangle {
    constructor(o, x, y, p1, p2, p3) {
        this.o = o; this.x = x; this.y = y;
        this.p1 = p1; this.p2 = p2; this.p3 = p3;
        this.normal = null;
        this.dist = null;
    }
    setNeighbours(t1, t2, t3) {
        this.t1 = t1; this.t2 = t2; this.t3 = t3;
    }
    draw(t) {
        const p1 = transform4([...this.p1.p,1], viewProjectionMatrix);
        const p2 = transform4([...this.p2.p,1], viewProjectionMatrix);
        const p3 = transform4([...this.p3.p,1], viewProjectionMatrix);

        const s = 50;
        const p = polygons.create();
        p.addPoints([p1[0]/p1[3]*s, -p1[1]/p1[3]*s],
                    [p2[0]/p2[3]*s, -p2[1]/p2[3]*s],
                    [p3[0]/p3[3]*s, -p3[1]/p3[3]*s]);
        
        const edge = (t1, t2) => {
            if (!t2) return false;
            const n1 = t1.getNormal();
            const n2 = t2.getNormal();
            sub3m(tempVec1, t2.p1.p, cameraPos);
            sub3m(tempVec2, t1.p1.p, cameraPos);
            return Math.sign(dot3(n2, tempVec1)) == Math.sign(dot3(n1, tempVec2));
        };
        
        if (this.o == 0) {
            if (this.x % 2 == 1 || edge(this, this.t2)) p.addSegments(p.cp[1], p.cp[2]);
            if (this.y % 2 == 0 || edge(this, this.t1)) p.addSegments(p.cp[0], p.cp[1]);
            if (edge(this, this.t3)) p.addSegments(p.cp[0], p.cp[2]);
        } else {
            if (this.x % 2 == 0 || edge(this, this.t1)) p.addSegments(p.cp[0], p.cp[1]);
            if (this.y % 2 == 1 || edge(this, this.t2)) p.addSegments(p.cp[1], p.cp[2]);
            if (edge(this, this.t3)) p.addSegments(p.cp[0], p.cp[2]);
        }
        
        if (this.x % 4 < 2) p.addHatching(Math.PI/4, 1);
        polygons.draw(t, p);
    }
    getDist() {
        if (this.dist !== null) return this.dist;
        // Use squared distances to avoid sqrt
        sub3m(tempVec1, this.p1.p, cameraPos);
        sub3m(tempVec2, this.p2.p, cameraPos);
        sub3m(tempVec3, this.p3.p, cameraPos);
        return this.dist = Math.max(
            len3sq(tempVec1),
            len3sq(tempVec2),
            len3sq(tempVec3)
        );
    }
    getNormal() {
        if (this.normal) return this.normal;
        sub3m(tempVec1, this.p1.p, this.p3.p);
        sub3m(tempVec2, this.p2.p, this.p3.p);
        this.normal = [0,0,0];
        cross3m(this.normal, tempVec1, tempVec2);
        return this.normal;
    }
}

const w = new world();
const r = new cloth(w, 77, 53, .1, [0,1,0]);

let seed = 1000;
w.points.sort((a,b) => random()-.5);

// The walk function will be called until it returns false.
function walk(i) {
    if (i < simulationSteps) {
        w.update(1/30, 8);
        
        turtle.jump(i*190/simulationSteps-95, 95);
        turtle.goto((i+1)*190/simulationSteps-95, 95);

        if (i == simulationSteps-1) {
            // Sort front to back using cached distances
            r.triangles.sort((a,b) => a.getDist() - b.getDist());
        }
    } else {
        r.triangles[i - simulationSteps].draw(turtle);
    }
    return i < r.triangles.length - 1 + simulationSteps;
}

function setupCamera(camPos, camLookat) {
    const viewMatrix = lookAt4m(camPos, camLookat, [0,1,0]);
    const projectionMatrix = perspective4m(0.25, 1);
    return multiply4m(projectionMatrix, viewMatrix);
}

// Original vec3 functions kept for compatibility
function scale3(a,b) { return [a[0]*b,a[1]*b,a[2]*b]; }
function dist3(a,b) { return Math.sqrt((a[0]-b[0])**2+(a[1]-b[1])**2+(a[2]-b[2])**2); }
function normalize3(a) { return scale3(a,1/len3(a)); }
function add3(a,b) { return [a[0]+b[0],a[1]+b[1],a[2]+b[2]]; }
function sub3(a,b) { return [a[0]-b[0],a[1]-b[1],a[2]-b[2]]; }
function cross3(a,b) { return [a[1]*b[2]-a[2]*b[1],a[2]*b[0]-a[0]*b[2],a[0]*b[1]-a[1]*b[0]]; }

// vec4 functions
function transform4(a,b) {
    const d=new Float32Array(4);
    for(let c=0;4>c;c++)d[c]=b[c]*a[0]+b[c+4]*a[1]+b[c+8]*a[2]+b[c+12];
    return d;
}
// mat4 functions
function lookAt4m(a,b,d) {
    const c=new Float32Array(16);
    b=normalize3(sub3(a,b));
    d=normalize3(cross3(d,b));
    const e=normalize3(cross3(b,d));
    c[0]=d[0];c[1]=e[0];c[2]=b[0];c[3]=0;
    c[4]=d[1];c[5]=e[1];c[6]=b[1];c[7]=0;
    c[8]=d[2];c[9]=e[2];c[10]=b[2];c[11]=0;
    c[12]=-(d[0]*a[0]+d[1]*a[1]+d[2]*a[2]);
    c[13]=-(e[0]*a[0]+e[1]*a[1]+e[2]*a[2]);
    c[14]=-(b[0]*a[0]+b[1]*a[1]+b[2]*a[2]);
    c[15]=1;
    return c;
}
function multiply4m(a,b) {
    const d=new Float32Array(16);
    for(let c=0;16>c;c+=4)
        for(let e=0;4>e;e++)
            d[c+e]=b[c+0]*a[0+e]+b[c+1]*a[4+e]+b[c+2]*a[8+e]+b[c+3]*a[12+e];
    return d;
}
function perspective4m(a,b) {
    const c=(new Float32Array(16)).fill(0,0);
    c[5]=1/Math.tan(a/2);
    c[0]=c[5]/b;
    c[10]=c[11]=-1;
    return c;
}

////////////////////////////////////////////////////////////////
// Polygon Clipping utility code - Created by Reinder Nijhoff 2019
////////////////////////////////////////////////////////////////

function Polygons() {
	const polygonList = [];
	const Polygon = class {
		constructor() {
			this.cp = [];
			this.dp = [];
			this.aabb = [];
		}
		addPoints(...points) {
		    let xmin = 1e5, xmax = -1e5, ymin = 1e5, ymax = -1e5;
			(this.cp = [...this.cp, ...points]).forEach( p => {
				xmin = Math.min(xmin, p[0]), xmax = Math.max(xmax, p[0]);
				ymin = Math.min(ymin, p[1]), ymax = Math.max(ymax, p[1]);
			});
		    this.aabb = [(xmin+xmax)/2, (ymin+ymax)/2, (xmax-xmin)/2, (ymax-ymin)/2];
		}
		addSegments(...points) {
		    points.forEach(p => this.dp.push(p));
		}
		addOutline() {
			for (let i = 0, l = this.cp.length; i < l; i++) {
				this.dp.push(this.cp[i], this.cp[(i + 1) % l]);
			}
		}
		draw(t) {
			for (let i = 0, l = this.dp.length; i < l; i+=2) {
				t.jump(this.dp[i]), t.goto(this.dp[i + 1]);
			}
		}
		addHatching(a, d) {
			const tp = new Polygon();
			tp.cp.push([-1e5,-1e5],[1e5,-1e5],[1e5,1e5],[-1e5,1e5]);
			const dx = Math.sin(a) * d,   dy = Math.cos(a) * d;
			const cx = Math.sin(a) * 200, cy = Math.cos(a) * 200;
			for (let i = 0.5; i < 150 / d; i++) {
				tp.dp.push([dx * i + cy,   dy * i - cx], [dx * i - cy,   dy * i + cx]);
				tp.dp.push([-dx * i + cy, -dy * i - cx], [-dx * i - cy, -dy * i + cx]);
			}
			tp.boolean(this, false);
			this.dp = [...this.dp, ...tp.dp];
		}
		inside(p) {
			let int = 0;
			for (let i = 0, l = this.cp.length; i < l; i++) {
				if (this.segment_intersect(p, [0.1, -1000], this.cp[i], this.cp[(i + 1) % l])) {
					int++;
				}
			}
			return int & 1;
		}
		boolean(p, diff = true) {
		    if (Math.abs(this.aabb[0] - p.aabb[0]) - (p.aabb[2] + this.aabb[2]) >= 0 &&
				Math.abs(this.aabb[1] - p.aabb[1]) - (p.aabb[3] + this.aabb[3]) >= 0) return this.dp.length > 0;
				
			const ndp = [];
			for (let i = 0, l = this.dp.length; i < l; i+=2) {
				const ls0 = this.dp[i];
				const ls1 = this.dp[i + 1];
				const int = [];
				for (let j = 0, cl = p.cp.length; j < cl; j++) {
					const pint = this.segment_intersect(ls0, ls1, p.cp[j], p.cp[(j + 1) % cl]);
					if (pint !== false) {
						int.push(pint);
					}
				}
				if (int.length === 0) {
					if (diff === !p.inside(ls0)) {
						ndp.push(ls0, ls1);
					}
				} else {
					int.push(ls0, ls1);
					const cmpx = ls1[0] - ls0[0];
					const cmpy = ls1[1] - ls0[1];
					int.sort( (a,b) =>  (a[0] - ls0[0]) * cmpx + (a[1] - ls0[1]) * cmpy - 
					                    (b[0] - ls0[0]) * cmpx - (b[1] - ls0[1]) * cmpy);
					 
					for (let j = 0; j < int.length - 1; j++) {
						if ((int[j][0] - int[j+1][0])**2 + (int[j][1] - int[j+1][1])**2 >= 0.001) {
							if (diff === !p.inside([(int[j][0]+int[j+1][0])/2,(int[j][1]+int[j+1][1])/2])) {
								ndp.push(int[j], int[j+1]);
							}
						}
					}
				}
			}
			return (this.dp = ndp).length > 0;
		}
		segment_intersect(l1p1, l1p2, l2p1, l2p2) {
			const d   = (l2p2[1] - l2p1[1]) * (l1p2[0] - l1p1[0]) - (l2p2[0] - l2p1[0]) * (l1p2[1] - l1p1[1]);
			if (d === 0) return false;
			const n_a = (l2p2[0] - l2p1[0]) * (l1p1[1] - l2p1[1]) - (l2p2[1] - l2p1[1]) * (l1p1[0] - l2p1[0]);
			const n_b = (l1p2[0] - l1p1[0]) * (l1p1[1] - l2p1[1]) - (l1p2[1] - l1p1[1]) * (l1p1[0] - l2p1[0]);
			const ua = n_a / d;
			const ub = n_b / d;
			if (ua >= 0 && ua <= 1 && ub >= 0 && ub <= 1) {
				return [l1p1[0] + ua * (l1p2[0] - l1p1[0]), l1p1[1] + ua * (l1p2[1] - l1p1[1])];
			}
			return false;
		}
	};
	return {
		list: () => polygonList,
		create: () => new Polygon(),
		draw: (turtle, p, addToVisList=true) => {
			for (let j = 0; j < polygonList.length && p.boolean(polygonList[j]); j++);
			p.draw(turtle);
			if (addToVisList) polygonList.push(p);
		}
	};
}

////////////////////////////////////////////////////////////////
// Slowpoke utility code
////////////////////////////////////////////////////////////////

function Slowpoke(x, y) {
    const linesDrawn = {};
    class Slowpoke extends Turtle {
        goto(x, y) {
            const p = Array.isArray(x) ? [...x] : [x, y];
            if (this.isdown()) {
                const o = [this.x(), this.y()];
                const h1 = o[0].toFixed(2)+'_'+p[0].toFixed(2)+o[1].toFixed(2)+'_'+p[1].toFixed(2);
                const h2 = p[0].toFixed(2)+'_'+o[0].toFixed(2)+p[1].toFixed(2)+'_'+o[1].toFixed(2);
                if (linesDrawn[h1] || linesDrawn[h2]) {
                    super.up();
                    super.goto(p);
                    super.down();
                    return;
                }
                linesDrawn[h1] = linesDrawn[h2] = true;
            } 
            super.goto(p);
        }
    }
    return new Slowpoke(x,y);
}

function random() {
    const bits = 20; 
    const mod = 1<<bits;
    let r = 1103515245 * (((seed+=12345) >> 1) ^ seed);
    r = 1103515245 * (r ^ (r >> 3));
    r = r ^ (r >> 16);
    return (r % mod) / mod;
}