{-# OPTIONS -Wall #-}

{- | 
Module      :  LPFPCore.Newton2
Copyright   :  (c) Scott N. Walck 2023
License     :  BSD3 (see LICENSE)
Maintainer  :  Scott N. Walck <walck@lvc.edu>
Stability   :  stable

Code from chapter 14 of the book Learn Physics with Functional Programming
-}

module LPFPCore.Newton2 where

velocityCF :: Mass
           -> Velocity          -- initial velocity
           -> [Force]           -- list of forces
           -> Time -> Velocity  -- velocity function

type R = Double

type Mass     = R
type Time     = R
type Position = R
type Velocity = R
type Force    = R

velocityCF :: Time -> Time -> [Time] -> Time -> Time
velocityCF Time
m Time
v0 [Time]
fs
    = let fNet :: Time
fNet = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Time]
fs       -- net force
          a0 :: Time
a0   = Time
fNet forall a. Fractional a => a -> a -> a
/ Time
m     -- Newton's second law
          v :: Time -> Time
v Time
t  = Time
v0 forall a. Num a => a -> a -> a
+ Time
a0 forall a. Num a => a -> a -> a
* Time
t  -- constant acceleration eqn
      in Time -> Time
v

positionCF :: Mass
           -> Position          -- initial position
           -> Velocity          -- initial velocity
           -> [Force]           -- list of forces
           -> Time -> Position  -- position function
positionCF :: Time -> Time -> Time -> [Time] -> Time -> Time
positionCF Time
m Time
x0 Time
v0 [Time]
fs
    = let fNet :: Time
fNet = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Time]
fs
          a0 :: Time
a0   = Time
fNet forall a. Fractional a => a -> a -> a
/ Time
m
          x :: Time -> Time
x Time
t  = Time
x0 forall a. Num a => a -> a -> a
+ Time
v0 forall a. Num a => a -> a -> a
* Time
t forall a. Num a => a -> a -> a
+ Time
a0forall a. Num a => a -> a -> a
*Time
tforall a. Floating a => a -> a -> a
**Time
2 forall a. Fractional a => a -> a -> a
/ Time
2
      in Time -> Time
x

velocityFt :: R                 -- dt for integral
           -> Mass
           -> Velocity          -- initial velocity
           -> [Time -> Force]   -- list of force functions
           -> Time -> Velocity  -- velocity function
velocityFt :: Time -> Time -> Time -> [Time -> Time] -> Time -> Time
velocityFt Time
dt Time
m Time
v0 [Time -> Time]
fs
    = let fNet :: Time -> Time
fNet Time
t = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Time -> Time
f Time
t | Time -> Time
f <- [Time -> Time]
fs]
          a :: Time -> Time
a Time
t = Time -> Time
fNet Time
t forall a. Fractional a => a -> a -> a
/ Time
m
      in Time -> Time -> (Time -> Time) -> Time -> Time
antiDerivative Time
dt Time
v0 Time -> Time
a

-- | Given a step size, a y-intercept, and a function, return a function
--   with the given y-intercept whose
--   derivative is the given function.
antiDerivative :: R -> R -> (R -> R) -> (R -> R)
antiDerivative :: Time -> Time -> (Time -> Time) -> Time -> Time
antiDerivative Time
dt Time
v0 Time -> Time
a Time
t = Time
v0 forall a. Num a => a -> a -> a
+ Time -> (Time -> Time) -> Time -> Time -> Time
integral Time
dt Time -> Time
a Time
0 Time
t

-- | Given a step size, a function, a lower limit, and an upper limit, return
--   the definite integral of the function.
integral :: R -> (R -> R) -> R -> R -> R
integral :: Time -> (Time -> Time) -> Time -> Time -> Time
integral Time
dt Time -> Time
f Time
a Time
b
    = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Time -> Time
f Time
t forall a. Num a => a -> a -> a
* Time
dt | Time
t <- [Time
aforall a. Num a => a -> a -> a
+Time
dtforall a. Fractional a => a -> a -> a
/Time
2, Time
aforall a. Num a => a -> a -> a
+Time
3forall a. Num a => a -> a -> a
*Time
dtforall a. Fractional a => a -> a -> a
/Time
2 .. Time
b forall a. Num a => a -> a -> a
- Time
dtforall a. Fractional a => a -> a -> a
/Time
2]]

positionFt :: R                 -- dt for integral
           -> Mass
           -> Position          -- initial position
           -> Velocity          -- initial velocity
           -> [Time -> Force]   -- list of force functions
           -> Time -> Position  -- position function
positionFt :: Time -> Time -> Time -> Time -> [Time -> Time] -> Time -> Time
positionFt Time
dt Time
m Time
x0 Time
v0 [Time -> Time]
fs
    = Time -> Time -> (Time -> Time) -> Time -> Time
antiDerivative Time
dt Time
x0 (Time -> Time -> Time -> [Time -> Time] -> Time -> Time
velocityFt Time
dt Time
m Time
v0 [Time -> Time]
fs)

