SDF Primitives

fn sdfBox(p: vec2<f32>, b: vec2<f32>) -> f32 {
  let d = abs(p) - b;
  return length(max(d, vec2(0.0))) + min(max(d.x, d.y), 0.0);
}

fn sdfBoxFrame(position: vec3<f32>, size: vec3<f32>, thickness: f32) -> f32 {
  let q = abs(position) - size;
  let w = abs(q + thickness) - thickness;
  return min(min(
    length(max(vec3<f32>(q.x, w.y, w.z), vec3<f32>(0.0))) + min(max(q.x, max(w.y, w.z)), 0.0),
    length(max(vec3<f32>(w.x, q.y, w.z), vec3<f32>(0.0))) + min(max(w.x, max(q.y, w.z)), 0.0)),
    length(max(vec3<f32>(w.x, w.y, q.z), vec3<f32>(0.0))) + min(max(w.x, max(w.y, q.z)), 0.0));
}

fn sdfCappedTorus(position: vec3<f32>, majorRadius: f32, minorRadius: f32, angle: f32) -> f32 {
  let sc = vec2<f32>(sin(angle), cos(angle));
  let q = vec3<f32>(abs(position.x), position.y, position.z);
  let k = select(
    length(q.xy), 
    dot(q.xy, sc),
    sc.y * q.x > sc.x * q.y
  );
  return sqrt(dot(q, q) + 
              majorRadius * majorRadius - 
              2.0 * majorRadius * k) - 
              minorRadius;
}

fn sdfCapsule(position: vec3<f32>, radius: f32, height: f32) -> f32 {
  let d = abs(length(position.xz)) - radius;
  let p = vec2<f32>(d, abs(position.y) - height * 0.5);
  return length(max(p, vec2<f32>(0.0))) + min(max(p.x, p.y), 0.0) - radius;
}

fn sdfCircle(p: vec2<f32>, r: f32) -> f32 {
  return length(p) - r;
}

fn sdfCone(position: vec3<f32>, radius: f32, height: f32) -> f32 {
  let q = vec2<f32>(length(position.xz), position.y);
  let h = height;
  let r = radius;
  
  // Calculate distance
  let d1 = -q.y - h;
  let d2 = max(q.x * h - q.y * r, q.y * h + q.x * r);
  
  return length(max(vec2<f32>(d1, d2), vec2<f32>(0.0))) + min(max(d1, d2), 0.0);
}

fn sdfCylinder(position: vec3<f32>, radius: f32, height: f32) -> f32 {
  let d = vec2<f32>(length(position.xz), abs(position.y)) - vec2<f32>(radius, height * 0.5);
  return min(max(d.x, d.y), 0.0) + length(max(d, vec2<f32>(0.0)));
}

fn sdfEllipsoid(position: vec3<f32>, radius: vec3<f32>) -> f32 {
  let k0 = length(position / radius);
  let k1 = length(position / (radius * radius));
  return k0 * (k0 - 1.0) / k1;
}

fn sdfGyroid(position: vec3<f32>, scale: f32, thickness: f32) -> f32 {
  let p = position * scale;
  return (abs(dot(sin(p), cos(p.zxy))) - thickness) / scale;
}

fn sdfHexagonalPrism(position: vec3<f32>, radius: f32, height: f32) -> f32 {
  // Project into 2D
  var p = abs(position);
  let k = vec3<f32>(-0.866025404, 0.5, 0.577350269);
  
  // Hexagon in xy-plane
  p = vec3<f32>(p.x + p.y * k.x, p.y * k.y, p.z);
  p = vec3<f32>(p.x - min(p.x, p.y), p.y, p.z);
  let d = vec2<f32>(length(vec2<f32>(p.x, p.y - radius * k.z)) - radius, abs(p.z) - height * 0.5);
  
  return min(max(d.x, d.y), 0.0) + length(max(d, vec2<f32>(0.0)));
}

