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

Declare(MDRConv);
# From IOPrunedMDRConv, signature (n,h,oblk,opat,iblk,ipat,isFreqData)
# but this one has oblk=1, iblk=1, opat=ipat=[[0..n(1)-1], ...]
# Left with n (was [1]), h, (was [2]), isFreqData (was [7]).
Class(MDRConv, TaggedNonTerminal, rec(
    a_n := self >> self.params[1],
    a_h := self >> self.params[2],
    a_isFreqData := self >> self.params[3],

    abbrevs := [
                (n, h) -> Checked(ForAll(n, IsPosIntSym),
                                  [_unwrap(n), h, false]),
                (n, h, isFreqData) -> Checked(ForAll(n, IsPosIntSym),
                                              [_unwrap(n), h, isFreqData])
               ],
    dims   := self >> let(nprod := Product(self.a_n()),
                          [nprod, nprod]),
    isReal := self >> true,
    normalizedArithCost := self >> let(n := self.a_n(),
                                       IntDouble(5 * n * d_log(n) / d_log(2))),
    TType := TReal,
    terminate := self >> When(self.a_isFreqData(),
               # case of isFreqData = true
               let(nlist := self.a_n(),
                   n1 := nlist[1],
                   nfreq := RClength(n1)/2, # WAS n1/2 + 1
                   nrest := Drop(nlist, 1),
                   idft := Tensor(IPRDFT(n1, -1), I(Product(nrest))) *
                                  Tensor(I(nfreq), L(2*Product(nrest), 2)) *
                                  Tensor(I(nfreq), RC(MDDFT(nrest, -1))),
                   tlist := MatSPL(idft) * self.a_h().list,
                   MDRConv(nlist,
                           FList(TReal, tlist),
                           false).terminate()
                   ),
               # case of isFreqData = false
               let(nlist := self.a_n(),
                   scat3d := Tensor(List(nlist, n->I(n))::[Mat([[1], [0]])]),
                   gath3d := Tensor(List(nlist, n->I(n))::[Mat([[1, 0]])]),
                   dft3dr := RC(MDDFT(nlist, -1)),
                   idft3dr := RC(MDDFT(nlist, 1)),
                   gfd := List((1/Product(nlist)) *
                               MatSPL(MDDFT(nlist, 1)) *
                               self.a_h().list,
                               ComplexAny),
                   gdiagr := RC(Diag(gfd)),
                   t := gath3d * idft3dr * gdiagr * dft3dr * scat3d,
                   t.terminate()
                  )
        ), # end of When

    hashAs := self >> ApplyFunc(ObjId(self),
                                [self.a_n(),
                                 fUnk(self.a_h().range(),
                                 self.a_h().domain())]::Drop(self.params, 2))
));

# used to represent convolution with real valued conjugate even symbol
Class(MDRConvR, MDRConv);


text:
Read("SchreierSims.gi");
LoadPackage( "AtlasRep" );

MicroSeconds := function()
  local t;
  t := IO_gettimeofday();
  return t.tv_sec * 1000000 + t.tv_usec;
end;

# choose which version of Schreier-Sims-Algorithm will be applied
# randomized = true will apply the randomized Schreier-Sims,
# randomited = false will apply the deterministic algorithm
randomized := true;

# group is the short cut for the group that will be investigated
# is used by AtlasGroups to load group information
group := "M11";

# define your group here:
G := AtlasGroup(group);
#G := SymmetricGroup(100);

# Rubik group
#G:= Group(
#( 1, 3, 8, 6)( 2, 5, 7, 4)( 9,33,25,17)(10,34,26,18)(11,35,27,19),
#( 9,11,16,14)(10,13,15,12)( 1,17,41,40)( 4,20,44,37)( 6,22,46,35),
#(17,19,24,22)(18,21,23,20)( 6,25,43,16)( 7,28,42,13)( 8,30,41,11),
#(25,27,32,30)(26,29,31,28)( 3,38,43,19)( 5,36,45,21)( 8,33,48,24),
#(33,35,40,38)(34,37,39,36)( 3, 9,46,32)( 2,12,47,29)( 1,14,48,27),
#(41,43,48,46)(42,45,47,44)(14,22,30,38)(15,23,31,39)(16,24,32,40) );

S := Set(GeneratorsOfGroup(G));
n := LargestMovedPoint(G);
B := [];

Print("*********************************************************\n");
if randomized then
    Print("Apply randomized Schreier-Sims-Algorithm to initial sets...\n\n");
else
   Print("Apply deterministic Schreier-Sims-Algorithm to initial sets...\n\n");
fi;
Print("Info:\n");
Print(" - Group: "); Print(group); Print(" on "); Print(n); Print(" points\n");
Print(" - Size: "); Print(Size(G)); Print("\n");

t1 := MicroSeconds();
if randomized then
    BSGS := RandomSchreierSims(B, S, 2, n);
else
    BSGS := SchreierSims(B, S, n);
fi;
Add(BSGS[2], [()]);
Print("\nProcedure done. Possible BSGS calculated. Time needed: "); Print((MicroSeconds()-t1) * 1.0 / 1000000); Print("\n");

Print("Base: "); Print(BSGS[1]); Print("\n");
Print("Size of orbits: "); Print(BSGS[3]); Print("\n");

if randomized then
    Print("Length of SGS: "); Print(Length(BSGS[2][1])); Print("\n");
    Print("Remove redundant generators from SGS ...\n");
    RemoveRedundantGens(BSGS, S, n);
    Print("SGS reduced. New length of SGS: "); Print(Length(BSGS[2][1])); Print("\n");
fi;

#dummy := [];
#for i in [1..Length(BSGS[1])] do
#    Add(dummy, OrbitSV(BSGS[1][i], BSGS[2][i], n));
#od;
#Add(BSGS, dummy);

#MakeImmutable(BSGS);

ext: gi

max_stars_repo_name: AlexanderKlemps/parallel-schreier-sims

max_issues_repo_name: AlexanderKlemps/parallel-schreier-sims

max_forks_repo_name: AlexanderKlemps/parallel-schreier-sims


text:
--/**
-- * <copyright>
-- *
-- * Copyright (c) 2010, 2009 IBM Corporation and others.
-- * All rights reserved.   This program and the accompanying materials
-- * are made available under the terms of the Eclipse Public License v2.0
-- * which accompanies this distribution, and is available at
-- * http://www.eclipse.org/legal/epl-v20.html
-- *
-- * Contributors:
-- *   See (or edit) Notice Declaration below
-- *
-- * </copyright>
-- */
--
-- The EssentialOCL Parser
--

%Define
    -- Redefinition of macros used in the parser template
    --
    $default_repair_count /.getDefaultRepairCount()./
	$super_parser_class /.AbstractOCLParser./
    $prs_stream_class /.DerivedPrsStream./

	-- Definition of new macros used by the grammar file
	-- which may be redefined by extended files.
    $copyright_contributions /.*./

	-- Definition of new macros used by the grammar file
	-- which are not intended to be extended.
	$lpg_ns /.lpg.runtime./ -- package namespace of the LPG Runtime API

	-- Some useful macros
    $NewCase
    /. $Header
                case $rule_number:./

    $EmptyListAction -- Deprecated, code inline with correct generic parameter type
    /. $Header
                case $rule_number:
                    setResult(new BasicEList<Object>());
                    break;./

    -- BeginJava and EndJava need to be reworked in order to be able to properly use $NewCase macro

    -- BeginJava does nothing
	-- block-actions should call BeginCode, instead
    $BeginJava /../

  	-- EndJava does nothing
	-- block-actions should call EndCode, instead
	$EndJava /../

	$BeginCode
	/.$BeginAction
					$symbol_declarations./

	$EndCode /.$EndAction./

%End

%Notice
    /./**
 * Essential OCL Grammar
 * <copyright>
 *
 * Copyright (c) 2010, 2010 IBM Corporation and others.
 * All rights reserved.   This program and the accompanying materials
 * are made available under the terms of the Eclipse Public License v2.0
 * which accompanies this distribution, and is available at
 * http://www.eclipse.org/legal/epl-v20.html
 *
 * Contributors:
 *   IBM - Initial API and implementation
 *   E.D.Willink - Elimination of some shift-reduce conflicts
 *   E.D.Willink - Remove unnecessary warning suppression
 *   E.D.Willink - Bugs 184048, 225493, 243976, 259818, 282882, 287993, 288040, 292112, 295166
 *   Borland - Bug 242880
 *   Adolfo Sanchez-Barbudo Herrera (Open Canarias):
 *        - 242153: LPG v 2.0.17 adoption.
 *        - 299396: Introducing new LPG templates
 *        - 300534: Removing the use of deprecated macros.
 *******************************************************************************/
    ./
%End

%Globals
    /.import org.eclipse.emf.common.util.BasicEList;
	import org.eclipse.emf.common.util.EList;
	import org.eclipse.ocl.cst.BooleanLiteralExpCS;
	import org.eclipse.ocl.cst.CSTNode;
	import org.eclipse.ocl.cst.CallExpCS;
	import org.eclipse.ocl.cst.CollectionLiteralExpCS;
	import org.eclipse.ocl.cst.CollectionLiteralPartCS;
	import org.eclipse.ocl.cst.CollectionTypeCS;
	import org.eclipse.ocl.cst.CollectionTypeIdentifierEnum;
	import org.eclipse.ocl.cst.FeatureCallExpCS;
	import org.eclipse.ocl.cst.IfExpCS;
	import org.eclipse.ocl.cst.IntegerLiteralExpCS;
	import org.eclipse.ocl.cst.InvalidLiteralExpCS;
	import org.eclipse.ocl.cst.IsMarkedPreCS;
	import org.eclipse.ocl.cst.IterateExpCS;
	import org.eclipse.ocl.cst.IteratorExpCS;
	import org.eclipse.ocl.cst.LetExpCS;
	import org.eclipse.ocl.cst.NullLiteralExpCS;
	import org.eclipse.ocl.cst.OCLExpressionCS;
	import org.eclipse.ocl.cst.OperationCallExpCS;
	import org.eclipse.ocl.cst.PathNameCS;
	import org.eclipse.ocl.cst.PrimitiveTypeCS;
	import org.eclipse.ocl.cst.RealLiteralExpCS;
	import org.eclipse.ocl.cst.SimpleNameCS;
	import org.eclipse.ocl.cst.SimpleTypeEnum;
	import org.eclipse.ocl.cst.StringLiteralExpCS;
	import org.eclipse.ocl.cst.TupleLiteralExpCS;
	import org.eclipse.ocl.cst.TupleTypeCS;
	import org.eclipse.ocl.cst.TypeCS;
	import org.eclipse.ocl.cst.UnlimitedNaturalLiteralExpCS;
	import org.eclipse.ocl.cst.VariableCS;
	import org.eclipse.ocl.cst.VariableExpCS;
	import org.eclipse.ocl.lpg.DerivedPrsStream;

	import $lpg_ns.BadParseException;
	import $lpg_ns.BadParseSymFileException;
	import $lpg_ns.DiagnoseParser;
	import $lpg_ns.ErrorToken;
	import $lpg_ns.IToken;
	import $lpg_ns.ILexStream;
	import $lpg_ns.Monitor;
	import $lpg_ns.NullExportedSymbolsException;
	import $lpg_ns.NullTerminalSymbolsException;
	import $lpg_ns.ParseTable;
	import $lpg_ns.RuleAction;
	import $lpg_ns.UndefinedEofSymbolException;
	import $lpg_ns.UnimplementedTerminalsException;
    ./
%End

%KeyWords
-- Reserved keywords
    and implies not or xor
    if then else endif
    let in
    false true
    null invalid
    self

-- Restricted keywords
    Bag Collection OrderedSet Sequence Set
    Tuple
    Boolean Integer Real String UnlimitedNatural
    OclAny OclInvalid OclVoid
%End

-- Terminals
%Identifier
    IDENTIFIER
%End

%Terminals

    QUOTED_IDENTIFIER INTEGER_LITERAL REAL_LITERAL STRING_LITERAL

    PLUS     ::= '+'
    MINUS    ::= '-'
    MULTIPLY ::= '*'
    DIVIDE   ::= '/'

    GREATER       ::= '>'
    LESS          ::= '<'
    EQUAL         ::= '='
    GREATER_EQUAL ::= '>='
    LESS_EQUAL    ::= '<='
    NOT_EQUAL     ::= '<>'

    LPAREN   ::= '('
    RPAREN   ::= ')'
    LBRACE   ::= '{'
    RBRACE   ::= '}'
    LBRACKET ::= '['
    RBRACKET ::= ']'

    ARROW      ::= '->'
    BAR        ::= '|'
    COMMA      ::= ','
    COLON      ::= ':'
    COLONCOLON ::= '::'
    SEMICOLON  ::= ';'
    DOT        ::= '.'
    DOTDOT     ::= '..'
%End

