Can I improve this sloppy code used with Manipulate?

Question

I'm not an advanced user of Mathematica, and yet I managed to get the Manipulate command work to do something interesting. The problem I have is that I have the feeling that the code I wrote is really sloppy and can be improved/simplified greatly (specially the output could be made much better). Say you have a triangle defined by points p1,p2,p3, and there is a line that's defined by points l1,l2. Given the locations of all this points, the line is used to find intersection points with the edges of the triangle and these new points are used to subdivide the original triangle in sub-triangles and I use area coordinates to do some computation. The code I currently have is as follows:

At[X_] := 1/2 Det[( {
     {1, X[[1, 1]], X[[1, 2]]},
     {1, X[[2, 1]], X[[2, 2]]},
     {1, X[[3, 1]], X[[3, 2]]}
    } )]
N1t[x_, y_, X_] := (1/(
 2 At[X]))(X[[2, 1]] X[[3, 2]] - 
   X[[3, 1]] X[[2, 
     2]] + (X[[2, 2]] - X[[3, 2]]) x + (X[[3, 1]] - X[[2, 1]]) y)
N2t[x_, y_, X_] := (1/(
 2 At[X]))(X[[3, 1]] X[[1, 2]] - 
   X[[1, 1]] X[[3, 
     2]]  + (X[[3, 2]] - X[[1, 2]]) x + (X[[1, 1]] - X[[3, 1]]) y)
N3t[x_, y_, X_] := (1/(
 2 At[X]))(X[[1, 1]] X[[2, 2]] - 
   X[[2, 1]] X[[1, 
     2]] + (X[[1, 2]] - X[[2, 2]]) x + (X[[2, 1]] - X[[1, 1]]) y)
Manipulate[
 {P1, P2, P3, L1, L2} = {p1, p2, p3, l1, l2};
 Module[{m = (l2[[2]] - l1[[2]])/(l2[[1]] - l1[[1]]), 
   mi = (p3[[2]] - p1[[2]])/(p3[[1]] - p1[[1]]), 
   mj = (p3[[2]] - p2[[2]])/(p3[[1]] - p2[[1]]), a, b = -1, c, ci, cj,
    xi, yi, xj, yj, pi, pj},
  a = m;
  c = l1[[2]] - m l1[[1]];
  ci = -mi p1[[1]];
  cj = -mj p2[[1]];
  xi = xi /. Solve[m xi + c == mi xi + ci, xi][[1]];
  yi = m xi + c;
  pi = {xi, yi};
  xj = xj /. Solve[m xj + c == mj xj + cj, xj][[1]];
  yj = m xj + c;
  pj = {xj, yj};
  eq1 = FullSimplify[
    wi1 N2t[x, y, {p3, pi, pj}] == -N3t[x, y, {p1, p2, p3}]];
  eq2 = FullSimplify[
    wi2 N3t[x, y, {p1, pj, pi}] ==  N1t[x, y, {p1, p2, p3}]];
  wi = FullSimplify[Solve[eq1 && eq2, {wi1, wi2}]] /. {x -> pi[[1]], 
     y -> pi[[2]]};
  eq3 = FullSimplify[
    wj1 N3t[x, y, {p3, pi, pj}] == -N3t[x, y, {p1, p2, p3}]];
  eq4 = FullSimplify[
    wj2 N3t[x, y, {p1, p2, pj}] ==  N2t[x, y, {p1, p2, p3}]];
  wj = FullSimplify[Solve[eq3 && eq4, {wj1, wj2}]] /. {x -> pj[[1]], 
     y -> pj[[2]]};
  {Plot[m (x - l1[[1]]) + l1[[2]], {x, 0, 1}, Axes -> None,
    Epilog -> {
      Thick, {Line[{p1, p2, p3, p1}], Line[{p1, p2, pj, p1}], 
       Line[{p3, pi, pj, p3}]}
      },
    PlotRange -> {{-0.1, 1.1}, {0, 1}},
    AspectRatio -> Automatic], Style["pi", Bold], pi, 
   Style["pj", Bold], pj, Style["wi (both)", Bold], wi, 
   Style["wi", Bold], -wi1 /. wi[[1, 1]], 
   Style["wj (both)", Bold], wj, 
   Style["wj", Bold], -wj1 /. wj[[1, 1]]}
  ],
 {{p1, {0, 0}}, Locator},
 {{p2, {1, 0}}, Locator},
 {{p3, {0.5, 1}}, Locator},
 {{l1, {-0.1, 1/2}}, Locator},
 {{l2, {1.1, 1/2}}, Locator}]

