\documentclass{article} \usepackage{open-axiom} \begin{document} \title{\$SPAD/src/algebra intclos.spad} \author{Victor Miller, Barry Trager, Clifton Williamson} \maketitle \begin{abstract} \end{abstract} \eject \tableofcontents \eject \section{package TRIMAT TriangularMatrixOperations} <>= )abbrev package TRIMAT TriangularMatrixOperations ++ Fraction free inverses of triangular matrices ++ Author: Victor Miller ++ Date Created: ++ Date Last Updated: 24 Jul 1990 ++ Keywords: ++ Examples: ++ References: ++ Description: ++ This package provides functions that compute "fraction-free" ++ inverses of upper and lower triangular matrices over a integral ++ domain. By "fraction-free inverses" we mean the following: ++ given a matrix B with entries in R and an element d of R such that ++ d * inv(B) also has entries in R, we return d * inv(B). Thus, ++ it is not necessary to pass to the quotient field in any of our ++ computations. TriangularMatrixOperations(R,Row,Col,M): Exports == Implementation where R : IntegralDomain Row : FiniteLinearAggregate R Col : FiniteLinearAggregate R M : MatrixCategory(R,Row,Col) Exports ==> with UpTriBddDenomInv: (M,R) -> M ++ UpTriBddDenomInv(B,d) returns M, where ++ B is a non-singular upper triangular matrix and d is an ++ element of R such that \spad{M = d * inv(B)} has entries in R. LowTriBddDenomInv:(M,R) -> M ++ LowTriBddDenomInv(B,d) returns M, where ++ B is a non-singular lower triangular matrix and d is an ++ element of R such that \spad{M = d * inv(B)} has entries in R. Implementation ==> add UpTriBddDenomInv(A,denom) == AI := zero(nrows A, nrows A)$M offset := minColIndex AI - minRowIndex AI for i in minRowIndex AI .. maxRowIndex AI for j in minColIndex AI .. maxColIndex AI repeat qsetelt!(AI,i,j,(denom exquo qelt(A,i,j))::R) for i in minRowIndex AI .. maxRowIndex AI repeat for j in offset + i + 1 .. maxColIndex AI repeat qsetelt!(AI,i,j, - (((+/[qelt(AI,i,k) * qelt(A,k-offset,j) for k in i+offset..(j-1)]) exquo qelt(A, j-offset, j))::R)) AI LowTriBddDenomInv(A, denom) == AI := zero(nrows A, nrows A)$M offset := minColIndex AI - minRowIndex AI for i in minRowIndex AI .. maxRowIndex AI for j in minColIndex AI .. maxColIndex AI repeat qsetelt!(AI,i,j,(denom exquo qelt(A,i,j))::R) for i in minColIndex AI .. maxColIndex AI repeat for j in i - offset + 1 .. maxRowIndex AI repeat qsetelt!(AI,j,i, - (((+/[qelt(A,j,k+offset) * qelt(AI,k,i) for k in i-offset..(j-1)]) exquo qelt(A, j, j+offset))::R)) AI @ \section{package IBATOOL IntegralBasisTools} <>= )abbrev package IBATOOL IntegralBasisTools ++ Functions common to both integral basis packages ++ Author: Victor Miller, Barry Trager, Clifton Williamson ++ Date Created: 11 April 1990 ++ Date Last Updated: 20 September 1994 ++ Keywords: integral basis, function field, number field ++ Examples: ++ References: ++ Description: ++ This package contains functions used in the packages ++ FunctionFieldIntegralBasis and NumberFieldIntegralBasis. IntegralBasisTools(R,UP,F): Exports == Implementation where R : EuclideanDomain with squareFree: $ -> Factored $ ++ squareFree(x) returns a square-free factorisation of x UP : UnivariatePolynomialCategory R F : FramedAlgebra(R,UP) Mat ==> Matrix R NNI ==> NonNegativeInteger Ans ==> Record(basis: Mat, basisDen: R, basisInv:Mat) Exports ==> with diagonalProduct: Mat -> R ++ diagonalProduct(m) returns the product of the elements on the ++ diagonal of the matrix m matrixGcd: (Mat,R,NNI) -> R ++ matrixGcd(mat,sing,n) is \spad{gcd(sing,g)} where \spad{g} is the ++ gcd of the entries of the \spad{n}-by-\spad{n} upper-triangular ++ matrix \spad{mat}. divideIfCan!: (Matrix R,Matrix R,R,Integer) -> R ++ divideIfCan!(matrix,matrixOut,prime,n) attempts to divide the ++ entries of \spad{matrix} by \spad{prime} and store the result in ++ \spad{matrixOut}. If it is successful, 1 is returned and if not, ++ \spad{prime} is returned. Here both \spad{matrix} and ++ \spad{matrixOut} are \spad{n}-by-\spad{n} upper triangular matrices. leastPower: (NNI,NNI) -> NNI ++ leastPower(p,n) returns e, where e is the smallest integer ++ such that \spad{p **e >= n} idealiser: (Mat,Mat) -> Mat ++ idealiser(m1,m2) computes the order of an ideal defined by m1 and m2 idealiser: (Mat,Mat,R) -> Mat ++ idealiser(m1,m2,d) computes the order of an ideal defined by m1 and m2 ++ where d is the known part of the denominator idealiserMatrix: (Mat, Mat) -> Mat ++ idealiserMatrix(m1, m2) returns the matrix representing the linear ++ conditions on the Ring associatied with an ideal defined by m1 and m2. moduleSum: (Ans,Ans) -> Ans ++ moduleSum(m1,m2) returns the sum of two modules in the framed ++ algebra \spad{F}. Each module \spad{mi} is represented as follows: ++ F is a framed algebra with R-module basis \spad{w1,w2,...,wn} and ++ \spad{mi} is a record \spad{[basis,basisDen,basisInv]}. If ++ \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ a basis \spad{v1,...,vn} for \spad{mi} is given by ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of 'basis' contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. Implementation ==> add import ModularHermitianRowReduction(R) import TriangularMatrixOperations(R, Vector R, Vector R, Matrix R) diagonalProduct m == ans : R := 1 for i in minRowIndex m .. maxRowIndex m for j in minColIndex m .. maxColIndex m repeat ans := ans * qelt(m, i, j) ans matrixGcd(mat,sing,n) == -- note: 'matrix' is upper triangular; -- no need to do anything below the diagonal d := sing for i in 1..n repeat for j in i..n repeat if not zero?(mij := qelt(mat,i,j)) then d := gcd(d,mij) one? d => return d d divideIfCan!(matrix,matrixOut,prime,n) == -- note: both 'matrix' and 'matrixOut' will be upper triangular; -- no need to do anything below the diagonal for i in 1..n repeat for j in i..n repeat (a := (qelt(matrix,i,j) exquo prime)) case "failed" => return prime qsetelt!(matrixOut,i,j,a :: R) 1 leastPower(p,n) == -- efficiency is not an issue here e : NNI := 1; q := p while q < n repeat (e := e + 1; q := q * p) e idealiserMatrix(ideal,idealinv) == -- computes the Order of the ideal n := rank()$F bigm := zero(n * n,n)$Mat mr := minRowIndex bigm; mc := minColIndex bigm v := basis()$F for i in 0..n-1 repeat r := regularRepresentation qelt(v,i + minIndex v) m := ideal * r * idealinv for j in 0..n-1 repeat for k in 0..n-1 repeat bigm(j * n + k + mr,i + mc) := qelt(m,j + mr,k + mc) bigm idealiser(ideal,idealinv) == bigm := idealiserMatrix(ideal, idealinv) transpose squareTop rowEch bigm idealiser(ideal,idealinv,denom) == bigm := (idealiserMatrix(ideal, idealinv) exquo denom)::Mat transpose squareTop rowEchelon(bigm,denom) moduleSum(mod1,mod2) == rb1 := mod1.basis; rbden1 := mod1.basisDen; rbinv1 := mod1.basisInv rb2 := mod2.basis; rbden2 := mod2.basisDen; rbinv2 := mod2.basisInv -- compatibility check: doesn't take much computation time (not square? rb1) or (not square? rbinv1) or (not square? rb2) _ or (not square? rbinv2) => error "moduleSum: matrices must be square" ((n := nrows rb1) ~= (nrows rbinv1)) or (n ~= (nrows rb2)) _ or (n ~= (nrows rbinv2)) => error "moduleSum: matrices of imcompatible dimensions" (zero? rbden1) or (zero? rbden2) => error "moduleSum: denominator must be non-zero" den := lcm(rbden1,rbden2); c1 := den quo rbden1; c2 := den quo rbden2 rb := squareTop rowEchelon(vertConcat(c1 * rb1,c2 * rb2),den) rbinv := UpTriBddDenomInv(rb,den) [rb,den,rbinv] @ \section{package FFINTBAS FunctionFieldIntegralBasis} <>= )abbrev package FFINTBAS FunctionFieldIntegralBasis ++ Integral bases for function fields of dimension one ++ Author: Victor Miller ++ Date Created: 9 April 1990 ++ Date Last Updated: 20 September 1994 ++ Keywords: ++ Examples: ++ References: ++ Description: ++ In this package R is a Euclidean domain and F is a framed algebra ++ over R. The package provides functions to compute the integral ++ closure of R in the quotient field of F. It is assumed that ++ \spad{char(R/P) = char(R)} for any prime P of R. A typical instance of ++ this is when \spad{R = K[x]} and F is a function field over R. FunctionFieldIntegralBasis(R,UP,F): Exports == Implementation where R : EuclideanDomain with squareFree: $ -> Factored $ ++ squareFree(x) returns a square-free factorisation of x UP : UnivariatePolynomialCategory R F : FramedAlgebra(R,UP) I ==> Integer Mat ==> Matrix R NNI ==> NonNegativeInteger Exports ==> with integralBasis : () -> Record(basis: Mat, basisDen: R, basisInv:Mat) ++ \spad{integralBasis()} returns a record ++ \spad{[basis,basisDen,basisInv]} containing information regarding ++ the integral closure of R in the quotient field of F, where ++ F is a framed algebra with R-module basis \spad{w1,w2,...,wn}. ++ If \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ the \spad{i}th element of the integral basis is ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of \spad{basis} contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. localIntegralBasis : R -> Record(basis: Mat, basisDen: R, basisInv:Mat) ++ \spad{integralBasis(p)} returns a record ++ \spad{[basis,basisDen,basisInv]} containing information regarding ++ the local integral closure of R at the prime \spad{p} in the quotient ++ field of F, where F is a framed algebra with R-module basis ++ \spad{w1,w2,...,wn}. ++ If \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ the \spad{i}th element of the local integral basis is ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of \spad{basis} contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. Implementation ==> add import IntegralBasisTools(R, UP, F) import ModularHermitianRowReduction(R) import TriangularMatrixOperations(R, Vector R, Vector R, Matrix R) squaredFactors: R -> R squaredFactors px == */[(if ffe.exponent > 1 then ffe.factor else 1$R) for ffe in factors squareFree px] iIntegralBasis: (Mat,R,R) -> Record(basis: Mat, basisDen: R, basisInv:Mat) iIntegralBasis(tfm,disc,sing) == -- tfm = trace matrix of current order n := rank()$F; tfm0 := copy tfm; disc0 := disc rb := scalarMatrix(n, 1); rbinv := scalarMatrix(n, 1) -- rb = basis matrix of current order -- rbinv = inverse basis matrix of current order -- these are wrt the original basis for F rbden : R := 1; index : R := 1; oldIndex : R := 1 -- rbden = denominator for current basis matrix -- index = index of original order in current order not sizeLess?(1, sing) => [rb, rbden, rbinv] repeat -- compute the p-radical idinv := transpose squareTop rowEchelon(tfm, sing) -- [u1,..,un] are the coordinates of an element of the p-radical -- iff [u1,..,un] * idinv is in sing * R^n id := rowEchelon LowTriBddDenomInv(idinv, sing) -- id = basis matrix of the p-radical idinv := UpTriBddDenomInv(id, sing) -- id * idinv = sing * identity -- no need to check for inseparability in this case rbinv := idealiser(id * rb, rbinv * idinv, sing * rbden) index := diagonalProduct rbinv rb := rowEchelon LowTriBddDenomInv(rbinv, rbden * sing) g := matrixGcd(rb,sing,n) if sizeLess?(1,g) then rb := (rb exquo g) :: Mat rbden := rbden * (sing quo g) rbinv := UpTriBddDenomInv(rb, rbden) disc := disc0 quo (index * index) indexChange := index quo oldIndex; oldIndex := index sing := gcd(indexChange, squaredFactors disc) not sizeLess?(1, sing) => return [rb, rbden, rbinv] tfm := ((rb * tfm0 * transpose rb) exquo (rbden * rbden)) :: Mat integralBasis() == n := rank()$F; p := characteristic$F (not zero? p) and (n >= p) => error "integralBasis: possible wild ramification" tfm := traceMatrix()$F; disc := determinant tfm sing := squaredFactors disc -- singularities of relative Spec iIntegralBasis(tfm,disc,sing) localIntegralBasis prime == n := rank()$F; p := characteristic$F (not zero? p) and (n >= p) => error "integralBasis: possible wild ramification" tfm := traceMatrix()$F; disc := determinant tfm (disc exquo (prime * prime)) case "failed" => [scalarMatrix(n,1),1,scalarMatrix(n,1)] iIntegralBasis(tfm,disc,prime) @ \section{package WFFINTBS WildFunctionFieldIntegralBasis} <>= )abbrev package WFFINTBS WildFunctionFieldIntegralBasis ++ Authors: Victor Miller, Clifton Williamson ++ Date Created: 24 July 1991 ++ Date Last Updated: 20 September 1994 ++ Basic Operations: integralBasis, localIntegralBasis ++ Related Domains: IntegralBasisTools(R,UP,F), ++ TriangularMatrixOperations(R,Vector R,Vector R,Matrix R) ++ Also See: FunctionFieldIntegralBasis, NumberFieldIntegralBasis ++ AMS Classifications: ++ Keywords: function field, integral basis ++ Examples: ++ References: ++ Description: ++ In this package K is a finite field, R is a ring of univariate ++ polynomials over K, and F is a framed algebra over R. The package ++ provides a function to compute the integral closure of R in the quotient ++ field of F as well as a function to compute a "local integral basis" ++ at a specific prime. WildFunctionFieldIntegralBasis(K,R,UP,F): Exports == Implementation where K : FiniteFieldCategory --K : Join(Field,Finite) R : UnivariatePolynomialCategory K UP : UnivariatePolynomialCategory R F : FramedAlgebra(R,UP) I ==> Integer Mat ==> Matrix R NNI ==> NonNegativeInteger SAE ==> SimpleAlgebraicExtension RResult ==> Record(basis: Mat, basisDen: R, basisInv:Mat) IResult ==> Record(basis: Mat, basisDen: R, basisInv:Mat,discr: R) MATSTOR ==> StorageEfficientMatrixOperations Exports ==> with integralBasis : () -> RResult ++ \spad{integralBasis()} returns a record ++ \spad{[basis,basisDen,basisInv]} containing information regarding ++ the integral closure of R in the quotient field of F, where ++ F is a framed algebra with R-module basis \spad{w1,w2,...,wn}. ++ If \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ the \spad{i}th element of the integral basis is ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of \spad{basis} contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. localIntegralBasis : R -> RResult ++ \spad{integralBasis(p)} returns a record ++ \spad{[basis,basisDen,basisInv]} containing information regarding ++ the local integral closure of R at the prime \spad{p} in the quotient ++ field of F, where F is a framed algebra with R-module basis ++ \spad{w1,w2,...,wn}. ++ If \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ the \spad{i}th element of the local integral basis is ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of \spad{basis} contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. Implementation ==> add import IntegralBasisTools(R, UP, F) import ModularHermitianRowReduction(R) import TriangularMatrixOperations(R, Vector R, Vector R, Matrix R) import DistinctDegreeFactorize(K,R) listSquaredFactors: R -> List R listSquaredFactors px == -- returns a list of the factors of px which occur with -- exponent > 1 ans : List R := empty() factored := factor(px)$DistinctDegreeFactorize(K,R) for f in factors(factored) repeat if f.exponent > 1 then ans := concat(f.factor,ans) ans iLocalIntegralBasis: (Vector F,Vector F,Matrix R,Matrix R,R,R) -> IResult iLocalIntegralBasis(bas,pows,tfm,matrixOut,disc,prime) == n := rank()$F; standardBasis := basis()$F -- 'standardBasis' is the basis for F as a FramedAlgebra; -- usually this is [1,y,y**2,...,y**(n-1)] p2 := prime * prime; sae := SAE(K,R,prime) p := characteristic$F; q := size()$sae lp := leastPower(q,n) rb := scalarMatrix(n,1); rbinv := scalarMatrix(n,1) -- rb = basis matrix of current order -- rbinv = inverse basis matrix of current order -- these are wrt the orginal basis for F rbden : R := 1; index : R := 1; oldIndex : R := 1 -- rbden = denominator for current basis matrix -- index = index of original order in current order repeat -- pows = [(w1 * rbden) ** q,...,(wn * rbden) ** q], where -- bas = [w1,...,wn] is 'rbden' times the basis for the order B = 'rb' for i in 1..n repeat bi : F := 0 for j in 1..n repeat bi := bi + qelt(rb,i,j) * qelt(standardBasis,j) qsetelt!(bas,i,bi) qsetelt!(pows,i,bi ** p) coor0 := transpose coordinates(pows,bas) denPow := rbden ** ((p - 1) :: NNI) (coMat0 := coor0 exquo denPow) case "failed" => error "can't happen" -- the jth column of coMat contains the coordinates of (wj/rbden)**q -- with respect to the basis [w1/rbden,...,wn/rbden] coMat := coMat0 :: Matrix R -- the ith column of 'pPows' contains the coordinates of the pth power -- of the ith basis element for B/prime.B over 'sae' = R/prime.R pPows := map(reduce,coMat)$MatrixCategoryFunctions2(R,Vector R, Vector R,Matrix R,sae,Vector sae,Vector sae,Matrix sae) -- 'frob' will eventually be the Frobenius matrix for B/prime.B over -- 'sae' = R/prime.R; at each stage of the loop the ith column will -- contain the coordinates of p^k-th powers of the ith basis element frob := copy pPows; tmpMat : Matrix sae := new(n,n,0) for r in 2..leastPower(p,q) repeat for i in 1..n repeat for j in 1..n repeat qsetelt!(tmpMat,i,j,qelt(frob,i,j) ** p) times!(frob,pPows,tmpMat)$MATSTOR(sae) frobPow := frob ** lp -- compute the p-radical ns := nullSpace frobPow for i in 1..n repeat for j in 1..n repeat qsetelt!(tfm,i,j,0) for vec in ns for i in 1.. repeat for j in 1..n repeat qsetelt!(tfm,i,j,lift qelt(vec,j)) id := squareTop rowEchelon(tfm,prime) -- id = basis matrix of the p-radical idinv := UpTriBddDenomInv(id, prime) -- id * idinv = prime * identity -- no need to check for inseparability in this case rbinv := idealiser(id * rb, rbinv * idinv, prime * rbden) index := diagonalProduct rbinv rb := rowEchelon LowTriBddDenomInv(rbinv,rbden * prime) if divideIfCan!(rb,matrixOut,prime,n) = 1 then rb := matrixOut else rbden := rbden * prime rbinv := UpTriBddDenomInv(rb,rbden) indexChange := index quo oldIndex oldIndex := index disc := disc quo (indexChange * indexChange) (not sizeLess?(1,indexChange)) or ((disc exquo p2) case "failed") => return [rb, rbden, rbinv, disc] integralBasis() == traceMat := traceMatrix()$F; n := rank()$F disc := determinant traceMat -- discriminant of current order zero? disc => error "integralBasis: polynomial must be separable" singList := listSquaredFactors disc -- singularities of relative Spec runningRb := scalarMatrix(n,1); runningRbinv := scalarMatrix(n,1) -- runningRb = basis matrix of current order -- runningRbinv = inverse basis matrix of current order -- these are wrt the original basis for F runningRbden : R := 1 -- runningRbden = denominator for current basis matrix empty? singList => [runningRb, runningRbden, runningRbinv] bas : Vector F := new(n,0); pows : Vector F := new(n,0) -- storage for basis elements and their powers tfm : Matrix R := new(n,n,0) -- 'tfm' will contain the coordinates of a lifting of the kernel -- of a power of Frobenius matrixOut : Matrix R := new(n,n,0) for prime in singList repeat lb := iLocalIntegralBasis(bas,pows,tfm,matrixOut,disc,prime) rb := lb.basis; rbinv := lb.basisInv; rbden := lb.basisDen disc := lb.discr -- update 'running integral basis' if newly computed -- local integral basis is non-trivial if sizeLess?(1,rbden) then mat := vertConcat(rbden * runningRb,runningRbden * rb) runningRbden := runningRbden * rbden runningRb := squareTop rowEchelon(mat,runningRbden) runningRbinv := UpTriBddDenomInv(runningRb,runningRbden) [runningRb, runningRbden, runningRbinv] localIntegralBasis prime == traceMat := traceMatrix()$F; n := rank()$F disc := determinant traceMat -- discriminant of current order zero? disc => error "localIntegralBasis: polynomial must be separable" (disc exquo (prime * prime)) case "failed" => [scalarMatrix(n,1), 1, scalarMatrix(n,1)] bas : Vector F := new(n,0); pows : Vector F := new(n,0) -- storage for basis elements and their powers tfm : Matrix R := new(n,n,0) -- 'tfm' will contain the coordinates of a lifting of the kernel -- of a power of Frobenius matrixOut : Matrix R := new(n,n,0) lb := iLocalIntegralBasis(bas,pows,tfm,matrixOut,disc,prime) [lb.basis, lb.basisDen, lb.basisInv] @ \section{package NFINTBAS NumberFieldIntegralBasis} <>= )abbrev package NFINTBAS NumberFieldIntegralBasis ++ Author: Victor Miller, Clifton Williamson ++ Date Created: 9 April 1990 ++ Date Last Updated: 20 September 1994 ++ Basic Operations: discriminant, integralBasis ++ Related Domains: IntegralBasisTools, TriangularMatrixOperations ++ Also See: FunctionFieldIntegralBasis, WildFunctionFieldIntegralBasis ++ AMS Classifications: ++ Keywords: number field, integral basis, discriminant ++ Examples: ++ References: ++ Description: ++ In this package F is a framed algebra over the integers (typically ++ \spad{F = Z[a]} for some algebraic integer a). The package provides ++ functions to compute the integral closure of Z in the quotient ++ quotient field of F. NumberFieldIntegralBasis(UP,F): Exports == Implementation where UP : UnivariatePolynomialCategory Integer F : FramedAlgebra(Integer,UP) FR ==> Factored Integer I ==> Integer Mat ==> Matrix I NNI ==> NonNegativeInteger Ans ==> Record(basis: Mat, basisDen: I, basisInv:Mat,discr: I) Exports ==> with discriminant: () -> Integer ++ \spad{discriminant()} returns the discriminant of the integral ++ closure of Z in the quotient field of the framed algebra F. integralBasis : () -> Record(basis: Mat, basisDen: I, basisInv:Mat) ++ \spad{integralBasis()} returns a record ++ \spad{[basis,basisDen,basisInv]} ++ containing information regarding the integral closure of Z in the ++ quotient field of F, where F is a framed algebra with Z-module ++ basis \spad{w1,w2,...,wn}. ++ If \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ the \spad{i}th element of the integral basis is ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of \spad{basis} contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. localIntegralBasis : I -> Record(basis: Mat, basisDen: I, basisInv:Mat) ++ \spad{integralBasis(p)} returns a record ++ \spad{[basis,basisDen,basisInv]} containing information regarding ++ the local integral closure of Z at the prime \spad{p} in the quotient ++ field of F, where F is a framed algebra with Z-module basis ++ \spad{w1,w2,...,wn}. ++ If \spad{basis} is the matrix \spad{(aij, i = 1..n, j = 1..n)}, then ++ the \spad{i}th element of the integral basis is ++ \spad{vi = (1/basisDen) * sum(aij * wj, j = 1..n)}, i.e. the ++ \spad{i}th row of \spad{basis} contains the coordinates of the ++ \spad{i}th basis vector. Similarly, the \spad{i}th row of the ++ matrix \spad{basisInv} contains the coordinates of \spad{wi} with ++ respect to the basis \spad{v1,...,vn}: if \spad{basisInv} is the ++ matrix \spad{(bij, i = 1..n, j = 1..n)}, then ++ \spad{wi = sum(bij * vj, j = 1..n)}. Implementation ==> add import IntegralBasisTools(I, UP, F) import ModularHermitianRowReduction(I) import TriangularMatrixOperations(I, Vector I, Vector I, Matrix I) frobMatrix : (Mat,Mat,I,NNI) -> Mat wildPrimes : (FR,I) -> List I tameProduct : (FR,I) -> I iTameLocalIntegralBasis : (Mat,I,I) -> Ans iWildLocalIntegralBasis : (Mat,I,I) -> Ans frobMatrix(rb,rbinv,rbden,p) == n := rank()$F; b := basis()$F v : Vector F := new(n,0) for i in minIndex(v)..maxIndex(v) for ii in minRowIndex(rb)..maxRowIndex(rb) repeat a : F := 0 for j in minIndex(b)..maxIndex(b) for jj in minColIndex(rb)..maxColIndex(rb) repeat a := a + qelt(rb,ii,jj) * qelt(b,j) qsetelt!(v,i,a**p) mat := transpose coordinates v ((transpose(rbinv) * mat) exquo (rbden ** p)) :: Mat wildPrimes(factoredDisc,n) == -- returns a list of the primes <=n which divide factoredDisc to a -- power greater than 1 ans : List I := empty() for f in factors(factoredDisc) repeat if f.exponent > 1 and f.factor <= n then ans := concat(f.factor,ans) ans tameProduct(factoredDisc,n) == -- returns the product of the primes > n which divide factoredDisc -- to a power greater than 1 ans : I := 1 for f in factors(factoredDisc) repeat if f.exponent > 1 and f.factor > n then ans := f.factor * ans ans integralBasis() == traceMat := traceMatrix()$F; n := rank()$F disc := determinant traceMat -- discriminant of current order disc0 := disc -- this is disc(F) factoredDisc := factor(disc0)$IntegerFactorizationPackage(Integer) wilds := wildPrimes(factoredDisc,n) sing := tameProduct(factoredDisc,n) runningRb := scalarMatrix(n, 1); runningRbinv := scalarMatrix(n, 1) -- runningRb = basis matrix of current order -- runningRbinv = inverse basis matrix of current order -- these are wrt the original basis for F runningRbden : I := 1 -- runningRbden = denominator for current basis matrix one? sing and empty? wilds => [runningRb, runningRbden, runningRbinv] -- id = basis matrix of the ideal (p-radical) wrt current basis matrixOut : Mat := scalarMatrix(n,0) for p in wilds repeat lb := iWildLocalIntegralBasis(matrixOut,disc,p) rb := lb.basis; rbinv := lb.basisInv; rbden := lb.basisDen disc := lb.discr -- update 'running integral basis' if newly computed -- local integral basis is non-trivial if sizeLess?(1,rbden) then mat := vertConcat(rbden * runningRb,runningRbden * rb) runningRbden := runningRbden * rbden runningRb := squareTop rowEchelon(mat,runningRbden) runningRbinv := UpTriBddDenomInv(runningRb,runningRbden) lb := iTameLocalIntegralBasis(traceMat,disc,sing) rb := lb.basis; rbinv := lb.basisInv; rbden := lb.basisDen disc := lb.discr -- update 'running integral basis' if newly computed -- local integral basis is non-trivial if sizeLess?(1,rbden) then mat := vertConcat(rbden * runningRb,runningRbden * rb) runningRbden := runningRbden * rbden runningRb := squareTop rowEchelon(mat,runningRbden) runningRbinv := UpTriBddDenomInv(runningRb,runningRbden) [runningRb,runningRbden,runningRbinv] localIntegralBasis p == traceMat := traceMatrix()$F; n := rank()$F disc := determinant traceMat -- discriminant of current order (disc exquo (p*p)) case "failed" => [scalarMatrix(n, 1), 1, scalarMatrix(n, 1)] lb := p > rank()$F => iTameLocalIntegralBasis(traceMat,disc,p) iWildLocalIntegralBasis(scalarMatrix(n,0),disc,p) [lb.basis,lb.basisDen,lb.basisInv] iTameLocalIntegralBasis(traceMat,disc,sing) == n := rank()$F; disc0 := disc rb := scalarMatrix(n, 1); rbinv := scalarMatrix(n, 1) -- rb = basis matrix of current order -- rbinv = inverse basis matrix of current order -- these are wrt the original basis for F rbden : I := 1; index : I := 1; oldIndex : I := 1 -- rbden = denominator for current basis matrix -- id = basis matrix of the ideal (p-radical) wrt current basis tfm := traceMat repeat -- compute the p-radical = p-trace-radical idinv := transpose squareTop rowEchelon(tfm,sing) -- [u1,..,un] are the coordinates of an element of the p-radical -- iff [u1,..,un] * idinv is in p * Z^n id := rowEchelon LowTriBddDenomInv(idinv, sing) -- id = basis matrix of the p-radical idinv := UpTriBddDenomInv(id, sing) -- id * idinv = sing * identity -- no need to check for inseparability in this case rbinv := idealiser(id * rb, rbinv * idinv, sing * rbden) index := diagonalProduct rbinv rb := rowEchelon LowTriBddDenomInv(rbinv, sing * rbden) g := matrixGcd(rb,sing,n) if sizeLess?(1,g) then rb := (rb exquo g) :: Mat rbden := rbden * (sing quo g) rbinv := UpTriBddDenomInv(rb, rbden) disc := disc0 quo (index * index) indexChange := index quo oldIndex; oldIndex := index one? indexChange => return [rb, rbden, rbinv, disc] tfm := ((rb * traceMat * transpose rb) exquo (rbden * rbden)) :: Mat iWildLocalIntegralBasis(matrixOut,disc,p) == n := rank()$F; disc0 := disc rb := scalarMatrix(n, 1); rbinv := scalarMatrix(n, 1) -- rb = basis matrix of current order -- rbinv = inverse basis matrix of current order -- these are wrt the original basis for F rbden : I := 1; index : I := 1; oldIndex : I := 1 -- rbden = denominator for current basis matrix -- id = basis matrix of the ideal (p-radical) wrt current basis p2 := p * p; lp := leastPower(p::NNI,n) repeat tfm := frobMatrix(rb,rbinv,rbden,p::NNI) ** lp -- compute Rp = p-radical idinv := transpose squareTop rowEchelon(tfm, p) -- [u1,..,un] are the coordinates of an element of Rp -- iff [u1,..,un] * idinv is in p * Z^n id := rowEchelon LowTriBddDenomInv(idinv,p) -- id = basis matrix of the p-radical idinv := UpTriBddDenomInv(id,p) -- id * idinv = p * identity -- no need to check for inseparability in this case rbinv := idealiser(id * rb, rbinv * idinv, p * rbden) index := diagonalProduct rbinv rb := rowEchelon LowTriBddDenomInv(rbinv, p * rbden) if divideIfCan!(rb,matrixOut,p,n) = 1 then rb := matrixOut else rbden := p * rbden rbinv := UpTriBddDenomInv(rb, rbden) indexChange := index quo oldIndex; oldIndex := index disc := disc quo (indexChange * indexChange) one? indexChange or gcd(p2,disc) ~= p2 => return [rb, rbden, rbinv, disc] discriminant() == disc := determinant traceMatrix()$F intBas := integralBasis() rb := intBas.basis; rbden := intBas.basisDen index := ((rbden ** rank()$F) exquo (determinant rb)) :: Integer (disc exquo (index * index)) :: Integer @ \section{License} <>= --Copyright (c) 1991-2002, The Numerical ALgorithms Group Ltd. --All rights reserved. -- --Redistribution and use in source and binary forms, with or without --modification, are permitted provided that the following conditions are --met: -- -- - Redistributions of source code must retain the above copyright -- notice, this list of conditions and the following disclaimer. -- -- - Redistributions in binary form must reproduce the above copyright -- notice, this list of conditions and the following disclaimer in -- the documentation and/or other materials provided with the -- distribution. -- -- - Neither the name of The Numerical ALgorithms Group Ltd. nor the -- names of its contributors may be used to endorse or promote products -- derived from this software without specific prior written permission. -- --THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS --IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED --TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A --PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER --OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, --EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, --PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR --PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF --LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING --NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS --SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. @ <<*>>= <> <> <> <> <> <> @ \eject \begin{thebibliography}{99} \bibitem{1} nothing \end{thebibliography} \end{document}