[VA] SRC #012a - Then and Now: Probability

[Valentin Albillo]

.
Hi, all,

After a 7-month hiatus here's my brand-new multi-part SRC #012 - Then an Now, where I'll convincingly demonstrate that some advanced vintage HP calcs which were great problem-solvers back THEN in the 80's (some 40 years ago !), are NOW still perfectly capable of solving recently-proposed non-trivial problems intended to be tackled using modern 2020-era personal computers, not ancient pocket calcs.

To that effect, in the following weeks I'll be proposing a number of such hard problems for you to try and solve using EXCLUSIVELY VINTAGE HP CALCULATORS (physical or virtual,) coding in either RPN, RPL or HP-71B language AND NOTHING ELSE: NO CAS/XCAS, MATHEMATICA, EXCEL, C/C++, PYTHON, etc. Besides, you must post actual code, not just LENGTHY THEORY SESSIONS/EXPOSITIONS (A. C., I'm looking at you ! ).

Finally, NO CODE PANELS at all, just post your RPN/RPL/71B code as-is or formatted however your prefer. Please consider that I`m taking the trouble to use a lot of time and effort to carefully format and solve these problems for your entertainment and potential benefit so please be fair to me and respect those simple rules: only vintage HP calcs, only RPN/RPL/71B code, no math sessions/expositions, no CODE panels. That's it !

'Nuff said, let's begin with one of the easiest problems from the lot, namely:

Problem 1: Probability

A man starts at the top of an equilateral triangular grid having R rows of points and then takes random steps from point to point. Write a program to compute the probability P that after S such steps he ends up in the bottom row, and run it for the case R = 30 rows and S = 60 steps, see the figure below


Once you've found that probability you'll find it very easy to answer any number of additional questions, e.g. What's the probability he ends up in any edge ? In any corner ? In the first 7 rows ? What's the point which has the highest probability ? The lowest ? As a quick check, adding up the probabilities for all the points should return 1 ... after all, he must end up on some point or another ! .

You should strive for 10-12 correct digits (give or take a few ulps) depending on whether you're using a 10- or 12-digit HP model, and of course the faster the running time, the better. If desired, you can check the correctness of your code by running the simpler R = 5 rows, S = 4 steps case, which should return a probability P = 23/288.

If I see enough interest, in a few days I'll post my own original solution for the HP-71B, which is a short program capable of quickly solving the generic problem for any number of rows and steps. I'll also comment on accurate results and possible optimizations, as well as on some other related probabilities and statistics, and once everything is said and done I'll post the next Problem 2.

Let's see your very own clever solutions AND remember the above rules, please.

V.

[ttw]

I am assuming that the probabilities are equal for all branches from any point but the bottom row. (The probabilities fall into various classes depending on the location of the points.) I also assume that the bottom row is an absorbing barrier; the red lines on the border have probability zero (or are non-existent.)

The point being that it's not a copy of a Quincunx.

[ijabbott]

(10-06-2022 03:18 AM)ttw Wrote:  I am assuming that the probabilities are equal for all branches from any point but the bottom row. (The probabilities fall into various classes depending on the location of the points.) I also assume that the bottom row is an absorbing barrier; the red lines on the border have probability zero (or are non-existent.)

The point being that it's not a copy of a Quincunx.

I do not think there is any significance to the red lines, in which case it would be have been better to only colour the points. I would assume that all edges emanating from a point have an equal probability of being traversed on the next step from that point.

[ijabbott]

I wonder if it would be easier to work it out starting at the top and reaching the bottom, or vice versa?

[Valentin Albillo]

.
Hi, ttw,

(10-06-2022 03:18 AM)ttw Wrote:  I am assuming that the probabilities are equal for all branches from any point [...]

Correct.

Quote:[...] but the bottom row.

Nope. The bottom row is like any other row, nothing special about it. I could've asked the probability for the penultimate row or any arbitrary row instead, even the top "row" (the single starting point at the top) for that matter.

Quote:I also assume that the bottom row is an absorbing barrier; the red lines on the border have probability zero (or are non-existent.)

Wrong assumption, the bottom row doesn't absorb anything and the red color is just cosmetic, to highlight the bottom row, it means nothing.

Thanks for your interest and let's see your RPN/RPL/71B code ! You can do it !

Regards.
V.

[C.Ret]

Hi there,

Just three small questions that I ask everyone to confirm that I have understood the problem and not read it too quickly - I have this bad habit -.

           7        
          / \       
         1 ─ 2      
        / \ / \     
       6 ─ A ─ 3    
      / \ / \ / \   
     . ─ 5 ─ 4 ─ .  
    / \ / \ / \ / \ 
   . ─ . ─ 8 ─ . ─ .

From point A, a step can only lead the man only to points 1 2 3 4 5 and 6. To go to points 7 and 8, you need at least two steps by passing through one intermediate points?


My program finds that on a triangle of R=5 rows, the probability that the man is on the last row after S=4 steps is P ~ 0.07961.
Is this a correct value?

On this same triangle R=5, my program finds a probability P ~ 0.1766 that the man is on the last line after S =7 steps.
Does this seem possible to you?

EDIT:
I upgraded my code to avoid rounding errors as much as possible.
Now I find P = 1656 / 20736 for R=5 and S=4 and P = 6329160 / 35831808 for R=5 and S=7.

[Valentin Albillo]

.
Hi, C.Ret,

(10-08-2022 04:57 PM)C.Ret Wrote:  My program finds that on a triangle of R=5 rows, the probability that the man is on the last row after S=4 steps is P ~ 0.07961. Is this a correct value?

Well, in my OP I said:

"If desired, you can check the correctness of your code by running the simpler R = 5 rows, S = 4 steps case, which should return a probability P = 23/288"

As 23/288 evaluates to 0.079861... then yes, your result is correct, except that you omitted the 8, a simple typo.

Quote:EDIT:
I upgraded my code to avoid rounding errors as much as possible. Now I find P = 1656 / 20736 for R=5 and S=4 [...]

But 1656/20736 immediately simplifies to 23/288, which is the correct result, as stated in my OP. Why would you give your result as a fraction not reduced to its lowest terms ?

Thanks for your interest and regards.
V.

[C.Ret]

(10-09-2022 01:29 AM)Valentin Albillo Wrote:  As 23/288 evaluates to 0.079861... then yes, your result is correct, except that you omitted the 8, a simple typo.

Quote:EDIT:
I upgraded my code to avoid rounding errors as much as possible. Now I find P = 1656 / 20736 for R=5 and S=4 [...]

But 1656/20736 immediately simplifies to 23/288, which is the correct result, as stated in my OP. Why would you give your result as a fraction not reduced to its lowest terms ?

Thanks Valentin.

I was very tired last night, which certainly explains the errors in copying, but also my lack of lucidity. I didn't even see that my fraction was reduced to 23/288.
My general condition explains this but does not excuse me.

Today, I wake up in great shape and I will go over my notes and resume my program. I will post it here and explain how it works. The major problem now is its efficiency because in the current version, on my poor HP-71B, it needs no less than 1'27" for (R,S)=(5, 4) and 2'32" for (R, S)=(5, 7).

And of all the authorized machines I own, this HP-71B is by far the fastest.

I expect it to take over 4 hours to calculate (R,S)=(50,60).

I recently got a MATH module, I hope to find the algorithmic way to take advantage of it. For the moment, the calculation is done using arrays and numerous nested loops. It looks very much like a matrix product.

I'm going to have a good breakfast and finish some work in the garden to prepare it for the arrival of winter while thinking about it.

I will post this evening the fruit of my developments, trying to do it before falling asleep and being unable to reread myself...

Best regards.
C.Ret

[PeterP]

Valentin, while it feels out of my reach, it is nonetheless very intriguing as your puzzles always are. Please apologize my question as a burden on your time: we are looking for the probability that we are in the last row after S steps, correct? Not the probability that we can reach the last row in at most S steps?

For S = R - 1 as in the example that is identical, but for S > R one can reach the last row and then step away from it again.

[PeterP]

I am afraid I am missing something in Valentin’s description and I was hoping someone from the community can help me find my error.

I can only find 144 possible 4-step paths in the triangle with 5 rows. Rather than, as the answer/hint in the OP points out 288 (probability of success = good outcomes / all possible outcomes)

Each dot can reach each neighboring dot in 1 step. Each dot has either 2 (the corners of the pyramid), 4 (the middle part of the edges), or 6 (the inner part of the pyramid) neighbors.

I have now gone through the manual process of writing out all possible paths by naming the dots in the Pyramide from 1 to 15 (ie the starting point is dot 1, the second row is dots 2 and 3, the third row is dots 4, 5, and 6, and so forth )

That manual process yields the same outcome as my code = 144 possible 4-step paths, of which 16 end up in the last row. (How often you hit each dot in the last row seems to be the binomial triangle row equivalent to the number of steps, so 1-4-6-4-1 or a total of 16 times)

Clearly I am missing something fundamentally here. Maybe some kind soul can give me a pointer on what I am missing?

[Valentin Albillo]

.
Hi, PeterP, long time no see !

(10-09-2022 04:45 PM)PeterP Wrote:  Valentin, while it feels out of my reach, [...]

Not at all, you're being far too humble, I've seen you solve challenges ten times as difficult. This one's pretty easy, and the algorithm is the same for the huge 30-row grid as for a much smaller one.

Quote:we are looking for the probability that we are in the last row after S steps, correct? Not the probability that we can reach the last row in at most S steps? For S = R - 1 as in the example that is identical, but for S > R one can reach the last row and then step away from it again.

Correct. You perform S random steps from the initial position, then check whether you're in the last row or not and that's it; the man can step away, step towards, or dance a cha-cha-cha for that matter.

Thanks for your interest, looking forward to your code and results, and best regards.
V.

[Valentin Albillo]

(10-09-2022 04:57 AM)C.Ret Wrote:  

(10-09-2022 01:29 AM)Valentin Albillo Wrote:  As 23/288 evaluates to 0.079861... then yes, your result is correct, except that you omitted the 8, a simple typo. But 1656/20736 immediately simplifies to 23/288, which is the correct result, as stated in my OP. Why would you give your result as a fraction not reduced to its lowest terms ?

Thanks Valentin.

My pleasure.

Quote:Today, I wake up in great shape and I will go over my notes and resume my program. I will post it here and explain how it works.

That's great ! Please do !

Quote:The major problem now is its efficiency because in the current version, on my poor HP-71B, [...] I expect it to take over 4 hours to calculate (R,S)=(50,60).

I estimate that my own, non-optimized original solution would solve the (30,60) case {not (50,60), yet another typo} in less than 30 min. when running on a physical HP-71B. Some pretty obvious optimization would make it run in half the time.

Quote:I will post this evening the fruit of my developments, trying to do it before falling asleep and being unable to reread myself...

Ok, good luck !

Thanks and best regards.
V.

[pier4r]

Nice Problems! Do you come across those on your own, due to your work or tinkering, or do you see those (at least in part) in other places? They are really "tasty"!

One observation:

(10-05-2022 08:38 PM)Valentin Albillo Wrote:  Finally, NO CODE PANELS at all, just post your RPN/RPL/71B code as-is or formatted however your prefer. Please consider that I`m taking the trouble to use a lot of time and effort to carefully format and solve these problems for your entertainment and potential benefit so please be fair to me and respect those simple rules: only vintage HP calcs, only RPN/RPL/71B code, no math sessions/expositions, no CODE panels. That's it !