I don't like the fact that I'm repeating code to compute intersection points between the line and the sides of the triangle. I also don't like the fact that I can move these points even to places where the computation no longer makes any sense (for example, have the two line points not intersect the triangle). How can I improve this script? Is there a way to constraint the location of those points in Manipulate?

UPDATE

After suggestions, this is my new version of the code.

In[1]:= At[X_] := 1/2 Det[( {
     {1, X[[1, 1]], X[[1, 2]]},
     {1, X[[2, 1]], X[[2, 2]]},
     {1, X[[3, 1]], X[[3, 2]]}
    } )]
N1t[x_, y_, X_] := 
 1/(2 At[X]) (X[[2, 1]] X[[3, 2]] - 
    X[[3, 1]] X[[2, 
      2]] + (X[[2, 2]] - X[[3, 2]]) x + (X[[3, 1]] - X[[2, 1]]) y)
N2t[x_, y_, X_] := 
 1/(2 At[X]) (X[[3, 1]] X[[1, 2]] - 
    X[[1, 1]] X[[3, 
      2]]  + (X[[3, 2]] - X[[1, 2]]) x + (X[[1, 1]] - X[[3, 1]]) y)
N3t[x_, y_, X_] := 
 1/(2 At[X]) (X[[1, 1]] X[[2, 2]] - 
    X[[2, 1]] X[[1, 
      2]] + (X[[1, 2]] - X[[2, 2]]) x + (X[[2, 1]] - X[[1, 1]]) y)

In[5]:= eq1[p1_, p2_, p3_, pi_, pj_] := 
  wi1 N2t[x, y, {p3, pi, pj}] == -N3t[x, y, {p1, p2, p3}];
eq2[p1_, p2_, p3_, pi_, pj_] := 
  wi2 N3t[x, y, {p1, pj, pi}] ==  N1t[x, y, {p1, p2, p3}];
Weight1[p1_, p2_, p3_, pi_, pj_] := 
  Solve[eq1[p1, p2, p3, pi, pj] && eq2[p1, p2, p3, pi, pj], {wi1, 
    wi2}];

In[8]:= eq3[p1_, p2_, p3_, pi_, pj_] := 
  wj1 N3t[x, y, {p3, pi, pj}] == -N3t[x, y, {p1, p2, p3}];
eq4[p1_, p2_, p3_, pi_, pj_] := 
  wj2 N3t[x, y, {p1, p2, pj}] ==  N2t[x, y, {p1, p2, p3}];
Weight2[p1_, p2_, p3_, pi_, pj_] := 
  Solve[eq3[p1, p2, p3, pi, pj] && eq4[p1, p2, p3, pi, pj], {wj1, 
    wj2}];

In[13]:= Manipulate[
 Module[{m = (l2[[2]] - l1[[2]])/(l2[[1]] - l1[[1]]), 
   mi = (p3[[2]] - p1[[2]])/(p3[[1]] - p1[[1]]), 
   mj = (p3[[2]] - p2[[2]])/(p3[[1]] - p2[[1]]), a, b = -1, c, ci, cj,
    xi, yi, xj, yj, pi, pj},
  a = m;
  c = l1[[2]] - m l1[[1]];
  ci = -mi p1[[1]];
  cj = -mj p2[[1]];
  xi = xi /. Solve[m xi + c == mi xi + ci, xi][[1]];
  yi = m xi + c;
  pi = {xi, yi};
  xj = xj /. Solve[m xj + c == mj xj + cj, xj][[1]];
  yj = m xj + c;
  pj = {xj, yj};
  wi = Weight1[p1, p2, p3, pi, pj] /. {x -> pi[[1]], y -> pi[[2]]};
  wj = Weight2[p1, p2, p3, pi, pj] /. {x -> pj[[1]], y -> pj[[2]]};
  {Plot[m (x - l1[[1]]) + l1[[2]], {x, 0, 1}, Axes -> None,
    Epilog -> {
      Thick, {Line[{p1, p2, p3, p1}], Line[{p1, p2, pj, p1}], 
       Line[{p3, pi, pj, p3}]}
      },
    PlotRange -> {{-0.1, 1.1}, {0, 1}},
    AspectRatio -> Automatic], Style["pi", Bold], pi, 
   Style["pj", Bold], pj, Style["wi (both)", Bold], wi, 
   Style["wi", Bold], -wi1 /. wi[[1, 1]], 
   Style["wi check (1-wi)", Bold], 1 + wi1 /. wi[[1, 1]], 
   Style["wj (both)", Bold], wj, 
   Style["wj", Bold], -wj1 /. wj[[1, 1]], 
   Style["wj check (1-wj)", Bold] , 1 + wj1 /. wj[[1, 1]]}
  ],
 {{p1, {0, 0}}, Locator},
 {{p2, {1, 0}}, Locator},
 {{p3, {0.5, 1}}, Locator},
 {{l1, {-0.1, 1/2}}, Locator},
 {{l2, {1.1, 1/2}}, Locator}]

