fp64 distance estimate
Mandelbrot boundary rendered by distance estimate (d = ½·|z|·ln|z|/|dz|) with precision split mid-formula.
// ═══ typeshade example — fp64 Mandelbrot distance estimate ═══//// MIXED precision on purpose: the ORBIT (z) iterates in f64 — its absolute// position is what deep zoom destroys — while the DERIVATIVE (dz ← 2·z·dz + 1)// iterates in plain f32 from a per-step narrowed z, because the distance// estimate d = ½·|z|·ln|z| / |dz| only ever needs |dz| to a few digits.// Precision goes where it pays: the classic df64 discipline (opt in per VALUE,// not per shader). The boundary distance, normalised by the view span, shades// glowing filaments that stay crisp at any depth on the f64 side — the f32// left half collapses past a ~1e-7 span. Camera on a needle filament of the// period-3 minibrot neighbourhood (CPU-verified: 8–10 distinct distance// buckets per 40² window from 1e-5 down to 1e-12 spans).//// The `_fp64` guard uniform is auto-injected by the lowering; the render// harnesses bind it to 1.0f by probing the program for the Fp64Guard block.
import { fn, module, vec3, vec4, f32, f64, pow, sqrt, log, exp, min, max, toF32, toF64, f32T, vec2fT, vec2f64T, If, Loop, Var, Let, u32, uniformStruct,} from '../src/index.js';import { VsOut, vs } from './_fullscreen.js';import type { ShaderExample } from './_shared.js';
// A needle filament beside the period-3 minibrot (x ≈ −1.749); y = 0 keeps// the never-escaping real axis (and so the set's boundary) mid-frame forever.const CENTER_X = -1.7489;const CENTER_Y = 0;const ITER = 160;const ESCAPE_M2 = 1e6; // large escape radius tightens the distance estimate
const U = uniformStruct( 'Uniforms', { group: 0, binding: 0, as: 'u' }, { center: vec2f64T, // one DF64Vec2 slot — host packs [hi.x, hi.y, lo.x, lo.y] resolution: vec2fT, zoom_exp: f32T, // view span = 10^-zoom_exp complex units fp64: f32T, // toggle: 1 = split-screen f32 | f64 (canonical), 0 = all-f32 },);
const fsDe = fn( 'fs_de', { vo: VsOut }, (p) => { const span = Let(pow(f32(10.0), U.field.zoom_exp.neg())); const half = Let(p.vo.uv.x.mul(2.0)); const sx = Let(half.sub(p.vo.uv.x.lt(0.5).select(0.0, 1.0))); const dx = Let(sx.sub(0.5).mul(span)); const dy = Let( p.vo.uv.y.sub(0.5).mul(span).mul(U.field.resolution.y.div(U.field.resolution.x).mul(2.0)), );
const m2f = Var(f32(0)); // |z|² when the orbit stopped (≤ ESCAPE_M2 = interior) const dm2 = Var(f32(1)); // |dz|² at that moment If(p.vo.uv.x.lt(0.5).or(U.field.fp64.lt(0.5)), () => { // f32 twin — orbit AND derivative in f32, center narrowed. const cx = Let(toF32(U.field.center.x).add(dx)); const cy = Let(toF32(U.field.center.y).add(dy)); const zx = Var(f32(0)); const zy = Var(f32(0)); const ux = Var(f32(0)); // dz const uy = Var(f32(0)); Loop( u32(0), (j) => j.lt(u32(ITER)), () => { If(zx.mul(zx).add(zy.mul(zy)).le(ESCAPE_M2), () => { const nux = Let(zx.mul(ux).sub(zy.mul(uy)).mul(2.0).add(1.0)); uy.assign(zx.mul(uy).add(zy.mul(ux)).mul(2.0)); ux.assign(nux); const nzx = Let(zx.mul(zx).sub(zy.mul(zy)).add(cx)); zy.assign(zx.mul(zy).mul(2.0).add(cy)); zx.assign(nzx); }); }, ); m2f.assign(zx.mul(zx).add(zy.mul(zy))); dm2.assign(ux.mul(ux).add(uy.mul(uy))); }).else(() => { // f64 orbit / f32 derivative — z narrowed once per step for the dz twin. const cx = Let(U.field.center.x.add(toF64(dx))); const cy = Let(U.field.center.y.add(toF64(dy))); const zx = Var(f64(0)); const zy = Var(f64(0)); const ux = Var(f32(0)); const uy = Var(f32(0)); Loop( u32(0), (j) => j.lt(u32(ITER)), () => { If(toF32(zx.mul(zx).add(zy.mul(zy))).le(ESCAPE_M2), () => { const zx32 = Let(toF32(zx)); const zy32 = Let(toF32(zy)); const nux = Let(zx32.mul(ux).sub(zy32.mul(uy)).mul(2.0).add(1.0)); uy.assign(zx32.mul(uy).add(zy32.mul(ux)).mul(2.0)); ux.assign(nux); const nzx = Let(zx.mul(zx).sub(zy.mul(zy)).add(cx)); zy.assign(zx.mul(zy).mul(2.0).add(cy)); zx.assign(nzx); }); }, ); m2f.assign(toF32(zx.mul(zx).add(zy.mul(zy)))); dm2.assign(ux.mul(ux).add(uy.mul(uy))); });
// d = ½·|z|·ln|z| / |dz|, normalised by the span → depth-invariant shading. const mz = Let(sqrt(max(m2f, 1.0))); const de = Let( mz .mul(log(mz)) .mul(0.5) .div(sqrt(max(dm2, 1e-30))), ); const t = Let(min(de.div(span.mul(0.012)), 40.0)); // Interior (never escaped) → 0 glow; filaments glow where d/span → 0. const escaped = Let(m2f.gt(ESCAPE_M2).select(1.0, 0.0)); const glow = Let(exp(t.neg().mul(1.2)).mul(escaped)); const body = Let(exp(t.neg().mul(0.25)).mul(escaped)); const rgb = vec3(0.02, 0.03, 0.08) .add(vec3(0.12, 0.2, 0.42).mul(body)) .add(vec3(1.0, 0.85, 0.45).mul(glow)); return vec4(rgb, f32(1)); }, { stage: 'fragment', retAttr: '@location(0)' },);
// `_fp64` guard lands at (group 0, binding 1) automatically.const fp64MandelbrotDeModule = module({ funcs: [vs, fsDe], uses: [U, VsOut],});
export const fp64MandelbrotDe: ShaderExample = { id: 'fp64-mandelbrot-de', title: 'fp64 distance estimate', blurb: 'Mandelbrot boundary rendered by DISTANCE ESTIMATE (d = ½·|z|·ln|z|/|dz|) with precision split mid-formula: the orbit z iterates in emulated f64 (its absolute position is what deep zoom destroys) while the derivative dz iterates in plain f32 (the estimate needs |dz| to a few digits only). Filaments keep glowing at any depth on the right; the all-f32 left half collapses past a ~1e-7 span. Drag to pan, wheel to zoom, flip the fp64 toggle to compare.', category: 'generic', file: 'fp64-mandelbrot-de.ts', module: fp64MandelbrotDeModule, renderable: true, splitLabels: ['f32', 'f64 orbit + f32 dz'], controls: { center: { kind: 'pan2d', value: [CENTER_X, CENTER_Y], zoomExpField: 'zoom_exp', unitsPerWidth: 2, }, resolution: { kind: 'resolution' }, zoom_exp: { kind: 'slider', label: 'Zoom 10^-x', min: 1, max: 13, step: 0.05, value: 5, wheel: true, }, fp64: { kind: 'toggle', label: 'fp64 emulation', value: true }, },};struct Uniforms { @align(16) center: DF64Vec2, resolution: vec2<f32>, zoom_exp: f32, fp64: f32,}
struct VsOut { @builtin(position) pos: vec4<f32>, @location(0) uv: vec2<f32>,}
struct DF64Vec2 { hi: vec2<f32>, lo: vec2<f32>,}
@group(0) @binding(0) var<uniform> u: Uniforms;@group(0) @binding(1) var _fp64: texture_2d<f32>;
@vertexfn vs(@builtin(vertex_index) vi: u32) -> VsOut { let _cse0 = ((f32((vi & 1u)) * 4.0) - 1.0); let _cse1 = ((f32((vi >> 1u)) * 4.0) - 1.0); return VsOut(vec4<f32>(_cse0, _cse1, 0.0, 1.0), vec2<f32>(((_cse0 * 0.5) + 0.5), ((_cse1 * 0.5) + 0.5)));}
@fragmentfn fs_de(vo: VsOut) -> @location(0) vec4<f32> { let _fp64_g = textureLoad(_fp64, vec2<i32>(0, 0), 0).x; let _v0 = pow(10.0, (-u.zoom_exp)); let _v1 = (vo.uv.x * 2.0); let _cse0 = (vo.uv.x < 0.5); let _v2 = (_v1 - select(1.0, 0.0, _cse0)); let _v3 = ((_v2 - 0.5) * _v0); let _v4 = (((vo.uv.y - 0.5) * _v0) * ((u.resolution.y / u.resolution.x) * 2.0)); var _v5: f32 = 0.0; var _v6: f32 = 1.0; let _cse1 = vec2<f32>(u.center.hi.x, u.center.lo.x); let _cse2 = vec2<f32>(u.center.hi.y, u.center.lo.y); if ((_cse0 || (u.fp64 < 0.5))) { let _v7 = (df64_narrow(_cse1) + _v3); let _v8 = (df64_narrow(_cse2) + _v4); var _v9: f32 = 0.0; var _v10: f32 = 0.0; var _v11: f32 = 0.0; var _v12: f32 = 0.0; for (var _v13: u32 = 0u; (_v13 < 160u); _v13 = (_v13 + 1u)) { let _gv1 = (_v9 * _v9); let _gv2 = (_v10 * _v10); if (((_gv1 + _gv2) <= 1000000.0)) { let _v14 = ((((_v9 * _v11) - (_v10 * _v12)) * 2.0) + 1.0); _v12 = (((_v9 * _v12) + (_v10 * _v11)) * 2.0); _v11 = _v14; let _v15 = ((_gv1 - _gv2) + _v7); _v10 = (((_v9 * _v10) * 2.0) + _v8); _v9 = _v15; } } _v5 = ((_v9 * _v9) + (_v10 * _v10)); _v6 = ((_v11 * _v11) + (_v12 * _v12)); } else { let _v16 = df64_add(_cse1, vec2<f32>(_v3, 0.0), _fp64_g); let _v17 = df64_add(_cse2, vec2<f32>(_v4, 0.0), _fp64_g); let _cse3 = vec2<f32>(0.0, 0.0); var _v18: vec2<f32> = _cse3; var _v19: vec2<f32> = _cse3; var _v20: f32 = 0.0; var _v21: f32 = 0.0; for (var _v22: u32 = 0u; (_v22 < 160u); _v22 = (_v22 + 1u)) { let _gv3 = df64_sqr(_v18, _fp64_g); let _gv4 = df64_sqr(_v19, _fp64_g); if ((df64_narrow(df64_add(_gv3, _gv4, _fp64_g)) <= 1000000.0)) { let _v23 = df64_narrow(_v18); let _v24 = df64_narrow(_v19); let _v25 = ((((_v23 * _v20) - (_v24 * _v21)) * 2.0) + 1.0); _v21 = (((_v23 * _v21) + (_v24 * _v20)) * 2.0); _v20 = _v25; let _v26 = df64_add(df64_sub(_gv3, _gv4, _fp64_g), _v16, _fp64_g); _v19 = df64_add((df64_mul(_v18, _v19, _fp64_g) * 2.0), _v17, _fp64_g); _v18 = _v26; } } _v5 = df64_narrow(df64_add(df64_sqr(_v18, _fp64_g), df64_sqr(_v19, _fp64_g), _fp64_g)); _v6 = ((_v20 * _v20) + (_v21 * _v21)); } let _v27 = sqrt(max(_v5, 1.0)); let _v28 = (((_v27 * log(_v27)) * 0.5) / sqrt(max(_v6, 1e-30))); let _v29 = min((_v28 / (_v0 * 0.012)), 40.0); let _v30 = select(0.0, 1.0, (_v5 > 1000000.0)); let _gv0 = (-_v29); let _v31 = (exp((_gv0 * 1.2)) * _v30); let _v32 = (exp((_gv0 * 0.25)) * _v30); return vec4<f32>(((vec3<f32>(0.02, 0.03, 0.08) + (vec3<f32>(0.12, 0.2, 0.42) * _v32)) + (vec3<f32>(1.0, 0.85, 0.45) * _v31)), 1.0);}
fn df64_twoSum(a: f32, b: f32, _fp64_g: f32) -> vec2<f32> { let _v0 = (a + b); let _v1 = (((_v0 * _fp64_g) - a) * _fp64_g); let _v2 = (((a - ((_v0 - _v1) * _fp64_g)) * _fp64_g) + (b - _v1)); return vec2<f32>(_v0, _v2);}
fn df64_quickTwoSum(a: f32, b: f32, _fp64_g: f32) -> vec2<f32> { let _v0 = ((a + b) * _fp64_g); let _v1 = (b - ((_v0 - a) * _fp64_g)); return vec2<f32>(_v0, _v1);}
fn df64_split(a: f32, _fp64_g: f32) -> vec2<f32> { let _v0 = (a * (_fp64_g * 4097.0)); let _v1 = ((_v0 * _fp64_g) - (_v0 - a)); let _v2 = ((a * _fp64_g) - _v1); return vec2<f32>(_v1, _v2);}
fn df64_twoProd(a: f32, b: f32, _fp64_g: f32) -> vec2<f32> { let _v0 = (a * b); let _v1 = df64_split(a, _fp64_g); let _v2 = df64_split(b, _fp64_g); let _v3 = (((((_v1.x * _v2.x) - _v0) + (_v1.x * _v2.y)) + (_v1.y * _v2.x)) + (_v1.y * _v2.y)); return vec2<f32>(_v0, _v3);}
fn df64_twoSqr(a: f32, _fp64_g: f32) -> vec2<f32> { let _v0 = (a * a); let _v1 = df64_split(a, _fp64_g); let _v2 = (((((_v1.x * _v1.x) - _v0) * _fp64_g) + (((_v1.x * _v1.y) * 2.0) * _fp64_g)) + ((_v1.y * _v1.y) * _fp64_g)); return vec2<f32>(_v0, _v2);}
fn df64_add(a: vec2<f32>, b: vec2<f32>, _fp64_g: f32) -> vec2<f32> { var _v0: vec2<f32> = df64_twoSum(a.x, b.x, _fp64_g); let _v1 = df64_twoSum(a.y, b.y, _fp64_g); _v0.y = (_v0.y + _v1.x); _v0 = df64_quickTwoSum(_v0.x, _v0.y, _fp64_g); _v0.y = (_v0.y + _v1.y); _v0 = df64_quickTwoSum(_v0.x, _v0.y, _fp64_g); return _v0;}
fn df64_sub(a: vec2<f32>, b: vec2<f32>, _fp64_g: f32) -> vec2<f32> { return df64_add(a, (-b), _fp64_g);}
fn df64_mul(a: vec2<f32>, b: vec2<f32>, _fp64_g: f32) -> vec2<f32> { var _v0: vec2<f32> = df64_twoProd(a.x, b.x, _fp64_g); _v0.y = (_v0.y + (a.x * b.y)); _v0 = df64_quickTwoSum(_v0.x, _v0.y, _fp64_g); _v0.y = (_v0.y + (a.y * b.x)); return df64_quickTwoSum(_v0.x, _v0.y, _fp64_g);}
fn df64_sqr(a: vec2<f32>, _fp64_g: f32) -> vec2<f32> { var _v0: vec2<f32> = df64_twoSqr(a.x, _fp64_g); _v0.y = (_v0.y + ((a.x * a.y) * 2.0)); return df64_quickTwoSum(_v0.x, _v0.y, _fp64_g);}
fn df64_narrow(a: vec2<f32>) -> f32 { return (a.x + a.y);}#version 300 esprecision highp float;precision highp int;
out vec2 uv;
void main() { uint vi = uint(gl_VertexID); float _cse0 = ((float((vi & 1u)) * 4.0) - 1.0); float _cse1 = ((float((vi >> 1u)) * 4.0) - 1.0); gl_Position = vec4(_cse0, _cse1, 0.0, 1.0); uv = vec2(((_cse0 * 0.5) + 0.5), ((_cse1 * 0.5) + 0.5));}#version 300 esprecision highp float;precision highp int;
struct DF64Vec2 { vec2 hi; vec2 lo;};layout(std140) uniform Uniforms { DF64Vec2 center; vec2 resolution; float zoom_exp; float fp64;} u;
uniform highp sampler2D _fp64;vec2 df64_twoSum(float a, float b, float _fp64_g) { float _v0 = (a + b); float _v1 = (((_v0 * _fp64_g) - a) * _fp64_g); float _v2 = (((a - ((_v0 - _v1) * _fp64_g)) * _fp64_g) + (b - _v1)); return vec2(_v0, _v2);}
vec2 df64_quickTwoSum(float a, float b, float _fp64_g) { float _v0 = ((a + b) * _fp64_g); float _v1 = (b - ((_v0 - a) * _fp64_g)); return vec2(_v0, _v1);}
vec2 df64_split(float a, float _fp64_g) { float _v0 = (a * (_fp64_g * 4097.0)); float _v1 = ((_v0 * _fp64_g) - (_v0 - a)); float _v2 = ((a * _fp64_g) - _v1); return vec2(_v1, _v2);}
vec2 df64_twoProd(float a, float b, float _fp64_g) { float _v0 = (a * b); vec2 _v1 = df64_split(a, _fp64_g); vec2 _v2 = df64_split(b, _fp64_g); float _v3 = (((((_v1.x * _v2.x) - _v0) + (_v1.x * _v2.y)) + (_v1.y * _v2.x)) + (_v1.y * _v2.y)); return vec2(_v0, _v3);}
vec2 df64_twoSqr(float a, float _fp64_g) { float _v0 = (a * a); vec2 _v1 = df64_split(a, _fp64_g); float _v2 = (((((_v1.x * _v1.x) - _v0) * _fp64_g) + (((_v1.x * _v1.y) * 2.0) * _fp64_g)) + ((_v1.y * _v1.y) * _fp64_g)); return vec2(_v0, _v2);}
vec2 df64_add(vec2 a, vec2 b, float _fp64_g) { vec2 _v0 = df64_twoSum(a.x, b.x, _fp64_g); vec2 _v1 = df64_twoSum(a.y, b.y, _fp64_g); _v0.y = (_v0.y + _v1.x); _v0 = df64_quickTwoSum(_v0.x, _v0.y, _fp64_g); _v0.y = (_v0.y + _v1.y); _v0 = df64_quickTwoSum(_v0.x, _v0.y, _fp64_g); return _v0;}
vec2 df64_sub(vec2 a, vec2 b, float _fp64_g) { return df64_add(a, (-b), _fp64_g);}
vec2 df64_mul(vec2 a, vec2 b, float _fp64_g) { vec2 _v0 = df64_twoProd(a.x, b.x, _fp64_g); _v0.y = (_v0.y + (a.x * b.y)); _v0 = df64_quickTwoSum(_v0.x, _v0.y, _fp64_g); _v0.y = (_v0.y + (a.y * b.x)); return df64_quickTwoSum(_v0.x, _v0.y, _fp64_g);}
vec2 df64_sqr(vec2 a, float _fp64_g) { vec2 _v0 = df64_twoSqr(a.x, _fp64_g); _v0.y = (_v0.y + ((a.x * a.y) * 2.0)); return df64_quickTwoSum(_v0.x, _v0.y, _fp64_g);}
float df64_narrow(vec2 a) { return (a.x + a.y);}in vec2 uv;layout(location = 0) out vec4 _ret;
void main() { float _fp64_g = texelFetch(_fp64, ivec2(0, 0), 0).x; float _v0 = pow(10.0, (-u.zoom_exp)); float _v1 = (uv.x * 2.0); bool _cse0 = (uv.x < 0.5); float _v2 = (_v1 - (_cse0 ? 0.0 : 1.0)); float _v3 = ((_v2 - 0.5) * _v0); float _v4 = (((uv.y - 0.5) * _v0) * ((u.resolution.y / u.resolution.x) * 2.0)); float _v5 = 0.0; float _v6 = 1.0; vec2 _cse1 = vec2(u.center.hi.x, u.center.lo.x); vec2 _cse2 = vec2(u.center.hi.y, u.center.lo.y); if ((_cse0 || (u.fp64 < 0.5))) { float _v7 = (df64_narrow(_cse1) + _v3); float _v8 = (df64_narrow(_cse2) + _v4); float _v9 = 0.0; float _v10 = 0.0; float _v11 = 0.0; float _v12 = 0.0; for (uint _v13 = 0u; (_v13 < 160u); _v13 = (_v13 + 1u)) { float _gv1 = (_v9 * _v9); float _gv2 = (_v10 * _v10); if (((_gv1 + _gv2) <= 1000000.0)) { float _v14 = ((((_v9 * _v11) - (_v10 * _v12)) * 2.0) + 1.0); _v12 = (((_v9 * _v12) + (_v10 * _v11)) * 2.0); _v11 = _v14; float _v15 = ((_gv1 - _gv2) + _v7); _v10 = (((_v9 * _v10) * 2.0) + _v8); _v9 = _v15; } } _v5 = ((_v9 * _v9) + (_v10 * _v10)); _v6 = ((_v11 * _v11) + (_v12 * _v12)); } else { vec2 _v16 = df64_add(_cse1, vec2(_v3, 0.0), _fp64_g); vec2 _v17 = df64_add(_cse2, vec2(_v4, 0.0), _fp64_g); vec2 _cse3 = vec2(0.0, 0.0); vec2 _v18 = _cse3; vec2 _v19 = _cse3; float _v20 = 0.0; float _v21 = 0.0; for (uint _v22 = 0u; (_v22 < 160u); _v22 = (_v22 + 1u)) { vec2 _gv3 = df64_sqr(_v18, _fp64_g); vec2 _gv4 = df64_sqr(_v19, _fp64_g); if ((df64_narrow(df64_add(_gv3, _gv4, _fp64_g)) <= 1000000.0)) { float _v23 = df64_narrow(_v18); float _v24 = df64_narrow(_v19); float _v25 = ((((_v23 * _v20) - (_v24 * _v21)) * 2.0) + 1.0); _v21 = (((_v23 * _v21) + (_v24 * _v20)) * 2.0); _v20 = _v25; vec2 _v26 = df64_add(df64_sub(_gv3, _gv4, _fp64_g), _v16, _fp64_g); _v19 = df64_add((df64_mul(_v18, _v19, _fp64_g) * 2.0), _v17, _fp64_g); _v18 = _v26; } } _v5 = df64_narrow(df64_add(df64_sqr(_v18, _fp64_g), df64_sqr(_v19, _fp64_g), _fp64_g)); _v6 = ((_v20 * _v20) + (_v21 * _v21)); } float _v27 = sqrt(max(_v5, 1.0)); float _v28 = (((_v27 * log(_v27)) * 0.5) / sqrt(max(_v6, 1e-30))); float _v29 = min((_v28 / (_v0 * 0.012)), 40.0); float _v30 = ((_v5 > 1000000.0) ? 1.0 : 0.0); float _gv0 = (-_v29); float _v31 = (exp((_gv0 * 1.2)) * _v30); float _v32 = (exp((_gv0 * 0.25)) * _v30); _ret = vec4(((vec3(0.02, 0.03, 0.08) + (vec3(0.12, 0.2, 0.42) * _v32)) + (vec3(1.0, 0.85, 0.45) * _v31)), 1.0);}This example is steered by a control the page has no value for, so the page shows no picture.
The WGSL and GLSL tabs are the compiler's own output at commit 26de7be8, read from the goldens its emit suite bakes (emit-goldens.test.ts).
This example is written against the fn() builder API, which the editor in the Playground does not take.