1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
|
\documentclass{article}
\usepackage{axiom}
\begin{document}
\title{\$SPAD/src/input herm.input}
\author{The Axiom Team}
\maketitle
\begin{abstract}
\end{abstract}
\eject
\tableofcontents
\eject
<<*>>=
-------------------------------- newlodo.input ---------------------------
-- LODO2(M,A) is the domain of linear ordinary differential operators over
-- an A-module M, where A is a differential ring. This includes the
-- cases of operators which are polynomials in D acting upon scalars or
-- vectors depending on a single variable. The coefficients of the
-- operator polynomials can be integers, rational functions, matrices
-- or elements of other domains.
------------------------------------------------------------------------
-- Differential operators with constant coefficients
------------------------------------------------------------------------
)clear all
RN:=FRAC INT
Dx: LODO2(RN, UP(x,RN))
Dx := D() -- definition of an operator
a := Dx + 1
b := a + 1/2*Dx**2 - 1/2
p: UP(x,RN) := 4*x**2 + 2/3 -- something to work on
a p -- application of an operator to a polynomial
(a*b) p = a b p -- multiplication is defined by this identity
c := (1/9)*b*(a + b)**2 -- exponentiation follows from multiplication
(a**2 - 3/4*b + c) (p + 1) -- general application of operator expressions
------------------------------------------------------------------------
-- Differential operators with rational function coefficients
------------------------------------------------------------------------
)clear all
RFZ := FRAC UP(x,INT)
(Dx, a, b): LODO1 RFZ
Dx := D()
b := 3*x**2*Dx**2 + 2*Dx + 1/x
a := b*(5*x*Dx + 7)
p: RFZ := x**2 + 1/x**2
(a*b - b*a) p -- operator multiplication is not commutative
-- When the coefficients of the operator polynomials come from a field
-- it is possible to define left and right division of the operators.
-- This allows the computation of left and right gcd's via remainder
-- sequences, and also the computation of left and right lcm's.
leftDivide(a,b) -- result is the quotient/remainder pair
a - (b * %.quotient + %.remainder)
rightDivide(a,b)
a - (%.quotient * b + %.remainder)
-- A GCD doesn't necessarily divide a and b on both sides.
e := leftGcd(a,b)
leftRemainder(a, e) -- remainder from left division
rightRemainder(a, e) -- remainder from right division
-- An LCM is not necessarily divisible from both sides.
f := rightLcm(a,b)
leftRemainder(f, b)
rightRemainder(f, b) -- the remainder is non-zero
------------------------------------------------------------------------
--
-- Problem: find the first few coefficients of exp(x)/x**i in
-- Dop phi
-- where
-- Dop := D**3 + G/x**2 * D + H/x**3 - 1
-- phi := sum(s[i]*exp(x)/x**i, i = 0..)
------------------------------------------------------------------------
)clear all
Dx: LODO(EXPR INT, f +-> D(f, x))
Dx := D()
Dop:= Dx**3 + G/x**2*Dx + H/x**3 - 1
n == 3
phi == reduce(+,[subscript(s,[i])*exp(x)/x**i for i in 0..n])
phi1 == Dop(phi) / exp x
phi2 == phi1 *x**(n+3)
phi3 == retract(phi2)@(POLY INT)
pans == phi3 ::UP(x,POLY INT)
pans1 == [coefficient(pans, (n+3-i) :: NNI) for i in 2..n+1]
leq == solve(pans1,[subscript(s,[i]) for i in 1..n])
leq
n==4
leq
n==7
leq
------------------------------------------------------------------------
-- Differential operators with matrix coefficients acting on vectors.
------------------------------------------------------------------------
)clear all
PZ := UP(x,INT); Vect := DPMM(3, PZ, SQMATRIX(3,PZ), PZ);
Modo := LODO2(SQMATRIX(3,PZ), Vect);
p := directProduct([3*x**2 + 1, 2*x, 7*x**3 + 2*x]::(VECTOR(PZ)))@Vect
m := [[x**2, 1, 0], [1, x**4, 0], [0, 0, 4*x**2]]::(SQMATRIX(3,PZ))
-- Vect is a left SM(3,PZ)-module
q: Vect := m * p
-- Operator combination and application
Dx: Modo := D()
a: Modo := 1*Dx + m
b: Modo := m*Dx + 1
a*b
a p
b p
(a+b) (p + q)
@
\eject
\begin{thebibliography}{99}
\bibitem{1} nothing
\end{thebibliography}
\end{document}
|
