{-# LANGUAGE TypeOperators, FlexibleContexts, TypeFamilies #-}
import Prelude ()
import Data.Manifold.TreeCover
import Data.Random
import Data.Random.Manifold
import Data.Manifold
import Data.Manifold.Web
import Data.Manifold.DifferentialEquation
import Math.LinearMap.Category
import Data.VectorSpace
import Data.Basis (Basis)
import Linear(V2(..), ex, ey)
import Data.Semigroup
import qualified Data.Foldable as Hask
import Control.Lens
:opt no-lint
import Control.Category.Constrained.Prelude
import Control.Arrow.Constrained
From diagrams:
import Diagrams.Prelude (p2, circle, (&), (^&), moveTo, opacity, fromVertices, Point(P))
type X = ℝ
type T = ℝ
type U = ℝ
type Ðx'U = ℝ
type x × y = ℝ²
et = ey :: Basis ℝ²
From dynamic-plot:
import Graphics.Dynamic.Plot.R2
import Data.Colour.Names
import Data.Colour.Manifold
colourscheme :: Shade' ℝ -> Shade (Colour ℝ)
colourscheme (Shade' u du) = interp (Shade u $ dualNorm du :: Shade ℝ)
where Just interp = rangeOnGeodesic darkblue orange
prettyWebPlot :: PointsWeb ℝ² y -> DynamicPlottable
prettyWebPlot w = plot [ diagramPlot . opacity 0.5 $ fromVertices [P r₁, P r₂]
| ((r₁@(V2 x₁ y₁),_),(r₂@(V2 x₂ y₂),_)) <- edg ]
where edg = webEdges w
colourscheme_heat :: Shade' (U, Ðx'U) -> Shade (Colour ℝ)
colourscheme_heat = cm . factoriseShade
where cm :: (Shade' U, Shade' Ðx'U) -> Shade (Colour ℝ)
cm (Shade' u du, _) = interp (Shade u $ dualNorm du :: Shade ℝ)
Just interp = rangeOnGeodesic darkblue orange
μ :: LocalLinear (X × T) (U, Ðx'U) +> (U, Ðx'U)
μ = arr.LinearFunction $
\(LinearMap (V2 (ðx'u, ðx'ðx'u)
(ðt'u, ðt'ðx'u))) -> (ðx'ðx'u - ðt'u, ðx'u)
heatEqn :: DifferentialEqn (X × T) (U, Ðx'U)
-- Shade (X×T, U,Ðx'U) -> Shade' (X×T +> (U,Ðx'U)) -> ?(Shade' (X×T +> (U,Ðx'U)))
heatEqn = constLinearPDE $ μ
--(`linMapAt` V2 0 2) <$>
heatEqn ((V2 0 0,(0,0)):±[(V2 1 0,(0,0)), (V2 0 1,(0,0)), (V2 0 0,(1,0)), (V2 0 0,(0,1))])
(Shade' (LinearMap $ V2 (0,1) (0,0)) . spanNorm $ Tensor <$>
[V2 (1,0) (0,0), V2 (0,1) (0,0), V2 (0,0) (1,0), V2 (0,0) (0,1)] )
Option {getOption = Just (Shade' {_shade'Ctr = (0.0,0.0), _shade'Narrowness = spanNorm [(0.0,1.0)]},Shade' {_shade'Ctr = Left ().<V2 0.0 0.3333333333333333 ^+^ Right ().<V2 0.6666666666666667 0.0, _shade'Narrowness = spanNorm [ex.⊗(1.3107943462579543,0.0) ^+^ ey.⊗(0.0,0.0),ex.⊗(0.0,1.3107943462579543) ^+^ ey.⊗(-0.6866065623255952,0.0),ex.⊗(0.0,0.0) ^+^ ey.⊗(1.116581052478165,0.0),ex.⊗(0.0,0.0) ^+^ ey.⊗(0.0,0.904534033733291)]})}
initState_heat :: X -> (U, Ðx'U)
initState_heat x = ( tanh (s * (1 - x^2))
, - 2 * s*x / cosh (s * (1 - x^2))^2 )
where s = 0.5
-- plotWindow $ continFnPlot <$> [fst . initState_heat, snd . initState_heat]
tf_heat :: Needle X -> Needle T -> PointsWeb (X × T) (Shade' (U, Ðx'U))
tf_heat δx₀ δt₀ = fromWebNodes euclideanMetric
$ [(V2 x 0, initState_heat x|±|[(0.1,0), (0,0.1)]) | x<-[-2, δx₀-2 .. 0] ]
++ [(V2 x t, zeroV|±|[(1, 0), (0,1) ]) | x<-[-2, δx₀-2 .. 0], t<-[δt₀, δt₀*2 .. 1] ]
GraphWindowSpecR2{lBound=-1.333333333333333, rBound=1.333333333333333, bBound=-1.0328396847046188, tBound=1.0328396847046166, xResolution=640, yResolution=480}
startSt_heat = tf_heat 0.13 0.13
forM_ [ iterateFilterDEqn_static (HighlightInconsistencies $ (-1, 0)|±|[(10, 0),(0, 10)]) heatEqn startSt_heat ]
$ \tfs ->
plotWindow
[ prettyWebPlot $ head tfs
, plotLatest [ plot (fmap colourscheme_heat tfi)
& legendName ("i = "++show i)
| (i,tfi) <- zip [0..] tfs ] ]