pedalCoast :: Time -> Force
pedalCoast :: Time -> Time
pedalCoast Time
t
    = let tCycle :: Time
tCycle = Time
20
          nComplete :: Int
          nComplete :: Int
nComplete = forall a b. (RealFrac a, Integral b) => a -> b
truncate (Time
t forall a. Fractional a => a -> a -> a
/ Time
tCycle)
          remainder :: Time
remainder = Time
t forall a. Num a => a -> a -> a
- forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
nComplete forall a. Num a => a -> a -> a
* Time
tCycle
      in if Time
remainder forall a. Ord a => a -> a -> Bool
< Time
10
         then Time
10
         else Time
0

fAir :: R  -- drag coefficient
     -> R  -- air density
     -> R  -- cross-sectional area of object
     -> Velocity
     -> Force
fAir :: Time -> Time -> Time -> Time -> Time
fAir Time
drag Time
rho Time
area Time
v = -Time
drag forall a. Num a => a -> a -> a
* Time
rho forall a. Num a => a -> a -> a
* Time
area forall a. Num a => a -> a -> a
* forall a. Num a => a -> a
abs Time
v forall a. Num a => a -> a -> a
* Time
v forall a. Fractional a => a -> a -> a
/ Time
2

newtonSecondV :: Mass
              -> [Velocity -> Force]  -- list of force functions
              -> Velocity             -- current velocity
              -> R                    -- derivative of velocity
newtonSecondV :: Time -> [Time -> Time] -> Time -> Time
newtonSecondV Time
m [Time -> Time]
fs Time
v0 = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [Time -> Time
f Time
v0 | Time -> Time
f <- [Time -> Time]
fs] forall a. Fractional a => a -> a -> a
/ Time
m

updateVelocity :: R                    -- time interval dt
               -> Mass
               -> [Velocity -> Force]  -- list of force functions
               -> Velocity             -- current velocity
               -> Velocity             -- new velocity
updateVelocity :: Time -> Time -> [Time -> Time] -> Time -> Time
updateVelocity Time
dt Time
m [Time -> Time]
fs Time
v0
    = Time
v0 forall a. Num a => a -> a -> a
+ (Time -> [Time -> Time] -> Time -> Time
newtonSecondV Time
m [Time -> Time]
fs Time
v0) forall a. Num a => a -> a -> a
* Time
dt

velocityFv :: R                    -- time step
           -> Mass
           -> Velocity             -- initial velocity v(0)
           -> [Velocity -> Force]  -- list of force functions
           -> Time -> Velocity     -- velocity function
velocityFv :: Time -> Time -> Time -> [Time -> Time] -> Time -> Time
velocityFv Time
dt Time
m Time
v0 [Time -> Time]
fs Time
t
    = let numSteps :: Int
numSteps = forall a. Num a => a -> a
abs forall a b. (a -> b) -> a -> b
$ forall a b. (RealFrac a, Integral b) => a -> b
round (Time
t forall a. Fractional a => a -> a -> a
/ Time
dt)
      in forall a. (a -> a) -> a -> [a]
iterate (Time -> Time -> [Time -> Time] -> Time -> Time
updateVelocity Time
dt Time
m [Time -> Time]
fs) Time
v0 forall a. [a] -> Int -> a
!! Int
numSteps

bikeVelocity :: Time -> Velocity
bikeVelocity :: Time -> Time
bikeVelocity = Time -> Time -> Time -> [Time -> Time] -> Time -> Time
velocityFv Time
1 Time
70 Time
0 [forall a b. a -> b -> a
const Time
100,Time -> Time -> Time -> Time -> Time
fAir Time
2 Time
1.225 Time
0.6]

newtonSecondTV :: Mass
               -> [(Time,Velocity) -> Force]  -- force funcs
               -> (Time,Velocity)             -- current state
               -> (R,R)                       -- deriv of state
newtonSecondTV :: Time -> [(Time, Time) -> Time] -> (Time, Time) -> (Time, Time)
newtonSecondTV Time
m [(Time, Time) -> Time]
fs (Time
t,Time
v0)
    = let fNet :: Time
fNet = forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [(Time, Time) -> Time
f (Time
t,Time
v0) | (Time, Time) -> Time
f <- [(Time, Time) -> Time]
fs]
          acc :: Time
acc = Time
fNet forall a. Fractional a => a -> a -> a
/ Time
m
      in (Time
1,Time
acc)

updateTV :: R                           -- time interval dt
         -> Mass
         -> [(Time,Velocity) -> Force]  -- list of force funcs
         -> (Time,Velocity)             -- current state
         -> (Time,Velocity)             -- new state
