Use automatic differentiation and distance estimate to make details less noisy

Cargo.toml

@@ -1,13 +1,15 @@
 [package]
 name = "multiman"
-version = "0.1.0"
+version = "0.2.0"
 authors = ["Jan Trefil <8711792+htrefil@users.noreply.github.com>"]
 edition = "2018"
 
-# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html
-
 [dependencies]
-structopt = "0.3.3"
-num-complex = "0.2.3"
-image = "0.22.3"
-num_cpus = "1.11.1"
+structopt = { version = "0.3.26", default_features = false }
+num-complex = "0.4.0"
+num_cpus = "1.13.1"
+
+[dependencies.image]
+version = "0.23.14"
+default_features = false
+features = ["png", "pnm", "webp", "bmp"]

README.md

@@ -1,33 +1,46 @@
 # multiman
 A parallel Mandelbrot and Julia set-style fractal renderer supporting custom user expressions.
-It leverages all available CPU cores to render images in parallel.
+It leverages automatic differentiation and a distance estimate to make the details less noisy.
+All available CPU cores are used to render images in parallel.
 
 ## Usage
 ```
-multiman 0.1.0
+multiman 0.2.0
 
 USAGE:
-    multiman <init> <iter> <width> <height> <output-path>
+    multiman <init> <first> <iter> <width> <height> <output-path>
 
 FLAGS:
     -h, --help       Prints help information
     -V, --version    Prints version information
 
 ARGS:
-    <init>           Initialization expression
+    <init>           Initialization expression (the value of the pixel)
+    <first>          Start of the iteration
     <iter>           Iteration expression
     <width>          Width of the image
     <height>         Height of the image
     <output-path>    Path of the resulting image
 ```
 
+In the expressions, the following variables are available:
+```
+w    The width of the image
+h    The height of the image
+x    The current x coordinate
+y    The current y coordinate
+c    The initial complex value of the current pixel
+z    The value of the previous iteration
+```
+They can be combined using arithmetic operations and exponentiation (using the symbol `^`).
+
 ## Examples
 ```
-multiman "(X / WIDTH * 2 - 1) + (Y / HEIGHT * 2 - 1.0) * i" "Z * Z + (0.1 + 0.65i)" 500 500 examples/1.png
+multiman "(x / w * 2 - 1) + (y / h * 2 - 1.0) * i" "c" "z * z + (0.1 + 0.65i)" 1000 1000 examples/1.png
 ```
 ![1](examples/1.png)
 
 ```
-multiman "0i" "Z * Z + (X * 2 / WIDTH - 1.5 + (Y * 2 / HEIGHT - 1) * i)" 500 500 examples/2.png
+multiman "x * 2 / w - 1.5 + (y * 2 / h - 1) * i" "0" "z * z + c" 1000 1000 examples/2.png
 ```
-![2](examples/2.png)
+![2](examples/2.png)

src/ast.rs

@@ -4,7 +4,7 @@ use std::fmt::{self, Debug, Formatter};
 #[derive(Clone, Debug, PartialEq)]
 pub enum ExprKind {
     Bin(BinOp, Box<Expr>, Box<Expr>),
-    Var(String),
+    Var(u8),
     Real(f64),
     Imag(f64),
 }