%Headers
	/.

	public $environment_class getOCLEnvironment() {
		return getLexer().getOCLEnvironment();
	}

	@Override
	public $super_lexer_class getLexer() {
		return ($super_lexer_class) super.getLexer();
	}

	// Some methods for backwards compatibility
	/**
	*

	* Before 3.0, this method was used with the now-deprecated  "dollar"getToken macro (which
	* provided token index in the prsStream) to obtain an IToken f a rule given the index of the
	* right hand side token in the said rule. In 3.0 a convenience method has been introduced
	* in order to directly return the IToken, given the index of the right hand side token in the rule.
	*

	*
	*

	* In an action-block of a rule, instead of doing getIToken("dollar"getToken(i))
	* you should do getRhsTokenText(i)
	*

	* @param i the right hand side token index
	* @return the correspondent IToken.
	*
	* @since 3.0
	*/
	@Deprecated
	protected IToken getIToken(int i) {
		return prsStream.getIToken(i);
	}

	/**
	*

	* Before 3.0, this method was used with the now-deprecated "dollar"getToken macro (which
	* provided token index in the prsStream) to obtain an IToken f a rule given the index of the
	* right hand side token in the said rule. In 3.0 a convenience method has been introduced
	* in order to directly return the IToken, given the index of the right hand side token in the rule.
	*

	*
	*

	* In an action-block of a rule, instead of doing getTokenText("dollar"getToken(i))
	* you should do getRhsTokenText(i)
	*

	* @param i the right hand side token index
	* @result the text of the correspondent right hand side IToken.
	*/
	@Deprecated
	protected String getTokenText(int i) {
		return prsStream.getTokenText(i);
	}

	/**
	* A convenience method to obtain the text of a right hand side IToken.
	*
	* @param i the right hand side token index
	* @result the text of the correspondent right hand side IToken.
	*
	* @since 3.0
	*/
	protected String getRhsTokenText(int i) {
		return prsStream.getTokenText(getRhsTokenIndex(i));
	}
	./
%End

%Rules

-----------------------------------------------------------------------
--  Names
-----------------------------------------------------------------------

    reservedKeyword -> and
    reservedKeyword -> else
    reservedKeyword -> endif
    reservedKeyword -> if
    reservedKeyword -> implies
    reservedKeyword -> in
    reservedKeyword -> let
    reservedKeyword -> not
    reservedKeyword -> or
    reservedKeyword -> then
    reservedKeyword -> xor

    tupleKeywordCS ::= Tuple
    reservedKeywordCS ::= reservedKeyword
    restrictedKeywordCS -> CollectionTypeIdentifier
--  restrictedKeywordCS -> BooleanLiteralExpCS
--  restrictedKeywordCS -> InvalidLiteralExpCS
--  restrictedKeywordCS -> NullLiteralExpCS
--  restrictedKeywordCS -> selfKeywordCS
    restrictedKeywordCS -> PrimitiveTypeIdentifier
    restrictedKeywordCS -> tupleKeywordCS

    SimpleName ::= IDENTIFIER
    SimpleName -> QuotedSimpleNameCS
    QuotedSimpleNameCS ::= QUOTED_IDENTIFIER
    QuotedSimpleNameCS ::= QuotedSimpleNameCS STRING_LITERAL

    unreservedSimpleNameCS -> SimpleName
    unreservedSimpleNameCS -> restrictedKeywordCS

	InfixOperator -> '*'
	InfixOperator -> '/'
	InfixOperator -> '+'
	InfixOperator -> '-'
	InfixOperator -> '>'
	InfixOperator -> '<'
	InfixOperator -> '>='
	InfixOperator -> '<='
	InfixOperator -> '='
	InfixOperator -> '<>'
	InfixOperator -> 'and'
	InfixOperator -> 'or'
	InfixOperator -> 'xor'
	InfixOperator -> 'implies'
	InfixOperator -> '.'
	InfixOperator -> '->'

	PrefixOperator -> '-'
	PrefixOperator -> 'not'

-----------------------------------------------------------------------
--  Types
-----------------------------------------------------------------------
    PrimitiveTypeIdentifier ::= Boolean
    PrimitiveTypeIdentifier ::= Integer
    PrimitiveTypeIdentifier ::= Real
    PrimitiveTypeIdentifier ::= String
    PrimitiveTypeIdentifier ::= UnlimitedNatural
    PrimitiveTypeIdentifier ::= OclAny
    PrimitiveTypeIdentifier ::= OclInvalid
    PrimitiveTypeIdentifier ::= OclVoid

    PrimitiveTypeCS ::= PrimitiveTypeIdentifier

    CollectionTypeIdentifier ::= Set
    CollectionTypeIdentifier ::= Bag
    CollectionTypeIdentifier ::= Sequence
    CollectionTypeIdentifier ::= Collection
    CollectionTypeIdentifier ::= OrderedSet

    CollectionTypeCS ::= CollectionTypeIdentifier
    CollectionTypeCS ::= CollectionTypeIdentifier '(' TypeExpCS ')'

    TupleTypeCS ::= Tuple
    TupleTypeCS ::= Tuple '(' tupleTypePartCSlistopt ')'

    tupleTypePartCSlistopt ::= %empty
    tupleTypePartCSlistopt -> tupleTypePartCSlist

    tupleTypePartCSlist ::= tuplePartCS
    tupleTypePartCSlist ::= tupleTypePartCSlist ',' tuplePartCS

    tuplePartCS ::= SimpleName ':' TypeExpCS

-----------------------------------------------------------------------
--  Declarations
-----------------------------------------------------------------------
    VariableDeclarationCS -> SimpleName
    VariableDeclarationCS -> SimpleName '=' OCLExpressionCS
    VariableDeclarationCS -> SimpleName ':' TypeExpCS
    VariableDeclarationCS -> SimpleName ':' TypeExpCS '=' OCLExpressionCS

-----------------------------------------------------------------------
--  Literals
-----------------------------------------------------------------------
--    LiteralExpCS -> CollectionLiteralExpCS
--    LiteralExpCS -> TupleLiteralExpCS
--    LiteralExpCS -> PrimitiveLiteralExpCS
--    LiteralExpCS -> TypeLiteralExpCS

    CollectionLiteralExpCS ::= CollectionTypeCS '{' CollectionLiteralPartCSlistopt '}'

    CollectionLiteralPartCSlistopt ::= %empty
    CollectionLiteralPartCSlistopt -> CollectionLiteralPartCSlist

    CollectionLiteralPartCSlist ::= CollectionLiteralPartCS
    CollectionLiteralPartCSlist ::= CollectionLiteralPartCSlist ',' CollectionLiteralPartCS

    CollectionLiteralPartCS ::= OCLExpressionCS
    CollectionLiteralPartCS ::= OCLExpressionCS '..' OCLExpressionCS

    PrimitiveLiteralExpCS -> NumberLiteralExpCS
    PrimitiveLiteralExpCS -> StringLiteralExpCS
    PrimitiveLiteralExpCS -> BooleanLiteralExpCS
    PrimitiveLiteralExpCS -> UnlimitedNaturalLiteralExpCS
    PrimitiveLiteralExpCS -> InvalidLiteralExpCS
    PrimitiveLiteralExpCS -> NullLiteralExpCS

    TupleLiteralExpCS ::= Tuple '{' TupleLiteralPartCSlist '}'

    TupleLiteralPartCSlist ::= TupleLiteralPartCS
    TupleLiteralPartCSlist ::= TupleLiteralPartCSlist ',' TupleLiteralPartCS

    TupleLiteralPartCS -> SimpleName '=' OCLExpressionCS
    TupleLiteralPartCS -> SimpleName ':' TypeExpCS '=' OCLExpressionCS

    NumberLiteralExpCS ::= INTEGER_LITERAL
    NumberLiteralExpCS ::= REAL_LITERAL

    StringLiteralExpCS ::= STRING_LITERAL
    StringLiteralExpCS ::= StringLiteralExpCS STRING_LITERAL

    BooleanLiteralExpCS ::= true
    BooleanLiteralExpCS ::= false

    UnlimitedNaturalLiteralExpCS ::= '*'

    InvalidLiteralExpCS ::= invalid

    NullLiteralExpCS ::= null

    TypeLiteralExpCS ::= PrimitiveTypeCS
    TypeLiteralExpCS ::= CollectionTypeCS
    TypeLiteralExpCS ::= TupleTypeCS

    TypeNameExpCS -> SimpleName
    TypeNameExpCS ::= SimpleName '(' ')'
    TypeNameExpCS ::= SimpleName '(' UntypedExpressionCSlist ')'
    TypeNameExpCS ::= TypeNameExpCS '::' unreservedSimpleNameCS

    TypeExpCS -> TypeNameExpCS
    TypeExpCS -> TypeLiteralExpCS

-----------------------------------------------------------------------
--  Expressions
-----------------------------------------------------------------------
    OCLExpressionCS -> InfixedExpCS
    OCLExpressionCS -> InfixedLetExpCS

	InfixedExpCS -> PrefixedExpCS
	InfixedExpCS ::= InfixedExpCS InfixOperator PrefixedExpCS

	InfixedLetExpCS -> PrefixedLetExpCS
	InfixedLetExpCS ::= InfixedExpCS InfixOperator PrefixedLetExpCS

	PrefixedExpCS -> PrimaryExpCS
	PrefixedExpCS ::= PrefixOperator PrefixedExpCS

	PrefixedLetExpCS ::= LetExpCS
	PrefixedLetExpCS ::= PrefixOperator PrefixedLetExpCS

	PrimaryExpCS -> NavigatingExpCS
	PrimaryExpCS -> SelfExpCS
	PrimaryExpCS -> PrimitiveLiteralExpCS
	PrimaryExpCS -> TupleLiteralExpCS
	PrimaryExpCS -> CollectionLiteralExpCS
--	PrimaryExpCS -> TypeExpCS
--    PrimaryExpCS -> NameExpCS
    PrimaryExpCS -> TypeLiteralExpCS
	PrimaryExpCS -> IfExpCS
	PrimaryExpCS -> NestedExpCS

    NameExpCS -> SimpleName
    NameExpCS ::= NameExpCS '::' unreservedSimpleNameCS

	IndexExpCS -> NameExpCS
	IndexExpCS ::= NameExpCS '[' ExpCSlist ']'
	IndexExpCS ::= NameExpCS '[' ExpCSlist ']' '[' ExpCSlist ']'

	ExpCSlist ::= OCLExpressionCS
	ExpCSlist ::= ExpCSlist ',' OCLExpressionCS

	NavigatingExpCSbase -> IndexExpCS

	NavigatingExpCS -> NavigatingExpCSbase
	NavigatingExpCS ::= NavigatingExpCSbase '(' ')'
	NavigatingExpCS ::= NavigatingExpCSbase '(' TypedExpressionCS ')'
	NavigatingExpCS ::= NavigatingExpCSbase '(' TypedExpressionCS ',' UntypedExpressionCSlist ')'
	NavigatingExpCS ::= NavigatingExpCSbase '(' TypedExpressionCS '|' UntypedExpressionCSlist ')'
	NavigatingExpCS ::= NavigatingExpCSbase '(' TypedExpressionCS ',' TypedExpressionCS '|' UntypedExpressionCSlist ')'
	NavigatingExpCS ::= NavigatingExpCSbase '(' TypedExpressionCS ';' iteratorAccumulatorCS '|' UntypedExpressionCSlist ')'

    UntypedExpressionCSlist ::= UntypedExpressionCS
    UntypedExpressionCSlist ::= UntypedExpressionCSlist ',' UntypedExpressionCS

	UntypedExpressionCS -> OCLExpressionCS

	TypedExpressionCS ::= OCLExpressionCS
	TypedExpressionCS ::= OCLExpressionCS ':' TypeExpCS

	iteratorAccumulatorCS ::= SimpleName ':' TypeExpCS '=' OCLExpressionCS

    IfExpCS ::= if OCLExpressionCS then OCLExpressionCS else OCLExpressionCS endif

    LetExpCS ::= let LetVariableCSlist in OCLExpressionCS

    LetVariableCSlist ::= LetVariableCS
    LetVariableCSlist ::= LetVariableCSlist ',' LetVariableCS

    LetVariableCS ::= SimpleName '=' OCLExpressionCS
    LetVariableCS ::= SimpleName ':' TypeExpCS '=' OCLExpressionCS

	NestedExpCS  ::= '(' OCLExpressionCS ')'

    SelfExpCS ::= self
%End


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

# tag to force GT_NthLoop expansions to be smaller, results in less
# ruletrees being generated, which is sometimes desirable.
#
# NOTE: it is used in the GT_NthLoop rule at currently needs no
#       parameters

Class(ALimitNthLoop, AGenericTag);

# Input/Output tag which describes the input/outputs of the given block.
# Designed with OL in mind, although written for transforms.
#
# You initialize the tag with a list of pairs which specify which input
# and outputs are connected.
#
# In the case of transforms, this means AIO([1,1]) is the only acceptable
# input. And it means the block is done inplace. Normal operation is
# out-of-place.
#
# for OL, this can be AIO([1,3],[2,1]), which means input1 -> output3 and
# input2 -> output1
#
Class(AIO, AGenericTag);


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

RewriteRules(RulesStrengthReduce, rec(
    re_vpack := Rule([re, vpack], x -> ApplyFunc(vpack, List(x.args[1].args, re))),
    im_vpack := Rule([im, vpack], x -> ApplyFunc(vpack, List(x.args[1].args, im))),

    velem_vpack := Rule([velem, @(1, vpack), @(2, Value)], x -> @(1).val.args[1+@(2).val.v]),

    velem_vdup := Rule([velem, @(1, vdup), @(2)], x -> @(1).val.args[1]),

    fStretch_Drop := Rule( @(1, fStretch, x -> x.num = x.den), e -> e.func),
));


text:
#
# JupyterKernel: Jupyter kernel using ZeroMQ
#
# Implementations
#

# This is a bit ugly: The global variable _KERNEL is assigned to
# the jupyter kernel object at so that we can use it from everywhere.
_KERNEL := "";

## This is plainly wrong, it just reessembles the behaviour of
## `UPDATE_STAT` in newer GAP version we need here.
if not IsBound( UPDATE_STAT ) then
    BindGlobal( "UPDATE_STAT",
      function( string, value )
        time := value;
    end );
fi;

BindConstant( "JUPYTER_KERNEL_MODE_CONTROL", 1 );
BindConstant( "JUPYTER_KERNEL_MODE_EXEC", 2 );

InstallGlobalFunction( JUPYTER_LogProtocol,
function(filename)
    _KERNEL!.ProtocolLog := OutputTextFile(filename, false);
    SetPrintFormattingStatus(_KERNEL!.ProtocolLog, false);
end);