UPDATE 2

After some work, this is my most current version (the code is compact, so it's not optimized for speed)

In[1]:= (* Shape functions *)
At[X_] := 1/2 Det[( {
     {1, X[[1, 1]], X[[1, 2]]},
     {1, X[[2, 1]], X[[2, 2]]},
     {1, X[[3, 1]], X[[3, 2]]}
    } )]
N1t[x_, y_, X_] := 
 1/(2 At[X]) (X[[2, 1]] X[[3, 2]] - 
    X[[3, 1]] X[[2, 
      2]] + (X[[2, 2]] - X[[3, 2]]) x + (X[[3, 1]] - X[[2, 1]]) y)
N2t[x_, y_, X_] := 
 1/(2 At[X]) (X[[3, 1]] X[[1, 2]] - 
    X[[1, 1]] X[[3, 
      2]]  + (X[[3, 2]] - X[[1, 2]]) x + (X[[1, 1]] - X[[3, 1]]) y)
N3t[x_, y_, X_] := 
 1/(2 At[X]) (X[[1, 1]] X[[2, 2]] - 
    X[[2, 1]] X[[1, 
      2]] + (X[[1, 2]] - X[[2, 2]]) x + (X[[2, 1]] - X[[1, 1]]) y)

(* Intersection between two lines *)    
intersect[{p1x_, p1y_}, {p2x_, p2y_}, {q1x_, q1y_}, {q2x_, q2y_}] := 
 With[{d1 = Det[( {
       {p1x, p1y},
       {p2x, p2y}
      } )], d2 = Det[( {
       {q1x, q1y},
       {q2x, q2y}
      } )], d3 = Det[( {
       {p1x - p2x, p1y - p2y},
       {q1x - q2x, q1y - q2y}
      } )]}, {Det[( {
      {d1, p1x - p2x},
      {d2, q1x - q2x}
     } )]/d3, Det[( {
      {d1, p1y - p2y},
      {d2, q1y - q2y}
     } )]/d3}]

Interactive environment

In[3]:= Manipulate[
 DynamicModule[{m = (l2[[2]] - l1[[2]])/(l2[[1]] - l1[[1]]), pi, pj, wi, wj},

  (* Find intersection points *)
  pi = {pix, piy} = intersect[p1, p3, l1, l2];
  pj = {pjx, pjy} = intersect[p2, p3, l1, l2];

  (* Compute weight equations *)

  Wi1eq = Simplify[
     wi1 N2t[pix, piy, {p3, pi, pj}]] == -Simplify[
      N3t[pix, piy, {p1, p2, p3}]];
  Wi2eq = Simplify[wi2 N3t[pix, piy, {p1, pj, pi}]] == 
    Simplify[N1t[pix, piy, {p1, p2, p3}]];
  Wj1eq = Simplify[
     wj1 N3t[pjx, pjy, {p3, pi, pj}]] == -Simplify[
      N3t[pjx, pjy, {p1, p2, p3}]];
  Wj2eq = Simplify[wj2 N3t[pjx, pjy, {p1, p2, pj}]] == 
    Simplify[N2t[pjx, pjy, {p1, p2, p3}]];

  (* Output *)
  Column[{
    Plot[m (x - l1[[1]]) + l1[[2]], {x, 0, 1}, Axes -> None,
     Epilog -> {Thick, PointSize[0.056],
       {Line[{p1, p2, p3, p1}],
        Line[{p1, p2, pj, p1}],
        Line[{p3, pi, pj, p3}],
        {Red, Point[pi]},
        {Red, Point[pj]},
        }
       },
     PlotRange -> {{-1/2, 3/2}, {-0.1, 1.1}},
     AspectRatio -> Automatic],

    (* Node i *)
    Style["\nNode i", Bold, FontSize -> 14],
    Row[{Style["\!\(\*SubscriptBox[\(p\), \(i\)]\) = ", Bold], pi}],
    Row[{Style[
       "\!\(\*SubscriptBox[\(w\), \(i1\)]\),\!\(\*SubscriptBox[\(w\), \(i2\)]\
\) = ", Bold], {Wi1eq[[2]], Wi2eq[[2]]}}],
    Row[{Style["\!\(\*SubscriptBox[\(w\), \(i\)]\) = ", Bold], -Wi1eq[[2]]}],
    Row[{Style[
       "check: 1 - \!\(\*SubscriptBox[\(w\), \(i\)]\) = \
\!\(\*SubscriptBox[\(w\), \(i2\)]\):", Bold], 1 + Wi1eq[[2]]}],

    (* Node j *)
    Style["\nNode j", Bold, FontSize -> 14],
    Row[{Style["\!\(\*SubscriptBox[\(p\), \(j\)]\) = ", Bold], pj}],
    Row[{Style[
       "\!\(\*SubscriptBox[\(w\), \(j1\)]\),\!\(\*SubscriptBox[\(w\), \(j2\)]\
\) = ", Bold], {Wj1eq[[2]], Wj2eq[[2]]}}],
    Row[{Style["\!\(\*SubscriptBox[\(w\), \(j\)]\) = ", Bold], -Wj1eq[[2]]}],
    Row[{Style[
       "check: 1 - \!\(\*SubscriptBox[\(w\), \(j\)]\) = \
\!\(\*SubscriptBox[\(w\), \(j2\)]\): ", Bold], 1 + Wj1eq[[2]]}]
    }]
  ],
 {{p1, {0, 0}}, Locator},
 {{p2, {1, 0}}, Locator},
 {{p3, {0.5, 1}}, Locator},
 {{l1, {-0.1, 1/2}}, Locator},
 {{l2, {1.1, 1/2}}, Locator}]

Answer 1

Analysis

Key step is to solve for wi1, wi2, wj1 and wj2 symbolically in terms of the input points l1, l2, p1, p2 and p3.

The symbolic solution will be used in the Manipulate which will dramatically speed up the process.

Below is the code, some are a direct copy, others have changes.

At[X_] := 1/2 Det[({{1, X[[1, 1]], X[[1, 2]]},
     {1, X[[2, 1]], X[[2, 2]]}, {1, X[[3, 1]], X[[3, 2]]}})]

N1t[x_, y_, X_] := (1/(2 At[X])) (X[[2, 1]] X[[3, 2]] - X[[3, 1]] X[[2, 2]] +
    (X[[2, 2]] - X[[3, 2]]) x + (X[[3, 1]] - X[[2, 1]]) y)

N2t[x_, y_, X_] := (1/(2 At[X])) (X[[3, 1]] X[[1, 2]] - X[[1, 1]] X[[3, 2]] +
    (X[[3, 2]] - X[[1, 2]]) x + (X[[1, 1]] - X[[3, 1]]) y)

N3t[x_, y_, X_] := (1/(2 At[X])) (X[[1, 1]] X[[2, 2]] - X[[2, 1]] X[[1, 2]] +
    (X[[1, 2]] - X[[2, 2]]) x + (X[[2, 1]] - X[[1, 1]]) y)

While deriving the solution we will define the points in terms of {x, y} pairs. Later the actual points in the Manipulate will be used.

l1 = {l1x, l1y};
l2 = {l2x, l2y};
p1 = {p1x, p1y};
p2 = {p2x, p2y};
p3 = {p3x, p3y};

Here are several intermediate results defined in terms of these points.

m = (l2[[2]] - l1[[2]])/(l2[[1]] - l1[[1]]);
mi = (p3[[2]] - p1[[2]])/(p3[[1]] - p1[[1]]);
mj = (p3[[2]] - p2[[2]])/(p3[[1]] - p2[[1]]);

c = l1[[2]] - m l1[[1]];
ci = -mi p1[[1]];
cj = -mj p2[[1]];

Some of the intermediate results used Solve. Instead it is recommended to directly use the solution.

xi = (ci - c)/(m - mi);
yi = m xi + c;
pi = {pix, piy} = {xi, yi};

xj = (cj - c)/(m - mj);
yj = m xj + c;
pj = {pjx, pjy} = {xj, yj};

Solution

Now we will solve the four equations. In the original code replacement is used in eq1 and eq2 for x -> pi[[1]] and y -> pi[[1]]. Here we will do that at the beginning. The net result is that the eq1 through eq4 results can be used directly for wi1 through wj2.

eq1 = Simplify[wi1 N2t[pix, piy, {p3, pi, pj}]] ==
      - Simplify[N3t[pix, piy, {p1, p2, p3}]]

The output is an symbolic result for wi1. Next we will define a function for wi1 that uses the symbolic results.

In the Manipulate the input arguments will be the points l1, l2, p1, p2 and p3. The arguments in the function definition are defined in terms of l1x through p3y.

wi1Fun[{l1x_, l1y_}, {l2x_, l2y_}, {p1x_, p1y_}, {p2x_, p2y_},
      {p3x_, p3y_}] := -(((l1x - l2x) p1y (p1x - p2x) (p1y - p3y) + 
      l2y (p1y p2x p3x + p1x^2 (p1y + p2y - p3y) + 
         l1x (p1y p2x - p1y p3x + p2y p3x - p2x p3y) + 
         p1x (-2 p1y p2x - l1x p2y - p2y p3x + l1x p3y + p2x p3y)) - 
      l1y (p1y p2x p3x + p1x^2 (p1y + p2y - p3y) + 
         l2x (p1y p2x - p1y p3x + p2y p3x - p2x p3y) + 
         p1x (-2 p1y p2x - l2x p2y - p2y p3x + l2x p3y + 
            p2x p3y)))/((l1y (p1x - p3x) + 
        l2y (-p1x + p3x) - (l1x - l2x) (p1y - p3y)) (p1y (p2x - p3x) +
         p2y p3x - p2x p3y + p1x (-p2y + p3y))))

We repeat in a similar fashion for wi2, wj1 and wj2.

eq2 = Simplify[wi2 N3t[pix, piy, {p1, pj, pi}]] == 
  Simplify[N1t[pix, piy, {p1, p2, p3}]]

wi2Fun[{l1x_, l1y_}, {l2x_, l2y_}, {p1x_, p1y_}, {p2x_, p2y_},
  {p3x_, p3y_}] := (-p2y p3x + ((l1y (l2x - p1x) + l2y (-l1x + p1x))
    (p2x - p3x) (p1y - p3y))/(
    l2y (p1x - p3x) + l1y (-p1x + p3x) + (l1x - l2x) (p1y - p3y)) + 
    p2x p3y + ((p2y - p3y) (l1y l2x (p1x - p3x) + 
         l2x p1x (p1y - p3y) + 
         l1x (-l2y p1x - p1x p1y + l2y p3x + p1x p3y)))/(l1y (p1x - 
          p3x) + l2y (-p1x + p3x) - (l1x - l2x) (p1y - 
          p3y)))/(-p2y p3x + p1y (-p2x + p3x) + p1x (p2y - p3y) + 
    p2x p3y)

eq3 = Simplify[wj1 N3t[pjx, pjy, {p3, pi, pj}]] ==
     -Simplify[N3t[pjx, pjy, {p1, p2, p3}]]

wj1Fun[{l1x_, l1y_}, {l2x_, l2y_}, {p1x_, p1y_}, {p2x_, p2y_},
  {p3x_, p3y_}] := -((-p1y p2x + 
      p1x p2y + ((p1x - p2x) (l1y (l2x - p2x) + 
         l2y (-l1x + p2x)) (p2y - p3y))/(
      l2y (p2x - p3x) + 
       l1y (-p2x + p3x) + (l1x - l2x) (p2y - 
          p3y)) + ((p1y - p2y) (l1y l2x (-p2x + p3x) + 
           l2x p2x (-p2y + p3y) + 
           l1x (l2y p2x + p2x p2y - l2y p3x - p2x p3y)))/(l2y (p2x - 
            p3x) + l1y (-p2x + p3x) + (l1x - l2x) (p2y - 
            p3y)))/(-p2y p3x + p1y (-p2x + p3x) + p1x (p2y - p3y) + 
      p2x p3y))

eq4 = Simplify[wj2 N3t[pjx, pjy, {p1, p2, pj}]] == 
      Simplify[N2t[pjx, pjy, {p1, p2, p3}]]

wj2Fun[{l1x_, l1y_}, {l2x_, l2y_}, {p1x_, p1y_}, {p2x_, p2y_},
   {p3x_, p3y_}] := (p1y p3x - ((l1y (l2x - p2x) + l2y (-l1x + p2x)) (p1x - 
       p3x) (p2y - p3y))/(
    l2y (p2x - p3x) + l1y (-p2x + p3x) + (l1x - l2x) (p2y - p3y)) - 
    p1x p3y + ((p1y - p3y) (l1y l2x (-p2x + p3x) + 
         l2x p2x (-p2y + p3y) + 
         l1x (l2y p2x + p2x p2y - l2y p3x - p2x p3y)))/(l1y (p2x - 
          p3x) + l2y (-p2x + p3x) - (l1x - l2x) (p2y - 
          p3y)))/(-p2y p3x + p1y (-p2x + p3x) + p1x (p2y - p3y) + 
    p2x p3y)

Manipulate

In the Manipulate we still have to define several intermediate terms, not because they are required to generate wi1 through wj2, but rather because they are used in the graphics.

It is recommended to use DynamicModule rather than Module inside of Manipulate. I also made a column with rows for the output (completely subjective and optional).

Manipulate[

 DynamicModule[
  {
   m = (l2[[2]] - l1[[2]])/(l2[[1]] - l1[[1]]),
   mi = (p3[[2]] - p1[[2]])/(p3[[1]] - p1[[1]]), 
   mj = (p3[[2]] - p2[[2]])/(p3[[1]] - p2[[1]]),

   c, ci, cj,
   xi, xj,
   yi, yj,
   pi, pj,

   wi1 = wi1Fun[l1, l2, p1, p2, p3],
   wi2 = wi2Fun[l1, l2, p1, p2, p3],
   wj1 = wj1Fun[l1, l2, p1, p2, p3],
   wj2 = wj2Fun[l1, l2, p1, p2, p3],
   wi,
   wj
   },

  c = l1[[2]] - m l1[[1]];
  ci = -mi p1[[1]];
  cj = -mj p2[[1]];
  xi = (ci - c)/(m - mi);
  yi = m xi + c;
  pi = {xi, yi};
  xj = (cj - c)/(m - mj);
  yj = m xj + c;
  pj = {xj, yj};

  wi = {wi1, wi2};
  wj = {wj1, wj2};

  Column[{
    Plot[m (x - l1[[1]]) + l1[[2]], {x, 0, 1},
     Axes -> None,
     Epilog -> {
       Thick,
       {
        Line[{p1, p2, p3, p1}],
        Line[{p1, p2, pj, p1}],
        Line[{p3, pi, pj, p3}]}
       },
     PlotRange -> {{-0.1, 1.1}, {0, 1}},
     AspectRatio -> Automatic],
    "",
    Row[{Style["  pi = ", Bold], pi}],
    Row[{Style["  pj = ", Bold], pj}],
    Row[{Style["  wi = ", Bold], wi}],
    Row[{Style["  wj = ", Bold], wj}]
    }]
  ],

 {{p1, {0, 0}}, Locator},
 {{p2, {1, 0}}, Locator},
 {{p3, {0.5, 1}}, Locator},
 {{l1, {-0.1, 1/2}}, Locator},
 {{l2, {1.1, 1/2}}, Locator}
 ]

The output from Manipulate looks like the figure below before altering any of the locator points.