text:
# Copyright (c) 2018-2020, Carnegie Mellon University
# See LICENSE for details

Declare(ScatH, VScat_L);

Class(GathH, SumsBase, BaseMat, rec(
    rChildren := self >> [self.N, self.n, self.base, self.stride],
    rSetChild := rSetChildFields("N", "n", "base", "stride"),
    new := (self, N, n, base, stride) >>
	SPL(WithBases(self, rec(
      	N:= N, n := n, base := base, stride := stride, dimensions := [n, N]))),
    dims := self >> self.dimensions,
    sums := self >> self,
    area := self >> self.n,
    isReal := self >> true,
    transpose := self >> ScatH(self.N, self.n, self.base, self.stride),
    conjTranspose := self >> self.transpose(),
    inverse := self >> self.transpose(),
    toAMat := self >> let(
	   n := self.n,
       N := self.N,
       base := self.base,
       stride := self.stride,
       AMatMat(List([0..n-1], row -> let(
         idx := base+row,
		 BasisVec(N, idx))))
    ),
    isIdentity := self >> (self.n=self.N) and self.base=0 and self.stride = 1,
));

Class(ScatH, SumsBase, BaseMat, rec(
    rChildren := self >> [self.N, self.n, self.base, self.stride],
    rSetChild := rSetChildFields("N", "n", "base", "stride"),
    new := (self, N, n, base, stride) >>
SPL(WithBases(self, rec(
      	N:= N, n := n, base := base, stride := stride, dimensions := [N, n]))),
    dims := self >> self.dimensions,
    sums := self >> self,
    area := self >> self.n,
    isReal := self >> true,
    transpose := self >> GathH(self.N, self.n, self.base, self.stride),
    conjTranspose := self >> self.transpose(),
    inverse := self >> self.transpose(),
    toAMat := self >> let(
	   n := self.n,
       N := self.N,
       base := self.base,
       stride := self.stride,
       TransposedAMat(AMatMat(List([0..n-1], row -> let(
         idx := base+row*stride,
		 BasisVec(N, idx))))
    )),
    isIdentity := self >> (self.n=self.N) and self.base=0 and self.stride = 1,
));

Class(VGath_L, SumsBase, BaseMat, rec(
    rChildren := self >> [self.N, self.n, self.base],
    rSetChild := rSetChildFields("N", "n", "base"),
    new := (self, N, n, base) >>
	SPL(WithBases(self, rec(
      	N:= N, n := n, base := base, dimensions := [2, N]))),
    dims := self >> self.dimensions,
    sums := self >> self,
    area := self >> self.n,
    isReal := self >> true,
    transpose := self >> VScat_L(self.N, self.n, self.base),
    conjTranspose := self >> self.transpose(),
    inverse := self >> self.transpose(),
    toAMat := self >> let(
	   n := self.n,
       N := self.N,
       base := self.base,
       AMatSPL(L(2*n, 2) * Tensor(SPLAMat(AMatMat(List([0..n-1], row -> let(
         idx := base+row,
		 BasisVec(N, idx))))), I(2)))
    ),
    isIdentity := False
));

Class(VScat_L, SumsBase, BaseMat, rec(
    rChildren := self >> [self.N, self.n, self.base],
    rSetChild := rSetChildFields("N", "n", "base"),
    new := (self, N, n, base) >>
SPL(WithBases(self, rec(
      	N:= N, n := n, base := base, dimensions := [N, 2]))),
    dims := self >> self.dimensions,
    sums := self >> self,
    area := self >> self.n,
    isReal := self >> true,
    transpose := self >> VGath_L(self.N, self.n, self.base),
    conjTranspose := self >> self.transpose(),
    inverse := self >> self.transpose(),
    toAMat := self >> TransposedAMat(self.transpose().toAMat()),
    isIdentity := False
));

max_stars_repo_name: spiral-software/spiral-package-ffte

max_issues_repo_name: spiral-software/spiral-package-ffte

max_forks_repo_name: spiral-software/spiral-package-ffte


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

Class(ScratchSumsGen, LegacySumsGen);

ext: gi


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

Declare(BaseOverlap);

# BaseOverlap
# subclasses should define these fields:
#    _transpose_class
#
Class(BaseOverlap, BaseOperation, rec(

    new := (self, isize, overlap, spls) >>
       self._new(isize,overlap,spls),

    _new := meth(self, isize, overlap, spls)
        local result;
    if not IsList(spls) then spls := [spls]; fi;
    #Constraint(IsInt(overlap));
    #Constraint(IsInt(isize));
        DoForAll(spls, s -> Constraint(IsSPL(s)));
    result := WithBases(self, rec(isize := isize,
                            overlap := overlap,
                    _children := spls));
    result.dimensions := result.dims();
    return SPL(result);
    end,
    #-----------------------------------------------------------------------
    rChildren := self >> Concat([self.isize, self.overlap, self._children], When(IsBound(self._setDims), [self._setDims],[])),
    rSetChild := rSetChildFields("isize", "overlap", "_children", "_setDims"),

    #-----------------------------------------------------------------------
    isPermutation := self >> false,
    #-----------------------------------------------------------------------
    transpose := self >> CopyFields(self,
	self._transpose_class._new(
            self.isize,
            self.overlap,
            List(self._children, c -> c.transpose())
	)
    ),
    conjTranspose := self >> CopyFields(self,
	self._transpose_class._new(
            self.isize,
            self.overlap,
            List(self._children, c -> c.conjTranspose())
	)
    ),
    #----------------------------- Helpers ----------------------------------

    # Create SparseSPL list for AMatSPL conversion of type "rowDirectSum" and
    # "rowtensor" that represent a row overlapped direct sum of identity
    # matrices of sizes given in list N with overlap v.
    _rowoverlap := function(N, v)
        local i,j,k,l,min,L;
	L:=[];
	i:=1;
	j:=1;
	for k in [1..Length(N)] do
            for l in [1..N[k]] do
                Add(L, [i, j, 1]);
		i := i+1;
		j := j+1;
	    od;
            j := j-v;
        od;
        # Make sure there are no negative indices
	min := Minimum(List([1..Length(L)],i->L[i][2]));
	if min < 1 then
            for k in [1..Length(L)] do
                L[k][2] := L[k][2] + 1 - min;
            od;
	fi;
	return L;
   end,

   # Create SparseSPL list for AMatSPL conversion of type "colDirectSum" and
   # "coltensor" that represent a column overlapped direct sum of identity
   # matrices of sizes given in list N with overlap v.
   _coloverlap := function(N, v)
       local i,j,k,L;
       L := BaseOverlap._rowoverlap(N, v);
       for k in [1..Length(L)] do
           L[k]{[1..2]}:=Reversed(L[k]{[1..2]});
       od;
       return L;
   end,

   ovDim := (ov, dims) ->
       FoldL(dims,
         (p, dim) -> let(dd := p[1] + dim - ov,
                   [ dd, Min2(p[2], dd-dim+1), Max2(p[3], dd) ]),
         [ov, 1, 1]),

   ovDimList := (ov, dims) -> Drop(
       ScanL(dims,
         (p, dim) -> let(dd := p[1] + dim - ov,
                   [ dd, Min2(p[2], dd-dim+1), Max2(p[3], dd) ]),
         [ov, 1, 1]),
       1),

   spans := (ov, dims) -> let(
       dl := BaseOverlap.ovDimList(ov,dims),
       ofs := 1 - dl[1][2],
       List([1..Length(dims)], i -> [ofs+dl[i][1]-dims[i]+1, ofs+dl[i][1]]))
));


text:
MANIFEST0000640005642700116100000000014612754727571011213 0ustar  cgibsonenginstaller_version: 4
image_type: unittest
version: gftv200-40.1
platforms: [ GFLT110 ]
multiloader: 1
loader.gflt110.bin0000640005642700116100000000001312746742433013163 0ustar  cgibsonenggflt110.binkernel.img0000640005642700116100000000001212677613710012021 0ustar  cgibsonengkernel.img

max_stars_repo_name: DentonGentry/gfiber-platform

max_issues_repo_name: DentonGentry/gfiber-platform

max_forks_repo_name: DentonGentry/gfiber-platform


text:
# our deficient list
L :=
[ "ABCD", "CABD", "ACDB", "DACB", "BCDA",
  "ACBD", "ADCB", "CDAB", "DABC", "BCAD",
  "CADB", "CDBA", "CBAD", "ABDC", "ADBC",
  "BDCA", "DCBA", "BACD", "BADC", "BDAC",
  "CBDA", "DBCA", "DCAB" ];

# convert L to permutations on 1..4
u := List(L, s -> List([1..4], i -> Position("ABCD", s[i])));

# set difference (with all permutations)
v := Difference(PermutationsList([1..4]), u);

# convert back to letters
s := "ABCD";
List(v, p -> List(p, i -> s[i]));

ext: gap

max_issues_repo_name: seanwallawalla-forks/RosettaCodeData

max_forks_repo_name: seanwallawalla-forks/RosettaCodeData


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

# size = scratchpad segment size
# nsgmts = number of segments, defaults to 1.
Class(ALStore, AGenericTag, rec(
    updateParams := meth(self)
        Checked(Length(self.params)=3);
        self.size := self.params[1];
        self.nsgmts := self.params[2];
        self.linesize := self.params[3];
    end,
    isRegCx := false
));

Class(ALStoreCx, ALStore, rec(
    updateParams := meth(self)
        Checked(Length(self.params)=3);
        self.size := self.params[1]/2;
        self.nsgmts := self.params[2];
        self.linesize := self.params[3];
    end,
    isRegCx := true
));

# APad tag.
#
# APad is a special buffering tag meant for things like scratchpads, hence
# the name. It takes four parameters. In order, they are:
# b - the block size in number of elements
# s - the segment size, in number of blocks
# u - the number of segments
# n - a string identifier
#
# The block size is the SMALLEST allowed transfer size when copying data.
# The breakdown rules which propagate the APad tag insure that sub-block
# numbers of contiguous elements are never moved when the tag is present.
#
# A segment is a discrete memory separate from other segments, in the case
# that the 'u' parameter is >1. In the case of DPA, each segment is
# connected to a different compute processor, and we do parallelization
# after the APad tag is dropped.
#
# The string identifier is for the platform writer. You can label things
# like "local memory" or "vector register file."
#
# Also, the software pipelining loop (rather than a standard ISum) is
# automatically used with the APad tag.

Class(APad, AGenericTag, rec(
    applied := false,

    b := (self) >> self.params[1],
    s := (self) >> self.params[2],
    u := (self) >> self.params[3],
    bs := (self) >> self.params[1] * self.params[2],
    bsu := (self) >> Product(DropLast(self.params, 1)),
    n := (self) >> self.params[4],

    apply := self >> CopyFields(self, rec(applied := true)),

    updateParams := meth(self)
        Checked(Length(self.params)=4);
    end,

    print := meth(self)
        Print(self.name, "(", PrintCS([self.b(), self.s(), self.u()]),
            ", \"", self.n(), "\")");
    end
));


text:
# Built-in function
DivisorsInt(Factorial(5));
# [ 1, 2, 3, 4, 5, 6, 8, 10, 12, 15, 20, 24, 30, 40, 60, 120 ]

# A possible implementation, not suitable to large n
div := n -> Filtered([1 .. n], k -> n mod k = 0);

div(Factorial(5));
# [ 1, 2, 3, 4, 5, 6, 8, 10, 12, 15, 20, 24, 30, 40, 60, 120 ]

# Another implementation, usable for large n (if n can be factored quickly)
div2 := function(n)
  local f, p;
  f := Collected(FactorsInt(n));
  p := List(f, v -> List([0 .. v[2]], k -> v[1]^k));
  return SortedList(List(Cartesian(p), Product));
end;

div2(Factorial(5));
# [ 1, 2, 3, 4, 5, 6, 8, 10, 12, 15, 20, 24, 30, 40, 60, 120 ]


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

#F ParStream(p, u)
#F p=# of procs
#F u = minimum packet size (specified)
Class(ParStream, AGenericTag, rec(isParStream := true));

#F MBufCell_spec(u1, u2, p)
#F u1= PkSize for AxI.
#F u2=PkSize for L(IxA) and (IxA)L
#F  p=# of procs (Optional: default=1)
#F  r=# of iters to leave for the children (experimental). Optional.
Class(MBufCell_spec, AGenericTag, rec(isMBufCell_spec := true));

#F MBufCell_its(b)
#F Performs a multibuffer loop with exactly b multibuffered iterations
#F Not to be called directly
Class(MBufCell_its, AGenericTag, rec(isMBufCell_its := true));

#F ==========================================================================
#F MemCell(<its>) - Cell from main mem tag
Class(MemCell, AGenericTag, rec(isMemCell := true));

#F MBuf wrapper tags
Class(MBuf_maxWrapper_WHT, AGenericTag, rec(isMaxWrapper_WHT := true));
Class(MBuf_maxWrapper_DFT, AGenericTag, rec(isMaxWrapper_DFT := true));
Class(MBuf_maxWrapper_DFT_vecrecur, AGenericTag, rec(isMaxWrapper_DFT_vecrecur := true));

#-----------------------------------------------------------------------------
# Deprecated
#-----------------------------------------------------------------------------

#F MBufCell(<its>) - Cell buffer tag
Class(MBufCell, AGenericTag, rec(isMBufCell := true));
Class(MBufCell_max, AGenericTag, rec(isMBufCell := true));
#Class(MBufCell_mbuf, AGenericTag, rec(isMBufCell_mbuf := true));

#F ==========================================================================
#F 2DBufVecRecur(<its>) - Cell from main mem tag
Class(2DBufVecRecur, AGenericTag, rec(is2DBufVecRecur := true));


