\documentclass{article} \usepackage{open-axiom} \begin{document} \title{\$SPAD/src/algebra solverad.spad} \author{Patrizia Gianni} \maketitle \begin{abstract} \end{abstract} \eject \tableofcontents \eject \section{package SOLVERAD RadicalSolvePackage} <>= )abbrev package SOLVERAD RadicalSolvePackage ++ Author: P.Gianni ++ Date Created: Summer 1990 ++ Date Last Updated: October 1991 ++ Basic Functions: ++ Related Constructors: SystemSolvePackage, FloatingRealPackage, ++ FloatingComplexPackage ++ Also See: ++ AMS Classifications: ++ Keywords: ++ References: ++ Description: ++ This package tries to find solutions ++ expressed in terms of radicals for systems of equations ++ of rational functions with coefficients in an integral domain R. RadicalSolvePackage(R): Cat == Capsule where R : Join(EuclideanDomain, CharacteristicZero) PI ==> PositiveInteger NNI==> NonNegativeInteger Z ==> Integer B ==> Boolean ST ==> String PR ==> Polynomial R UP ==> SparseUnivariatePolynomial PR LA ==> LocalAlgebra(PR, Z, Z) RF ==> Fraction PR RE ==> Expression R EQ ==> Equation SY ==> Symbol SU ==> SuchThat(List RE, List Equation RE) SUP==> SparseUnivariatePolynomial L ==> List P ==> Polynomial SOLVEFOR ==> PolynomialSolveByFormulas(SUP RE, RE) UPF2 ==> SparseUnivariatePolynomialFunctions2(PR,RE) Cat ==> with radicalSolve : (RF,SY) -> L EQ RE ++ radicalSolve(rf,x) finds the solutions expressed in terms of ++ radicals of the equation rf = 0 with respect to the symbol x, ++ where rf is a rational function. radicalSolve : RF -> L EQ RE ++ radicalSolve(rf) finds the solutions expressed in terms of ++ radicals of the equation rf = 0, where rf is a ++ univariate rational function. radicalSolve : (EQ RF,SY) -> L EQ RE ++ radicalSolve(eq,x) finds the solutions expressed in terms of ++ radicals of the equation of rational functions eq ++ with respect to the symbol x. radicalSolve : EQ RF -> L EQ RE ++ radicalSolve(eq) finds the solutions expressed in terms of ++ radicals of the equation of rational functions eq ++ with respect to the unique symbol x appearing in eq. radicalSolve : (L RF,L SY) -> L L EQ RE ++ radicalSolve(lrf,lvar) finds the solutions expressed in terms of ++ radicals of the system of equations lrf = 0 with ++ respect to the list of symbols lvar, ++ where lrf is a list of rational functions. radicalSolve : L RF -> L L EQ RE ++ radicalSolve(lrf) finds the solutions expressed in terms of ++ radicals of the system of equations lrf = 0, where lrf is a ++ system of univariate rational functions. radicalSolve : (L EQ RF,L SY) -> L L EQ RE ++ radicalSolve(leq,lvar) finds the solutions expressed in terms of ++ radicals of the system of equations of rational functions leq ++ with respect to the list of symbols lvar. radicalSolve : L EQ RF -> L L EQ RE ++ radicalSolve(leq) finds the solutions expressed in terms of ++ radicals of the system of equations of rational functions leq ++ with respect to the unique symbol x appearing in leq. radicalRoots : (RF,SY) -> L RE ++ radicalRoots(rf,x) finds the roots expressed in terms of radicals ++ of the rational function rf with respect to the symbol x. radicalRoots : (L RF,L SY) -> L L RE ++ radicalRoots(lrf,lvar) finds the roots expressed in terms of ++ radicals of the list of rational functions lrf ++ with respect to the list of symbols lvar. contractSolve: (EQ RF,SY) -> SU ++ contractSolve(eq,x) finds the solutions expressed in terms of ++ radicals of the equation of rational functions eq ++ with respect to the symbol x. The result contains new ++ symbols for common subexpressions in order to reduce the ++ size of the output. contractSolve: (RF,SY) -> SU ++ contractSolve(rf,x) finds the solutions expressed in terms of ++ radicals of the equation rf = 0 with respect to the symbol x, ++ where rf is a rational function. The result contains new ++ symbols for common subexpressions in order to reduce the ++ size of the output. Capsule ==> add import DegreeReductionPackage(PR, R) import SOLVEFOR SideEquations: List EQ RE := [] ContractSoln: B := false ---- Local Function Declarations ---- solveInner:(PR, SY, B) -> SU linear: UP -> List RE quadratic: UP -> List RE cubic: UP -> List RE quartic: UP -> List RE rad: PI -> RE wrap: RE -> RE New: RE -> RE makeEq : (List RE,L SY) -> L EQ RE select : L L RE -> L L RE isGeneric? : (L PR,L SY) -> Boolean findGenZeros : (L PR,L SY) -> L L RE findZeros : (L PR,L SY) -> L L RE New s == s = 0 => 0 S := new()$Symbol ::PR::RF::RE SideEquations := append([S = s], SideEquations) S linear u == [(-coefficient(u,0))::RE /(coefficient(u,1))::RE] quadratic u == quadratic(map(coerce,u)$UPF2)$SOLVEFOR cubic u == cubic(map(coerce,u)$UPF2)$SOLVEFOR quartic u == quartic(map(coerce,u)$UPF2)$SOLVEFOR rad n == n::Z::RE wrap s == (ContractSoln => New s; s) ---- Exported Functions ---- -- find the zeros of components in "generic" position -- findGenZeros(rlp:L PR,rlv:L SY) : L L RE == pp:=rlp.first v:=first rlv rlv:=rest rlv res:L L RE:=[] res:=append([reverse cons(r,[eval( (-coefficient(univariate(p,vv),0)::RE)/(leadingCoefficient univariate(p,vv))::RE, kernel(v)@Kernel(RE),r) for vv in rlv for p in rlp.rest]) for r in radicalRoots(pp::RF,v)],res) res findZeros(rlp:L PR,rlv:L SY) : L L RE == parRes:=[radicalRoots(p::RF,v) for p in rlp for v in rlv] parRes:=select parRes res:L L RE :=[] res1:L RE for par in parRes repeat res1:=[par.first] lv1:L Kernel(RE):=[kernel rlv.first] rlv1:=rlv.rest p1:=par.rest while p1~=[] repeat res1:=cons(eval(p1.first,lv1,res1),res1) p1:=p1.rest lv1:=cons(kernel rlv1.first,lv1) rlv1:=rlv1.rest res:=cons(res1,res) res radicalSolve(pol:RF,v:SY) == [equation(v::RE,r) for r in radicalRoots(pol,v)] radicalSolve(p:RF) == zero? p => error "equation is always satisfied" lv:=removeDuplicates concat(variables numer p, variables denom p) empty? lv => error "inconsistent equation" #lv>1 => error "too many variables" radicalSolve(p,lv.first) radicalSolve(eq: EQ RF) == radicalSolve(lhs eq -rhs eq) radicalSolve(eq: EQ RF,v:SY) == radicalSolve(lhs eq - rhs eq,v) radicalRoots(lp: L RF,lv: L SY) == parRes:=triangularSystems(lp,lv)$SystemSolvePackage(R) parRes= list [] => [] -- select the components in "generic" form rlv:=reverse lv rpRes:=[reverse res for res in parRes] listGen:= [res for res in rpRes|isGeneric?(res,rlv)] result:L L RE:=[] if listGen~=[] then result:="append"/[findGenZeros(res,rlv) for res in listGen] for res in listGen repeat rpRes:=delete(rpRes,position(res,rpRes)) -- non-generic components rpRes = [] => result append("append"/[findZeros(res,rlv) for res in rpRes], result) radicalSolve(lp:L RF,lv:L SY) == [makeEq(lres,lv) for lres in radicalRoots(lp,lv)] radicalSolve(lp: L RF) == lv:="setUnion"/[setUnion(variables numer p,variables denom p) for p in lp] [makeEq(lres,lv) for lres in radicalRoots(lp,lv)] radicalSolve(le:L EQ RF,lv:L SY) == lp:=[rhs p -lhs p for p in le] [makeEq(lres,lv) for lres in radicalRoots(lp,lv)] radicalSolve(le: L EQ RF) == lp:=[rhs p -lhs p for p in le] lv:="setUnion"/[setUnion(variables numer p,variables denom p) for p in lp] [makeEq(lres,lv) for lres in radicalRoots(lp,lv)] contractSolve(eq:EQ RF, v:SY)== solveInner(numer(lhs eq - rhs eq), v, true) contractSolve(pq:RF, v:SY) == solveInner(numer pq, v, true) radicalRoots(pq:RF, v:SY) == lhs solveInner(numer pq, v, false) -- test if the ideal is radical in generic position -- isGeneric?(rlp:L PR,rlv:L SY) : Boolean == "and"/[degree(f,x)=1 for f in rest rlp for x in rest rlv] ---- select the univariate factors select(lp:L L RE) : L L RE == lp=[] => list [] [:[cons(f,lsel) for lsel in select lp.rest] for f in lp.first] ---- Local Functions ---- -- construct the equation makeEq(nres:L RE,lv:L SY) : L EQ RE == [equation(x :: RE,r) for x in lv for r in nres] solveInner(pq:PR,v:SY,contractFlag:B) == SideEquations := [] ContractSoln := contractFlag factors:= factors (factor pq)$MultivariateFactorize(SY,IndexedExponents SY,R,PR) constants: List PR := [] unsolved: List PR := [] solutions: List RE := [] for f in factors repeat ff:=f.factor not member?(v, variables (ff)) => constants := cons(ff, constants) u := univariate(ff, v) t := reduce u u := t.pol n := degree u l: List RE := n = 1 => linear u n = 2 => quadratic u n = 3 => cubic u n = 4 => quartic u unsolved := cons(ff, unsolved) [] for s in l repeat if t.deg > 1 then s := wrap s T0 := expand(s, t.deg) for i in 1..f.exponent repeat solutions := append(T0, solutions) re := SideEquations [solutions, SideEquations]$SU @ \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}