Skip to content

Added Forward Euler to Haskell #215

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 6 commits into from
Jul 8, 2018
Merged

Added Forward Euler to Haskell #215

Show file tree
Hide file tree
Changes from 2 commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
da1a022
Added Forward Euler in Haskell
Thurii Jul 4, 2018
eb1cd8d
Add explicit types and apply hindent
Thurii Jul 4, 2018
537be2c
Generalize solveEuler and checkResults
Thurii Jul 5, 2018
78ed76a
Further generalize checkResult and slim down solveEuler
Thurii Jul 6, 2018
0dd8148
Move timestep over to solveEuler
Thurii Jul 7, 2018
1a39112
Fix Typo
Thurii Jul 8, 2018
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 20 additions & 0 deletions chapters/differential_equations/euler/code/haskell/euler.hs
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,20 @@
solveEuler :: Double -> Int -> [Double]
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For implementing a general acceleration function it could be:

solveEuler :: (Double -> Double) -> Int -> [Double]
solveEuler accFun n = take n $ iterate accFun 1

(note that we don't need the 'timestep' anymore because we can infer it will be implemented in the acceleration function)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm maybe overdoing it, but if you want to be more general, something like
solveEuler :: (Particle -> Acceleration) -> Time -> Particle -> [Particle]
would be good, where Particle represents a particle in space-time, something like like
type Particle = (Position, Speed, Acceleration, Time)
the first parameter is the physical model, the second parameter is the time step, because it shouldn't be part of the physical model, and the third parameter is the initial condition.
I'll admit that this kind of general thinking it not represented in any of the other implementations.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think using Doubles for now is good, but how would you implement the timestep given to the solveEuler function? I think it is unnecessary if we implement it in the acceleration function

Copy link
Member

@jiegillet jiegillet Jul 5, 2018 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Well, in my Verlet implementation (see link in my other comment), everything is either double or list of doubles (to use in more than 1 dimension). Giving them a name only makes it easy to reason about them.
We can't use the tilmestep in the acceleration function, the acceleration function is the f in: dy(t)/dt=f(t,y(t)). You give f any time, place and speed, it tells you what the acceleration is. Once you know the acceleration, you calculate the new position using the tilmestep and whichever method you like, in this case Euler. So these are used in different places.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this is what you meant:

solveEuler :: (Double -> Double) -> Double -> Int -> [Double]
solveEuler accFun timestep n =
  take n $ iterate (\x -> x + accFun x * timestep) 1

as now you can pass the timestep and the acceleration function and it calculates the solution. the function (\x -> x + accFun x * timestep) should calculate the change in x for a given time and acceleration generated by the acceleration function.
Is this correct or I misunderstood what you are trying to say with the acceleration function?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is what I have right now:

solveEuler :: Num a => (a -> a) -> a -> Int -> [a]
solveEuler func initVal n = take n $ iterate func initVal

Just pure forward euler and func can be whatever you need.

solveEuler timestep n = take n $ iterate f 1
where
f x = x - 3 * x * timestep
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am not sure whats more readable,

solveEuler :: Double -> Int -> [Double]
solveEuler timestep n = take n $ iterate f 1
    where
        f x = x - 3 * x * timestep

or

solveEuler :: Double -> Int -> [Double]
solveEuler timestep n = take n $ iterate (\x -> x - 3 * x * timestep) 1

What do you think?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The second one


checkResult :: [Double] -> Double -> Double -> Bool
checkResult results threshold timestep =
and $ zipWith check' results [exp $ -3 * i * timestep | i <- [0 ..]]
where
check' result solution = abs (result - solution) < threshold

main :: IO ()
main =
let timestep = 0.01
n = 1
threshold = 0.01
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

adding a value for accFun:

accFun x = x - 3 * x * timestep

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sould the acceleration function be:

accFun x = -3 * x

in putStrLn $
if checkResult (solveEuler timestep n) threshold timestep
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And then adding the function as a parameter (we do not need timestep anymore):

if checkResult (solveEuler accFun n) threshold timestep

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So now with the timestep:

if checkResult (solveEuler accFun timestep n) threshold timestep

Copy link
Contributor Author

@Thurii Thurii Jul 5, 2018 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you think it's better to leave timestep as a parameter in solveEuler, or have it be a part of accFun?
I'm thinking timestep is more closely related to solveEuler and should be a parameter so that you can try different timesteps without having to modify accFunc.

then "All values within threshold"
else "Value(s) not in threshold"
4 changes: 3 additions & 1 deletion chapters/differential_equations/euler/euler.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -121,6 +121,9 @@ Full code for the visualization follows:
{% sample lang="py" %}
### Python
[import, lang:"python"](code/python/euler.py)
{% sample lang="hs" %}
### Haskell
[import, lang:"haskell"](code/haskell/euler.hs)
{% endmethod %}

<script>
Expand All @@ -145,4 +148,3 @@ $$
\newcommand{\bfomega}{\boldsymbol{\omega}}
\newcommand{\bftau}{\boldsymbol{\tau}}
$$