text:
#############################################################################
##
#W  rws.gi                  automgrp package                   Yevgen Muntyan
#W                                                             Dmytro Savchuk
##
#Y  Copyright (C) 2003 - 2018 Yevgen Muntyan, Dmytro Savchuk
##

DeclareGlobalFunction("__AG_FindRelators");

# It's not an IsRewritingSytem object
DeclareCategory("IsAGRewritingSystem", IsObject);

DeclareRepresentation("IsAGRewritingSystemRep",
                      IsComponentObjectRep,
                      ["rels",    # list of relators, they are words in the family's free group
                       "kb",      # Knuth-Bendix rewriting system
                       "fam",     # the automata family
                       "fpg_fam", # family of the elements from FP group <freegroup / rels>
                       "fpm_fam", # family of the elements of FP monoid obtained from the above group
                       "mhom",    # isomorphism from the FP group to the FP monoid
                      ]);

__AG_CreateRewritingSystem := function(fam, rels)
  local grp, fr_grp, fp_grp, fp_mon, fmon,
        ord, m_gens, rws, rws_data;

  rws_data := rec(fam := fam);

  if rels = fail then
    rels := __AG_FindRelators(fam, AG_Globals.max_rws_relator_len);
  fi;
  rws_data.rels := Difference(rels, [Word(One(fam))]);

  if IsAutomFamily(fam) then
    grp := GroupOfAutomFamily(fam);
  else
    grp := GroupOfSelfSimFamily(fam);
  fi;

  fr_grp := UnderlyingFreeGroup(grp);
  fp_grp := fr_grp / rels;
  rws_data.mhom := IsomorphismFpMonoid(fp_grp);
  fp_mon := Image(rws_data.mhom);
  fmon := FreeMonoidOfFpMonoid(fp_mon);

  rws_data.fpg_fam := FamilyObj(One(fp_grp));
  rws_data.fpm_fam := FamilyObj(One(fp_mon));

  m_gens := GeneratorsOfMonoid(fmon);
#  Print("mgens1=",m_gens,"\n");
#  m_gens := Permuted(m_gens, Product(List([1..Length(m_gens)/2], k->(2*k-1, 2*k))));
#  Print("mgens2=",m_gens,"\n");
  ord := ShortLexOrdering(FamilyObj(One(fmon)), m_gens);

  rws_data.kb := KnuthBendixRewritingSystem(fp_mon, ord);
  MakeConfluent(rws_data.kb);

  rws := Objectify(NewType(NewFamily("AGRewritingSystem"),
                           IsAGRewritingSystem and IsAGRewritingSystemRep),
                   rws_data);

  return rws;
end;

__AG_ReducedForm := function(rws, word)
  local reduced, word_mon;
  word_mon := ImageElm(rws!.mhom, ElementOfFpGroup(rws!.fpg_fam, word));
  reduced := ReducedForm(rws!.kb, UnderlyingElement(word_mon));
  word_mon := ElementOfFpMonoid(rws!.fpm_fam, reduced);
  return UnderlyingElement(PreImagesRepresentative(rws!.mhom, word_mon));
end;

InstallGlobalFunction(__AG_FindRelators,
function(fam, max_len)
  local fr_grp, w, elm, rels;

  if fam!.rws = fail then
    rels := [];
  else
    rels := List(fam!.rws!.rels);
  fi;

  fr_grp := UnderlyingFreeGroup(fam);

  # XXX FindRelations?
  for w in fr_grp do
    if Length(w) > max_len then
      break;
    fi;
    if IsAutomFamily(fam) then
      elm := Autom(w, fam);
    else
      elm := SelfSim(w, fam);
    fi;
    if IsOne(elm) then
      Add(rels, w);
    fi;
  od;

  return Difference(rels, [One(fr_grp)]);
end);

__AG_UpdateRewritingSystem := function(arg)
  local fam, max_len, new_rels;

  fam := arg[1];
  if Length(arg) > 1 then
    max_len := arg[2];
  else
    max_len := AG_Globals.max_rws_relator_len;
  fi;

  new_rels := __AG_FindRelators(fam, max_len);
  if fam!.rws = fail or new_rels <> fam!.rws!.rels then
    fam!.rws := __AG_CreateRewritingSystem(fam, new_rels);
  fi;

  return fail;
end;

__AG_UseRewritingSystem := function(fam, use)
  if fam!.use_rws <> use then
    if use then
      if fam!.rws = fail then
        fam!.rws := __AG_CreateRewritingSystem(fam, fail);
      fi;
    fi;
    fam!.use_rws := use;
  fi;
end;

__AG_RewritingSystem := function(fam)
  return fam!.rws;
end;

__AG_AddRelators := function(fam, rels)
  local old_rels;
  if fam!.rws = fail then
    old_rels := [];
  else
    old_rels := fam!.rws!.rels;
  fi;
  rels := Difference(Union(old_rels, rels), [One(UnderlyingFreeGroup(fam))]);
  if rels <> old_rels then
    fam!.rws := __AG_CreateRewritingSystem(fam, rels);
  fi;
end;

InstallMethod(AG_ReducedForm, [IsAGRewritingSystem, IsAssocWord],
function(rws, w)
  return __AG_ReducedForm(rws, w);
end);

InstallMethod(AG_ReducedForm, [IsAGRewritingSystem, IsList and IsAssocWordCollection],
function(rws, words)
  return Difference(List(words, w -> __AG_ReducedForm(rws, w)), [One(words[1])]);
end);

__AG_ReducedFormOfGroup := function(group, fam, construct)
  local gens, rws;

  if not fam!.use_rws then
    return group;
  fi;

  gens := GeneratorsOfGroup(group);
  gens := List(AG_ReducedForm(fam!.rws, List(gens, a -> Word(a))), w -> construct(w, fam));
  gens := Difference(gens, [One(group)]);
  if IsEmpty(gens) then
    return Group(One(group));
  else
    return GroupWithGenerators(gens);
  fi;
end;

InstallMethod(AG_ReducedForm, [IsAutomGroup],
function(group)
  return __AG_ReducedFormOfGroup(group, UnderlyingAutomFamily(group), Autom);
end);

InstallMethod(AG_ReducedForm, [IsSelfSimGroup],
function(group)
  return __AG_ReducedFormOfGroup(group, UnderlyingSelfSimFamily(group), SelfSim);
end);

__AG_ReducedFormOfElm := function(g, construct)
  local fam;

  fam := FamilyObj(g);

  if not fam!.use_rws then
    return g;
  else
    return construct(AG_ReducedForm(fam!.rws, Word(g)), fam);
  fi;
end;

InstallMethod(AG_ReducedForm, [IsAutom],
function(g)
  return __AG_ReducedFormOfElm(g, Autom);
end);

InstallMethod(AG_ReducedForm, [IsSelfSim],
function(g)
  return __AG_ReducedFormOfElm(g, SelfSim);
end);

InstallMethod(AG_UseRewritingSystem, [IsObject],
function(obj)
  AG_UseRewritingSystem(obj, true);
end);

InstallMethod(AG_UseRewritingSystem, [IsAutomFamily, IsBool],
function(fam, use)
  __AG_UseRewritingSystem(fam, use);
end);

InstallMethod(AG_UseRewritingSystem, [IsAutomGroup, IsBool],
function(grp, use)
  AG_UseRewritingSystem(UnderlyingAutomFamily(grp), use);
end);

InstallMethod(AG_UseRewritingSystem, [IsSelfSimFamily, IsBool],
function(fam, use)
  __AG_UseRewritingSystem(fam, use);
end);

InstallMethod(AG_UseRewritingSystem, [IsSelfSimGroup, IsBool],
function(grp, use)
  AG_UseRewritingSystem(UnderlyingSelfSimFamily(grp), use);
end);

InstallMethod(AG_UpdateRewritingSystem, [IsObject],
function(obj)
  AG_UpdateRewritingSystem(obj, AG_Globals.max_rws_relator_len);
end);

InstallMethod(AG_UpdateRewritingSystem, [IsAutomFamily, IsPosInt],
function(fam, max_len)
  __AG_UpdateRewritingSystem(fam, max_len);
end);

InstallMethod(AG_UpdateRewritingSystem, [IsAutomGroup, IsPosInt],
function(grp, max_len)
  AG_UpdateRewritingSystem(UnderlyingAutomFamily(grp), max_len);
end);

InstallMethod(AG_UpdateRewritingSystem, [IsSelfSimFamily, IsPosInt],
function(fam, max_len)
  __AG_UpdateRewritingSystem(fam, max_len);
end);

InstallMethod(AG_UpdateRewritingSystem, [IsSelfSimGroup, IsPosInt],
function(grp, max_len)
  AG_UpdateRewritingSystem(UnderlyingSelfSimFamily(grp), max_len);
end);

InstallMethod(AG_RewritingSystem, [IsAutomFamily],
function(fam)
  return __AG_RewritingSystem(fam);
end);

InstallMethod(AG_RewritingSystem, [IsAutomGroup],
function(grp)
  return AG_RewritingSystem(UnderlyingAutomFamily(grp));
end);

InstallMethod(AG_RewritingSystem, [IsSelfSimFamily],
function(fam)
  return __AG_RewritingSystem(fam);
end);

InstallMethod(AG_RewritingSystem, [IsSelfSimGroup],
function(grp)
  return AG_RewritingSystem(UnderlyingSelfSimFamily(grp));
end);

InstallMethod(AG_RewritingSystemRules, [IsObject],
function(obj)
  return Rules(AG_RewritingSystem(obj)!.kb);
end);

InstallMethod(AG_AddRelators, [IsAutomGroup, IsList and IsAssocWordCollection],
function(obj, rels) __AG_AddRelators(UnderlyingAutomFamily(obj), rels); end);

InstallMethod(AG_AddRelators, [IsAutomFamily, IsList and IsAssocWordCollection],
function(obj, rels) __AG_AddRelators(obj, rels); end);

InstallMethod(AG_AddRelators, [IsSelfSimGroup, IsList and IsAssocWordCollection],
function(obj, rels) __AG_AddRelators(UnderlyingSelfSimFamily(obj), rels); end);

InstallMethod(AG_AddRelators, [IsSelfSimFamily, IsList and IsAssocWordCollection],
function(obj, rels) __AG_AddRelators(obj, rels); end);

InstallMethod(AG_AddRelators, [IsObject, IsList and IsAutomCollection],
function(obj, list)
  AG_AddRelators(obj, List(list, g -> Word(g)));
end);

InstallMethod(AG_AddRelators, [IsObject, IsList and IsSelfSimCollection],
function(obj, list)
  AG_AddRelators(obj, List(list, g -> Word(g)));
end);

max_stars_repo_name: gap-packages/automgrp

max_issues_repo_name: gap-packages/automgrp

max_forks_repo_name: gap-packages/automgrp


text:
# SIMT-operators transfer their dim tag info to SIMT-SUM/ISum operators when SPL is lowered to SigmaSPL

