2

I am trying to reason how to convert an imperative style program into a functional one like Haskell.

The function is:

void calcPerim(point polycake[], int v, int y, double *perim1, double *perim2){
    int next = 0;
    int index = 0;
    point points[2];
    *perim1 = 0.0;
    *perim2 = 0.0;

    for(int i = 0; i < v; i++)
    {
        next = (i + 1) % v;
        if(polycake[i].y < y && polycake[next].y < y)
            (*perim1) += distance(polycake[i], polycake[next]);
        else if(polycake[i].y > y && polycake[next].y > y)
            (*perim2) += distance(polycake[i], polycake[next]);
        else
        {
            points[index] = intersectPoint(polycake[i], polycake[next], y);
            if(polycake[i].y < y)
            {
                (*perim1) += distance(polycake[i], points[index]);
                (*perim2) += distance(polycake[next],points[index]);
            }
            else
            {
                (*perim2) += distance(polycake[i], points[index]);
                (*perim1) += distance(polycake[next],points[index]);
            }
            index++;
        }
    }

    (*perim1) += distance(points[0], points[1]);
    (*perim2) += distance(points[0], points[1]);
}

I am finding it difficult to understand how I can turn this into a functional approach when it is updating two variables at the same time in some cases. Would it make sense when translating this into recursion to pass in a tuple (perim1, perim2)?

6
  • Yes, using recursion to replace a for loop is a right direction to move into. Commented Mar 12, 2018 at 15:11
  • 4
    This function doesn't take two arguments, it returns two results. Any value passed in to *perim1 or *perim2 is obliderated when they are set to 0.0. Commented Mar 12, 2018 at 15:12
  • 1
    By passing and updating a tuple in the recursion. Commented Mar 12, 2018 at 15:25
  • 3
    I'd start from something like zip polycake (drop 1 polycake) and loop over that using some fold. Perhaps also use cycle polycake to make the access circular, and take v to stop the iteration. Commented Mar 12, 2018 at 16:02
  • 1
    What is this supposed to calculate? The perimeter of some polygon? Where did you got this from? Commented Mar 12, 2018 at 18:39

3 Answers 3

9

It might be a good idea to not translate it straight to Haskell but rather first to C++, which already allows to you structure it in a much more functional way.

First thing, as Cirdec commented, this function doesn't really take perim1 as arguments – those are “output arguments” as Fortran people would say, i.e. they're really results. Also, the v parameter seems to be basically just length of the input array. So in C++ you can reduce it to:

std::pair<double, double> calcPerim(std::vector <point> polycake, int y){
   double perim1 = 0, perim2 = 0;
   ...
   return std::make_pair(perim1, perim2);
}