updateTV :: Time
-> Time -> [(Time, Time) -> Time] -> (Time, Time) -> (Time, Time)
updateTV Time
dt Time
m [(Time, Time) -> Time]
fs (Time
t,Time
v0)
    = let (Time
dtdt, Time
dvdt) = Time -> [(Time, Time) -> Time] -> (Time, Time) -> (Time, Time)
newtonSecondTV Time
m [(Time, Time) -> Time]
fs (Time
t,Time
v0)
      in (Time
t  forall a. Num a => a -> a -> a
+ Time
dtdt forall a. Num a => a -> a -> a
* Time
dt
         ,Time
v0 forall a. Num a => a -> a -> a
+ Time
dvdt forall a. Num a => a -> a -> a
* Time
dt)

statesTV :: R                           -- time step
         -> Mass
         -> (Time,Velocity)             -- initial state
         -> [(Time,Velocity) -> Force]  -- list of force funcs
         -> [(Time,Velocity)]           -- infinite list of states
statesTV :: Time
-> Time -> (Time, Time) -> [(Time, Time) -> Time] -> [(Time, Time)]
statesTV Time
dt Time
m (Time, Time)
tv0 [(Time, Time) -> Time]
fs
    = forall a. (a -> a) -> a -> [a]
iterate (Time
-> Time -> [(Time, Time) -> Time] -> (Time, Time) -> (Time, Time)
updateTV Time
dt Time
m [(Time, Time) -> Time]
fs) (Time, Time)
tv0

velocityFtv :: R                           -- time step
            -> Mass
            -> (Time,Velocity)             -- initial state
            -> [(Time,Velocity) -> Force]  -- list of force funcs
            -> Time -> Velocity            -- velocity function
velocityFtv :: Time
-> Time -> (Time, Time) -> [(Time, Time) -> Time] -> Time -> Time
velocityFtv Time
dt Time
m (Time, Time)
tv0 [(Time, Time) -> Time]
fs Time
t
    = let numSteps :: Int
numSteps = forall a. Num a => a -> a
abs forall a b. (a -> b) -> a -> b
$ forall a b. (RealFrac a, Integral b) => a -> b
round (Time
t forall a. Fractional a => a -> a -> a
/ Time
dt)
      in forall a b. (a, b) -> b
snd forall a b. (a -> b) -> a -> b
$ Time
-> Time -> (Time, Time) -> [(Time, Time) -> Time] -> [(Time, Time)]
statesTV Time
dt Time
m (Time, Time)
tv0 [(Time, Time) -> Time]
fs forall a. [a] -> Int -> a
!! Int
numSteps

pedalCoastAir :: [(Time,Velocity)]
pedalCoastAir :: [(Time, Time)]
pedalCoastAir = Time
-> Time -> (Time, Time) -> [(Time, Time) -> Time] -> [(Time, Time)]
statesTV Time
0.1 Time
20 (Time
0,Time
0)
                [\(Time
t,Time
_) -> Time -> Time
pedalCoast Time
t
                ,\(Time
_,Time
v) -> Time -> Time -> Time -> Time -> Time
fAir Time
2 Time
1.225 Time
0.5 Time
v]

pedalCoastAir2 :: Time -> Velocity
pedalCoastAir2 :: Time -> Time
pedalCoastAir2 = Time
-> Time -> (Time, Time) -> [(Time, Time) -> Time] -> Time -> Time
velocityFtv Time
0.1 Time
20 (Time
0,Time
0)
                 [\( Time
t,Time
_v) -> Time -> Time
pedalCoast Time
t
                 ,\(Time
_t, Time
v) -> Time -> Time -> Time -> Time -> Time
fAir Time
1 Time
1.225 Time
0.5 Time
v]

velocityCF' :: Mass
            -> Velocity          -- initial velocity
            -> [Force]           -- list of forces
            -> Time -> Velocity  -- velocity function
velocityCF' :: Time -> Time -> [Time] -> Time -> Time
velocityCF' Time
m Time
v0 [Time]
fs Time
t = forall a. HasCallStack => a
undefined Time
m Time
v0 [Time]
fs Time
t

sumF :: [R -> R] -> R -> R
sumF :: [Time -> Time] -> Time -> Time
sumF = forall a. HasCallStack => a
undefined

positionFv :: R                    -- time step
           -> Mass
           -> Position             -- initial position x(0)
           -> Velocity             -- initial velocity v(0)
           -> [Velocity -> Force]  -- list of force functions
           -> Time -> Position     -- position function
positionFv :: Time -> Time -> Time -> Time -> [Time -> Time] -> Time -> Time
positionFv = forall a. HasCallStack => a
undefined

positionFtv :: R                    -- time step
            -> Mass
            -> Position             -- initial position x(0)
            -> Velocity             -- initial velocity v(0)
            -> [(Time,Velocity) -> Force]  -- force functions
            -> Time -> Position     -- position function
positionFtv :: Time
-> Time -> Time -> Time -> [(Time, Time) -> Time] -> Time -> Time
positionFtv = forall a. HasCallStack => a
undefined

updateExample :: (Time,Velocity)  -- starting state
              -> (Time,Velocity)  -- ending state
updateExample :: (Time, Time) -> (Time, Time)
updateExample = forall a. HasCallStack => a
undefined