InstallGlobalFunction( JUPYTER_UnlogProtocol,
function()
    local tmp;
    # In case `CloseStream` causes messages to be printed
    tmp := _KERNEL!.ProtocolLog;
    Unbind(_KERNEL!.ProtocolLog);
    CloseStream(tmp);
end);

InstallGlobalFunction( NewJupyterKernel,
function(conf)
    local pid, address, kernel, poll, msg, status, res;

    address := Concatenation(conf.transport, "://", conf.ip, ":");
    kernel := rec( config := Immutable(conf)
                 , Username := "username"
                 , ProtocolVersion := "5.3"
                 , ZmqIdentity := HexStringUUID( RandomUUID() )
                 , SessionKey := conf.key
                 , SessionID := ""
                 , ExecutionCount := 0);

    kernel.MsgHandlers := rec( kernel_info_request := function(msg)
                                 kernel!.SessionID := msg.header.session;
                                 return JupyterMsg( kernel
                                                  , "kernel_info_reply"
                                                  , msg.header
                                                  , rec( protocol_version := kernel!.ProtocolVersion
                                                       , implementation := "GAP"
                                                       , implementation_version := GAPInfo.PackagesInfo.jupyterkernel[1].Version
                                                       , language_info := rec( name := "GAP 4"
                                                                             , version := GAPInfo.Version
                                                                             , mimetype := "text/x-gap"
                                                                             , file_extension := ".g"
                                                                             , pygments_lexer := "gap"
                                                                             , codemirror_mode := "gap"
                                                                             , nbconvert_exporter := "" )
                                                       , banner := Concatenation( "GAP Jupyter kernel ", GAPInfo.PackagesInfo.jupyterkernel[1].Version, "\n",
                                                                                "Running on GAP ", GAPInfo.BuildVersion, "\n")
                                                       , help_links := [ rec( text := "GAP website", url := "https://www.gap-system.org/")
                                                                       , rec( text := "GAP documentation", url := "https://www.gap-system.org/Doc/doc.html")
                                                                       , rec( text := "GAP tutorial", url := "https://www.gap-system.org/Manuals/doc/chap0.html")
                                                                       , rec( text := "GAP reference", url := "https://www.gap-system.org/Manuals/doc/ref/chap0.html") ]
                                                       , status := "ok" )
                                                  , rec() );
                               end,

                               history_request := function(msg)
                                   msg.header.msg_type := "history_reply";
                                   msg.content := rec( history := [] );
                               end,

                               execute_request := function(msg)
                                   local publ, res, rep, r, str, data, metadata, t;

                                   JupyterMsgSend(kernel, kernel!.IOPub, JupyterMsg( kernel
                                                                       , "execute_input"
                                                                       , msg.header
                                                                       , rec( code := msg.content.code
                                                                            , execution_count := kernel!.ExecutionCount )
                                                                       , rec() ) );
                                   str := InputTextString(msg.content.code);

                                   # READ_ALL_COMMANDS was changed from 4.10. We make
                                   # JupyterKernel compatible for the time being (until
                                   # 4.10 is released at least)
                                   t := NanosecondsSinceEpoch();
                                   if CompareVersionNumbers(GAPInfo.Version, "4.10") then
                                       res := READ_ALL_COMMANDS(str, false, false, IdFunc);
                                   else
                                       res := READ_ALL_COMMANDS(str, false);
                                   fi;
                                   # This is probably supremely naughty; we overwrite GAP's
                                   # global time variable
                                   UPDATE_STAT( "time", QuoInt((NanosecondsSinceEpoch() - t), 1000000) );

                                   # Flush StdOut...
                                   Print("\c");
                                   for r in res do
                                       if r[1] = true then
                                           kernel!.ExecutionCount := kernel!.ExecutionCount + 1;

                                           # r[2] contains the result, r[3] is true if a dual semicolon was parsed
                                           if IsBound(r[2]) and r[3] = false then
                                               # FIXME: This is probably doable slightly more nicely
                                               rep := JupyterRender(r[2]);
                                               metadata := JupyterRenderableMetadata(rep);
                                               data := JupyterRenderableData(rep);
                                               # Only send a result message when there is a result
                                               # value
                                               # publ.execution_count := kernel!.ExecutionCount;
                                               JupyterMsgSend(kernel, kernel!.IOPub, JupyterMsg( kernel
                                                                                , "execute_result"
                                                                                , msg.header
                                                                                , rec( transient := rec()
                                                                                , data := data
                                                                                , metadata := metadata
                                                                                , execution_count := kernel!.ExecutionCount )
                                                                                , rec() ) );
                                           fi;
                                       fi;
                                   od;
                                   publ := JupyterMsg( kernel
                                                     , "execute_reply"
                                                     , msg.header
                                                     , rec( status := "ok"
                                                          , execution_count := kernel!.ExecutionCount )
                                                     , rec() );
                                   return publ;
                               end,

                               inspect_request := function(msg)
                                   return JupyterMsg( kernel
                                                    , "inspect_reply"
                                                    , msg.header
                                                    , JUPYTER_Inspect( msg.content.code
                                                                     , msg.content.cursor_pos )
                                                    , rec() );
                               end,

                               complete_request := function(msg)
                                   return JupyterMsg( kernel
                                                    , "complete_reply"
                                                    , msg.header
                                                    , JUPYTER_Complete( msg.content.code
                                                                      , msg.content.cursor_pos )
                                                    , rec() );
                               end,

                               history_request := function(msg)
                                   return JupyterMsg( kernel
                                                    , "history_reply"
                                                    , msg.header
                                                    , rec( history := [] )
                                                    , rec() );
                               end,

                               is_complete_request := function(msg)
                                   return JupyterMsg( kernel
                                                    , "is_complete_reply"
                                                    , msg.header
                                                    , rec( status := "complete" )
                                                    , rec() );
                               end,

                               comm_open := function(msg)
                                   return JupyterMsg( kernel
                                                    , "comm_open_reply"
                                                    , msg.header
                                                    , rec( status := "ok" )
                                                    , rec() );
                               end,

                               comm_info_request := function(msg)
                                   return JupyterMsg( kernel
                                                    , "comm_info_reply"
                                                    , msg.header
                                                    , rec( comms := rec(), status := "ok" )
                                                    , rec() );
                               end,

                               interrupt_request := function(msg)
                                   local status;
                                   # This is SIGINT
                                   status := IO_kill(pid, 2);
                                   return JupyterMsg( kernel
                                                    , "interrupt_reply"
                                                    , msg.header
                                                    , rec()
                                                    , rec() );

                               end,
                               shutdown_request := function(msg)
                                   kernel!.quitting := true;
                                   return JupyterMsg( kernel
                                                 , "shutdown_reply"
                                                 , msg.header
                                                 , rec( restart := msg.content.restart )
                                                 , rec() );
                               end );

    kernel.SignalBusy := function()
        JupyterMsgSend( kernel, kernel!.IOPub
                      , JupyterMsg( kernel
                                  , "status"
                                  , kernel!.CurrentMsg
                                  , rec( execution_state := "busy" )
                                  , rec() ) );
    end;
    kernel.SignalIdle := function()
        JupyterMsgSend( kernel, kernel!.IOPub
                      , JupyterMsg( kernel
                                  , "status"
                                  , kernel!.CurrentMsg
                                  , rec( execution_state := "idle" )
                                  , rec() ) );
    end;

    kernel.HandleShellMsg := function(msg)
        local hdl_dict, f, t, reply;

        # We store the currently processed
        # message header, because we need it
        # for replies
        kernel!.CurrentMsg := msg.header;

        kernel!.SignalBusy();
        t := msg.header.msg_type;
        if IsBound(kernel!.MsgHandlers.(t)) then
            # Currently we send the "reply" to each "request" on the Shell socket
            # here. We might opt to move the sending into the handler functions,
            # since at least "execute" has to send more than one message anyway
            JupyterMsgSend(kernel, kernel!.Shell, kernel!.MsgHandlers.(t)(msg) );

            kernel!.SignalIdle();
            return true;
        else
            Print("unhandled message type: ", msg.header.msg_type, "\n");
            kernel!.SignalIdle();
            return fail;
        fi;

    end;

    kernel.HandleControlMsg := function(msg)
        local hdl_dict, f, t, reply;

        kernel!.CurrentMsg := msg.header;

        t := msg.header.msg_type;
        if IsBound(kernel!.MsgHandlers.(t)) then
            if t in [ "interrupt_request", "shutdown_request" ] then
                JupyterMsgSend(kernel, kernel!.Control, kernel!.MsgHandlers.(t)(msg) );
            fi;
            return true;
        fi;

    end;

    _KERNEL := kernel;

    # This should happen in "Run" somehow, as currently the creation
    # of a Jupyter Kernel breaks the running GAP session, taking
    # every hope of debugging the kernel
    pid := IO_fork();
    if pid = fail then
        return fail;
    elif pid > 0 then # we are the parent and do heartbeat and control messages
        kernel.mode := JUPYTER_KERNEL_MODE_CONTROL;
        kernel.HB := ZmqRouterSocket( Concatenation(address, String(conf.hb_port) ) );
        kernel.Control := ZmqRouterSocket( Concatenation(address, String(conf.control_port) )
                                         , kernel!.ZmqIdentity);
        kernel.quitting := false;
        kernel.Loop := function()
            local topoll, poll, i, msg, res;
            topoll := [ kernel!.HB, kernel!.Control ];
            while true do
                poll := ZmqPoll( topoll, [], 5000 );
                if 1 in poll then
                    msg := ZmqReceiveList(kernel!.HB);
                    ZmqSend(kernel!.HB, msg);
                fi;
                if 2 in poll then
                    msg := JupyterMsgRecv(kernel, kernel!.Control);
                    res := kernel!.HandleControlMsg(msg);
                    if res = fail then
                        Print("failed to handle message\n");
                    fi;
                fi;
                if kernel!.quitting then
                    IO_kill(pid, 3);
                    status := IO_WaitPid(pid, true);
                    QUIT_GAP(0);
                fi;
                # Check whether child has gone away
                status := IO_WaitPid(pid, false);
                if IsRecord(status) then
                    # TODO find out what these statuses mean
                    if status.pid = pid and status.status in [ 3, 9, 15, 131 ] then
                        QUIT_GAP(0);
                    fi;
                fi;
            od;
        end;
    else
        kernel.mode  := JUPYTER_KERNEL_MODE_EXEC; # Handler
        kernel.IOPub := ZmqPublisherSocket( Concatenation(address, String(conf.iopub_port))
                                          , kernel!.ZmqIdentity);
        kernel.Shell := ZmqDealerSocket( Concatenation(address, String(conf.shell_port))
                                       , kernel!.ZmqIdentity);
        kernel.StdIn := ZmqRouterSocket( Concatenation(address, String(conf.stdin_port))
                                       , kernel!.ZmqIdentity);

        # TODO: This is of course still hacky, but better than before
        kernel!.StdOut := OutputStreamZmq(kernel, kernel!.IOPub);
        kernel!.StdErr := OutputStreamZmq(kernel, kernel!.IOPub, "stderr");
        # TODO: Hack to be able to change ERROR_OUTPUT.
        MakeReadWriteGlobal("ERROR_OUTPUT");
        ERROR_OUTPUT := kernel!.StdErr;
        MakeReadOnlyGlobal("ERROR_OUTPUT");
        OutputLogTo(kernel!.StdOut);

        # Jupyter Heartbeat and Control channel is handled by a fork'ed GAP process
        # (yes, really, its better than starting a separate thread, because it
        # doesn't need special pthread code downside is that it doesn't work on
        # cygwin, of course, but maybe we could just ExecuteProcess on windows, or
        # wait for bash on windows to become popular enoug.
        kernel.Loop := function()
            # To catch SIGINT when the kernel is idle
            while true do
                CALL_WITH_CATCH(function()
                                   local topoll, poll, i, msg, res;

                                   topoll := [ kernel!.Shell, kernel!.StdIn ];
                                   while true do
                                       poll := ZmqPoll(topoll, [], 5000);
                                       if 1 in poll then
                                           msg := JupyterMsgRecv(kernel, topoll[1]);
                                           res := kernel!.HandleShellMsg(msg);
                                           if res = fail then
                                               Print("failed to handle message\n");
                                           fi;
                                       fi;
                                       if 2 in poll then
                                           msg := ZmqReceiveList(topoll[2]);
                                       fi;
                                   od;
                               end, []);
            od;
        end;
    fi;

    Objectify(GAPJupyterKernelType, kernel);
    return kernel;
end);

InstallMethod( ViewString
             , "for Jupyter kernels"
             , [ IsGAPJupyterKernel ]
             , x -> "<GAP Jupyter Kernel>" );

InstallMethod( Run
             , "for Jupyter kernel"
             , [ IsGAPJupyterKernel ]
             , function(x)
                 # TODO: we should really not be doing this.
                 MakeReadWriteGlobal("HELP_SHOW_MATCHES");
                 UnbindGlobal("HELP_SHOW_MATCHES");
                 DeclareSynonym("HELP_SHOW_MATCHES", JUPYTER_HELP_SHOW_MATCHES);

                 MakeReadWriteGlobal("HELP");
                 UnbindGlobal("HELP");
                 DeclareSynonym("HELP", JUPYTER_HELP);

                 SetUserPreference("browse", "SelectHelpMatches", false);
                 SetUserPreference("Pager", "tail");
                 SetUserPreference("PagerOptions", "");
                 # This is of course complete nonsense if you're running the jupyter notebook
                 # on your local machine.
                 SetHelpViewer("jupyter_online");
                 x!.Loop();
             end);

