Package $ext      -- Ext Modules

Alias Ext := $ext;
  
Define About()
  PrintLn "
KeyWords : Ext, Depth
Author   : A.Damiano
E-mail   : alberto@tlc185.com
Version  : CoCoA 4.6
Date     : May 2006
Thanks To: The author is thankful to M.Kreuzer and L. Robbiano for their
extreme patience.  The author is also thankful to A.Bigatti for teaching him
how to transform his raw code into a package and for her constant support.
"
End; -- About


  -----------------[   Manual   ]----------------------

Define Man()
  PrintLn "
------------------[ PRELIMINARIES ]-------------------

This is a new version with respect to the ext.pckg that was distributed with CoCoA 4.3 and 4.4

VERY IMPORTANT: only homogeneous modules are supported. Use non-homogeneous objects at your own risk!

Main differences with the previous version include:

1) SHIFTS have been removed, consequently only standard homogeneous modules and quotients are supported

2) as a consequence of 1), the type Tagged(""shifted"") has been removed. Ext will just be a Tagged(""Quotient"")

3) The former functions Presentation(), HomPresentation() and Ker(Presentation) have been removed

4) The algorithm uses Res() to compute the maps needed, and not SyzOfGens any longer, believed to cause troubles

5) The function Ext always has THREE variables, see syntax...

Finally, please note that this package is still being tested so please report bugs to alberto@tlc185.com




---------------------[ SYNTAX of EXT ]----------------------

Some directions on how to use these functions:

Ext.Ext(I: INT,  M: Tagged(""Quotient""), N: Tagged(""Quotient"")):Tagged(""Quotient"");

Ext.Ext(I: LIST, M: Tagged(""Quotient""), N: Tagged(""Quotient"")):Tagged(""ExtList"");

NOTE:in the second form I is a LIST of non-negative integers of type INT and the output is a LIST of 
pairs (lists of two elements) of the type [J,E] where J is an element of I and E the corresponding
Ext module calculated using Ext(J,M,N). In this case the modules are printed in a ""nicer way"".

IMPORTANT:the only exception to M or N (or even the output) is when they are either 
a zero module or a free module. In these cases their type is MOD.  


This function calculates a presentation (i.e. a quotient of a free module) for the module

Ext^I_R(M,N)     , where R is the current R

using the maps Phi*_n  in the dual of a minimal free resolution and presenting 

Ker(Phi*_I)/Im(Phi*_{I-1})   

as a quotient using syzygies. If the Ext module required is zero, it will simply return

Module([0]);
--------------------------

The variable I could also be a LIST of non-negative integers. In this case 
the function Ext(...) prints all the Ext modules correpsonding to the integers in I.
The output is of special type Tagged(""ExtList"") which is basically just a list so that
it can be assigned to a variable and each component of the list can be utilized individually.


NOTE:The input is pretty flexible in terms of what you can use for M and N. 
Definitely if you need to use a zero module (do you really?) you can define it in as is, for example as

M:=Module([0]);

or

N:=Ideal(0);

This is in fact not the only way to define a zero module: one could use

M:=R^2/Module([1,0],[0,1]);

or

N:=R/Ideal(1);

or other things like that. Ext(...) should recognize that and return the zero module right away. 

It is also possible to just use a free module as M or N. Just define something like

M:=R^3;

or 

N:=Module([1,0],[0,1]);

etc.

VERY IMPORTANT: CoCoA cannot accept the ring R as one of the inputs, so if you want
to calculate the module Ext^I_R(M,R) you need to type something like

Ext.Ext(I,M,Ideal(1));

or

Ext.Ext(I,M,R^1);


or

Ext.Ext(I,M,R/Ideal(0));

or using any other way to define the free module R^1 that may come to your mind,
-hopefully- this will be recognized as R^1.
Let me know if you find a sintax that is not supported,
those are always the most tedious cases to handle. 
Note that in this case there is a special function called to calculate Ext
that should be a little faster than the original version of the Ext package.
I would consider in future updates having a specific function Ext_R(I,M) for this particular 
(and most used) case, like Singular does.



------------------[ EXAMPLES using EXT ]------------------

Use R::=Q[x,y,z];

I:=Ideal(x^5,y^3,z^2);

Ideal(0):(I);

Ideal(0)
-------------------------------
Hom.Hom(R^1/Module(I),R^1);   -- from Hom package, see $contrib/hom