Now, you have this mutating for loop. In a functional language, the general approach would be to replace that with recursion. For this, you need to make all mutable-state variables function parameters. That includes i, index, points and the perim accumulators (so they're back, in a way... but now as input arguments). You don't need next (which is anyways re-computed from scratch in each iteration).

std::pair<double, double> calcPerim_rec
                   ( std::vector<point> polycake, int y
                   , int i, int index, std::array<point,2> points
                   , double perim1Acc, double perim2Acc ){
   ...
}

...to be used by

std::pair<double, double> calcPerim(std::vector<point> polycake, int y){
   return calcPerim_rec(polycake, y, 0, 0, {}, 0, 0);
}

The recursive function looks very similar to your original loop body; you just need to phrase the end condition:

std::pair<double, double> calcPerim_rec
                   ( std::vector<point> polycake, int y
                   , int i, int index, std::array<point,2> points
                   , double perim1Acc, double perim2Acc ){
    if (i < polycake.length()) {
        int next = (i + 1) % polycake.length();
        if(polycake[i].y < y && polycake[next].y < y)
            perim1Acc += distance(polycake[i], polycake[next]);
        else if(polycake[i].y > y && polycake[next].y > y)
            perim2Acc += distance(polycake[i], polycake[next]);
        else
        {
            points[index] = intersectPoint(polycake[i], polycake[next], y);
            if(polycake[i].y < y)
            {
                perim1Acc += distance(polycake[i], points[index]);
                perim2Acc += distance(polycake[next],points[index]);
            }
            else
            {
                perim2Acc += distance(polycake[i], points[index]);
                perim1Acc += distance(polycake[next],points[index]);
            }
            ++index;
        }
        ++i;
        return calcPerim_rec
                 ( polycake, y, i, index, points, perim1Acc, perim2Acc );
    } else {
       perim1Acc += distance(points[0], points[1]);
       perim2Acc += distance(points[0], points[1]);
       return std::make_pair(perim1Acc, perim2Acc);
    }
}

There's still quite a bit of mutability going on, but we've already encapsulated it to happen all on local variables of the recursion function call, instead of variables lying around during the loop execution. And each of these variables is only updated once, followed by the recursive call, so you can just skip the mutation and simply pass a value plus update to the recursive call:

std::pair<double, double> calcPerim_rec
                   ( std::vector<point> polycake, int y
                   , int i, int index, std::array<point,2> points
                   , double perim1Acc, double perim2Acc ){
    if (i < polycake.length()) {
        int next = (i + 1) % polycake.length();
        if(polycake[i].y < y && polycake[next].y < y)
            return calcPerim_rec
             ( polycake, y, i+1, index, points
             , perim1Acc + distance(polycake[i], polycake[next])
             , perim2Acc
             );
        else if(polycake[i].y > y && polycake[next].y > y)
            return calcPerim_rec
             ( polycake, y, i+1, index, points
             , perim1Acc
             , perim2Acc + distance(polycake[i], polycake[next])
             );
        else
        {
            points[index] = intersectPoint(polycake[i], polycake[next], y);
            if(polycake[i].y < y)
            {
                return calcPerim_rec
                 ( polycake, y, i+1, index+1
                 , points
                 , perim1Acc + distance(polycake[i], points[index])
                 , perim2Acc + distance(polycake[next],points[index])
                 );
            }
            else
            {
                return calcPerim_rec
                 ( polycake, y, i+1, index+1
                 , points
                 , perim1Acc + distance(polycake[i], points[index])
                 , perim2Acc + distance(polycake[next],points[index])
                 );
            }
        }
    } else {
       return std::make_pair( perim1Acc + distance(points[0], points[1])
                            , perim2Acc + distance(points[0], points[1]) );
    }
}

Well, quite a bit of awkward passing-on of parameters, and we still have a mutation of points – but essentially, the code can now be translated to Haskell.

import Data.Vector (Vector, (!), length) as V

calcPerim_rec :: Vector Point -> Int -> Int -> Int -> Int -> [Point] -> (Double, Double) -> (Double, Double)
calcPerim_rec polycake y i index points (perim1Acc, perim2Acc)
 = if i < V.length polycake
    then
     let next = (i + 1) `mod` V.length polycake
     in if yCoord (polycake!i) < y && yCoord (polycake!next) < y
         then calcPerim_rec polycake v y (i+1) index points
                (perim1Acc + distance (polycake!i) (polycake!next)
                perim2Acc
         else
          if yCoord (polycake!i) > y && yCoord (polycake!next) > y)
           then calcPerim_rec polycake v y (i+1) index points
                     perim1Acc
                     (perim2Acc + distance (polycake!i) (polycake!next))
           else
            let points' = points ++ [intersectPoint (polycake!i) (polycake!next) y]
            in if yCoord (polycake!i) < y
                then calcPerim_rec polycake v y (i+1) (index+1)
                       points'
                       (perim1Acc + distance (polycake!i) (points!!index))
                       (perim2Acc + distance (polycake!next) (points!!index))
                else calcPerim_rec polycake v y (i+1) (index+1)
                       points'
                       (perim1Acc + distance (polycake!i) points!!index))
                       (perim2Acc + distance (polycake!next) points!!index))
    else ( perim1Acc + distance (points!!0) (points!!1)
         , perim2Acc + distance (points!!0) (points!!1) )

There's a lot here that could be stylistically improved, but it should in essence work.

A good first thing to actually make it idiomatic is to try and get rid of indices. Indices are strongly eschewed in Haskell, and can often be avoided when you properly work with lists instead of arrays.

Sign up to request clarification or add additional context in comments.

1 Comment

thanks for the detailed answer! It's been tough trying to switch my mode of thinking from the mutable world of C based languages to functional thinking, this helps a lot!
3

It's rarely a good idea to first write a C version and then try to translate it to Haskell.

Instead, consider what you're trying to do, rather than how you're trying to do it.

It appears that given a series of point representing a polygon and a horizontal line at height y, you want to divide it into two polygons at line y and return the perimeter of both. The algorithm assumes the polygon is convex on the vertical axis:

A polygon split by a horizontal line

You're doing this by:

  1. Dividing the segments into those entirely over and entirely under y
  2. Segments that cross y are split into two parts, the one above and the one below y, indicated by red dots.
  3. Adding the intersection line between the two split points (cyan) to both polygons.

We can just implement that logic directly, rather than trying to emulate the iterative approach. Here's an example:

type Length = Double
type Point = (Double, Double)
type Segment = (Point, Point)

-- Check whether a segment is over, under or on the line given by y
segmentCompare :: Double -> Segment -> Ordering
segmentCompare y (p,q) =                                                                          
    case () of                                                                                        
        _ | all (`isUnder` y) [p,q] -> LT
        _ | all (`isOver` y)  [p,q] -> GT
        _ -> EQ