While I agree that given the effort you took it would be only fair to follow your requests, I was thinking that limiting explorations may take away some fun for some people (you mentioned a couple of users too). Further one thing leads to another if the discussion is lively and thus things can be interesting also with "less pure" discussions. Therefore - asking mostly the mods here - would it be possible if the community (not necessarily Valentin, as he already put a lot of effort in the #1 post) opens an extra thread to put there all the math discussion and the non-HP-calc code ?

In other words having two threads, one "pure" and the other for all the other discussions. This to compromise on limiting the explorations and the sharing.

[rprosperi]

(10-10-2022 09:50 AM)pier4r Wrote:  [snip]
Therefore - asking mostly the mods here - would it be possible if the community (not necessarily Valentin, as he already put a lot of effort in the #1 post) opens an extra thread to put there all the math discussion and the non-HP-calc code ?

In other words having two threads, one "pure" and the other for all the other discussions. This to compromise on limiting the explorations and the sharing.

Sure, go ahead, no harm in encouraging related discussions in a parallel thread, it lets folks participate in the problem without violating Valentin's requested 'rules', which I believe were put in place to allow easy capture into PDF documents, which remain problem-focused in the style he prefers, to preserve for subsequent publication, likely on his excellent site. As you say, it's only reasonable that folks respect the rules given the significant effort Valentin puts into creating, documenting and monitoring these great articles.

[Vincent Weber]

Hi Valentin,

Very interesting problem, many thanks!

I have a solution that is somewhat working (albeit not optimised at all), using a HP vintage calculator that is... The HP 27S, with its powerful formula-based language (using LET and GET functions for intermediate results). Before I post it, does it qualify? The calculator is vintage, but the formula language is not explicitely stated in your RPN/RPL/Basic list, so I'm in doubt...

Best regards,

Vincent

[Valentin Albillo]

.
Hi, Vincent, glad to hear from you,

Vincent Weber Wrote:Very interesting problem, many thanks!

You're welcome, thanks to you for your kind words.

Quote:I have a solution that is somewhat working (albeit not optimised at all), using a HP vintage calculator that is... The HP 27S, with its powerful formula-based language (using LET and GET functions for intermediate results). Before I post it, does it qualify?

Of course it does. I didn't consider formula-based languages (do they have a name ?) because I didn't think they could tackle this kind of challenge but if you have a working solution, no matter how raw, I think it's pretty interesting so please post it and, if possible, explain in detail its workings.

Thanks and best regards.
V.

[Vincent Weber]

Many thanks, Valentin.

Here is my formula. Although I used Plus42, I made sure it does not use any of the Plus42 extensions, so that it can work on a vanilla 17B,27S or 19BII.

P=Σ(T:1:N:1:0×(L(J:1)+L(K:1)+Σ(I:1: S:1:L(Q:IF(J=1:2:IF(K=1 OR K=J:4:6)))+L(Z:INT(RAN#×Q))+IF(Z=0:L(J:J+1):IF(Z=1:L(J:J+1)+L(K:K+1):IF(Z=2:IF(K​=1:L(K:2):L(K:K-1)):IF(Z=3:IF(K=1:L(J:J-1):IF(K=J:L(J:J-1)+L(K:K-1):L(K:K+1))):IF(Z=4:L(J:J-1)+L(K:K-1):L(J:J-1))))))))+IF(J=R:1:0))÷N

This works by taking enough samples (N, 10.000 for instance) of scenarios defined by random numbers, like a Monte-Carlo simulation. I initially considered brute force, e. g. trying every single possible scenario, but if you ignore the edge cases you have something like 6 possible moves ^ S possibilities, which is astronomical if S=60, not feasible...

Now the weakness of my code is the random numbers generation. The 27S language only has RAN# (pseudo-random numbers between 0 and 1), has no SEED function, and I am well aware that using INT(RAND#*6) as a proxy for dice rolling is grossly wrong, biased towards lower numbers... I just don't see what else to do.

With S=4 and R=5, with N big enough (e.g. 1 million, which takes literally minutes, even on high performance Plus42 - I don't dare to imagine how much time it would take on a real machine), I get mediocre results that tend to be somewhat close to 23/288.

With S=30 and R=60 I get almost every time a probability of 0. Strange! It does not go very deep down the rows, numbers in the range of 10-20.

The algorithm is simple: J is the row number (from 1 to R), K is the position within the row (from 1 to J, as row J has J positions). If we are at the edges (K=1 or K=J), special treatment arises: only 4 possible moves (even only 2 at the top of the pyramid), so only 4 integer choices for the random numbers: 0 for down-left, 1 for down-right, 2 for horizontal-left (or right if we are at the left edge). The general case adds 3 more choices: 3 for horizontal-right, 4 for up-left, 5 for up-right. J and K are updated accordingly, until S moves are done. Then if J=R 1 is sumed, otherwise 0. In the end P is the number of successful scenarios, just divide by N to get the probability.

This is still an alpha version, I will try to improve it, especially on the random number generation...

Cheers,

Vincent

[Fernando del Rey]

Thanks Valentin, as always, nice and interesting challenge!

I have written a still very raw little program for the plain vanilla 71B (about 25 lines of code), which is giving me the following result for the case of a 30-row grid and after 60 steps.

9.55109846817E-6

As the number of arithmetic operations (divisions and sums) is quite large, I fear that rounding errors may have crept up enough to significantly affect the result.

Using Emu71/Win in my PC, execution time is about 6 seconds for the 30-rowns/60-steps case. I have not yet tried it in the real 71, but I’m afraid execution time may be too long to be practical. I’ll give it a try next.

I do have an idea how to improve the program for execution speed, but for the 30/60 case at hand the reduction in execution time would be about 25% only. Not a great deal.

Before posting my code or going further, please let me know if the result I am obtaining is not too far away from the correct figure.

Thanks again for your efforts to keep us all entertained!

[Fernando del Rey]

Ups!

I was a bit too quick posting my result!

I found a silly error in my code. After correcting the error, the result I get for the 30/60 case is:

9.51234350207E-6

[Valentin Albillo]

.
Hi, Fernando,

Can't reply at length right now but I'll answer your question, namely:

(10-11-2022 10:22 AM)Fernando del Rey Wrote:  I found a silly error in my code. After correcting the error, the result I get for the 30/60 case is:
9.51234350207E-6

Fully correct to the 12 digits you provide. Now go on and eventually post your code and, if at all possible, some comments on its inner workings.

Thanks for your interest, kind words and above all, your results and forthcoming code !

Best regards.
V.