Class(SIMTSumsGen, DefaultSumsGen, rec(

    _directSumSIMTSums := (self, obj, roverlap, coverlap, opts) >>
        let(cspans  :=  BaseOverlap.spans(roverlap, List(obj.children(), Cols)),
            rspans  :=  BaseOverlap.spans(coverlap, List(obj.children(), Rows)),
            nblocks  :=  obj.numChildren(),
            C  :=  Cols(obj),
            R  :=  Rows(obj),
            OP  :=  When(coverlap > 0, SUMAcc, SIMTSUM),
            OP( When(OP=SIMTSUM, [obj.simt_dim], []) :: Map([1..nblocks],
                 i -> Compose(
                      WID(R, 1 + rspans[i][2] - rspans[i][1], rspans[i][1]-1),
                      self(obj.child(i), opts),
                      RID(C, 1 + cspans[i][2] - cspans[i][1], cspans[i][1]-1) )))
            ),

    SIMTDirectSum := (self, o, opts) >> Cond(o.isPermutation(),
    Prm(ApplyFunc(fDirsum, List(o.children(), c->self(c, opts).func))),
        self._directSumSIMTSums(o, 0, 0, opts)
    ),

    IDirSum := (self, o, opts) >> let(
        bkcols := Cols(o.child(1)),
        bkrows := Rows(o.child(1)),
        nblocks := o.domain,
        cols := Cols(o),
        rows := Rows(o),

        SIMTISum(o.simt_dim, o.var, o.domain,
            Compose(Scat(fTensor(fBase(nblocks, o.var), fId(bkrows))),
                self(o.child(1), opts),
                Gath(fTensor(fBase(nblocks, o.var), fId(bkcols)))))
    ),

    SIMTTensor := meth(self, o, opts)
        local ch, col_prods, row_prods, col_prods_rev, row_prods_rev, i, j1, j2, j1b, j2b, bkcols, bkrows, prod, term;

        if o.isPermutation() then
            return Prm(ApplyFunc(o.fTensor, List(o.children(), c->self(c, opts).func)));
        elif ForAll(o.children(), x->ObjId(x)=Diag) then
            return Diag(ApplyFunc(diagTensor, List(o.children(), c->c.element)));
        fi;

        ch := o.children();
        col_prods := ScanL(ch, (x,y)->x*Cols(y), 1);
        col_prods_rev := Drop(ScanR(ch, (x,y)->x*Cols(y), 1), 1);

        row_prods := ScanL(ch, (x,y)->x*Rows(y), 1);
        row_prods_rev := Drop(ScanR(ch, (x,y)->x*Rows(y), 1), 1);

        prod := [];
        for i in [1..Length(ch)] do
            if not IsIdentitySPL(ch[i]) then
                bkcols := Cols(ch[i]);
                bkrows := Rows(ch[i]);
                j1 := Ind(col_prods[i]);
                j2 := Ind(row_prods_rev[i]);
                j1b := When(j1.range = 1, 0, fBase(j1));
                j2b := When(j2.range = 1, 0, fBase(j2));

                term := Scat(o.scatTensor(i)(Filtered([j1b, fId(bkrows), j2b], x->x<>0))) *
                        self(ch[i], opts) *
                    Gath(o.gathTensor(i)(Filtered([j1b, fId(bkcols), j2b], x->x<>0)));

#                if j2.range <> 1 then
                    term := SIMTISum(o.simt_dim, j2, j2.range, term);
#                fi;

#                if j1.range <> 1 then
                    term := SIMTISum(o.simt_dim, j1, j1.range, term);
#                fi;

                Add(prod, term);
            fi;
        od;
        return Compose(prod);
    end,

    SIMTISum := (self, o, opts) >> CopyFields(o, rec(_children := [self(o._children[1], opts)])),

    SIMTIDirSum := (self, o, opts) >> let(
        bkcols := Cols(o.child(1)),
        bkrows := Rows(o.child(1)),
        nblocks := o.domain,
        cols := Cols(o),
        rows := Rows(o),

        SIMTISum(o.simt_dim, o.var, o.domain,
            Compose(Scat(fTensor(fBase(nblocks, o.var), fId(bkrows))),
                self(o.child(1), opts),
                Gath(fTensor(fBase(nblocks, o.var), fId(bkcols)))))
    ),

    SIMTVStack := (self, o, opts) >> let(
        ch   := List(o.children(), e -> self(e, opts)),
        l    := Length(ch),
        ofs  := ScanL(ch, (p,x)->p+Rows(x), 0),
        SIMTSUM(o.simt_dim, List([1..l], i ->
            Scat(H(Rows(o), Rows(ch[i]), ofs[i], 1)) * ch[i] * Gath(fId(Cols(o)))))
    ),

    SIMTHStack1 := (self, o, opts) >> let(
        ch   := List(o.children(), e -> self(e, opts)),
        l    := Length(ch),
        it0  := Scat(fId(Rows(o))) * ch[1] * Gath(H(Cols(o), Cols(ch[1]), 0, 1)),
        ofs  := ScanL(ch, (p,x)->p+Cols(x), 0),
        its  := List([2..l], i ->
                   Scat(fId(Rows(o))) * ch[i] * Gath(H(Cols(o), Cols(ch[i]), ofs[i], 1))),
        SIMTSUM(o.simt_dim, [it0] :: its )
    ),

    SIMTIterVStack := (self, o, opts) >> let(
        bkcols := Cols(o.child(1)),
        bkrows := Rows(o.child(1)),
        nblocks := o.domain,
        cols := Cols(o), rows := Rows(o),
        SIMTISum(o.simt_dim, o.var, o.domain,
            Scat(fTensor(fBase(nblocks, o.var), fId(bkrows))) *
            self(o.child(1), opts) *
            Gath(fId(bkcols)))
    ),

    SIMTIterHStack1 := (self, o, opts) >> let(
        bkcols := Cols(o.child(1)),
        bkrows := Rows(o.child(1)),
        nblocks := o.domain,
        cols := Cols(o), rows := Rows(o),
        j := Ind(nblocks),
        ch := self(o.child(1), opts),

        SIMTISum(o.simt_dim, j, j.range,
        Scat(fId(bkrows)) *
        SubstVars(Copy(ch), rec((o.var.id) := j)) *
        Gath(fTensor(fBase(nblocks, j), fId(bkcols)))
        )
    ),

    ));

max_stars_repo_name: broderickpt/spiral-package-simt

max_issues_repo_name: broderickpt/spiral-package-simt

max_forks_repo_name: broderickpt/spiral-package-simt


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

Declare(ImageData, ImageVar, ImageUnk);

Class(ImageUnk, SumsBase, BaseContainer, rec(
    new := (self, dims) >> Checked(IsList(dims) and Length(dims) = 2 and ForAll(dims, IsPosInt),
        SPL(WithBases(self, rec(dimensions := dims, _children := [], dims := self >> self.dimensions)))),
    print := (self, i, is) >> Print("ImageUnk(", self.dims(), ")"),
    codeletShape := self >> let(d := self.dims(), Concat("ImageUnk_", StringInt(d[1]),"x",StringInt(d[2]))),
    fftImage := self >> ApplyFunc(O, TRDFT2D(self.dims()).dims()/self.dims()[2]),
    freqImage := self >> let(lindims := Product(TRDFT2D(self.dims()).dims())/(2*Product(self.dims())),
        fPrecompute(fConst(TComplex, lindims, 1+E(4)))),
    timeImage := self >> O(self.dims()[1], self.dims()[2]),
    randomImage := (self, m, n) >> ImageUnk(m, n),
    toAMat := self >> NullAMat(self.dimensions),
    unknownImage := self >> self
));

Class(ImageData, Mat, rec(
    print := (self, i, is) >> Print("ImageData(", Mat(self.element).dims(), ")"),
    codeletShape := self >> let(d := self.dims(), Concat("ImageData_", StringInt(d[1]),"x",StringInt(d[2]))),
    fftImage := self >> TRDFT2D(self.dims()).compute(self.element),
    freqImage := self >> let(fftRData := 1/Product(self.dims()) * Flat(self.fftImage()),
        fftCxData := List([1..Length(fftRData)/2], i->Complex(fftRData[2*i-1], fftRData[2*i])),
        FData(List(fftCxData, i->TComplex.value(i)))),
    timeImage := self >> self.element,
    randomImage := (self, m, n) >> let(data := List([1..m], i->List([1..n], j->FloatRat(Random([1..1000]/1000)))),
        ImageData(data)),
    unknownImage := self >> ImageUnk(self.dims())
));

Class(ImageVar, SumsBase, BaseContainer, rec(
    new := (self, dims) >> Checked(IsList(dims) and Length(dims) = 2 and ForAll(dims, IsPosInt),
        SPL(WithBases(self, rec(dimensions := dims, _children := [], dims := self >> self.dimensions)))),
    print := (self, i, is) >> Print("ImageVar(", self.dims(), ")"),
    codeletShape := self >> let(d := self.dims(), Concat("ImageVar_", StringInt(d[1]),"x",StringInt(d[2]))),
    fftImage := self >> ApplyFunc(O, TRDFT2D(self.dims()).dims()/self.dims()[2]),
    freqImage := self >> let(lindims := Product(TRDFT2D(self.dims()).dims())/(2*Product(self.dims())),
        type := TArray(TComplex, lindims), d := param(type, "freqImage"), FDataOfs(d, lindims, 0)),
    timeImage := self >> O(self.dims()[1], self.dims()[2]),
    randomImage := (self, m, n) >> ImageVar(m, n),
    toAMat := self >> NullAMat(self.dimensions),
    unknownImage := self >> ImageUnk(self.dims())
));

#F TRConv2D(img)
# For an image
# img := [ [ 1, 2, 3, 4 ],
#          [ 5, 6, 7, 8 ],
#          [ 9, 10, 11, 12 ],
#          [ 13, 14, 15, 16 ] ];
#
# MatSPL(TRConv2D(RealImage(img))) =
#
#[ [ 1, 4, 3, 2, 13, 16, 15, 14, 9, 12, 11, 10, 5, 8, 7, 6 ],
#  [ 2, 1, 4, 3, 14, 13, 16, 15, 10, 9, 12, 11, 6, 5, 8, 7 ],
#  [ 3, 2, 1, 4, 15, 14, 13, 16, 11, 10, 9, 12, 7, 6, 5, 8 ],
#  [ 4, 3, 2, 1, 16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5 ],
#  [ 5, 8, 7, 6, 1, 4, 3, 2, 13, 16, 15, 14, 9, 12, 11, 10 ],
#  [ 6, 5, 8, 7, 2, 1, 4, 3, 14, 13, 16, 15, 10, 9, 12, 11 ],
#  [ 7, 6, 5, 8, 3, 2, 1, 4, 15, 14, 13, 16, 11, 10, 9, 12 ],
#  [ 8, 7, 6, 5, 4, 3, 2, 1, 16, 15, 14, 13, 12, 11, 10, 9 ],
#  [ 9, 12, 11, 10, 5, 8, 7, 6, 1, 4, 3, 2, 13, 16, 15, 14 ],
#  [ 10, 9, 12, 11, 6, 5, 8, 7, 2, 1, 4, 3, 14, 13, 16, 15 ],
#  [ 11, 10, 9, 12, 7, 6, 5, 8, 3, 2, 1, 4, 15, 14, 13, 16 ],
#  [ 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 16, 15, 14, 13 ],
#  [ 13, 16, 15, 14, 9, 12, 11, 10, 5, 8, 7, 6, 1, 4, 3, 2 ],
#  [ 14, 13, 16, 15, 10, 9, 12, 11, 6, 5, 8, 7, 2, 1, 4, 3 ],
#  [ 15, 14, 13, 16, 11, 10, 9, 12, 7, 6, 5, 8, 3, 2, 1, 4 ],
#  [ 16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1 ] ]
#
Class(TRConv2D, TaggedNonTerminal, rec(
    abbrevs := [ img -> [img]],
    dims := self >> self.params[1].dims(),
    isReal := True,
    terminate := self >>
        TIRDFT2D(self.params[1].dims()).terminate() *
        RCDiag(FList(TReal, 1/Product(self.params[1].dims())*Flat(self.params[1].fftImage()))) *
        TRDFT2D(self.params[1].dims()).terminate(),
    normalizedArithCost := (self) >>
        TRDFT2D(self.params[1].dims()).normalizedArithCost() +
        TIRDFT2D(self.params[1].dims()).normalizedArithCost() +
        3 * Rows(TRDFT2D(self.params[1].dims())),
    symTerminate := self >> let(N := self.params[1].dims(),
        imgv := Flat(self.params[1].timeImage()),
        tf := MDDFT(N),
        ti := MDDFT(N, -1),
        d := MatSPL(tf) * imgv,
        t := ti * (1/Product(N))*Diag(d) * tf,
        MatSPL(t)),
    forwardTransform := self >> let(t := TRDFT2D(self.params[1].dims()), tags := self.getTags(),
        When(
            Length(tags) >= 1 and ObjId(tags[1]) = paradigms.vector.AVecReg, let(n := t.dims()[1], v := tags[1].v,
                TCompose([TPrm(fTensor(fId(n/(2*v)), L(2*v, 2))), t]).withTags(tags)),
            t.withTags(tags))),

    HashId := self >> let(conv := self.params[1].dims(), When(IsBound(self.tags), Concatenation(conv, self.tags), conv)),

    hashAs := meth(self)
        local conv;
        conv := Copy(self);
        conv.params[1] := self.params[1].unknownImage();
        return conv;
    end
));

#F 2D RConv Rule
NewRulesFor(TRConv2D, rec(
    TRConv2D_TRDFT2D_tSPL := rec(
        switch := true,
        applicable := (self, t) >> true,
        children := (self, t) >>
            [[
                TCompose([
                    TIRDFT2D(t.params[1].dims()),
                    TRCDiag(fPrecompute(t.params[1].freqImage())),
                    TRDFT2D(t.params[1].dims())
                ]).withTags(t.getTags())
            ]],
        apply := (self, t, C, Nonterms) >> C[1]
    )
));

#Class(TRCorr2D, TaggedNonTerminal, rec(
#    abbrevs := [ img -> [img]],
#    dims := self >> self.params[1].dims(),
#    isReal := True,
#    terminate := self >> let(tf := TRDFT2D(self.params[1].dims()), d := 1/Product(self.params[1].dims()) * Flat(tf.compute(self.params[1].element)),
#        TIRDFT2D(self.params[1].dims()).terminate() *
#        RCDiag(FList(TReal, List([1..Length(d)], j->When(IsOddInt(j), d[j], -d[j])))) *
#        tf.terminate(),
#    normalizedArithCost := (self) >>
#        TRDFT2D(self.params[1].dims()).normalizedArithCost() +
#        TIRDFT2D(self.params[1].dims()).normalizedArithCost() +
#        3 * Rows(TRDFT2D(self.params[1].dims()))),
#    symTerminate := self >> let(N := self.params[1].dims(),
#        img := self.params[1].element,
#        rimg := Reversed(List(img, Reversed)),
#        imgv := Flat(rimg),
#        tf := MDDFT(N),
#        ti := MDDFT(N, -1),
#        d := MatSPL(tf) * imgv,
#        t := ti * (1/Product(N))*Diag(d) * tf,
#        MatSPL(t))
#));


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

Declare(TTensorI_OL);
Declare(TTensorI_OL_Vectorize_AVecLast);

#Class(RewriteBarrier,RewritableObject);
#Class(RewriteBarrier,DFT);
#RewriteBarrier.numops:=self>>1;

DropVectorTag:=function(nt)
    return nt.withoutTag(AVecReg);
#D    return SetTag(nt,Filtered(GetTags(nt),t->ObjId(t)<>AVecReg));
end;

DropParTag:=function(nt)
   return nt.withoutTag(AParSMP);
#D    return SetTag(nt,Filtered(GetTags(nt),t->ObjId(t)<>AParSMP));
end;

Devectorize:=function(l,vlen)
  return List(l,e->
      Cond(
          ObjId(e)=TArray and e.t=TUnknown, TArray(e.t,e.size*vlen),
          ObjId(e)=TArray and ObjId(e.t)=TVect, TArray(e.t.t,e.size*e.t.size),
          ObjId(e)=TVect, TArray(e.t,e.size),
          e)
  );
end;

