formula-manipulation.a68

     
   1  COMMENT
   2  
   3  @section Synopsis
   4  
   5  Symbolic computing in Algol 68.
   6  
   7  Based on example 11.10 in the Revised Report on Algol 68.
   8  
   9  COMMENT
  10  
  11  # Data structure #
  12  
  13  MODE FORMULA = UNION (REF CONST, REF VAR, REF DYADIC, REF MONADIC, REF CALL),
  14       DYADIC = STRUCT (FORMULA lhs, INT operator, FORMULA rhs),
  15       MONADIC = STRUCT (INT operator, FORMULA rhs),
  16       CALL = STRUCT (REF FUNCTION name, FORMULA parameter),
  17       FUNCTION = STRUCT (REF VAR bound var, FORMULA body),
  18       VAR = STRUCT (STRING name, NUMBER value),
  19       CONST = STRUCT (NUMBER value),
  20       NUMBER = LONG LONG REAL;
  21  
  22  # Access operators #
  23  
  24  OP VALUE = (REF CONST c) REF NUMBER: value OF c,
  25     VALUE = (REF VAR v) REF NUMBER: value OF v,
  26     NAME = (REF VAR v) REF STRING: name OF v,
  27     NAME = (REF CALL c) REF FUNCTION: name OF c,
  28     PARAMETER = (REF CALL c) REF FORMULA: parameter OF c,
  29     LHS = (REF DYADIC t) REF FORMULA: lhs OF t,
  30     RHS = (REF DYADIC t) REF FORMULA: rhs OF t,
  31     RHS = (REF MONADIC t) REF FORMULA: rhs OF t,
  32     OPERATOR = (REF DYADIC t) REF INT: operator OF t,
  33     OPERATOR = (REF MONADIC t) REF INT: operator OF t,
  34     BOUND = (REF FUNCTION f) REF REF VAR: bound var OF f,
  35     BODY = (REF FUNCTION f) REF FORMULA: body OF f;
  36  
  37  # Generate objects #
  38  
  39  PROC make dyadic = (FORMULA u, INT op, FORMULA v) REF DYADIC:
  40       NEW DYADIC := (u, op, v);
  41  
  42  PROC make monadic = (INT op, FORMULA v) REF MONADIC:
  43       NEW MONADIC := (op, v); 
  44  
  45  PROC make call = (REF FUNCTION name, FORMULA parameter) REF CALL:
  46       NEW CALL := (name, parameter);
  47  
  48  PROC make function = (REF VAR bound var, FORMULA body) REF FUNCTION:
  49       NEW FUNCTION := (bound var, body);
  50  
  51  PROC make var = (STRING name, NUMBER value) REF VAR: NEW VAR := (name, value);
  52  
  53  PROC make const = (NUMBER x) REF CONST: (NEW CONST c; VALUE c := x; c); 
  54  
  55  PROC zero = REF CONST: make const (0), 
  56       one = REF CONST: make const (1),
  57       two = REF CONST: make const (2),
  58       three = REF CONST: make const (3);
  59  
  60  # Basic routines and operators #
  61  
  62  PROC is var = (FORMULA f) BOOL: (f | (REF VAR): TRUE | FALSE);
  63  
  64  PROC is const = (FORMULA f) BOOL: (f | (REF CONST): TRUE | FALSE);
  65  
  66  PROC const value = (FORMULA f) NUMBER: (f | (REF CONST v): VALUE v);
  67  
  68  PROC var name = (FORMULA f) STRING: (f | (REF VAR v): NAME v);
  69  
  70  INT plus = 1, minus = 2, times = 3, divide = 4, up = 5;
  71  
  72  OP COPY = (FORMULA u) FORMULA:
  73     CASE u
  74     IN (REF CONST v): make const (VALUE v),
  75        (REF VAR v): make var (NAME v, VALUE v),
  76        (REF DYADIC v): make dyadic (COPY LHS v, OPERATOR v, COPY RHS v),
  77        (REF MONADIC v): make monadic (OPERATOR v, COPY RHS v),
  78        (REF CALL v): make call (NAME v, PARAMETER v) 
  79     ESAC;
  80  
  81  OP = = (REF FUNCTION u, v) BOOL: BOUND u = BOUND v AND BODY u = BODY v;
  82  
  83  OP = = (FORMULA a, b) BOOL:
  84     CASE a
  85     IN (REF CONST u): (b | (REF CONST v): VALUE u = VALUE v | FALSE),
  86        (REF VAR u): (b | (REF VAR v): NAME u = NAME v | FALSE), 
  87        (REF DYADIC u): 
  88           (b | (REF DYADIC v): LHS u = LHS v AND RHS u = RHS v AND OPERATOR u = OPERATOR v | FALSE),
  89        (REF MONADIC u): 
  90           (b | (REF MONADIC v): RHS u = RHS v AND OPERATOR u = OPERATOR v | FALSE),
  91        (REF CALL u):
  92           (b | (REF CALL v): NAME u = NAME v AND PARAMETER u = PARAMETER v | FALSE) 
  93     OUT FALSE
  94     ESAC;
  95  
  96  OP /= = (FORMULA a, b) BOOL: NOT (a = b);
  97  
  98  # Basic math #
  99  
 100  OP + = (FORMULA a) FORMULA: a;
 101  
 102  OP - = (FORMULA a) FORMULA: 
 103     (a = zero | a | make monadic (minus, a));
 104  
 105  OP + = (FORMULA a, b) FORMULA: 
 106     (a = zero | b |: b = zero | a | make dyadic (a, plus, b));
 107  
 108  OP - = (FORMULA a, b) FORMULA: 
 109     (b = zero | a | make dyadic (a, minus, b));
 110  
 111  OP * = (FORMULA a, b) FORMULA:
 112     IF a = zero OR b = zero
 113     THEN zero
 114     ELSE (a = one | b |: b = one | a | make dyadic (a, times, b))
 115     FI;
 116  
 117  OP / = (FORMULA a, b) FORMULA:
 118     IF a = zero AND NOT (b = zero)
 119     THEN zero
 120     ELSE (b = one | a | make dyadic (a, divide, b))
 121     FI;
 122  
 123  OP ^ = (FORMULA a, REF CONST b) FORMULA:
 124     IF a = one OR (b IS zero)
 125     THEN one
 126     ELSE (b IS one | a | make dyadic (a, up, b))
 127     FI;
 128  
 129  # Applications of above definitions: derivative, evaluation, simplification #
 130  
 131  PROC derivative = (FORMULA e, # with respect to # REF VAR x) FORMULA:
 132       # derivative a formula #
 133       CASE e
 134       IN (REF CONST): zero,
 135          (REF VAR v): (v = x | one | zero),
 136          (REF DYADIC f):
 137             CASE FORMULA u = LHS f, v = RHS f;
 138                  FORMULA deriv u = derivative (u, x), 
 139                          deriv v = derivative (v, x);
 140                  OPERATOR f
 141             IN deriv u + deriv v,
 142                deriv u - deriv v,
 143                u * deriv v + deriv u * v,
 144                (deriv u - f * deriv v) / v,
 145                (v | (REF CONST c): v * u ^ make const (VALUE c - 1) * deriv u)
 146             ESAC,
 147          (REF MONADIC f):
 148             CASE FORMULA v = RHS f;
 149                  FORMULA deriv v = derivative (v, x);
 150                  OPERATOR f
 151             IN + deriv v,
 152                - deriv v
 153             ESAC,
 154          (REF CALL c):
 155             BEGIN
 156                REF FUNCTION f = NAME c;
 157                FORMULA g = PARAMETER c;
 158                REF VAR y = BOUND f;
 159                REF FUNCTION deriv f = make function (y, derivative (BODY f, y));
 160                (make call (deriv f, g)) * derivative (g, x)
 161             END
 162       ESAC;
 163  
 164  PROC evaluate = (FORMULA e) NUMBER:
 165       # Value of a formula #
 166       CASE e
 167       IN (REF CONST c): VALUE c,
 168          (REF VAR v): VALUE v,
 169          (REF DYADIC f):
 170             CASE NUMBER u = evaluate (LHS f), v = evaluate (RHS f);
 171                  OPERATOR f
 172             IN u + v, u - v, u * v, u / v, u ^ SHORTEN SHORTEN ROUND v
 173             ESAC,
 174          (REF MONADIC f):
 175             CASE NUMBER v = evaluate (RHS f);
 176                  OPERATOR f
 177             IN v, - v
 178             ESAC,
 179          (REF CALL c):
 180             BEGIN
 181                REF FUNCTION f = NAME c;
 182                VALUE BOUND f := evaluate (PARAMETER c);
 183                evaluate (BODY f)
 184             END
 185       ESAC;
 186  
 187  OP SIMPLIFY = (FORMULA u) FORMULA:
 188     # Example simplifications - extend as you see fit #
 189     CASE u
 190     IN (REF CONST v): make const (VALUE v),
 191        (REF VAR v): make var (NAME v, VALUE v),
 192        (REF DYADIC v):
 193           IF FORMULA f = SIMPLIFY LHS v, g = SIMPLIFY RHS v;
 194              is const (f) THEF is const (g)
 195           THEN make const (evaluate (make dyadic (f, OPERATOR v, g)))
 196           ELIF OPERATOR v = plus
 197           THEN (f = g | make dyadic (two, times, f) | make dyadic (f, plus, g))
 198           ELIF OPERATOR v = minus
 199           THEN IF is const (f) THEF const value (f) = 0
 200                THEN make monadic (minus, g)
 201                ELSE (f = g | zero | make dyadic (f, minus, g))
 202                FI
 203           ELIF OPERATOR v = times
 204           THEN (is const (g) | make dyadic (g, times, f) |: f = g | make dyadic (f, up, two) | make dyadic (f, times, g))
 205           ELSE make dyadic (f, OPERATOR v, g)
 206           FI,
 207        (REF MONADIC v):
 208           IF FORMULA g = SIMPLIFY RHS v;
 209              is const (g)
 210           THEN make const (evaluate (make monadic (OPERATOR v, g)))
 211           ELSE make monadic (OPERATOR v, g)
 212           FI,
 213        (REF CALL v): make call (NAME v, PARAMETER v) 
 214     ESAC;
 215  
 216  # A small demonstration #
 217  
 218  OP FMT = (NUMBER x) STRING:
 219     (x = ENTIER x | whole (x, 0) | fixed (x, 0, 4));
 220  
 221  PROC write = (FORMULA e) VOID:
 222       CASE e
 223       IN (REF CONST c): print (FMT VALUE c),
 224          (REF VAR v): print (NAME v),
 225          (REF DYADIC f):
 226             BEGIN
 227                print ("(");
 228                write (LHS f);
 229                print ((OPERATOR f | " + ", " - ", " * ", " / ", " ^ "));
 230                write (RHS f);
 231                print (")")
 232             END,
 233          (REF MONADIC f):  
 234             BEGIN
 235                print (("(", (OPERATOR f | "+ ", "- ")));
 236                write (RHS f);
 237                print (")")
 238             END
 239       ESAC;
 240  
 241  PROC print and simplify = (FORMULA f) VOID:
 242       IF write (f);
 243          FORMULA g = SIMPLIFY f;
 244          f /= g
 245       THEN print (" = ");
 246            print and simplify (g)
 247       ELIF ~ is const (f)
 248       THEN print (" = ");
 249            print (FMT evaluate (f)) 
 250       FI;
 251  
 252  PROC demo = (FORMULA f, REF VAR z) VOID:
 253       BEGIN print ((new line, "f(x, y) = "));
 254             print and simplify (f);
 255             FORMULA df := derivative (f, z);
 256             print ((new line, "  df/d", NAME z, " = "));
 257             print and simplify (df)
 258       END;
 259  
 260  REF VAR x = make var ("x", -1), y = make var ("y", 1);
 261  
 262  printf (($"x = ", g, ", y = ", g$, FMT VALUE x, FMT VALUE y));
 263  demo (one + two * three, y);
 264  demo (x + x + zero * y, x);
 265  demo (x * two, x);
 266  demo (x * x * x + y * y, x);
 267  demo (x * x + two * x * y + y * y, x);
 268  demo (x + y / x, x);
 269  demo (x + y / x, y);
 270  new line (standout)
 271