InstallGlobalFunction( JUPYTER_KernelStart_HPC,
function(conf)
    Error("HPC-GAP is not supported with this code.");
    QUIT_GAP(0);
end);

InstallGlobalFunction( JUPYTER_KernelStart_GAP,
function(configfile)
    local instream, conf, address, kernel, s;

    instream := InputTextFile(configfile);
    conf := JsonStreamToGap(instream);

    kernel := NewJupyterKernel(conf);
    Run(kernel);
end);

max_stars_repo_name: ZachNewbery/JupyterKernel

max_issues_repo_name: ZachNewbery/JupyterKernel

max_forks_repo_name: ZachNewbery/JupyterKernel


text:
#############################################################################
##
##                                ActionsForCAP package
##
##  Copyright 2015, Sebastian Gutsche, TU Kaiserslautern
##                  Sebastian Posur,   RWTH Aachen
##
##
#############################################################################

DeclareRepresentation( "IsLeftCoactionObjectRep",
                       IsLeftCoactionObject and IsAttributeStoringRep,
                       [ ] );

BindGlobal( "TheFamilyOfLeftCoactionObjects",
            NewFamily( "TheFamilyOfLeftCoactionObjects" ) );

BindGlobal( "TheTypeOfLeftCoactionObjects",
            NewType( TheFamilyOfLeftCoactionObjects,
                     IsLeftCoactionObjectRep ) );

DeclareRepresentation( "IsRightCoactionObjectRep",
                       IsRightCoactionObject and IsAttributeStoringRep,
                       [ ] );

BindGlobal( "TheFamilyOfRightCoactionObjects",
            NewFamily( "TheFamilyOfRightCoactionObjects" ) );

BindGlobal( "TheTypeOfRightCoactionObjects",
            NewType( TheFamilyOfRightCoactionObjects,
                     IsRightCoactionObjectRep ) );

#############################
##
## Constructors
##
#############################

##
InstallMethod( LeftCoactionObject,
               [ IsCapCategoryMorphism, IsCapCategoryObject ],

  function( structure_morphism, coacting_object )

    return LeftCoactionObject( structure_morphism, LeftCoactionsCategory( coacting_object ) );

end );

##
InstallMethod( RightCoactionObject,
               [ IsCapCategoryMorphism, IsCapCategoryObject ],

  function( structure_morphism, coacting_object )

    return RightCoactionObject( structure_morphism, RightCoactionsCategory( coacting_object ) );

end );

max_stars_repo_name: ajnavarro/language-dataset

max_issues_repo_name: ajnavarro/language-dataset

max_forks_repo_name: ajnavarro/language-dataset


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

Declare(TPrmMulti,CodeBlock);
IsTwoPower := i >> 2 ^ Log2Int(i) = i;

get_l_power := function(t,exp)
  local p;
  if (exp=0) then
    p := L(2^t,1);
  else
    p:= L(2^t,2^(t-exp));
  fi;
  return p;
end;

# Declaration of SortBase, a 2x2 sorter.
Class(SortBase, BaseMat, rec(
   abbrevs   := [()-> []],
   new       := (self) >> SPL( WithBases(self, rec()) ).setDims(),
   dims      := self >> [ 2, 2 ],
   isReal    := True,
   sums     := self >> self,
   rChildren := self >> [],
   rSetChild := rSetChildFields(),
   toAMat := self >> Error("not supported"),
   transpose := self >> self,
));

# Declaration of SortBase, a 1x1 sorter.
Class(SortBase_w1, BaseMat, rec(
   abbrevs   := [(a)-> [a]],
   new       := (self, a) >> SPL( WithBases(self, rec(dimensions:=[1,1], a := a))),
   print := (self, i, is) >> Print(self.name, "(", self.a, ")"),
   dims      := self >> [ 1, 1 ],
   isReal    := True,
   sums     := self >> self,
   rChildren := self >> [self.a],
   rSetChild := rSetChildFields("a"),
   toAMat := self >> Error("not supported"),
   transpose := self >> self,
));

Class(SortConfigBase_w1, BaseMat, rec(
   abbrevs   := [(a)-> [a] , (b)-> [b]],
   new       := (self, a, b) >> SPL( WithBases(self, rec(dimensions:=[1,1], a := a, b:=b))),
   print := (self, i, is) >> Print(self.name, "(", self.a, ",", self.b, ")"),
   dims      := self >> [ 1, 1 ],
   isReal    := True,
   sums     := self >> self,
   rChildren := self >> [self.a, self.b],
   rSetChild := rSetChildFields("a","b"),
   toAMat := self >> Error("not supported"),
   transpose := self >> self,
));

# Declaration of SortConfigBase, a 2x2 Configurable sorter.
Class(SortConfigBase, BaseMat, rec(
   abbrevs   := [(a)-> [a]],
   new       := (self, a) >> SPL( WithBases(self, rec(dimensions:=[2,2], a := a))),
   print := (self, i, is) >> Print(self.name, "(", self.a, ")"),
   dims      := self >> [ 2, 2 ],
   isReal    := True,
   sums     := self >> self,
   rChildren := self >> [self.a],
   rSetChild := rSetChildFields("a"),
   toAMat := self >> Error("not supported"),
   transpose := self >> self,
));

# Declaration of SortLinearBase, used for linear sorter with w = 1
Class(LinearSortBase, BaseMat, rec(
   abbrevs   := [()-> []],
   new       := (self) >> SPL( WithBases(self, rec()) ).setDims(),
   dims      := self >> [ 2, 2 ],
   isReal    := True,
   sums     := self >> self,
   rChildren := self >> [],
   rSetChild := rSetChildFields(),
   toAMat := self >> Error("not supported"),
   transpose := self >> self,
));

#taken from: .../spiral/compiler/dag.gi
regassign.op_in  := self >> ConcatList(self.loc.rChildren(), ArgsExp) :: ArgsExp(self.exp);
regassign.op_out := self >> [self.loc];
regassign.op_inout := self >> [];
regassign.getNoScalar := self >> When(IsBound(self.exp.getNoScalar), self.exp.getNoScalar(), []);

# This tells Spiral how to translate 'SortBase' into code.
HDLCodegen.SortBase := (self, o, y, x, opts) >>
    chain(
	assign(nth(y,0), cond(leq(nth(x,0), nth(x,1)), nth(x,0), nth(x,1))),
	assign(nth(y,1), cond(leq(nth(x,0), nth(x,1)), nth(x,1), nth(x,0)))
    );

HDLCodegen.SortBase_w1 := (self, o, y, x, opts) >>
   let(
     t0 := TempVar(x.t.t),
     t1 := TempVar(x.t.t),
     t3 := TempVar(x.t.t),
     t5 := TempVar(x.t.t),
     t6 := TempVar(x.t.t),

     chain(
	assign(t3, imod(o.a, 2)),
	regassign(t0, cond(t3, t0, nth(x,0))),
	assign(t5, cond(leq(t0, nth(x,0)), nth(x,0), t0)),
	assign(t6, cond(leq(t0, nth(x,0)), t0, nth(x,0))),
	regassign(t1, cond(t3, t5, t1)),
	assign(nth(y,0), cond(t3, t6, t1))
     )
    );

HDLCodegen.SortConfigBase_w1 := (self, o, y, x, opts) >>
   let(
     t0 := TempVar(x.t.t),
     t1 := TempVar(x.t.t),
     t3 := TempVar(x.t.t),
     t5 := TempVar(x.t.t),
     t6 := TempVar(x.t.t),
     t7 := TempVar(x.t.t),
     t8 := TempVar(x.t.t),
     t9 := TempVar(x.t.t),
     t10 := TempVar(x.t.t),
     t2 := TempVar(x.t.t),

     chain(
	assign(t3, imod(o.a, 2)),
	regassign(t0, cond(t3, t0, nth(x,0))),
	assign(t7, eq(o.b,0)),
	assign(t8, eq(o.b,1)),

	assign(t2, leq(t0, nth(x,0))),
	assign(t9, cond(t2, t0, nth(x,0))),
	assign(t10, cond(t2, nth(x,0), t0)),
	assign(t5, cond(t7, nth(x,0) , t8, t9, t10)),
	assign(t6, cond(t7, t0 , t8, t10, t9)),

	regassign(t1, cond(t3, t5, t1)),
	assign(nth(y,0), cond(t3, t6, t1))
     )
    );

HDLCodegen.LinearSortBase := (self, o, y, x, opts) >>
   let(
     #x[1] indicates the start of a new list. When x[1] equals 1, y[0]
     #is set to the value of the register t0, and t0 is set to the value of x[0]

     t0 := TempVar(x.t.t),
     t1 := TempVar(x.t.t),
     t2 := TempVar(x.t.t),
     t3 := TempVar(x.t.t),
     t4 := TempVar(x.t.t),
     t5 := TempVar(x.t.t),
     r0 := TempVar(x.t.t),
     x0 := TempVar(x.t.t),
     x1 := TempVar(x.t.t),

     chain(

	 assign(x0, nth(x,0)),
	 assign(x1, nth(x,1)),
	 assign(t3, eq(x1, 1)),
	 regassign(r0, cond(t3, x0, t1)),
	 #regassign(r0, t0),
	 assign(t4, leq(r0, x0)),
	 assign(t1, cond(t4, x0, r0)),
	 assign(t2, cond(t4, r0, x0)),
	 assign(nth(y, 0), cond(t3, r0, t2)),
	 #regassign(t5, x1),
	 assign(nth(y, 1), x1)
     )
    );

# HDLCodegen.LinearSortBase := (self, o, y, x, opts) >>
#    let(
#      t0 := TempVar(x.t.t),
#      #x[1] indicates the start of a new list. When x[1] equals 1, y[0]
#      #is set to the value of the register t0, and t0 is set to the value of x[0]
#      chain(
# 	assign(nth(y,0),cond(eq(nth(x,1),1),t0,leq(t0, nth(x,0)),t0,nth(x,0))),
# 	regassign(t0,cond(eq(nth(x,1),1),nth(x,0),leq(t0, nth(x,0)),nth(x,0),t0)),
# 	assign(nth(y,1),nth(x,1))
#      )
#     );

HDLCodegen.SortConfigBase := (self, o, y, x, opts) >>
    let(
	t0 := TempVar(x.t.t),
	t1 := TempVar(x.t.t),
	t2 := TempVar(x.t.t),
	t3 := TempVar(x.t.t),
	chain(
	    assign(t2, nth(x,0)),
	    assign(t3, nth(x,1)),
	    assign(t0, cond(leq(t2, t3), t2, t3)),
	    assign(t1, cond(leq(t2, t3), t3, t2)),
	    assign(nth(y,0), cond(eq(o.a,0), t2, eq(o.a,1), t1, t0)),
	    assign(nth(y,1), cond(eq(o.a,0), t3, eq(o.a,1), t0, t1))
	)
    );

# Declaration of a sorter of size n
Class(Sort, TaggedNonTerminal, rec(
    abbrevs := [
    (n)       -> Checked(IsPosIntSym(n), [_unwrap(n)]),
    ],

    hashAs := self >> ObjId(self)(self.params[1]).withTags(self.getTags()),

    dims := self >> [ self.params[1], self.params[1] ],

    terminate := self >> Error("not supported"), # we could probably support this
));

# SortIJPerm(n): DirectSum(I(n/2), J(n/2))
# We do this so we can define a .permBits() function.
# This will let us easily generate large-size hardware implementations.

# (If we do not do this, the perm tool will have to manually compute the
# bit matrix representation, which will include making an n-times-n matrix.
Class(SortIJPerm, PermClass, rec(
    def := (n) -> Checked(
        IsPosIntSym(n),
        rec(size := n)),

    lambda := self >> let(
        n := self.params[1],
	fDirsum(fId(n/2), J(n/2)).lambda()
    ),

    transpose := self >> self,
    isSymmetric := self >> true,

    permBits := meth(self)
        local n, k, a, b, tmp, i;
	n := self.params[1];
	k := LogInt(n, 2);
	a := [List([1..k], i->0)];
	for i in [1 .. k-1] do
	    tmp := Concatenation([1], List([1..k-1], i->0));
	    Append(a, [tmp]);
        od;
	a := a * GF(2).one;
	b := MatSPL(I(k))*GF(2).one;

	return (a+b);
    end,
));