Module([0])
-------------------------------
Ext.Ext(0,R/I,R^1);   --- all those things should be isomorphic

Module([0])
-------------------------------
Ext.Ext(0..4,R/I,R/Ideal(0));


Ext^0 = Module([0])

Ext^1 = Module([0])

Ext^2 = Module([0])

Ext^3 = R^1/Module([x^5], [-y^3], [z^2])

Ext^4 = Module([0])

-------------------------------
N:=Module([x^2,y],[x+z,0]);

Ext.Ext(0..4,R/I,R^2/N);


Ext^0 = Module([0])

Ext^1 = Module([0])

Ext^2 = R^8/Module([0, 0, 0, -1, 0, 0, 0, 0], [0, 1, 0, 0, 0, 0, 0, 0], [1, x - z, 0, 0, 0, 0, 0, 0], [0, 0, y, x - z, 1, 0, 0, 0], [0, 0, y, x - z, 0, 0, 0, 0], [0, x^3, 0, 0, 0, x + z, 0, 0], [y^2, 0, -z^2, 0, 0, 0, 0, 0], [0, x^4 - x^3z, 0, 0, 0, -z^2, -1, 0], [-y^3z, x^2y^3 - xy^3z + y^3z^2, 0, -x^4 + x^3z - x^2z^2, 0, y^3, 0, 0], [0, x^4 - x^3z, 0, 0, 0, -z^2, 0, 0], [0, 0, 0, 0, 0, 0, 0, -1], [-y^2z^3, -x^4y^2 + x^3y^2z, x^5 + z^5, 0, 0, y^2z^2, 0, 0])

Ext^3 = R^2/Module([x^5, 0], [0, x^5], [-y^3, 0], [0, -y^3], [z^2, 0], [0, z^2], [x^2, y], [x + z, 0])

Ext^4 = Module([0])

-------------------------------

As you can see the Ext^2 is not  presented minimally. This is an issue I will try to resolve for the
next version. 


------------------[ SYNTAXT of DEPTH ]------------------- 

The function Depth(...) has a double sintax:

Ext.Depth(I: IDEAL, M: Tagged(""Quotient"")): INT;

or