fn sdfIcosahedron(position: vec3<f32>, size: f32) -> f32 {
  var p = position;
  let s = size;
  
  // Constants for icosahedron
  let phi = 1.618033988749895;
  let a = s;
  let b = s * phi;
  
  // Compute distance to icosahedron
  p = abs(p / s);
  let d = p.x * p.y * p.z;
  let m = max(max(p.x, p.y), p.z);
  let n = min(min(p.x, p.y), p.z);
  let mid = p.x + p.y + p.z - m - n;
  
  // Calculate the signed distance
  let q = select(mid, d, m < phi * n);
  return (length(p) - phi) * s;
}

fn sdfIntersection(d1: f32, d2: f32) -> f32 {
  return max(d1, d2);
}

fn sdfJulia(position: vec3<f32>, c: vec4<f32>, iterations: f32, bailout: f32) -> f32 {
  var z = vec4<f32>(position, 0.0);
  var dz = vec4<f32>(1.0, 0.0, 0.0, 0.0);
  var m = dot(z, z);
  var i = 0;
  
  // Julia set iteration
  for (i = 0; i < i32(iterations) && m < bailout * bailout; i += 1) {
    // Quaternion multiplication for dz = 2.0 * z * dz
    dz = 2.0 * vec4<f32>(
      z.x * dz.x - z.y * dz.y - z.z * dz.z - z.w * dz.w,
      z.x * dz.y + z.y * dz.x + z.z * dz.w - z.w * dz.z,
      z.x * dz.z - z.y * dz.w + z.z * dz.x + z.w * dz.y,
      z.x * dz.w + z.y * dz.z - z.z * dz.y + z.w * dz.x
    );
    
    // Quaternion multiplication for z = z * z + c
    z = vec4<f32>(
      z.x * z.x - z.y * z.y - z.z * z.z - z.w * z.w,
      z.x * z.y + z.y * z.x + z.z * z.w - z.w * z.z,
      z.x * z.z - z.y * z.w + z.z * z.x + z.w * z.y,
      z.x * z.w + z.y * z.z - z.z * z.y + z.w * z.x
    ) + c;
    
    m = dot(z, z);
  }
  
  // Compute the distance
  let dist = 0.5 * log(m) * sqrt(m) / length(dz);
  return dist;
}

fn sdfOctahedron(position: vec3<f32>, size: f32) -> f32 {
  let p = abs(position);
  let m = p.x + p.y + p.z - size;
  
  // Calculate the distance
  var q: vec3<f32>;
  if (3.0 * p.x < m) {
    q = p;
  } else if (3.0 * p.y < m) {
    q = vec3<f32>(p.x, p.z, p.y);
  } else if (3.0 * p.z < m) {
    q = vec3<f32>(p.x, p.y, p.z);
  } else {
    q = p;
  }
  
  let k = clamp(0.5 * (q.z - q.y + size), 0.0, size);
  return length(vec3<f32>(q.x, q.y - size + k, q.z - k));
}

fn sdfPlane(position: vec3<f32>, normal: vec3<f32>) -> f32 {
  let n = normalize(normal);
  return dot(position, n);
}

fn sdfPyramid(position: vec3<f32>, size: f32, height: f32) -> f32 {
  // Normalize position
  var p = position;
  let h = height;
  let m2 = h * h + size * size;
  
  // Project into 2D
  let q = abs(p);
  p.y -= h;
  p.y = max(p.y, 0.0);
  
  // Distance calculation
  var d: f32;
  if (max(q.x, q.z) < size) {
    d = length(vec2<f32>(length(p.xz), p.y)) - sqrt(m2);
  } else {
    d = length(vec2<f32>(length(max(abs(p.xz) - vec2<f32>(size), vec2<f32>(0.0))), p.y));
  }
  
  // Account for position below base
  d = select(d, length(p) - sqrt(m2), p.y < 0.0);
  
  return d;
}

fn sdfRhombus(position: vec3<f32>, dimensions: vec3<f32>, sharpness: f32) -> f32 {
  var p = abs(position);
  let b = dimensions;
  let e = sharpness;
  
  // Calculate distance to rhombus
  p = p - b;
  let q = abs(p.x + p.y + p.z) + e;
  let h = max(vec3<f32>(q) - vec3<f32>(e), vec3<f32>(0.0));
  
  return min(max(p.x, max(p.y, p.z)), 0.0) + length(h);
}