# Rule to breakdown Sort(n), where n is a power of two.
NewRulesFor(Sort, rec(
    Sort_Stream := rec(
        info         := "Streaming sorting network",

        applicable   := nt -> Length(nt.params) = 1 and IsTwoPower(nt.params[1]),

        children := (self, nt) >> let(

	    tag_w_tmp := nt.tags[1],
	    tag_w := tag_w_tmp.bs,
	    t := Log2Int(nt.params[1]),
	    p := Ind(2^t),
	    get_bb := w -> Cond(w=1, TTensorInd(SortBase_w1(p), p, APar, APar),TTensorI(SortBase(), 2^(t-1), APar, APar)),

	    [[ TCompose(
	           [TCompose(List([1..t-1], i ->
	               TCompose([
		           #TTensorI(SortBase(), 2^(t-1), APar, APar),
		           get_bb(tag_w),
		           TCompose(List([2..(t-i+1)], j ->
		               TCompose([
			           TTensorI(TPrm(Tensor(I(2), L(2^(j-1), 2^(j-2))) * L(2^j,2)), 2^(t-j), APar, APar),
			           get_bb(tag_w)
			           #TTensorI(SortBase(), 2^(t-1), APar, APar)
			       ])
		           )),
		        #   TTensorI(TPrm(L(2^(t-i+1), 2^(t-i)) * DirectSum(I(2^(t-i)), J(2^(t-i)))  ), 2^(i-1), APar, APar)
			   TTensorI(TPrm(L(2^(t-i+1), 2^(t-i)) * SortIJPerm(2^(t-i+1))), 2^(i-1), APar, APar)
		       ])
	            )),
		    get_bb(tag_w)]
		    #TTensorI(SortBase(), 2^(t-1), APar, APar)]
                ).withTags(nt.getTags())
            ]]
	),

        apply        := (nt, c, cnt) -> c[1],

    ),

   Sort_Stream6 := rec(
        info         := "Version of Sort_Stream (sortAlg1) removing J permutations and with configurable 2-input sorters",

        applicable   := nt -> Length(nt.params) = 1 and IsTwoPower(nt.params[1]),

        children := (self, nt) >> let(
	    t := Log2Int(nt.params[1]),
            k := Ind(2^(t-1)),
            z := (s) >> t-s,

            c1 := (s) >> logic_and(eq(bit_sel(k, z(s)), 1), neq(s, 1)),
            access_f := (s) >> cond(c1(s), 1, 2),

	    [[ TCompose(
	           [TCompose(List([1..t-1], i ->
	               TCompose([
			   TTensorInd(SortConfigBase(access_f(i)), k, APar, APar),
		           TCompose(List([2..(t-i+1)], j ->
		               TCompose([
			           TTensorI(TPrm(Tensor(I(2), L(2^(j-1), 2^(j-2))) * L(2^j,2)), 2^(t-j), APar, APar),
			           TTensorInd(SortConfigBase(access_f(i)), k, APar, APar)
			       ])
		           )),
		           TTensorI(TPrm(L(2^(t-i+1), 2^(t-i))), 2^(i-1), APar, APar)
		       ])
	            )),
		    TTensorInd(SortConfigBase(access_f(t)), k, APar, APar)]
		    #TTensorI(SortBase(), 2^(t-1), APar, APar)] #what is this line?
                ).withTags(nt.getTags())
            ]]
	),

        apply        := (nt, c, cnt) -> c[1],

    ),

   Sort_Stream_Iter := rec(
        info         := "Stream/Iter sorting network",

	depth_params := [],

        applicable   := (self, nt) >> Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and
	                              Length(self.depth_params) = Log2Int(nt.params[1])-1 and
				      let(
					  t := Log2Int(nt.params[1]),
					  d := self.depth_params,
				      ForAll(List([1..t-1], i -> IsInt((t-i+1)/d[i])), j -> j)
				      ),

        children := (self, nt) >> let(
   	    t := Log2Int(nt.params[1]),

	    tag_w_tmp := nt.tags[1],
	    tag_w := tag_w_tmp.bs,
	    p := Ind(2^t),

	    get_bb := w -> Cond(w=1, TTensorInd(SortBase_w1(p), p, APar, APar),TTensorI(SortBase(), 2^(t-1), APar, APar)),

	    stage := i >> TCompose([
			     get_bb(tag_w),
			     #TTensorI(SortBase(), 2^(t-1), APar, APar),
			     TTensorI(TPrm(L(2^(t-i+1), 2^(t-i))) , 2^(i-1), APar, APar)
			  ]),

	    full_stage := i >> TCompose(List([1..self.depth_params[i]], j2 -> stage(i))),

  	    [[ TCompose(
	           [TCompose(List([1..t-1], i -> let(
		       j := Ind(t-i+1),
		       d := self.depth_params[i],
		       j1 := Ind(d),
		       j2 := Ind((t-i+1)/d),
	               TCompose([
			   Cond(d = ((t-i+1)),
                               full_stage(i),
				TCompose(List([1..d],j1 ->
				TICompose(j2, (t-i+1)/d, stage(i))
				))
			   ),
		       	   #TTensorI(TPrm(DirectSum(I(2^(t-i)), J(2^(t-i)))), 2^(i-1), APar, APar)
		       	   TTensorI(TPrm(SortIJPerm(2^(t-i+1))), 2^(i-1), APar, APar)
		       ])
	            ))),
		    get_bb(tag_w)] #last iterations left outside of TCompose
		    #TTensorI(SortBase(), 2^(t-1), APar, APar)] #last iterations left outside of TCompose
               ).withTags(nt.getTags())
            ]]
	),

        apply        := (nt, c, cnt) -> c[1],

   ),

   #old version -works-
#      Sort_Stream_Iter := rec(
#        info         := "Stream/Iter sorting network",
#
#	depth_params := [],
#
#        applicable   := (self, nt) >> Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and
#	                              Length(self.depth_params) = Log2Int(nt.params[1])-1 and
#				      let(
#					  t := Log2Int(nt.params[1]),
#					  d := self.depth_params,
#				      ForAll(List([1..t-1], i -> IsInt((t-i+1)/d[i])), j -> j)
#				      ),
#
#
#        children := (self, nt) >> let(
#   	    t := Log2Int(nt.params[1]),
#
#	    stage := i >> TCompose([
#			     TTensorI(SortBase(), 2^(t-1), APar, APar),
#			     TTensorI(TPrm(L(2^(t-i+1), 2^(t-i))) , 2^(i-1), APar, APar)
#			  ]),
#
#	    full_stage := i >> TCompose(List([1..self.depth_params[i]], j2 -> stage(i))),
#
#
#  	    [[ TCompose(
#	           [TCompose(List([1..t-1], i -> let(
#		       j := Ind(t-i+1),
#		       d := self.depth_params[i],
#		       j1 := Ind((t-i+1)/d),
#	               TCompose([
#			   Cond(j1.range = 1,
#			       full_stage(i),
#			       TICompose(j1, (t-i+1)/d, full_stage(i))
#			   ),
#		       	   TTensorI(TPrm(DirectSum(I(2^(t-i)), J(2^(t-i)))), 2^(i-1), APar, APar)
#		       ])
#	            ))),
#		    TTensorI(SortBase(), 2^(t-1), APar, APar)]
#               ).withTags(nt.getTags())
#            ]]
#	),
#
#        apply        := (nt, c, cnt) -> c[1],
#
#   ),

   Sort_Stream5 := rec(
        info         := "Version of Sort_Stream_Iter (sortAlg2) removing J permutations and with configurable 2-input sorters",

	depth_params := [],

        applicable   := (self, nt) >> Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and
	                              Length(self.depth_params) = Log2Int(nt.params[1])-1 and
				      let(
					  t := Log2Int(nt.params[1]),
					  d := self.depth_params,
				      ForAll(List([1..t-1], i -> IsInt((t-i+1)/d[i])), j -> j)
				      ),

        children := (self, nt) >> let(
   	    t := Log2Int(nt.params[1]),
	    k := Ind(2^(t-1)),
	    z := (s) >> t-s,

	    c1 := (s) >> logic_and(eq(bit_sel(k, z(s)), 1), neq(s, 1)),
	    access_f := (s) >> cond(c1(s), 1, 2),

	    stage := i >> TCompose([
			     TTensorInd(SortConfigBase(access_f(i)), k, APar, APar),
			     TTensorI(TPrm(L(2^(t-i+1), 2^(t-i))) , 2^(i-1), APar, APar)
			  ]),

	    full_stage := i >> TCompose(List([1..self.depth_params[i]], j2 -> stage(i))),

  	    [[ TCompose(
	           [TCompose(List([1..t-1], i -> let(
		       j := Ind(t-i+1),
		       d := self.depth_params[i],
		       j1 := Ind(d),
		       j2 := Ind((t-i+1)/d),
		       Cond(d = ((t-i+1)),
                                full_stage(i),
				TCompose(List([1..d],j1 ->
				TICompose(j2, (t-i+1)/d, stage(i))
				))
			       #unroll this change to TCompose
			)
	            ))),
		    TTensorInd(SortConfigBase(access_f(t)), k, APar, APar)] #last iterations left outside of TCompose
               ).withTags(nt.getTags())
            ]]
	),

        apply        := (nt, c, cnt) -> c[1],

   ),

  Sort_Stream3 := rec(
        info         := "",

	depth_out := 1,
        depth_in := [],

        applicable   := (self, nt) >> Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and
			     	      Length(self.depth_in) = Log2Int(nt.params[1])-1  and
				      let(
                                          t := Log2Int(nt.params[1]), d_out := self.depth_out,
                                          d_in := self.depth_in, IsInt(t/d_out) and ((d_out=1) or ((d_out=t) and
                                      	  (ForAll(List([1..t-1], i -> IsInt((t-i+1)/d_in[i])), j -> j))))
                                      ),

        children := (self, nt) >> let(
	    d_out := self.depth_out,
	    d_in := self.depth_in,
            t := Log2Int(nt.params[1]),
	    k := Ind(2^(t-1)),
	    iter_d1 := Ind(t),

            z := (lp, jp) >> (t-1)-(lp+jp),
            c2 := (lp, jp) >> logic_and(eq(bit_sel(k, z(lp, jp)), 1), neq(lp, 0)),
            access_f := (lp, jp) >> cond(c2(lp, jp), 1, 2),
            stage := (lp, jp) >> TCompose([
                       TTensorInd(SortConfigBase(access_f(lp, jp)), k, APar, APar),
                       TPrm(L(2^t, 2^(t-1)))
                   ]),
	    full_stage := i >> TCompose(List([0..d_in[i+1]-1], j -> stage(i,j))),
	    [[ Cond(d_out=1,
			TICompose(iter_d1,t,let(
				iter_d2 := Ind(t-iter_d1),
	    			p_list := List([0..t-1], i-> get_l_power(t,i)),
				TCompose([
					TICompose(iter_d2,t-iter_d1, stage(iter_d1,iter_d2)),
		  	 		TPrmMulti(p_list,iter_d1)
				])
			)).withTags(nt.getTags()),
		    d_out=t,
			TCompose([
				TCompose(List([0..t-2], iter_d3 -> let(
				  d_stage := d_in[iter_d3+1],
            			  iter_d5 := Ind((t-iter_d3)/d_stage),
				  TCompose([
					Cond(d_stage=(t-iter_d3),
						full_stage(iter_d3),
						TCompose(List([0..d_stage-1],iter_d4 ->
						  TICompose(iter_d5,(t-iter_d3)/d_stage,stage(iter_d3,iter_d4*iter_d5))))
					),
					get_l_power(t,iter_d3)
				  ])
				))),
				stage(t-1,0),
				get_l_power(t,t-1)
			]).withTags(nt.getTags())
	    )]]

	),
        apply        := (nt, c, cnt) -> c[1],
   ),

   Sort_Stream4 := rec(
        info         := "",

	depth := 1,

        applicable   := (self, nt) >> Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and
	                              let (d := self.depth, t := Log2Int(nt.params[1]), IsInt(t*t/self.depth) and
				      (IsInt(d/t) or IsInt(t/d))),

        children := (self, nt) >> let(
   	    t := Log2Int(nt.params[1]),
	    d := self.depth,
	    d1 := cond(leq(d,t), d, t).ev(),
	    d2 := cond(leq(d,t), 1, d/t).ev(),

	    k := Ind(2^(t-1)),
	    j := Ind(t),
	    l := Ind(t),
	    n := Ind(t/d2),
	    v := Ind(t/d1),

            d_tmp := cond(leq(d,t), t/d, 0).ev(),
            d2_tmp := cond(leq(t,d), d/t, 0).ev(),
	    s_tmp := ((t*t)/d),
            m2 := Ind(d_tmp),
	    s2 := Ind(s_tmp),
	    v2 := Ind(d2_tmp),
	    n2 := Ind(d),
	    n3 := Ind(d),
	    l2 := Ind(t),

	    tag_w_tmp := nt.tags[1],
	    tag_w := tag_w_tmp.bs,
	    p := Ind(2^t),

	    c1 := (lp, jp) >> lt((t-1), (lp+jp)),
	    z := (lp, jp) >> (t-1)-(lp+jp),
	    z_w1 := (lp, jp) >> (t-1)-(lp+jp)+1,
	    c2 := (lp, jp) >> logic_and(eq(bit_sel(k, z(lp, jp)), 1), neq(lp, 0)),
	    c2_w1 := (lp, jp) >> logic_and(eq(bit_sel(p, z_w1(lp, jp)), 1), neq(lp, 0)),

	    access_f := (lp, jp) >> cond(c1(lp, jp), 0, c2(lp, jp), 1, 2),
	    access_f_w1 := (lp, jp) >> cond(c1(lp, jp), 0, c2_w1(lp, jp), 1, 2),

	    get_bb := (lp,jp) -> Cond(tag_w=1, TTensorInd(SortConfigBase_w1(p,access_f_w1(lp,jp)), p, APar, APar),TTensorInd(SortConfigBase(access_f(lp, jp)), k, APar, APar)),

	    stage := (lp, jp) >> TCompose([
		       get_bb(lp, jp),
		       #TTensorInd(SortConfigBase(access_f(lp, jp)), k, APar, APar),
		       TPrm(L(2^t, 2^(t-1)))
	           ]),

	    full_stage := np_1 >> TCompose(List([0..t-1], m_1 -> TCompose(List([0..t-1], j_1 -> stage(m_1, j_1))))),
	    full_stage1 := np >> TCompose(List([0..d2-1], m -> TCompose(List([0..t-1], j -> stage(d2*np+m, j))))),
	    full_stage2 := vp >> TCompose(List([0..d1-1], s -> stage(l, vp+s))),
            full_stage1b := np2 >> TICompose(m2,d_tmp, TICompose(j,t, stage((t/d)*np2+m2, j))),

             # Old: problem is that it's assuming the l2 above, which is an unassigned iterator.
	     # There is also a problem with the vp2+s2 parameter: you need to multiply vp2 by the number of iterations.
             # full_stage2b := (vp2) >> TICompose(s2,s_tmp, stage(l2, vp2+s2)),
            full_stage2b := (vp2, l3) >> TICompose(s2,s_tmp, stage(l3, vp2*s_tmp+s2)),

  	    [[ Cond(d=t*t, full_stage(0).withTags(nt.getTags()),
		    d<t, TCompose(List([0..d-1], n2 -> full_stage1b(n2))).withTags(nt.getTags()),
                    d=t, TCompose(List([0..d-1], n3 -> full_stage1b(n3))).withTags(nt.getTags()),
                    d>t, TCompose(List([0..t-1], l3 -> TCompose(List([0..d2_tmp-1], v2 -> full_stage2b(v2, l3))))).withTags(nt.getTags())) #the problem seems to be the outer most TCompose works if it was TICompose
		    #d>t, TICompose(n, t/d2, full_stage1(n)).withTags(nt.getTags())) #old one but will leave it as new one does not work yet
	        ]]
	),

#	Old version -works-
#	   Sort_Stream4 := rec(
#        info         := "",
#
#	depth := 1,
#
#        applicable   := (self, nt) >> Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and
#	                              let (d := self.depth, t := Log2Int(nt.params[1]), IsInt(t*t/self.depth) and
#				      (IsInt(d/t) or IsInt(t/d))),
#
#
#        children := (self, nt) >> let(
#   	    t := Log2Int(nt.params[1]),
#	    d := self.depth,
#	    d1 := cond(leq(d,t), d, t).ev(),
#	    d2 := cond(leq(d,t), 1, d/t).ev(),
#
#	    k := Ind(2^(t-1)),
#	    j := Ind(t),
#	    l := Ind(t),
#	    n := Ind(t/d2),
#	    v := Ind(t/d1),
#
#
#	    c1 := (lp, jp) >> lt((t-1), (lp+jp)),
#	    z := (lp, jp) >> (t-1)-(lp+jp),
#	    c2 := (lp, jp) >> logic_and(eq(bit_sel(k, z(lp, jp)), 1), neq(lp, 0)),
#
#	    access_f := (lp, jp) >> cond(c1(lp, jp), 0, c2(lp, jp), 1, 2),
#
#	    stage := (lp, jp) >> TCompose([
#		       TTensorInd(SortConfigBase(access_f(lp, jp)), k, APar, APar),
#		       TPrm(L(2^t, 2^(t-1)))
#	           ]),
#
#	    full_stage1 := np >> TCompose(List([0..d2-1], m -> TCompose(List([0..t-1], j -> stage(d2*np+m, j))))),
#	    full_stage2 := vp >> TCompose(List([0..d1-1], s -> stage(l, vp+s))),
#
#  	    [[ Cond(d=t*t, full_stage1(0).withTags(nt.getTags()),
#		    d>t, TICompose(n, t/d2, full_stage1(n)).withTags(nt.getTags()),
#		    d<t, TICompose(l, t, TICompose(v, t/d1, full_stage2(d1*v))).withTags(nt.getTags()),
#		    d=t, TICompose(l, t, full_stage2(0)).withTags(nt.getTags()))
#	        ]]
#	),
#
#        apply        := (nt, c, cnt) -> c[1],
#   ),

	# this works, but only for d=1.
#         children := (self, nt) >> let(
#    	    t := Log2Int(nt.params[1]),
# 	    k := Ind(2^(t-1)),
# 	    l := Ind(t),
# 	    j := Ind(t),

# 	    c1 := lt((t-1), (l+j)),
# 	    z := (t-1)-(l+j),
# 	    c2 := logic_and(eq(bit_sel(k, z), 1), neq(l, 0)),

# 	    access_f := cond(c1, 0, c2, 1, 2),

# 	    d1 := cond(leq(d,t), d, t),
# 	    d2 := cond(leq(d,t), 1, d/t),

# 	    c1 := lt((t-1), (l+j)),
# 	    z := (t-1)-(l+j),
# 	    c2 := logic_and(eq(bit_sel(k, z), 1), neq(l, 0)),

# 	    access_f := cond(c1, 0, c2, 1, 2),

# 	    d1 := cond(leq(d,t), d, t),
# 	    d2 := cond(leq(d,t), 1, d/t),

#   	    [[ TICompose(l, t,
# 		   TICompose(j, t, TCompose([
# 		       TTensorInd(SortConfigBase(access_f), k, APar, APar),
# 		       TPrm(TL(2^t, 2^(t-1), 1, 1))
# 	           ]))).withTags(nt.getTags())
#             ]]
# 	),

        apply        := (nt, c, cnt) -> c[1],

   ),

  Linear_Sort := rec(
        info         := "",

        applicable   := (self, nt) >> IsTwoPower(nt.params[1]),

        children := (self, nt) >> let(
            t := Log2Int(nt.params[1]),
	    [[ TCompose(List([1..2^t], i ->
		TTensorI(LinearSortBase(), nt.params[1], APar, APar)
			)).withTags(nt.getTags()),
	    ]]

	),
        apply        := (nt, c, cnt) -> c[1],
   ),

));

