1 /* MetaEthical.AI
 2    ==============
 3    Read the Introduction at:
 4    https://medium.com/@june_ku/formal-metaethics-and-metasemantics-for-ai-alignment-2e72533cad6d
 5 
 6    Developed by June Ku in the set theoretic programming language setlX v2.2.0
 7    and interspersed with considerable non-technical philosophical commentary.
 8      https://randoom.org/Software/SetlX/ 
 9      https://download.randoom.org/setlX/pc/archive/setlX_v2-2-0.binary_only.zip
10        - contains among other things, tutorial.pdf introduction and reference
11 
12    I am very grateful to Howard Nye for our work together on metaethics as well 
13    as many helpful comments and conversations during much of this development.
14 */
15 load('lib.stlx');
16 load('causal_model.stlx');
17 load('causal_markov_model.stlx');
18 load('decision_algorithm.stlx');
19 
20 /* MetaEthical.AI Utility Function
21    ===============================
22    Takes a causal markov model of the world w and set of brains bs in w and 
23    returns a metaethical.ai utility function. 
24 
25    Abstract: We construct a fully technical ethical goal function for AI by 
26    directly tackling the philosophical problems of metaethics and mental 
27    content. To simplify our reduction of these philosophical challenges into 
28    "merely" engineering ones, we suppose that unlimited computation and a 
29    complete low-level causal model of the world and the adult human brains in 
30    it are available.
31 
32    Given such a model, the AI attributes beliefs and values to a brain in two 
33    stages. First, it identifies the syntax of a brainís mental content by 
34    selecting a decision algorithm which is i) isomorphic to the brainís causal 
35    processes and ii) best compresses its behavior while iii) maximizing charity.
36    The semantics of that content then consists first in sense data that 
37    primitively refer to their own occurrence and then in logical and causal 
38    structural combinations of such content.
39 
40    The resulting decision algorithm can capture how we decide what to *do*, but 
41    it can also identify the ethical factors that we seek to determine when we 
42    decide what to *value* or even how to *decide*. Unfolding the implications 
43    of those factors, we arrive at what we should do. All together, this allows 
44    us to imbue the AI with the necessary concepts to determine and do what we 
45    should want it to do.  
46 */
47 metaethical_ai_u := procedure(w, bs) {
48   /* Given a set of brain models, associate them with the decision algorithms
49      they implement.  */
50   d := { [b, decision_algorithm.implemented_by(b, {})] : b in bs };
51   /* Then map each brain to its rational self's values (understood extensionally
52      i.e. cashing out the meaning of their mental concepts in terms of the 
53      world events they refer to).  */
54   ext_rufs := { [b, d[b].ext_ruf(w,b)] : b in bs };
55 
56   state_space := w.cm().poss_states();
57   // Possible social choice utility functions
58   psc_us := choices({ {[state, r_i] : r_i in [1..#state_space] }
59                       : state in state_space });
60   sc_u := sort_set(psc_us, procedure(psc_u1, psc_u2) {
61                              return meai_psc_dist(psc_u1, ext_rufs) <
62                                     meai_psc_dist(psc_u2, ext_rufs);
63                            })[1];
64   return sc_u;
65 };
66 
67   meai_psc_dist := procedure(psc_u, ext_rufs) {
68     return +/ { ord_u_dist(card_u_to_ord(ext_ruf), card_u_to_ord(psc_u)) ** 2
69                 : ext_ruf in ext_rufs };
70   };
71 
72 /* MIT License
73    -----------
74    Copyright 2019 June Ku 
75 
76    Permission is hereby granted, free of charge, to any person obtaining a copy 
77    of this software and associated documentation files (the "Software"), to deal
78    in the Software without restriction, including without limitation the rights 
79    to use, copy, modify, merge, publish, distribute, sublicense, and/or sell 
80    copies of the Software, and to permit persons to whom the Software is 
81    furnished to do so, subject to the following conditions:
82 
83    The above copyright notice and this permission notice shall be included in 
84    all copies or substantial portions of the Software.
85 
86    THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR 
87    IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 
88    FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE 
89    AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 
90    LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING 
91    FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS 
92    IN THE SOFTWARE.
93 */
   1 class decision_algorithm(i, o, p, cp, u, m, m2, n, r, t_s, t_o,
   2                          cached_f, cached_lf, cached_ref) {
   3 // Definition and Structure of Decision Algorithms 
   4 // ===============================================
   5 
   6   i := i; /* Set of input variables (or technically names, ie strings)
   7                which are a subset of range(p)
   8              e.g. { 'P(I_1)', 'P(I_2)', 'P(I_3)' } 
   9              i.e. the agent's sense data. 
  10              Generally, p should include it, e.g. { ['I_1', 'P(I_1)'], ... }
  11   Technical note: There's some abuse of notation here. Each string like 'I_1' 
  12   can depending on context be treated as:
  13     i) the name of a decision algorithm variable, e.g. as parents() of other 
  14        variables. Assignments of value to this variable would be a subjective 
  15        probability 0-1. That would form a local/partial decision state, which 
  16        could get mapped by an implementation functor f() from a brain state.
  17     ii) an individual constant term. A predicate expression in the domain of p 
  18         that has this string as an argument would be the prototypical context 
  19         illustrating this usage. Semantically, these would refer to the brain 
  20         events that implement them.  
  21     iii) a couple possibilities for this last option: 
  22       a) a 0-ary predicate or sentential expression. It takes no inputs and 
  23          forms a proposition on its own. 
  24       b) a 1-place predicate that takes the constant in ii and simply asserts 
  25          its existence, perhaps loosely inspired by Free Logic's E!. 
  26       In either case, its truth-maker/referent would be the same as option ii's.
  27       Treating the string as a sentential expression or forming an expression 
  28       with it as the argument of a connective would indicate this usage.  */
  29 
  30   o := o; /* Set of output variables
  31              e.g. { 'O_1', 'O_2', 'O_3' }
  32              i.e. the agent's motor output. 
  33              See also the technical note for i. */
  34 
  35   /* Function from logical formulas (setlx functors/strings/terms) to 
  36        probability variables
  37      e.g. { ['A', 'P(A)'], ['B', 'P(B)'], [And('A', 'B'), 'P(A & B)'], 
  38             [BoxArrow('A', 'B'), 'P(A []-> B)'] } */
  39   p := p;
  40 
  41   /* Function from ordered pairs [x, y] to conditional probability variables
  42      e.g. { [[x, y], 'P(x|y)'], ... } 
  43      or more concretely: { [['C', And('A', 'B')], 'P(C|A,B)'], ... } */
  44   cp := cp;
  45 
  46   /* Function from logical formulas to utility variables 
  47      e.g. { ['A', 'U(A)'], ... } 
  48      u has the same structure as p above */
  49   u := u;
  50 
  51   /* Set of memory variables 
  52      This is a catch-all for other cognitive states that don't fit into the 
  53      other sets of variables yet may be important to model, e.g. intermediate 
  54      states in a calculation performing inference or decision. One can think of
  55      it as roughly analogous to a computer's RAM.  */
  56   m := m;
  57 
  58   /* Set of memory variables for t_o 
  59      Same idea as m above but these only mediate between p + cp + u and o. 
  60      This helps ensure it is a genuinely decision algorithmic explanation.  */
  61   m2 := m2;
  62 
  63   /* List of accepted norms / higher-order utilities
  64      [ n_1, n_2, ... ]
  65      each n_i == { ['e', e_n_i], ['m2', m2_n_i], ['o', o_n_i] }
  66        e_n_i is a function from logical formulas to evaluation variables
  67          e.g. { ['A', 'E_n_i(A)'], ... }
  68          They are basically utility function variables but since they are
  69          higher-order, i.e. they influence normative judgments or prescriptions
  70          rather than motor output, we call them evaluation variables instead
  71          when we want to emphasize that they don't directly govern action.
  72        m2_n_i is a set of memory variables that only moderate p + e_n_i to 
  73          o_n_i
  74        o_n_i is a set of output variables that causally influence various u or
  75          e_n_j. In simplified models, we can make o_n_i trivially cause u or 
  76          e_n_i states by just being identical to them 
  77 
  78      These are essentially decision criteria that govern how we change our
  79      preferences or other higher-order norms. The way in which they govern is 
  80      analogous and mathematically isomorphic to how first-order utility 
  81      functions or preferences govern action. There, beliefs that some motor
  82      output would cause some preferred state to be achieved will tend to cause
  83      such output and therefore the desired state. In the case of these accepted
  84      norms, beliefs that the output variables (which we can also understand as 
  85      normative judgments or prescriptions) will cause or constitute a change 
  86      in values (which is what the exprs of the domain of e_n_i would refer to) 
  87      will tend to cause such outputs and therefore changes in value. 
  88 
  89      It is very important to note that what makes these preferences / utilities 
  90      higher-order is their causal or functional roles and not the content of 
  91      the expressions associated with the evaluation variables. One can have, in 
  92      our sense, first-order preferences which prefer states involving other
  93      preferences. Those might be referentially "higher-order" in some sense but 
  94      they would still be first-order in our sense so long as they govern
  95      actions rather than normative judgments or prescriptions.
  96 
  97      See http://www.metaethical.ai/norm_descriptivism.pdf for more.
  98   */
  99   n := n;
 100 
 101   // Internal states
 102   s := procedure() {
 103     if (n == {} || n == [] || n == om) {
 104       n_vs := {};
 105     } else {
 106       n_vs := +/ { range(n_i['e']) + n_i['m2'] + n_i['o'] : n_i in n };
 107     }
 108     return range(p) + range(cp) + range(u) + m + m2 + n_vs - i;
 109   };
 110 
 111   // Function from i + s() + o to ranges of values for those variables
 112   r := r;
 113 
 114   /* Higher order function from a function assigning values to all input and 
 115      internal states to a function assigning values to all next internal states 
 116        is -> s
 117   */
 118   t_s := t_s;
 119 
 120   /* Higher order function from a function assigning values to all probability 
 121      and utility variables as well as m2 to a function assigning values to 
 122      outputs 
 123        p u m2 -> o
 124   */
 125   t_o := t_o;
 126 
 127   /* State transitions
 128      Returns a higher-order function from functions assigning values to i + s 
 129      to functions assigning values to s + o. This is the intuitive combination
 130      of t_s and t_o. */
 131   t := procedure() {
 132     return { [    // a higher-order function which takes 
 133               is, // a function assigning values to i + s and returns a function
 134               t_s[is] + t_o[just_p_u_m2(is)]
 135                   // composed of the union of a func assigning values to s and 
 136                   // a func assigning values to o
 137              ] : is in domain(t_s) };
 138   };
 139 
 140   /* We explicitly cache f after computing it. This allows us to directly set
 141      its value when testing or when performing the brute force search would be
 142      prohibitively costly.  */
 143   cached_f := cached_f;
 144   cached_lf := cached_lf;
 145 
 146   /* ns as in the plural of n. 
 147      Takes 'e', 'm2' or 'o' and returns the corresponding variables from all n.
 148      e.g. ns('e') would return all e_n_i */
 149   ns := procedure(str) {
 150     if (str == 'e') {
 151       return +/ { range(n_i[str]) : n_i in n };
 152     } else {
 153       return +/ { n_i[str] : n_i in n };
 154     }
 155   };
 156 
 157   p_vars := procedure() {
 158     return range(p) + range(cp);
 159   };
 160 
 161   p_u_vars := procedure() {
 162     return p_vars() + range(u);
 163   };
 164 
 165   /* Takes a function (is) assigning values to input and internal state vars
 166      and returns its partial function with only vars in p, u and m2  */
 167   just_p_u_m2 := procedure(is) {
 168     return { [var, val] : [var, val] in is | var in p_u_vars() + m2 };
 169   };
 170 
 171 // Variable Relationships 
 172 // ----------------------
 173 
 174   /* Takes a variable and returns its parents, i.e. the variables that 
 175      are necessary to determine its next state in the state transition t(). */
 176   parents := procedure(y) {
 177     result := { z : z in (i + s()) | is_necessary_parent_of(z, y) };
 178     return result;
 179   };
 180 
 181     is_necessary_parent_of := procedure(z, y) {
 182       vs := {v : v in i + s() | v != z}; // all other variables but z
 183       f_vs := possible_states(vs, r);
 184       bool := exists(f_v in f_vs | vals_diff(f_v, y, z));
 185       return bool;
 186     };
 187 
 188     /* Check whether there exists an assignment of values to the other 
 189        variables such that there are still different values of y depending 
 190        on values of z.  */
 191     vals_diff := procedure(f_v, y, z) {
 192       vals := { t_s[ f_v + {[z,r_i]} ][y] : r_i in r[z] };
 193       return #vals > 1;
 194     };
 195 
 196   // Takes a variable and returns its children
 197   children := procedure(y) {
 198     return { v : v in s() + o | y in parents(v) };
 199   };
 200 
 201   /* Takes a variable v and a set of variables vs (with v usually in vs)
 202      and returns the first ancestors within any given branch of v's ancestors 
 203      that are not in vs */
 204   non_shared_ancestors := procedure(v, vs) {
 205     result := {};
 206     for (parent in parents(v)) {
 207       if (parent in vs) {
 208         result := result + non_shared_ancestors(parent, vs);
 209       } else {
 210         result := result + { parent };
 211       }
 212     }
 213     return result;
 214   };
 215 
 216   /* Takes a variable v and a set of variables vs (with v usually in vs)
 217      and returns the first descendants within any given branch of v's 
 218      descendants that are not in vs */
 219   non_shared_descendants := procedure(v, vs) {
 220     result := {};
 221     for (child in children(v)) {
 222       if (child in vs) {
 223         result := result + non_shared_descendants(child, vs);
 224       } else {
 225         result := result + { child };
 226       }
 227     }
 228     return result;
 229   };
 230 
 231 
 232 // Syntax and Implementation
 233 // =========================
 234 // See also better_explanation and the section starting with complexity.
 235 
 236   /* Function from a causal_markov_moodel of a brain to an implementation 
 237      mapping f which is a higher-order function from a function assigning 
 238      states to the brain's variables to a function assigning states to internal 
 239      state variables s() meeting certain criteria. 
 240   */
 241   f := procedure(b) {
 242     if (cached_f == om) {
 243       cached_f := compute_f(b);
 244     }
 245     return cached_f;
 246   };
 247 
 248   compute_f := cachedProcedure(b) {
 249     fs := { pf : pf in possible_fs(b) | commutes(pf, b) };
 250     /* Todo: Sort by complexity and return the simplest? 
 251        Or also factor in ambitiousness? Coverage of more states?  */
 252     if (#fs > 0) {
 253       return first(fs);
 254     } else {
 255       return {};
 256     }
 257   };
 258 
 259   /* Function from a causal_markov_model of a brain to the possible 
 260      implementation mappings f which map from a function assigning values to 
 261      the exogenous and endogenous variables of the causal_markov_model to 
 262      functions assigning values to the input, internal state, and output 
 263      variables. */
 264   possible_fs := procedure(b) {
 265     return function_space( b.poss_states(), poss_states() );
 266   };
 267 
 268   poss_states := cachedProcedure() {
 269     return possible_states(i + s() + o, r);
 270   };
 271 
 272   is_implemented_by := procedure(b) {
 273     return !(f(b) == {});
 274   };
 275 
 276   /* Implementation Functor
 277      ----------------------
 278      Check whether a possible implementation functor pf commutes. See f.
 279 
 280      This basically implements Chalmers' criterion for when some causal system 
 281      implements a computation: http://consc.net/papers/rock.html
 282 
 283      In category theory terms, we can say that f is an implementation functor.
 284      Given a category of brain states and the causal relations between them
 285      and another of decision states and cognitive state transitions between
 286      them, f would be a functor mapping i) brain states to decision states and
 287      ii) causal transitions to decision state transitions such that it doesn't
 288      matter whether one goes from a) a brain state to the brain state it causes 
 289      and then to the decision state f maps that brain state to or b) first from 
 290      the same original brain state to the decision state f says it implements 
 291      and then to the next decision state this decision algorithm transitions to. 
 292      Either way, you would end up at the same final decision state. 
 293   */
 294   commutes := procedure(pf, b) {
 295     /*                 pf
 296       bs_{b.uv_t_i-1} ------> ds_iso_t_i-1
 297        |                       |
 298        | b.response            | restrict_domain
 299        |                       | t()
 300        v                       v
 301       bs_{b.v_ti} ----------> ds_so_t_i
 302            +u  pf  restrict_domain
 303 
 304       Roughly check that t( pf ) == pf( b.response ) 
 305     */
 306     return forall(bs in domain(pf) |
 307        t()[restrict_domain(pf[bs], i + s())] ==
 308        restrict_domain(pf[
 309          b.response(b.v, bs) + { [u_var, first(b.r[u_var])] : u_var in b.u }
 310        ], s() + o)
 311      );
 312   };
 313 
 314   /* Takes a decision state assigning values to i + s() + o and returns a
 315      function assigning values to just s() + o  */
 316   drop_i := procedure(ds) {
 317     return { [var, val] : [var, val] in ds | !(var in i) };
 318   };
 319 
 320   /* Local / partial f
 321      Returns a function mapping partial causal states to partial decision states
 322      which commutes, i.e. it maps partial causal states to the partial decision 
 323      state which is forced by the complete mapping f.  */
 324   lf := procedure(b) {
 325     if (cached_lf == om) {
 326       cached_lf := compute_lf(b);
 327     }
 328     return cached_lf;
 329   };
 330 
 331   compute_lf := cachedProcedure(b) {
 332     /* Possible local fs 
 333        Function space from union of partial causal states to the union of 
 334        partial decision states.  */
 335     plfs := function_space(
 336       +/ { partial_functions(cs) : cs in b.poss_states() }, // causal states
 337       +/ { partial_functions(ds) : ds in poss_states() } // decision states
 338     );
 339     lfs := [ plf : plf in plfs | lf_commutes(plf, b) ];
 340     lfs := sort_list(lfs, procedure(plf1, plf2) { return #plf1 > #plf2; });
 341     return lfs[1];
 342   };
 343 
 344   /* Function from possible local function plf and a causal model of brain b 
 345      to boolean value whether it commutes. */
 346   lf_commutes := procedure(plf, b) {
 347     return forall(pcs in domain(plf) |  // pcs = partial causal state
 348       /* The set of complete decision states that is mapped by f from each 
 349          complete causal state compatible with the partial causal state is 
 350          equal to the set of complete decision states compatible with 
 351          the partial decision states mapped by the possible lf from the partial
 352          causal state. */
 353       { f(b)[ccs] : ccs in compatible_complete_states(pcs, b.poss_states()) } ==
 354       compatible_complete_states(plf[pcs], poss_states())
 355     );
 356   };
 357 
 358   /* Local / partial f inverse
 359   /* Returns a "function" from local / partial decision states to a set of 
 360      local / partial causal states, namely the inverse of lf.  */
 361   lf_inv := cachedProcedure(b) {
 362     return { [ds, { cs : cs in domain(lf(b)) | lf(b)[cs] == ds} ] :
 363              ds in range(lf(b)) };
 364   };
 365 
 366   // Helper function. Same as lf_inv but time subscripts are added to the b vars
 367   lf_inv_at := procedure(b, t_i) {
 368     return { [ds, { assign_time(cs, t_i) : cs in css }]
 369              : [ds, css] in lf_inv(b) };
 370   };
 371 
 372   /* Takes an event in b.cm() and returns the local decision state that the 
 373      brain state implements.  */
 374   lf_cm := procedure(b, b_event) {
 375     return lf(b)[drop_time(b_event)];
 376   };
 377 
 378   // Returns the decision state variables implemented by a local brain state.
 379   lf_cm_vs := procedure(b, b_event) {
 380     d_event := lf_cm(b, b_event);
 381     if (d_event == om) {
 382       return {};
 383     } else {
 384       return domain(d_event);
 385     }
 386   };
 387 
 388   /* Similar to above but returns the formula in the domain(p) that gets 
 389      mapped to the p var.  */
 390   lf_cm_p_inv := procedure(b, b_event) {
 391     vs := lf_cm_vs(b, b_event);
 392     if (#vs == 1) {
 393       return p_inv(first(vs));
 394     } else {
 395       return om;
 396     }
 397   };
 398 
 399   plf_to_pf := procedure(b, plf) {
 400     /* For each decision state var, construct the function from complete 
 401        causal states to the state of that variable. */
 402     pf_var := {
 403       [var, +/ { { [ccs, pds]
 404                    : ccs in compatible_complete_states(pcs, b.poss_states())
 405                  } : [pcs, pds] in plf | first(first(pds)) == var }
 406       ] : var in i + s() + o
 407     };
 408     return { [bs, +/ { (pf_var[var])[bs] : var in i + s() + o }]
 409              : bs in b.poss_states() };
 410   };
 411 
 412 
 413   // This is mostly just checking for syntactic validity or well-formedness.
 414   // See also all_lte.
 415   is_valid := procedure() {
 416     return m2s_are_valid() && no_e_in_u();
 417   };
 418 
 419   /* Ensures that m2s merely moderate between p + u and o.
 420      I'm no longer sure if this is necessary since we measure and punish 
 421      attributing irrationality but it also doesn't seem too bad of a 
 422      simplifying assumption for now.  */
 423   m2s_are_valid := procedure() {
 424     return forall (m2_i in m2 |
 425       forall (v in non_shared_ancestors(m2_i, m2) | v in p_u_vars()) &&
 426       forall (v in non_shared_descendants(m2_i, m2) | v in o)
 427     ) && forall (n_i in n |
 428       forall (m2_j_n_i in n_i['m2'] |
 429         forall (v in non_shared_ancestors(m2_j_n_i, n_i['m2']) |
 430                 v in p_vars() + range(n_i['e']) ) &&
 431         forall (v in non_shared_descendants(m2_j_n_i, n_i['m2']) |
 432                 v in n_i['o'])
 433       )
 434     );
 435   };
 436 
 437   /* This helps ensure higher-order ambitiousness by preventing the smuggling
 438      of evaluation functions within u. For each p expr of the form x []-> y,
 439      make sure it does not take the form that x-ing would bring about y which 
 440      the agent values (but x is not in the normal outputs o) and this 
 441      recognition causes x-ing. If it is of that form, we want a decision 
 442      algorithm that includes that explicitly within n.  
 443      Todo: Generalize to no_e_in_e and cover more cases, e.g. indirect 
 444      reference or equivalent expressions or influencing as grandparents.  */
 445   no_e_in_u := procedure() {
 446     return forall(expr in domain(p) | fct(expr) != "BoxArrow" || !e_in_u(expr));
 447   };
 448 
 449     e_in_u := procedure(expr) {
 450       x := args(expr)[1]; // like an n[i]['o']
 451       y := args(expr)[2];
 452       return y in domain(u) && !(x in o) && p[expr] in parents(x);
 453     };
 454 
 455 // Semantics of Mental Content 
 456 // ===========================
 457 
 458   /* Takes a set of pairs [w, b] and returns a function from possible worlds 
 459      to extensions / referents. 
 460        See: https://plato.stanford.edu/entries/two-dimensional-semantics/ 
 461   */
 462   intension := procedure(poss_wbs) {
 463     return { [wb[1], ref(wb[1], wb[2])] : wb in poss_wbs };
 464   };
 465 
 466   /* Reference
 467      ---------
 468      Reference function from causal markov model of world w and brain b 
 469      to a function which takes a brain and returns its set of truth-makers if 
 470      any, i.e. a set of poss events in the world or equivalently, functions
 471      assigning values to vars in the world. 
 472 
 473      There are a few ways that reference is grounded.
 474 
 475        1. Input variable states refer to the brain states implementing them.
 476           You can think of their content as something like "this sense datum
 477           is occurring" or perhaps even as just the pure demonstrative "this".
 478           It seems intuitive that these states can represent themselves since
 479           everything trivially carries information about itself. 
 480 
 481        2. The meanings of logical connectives (e.g. And, Or, Implies) and the 
 482           causal connective BoxArrow are grounded by the axioms that define 
 483           them. These meanings are enforced by punishing any attributions of 
 484           these connectives into credences that violate their axioms. See 
 485           incoherence. Cf inferential / conceptual role semantics:
 486             https://plato.stanford.edu/entries/content-narrow/#rol  
 487 
 488        3. The rest can be understood as a kind of ramsification. Other posited 
 489           objects/events, predicates and relations refer to those things
 490           (if any) which play the roles posited for them. 
 491              https://en.wikipedia.org/wiki/Ramsey%E2%80%93Lewis_method
 492              http://www.jimpryor.net/teaching/courses/mind/notes/ramseylewis.html
 493           We soften the requirement of strict truth by scoring by squared error.
 494           This amounts to a principle of charity favoring interpretations
 495           that make more beliefs turn out to be true (while also adhering to
 496           the constraints above).  
 497   */
 498   ref := procedure(w, b) {
 499     if (cached_ref == om) {
 500       cached_ref := compute_ref(w, b);
 501     }
 502     return cached_ref;
 503   };
 504 
 505   compute_ref := cachedProcedure(w, b) {
 506     /* Base expressions are strings or predicates. We'll build composite 
 507        expressions together with corresponding references out of combinations of
 508        these and logical / causal connectives. */
 509     base_exprs := { [expr, e_arity] : [expr, e_arity] in str_exprs() + preds()
 510                                     | !(p[expr] in i) };
 511     pbers := poss_base_expr_refs(w, base_exprs);
 512     poss_atom_refs := { i_ref(w, b) + base_ref
 513                         : base_ref in poss_base_refs(w, b, pbers) };
 514     poss_refs := [ full_ref(w, b, poss_atom_ref)
 515                    : poss_atom_ref in poss_atom_refs ];
 516     return sort_list(poss_refs, less_sq_err)[1];
 517   };
 518 
 519     /* Input vars refer to their own brain-state implementors. */
 520     i_ref := procedure(w, b) {
 521       return {
 522         [b_event, compatible_complete_states(b_event, w.cm().poss_states())]
 523         : b_event in b.actual_sync_events()
 524         | #lf_cm(b, b_event) == 1 && lf_cm(b, b_event) in poss_i_events()
 525       };
 526     };
 527 
 528     str_exprs := procedure() {
 529       return { [expr, 0] : expr in domain(p)
 530                          | !(p[expr] in i) && isString(expr) };
 531     };
 532 
 533     preds := procedure() {
 534       return { [fct(expr), arity(expr)]
 535                : expr in domain(p)
 536                | isTerm(expr) && !(fct(expr) in connectives()) };
 537     };
 538 
 539     // Set of choices of poss assignments of world w states/sets/tuples to expr
 540     poss_base_expr_refs := procedure(w, base_exprs) {
 541       return choices({ {expr} >< range_for_arity(e_arity, w.cm().poss_states())
 542                        : [expr, e_arity] in base_exprs });
 543     };
 544 
 545     /* Use poss_base_expr_refs to construct poss reference functions assigning
 546        reference to brain events rather than exprs.  */
 547     poss_base_refs := procedure(w, b, pbers) {
 548       return {
 549         +/ {
 550           // For each brain event that represents the given expr
 551           // Assign to it the reference of that expr according to pber
 552           { [b_event, pber[expr]] : b_event in b.actual_sync_events()
 553                                   | lf_cm_p_inv(b, b_event) == expr }
 554           // Do that for each base expr and collect them into a single function
 555           : expr in domain(pber) }
 556         // for each poss base expr reference func
 557         : pber in pbers
 558       };
 559     };
 560 
 561     // Full Reference
 562     // Much of the below builds on fairly standard first-order structures or 
 563     // models, e.g., https://en.wikipedia.org/wiki/First-order_logic#Semantics 
 564     full_ref := procedure(w, b, atom_ref) {
 565       /* Full ref maps each brain event which implements a credence to the 
 566          result of reducing the reference of the corresponding expression. */
 567       return { [b_event, rr(
 568                             lf_cm_p_inv(b, b_event), // expr it is credence in 
 569                             [ w, b, atom_ref,
 570                               // vars are all the same time so just take 1st
 571                               var_time(first(first(b_event))),
 572                               {} ] // empty context
 573                          ) ]
 574                : b_event in actual_p_events(b) };
 575     }; // end full_ref
 576 
 577       /* Reduce reference, takes an expression and an array of further
 578          parameters consisting of a world, brain, atomic reference function, 
 579          time and context and returns a set of truth-making w events or the set 
 580          of w events the expression is true in. The set of truthmakers is our 
 581          way of enforcing a canonical boolean algebraic form. The different 
 582          elements act as disjuncts of events and combining assignments of vals 
 583          to vars into a map acts as conjunction. It may seem counterintuitive 
 584          and take some getting used to but we use this canonical form even when
 585          there there is only one disjunct. */
 586       rr := procedure(expr, arr) {
 587         w := arr[1]; b := arr[2]; atom_ref := arr[3]; ti := arr[4];
 588         context := arr[5];  // context is for binding of variables, see Exists
 589         if (isTerm(expr)) {
 590           if (#args(expr) >= 1) { a1 := args(expr)[1]; }
 591           if (#args(expr) >= 2) { a2 := args(expr)[2]; }
 592         }
 593         switch {
 594           case expr in domain(context) :
 595             return context[expr];
 596           case p[expr] in i :
 597             return compatible_complete_states(
 598                      first({ b_event : b_event in pow(b.a[ti])
 599                                      | lf_cm_p_inv(b, b_event) == expr }),
 600                      w.cm().poss_states()
 601                    );
 602           case isString(expr) :
 603             return atom_ref[expr];
 604           case isSet(expr) :
 605             /* This shouldn't happen but just in case, return the set.
 606                Maybe throw an error instead? Or it's probably just rr(rr(.. */
 607             print('Warning: Expected an expr but received a set: ' + expr);
 608             return expr;
 609           case fct(expr) == 'At_t' :
 610             if (ti + a2 >= 0) { // a2 should already be negative
 611               return rr(a1, [w, b, atom_ref, ti + a2, context]);
 612             } else {
 613               return {};
 614             }
 615           case fct(expr) == 'And' :
 616             return rr(a1, arr) * rr(a2, arr);
 617           case fct(expr) == 'Or' :
 618             return rr(a1, arr) + rr(a2, arr);
 619           case fct(expr) == 'Not' :
 620             if (isSet(a1)) { // Set of states
 621               return w.cm().poss_states() - a1;
 622             } else {
 623               // String or functor, i.e. predicate or relation application
 624               // [~A] = [~[A]]
 625               return rr(Not(rr(a1, arr)), arr);
 626             }
 627           case fct(expr) == 'Implies' :
 628             return rr( Or(Not(a1), a2), arr );
 629           case fct(expr) == 'Forall' :
 630             // Vx Px = ~Ex ~Px
 631             return rr(Not(Exists(a1, Not(a2))), arr);
 632           case fct(expr) == 'Exists' :
 633             // Set of truth-makers for a2 when any {w_event} is subbed for a1
 634             return +/ { rr(a2, [ w, b, atom_ref, ti,
 635                                  context + { [a1, {w_event}] } ])
 636                         : w_event in range_for_arity(0, w.cm().poss_states()) };
 637           case fct(expr) == 'BoxArrow' : // A []-> B, subjunctive conditional
 638             // Antecedent
 639             if (isSet(a1)) {
 640               ant := a1;
 641             } else {
 642               ant := rr(a1, arr);
 643             }
 644             // Consequent
 645             if (isSet(a2)) {
 646               cons := a2;
 647             } else {
 648               cons := rr(a2, arr);
 649             }
 650             if (#ant > 1) { // Antecedent has disjuncts.
 651               // This or that causes the effect if this does or that does.
 652               return +/ { rr(BoxArrow(disjunct, cons), arr) : disjunct in ant };
 653             } else if (#ant == 1) {
 654               if (#cons > 1) {
 655                 // A thing causes this or that if it causes this or causes that.
 656                 return +/ { rr(BoxArrow(ant, disjunct), arr) : disjunct in cons };
 657               } else if (#cons == 1) {
 658                 if ( w.cm().response( domain(first(cons)), first(ant) ) ==
 659                      first(cons) ) {
 660                   /* If this aspect of the causal model is true, it's true in
 661                      all states. Fire causes smoke can be true even in a     
 662                      possible world where all fire is prevented. A more 
 663                      expressive language could substitute a more meaningful 
 664                      truthmaker here. */
 665                   return w.cm().poss_states();
 666                 } else {
 667                   return {};
 668                 }
 669               } else {
 670                 // Trivially true because everything causes the empty effect? 
 671                 return w.cm().poss_states();
 672               }
 673             } else {
 674               return {}; // Nothing is caused by nothingness?
 675             }
 676             return om;
 677           default : // Predicate or Relation
 678             tuple := [rr(arg, arr) : arg in args(expr)];
 679             if (tuple in atom_ref[fct(expr)]) {
 680               // If a predication is true, it's true in all states.
 681               return w.cm().poss_states();
 682             } else {
 683               return {};
 684             }
 685         } // end switch
 686       }; // end rr
 687 
 688     less_sq_err := procedure(w, b, ref_1, ref_2) {
 689       return sq_err(w, b, ref_1) < sq_err(w, b, ref_2);
 690     };
 691 
 692     sq_err := procedure(w, b, ref_i) {
 693       return +/ {
 694         ( true_p(w, ref_i, b_event) - first(range(lf_cm(b, b_event))) ) ** 2
 695         : b_event in actual_p_events(b)
 696       };
 697     };
 698 
 699     true_p := procedure(w, ref_i, b_event) {
 700       if (w.actual_state() in ref_i[b_event]) {
 701         return 1;
 702       } else {
 703         return 0;
 704       }
 705     };
 706 
 707   // Takes an expression (in a world, brain and time) and returns its referent.
 708   ref_expr := procedure(w, b, t_i, expr) {
 709     p_var := p[expr];
 710     ds := f(b)[drop_time(b.a[t_i])]; // decision state
 711     p_event := { [p_var, ds[p_var]] };
 712     b_event := first( lf_inv_at(b, t_i)[p_event] );
 713     return ref(w, b)[b_event];
 714   };
 715 
 716   // Returns the brain events that implement a placement of probabibility.
 717   actual_p_events := procedure(b) {
 718     return { b_event : b_event in b.actual_sync_events()
 719              |      #lf_cm_vs(b, b_event) == 1 &&         // just one var
 720                first(lf_cm_vs(b, b_event)) in range(p) }; // and it's a p_var
 721   };
 722 
 723   /* Function from causal_markov_model b to a function from times to 
 724      probability states, i.e. assignments of values to range(p) + range(cp).
 725      At each time, use f(b) to see what decision state it maps the 
 726      brain state to and collect the states for range(p) + range(cp). 
 727 
 728      Actual world b.a is a function from times to assignments of values to 
 729      variables at that time.
 730         e.g. b.a := { [1, assign_time({ [y, true] : y in b.uv() }, 1)], ...} */
 731   prob_states := procedure(b) {
 732     return { [ti, prob_states_t(b, ti) ] : ti in [1..b.n] };
 733   };
 734 
 735   /* Probability States at Time ti 
 736      Takes a brain b at a time ti and returns the decision state that is 
 737      implemented by that brain state but restricted to just the probability 
 738      variables.  */
 739   prob_states_t := procedure(b, ti) {
 740     bs := { [var, val] : [var, val] in drop_time(b.a[ti]) | var in b.uv() };
 741     ds := f(b)[bs];
 742     return { [var, val] : [var, val] in ds | var in p_vars() };
 743   };
 744 
 745 
 746 // Utility / Evaluation Functions
 747 // ==============================
 748 
 749   // Collect all utility/evaluation function vars with their output vars
 750   eo := procedure() {
 751     return {[u, o]} + { [n_i['e'], n_i['o']] : n_i in n };
 752   };
 753 
 754   /* Instrumental irrationality, including of evaluation as well as utility 
 755      functions. */
 756   instr_irrat := procedure(w, b) {
 757     return +/ { e_dist(w, b, ti) : ti in [1..(b.n - 1)] } / (b.n - 1);
 758   };
 759 
 760     /* Takes a brain b and time ti and returns the distance between the actual 
 761        outputs and the optimal outputs according to the utility / evaluation 
 762        functions and []-> beliefs.  */
 763     e_dist := procedure(w, b, ti) {
 764       /* These are all weighted equally right now, but we may want to 
 765          consider different weights. Should first-order utility be treated 
 766          differently? What about highest order? Does complexity matter?  */
 767       return +/ { e_i_dist(w, b, ti, e_i, o_i) : [e_i, o_i] in eo() } / #eo();
 768     }; // end e_dist
 769 
 770     // Similar to above but for a specific util/eval function and output vars
 771     e_i_dist := procedure(w, b, ti, e_i, o_i) {
 772       poss_o_i_states := possible_states(o_i, r);
 773       actual_d_st0 := f(b)[ drop_time(b.a[ti]  ) ];
 774       actual_d_st1 := f(b)[ drop_time(b.a[ti+1]) ];
 775       actual_o_i_state := { [o_i_j, actual_d_st1[o_i_j]] : o_i_j in o_i };
 776       sorted_o_i_states := sort_list( [x : x in poss_o_i_states],
 777         procedure(x, y) {
 778           return exp_u(w, b, ti, actual_d_st0, e_i, x) >
 779                  exp_u(w, b, ti, actual_d_st0, e_i, y);
 780         } );
 781       ix := index_of(actual_o_i_state, sorted_o_i_states);
 782       /* Cardinal measures seem too easy to abuse so we convert them to an 
 783          ordinal measure. */
 784       return ix / #sorted_o_i_states;
 785     }; // end e_i_dist
 786 
 787     // Expected Utility
 788     exp_u := procedure(w, b, ti, actual_d_st, e_i, o_i_state) {
 789       /* Using our reference function, find an expression that refers to 
 790          the given output state. */
 791       o_expr_b_events := {
 792         b_event : b_event in actual_p_events(b)
 793                 | var_time(first(b_event)) == ti &&
 794                   ref(w, b)[b_event] != om &&
 795                   exists(w_state in ref(w, b)[b_event] |
 796                     exists(b_ev in lf_inv(b)[o_i_state] |
 797                       restricted_eq(w_state, assign_time(b_ev, ti+1))
 798                     )
 799                   )
 800       };
 801       assert(o_expr_b_events != {}, "No o_expr_b_events found");
 802       o_expr := lf_cm_p_inv(b, first(o_expr_b_events));
 803       // Todo: Handle multiple exprs referring to o_i_state
 804       /* Is the above too indirect and complicated? There may be worries of 
 805          not capturing the correct mode of presentation? Could we "hardcode"
 806          something into the language instead? Like with input vars or
 807          memories using At_t. Or the simplest may be to allow it to overlap
 808          with input vars. That is actually still allowed by the above. 
 809            Is this not so easy though since it needs to refer to a *future* o 
 810          state?
 811 
 812          Add up the expectations of utility. Here we are interpreting our   
 813          utility vars as additive values. U(x & y) = U(x) + U(y). If, 
 814          however, there are exceptions, they can still be encoded by 
 815          having a util var associated with the conjunction. Then, if U(x)
 816          and U(y) are still defined, you may need to be careful to encode 
 817          only the intended difference from additive in the conjunction 
 818          rather than doubly counting its full value. 
 819            Here, we are assuming the P(o []-> u) beliefs reflect the agent's
 820          judgment of overall utility. If the agent is poor at aggregating 
 821          utilities, their final judgment might come apart from these 
 822          credences. There may also be issues with the timing of these
 823          and any aggregated judgments. While this may require much 
 824          technical work to thoroughly address, we hope it is at least 
 825          conceptually and philosophical clear what is to be done. */
 826       return +/ {  actual_d_st[ p[BoxArrow(o_expr, u_expr)] ] // P(o []-> u)
 827                  * actual_d_st[e_i_j] : [u_expr, e_i_j] in e_i }; // * u
 828     }; // end exp_u
 829 
 830 // Probabilistic Coherence 
 831 // =======================
 832 
 833   incoherence := procedure(b) {
 834     /* Synchronic irrationality
 835        Avg of distances to the closest coherent probability distribution 
 836        at each time normalized by max distance. */
 837     sync_irrat := avg({ p_dist(b, ti) / max_p_dist(b, ti) : ti in [1..b.n] });
 838     /* Diachronic irrationality
 839        Probabilities either stay the same or get updated to be close to what
 840        a (locally) ideally rational inference would result in.  */
 841     dia_irrat := 1 - p_conn_continuity([ f(b)[drop_time(b.a[ti])]
 842                                          : ti in [1..b.n] ]);
 843     return avg([sync_irrat, dia_irrat]);
 844   };
 845 
 846   /* Maximum Probability Distance 
 847      Find the probability distribution that is the farthest distance from the 
 848      given brain's probability distribution and return that distance.
 849      This is used to normalize probability distance scores.  */
 850   max_p_dist := procedure(b, ti) {
 851     poss_p_states := [poss_p : poss_p in possible_states(p_vars(), r)];
 852     p_states := prob_states_t(b, ti);
 853     farthest_pp := sort_list(poss_p_states, procedure(pp1,pp2) {
 854                                return prob_distance(p, cp, p_states, pp1) >
 855                                       prob_distance(p, cp, p_states, pp2);
 856                              })[1];
 857     return prob_distance(p, cp, p_states, farthest_pp);
 858   };
 859 
 860   /* Returns the set of possible causal models formed from the base variables
 861      in the domain of p. */
 862   poss_cms := procedure() {
 863     base_vars := [ p[expr] : expr in domain(p) |
 864                    isString(expr) || (isTerm(expr) && fct(expr) == 'At_t') ];
 865     vs := { v : v in base_vars };
 866     orderings := permutations(base_vars);
 867     return +/ { poss_cms_from_l(l, vs) : l in orderings };
 868   };
 869 
 870     poss_cms_from_l := procedure(l, vs) { // l is an ordering (list)
 871       // Each var's poss parents are sets containing the preceding vars in l.
 872       poss_v_parents := { [l[j], pow({ v : v in l[1..j-1] })] : j in [2..#l] };
 873       poss_v_parents[l[1]] := {{}};
 874       poss_parents := possible_states(vs, poss_v_parents);
 875 
 876       return +/ { poss_cms_from_p_pars(vs, poss_parents, p_pars)
 877                   : p_pars in poss_parents };
 878     };
 879 
 880     /* Given a set of parenthood relations, construct its set of possible 
 881        causal models. For each variable, construct the set of possible 
 882        functions from its parents. */
 883     poss_cms_from_p_pars := procedure(vs, poss_parents, p_pars) {
 884       cm_u := { v : v in vs | p_pars[v] == {} };
 885       cm_v := vs - cm_u;
 886       cm_r := { [v, {true, false}] : v in vs };
 887       poss_v_cm_fs := {
 888         [v, function_space(possible_states(p_pars[v], cm_r), {true, false})]
 889         : v in cm_v };
 890       poss_cm_fs := possible_states(cm_v, poss_v_cm_fs);
 891       return { causal_model(cm_u, cm_v, cm_r, cm_f, om) : cm_f in poss_cm_fs };
 892     };
 893 
 894   /* Takes a causal_model and returns a set of coherent probability 
 895      distributions compatible with the cm. Each distribution is given by
 896      a pair of assignments of values to p and cp.  */
 897   coh_pcms := procedure(cm) {
 898     // Coherent probability distributions for exogenous variables u
 899     coh_u_probs := { p_u : p_u in poss_p_us(cm.u)
 900                          | p_u_is_coherent(p_u, cm.u) };
 901     return { [probs(cm, p_u) + subj_conds(cm, p_u), cond_probs(cm, p_u)]
 902              : p_u in coh_u_probs };
 903   };
 904 
 905     poss_p_us := procedure(us) {
 906       return choices({ {state} >< p_range() : state in state_space(us) });
 907     };
 908 
 909       p_range := procedure() {
 910         return r[first(range(p))]; // assume all probs share the same range
 911       };
 912 
 913       /* For each prob variable, choose either it being true or false */
 914       state_space := procedure(vs) {
 915         return choices({ {[v, true], [v, false]} : v in vs });
 916       };
 917 
 918     p_u_is_coherent := procedure(p_u, base_vars) {
 919       poss_exprs := da().add_up_to(base_vars, 2 ** #base_vars);
 920       poss_probs := possible_states(poss_exprs,
 921                                     { [expr, p_range()] : expr in poss_exprs });
 922       poss_ord_pairs := poss_exprs >< poss_exprs;
 923       poss_cond_probs := possible_states(poss_ord_pairs,
 924                                 { [pair, p_range()] : pair in poss_ord_pairs });
 925       poss_p_cps := poss_probs >< poss_cond_probs;
 926       return exists(p_cp in poss_p_cps
 927                     | is_superset_of(p_cp, p_u) && superset_is_coherent(p_cp) );
 928     };
 929 
 930       is_superset_of := procedure(p_cp, p_u) {
 931         return p_u < first(p_cp);
 932       };
 933 
 934       superset_is_coherent := procedure(p_cp) {
 935         d := decision_algorithm.new();
 936         d.p := { [expr, 'P(' + str(expr) + ')']
 937                  : [expr, p_val] in first(p_cp) };
 938         d.cp := { [l, 'P(' + str(l[1]) + '|' + str(l[2]) + ')']
 939                   : [l, p_val] in last(p_cp) };
 940         pp := { ['P(' + str(expr) + ')',                   p_val]
 941                 : [expr, p_val] in first(p_cp) } +
 942               { ['P(' + str(l[1]) + '|' + str(l[2]) + ')', p_val]
 943                 : [l, p_val] in last(p_cp) };
 944         return d.is_coherent(pp);
 945       };
 946 
 947     /* Given p_u, a prob distr over states of the exogenous variables u, use 
 948        the causal_model to calculate the probabilities of endogenous vars vs. */
 949     probs := procedure(cm, p_u) {
 950       // For each state of exogenous vars u, use the causal model to fill in 
 951       // what the state of endogenous vars v would be. 
 952       p_uv := { [u_state + cm.response(cm.v, u_state), p_val]
 953                 : [u_state, p_val] in p_u };
 954       // For each poss event of the causal model, turn the event into an expr
 955       // and assign it probability equal to the sum of the probabilities of 
 956       // uv states in which the event is realized.
 957       return { [ func_to_expr(g_vs),
 958                  +/ { p_val : [uv_state, p_val] in p_uv | g_vs < uv_state }]
 959                : g_vs in cm.poss_events() };
 960     };
 961 
 962     /* Next we calculate the conditional probabilities for a given p_u */
 963     cond_probs := procedure(cm, p_u) {
 964       return +/ { calc_cond_prob(g, h) : [g,h] in cm.poss_events() ** 2 };
 965     };
 966 
 967       calc_cond_prob := procedure(cm, p_u, g, h) {
 968         if (compatible_values(g, h)) {
 969           a := func_to_expr(g);
 970           b := func_to_expr(h);
 971           c := func_to_expr(g + h);
 972           /* Using P(A|B) == P(And(A,B)) / P(B)  (Kolmogorov definition) */
 973           return { [ cp[[a,b]], probs(cm, p_u)[c] / probs(cm, p_u)[b] ] };
 974         } else {
 975           return 0;
 976         }
 977       };
 978 
 979     // Just set to 100% if the subjunctive conditional is true in cm
 980     subj_conds := procedure(cm, p_u) {
 981       return { [p[BoxArrow(func_to_expr(g), func_to_expr(h))], 1.0]
 982                 : [g,h] in cm.poss_events() ** 2
 983                 | cm.response(domain(h), g) == h };
 984     };
 985 
 986   // Note that the function it returns assigns true/false to prob vars and not
 987   // atomic expressions. You may want expr_to_events instead.
 988   expr_to_p_funcs := procedure(expr, dom) {
 989     return { { [p[term], bool] : [term, bool] in event }
 990              : event in expr_to_events(expr, dom) };
 991   };
 992 
 993   /* Coherent Causal Model Probabilities
 994      Coherent probability distributions over p and cp, including the 
 995      subjunctive conditionals, which can be seen as pieces of causal models. 
 996      We start with probability distributions over all causal models (from
 997      the set of base variables in p). Collecting the coherent distributions 
 998      over p and cp within each model, we then distribute the probability of
 999      the causal model to those distributions and sum them all up again. */
1000   coh_cm_ps := procedure() {
1001     // Possible probabilities of causal models
1002     poss_p_cms := choices({ {poss_cm} >< p_range : poss_cm in poss_cms() });
1003     /* The probabilities within any distribution must sum to 1. 
1004        (Each causal model considered is global and includes all base vars so 
1005        there is no worries of overlapping by one cm being a "superset" of 
1006        another.) */
1007     coh_p_cms := { p_cms : p_cms in poss_p_cms | +/ range(p_cms) == 1.0 };
1008 
1009     // Collect from each coherent probability distribution over causal models
1010     // the coherent probability distributions over p and cp that are allowed.
1011     return +/ { sum_p_cms(p_cms) : p_cms in coh_p_cms };
1012 
1013     /* Todo: Double-check whether it's necessary to expand the probability 
1014        distributions to non-canonical expressions in p? Now that prob_distance
1015        handles logical equivalences, it seems unnecessary.  */
1016   };
1017 
1018     /* Given a coherent probability distribution over causal models... */
1019     sum_p_cms := procedure(p_cms) {
1020       /* For each causal model and its probability, collect the possible
1021          choices of a coherent probability distribution from each, where
1022          those probability distributions have had the probability of the 
1023          causal model distributed over them.  */
1024       sets_of_p_cps := choices({
1025         { distr_p_cm(p_cp, p_cm) : p_cp in coh_pcms(cm) }
1026         : [cm, p_cm] in p_cms
1027       });
1028 
1029       /* ... For each choice of p_cps from each causal model, sum them up. */
1030       return { foldr(add_p_distr, {}, [p_cp : p_cp in p_cps])
1031                : p_cps in sets_of_p_cps };
1032     };
1033 
1034       /* Distribute probability of causal model.
1035          Takes p_cp, a set of pairs of coherent probability distributions over p
1036          and cp, and a probability of the causal model p_cm and then 
1037          distributes the probability of the causal model so that the probability
1038          of each expression is reduced by multiplying the p_cm. */
1039       distr_p_cm := procedure(p_cp, p_cm) {
1040         // p_cp == [cm_p, cm_cp]
1041         cm_p := p_cp[1];
1042         cm_cp := p_cp[2];
1043         return [{ [p_var, p_val * p_cm] : [p_var, p_val] in cm_p },
1044                 { [cp_l, cp_val * p_cm] : [cp_l, cp_val] in cm_cp }];
1045       };
1046 
1047       /* Given two pairs of probability distributions over p_vars
1048          add up the probabilities of expressions. */
1049       add_p_distr := procedure(p_cp1, p_cp2) {
1050         p1  := p_cp1[1];
1051         cp1 := p_cp1[2];
1052         p2  := p_cp2[1];
1053         cp2 := p_cp2[2];
1054         return [sum_p(p1, p2), sum_p(cp1, cp2)];
1055       };
1056 
1057         sum_p := procedure(p1, p2) {
1058           /* If one has a variable where the other does not have any variable
1059              with an equivalent expression, just go with the probability 
1060              defined by the one. */
1061           return { [v1, p1[v1]] : v1 in domain(p1)
1062                    | !exists(v2 in domain(p2) | equiv_p_expr(p, cp, v1, v2)) } +
1063                  { [v2, p2[v2]] : v2 in domain(p2)
1064                    | !exists(v1 in domain(p1) | equiv_p_expr(p, cp, v1, v2)) } +
1065           /* Otherwise, for any variables that have equivalent expressions,
1066              add up the probabilities for those variables. (Since they've 
1067              already been reduced in distr_p_cm, we don't need to average 
1068              the sum.)  */
1069                  +/ { { [v1, p1[v1] + p2[v2]],[v2, p1[v1] + p2[v2]]  }
1070                       : [v1, v2] in domain(p1) >< domain(p2)
1071                       | equiv_p_expr(p, cp, v1, v2) };
1072         };
1073 
1074           equiv_p_expr := procedure(p, cp, v1, v2) {
1075             e1 := p_cp_inv(p, cp, v1);
1076             e2 := p_cp_inv(p, cp, v2);
1077             return equiv_expr(e1, e2, all_atomic_formulas(domain(p)));
1078           };
1079 
1080   /* Takes a brain and time and finds the closest perfectly coherent 
1081      probability distribution and returns the distance to it. */
1082   p_dist := procedure(b, ti) {
1083     p_states := prob_states_t(b, ti);
1084     return prob_distance(p, cp, sorted_rat(p_states)[1], p_states);
1085   };
1086 
1087   /* Takes a prob state and returns a list of coherent probability 
1088      distributions sorted by how close they are p_states, closer ones first. */
1089   sorted_rat := procedure(p_states) {
1090     return sort_list( [ coh_p + coh_cp : [coh_p, coh_cp] in coh_cm_ps() ],
1091       /* A procedure measuring which of two probability states is closer to 
1092          the original one given.  */
1093       procedure(pp1, pp2) {
1094         return prob_distance(p, cp, pp1, p_states) <
1095                prob_distance(p, cp, pp2, p_states);
1096       }
1097     );
1098   };
1099 
1100   /* Takes a p or cp var in range(p) + range(cp) and returns its preimage 
1101      formula(s) */
1102   p_inv := procedure(v) {
1103     if (v in range(p)) {
1104       return inv(p)[v];
1105     } else {
1106       return inv(cp)[v];
1107     }
1108   };
1109 
1110   // Same as p_inv but with p and cp explicitly passed in due to setlx bug
1111   p_cp_inv := procedure(p, cp, v) {
1112     if (v in range(p)) {
1113       return inv(p)[v];
1114     } else {
1115       return inv(cp)[v];
1116     }
1117   };
1118 
1119   // Same as p_inv but if in cp, it returns the first expr instead of the pair
1120   p_var_inv := procedure(v) {
1121     if (v in range(p)) {
1122       return inv(p)[v];
1123     } else {
1124       return first(inv(cp)[v]);
1125     }
1126   };
1127 
1128 
1129   /* Obeys the axioms of probability.
1130        See especially p 45 (or 64 in pdf) of E.T. Jaynes, Probability Theory: 
1131        The Logic of Science. https://bayes.wustl.edu/etj/prob/book.pdf  */
1132   old_is_coherent := procedure(pp) {
1133     /* For easier handling, turn domain(p) into a datatype like domain(cp) 
1134        so that e.g. P(A) looks more like P(A|om) */
1135     force_list := procedure(v) {
1136       if (isList(p_inv(v))) {
1137         return p_inv(v);
1138       } else {
1139         return [p_inv(v), om];
1140       }
1141     };
1142     pp := { [force_list(v), num] : [v,num] in pp };
1143     obeys_axiom_3 := procedure(l) {
1144       if (isTerm(l[1]) && fct(l[1]) == "And") { // [And(a,b), c]
1145         [a,b] := args(l[1]);
1146         c := l[2];
1147         // P(a,b|c) = P(a|b,c) * P(b|c)
1148         return pp[l] == pp[ [a, And(b,c)] ] * pp[ [b,c] ];
1149       } else {
1150         return true;
1151       }
1152     };
1153     return forall(l in domain(pp) | pp[l] >= 0) &&
1154            // P(A|x) + P(~A|x) = 1
1155            // Todo: Simplify double negation to avoid infinite NOT() regression?
1156            forall(l in domain(pp) | pp[l] == 1 - pp[ [Not(l[1]), l[2]] ]) &&
1157            forall(l in domain(pp) | obeys_axiom_3(l));
1158   }; // end old_is_coherent
1159 
1160   /* From Lasky, Kathryn. "Axiomatic First-Order Probability" 
1161      http://ceur-ws.org/Vol-527/paper5.pdf 
1162 
1163       Laskey axioms
1164         P1. probabilities are between 0 and 1
1165         P2. probabilities in tautologies are 1
1166               or just logical truths
1167         P3. if probability in two statements both being true are 0, then 
1168               probability of one or the other is equal to adding probability 
1169               in each
1170         P4. Bayesian conditioning
1171         P5. interchangeability of logically equivalent statements
1172 
1173      Todo: Double-check if it's necessary to handle undefined/om probability 
1174            values. It seems p_u_is_coherent checks its coherence by seeing if
1175            a superset is perfectly coherent so I think that strategy should 
1176            work. But there may be some slight technical issues making sure 
1177            we're allowing or disallowing om in the right places there.
1178   */
1179   is_coherent := procedure(pp) {
1180     pp := { [force_list(v), num] : [v,num] in pp };
1181     return forall(l in domain(pp) |
1182       // P1: Probabilities are between 0 and 1.
1183       0 <= pp[l] && pp[l] <= 1 &&
1184       // P2: Logical truths have probability 1.
1185       (!is_logical_truth(l[1]) || pp[l] == 1) &&
1186       forall(l2 in domain(pp) |
1187         ( l[2] != l2[2] || (
1188           // z is shared in P(_|z) && P(_|z)
1189 
1190           // P3: If P(x^y|z) = 0, then P(x v y|z) = P(x|z) + P(y|z)
1191           ( pp[[And(l[1], l2[1]), l[2]]] != 0 ||
1192             pp[[ Or(l[1], l2[1]), l[2]]] == pp[l] + pp[l2] ) &&
1193 
1194           // P4: Bayesian conditioning. P(x^y|z) = P(y|z) * P(x|y,z)  
1195           ( fct(l[1]) != "And" || !(l2[1] in args(l[1])) ||
1196             forall(l3 in domain(p) |
1197               !(l3[1] in args(l[1])) || l3[1] == l2[1] ||
1198               l3[2] != And(l2[1], l[2]) || pp[l] == pp[l2] * pp[l3]
1199             ) ) &&
1200 
1201           // P5: Interchangeability of equivalent expressions. Part 1
1202           //     If x<->y, then Vz P(x|z) = P(y|z) 
1203           ( !equiv_expr(l[1], l2[1], atomic_formulas(domain(p))) ||
1204             pp[l] == pp[l2] )
1205         ) ) /* end || */ && (
1206           // P5 Part 2: If x <-> y, then Vz P(z|x) = P(z|y) 
1207           !equiv_expr(l[2], l2[2], atomic_formulas(domain(p))) ||
1208           l[1] != l2[1] || pp[l] == pp[l2]
1209         )
1210       ) // end forall l2
1211     ); // end forall l
1212   }; // end is_coherent
1213 
1214     /* For easier handling, turn domain(p) into a datatype like domain(cp) 
1215        so that e.g. P(A) looks more like P(A|om) */
1216     force_list := procedure(v) {
1217       if (isList(p_inv(v))) {
1218         return p_inv(v);
1219       } else {
1220         return [p_inv(v), om];
1221       }
1222     };
1223 
1224   /* Takes two functions mapping p & cp vars to probabilities and computes
1225      the difference. See Staffel, Julia. Measuring the Overall Incoherence of 
1226      Credence Functions. Synthese, 2015. DOI 10.1007/s11229-014-0640-x 
1227 
1228      Technical Note: Takes p and cp (which should just be the usual this.p and 
1229      this.cp) but needs to be passed explicitly because of a setlx bug: The 
1230      second and subsequent calls to an instance procedure within a nested 
1231      procedure is not able to access instance variables / procedures. 
1232 
1233      See ambitiousness for a more general workaround. There and elsewhere we
1234      assign *this* (the decision algorithm instance) to *self* and then pass
1235      in self as the first parameter. The receiving instance procedure then
1236      makes sure to call any instance variables and procedures on self rather 
1237      than their implicit and buggy this. The result may be likened to Python's
1238      way of handling instance methods.
1239   */
1240   prob_distance := procedure(p, cp, pp1, pp2) {
1241     /* First pass: (but doesn't handle undefineds and logical equivalents)
1242     return +/ { abs(pp1[v] - pp2[v]) : v in domain(pp1) | v in domain(pp2) };
1243     // Second pass: (but doesn't handle double-counting)
1244     return +/ { abs(pp1[v] - pp2[v2]) : [v, v2] in domain(pp1) >< domain(pp2)
1245                                       | v == v2 || equiv_p_expr(v, v2) };     
1246        Consider three people where the "true" probability is 90%
1247          Person 1: A & B 80%, B & A 70%
1248          Person 2: A & B 80%
1249          Person 3: A & B 80%, B & A 80%
1250        Averaging probabilities in logical equivalents favors person 2 and 3 
1251        over 1 but not 2 over 3.
1252     */
1253     return +/ { abs( pp_equiv(p, cp, pp1)[p_class] -
1254                      pp_equiv(p, cp, pp2)[p_class] )
1255                 : p_class in p_equiv_classes(p, cp)
1256                 | pp_equiv(p, cp, pp1)[p_class] != om &&
1257                   pp_equiv(p, cp, pp2)[p_class] != om };
1258   };
1259 
1260     // Equivalence classes of probability variables based on logical equivalence
1261     p_equiv_classes := procedure(p, cp) {
1262       return { { x : x in range(p) + range(cp)
1263                    | equiv_p_expr(p, cp, p_var, x) }
1264                : p_var in range(p) + range(cp) };
1265     };
1266 
1267     // Converts a possibility probability distribution pp to one that ranges
1268     // over equivalence classes of logically equivalent p vars instead of 
1269     // directly over the p vars.
1270     pp_equiv := procedure(p, cp, pp) {
1271       return { [p_class, avg_p_class(pp, p_class)]
1272                : p_class in p_equiv_classes(p, cp)
1273                | avg_p_class(pp, p_class) != om };
1274     };
1275 
1276       // Return the average of the values pp assigns to the p vars in p_subset
1277       avg_p_class := procedure(pp, p_subset) {
1278         p_vals := [ pp[p_var] : p_var in domain(pp) | p_var in p_subset ];
1279         if (#p_vals == 0) {
1280           return om;
1281         } else {
1282           return avg(p_vals);
1283         }
1284       };
1285 
1286 // Helpers 
1287 // -------
1288 
1289   poss_io_states := cachedProcedure() {
1290     return possible_states(i + o, r);
1291   };
1292 
1293   poss_s_states := cachedProcedure() {
1294     return possible_states(s(), r);
1295   };
1296 
1297   poss_m_states := cachedProcedure() {
1298     return possible_states(m, r);
1299   };
1300 
1301   poss_m2_states := cachedProcedure() {
1302     return possible_states(m2, r);
1303   };
1304 
1305   poss_ms_states := cachedProcedure() {
1306     return possible_states(m + m2, r);
1307   };
1308 
1309   poss_i_states := cachedProcedure() {
1310     return possible_states(i, r);
1311   };
1312 
1313   poss_is_states := cachedProcedure() {
1314     return possible_states(i + s(), r);
1315   };
1316 
1317   poss_o_states := cachedProcedure() {
1318     return possible_states(o, r);
1319   };
1320 
1321   poss_u_states := cachedProcedure() {
1322     return possible_states(range(u), r);
1323   };
1324 
1325 
1326 // Best Compression
1327 // ================
1328 
1329   /* Complexity of Intentional Explanation
1330      -------------------------------------
1331      This can be seen as a more realist and more formal version of Dennett's
1332      Intentional Stance. 
1333        See e.g. https://sites.google.com/site/minddict/intentional-stance-the
1334 
1335      On a first pass, the best intentional explanation can be taken to be the
1336      one which best compresses the brain's behavior.
1337        See better_explanation for further considerations.
1338 
1339      Also see https://en.wikipedia.org/wiki/Minimum_description_length:
1340        "The minimum description length (MDL) principle is a formalization of 
1341         Occam's razor in which the best hypothesis (a model and its parameters) 
1342         for a given set of data is the one that leads to the best compression 
1343         of the data." 
1344 
1345      Since b and d already contain information on their transitions from any 
1346      possible state, we don't need to run this for each time step.
1347   */
1348   complexity := cachedProcedure(b) {
1349     // Complexity of b given d and f + Complexity of d and f
1350     raw_score := k_given(b.t_to_str(), df_to_str(b)) + k(df_to_str(b));
1351     /* Now normalize. Worst case is it doesn't help encode b any better and 
1352        it's just as complex as b. Best case is if it could take 0 bits.  */
1353     worst_score := 2 * k(b.t_to_str());
1354     return 1 - raw_score / worst_score;
1355   };
1356 
1357   to_s := cachedProcedure() {
1358     return str([i, p, cp, u, m, m2, n, o, r, t()]);
1359   };
1360 
1361   f_to_s := procedure(b) {
1362     return str(f(b));
1363   };
1364 
1365   df_to_str := procedure(b) {
1366     return to_s() + "@" + f_to_s(b);
1367     // '@' is just a separator, an otherwise unused character
1368   };
1369 
1370   /* Ambitiousness
1371      -------------
1372      If we only optimized for the best compression of the brain, we face the 
1373      following problem. Suppose the syntax we use is somewhat inefficient. 
1374      Then, the smallest K(b|df) + K(df) might be produced by an empty d with 
1375      all the work being done by a turing machine which generates b using 
1376      essentially the same decision algorithm we would have wanted to use as d 
1377      but with a more efficient syntax. Or it may only put some into the 
1378      official d and sneak in the rest in the turing machine. It seems what
1379      we'd need is something like a not-a-decision-algorithmic-explanation
1380      predicate to apply to the turing machine. Or we want d to be ambitious, 
1381      ideally packing as much of a decision-algorithmic explanation as possible
1382      or at least positively weighing its ambitiousness against potential 
1383      increases in complexity.
1384      
1385      For a given d, we consider the auxiliary hypothesis h, which together w/
1386      d best explains b. We then look for a decomposition of h into [h_i, h_j] 
1387      where [h_i, h_j] is similar to h, and [h_i, h_j] makes d least ambitious, 
1388      meaning that the d* which is a superset of d and maximizes coherence and 
1389      similarity to d + h_i, improves coherence the least with the greatest 
1390      complexity cost. d's ambitiousness in the worst case decomposition is d's 
1391      ambitiousness. In other words, the more there is a d* which is highly 
1392      similar to d + some h_i and very coherent, the less ambitious d is. It 
1393      would leave out too much that could and should be explained using a 
1394      decision algorithmic explanation. If d is very ambitious, there's more 
1395      like just noise left for the auxiliary hypothesis or you'd have to stray 
1396      farther from the auxiliary hypothesis h to subsume it under a decision 
1397      algorithmic explanation.
1398   */
1399   ambitiousness := procedure(b) {
1400     poss_hs := all_strs_lte(#b.t_to_str()) >< all_strs_lte(#b.t_to_str());
1401     self := this; // workaround for setlx bug described in prob_distance
1402     poss_hs := sort_set(poss_hs,
1403       /* Procedure for worse decomposition. Calling it worse rather than
1404          better because we want the [h1, h2] which represents the worst case
1405          scenario as far as ambitiousness goes.  
1406       */
1407       procedure(hs_i, hs_j) {
1408         return hs_score(self, b, hs_i) > hs_score(self, b, hs_j);
1409       });
1410     hs := poss_hs[1];
1411     return ambitiousness_with_hs(b, hs);
1412   };
1413 
1414     hs_score := procedure(self, b, hs_k) {
1415       return self.hs_similarity(b, hs_k) - self.ambitiousness_with_hs(b, hs_k);
1416     };
1417 
1418     hs_similarity := procedure(b, hs_k) { // hs_k == [h1, h2]
1419       return ( k_given(h(b), str(hs_k)) + k_given(str(hs_k), h(b)) ) * -1;
1420       // or ** -1 ?
1421     };
1422 
1423   // Stringified auxiliary hypothesis which together with d explains b.t_to_str
1424   h := procedure(b) {
1425     poss_hs := all_strs_lte(#b.t_to_str());
1426     self := this; // workaround for setlx bug described in prob_distance
1427     poss_hs := sort_set(poss_hs, procedure(h_i, h_j) { // better h explanation
1428                  return h_complexity(self, b, h_i) < h_complexity(self, b, h_j);
1429                });
1430     return poss_hs[1];
1431   };
1432 
1433     h_complexity := procedure(self, b, h_k) {
1434       return k_given(b.t_to_str(), self.df_to_str(b) + "|" + h_k);
1435     };
1436 
1437   ambitiousness_with_hs := procedure(b, hs_k) {
1438     poss_d_stars := { poss_d_star : poss_d_star in supersets_lte(b)
1439                       | poss_d_star.is_valid() &&
1440                         poss_d_star.is_implemented_by(b) };
1441     self := this; // workaround for setlx bug described in prob_distance
1442     poss_d_stars := sort_set(poss_d_stars,
1443       procedure(d_star_i, d_star_j) { // better_d_star
1444         return d_star_sim(self, b, hs_k, d_star_i) - d_star_i.incoherence(b) >
1445                d_star_sim(self, b, hs_k, d_star_j) - d_star_j.incoherence(b);
1446       });
1447     d_star := poss_d_stars[1];
1448 
1449     delta_incoh := ( d_star.incoherence(b) - incoherence(b) ) / incoherence(b);
1450     delta_compl := (k(d_star.df_to_str(b)) - k(df_to_str(b))) / k(df_to_str(b));
1451     /* More incoherent -> more ambitious
1452          because it's like less of real Decision Algorithmic stuff left in h1
1453        Larger complexity -> less ambitious 
1454          because there's more Decision Algorithmic stuff left in hs_1 & not in d
1455          but maybe replace complexity with more like raw size? #states?  */
1456     return delta_incoh - delta_compl;
1457   };
1458 
1459     d_star_sim := procedure(self, b, hs_k, d_star_k) {
1460       raw_sim := (  // roughly, K(d*|d,h1) + K(d,h1|d*_s)
1461         k_given( d_star_k.df_to_str(b), self.df_to_str(b) + "|" + hs_k[1]) +
1462         k_given( self.df_to_str(b) + "|" + hs_k[1],
1463                  d_star_k.simplest_version(b).df_to_str(b) )
1464       );
1465       /* Normalize. Best similarity is if 0 extra bits are required. 
1466          Worst is K(d*) + K(d); the given bits didn't help at all. */
1467       worst_sim := k(d_star_k.df_to_str(b)) + k(self.df_to_str(b));
1468       return 1 - raw_sim / worst_sim;
1469     };
1470 
1471   /* Returns decision algorithms with # of possible states <= #b.poss_states
1472      which are "supersets" of self. It has additional variables and therefore
1473      states, but when considering just the variables shared with self, it is 
1474      isomorphic to self.  */
1475   supersets_lte := procedure(b) {
1476     return { d : d in decision_algorithm.all_lte(#b.poss_states())
1477                | is_subset_of(d) };
1478   };
1479 
1480   // Returns whether this decision algorithm is a "subset" of the given d
1481   is_subset_of := procedure(d) {
1482     // i, o, p, cp, u, m, m2, n, r, t_s, t_o
1483     var_pairs := { [i, d.i], [o, d.o], [p, d.p], [cp, d.cp],
1484                    [u, d.u], [m, d.m], [m2, d.m2] };
1485     vars_are_subsets := forall([vs, d_vs] in var_pairs | vs < d_vs);
1486     if (!vars_are_subsets || #n > #d.n) { return false; }
1487 
1488     n_are_subsets := forall(j in [1..#n] |
1489                        forall(ltr in {'e','m2','o'} |
1490                          n[j][ltr] < d.n[j][ltr] ));
1491     if (!n_are_subsets) { return false; }
1492 
1493     /* For all transitions from a i+s state to a s+o state in self, it is the
1494        case that any transition in d starting from a superset of i+s ends in
1495        a d state that is a superset of the s+o state.  */
1496     return forall([is, so] in t() |
1497              forall(ccs in compatible_complete_states(is, d.poss_is_states()) |
1498                compat(so, d.t()[ccs])
1499              )
1500            );
1501   };
1502 
1503   is_simplest_version := procedure(b) {
1504     return is_equal_to(simplest_version(b));
1505   };
1506 
1507   /* Using simplest_version handles worries with using m as arbitrary data
1508      to smuggle in unlabelled decision algorithmic information. Only the 
1509      decision algorithmic aspects of a decision algorithm should do the 
1510      explaining. This should also handle worries of higher-order utility
1511      information being smuggled in within the memory variables.  */
1512   simplest_version := procedure(b) {
1513     alt_ds := [ d : d in decision_algorithm.all(b)
1514                   | is_isomorphic_to(d) && d.is_implemented_by(b) ];
1515     alt_ds := sort_list(alt_ds, decision_algorithm.simpler);
1516     return alt_ds[1];
1517   };
1518 
1519   is_equal_to := procedure(d) {
1520     l1 := [  i,   o,   p,   cp,   u,   m,   m2,   n,   r,   t_s,   t_o];
1521     l2 := [d.i, d.o, d.p, d.cp, d.u, d.m, d.m2, d.n, d.r, d.t_s, d.t_o];
1522     return l1 == l2;
1523   };
1524 
1525   is_isomorphic_to := procedure(d) {
1526     /* Everything but the memory variables is the same 
1527        And there is an isomorphism between m+m2 states and d.m+d.m2 states
1528        Todo: Also add n minus m2s in n to simple equality comparison.
1529     */
1530     l1 := [  i,   o,   x,   p,   cp,   u,   r,   t_s,   t_o];
1531     l2 := [d.i, d.o, d.x, d.p, d.cp, d.u, d.r, d.t_s, d.t_o];
1532     return l1 == l2 && exists(fm in possible_fms(d) | fm_commutes(fm, d));
1533   };
1534 
1535     // Todo: Also do higher-order m2s
1536     possible_fms := procedure(d) {
1537       /* Construct the space of functions from states of m+m2 to d.m+d.m2,
1538          respecting the m / m2 boundary */
1539       fm_ms  := {
1540         { [ m_st + m2_st,  fm_m[m_st] + fm_m2[m2_st] ]
1541           : [m_st, m2_st] in poss_m_states() >< poss_m2_states() }
1542         : [fm_m, fm_m2] in
1543           function_space(possible_states(m,  r), possible_states(d.m,  d.r)) ><
1544           function_space(possible_states(m2, r), possible_states(d.m2, d.r))
1545       };
1546 
1547       // For convenience map everything else onto its counterpart in d
1548       fms := {
1549         { [iso, apply_fm_m(fm_m, iso) ] : iso in possible_states(i+s()+o, r) }
1550         : fm_m in fm_ms
1551       };
1552       return fms;
1553     };
1554 
1555       // Todo: Also do higher-order m2s.
1556       apply_fm_m := procedure(fm_m, iso) {
1557         iso_m := { [var, val] : [var, val] in iso | var in m || var in m2 };
1558         return { [v, iso[v]] : v in i + s() + o - m - m2 } + fm_m[iso_m];
1559       };
1560 
1561     fm_commutes := procedure(fm, d) {
1562       // Restrict fm to just i + s() variables
1563       fm_is := { [{ [var, val] : [var, val] in iso1 | var in i + s() },
1564                   { [var, val] : [var, val] in iso2 | var in d.i + d.s() }]
1565                     : [iso1, iso2] in fm };
1566       // Restrict fm to just s() + o variables
1567       fm_so := { [{ [var, val] : [var, val] in iso1 | var in s() + o },
1568                   { [var, val] : [var, val] in iso2 | var in d.s() + d.o }]
1569                     : [iso1, iso2] in fm };
1570       return { fm_so[ t()[g_is] ] == t()[ fm_is[g_is] ] : [g_is, h_so] in t()
1571              } == { true };
1572     };
1573 
1574 // Rationality
1575 // ===========
1576 
1577   /* Rational Utility Function
1578      -------------------------
1579      Takes the world w and brain b and collects all possible continuations of
1580      brain/decision states using bs_msr_paths. Then it weights them by 
1581      agential_identity, optimal_rxs, and subjective likelihood and finds the 
1582      center-of-mass / social-choice of the resulting first order utility 
1583      functions. 
1584 
1585      This comes out of the metaethical view developed by Howard Nye and June
1586      Ku and written up at http://www.metaethical.ai/norm_descriptivism.pdf. 
1587      There we focused on the simple case in which the norms we accept or our
1588      higher-order preferences formed a hierarchy with a single fundamental norm 
1589      or highest-order utility function governing our lower-order norms /
1590      functions. (Also see the comments above n := n towards the top.) In that 
1591      case, what we have most reason to do would be what that highest order 
1592      utility function would prescribe, perhaps through the prescription of 
1593      certain lower-level utilty functions.
1594 
1595      Elsewhere, we have acknowledged that things might be more complicated. We 
1596      might have a network of accepted norms / higher-order preferences in which 
1597      they all influence one another rather than a strict top-down hierarchy.
1598      We have suggested that in that case, the weights of each accepted norm and
1599      the parameters governing the network behavior could be treated as the 
1600      fundamental norm. Developing this further, we might model it as a markov 
1601      model with states of norm acceptances influencing other acceptances and 
1602      look for a stationary distribution, i.e. an equilibrium state in which the 
1603      state it transitions to is the same state it started with. 
1604 
1605      However, even that makes a simplifying assumption that the network 
1606      topology of which norms influence which others remains constant. But if 
1607      we imagine a brain's behavior in response to varying inputs, it seems
1608      likely that sometimes, which higher-order preferences influence which 
1609      others would also vary. To handle this, we implement a parliamentary model
1610      to aggregate the various influences across different possible 
1611      continuations of the agent. 
1612        This is inspired by Nick Bostrom and Toby Ord on Parliamentary Models 
1613      for Moral Uncertainty. Although it's applied in a different context, it's
1614      similar in that both are aggregating many utility functions into one.
1615        http://www.overcomingbias.com/2009/01/moral-uncertainty-towards-a-solution.html
1616 
1617      In the end, this does bring us closer to the ideal response / advisor
1618      family of metaethical theories but I think our focus on the influence of
1619      optimal normative judgments remains distinctive and arguably avoids some
1620      worries. We are reading off the first-order utility function from ideal
1621      selves rather than the actions so we have less distortion from the 
1622      idealization process to worry about since there's less to change. We
1623      address remaining worries about distortion with agential identity to 
1624      measure how much one remains the same agent. 
1625 
1626      The idealization is mostly done by weighting more heavily selves who are
1627      the result of normative judgments that best satisfy the decision criteria
1628      governing them. We also incorporate subjective likelihood of the inputs
1629      determining possible continuations to factor in normal operating 
1630      conditions without letting in external factors and violating supervenience
1631      on just the brain.
1632   */
1633   ruf := procedure(w, b) {
1634     t_z := 10 ** 9; // Number of time steps of continuations to consider
1635     /* Todo: Find a more principled way
1636          * "Ask" each branch how much longer to continue?
1637          * Compute infinitely but in a way that asymptotes? 
1638          * Keep computing until it converges, with minimum time? 
1639              If it never converges? */
1640     voter_wts := {};
1641     t_i := 1;
1642     bs1 := drop_time(b.a[b.n]);
1643     while (t_i <= t_z) {
1644       bmps := bs_msr_paths(w, b, bs1, t_i);
1645       voter_wts := voter_wts + {
1646         [ additive_to_ordinal(f(b)[first(last(bmp))], 0), // voter util func
1647           time_wt(t_i) * voter_weight(w, b, bmp) ]        // voter weight
1648         : bmp in last(bmps)
1649       };
1650       t_i := t_i + 1;
1651     }
1652     // Poss social choice utility functions
1653     psc_u_dom := domain(p) + domain(u);
1654     psc_u_rng := r[first(range(u))] + { om };
1655     psc_us := choices({ { [expr, r_i] : r_i in psc_u_rng }
1656                         : expr in psc_u_dom });
1657     psc_us := { { [expr, r_i] : [expr, r_i] in psc_u | r_i != om }
1658                 : psc_u in psc_us }; // filter out om
1659     return sort_set(psc_us,
1660                     procedure(psc_u1, psc_u2) {
1661                       return rat_psc_dist(voter_wts, psc_u1) <
1662                              rat_psc_dist(voter_wts, psc_u2);
1663                     })[1];
1664   };
1665 
1666     time_wt := procedure(t_i) { return 1; }; // Placeholder constant function
1667 
1668     voter_weight := procedure(w, b, bmp) {
1669       msr := last(last(bmp));
1670       bss := [ bs : [bs, msr] in bmp ]; // Brain StateS
1671 
1672       cached_ref := om;
1673       wb := new_wb(w, b, bss);
1674       w2 := first(wb);
1675       b2 := last(wb);
1676 
1677       dss := [ f(b2)[bs] : bs in bss ];
1678       agential_id := agential_identity(dss, w2, b2);
1679       return msr * agential_id * optimal_rxs(w2, b2, dss);
1680     };
1681 
1682     // Possible Social Choice for Rationality (as opposed to in optimal_rxs)
1683     rat_psc_dist := procedure(voter_wts, psc_u) {
1684       ord_u := card_u_to_ord_u( add_u_to_card_u(psc_u) );
1685       return +/ { voter_wts[voter] * ord_u_dist(voter, ord_u) ** 2
1686                   : voter in domain(voter_wts) };
1687     };
1688 
1689   new_wb := procedure(w, b, bss) {
1690     b2   := b;
1691     b2.a := b2.a + [assign_time(bss[t_j], b.n + t_j) : t_j in [1..#bss]];
1692     b2.n := #(b2.a);
1693 
1694     w2   := w;
1695     w2.a := w2.a + [assign_time(bss[t_j], b.n + t_j) : t_j in [1..#bss]];
1696     w2.n := #(w2.a);
1697 
1698     return [w, b];
1699   };
1700 
1701   /* Optimal Normative Judgments / Prescriptions
1702      -------------------------------------------
1703      Returns a score for how much the list of decision states has had outputs
1704      of higher-order utility functions (i.e. normative judgments / 
1705      prescriptions) that maximized higher-order utility. If the state space
1706      of prescriptions (n[j]['o']) overlap, we favor those states voted for by 
1707      more accepted norms through a synchronic social choice function.
1708   */
1709   optimal_rxs := procedure(w, b, dss) {
1710     opts := [ optimality(w, b, dss[j], dss[j+1], b.n - (#dss - j) )
1711               : j in [1..(#dss - 1)] ];
1712     return 0.5 * min(opts) + 0.5 * avg(opts);
1713   };
1714 
1715     optimality := procedure(w, b, ds, ds2, t_i) {
1716       ord_us := { additive_to_ordinal(ds, j) : j in [1..#n] };
1717 
1718       /* Social choice to balance multiple synchronic utility functions if
1719          they conflict.  */
1720       psc_base_exprs := +/ { ue_base_exprs(j) : j in [1..#n] };
1721       psc_state_space := choices({ {[subform, true], [subform, false]}
1722                                    : subform in psc_base_exprs });
1723       sync_psc_us := choices({ { [state, r_i] : r_i in r[first(range(u))] }
1724                                : state in psc_state_space });
1725       sc_u := sort_set(psc_us,
1726                        procedure(psc_u1, psc_u2) {
1727                          return sync_psc_dist(ord_us, psc_u1) <
1728                                 sync_psc_dist(ord_us, psc_u2);
1729                        })[1];
1730 
1731       n_o_vars := +/ { n_i['o'] : n_i in n }; // Prescription vars
1732       // Observed state of those vars in ds2
1733       obs_rxs := { [n_o_var, ds2[n_o_var]] : n_o_var in n_o_vars };
1734 
1735       // Collect all possible prescriptions along with their (higher) utilities.
1736       util_rxs := {
1737         [rxs, util_of_rxs(w, b, t_i, ds, ds2, psc_base_exprs, sc_u, rxs)]
1738         : rxs in poss_rx_states()
1739       };
1740       // Sort them by optimality.
1741       util_rxs := sort_set(util_rxs, procedure(rxs_u1, rxs_u2) {
1742                                        return last(rxs_u1) > last(rxs_u2);
1743                                      });
1744       // Return a normalized score for the obs_rxs by their percentile rank.
1745       return 1 - index_of(
1746         [ obs_rxs,
1747           util_of_rxs(w, b, t_i, ds, ds2, psc_base_exprs, sc_u, obs_rxs) ],
1748         util_rxs
1749       ) / #util_rxs;
1750     }; // end optimality
1751 
1752       sync_psc_dist := procedure(ord_us, sync_psc_u) {
1753         ord_psc_u := card_u_to_ord_u(sync_psc_u);
1754         return +/ { ord_u_dist(ord_u, ord_psc_u) ** 2 : ord_u in ord_us };
1755       };
1756 
1757       util_of_rxs := procedure(w, b, t_i, ds, ds2, psc_base_exprs, sc_u, rxs) {
1758         /* For the sake of simplicity, we assume that higher-order utilities 
1759            depend on only intrinsic properties of brains / decision algorithms.
1760            This means that being given the future brain states bss is 
1761            sufficient to tell us the utility of various prescriptions under 
1762            this assumed future. This assumption can be quite plausible for 
1763            higher-order utilities that enforce consistency, coherence or some
1764            form of reflective equilibrium. If we abandon this assumption, we
1765            may need to model a set of possible worlds and a probability 
1766            distribution or measure over them.
1767 
1768            We might typically expect higher-order utilities to be placed in 
1769            some change in lower-order utility and then a recognition that some
1770            prescription would more-or-less immediately cause that change would 
1771            cause the prescription and in turn, the intended change. But since
1772            we can neither assume that the higher-order utility is placed in 
1773            a specific, concrete lower-order utility change in the near future as
1774            opposed to more open-ended properties nor assume the immediacy of the
1775            prescription's causation of a lower-order utility change, we just 
1776            look very far into the future and sum up the utility of any events 
1777            in that 4-dimensional space-time block that would be included in the
1778            reference of higher-order utility expressions. Since the network
1779            topology of the higher-order utility functions may be path dependent,
1780            we also use the bs_msr_paths to aggregate over the different futures
1781            made possible by various input states weighted by their subjective
1782            likelihood.
1783         */
1784         new_ds2 := rxs + { [var, val] : [var, val] in ds2
1785                                       | !(var in domain(rxs)) };
1786         bs := inv(f(b))[new_ds2];
1787         bmps := bs_msr_paths(w, b, bs, 10 ** 9);
1788         future_wts := {                      // collect a ... 
1789           [ [first(bs_msr) : bs_msr in bmp], // list of brain states, bss,
1790             last(last(bmp)) ]                // & measure from the final state
1791           : bmp in last(bmps)                // for each bs_msr_path
1792         }; // to get { [bss, msr], ... }
1793 
1794         return +/ {
1795           msr * util_of_future(w, b, t_i, ds, psc_base_exprs, sc_u, bss)
1796           : [bss, msr] in future_wts
1797         };
1798       }; // end util_of_rxs
1799 
1800         util_of_future := procedure(w, b, t_i, ds, psc_base_exprs, sc_u, bss) {
1801           wb := new_wb(w, b, bss);
1802           w2 := first(wb);
1803           b2 := last(wb);
1804 
1805           future_state := { [expr, eval_expr(w2, b2, t_i, ds, bss, expr)]
1806                             : expr in psc_base_exprs };
1807           return sc_u[future_state];
1808         };
1809 
1810         // Returns a boolean for whether the expr is true in that future
1811         eval_expr := procedure(w, b, t_i, ds, bss, expr) {
1812           p_var := p[expr];
1813           p_event := { [p_var, ds[p_var]] };
1814           b_event := first( lf_inv_at(b,t_i)[p_event] );
1815           ref_b_event := ref(w,b)[b_event];
1816           disjuncts := { restricted_eq(w.actual_state(), disjunct)
1817                          : disjunct in ref_b_event };
1818           return #{ disj : disj in disjuncts | disj == true } > 0;
1819         };
1820 
1821   /* Extensional rational utility */
1822   ext_ruf := procedure(w, b) {
1823     uf := ruf(w,b);
1824     // Map expressions to their referents.
1825     expr_ref := { [expr, ref_expr(w, b, b.n, expr)] : [expr, r_i] in uf };
1826     // Map every possible world state to its utility.
1827     return { [w_state, w_util(uf, expr_ref, w_state)]
1828              : w_state in w.cm().poss_states() };
1829   };
1830 
1831     w_util := procedure(uf, expr_ref, w_state) {
1832       /* Truths, or utility relevant truths, I suppose. Collect the expressions 
1833          in the rational utility function that had this w_state among its set of
1834          truthmakers / referents.  */
1835       truths := { expr : expr in domain(expr_ref) | w_state in expr_ref[expr] };
1836       // Sum the utilities 
1837       return +/ { uf[expr] : expr in truths };
1838     };
1839 
1840   /* Agential Identity
1841      -----------------
1842      Agential Identity is much like Personal Identity a la Derek Parfit's 
1843      Reasons and Persons. But where personal identity roughly concerns when 
1844      someone is the same person for the purposes of care, agential identity
1845      concerns when someone is the same agent for purposes of trusting their
1846      normative judgments. 
1847 
1848      If we are thinking in terms of ideal advisor / observer theories, we might
1849      be considering hypothetical versions of ourselves with more information 
1850      or other qualities. But consider the full knowledge of what heroin feels 
1851      like. If we had that information, we might become addicts and change our 
1852      values. But this probably isn't a change that would make us want to defer
1853      to that self's values. We may protest that something in the process of
1854      adding this knowledge has changed our values, not illuminated them.
1855      
1856      The concept of agential identity is meant to measure how much other 
1857      possible selves truly reflect *our* values and agency. Generally, this 
1858      will be measured either by similarity to our current selves or by
1859      identifying the changes as cognitive changes made through inference or
1860      applying higher-order decision criteria to values. 
1861   */
1862   agential_identity := procedure(dss, w, b) { // dss : list of decision states
1863     return 0.5 * agential_continuity(dss, w, b) +
1864            // kinda like: 0.5 * agential_connection(dss[1], {}, dss[#dss])
1865            0.5 * (0.5 * p_conn(dss[1], last(dss)) +
1866                   0.5 * u_conn(dss[1], {}, last(#dss)));
1867   };
1868   agential_continuity := procedure(dss, w, b) {
1869     conns := [agential_connection(dss[j], dss[j+1], dss[j+2], w, b, j)
1870               : j in [1..(#dss - 2)] ];
1871     return 0.5 * min(conns) + 0.5 * avg(conns);
1872   };
1873   agential_connection := procedure(ds1, ds2, ds3, w, b, t_i) {
1874     return 0.5 * avg([p_conn(ds1, ds2), p_conn(ds2, ds3)]) +
1875            0.5 * u_conn(ds1, ds2, ds3, w, b, t_i);
1876   };
1877   /* Agential continuity restricted to p_conn.
1878      Used for measuring diachronic rationality in incoherence. */
1879   p_conn_continuity := procedure(dss) {
1880     return avg({ p_conn(dss[j], dss[j+1]) : j in [1..(#dss - 1)] });
1881   };
1882   /* Measures connectedness of credences between two decision states. They
1883      are connected by either having similar probability values or by reflecting
1884      cognitive changes, i.e. updating using conditional probabilities and
1885      resulting in probability values that are close to (locally) ideally 
1886      rational inferences. */
1887   p_conn := procedure(ds1, ds2) {
1888     return avg({
1889       // Gain full credit of the higher of similarity or cognitive score and
1890       max( sim_cog_scores(ds1, ds2, p_i) ) +
1891         // Close the remaining distance proportional to
1892         (1 - max( sim_cog_scores(ds1, ds2, p_i) )) *
1893         // the other of the similarity or cognitive score
1894         min( sim_cog_scores(ds1, ds2, p_i) )
1895       // For every p var
1896       : p_i in p_vars()
1897     });
1898   };
1899 
1900     sim_cog_scores := procedure(ds1, ds2, p_i) {
1901       return [p_i_sim(ds1, ds2, p_i), cog_change(ds1, ds2, p_i)];
1902     };
1903 
1904     // Similarity / closeness. Returns 0-1.
1905     p_i_sim := procedure(ds1, ds2, p_i) {
1906       max_dist := max({ abs(ds1[p_i] - r_i) : r_i in r[p_i] });
1907       return 1 - abs(ds1[p_i] - ds2[id_var(p_i)]) / max_dist;
1908     };
1909 
1910     // Cognitive change. Returns 0-1.
1911     cog_change := procedure(ds1, ds2, p_i) {
1912       if (p_i in range(p)) { // rather than range(cp)
1913         return 0;
1914       } else {
1915         arr := p_inv(p_i);
1916         hyp := arr[1];
1917         evd := arr[2];
1918         ancs := non_shared_ancestors(p_i, m);
1919         if (p[hyp] in ancs && p[evd] in ancs && cp[[evd, hyp]] in ancs) {
1920           /* For simplicity we just use the p vals in the immediate parent 
1921              (ds1). Really we should construct and use 
1922              non_shared_ancestors_with_depth and pull the p values from the 
1923              correct time slice. Similarly for id_expr. */
1924 
1925           // Rational P(H|E) = P(E|H) * P(H) / P(E)
1926           rat_h_e := ds1[cp[[evd, hyp]]] * ds1[p[hyp]] / ds1[p[evd]];
1927           // Distance from rational
1928           dist_rat := abs(rat_h_e - ds2[ cp[ [id_expr(hyp), id_expr(evd)] ] ]);
1929           // Maximum distance possible
1930           max_dist := max({ abs(rat_h_e - r_i) : r_i in r[p_i] });
1931           return 1 - dist_rat / max_dist;
1932         } else {
1933           return 0;
1934         }
1935       } // end else 
1936     }; // end cog_change
1937 
1938   /* Takes an expression e and returns the expression at the next time step
1939      that would be identical. In many cases, this is just the same expression e,
1940      but we define input variables as indexicals so the same expression next
1941      time step does not refer to the expression in this time step. So, we
1942      introduce a new term At(expr, -1) which at that time step, would refer to 
1943      the same as this time's expr.  
1944   */
1945   id_expr := procedure(expr) {
1946     if (isList(expr)) { // in case it's called on p_inv of cp var
1947       return [ id_expr(expr_i) : expr_i in expr ];
1948     } else {
1949       if (p[expr] in i) { // current input
1950         return At_t(expr, -1);
1951       } else if (isTerm(expr) && fct(expr) == "At_t") { // past, add 1 more step
1952         return At_t(args(expr)[1], args(expr)[2] - 1);
1953       } else {
1954         return expr;
1955       }
1956     }
1957   };
1958 
1959   /* Takes a prob variable and finds its expr preimage, then finds the expr that
1960      in the next time would be identical to this expr. Returns that expr var. */
1961   id_var := procedure(p_var) {
1962     return p[id_expr(p_inv(p_var))];
1963   };
1964 
1965   poss_i_events := procedure() {
1966     return +/ { possible_states(isub, r) : isub in ne_pow(i) };
1967   };
1968 
1969   e_vars := procedure() {
1970     return +/ { range(n_i['e']) : n_i in n };
1971   };
1972 
1973   poss_ue_events := procedure() {
1974     return +/ { possible_states(ue_sub, r)
1975                 : ue_sub in ne_pow(range(u) + e_vars()) };
1976   };
1977 
1978   poss_p_i_events := procedure() {
1979     return +/ { { [p_i, r_i] : r_i in r[p_i] } : p_i in range(p) };
1980   };
1981 
1982   poss_p_var_events := procedure() {
1983     return +/ { { [p_i, r_i] : r_i in r[p_i] } : p_i in p_vars() };
1984   };
1985 
1986   // Possible states of all prescription vars
1987   poss_rx_states := procedure() {
1988     n_o_vars := +/ { n_i['o'] : n_i in n };
1989     return possible_states(n_o_vars, r);
1990   };
1991 
1992   // Possible single prescriptions for some higher-order utility
1993   poss_rxs := procedure() {
1994     n_os := { n_i['o'] : n_i in n };
1995     return +/ { possible_states(n_o, r) : n_o in n_os };
1996   };
1997 
1998   /* Utility Function Connectedness
1999      ------------------------------
2000      Measures how much the utility function has either remained similar or 
2001      changed as the result of a cognitive process involving i) recognition that 
2002      some prescription (rx) would cause a utility change, ii) that recognition
2003      in fact causes the prescription and iii) that prescription causes the 
2004      utility change. These are to be contrasted with changes to utilities from
2005      mere noise, associations or biases.  */
2006   u_conn := procedure(w, b, t_i, ds1, ds2, ds3) {
2007     actual_rxs := { rx : rx in poss_rxs()
2008                        | rx_to_tuple(w, b, t_i, ds1, ds2, ds3, rx) != {} };
2009     rx_subsets := pow(actual_rxs);
2010 
2011     dists := { dist_to(w, b, t_i, ds1, ds2, ds3, rx_subset)
2012                : rx_subset in rx_subsets };
2013     return 1 - min(dists);
2014   }; // end u_conn
2015 
2016     /* Takes a prescription and returns the tuples if any which realize i, ii
2017        and iii above.  */
2018     rx_to_tuple := procedure(w, b, t_i, ds1, ds2, ds3, rx) {
2019       poss_tuples := poss_p_var_events() >< poss_i_events()  ><
2020                                p_range() >< poss_ue_events();
2021       return { tuple : tuple in poss_tuples
2022                      | observed(w, b, t_i, rx, ds1, ds2, ds3, tuple) };
2023     };
2024 
2025       observed := procedure(w, b, t_i, rx, ds1, ds2, ds3, tuple) {
2026         f_p_i  := tuple[1]; // Assignment of value to p_i, rx []-> ue
2027         g_i_rx := tuple[2]; /* Assignment of value to a prescription which is
2028                                introspectively accessible as i vars */
2029         p_ue_val := tuple[3]; // Credence placed in ue (the consequent of []->)
2030         h_ue   := tuple[4]; // ue event
2031         return f_p_i != {} && g_i_rx != {} && h_ue != {} &&
2032           // g_i_rx and rx are implemented by the same brain event and 
2033           // the decision events occur at the right times
2034           lf_inv(g_i_rx) == lf_inv(rx) && restricted_eq(ds2, g_i_rx) &&
2035           restricted_eq(ds2, rx) && restricted_eq(ds1, f_p_i) &&
2036           restricted_eq(ds3, h_ue) &&
2037           // first(domain(f_p_i)) is p_i and p_i is _ []-> _
2038           fct(p_var_inv(first(domain(f_p_i)))) == "BoxArrow" &&
2039           /* The ref of the brain event of ds1 implementing the placement of 
2040              prob in the p_i expr == lf_inv h_ue in two time steps
2041              i.e. the consequent of []-> refers to the intended ue change */
2042           /* Todo: Allow for changes that had i, ii, and iii but only partially 
2043              updated the utility / eval functions.  */
2044           ref(w,b)[lf_inv_at(b, t_i)[
2045             {[ p[args(p_var_inv(first(domain(f_p_i))))[2]], p_ue_val ]}
2046           ]] == lf_inv_at(b, t_i + 2)[h_ue] &&
2047           // p_i influences each var in rx
2048           first(domain(f_p_i)) in { parents(v) : v in domain(rx) } &&
2049           // each var in rx influences each var in the utility change 
2050           forall(v in g_i_rx | v in { parents(ue) : ue in domain(h_ue) });
2051       };
2052 
2053     /* Given a subset of observed prescriptions, we compare the resulting 
2054        utility function to what was prescribed within that subset. The overall 
2055        u_conn score will take the minimum distance seen this way.
2056          Todo: Allow for changes that had i, ii, and iii but only partially 
2057        updated the utility / eval functions.  */
2058     dist_to := procedure(w, b, t_i, ds1, ds2, ds3, rx_subset) {
2059       u1 := { [u_i, ds1[u_i]] : u_i in u };
2060       e1 := { u1 } + { { [e_i, ds1[e_i]] : e_i in n_i['e'] } : n_i in n };
2061 
2062       rxd := { [e_i, rx_e_i_val(w, b, t_i, ds1, ds2, ds3, rx_subset, e_i)]
2063                : e_i in domain(e1) };
2064       rxd_u := { [u_i, rxd[u_i]] : u_i in u };
2065       rxd_n := [ { [e_i, rxd[e_i]] : e_i in n_i['e'] } : n_i in n ];
2066       rxd_e := [ rxd_u ] + rxd_n;
2067       // The prescribed util / eval change
2068       ord_rxd_e := [ additive_to_ordinal(rxd_e[j], j) : j in [1..#rxd_e] ];
2069       // The actual change in ds3
2070       ord_e := [ additive_to_ordinal(e_func(ds3, j), j) : j in [1..#rxd_e] ];
2071 
2072       e_dists := [ ord_u_dist(ord_e[j], ord_rxd_e[j]) : j in [1..#ord_e] ];
2073 
2074       // Todo: Add weights? 
2075       return avg(e_dists);
2076     }; // end dist_to
2077 
2078       /* Given a subset of observed prescriptions and an e_i var, return the
2079          value prescribed for the e_i var.  */
2080       rx_e_i_val := procedure(w, b, t_i, ds1, ds2, ds3, rx_subset, e_i) {
2081         /* The prescribed values (rxds) are the values of the e_i var in the 
2082            ue event from rx_to_tuple (for any rxs that prescribe for e_i). */
2083         rxds := { rx_to_tuple(w, b, t_i, ds1, ds2, ds3, rx_i)[4][e_i]
2084           : rx_i in rx_subset  // h_ue[e_i]
2085           | e_i in domain(rx_to_tuple(w, b, t_i, ds1, ds2, ds3, rx_i)[4]) };
2086             // e_i in domain(h_ue)
2087         if (#rxds > 0) {
2088           // If #rxds > 1, they should all be the same anyway bc it was observed
2089           return first(rxds);
2090           /* Todo: In later versions, we probably want to allow for non-exact 
2091              changes and then choose the one closest to the observed.  */
2092         } else {
2093           return e1[e_i]; // No change
2094         }
2095       };
2096 
2097   // Utility / Evaluation Function implemented in a given decision state ds
2098   e_func := procedure(ds, j) {
2099     if (j == 0) {
2100       return { [u_i, ds[u_i]] :  u_i in range(u) };
2101     } else {
2102       return { [e_i, ds[e_i]] : e_i in range(n[j]['e']) };
2103     }
2104   };
2105 
2106   ue := procedure(j) {
2107     if (j == 0) {
2108       return u;
2109     } else {
2110       return n[j]['e'];
2111     }
2112   };
2113 
2114   ue_base_exprs := procedure(j) {
2115     return +/ { subformulas(x_i) : x_i in domain(ue(j)) };
2116   };
2117 
2118   /* Takes a decision state, i.e. a function from (among other things) utility
2119      variables to numbers (utilities), and an index j and returns a set of 
2120      SuccEq (succeeds or equals) setlx functors. If j is 0, then u is converted
2121      to ordinal, otherwise, n[j]['e'] is. */
2122   additive_to_ordinal := procedure(ds, j) {
2123     return add_u_to_card_u(e_func(ds,j));
2124   };
2125 
2126   /* Brain-State Measure Paths 
2127      -------------------------
2128      Returns a list where for each future time, there is a set of 
2129      bs_msr_paths, each of which in turn is a list of bs_msrs for each time 
2130      from the first to the given time. A bs_msr is a pair [bs, msr]. 
2131     
2132      One can think in terms of a tree. At the root is the given starting brain 
2133      state. Each branch is a continuation with different inputs and a measure 
2134      for that branch based on subjective likelihood. A path is a full branch 
2135      line from root to leaf collecting pairs of brain states and measures along
2136      the way. 
2137 
2138      e.g. [bs_msr_paths_t1, ..., bs_msr_paths_tz]
2139        bs_msr_paths_tj := { bs_msrs_tj_1, bs_msr_tj_2, ... }
2140          bs_msrs_tj_k := [ bs_msr_t1, ..., bs_msr_tj ]
2141            bs_msr_tj := [bs_tj, msr]
2142 
2143                      [ bs1, 1 ]
2144                     /   |      \
2145                is1 /    |is2    \  is3   (input states leading to continuations)
2146                   /     |        \  
2147                  /      |         \
2148      [bs_t2a, 1/6]  [bs_t2b, 1/2]  [bs_t2c, 1/3]  
2149       /    |    \     /   |   \      /   |   \
2150      ...
2151 
2152      There's some duplication and space inefficiency in having the previous 
2153      times' brain states and measures present both in the earlier elements of
2154      the top level list as well as within the set of paths for each element.
2155      We probably could do without it but it makes it a little more convenient
2156      when calculating the voter weights for each time's voters in computing
2157      the rational utility function. Some aspects of those weights like the 
2158      agential identity depend on the entire path taken and not just the final
2159      brain state.
2160   */
2161   bs_msr_paths := procedure(w, b, bs1, t_j) {
2162     if (t_j == 1) {
2163       return [{[[bs1, 1]]}];
2164     } else {
2165       return bs_msr_paths(w, b, bs1, t_j - 1) +
2166         [+/ {
2167            { bmp + [bs_msr]
2168              : bs_msr in next_bs_msrs(w, b, last(bmp)[1], last(bmp)[2], t_j-1) }
2169            : bmp in last(bs_msr_paths(w, b, bs1, t_j - 1))
2170          }];
2171     }
2172   };
2173 
2174   /* Takes a brain, brain state, measure and prev t_i and returns the possible
2175      continuations in the next time step together with their measure 
2176      roughly proportional to the subjective likelihood.  */
2177   next_bs_msrs := procedure(w, b, bs_i, msr, t_j) {
2178     bs_j := b.response(bs_i);
2179     // Possible brain input states from each di possible decision input states.
2180     pbis := { [di, lf_inv(b)[di]] : di in poss_i_states() };
2181     /* Possibly empty brain input states since we now filter by compatibility
2182        with the brain response bs_j. (Recall that input vars may overlap with
2183        other vars that bs_j may dictate values for.) */
2184     pe_bis := { [di, { pbi : pbi in pbi_s | compat(pbi, bs_j) }]
2185                 : [di, pbi_s] in pbis };
2186     // pbi_s = possible brain input sets
2187     bis := { [di, pbi_s] : [di, pbi_s] in pe_bis | pbi_s != {} };
2188 
2189     /* For each brain input state of each decision input state, map the total
2190        brain state including the response bs_j to its new measure. 
2191 
2192        Ideally, the new measure would be proportional to the brain input 
2193        state's subjective probability. This was under development but doing it 
2194        in a way that is biologically plausible presents a challenge. To avoid 
2195        an implausible number of explicit probabilities, we'd probably have to 
2196        look at the closest coherent extensions of the current probability 
2197        distribution as well as allow for indirect reference, which introduces 
2198        its own challenges in translating probabilities from direct to indirect.
2199        
2200        Although I by no means think these challenges insurmountable, the 
2201        technical complications involved seemed too great for a first version
2202        and not worth delaying the initial release over. As a placeholder here,
2203        we simply assume each possible decision input state is equally probable
2204        and each brain state implementing a given decision input state is
2205        similarly equally probable to other such implementations.
2206 
2207        While incorporating the subjective probability seems preferable, it's
2208        also worth mentioning that it may not actually be necessary. The main
2209        motivation is to make continuations of the agent under very weird 
2210        scenarios matter less. But if the worry there is that those 
2211        circumstances are distorting the agent's values, we already have the 
2212        agential_identity measure to discount them. In other words, agential
2213        identity may screen off a circumstance's weirdness from its relevance.
2214     */
2215     return +/ {
2216       { [ bi + bs_j, msr * (1 / #domain(bis)) * (1 / #bis[di]) ]
2217         : bi in bis[di] }
2218       : di in domain(bis)
2219     };
2220   };
2221 
2222 
2223 // Class Methods
2224 // =============
2225 
2226   static {
2227     /* Function from a causal_model of a brain b and a set of decision 
2228        algorithm candidates to the decision algorithm that is the best 
2229        compression of the brain.  */
2230     implemented_by := cachedProcedure(b, ds) {
2231       if (ds == {}) {
2232         ds := decision_algorithm.all(b);
2233       }
2234       /* Requiring simplest version helps ensure higher-order ambitiousness,
2235          i.e. prevents smuggling in higher-order norms within the n[j]['m2'] 
2236          variables. In future iterations, this might be softened.  */
2237       ds := [d : d in ds | d.is_implemented_by(b) && d.is_simplest_version(b)];
2238       ds := sort_list(ds, decision_algorithm.better_explanation(b));
2239       return ds[1];
2240     };
2241 
2242     /* Given a causal model of a brain, return the set of possible decision 
2243        algorithms which could be implemented by it based simply on number of 
2244        states. A decision algorithm can't be more fine-grained than its 
2245        substrate.  */
2246     all := cachedProcedure(b) {
2247       return { d : d in decision_algorithm.all_lte(#b.poss_states()) };
2248     };
2249 
2250     // All with #poss_states <= z
2251     // Todo: Add possible At_t exprs and predicates and quantifiers
2252     all_lte := cachedProcedure(z) {
2253       poss_i := da().ltr_up_to('I', z/2);
2254       var_ltrs := { 'A', 'B', 'C' };
2255       poss_base_vars := da().vars_up_to(var_ltrs, n);
2256       poss_exprs := { da().add_up_to(pbe, z/2)
2257                       : pbe in poss_base_vars };
2258       poss_p := da().poss_expr_to_ltr(poss_exprs, 'P');
2259       // e.g. { [[x, y], 'P(x|y)'], ... } 
2260       poss_cp := pow({ [[x, y], 'P(' + str(x) + '|' + str(y) + ')']
2261                        : [x, y] in poss_exprs >< poss_exprs });
2262       poss_o := da().ltr_up_to('O', z/2);
2263       poss_u := da().poss_expr_to_ltr(poss_exprs, 'U');
2264       poss_m := da().ltr_up_to('M', z/2);
2265       poss_m2 := da().ltr_up_to('M2', z/2);
2266       poss_n := +/ { choices({ da().poss_nj(poss_exprs, z, j) : j in [1..x] })
2267                      : x in [1..z] };
2268       return +/ {
2269         +/ { da().poss_ds(d_i, d_o, d_p, d_cp, d_u, d_m, d_m2, d_n, poss_r)
2270              : poss_r in da().poss_rs(z, d_i, d_o, d_p, d_cp,
2271                                          d_u, d_m, d_m2, d_n) }
2272         : [d_i, d_o, d_p, d_cp, d_u, d_m, d_m2, d_n] in
2273           da().cart_prod(poss_i, poss_o, poss_p, poss_cp,
2274                          poss_u, poss_m, poss_m2, poss_n)
2275       };
2276     }; // end all_lte
2277 
2278       ltr_up_to := procedure(ltr, j) {
2279         return { ltr + '_' + y : y in [1..j] };
2280       };
2281 
2282       vars_up_to := procedure(ltrs, j) {
2283         return +/ { da().ltr_up_to(ltr, j) : ltr in ltrs };
2284       };
2285 
2286       add_one := procedure(exprs) {
2287         unary := { 'Not' };
2288         binary := { 'And', 'Or', 'Implies', 'BoxArrow' };
2289         add_un := { makeTerm(op, [expr]) : [op, expr] in unary >< exprs };
2290         add_bi := { makeTerm(op, [expr1, expr2])
2291                     : [op, [expr1, expr2]] in binary >< (exprs >< exprs) };
2292         return exprs + add_un + add_bi;
2293       };
2294 
2295       add_up_to := procedure(base_vars, j) {
2296         if (j == 0) {
2297           return base_vars;
2298         } else {
2299           return da().add_one( da().add_up_to(base_vars, j - 1) );
2300         }
2301       };
2302 
2303       expr_to_ltr_var := procedure(expr, ltr) {
2304         return ltr + '(' + replace(str(expr), '"', '') + ')';
2305       };
2306 
2307       poss_expr_to_ltr := procedure(poss_exprs, ltr) {
2308         return { { [expr, da().expr_to_ltr_var(expr, ltr)] : expr in subset }
2309                  : subset in pow(poss_exprs) };
2310       };
2311 
2312       poss_nj := procedure(poss_exprs, z, j) {
2313         poss_e_nj := da().poss_expr_to_ltr(poss_exprs, 'E_n' + j);
2314         poss_m2_nj := da().ltr_up_to('M2_n' + j, z/2);
2315         poss_o_nj := da().ltr_up_to('O_n' + j, z/2);
2316         return { { ['e', e_n], ['m', m_n], ['o', o_n] }
2317           : [e_n, m_n, o_n] in cart_prod([poss_e_nj, poss_m2_nj, poss_o_nj]) };
2318       };
2319 
2320       poss_rs := procedure(z, d_i, d_o, d_p, d_cp, d_u, d_m, d_m2, d_n) {
2321         poss_p_vals := { { j / x : j in [0..x] } : x in [1..z] };
2322         poss_ue_vals := ne_pow({ j : j in [1..2**z] });
2323         poss_z_vals := { { j : j in [1..x] } : x in [1..z] };
2324         poss_o_vals := poss_z_vals;
2325 
2326         poss_r_ps := {
2327           { [i_var,  p_vals] : i_var  in d_i         } +
2328           { [p_var,  p_vals] : p_var  in range(d_p)  } +
2329           { [cp_var, p_vals] : cp_var in range(d_cp) }
2330           : p_vals in poss_p_vals
2331         };
2332 
2333         return +/ choices(
2334           { poss_r_ps, da().poss_r_ue(poss_ue_vals, d_u) } +
2335           { da().poss_r_ue(poss_ue_vals, n_i['e']) : n_i in d_n } +
2336           { da().poss_r_z(poss_z_vals, d_z) : d_z in { d_o, d_m, d_m2 } }
2337         );
2338       }; // end poss_rs
2339 
2340         poss_r_ue := procedure(poss_ue_vals, ue) {
2341           return { { [ue_var, ue_vals] : ue_var in range(ue) }
2342                    : ue_vals in poss_ue_vals };
2343         };
2344 
2345         poss_r_z := procedure(poss_z_vals, d_z) {
2346           return choices({ {z_var} >< poss_z_vals : z_var in d_z });
2347         };
2348 
2349       poss_ds := procedure(z, d_i, d_o, d_p, d_cp, d_u, d_m, d_m2, d_n, d_r) {
2350         poss_is_state := possible_states(d_i + d_p + d_cp +
2351                                          d_u + d_m + d_m2 + d_n, d_r);
2352         poss_s_state  := possible_states(d_p + d_cp + d_u +
2353                                          d_m + d_m2 + d_n, d_r);
2354         poss_t_s := function_space(poss_is_state, poss_s_state);
2355 
2356         poss_p_u_m2 := possible_states(d_p + d_u + d_m2, d_r);
2357         poss_o_state := possible_states(d_o, d_r);
2358         poss_t_o := function_space(poss_p_u_m2, poss_o_state);
2359 
2360         das := {
2361           decision_algorithm(d_i, d_o, d_p, d_cp, d_u, d_m, d_m2, d_n,
2362                              d_r, d_t_s, d_t_o, om, om, om)
2363           : [d_t_s, d_t_o] in poss_t_s >< poss_t_o
2364         };
2365         return { da : da in das | da.is_valid() && #da.poss_states() <= z};
2366       };
2367 
2368     better_explanation := cachedProcedure(b) {
2369       return procedure(d1, d2) {
2370                return da().expl_score(b, d1) > da().expl_score(b, d2);
2371              };
2372     };
2373 
2374       expl_score := procedure(b, d_i) {
2375         return d_i.ambitiousness(b) - d_i.complexity(b)
2376                - d_i.instr_irrat(b) - d_i.incoherence(b);
2377       };
2378 
2379     simpler := cachedProcedure() {
2380       return procedure(d1, d2) {
2381                return k(d1.to_s) < k(d2.to_s);
2382              };
2383     };
2384 
2385     new := procedure() {
2386       return decision_algorithm({}, {}, {}, {}, {}, {}, {}, [], {}, {}, {},
2387                                 om, om, om);
2388     };
2389 
2390   }
2391 }
2392 
2393 // alias 
2394 da := procedure() { return decision_algorithm; };
   1 print("\nTesting Decision Algorithm");
   2 print(  "==========================");
   3 
   4 test_s := procedure() {
   5   d := decision_algorithm({'I'}, {'O'}, { ['A', 'P(A)'] }, // i, o, p
   6     { [['A', 'A'], 'P(A|A)'] }, {['A','U(A)']}, {'M'}, {'M2'}, // cp, u, m, m2
   7     [{ ['e', {['B','E_n1(B)']}], ['m2', {'M2_n1'}], ['o', {'O_n1'}] }], // n
   8     {}, {}, {}, om, om, om); // r, t_s, t_o, cacheds
   9   assert(d.s() == { 'P(A)', 'P(A|A)', 'U(A)', 'M', 'M2',
  10                     'E_n1(B)', 'M2_n1', 'O_n1' },
  11          "Error: s is incorrect.");
  12   print('[ok] s');
  13 };
  14 test_s();
  15 
  16 test_t := procedure() {
  17   d := decision_algorithm.new();
  18   d.i := { 'I' };
  19   d.o := { 'O' };
  20   d.p := { ['A', 'P'] };
  21   d.t_s := { // is                         -> s
  22              [{ ['I', false], ['P', false] }, { ['P', false] }],
  23              [{ ['I', false], ['P', true ] }, { ['P', true ] }],
  24              [{ ['I', true ], ['P', false] }, { ['P', true ] }],
  25              [{ ['I', true ], ['P', true ] }, { ['P', false] }]
  26            }; // xor
  27   d.t_o := { // p u m2       -> o
  28              [{ ['P', false] }, { ['O', false] }],
  29              [{ ['P', true ] }, { ['O', true ] }]
  30            }; // identity
  31   assert(d.t() == {
  32            // is                         ->  so
  33            [{["I", false], ["P", false]}, {["O", false], ["P", false]}],
  34            [{["I", false], ["P", true ]}, {["O", true ], ["P", true ]}],
  35            [{["I", true ], ["P", false]}, {["O", false], ["P", true ]}],
  36            [{["I", true ], ["P", true ]}, {["O", true ], ["P", false]}]
  37          },
  38          "Error: t is incorrect.");
  39   print('[ok] t');
  40 };
  41 test_t();
  42 
  43 test_ns := procedure() {
  44   d := decision_algorithm.new();
  45   d.n := [{ ['e', {['A','E_n1(A)']}], ['m2', {'M2_n1'}], ['o', {'O_n1'}] },
  46           { ['e', {['B','E_n2(B)']}], ['m2', {'M2_n2'}], ['o', {'O_n2'}] } ];
  47   assert(d.ns( 'e') == { 'E_n1(A)', 'E_n2(B)' } &&
  48          d.ns('m2') == {   'M2_n1', 'M2_n2'   } &&
  49          d.ns( 'o') == {    'O_n1', 'O_n2'    },
  50          "Error: ns is incorrect.");
  51   print('[ok] ns');
  52 };
  53 test_ns();
  54 
  55 test_new := procedure() {
  56   assert(decision_algorithm.new().t_o == {},
  57          "Error: new is incorrect.");
  58   print('[ok] new');
  59 };
  60 test_new();
  61 
  62 test_p_u_vars := procedure() {
  63   d := decision_algorithm.new();
  64   d.p  := { [       'A',   'P(A)'] };
  65   d.cp := { [['A', 'A'], 'P(A|A)'] };
  66   d.u  := { [       'A',   'U(A)'] };
  67   assert(d.p_u_vars() == {'P(A)', 'P(A|A)', 'U(A)'},
  68          "Error: is incorrect.");
  69   print('[ok] p_u_vars');
  70 };
  71 test_p_u_vars();
  72 test_p_vars := procedure() { test_p_u_vars(); };
  73 
  74 test_just_p_u_m2 := procedure() {
  75   d := decision_algorithm.new();
  76   d.i  := { 'I' };
  77   d.p  := { ['A', 'P(A)'] };
  78   d.u  := { ['A', 'U(A)'] };
  79   d.m2 := { 'M2' };
  80   assert(d.just_p_u_m2({ ['I', true],
  81                          ['P(A)', 0.3], ['U(A)', 4], ['M2', true] }) ==
  82                        { ['P(A)', 0.3], ['U(A)', 4], ['M2', true] },
  83          "Error: just_p_u_m2 is incorrect.");
  84   print('[ok] just_p_u_m2');
  85 };
  86 test_just_p_u_m2();
  87 
  88 test_parents := procedure() {
  89   d := decision_algorithm.new();
  90   d.m := { 'M1', 'M2' };
  91   d.r := { ['M1', {true, false}], ['M2', {true,false}] };
  92   d.t_s := { // is                  -> s
  93     [{ ['M1', false], ['M2', false] }, { ['M1', true ], ['M2', true ] }],
  94     [{ ['M1', false], ['M2', true ] }, { ['M1', false], ['M2', true ] }],
  95     [{ ['M1', true ], ['M2', false] }, { ['M1', false], ['M2', false] }],
  96     [{ ['M1', true ], ['M2', true ] }, { ['M1', true ], ['M2', false] }]
  97   }; // M1 = ~xor; M2 = ~M1
  98   assert(d.parents('M1') == {'M1', 'M2'} && d.parents('M2') == {'M1'},
  99          "Error: parents is incorrect.");
 100   print('[ok] parents');
 101   assert(d.children('M1') == { 'M1', 'M2' } && d.children('M2') == {'M1'},
 102         "Error: children is incorrect.");
 103   print('[ok] children');
 104 };
 105 test_children               := procedure() { test_parents(); };
 106 test_is_necessary_parent_of := procedure() { test_parents(); };
 107 test_vals_diff              := procedure() { test_parents(); };
 108 test_parents();
 109 
 110 /*   ,------.
 111     I--P--Y--X
 112      `--U-'
 113          `--Z
 114 */
 115 test_non_shared_ancestors := procedure() {
 116   d := decision_algorithm.new();
 117   d.i := {'I'};
 118   d.p := {['A', 'P']};
 119   d.u := {['B', 'U']};
 120   d.m2 := {'X', 'Y', 'Z'};
 121   d.r := { ['I', tf()], ['P', {0,1}], ['U', {0,1}],
 122            ['X', tf()], ['Y', tf()], ['Z', tf()] };
 123   ts_proc := procedure(is_state) {
 124     return { ['P', !is_state['I']], ['U', is_state['I']],
 125              ['X', !is_state['Y'] || is_state['I']],
 126              ['Y', is_state['P'] + is_state['U'] % 2],
 127              ['Z', is_state['U'] == 0] };
 128   };
 129   d.t_s := proc2func(ts_proc, d.poss_is_states());
 130   assert(d.non_shared_ancestors('X', d.m2) == { 'I', 'P', 'U' } &&
 131          d.non_shared_ancestors('Y', d.m2) == { 'P', 'U' } &&
 132          d.non_shared_ancestors('Z', d.m2) == { 'U' },
 133          "Error: non_shared_ancestors is incorrect.");
 134   print('[ok] non_shared_ancestors');
 135 
 136   assert(d.poss_i_states() == { {['I', true]}, {['I', false]} },
 137          "Error: poss_i_states is incorrect.");
 138   print('[ok] poss_i_states');
 139 
 140   assert({ {["P", 0], ["U", 0], ["X",  true], ["Y",  true], ["Z",  true]},
 141            {["P", 1], ["U", 1], ["X", false], ["Y", false], ["Z", false]}
 142          } < d.poss_s_states(),
 143          "Error: poss_s_states is incorrect");
 144   print('[ok] poss_s_states');
 145 
 146   assert(d.poss_u_states() == { {['U', 0]}, {['U', 1]} },
 147          "Error: poss_u_states is incorrect.");
 148   print('[ok] poss_u_states');
 149 
 150   assert(d.is_equal_to(d),
 151          "Error: is_equal_to is incorrect");
 152   print('[ok] is_equal_to true');
 153 
 154   d2 := d;
 155   d2.i := d.i + { 'I2' };
 156   assert(!d2.is_equal_to(d),
 157          "Error: is_equal_to is incorrect");
 158   print('[ok] is_equal_to false');
 159 };
 160 test_non_shared_ancestors();
 161 test_poss_is_states := procedure() { test_non_shared_ancestors(); };
 162 test_poss_i_states  := procedure() { test_non_shared_ancestors(); };
 163 test_poss_s_states  := procedure() { test_non_shared_ancestors(); };
 164 test_poss_u_states  := procedure() { test_non_shared_ancestors(); };
 165 test_is_equal_to    := procedure() { test_non_shared_ancestors(); };
 166 
 167 
 168 
 169 /* X--Y--O
 170     `--M--P
 171 */
 172 test_non_shared_descendants := procedure() {
 173   d := decision_algorithm.new();
 174   d.m2 := {'X', 'Y'};
 175   d.o := {'O'};
 176   d.m := {'M'};
 177   d.p := {['A','P']};
 178   d.r := { ['X', tf()], ['Y', tf()], ['O', tf()], ['M', tf()], ['P', {0,1}] };
 179   ts_proc := procedure(is_state) {
 180     return { ['X', true], ['Y', !is_state['X']], ['O', is_state['Y']],
 181              ['M', is_state['X']], ['P', #{ j : j in {1} | is_state['M']}] };
 182   };
 183   d.t_s := proc2func(ts_proc, d.poss_is_states());
 184   to_proc := procedure(p_u_m2) {
 185     return {[ 'O', p_u_m2['Y'] ]};
 186   };
 187   d.t_o := proc2func(to_proc, { d.just_p_u_m2(is) : is in d.poss_is_states() });
 188   assert(d.non_shared_descendants('X', d.m2) == {'M', 'O'} &&
 189          d.non_shared_descendants('Y', d.m2) == {'O'},
 190          "Error: non_shared_descendants is incorrect.");
 191   print('[ok] non_shared_descendants');
 192 };
 193 test_non_shared_descendants();
 194 
 195 d1 := cachedProcedure() {
 196   d   := decision_algorithm.new();
 197   d.i := { 'I' };
 198   d.p := { ['A', 'P'] };
 199   d.o := { 'O' };
 200   d.r := { [var, tf()] : var in {'I','O'} } + { ['P', {0,1}] };
 201   d.t_s := { [{['I', false], ['P', 0]}, {['P', 1]}], // eg P_t2 = ~I xor P_t1
 202              [{['I', false], ['P', 1]}, {['P', 0]}],
 203              [{['I',  true], ['P', 0]}, {['P', 0]}],
 204              [{['I',  true], ['P', 1]}, {['P', 1]}]
 205            };
 206   d.t_o := { [ {['P', 0]}, {['O',  true]} ],  // O = ~P
 207              [ {['P', 1]}, {['O', false]} ] };
 208   return d;
 209 };
 210 
 211 b1 := cachedProcedure() {
 212   b   := causal_markov_model.new();
 213   b.u := { 'U' };
 214   b.v := { 'V_P1', 'V_P2', 'O' };
 215   b.r := { [var, tf()] : var in {'V_P1', 'V_P2', 'O'} } + {['U',{1,2,3}]};
 216   /* P = V_P1 <-> V_P2
 217      O = ~ P
 218      O = ~ (V_P1 <-> V_P2) 
 219   */
 220   b.f['O'] := proc2func(procedure(cmm_st) {
 221                           return !(cmm_st['V_P1'] == cmm_st['V_P2']);
 222                         }, b.poss_states());
 223   /* P_t2 = I xor P_t1
 224           = U xor (V_P1 <-> V_P2)     
 225     t1  U V_P1 V_P2 P  ~I  t2  P  V_P1  V_P2
 226         1   F   F   T   T      F   T    F
 227         1   F   T   F   T      T   T    T
 228         1   T   F   F   T      T   F    F
 229         1   T   T   T   T      F   F    T
 230         2   F   F   T   T      F    
 231         2   F   T   F   T      T
 232         2   T   F   F   T      T
 233         2   T   T   T   T      F
 234         3   F   F   T   F      T   F    F    
 235         3   F   T   F   F      F   T    F
 236         3   T   F   F   F      F   F    T
 237         3   T   T   T   F      T   T    T
 238 
 239     two ways for U to be false {1,2}. {3} is truth
 240   */
 241   b.f['V_P1'] := proc2func(procedure(cmm_st) {
 242                    if (cmm_st['U'] == 3) {
 243                      return cmm_st['V_P2'];
 244                    } else {
 245                      return !cmm_st['V_P1'];
 246                    }
 247                  }, b.poss_states());
 248   b.f['V_P2'] := proc2func(procedure(cmm_st) {
 249                    if (cmm_st['U'] == 3) {
 250                      return cmm_st['V_P1'];
 251                    } else {
 252                      return cmm_st['V_P2'];
 253                    }
 254                  }, b.poss_states());
 255   b.n := 3;
 256   return b;
 257 };
 258 
 259 plf1 := cachedProcedure() {
 260   return { [{ ['U', 1]                         }, { ['I', false] }],
 261            [{ ['U', 2]                         }, { ['I', false] }],
 262            [{ ['U', 3]                         }, { ['I',  true] }],
 263            [{ ['V_P1', false], ['V_P2', false] }, { ['P',     1] }],
 264            [{ ['V_P1', false], ['V_P2',  true] }, { ['P',     0] }],
 265            [{ ['V_P1',  true], ['V_P2', false] }, { ['P',     0] }],
 266            [{ ['V_P1',  true], ['V_P2',  true] }, { ['P',     1] }],
 267            [{ ['O', false]                     }, { ['O', false] }],
 268            [{ ['O',  true]                     }, { ['O',  true] }]
 269          };
 270 };
 271 
 272 pf1 := cachedProcedure() {
 273   return d1().plf_to_pf(b1(), plf1());
 274 };
 275 
 276 test_poss_states := procedure() {
 277   assert(d1().poss_states() == {
 278            {["I", false], ["O", false], ["P", 0]},
 279            {["I", false], ["O", false], ["P", 1]},
 280            {["I", false], ["O",  true], ["P", 0]},
 281            {["I", false], ["O",  true], ["P", 1]},
 282            {["I",  true], ["O", false], ["P", 0]},
 283            {["I",  true], ["O", false], ["P", 1]},
 284            {["I",  true], ["O",  true], ["P", 0]},
 285            {["I",  true], ["O",  true], ["P", 1]}
 286          }, "Error: poss_states is incorrect.");
 287   print('[ok] poss_states');
 288 };
 289 test_poss_states();
 290 
 291 test_commutes := procedure() {
 292   // print('pf ' + pf1());  
 293 
 294   // This may take a while so it's commented out but it worked last time.
 295   // assert(d1().commutes(pf1(), b1()) == true, "Error: commutes true is incorrect.");
 296   // print('[ok] commutes true');
 297   print('[..] commutes true');
 298 };
 299 test_commutes();
 300 
 301 test_lf_commutes := procedure() {
 302   d := d1();
 303   d.cached_f := pf1();
 304   assert(d.lf_commutes(plf1(), b1()) == true,
 305          "Error: lf_commutes true is incorrect");
 306   print('[ok] lf_commutes true');
 307 
 308   plf := plf1();
 309   plf[ {['O', false]} ] := { ['O', true] };
 310   assert(d.lf_commutes(plf, b1()) == false,
 311          "Error: lf_commutes false is incorrect");
 312   print('[ok] lf_commutes false');
 313 };
 314 test_lf_commutes();
 315 
 316 d2 := procedure() {
 317   d := d1();
 318   d.cached_f  := pf1();
 319   d.cached_lf := plf1();
 320   return d;
 321 };
 322 
 323 test_is_implemented_by := procedure() {
 324   assert(d2().is_implemented_by(b1()) == true,
 325          "Error: is_implemented_by is incorrect.");
 326   print('[ok] is_implemented_by');
 327 };
 328 test_is_implemented_by();
 329 
 330 test_lf_inv := procedure() {
 331   assert(d2().lf_inv(b1())[ {['P',1]} ] == {
 332            { ['V_P1', false], ['V_P2', false] },
 333            { ['V_P1',  true], ['V_P2',  true] }
 334          }, "Error: lf_inv is incorrect");
 335   print('[ok] lf_inv');
 336 };
 337 test_lf_inv();
 338 test_plf_to_pf := procedure() { test_lf_inv(); };
 339 
 340 test_lf_inv_at := procedure() {
 341   assert(d2().lf_inv_at(b1(), 2)[ {['P',1]} ] == {
 342            { ['V_P1_t2', false], ['V_P2_t2', false] },
 343            { ['V_P1_t2',  true], ['V_P2_t2',  true] }
 344          }, "Error: lf_inv_at is incorrect");
 345   print('[ok] lf_inv_at');
 346 };
 347 test_lf_inv_at();
 348 
 349 test_lf_cm := procedure() {
 350   assert(d2().lf_cm(b1(), {['V_P1_t2',false],['V_P2_t2',false]}) == {['P',1]},
 351          "Error: lf_cm is incorrect");
 352   print('[ok] lf_cm');
 353 };
 354 test_lf_cm();
 355 
 356 test_lf_cm_vs := procedure() {
 357   assert(d2().lf_cm_vs(b1(), {['V_P1_t2',false],['V_P2_t2',false] }) == {'P'},
 358          "Error: lf_cm_vs is incorrect");
 359   print('[ok] lf_cm_vs');
 360 };
 361 test_lf_cm_vs();
 362 
 363 test_lf_cm_p_inv := procedure() {
 364   assert(d2().lf_cm_p_inv(b, { ['V_P1_t2', false], ['V_P2_t2', false] }) == 'A',
 365          "Error: lf_cm_p_inv is incorrect");
 366   print('[ok] p_inv');
 367   print('[ok] lf_cm_p_inv');
 368 };
 369 test_lf_cm_p_inv();
 370 test_p_inv := procedure() { test_lf_cm_p_inv(); };
 371 
 372 test_commutes_f := procedure() {
 373   // print('b.f[O] ' + b.f['O']);
 374   b := b1();
 375   b.f['O'] := proc2func(procedure(cmm_st) {
 376                           return cmm_st['V_P1'] && cmm_st['V_P2'];
 377                         }, b.poss_states());
 378   // print('b.f[O] ' + b.f['O']);
 379   // bs := {["O", true], ["U", 3], ["V_P1", true], ["V_P2", true]};
 380   // print('pf[bs] ' + pf[bs]);
 381   // print('t()[pf[bs]] ' + d.t()[pf[bs]]);
 382   // print('b.resp ' + b.response(b.v, bs));
 383   // print('pf[b.resp] ' + pf[b.response(b.v, bs)]);
 384 
 385   // assert(d1().commutes(pf1(), b) == false, "Error: commutes false is incorrect.");
 386   // print('[ok] commutes false');
 387   print('[..] commutes false');
 388 };
 389 test_commutes_f();
 390 
 391 test_drop_i := procedure() {
 392   d := decision_algorithm.new();
 393   d.i := { 'P(I_1)', 'P(I_2)' };
 394   ds  := { ['P(A)', 1], ['P(I_1)', true], ['P(I_2)', false], ['M1', true] };
 395   assert(d.drop_i(ds) == { ['P(A)', 1], ['M1', true] },
 396          "Error: drop_i is incorrect.");
 397   print('[ok] drop_i');
 398 };
 399 test_drop_i();
 400 
 401 /*   
 402     I--P--Y--X--O
 403      `--U-'
 404          `--Z--W
 405 */
 406 test_m2s_are_valid := procedure() {
 407   d    := decision_algorithm.new();
 408   d.i  := { 'I' };
 409   d.p  := { ['A', 'P'] };
 410   d.u  := { ['B', 'U'] };
 411   d.m2 := { 'X', 'Y', 'Z' };
 412   d.o  := { 'W', 'O' };
 413   d.r  := { [var, tf()] : var in d.i + d.s() };
 414   ts_proc := procedure(is_state) {
 415     return { ['P', !is_state['I']], ['U', is_state['I']],
 416              ['X', !is_state['Y']],
 417              ['Y', is_state['P'] && is_state['U']],
 418              ['Z', is_state['U'] == 0] };
 419   };
 420   d.t_s := proc2func(ts_proc, d.poss_is_states());
 421   to_proc := procedure(p_u_m2) {
 422     return { ['O', p_u_m2['X']], ['W', p_u_m2['W']] };
 423   };
 424   d.t_o := proc2func(to_proc, { d.just_p_u_m2(is) : is in d.poss_is_states() });
 425   assert(d.m2s_are_valid() == true, "Error: is incorrect.");
 426   print('[ok] m2s_are_valid true');
 427 
 428   // I------X
 429   ts_proc2 := procedure(is_state) {
 430     return { ['P', !is_state['I']], ['U', is_state['I']],
 431              ['X', !is_state['Y'] || is_state['I']],
 432              ['Y', is_state['P'] && is_state['U']],
 433              ['Z', is_state['U'] == 0] };
 434   };
 435   d.t_s := proc2func(ts_proc2, d.poss_is_states());
 436   assert(d.m2s_are_valid() == false, "Error: is incorrect.");
 437   print('[ok] m2s_are_valid false');
 438 
 439   // W-----U
 440   ts_proc3 := procedure(is_state) {
 441     return { ['P', !is_state['I']], ['U', is_state['I'] && is_state['W']],
 442              ['X', !is_state['Y']],
 443              ['Y', is_state['P'] && is_state['U']],
 444              ['Z', is_state['U'] == 0] };
 445   };
 446   d.t_s := proc2func(ts_proc2, d.poss_is_states());
 447   assert(d.m2s_are_valid() == false, "Error: is incorrect.");
 448   print('[ok] m2s_are_valid false 2');
 449 };
 450 test_m2s_are_valid();
 451 
 452 /*
 453 test_no_e_in_u := procedure() {
 454   d := decision_algorithm.new();
 455   assert(d.no_e_in_u() == ,
 456          "Error: no_e_in_u is incorrect.");
 457   print('[ok] no_e_in_u false');
 458 
 459   assert(d.no_e_in_u() == ,
 460          "Error: no_e_in_u is incorrect.");
 461   print('[ok] no_e_in_u true');
 462 };
 463 test_no_e_in_u();
 464 */
 465 
 466 /* test_i_ref := procedure() {
 467   print(d1().i_ref(b2(), b2())); // out of memory
 468 //  assert(d1().i_ref(b1(), b1()) == 
 469 //         "Error: i_ref is incorrect.");
 470   print('[  ] i_ref');
 471 };
 472 test_i_ref();
 473 */
 474 
 475 test_str_exprs := procedure() {
 476   d := decision_algorithm.new();
 477   d.i := { 'P(I_1)' };
 478   d.p := { ['A', 'P(A)'], ['B', 'P(B)'], ['I_1', 'P(I_1)'],
 479            [And('A', 'B'), 'P(A^B)'],
 480            [BoxArrow('A', 'B'), 'P(A[]->B)'] };
 481   assert(d.str_exprs() == { ['A',0] , ['B',0] },
 482          "Error: is incorrect.");
 483   print('[ok] str_exprs');
 484 };
 485 test_str_exprs();
 486 
 487 test_preds := procedure() {
 488   d := decision_algorithm.new();
 489   d.i := { 'P(I_1)' };
 490   d.p := { ['A', 'P(A)'], ['B', 'P(B)'], ['I_1', 'P(I_1)'],
 491            [And('A', 'B'), 'P(A^B)'], [Q('A'), 'P(Q(A))'],
 492            [R('A','B'), 'P(R(A,B))'],
 493            [BoxArrow('A', 'B'), 'P(A[]->B)'] };
 494   assert(d.preds() == { ['Q', 1], ['R', 2] },
 495          "Error: is incorrect.");
 496   print('[ok] preds arity');
 497 };
 498 test_preds();
 499 
 500 w1 := procedure() {
 501   w := causal_markov_model.new();
 502   w.u := {};
 503   w.v := { 'X', 'Y' };
 504   w.n := 2;
 505   w.r := { ['X', tf()], ['Y', tf()] };
 506   w.f := { ['X', { [{['X', false], ['Y', false]}, false],
 507                    [{['X', false], ['Y',  true]},  true],
 508                    [{['X',  true], ['Y', false]},  true],
 509                    [{['X',  true], ['Y',  true]},  true]
 510                  }],
 511            ['Y', { [{['X', false], ['Y', false]},  true],
 512                    [{['X', false], ['Y',  true]},  true],
 513                    [{['X',  true], ['Y', false]},  true],
 514                    [{['X',  true], ['Y',  true]}, false]
 515                  }]
 516          };
 517 };
 518 
 519 test_poss_base_expr_refs := procedure() {
 520   d := decision_algorithm.new();
 521   w := w1();
 522   // print(d.poss_base_expr_refs(w, { ['A', 0] }));
 523   // print(d.poss_base_expr_refs(w, { ['A', 0], [Q('A'), 1], [R('A','B'), 2] }));
 524   // print(d.poss_base_expr_refs(w, { ['A', 0], [Q('A'), 1] }));
 525   // 2 ** 81 subsets
 526   // assert(d.poss_base_expr_refs(w, { 'A', Q('A'), R('A','B') }),
 527   //        "Error: is incorrect.");
 528   // print('[ok] poss_base_expr_refs');
 529 };
 530 test_poss_base_expr_refs();
 531 
 532 /*
 533 test_poss_base_refs := procedure() {
 534   d := decision_algorithm.new();
 535   d.cached_lf := { };
 536   w := w1();
 537   w_sts := w1.cm().poss_states();
 538   pbers := { 
 539     { ['A', { w_sts[1], w_sts[2] }],
 540       ['Q', { [{ w_sts[1], w_sts[2] }], 
 541               [{ w_sts[4] }] }] },
 542     { ['A', { w_sts[3], w_sts[4] }],
 543       ['Q', { [{ w_sts[2], w_sts[3] }], 
 544               [{ w_sts[5] }] }] }
 545   };  
 546   // assert(d.poss_base_refs(w, b, pbers) == ,
 547   //         "Error: poss_base_refs is incorrect.");
 548   print('[ok] poss_base_refs');
 549 };
 550 // test_poss_base_refs();
 551 */
 552 
 553 test_rr := procedure() {
 554   d := decision_algorithm.new();
 555   b := causal_markov_model.new();
 556   b.v := { 'V_1' };
 557   b.r := { ['V_1', tf()] };
 558   b.f := {
 559     ['V_1', { [{ ['V_1', false] },  true],
 560               [{ ['V_1',  true] }, false]
 561             }]
 562   };
 563   b.n := 3; // bug when set to 2?
 564   // rr(expr, [w, b, atom_ref, ti, context])
 565   assert(d.rr('A', [b, b, {}, 1, { ['A', { {['V_1_t1', true]} }] }]) ==
 566                                          { {['V_1_t1', true]} },
 567          "Error: rr context is incorrect.");
 568   print('[ok] rr context');
 569 
 570   d.i := { 'P(I_1)' };
 571   d.p := { ['I_1','P(I_1)'] };
 572   b.a := [ {['V_1_t1',  true]},
 573            {['V_1_t2', false]},
 574            {['V_1_t3',  true]} ];
 575   d.cached_lf := {
 576     [ {['V_1', false]}, {['P(I_1)', false]} ],
 577     [ {['V_1',  true]}, {['P(I_1)',  true]} ]
 578   };
 579   all_sts := b.cm().poss_states();
 580   assert(d.rr('I_1', [b, b, {}, 1, {}]) ==
 581            compat_compl_sts( {['V_1_t1', true]}, all_sts ),
 582          "Error: rr i is incorrect.");
 583   print('[ok] rr i');
 584 
 585   assert(d.rr(At_t('I_1', -1), [b, b, {}, 2, {}]) ==
 586            compat_compl_sts( {['V_1_t1', true]}, all_sts ),
 587          "Error: rr At_t is incorrect.");
 588   print('[ok] rr At_t');
 589 
 590   atom_ref := { ['B', compat_compl_sts( {['V_1_t1', true]}, all_sts )],
 591                 ['IsAorB', { [compat_compl_sts({['V_1_t1', true]}, all_sts)],
 592                              [compat_compl_sts({['V_1_t2', true]}, all_sts)] }]
 593               };
 594   assert(d.rr('B', [b, b, atom_ref, 1, {}]) ==
 595            compat_compl_sts({['V_1_t1', true]}, all_sts) &&
 596          d.rr('IsAorB', [b, b, atom_ref, 1, {}]) == {
 597            [compat_compl_sts({['V_1_t1', true]}, all_sts)],
 598            [compat_compl_sts({['V_1_t2', true]}, all_sts)]
 599          },
 600          "Error: rr atom_ref is incorrect.");
 601   print('[ok] rr atom_ref');
 602 
 603   atom_ref := { ['C', compat_compl_sts( {['V_1_t1', true]}, all_sts )],
 604                 ['D', compat_compl_sts( {['V_1_t2', true]}, all_sts )] };
 605   assert(d.rr(And('C', 'D'), [b, b, atom_ref, 1, {}]) ==
 606            compat_compl_sts({['V_1_t1', true], ['V_1_t2', true]}, all_sts),
 607          "Error: rr And is incorrect");
 608   print('[ok] rr And');
 609 
 610   assert(d.rr(Or('C', 'D'), [b, b, atom_ref, 1, {}]) ==
 611            compat_compl_sts({['V_1_t1', true]}, all_sts) +
 612            compat_compl_sts({['V_1_t2', true]}, all_sts),
 613          "Error: rr Or is incorrect");
 614   print('[ok] rr Or');
 615 
 616   assert(d.rr(Not('C'), [b, b, atom_ref, 1, {}]) ==
 617            compat_compl_sts({['V_1_t1', false]}, all_sts),
 618          "Error: rr Not is incorrect");
 619   print('[ok] rr Not');
 620 
 621   assert(d.rr(Implies('C', 'D'), [b, b, atom_ref, 1, {}]) ==
 622            compat_compl_sts({['V_1_t1', false]}, all_sts) +
 623            compat_compl_sts({['V_1_t2', true]}, all_sts),
 624          "Error: rr Implies is incorrect");
 625   print('[ok] rr Implies');
 626 
 627   print('[  ] rr');
 628 };
 629 test_rr();
 630 
 631 test_true_p := procedure() {
 632   d := decision_algorithm.new();
 633   b_event := {['V', true]};
 634   ref_i := { [ b_event, { {['X', true]}, {['Y', true]} }] };
 635   w := causal_markov_model.new();
 636   // TODO: rewrite to use states rather than events
 637 // or just fold into sq_err
 638   w.a := [{ ['Y', true] }];
 639   assert(d.true_p(w, ref_i, b_event) == 1, "Error: true_p true is incorrect.");
 640   print('[ok] true_p true');
 641 
 642   w.a := [{ ['X', true] }];
 643   assert(d.true_p(w, ref_i, b_event) == 1, "Error: true_p true is incorrect.");
 644   print('[ok] true_p true 2');
 645 
 646   w.a := [{ ['X', false] }];
 647   assert(d.true_p(w, ref_i, b_event) == 0, "Error: true_p false is incorrect.");
 648   print('[ok] true_p false');
 649 };
 650 //test_true_p();
 651 
 652 test_sq_err := procedure() {
 653   d := decision_algorithm.new();
 654   d.p := { ['A', 'P(A)'] };
 655   d.cached_lf := { [ {['V',  true]}, {['P(A)', 0.75]} ],
 656                    [ {['V', false]}, {['P(A)', 0.33]} ] };
 657   b := causal_markov_model.new();
 658   b.a := [{['V_t1', true]}, {['V_t2', false]}];
 659   ref_i := { [{['V_t1',  true]}, { {['V_t1', true], ['V_t2', false]},
 660                                    {['V_t1', true], ['V_t2',  true]} }],
 661              [{['V_t2', false]}, { {['V_t1', true], ['V_t2', false]},
 662                                    {['V_t1', true], ['V_t2',  true]} }] };
 663   assert(abs(d.sq_err(b, b, ref_i) - 0.5114) < 0.01,
 664          "Error: sq_err is incorrect.");
 665   print('[ok] true_p');
 666   print('[ok] sq_err');
 667 };
 668 test_sq_err();
 669 test_true_p := procedure() { test_sq_err(); };
 670 
 671 b2 := procedure() {
 672   b := b1();
 673   b.n := 2;
 674   b.a := [
 675     { ['U_t1', 1], ['V_P1_t1', false], ['V_P2_t1', false], ['O_t1',  true] },
 676     { ['U_t2', 2], ['V_P1_t2',  true], ['V_P2_t2', false], ['O_t2', false] }
 677   ];
 678   return b;
 679 };
 680 
 681 test_ref_expr := procedure() {
 682   d := d2();
 683   b := b2();
 684   bs1 := {{['U', 1], ['V_P1_t1', false], ['V_P2_t1', false], ['O', false]}};
 685   d.cached_ref := {
 686     [ {['V_P1_t1', false], ['V_P2_t1', false]}, bs1 ]
 687   };
 688   assert(d.ref_expr(b, b, 1, 'A') == bs1,
 689          "Error: ref_expr is incorrect.");
 690   print('[ok] ref_expr');
 691 };
 692 test_ref_expr();
 693 
 694 test_actual_p_events := procedure() {
 695   assert(d2().actual_p_events(b2()) == {
 696            { ['V_P1_t1', false], ['V_P2_t1', false] },
 697            { ['V_P1_t2',  true], ['V_P2_t2', false] }
 698          },
 699          "Error: actual_p_events is incorrect.");
 700   print('[ok] actual_p_events');
 701 };
 702 test_actual_p_events();
 703 
 704 test_prob_states := procedure() {
 705   assert(d2().prob_states(b2()) == { [1, {['P', 1]}], [2, {['P', 0]}] },
 706          "Error: prob_states is incorrect.");
 707   print('[ok] prob_states');
 708   print('[ok] prob_states_t');
 709 };
 710 test_prob_states();
 711 test_prob_states_t := procedure() { test_prob_states(); };
 712 
 713 d3 := procedure() {
 714   d := decision_algorithm.new();
 715   d.u := { ['A', 'U(A)'] };
 716   d.o := { 'O' };
 717   d.n := [{ ['e', {['B','E1(B)']}], ['m2',{'M2_1'}], ['o',{'O_1'}] },
 718           { ['e', {['C','E2(C)']}], ['m2',{'M2_2'}], ['o',{'O_2'}] }];
 719   return d;
 720 };
 721 
 722 test_eo := procedure() {
 723   assert(d3().eo() == { [{ ['A',  'U(A)'] }, { 'O' }],
 724                         [{ ['B', 'E1(B)'] }, {'O_1'}],
 725                         [{ ['C', 'E2(C)'] }, {'O_2'}] },
 726          "Error: eo is incorrect.");
 727   print('[ok] eo');
 728 };
 729 test_eo();
 730 
 731 test_instr_irrat := procedure() {
 732   print('test instr irrat');
 733   d := decision_algorithm.new();
 734   d.o := { 'O' };
 735   d.p := { ['A', 'P(A)'],
 736            [Not('A'), 'P(~A)'],
 737            [BoxArrow('A', 'B'), 'P(A []-> B)'],
 738            [BoxArrow(Not('A'), 'B'), 'P(~A []-> B)'] };
 739   d.u := { ['B', 'U(B)'] };
 740   d.r := { ['P(A)', {0.3, 0.7}],
 741            ['P(~A)', {0.1, 0.9}],
 742            ['P(A []-> B)', {0.2, 0.8}],
 743            ['P(~A []-> B)', {0.25, 0.75}],
 744            ['U(B)', { 5 }], ['O', tf()] };
 745   b := causal_markov_model.new();
 746   b.u := {};
 747   b.v := { 'V_A', 'V_~A', 'V_A_[]->_B', 'V_~A_[]->_B', 'V_O', 'V_U' };
 748   b.r := { [var, tf()] : var in b.v };
 749   false_func := proc2func(procedure(bs) { return false; }, b.poss_states());
 750   print('false func');
 751   b.f := {
 752     ['V_A', false_func],
 753     ['V_~A', false_func],
 754     ['V_A_[]->_B', false_func],
 755     ['V_~A_[]->_B', false_func],
 756     ['V_O', false_func],
 757     ['V_U', false_func]
 758   };
 759   b.n := 3;
 760   b.a := [
 761     { ['V_A_t1',  true], ['V_~A_t1', false],
 762       ['V_~A_[]->_B_t1', false], ['V_~A_[]->_B_t1', true],
 763       ['V_O_t1', false], ['V_U_t1', false] },
 764     { ['V_A_t2',  true], ['V_~A_t2', false],
 765       ['V_~A_[]->_B_t2', false], ['V_~A_[]->_B_t2', true],
 766       ['V_O_t2', false], ['V_U_t2', false] },
 767     { ['V_A_t3',  true], ['V_~A_t3', false],
 768       ['V_~A_[]->_B_t3', false], ['V_~A_[]->_B_t3', true],
 769       ['V_O_t3', false], ['V_U_t3', false] }
 770   ];
 771   d.cached_lf := {
 772     [{['V_A', false]},         {['P(A)', 0.3]}        ],
 773     [{['V_A',  true]},         {['P(A)', 0.7]}        ],
 774     [{['V_~A', false]},        {['P(~A)', 0.1]}       ],
 775     [{['V_~A',  true]},        {['P(~A)', 0.9]}       ],
 776     [{['V_A_[]->_B', false]},  {['P(A []-> B)', 0.2]} ],
 777     [{['V_A_[]->_B',  true]},  {['P(A []-> B)', 0.8]} ],
 778     [{['V_~A_[]->_B', false]}, {['P(~A []-> B)', 0.2]}],
 779     [{['V_~A_[]->_B',  true]}, {['P(~A []-> B)', 0.8]}],
 780     [{['V_O', false]},         {['O', false]}         ],
 781     [{['V_O',  true]},         {['O',  true]}         ],
 782     [{['V_U', false]},         {['U(B)', false]}      ],
 783     [{['V_U',  true]},         {['U(B)',  true]}      ]
 784   };
 785   print('cached lf');
 786   d.cached_f := d.plf_to_pf(b, d.cached_lf);
 787   print('cached f');
 788   b_cm := b.cm();
 789   print('b cm');
 790   b_cm_poss_states := b_cm.poss_states();
 791   print(#b_cm_poss_states);
 792   o_t2 := compat_compl_sts({['V_O_t2', true]}, b_cm_poss_states);
 793   print('o_t2');
 794   o_t3 := compat_compl_sts({['V_O_t3', true]}, b_cm_poss_states);
 795   no_t2 := compat_compl_sts({['V_O_t2', false]}, b_cm_poss_states);
 796   no_t3 := compat_compl_sts({['V_O_t3', false]}, b_cm_poss_states);
 797   print('no_t3');
 798   d.cached_ref := {
 799     [{['V_A_t1', true ]}, o_t2],
 800     [{['V_A_t1', false]}, o_t2],
 801 
 802     [{['V_~A_t1', true ]}, no_t2],
 803     [{['V_~A_t1', false]}, no_t2],
 804 
 805     [{['V_A_t2', true ]}, o_t3],
 806     [{['V_A_t2', false]}, o_t3],
 807 
 808     [{['V_~A_t2', true ]}, no_t3],
 809     [{['V_~A_t2', false]}, no_t3]
 810 
 811   };
 812   // print(d.instr_irrat(b, b));
 813   // assert(d.instr_irrat() == ,
 814   //       "Error: instr_irrat is incorrect.");
 815   print('[.?] instr_irrat');
 816 };
 817 // test_instr_irrat();
 818 
 819 test_max_p_dist := procedure() {
 820   d := decision_algorithm.new();
 821   d.p := { ['A', 'P(A)'], ['B', 'P(B)'] };
 822   thirds := { 0, 0.33, 0.66, 1 };
 823   d.r := { ['P(A)', thirds], ['P(B)', thirds] };
 824   b := causal_markov_model.new();
 825   b.v := { 'V_A', 'V_B' };
 826   four := { 1,2,3,4 };
 827   b.r := { ['V_A', four], ['V_B', four] };
 828   d.cached_f := {
 829     [{['V_A', 1], ['V_B', 1]}, {['P(A)',    0], ['P(B)',    0]}],
 830     [{['V_A', 1], ['V_B', 2]}, {['P(A)',    0], ['P(B)', 0.33]}],
 831     [{['V_A', 1], ['V_B', 3]}, {['P(A)',    0], ['P(B)', 0.66]}],
 832     [{['V_A', 1], ['V_B', 4]}, {['P(A)',    0], ['P(B)',    1]}],
 833     [{['V_A', 2], ['V_B', 1]}, {['P(A)', 0.33], ['P(B)',    0]}],
 834     [{['V_A', 2], ['V_B', 2]}, {['P(A)', 0.33], ['P(B)', 0.33]}],
 835     [{['V_A', 2], ['V_B', 3]}, {['P(A)', 0.33], ['P(B)', 0.66]}],
 836     [{['V_A', 2], ['V_B', 4]}, {['P(A)', 0.33], ['P(B)',    1]}],
 837     [{['V_A', 3], ['V_B', 1]}, {['P(A)', 0.66], ['P(B)',    0]}],
 838     [{['V_A', 3], ['V_B', 2]}, {['P(A)', 0.66], ['P(B)', 0.33]}],
 839     [{['V_A', 3], ['V_B', 3]}, {['P(A)', 0.66], ['P(B)', 0.66]}],
 840     [{['V_A', 3], ['V_B', 4]}, {['P(A)', 0.66], ['P(B)',    1]}],
 841     [{['V_A', 4], ['V_B', 1]}, {['P(A)',    1], ['P(B)',    0]}],
 842     [{['V_A', 4], ['V_B', 2]}, {['P(A)',    1], ['P(B)', 0.33]}],
 843     [{['V_A', 4], ['V_B', 3]}, {['P(A)',    1], ['P(B)', 0.66]}],
 844     [{['V_A', 4], ['V_B', 4]}, {['P(A)',    1], ['P(B)',    1]}]
 845   };
 846   b.a := [{['V_A_t1', 2], ['V_B_t1', 3]}];
 847   assert(d.max_p_dist(b, 1) == 1.33,
 848          "Error: max_p_dist is incorrect.");
 849   print('[ok] max_p_dist');
 850 };
 851 test_max_p_dist();
 852 
 853 test_poss_cms := procedure() {
 854   d := decision_algorithm.new();
 855   d.p := { ['A', 'P(A)'], ['B', 'P(B)'], ['C', 'P(C)'] };
 856   d.r := { [var, tf()] : var in range(d.p) };
 857   assert(exists(cm in d.poss_cms() |
 858            cm.f == {
 859              ["P(A)", { [{["P(B)", false]},  true],
 860                         [{["P(B)",  true]}, false] }],
 861              ["P(B)", { [{["P(C)", false]}, false],
 862                         [{["P(C)",  true]},  true] }]
 863            }) &&
 864          exists(cm in d.poss_cms() |
 865            cm.f == {
 866              ["P(A)", { [{["P(B)", false]}, false],
 867                         [{["P(B)",  true]}, false]} ],
 868              ["P(C)", { [{["P(A)", false]},  true],
 869                         [{["P(A)",  true]},  true] }]
 870            }),
 871          "Error: poss_cms is incorrect.");
 872   print('[ok] poss_cms');
 873 };
 874 test_poss_cms();
 875 test_poss_cms_from_l      := procedure() { test_poss_cms(); };
 876 test_poss_cms_from_p_pars := procedure() { test_poss_cms(); };
 877 
 878 test_e_func := procedure() {
 879   d := d3();
 880   ds := { ['U(A)', 5], ['E1(B)', 3], ['E2(C)', 8] };
 881   assert(d.e_func(ds, 0) == { ['U(A)',  5] } &&
 882          d.e_func(ds, 1) == { ['E1(B)', 3] } &&
 883          d.e_func(ds, 2) == { ['E2(C)', 8] },
 884          "Error: e_func is incorrect.");
 885   print('[ok] e_func');
 886 };
 887 test_e_func();
 888 
 889 test_ue := procedure() {
 890   d := d3();
 891   assert(d.ue(0) == d.u && d.ue(1) == d.n[1]['e'] && d.ue(2) == d.n[2]['e'],
 892          "Error: ue is incorrect.");
 893   print('[ok] ue');
 894 };
 895 test_ue();
 896 
 897 test_ue_base_exprs := procedure() {
 898   d := d3();
 899   assert(d.ue_base_exprs(0) == { 'A' } &&
 900          d.ue_base_exprs(1) == { 'B' } &&
 901          d.ue_base_exprs(2) == { 'C' },
 902          "Error: ue_base_exprs is incorrect.");
 903   print('[ok] ue_base_exprs');
 904 };
 905 test_ue_base_exprs();
 906 
 907 /*
 908 test_state_space := procedure() {
 909   d := decision_algorithm.new();
 910   d.p := { ['A', 'P(A)'], ['B', 'P(B)'] };
 911   assert(d.state_space(range(d.p)) == { 
 912            {Not('A'), Not('B')},
 913            {Not('A'),     'B' },
 914            {    'A' , Not('B')},
 915            {    'A' ,     'B' }
 916          },
 917          "Error: state_space is incorrect.");
 918   print('[ok] state_space');
 919   assert(d.event_space(range(d.p)) == d.state_space(range(d.p)) + {
 920            { Not('A') }, { 'A' }, { Not('B') }, { 'B' }, {}
 921          },
 922          "Error: event_space is incorrect.");
 923   print('[ok] event_space');
 924   // print(d.expr_space(d.state_space(range(d.p))));
 925   assert(d.expr_space( d.state_space(range(d.p)) ) == {
 926            And("A",      "B"), And(    "A" , Not("B")), 
 927            And("B", Not("A")), And(Not("A"), Not("B"))
 928          },
 929          "Error: expr_space is incorrect.");
 930   print('[ok] expr_space');
 931 
 932   d.p[And("A", "B")] := 'P(A^B)';
 933   d.p[And("A", Not("B"))] := 'P(A^~B)';
 934   d.p[And("B", Not("A"))] := 'P(B^~A)';
 935   d.p[And(Not("A"), Not("B"))] := 'P(~A^~B)';
 936   d.r := { ['P(A^B)', {0,1}], ['P(A^~B)', {0,1}], 
 937            ['P(B^~A)', {0,1}], ['P(~A^~B)', {0,1}] };
 938   assert({["P(A^B)", 1], ["P(A^~B)", 0], ["P(B^~A)", 0], ["P(~A^~B)", 0]}
 939            in d.poss_p_vs({'P(A)', 'P(B)'}),
 940          "Error: poss_p_vs is incorrect");
 941   print('[ok] poss_p_vs');
 942 };
 943 test_state_space();
 944 test_event_space := procedure() { test_state_space(); };
 945 test_expr_space := procedure() { test_expr_space(); };
 946 test_poss_p_vs := procedure() { test_poss_p_vs(); };
 947 */
 948 
 949 // deprecated?
 950 /*
 951 test_p_of_vs := procedure() {
 952   d := decision_algorithm.new();
 953   d.p := { ['A', 'P(A)'], ['B', 'P(B)'], ['C', 'P(C)'] };
 954   cm := causal_model(u, v, r, f, om);
 955   cm.u := { 'P(A)', 'P(B)' };
 956   cm.v := { 'P(C)' };
 957   cm.r := { [var, tf()] : var in cm.u + cm.v };
 958   cm.f := { 
 959     ['P(C)', { [ {['P(A)', false], ['P(B)', false]}, false],
 960                [ {['P(A)', false], ['P(B)',  true]},  true],
 961                [ {['P(A)',  true], ['P(B)', false]},  true],
 962                [ {['P(A)',  true], ['P(B)',  true]},  true]
 963              }]
 964   };
 965   assert(d.p_of_vs(cm, { ['P(A)', 0.33], ['P(B)', 0.5] }) == ,
 966          "Error: p_of_vs is incorrect.");
 967   print('[ok] p_of_vs');
 968 };
 969 test_p_of_vs();
 970 */
 971 
 972 test_coh_pcms := procedure() {
 973   cm := causal_model.new();
 974   cm.u := { 'A' };
 975   cm.v := { 'B' };
 976   cm.r := { ['A', tf()], ['B', tf()] };
 977   cm.f := { ['B', { [{ ['A', false] },  true],
 978                     [{ ['A',  true] }, false] }] };
 979   d := decision_algorithm.new();
 980   d.p := { ['A', 'P(A)'] };
 981   d.r := { ['P(A)', { 0, 1 }] };
 982 //  print(d.coh_pcms(cm));
 983 //  assert(d.coh_pcms(cm) == ,
 984 //         "Error: coh_pcms is incorrect.");
 985   print('[.?] coh_pcms');
 986 };
 987 // test_coh_pcms();
 988 
 989 /* Template
 990 test_ := procedure() {
 991   d := decision_algorithm.new();
 992   assert(d.() == ,
 993          "Error:  is incorrect.");
 994   print('[ok] ');
 995 };
 996 test_();
 997 */
 998 
 999 print(' ');
1000 
1001 test_expr_to_p_funcs := procedure() {
1002   d := decision_algorithm.new();
1003   d.p := { ['A', 'P(A)'], ['B', 'P(B)'] };
1004   assert(d.expr_to_p_funcs(Implies('A', 'B'), {'A', 'B' }) == {
1005            {["P(A)", false], ["P(B)", false]},
1006            {["P(A)", false], ["P(B)",  true]},
1007            {["P(A)",  true], ["P(B)",  true]}
1008          }, "Error: expr_to_p_funcs is incorrect.");
1009   print('[ok] expr_to_p_funcs');
1010 
1011   d.p[Implies('A', 'B')] := 'P(A->B)';
1012   d.p[Not(And('A', Not('B')))] := 'P(~(A^~B))';
1013   assert(d.equiv_p_expr(d.p, d.cp, 'P(A->B)', 'P(~(A^~B))') == true,
1014          "Error: equiv_p_expr is incorrect");
1015   print('[ok] equiv_p_expr true');
1016 };
1017 test_expr_to_p_funcs();
1018 
1019 test_sum_p := procedure() {
1020   d := decision_algorithm.new();
1021   d.p := { ['A', 'P(A)'], ['B', 'P(B)'], ['C', 'P(C)'] };
1022   p1 := { ['P(A)', 0.6], ['P(C)', 0.3] };
1023   p2 := { ['P(B)', 0.2], ['P(C)', 0.5] };
1024   assert(d.sum_p(p1, p2) == { ['P(A)', 0.6], ['P(B)', 0.2], ['P(C)', 0.8] },
1025          "Error: sum_p is incorrect.");
1026   print('[ok] sum_p');
1027 };
1028 test_sum_p();
1029 
1030 test_add_p_distr := procedure() {
1031   d := decision_algorithm.new();
1032   d.p := { ['A', 'P(A)'], ['B', 'P(B)'], ['C', 'P(C)'] };
1033   d.cp := { [['C', 'A'], 'P(C|A)'], [['C', 'B'], 'P(C|B)'] };
1034   p_cp1 := [{ ['P(A)', 0.2] }, { ['P(C|A)', 0.3], ['P(C|B)', 0.4] }];
1035   p_cp2 := [{ ['P(A)', 0.4] }, { ['P(C|A)', 0.1], ['P(C|B)', 0.3] }];
1036   p_cp := d.add_p_distr(p_cp1, p_cp2);
1037   assert((p_cp[2])['P(C|A)'] == 0.4 && (p_cp[2])['P(C|B)'] == 0.7 &&
1038          (p_cp[1])['P(A)'] - 0.6 < 0.01,
1039          "Error: add_p_distr is incorrect.");
1040   print('[ok] add_p_distr');
1041 };
1042 test_add_p_distr();
1043 
1044 test_prob_distance := procedure() {
1045   d := decision_algorithm.new();
1046   d.p := { ['A', 'P(A)'], ['B', 'P(B)'] };
1047   d.cp := { [['B', 'A'], 'P(B|A)'] };
1048   pp1 := { ['P(A)', 0.4], ['P(B)', 0.5], ['P(B|A)', 0.6] };
1049   pp2 := { ['P(A)', 0.1], ['P(B)', 0.9], ['P(B|A)', 0.2] };
1050   assert(d.prob_distance(d.p, d.cp, pp1, pp2) == 1.1,
1051          "Error: prob_distance is incorrect.");
1052   print('[ok] prob_distance');
1053 };
1054 test_prob_distance();
1055 
1056 test_poss_io_states := procedure() {
1057   d := decision_algorithm.new();
1058   d.i := { 'P(I)' };
1059   d.o := { 'O' };
1060   d.r := { ['P(I)', {0,1}], ['O', tf()]};
1061   assert(d.poss_io_states() == {
1062           { ['P(I)', 0], ['O', false] },
1063           { ['P(I)', 0], ['O',  true] },
1064           { ['P(I)', 1], ['O', false] },
1065           { ['P(I)', 1], ['O',  true] }
1066          },
1067          "Error: poss_io_states is incorrect.");
1068   print('[ok] poss_io_states');
1069 };
1070 test_poss_io_states();
1071 
1072 test_poss_m_states := procedure() {
1073   d := decision_algorithm.new();
1074   d.m := { 'M1', 'M2' };
1075   d.r := { ['M1', tf()], ['M2', tf()] };
1076   assert(d.poss_m_states() == {
1077            { ['M1', false], ['M2', false] },
1078            { ['M1', false], ['M2',  true] },
1079            { ['M1',  true], ['M2', false] },
1080            { ['M1',  true], ['M2',  true] }
1081          },
1082          "Error: poss_m_states is incorrect.");
1083   print('[ok] poss_m_states');
1084 };
1085 test_poss_m_states();
1086 
1087 test_poss_o_states := procedure() {
1088   d := decision_algorithm.new();
1089   d.o := { 'O1', 'O2' };
1090   d.r := { ['O1', tf()], ['O2', tf()] };
1091   assert(d.poss_o_states() == {
1092            { ['O1', false], ['O2', false] },
1093            { ['O1', false], ['O2',  true] },
1094            { ['O1',  true], ['O2', false] },
1095            { ['O1',  true], ['O2',  true] }
1096          },
1097          "Error: poss_o_states is incorrect.");
1098   print('[ok] poss_o_states');
1099 };
1100 test_poss_o_states();
1101 
1102 test_possible_fms := procedure() {
1103   d := decision_algorithm.new();
1104   d.o := { 'O' };
1105   d.m := { 'M_1', 'M_2' };
1106   d.m2 := { 'M2_1', 'M2_2' };
1107   d.r := { [var, tf()] : var in d.o + d.m + d.m2 };
1108   // d.r := { [var, tf()] : var in d.o + d.m };
1109   ds := decision_algorithm.new();
1110   ds.o := { 'O' };
1111   ds.m := { 'ds.M_1', 'ds.M_2' };
1112   ds.m2 := { 'ds.M2_1', 'ds.M2_2' };
1113   ds.r := { [var, tf()] : var in ds.o + ds.m + ds.m2 };
1114   // ds.r := { [var, tf()] : var in ds.o + ds.m };
1115   p_fms := [fm : fm in d.possible_fms(ds)];
1116   print(p_fms);
1117   // print('[.?] possible_fms');
1118 };
1119 // test_possible_fms();
1120 test_apply_fm_m := procedure() { test_possible_fms(); };
1121 test_fm_commutes := procedure() { test_possible_fms(); };
1122 
1123 test_force_list := procedure() {
1124   d := decision_algorithm.new();
1125   d.p := { ['A', 'P(A)'] };
1126   d.cp := { [['B', 'A'], 'P(B|A)'] };
1127   assert(d.force_list('P(A)'  ) == ['A',  om] &&
1128          d.force_list('P(B|A)') == ['B', 'A'],
1129          "Error: force_list  is incorrect.");
1130   print('[ok] force_list ');
1131 };
1132 test_force_list();
1133 
1134 test_p_i_sim := procedure() {
1135   d := decision_algorithm.new();
1136   d.p := { ['A', 'P(A)'] };
1137   d.r := { ['P(A)', {0, 0.25, 0.75, 1}] };
1138   ds1 := { ['P(A)', 0.25] };
1139   ds2 := { ['P(A)', 0] };
1140   assert(abs(d.p_i_sim(ds1, ds2, 'P(A)') - 2 / 3) < 0.01,
1141          "Error: p_i_sim is incorrect.");
1142   print('[ok] p_i_sim');
1143 };
1144 test_p_i_sim();
1145 
1146 test_id_expr := procedure() {
1147   d := decision_algorithm.new();
1148   d.i := {'P(I)'};
1149   d.p := {['A', 'P(A)'], ['I', 'P(I)'] };
1150   assert(d.id_expr('A') == 'A' &&
1151          d.id_expr('I') == At_t('I', -1) &&
1152          d.id_expr(At_t('I', -1)) == At_t('I', -2),
1153          "Error: id_expr is incorrect.");
1154   print('[ok] id_expr');
1155 
1156   d.p := d.p + { [At_t('I', -1), 'P(I_-1)'],
1157                  [At_t('I', -2), 'P(I_-2)'] };
1158   assert(d.id_var('P(I)')    == 'P(I_-1)' &&
1159          d.id_var('P(I_-1)') == 'P(I_-2)',
1160          "Error: id_var is incorrect.");
1161   print('[ok] id_var');
1162 };
1163 test_id_expr();
1164 test_id_var := procedure() { test_id_expr(); };
1165 
1166 test_poss_i_events := procedure() {
1167   d := decision_algorithm.new();
1168   d.i := { 'P(I)', 'P(I2)' };
1169   d.r := { ['P(I)', {0,1}], ['P(I2)', {0,1}] };
1170   assert(d.poss_i_events() == {
1171            {["P(I)", 0]}, {["P(I)", 1]},
1172            {["P(I2)", 0]}, {["P(I2)", 1]},
1173            {["P(I)", 0], ["P(I2)", 0]},
1174            {["P(I)", 0], ["P(I2)", 1]},
1175            {["P(I)", 1], ["P(I2)", 0]},
1176            {["P(I)", 1], ["P(I2)", 1]}
1177          }, "Error: poss_i_events is incorrect.");
1178   print('[ok] poss_i_events');
1179 };
1180 test_poss_i_events();
1181 
1182 test_e_vars := procedure() {
1183   d := decision_algorithm.new();
1184   d.n := [{ ['e', {['A','E_n1(A)']}], ['m2', {'M2_n1'}], ['o', {'O_n1'}] },
1185           { ['e', {['B','E_n2(B)']}], ['m2', {'M2_n2'}], ['o', {'O_n2'}] } ];
1186   assert(d.e_vars() == { 'E_n1(A)', 'E_n2(B)' },
1187          "Error: e_vars  is incorrect.");
1188   print('[ok] e_vars ');
1189 };
1190 test_e_vars();
1191 
1192 test_poss_ue_events := procedure() {
1193   d := decision_algorithm.new();
1194   d.n := [{ ['e', {['A','E']}], ['m2', {'M2'}], ['o', {'O'}] }];
1195   d.u := { ['B', 'U'] };
1196   d.r := { ['E', {0,1}], ['U', {0,1}] };
1197   assert(d.poss_ue_events() == {
1198            {["E", 0]}, {["E", 0], ["U", 0]}, {["E", 0], ["U", 1]}, {["E", 1]},
1199            {["E", 1], ["U", 0]}, {["E", 1], ["U", 1]}, {["U", 0]}, {["U", 1]}
1200          }, "Error: poss_ue_events is incorrect.");
1201   print('[ok] poss_ue_events');
1202 };
1203 test_poss_ue_events();
1204 
1205 test_poss_p_i_events := procedure() {
1206   d := decision_algorithm.new();
1207   d.p := { ['A', 'P(A)'], ['B', 'P(B)'] };
1208   d.r := { ['P(A)', {0,1}], ['P(B)', {0,1}] };
1209   assert(d.poss_p_i_events() == { ["P(A)", 0], ["P(A)", 1], ["P(B)", 0], ["P(B)", 1] },
1210          "Error: poss_p_i_events is incorrect.");
1211   print('[ok] poss_p_i_events');
1212 };
1213 test_poss_p_i_events();
1214 
1215 test_poss_rx_states := procedure() {
1216   d := decision_algorithm.new();
1217   d.n := [{ ['o', {'O_n1'}] }, { ['o', {'O_n2'}] }];
1218   d.r := { ['O_n1', tf()], ['O_n2', tf()] };
1219   assert(d.poss_rx_states() == {
1220            {["O_n1", false], ["O_n2", false]},
1221            {["O_n1", false], ["O_n2",  true]},
1222            {["O_n1",  true], ["O_n2", false]},
1223            {["O_n1",  true], ["O_n2",  true]}
1224          }, "Error: poss_rx_states is incorrect.");
1225   print('[ok] poss_rx_states');
1226 };
1227 test_poss_rx_states();
1228 
1229 test_poss_rxs := procedure() {
1230   d := decision_algorithm.new();
1231   d.n := [{ ['o', {'O_n1_1', 'O_n1_2'}] }, { ['o', {'O_n2'}] } ];
1232   d.r := { [no_var, tf()] : no_var in {'O_n1_1', 'O_n1_2', 'O_n2'} };
1233   assert(d.poss_rxs() == {
1234            {["O_n1_1", false], ["O_n1_2", false]},
1235            {["O_n1_1", false], ["O_n1_2",  true]},
1236            {["O_n1_1",  true], ["O_n1_2", false]},
1237            {["O_n1_1",  true], ["O_n1_2",  true]},
1238            {["O_n2", false]},
1239            {["O_n2",  true]}
1240          }, "Error: poss_rxs is incorrect.");
1241   print('[ok] poss_rxs');
1242 };
1243 test_poss_rxs();
1244 
1245 /*
1246 test_n_disj := procedure() {
1247   d := decision_algorithm.new();
1248   b := causal_markov_model.new();
1249   b.u := { 'U1' };
1250   b.r := { ['U1', tf()] };
1251   be := { ['U1', true] };
1252   not_be := { ['U1', false] };
1253   w := b;
1254   d.cached_ref := { [be, { be, not_be }] };
1255   assert(d.n_disj(w, b, be) == 2,
1256          "Error: n_disj is incorrect.");
1257   print('[ok] n_disj');
1258 
1259   assert(d.n_ev_disj(w, b, not_be, be) == 1,
1260          "Error: n_ev_disj is incorrect.");
1261   print('[ok] n_ev_disj');
1262 
1263   assert(d.spec(w, b, not_be, be) == 1/2,
1264          "Error: spec is incorrect.");
1265   print('[ok] spec');
1266 };
1267 test_n_disj();
1268 test_n_ev_disj := procedure() { test_n_disj(); };
1269 test_spec := procedure() { test_spec(); };
1270 */
1271 
1272 test_ltr_up_to := procedure() {
1273   assert(da().ltr_up_to('A', 3) == { 'A_1', 'A_2', 'A_3' },
1274          "Error: ltr_up_to is incorrect.");
1275   print('[ok] ltr_up_to');
1276 };
1277 test_ltr_up_to();
1278 
1279 test_vars_up_to := procedure() {
1280   assert(da().vars_up_to({ 'A', 'B' }, 2) ==
1281            { 'A_1', 'A_2', 'B_1', 'B_2' },
1282          "Error: vars_up_to is incorrect.");
1283   print('[ok] vars_up_to');
1284 };
1285 test_vars_up_to();
1286 
1287 test_expr_to_ltr_var := procedure() {
1288   assert(da().expr_to_ltr_var(Implies('A', 'B'), 'U') == 'U(Implies(A, B))',
1289          "Error: expr_to_ltr_var is incorrect.");
1290   print('[ok] expr_to_ltr_var');
1291 };
1292 test_expr_to_ltr_var();
1293 
1294 test_poss_r_z  := procedure() {
1295   poss_z_vals := { { j : j in [1..x] } : x in [1..3] };
1296   assert(da().poss_r_z(poss_z_vals, da().ltr_up_to('A', 2)) == {
1297            {["A_1", {1}],       ["A_2", {1}]},
1298            {["A_1", {1}],       ["A_2", {1, 2}]},
1299            {["A_1", {1}],       ["A_2", {1, 2, 3}]},
1300            {["A_1", {1, 2}],    ["A_2", {1}]},
1301            {["A_1", {1, 2}],    ["A_2", {1, 2}]},
1302            {["A_1", {1, 2}],    ["A_2", {1, 2, 3}]},
1303            {["A_1", {1, 2, 3}], ["A_2", {1}]},
1304            {["A_1", {1, 2, 3}], ["A_2", {1, 2}]},
1305            {["A_1", {1, 2, 3}], ["A_2", {1, 2, 3}]}},
1306          "Error: poss_r_z is incorrect.");
1307   print('[ok] poss_r_z ');
1308 };
1309 test_poss_r_z();
1310 
1311 test_poss_r_ue := procedure() {
1312   u := { ['A', 'U(A)'], ['B', 'U(B)'] };
1313   poss_ue_vals := pow({ j : j in [1..2**2] }) - {{}};
1314   assert({ {["U(A)", {1, 2, 4}], ["U(B)", {1, 2, 4}]},
1315            {["U(A)", {2, 3, 4}], ["U(B)", {2, 3, 4}]}
1316          } < da().poss_r_ue(poss_ue_vals, u),
1317          "Error: poss_r_ue is incorrect.");
1318   print('[ok] poss_r_ue');
1319 };
1320 test_poss_r_ue();
  1 class causal_markov_model(n, u, v, r, f, a, cached_cm) {
  2   n := n; // integer number of time slices
  3   u := u; // Set of exogenous variables
  4           //   e.g. u := { 'O' };
  5   v := v; // Set of endogenous variables
  6           //   e.g. v := { 'F', 'L', 'ML' };
  7   r := r; // Function assigning a range of values for each variable
  8           //   e.g. r := { [x, { false, true }] : x in u + v };
  9   /* Higher order function from y in v to a higher order function from a
 10      function assigning values to exogenous and endogenous vars to a value 
 11      for y.
 12        eg f := { [y, f_y], [z, f_z], ... }
 13          f_y := { [f_y_1, true], [f_y_2, false], ... } 
 14           f_y_1 := { [y, false], [z, false], ... }       (row in truth table)
 15      A variable (technically, its past time slice) can also be its own parent
 16   */
 17   f := f;
 18 
 19   /* A list of actual states, one for each time slice. Each actual state is a 
 20      function assigning values to cm() vars of the given time.
 21        e.g. [ {['y_t1',true],...}, {['y_t2',false],...}, ... ] */
 22   a := a;
 23 
 24   uv := procedure() {
 25     return u + v;
 26   };
 27 
 28   u_t := procedure(t) {
 29     return { x + '_t' + t : x in u };
 30   };
 31   v_t := procedure(t) {
 32     return { x + '_t' + t : x in v };
 33   };
 34   uv_t := procedure(t) {
 35     return u_t(t) + v_t(t);
 36   };
 37 
 38   cm := procedure() {
 39     if (cached_cm == om) {
 40       cached_cm := compute_cm();
 41     }
 42     return cached_cm;
 43   };
 44 
 45     compute_cm := procedure() {
 46       return causal_model(cm_u(), cm_v(), cm_r(), cm_f(), om);
 47     };
 48 
 49     cm_u := procedure() {
 50       u_ts := { u_t(t) : t in [1..n] };
 51       if (u_ts == {}) {
 52         u_ts := {{}};
 53       }
 54       return +/ u_ts + v_t(1);
 55       // Includes the first v_t since there's no previous time parents
 56     };
 57 
 58     cm_v := procedure() {
 59       return +/ { v_t(t) : t in [2..n] };
 60       // Starts at 2 because the first, parentless v_t is counted in cm_u
 61     };
 62 
 63     cm_r := procedure() {
 64       return { [x, r[orig_var(x)]] : x in cm_u() + cm_v() };
 65     };
 66 
 67     cm_f := procedure() {
 68       return { [x, fix_time(x, f[orig_var(x)], var_time(x) - 1)] : x in cm_v() };
 69     };
 70 
 71       fix_time := procedure(x, f_x, t) {
 72         // Add time subscript to f_x
 73         f_prev_t := { [assign_time(g, t), f_x[g]] : g in domain(f_x) };
 74         // Extend from one timeslice to states over all vars but x
 75         ext_xs := (cm_u() + cm_v()) - uv_t(t) - {x}; // extended vars
 76         poss_ext_states := possible_states(ext_xs, cm_r()); // extended states
 77         // Output remains f_prev_t[prev_t_st] since it's a markov model
 78         return +/ { { [prev_t_st + ext_state, f_prev_t[prev_t_st] ]
 79                       : ext_state in poss_ext_states }
 80                     : prev_t_st in domain(f_prev_t)
 81                   };
 82       };
 83 
 84   poss_states := procedure() {
 85     return possible_states(u + v, r);
 86   };
 87 
 88   poss_events := procedure() {
 89     return +/ { 2 ** ps : ps in poss_states() };
 90   };
 91 
 92   actual_events := procedure() {
 93     return pow(actual_state()) - {{}};
 94   };
 95 
 96   actual_state := procedure() {
 97     return +/ { actual_state : actual_state in a };
 98   };
 99 
100   /* Set of all actual synchronic events, 
101      i.e. those taking place at a single time. */
102   actual_sync_events := procedure() {
103     return +/ { pow(actual_state) : actual_state in a } - {{}};
104   };
105 
106   /* Use the causal model on the first two times */
107   // g needs to be full state?
108   response :=  procedure(ys, g) {
109     g1 := assign_time(g, 1);
110     ys2 := { y + '_t' + 2 : y in ys };
111     return drop_time(cm().response(ys2, g1));
112   };
113 
114   /* Form the full state transition table and convert it to a string along with
115      the vars and ranges. */
116   t_to_str := procedure() {
117     tt := { [ps, response(v, ps)] : ps in poss_states() };
118     return str([tt, n, u, v, r]);
119   };
120 
121   static {
122     new := procedure() {
123       return causal_markov_model(0, {}, {}, {}, {}, [], om);
124     };
125   }
126 }
  1 print("\nTesting Causal Markov Model");
  2 print(  "===========================");
  3 
  4 cmm1 := cachedProcedure() {
  5   cmm := causal_markov_model.new();
  6   cmm.n := 3;
  7   cmm.u := { 'U' };
  8   cmm.v := { 'V' };
  9   cmm.r := { ['U', tf()], ['V', tf()] };
 10   cmm.f := { // eg v_t2 = U_t1 xor V_t1
 11     ['V', { [{ ['U', false], ['V', false] }, false],
 12             [{ ['U', false], ['V',  true] },  true],
 13             [{ ['U',  true], ['V', false] },  true],
 14             [{ ['U',  true], ['V',  true] }, false]
 15           }]
 16   };
 17   cmm.a := {
 18     {['U_t1',  true], ['V_t1', false]},
 19     {['U_t2',  true], ['V_t2',  true]},
 20     {['U_t3', false], ['V_t3',  true]}
 21   };
 22   return cmm;
 23 };
 24 
 25 test_cm_u := procedure() {
 26   assert(cmm1().cm_u() == { 'U_t1', 'U_t2', 'U_t3', 'V_t1' },
 27          "Error: cm_u is incorrect.");
 28   print('[ok] cm_u');
 29 };
 30 test_cm_u();
 31 
 32 test_cm_v := procedure() {
 33   assert(cmm1().cm_v() == { 'V_t2', 'V_t3' },
 34          "Error: cm_v is incorrect.");
 35   print('[ok] cm_v');
 36 };
 37 test_cm_v();
 38 
 39 test_cm_r := procedure() {
 40   assert(cmm1().cm_r() == { ['U_t1', tf()], ['U_t2', tf()], ['U_t3', tf()],
 41                             ['V_t1', tf()], ['V_t2', tf()], ['V_t3', tf()] },
 42          "Error: is incorrect.");
 43   print('[ok] cm_r');
 44 };
 45 test_cm_r();
 46 
 47 test_cm_f := procedure() {
 48   assert( { [ { ["U_t1", false], ["U_t2", false], ["U_t3", false],
 49                 ["V_t1", false],                  ["V_t3", false]  }, false],
 50             [ { ["U_t1", false], ["U_t2", false], ["U_t3", false],
 51                 ["V_t1", false],                  ["V_t3",  true]  }, false]
 52           } < cmm1().cm_f()['V_t2'],
 53           "Error: cm_f is incorrect.");
 54   print('[ok] cm_f');
 55 };
 56 test_cm_f();
 57 
 58 test_poss_states := procedure() {
 59   assert(cmm1().poss_states() == {
 60            { ["U", false], ["V", false] },
 61            { ["U", false], ["V",  true] },
 62            { ["U",  true], ["V", false] },
 63            { ["U",  true], ["V",  true] }
 64          },
 65          "Error: poss_states is incorrect.");
 66   print('[ok] poss_states');
 67 };
 68 test_poss_states();
 69 
 70 test_poss_events := procedure() {
 71   assert(cmm1().poss_events() == { {},
 72            {["U", false]}, {["U", true]}, {["V", false]}, {["V", true]},
 73            {["U", false], ["V", false]},
 74            {["U", false], ["V",  true]},
 75            {["U",  true], ["V", false]},
 76            {["U",  true], ["V",  true]}
 77          },
 78          "Error: poss_events is incorrect.");
 79   print('[ok] poss_events');
 80 };
 81 test_poss_events();
 82 
 83 test_actual_events := procedure() {
 84   assert(// test a random event
 85          { ["U_t2", true], ["U_t3", false], ["V_t1", false], ["V_t2", true]
 86          } in cmm1().actual_events(),
 87          "Error: actual_events is incorrect.");
 88   print('[ok] actual_events');
 89   print('[ok] actual_state');
 90 };
 91 test_actual_events();
 92 test_actual_state := procedure() { test_actual_state(); };
 93 
 94 test_actual_sync_events := procedure() {
 95   assert(// test some random events
 96          !({ ["U_t2", true], ["U_t3", false], ["V_t1", false], ["V_t2", true]
 97            } in cmm1().actual_sync_events()) &&
 98          { {['U_t2', true]}, {['U_t3',false], ['V_t3',true]}
 99            } < cmm1().actual_sync_events() &&
100          !({['U_t1', false]} in cmm1().actual_sync_events()),
101          "Error: actual_sync_events is incorrect.");
102   print('[ok] actual_sync_events');
103 };
104 test_actual_sync_events();
105 
106 test_response := procedure() {
107   assert(cmm1().response({'V'}, {['U', true], ['V', false]})['V'] == true &&
108          cmm1().response({'V'}, {['U', true], ['V',  true]})['V'] == false,
109          "Error: response is incorrect.");
110   print('[ok] response');
111 };
112 test_response();
113 
114 test_t_to_str := procedure() {
115   assert(eval(cmm1().t_to_str()) == [
116            // tt
117            { [{["U", false], ["V", false]}, {["V", false]}],
118              [{["U", false], ["V",  true]}, {["V",  true]}],
119              [{["U",  true], ["V", false]}, {["V",  true]}],
120              [{["U",  true], ["V",  true]}, {["V", false]}] },
121            3, // n
122            {"U"},
123            {"V"},
124            {["U", {false, true}], ["V", {false, true}]} // r
125          ],
126          "Error: t_to_str is incorrect.");
127   print('[ok] t_to_str');
128 };
129 test_t_to_str();
  1 /* Adapted from Judea Pearl's Causal Models
  2      e.g. Pearl, Judea. Causes and Explanations: A Structural-Model Approach
  3           http://ftp.cs.ucla.edu/pub/stat_ser/R266-part1.pdf
  4 */
  5 class causal_model(u, v, r, f, cached_poss_states) {
  6   u := u; // Set of exogenous variables
  7           //   e.g. u := { 'O' };
  8   v := v; // Set of endogenous variables
  9           //   e.g. v := { 'F', 'L', 'ML' };
 10   r := r; // Function assigning a range of values for each variable
 11           //   e.g. r := { [x, { false, true }] : x in u + v };
 12   /* Higher order function from y in v to a higher order function from a
 13      function assigning values to exogenous and other endogenous vars to a 
 14      value for y.
 15        eg f := { [y, f_y], [z, f_z], ... }
 16             f_y := { [f_y_1, true], [f_y_2, false], ... } // rows in truth table
 17               f_y_1 := { [x, false], [z, false], ... }
 18   */
 19   f := f;
 20 
 21   uv := procedure() {
 22     return u + v;
 23   };
 24 
 25   /* Takes an assignment of values to endogenous variables and returns a 
 26      causal submodel. See 7.12 Pearl's Causality book.  */
 27   effect_of_do := procedure(g) {
 28     h := { [y, f[y]] : y in v | !(y in domain(g)) } +
 29          { [x, i_proc_to_func(x, procedure(h) {
 30                                    return g[x];
 31                                  })] : x in domain(g) };
 32     return causal_model(u, v, r, h, om);
 33   };
 34 
 35   /* Takes a set of variables ys and an assignment of values g to endogenous 
 36      and exogenous variables and returns an assignment of values to ys */
 37   response := procedure(ys, g) {
 38     sm := effect_of_do(g);
 39     result := { [y, sm.solve_for(y, g)] : y in ys };
 40     if (dev) { print(result); }
 41     return result;
 42   };
 43 
 44   /* Takes a variable y and an assignment of values g to exogenous variables 
 45      and returns a value for y */
 46   solve_for := procedure(y, g) {
 47     g := { [x, g[x]] : x in domain(g) | x in u };
 48     if (y in u) {
 49       return g[u];
 50     } else if (#parents(y) == 0) {
 51       return last(first(f[y]));
 52     } else { // y in v
 53       // todo: sometimes a conflict? 
 54       h := g + { procedure() {
 55                    if (z in parents(y)) {
 56                      return [z, solve_for(z, g)];
 57                    } else { // doesn't matter, pick first
 58                      return [z, first(r[z])];
 59                    }
 60                  }() : z in (u + v) - domain(g) | z != y };
 61       return f[y][h];
 62     }
 63   };
 64 
 65   // Takes a variable and returns its parents
 66   parents := procedure(y) {
 67     return { z : z in u + v | z != y && is_nec_parent_of(z, y) };
 68   };
 69 
 70     /* Check whether there exists an assignment of values to the other 
 71        variables such that there are still different values of y depending 
 72        on values of z.  */
 73     is_nec_parent_of := procedure(z, y) {
 74       xs := {x : x in u + v | x != y && x != z}; // all other variables but z
 75       g_xs := possible_states(xs, r);
 76       return exists(g_x in g_xs | values_diff(g_x, y, z));
 77     };
 78 
 79     values_diff := procedure(g_x, y, z) {
 80       vals := { f[y][ g_x + {[z,r_i]} ] : r_i in r[z] };
 81       return #vals > 1;
 82     };
 83 
 84   /* Instance method shorthands for call to static version */
 85   i_proc_to_func := procedure(y, p) {
 86     return causal_model.proc_to_func(y, u, v, r, p);
 87   };
 88 
 89   poss_states := procedure() {
 90     if (cached_poss_states == om) {
 91       cached_poss_states := compute_poss_states();
 92     }
 93     return cached_poss_states;
 94   };
 95 
 96     compute_poss_states := procedure() {
 97       return possible_states(u + v, r);
 98     };
 99 
100   poss_events := cachedProcedure() {
101     return +/ { pow(ps) : ps in poss_states() };
102   };
103 
104   sorted_v := procedure() {
105     return sort_list( [y : y in v],
106                       procedure(y, z) { return depth(y) < depth(z); } );
107   };
108 
109     depth := procedure(y) {
110       if (y in u) {
111         return 0;
112       } else {
113         return 1 + max({ depth(z) : z in parents(y) });
114       }
115     };
116 
117   static {
118     /* Helper method that takes an endogenous variable y, & u, v, r and a 
119        procedure and returns a function. The procedure should take an 
120        assignment of values to exogenous variables and other endogenous 
121        variables and return y's value. */
122     proc_to_func := procedure(y, u, v, r, p) {
123       dom := possible_states({ x : x in u + v | x != y }, r);
124       return { [g, p(g)] : g in dom };
125     };
126 
127     new := procedure() {
128       return causal_model({}, {}, {}, {}, om);
129     };
130   }
131 }
 1 print("\nTesting Causal Model");
 2 print(  "====================");
 3 
 4 cm1 := procedure() {
 5   u := { 'O' };
 6   v := { 'F', 'L', 'ML' };
 7   r := { [x, { false, true }] : x in u + v };
 8   f_f := causal_model.proc_to_func('F', u, v, r, procedure(h) {
 9            if (h['L'] || h['ML']) {
10              return true;
11            } else {
12              return false;
13            }
14          });
15   f := { ['F', { [{ ['O', false], ['L', false], ['ML', false] }, false],
16                  [{ ['O', false], ['L', false], ['ML', true ] }, true ],
17                  [{ ['O', false], ['L', true ], ['ML', false] }, true ],
18                  [{ ['O', false], ['L', true ], ['ML', true ] }, true ],
19                  [{ ['O', true ], ['L', false], ['ML', false] }, false],
20                  [{ ['O', true ], ['L', false], ['ML', true ] }, true ],
21                  [{ ['O', true ], ['L', true ], ['ML', false] }, true ],
22                  [{ ['O', true ], ['L', true ], ['ML', true ] }, true ] }],
23          ['L',  causal_model.proc_to_func('L', u, v, r, procedure(g) {
24                   return false;
25                 })],
26          ['ML', causal_model.proc_to_func('ML', u, v, r, procedure(g) {
27                   return false;
28                 })] };
29   assert(f_f == f['F'], 'proc_to_func helper method is incorrect');
30   print('[ok] proc_to_func');
31 
32   cm := causal_model(u,v,r,f);
33   assert(cm.f['F'] == f_f, 'causal model class is incorrect');
34 
35   return cm;
36 };
37 test_proc_to_func := procedure() { cm1(); };
38 
39 test_effect_of_do := procedure() {
40   sm := cm1().effect_of_do({ ['F', true], ['L', true] });
41   assert((sm.f['F' ])[{ ['O', false], ['L', false], ['ML', false] }] <==> true,
42            'effect_of_do F is incorrect');
43   assert((sm.f['L' ])[{ ['O', false], ['F', false], ['ML', false] }] <==> true,
44            'effect_of_do L is incorrect');
45   assert((sm.f['ML' ])[{ ['O', false], ['F', false], ['L', false] }] <==> false,
46            'effect_of_do ML is incorrect');
47   print('[ok] effect_of_do');
48 };
49 
50 test_parents := procedure() {
51   assert(cm1().parents('F') == { 'L', 'ML' }, 'parents is incorrect');
52   print('[ok] parents');
53 };
54 
55 test_response := procedure() {
56   assert(cm1().response({'F'}, { ['L', false], ['O', true] })['F'] == false,
57          'response is incorrect');
58   assert(cm1().response({'F'}, { ['L', true], ['O', true] })['F'] == true,
59          'response is incorrect');
60   print('[ok] response');
61 };
  1 dev := false;
  2 
  3 tf := procedure() { return { true, false }; };
  4 
  5 // Non-empty powerset. Powerset without the empty set.
  6 ne_pow := procedure(s) {
  7   return pow(s) - {{}};
  8 };
  9 
 10 min := procedure(l) {
 11   mins := { e : e in l | forall(e2 in [x :  x in l | x != e ] | e2 >= e) };
 12   return first(mins);
 13 };
 14 
 15 max := procedure(l) {
 16   maxs := { e : e in l | forall(e2 in [x :  x in l | x != e ] | e >= e2) };
 17   return first(maxs);
 18 };
 19 
 20 avg := procedure(l) {
 21   return (+/ l)  / #l;
 22 };
 23 
 24 index_of := procedure(e, xs) {
 25   i := 1;
 26   while (xs[i] != e) {
 27     i := i + 1;
 28   }
 29   return i;
 30 };
 31 
 32 equiv_expr := procedure(e1, e2, dom) {
 33   if (dev) { print('equiv_expr dom ' + dom); }
 34   if (isList(e1) && !isList(e2) ||
 35       isList(e2) && !isList(e1)   ) {
 36     return false;
 37   } else if (isList(e1)) {
 38     return forall(j in [1..#e1] | equiv_expr(e1[j], e2[j], dom) );
 39   } else {
 40     // todo: handle quantifiers, []->?
 41     return expr_to_events(e1, dom) == expr_to_events(e2, dom);
 42   }
 43 };
 44 
 45 // Possible references for arity
 46 old_range_for_arity := cachedProcedure(e_arity, poss_events) {
 47   /* Singleton events, each being a set of just one event, a trivial 
 48      disjunction of truth-making events given our convention here of using 
 49      sets of events to express the disjunction of them. */
 50   sngl_events := { {event} : event in poss_events };
 51   if (e_arity == 0) {
 52     return sngl_events;
 53   /* Wrapping the sngl event below as a list/1-tuple might be a little 
 54      counterintuitive but it allows it to be treated consistently with
 55      higher arities */
 56   } else if (e_arity == 1) {
 57       return 2 ** { [se] : se in sngl_events };
 58   } else {
 59       return 2 ** iter_prod(sngl_events, e_arity);
 60   }
 61 };
 62 
 63 // Possible references for arity of expr
 64 range_for_arity := cachedProcedure(e_arity, poss_sts) {
 65   /* Propositions are sets of possible states. They can also be treated as 
 66      objects. You might think of them as the essential realization conditions
 67      for objects. */
 68   poss_objects := pow(poss_sts);
 69   // handling empty sets?
 70   if (e_arity == 0) {
 71     return poss_objects;
 72   /* Wrapping the sngl event below as a list/1-tuple might be a little 
 73      counterintuitive but it allows it to be treated consistently with
 74      higher arities */
 75   } else if (e_arity == 1) {
 76       return pow({ [x] : x in poss_objects - {{}} });
 77   } else {
 78       return pow(iter_prod(poss_objects - {{}}, e_arity));
 79       // return 2 ** iter_prod(sngl_events, e_arity);
 80   }
 81 };
 82 
 83 connectives := procedure() {
 84   return {'And','Not','Or','Implies','BoxArrow','At_t','Forall','Exists'};
 85 };
 86 
 87 /* Unused?
 88 fcts_n_strs := procedure(expr) {
 89   if (isTerm(expr)) {
 90     return ({ fct(expr) } - connectives()) + 
 91            (+/ { fcts_n_strs(expr) : arg_i in args(expr) });
 92   } else {
 93     return { expr };
 94   }
 95 };
 96 */
 97 
 98 maybe_not := procedure(var, val) {
 99   if (val == true) {
100     return var;
101   } else {
102     return Not(var);
103   }
104 };
105 
106 inv_maybe_not := procedure(expr) {
107   if (isTerm(expr) && fct(expr) == "Not") {
108     return { [expr, false] };
109   } else {
110     return { [expr,  true] };
111   }
112 };
113 
114 /* Takes an expr and returns the set of vars, atomic propositions and objects
115    that occur in it.  */
116 atomic_exprs := procedure(expr) {
117   if (isTerm(expr)) {
118     return +/ { atomic_exprs(arg_i) : arg_i in args(expr) };
119   } else {
120     return { expr };
121   }
122 };
123 
124 /* Takes an expr and returns the set of atomic propositions and objects
125    that occur in it (but not its variables).  */
126 atomic_terms := procedure(expr) {
127   return atomic_exprs(expr) - parse_vars(expr);
128 };
129 
130 /* Takes a function assigning true or false to vars and turns it into an 
131    expression that is a conjunction of either var or Not(var) depending 
132    on each variable's truth value.  */
133 func_to_expr := procedure(g) {
134   sorted_vs := sort([ v : v in domain(g) ]);
135   return iter_fct("And", [ maybe_not(v, g[v]) : v in sorted_vs ]);
136 };
137 
138 all_atomic_formulas := procedure(exprs) {
139   return +/ { atomic_formulas(expr) : expr in exprs };
140 };
141 
142 atomic_formulas := procedure(expr) {
143   return atomic_forms(expr) - parse_vars(expr);
144 };
145 
146 atomic_forms := procedure(expr) {
147   if (dev) { print('atomic forms ' + expr); }
148   if (!isTerm(expr) || !(fct(expr) in connectives())) {
149     return { expr };
150   } else {
151     return +/ { atomic_forms(arg_i) : arg_i in args(expr) };
152   }
153 };
154 
155 expr_to_events := procedure(expr, dom) {
156   state_space := exprs_to_state_space(dom);
157   if (dev) { print('st sp: ' + state_space); }
158   expand := procedure(events) {
159     if (dev) { print('expand ' + events); }
160     return { state : state in state_space
161              | exists(event in events | restricted_eq(state, event)) };
162   };
163   ete := procedure(expr) {
164     if (dev) { print('ete(' + expr + ')'); }
165     if (isString(expr)) {
166       return { state : state in state_space | state[expr] == true };
167     }
168     if (isTerm(expr)) {
169       if (#args(expr) >= 1) { a1 := args(expr)[1]; }
170       if (#args(expr) >= 2) { a2 := args(expr)[2]; }
171 
172       if (!(fct(expr) in connectives())) {  // applied predicate / relation
173         return { [expr, true] };
174       } else if (fct(expr) == "And") {
175         return expand(ete(a1)) * expand(ete(a2));
176         return expand(ete(a1)) * expand(ete(a2));
177         // return { g + h : [g,h] in ete(a1) >< ete(a2) | isMap(g + h) };
178       } else if (fct(expr) == "Or") {
179         return expand(ete(a1)) + expand(ete(a2));
180       // a -> b  ==  ~a v b
181       } else if (fct(expr) == "Implies") {
182         return expand(ete( Or(Not(a1), a2) ));
183       } else if (fct(expr) == "Exists") {
184         // Ex Qx <-> ~ Vx ~Qx
185         return expand(ete( Not( Forall(a1, Not(a2)) ) ));
186       } else if (fct(expr) == "Forall") {
187         terms := { var : var in dom | isString(var) };
188         return { state : state in state_space |
189           forall(term in terms | is_true_in(sub(a2, a1, term), state, dom)) };
190       } else if (fct(expr) == "Not") {
191         if (dev) { print('Not'); }
192         switch {
193           case isString(a1) :
194             if (dev) { print('Not(str)'); }
195             return { {[a1, false]} };
196           // Pair, i.e. { [var, val] }
197           case isSet(a1) && #a1 == 1 && isMap(a1) :
198             return { { [first(a1)[1], !first(a1)[2]] } };
199           case isSet(a1) && isMap(a1) :
200             if (dev) { print('~(a ^ b) = ~a v ~b'); }
201             // ~(a ^ b) = ~a v ~b
202             return +/ { expand(ete(Not({pair}))) : pair in a1 };
203           // Singleton map, eg { g } or { {[x,y],[w,z],...} }
204           case isSet(a1) && #a1 == 1 :
205             if (dev) { print('~ { g } = ~g'); }
206             // ~ { g } = ~g
207             return expand(ete(Not(first(a1))));
208           // Implicit disjunction of maps, eg { g, h, ... }
209           case isSet(a1) :
210             if (dev) { print('~(a v b) = ~a ^ ~b'); }
211             // ~(a v b) = ~a ^ ~b
212             head := first(a1);
213             return expand(ete(  And( Not(head), Not(a1 - {head}) )  ));
214           // Do we allow Not(a[]->b)? How?
215           // case isTerm(a1) && fct(a1) == 'BoxArrow' :
216           default : // e.g. functor term, i.e. predicate or other connective?
217             if (dev) { print('[~a] = [~[a]]'); }
218             // [~a] = [~[a]]
219             return expand(ete(Not( expand(ete(a1)))));
220         }
221       }
222     } else { // not term, atomic expr
223       return ({ [expr, true] });
224     }
225   }; // end ete
226   return expand(ete(expr));
227 };
228 
229 /* https://en.wikipedia.org/wiki/Fold_(higher-order_function)  */
230 foldr := procedure(proc, z, l) {
231   if (#l == 0) {
232     return z;
233   } else if (#l == 1) {
234     return proc(l[1], z);
235   } else {
236     return proc(l[1], foldr(proc, z, l[2..]));
237   }
238 };
239 
240 /* Takes a setlx functor and returns a procedure which takes two expressions
241    and if either is om, returns the other, or else constructs a term out of
242    the functor and two expressions.  */
243 fctize := procedure(fct) {
244   return procedure(e1, e2) {
245     if (e1 == om) { return e2; } else if (e2 == om) { return e1; } else {
246       return makeTerm(fct, [e1,e2]);
247     }
248   };
249 };
250 
251 /* Takes a setlx functor and list of expressions and constructs an expression
252    by repeatedly applying the functor.  */
253 iter_fct := procedure(fct, l) {
254   return foldr(fctize(fct), om, l);
255 };
256 
257 /* Takes a function f and returns its inverse. To avoid having to return a set 
258    or list, it just returns the first preimage that f maps to the img. */
259 inv := procedure(f) {
260   return { [img, ([ pre : pre in domain(f) | f[pre] == img])[1] ] : img in range(f) };
261 };
262 
263 /* Takes a setlx functor and returns its arity */
264 arity := procedure(fct) {
265   if (isTerm(fct)) {
266     return #args(fct);
267   } else {
268     return 0;
269   }
270 };
271 
272 /* Takes a set like s >< s >< s and returns tuples instead of nested pairs */
273 flatten := procedure(set) {
274   return { flatten_list(e) : e in set | isList(e)} +
275          { e : e in set | !isList(e) };
276 };
277 
278 flatten_list := procedure(list) {
279   if (isList(list)) {
280     return +/ [ flatten_list(l) : l in list ];
281   } else {
282     return [list];
283   }
284 };
285 
286 /* Takes a list of sets and returns the set of tuples taking a value from 
287    each set.  */
288 cart_prod := procedure(list) {
289   if (#list == 1) {
290     return first(list);
291   } else {
292     return flatten(first(list) >< cart_prod(list[2..]));
293   }
294 };
295 
296 /* Iterated Cartesian Product 
297    Yields the cartesian product of a set with itself n times.  */
298 iter_prod := procedure(set, n) {
299   if (n == 1) {
300     return set;
301   } else {
302     return flatten(iter_prod(set, n - 1) >< set);
303   }
304 };
305 
306 iter_mult := procedure(set) {
307   if (#set == 1) {
308     return first(set);
309   } else {
310     return first(set) * iter_mult(set - { first(set) });
311   }
312 };
313 
314 /* Takes a set of sets and returns a set of all possible choices taking an 
315    element from each set. */
316 choices := procedure(set) {
317   list := [ [elem : elem in subset] : subset in set ];
318   return set_choices(list);
319 };
320 
321 // Same as list_choices but returns a set rather than list of the choices.
322 set_choices := procedure(r) {
323   return { { o : o in p } : p in list_choices(r) };
324 };
325 
326 /* Takes a list of lists and returns a list of all possible choices taking
327    an element from each list  */
328 list_choices := procedure(r) {
329   if (#r == 1) {
330     return [[x] : x in r[1]];
331   } else {
332     return +/ [  [[x] + row : row in list_choices(r[2..])]  : x in r[1] ];
333   }
334 };
335 
336 /* Generate the set of all functions from dom to rng. 
337    Generate atomic maps. Collect the choices of those atomic maps. */
338 function_space := procedure(dom, rng) {
339   // Collect into a set the elements from a list of choices of atomic maps
340   // p for preimage and permutation? o is ordered pair?
341   return { { o : o in p } : p in list_choices([ [ [p, i] : i in rng ] : p in dom]) };
342 };
343 
344 /* Takes a set of variables s and a function r assigning a set which is the 
345    range of possible values to each variable and returns the set of possible 
346    functions assigning values to those variables.  */
347 possible_states := cachedProcedure(s, r) {
348   if (s == {} || r == {}) {
349     return {{}};
350   }
351   return choices({ {v} >< r[v] : v in s });
352 };
353 
354 /* Takes a function and returns all partial functions, including empty.  */
355 partial_functions := cachedProcedure(f) {
356   return 2 ** f;
357 };
358 
359 /* Takes two functions assigning values to vars and tests whether they are
360    compatible, i.e. do not assign different values to the same var.  */
361 compatible_values := procedure(f, g) {
362   h := f + g;
363   return forall(var in domain(f) + domain(g) | h[var] != om);
364 };
365 compat := procedure(f, g) { return compatible_values(f, g); }; // alias
366 
367 /* Takes partial state ps, which is a function assigning values to a set of
368    variables, and set of complete states cs, a set of functions assigning 
369    values to a larger set of variables and returns the subset of cs compatible 
370    with ps.  */
371 compatible_complete_states := cachedProcedure(ps, cs) {
372   return { c : c in cs
373            | forall(d in domain(ps)
374                     | c[d] == ps[d] )
375          };
376 };
377 compat_compl_sts := procedure(ps, cs) { // alias 
378   return compatible_complete_states(ps, cs);
379 };
380 
381 proc2func := procedure(p, dom) {
382   return { [d, p(d)] : d in dom };
383 };
384 
385 // Breaks apart a var label of form 'V_ti' into ['V', i]
386 parse_var := procedure(label) {
387   j := #label;
388   while (int(label[j]) != om && j != 1) {
389     j -= 1;
390   }
391   t := int(label[(j + 1)..]);
392   assert(t != om, "Error: Expected var label to have time subscript.");
393   var := label[1..(j - 2)];
394   return [var, t];
395 };
396 
397 orig_var := procedure(label) {
398   return parse_var(label)[1];
399 };
400 
401 var_time := procedure(label) {
402   return parse_var(label)[2];
403 };
404 
405 /* Takes an assignment of values to variables and adds a time subscript to the 
406    variables' labels. */
407 assign_time := procedure(g, t) {
408   return { [ y+'_t'+t, g[y] ] : y in domain(g) };
409 };
410 
411 /* Takes an assignment of values to variables and removes time subscripts 
412    from the variables' labels. */
413 drop_time := procedure(g) {
414   return { [ orig_var(y), g[y] ] : y in domain(g) };
415 };
416 
417 /* From Hilbert System, Wikipedia 
418      https://en.wikipedia.org/w/index.php?title=Hilbert_system&oldid=873786806
419 */
420 is_logical_truth := procedure(expr) {
421   vars := parse_vars(expr);
422 
423   /* Axiom schemes
424      Redundant: 
425        P1: phi -> phi
426        as_p1 := procedure(phi) { return Implies(phi, phi); };  
427       
428      P2: phi -> (psi -> phi)  */
429   as_p2 := procedure(phi, psi) { return Implies(phi, Implies(psi, phi)); };
430   // P3: (phi -> (psi -> xi)) -> ((phi -> psi) -> (phi -> xi))
431   as_p3 := procedure(phi, psi, xi) {
432     return Implies(
433              Implies(phi, Implies(psi, xi)),
434              Implies(Implies(phi, psi), Implies(phi, xi))
435            );
436   };
437   // P4: (~phi -> ~psi) -> (psi -> phi)
438   as_p4 := procedure(phi, psi) {
439     return Implies(Implies(Not(phi), Not(psi)), Implies(psi, phi));
440   };
441   // Q5: Vx(phi) -> phi[x := t] where t may be substituted for x in phi
442   as_q5 := procedure(var, phi) {
443     return { Implies( Forall(var, phi), sub(phi, var, term))
444              : term in atomic_terms(expr) };
445   };
446   // Q6: Vx(phi -> psi) -> (Vx(phi) -> Vx(psi))
447   as_q6 := procedure(var, phi, psi) {
448     return Implies( Forall(var, Implies(phi, psi)),
449                     Implies(Forall(var, phi), Forall(var, psi)) );
450   };
451   // Q7: phi -> Vx phi where x not in phi
452   as_q7 := procedure(phi, var) {
453     if (var in atomic_exprs(phi)) {
454       return {};
455     } else {
456       return { Implies(phi, Forall(var, phi)) };
457     }
458   };
459   /* Conservative Extensions
460      Conjunction:
461        Introduction: phi -> (psi -> (phi ^ psi))  */
462   conj_intr := procedure(phi, psi) {
463     return Implies(phi, Implies(psi, And(phi, psi)));
464   };
465   //   Elimination: (phi ^ psi) -> phi   and   (phi ^ psi) -> psi
466   conj_elim := procedure(phi, psi) {
467     return { Implies(And(phi, psi), phi),
468              Implies(And(phi, psi), psi) };
469   };
470   /* Disjunction:
471        Introduction: phi -> (phi v psi)   and   psi -> (phi v psi)  */
472   disj_intr := procedure(phi, psi) {
473     return { Implies(phi, Or(phi, psi)),
474              Implies(psi, Or(phi, psi)) };
475   };
476   //   Elimination: (phi -> xi) -> ((psi -> xi) -> ((phi v psi) -> xi))
477   disj_elim := procedure(phi, psi, xi) {
478     return Implies( Implies(phi, xi),
479                     Implies( Implies(psi, xi),
480                              Implies( Or(phi, psi), xi)));
481   };
482   /* Existential Quantification 
483      Introduction: Vx(phi -> Ey(phi[x := y]))  */
484   ex_intr := procedure(phi, x_var, y_var) {
485     return Forall(x_var, Implies(phi, Exists(y_var, sub(phi, x_var, y_var))));
486   };
487   // Elimination: Vx(phi -> psi) -> (Ex(phi) -> psi) where !(x in psi)
488   ex_elim := procedure(phi, psi, var) {
489     if (var in atomic_exprs(psi)) {
490       return {};
491     } else {
492       return { Implies( Forall(var, Implies(phi, psi)),
493                         Implies( Exists(var, phi), psi)) };
494     }
495   };
496   // Modus Ponens Rule of Inference
497   // phi, phi -> psi |- psi
498   modus_ponens := procedure(phi, psi) {
499     if (isTerm(psi) && fct(psi) == 'Implies' && args(psi)[1] == phi) {
500       return { args(psi)[2] };
501     } else {
502       return {};
503     }
504   };
505 
506   derive_more := procedure(exprs) {
507     return
508       (+/ { as_q7(phi, var) : [phi, var] in exprs >< vars }) +
509       (+/ { { as_p2(phi, psi), as_p4(phi, psi), conj_intr(phi, psi) } +
510               { as_q6(var, phi, psi) : var in vars } +
511               conj_elim(phi, psi) + disj_intr(phi, psi) + modus_ponens(phi, psi)
512             : [phi, psi] in exprs >< exprs }) +
513       { as_p3(phi, psi, xi) : [phi, psi, xi] in iter_prod(exprs, 3) } +
514       (+/ { as_q5(var, phi, vars) : [var, phi] in vars >< exprs }) +
515       (+/ { ex_elim(phi, psi, var)
516             : [phi, psi, var] in cart_prod([exprs, exprs, vars]) } ) +
517       { ex_intr(phi, x_var, y_var)
518             : [phi, x_var, y_var] in cart_prod([exprs,vars,vars]) };
519   };
520   derive := procedure(exprs, j) {
521     if (j == 0) {
522       return exprs;
523     } else {
524       return derive_more(derive(exprs, j - 1));
525     }
526   };
527   return expr in derive(subformulas(expr), 10 ** 6);
528 };
529 
530 parse_vars := procedure(expr) {
531   if (isTerm(expr)) {
532     parse_more := +/ { parse_vars(arg) : arg in args(expr) };
533     if (fct(expr) in {'Forall', 'Exists'}) {
534       return { args(expr)[1] } + parse_more;
535     } else {
536       return parse_more;
537     }
538   } else {
539     return {};
540   }
541 };
542 
543 /* Substitute any occurrences of var in expr with sub_var instead. */
544 sub := procedure(expr, var, sub_var) {
545   if (isTerm(expr)) {
546     return makeTerm(fct(expr), [sub(arg, var, sub_var) : arg in args(expr)]);
547   } else {
548     if (expr == var) {
549       return sub_var;
550     } else {
551       return expr;
552     }
553   }
554 };
555 
556 subformulas := procedure(expr) {
557   if (isTerm(expr) && fct(expr) in connectives()) {
558     return { expr } + (+/ { subformulas(arg) : arg in args(expr) })
559                     - parse_vars(expr);
560   } else {
561     return { expr };
562   }
563 };
564 
565 /* Hausdorff Distance
566    "It is the greatest of all the distances from a point in one set to the 
567     closest point in the other set."   
568       -- https://en.wikipedia.org/wiki/Hausdorff_distance 
569 */
570 hausdorff_dist := procedure(s1, s2, dist_f) {
571   dist_to_closest := procedure(e1) {
572     return min([ dist_f(e1, e2) : e2 in s2 ]);
573   };
574   return max([ dist_to_closest(e1) : e1 in s1 ]);
575 };
576 
577 /* Kendalls Tau Distance
578    Takes two strict ordinal utility functions with the same domain and counts 
579    the number of discordant pairs. Returns a normalized distance from 0 - 1.
580    Each u1, u2 is a set of setlx functors Succ(state1, state2). (Succ stands 
581    for succeeds.) States are logical expressions, i.e. setlx functors of 
582    combinators on string base propositions. 
583 
584    "Kendall tau distance is also called bubble-sort distance since it is 
585     equivalent to the number of swaps that the bubble sort algorithm would 
586     take to place one list in the same order as the other list."
587       -- https://en.wikipedia.org/wiki/Kendall_tau_distance  
588 */
589 kendalls_tau_dist := procedure(u1, u2) {
590   pairs := { [x,y] : [x,y] in states_u(u1) >< states_u(u1) | x != y };
591   // Discordant (unordered) pairs
592   disc_un_pairs := #{ [x,y] : [x,y] in pairs
593                       | (Succ(x,y) in u1 && Succ(y,x) in u2) ||
594                         (Succ(x,y) in u2 && Succ(y,x) in u1)
595                     } / 2;
596   n := #states_u(u1);
597   return disc_un_pairs / (n * (n - 1) / 2);
598 };
599 
600 /* Takes an ordinal utility function expressed as a set of SuccEq setlx
601    functors comparing events or states and returns the set of atomic variables
602    those events or states range over.  */
603 atom_u := procedure(u) {
604   states := +/ { { args(pref)[1], args(pref)[2] } : pref in u };
605   return +/ { domain(state) : state in states };
606 };
607 
608 /* Takes an ordinal utility function and returns the states that are compared
609    in them.  */
610 states_u := procedure(u) {
611   return +/ { { args(pref)[1], args(pref)[2] } : pref in u };
612 };
613 
614 /* Cardinal Utility Function to Ordinal Utility Function */
615 card_u_to_ord_u := procedure(card_u, state_space) {
616   compare_sts := procedure(st1, st2) {
617     if (card_u[st1] >= card_u[st2]) {
618       return SuccEq(st1, st2); // Succeeds or equals
619     } else {
620       return SuccEq(st2, st1);
621     }
622   };
623   return { compare_sts(st1, st2) : [st1,st2] in state_space >< state_space };
624 };
625 
626 /* Additive Utility Function to Cardinal Utility Function: 
627    Takes a utility function expressed as { [expr, val], ... } where the 
628    utility of a state is the sum of all its true exprs' val.
629    Returns a cardinal utility function assigning values to states rather than
630    expressions.
631 */
632 add_u_to_card_u := procedure(add_u) {
633   atoms := +/ { atomic_formulas(expr) : expr in domain(add_u) };
634   state_space := exprs_to_state_space(domain(add_u));
635   util := procedure(state) {
636     sum := +/ { val : [expr, val] in add_u | is_true_in(expr, state, atoms) };
637     if (sum == om) { sum := 0; }
638     return sum;
639   };
640   return { [state, util(state)] : state in state_space };
641 };
642 
643 /* Expressions to State Space:
644    Return the possible ways of assigning true or false to each expr in exprs. */
645 exprs_to_state_space := procedure(exprs) {
646   if (dev) { print('exprs_to_state_space ' + exprs); }
647   atoms := +/ { atomic_formulas(expr) : expr in exprs };
648   return choices({ {[atom,true], [atom,false]} : atom in atoms });
649 };
650 
651 /* Expression is True in State:
652    Takes an expression and a state, a function assigning true/false to the 
653    atomic terms in expr, and returns true or false for whether the expr
654    is true in the state.
655 */
656 is_true_in := procedure(expr, state, dom) {
657   // return state in expr_to_events(expr, dom); // assumes event is state
658   return exists(event in expr_to_events(expr, dom)
659                 | restricted_eq(state, event) );
660 };
661 
662 /* Ordinality Utility Function Distance:
663    Takes two weak ordinal utility functions expressed as a set of SuccEq setlx 
664    functors and returns the distance between them.
665 
666    See hausdorff_dist and kendalls_tau_dist.
667    See also: https://en.wikipedia.org/wiki/Kemeny-Young_method 
668      Interestly, since we read the utility functions off the brains rather
669      than elicit them through votes, we don't need to be concerned with 
670      strategic manipulation through insincere votes.
671 */
672 ord_u_dist := procedure(u1, u2) {
673   atom_us := atom_u(u1) + atom_u(u2);
674   state_space := exprs_to_state_space(atom_us);
675   sharpenings := procedure(u) {
676     poss_total_orders := permutations([s : s in state_space]);
677     list_to_ord := procedure(l) {
678       if (l == []) {
679         return {};
680       } else {
681         return { Succ(l[1], s) : s in l[2..] } + list_to_ord(l[2..]);
682       }
683     };
684     ord_us := { list_to_ord(pto) : pto in poss_total_orders };
685     ext_u := atom_us - atom_u(u); // extending u vars
686     if (ext_u == {}) {
687       // return { u };
688       exts := {{}};
689     } else {
690       exts := choices({ {[atom,true], [atom,false]} : atom in ext_u });
691     }
692     return {
693       ord_u : ord_u in ord_us
694       | forall(pref in u // u: a1 >= a2 
695           | Succ(args(pref)[1], args(pref)[2]) in ord_u || // ord_u: a1 > a2
696             SuccEq(args(pref)[2], args(pref)[1]) in u || // u: a2 >= a1, ie ==
697             forall(ext in exts
698               | Succ(args(pref)[1] + ext, args(pref)[2] + ext) in ord_u
699                   // a1 + z > a2 + z
700             ) // end forall ext
701         ) // end forall pref
702     }; // end return
703   }; // end sharpenings  
704   return hausdorff_dist( sharpenings(u1),
705                          sharpenings(u2),
706                          kendalls_tau_dist );
707 };
708 
709 /* p 56 of setlx tutorial */
710 sort_list := procedure(l, cmp) {
711   if (#l < 2) { return l; }
712   m := #l \ 2;
713   [l1, l2] := [l[.. m], l[m+1 ..]];
714   [s1, s2] := [sort_list(l1, cmp), sort_list(l2, cmp)];
715   return merge(s1, s2, cmp);
716 };
717 
718 sort_set := procedure(s, cmp) {
719   return sort_list([x : x in s], cmp);
720 };
721 
722 merge := procedure(l1, l2, cmp) {
723   match ([l1, l2]) {
724     case [[], _] : return l2;
725     case [_, []] : return l1;
726     case [[x1|r1], [x2|r2]] :
727     if (cmp(x1, x2)) {
728       return [x1 | merge(r1, l2, cmp)];
729     } else {
730       return [x2 | merge(l1, r2, cmp)];
731     }
732   }
733 };
734 
735 setlx_chars := procedure() {
736   str := " !\"#%&'()*+,-./0123456789:;<=>?" +
737          "ABCDEFGHIJKLMNOPQRSTUVWXYZ[\\]_" +
738          "abcdefghijklmnopqrstuvwxyz{|}~";
739   return { str[j] : j in [1..#str] };
740 };
741 
742 all_strs_lte := procedure(s) {
743   poss_strs := {''};
744   length := 0;
745   while (length < #s) {
746     length := length + 1;
747     poss_strs := poss_strs + { poss_str + char : [poss_str, char] in
748                                                  poss_strs >< setlx_chars() };
749   }
750   return poss_strs;
751 };
752 
753 restrict_domain := procedure(event, vars) {
754   return { [var, val] : [var, val] in event | var in vars };
755 };
756 
757 restricted_eq := procedure(state, event) {
758   return event == restrict_domain(state, domain(event));
759 };
760 
761 /* Kolmogorov Complexity of a string
762    https://en.wikipedia.org/wiki/Kolmogorov_complexity  
763 
764    k is an incomputable function but we include a naive attempt to define it
765    anyway. Although it may seem possible this way, it will not work because
766    eval(poss_str) may not terminate. In practice, some finite approximation
767    could be used. For instance, if setlx allowed us to catch a timeout error
768    as many other languages do, it would be particularly easy by just trying 
769    each program string up to a fixed length of time.
770 */
771 k := procedure(s) {
772   poss_strs := {''};
773   success := false;
774   length := 0;
775   while (!success) {
776     length := length + 1;
777     poss_strs := { poss_str + char
778                    : [poss_str, char] in poss_strs >< setlx_chars() };
779     if ({ eval(poss_str) == s : poss_str in poss_strs } ||
780         length == #s + 2 /* '"'+s+'"' */                  )  {
781       success := true;
782     }
783   }
784   return length;
785 };
786 
787 /* Conditional Kolmogorov Complexity  K(s|i)  
788    The length of the shortest program which outputs a string s when given the 
789    string i as input. Again, this is incomputable but can be finitely 
790    approximated.
791 */
792 k_given := procedure(s, i) {
793   eval_str := procedure(str) {
794     try {
795       return eval(str + '("' + i + '");');
796     } catch(e) {
797       return '';
798     }
799   };
800   strs := { str : str in all_str_lte(#s + 2) | eval_str(str) == s };
801   strs := sort(strs, procedure(s1, s2) { return #s1 < #s2; });
802   return #strs[1];
803 };
  1 print('Testing Lib');
  2 print('===========');
  3 
  4 test_ne_pow := procedure() {
  5   assert(ne_pow({1,2}) == {{1}, {2}, {1,2}},
  6          "Error: ne_pow is incorrect.");
  7   print('[ok] ne_pow');
  8 };
  9 test_ne_pow();
 10 
 11 test_min := procedure() {
 12   assert(min([3,4,5]) == 3, "Error: min is incorrect.");
 13   print('[ok] min');
 14 };
 15 test_min();
 16 
 17 test_max := procedure() {
 18   assert(max([3,4,5]) == 5, "Error: max is incorrect.");
 19   print('[ok] max');
 20 };
 21 test_max();
 22 
 23 test_avg := procedure() {
 24   assert(avg([3,4,5]) == 4, "Error: avg is incorrect.");
 25   print('[ok] avg');
 26 };
 27 test_avg();
 28 
 29 test_index_of := procedure() {
 30   assert(index_of(5, [3,4,5]) == 3, "Error: index_of is incorrect.");
 31   print('[ok] index_of');
 32 };
 33 test_index_of();
 34 
 35 test_equiv_expr := procedure() {
 36   e1 := And('A', 'B');
 37   e2 := And('B', 'A');
 38   assert(equiv_expr(e1, e2, { 'A', 'B' }), "Error: equiv_expr is incorrect.");
 39   print('[ok] equiv_expr: A ^ B == B ^ A');
 40 
 41   assert(equiv_expr(Exists('x', Odd('x')),
 42                     Not(Forall('y', Not(Odd('y')))),
 43                     { Odd('a'), Odd('b'), Odd('c') }),
 44          "Error: equiv_expr is incorrect.");
 45   print('[ok] equiv_expr: Ex Odd(x) == ~( Vy ~Odd(y) )');
 46   // print('[ok] expr_to_events');
 47 };
 48 test_equiv_expr();
 49 
 50 test_maybe_not := procedure() {
 51   assert(maybe_not('A', true) == 'A',
 52          "Error: maybe_not true is incorrect.");
 53   print('[ok] maybe_not true');
 54   assert(maybe_not('A', false) == Not('A'),
 55          "Error: maybe_not false is incorrect.");
 56   print('[ok] maybe_not false');
 57 };
 58 test_maybe_not();
 59 
 60 test_atomic_exprs := procedure() {
 61   assert(atomic_exprs(Exists('x', And('B', Odd('a')))) == { 'x', 'B', 'a' },
 62          "Error: atomic_exprs is incorrect.");
 63   print('[ok] atomic_exprs');
 64 };
 65 test_atomic_exprs();
 66 
 67 test_atomic_terms := procedure() {
 68   assert(atomic_terms(Exists('x', And('B', Odd('a')))) == { 'B', 'a' },
 69          "Error: atomic_terms is incorrect.");
 70   print('[ok] atomic_terms');
 71 };
 72 test_atomic_terms();
 73 
 74 test_func_to_expr := procedure() {
 75   assert(func_to_expr({ ['A', true], ['B', false], ['C', false] }) ==
 76          And("A", And(Not("B"), Not("C"))),
 77          "Error: func_to_expr is incorrect.");
 78   print('[ok] func_to_expr');
 79   print('[ok] iter_fct');
 80   print('[ok] foldr');
 81   print('[ok] fctize');
 82 };
 83 test_func_to_expr();
 84 
 85 test_atomic_formulas := procedure() {
 86   assert(atomic_formulas(Forall('x', Implies(Odd('a'),'B') )) ==
 87          { Odd('a'), 'B' },
 88          "Error: atomic_formulas is incorrect.");
 89   print('[ok] atomic_formulas');
 90   print('[ok] atomic_forms');
 91 };
 92 test_atomic_formulas();
 93 test_atomic_forms := procedure() { test_atomic_formulas(); };
 94 
 95 test_expr_to_events := procedure() {
 96   assert(expr_to_events(Not('A'), {'A'}) == { { ['A', false] } },
 97          "Error: expr_to_events ~A is incorrect.");
 98   assert(expr_to_events('A', {'A', 'B' }) ==
 99          { { ['A', true], ['B', false] },
100            { ['A', true], ['B', true ] } },
101          "Error: expr_to_events A is incorrect.");
102   assert(expr_to_events(Implies('A', 'B'), {'A', 'B' }) ==
103          { { ['A', true ], ['B', true ] },
104            { ['A', false], ['B', true ] },
105            { ['A', false], ['B', false] } },
106          "Error: expr_to_events A->B is incorrect.");
107   dom := { 'A', 'B', 'C' };
108   assert(expr_to_events(And('A', 'B'), dom) ==
109          { { ['A', true], ['B', true], ['C', true ] },
110            { ['A', true], ['B', true], ['C', false] }  },
111          "Error: expr_to_events A^B is incorrect.");
112   print('[ok] expr_to_events');
113   print('[  ] more expr_to_events');
114 };
115 test_expr_to_events();
116 
117 test_inv := procedure() {
118   assert(inv({ ['a', 1], ['b', 2] }) == { [1, 'a'], [2, 'b'] },
119          "Error: inv is incorrect.");
120   print('[ok] inv');
121 };
122 test_inv();
123 
124 test_arity := procedure() {
125   assert(arity('A') == 0 && arity(B()) == 0 && arity(Odd('x')) == 1 &&
126          arity(GreaterThan('x', 'y')) == 2 &&
127          arity(Relation('x', 'y', 'z')) == 3,
128          "Error: arity is incorrect.");
129   print('[ok] arity');
130 };
131 test_arity();
132 
133 test_flatten := procedure() {
134   assert([1,'a','C'] in flatten({1, 2, 3} >< {'a','b','c'} >< {'A', 'B', 'C'}),
135          "Error: flatten is incorrect.");
136   print('[ok] flatten');
137   print('[ok] flatten_list');
138 };
139 test_flatten();
140 
141 test_cart_prod := procedure() {
142   assert([2,'a','B'] in cart_prod([{1, 2, 3}, {'a','b','c'}, {'A', 'B', 'C'}]),
143          "Error: cart_prod is incorrect.");
144   print('[ok] cart_prod');
145 };
146 test_cart_prod();
147 
148 test_iter_prod := procedure() {
149   assert([3,2,1] in iter_prod({1,2,3}, 3),
150          "Error: is incorrect.");
151   print('[ok] iter_prod ');
152 };
153 test_iter_prod();
154 
155 test_iter_mult := procedure() {
156   assert(iter_mult({ 3, 5, 7 }) == 105,
157          "Error: iter_mult is incorrect.");
158   print('[ok] iter_mult');
159 };
160 test_iter_mult();
161 
162 test_choices := procedure() {
163   assert({3,'b','A'} in choices({ {1,2,3}, {'a','b','c'}, {'A','B','C'} }),
164          "Error: is incorrect.");
165   print('[ok] choices ');
166   print('[ok] set_choices ');
167   print('[ok] list_choices ');
168 };
169 test_choices();
170 
171 test_function_space := procedure() {
172   dom := { 'a', 'b' };
173   rng := { 1, 2 };
174   assert(function_space(dom, rng) == { {["a", 1], ["b", 1]},
175                                        {["a", 1], ["b", 2]},
176                                        {["a", 2], ["b", 1]},
177                                        {["a", 2], ["b", 2]}
178                                      },
179          'Error: function_space is incorrect');
180   print('[ok] function_space');
181 };
182 test_function_space();
183 
184 test_possible_states := procedure() {
185   s := { 'a', 'b' };
186   r := { ['a', {0,1}], ['b', {3,4}] };
187   assert(possible_states(s,r) == { {["a", 0], ["b", 3]},
188                                    {["a", 0], ["b", 4]},
189                                    {["a", 1], ["b", 3]},
190                                    {["a", 1], ["b", 4]}
191                                  },
192          'Error: possible_states is incorrect');
193   print('[ok] possible_states');
194 };
195 test_possible_states();
196 
197 test_partial_functions := procedure() {
198   f := { [1,2], [2,4] };
199   assert(partial_functions(f) == {{}, {[1, 2]}, {[1, 2], [2, 4]}, {[2, 4]}},
200          'Error: partial_functions is incorrect');
201 };
202 test_partial_functions();
203 
204 test_compatible_values := procedure() {
205   assert(compatible_values(
206            { ['A',true], ['B',false] },
207            { ['C',false], ['B',false] }
208          ) == true,
209          "Error: compatible_values is incorrect.");
210   assert(compatible_values(
211            { ['A',true], ['B',false] },
212            { ['C',false], ['B',true] }
213          ) == false,
214          "Error: compatible_values is incorrect.");
215   print('[ok] compatible_values');
216 };
217 test_compatible_values();
218 
219 test_compatible_complete_states := procedure() {
220   ps := { ['X', true] };
221   cs := possible_states( { 'X', 'Y' },
222                          { ['X', {true, false}],
223                            ['Y', {true, false}]  });
224   assert(compatible_complete_states(ps, cs) == {
225          {["X", true], ["Y", false]},
226          {["X", true], ["Y", true]}
227        }, 'Error: compatible_complete_states is incorrect');
228   print('[ok] compatible_complete_states');
229 };
230 test_compatible_complete_states();
231 
232 test_proc2func := procedure() {
233   assert(proc2func(procedure(x) { return 2 * x; }, { 1,2,3 }) ==
234          { [1,2], [2,4], [3,6] },
235          "Error: is incorrect.");
236   print('[ok] proc2func');
237 };
238 test_proc2func();
239 
240 test_parse_var := procedure() {
241   assert(parse_var('V_t4') == ['V', 4],
242          "Error: parse_var is incorrect.");
243   print('[ok] parse_var');
244 };
245 test_parse_var();
246 
247 test_orig_var := procedure() {
248   assert(orig_var('V_t3') == 'V',
249          "Error: orig_var is incorrect.");
250   print('[ok] orig_var ');
251 };
252 test_orig_var();
253 
254 test_var_time := procedure() {
255   assert(var_time('V_t2') == 2,
256          "Error: is incorrect.");
257   print('[ok] var_time');
258 };
259 test_var_time();
260 
261 test_assign_time := procedure() {
262   assert(assign_time({ ['A',    true], ['B',     false] }, 4) ==
263                      { ['A_t4', true], ['B_t4', false] },
264          "Error: assign_time is incorrect.");
265   print('[ok] assign_time');
266 };
267 test_assign_time();
268 
269 test_drop_time := procedure() {
270   assert(drop_time({ ['A_t3', true], ['B_t3', false] }) ==
271                    { ['A',     true], ['B',    false] },
272          "Error: drop_time is incorrect.");
273   print('[ok] drop_time');
274 };
275 test_drop_time();
276 
277 /*
278 test_is_logical_truth := procedure() {
279   assert(true,
280          "Error: is_logical_truth is incorrect.");
281   print('[  ] is_logical_truth');
282 };
283 test_is_logical_truth();
284 */
285 
286 test_parse_vars := procedure() {
287   assert(parse_vars(Forall('x', Exists('y', And(Odd('b'), LessThan('y','a')))))
288          == { 'x', 'y' },
289          "Error: is incorrect.");
290   print('[ok] parse_vars');
291 };
292 test_parse_vars();
293 
294 test_sub := procedure() {
295   assert(sub(Forall('x', Exists('y', And(Odd('b'), LessThan('y','x')))),'x','z')
296          ==  Forall('z', Exists('y', And(Odd('b'), LessThan('y','z')))),
297          "Error: sub is incorrect.");
298   print('[ok] sub');
299 };
300 test_sub();
301 
302 test_subformulas := procedure() {
303   assert(subformulas(Forall('x', Exists('y', And(Odd('b'), LessThan('y','x')))))
304          == { Odd('b'), LessThan('y','x'), And(Odd('b'), LessThan('y','x')),
305               Exists('y', And(Odd('b'), LessThan('y','x'))),
306               Forall('x', Exists('y', And(Odd('b'), LessThan('y','x'))))
307             },
308          "Error: subformulas is incorrect.");
309   print('[ok] subformulas');
310 };
311 test_subformulas();
312 
313 test_hausdorff_dist := procedure() {
314   assert(hausdorff_dist({1, 8, 17}, {4, 7, 15},
315                         procedure(x, y) { return (x - y) ** 2; }) == 9,
316          "Error: hausdorff_dist is incorrect.");
317   print('[ok] hausdorff_dist');
318 };
319 test_hausdorff_dist();
320 
321 test_kendalls_tau_dist := procedure() {
322   // In actual usage, instead of letters, states/events would be compared.
323   // A > B > C > D
324   u1 := { Succ('A', 'B'), Succ('A', 'C'), Succ('A', 'D'),
325           Succ('B', 'C'), Succ('B', 'D'),
326           Succ('C', 'D') };
327   // C > A > D > B
328   // Bubble sort: abcd => aCBd => CAbd => caDB 
329   u2 := { Succ('C', 'A'), Succ('C', 'D'), Succ('C', 'B'),
330           Succ('A', 'D'), Succ('A', 'B'),
331           Succ('D', 'B') };
332   assert(kendalls_tau_dist(u1, u2) == 1/2, // = 3 / (4 * (4-1) / 2)
333          "Error: is incorrect.");
334   print('[ok] kendalls_tau_dist');
335   print('[ok] states_u');
336 };
337 test_kendalls_tau_dist();
338 
339 test_atom_u := procedure() {
340   assert(atom_u({ SuccEq({ ['A', true ], ['B', false] },
341                          { ['A', false], ['B', false] }) }) == { 'A', 'B' },
342          "Error: atom_u is incorrect.");
343   print('[ok] atom_u');
344 };
345 test_atom_u();
346 
347 test_card_u_to_ord_u := procedure() {
348   // In actual usage, instead of letters, states/events would be compared.
349   assert(card_u_to_ord_u({ ['A', 5], ['B', 3], ['C', 3], ['D', 1] },
350                          { 'A', 'B', 'C', 'D' }) ==
351          { SuccEq('A', 'B'), SuccEq('A', 'C'), SuccEq('A', 'D'),
352            SuccEq('B', 'C'), SuccEq('B', 'D'),
353            SuccEq('C', 'B'), SuccEq('C', 'D'),
354            // reflexive
355            SuccEq('A', 'A'), SuccEq('B', 'B'),
356            SuccEq('C', 'C'), SuccEq('D', 'D')
357          },
358          "Error: card_u_to_ord_u is incorrect.");
359   print('[ok] card_u_to_ord_u');
360 };
361 test_card_u_to_ord_u();
362 
363 test_add_u_to_card_u := procedure() {
364   assert(add_u_to_card_u({ ['A', 3], ['B', 5], ['C', 9],
365                            [And('A','B'), -2] }) ==
366          { [{ ['A', false], ['B', false], ['C', false] }, 0],
367            [{ ['A', false], ['B', false], ['C', true ] }, 9],
368            [{ ['A', false], ['B', true ], ['C', false] }, 5],
369            [{ ['A', false], ['B', true ], ['C', true ] }, 14],
370            [{ ['A', true ], ['B', false], ['C', false] }, 3],
371            [{ ['A', true ], ['B', false], ['C', true ] }, 12],
372            [{ ['A', true ], ['B', true ], ['C', false] }, 6],
373            [{ ['A', true ], ['B', true ], ['C', true ] }, 15]
374          },
375          "Error: add_u_to_card_u is incorrect.");
376   print('[ok] add_u_to_card_u');
377 };
378 test_add_u_to_card_u();
379 
380 test_ord_u_dist := procedure() {
381   dom := exprs_to_state_space({ 'A', 'B' });
382   ldom := [ x : x in dom ];
383   u1 := card_u_to_ord_u({
384     [ldom[1], 10], [ldom[2], 7], [ldom[3], 3 ], [ldom[4], 0]
385   }, dom);
386   u2 := card_u_to_ord_u({
387     [ldom[1], 9], [ldom[2], 1], [ldom[3], 11 ], [ldom[4], 4]
388   }, dom);
389   assert(ord_u_dist(u1, u2) == 1/2,
390          "Error: ord_u_dist is incorrect.");
391   print('[ok] ord_u_dist');
392   print('[  ] ord_u_dist w/ sharpenings');
393 };
394 test_ord_u_dist();
395 
396 test_sort_list := procedure() {
397   less := procedure(x, y) { return x < y; };
398   greater := procedure(x, y) { return y < x; };
399 
400   l := [1,3,5,4,2];
401   s1 := sort_list(l, less );
402   s2 := sort_list(l, greater);
403   assert(s1 == [1,2,3,4,5], 'Error: sort_list is incorrect.');
404   assert(s2 == [5,4,3,2,1], 'Error: sort_list is incorrect.');
405   print('[ok] sort_list');
406 };
407 test_sort_list();
408 
409 test_sort_list := procedure() {
410   less := procedure(x, y) { return x < y; };
411   s := {1,3,5,4,2};
412   assert(sort_set(s, less) == [1,2,3,4,5], 'Error: sort_list is incorrect.');
413   print('[ok] sort_set');
414 };
415 test_sort_list();
416 
417 test_restrict_domain := procedure() {
418   assert(restrict_domain({ ['a', 1], ['b', 3], ['c', 5] }, {'b', 'c'})
419            == { ['b', 3], ['c', 5] },
420          "Error: is incorrect.");
421   print('[ok] restrict_doman');
422 };
423 test_restrict_domain();
424 
425 test_restricted_eq := procedure() {
426   assert(restricted_eq({['a', 1], ['b', 2], ['c', 3]},
427                                  {['b', 2], ['c', 3]} ) == true,
428          "Error: restricted_eq true is incorrect.");
429   print('[ok] restricted_eq true');
430   assert(restricted_eq({['a', 1], ['b', 2], ['c', 3]},
431                                  {['b', 2], ['c', 3], ['d', 4]} ) == false,
432          "Error: restricted_eq false is incorrect.");
433   print('[ok] restricted_eq false');
434 };
435 test_restricted_eq();
436 
437 test_range_for_arity := procedure() {
438   poss_sts := possible_states({ 'A', 'B' }, { ['A', tf()], ['B', tf()] });
439   sts := { first(poss_sts), { ['A', true], ['B', false]} };
440   assert(sts in range_for_arity(0, poss_sts),
441          "Error: range_for_arity 0 is incorrect.");
442   print('[ok] range_for_arity 0');
443   sts2 := { { ['A', false], ['B', false]}, { ['A', false], ['B', true]} };
444   // assert({ [sts], [sts2] } in range_for_arity(1, poss_sts), 
445   //       "Error: range_for_arity 1 is incorrect.");
446   // print('[ok] range_for_arity 1');
447   print('[..] range_for_arity 1');
448   // assert({ [sts, sts2] } in range_for_arity(2, poss_sts),
449   //       "Error: range_for_arity 2 is incorrect.");
450   // print('[ok] range_for_arity 2');
451 };
452 test_range_for_arity();
453 
454 test_exprs_to_state_space := procedure() {
455   assert(exprs_to_state_space({ "A", "B", Implies("A", "B"),
456                                 Not(And("A", Not("B"))) }) == {
457            {['A', false], ['B', false]}, {['A', false], ['B', true]},
458            {['A',  true], ['B', false]}, {['A',  true], ['B', true]}
459          }, "Error: exprs_to_state_space is incorrect.");
460   print('[ok] exprs_to_state_space');
461 };
462 test_exprs_to_state_space();
1 /* SetlX 2.2.0 Cheat Sheet 2 Here I've compiled a short list of examples illustrating the setlX language. 3 4 There is also a slightly expanded quick_ref.stlx version which demonstrates 5 nearly all the setlX needed in order to understand MetaEthical.AI. The main 6 difference is that this cheat sheet focuses on aspects that are more unique 7 to the setlX language and omits many commands that are more obvious or 8 self-explanatory, especially among programmers. 9 10 A much more thorough guide can be found in the official setlX tutorial.pdf. 11 */ 12 13 // Standard Arithmetic Operations 14 x := (1 + 2) * 3 / 4; 15 // ~< Result: 9/4 >~ 16 /* Note that lines always end with a semi-colon. 17 Assignment of a value to a variable is done with ":=", not "=". 18 Variable names start with a lowercase letter because capitalized names are 19 reserved for setlx functors. */ 20 21 // Exponents 22 2 ** 5; 23 // ~< Result: 32 >~ 24 25 // String Length 26 #"abc"; 27 // ~< Result: 3 >~ 28 29 30 // Sets 31 // ==== 32 a := { 1, 2, 3 }; 33 b := { 3, 4, 5 }; 34 35 // Cardinality 36 #b; 37 // ~< Result: 3 >~ 38 39 // Union 40 a + b; 41 // ~< Result: {1, 2, 3, 4, 5} >~ 42 43 // Intersection 44 a * b; 45 // ~< Result: {3} >~ 46 47 // Difference 48 a - b; 49 // ~< Result: {1, 2} >~ 50 51 // Check Membership 52 1 in a; 53 // ~< Result: true >~ 54 55 // Check if Subset 56 {1, 2} < {1, 2, 3}; 57 // ~< Result: true >~ 58 59 // Powerset 60 pow({ 1, 2 }); 61 // ~< Result: {{}, {1}, {1, 2}, {2}} >~ 62 2 ** { 1, 2 }; // alternate syntax 63 // ~< Result: {{}, {1}, {1, 2}, {2}} >~ 64 65 // Cartesian Product 66 pairs := { 'a', 'b' } >< { 1, 2 }; 67 // ~< Result: {["a", 1], ["a", 2], ["b", 1], ["b", 2]} >~ 68 69 // Set comprehension 70 { x * 2 : x in a }; 71 // ~< Result: {2, 4, 6} >~ 72 73 // Set comprehension with filter 74 { x : x in a | x < 3 }; 75 // ~< Result: {1, 2} >~ 76 77 // Set comprehension with list pattern matching 78 // Also, since we don't use the first element of the ordered pair, we can 79 // discard it with _ instead of giving it a variable name. 80 { y : [_, y] in pairs }; 81 // ~< Result: {1, 2} >~ 82 83 // Quantification 84 exists(x in a | x == 2); 85 // ~< Result: true >~ 86 forall(x in a | x < 5); 87 // ~< Result: true >~ 88 89 // Returns the first element according to the language's internal order. 90 // If the set is not all numbers or lists of them, the ordering should not be 91 // relied upon. 92 first(a); 93 // ~< Result: 1 >~ 94 95 // Sum of a Set of... 96 // Integers 97 +/ {1, 2, 3}; 98 // ~< Result: 6 >~ 99 100 // Strings 101 +/ {'a', 'b', 'c'}; 102 // ~< Result: "abc" >~ 103 104 // Sets 105 +/ { {1, 2}, {3, 4} }; 106 // ~< Result: {1, 2, 3, 4} >~ 107 108 // Lists 109 +/ { [1, 2], [3, 4] }; 110 // ~< Result: [1, 2, 3, 4] >~ 111 112 113 // Lists 114 // ===== 115 l := [ 1, 2, 3 ]; 116 m := [ 3, 4, 5 ]; 117 118 // Note that list indexes start at 1! 119 // Not 0 as with most other programming languages. 120 m[1]; 121 // ~< Result: 3 >~ 122 123 // Concatenation 124 l + m; 125 // ~< Result: [ 1, 2, 3, 3, 4, 5] >~ 126 // Note that unlike with sets, the 3 is duplicated 127 128 // List comprehension 129 [ x * 2 : x in l ]; 130 // ~< Result: [2, 4, 6] >~ 131 132 // List comprehension with filter 133 [ x : x in l | x < 3 ]; 134 // ~< Result: [1, 2] >~ 135 136 // Check Inclusion 137 1 in l; 138 // ~< Result: true >~ 139 140 // Quantification 141 exists(x in l | x == 2); 142 // ~< Result: true >~ 143 forall(x in l | x < 5); 144 // ~< Result: true >~ 145 146 // Range 147 m[1..2]; 148 // ~< Result: [3, 4] >~ 149 150 m[..2]; 151 // ~< Result: [3, 4] >~ 152 153 m[2..]; 154 // ~< Result: [4, 5] >~ 155 156 first([3,2,1]); 157 // ~< Result: 3 >~ 158 159 last([3,2,1]); 160 // ~< Result: 1 >~ 161 162 permutations([1, 2, 3]); 163 // ~< Result: {[1, 2, 3], [1, 3, 2], [2, 1, 3], [2, 3, 1], [3, 1, 2], [3, 2, 1]} >~ 164 165 // Sum of a List of: 166 // Integers 167 +/ [1, 2, 3]; 168 // ~< Result: 6 >~ 169 170 // Strings 171 +/ ['a', 'b', 'c']; 172 // ~< Result: "abc" >~ 173 174 // Sets 175 +/ [ {1, 2}, {3, 4} ]; 176 // ~< Result: {1, 2, 3, 4} >~ 177 178 // Lists 179 +/ [ [1, 2], [3, 4] ]; 180 // ~< Result: [1, 2, 3, 4] >~ 181 182 183 // Functions 184 // ========= 185 // They are just sets of ordered pairs where each element of the domain is 186 // mapped to one element of the range. 187 f := { [1,2], [2,4], [3,6] }; 188 189 // You can find the value given an input 190 f[2]; 191 // ~< Result: 4 >~ 192 193 // But if it's not a function (an input is mapped to more than one output) 194 // then it's om (which is like javascript's null). 195 g := { [1,2], [1,3] }; 196 g[1]; 197 // ~< Result: om >~ 198 199 domain(f); 200 // ~< Result: {1, 2, 3} >~ 201 202 range(f); 203 // ~< Result: {2, 4, 6} >~ 204 205 isMap({ [1,2], [2,4] }); 206 // ~< Result: true >~ 207 isMap({ [1,2], [1,4] }); 208 // ~< Result: false >~ 209 210 211 // Control Structures 212 // ================== 213 switch { 214 case false : 215 print('hi'); 216 case false : 217 print('hey'); 218 default : 219 print('bye'); 220 } 221 // bye 222 223 224 // Procedures 225 sum := procedure(x, y) { 226 return x + y; 227 }; 228 sum(3, 4); 229 // ~< Result: 7 >~ 230 231 // It can also be cached to avoid recomputing when given the same inputs 232 sum := cachedProcedure(x, y) { 233 return x + y; 234 }; 235 sum(3, 4); 236 // ~< Result: 7 >~ 237 238 239 // Object-oriented Programming 240 class point(x, y) { 241 x := x; 242 y := y; 243 244 // Instance procedure 245 dist_to := procedure(pt) { 246 return sqrt( (pt.x - x) ** 2 + (pt.y - get_y()) ** 2 ); 247 }; 248 249 // Not really necessary but just demonstrating calling another instance proc 250 get_y := procedure() { 251 // Also just demonstrating 'this' keyword, usually unneeded and implicit 252 return this.y; 253 }; 254 255 // Class procedures 256 static { 257 orig := procedure() { 258 return point(0, 0); 259 }; 260 } 261 } 262 pt := point(1, 1); 263 pt2 := point(4, 5); 264 pt.dist_to(pt2); 265 // ~< Result: 5.0 >~ 266 pt3 := point(3, 4); 267 pt3.dist_to(point.orig()); 268 // ~< Result: 5.0 >~ 269 270 271 /* SetlX Functors 272 ============== 273 SetlX allows for the creation of data structures composed from what they call 274 functors and arguments. Functor names always start with a capital letter and 275 then are followed by parentheses enclosing its arguments, which can be of 276 more or less any type. And that's basically all they are, like strings that 277 can take arguments to enable tree-like structures. 278 279 Here we demonstrate its usage to represent logical expressions where 280 arguments are either strings or other functor terms. (If you're familiar 281 with functors from category theory, the usage of that term here doesn't 282 really have anything to do with that notion.) 283 */ 284 And('A', 'B'); 285 // ~< Result: And("A", "B") >~ 286 287 Implies(And('A', 'B'), 'C'); 288 // ~< Result: Implies(And("A", "B"), "C") >~ 289 290 fct(And('A', 'B')); 291 // ~< Result: "And" >~ 292 293 args(And('A', 'B')); 294 // ~< Result: ["A", "B"] >~ 295 296 makeTerm('And', ['A', 'B']); 297 // ~< Result: And("A", "B") >~ 298 299 300 // Meta-programming 301 eval("{1} + {2, 3}"); 302 // ~< Result: {1, 2, 3} >~
1 /* SetlX 2.2.0 Quick Reference 2 This list of examples introduces much of the setlX language. 3 4 Nearly all of MetaEthical.AI is just a combination of the basic components 5 demonstrated here. There is also a slightly more compact cheat_sheet.stlx 6 that focuses on just the aspects that are more unique to setlX and omits 7 many commands that are more obvious or self-explanatory, especially among 8 programmers. 9 10 A much more thorough guide can be found in the official setlX tutorial.pdf. 11 */ 12 13 14 // This is a comment. 15 /* This is a 16 multi-line comment */ 17 18 19 // Standard Arithmetic Operations 20 x := (1 + 2) * 3 / 4; 21 // ~< Result: 9/4 >~ 22 /* Note that lines always end with a semi-colon. 23 Assignment of a value to a variable is done with ":=", not "=". 24 Variable names start with a lowercase letter because capitalized names are 25 reserved for setlx functors. */ 26 27 // Exponents 28 2 ** 5; 29 // ~< Result: 32 >~ 30 31 sqrt(81); 32 // ~< Result: 9.0 >~ 33 34 // Absolute value 35 abs(-7); 36 // ~< Result: 7 >~ 37 38 39 // Strings 40 // ======= 41 str1 := "abc"; 42 str2 := 'def'; // single quotes work too 43 44 // String concatenation 45 str3 := str1 + str2; 46 // ~< Result: "abcdef" >~ 47 48 // Length 49 #str3; 50 // ~< Result: 6 >~ 51 52 53 // Boolean values 54 // ============== 55 p := true; 56 q := false; 57 58 // Disjunction 59 p || q; 60 // ~< Result: true >~ 61 62 // Conjunction 63 p && q; 64 // ~< Result: false >~ 65 66 // Equality / Biconditional 67 p == q; 68 // ~< Result: false >~ 69 70 p != q; 71 // ~< Result: true >~ 72 73 // Negation 74 !p; 75 // ~< Result: false >~ 76 77 78 // Sets 79 // ==== 80 a := { 1, 2, 3 }; 81 b := { 3, 4, 5 }; 82 83 // Cardinality 84 #b; 85 // ~< Result: 3 >~ 86 87 // Union 88 a + b; 89 // ~< Result: {1, 2, 3, 4, 5} >~ 90 91 // Intersection 92 a * b; 93 // ~< Result: {3} >~ 94 95 // Difference 96 a - b; 97 // ~< Result: {1, 2} >~ 98 99 // Check Membership 100 1 in a; 101 // ~< Result: true >~ 102 103 // Check if Subset 104 {1, 2} < {1, 2, 3}; 105 // ~< Result: true >~ 106 {3} < {1, 2}; 107 // ~< Result: false >~ 108 109 // Powerset 110 pow({ 1, 2 }); 111 // ~< Result: {{}, {1}, {1, 2}, {2}} >~ 112 2 ** { 1, 2 }; // alternate syntax 113 // ~< Result: {{}, {1}, {1, 2}, {2}} >~ 114 115 // Cartesian Product 116 pairs := { 'a', 'b' } >< { 1, 2 }; 117 // ~< Result: {["a", 1], ["a", 2], ["b", 1], ["b", 2]} >~ 118 119 // Set comprehension 120 { x * 2 : x in a }; 121 // ~< Result: {2, 4, 6} >~ 122 123 // Set comprehension with filter 124 { x : x in a | x < 3 }; 125 // ~< Result: {1, 2} >~ 126 127 // Set comprehension with list pattern matching 128 // Also, since we don't use the first element of the ordered pair, we can 129 // discard it with _ instead of giving it a variable name. 130 { y : [_, y] in pairs }; 131 // ~< Result: {1, 2} >~ 132 133 // Quantification 134 exists(x in a | x == 2); 135 // ~< Result: true >~ 136 forall(x in a | x < 5); 137 // ~< Result: true >~ 138 139 // Returns the first element according to the language's internal order. 140 // If the set is not all numbers or lists of them, the ordering should not be 141 // relied upon. 142 first(a); 143 // ~< Result: 1 >~ 144 145 // Sum of a Set of... 146 // Integers 147 +/ {1, 2, 3}; 148 // ~< Result: 6 >~ 149 150 // Strings 151 +/ {'a', 'b', 'c'}; 152 // ~< Result: "abc" >~ 153 154 // Sets 155 +/ { {1, 2}, {3, 4} }; 156 // ~< Result: {1, 2, 3, 4} >~ 157 158 // Lists 159 +/ { [1, 2], [3, 4] }; 160 // ~< Result: [1, 2, 3, 4] >~ 161 162 163 // Lists 164 // ===== 165 l := [ 1, 2, 3 ]; 166 m := [ 3, 4, 5 ]; 167 168 // Note that list indexes start at 1! 169 // Not 0 as with most other programming languages. 170 m[1]; 171 // ~< Result: 3 >~ 172 173 // Concatenation 174 l + m; 175 // ~< Result: [ 1, 2, 3, 3, 4, 5] >~ 176 // Note that unlike with sets, the 3 is duplicated 177 178 // List comprehension 179 [ x * 2 : x in l ]; 180 // ~< Result: [2, 4, 6] >~ 181 182 // List comprehension with filter 183 [ x : x in l | x < 3 ]; 184 // ~< Result: [1, 2] >~ 185 186 // Check Inclusion 187 1 in l; 188 // ~< Result: true >~ 189 190 // Quantification 191 exists(x in l | x == 2); 192 // ~< Result: true >~ 193 forall(x in l | x < 5); 194 // ~< Result: true >~ 195 196 // Range 197 m[1..2]; 198 // ~< Result: [3, 4] >~ 199 200 m[..2]; 201 // ~< Result: [3, 4] >~ 202 203 m[2..]; 204 // ~< Result: [4, 5] >~ 205 206 first([3,2,1]); 207 // ~< Result: 3 >~ 208 209 last([3,2,1]); 210 // ~< Result: 1 >~ 211 212 permutations([1, 2, 3]); 213 // ~< Result: {[1, 2, 3], [1, 3, 2], [2, 1, 3], [2, 3, 1], [3, 1, 2], [3, 2, 1]} >~ 214 215 // Sum of a List of: 216 // Integers 217 +/ [1, 2, 3]; 218 // ~< Result: 6 >~ 219 220 // Strings 221 +/ ['a', 'b', 'c']; 222 // ~< Result: "abc" >~ 223 224 // Sets 225 +/ [ {1, 2}, {3, 4} ]; 226 // ~< Result: {1, 2, 3, 4} >~ 227 228 // Lists 229 +/ [ [1, 2], [3, 4] ]; 230 // ~< Result: [1, 2, 3, 4] >~ 231 232 233 // Functions 234 // ========= 235 // They are just sets of ordered pairs where each element of the domain is 236 // mapped to one element of the range. 237 f := { [1,2], [2,4], [3,6] }; 238 239 // You can find the value given an input 240 f[2]; 241 // ~< Result: 4 >~ 242 243 // But if it's not a function (an input is mapped to more than one output) 244 // then it's om (which is like javascript's null). 245 g := { [1,2], [1,3] }; 246 g[1]; 247 // ~< Result: om >~ 248 249 domain(f); 250 // ~< Result: {1, 2, 3} >~ 251 252 range(f); 253 // ~< Result: {2, 4, 6} >~ 254 255 256 // Control Structures 257 // ================== 258 if (true) { 259 print('hi'); 260 } 261 // hi 262 263 if (false) { 264 print('hi'); 265 } else { 266 print('bye'); 267 } 268 // bye 269 270 if (false) { 271 print('hi'); 272 } else if (true) { 273 print('hey'); 274 } else { 275 print('bye'); 276 } 277 // hey 278 279 switch { 280 case false : 281 print('hi'); 282 case false : 283 print('hey'); 284 default : 285 print('bye'); 286 } 287 // bye 288 289 s := {}; 290 for (i in [1, 2, 3]) { 291 s += { i*2 }; 292 } 293 s; 294 // ~< Result: {2, 4, 6} >~ 295 296 x := 0; 297 while (x < 5) { 298 x += 1; 299 } 300 x; 301 // ~< Result: 5 >~ 302 303 // Todo: match, _ 304 305 306 // Procedures 307 sum := procedure(x, y) { 308 return x + y; 309 }; 310 sum(3, 4); 311 // ~< Result: 7 >~ 312 313 // It can also be cached to avoid recomputing when given the same inputs 314 sum := cachedProcedure(x, y) { 315 return x + y; 316 }; 317 sum(3, 4); 318 // ~< Result: 7 >~ 319 320 321 // Object-oriented Programming 322 class point(x, y) { 323 x := x; 324 y := y; 325 326 // Instance procedure 327 dist_to := procedure(pt) { 328 return sqrt( (pt.x - x) ** 2 + (pt.y - get_y()) ** 2 ); 329 }; 330 331 // Not really necessary but just demonstrating calling another instance proc 332 get_y := procedure() { 333 // Also just demonstrating 'this' keyword, usually unneeded and implicit 334 return this.y; 335 }; 336 337 // Class procedures 338 static { 339 orig := procedure() { 340 return point(0, 0); 341 }; 342 } 343 } 344 pt := point(1, 1); 345 pt2 := point(4, 5); 346 pt.dist_to(pt2); 347 // ~< Result: 5.0 >~ 348 pt3 := point(3, 4); 349 pt3.dist_to(point.orig()); 350 // ~< Result: 5.0 >~ 351 352 353 /* SetlX Functors 354 ============== 355 SetlX allows for the creation of data structures composed from what they call 356 functors and arguments. Functor names always start with a capital letter and 357 then are followed by parentheses enclosing its arguments, which can be of 358 more or less any type. And that's basically all they are, like strings that 359 can take arguments to enable tree-like structures. 360 361 Here we demonstrate its usage to represent logical expressions where 362 arguments are either strings or other functor terms. (If you're familiar 363 with functors from category theory, the usage of that term here doesn't 364 really have anything to do with that notion.) 365 */ 366 And('A', 'B'); 367 // ~< Result: And("A", "B") >~ 368 369 Implies(And('A', 'B'), 'C'); 370 // ~< Result: Implies(And("A", "B"), "C") >~ 371 372 fct(And('A', 'B')); 373 // ~< Result: "And" >~ 374 375 args(And('A', 'B')); 376 // ~< Result: ["A", "B"] >~ 377 378 makeTerm('And', ['A', 'B']); 379 // ~< Result: And("A", "B") >~ 380 381 382 // Type checking 383 // ============= 384 isTerm(And('A', 'B')); 385 // ~< Result: true >~ 386 387 isString('abc'); 388 // ~< Result: true >~ 389 390 isSet({}); 391 // ~< Result: true >~ 392 393 isList([]); 394 // ~< Result: true >~ 395 396 isMap({ [1,2], [2,4] }); 397 // ~< Result: true >~ 398 isMap({ [1,2], [1,4] }); 399 // ~< Result: false >~ 400 401 402 int('5'); 403 // ~< Result: 5 >~ 404 405 str(5); 406 // ~< Result: "5" >~ 407 408 str({}); 409 // ~< Result: "{}" >~ 410 411 412 // Meta-programming 413 eval("{1} + {2, 3}"); 414 // ~< Result: {1, 2, 3} >~