-
-
Notifications
You must be signed in to change notification settings - Fork 68
/
coordinate.cljc
139 lines (126 loc) · 5.64 KB
/
coordinate.cljc
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
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
;
; Copyright © 2017 Colin Smith.
; This work is based on the Scmutils system of MIT/GNU Scheme:
; Copyright © 2002 Massachusetts Institute of Technology
;
; This is free software; you can redistribute it and/or modify
; it under the terms of the GNU General Public License as published by
; the Free Software Foundation; either version 3 of the License, or (at
; your option) any later version.
;
; This software is distributed in the hope that it will be useful, but
; WITHOUT ANY WARRANTY; without even the implied warranty of
; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
; General Public License for more details.
;
; You should have received a copy of the GNU General Public License
; along with this code; if not, see <http://www.gnu.org/licenses/>.
;
(ns sicmutils.calculus.coordinate
(:require [sicmutils.calculus.basis :as b]
[sicmutils.calculus.manifold :as m]
[sicmutils.calculus.vector-field :as vf]
[sicmutils.calculus.form-field :as ff]
[sicmutils.matrix :as matrix]
[sicmutils.structure :as s]
[sicmutils.util :as u]))
(defn coordinate-functions
[coordinate-system]
(let [prototype (m/coordinate-prototype coordinate-system)]
(s/mapr (fn [access-chain]
(comp (apply s/component access-chain)
#(m/point->coords coordinate-system %)))
(s/structure->access-chains prototype))))
(defn quotify-coordinate-prototype
"Scmutils wants to allow forms like this:
(using-coordinates (up x y) R2-rect ...)
Note that x, y are unquoted. This function converts such an unquoted for
into a quoted one that could be evaluated to return an up-tuple of the symbols:
(up 'x 'y)
Such an object is useful for s/mapr. The function xf is applied before quoting."
[xf p]
(let [q (fn q [p]
(cond (and (sequential? p)
('#{up down} (first p))) `(~(first p) ~@(map q (rest p)))
(vector? p) (mapv q p)
(symbol? p) `'~(xf p)
:else (u/illegal "Invalid coordinate prototype")))]
(q p)))
(defn ^:private symbols-from-prototype
[p]
(cond (and (sequential? p)
('#{up down} (first p))) (mapcat symbols-from-prototype (rest p))
(vector? p) (mapcat symbols-from-prototype p)
(symbol? p) `(~p)
:else (u/illegal (str "Invalid coordinate prototype: " p))))
(defmacro let-coordinates
"Example:
(let-coordinates [[x y] R2-rect
[r theta] R2-polar]
body...)"
[bindings & body]
(when-not (even? (count bindings))
(u/illegal "let-coordinates requires an even number of bindings"))
(let [pairs (partition 2 bindings)
prototypes (map first pairs)
c-systems (mapv second pairs)
coordinate-names (mapcat symbols-from-prototype prototypes)
coordinate-vector-field-names (map vf/coordinate-name->vf-name coordinate-names)
coordinate-form-field-names (map ff/coordinate-name->ff-name coordinate-names)]
`(let [[~@c-systems :as c-systems#]
(mapv m/with-coordinate-prototype
~c-systems
~(mapv #(quotify-coordinate-prototype identity %) prototypes))
c-fns# (map coordinate-functions c-systems#)
c-vfs# (map vf/coordinate-basis-vector-fields c-systems#)
c-ffs# (map ff/coordinate-basis-oneform-fields c-systems#)
~(vec coordinate-names) (flatten c-fns#)
~(vec coordinate-vector-field-names) (flatten c-vfs#)
~(vec coordinate-form-field-names) (flatten c-ffs#)]
~@body)))
(defmacro using-coordinates
"Example:
(using-coordinates (up x y) R2-rect
body...)
Note: this is just a macro wrapping let-coordinates, the use of which is
preferred."
[coordinate-prototype coordinate-system & body]
`(let-coordinates [~coordinate-prototype ~coordinate-system] ~@body))
(defn coordinate-system->vector-basis
[coordinate-system]
(vf/coordinate-basis-vector-fields coordinate-system))
(defn coordinate-system->oneform-basis
[coordinate-system]
(ff/coordinate-basis-oneform-fields coordinate-system))
(defn ^:private c:generate
[n orientation f]
(if (= n 1)
(f 0)
(s/generate n orientation f)))
(defn vector-basis->dual
[vector-basis coordinate-system]
(let [prototype (m/coordinate-prototype coordinate-system)
vector-basis-coefficient-functions (s/mapr #(vf/vector-field->components % coordinate-system) vector-basis)
guts (fn [coords]
(matrix/s:transpose (s/compatible-shape prototype)
(matrix/s:inverse
(s/compatible-shape prototype)
(s/mapr #(% coords) vector-basis-coefficient-functions)
prototype)
prototype))
oneform-basis-coefficient-functions (c:generate (:dimension (m/manifold coordinate-system))
::s/up
#(comp (s/component %) guts))
oneform-basis (s/mapr #(ff/components->oneform-field % coordinate-system) oneform-basis-coefficient-functions)]
oneform-basis))
(defn coordinate-system->basis
"Return the standard basis object for the coordinate system."
[coordinate-system]
(b/make-basis
(vf/coordinate-basis-vector-fields coordinate-system)
(ff/coordinate-basis-oneform-fields coordinate-system)))
(defn Jacobian
"Compute the Jacobian of transition from from-basis to to-basis."
[to-basis from-basis]
(s/mapr (b/basis->oneform-basis to-basis)
(b/basis->vector-basis from-basis)))