#---------------------------------------------------------------------
#---------------------------------------------------------------------
#---------------------------------------------------------------------
# Old stuff, not used now

#     Sort_Stream_old := rec(
#         info         := "Streaming sorting network",

#         switch       := true,
#         applicable   := nt ->
#             Length(nt.params) = 1 and IsTwoPower(nt.params[1]),

#         children := (self, nt) >> let(
# 	    t := Log2Int(nt.params[1]),
# 	    [[ TCompose(List([1..t], i ->
# 	 	TCompose([
# 		    TTensorI(BitonicSort(2^(t-i+1)), 2^(i-1), APar, APar),
# 		    TTensorI(TPrm(DirectSum(I(2^(t-i)), J(2^(t-i)))), 2^(i-1), APar, APar)
# 		])
#             )).withTags(nt.getTags())]]
#         ),

#         apply        := (nt, c, cnt) -> c[1],

#     ),

# Declaration of a bitonic sorter of size n
# Class(BitonicSort, TaggedNonTerminal, rec(
#     abbrevs := [
#     (n)       -> Checked(IsPosIntSym(n), [_unwrap(n)]),
#     ],

#     hashAs := self >> ObjId(self)(self.params[1]).withTags(self.getTags()),

#     dims := self >> [ self.params[1], self.params[1] ],

#     terminate := self >> Error("not supported"), # we could probably support this
# ));

# Old stuff, not used now.
# # Rules to break-down BitonicSort()
# NewRulesFor(BitonicSort, rec(
#     BitonicSort_Stream := rec(
#         info         := "Streaming Bitonic sorting network",

#         switch       := true,
#         applicable   := nt ->
#             Length(nt.params) = 1 and IsTwoPower(nt.params[1]) and nt.params[1] > 2,

#         children := (self, nt) >> let(
# 	    k := Log2Int(nt.params[1]),
# 	    [[ TCompose([ TTensorI(SortBase(), 2^(k-1), APar, APar),
# 		  TCompose(List([2..k], j -> TCompose([
# 		      TTensorI(TPrm(Tensor(I(2), L(2^(j-1), 2^(j-2))) * L(2^j, 2)), 2^(k-j), APar, APar),
# 		      TTensorI(SortBase(), 2^(k-1), APar, APar)
# 		  ]))),
# 		  TPrm(L(2^k, 2^(k-1)))
# 	       ]).withTags(nt.getTags())
# 	    ]]
#         ),

#         apply        := (nt, c, cnt) -> c[1],
#     ),

#     BitonicSort_Base := rec(
# 	info := "Bitonic soring network base rule",
# 	applicable := nt -> nt.params[1] = 2,
# 	chldren := (self, nt) >> [[ ]],
# 	apply := (nt, c, cnt) -> SortBase()
#     )

# ));


text:
#
# Read("~/Workspace/groupsSB/epi/B2/generic.UU.gi");
#

type:="B";
rank:=2;
nr_pos_roots:=4;

Read("~/Workspace/groupsSB/epi/group.gi");

handleUaUbUc:=function()
	local rels, vals,i,v;
	# UaUb=Uc = U(a1+b1,a2+b2,a3+b3+a2*b1)
	rels:=Set(Concatenation(Ua*Ub-Uc));

	vals:=[];
	Append(vals,[[[vars[22]],[vars[2]+vars[12]]]]);
	Append(vals,[[[vars[21]],[vars[1]+vars[11]]]]);
	Append(vals,[[[vars[23]],[vars[3]+vars[13]-vars[2]*vars[11]]]]);
#a_{2}^2*b_{1}+2*a_{2}*b_{1}*b_{2}-2*a_{3}*b_{2}-a_{4}-b_{4}
	Append(vals,[[[vars[24]],[-vars[2]^2*vars[11]-2*vars[2]*vars[11]*vars[12]+2*vars[3]*vars[12]+vars[4]+vars[14]]]]);

	for i in [1..Length(rels)] do
		for v in vals do
			Print(rels[i],"\n");
			rels[i]:=One(APR)*Value(rels[i],v[1],v[2]);
		od;
	od;
	return [vals,rels];
end;

max_stars_repo_name: iuliansimion/groupsSB

max_issues_repo_name: iuliansimion/groupsSB

max_forks_repo_name: iuliansimion/groupsSB


text:
ZapGremlins := function(s)
  local upper, lower, c, i, n, t;
  upper := "ABCDEFGHIJKLMNOPQRSTUVWXYZ";
  lower := "abcdefghijklmnopqrstuvwxyz";
  t := [ ];
  i := 1;
  for c in s do
    n := Position(upper, c);
    if n <> fail then
      t[i] := lower[n];
      i := i + 1;
    else
      n := Position(lower, c);
      if n <> fail then
        t[i] := c;
        i := i + 1;
      fi;
    fi;
  od;
  return t;
end;

IsPalindrome := function(s)
  local t;
  t := ZapGremlins(s);
  return t = Reversed(t);
end;

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

Declare(requiredFreeVars);
Declare(buildMap);
Declare(collectLoopVars);
Declare(buildGlobalList);
Declare(compareMatrices);
Declare(doAllPromotes);
Declare(permToData);
Declare(PadLeft);
Declare(PadRight);
Declare(isContainedInComposeDistsG);
Declare(isContainedInComposeDistsS);

################################## CANDIDATES FOR REMOVAL -- BEGIN #########################3
Class(ProducOutputSDGR, RuleSet);
RewriteRules(ProducOutputSDGR, rec(

));

Class(RulesNoPull_Dist, RuleSet);
RewriteRules(RulesNoPull_Dist, rec(

    vtensor_nopulldist := Rule([@(1, VTensor), @(2, NoPull_Dist)],
        e->NoPull_Dist(VTensor(@(2).val.children(), @(1).val.vlen))
    )
));

#F Removes NoPull_Dist
Class(RemoveNoPull_Dist, RuleSet);
RewriteRules(RemoveNoPull_Dist, rec(
    remove_sizedkernel := Rule(@(1, NoPull_Dist), e->@(1).val.child(1))
));
################################## CANDIDATES FOR REMOVAL -- END #########################3

#F To convert DistSum(A)*DistSum(B)*... -> DistContainer(Seq(A,B,..))
Class(RulesComposeDists, RuleSet);
RewriteRules(RulesComposeDists, rec(
    dist_to_seq := Rule(@(1, Compose, e->ForAll(e.children(), c->(ObjId(c)=DistSum or ObjId(c)=Comm_Cell))),
        e -> ComposeDists(@(1).val.children()) )
));

#F DistSum|*|DistSum -> DistSum
Class(RulesDistMerge, RuleSet);
RewriteRules(RulesDistMerge, rec(
    #NOTE: If these DistSums have ScatSend or GathRecv, this rule will do the wrong thing!
    # Assuming this rule will only match DMP type DFT algos.
    dist_merge := ARule(ComposeDists, [@(1,DistSum), @(2,DistSum)],
      e -> let(a := @(1).val, b := @(2).val,
          [ DistSum(a.P, a.var, a.domain, a.child(1) * b.child(1)) ])
      ),

    # Above rule produces DistSum(S * child1 * G * S * child2 * G).
    # This rule gets rid of the G * S *

    gathdist_scat_dist := ARule(Compose, [GathDist, ScatDist],
      e -> [])
));

#F Removes Buf() so other rules can work!
Class(RemoveBuf, RuleSet);
RewriteRules(RemoveBuf, rec(
    remove_buf := Rule(@(1, Buf), e->@(1).val.child(1))
));

#F Convert initial and final GathRecv/ScatSend to Dist (null)
Class(CellDFTBlockCyclicLayoutHack, RuleSet);
RewriteRules(CellDFTBlockCyclicLayoutHack, rec(

    GathRecvDataLayout := Rule([@(1, GathRecv), @(2, [H, fTensor, fCompose, fId]), ...],
        e->GathDist(@1.val.func.range(), @1.val.pkSize, @1.val.P, @1.val.i)),

    ScatSendDataLayout := Rule([@(1, ScatSend), @(2, [H, fTensor, fCompose, fId]), ...],
        e->ScatDist(@1.val.func.range(), @1.val.pkSize, @1.val.P, @1.val.i))

));

#F Blockcyclic layout hack
CellDFTBlockCyclicLayoutHackWrap := function(sums, opts)
    if not (IsBound(opts.doNotUseBlockCyclic) and opts.doNotUseBlockCyclic = true) then
      return(CellDFTBlockCyclicLayoutHack(sums));
    else
      return(sums);
    fi;
end;