Ext.Depth(M: Tagged(""Quotient""): INT;

It calculates the I-depth of the module M.

I is a homogeneous IDEAL and M is a nonzero module,
expressed in the same way you would use for the function Ext, so it has to be either free or 
a quotient of the type R^s/W. However, M does not have to be necessarily homogeneous.

If I is not specified like in the second form, then we assume that I is the maximal ideal
generated by the indeterminates. 

If M is homogeneous and I is the ideal generated by the indeterminates (or if it is not specified),
then the function Depth uses the Auslander-Buchsbaum formula: 

Depth_I(M)=N-pd(M)

where N is the number of intederminates and pd(M) is the lenght of a minimal free resolution of M.

Otherwise, it calculates the modules Ext^I_R(R/I,M) and return the first index I such that Ext^I is
not zero.

-------------------[ EXAMPLES using DEPTH ]-----------------------

Use R::=Q[x,y,z];

Ext.Depth(Ideal(1)); -- it's the (x,y,z)-depth of the entire ring so it's 3

3
-------------------------------


I:=Ideal(x^5,y^3,z^2);

--- let's check that it is zerodimensional and CM

Dim(R/I);
Ext.Depth(R/I);

0
-------------------------------
0
-------------------------------


N:=Module([x^2,y],[x+z,0]);

Ext.Depth(I,R^2/N);  --- a max reg sequence would be (z^2,y^3)

2
-------------------------------

Use R::=Q[x,y,z,t,u,v];

N:=Module([x,y],[-y,x],[z,t],[-t,z],[u,v],[-v,u]);  -- Cauchy-Riemann system in three complex vars!

--- is it CM?

Depth(R^2/N);
Dim(R^2/N);

3
-------------------------------
3
-------------------------------

--- yes!


M:=Module([x,y,z],[t,v,u]);

Res(R^3/M);

0 --> R^2(-1) --> R^3
-------------------------------
Ext.Depth(R^3/M); -- using Auslander Buchsbaum 6-1=5

5
-------------------------------

Dim(R^3/M);  -- not CM

6
-------------------------------

Depth(Ideal(x,y,z,t),R^2/N); -- can you find a reg sequence here?

2
-------------------------------


Please send your (always very appreciated) comments, questions, bugs reports and ideas to

alberto@tlc185.com
"
End; --- Man



  ----------------------[ MAIN  FUNCTIONS ]-----------------------------

Define Ext(I,RM,RN);
  -- I:= a positive INTeger I, or a LIST of positive INTegers, and two QUOTIENT modules (may be zeromodules)
  -- O:= a presentation of the module Ext^I_R(RM,RN) if I is an Integer 
  --
  --     if I is a list it PRINTS (only prints!!) all modules Ext^J_R(RM,RN) for each J in I
  --
  If Type(I)=LIST Then
    Exts:=[[J,Ext.Ext(J,RM,RN)]|J In I];
    Return Ext.Tagged(Exts,"ExtList");
  Elsif Type(I)=INT Then
    If I<0 Then Error("Ext: Index must be positive");EndIf;
    If Ext.IsZeroModule(RM) Or Ext.IsZeroModule(RN) Then
      Return Module([0]);
    Elsif Ext.IsFreeModule(RM) Then
      If IsZero(I) Then
	S:=NumComps(RM);
	T:=NumComps(RN);
	If Ext.IsFreeModule(RN) Then Return R^(S*T);EndIf;
	N:=@RN;
	U:=Ext.SpecialHom(S,N);
	Return Ext.MinimalPresentation(CurrentRing()^(S*T)/U);
      Else 
	Return Module([0]);
      EndIf;
    Else ---- RM is not a free module so it is of the type R^S/M
      M:=@RM;
      ResM:=Res(M);
      If Ext.IsFreeModule(RN) And NumComps(RN)=1 Then  --- Then it is Ext^I(R^S/M,R)
	Return Ext.MinimalPresentation(Ext.NoShiftExt_FromRes(I,ResM));
      Elsif Ext.IsFreeModule(RN) Then
	N:=Ext.ZeroModule(RN);
	Return Ext.MinimalPresentation(Ext.NoShiftExt_RM_RN_FromRes(I,ResM,N));
      Else
	N:=@RN;
	Return Ext.MinimalPresentation(Ext.NoShiftExt_RM_RN_FromRes(I,ResM,N));
      EndIf;
    EndIf;
  Else
    Error("Ext: first argument must be INT or LIST of INT");
  EndIf;
EndDefine;


  ---------------------- DEPTH

Define Depth(...);
  If Len(ARGV)=2 Then 
    I:=ARGV[1];
    M:=ARGV[2];
  Elsif Len(ARGV)=1 Then
    M:=ARGV[1];
    I:=Ideal(Indets());
  Else Error("Bad arguments: expected a module and (optional) an ideal");
  EndIf;
  If Type(M)<>TAGGED("Quotient") Then Error("Depth: quotient ideal or module required");EndIf;
  M:=Tagged(Module(@M),"Quotient");
  If Ext.IsZeroModule(M) Then Error("Depth: nonzero module required");EndIf;
  If I=Ideal(Indets()) And IsHomog(@M) Then 
    ResMod:=Res(M);
    ResLen:=Len(ResMod[2]);
    Return NumIndets()-ResLen+1;
  Else
    N := NumIndets();
    D := 0;
    ResIdeal:=Res(I);
    P := Ext.NoShiftExt_RM_RN_FromRes(D,ResIdeal,@M);
    While D<N And P=Module([0]) Do
      D := D+1;
      P := Ext.NoShiftExt_RM_RN_FromRes(D,ResIdeal,@M);
    EndWhile;
    Return D;
  Endif;
EndDefine;




  -----------------[   pretty printing   ]-----------------------

Define Tag() Return  Ext.PkgName()+"."; EndDefine;
Define Tagged(X, T) Return  Tagged(X, Ext.Tag()+T); EndDefine;

Define Print_ExtList(X)
  ForEach M In X Do
    Print "Ext^";
    Print M[1];
    Print " = ";
    PrintLn M[2];
    PrintLn;
  EndForEach;
EndDefine;


  ----------------------[ MISC ]---------------------------------

Define PhiToMat(Phi);
  -- I: a LIST of lists or vectors defining the image of the map Phi
  -- O: the MATRIX whose columns are the elements of the input
  --
  Return Transposed(Mat(Phi));
EndDefine;

Define MatToPhi(A);
  -- I: a matrix representing a map (image generated by columns);
  -- O: a LIST of vectors, the columns of A
  --
  Return Gens(Ext.MatToModule(A));
EndDefine;


Define MatToModule(A);
  -- I: a MATRIX
  -- O: the MODULE generated by the columns of the matrix A
  --
  Return Module(Transposed(A));
EndDefine;

Define PhiToModule(Phi);
  -- I: a LIST of lists or vectors defining the image of the map Phi
  -- O: the MODULE generated by the elements of Phi;
  --
  Return Module(Phi); -- oh well
EndDefine;

Define ZeroModule(M);
  -- I: a Module (or ideal) in whichever form
  -- J: the zero module with the same number of components of M
  S:=NumComps(M);
  V:=NewVector(S,0);
  Return Module(V);
EndDefine;


Define PunctualProduct(V,W);
  -- I: two vectors or lists of the same length
  --
  -- O: the vector or list given by the componentwise product
  --
  If Len(V)<>Len(W) Then Error("PunctualProduct: expected two LIST/VECTORS of the same length");EndIf;
  Return Vector([V[J]*W[J]|J In 1..Len(V)]);
EndDefine;

Define Grado(...);
  -- I: pretty much anything for the first argument, the second is optional (an indet)
  -- O: the degree of the first input (wrt the indet given by the second arg if any)
  --   
  -- NOTE: the degree of zero is set to be zero!
  --
  If Len(ARGV)=1 Then
    W := ARGV[1];
    If IsZero(W) Then
      Return 0;
    Else
      Return Deg(W);
    EndIf;
  Elsif Len(ARGV)=2 Then
    W := ARGV[1];
    IND := ARGV[2];
    If IsZero(W) Then Return 0;
    Else
      Return Deg(W,IND);
    EndIf;
  EndIf;
EndDefine;

Define TensorMat(A,B);
  -- I: two MATs of any dimension
  -- O: the tensor product matrix as in def 3.3.12
  RowNum := Len(A);
  ColNum := Len(A[1]);
  T1 := [[A[I,J]*B|J In 1..ColNum]|I In 1..RowNum];
  Return BlockMatrix(T1);
EndDefine;

Define Truncate(E,ToKeep);
  -- I: E is a LIST of polys, a VECTOR, a LIST of vectors or lists of polys, a MODULE 
  --	  ToKeep is an INT
  --  
  -- O: for a LIST of polys or vectors, it returns the truncation up the the first ToKeep components of each element
  --    for a MODULE it truncates each generator and returns a module
  --
  If Type(E)=LIST And Type(E[1])<>LIST And Type(E[1])<>VECTOR Then
    If Len(E)<ToKeep Then Error("Truncation index is too big");EndIf;
    Return [E[I]|I In 1..ToKeep];
  Elsif Type(E)=VECTOR Then
    If Len(E)<ToKeep Then Error("Truncation index is too big");EndIf;
    Return Vector([E[I]|I In 1..ToKeep]);
  Elsif Type(E)=LIST And (Type(E[1])=LIST Or Type(E[1])=VECTOR) Then 
    Return [Truncate(V,ToKeep)|V In E]
  Elsif Type(E)=MODULE Then
    Return Module(Truncate(Gens(E),ToKeep));
  Elsif Type(E)=TAGGED("Quotient") Then
    M:=@E;
    NewM:=Ext.Truncate(M);
    Return Tagged(Module(NewM),"Quotient");
  Else
    Error("Truncate: expected a LIST of POLY, a VECTOR, a LIST of such, or a MODULE");
  EndIf;
EndDefine;


  ---------------- RES -----------------------

Define ShiftsFromRes(RES);
  -- I: a tagged resolution of a module given with generators (not a quotient)
  --
  -- O: a list of shifts, i.e. a list of lists of INT, each representing the positive shifts at a step of RES
  RealRes:=@RES;
  RealRes:=RealRes[2]; --- this is a list of various lists each containing a syzygy mod, the betti nrs and the shifts
  Return [ShiftVector(M)|M In RealRes];
EndDefine;

Define ShiftVector(M);
  -- I: one of the elements in the untagged resolution, i.e. a list of a (shifted) module, 
  --		and pairs [E,S] where E is the corresponding Betti in degree -S
  -- O: a list of shifts, each being a list of positive INT
  V := [];
  For I := 2 To Len(M) Do -- the first element of M is a module
    T := M[I];
    Expon := T[1];
    Shif := -T[2];
    S := NewList(Expon,Shif);
    V := Concat(V,S);
  EndFor;
  Return V;
EndDefine;

Define MapsFromRes(RES);
  -- I: a tagged resolution of a module given with generators (not a quotient)
  --
  -- O: the LIST of the maps in the resolution, starting from the right. 
  --    each map is a list of vectors
  RealRes:=@RES;
  RealRes:=RealRes[2]; --- this is a list of various lists each containing a syzygy mod, the betti nrm and the shifts
  RealMod:=Minimalized(RealRes[1,1]);
  Phis:=[Gens(RealMod)];
  For I:=2 To Len(RealRes) Do -- if the module has no syzygies the loop is not executed
    Phi:=Gens(RealRes[I,1]);
    Append(Phis,Phi);
  EndFor;
  Return Phis;
EndDefine;

Define DualComplex(Phis,Shifts);
  -- I: a LIST of maps Phis like the one coming from MapsFromRes 
  --    a LIST of lists of integers representing shifts, like the output of ShiftsFromRes
  --
  -- O: a pair of LISTs, the first is the list of the transposed maps and the second is the list of shifts for the dual complex
  NewPhis:=[];
  ForEach Phi In Phis Do
    A:=Ext.PhiToMat(Phi);
    NewPhi:=Ext.MatToPhi(Transposed(A));
    Append(NewPhis,NewPhi);
  EndForEach;
  Return [NewPhis,Shifts];
EndDefine;

  --------------------- KER ---------------------

Define NoShiftKer_Free(Phi);
  -- I: A map Phi:R^r-->R^s  as a LIST of generators of Im(Phi)
  -- O: LIST of Generators For Ker(Phi)as constructed In prop 3.3.1
  W:=Ext.PhiToModule(Phi);
  S:=SyzOfGens(W);
  ToKeep:=Len(Phi);
  Return [Ext.Truncate(V,ToKeep)|V In Gens(S)];
EndDefine;

Define NoShiftKer(Phi,N);
  -- I: A map Phi:R^r-->R^s/N  as a LIST of generators of Im(Phi)
  -- O: LIST of Generators For Ker(Phi) as constructed In prop 3.3.1
  W:=Ext.PhiToModule(Phi);
  V:=Module(Concat(W.Gens,N.Gens));
  S:=SyzOfGens(V);
  ToKeep:=Len(Phi);-- r;
  Return [Ext.Truncate(V,ToKeep)|V In Gens(S)];
EndDefine;


Define NoShiftQuotient(M,N);
  -- I: two modules M And N
  -- O: a presentation of the quotient M/(M intersection N), tagged "quotient"
  --    Using ???
  If M=N Then Return Module([0]);EndIf;
  G:=Gens(M);
  H:=Gens(N);
  All:=Concat(G,H);
  S:=Syz(All);
  ToKeep:=Len(G);
  Denominator:=Module([Ext.Truncate(V,ToKeep)|V In Gens(S)]);
  Minimalize(Denominator);
  Return (Tagged(Denominator,"Quotient"));
EndDefine; -- NoShiftQuotient

  -------------------- EXT ---------------


Define NoShiftExt_FromRes(I,RES);
  -- I: a free resolution For a module M
  -- O: a presentation of Ext^I_R(M,R);
  Maps:=Ext.MapsFromRes(RES);
  Dual:=Ext.DualComplex(Maps,[0]);
  DualMaps:=Dual[1];
  If IsZero(I) Then
    RealMod:=Ext.PhiToModule(Maps[1]);----------------------RealMod;
    --   Return Tagged(RealMod,"Quotient");
    Return $hom.HomR_MR_N(RealMod,Module([1]));
  ElsIf I<Len(Maps) Then
    PhiUp:=DualMaps[I+1];--PhiUp;
    PhiDown:=DualMaps[I];--PhiDown;
    K:=Ext.NoShiftKer_Free(PhiUp);
    Kernel:=Ext.PhiToModule(K);--Kernel;
    Im:=Ext.PhiToModule(PhiDown);--Im;
    Return Ext.NoShiftQuotient(Kernel,Im);
  ElsIf I=Len(Maps) Then
    PhiDown:=Last(DualMaps);
    Im:=Ext.PhiToModule(PhiDown);
    Return Tagged(Im,"Quotient");
  Else
    Return Module([0]);
  EndIf;
EndDefine;

Define IsZeroExt_FromRes(I,RES);
  -- I: a free resolution For a module M
  -- O: TRUE If  Ext^I_R(M,R)=0;
  Maps:=Ext.MapsFromRes(RES);
  Dual:=Ext.DualComplex(Maps,[0]);
  DualMaps:=Dual[1];
  If IsZero(I) Then
    RealMod:=Ext.PhiToModule(Maps[1]);----------------------RealMod;
    --   Return Tagged(RealMod,"Quotient");
    Return $hom.HomR_MR_N(RealMod,Module([1]));
    
    RealMod:=Ext.PhiToModule(Maps[1]);
    N:=NumComps(RealMod);
    Return RealMod=CurrentRing()^N;
  ElsIf I<Len(Maps) Then
    PhiUp:=DualMaps[I+1];--PhiUp;
    PhiDown:=DualMaps[I];--PhiDown;
    K:=Ext.NoShiftKer_Free(PhiUp);
    Kernel:=Ext.PhiToModule(K);--Kernel;
    Im:=Ext.PhiToModule(PhiDown);--Im;
    Return Kernel=Im;
  ElsIf I=Len(Maps) Then
    Return FALSE;
  Else
    Return TRUE;
  EndIf;
EndDefine;

Define NoShiftExt_RM_RN_FromRes(I,RES,N)
  -- I: an INT I being the index of the Ext module
  --    a free RESolution of the module M (must be homog)
  --    a MODule N
  -- O: a the I-th Ext_R module of (R^s/M,R^t/N);
  Maps:=Ext.MapsFromRes(RES);
  T:=NumComps(N);
  N:=Module(@N);
  H:=Ext.PhiToMat(N.Gens);
  Dual:=Ext.DualComplex(Maps,[0]); --- Only regarding the frist module, so all phis are transposed
  DualMaps:=Dual[1];
  --DualMaps;
  If IsZero(I) Then
    RealMod:=Ext.PhiToModule(Maps[1]);----------------------RealMod;
    Return $hom.HomR_MR_N(RealMod,N);
  ElsIf I<Len(Maps) Then
    PhiUp:=DualMaps[I+1];
    --PhiUp;
    NumeratorMap:=Ext.TensorMat(PhiUp,Identity(T)); -- this is the map we need the Ker of 
    PhiDown:=DualMaps[I];
    --PhiDown;
    DenominatorMap:=Ext.TensorMat(PhiDown,Identity(T));
    --Len(PhiUp)=Len(PhiDown[1]);
    S:=[Len(PhiDown),Len(PhiUp),Len(PhiUp[1])]; --- S_i-1, S_i and S_i+1
    --S;
    U:=[Ext.MatToModule(Ext.TensorMat(Identity(S[J]),H))|J In 1..3];
    --U[2];
    K:=Ext.NoShiftKer(NumeratorMap,U[3]);
    Kernel:=Ext.PhiToModule(K);
    --Kernel;
    Im:=Ext.PhiToModule(DenominatorMap);
    --Im;
    Return Ext.NoShiftQuotient(Kernel,(Im+U[2]));
  ElsIf I=Len(Maps) Then
    PhiDown:=Last(DualMaps);
    --PhiDown;
    SLast:=Len(PhiDown[1]);
    --SLast;
    Im:=Ext.PhiToModule(PhiDown);
    DenominatorMap:=Ext.TensorMat(PhiDown,Identity(T));
    --DenominatorMap;
    U:=Ext.MatToModule(Ext.TensorMat(Identity(SLast),H));
    --U;
    Return Tagged(Ext.PhiToModule(DenominatorMap)+U,"Quotient");
  Else
    Return Module([0]);
  EndIf;
EndDefine;
  
  ---------------------- IS IT ZERO? IS IT FREE?? --------------

Define IsZeroModule(Q);
  -- I: a MODule Q in whichever form
  -- O: TRUE if it is isomorphic to the zeromodule
  If Type(Q)=MODULE Or Type(Q)=IDEAL Then
    Q:=Module(Q);
    ZQ:=Ext.ZeroModule(Q);
    Return ZQ=Q;
  Elsif Type(Q)=TAGGED("Quotient") Then
    M:=@Q;
    S:=NumComps(M);  --- Q=R^S/M
    Return M=CurrentRing()^S;
  Else 
    Error("IsZeroModule: Bad Parameters");
  EndIf;
EndDefine;

Define IsFreeModule(Q);
  -- I: a MODule Q in whichever form
  -- O: TRUE if it is THE free module R^s with s=numComps(Q)
  --
  -- NOTE: IsFreeModule(Ideal(x))=FALSE!!!
  -- 
  If Type(Q)=TAGGED("Quotient") Then
    M:=@Q;
    ZM:=Ext.ZeroModule(M);
    Return ZM=M;
  Elsif Type(Q)=MODULE Or Type(Q)=IDEAL Then
    Q:=Module(Q);
    S:=NumComps(Q);
    Return Q=CurrentRing()^S;
  Else 
    Error("IsFreeModule: Bad Parameters");
  EndIf;
EndDefine;



Define SpecialHom(S,N);
  -- O: The module U such that Hom(R^S,R^T/N)=R^(S*T)/U
  --
  A:=Ext.PhiToMat(N.Gens);
  Return Ext.TensorMat(Identity(S),A);
EndDefine;

-------------- FINALLY, Minimal presentation of Ext!  Added May 2006 --------------


Define MinimalPresentation(Q); 
  -- I: A QUOTIENT Module of the type R^s/M, or a Zero module
  -- O: An Isomorphic quotient R^t/N, minimally presented Using [Kreuzer Robbiano II]  
  -- 
  -- Basically it's the gauss elimination On the presentation matrix Mat(M) where the pivots are the nonzero constants 
  -- 
  If IsZero(Q) Then Return Q; EndIf;
  M:=@Q; 
  A:=Mat(M); -- A; -- the matrix will be erre times esse 
  Erre:=Len(A); -- Erre;
  Esse:=Len(A[1]); -- Esse;
  N:=NumIndets(); 
  Zeros:=NewList(N,0); -- Zeros;
  EvalA:=Eval(A,Zeros); -- EvalA; 
  MinimalityCondition:=IsZero(EvalA); -- TRUE iif the matrix does Not have nonzero constants; 
  While Not(MinimalityCondition) Do 
    Found:=FALSE; 
    For I:=1 To Erre Do 
      For J:=1 To Esse Do 
		If Not(IsZero(A[I,J])) And Deg(A[I,J])=0 Then 
		  NewA:=Ext.GaussReduction(A,I,J); -- PrintLn(NewA); 
		  Found:=TRUE; 
		EndIf; 
		If Found Then Break; EndIf; -- out of the column loop 
      EndFor; 
      If Found Then Break; EndIf; -- out of the row loop 
    EndFor; 
    If Not(Found) Then
		Return Tagged(Module(A),"Quotient");
    EndIf;
    A:=NewA; -- A;
    If A=Mat[[]] Then 
      Return Module([0]); 
    EndIf; 
    If IsZero(A) Then 
      Esse:=Len(A[1]); 
      Return CurrentRing()^Esse; 
    EndIf; 
    ---- In any other case A is a nonzero rectangluar matrix, so we keep pivoting 
    Erre:=Len(A); -- Erre; 
    Esse:=Len(A[1]); -- Esse;
    EvalA:=Eval(A,Zeros); -- EvalA; 
    MinimalityCondition:=IsZero(EvalA); 
  EndWhile; 
  NewM:=Module(ReducedGBasis(Module(A))); 
  Return Tagged(NewM,"Quotient"); 
EndDefine; 
 
Define GaussReduction(A,I,J); 
  -- I: a Matrix A whose element A[I,J] is a nonzero constant 
  -- O: performs a Gauss reduction On A with pivot A[I,J], Then erase row I And column J from A 
  Pivot:=A[I,J];-- Pivot;
  Erre:=Len(A); 
  Esse:=Len(A[1]); 
  -- If the matrix is 1 by 1, then the outcome will be the empty matirx!
  If Erre=1 And Esse=1 Then 
	Return Mat[[]];
  EndIf;
  NewA:=[]; 
  For K:=1 To Erre Do 
    NewRow:=A[K]-(A[K,J]/Pivot)*A[I];-- NewRow;  -- If K=I it's the zero row Vector 
    Append(NewA,NewRow); -- PrintLn(NewA);
  EndFor; 
  Rows:=Diff(1..Erre,[I]);-- Rows;
  Cols:=Diff(1..Esse,[J]);-- Cols;
  Return Submat(Mat(NewA),Rows,Cols); 
EndDefine; -- GaussReduction 

  
  ------------------------------------------------------------------------- 


End;   -- Package