Class(VTensor_OL, Tensor, rec(
    new := (self, L) >> SPL(WithBases(self, rec(
        _children := [L[1]],
        vlen := L[2]))),
    print := (self,i,is) >> Print(self.name, "(",
        self.child(1).print(i+is,is), ", ", self.vlen,")"),
    sums := self >> Inherit(self, rec(_children := [self.child(1).sums()])),
    isPermutation := False,
    rng := meth(self)         #FULL HACK
       return Devectorize(self.child(1).rng(),self.vlen);
    end,
    dmn := meth(self)         #FULL HACK
       return Devectorize(self.child(1).dmn(),self.vlen);
    end,
    dims := meth(self)
       if (IsBound(self.rng)and IsBound(self.dmn)) then
          return [StripList(List(self.rng(),l->l.size)),
                  StripList(List(self.dmn(),l->l.size))];
       fi;
    end,
));

#VTensor.dmn:=  meth(self)         #FULL HACK
#  local x;
#  x:=self.child(1).dmn()[1];
#  if ObjId(x)=TArray and x.t=TUnknown then
#      return [TArray(TUnknown,x.size*self.vlen)];
#  else
#      return Devectorize(self.child(1).dmn());
#  fi;
#end;
#VTensor.rng:=  meth(self)         #FULL HACK
#  local x;
#  x:=self.child(1).rng()[1];
#  if ObjId(x)=TArray and x.t=TUnknown then
#      return [TArray(TUnknown,x.size*self.vlen)];
#  else
#      return Devectorize(self.child(1).rng());
#  fi;
#end;

BlockVPerm.dmn := meth(self)
       return [TArray(self.child(1).dmn()[1].t,self.child(1).dmn()[1].size*self.n)];
end;
BlockVPerm.rng := meth(self)
       return [TArray(self.child(1).rng()[1].t,self.child(1).rng()[1].size*self.n)];
end;

Class(ScatQuestionMark, Scat, rec());
Class(ICScatAcc,Scat,rec(codeletName:="ICSA"));

Class(VScatQuestionMark, VScat, rec());

Class(VScat_svQuestionMark, VScat_sv, rec());

Class(ScatInit, BaseMat, rec(
   new := meth(self, f,con,c)
        local res;
        res := SPL(WithBases(self, rec(func:=f,cond:=con, _children:=c)));
        return res;
   end,
   rng:=self>>self._children.rng(),
   dmn:=self>>self._children.dmn(),
   print := meth(self,i,is)
      Print(self.name, "(", self.func, ", ",self.cond,", ");
      Print("\n", Blanks(i+is));
      SPLOps.Print(self._children, i + is, is);
      Print("\n", Blanks(i),")");
      return;
   end,
   rChildren := self >> [self._children,self.func],
   rSetChild := meth ( self, n, newChild )
     if n= 1  then
         self._children :=newChild;
     elif n=2 then
         self.func := newChild;
else
        Error("<n> must be 1");
     fi;
   end

));

Class(ScatInitProbe,ScatInit,rec());
Class(ScatInitFixed,ScatInit,rec());

Class(KroneckerSymbol, BaseMat, rec(
   abbrevs := [ arg -> [Flat(arg)] ],
   new := meth(self, l)
        local res;
        res := SPL(WithBases(self, rec(element:=l)));
        return res;
   end,
   isExp:=true,
   dims:=self>>[1,1]  #avoid recursive definition for the old-school Tensor
));

## Base Vector
Class(BV, BaseMat, rec(
   new := (self, i) >> #Checked(IsVar(i),SPL(WithBases(self, rec(element:=i))))
       SPL(WithBases(self, rec(element:=i))),
   dims:=self>>self.element.dimensions #[0,0]  #avoid recursive definition for the old-school Tensor
));

## Multiplication operator
## Multiplication(1,n) is I(n)
## Multiplication(2,n) is a point-wise multiplication of two vectors of size n
 Class(Multiplication, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
        local res;
        if l[1]=1 then return Prm(fId(l[2])); fi;
        res := SPL(WithBases(self, rec(element:=l,TType:=Replicate(l[1],TUnknown))));
        return res;
    end,
    isPermutation := self >> false,
#    dmn:=self>>Replicate(self.element[1],TArray(self.TType,self.element[2])),
#HACK
#    dmn:=meth(self) local a; a:=Replicate(self.element[1],TArray(self.TType,self.element[2]));a[1]:=TArray(TReal,1);return a; end,
    dmn:=self >>List(self.TType,x->TArray(x,self.element[2])),
    rng:=self>>let(a:=Try(First(self.TType,x->ObjId(x)=TVect)),t:=Cond(a[1],a[2],self.TType[1]),[TArray(t,self.element[2])]),
    sums:= self>>self,
    numops:=self>>self.element[2]*(self.element[1]-1),
    transpose := self >>self,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ", ",self.element[2],")"); self.printA();
      return;
     end
 ));
 Class(ICMultiplication,Multiplication, rec(

 ));

  _mk_advdim_r := (e) -> List(e, x -> When(IsList(x), _mk_advdim_r(x), [x]));
  _mk_advdim := (d) -> _mk_advdim_r(When(IsList(d), d, [d]));

  Class(Glue, BaseMat, rec(
    abbrevs := [(n,size) ->[n,size]],
    new := meth(self,n,size)
        local res;
        res := SPL(WithBases(self, rec(element:=[n,size],dimensions:=[size*n,Replicate(n,size)],TType:=Replicate(n,TUnknown))));
        return res;
    end,
    isPermutation := self >> false,
    dmn:=self >>List([1..(self.element[1])],x->TArray(self.TType[x],self.element[2])),
    rng:=self>>[TArray(self.TType[1],self.element[1]*self.element[2])],
    sums:= self>>self,
    transpose:=self>>Copy(self),
    numops:=self>>self.element[2]*(self.element[1]-1),
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ", ",self.element[2],")"); self.printA();
      return;
     end,
    rChildren := self >> self.element,
    rSetChild := meth( self, n, newChild ) self.element[n] := Checked(n=1 or n=2, newChild); end,
    advdims := self >> let( d := self.dims(), [ _mk_advdim(d[1]), _mk_advdim(d[2]) ]),
    normalizedArithCost := (self) >> 0,
    isReal := self >> true, # makes no sense
 ));

  Class(Split, BaseMat, rec(
    abbrevs := [(size,n) ->[size,n]],
    new := meth(self,size,n)
        local res;
        res := SPL(WithBases(self, rec(element:=[size,n],dimensions:=[Replicate(n,size / n),size],TType:=[TUnknown])));
        return res;
    end,
    isPermutation := self >> false,
    rng:=self >>List([1..self.element[2]],x->TArray(self.TType[1],self.element[1] / self.element[2])),
    dmn:=self>>[TArray(self.TType[1],self.element[1])],
    sums:= self>>self,
    transpose:= self>>Copy(self),
    numops:=self>>(self.element[1]),
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ", ",self.element[2],")"); self.printA();
      return;
     end,
    rChildren := self >> self.element,
    rSetChild := meth( self, n, newChild ) self.element[n] := Checked(n=1 or n=2, newChild); end,
    advdims := self >> let( d := self.dims(), [ _mk_advdim(d[1]), _mk_advdim(d[2]) ]),
    normalizedArithCost := (self) >> 0,
    isReal := self >> true, # makes no sense
 ));

## NoOp, used with in COND_OL()
Class(NoOp, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(1, TUnknown))));
        fi;
    end,
    isPermutation := self >> false,
    dmn:=self >> [TArray(self.TType[1], self.element[1])],
    rng:=self >> [TArray(self.TType[1], self.element[1])],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