# Terminate VRCs (VRC(Scat/Gath) -> Scat/Gath(pkSize*2..)
# (Removed VRCL, VRCR from this mix Tue 09 Sep 2008 09:55:36 PM EDT)
Class(CellVRCTerm, RuleSet);
RewriteRules(CellVRCTerm, rec(
    VRC_ScatDist_Term := Rule([@(1, [RC, VRC,VRCLR]), @(2, [ScatDist])],
    e->ScatDist(@(2).val.N, @(2).val.pkSize*2, @(2).val.P, @(2).val.i)),

    VRC_GathDist_Term := Rule([@(1, [RC, VRC,VRCLR]), @(2, [GathDist])],
    e->GathDist(@(2).val.N, @(2).val.pkSize*2, @(2).val.P, @(2).val.i)),

    VRC_ScatSend_Term := Rule([@(1, [RC, VRC,VRCLR]), @(2, [ScatSend])],
    e->ScatSend(@(2).val.func, @(2).val.pkSize*2, @(2).val.P, @(2).val.i)),

    VRC_GathRecv_Term := Rule([@(1, [RC, VRC,VRCLR]), @(2, [GathRecv])],
    e->GathRecv(@(2).val.func, @(2).val.pkSize*2, @(2).val.P, @(2).val.i)),

    VRC_Comm_Cell := Rule([@(1,[RC, VRC,VRCLR]), @(2,Comm_Cell)],
    e->Comm_Cell(@(2).val.P, @(2).val.pkSize*2))

));

# VRC(DistSum.. -> DistSum(VRC...
RewriteRules(RulesVRC, rec(
    VRC_DistSum := Rule([@(1, [VRC,VRCL,VRCR,VRCLR]), @(2, [DistSum, DistSumLoop])],
       e->let(s := @(2).val,
       CopyFields(s, rec(_children := List(s.children(), c->ObjId(@(1).val)(c, @(1).val.v)),
                  dimensions := @(1).val.dimensions))
       ))
));

RewriteRules(RulesRC, rec(
    RC_DistSum := Rule([RC, @(1, DistSum)],
        e -> let(s:=@(1).val, DistSum(s.P, s.var, s.domain, RC(s.child(1))))),
));

# Pull Diag into DistSum (but not into Gath/Scat)
RewriteRules(RulesDiag, rec(
    CellPullInDiagRight := ARule(Compose,  [ @(1, [RCDiag, Diag, Prm, Scat]), @(2, DistSum) ],
     e -> let(s:=@(2).val, [ DistSum(s.P, s.var, s.domain, @(1).val * s.child(1)) ])),

    CellPullInDiagLeft := ARule(Compose, [ @(1, DistSum), @(2, [Prm, Gath, Diag, RCDiag]) ],
     e -> let(s:=@(1).val, [ DistSum(s.P, s.var, s.domain, s.child(1) * @(2).val) ])),

));

# Pull Diag into SAG (D*SAG -> SDAG and SAG*D -> SADG)
RewriteRules(RulesDiagStandalone, rec(
 # Gath * Diag
 CellCommuteGathDiag := ARule( Compose,
       [ @(1, [ GathDist, GathRecv ]), @(2, Diag) ], # o 1-> 2->
  e -> [ Diag(fCompose(@2.val.element, fTensor(@1.val.func, fId(@1.val.pkSize)))).attrs(@(2).val), @1.val ]),

 # Diag * Scat
 CellCommuteDiagScat := ARule( Compose,
       [ @(1, Diag), @(2, [ScatDist, ScatSend]) ], # <-1 <-2 o
  e -> [ @2.val, Diag(fCompose(@1.val.element, fTensor(@2.val.func, fId(@2.val.pkSize))  )).attrs(@(1).val) ]),

 # Gath * RCDiag
 CellCommuteGathRCDiag := ARule( Compose,
       [ [@(1, [ GathDist, GathRecv ]), [@(0,fTensor), ..., [fId,@(2).cond(IsEvenInt)]]],
      @(4, RCDiag) ],
  e -> [ RCDiag(fCompose(@(4).val.element, @(0).val), @(4).val.post),
         @(1).val ]),

 # RCDiag * Scat
 CellCommuteRCDiagScat := ARule( Compose,
       [ @(4, RCDiag),
     [@(1, [Scat, ScatDist, ScatSend]), [@(0,fTensor), ..., [fId,@(2).cond(IsEvenInt)]]] ],
  e -> [ @(1).val,
         RCDiag(fCompose(@(4).val.element, @(0).val), @(4).val.post) ]),

));

# Cell DMP: Fuse and remove PTensors that are next to Sigmas
Class(PTensorRules, RuleSet);
RewriteRules(PTensorRules, rec(
PTensorFuseRight := ARule(Compose, [ @(1,DistSum), @(2,PTensor) ],
    e -> [ DistSum( @(1).val.P, @(1).val.var, @(1).val.domain, Compose( @(1).val.child(1), @(2).val ) ) ] ),

PTensorFuseLeft := ARule(Compose, [  @(2,PTensor), @(1,DistSum) ],
    e -> [ DistSum( @(1).val.P, @(1).val.var, @(1).val.domain, Compose( @(2).val,  @(1).val.child(1) ) ) ] ),

PTensorFlipRight := ARule(Compose, [ @(1,GathDist), @(2,PTensor) ],
    e -> [ @(2).val.L, @(1).val ] ),

PTensorFlipLeft := ARule(Compose, [ @(2,PTensor), @(1,ScatDist) ],
    e -> [ @(1).val, @(2).val.L ] )
));

#i := Ind(P),
# Cell DMP: Convert standalone PTensors to Sigmas
Class(PTensorConvertRules, RuleSet);
RewriteRules(PTensorConvertRules, rec(
PTensorConvert := Rule( @(1,PTensor),
    e -> let(P := @(1).val.P,
             i := var("spuid", TInt, P),
             N := @(1).val.dims()[1],
            DistSum(P, i, P, ScatDist(N, 1, P, i) * @(1).val.L * GathDist(N, 1, P, i))
          )
    )
));

#F Gets all the DistSums under a single Compose so we can deal with them effectively
#F Mainly there because of VContainer
Class(DistSumChains, RuleSet);
RewriteRules(DistSumChains, rec(
    ChainDistSumsLeft := ARule(Compose,  [ @(1,DistSum), @(2,VContainer, e->ForAny(e.rChildren()[1].rChildren(), i->ObjId(i)=DistSum)) ],
     e -> let(ds:=@(1).val, vc:=@(2).val, [ VContainer(Compose(ds, vc.rChildren()), vc.isa)  ])),

    ChainDistSumsRight := ARule(Compose,  [ @(2,VContainer,  e->ForAny(e.rChildren()[1].rChildren(), i->ObjId(i)=DistSum)), @(1,DistSum) ],
     e -> let(ds:=@(1).val, vc:=@(2).val, [ VContainer(Compose(vc.rChildren(), ds), vc.isa)  ]))
));

#F Fuse output (SDGR*SS...) to just an SS
Class(FixBorder, RuleSet);
RewriteRules(FixBorder, rec(
    FixLeftBorder := Rule( [@(1,ComposeDists), [DistSum, ..., [Compose, ScatDist, GathRecv]], ...],
     e ->  let(c := @(1).val.rChildren(), ComposeDists(ListWithout(c,1))) )
));

#F Mark Cell Scat/Gaths as Inplace (ultimately, no-ops)
#F GDs, SDs, and GRs are no-ops for the Cell (SS becomes explicit DMA - SCATSEND)
applyCellInplace := function(sums, opts)
    SubstBottomUp(sums, GathDist, e->Inplace(e));
    SubstBottomUp(sums, GathRecv, e->Inplace(e));
    SubstBottomUp(sums, ScatDist, e->Inplace(e));

    # BB(Inplace()) -> Inplace() (So that a copy operation need not take place)
    SubstBottomUp(sums, [BB, Inplace], e->e.rChildren()[1]);
    return(sums);
end;

#NOTE: In the following padding functions, we previously produced pads with
#maximum packet size so that they merged nicely with multibuffered scats and
#gaths. But now, we don't do any merging between scats and gaths of the
#distributed and multibuffer paradigms. So it's better to get these padding
#functions to produce packet sizes that match the adjacent sigmaspl expressions.

#F This produces an SDGR on the left
PadLeft := function(e)
    local ss, N, pkSize, P,i, GR, SD;
    ss      := e.leftMostParScat();
    N       := ss.func.range();
    pkSize  := ss.pkSize;
    P       := ss.P;
    i       := var("spuid", TInt, P);
    GR      := GathRecv(fTensor(fBase(P, i), fId(N/P)), pkSize, P, i);
    SD      := ScatDist(N, pkSize, P, i);
    #SD     := ScatDist(P, (N*pkSize)/P, P, i); # Produce ScatDist with maximum packet size
    return(ComposeDists(DistSum(P, i, P, SD*GR),   e.rChildren()));
end;

#F This produces an SDGR on the left
PadLeftSep := function(e)
    local ss, N, pkSize, P,i, GR, SD;
    ss      := e.leftMostParScat();
    N       := ss.func.range();
    pkSize  := ss.pkSize;
    P       := ss.P;
    i       := var("spuid", TInt, P);
    GR      := GathRecv(fTensor(fBase(P, i), fId(N/P)), pkSize, P, i);
    SD      := ScatDist(N, pkSize, P, i);
    #SD     := ScatDist(P, (N*pkSize)/P, P, i); # Produce ScatDist with maximum packet size
    return(ComposeDists(DistSum(P, i, P, SD*GR), e));
end;

#F This produces a SSGD on the right
PadRight := function(e)
   local gr, N, P, i, pkSize, SS, GD, di;
   gr     := e.rightMostParGath();
   N      := gr.func.range();
   pkSize := gr.pkSize;
   P      := gr.P;
   i      := var("spuid", TInt, P);
   SS     := ScatSend(fAdd(N, N/P, i*(N/P)), pkSize, P, i);
   GD     := GathDist(N, pkSize, P, i);
   #GD     := GathDist(P, (N*pkSize)/P, P, i); # Produce GathDist with maximum packet size
   return(ComposeDists(e.rChildren(), DistSum(P, i, P, SS*GD)));
end;

#F This produces a SSGD on the right
PadRightSep := function(e)
   local gr, N, P, i, pkSize, SS, GD, di;
   gr     := e.rightMostParGath();
   N      := gr.func.range();
   pkSize := gr.pkSize;
   P      := gr.P;
   i      := var("spuid", TInt, P);
   SS     := ScatSend(fAdd(N, N/P, i*(N/P)), pkSize, P, i);
   GD     := GathDist(N, pkSize, P, i);
   #GD     := GathDist(P, (N*pkSize)/P, P, i); # Produce GathDist with maximum packet size
   return(ComposeDists(e, DistSum(P, i, P, SS*GD)));
end;

# These are required only for the CellSMP model, but shouldn't affect DMP, since DMP doesn't inject any SSs or GRs.
#F applyCellRules:
#F Assumes:
#F - There are no Buf()s in code
#F - Composes within a parallel region have already been converted to ComposeDistss
applyCellRules := function(sums, opts)
    local sumsorig, ss, sdgr, gr, ssgd, border;

#applyCellRules should not touch anything within a GT_ParStream. However, a
#sums expression can contain a mixture of ParStream and StreamChip generated
#expressions. We want to leave the former untouched, but want to act on the
#latter. We do this by marking things to be left untouched, and unmarking them at the end. This has 2 PRORLBEMS:
#
# 1) We shouldn't run ComposeDists before applyCellRules, though we also need to. Fix this.
# 2) If there's a GT_StreamPar rule, this whole thing falls apart.
#
# Essentially, we need a way of figuring out which were generated by StreamCore vs. StreamChip.
# Another approach: we could fuse the par and mbuf loops together via StreamCore's breakdown rule.
# In general, we could take care of StreamCore structures before we execute
# applyCellRules -- we could run DistMBuf in sigmaspl earlier, for instance.

    # Promote {Scat,Gath}Dist to ScatSend/GathRecv as needed
    SubstTopDown(sums, [ComposeDists, DistSum, DistSum], doAllPromotes);

    # Right-Pad DistSums that are standing by themselves (NOTE: cleanup this and the functions it calls)
    # Left-padding will be taken care of by statement below
    # NOTE: assumption here is, this is only for DFTxI.
    if IsBound(opts.doNotUseBlockCyclic) and opts.doNotUseBlockCyclic = true then
      SubstBottomUp(sums, @(1, DistSum, e->ObjId(e.rightMostParGath())=GathRecv), isContainedInComposeDistsG);
      SubstBottomUp(sums, @(1, DistSum, e->ObjId(e.leftMostParScat())=ScatSend),  isContainedInComposeDistsS);
    fi;

    # I/O padding: Don't do any input/output padding if these might get cancelled as a part of the data format change
    if IsBound(opts.doNotUseBlockCyclic) and opts.doNotUseBlockCyclic = true then
        SubstTopDown(sums, @(1, ComposeDists, e->ObjId(e.leftMostParScat())=ScatSend),  e->PadLeft(e));
        SubstTopDown(sums, @(1, ComposeDists, e->ObjId(e.rightMostParGath())=GathRecv),  e->PadRight(e));
    fi;

    #Error("BP");

    # NOTE: We shouldn't be doing this if we have a compose of multiple ParStream structures, for instance
    # Create table to map ScatSend/GathRecv pairs, so DMA can be done

    # NOTE: Why is this an else if structure? For StreamParChip, where we have
    # multiple stream stages each of which contain parallelism, all these
    # patterns could exist within the same sums expression.

    SubstTopDown(sums, [ComposeDists, DistSum, DistSum], permToData);
    SubstTopDown(sums, [ComposeDists, DistSum, DistSum, DistSum], permToData);
    SubstTopDown(sums, [ComposeDists, DistSum, DistSum, DistSum, DistSum], permToData);

    #if Length(Collect(sums, [ComposeDists, DistSum, DistSum])) = 1 then
    #     SubstTopDown(sums, [ComposeDists, DistSum, DistSum], permToData);
    #     else if Length(Collect(sums, [ComposeDists, DistSum, DistSum, DistSum])) = 1 then
    #               SubstTopDown(sums, [ComposeDists, DistSum, DistSum, DistSum], permToData);
    #     else if Length(Collect(sums, [ComposeDists, DistSum, DistSum, DistSum, DistSum])) = 1 then
    #               SubstTopDown(sums, [ComposeDists, DistSum, DistSum, DistSum, DistSum], permToData);
    #          fi;
    #          fi;
    #fi;

    # Fuse output (SDGR*SS...) to just an SS
    FixBorder(sums);

    return(sums);
