fp64 distance estimate
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.
One source, compiled every way
The DSL source below, compiled to WebGPU (WGSL), WebGL2 (GLSL ES 3.00 — two stages), and a reflection JSON the host binds from. Nothing but the DSL tab is hand-written.
// ═══ @xgis/shader-dsl 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.ts'import { VsOut, vs } from './_fullscreen.ts'import type { ShaderExample } from './_shared.ts'
// 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.7489const CENTER_Y = 0const ITER = 160const 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 { 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 _licm0 = vec2<f32>(2.0, 0.0); let _cse0 = (vo.uv.x < 0.5); let _cse1 = vec2<f32>(u.center.hi.x, u.center.lo.x); let _cse2 = vec2<f32>(u.center.hi.y, u.center.lo.y); let _cse3 = vec2<f32>(0.0, 0.0); let _v0 = pow(10.0, (-u.zoom_exp)); let _v1 = (vo.uv.x * 2.0); 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; 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)) { if ((((_v9 * _v9) + (_v10 * _v10)) <= 1000000.0)) { let _v14 = ((((_v9 * _v11) - (_v10 * _v12)) * 2.0) + 1.0); _v12 = (((_v9 * _v12) + (_v10 * _v11)) * 2.0); _v11 = _v14; let _v15 = (((_v9 * _v9) - (_v10 * _v10)) + _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)); let _v17 = df64_add(_cse2, vec2<f32>(_v4, 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)) { if ((df64_narrow(df64_add(df64_mul(_v18, _v18), df64_mul(_v19, _v19))) <= 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(df64_mul(_v18, _v18), df64_mul(_v19, _v19)), _v16); _v19 = df64_add(df64_mul(df64_mul(_v18, _v19), _licm0), _v17); _v18 = _v26; } } _v5 = df64_narrow(df64_add(df64_mul(_v18, _v18), df64_mul(_v19, _v19))); _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 _v31 = (exp(((-_v29) * 1.2)) * _v30); let _v32 = (exp(((-_v29) * 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) -> vec2<f32> { let _cse0 = textureLoad(_fp64, vec2<i32>(0, 0), 0).x; let _v0 = (a + b); let _v1 = (((_v0 * _cse0) - a) * _cse0); let _v2 = (((((a - ((_v0 - _v1) * _cse0)) * _cse0) * _cse0) * _cse0) + (b - _v1)); return vec2<f32>(_v0, _v2);}
fn df64_quickTwoSum(a: f32, b: f32) -> vec2<f32> { let _cse0 = textureLoad(_fp64, vec2<i32>(0, 0), 0).x; let _v0 = ((a + b) * _cse0); let _v1 = (b - ((_v0 - a) * _cse0)); return vec2<f32>(_v0, _v1);}
fn df64_split(a: f32) -> vec2<f32> { let _cse0 = textureLoad(_fp64, vec2<i32>(0, 0), 0).x; let _v0 = (a * (_cse0 * 4097.0)); let _v1 = ((_v0 * _cse0) - (_v0 - a)); let _v2 = ((a * _cse0) - _v1); return vec2<f32>(_v1, _v2);}
fn df64_twoProd(a: f32, b: f32) -> vec2<f32> { let _v0 = (a * b); let _v1 = df64_split(a); let _v2 = df64_split(b); 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_add(a: vec2<f32>, b: vec2<f32>) -> vec2<f32> { var _v0: vec2<f32> = df64_twoSum(a.x, b.x); let _v1 = df64_twoSum(a.y, b.y); _v0.y = (_v0.y + _v1.x); _v0 = df64_quickTwoSum(_v0.x, _v0.y); _v0.y = (_v0.y + _v1.y); _v0 = df64_quickTwoSum(_v0.x, _v0.y); return _v0;}
fn df64_sub(a: vec2<f32>, b: vec2<f32>) -> vec2<f32> { return df64_add(a, (-b));}
fn df64_mul(a: vec2<f32>, b: vec2<f32>) -> vec2<f32> { var _v0: vec2<f32> = df64_twoProd(a.x, b.x); _v0.y = (_v0.y + (a.x * b.y)); _v0 = df64_quickTwoSum(_v0.x, _v0.y); _v0.y = (_v0.y + (a.y * b.x)); return df64_quickTwoSum(_v0.x, _v0.y);}
fn df64_narrow(a: vec2<f32>) -> f32 { return (a.x + a.y);}#version 300 esprecision highp float;precision highp int;
struct VsOut { vec4 pos; vec2 uv;};out vec2 uv;
VsOut vs_impl(uint vi) { float _cse0 = ((float((vi & 1u)) * 4.0) - 1.0); float _cse1 = ((float((vi >> 1u)) * 4.0) - 1.0); return VsOut(vec4(_cse0, _cse1, 0.0, 1.0), vec2(((_cse0 * 0.5) + 0.5), ((_cse1 * 0.5) + 0.5)));}
void main() { VsOut _out = vs_impl(uint(gl_VertexID)); gl_Position = _out.pos; uv = _out.uv;}#version 300 esprecision highp float;precision highp int;
struct VsOut { vec4 pos; vec2 uv;};
struct DF64Vec2 { vec2 hi; vec2 lo;};layout(std140) uniform Uniforms { DF64Vec2 center; vec2 resolution; float zoom_exp; float fp64;} u;
uniform sampler2D _fp64;vec2 df64_twoSum(float a, float b);vec2 df64_quickTwoSum(float a, float b);vec2 df64_split(float a);vec2 df64_twoProd(float a, float b);vec2 df64_add(vec2 a, vec2 b);vec2 df64_sub(vec2 a, vec2 b);vec2 df64_mul(vec2 a, vec2 b);float df64_narrow(vec2 a);vec2 df64_twoSum(float a, float b) { float _cse0 = texelFetch(_fp64, ivec2(0, 0), 0).x; float _v0 = (a + b); float _v1 = (((_v0 * _cse0) - a) * _cse0); float _v2 = (((((a - ((_v0 - _v1) * _cse0)) * _cse0) * _cse0) * _cse0) + (b - _v1)); return vec2(_v0, _v2);}
vec2 df64_quickTwoSum(float a, float b) { float _cse0 = texelFetch(_fp64, ivec2(0, 0), 0).x; float _v0 = ((a + b) * _cse0); float _v1 = (b - ((_v0 - a) * _cse0)); return vec2(_v0, _v1);}
vec2 df64_split(float a) { float _cse0 = texelFetch(_fp64, ivec2(0, 0), 0).x; float _v0 = (a * (_cse0 * 4097.0)); float _v1 = ((_v0 * _cse0) - (_v0 - a)); float _v2 = ((a * _cse0) - _v1); return vec2(_v1, _v2);}
vec2 df64_twoProd(float a, float b) { float _v0 = (a * b); vec2 _v1 = df64_split(a); vec2 _v2 = df64_split(b); 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_add(vec2 a, vec2 b) { vec2 _v0 = df64_twoSum(a.x, b.x); vec2 _v1 = df64_twoSum(a.y, b.y); _v0.y = (_v0.y + _v1.x); _v0 = df64_quickTwoSum(_v0.x, _v0.y); _v0.y = (_v0.y + _v1.y); _v0 = df64_quickTwoSum(_v0.x, _v0.y); return _v0;}
vec2 df64_sub(vec2 a, vec2 b) { return df64_add(a, (-b));}
vec2 df64_mul(vec2 a, vec2 b) { vec2 _v0 = df64_twoProd(a.x, b.x); _v0.y = (_v0.y + (a.x * b.y)); _v0 = df64_quickTwoSum(_v0.x, _v0.y); _v0.y = (_v0.y + (a.y * b.x)); return df64_quickTwoSum(_v0.x, _v0.y);}
float df64_narrow(vec2 a) { return (a.x + a.y);}in vec2 uv;layout(location = 0) out vec4 _ret;
vec4 fs_de_impl(VsOut vo) { vec2 _licm0 = vec2(2.0, 0.0); bool _cse0 = (vo.uv.x < 0.5); vec2 _cse1 = vec2(u.center.hi.x, u.center.lo.x); vec2 _cse2 = vec2(u.center.hi.y, u.center.lo.y); vec2 _cse3 = vec2(0.0, 0.0); float _v0 = pow(10.0, (-u.zoom_exp)); float _v1 = (vo.uv.x * 2.0); float _v2 = (_v1 - (_cse0 ? 0.0 : 1.0)); float _v3 = ((_v2 - 0.5) * _v0); float _v4 = (((vo.uv.y - 0.5) * _v0) * ((u.resolution.y / u.resolution.x) * 2.0)); float _v5 = 0.0; float _v6 = 1.0; 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)) { if ((((_v9 * _v9) + (_v10 * _v10)) <= 1000000.0)) { float _v14 = ((((_v9 * _v11) - (_v10 * _v12)) * 2.0) + 1.0); _v12 = (((_v9 * _v12) + (_v10 * _v11)) * 2.0); _v11 = _v14; float _v15 = (((_v9 * _v9) - (_v10 * _v10)) + _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)); vec2 _v17 = df64_add(_cse2, vec2(_v4, 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)) { if ((df64_narrow(df64_add(df64_mul(_v18, _v18), df64_mul(_v19, _v19))) <= 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(df64_mul(_v18, _v18), df64_mul(_v19, _v19)), _v16); _v19 = df64_add(df64_mul(df64_mul(_v18, _v19), _licm0), _v17); _v18 = _v26; } } _v5 = df64_narrow(df64_add(df64_mul(_v18, _v18), df64_mul(_v19, _v19))); _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 _v31 = (exp(((-_v29) * 1.2)) * _v30); float _v32 = (exp(((-_v29) * 0.25)) * _v30); return 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);}
void main() { VsOut vo; vo.pos = gl_FragCoord; vo.uv = uv; _ret = fs_de_impl(vo);}{ "bindGroups": [ { "group": 0, "entries": [ { "group": 0, "binding": 0, "name": "u", "space": "uniform", "resourceKind": "uniform-buffer", "structName": "Uniforms" } ] } ], "uniforms": [ { "name": "Uniforms", "size": 32, "align": 16, "fields": [ { "name": "center", "type": "vec2<f64>", "offset": 0, "align": 16, "size": 16 }, { "name": "resolution", "type": "vec2<f32>", "offset": 16, "align": 8, "size": 8 }, { "name": "zoom_exp", "type": "f32", "offset": 24, "align": 4, "size": 4 }, { "name": "fp64", "type": "f32", "offset": 28, "align": 4, "size": 4 } ] } ], "storage": [], "entries": [ { "name": "vs", "stage": "vertex", "inputs": [ "u32" ], "output": "struct:VsOut" }, { "name": "fs_de", "stage": "fragment", "inputs": [ "struct:VsOut" ], "output": "vec4<f32>" } ], "overrides": []}
Run it locally:
npx tsx examples/print.ts fp64-mandelbrot-de
— prints the WGSL, GLSL, reflection. Source:
GitHub.
New to the IR? See Concepts; the full authoring + emit surface is in the API reference.