My function output is not what I was expecting, can someone have a look and suggest why this might be the case. Thanks

The assignment asks to write a function to find unique isosceles triangles of integer sides and area.

A=sqrt(m*(m-a)(m-b)(m-a)

where

m = (a+a+b)/2 or P/2.

for example P=16 produces a=5, b=6 A=12

    clear;
    for P = 1:100
      [a,b,A] = IsoTrig(P);
      fprintf('P: %d, a: %d, b: %d, A: %d\n', P, a, b, A);
    end;

function [a,b,A] = IsoTrig(P)
  if P<0 || fix(P) ~=P;
    disp('error, input is not a non negative integer')
    return
  end
  n = 0 % variable to sum loop passes
  m=P/2;
  for ii=1:P;
    a=ii/2;
    b=P-(2*a);
    if b>= 2*a || b<=0;
      continue
    end
    for A=sqrt((m*(m-a)*(m-b)*(m-a)));
      if fix(A)== A && A ~= 0 && fix(a)==a && fix(b) == b && n == 0
        a=a; 
        b=b; 
        A=A; 
        n = n+1;
      else fix(A) ~= A || fix(a) ~= a || fix (b) ~= b || A == 0 || n == 1;   
        a = 0;
        b = 0;
        A = 0;
        continue
      end
      [P,a,b,A]
    end
  end
end

• So we are talking about unique Heronian isosceles triangle.

*My code doesn't output what I want, the [P,a,b,A] gives me Heronian isosceles results but the function output itself is different.

*As can be seen my [P,a,b,A] call results is different to the fprintf results. The current program doesn't find any results to fprintf.

*I don't know how to isolate the unique triangles (where 1 heronian isosceles exists for given P) My code if working only reduces the output to 1 not per given P not find unique.

OUTPUT as stands
P: 101, a: 5.050000e+01, b: 0, A: 0
P: 102, a: 51, b: 0, A: 0
P: 103, a: 5.150000e+01, b: 0, A: 0
P: 104, a: 52, b: 0, A: 0
P: 105, a: 5.250000e+01, b: 0, A: 0
P: 106, a: 53, b: 0, A: 0
P: 107, a: 5.350000e+01, b: 0, A: 0
ans =
     108    30    48   432     <- from [P,a,b,A] call
P: 108, a: 54, b: 0, A: 0
P: 109, a: 5.450000e+01, b: 0, A: 0
P: 110, a: 55, b: 0, A: 0
P: 111, a: 5.550000e+01, b: 0, A: 0
ans =
     112    35    42   588
P: 112, a: 56, b: 0, A: 0
P: 113, a: 5.650000e+01, b: 0, A: 0
P: 114, a: 57, b: 0, A: 0
P: 115, a: 5.750000e+01, b: 0, A: 0
P: 116, a: 58, b: 0, A: 0
P: 117, a: 5.850000e+01, b: 0, A: 0
P: 118, a: 59, b: 0, A: 0
P: 119, a: 5.950000e+01, b: 0, A: 0
P: 120, a: 60, b: 0, A: 0
P: 121, a: 6.050000e+01, b: 0, A: 0
P: 122, a: 61, b: 0, A: 0
P: 123, a: 6.150000e+01, b: 0, A: 0
P: 124, a: 62, b: 0, A: 0
P: 125, a: 6.250000e+01, b: 0, A: 0
ans =
     126    35    56   588
P: 126, a: 63, b: 0, A: 0
P: 127, a: 6.350000e+01, b: 0, A: 0
ans =
     128    34    60   480
P: 128, a: 64, b: 0, A: 0
P: 129, a: 6.450000e+01, b: 0, A: 0
P: 130, a: 65, b: 0, A: 0
P: 131, a: 6.550000e+01, b: 0, A: 0
P: 132, a: 66, b: 0, A: 0
P: 133, a: 6.650000e+01, b: 0, A: 0
P: 134, a: 67, b: 0, A: 0
P: 135, a: 6.750000e+01, b: 0, A: 0
P: 136, a: 68, b: 0, A: 0
P: 137, a: 6.850000e+01, b: 0, A: 0
P: 138, a: 69, b: 0, A: 0
P: 139, a: 6.950000e+01, b: 0, A: 0
P: 140, a: 70, b: 0, A: 0
P: 141, a: 7.050000e+01, b: 0, A: 0
P: 142, a: 71, b: 0, A: 0
P: 143, a: 7.150000e+01, b: 0, A: 0
ans =
     144    37    70   420
P: 144, a: 72, b: 0, A: 0
P: 145, a: 7.250000e+01, b: 0, A: 0
P: 146, a: 73, b: 0, A: 0
P: 147, a: 7.350000e+01, b: 0, A: 0
P: 148, a: 74, b: 0, A: 0
P: 149, a: 7.450000e+01, b: 0, A: 0

Answer 1

[John D'Errico]

Without looking at your code, I first did some pencil and paper work.

Heron's formula for an isosceles triangle, with sides of {a,a,b}, would reduce as:

m = (a + a + b)/2 = a + b/2
A = sqrt(m*(m-a)*(m-a)*)(m-b))

but m-a is just b/2. And m-b is a-b/2. Substitute, and do some quick algebra...

A = sqrt(a^2 - b^2/4) * b/2

This may be a bit simpler to work with than the original formula. It does point out a restriction on a, that a>=b/2. Of course, that should be obvious, since an isosceles triangle MUST have that property.

Next, if you condition your search on the triangle perimeter, then we see another piece of information to reduce the search. If P is an odd integer, then so must be b. Likewise, if P is even, then b is also.

Better, is to look at it from the perspective of a. a must be an integer, at least as large as P/4, and no larger than P/2. (Think about it! Check my logic.) So for ANY given integer perimeter P, the possible side lengths are:

a = ceil(P/4):ceil(P/2 - 1);
b = P - 2*a;

It is always a good idea to test things out, just to be confident you have it right. In fact, I'll write a couple of tiny helper functions to compute a and b for any perimeter.

afun = @( P) ceil(P/4):ceil(P/2 - 1);
bfun = @( P) P - 2*afun( P);

First, try an even perimeter.

P = 12;
[afun( P);bfun( P)]
ans =
     3     4     5
     6     4     2

The first row of that array is a, the second row is b.

Now try an odd value for the perimeter.

P = 13;
[afun( P);bfun( P)]
ans =
     4     5     6
     5     3     1

Next, we require that the area is also an integer.

Areafun = @( P) sqrt(afun( P).^2 - bfun( P).^2/4) .* bfun( P)/2;

But can we do just a bit more before we go into actual numbers? We said before that if P is an odd number, then so must be b. What was the area again?

A = sqrt(a^2 - b^2/4) * b/2

If b is odd, then a^2-b^2/4 must ALWAYS be non-integer, and so must be the area. So we can conclude that to have an integer area, an integral sided isosceles triangle must always have an even perimeter.

So a simple code might just loop over the even integers, starting with P=4, searching for integer values of Areafun.

afun = @( P) ceil(P/4):ceil(P/2 - 1);
bfun = @( P) P - 2*afun( P);
Areafun = @( P) sqrt(afun( P).^2 - bfun( P).^2/4) .* bfun( P)/2;
for P = 4:2:100
  A = Areafun( P);
  % identify triangles with an integer, non-zero area
  ind = find((A == round(A)) & (A ~= 0));
  if ~isempty(ind)
    a = afun( P);
    b = bfun( P);
    disp('Solution(s) found')
    P
    [a(ind);b(ind);A(ind)]
  end
end

Could we do even better? Perhaps by thinking more deeply about the problem, perhaps we might.

[Scott]

Thanks for the suggestion, I tried working with Arrays myself but never succeeded. My method was using for loops with exclusions. The calculations work but I couldn't conceive of how to isolate unique Perimeter Heronian triangles where with the array approach it looks very achievable, exclude arrays that have more than one result.

[Scott]

We have stupid restrictions on us as well, so we cannot use round we must use fix() which is fine, isempty() is not a function we are permitted to use either. Just fix() and sqrt()

[John D'Errico]

Yes, homework is as it is. Sigh. :)

You can get around the lack of round by using fix. Or, for example, (rem(A,1)==0) is another way I often use to test for something being an integer.

You can test for an empty array by things like (numel(ind)==0), or use the length function, if that is allowed.

And of course, if you wanted to do so, you could unroll those vector computations into explicit loops. This is always possible. The code just gets much more messy.

The nice thing about the approach I took is it avoids messy loops as much as possible. Once you put things into simple computations on entire vectors or arrays, your code becomes MUCH cleaner. It is also arguably easier to debug.

Finally, my gut tells me that there are even better solutions to be found. Some careful examination of the solutions found and the equations as formulated would probably give a hint as to how one might solve this Diophantine problem in a far more elegant way, but that is mathematics for you. Or, you could always do it the easy way and use Google, searching for isosceles Heronian triangles. In fact, this is an interesting read on the subject. Or you can look here.

In those links you can find a simple formula that generates all Heronian isosceles triangles. More fun is to derive the formula. Of course, fun is a subjective word. :)

[John D'Errico]

So, as a followup, how might one derive the formula for all isosceles Heronian triangles?

If the base of the triangle has length b, with the other two sides having length a, how might one derive a formula for those side lengths?

We already know that the triangle base must be even, since I showed in my answer the perimeter must be an even integer. So the base length b must also be an even integer.

So break the isosceles triangle into two right triangles, back to back. The height of the triangle must also be an integer, if the area is integral, since the base of that right triangle is b/2 and b was even. Remember that the area of the Heronian triangle is simply h*(b/2).

So effectively, we can create any Heronian isosceles triangle by fabricating it from two pythagorean right triangles, back to back.

Since there is a well known formula to generate all pythagorean triangles, we are very quickly done... (Remainder of proof left to the student.)