## Constant
Class(Const, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(1, TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> [],
# [TArray(self.TType[1], 0)],
    rng:=self >> [TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

Class(Or, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

Class(And, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

Class(NotEqual, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

Class(Equal, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

Class(ExclusiveOr, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

## Minimums
## Minimums(n) -> Input is n-way vector
Class(Minimums, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1)],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

## Maximums
## Maximums(n) -> Input is n-way vector
Class(Maximums, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
       local res;
        if (Length(l) = 1) then
          return SPL(WithBases(self, rec(element:=l, TType:=Replicate(l[1], TUnknown))));
        fi;
	# Here
    end,
    isPermutation := self >> false,
    dmn:=self >> List(self.TType, x->TArray(x,1)),
    rng:=self >> [ TArray(self.TType[1], 1) ],
    sums:= self>>self,
    numops:=self>> 0,
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ")"); self.printA();
      return;
    end
));

#Addition
Class(Addition, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
        local res;
        if l[1]=1 then return Prm(fId(l[2])); fi;
        res := SPL(WithBases(self, rec(element:=l,TType:=Replicate(l[1],TUnknown))));
        return res.setDims();
    end,
    isPermutation := self >> false,

    dmn:=self >> List( [1..self.element[1]], i -> TArray(TUnknown, self.element[2])),
    rng:=self >> [TArray(TUnknown, self.element[2])],
    sums:= self>>self,
    numops:=self>>self.element[2]*(self.element[1]-1),
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ", ",self.element[2],")"); self.printA();
      return;
     end,
    transpose := self >> InertTranspose(self),
    area := self >> self.element[2]*(self.element[1]-1),
));

Class(Subtraction, BaseMat, rec(
    abbrevs := [ arg -> [Flat(arg)] ],
    new := meth(self, l)
        local res;
        if l[1]=1 then return Prm(fId(l[2])); fi;
        res := SPL(WithBases(self, rec(element:=l,TType:=Replicate(l[1],TUnknown))));
        return res;
    end,
    isPermutation := self >> false,
    dmn:=self >>List(self.TType,x->TArray(x,self.element[2])),
    rng:=self>>let(a:=Try(First(self.TType,x->ObjId(x)=TVect)),t:=Cond(a[1],a[2],self.TType[1]),[TArray(t,self.element[2])]),
    sums:= self>>self,
    numops:=self>>self.element[2]*(self.element[1]-1),
    print := meth(self,i,is)
      Print(self.name, "(", self.element[1], ", ",self.element[2],")"); self.printA();
      return;
     end
));

# Declare(Cross);
# Class(Cross, BaseOperation, rec(
#    abbrevs := [ arg -> [Flat(arg)] ],

#    new := meth(self, L)
#         if Length(L)=1 then return L[1]; fi;
#         return SPL(WithBases(self, rec(_children:=L)));
#    end,

#    codeletName:="C",
#    isPermutation := self >> false,
#    dmn:=self>>let(li:=[],Flat(List(self._children, x-> x.dmn()))),
#    rng:=self>>let(li:=[],Flat(List(self._children, x-> x.rng()))),
#    sums:= meth(self)
#    local i;
#    for i in [1..Length(self._children)] do
#       self._children[i]:=self._children[i].sums();
#    od;
#    return self;
#    end,
#    hasLeftGath:=meth(self)
#       local i;
#       for i in self._children do
# #          if (ObjId(i)=Compose and ObjId(i.rChildren()[1])=Gath) or
# #              (ObjId(i)=Gath) then
# #              return true;
#            if (ObjId(i)=Compose and IsBound(i.rChildren()[1].func) and i.rChildren()[1].func.free()<>Set([])) or (IsBound(i.func) and (not ObjId(i)=VGath_dup) and i.func.free()<>Set([])) and not(ObjId(i)=ISum) then
#                return true;
#           fi;
#       od;
#       return false;
#    end,
#    hasRightScat:=meth(self)
#       local i;
#       for i in self._children do
# #          if (ObjId(i)=Compose and ObjId(Last(i.rChildren()))=Scat) or
# #              (ObjId(i)=Scat) then
# #              return true;
# #          fi;
#            if (ObjId(i)=Compose and ((IsBound(Last(i.rChildren()).func) and Last(i.rChildren()).func.free()=Set([])) or (ObjId(Last(i.rChildren()))=ISum) or (ObjId(Last(i.rChildren()))=VGath_dup and Length(Last(i.rChildren()).func.free())=1))) or (IsBound(i.func) and i.func.free()=Set([])) or (ObjId(i)=VGath_dup and Length(i.func.free())=1) or ObjId(i)=ISum  then
#                return true;
#           fi;
#       od;
#       return false;
#    end,
#    splitCross:=meth(self)
#       local l,r;
#       l:=Copy(self);
#       r:=[];
#       for i in [1..Length(l._children)] do
# #          if (ObjId(l._children[i])=Compose and ObjId(Last(l._children[i].rChildren()))=Scat) then
#            if (ObjId(l._children[i])=Compose and ((IsBound(Last(l._children[i].rChildren()).func) and Last(l._children[i].rChildren()).func.free()=Set([])) or (ObjId(Last(l._children[i].rChildren()))=ISum) or (ObjId(Last(l._children[i].rChildren()))=VGath_dup and Length(Last(l._children[i].rChildren()).func.free())=1))) or (ObjId(l._children[i])=ISum) then
#               Add(r,Last(l._children[i].rChildren()));
#               l._children[i]._children:=DropLast(l._children[i]._children,1);
# #          elif (ObjId(l._children[i])=Scat) then
#            elif (IsBound(l._children[i].func) and l._children[i].func.free()=Set([])) or (ObjId(i)=VGath_dup and Length(i.func.free())=1) then
#               Add(r,l._children[i]);
#               l._children[i]:=Prm(fId(l._children[i].dims()[1]));
#           else
#               Add(r,Prm(fId(l._children[i].dims()[2])));
#           fi;
#       od;
#       return [l,Cross(r)];
#    end
#  )
# );

Class(CrossBlockTop,Cross);

Identities:=function(l)
   return Cross(List(l,x->Prm(fId(x.size))));
end;

BB.numops:= self>>self.child(1).numops();
NoPull.numops:= self>>self.child(1).numops();
NoPullRight.numops:= self>>self.child(1).numops();
NoPullLeft.numops:= self>>self.child(1).numops();

PushR.numops:= self>>self.child(1).numops();

ISum.numops:=self>>self.child(1).numops()*self.domain;

SUM.numops:=self>>Sum(List(self._children,c->c.numops()));

Cross.numops:=self>>Sum(List(self._children,c->c.numops()));

Compose.numops:=self>>Sum(List(self._children,c->c.numops()));

Scat.numops:=self>>self.dims()[2];

Gath.numops:=self>>self.dims()[1];

VGath.numops:=self>>self.dims()[1];
VGath_dup.numops:=self>>self.dims()[1];
VReplicate.numops:=self>>self.v;
VHAdd.numops:=self>>self.v*self.v;
VPerm.numops:=self>>self.dimensions[1];
Tensor.numops:=self>>self.dimensions[1];

ScatAcc.numops:=self>>self.dimensions[2]*2; #an add and a store
VScatAcc.numops:=self>>self.dimensions[2]*2; #an add and a store
VScat.numops:=self>>self.dimensions[2];

fBase.numops:=self>>0;
fTensor.numops:=self>>0;
fId.numops:=self>>0;
ScatInit.numops:=self>>self._children.numops();

Prm.numops:=self>>0;
I.numops:=self>>0;

VTensor_OL.numops:=self>>self._children[1].numops()*self.vlen;

Class(AOne,AGenericTag);
Class(AMul,AGenericTag, rec(
    __call__ := meth ( arg )
      local  result, self, params;
      self := arg[1];
      params := arg{[ 2 .. Length(arg) ]};
      result := WithBases(self, rec(
              params := params,
              operations := PrintOps));
      return result;
    end,
    print := self >> Print(self.__name__, "(", PrintCS(self.params), ")")
   ));

Class(VOLWrap, VWrapBase, rec(
    __call__ := (self,isa) >> Checked(IsSIMD_ISA(isa),
        WithBases(self, rec(operations:=PrintOps, isa:=isa))),

    wrap := (self,r,t) >> let(isa := self.isa, v := isa.v,
#This is OBVIOUSLY a hack. only deals with some kind of vectorization
            nontransforms.ol.TTensorI_OL_Vectorize_AVecLast(TTensorI_OL(t, [ AOne, AVec ],[ [ 0, 2 ] ], [ 1, v ], [ AVecReg(isa) ]), r)),

    twrap := (self, t) >> let(isa := self.isa, v := isa.v,
#This is OBVIOUSLY a hack. only deals with some kind of vectorization
           TTensorI_OL(t, [ AOne, AVec ],[ [ 0, 2 ] ], [ 1, v ], [ AVecReg(isa) ])
        ),

    print := self >> Print(self.name, "(", self.isa, ")")
));

Class(TTensorI_OL, Tagged_tSPL, rec(
    abbrevs := [
        (nt,g,s,v)  -> [nt,List([1..Length(g)],i->When(v[i]=1,AOne,g[i])),s,v]
    ],

    dmn := self >> let(
        nt := self.params[1],
        sizes := self.params[4],
        List([1..Length(sizes)], i->
            TArray(nt.dmn()[i].t,sizes[i] * nt.dmn()[i].size)
        )),

    rng := self >> let(
        nt := self.params[1],
        s := self.params[3],
        v := self.params[4],
        List([1..Length(s)], i->
            TArray(nt.rng()[i].t, Product(List(s[i], a ->
                When(a=0,
                    nt.rng()[i].size,
                    v[a]
                )
            )))
        )
    ),

    isReal := self >> self.params[1].isReal(),
    doNotMeasure:=true,
    doNotSaveInHashtable:=true,
    decomposePerformance:=true,
    transpose := self >> Copy(self),
#D    tagpos :=5,
));

#Changed by Marek, looks OK
NewRulesFor(TTensorI_OL, rec(
    TTensorI_OL_Base :=rec(
        applicable := (self,t) >> not(t.hasTag(AParSMP)) and t.getTags() = t.params[1].getTags(),
        freedoms := nt -> let(
            g := nt.params[2],
            nbvars := Length(Filtered(g,x->x=APar or x=AVec)),
            [Arrangements([1..nbvars],nbvars)]
        ),
        child := (nt,freedoms) -> [nt.params[1],InfoNt(freedoms)],
        recompose := (nt,cnt,cperf) -> cperf[1]*Product(List(Filtered(Zip2(nt.params[2],nt.params[4]),l->l[1] in [APar,AVec]),x->x[2])),
        apply := function(nt,c,cnt)
            local g,s,v,perm,ind1,ind,ind2,gathers,scatters,result,z;
            g := nt.params[2];
            s := nt.params[3];
            v := nt.params[4];
            perm := cnt[2].params[1][1];
            ind1 := List([1..Length(g)], i ->
                Cond(g[i]=APar or g[i]=AVec, fBase(Ind(v[i])),
                    g[i]=AOne, fId(1),
                    ObjId(g[i])=AMul, g[i],
                    Error("PV not known!")
                )
            );
            ind := List(ind1, x -> When(ObjId(x)=AMul,ind1[x.params[1]],x));
            gathers := Cross(List([1..Length(g)], i -> let(
                kernelsize := fId(cnt[1].dmn()[i].size),
                When(g[i]=APar or (ObjId(g[i])=AMul and g[i].params[2]=APar),
                    Gath(fTensor(ind[i],kernelsize)),
                    Gath(fTensor(kernelsize,ind[i]))
                )
            )));
            scatters := Cross(List([1..Length(s)], i ->
                ScatQuestionMark( fTensor(
                    List(s[i], y ->
                        When(y=0,
                            fId(cnt[1].rng()[i].size),
                            ind[y]
                        )
                    )
                ))
            ));
            result := scatters*c[1]*gathers;
            ind2:=[];
            for i in [1..Length(g)] do
                if g[i]=APar or g[i]=AVec then
                    Add(ind2,ind[i]);
                fi;
            od;
            for i in [1..Length(ind2)] do
                result:=ISum(ind2[perm[i]].params[2], ind2[perm[i]].params[1], result);
            od;

            return result;
        end
    ),
#D        applicable :=(self,t) >> Length(Filtered(GetTags(t),x->ObjId(x)=AParSMP))=0 and GetTags(t)=GetTags(t.params[1])
#D,
#D            freedoms := nt -> let(g:=nt.params[2],
#D                nbvars:=Length(Filtered(g,x->x=APar or x=AVec)),
#D                [Arrangements([1..nbvars],nbvars)]),
#D            child := (nt,freedoms) -> [nt.params[1],InfoNt(freedoms)],
#D            recompose := (nt,cnt,cperf) -> cperf[1]*Product(List(Filtered(Zip2(nt.params[2],nt.params[4]),l->l[1] in [APar,AVec]),x->x[2])),
#D            apply := function(nt,c,cnt)
#D                local g,s,v,perm,ind1,ind,ind2,gathers,scatters,result,z;
#D                g:=nt.params[2];
#D                s:=nt.params[3];
#D                v:=nt.params[4];
#D                perm:=cnt[2].params[1][1];
#D                ind1:=List([1..Length(g)],
#D                    i->Cond(g[i]=APar or g[i]=AVec,fBase(Ind(v[i])),
#D                    g[i]=AOne,fId(1),
#D                    ObjId(g[i])=AMul,g[i],
#D                    Error("PV not known!")));
#D                ind:=List(ind1,x->When(ObjId(x)=AMul,ind1[x.params[1]],x));
#D                gathers:= Cross(List([1..Length(g)],
#D                        i->let(kernelsize:=fId(cnt[1].dmn()[i].size),
#D                            When(g[i]=APar or (ObjId(g[i])=AMul and g[i].params[2]=APar),
#D                                Gath(fTensor(ind[i],kernelsize)),
#D                                Gath(fTensor(kernelsize,ind[i]))))));
#D                scatters:=Cross(List([1..Length(s)],
#D                            i->ScatQuestionMark(fTensor(
#D                                    List(s[i],y->
#D                                        When(y=0,fId(cnt[1].rng()[i].size),
#D                                            ind[y]))))));
#D                result:=scatters*c[1]*gathers;
#D                ind2:=[];
#D                for i in [1..Length(g)] do
#D                    if g[i]=APar or g[i]=AVec then
#D                        Add(ind2,ind[i]);
#D                    fi;
#D                od;
#D                for i in [1..Length(ind2)] do
#D                        result:=ISum(ind2[perm[i]].params[2],ind2[perm[i]].params[1],result);
#D                od;
#D
#D                return result;
#D            end
#D            ),

#Changed by Marek, looks OK
    TTensorI_OL_Parrallelize_AParFirst := rec(
        applicable :=(self,t) >>
            t.isTag(1, AParSMP)
            and 0 <> t.params[3][1][1]
            and APar = t.params[2][t.params[3][1][1]]
            and t.firstTag().params[1] = t.params[4][t.params[3][1][1]]
            and 1 = t.params[3][1][1], #t.params[3][1][1]=1 is a trick to prevent // the 2nd input because TensorGeneral breaks at the moment

        freedoms := nt -> [[1]],

        child := function (nt,freedoms)
            local PV,sizes;
            sizes := Copy(nt.params[4]);
            sizes[nt.params[3][1][1]] := 1;
            if Length(Filtered(sizes, t -> t<>1)) > 0 then
                PV:=Copy(nt.params[2]);
                PV[nt.params[3][1][1]] := AOne;
                return [TTensorI_OL(
                    DropParTag(nt.params[1]),
                    PV,
                    nt.params[3],
                    sizes,
                    Drop(nt.params[5],1)
                )];
            else
                return [ DropParTag(Copy(nt.params[1])) ];
            fi;
        end,

        apply := function(nt,c,cnt)
            local myCross,a,b,index;

            index := Ind(nt.firstTag().params[1]);

            a := List(c[1].dims()[2], x -> fId(x));
            a[nt.params[3][1][1]] := fTensor(fBase(index), a[nt.params[3][1][1]]);

            b:=List(a, x -> Gath(x));

            myCross:=Cross(b);

#              return SMPSum(GetFirstTag(nt).params[1],index, GetFirstTag(nt).params[1],ScatQuestionMark(fTensor(fBase(index), fId(c[1].dims()[1]))) *c[1]* Cross(Gath(fTensor(fBase(index), fId(c[1].dims()[2][1]))),Gath(fId(c[1].dims()[2][2]))));
            return SMPSum(
                nt.firstTag().params[1],
                index,
                nt.firstTag().params[1],
                ScatQuestionMark(fTensor(
                    fBase(index),
                    fId(c[1].dims()[1])
                ))
                * c[1] * myCross
            );
        end
    ),
#D            applicable :=(self,t) >> FirstTagEq(t, AParSMP) and t.params[3][1][1]<>0 and
#D               t.params[2][t.params[3][1][1]]=APar and
#D               t.params[4][t.params[3][1][1]]=GetFirstTag(t).params[1] and
#D               t.params[3][1][1]=1, #t.params[3][1][1]=1 is a trick to prevent // the 2nd input because TensorGeneral breaks at the moment
#D            freedoms := nt -> [[1]],
#D            child := function (nt,freedoms)
#D              local PV,sizes;
#D              sizes:=Copy(nt.params[4]);
#D              sizes[nt.params[3][1][1]]:=1;
#D              if Length(Filtered(sizes, t->t<>1))>0 then
#D                  PV:=Copy(nt.params[2]);
#D                  PV[nt.params[3][1][1]]:=AOne;
#D              return [TTensorI_OL(DropParTag(nt.params[1]),
#D                          PV,nt.params[3],sizes,Drop(nt.params[5],1))];
#D              else
#D                  return [DropParTag(Copy(nt.params[1]))];
#D              fi;
#D            end,
#D            apply := function(nt,c,cnt)
#D              local myCross,a,b,index;
#D              index:=Ind(GetFirstTag(nt).params[1]);
#D              a:=List(c[1].dims()[2],x->fId(x));
#D              a[nt.params[3][1][1]]:=fTensor(fBase(index),a[nt.params[3][1][1]]);
#D              b:=List(a,x->Gath(x));
#D              myCross:=Cross(b);
#D#              return SMPSum(GetFirstTag(nt).params[1],index, GetFirstTag(nt).params[1],ScatQuestionMark(fTensor(fBase(index), fId(c[1].dims()[1]))) *c[1]* Cross(Gath(fTensor(fBase(index), fId(c[1].dims()[2][1]))),Gath(fId(c[1].dims()[2][2]))));
#D              return SMPSum(GetFirstTag(nt).params[1],index, GetFirstTag(nt).params[1],ScatQuestionMark(fTensor(fBase(index), fId(c[1].dims()[1]))) *c[1]* myCross);
#D            end),

    #That code vectorizes the last guy if it is a AVec of the vector size
    #Hack, only works with one output
#Changed by Marek, looks OK
    TTensorI_OL_Vectorize_AVecLast :=rec(
        applicable := (self,t) >>
            t.isTag(1, AVecReg)
            and Last(t.params[3][1])<>0
            and AVec = t.params[2][Last(t.params[3][1])]
            and t.params[4][Last(t.params[3][1])] = t.firstTag().v,

        freedoms := nt -> [[1]],

        child := function(nt, freedoms)
            local PV, sizes;

            sizes := Copy(nt.params[4]);
            sizes[Last(nt.params[3][1])] := 1;

            if Length(Filtered(sizes, t->t<>1))>0 then
                PV:=Copy(nt.params[2]);
                PV[Last(nt.params[3][1])]:=AOne;
                return [TTensorI_OL(
                    DropVectorTag(nt.params[1]).setWrap(VOLWrap(nt.firstTag().isa)),
                    PV,
                    nt.params[3],
                    sizes,
                    Drop(nt.params[5],1)
                )];
            else
                return [ DropVectorTag(Copy(nt.params[1])).setWrap(VOLWrap(nt.firstTag().isa)) ];
          fi;
        end,

        apply := function(nt,c,cnt)
            local myCross, myScat, mydims, v;

            v := nt.firstTag().v;

            #little hack for KernelDup
            if ObjId(nt.params[1]).name <> "KernelMMMDuped" then
                myCross := Cross(List(nt.dims()[2], t ->
                    VGath_dup(fId(t), v)
                ));
                myCross._children[Last(nt.params[3][1])] :=
                    VPrm_x_I(
                        fId( nt.dims()[2][Last(nt.params[3][1])] / v ),
                        nt.firstTag().v
                    );
            else
                mydims := nt.dims()[2];
                mydims[1] := mydims[1] * v;
                myCross := Cross(List(mydims, t ->
                    VPrm_x_I( fId(t/v), v )
                ));
            fi;

            myScat := VScat(fId( nt.dims()[1]/v ), v );

            return myScat * VTensor_OL(c[1], v) * myCross;
        end
    ),
#D            applicable :=(self,t) >> FirstTagEq(t, AVecReg) and Last(t.params[3][1])<>0 and t.params[2][Last(t.params[3][1])]=AVec and t.params[4][Last(t.params[3][1])]=GetFirstTag(t).v,
#D            freedoms := nt -> [[1]],
#D            child := function (nt,freedoms)
#D              local PV,sizes;
#D              sizes:=Copy(nt.params[4]);
#D              sizes[Last(nt.params[3][1])]:=1;
#D              if Length(Filtered(sizes, t->t<>1))>0 then
#D                  PV:=Copy(nt.params[2]);
#D                  PV[Last(nt.params[3][1])]:=AOne;
#D              return [TTensorI_OL(DropVectorTag(nt.params[1]).setWrap(VOLWrap(GetFirstTag(nt).isa)),
#D                          PV,nt.params[3],sizes,Drop(nt.params[5],1))];
#D              else
#D                  return [DropVectorTag(Copy(nt.params[1])).setWrap(VOLWrap(GetFirstTag(nt).isa))];
#D              fi;
#D            end,
#D            apply := function(nt,c,cnt)
#D              local myCross,myScat,mydims;
#D#little hack for KernelDup
#D              if (ObjId(nt.params[1]).name<>"KernelMMMDuped") then
#D                  myCross:=Cross(List(nt.dims()[2],t->VGath_dup(fId(t),GetFirstTag(nt).v)));
#D                  myCross._children[Last(nt.params[3][1])]:=
#D                  VPrm_x_I(fId(nt.dims()[2][Last(nt.params[3][1])]/GetFirstTag(nt).v),GetFirstTag(nt).v);
#D              else
#D                  mydims:=nt.dims()[2];
#D                  mydims[1]:=mydims[1]*GetFirstTag(nt).v;
#D                  myCross:=Cross(List(mydims,t->VPrm_x_I(fId(t/GetFirstTag(nt).v),GetFirstTag(nt).v)));
#D              fi;
#D              myScat:=VScat(fId(nt.dims()[1]/GetFirstTag(nt).v),GetFirstTag(nt).v);
#D              return myScat*VTensor_OL(c[1],GetFirstTag(nt).v)*myCross;
#D            end),

#Changed by Marek, looks OK
    TTensorI_OL_Vectorize_AParFirst :=rec(
        switch:=false,

        applicable := (self,t) >>
            t.isTag(1, AVecReg)
            and 0 <> First(t.params[3][1], x -> true)
            and APar = t.params[2][First(t.params[3][1],x -> true)]
            and t.firstTag().v = t.params[4][First(t.params[3][1], x -> true)],

        freedoms := nt -> [[1]],

        child := function (nt,freedoms)
            local PV,outorder,theguy;

            PV := Copy(nt.params[2]);
            PV[First(nt.params[3][1],x->true)] := AVec;

            theguy := Copy(First(nt.params[3][1], x->true));
            outorder := Drop(Copy(nt.params[3][1]),1);
            Add(outorder,theguy);

            return [
                TTensorI_OL( nt.params[1], PV, [outorder], nt.params[4], nt.params[5]),
                TL(nt.rng()[1].size, nt.firstTag().v, 1, 1, nt.params[5]),
                TL(
                    nt.dmn()[First(nt.params[3][1], x -> true)].size,
                    nt.dmn()[First(nt.params[3][1], x -> true)].size / nt.firstTag().v,
                    1,
                    1,
                    nt.params[5]
                )
            ];
        end,

        apply := function(nt,c,cnt)
            local i;

            i := Identities(nt.dmn());
            i._children[First(nt.params[3][1], x -> true)] := c[3];

            return c[2] * c[1] * i;
        end
    ),

#D            applicable :=(self,t) >> FirstTagEq(t, AVecReg) and First(t.params[3][1],x->true)<>0 and t.params[2][First(t.params[3][1],x->true)]=APar and t.params[4][First(t.params[3][1],x->true)]=GetFirstTag(t).v,
#D            freedoms := nt -> [[1]],
#D            child := function (nt,freedoms)
#D              local PV,outorder,theguy;
#D              PV:=Copy(nt.params[2]);
#D              PV[First(nt.params[3][1],x->true)]:=AVec;
#D              theguy:=Copy(First(nt.params[3][1],x->true));
#D              outorder:=Drop(Copy(nt.params[3][1]),1);
#D              Add(outorder,theguy);
#D              return [TTensorI_OL(nt.params[1],PV,[outorder],nt.params[4],nt.params[5]),TL(nt.rng()[1].size,GetFirstTag(nt).v,1,1,nt.params[5]),TL(nt.dmn()[First(nt.params[3][1],x->true)].size,nt.dmn()[First(nt.params[3][1],x->true)].size/GetFirstTag(nt).v,1,1,nt.params[5])];
#D            end,
#D            apply := function(nt,c,cnt)
#D              local i;
#D              i:=Identities(nt.dmn());
#D              i._children[First(nt.params[3][1],x->true)]:=c[3];
#D              return c[2]*c[1]*i;
#D            end),
));


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

# Here I define the new layer of "index-free" Sigma-SPL.  I will use
# this for library generation. The goal of index-free formulas is to
# eliminate loop indices completely. In such formulas objects may
# still depend on the loop indices, because gather/scatter/diagonal
# functions become multivariate (vs. univariate) functions, taking
# loop indices + point index as arguments.
#
# See "Generalized Problem Spec" write-up.
#

Class(GTIndexFunction, FuncClass, rec(
    isGTIndexFunction := true,

    # .upRank (=lift in spl.maude) is is the inverse of reduction .downRank,
    # it adds implicit dependency on the (new) inner loop, this happens when, eg.
    # gather is pulled into a loop:  Gath(f) ISum(s) => ISum(Gath(f.upRank()) * s)
    upRank := self >> Error("Not implemented"),

    upRankBy := (self, n) >> Checked(IsPosInt0(n),
        FoldL([1..n], (f, i) -> f.upRank(), self)),

    # .downRank (=red in spl.maude) is the "reduction" operator for
    # index mapping functions it eliminates implicit deendencies on loops
    # by introducing explicit references to loop variables ind(..)
    #
    # Eliminates dependency on loop <loopid>, 1 denotes innermost loop
    # <ind> is a loop index
    downRank := (self, loopid, ind) >> Error("Not implemented"),

    # returns the rank, i.e., the number of implicit nested loops
    rank := self >> Error("Not implemented"),

    # eliminated all loop dependencies
    downRankFull := (self, inds) >>
        FoldL(Reversed([1..Length(inds)]), (f, i) -> f.downRank(i, inds[i]), self),

    # split a loop, increase rank by 1,
    # if a loop had n iterations, the new outer loop will have <outer_its> iterations,
    # and next inner <inner_its> iterations
    split := (self, loopid, inner_its, outer_its) >> Error("not implemented"),

    # rotate the ranks (see Doc(frotate) for explanation)
    rotate := (self, n) >> Error("not implemented"),

    isUnitStride := self >> false,
    vecStride := (self, rank) >> self.params[4][rank+1],

));

Class(GTOps, rec(
    Print := PrintOps.Print,
    \= := (a,b) -> When(IsBound(a.equals), a.equals(b), b=a)
));

IsGTIndexFunction := o -> IsRec(o) and IsBound(o.isGTIndexFunction) and o.isGTIndexFunction;

# slow implementation, Haskell would love this
_dropTail := (lst, arg) -> Cond(Length(lst)=0,lst,Cond(Last(lst)=arg, _dropTail(DropLast(lst, 1),arg), lst));

Class(HHBase, GTIndexFunction, rec(
    abbrevs := [ (N,n,b,strides) -> [toExpArg(N), toExpArg(n), toExpArg(b),
                                     List(_dropTail(strides, 0), toExpArg)] ],
    def := (N, n, b, strides) -> Checked(IsList(strides), ForAll(strides, IsPosInt0Sym),
	rec(N := N, n := n)),

    nloops := self >> Length(self.params[4]),
    free := self >> FreeVars(self.params),

    toSpl := (self, inds, kernel_size) >> self.downRankFull(inds),
    rank := self >> Maximum(0, Length(self.params[4])-1),  # rank = number of implicit loop dependencies

    domain := self >> When(IsValue(self.params[2]), self.params[2].v, self.params[2]),
    range := self >> When(IsValue(self.params[1]), self.params[1].v, self.params[1]),

    base := self >> self.params[3],
    strides := self >> Cond(self.params[4] = [], [0], self.params[4]),

    # rotate(<n>)
    #   switch the order of loops (represented by ranks), by making <n>-th loop innermost
    rotate := meth(self, n)
       local params, strides, rank;
       rank := self.rank();
       if rank<=1 and n<=1 then return self;
       elif n > rank then return self.upRank();
       fi;

       params := ShallowCopy(self.params);
       strides := ShallowCopy(params[4]);
       strides := [strides[1], strides[n+1]] :: strides{[2..n]} :: strides{[n+2 .. rank+1]};
       params[4] := strides;
       return self.from_rChildren(params);
    end,

    isUnitStride := self >> false, # redefine if unit stride is possible
));

Declare(KH, BHH, HH, HHZ);

#F HH(<N>, <n>, <b>, <strides>) - generalized stride index mapping function
#F   N - range
#F   n - domain
#F
#F   HH represents a multivariate function (j1 .. jn are indices of enclosing loops)
#F       (i, j1, ..., jn) -> b + strides[1]*i + strides[2]*j1 + ... + strides[n+1]*jn
#F
#F   Observe that dependence on the loop indices exists, even though no explicit
#F   loop index is present. This is how index free representation works.
#F
Class(HH, HHBase, rec(
    upRank := self >> let(s:=self.strides(),
	HH(self.params[1], self.params[2], self.params[3], Concatenation([s[1], 0], Drop(s, 1)))),

    downRank := (self, loopid, ind) >> Cond(loopid > self.rank(), self, let(s:=self.strides(),
	HH(self.params[1], self.params[2], s[loopid+1] * ind + self.params[3], ListWithout(s, loopid+1)))),

    split := (self, loopid, inner_its, outer_its) >> Cond(loopid > self.rank(), self, let(
        strides := self.params[4],
        vs := self.params[4][1+loopid],
	HH(self.params[1],
           self.params[2],
           self.params[3],
           Concatenation(strides{[1..loopid+1]}, [inner_its * vs], strides{[loopid+2..Length(strides)]})))),

    lambda := self >> let(
        n := self.params[2],
        b := self.params[3],
        s := self.strides(),
        ind := List([1..Length(s)], i -> When(i=1, Ind(n), var.fresh_t("q", TInt))),
        Lambda(Reversed(ind), b + Sum([1..Length(s)], i -> s[i] * ind[i]))),

    isUnitStride := self >> self.strides()[1] = 1
));

Class(BHH, HHBase, rec(
    abbrevs := [ (N,n,b,strides,refl) -> [toExpArg(N), toExpArg(n), toExpArg(b), List(_dropTail(strides,0), toExpArg), toExpArg(refl)] ],
    def := (N, n, b, strides, refl) -> Checked(
        IsList(strides), ForAll(strides, IsPosInt0Sym), IsPosIntSym(refl),
	rec(N := N, n := n)),

    upRank := self >> let(s:=self.strides(),
	BHH(self.params[1], self.params[2], self.params[3], Concatenation([s[1], 0], Drop(s, 1)), self.params[5])),

    downRank := (self, loopid, ind) >> Cond(loopid > self.rank(), self, let(s:=self.strides(),
	BHH(self.params[1], self.params[2], self.params[3] + s[loopid+1] * ind,
           ListWithout(s, loopid+1), self.params[5]))),

    split := (self, loopid, inner_its, outer_its) >> Cond(loopid > self.rank(), self, let(
        strides := self.strides(),
        vs := self.params[4][1+loopid],
	BHH(self.params[1],
           self.params[2],
           self.params[3],
           Concatenation(strides{[1..loopid+1]}, [inner_its * vs], strides{[loopid+2..Length(strides)]}),
           self.params[5]))),

    lambda := self >> let(
        n := self.params[2],
        b := self.params[3],
        s := self.strides(),
        refl := self.params[5],
        inds := List([1..Length(s)], i -> When(i=1, Ind(n), var.fresh_t("q", TInt))),
        res := b + Sum([1..Length(s)], i -> s[i] * inds[i]),
        Lambda(Reversed(inds), cond(leq(inds[1], idiv(n-1,2)), res, refl-res))
    ),
));

#F KH(<N>, <n>, <b>, <strides>, <corr>) - generalized (I_m dirsum J_m dirsum I_m ...) L^mn_n thingy
#F   N - range
#F   n - domain
#F
#F   KH represents a multivariate function (j1 .. jn are indices of enclosing loops)
#F       (i, j1, ..., jn) -> b + strides[1]*i + strides[2]*j1 + ... + strides[n+1]*jn
#F
#F   Observe that dependence on the loop indices exists, even though no explicit
#F   loop index is present. This is how index free representation works.
#F
Class(KH, HHBase, rec(
    abbrevs := [ (N,n,b,strides,corr) -> [toExpArg(N), toExpArg(n), toExpArg(b), List(_dropTail(strides,0), toExpArg), List(corr, toExpArg)] ],
    def := (N, n, b, strides, corr) -> Checked(
        IsList(strides), ForAll(strides, IsPosInt0Sym),
        IsList(corr),    ForAll(corr, IsIntSym), Length(corr) = 2,
	rec(N := N, n := n)),

    upRank := self >> let(s:=self.strides(),
	KH(self.params[1], self.params[2], self.params[3],
           Concatenation([s[1], 0], Drop(s, 1)),
           self.params[5])),

    downRank := (self, loopid, ind) >> Checked(loopid=1, self.rank()=1, let(
        s:=self.strides(), corr:=self.params[5], b:=self.params[3],
	KH(self.params[1],
           self.params[2],
           b,
           [s[1]],
           [ s[2]*ind      + cond(neq(imod(ind,2),0), corr[2], corr[1]),
             s[1]-s[2]-ind + cond(neq(imod(ind,2),0), corr[1], corr[2]) ]))),

    lambda := self >> Checked(Length(self.strides())=1, let(
        n := self.params[2],
        b := self.params[3],
        s := self.strides()[1],
        i := Ind(n),
        corr_even := self.params[5][1], corr_odd := self.params[5][2],
        Lambda(i, b + s*i + imod(i+1, 2)*corr_even + imod(i, 2)*corr_odd)))
));

#F HHZ(<N>, <n>, <b>, <strides>) - generalized stride mod N index mapping (prime-factor FFT)
#F   N - range
#F   n - domain
#F
#F   HHZ represents a multivariate function (j1 .. jn are indices of enclosing loops)
#F       (i, j1, ..., jn) -> (b + strides[1]*i + strides[2]*j1 + ... + strides[n+1]*jn) mod N
#F
#F   Observe that dependence on the loop indices exists, even though no explicit
#F   loop index is present. This is same as in HH.
#F
Class(HHZ, HHBase, rec(
    split := (self, loopid, inner_its, outer_its) >> Cond(loopid > self.rank(), self, let(
        strides := self.params[4],
        vs := self.params[4][1+loopid],
	HHZ(self.params[1],
           self.params[2],
           self.params[3],
           Concatenation(strides{[1..loopid+1]}, [inner_its * vs], strides{[loopid+2..Length(strides)]})))),

    upRank := self >> let(s:=self.strides(),
	HHZ(self.params[1], self.params[2], self.params[3], Concatenation([s[1], 0], Drop(s, 1)))),

    downRank := (self, loopid, ind) >> Cond(loopid > self.rank(), self, let(s:=self.strides(),
	HHZ(self.params[1], self.params[2], self.params[3] + s[loopid+1]*ind, ListWithout(s, loopid+1)))),

    lambda := self >> let(hlambda := ApplyFunc(HH, self.params).lambda(),
        CopyFields(hlambda, rec(expr := imod(hlambda.expr, self.params[1]))))
));

Declare(UU);

##    UU(<N>,<n>,<base>,<stride>,<columns>,<leading_dim>,
##                                 < <vstride_cols, vstride_rows>, ... > )
##
##    index mapping:
##       (i, j1, ..., jn) ->
##          base + stride * (i / columns) + leading_dim * (i mod columns)
##               + vstride_cols[1]*j1 + ... + vstride_cols[n]*jn
##               + (vstride_rows[1]*j1 + ... + vstride_rows[n]*jn)* leading_dim
##
##    This is roughly the equivalent of an HH function for 2D inputs.

Class(UU, GTIndexFunction, rec(
    abbrevs := [ (N, n, bX, bY, s, c, ld, strides) -> [toExpArg(N), toExpArg(n), toExpArg(bX), toExpArg(bY),
                                            toExpArg(s), toExpArg(c), toExpArg(ld),
                                     List(_dropTail(strides,[0,0]), x->List(x,toExpArg))] ],
    def := (N, n, bX, bY, s, c, ld, strides) -> Checked(IsList(strides), ForAll(strides, x->
                                              IsList(x) and Length(x)=2 and ForAll(x,IsPosInt0Sym)),
                                                rec(N := N, n := n)),

    nloops := self >> Length(self.params[8]),
    free := self >> FreeVars(self.params),

    toSpl := (self, inds, kernel_size) >> self.downRankFull(inds),
    rank := self >> Length(self.params[8]),   # rank = number of implicit loop dependencies

    domain := self >> When(IsValue(self.params[2]), self.params[2].v, self.params[2]),
    range := self >> When(IsValue(self.params[1]), self.params[1].v, self.params[1]),

    upRank := self >> let(ss:=self.params[8],
	UU(self.params[1], self.params[2], self.params[3], self.params[4], self.params[5], self.params[6], self.params[7],
            Concatenation([[0,0]], ss))),

    downRank := (self, loopid, ind) >> Cond(loopid > self.rank(), self, let(ss:=self.params[8],
	UU(self.params[1], self.params[2], self.params[3] + ss[loopid][1] * ind, self.params[4] + ss[loopid][2] * ind * self.params[7],
            self.params[5], self.params[6], self.params[7], ListWithout(ss, loopid)))),

    split := (self, loopid, inner_its, outer_its) >> Cond(loopid > self.rank(), self, let(
        strides := self.params[8],
        vs := self.params[8][loopid],
        UU(self.params[1],
           self.params[2],
           self.params[3],
           self.params[4],
           self.params[5],
           self.params[6],
           self.params[7],
           Concatenation(strides{[1..loopid]}, [inner_its * vs], strides{[loopid+1..Length(strides)]})))),

    lambda := self >> let(
        n := self.params[2],
        bX := self.params[3],
        bY := self.params[4],
        s := self.params[5],
        c := self.params[6],
        ld := self.params[7],
        ss := self.params[8],
        ind := List([1..Length(ss)], i -> var.fresh_t("q", TInt)),
        index := Ind(n),
        Lambda(Concatenation(Reversed(ind),[index]), (bX  + imod(index,c) * s) + (bY+idiv(index,c) * ld) + Sum([1..Length(ss)], i -> ss[i][1] * ind[i]) + Sum([1..Length(ss)], i -> ss[i][2] * ind[i] * ld))),

    isUnitStride := self >> Error("isUnitStride is not supported for now with UU"),
));

Declare(XChain);

#F XChain( <permuted [1..n]> ) - fTensor chain of fBase's and single fId.
#F
#F   This represents a subset of functions captured by HH. Namely those
#F   functions that described tightly packed index space, i.e.,
#F   if f is an XChain in k-nested loop with dimensions in loop_dims, then
#F
#F   domain(f) * Product(loop_dims) = range(f)
#F
Class(XChain, GTIndexFunction, rec(
    def := perm -> Checked(IsList(perm), ForAll(perm, IsPosInt0),
	Set(Copy(perm))=[0..Length(perm)-1],
        rec()),

    range := self >> 0,
    domain := self >> 0,

    equals := (self, o) >> ObjId(self)=ObjId(o) and self.params[1]=o.params[1],

    printFull := self >> Print(self.name, "(", self.params[1], ")"),
    printShort := self >> Print("x(", PrintCS(self.params[1]), ")"),

    toSpl := (self, inds, kernel_size) >> let(fbases := List(inds, fBase),
	ApplyFunc(fTensor, List(self.params[1], i -> When(i=0, fId(kernel_size), fbases[i])))),

    toDiag := (self, loop_dims, kernel) >> let(fconst := List(loop_dims, d->fConst(d,1)),
	ApplyFunc(diagTensor, List(self.params[1], i -> When(i=0, kernel, fconst[i])))),

    without := (self, loopid) >> Checked(loopid >= 1, loopid <= Length(self.params[1])-1, let(
	lst := ListWithout(self.params[1], Position(self.params[1], loopid)),
	XChain(List(lst, x -> When(x < loopid, x, x-1))))),

    part := (self, loopid, ind, kernel_size, loop_dims) >>
        ApplyFunc(fTensor, List(self.params[1], i ->
	    Cond(i=0,      fId(kernel_size),
		 i=loopid, fBase(ind),
		 fId(loop_dims[i])))),

    composeWith := (self, f) >> When(ObjId(f)<>XChain, Error("Can only compose with another XChain"),
	let(fperm := f.params[1],
	    myperm := self.params[1],
	    pos := Position(myperm, 0),
	    XChain(Concatenation(1+myperm{[1..pos-1]}, fperm, 1+myperm{[pos+1..Length(myperm)]})))),
));

# NOTE: toSpl -> SubstTopDown( ind -> inds[i] )
# Normal functions retrofitted for GT
fId.toSpl  := (self, inds, ksize) >> self;
fId.without:= (self, loopid) >> self;
fId.part   := (self, loopid, loopvar, kernel_size, loop_dims) >> self;
#fId.toDiag := (self, loop_dims) >> fConst(TReal, self.params[1], 1);

fBase.toSpl  := (self, inds, ksize) >> Cond(
    ObjId(self.params[2])=ind, fBase(inds[self.params[2].n]),
    self);

fBase.without := (self, loopid) >> Cond(
    self.params[2]=ind(self.params[1], loopid), fBase(self.params[1], 0),
    self);

fBase.part := (self, loopid, loopvar, ksize, loop_dims) >> Cond(
    self.params[2]=ind(self.params[1], loopid), fAdd(self.params[1], self.params[1], loopvar),
    fId(self.params[1]));

fTensor.part := (self, loopid, loopvar, ksize, loop_dims) >>
     ApplyFunc(fTensor, List(self.children(), c -> c.part(loopid, loopvar, ksize, loop_dims)));

fTensor.without := (self, loopid) >>
     ApplyFunc(fTensor, List(self.children(), c -> c.without(loopid)));

fTensor.toSpl := (self, inds, ksize) >>
     ApplyFunc(ObjId(self), List(self.children(), c -> c.toSpl(inds, ksize)));

fCompose.toSpl := (self, inds, ksize) >>
     ApplyFunc(fCompose, List(self.children(), c -> c.toSpl(inds, ksize)));

FuncClass.toSpl := (self, inds, ksize) >> self;

#FuncClassOper.rank := self >> Maximum(List(self.children(), x->x.rank()));
#FuncClass.rank := self >> 0;
#DiagFunc.rank := self >> 0;

FuncClassOper.isUnitStride := self >> false;
FuncClass.isUnitStride     := self >> false;
fId.isUnitStride := self >> true;
fTensor.isUnitStride := self >> ForAll(self._children, x->x.isUnitStride());
#fBase.toDiag := (self, loop_dims) >> fConst(TReal, self.params[1], 1);
#fId.part   := (self, loopid, ind, kernel_size, loop_dims) >> self;

ListAddP := function(l1, l2, c)
    local i, lenl1, lenl2, res, maxlen, minlen;
    lenl1 := Length(l1);
    lenl2 := Length(l2);
    minlen := Minimum(lenl1, lenl2);
    maxlen := Maximum(lenl1, lenl2);
    res := Concatenation(l1, Replicate(maxlen-lenl1, c));
    for i in [1..lenl2] do
        res[i] := res[i] + l2[i];
    od;
    return res;
end;

ListAddZP := (l1, l2) -> ListAddP(l1,l2,0);
ListAddLP := (l1, l2) -> ListAddP(l1,l2,[0,0]);


text:
# Copyright (c) 2018-2021, Carnegie Mellon University
# See LICENSE for details

# This rule seems to be unnecessary
#

RulesFor(WHT, rec(
    #F WHT_Dirsum: switched off...
    #F
    #F see tensor products as direct sums so that they can be
    #F split separately
    #F
    #F   WHT_(2^k) =
    #F     prod_(i = 1)^r
    #F       dirsum_(a = 1)^2^(k_1 + .. + k_(i-1))
    #F         ( dirsum_(b = 1)^2^(k_(i+1) + .. + k_r) WHT_(2^(k_i))) ^
    #F           L^(2^(k_i + .. + k_r))_(2^(k_(i+1) + .. + k_r))
    #F         )
    #F
    WHT_Dirsum := rec (
	info             := "WHT_(2^k) -> prod dirsum (dirsum WHT_(2^ki))^L",
	forTransposition := false,
	switch           := false,
	isApplicable     := P -> P <> 1,
	isDerivable := meth ( self, S, C )
	    local C1, i, c, j;
	    if not (IsApplicableRule(self, S) and
		    ForAll(C, c -> ObjId(c) = WHT)) then return false; fi;
	    i  := 1;
	    C1 := [ ];
	    while i <= Length(C) do
	        c := C[i];
		for j in [i..i + 2^(S.params - c.params) - 1] do
		    if c.params <> C[j].params then
			return false;
		    fi;
		od;
		Add(C1, c);
		i := i + 2^(S.params - c.params);
	    od;
	    return Sum(C1, c -> c.params) = S.params;
	end,

	allChildren := function ( P )
	    local C, C1, p, c, i, k;
	    C := OrderedPartitions(P);
	    Unbind(C[Length(C)]);
	    C1 := [ ];
	    for p in C do
	         c := [ ];
		 for i in p do
		     for k in [1..2^(P - i)] do Add(c, WHT(i)); od;
		 od;
		 Add(C1, c);
	    od;
	    return C1;
	end,

	randomChildren := function ( P )
	    local C, sum, rand;
	    C   := [ ];
	    sum := 0;
	    while sum < P do
	        if sum = 0 then	rand := RandomList( [1..(P-1)] );
		else		rand := RandomList( [1..(P-sum)] );
		fi;
		Add( C, rand );
		sum := sum + rand;
	    od;
	    return Concatenation(List(C, c -> List([1..2^(P - c)], i -> WHT(c))));
	end,

	rule := function ( P, C )
	    local C1, c, left, right, factors, i;
	    i  := 1;
	    C1 := [ ];
	    while i <= Length(C) do
	        c := C[i];
		Add(C1, c);
	        i := i + 2^P/c.dimensions[1];
	    od;

	    left    := 1;
	    right   := 2^P / C1[1].dimensions[1];
	    factors := [ ];

            # first factor
	    Add(factors,
		Conjugate(
		    DirectSum(List([1..right], l -> C1[1])),
		    L(2^P, right))
	    );

	    left  := left * C1[1].dimensions[1];
	    right := right / C1[2].dimensions[1];

            # middle factors
	    for i in [2..Length(C1) - 1] do
  	        Add(factors, DirectSum(
		    List([1..left],
			k -> Conjugate(
			    DirectSum(List([1..right], l -> C1[i])),
			    L(right * C1[i].dimensions[1], right))
		    )
		));
		left  := left * C1[i].dimensions[1];
		right := right / C1[i+1].dimensions[1];
	    od;

            # last factor
	    Add(factors, DirectSum(List([1..left], k -> C1[Length(C1)])));
	    return Compose(factors);
	 end
    )
));


text:
method myform.resload()
{
   //this.btn1.caption = " textres( 143 )"
   //image = imageres(10)

}

method myform myform.compinit()
{

   this.width = 400
   this.height = 400
   //this.p_typeid = myform
   this.caption = "Gentee Forms Designer (Demo version 0.1)".ustr()
   //this.onresize = getid( "@m_resize", %{myform} )

   panleft = &newcomp( vPanel, (&this)->vcomp )
   uint xpanleft = panleft
   xpanleft as vPanel
   xpanleft.width = 203
   xpanleft.x = 0
   xpanleft.Border = $brdNone
   xpanleft.vertalign = $alvClient
   xpanleft.horzalign = $alhLeft//$ALIGN_VTOPBOTTOMOFF//$ALIGN_VCLIENT

   //uint xpanright = &newcomp( vPanel, this )//(&this)->vcomp )
   this.panright.owner = this
   uint xpanright as this.panright
   xpanright.x = xpanleft.x + xpanleft.width + 4
   xpanright.vertalign = $alvClient
   xpanright.horzalign = $alhLeftRight
   xpanright.right = 0
   xpanright.name = "right"

   uint t as newcomp( vTab, panleft->vcomp )->vTab
   t.y = 0
   t.height = 85
   t.horzalign = $alhClient
   t.vertalign = $alvTop
   /*uint tp as newcomp( vTabPage, t )->vTabPage
   tp.caption = "Controls"
   this.bt_arrow.owner = tp
   uint bt_bt as this.bt_arrow
   bt_bt.caption = " "
   bt_bt.btnstyle = $bsAsRadioBtn
   bt_bt.y =0
   bt_bt.x=0
   bt_bt.height = 20
   bt_bt.width = 65
   bt_bt.onclick.set( this, "bt_arrowc" )

   //bt_bt as newcomp( vBtn, panleft->vcomp )
   this.bt_btn.owner = tp
   bt_bt as this.bt_btn
   bt_bt.x=65
   bt_bt.y = 0
   bt_bt.height =20
   bt_bt.width =65
   bt_bt.caption = "Btn"
   bt_bt.onclick.set( this, "bt_btnc" )
   //this.bt_arrow.checked = 1
   bt_bt.btnstyle = $bsAsRadioBtn

   this.bt_panel.owner = tp
   bt_bt as this.bt_panel
   bt_bt.x=130
   bt_bt.y = 0
   bt_bt.height =20
   bt_bt.width =65
   bt_bt.caption = "Panel"
   bt_bt.onclick.set( this, "bt_panelc" )
   bt_bt.btnstyle = $bsAsRadioBtn

   this.bt_edit.owner = tp
   bt_bt as this.bt_edit
   bt_bt.x=0
   bt_bt.y = 20
   bt_bt.height =20
   bt_bt.width =65
   bt_bt.caption = "Edit"
   bt_bt.onclick.set( this, "bt_editc" )
   bt_bt.btnstyle = $bsAsRadioBtn

   uint h = 20
   uint y = 90
   bt_bt as newcomp( vBtn, panleft->vcomp )->vBtn
   bt_bt as vBtn
   bt_bt.y = y
   bt_bt.height =h
   bt_bt.caption = "Build"
   bt_bt.vertalign = $alvTop
   bt_bt.horzalign = $alhLeftRight
   bt_bt.x=10
   bt_bt.right =10
   bt_bt.onclick.set( this, "build" )//getid( "@build", %{myform} )

   y += h
   bt_bt as this.bsave
   bt_bt.owner = panleft->vcomp
   bt_bt.y = y
   bt_bt.height =h
   bt_bt.caption = "Save"
   bt_bt.vertalign = $alvTop
   bt_bt.horzalign = $alhLeftRight
   bt_bt.x=10
   bt_bt.right =10
   bt_bt.onclick.set( this, "save" )// =getid( "@save", %{myform} )

   y += h
   bt_bt as this.bopen
   bt_bt.owner = panleft->vcomp
   bt_bt.y = y
   bt_bt.height =h
   bt_bt.caption = "New"
   bt_bt.vertalign = $alvTop
   bt_bt.horzalign = $alhLeftRight
   bt_bt.x=10
   bt_bt.right =10
   bt_bt.onclick.set( this, "new" )// =getid( "@open", %{myform} )

   y += h
   bt_bt as this.babout
   bt_bt.owner = panleft->vcomp
   bt_bt.y = y
   bt_bt.height =h
   bt_bt.caption = "About"
   bt_bt.vertalign = $alvTop
   bt_bt.horzalign = $alhLeftRight
   bt_bt.x=10
   bt_bt.right =10
   bt_bt.onclick.set( this, "about" )

   this.edcur.owner = panleft->vcomp
   this.edcur.horzalign = $alhClient
   this.edcur.readonly = 1
   this.edcur.y = 175
   this.edcur.height = 24

   t as newcomp( vtab, panleft->vcomp )->vtab

   t.y = 200
   t.height = 200
   t.horzalign = $alhClient
   t.vertalign = $alvTopBottom

   tp as newcomp( vtabpage, t )->vtabpage

   //print( "\(&tp)\n" )
   tp.caption = "Properties"

   uint x as newcomp(vproplist, tp)->vproplist
   this.prl = &x
   x.horzalign = $alhClient
   x.vertalign = $alvClient
   plist = &x
   x.onprop.set( this, "propset" )
   x.ongetlist.set( this, "getlist" )

   uint tpe as newcomp( vtabpage, t )->vtabpage
   tpe.caption = "Events"
   //tpe.pageidx = 0
   x as newcomp(vproplist, tpe)->vproplist
   this.evl = &x
   x.horzalign = $alhClient
   x.vertalign = $alvClient//TopBottom
   x.onprop.set( this, "eventset" )
   x.ondblclick.set( this, "eventdblclick" )
   //x.ongetlist.set( this, "getlist" )

   cm.addcomp( "vctrl", vctrl )
   cm.addcomp( "vBtn",  vBtn )
   cm.addcomp( "vPanel", vPanel )
   cm.addcomp( "vEdit", vEdit )
   //cm.addcomp( "vlistbox", vlistbox )
   cm.addcomp( "vform", vform )
   */
   /*uint xwin as newcompdes( vform, xpanright )->vctrl
   edform = &xwin
   xwin.name = this.getnewname( xwin->vform.typename )
   xwin.p_designing = 1
   xwin.x = 0
   xwin.y = 0*/
   /*wined = &newcomp( vwined, xpanright )
   wined->vwined.onselect.set( this, "ctrlselect" )// =getid( "@ctrlselect", %{myform} )
   wined->vwined.ondelete.set( this, "ctrldelete" )
   wined->vwined.onnew.set( this, "wineditnew" )
   this.srcfile = "example"
   wined->vwined.select( edform->vctrl )*/
//print( "xxx\n" )
//   this.load()
//print( "xxx 2\n" )
   this.bt_arrow.Checked = 1
   return this
}

max_stars_repo_name: gentee/c-gentee

max_issues_repo_name: gentee/c-gentee

max_forks_repo_name: gentee/c-gentee


text:
--
-- An instance of this template must have a $Export section and the export_terminals option
--
-- Macros that may be redefined in an instance of this template
--
--     $eof_token
--     $additional_interfaces
--     $super_stream_class -- subclass com.ibm.lpg.Utf8LpgLexStream for getKind
--     $prs_stream_class -- use /.PrsStream./ if not subclassing
--
-- B E G I N N I N G   O F   T E M P L A T E   LexerTemplateD
--
%Options programming_language=csharp,margin=4
%Options table
%options action-block=("*.cs", "/.", "./")
%options ParseTable=LPG2.Runtime.ParseTable
%Options prefix=Char_

--
-- This template requires that the name of the EOF token be set
-- to EOF and that the prefix be "Char_" to be consistent with
-- KeywordTemplateD.
--
%Eof
    EOF
%End

--
-- This template also requires that the name of the parser EOF
-- Token to be exported be set to EOF_TOKEN
--
%Export
    EOF_TOKEN
%End

%Define
    --
    -- Macros that are be needed in an instance of this template
    --
    $eof_token /.$_EOF_TOKEN./

    $additional_interfaces /../
    $super_stream_class /.$file_prefix$Utf8LpgLexStream./
    $prs_stream_class /.IPrsStream./
    $super_class /.object./

    --
    -- Macros useful for specifying actions
    --
    $Header
    /.
                //
                // Rule $rule_number:  $rule_text
                //
                ./

    $DefaultAction
    /.$Header$case $rule_number: { ./

    $BeginAction /.$DefaultAction./

    $EndAction
    /.          break;
                }./

    $BeginJava
    /.$BeginAction
                $symbol_declarations./

    $EndJava /.$EndAction./

    $NoAction
    /.$Header$case $rule_number:
                    break; ./

    $BeginActions
    /.
        public void ruleAction( int ruleNumber)
        {
            switch(ruleNumber)
            {./

    $SplitActions
    /.
	            default:
	                ruleAction$rule_number(ruleNumber);
	                break;
	        }
	        return;
	    }

	    public void ruleAction$rule_number(int ruleNumber)
	    {
	        switch (ruleNumber)
	        {./

    $EndActions
    /.
                default:
                    break;
            }
            return;
        }./
%End

%Globals
    /.
    using LPG2.Runtime;
    using System;
    ./
%End

%Headers
    /.
    public class $action_type : $super_class , RuleAction$additional_interfaces
    {
        private $super_stream_class utf8LexStream;

        private static ParseTable prs = new $prs_type();
        public ParseTable getParseTable() { return prs; }

        private LexParser lexParser = new LexParser();
        public LexParser getParser() { return lexParser; }

        public int getToken(int i) { return lexParser.getToken(i); }
        public int getRhsFirstTokenIndex(int i) { return lexParser.getFirstToken(i); }
        public int getRhsLastTokenIndex(int i) { return lexParser.getLastToken(i); }

        public int getLeftSpan() { return lexParser.getToken(1); }
        public int getRightSpan() { return lexParser.getLastToken(); }

        public void resetKeywordLexer()
        {
            if (kwLexer == null)
                  this.kwLexer = new $kw_lexer_class(utf8LexStream.getInputBytes(), $_IDENTIFIER);
            else this.kwLexer.setInputBytes(utf8LexStream.getInputBytes());
        }

        public void reset(string filename, int tab)
        {
            utf8LexStream = new $super_stream_class(filename, tab);
            lexParser.reset((ILexStream) utf8LexStream, prs, (RuleAction) this);
            resetKeywordLexer();
        }

        public void reset(byte[] input_bytes, string filename)
        {
            reset(input_bytes, filename, 1);
        }

        public void reset(byte[] input_bytes, string filename, int tab)
        {
            utf8LexStream = new $super_stream_class(input_bytes, filename, tab);
            lexParser.reset((ILexStream) utf8LexStream, prs, (RuleAction) this);
            resetKeywordLexer();
        }

        public $action_type(string filename, int tab)
        {
            reset(filename, tab);
        }

        public $action_type(byte[] input_bytes, string filename, int tab)
        {
            reset(input_bytes, filename, tab);
        }

        public $action_type(byte[] input_bytes, string filename)
        {
            reset(input_bytes, filename, 1);
        }

        public $action_type() {}

        public ILexStream getILexStream() { return utf8LexStream; }

        /**
         * @deprecated replaced by {@link #getILexStream()}
         */
        public ILexStream getLexStream() { return utf8LexStream; }

        private void initializeLexer($prs_stream_class prsStream, int start_offset, int end_offset)
        {
            if (utf8LexStream.getInputBytes() == null)
                throw new NullReferenceException("LexStream was not initialized");
            utf8LexStream.setPrsStream(prsStream);
            prsStream.makeToken(start_offset, end_offset, 0); // Token list must start with a bad token
        }

        public void lexer($prs_stream_class prsStream)
        {
            lexer(null, prsStream);
        }

        public void lexer(Monitor monitor, $prs_stream_class prsStream)
        {
            if (utf8LexStream.getInputBytes() == null)
                throw new NullReferenceException("Utf8LexStream was not initialized");

            utf8LexStream.setPrsStream(prsStream);

            prsStream.makeToken(0, 0, 0); // Token list must start with a bad token

            lexParser.parseCharacters(monitor);  // Lex the input characters

            int i = utf8LexStream.getStreamIndex();
            prsStream.makeToken(i, i, $eof_token); // and end with the end of file token
            prsStream.setStreamLength(prsStream.getSize());

            return;
        }
    ./
%End

%Rules
    /.$BeginActions./
%End

%Trailers
    /.
        $EndActions
    }
    ./
%End

--
-- E N D   O F   T E M P L A T E
--

max_stars_repo_name: kuafuwang/LPG-VScode

max_issues_repo_name: kuafuwang/LPG-VScode

max_forks_repo_name: kuafuwang/LPG-VScode