-- Partition a list into (lt, eq, gt) based on f
partition3 :: (Segment -> Ordering) -> [Segment] -> ([Segment], [Segment], [Segment])
partition3 f = p' ([], [], [])
  where
    p' (lt, eq, gt) (x:xs) =
        case f x of
            LT -> p' (x:lt, eq, gt) xs
            EQ -> p' (lt, x:eq, gt) xs
            GT -> p' (lt, eq, x:gt) xs
    p' result [] = result

-- Split a crossing segment into an under part and over part, and return middle
divvy :: Double -> Segment -> (Segment, Segment, Point)
divvy y (start, end) =
    if start `isUnder` y
    then ((start, middle), (middle, end), middle)
    else ((middle, end), (start, middle), middle)
  where
    middle = intersectPoint y (start, end)

-- Split a polygon in two, or Nothing if it's not convex enough
splitPolygon :: Double -> [Point] -> Maybe ([Segment], [Segment])
splitPolygon y list = do
    let (under, crossing, over) = partition3 (segmentCompare y) pairs
    case crossing of
        -- No lines cross. Simple.
        [] -> return (under, over)
        -- Two segments cross. Divide them up.
        [(p1,p2),(q1,q2)] ->
            let (u1, o1, mid1) = divvy y (p1,p2)
                (u2, o2, mid2) = divvy y (q1, q2)
                split = (mid1, mid2) :: Segment
            in return (split:u1:u2:under, split:o1:o2:over)
        -- More segments cross. Algorithm doesn't work.
        rest -> fail "Can't split polygons concave at y"
  where
    pairs = zip list (drop 1 $ cycle list) :: [Segment]

-- Your original function that sums the perimeter of both polygons
calcPerim :: Double -> [Point] -> Maybe (Length, Length)
calcPerim y list = do
    (under, over) <- (splitPolygon y list :: Maybe ([Segment], [Segment]))
    return (sumSegments under, sumSegments over)

-- Self explanatory helpers
distance :: Segment -> Length
distance ((ax, ay), (bx, by)) = sqrt $ (bx-ax)^2 + (by-ay)^2

intersectPoint :: Double -> Segment -> Point
intersectPoint y ((px, py), (qx, qy)) =
    let slope     = (qx-px)/(qy-py)
        intercept = qy - slope*qx
        x = (y - intercept)/slope
    in
        if slope /= 0
        then (x,y)
        else (px, y)

sumSegments :: [Segment] -> Length
sumSegments = sum . map distance

isUnder :: Point -> Double -> Bool
isUnder (_, py) y = py < y
isOver (_, py) y = py > y

1 Comment

Assuming your interpretation is correct, can't you get around the two-crossing maximum by sorting the crossings and alternating inside/outside?
3

You can give this a try, it is a direct translation of your C algorithm to Haskell

data Point = Point {x :: Float, y :: Float}


calcPerim :: [Point] -> Int -> Int -> (Float, Float)
calcPerim ls v some_y = 
    let (x:xs)       = take v ls
        r            = zip (x:xs) (xs ++ [x])
        (u, c, o, _) = foldl someFunction (0, 0, [], fromIntegral some_y :: Float) r
        points_0     = last o
        points_1     = o !! ((length o) - 2)
        answer       = (u + (distance points_0 points_1), c + (distance points_0 points_1))
    in  answer


someFunction :: (Float, Float, [Point], Float) -> (Point, Point) -> (Float, Float, [Point], Float)
someFunction (perim_1, perim_2, points, some_y) (i, nxt)
    | y i < some_y && y nxt < some_y    = (perim_1 + (distance i nxt), perim_2, points, some_y)
    | y i > some_y && y nxt > some_y    = (perim_1, perim_2 + (distance i nxt), points, some_y)
    | y i < some_y                      = (perim_1 + (distance i temp_pt), perim_2 + (distance nxt temp_pt), temp_pt:points, some_y)
    | otherwise                         = (perim_1 + (distance nxt temp_pt), perim_2 + (distance i temp_pt), temp_pt:points, some_y)
    where temp_pt = intersection i nxt some_y


distance :: Point -> Point -> Float
distance p q = undefined


intersection :: Point -> Point -> Float -> Point
intersection p q f = undefined

I didn't run it. Not sure if I used the right fold.

2 Comments

I'm not sure, but perhaps v is meant to be (possibly) greater than length ls, in which case the original code cycles. So, maybe we should use take v (cycle ls) instead. (or v+1 ?)
@thatotherguy Ah, right. I was fooled by next = (i+1) % v, but this is not the only index which is being used.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.