diff options
Diffstat (limited to 'work/funcdiff_play.scm')
| -rw-r--r-- | work/funcdiff_play.scm | 212 |
1 files changed, 212 insertions, 0 deletions
diff --git a/work/funcdiff_play.scm b/work/funcdiff_play.scm new file mode 100644 index 0000000..429933c --- /dev/null +++ b/work/funcdiff_play.scm @@ -0,0 +1,212 @@ + +;;; Jan 29th + +(load "~/thesis/scmutils/src/calculus/load.scm") + +;(define R2 (make-manifold R^n-type 2)) ; doesn't work, so... +(define R2 (rectangular 2)) +;(define U (patch 'origin R2)) ; also undefined so going to AIM-2005-003.pdf + +(define R2 (rectangular 2)) +(define P2 (polar/cylindrical 2)) + + +(define R2-chi-inverse (R2 '->point)) +(define R2-chi (R2 '->coords)) +(define P2-chi (P2 '->coords)) +(define P2-chi-inverse (P2 '->point)) + +;;; Feb 1st + +(print-expression + ((compose R2-chi-inverse P2-chi) + (up 'x0 'y0))) +;= (up (sqrt (+ (expt x0 2) (expt y0 2))) (atan y0 x0)) + +(print-expression + ((compose R2-chi P2-chi-inverse) + (up 'r0 'theta0))) + +; (up (* r0 (cos theta0)) (* r0 (sin theta0))) + +(define R2->R (-> (UP Real Real) Real)) + +(define f + (compose (literal-function 'f-rect R2->R) + R2-chi)) +; there is a simpler syntax... + +(define R2-point (R2-chi-inverse (up 'x0 'y0))) +(define P2-point (P2-chi-inverse (up 'r0 'theta0))) + +(print-expression (f R2-point)) +; (f-rect (up x0 y0)) +(print-expression (f P2-point)) +; (f-rect (up (* r0 (cos theta0)) (* r0 (sin theta0)))) + +(define g (literal-manifold-function 'g-polar P2)) + +(print-expression (g R2-point)) +; (g-polar (up (sqrt (+ (expt x0 2) (expt y0 2))) (atan y0 x0))) + +(instantiate-coordinates R2 '(x y)) +(instantiate-coordinates P2 '(r theta)) + +(print-expression (x R2-point)) +;x0 + +(print-expression (theta R2-point)) +;(atan y0 x0) + +(define h (+ (* x (square r)) (cube y))) + +(print-expression (h P2-point)) +; (+ (* (expt r0 3) (expt (sin theta0) 3)) (* (expt r0 3) (cos theta0))) + +(print-expression (h R2-point)) +; (+ (expt x0 3) (* x0 (expt y0 2)) (expt y0 3)) + +(print-expression + (D h)) +; a-euclidean-derivative + +; pe is print-expression + +(pe ((D h) R2-point)) +; (down (+ (* 3 (expt x0 2)) (expt y0 2)) (+ (* 2 x0 y0) (* 3 (expt y0 2)))) + +;(pe (D (h R2-point))) +; ERROR + +(pe ((D (compose h R2-chi-inverse)) R2-point)) +; (down (+ (* 3 (expt x0 2)) (expt y0 2)) (+ (* 2 x0 y0) (* 3 (expt y0 2)))) + +; TODO: formal definition of operator vs. function? + +;;; Feb 2 + +(define (vector-field-procedure components coordinate-system) + (define (the-procedure f) + (compose (* (D (compose f (coordinate-system '->point))) + components) + (coordinate-system '->coords))) + the-procedure) + +(define (components->vecotr-field components coordinate-system) + (procedure->vector-field + (vector-field-procedure components coordinate-system))) + +(define v + (components->vector-field + (up (literal-function 'vx (-> (UP Real Real) Real)) + (literal-function 'vy (-> (UP Real Real) Real))) + R2)) + +; shorthand: +; (define w (literal-vector-field 'v R2)) + +;(pe ((v (literal-manifold-function 'f R2)) R2-point)) + +(define (coordinatize vector-field coordsys) + (define ((v f) x) + (let ((b (compose (vector-field (coordsys '->coords)) + (coordsys '->point)))) + (* ((D f) x) (b x)))) + (make-operator v)) + +#| +(pe (((coordinatize v R2) + (literal-function 'f (-> (UP Real Real) Real))) + (up 'x0 'y0))) +|# + +; Note: nice summary of vector field properties on p12 + +;(pe ((d/dx (square r)) R2-point)) +;(pe ((d/dx (square r)) P2-point)) + +(define J (- (* x d/dy) (* y d/dx))) + +(series:for-each pe + (((exp (* 'a J)) R2-chi) + ((R2 '->point) (up 1 0))) + 6) +;(up 1 0) +;(up 0 a) +;(up (* -1/2 (expt a 2)) 0) +;(up 0 (* -1/6 (expt a 3))) +;(up (* 1/24 (expt a 4)) 0) +;(up 0 (* 1/120 (expt a 5))) + +; now do evolution on coordinates + +(define ((((evolution order) + delta-t vector-field) + manifold-function) + manifold-point) + (series:sum + (((exp (* delta-t vector-field)) + manifold-function) + manifold-point) + order)) + +(pe ((((evolution 6) 'a J) R2-chi) + ((R2 '->point) (up 1 0)))) + +(pe ((((evolution 6) 2. J) R2-chi) + ((R2 '->point) (up 1 0)))) + +#| + +(define mywindow (frame -4 4 -4 4)) + +(define (plot-evolution win order initial a step) + (letrec + ((dostep + (lambda (val) + (cond + ((< val a) (let ((this-point ((((evolution order) val J) R2-chi) + ((R2 '->point) initial)))) + (plot-point win + (time this-point) + (coordinate this-point)) + (dostep (+ val step)))))))) + (dostep 0.))) + +(plot-evolution mywindow 6 (up 1. 0.) 6. .01) + + + +(define (explore-evolution window order length) + (define (iterate-mything i x y) + (if (< i length) + (let ((this-point ((((evolution order) i J) R2-chi) + ((R2 '->point) (up x y))))) + (plot-point window (time this-point) (coordinate this-point)) + (iterate-mything (+ .01 i) x y)) + (button-loop x y))) + (define (button-loop ox oy) + (pointer-coordinates + window + (lambda (x y button) + (let ((temp button)) + (cond ((eq? temp 0) (write-line (cons x (cons y (quote ())))) + (display " started.") + (iterate-mything 0 x y)) + ((eq? temp 1) (write-line (cons ox (cons oy (quote ())))) + (display " continued.") + (iterate-mything 0 ox oy)) + ((eq? temp 2) (write-line (cons x (cons y (quote ())))) + (display " hit.") + (button-loop ox oy))))))) + (newline) + (display "Left button starts a trajectory.") + (newline) + (display "Middle button continues a trajectory.") + (newline) + (display "Right button interrogates coordinates.") + (button-loop 0. 0.)) + +(explore-evolution mywindow 5 .4) + +|# |