fn sdfRoundBox(position: vec3<f32>, size: vec3<f32>, radius: f32) -> f32 {
  let q = abs(position) - size;
  return length(max(q, vec3<f32>(0.0))) + 
         min(max(q.x, max(q.y, q.z)), 0.0) - 
         radius;
}

fn sdfRoundedCone(position: vec3<f32>, radius1: f32, radius2: f32, height: f32, roundness: f32) -> f32 {
  // Calculate distances
  let p = position;
  let r1 = radius1 - roundness;
  let r2 = radius2 - roundness;
  let h = height;
  
  // Squared distance from axis
  let q = length(p.xz);
  
  // Project into 2D space
  let k1 = (r2 - r1) / h;
  let k2 = h / (r1 - r2);
  let projected = vec2<f32>(q - r1 + r1 * (p.y / h) * (r1 - r2) / r1, p.y - h);
  let ca = p.y * k1 - q;
  let cb = p.y - r1 * k2 + q * k2;
  
  var s: f32;
  if (ca < 0.0 && projected.y < 0.0) {
    s = length(projected) - roundness;
  } else if (ca > 0.0 && cb < 0.0) {
    s = -ca - roundness;
  } else {
    s = length(vec2<f32>(max(ca, 0.0), max(projected.y, 0.0))) - roundness;
  }
  
  return s;
}

fn sdfRoundedCylinder(position: vec3<f32>, radius: f32, height: f32, roundness: f32) -> f32 {
  // Calculate distances
  let radiusOffset = radius - roundness;
  let heightOffset = height * 0.5 - roundness;
  
  // Generate rounded cylinder
  let d = vec2<f32>(length(position.xz) - radiusOffset, abs(position.y) - heightOffset);
  return min(max(d.x, d.y), 0.0) + length(max(d, vec2<f32>(0.0))) - roundness;
}

fn sdfSphere(position: vec3<f32>, radius: f32) -> f32 {
  return length(position) - radius;
}

fn sdfSubtraction(d1: f32, d2: f32) -> f32 {
  return max(-d1, d2);
}

fn sdfTetrahedron(position: vec3<f32>, size: f32) -> f32 {
  var p = position;
  let s = size;
  
  // Set initial values
  let signVal = sign(p.x + p.y + p.z);
  p.x = abs(p.x);
  p.y = abs(p.y);
  p.z = abs(p.z);
  
  // Calculate the distance
  if (p.x < p.y) {
    let t = p.x;
    p.x = p.y;
    p.y = t;
  }
  if (p.x < p.z) {
    let t = p.x;
    p.x = p.z;
    p.z = t;
  }
  if (p.y < p.z) {
    let t = p.y;
    p.y = p.z;
    p.z = t;
  }
  
  let k = clamp((p.x + p.z - p.y) * 0.5, 0.0, p.z);
  return signVal * (length(vec3<f32>(p.x, p.y - s, p.z - k)) - s);
}

fn sdfTorus(position: vec3<f32>, majorRadius: f32, minorRadius: f32) -> f32 {
  let q = vec2<f32>(length(position.xz) - majorRadius, position.y);
  return length(q) - minorRadius;
}

fn sdfTriangularPrism(position: vec3<f32>, radius: f32, height: f32) -> f32 {
  var q = abs(position);
  
  // Triangle distance in xy-plane
  let k = sqrt(3.0);
  q.x = abs(q.x - q.y * k * 0.5);
  q.y = q.y * 0.866025404 + q.x * 0.5;
  
  // Combine with z distance
  let d1 = vec2<f32>(q.x - radius, q.y);
  let d2 = vec2<f32>(q.y - radius, q.x);
  let d = min(d1, d2);
  
  // Account for height
  let h = height * 0.5;
  let dz = q.z - h;
  let dz2 = max(dz, 0.0);
  
  return length(max(vec2<f32>(max(d.x, 0.0), dz2), vec2<f32>(0.0))) + min(max(d.x, dz), 0.0);
}

fn sdfUnion(d1: f32, d2: f32) -> f32 {
  return min(d1, d2);
}