You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('') and can be up to 35 characters long.
143 lines
5.8 KiB
143 lines
5.8 KiB
; Copyright (c) 2022 Thomas Baruchel


;


; Permission is hereby granted, free of charge, to any person obtaining a copy


; of this software and associated documentation files (the "Software"), to deal


; in the Software without restriction, including without limitation the rights


; to use, copy, modify, merge, publish, distribute, sublicense, and/or sell


; copies of the Software, and to permit persons to whom the Software is


; furnished to do so, subject to the following conditions:


;


; The above copyright notice and this permission notice shall be included in


; all copies or substantial portions of the Software.


;


; THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR


; IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,


; FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE


; AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER


; LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,


; OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE


; SOFTWARE.


;


; Installation


; 


; The functions can be used with or without compilation:


; * without compilation:


; load("carleman.lisp")$


; * with compilation (must be compiled only once):


; :lisp (compilefile "carleman.lisp");


; look for the compiled file like "convolution.o" or "carleman.fasl"


; and from now on:


; load("carleman.fasl")$




; Compute the Carleman matrix for series whose coefficients are in v.


; Return a list of lists.


; The vector v contains Maxima objects, but the car '(mlist simp) should


; be removed before calling the function.


(defun carleman (v)


(let ((n (listlength v)))


(loop repeat n


with w = (reverse v)


for u = (cons 1 (makelist (1 n) :initialelement 0))


then (loop for x on w


with s = NIL


do (setf s


(cons


(addn


(loop for a in x


for b in u


collect (mul a b)) NIL) s))


finally (return s))


collect u)))






; Formula (4.17)  but there seems to be some misprint in the PDF and


; the sign "" has been added here


(defun carlemandiagleft (m d)


(apply #'mapcar #'list


(loop for w in m


for i from 1


for y = (list (cdr w)) then (cons (nthcdr i w) y)


for d1 = 1 then (mul d1 d)


for r = (cdar m) then (cdr r) ; exactly the required number of 0's !!!


collect (loop with z = (list 1)


for u in y


for x = z


then (cons


(div


(addn


(loop for e in x


for f in u


collect (mul e f)) NIL)


(sub d1 d2)) x)


for d2 = (div d1 d) then (div d2 d)


finally (setf (cdr z) r)


(return x)))))






(defun carlemandiagmiddle (m d)


(loop for NIL in m


for z = (cdar m) then (cdr z)


for q = NIL then (cons 0 q)


for d1 = 1 then (mul d1 d)


collect (append q (cons d1 z))))






; Formula (4.16) in "Continuous time evolution form iterated maps and


; Carleman linearization" (Gralewicz and Kowalski)


(defun carlemandiagright (m d)


(loop for NIL in m


for d1 = 1 then (mul d1 d)


; transpose matrix m to z and iterate on rows of z


for z = (cdr (apply #'mapcar #'list m)) then (mapcar #'cdr (cdr z))


for q = NIL then (cons 0 q)


collect (loop with x = (list 1)


for y = x then (cdr y)


for u in z


for d2 = (mul d1 d) then (mul d2 d)


do (push (div


(addn


(loop for e in x


for f in u


collect (mul e f)) NIL)


(sub d1 d2)) (cdr y))


finally (return (append q x)))))






;;; MAXIMA interface


;;; ================




; Return the Carleman matrix of a function whose Taylor expansion is


; given as a list of coefficients.


(defun $carleman (v)


(simplifya (cons '($matrix)


(mapcar #'(lambda (x) (simplify (cons '(mlist) x)))


(carleman (cdr v)))) NIL))




; Let M be the Carleman matrix of a function having 0 as a fixed point


; (ie. f(0)=0) and f'(0) not in {0, 1} ; now, V(M) is such


; that M = V^(1) . L . V with L a diagonal matrix of eigenvalues.


;


; Return the diagonalized Carleman matrix of a function whose Taylor


; expansion is given as a list of coefficients.


; The first coefficient MUST be 0 (since f(0)=0 is a fixed point).


; The second coefficient MUST be some positive value different from 1.


;


; The function returns a list of three matrices whose product is the


; Carleman matrix.


(defun $carleman_diag (v)


(let ((m (carleman (cdr v)))


(d (caddr v)))


(simplifya (list '(mlist)


; left part


(simplifya


(cons '($matrix) (mapcar #'(lambda (x) (simplify (cons '(mlist) x)))


(carlemandiagleft m d))) NIL)


; diagonal matrix


(simplifya


(cons '($matrix) (mapcar #'(lambda (x) (simplify (cons '(mlist) x)))


(carlemandiagmiddle m d))) NIL)


; right part


(simplifya


(cons '($matrix) (mapcar #'(lambda (x) (simplify (cons '(mlist) x)))


(carlemandiagright m d))) NIL)) NIL)))