@@ -27,6 +27,7 @@ pub enum BinOp {
     Sub,
     Mul,
     Div,
+    Pow,
 }
 
 impl BinOp {
@@ -36,6 +37,7 @@ impl BinOp {
             TokenKind::Sub => BinOp::Sub,
             TokenKind::Mul => BinOp::Mul,
             TokenKind::Div => BinOp::Div,
+            TokenKind::Pow => BinOp::Pow,
             _ => return None,
         };
 
@@ -45,7 +47,7 @@ impl BinOp {
     pub fn precedence(&self) -> u32 {
         match self {
             BinOp::Add | BinOp::Sub => 1,
-            BinOp::Mul | BinOp::Div => 2,
+            BinOp::Mul | BinOp::Div | BinOp::Pow => 2,
         }
     }
 }

src/eval.rs

@@ -1,26 +1,34 @@
 use crate::ast::{BinOp, Expr, ExprKind};
 use crate::error::Error;
 use num_complex::Complex;
-use std::collections::HashMap;
 use std::ops::{Add, Div, Mul, Sub};
 
+#[derive(Default)]
 pub struct Context {
-    vars: HashMap<&'static str, Value>,
+    pub width: f64,
+    pub height: f64,
+    pub x: f64,
+    pub y: f64,
+    pub c: Value,
+    pub z: Value,
 }
 
 impl Context {
-    pub fn new() -> Context {
-        Context {
-            vars: HashMap::new(),
-        }
-    }
-
-    pub fn set(&mut self, name: &'static str, value: impl Into<Value>) {
-        self.vars.insert(name, value.into());
+    pub fn new() -> Self {
+        Default::default()
     }
 
-    pub fn update(&mut self, name: &str, value: impl Into<Value>) {
-        *self.vars.get_mut(name).unwrap() = value.into();
+    #[inline]
+    fn get(&self, name: u8) -> Option<Value> {
+        match name {
+            b'w' => Some(self.width.into()),
+            b'h' => Some(self.height.into()),
+            b'x' => Some(self.x.into()),
+            b'y' => Some(self.y.into()),
+            b'c' => Some(self.c.into()),
+            b'z' => Some(self.z.into()),
+            _ => None,
+        }
     }
 
     pub fn eval(&self, expr: &Expr) -> Result<Value, Error> {
@@ -34,12 +42,7 @@ impl Context {
                     BinOp::Sub => Ok(left - right),
                     BinOp::Mul => Ok(left * right),
                     BinOp::Div => {
-                        let zero = match right {
-                            Value::Real(num) => num == 0.0,
-                            Value::Complex(num) => num.re == 0.0 && num.im == 0.0,
-                        };
-
-                        if zero {
+                        if right.v == 0.0.into() {
                             return Err(Error {
                                 message: "Divide by zero",
                                 position: expr.position,
@@ -48,46 +51,71 @@ impl Context {
 
                         Ok(left / right)
                     }
+                    BinOp::Pow => {
+                        if right.v == 0.0.into() {
+                            return Err(Error {
+                                message: "Power by zero",
+                                position: expr.position,
+                            });
+                        }
+
+                        Ok(left.pow(right))
+                    }
                 }
             }
-            ExprKind::Var(ref var) => self.vars.get(var.as_str()).cloned().ok_or_else(|| Error {
+            ExprKind::Var(var) => self.get(var).ok_or_else(|| Error {
                 message: "Undefined variable",
                 position: expr.position,
             }),
-            ExprKind::Real(num) => Ok(Value::Real(num)),
-            ExprKind::Imag(num) => Ok(Value::Complex(Complex { re: 0.0, im: num })),
+            ExprKind::Real(num) => Ok(num.into()),
+            ExprKind::Imag(num) => Ok(Complex { re: 0.0, im: num }.into()),
         }
     }
 }
 
-#[derive(Clone, Copy, Debug)]
-pub enum Value {
-    Real(f64),
-    Complex(Complex<f64>),
+#[derive(Clone, Copy, Debug, Default)]
+pub struct Value {
+    /// Value of the expression
+    pub v: Complex<f64>,
+    /// Derivative of the expression
+    pub d: Complex<f64>,
+}
+
+impl Value {
+    fn pow(self, other: Value) -> Self {
+        let v = self.v.powc(other.v);
+        Value {
+            v,
+            d: (other.v * self.d / self.v + other.d * self.v.ln()) * v,
+        }
+    }
 }
 
 impl From<f64> for Value {
     fn from(value: f64) -> Value {
-        Value::Real(value)
+        Value {
+            v: value.into(),
+            d: 0.0.into(),
+        }
     }
 }
 
 impl From<Complex<f64>> for Value {
     fn from(value: Complex<f64>) -> Value {
-        Value::Complex(value)
+        Value {
+            v: value,
+            d: 0.0.into(),
+        }
     }
 }
 
 impl Add for Value {
     type Output = Value;
 
     fn add(self, other: Value) -> Value {
-        match (self, other) {
-            (Value::Real(a), Value::Real(b)) => Value::Real(a + b),
-            (Value::Complex(a), Value::Complex(b)) => Value::Complex(a + b),
-            (Value::Real(a), Value::Complex(b)) | (Value::Complex(b), Value::Real(a)) => {
-                Value::Complex(a + b)
-            }
+        Value {
+            v: self.v + other.v,
+            d: self.d + other.d,
         }
     }
 }
@@ -96,12 +124,9 @@ impl Sub for Value {
     type Output = Value;
 
     fn sub(self, other: Value) -> Value {
-        match (self, other) {
-            (Value::Real(a), Value::Real(b)) => Value::Real(a - b),
-            (Value::Complex(a), Value::Complex(b)) => Value::Complex(a - b),
-            (Value::Real(a), Value::Complex(b)) | (Value::Complex(b), Value::Real(a)) => {
-                Value::Complex(a - b)
-            }
+        Value {
+            v: self.v - other.v,
+            d: self.d - other.d,
         }
     }
 }
@@ -110,12 +135,9 @@ impl Mul for Value {
     type Output = Value;
 
     fn mul(self, other: Value) -> Value {
-        match (self, other) {
-            (Value::Real(a), Value::Real(b)) => Value::Real(a * b),
-            (Value::Complex(a), Value::Complex(b)) => Value::Complex(a * b),
-            (Value::Real(a), Value::Complex(b)) | (Value::Complex(b), Value::Real(a)) => {
-                Value::Complex(a * b)
-            }
+        Value {
+            v: self.v * other.v,
+            d: self.v * other.d + other.v * self.d,
         }
     }
 }
@@ -124,12 +146,9 @@ impl Div for Value {
     type Output = Value;
 
     fn div(self, other: Value) -> Value {
-        match (self, other) {
-            (Value::Real(a), Value::Real(b)) => Value::Real(a / b),
-            (Value::Complex(a), Value::Complex(b)) => Value::Complex(a / b),
-            (Value::Real(a), Value::Complex(b)) | (Value::Complex(b), Value::Real(a)) => {
-                Value::Complex(a / b)
-            }
+        Value {
+            v: self.v / other.v,
+            d: (self.d * other.v - other.d * self.v) / (other.d * other.d),
         }
     }
 }