Math & Utility

fn elasticWave(x: f32, amplitude: f32, frequency: f32, decay: f32, phase: f32) -> f32 {
  let d = max(0.001, decay);
  let decayTerm = exp(-d * x);
  let oscTerm = sin(frequency * x * 6.28318 + phase);
  return amplitude * decayTerm * oscTerm;
}

fn exponentialRamp(x: f32, base: f32, scale: f32, offset: f32) -> vec2<f32> {
  // Ensure base is positive and not 1 (which would make it linear)
  let b = select(base, 2.71828, abs(base - 1.0) < 0.001);
  
  // Calculate the exponential function
  let result = scale * pow(b, x) + offset;
  
  // Calculate the derivative
  let derivative = scale * pow(b, x) * log(b);
  
  return vec2<f32>(result, derivative);
}

fn logisticCurve(x: f32, midpoint: f32, steepness: f32, minValue: f32, maxValue: f32) -> vec2<f32> {
  // Scale factor for steepness
  let k = max(0.001, steepness);
  
  // Shift x relative to midpoint
  let z = -k * (x - midpoint);
  
  // Calculate the exponent
  let expTerm = exp(z);
  
  // Calculate the logistic function value
  let logistic = 1.0 / (1.0 + expTerm);
  
  // Scale to min-max range
  let range = maxValue - minValue;
  let value = minValue + range * logistic;
  
  // Calculate the derivative
  let derivative = range * k * expTerm / ((1.0 + expTerm) * (1.0 + expTerm));
  
  return vec2<f32>(value, derivative);
}

fn rotate2D(v: vec2<f32>, angle: f32) -> vec2<f32> {
  let c = cos(angle);
  let s = sin(angle);
  return vec2(v.x * c - v.y * s, v.x * s + v.y * c);
}

fn smoothStep(edge0: f32, edge1: f32, x: f32) -> f32 {
  let t = clamp((x - edge0) / (edge1 - edge0), 0.0, 1.0);
  return t * t * (3.0 - 2.0 * t);
}

fn smoothStepVec2(edge0: vec2f, edge1: vec2f, x: vec2f) -> vec2f {
  let t = clamp((x - edge0) / (edge1 - edge0), vec2f(0.0), vec2f(1.0));
  return t * t * (3.0 - 2.0 * t);
}

fn stepSequence(x: f32, steps: f32, smoothing: f32, minValue: f32, maxValue: f32) -> vec2<f32> {
  // Ensure at least 1 step and positive smoothing
  let numSteps = max(1.0, floor(steps));
  let smoothFactor = max(0.0, smoothing);
  
  // Normalize x to 0-1 range
  let normalizedX = fract(x);
  
  // Calculate the size of each step
  let stepSize = 1.0 / numSteps;
  
  // Calculate the current step (0 to numSteps-1)
  let currentStep = floor(normalizedX * numSteps);
  let nextStep = fract(currentStep + 1.0);
  
  // Calculate progress within the current step
  let stepProgress = fract(normalizedX * numSteps);
  
  // Calculate the progress values for current and next steps
  let currentStepValue = currentStep / (numSteps - 1.0);
  
  // Prepare next step value, handle the last step case
  var nextStepValue: f32 = 0.0;
  if (currentStep >= numSteps - 1.0) {
    nextStepValue = 1.0;
  } else {
    nextStepValue = nextStep / (numSteps - 1.0);
  }
  
  // Apply smoothing between steps if needed
  var result: f32 = 0.0;
  
  if (smoothFactor > 0.0 && stepProgress > (1.0 - smoothFactor) && numSteps > 1.0) {
    // Calculate smoothing factor
    let t = (stepProgress - (1.0 - smoothFactor)) / smoothFactor;
    
    // Smoothstep for better transition
    let smoothT = t * t * (3.0 - 2.0 * t);
    
    // Interpolate between current and next step
    result = mix(currentStepValue, nextStepValue, smoothT);
  } else {
    result = currentStepValue;
  }
  
  // Scale to min-max range
  let range = maxValue - minValue;
  let finalResult = minValue + result * range;
  
  return vec2<f32>(finalResult, currentStep);
}

fn taylorInvSqrt4(r: vec4f) -> vec4f { 
    return 1.79284291400159 - 0.85373472095314 * r; 
}