Mathematica

(* ::Package:: *)

(************************************************************************)\012(* This file was generated automatically by the Mathematica front end.  *)\012(* It contains Initialization cells from a Notebook file, which         *)\012(* typically will have the same name as this file except ending in      *)\012(* ".nb" instead of ".m".                                               *)\012(*                                                                      *)\012(* This file is intended to be loaded into the Mathematica kernel using *)\012(* the package loading commands Get or Needs.  Doing so is equivalent   *)\012(* to using the Evaluate Initialization Cells menu command in the front *)\012(* end.                                                                 *)\012(*                                                                      *)\012(* DO NOT EDIT THIS FILE.  This entire file is regenerated              *)\012(* automatically each time the parent Notebook file is saved in the     *)\012(* Mathematica front end.  Any changes you make to this file will be    *)\012(* overwritten.                                                         *)\012(************************************************************************)\012\012

(*
Copyright (c) 2026 Peter Overmann
 
SPDX-License-Identifier: MIT
 
This file is part of the "Creating Intelligence" project. It is licensed 
under the MIT License. You may obtain a copy of the License in the LICENSE 
file in the root directory of this repository.
 
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR 
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
*)


(*
This is the native Mathematica reference implementation 
of topological associative memory.
See "Memory Native.nb" for a performance-optimized version.
*)


Memory[config_Association] :=

    Module[ {mem = <||>, NA, PA, NB, PB, T = 1, residual, matching = {},
                capacity, zero, store, retrieve, clear, memorycount},
            
        {NA, PA} = config["A_parameters"];
        {NB, PB} = config["B_parameters"];
        
        If[ ! MatchQ[{NA, PA, NB, PB}, {__Integer}], 
                    Message[Memory::params, config]];   

        (* Memory capacity, needed for default threshold calculation. *)
        capacity = If[NA == NB && PA == PB, 
                Log[2.0] NA (NA - 1) (NA - 2) / (PA (PA - 1) (PA - 2)),
                Log[2.0] NA (NA - 1) NB / (PA (PA - 1) PB)];
            
        (* Default pattern matching threshold. *)
        While[T < PA && capacity^2 * Binomial[PA, T] * Binomial[NA-PA, PA-T] / 
                Binomial[NA, PA] >= 1, T++];

        (* User-defined (scaled) threshold . *)
        If[KeyExistsQ[config, "threshold"] && NumberQ[config["threshold"]],
            T = Round[config["threshold"] * PA]];
            
        If[T < 2, T = 2];

        (* Optional residual: the number of non-matching query elements 
            that may be absorbed as noise, in multiples of PA. 
            NA means unlimited (the default). *)
        residual = If[KeyExistsQ[config, "residual"] && NumberQ[config["residual"]] 
                && config["residual"] >= 0, 
                Floor[config["residual"] * PA + 10^-9], NA];        
        
        zero = Developer`ToPackedArray[ConstantArray[0, NB]];

        (* Store hetero-association A -> B in memory. *)
        store[A_List] := store[A, A];
        
        (* Store A -> B in memory. *)
        store[A_List, B_List] := Module[{v, sub2},
            (* Step 1: Expansion coding. *)
            sub2 = Subsets[A, {2}];
        
            (* Step 2: Memory update. *)
            v = zero; v[[B]] = 1; 
            (mem[#] = BitOr[Lookup[mem, Key[#], zero], v])& /@ sub2; 
            ];

        (* Memory retrieval. *)
        retrieve[A_List] := 
            Module[ {X, P, Y, Tquery, R, Ri, sub2, v, t, w, ws, h, cutoff},
            (* Step 1: Initialize. *)
            matching = {};
            X = A; 
            Tquery = Max[T, Length[A] - residual];
            
            While[ True,
                (* Step 2: Threshold check.  *)
                If[ (P = Length[X]) < Tquery, Return[{}]]; 
            
                (* Step 3: Expansion coding. *)
                sub2 = Subsets[Range[P], {2}];
            
                (* Step 4: Aggregate. *)
                Ri = ConstantArray[0, {P, NB}];
                (v = Lookup[mem, Key[{X[[#1]], X[[#2]]}], zero]; 
                        Ri[[#1]] += v; Ri[[#2]] += v) & @@@ sub2;
                R = (Plus @@ Ri) / 2; 
            
                (* Step 5: Select (kWTA with threshold Tquery(Tquery-1)/2) *)
                (* The first element of a negative Ordering is the index 
                    of the k-th largest value. Use this for the threshold check. *)
                t = Max[1, R[[First @ Ordering[R, -PB]]]];
                
                (* Step 6: Threshold check. *)
                If[ t < Tquery (Tquery - 1) / 2, Return[{}]]; 
                Y = Sort[Pick[Range[NB], UnitStep[R - t], 1]];
            
                (* Step 7: Per-element weights. *)
                w = Total /@ Ri[[All, Y]];
                
                (* Step 8: Convergence test.  *)
                ws = Sort[w]; 
                h = P;
                While[ (cutoff = ws[[-h]]) < Length[Y] * (h - 1), --h];

                If[ h < Tquery, Return[{}]];        
                
                (* Step 9: Refine. *)
                X = Pick[X, # >= cutoff & /@ w];
                
                (* Handle edge case. *)
                If[h === Tquery &&  h < Length[X], X = {}; Return[{}]];
    
                (* Finished if X has converged. *)
                If[Length[X] == P, matching = X; Return[Y]];
                
                (* Step 10: Iterate.  *)
                ]
            ];
            
        clear := (mem = <||>;);
        memorycount := Total[Values[mem], 2];
                
        <|
            "A_parameters" -> {NA, PA},
            "B_parameters" -> {NB, PB},
            "T" -> T, (* Absolute pattern matching threshold. *)
            "residual" -> residual,
            "store" -> store,
            "retrieve" -> retrieve,
            "matching" :> matching,
            "clear" :> clear,
            "memorycount" :> memorycount,
            "backend" -> "Reference"
        |>
        ]



Memory::params  = "Invalid hyperparameters: `1`";