end;

doAllPromotes := function(e, cx)
   local gathrecv, gathdist, scatsend, scatdist, i;

   # Assume: only one ScatSend and one GathRecv per DistSum.
   i := 1;

   # Check for GR/SD pair
   gathrecv := Collect(e._children[i],   GathRecv);
   scatdist := Collect(e._children[i+1], ScatDist);

   #Error("BP-doAllPromotes");

   if (Length(gathrecv) = 1 and Length(scatdist) = 1) then
     SubstTopDown(e, ScatDist,
        gs-> let(pkf := gs.pkSize / gathrecv[1].pkSize,
            ScatSend(fAdd(gs.N*pkf, gs.N*pkf/gs.P, gs.i*(gs.N*pkf/gs.P)), gs.pkSize/pkf, gs.P, gs.i))
     );
     return(e);
   fi;

   # Check for GD/SS pair
   gathdist := Collect(e._children[i],   GathDist);
   scatsend := Collect(e._children[i+1], ScatSend);

   if (Length(gathdist) = 1 and Length(scatsend) = 1) then
     SubstTopDown(e, GathDist,
        gs-> let(pkf := gs.pkSize / scatsend[1].pkSize,
            GathRecv(fAdd(gs.N, gs.N/gs.P, gs.i*(gs.N*pkf/gs.P)), gs.pkSize, gs.P, gs.i))
            );
     return(e);
   fi;

   return(e);

end;

#F permToData(e, cx) (e=expression, cx=context)
#F converts an expressions's function to an FList/FData
permToData := function(e, cx)
#permToData := function(d1, d2, cx)
   local gathrecv, scatsend, readMap, writeMap, loopnest, i;

   #gathrecv := e.rChildren()[1].rChildren()[2].rChildren()[Length(e.rChildren()[1].rChildren()[2].rChildren())];
   #scatsend := e.rChildren()[2].rChildren()[2].rChildren()[1];

   # Assume: only one ScatSend and one GathRecv per DistSum (no nested parallelism)
   for i in [1..(Length(Collect(e, DistSum))-1)] do
      gathrecv := Collect(e._children[i], GathRecv)[1];
      scatsend := Collect(e._children[i+1], ScatSend)[1];

      #PrintLine("permToData working on: ", gathrecv, " / ", scatsend);

      if gathrecv.func.__name__ = "FData" or scatsend.func.__name__ = "FData" then
         if gathrecv.func.__name__ = "FData" and scatsend.func.__name__ = "FData" then
            Error("Both are already FData. Why?\n");
         else
            Error("One of these funcs is an FData, the other is not.\n");
            #NOTE: Handle things correctly if both are already FDatas/FDataLists
         fi;
      fi;

      #loopnest := Concat(collectLoopVars(e, cx), [e.rChildren()[i].var, e.rChildren()[i+1].var]);
      loopnest := Concat(collectLoopVars(e, cx), [e._children[i].var, e._children[i+1].var]);

      [readMap, writeMap] := buildGlobalList(gathrecv, scatsend, loopnest);

      #NOTE: Use rSetChild or use .func directly here?
      SubstTopDown(e._children[i], GathRecv, e->GathRecv(FData([0..Length(readMap)-1]), e.pkSize, e.P, e.i));
      SubstTopDown(e._children[i+1], ScatSend, e->ScatSend(FData(readMap), e.pkSize, e.P, e.i));

      # Since we're doing a Gather-side normalize for the cell (the scatter
      # knows where to put an element, the gather is dumb), we write the
      # normalized FData instead of the actual readMap here.

   od;

   #PrintLine(gathrecv, "\n", scatsend, "\n", readMap, writeMap, "\n--\n", scatsend, "\n------------------\n\n");
   return e;
end;

#F Goes through all parents of e as given in cx.parents, and creates a list of
#F all parent loop variables in loop nest order.
collectLoopVars := function(e, cx)
   local parent, retList;
   retList := [];

   for parent in cx.parents do
      if (IsBound(parent.isSums) and parent.isSums = true) then
         if IsBound(parent.var) then
            retList := Concat(retList, [parent.var]);
         fi;
      fi;
   od;

   return retList;
end;

#F buildGlobalList(gathrecv, scatsend, loopnest)
#F Returns the read and write maps (as lists) for a given gath/scat pair
#F loopnest is a list of loop vars ordered by the gath/scat nesting
buildGlobalList := function(gathrecv, scatsend, loopnest)
   local i, j, gathMap, scatMap, t, readMap, writeMap;

   #NOTE: Need to add some error checking.

   gathMap := buildMap(gathrecv, loopnest);
   scatMap := buildMap(scatsend, loopnest);

   if Length(gathMap) <> Length(scatMap) then
      Error("buildGlobalList: packet sizes for scatter and gather don't match.");
   fi;

   t       := List([1..Length(gathMap)], i->0);
   readMap := List([1..Length(gathMap)], i->0);
   writeMap:= List([1..Length(gathMap)], i->0);

   for i in [1..Length(gathMap)] do
      t[gathMap[i]] := i;
   od;

   for i in [1..Length(gathMap)] do
      readMap[i] := t[scatMap[i]];
   od;

   for i in [1..Length(gathMap)] do
      for j in [1..Length(gathMap)] do
         if readMap[j] = i then
            writeMap[i] := j;
         fi;
      od;
   od;

   #NOTE: assuming that element values will range from 0..(n-1).
   readMap := List([1..Length(readMap)], i->readMap[i]-1);
   writeMap := List([1..Length(writeMap)], i->writeMap[i]-1);

   return [readMap, writeMap];
end;

#F buildMap(gathscat, loopnest)
#F Fully unrolls a gath/scat's function to a list based on the loop nest ordering given by loopnest
#F If loopnest has vars that the gathscat doesn't depend on, they will be ignored
buildMap := function(gathscat, loopnest)
   local map, reqFreeVars, loopVar;

   reqFreeVars := requiredFreeVars(loopnest, Flat(gathscat.free()));
   map := gathscat.func.lambda().tolist();
   for loopVar in reqFreeVars do
       map := Flat( List([0..(loopVar.range)-1], i->SubstVars(Copy(map), rec((loopVar.id) := V(i)))) );
   od;
   map := List([1..Length(map)], i->(map[i].ev())+1);

   return map;
end;

#F requiredFreeVars(loopnest, freevars)
#F
requiredFreeVars := function(loopnest, freevars)
   local returnList, i;
   returnList := [];
   for i in loopnest do
      if i in freevars then
         if not i in returnList then
            returnList := Concat(returnList, [i]);
         fi;
      fi;
   od;
   return returnList;
end;

isContainedInComposeDistsG := function(e,cx)
  local parent, foundComposeDist;
  foundComposeDist := false;
  for parent in cx.parents do
     if ObjId(parent)=ComposeDists then foundComposeDist := true; fi;
  od;

  if foundComposeDist=true then
     return(e);
  else
     return(PadRightSep(e));
  fi;
end;

isContainedInComposeDistsS := function(e,cx)
  local parent, foundComposeDist;
  foundComposeDist := false;
  for parent in cx.parents do
     if ObjId(parent)=ComposeDists then foundComposeDist := true; fi;
  od;

  if foundComposeDist=true then
     return(e);
  else
     return(PadLeftSep(e));
  fi;
end;


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

Declare(sve_gath, sve_scat, sve_ld2, sve_st2);

Class(TVectSVE, TReal);
TVectSVEx2 := TVect(TVectSVE, 2);

Class(TBoolSVE, TBool);
Class(TInt64SVE, TInt);

# SVE vector loop
Class(sve_loopn, loopn);

Class(sve_svcntd, Exp, rec(
    computeType := self >> TInt
));

# svfloat64_t svld1_gather_[u64]offset[_f64](svbool_t pg, const float64_t *base, svuint64_t offsets)

Class(sve_gath, Loc, rec(
    __call__ := (self, loc, n, pg, stride) >> WithBases(self,
        rec(operations := NthOps,
            loc := toExpArg(loc),
            n := toExpArg(n),
            pg := toExpArg(pg),
            stride := toExpArg(stride))).setType(),

    can_fold := False,

    rChildren := self >> [self.loc, self.n, self.pg, self.stride],
    rSetChild := rSetChildFields("loc", "n", "pg", "stride"),

    ev := self >> let(e := self.eval(), Cond(IsBound(e.v), e.v, e)),

    eval := self >> sve_gath(self.loc.eval(), self.n.eval(), self.pg.eval(), self.stride.eval()),

    computeType := self >> TVectSVE,

    isExpComposite := true
));

Class(sve_scat, Command, rec(
   __call__ := (self, loc, n, pg, stride, exp) >> WithBases(self,
       rec(operations := CmdOps,
       loc := toAssignTarget(loc),
       n := toExpArg(n),
       pg := toExpArg(pg),
       stride := toExpArg(stride),
       exp := toExpArg(exp)
       )),

   rChildren := self >> [self.loc, self.n, self.pg, self.stride, self.exp],
   rSetChild := rSetChildFields("loc", "n", "pg", "stride", "exp"),
   unroll := self >> self,

   print := (self,i,si) >> let(name := Cond(IsBound(self.isCompute) and self.isCompute,
                                            gap.colors.DarkYellow(self.__name__),
                                            IsBound(self.isLoad) and self.isLoad,
                                            gap.colors.DarkRed(self.__name__),
                                            IsBound(self.isStore) and self.isStore,
                                            gap.colors.DarkGreen(self.__name__),
                                            self.__name__),
                                Print(name, "(", self.loc, ", ", self.n, ", ", self.pg, ", ", self.stride, ", ", self.exp, ")"))
));

# svfloat64x2_t svld2[_f64](svbool_t pg, const float64_t *base)

Class(sve_ld2, Loc, rec(
    __call__ := (self, loc, n, pg) >> WithBases(self,
        rec(operations := NthOps,
            loc := toExpArg(loc),
            n := toExpArg(n),
            pg := toExpArg(pg))).setType(),

    can_fold := False,

    rChildren := self >> [self.loc, self.n, self.pg],
    rSetChild := rSetChildFields("loc", "n", "pg"),

    ev := self >> let(e := self.eval(), Cond(IsBound(e.v), e.v, e)),

    eval := self >> sve_ld2(self.loc.eval(), self.n.eval(), self.pg.eval()),

    computeType := self >> TVectSVEx2,

    isExpComposite := true
));

Class(sve_st2, Command, rec(
   __call__ := (self, loc, n, pg, exp) >> WithBases(self,
       rec(operations := CmdOps,
       loc := toAssignTarget(loc),
       n := toExpArg(n),
       pg := toExpArg(pg),
       exp := toExpArg(exp)
       )),

   rChildren := self >> [self.loc, self.n, self.pg, self.exp],
   rSetChild := rSetChildFields("loc", "n", "pg", "exp"),
   unroll := self >> self,

   print := (self,i,si) >> let(name := Cond(IsBound(self.isCompute) and self.isCompute,
                                            gap.colors.DarkYellow(self.__name__),
                                            IsBound(self.isLoad) and self.isLoad,
                                            gap.colors.DarkRed(self.__name__),
                                            IsBound(self.isStore) and self.isStore,
                                            gap.colors.DarkGreen(self.__name__),
                                            self.__name__),
                                Print(name, "(", self.loc, ", ", self.n, ", ", self.pg, ", ", self.exp, ")"))
));

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:
# Computes an isomorphic permutation group of the reduced degree.
# INPUT:
# > G - a finite group
# OUTPUT:
# > A group isomorphic to G which is a permutation group (and saved as such by GAP).
#   The degree is also reduced with SmallerDegreePermutationRepresentation (not guaranteed to be optimal, though).
InstallGlobalFunction( GSPermutationGroup, function( G )
    return Image( SmallerDegreePermutationRepresentation( Image( IsomorphismPermGroup( G ) ) ) );
end );

max_stars_repo_name: piotrmizerka/math-work-tools

max_issues_repo_name: piotrmizerka/math-work-tools

max_forks_repo_name: piotrmizerka/math-work-tools


text:
a := [85, 88, 75, 66, 25, 29, 83, 39, 97];
b := [68, 41, 10, 49, 16, 65, 32, 92, 28, 98];

# Compute a decimal approximation of a rational
Approx := function(x, d)
	local neg, a, b, n, m, s;
	if x < 0 then
		x := -x;
		neg := true;
	else
		neg := false;
	fi;
	a := NumeratorRat(x);
	b := DenominatorRat(x);
	n := QuoInt(a, b);
	a := RemInt(a, b);
	m := 10^d;
	s := "";
	if neg then
		Append(s, "-");
	fi;
	Append(s, String(n));
	n := Size(s) + 1;
	Append(s, String(m + QuoInt(a*m, b)));
	s[n] := '.';
	return s;
end;

PermTest := function(a, b)
	local c, d, p, q, u, v, m, n, k, diff, all;
	p := Size(a);
	q := Size(b);
	v := Concatenation(a, b);
	n := p + q;
	m := Binomial(n, p);
	diff := Sum(a)/p - Sum(b)/q;
	all := [1 .. n];
	k := 0;
	for u in Combinations(all, p) do
		c := List(u, i -> v[i]);
		d := List(Difference(all, u), i -> v[i]);
		if Sum(c)/p - Sum(d)/q > diff then
			k := k + 1;
		fi;
	od;
	return [Approx((1 - k/m)*100, 3), Approx(k/m*100, 3)];
end;

# in order, % less or greater than original diff
PermTest(a, b);
[ "87.197", "12.802" ]