Saturday, September 30, 2023

HP 15C: Quadratic Regression

HP 15C:  Quadratic Regression




Setup, Normal Equations, and Registers Used



This program fits bivariate date (x,y) to the quadratic polynomial:


y = c + b * x + a * x^2


This program uses the matrix feature of the HP 15C.  




The normal equations that are solved are:


n * c + Σx * b + Σx^2 * a = Σy

Σx *c + Σx^2 * b + Σx^3 *a = Σxy

Σx^2 * c + Σx^3 * b + Σx^4 * a = Σ(x^2 y)


Matrix A = 

[ [ n, Σx, Σx^2 ]

[ Σx, Σx^2, Σx^3 ]

[ Σx^2, Σx^3, Σx^4 ] ]


Matrix B = 

[ [ Σy ]

[ Σxy ]

[ Σ(x^2 y) ] ]


Matrix C = (Matrix B)^-1 Matrix A

[ [ c ] 

[ b ] 

[ a ] ]


Registers Used:


R0:  x data point, row pointer

R1:  y data point, column pointer


Default Statistics Registers:

R2:  n

R3:  Σx

R4:  Σx^2

R5:  Σy

R6:  Σy^2

R7:  Σxy


Additional Statistics Registers:

R8:  Σx^3

R9:  Σx^4

R.0:  Σ(x^2 y)   ("register point zero":  press the decimal point before the 0)



Instructions


1.  Run label A to clear the matrices and registers.  This needs to be done in order to get the most accurate results.  Zero is displayed to indicate when the calculator is ready.  


2.  For each point, enter y data point, press [ ENTER ], x data point, and run label B.  The number of data points (n) will be displayed.  


3.  To calculate the coefficients, run label C.   The coefficients c, b, and a are displayed in order.



Matrix Operations Used


MATRIX 0:  clear all the matrices 


MATRIX 1:  sets the row and counter pointer to 1,1.   


User Mode Set:  automatically advances the pointer to the right row by row.   In programs, turning on and off User Mode is not a step.  However, storage and recall operations in User Mode are marked with a "u" after the step number.  Be careful because if a matrix's pointers (row in R0 and column in R1) return to (1,1), while in USER mode, the next program is skipped.  A special acknowledgement to Torsten for pointing this to me.  



HP 15C Program Code: Quadratic Regression 



Program Memory:  67 steps, 90 bytes

Needs 15 additional memory registers to store the three matrices.


Comments begin with a hash symbol.  #



#  Label A:   Initialization

001 :  42,21,11 :   LBL A

002 :  42,34 :  CLEAR REG

003 :  42,16, 0 :  MATRIX 0

004 :  0 :   0

005 :  43,32  : RTN


# Label B:  Data Entry and Processing

006 : 42,21,12 :  LBL B

007 : 44, 0 :   STO 0

008 :  43,11 : x^2

009 :  34  :  x<>y

010 :  44, 1  : STO 1

011 :  20  :  ×

012 :  44,40,.0 :  STO+ .0  (# store-add to register point zero)

013 :  45, 0  : RCL 0  

014 :  3 :  3

015 :  14 : y^x

016 :  44,40, 8 :  STO+ 8

017 :  45,20, 0 :  RCL× 0

018 :  44,40, 9 : STO+ 9

019 :  45, 1 : RCL 1

020 :  45, 0 : RCL 0

021 :  49  : Σ+

022 :  43,32 : RTN


# Label C:  Calculation

023 : 42,21,13 : LBL C

024 : 3 :  3

025 : 36 : ENTER

026 : 42,23,11 : DIM A

027 : 42,16, 1 : MATRIX A


# Matrix A - Row 1

# Turn on USER Mode ( [ f ] [ RCL ] (USER))

028 : 45, 2 : RCL 2

029 u 44,11 :  STO A

030 : 45, 3 : RCL 3

031 u 44,11 : STO A

032 :  45, 4 :  RCL 4

033 u 44,11 : STO A


# Matrix A - Row 2

034 : 45, 3 : RCL 3

035 u 44,11 : STO A

036 : 45, 4 :  RCL 4

037 u 44,11 : STO A

038 :  45, 8 : RCL 8

039 u 44,11 : STO A


# Matrix A - Row 3

040 : 45, 4 : RCL 4

041 u 44,11 : STO A

042 : 45, 8:  RCL 8

043 u 44,11 : STO A

044 : 45, 9 : RCL 9

045 u 44,11 : STO A

# Turn off USER Mode (unless the next step would be skipped unnecessarily)



# Matrix B

046 : 42,16, 1 : MATRIX 1

047 : 3 :  3

048 : 36 : ENTER

049 : 1 : 1

# Turn on USER Mode

050 : 42,23,12 : DIM B

051 : 45, 5 : RCL 5

052 u 44,12 : STO B

053 : 45, 7 : RCL 7

054 u 44,12 : STO B

055 : 45,.0 : RCL .0   (# recall registers point-zero)

# Turn off USER Mode (unless the next step would be skipped unnecessarily and in this case, an Error 11 would occur)

056 :  44,12 : STO B



# Matrix C - Results

057 : 42,26,13 : RESULT C

058 : 45,16,12 : RCL MATRIX B

059 : 45,16,11 : RCL MATRIX A

060 : 10 : ÷

061 : 42,16, 1 : MATRIX 1

062 u 45,13 : RCL C

063 : 31 : R/S

064 u 45,13 : RCL C

065 : 31 : R/S

066 u 45,13 : RCL C

067 : 43,32 : RTN


Examples


Example 1:

(3, 1.3)

(4, 1.6)

(5, 1.5)

(6, 1.4)


c = -0.54

b = 0.92

a = -0.1


y = -0.54 + 0.92 x - 0.1 x^2


Example 2:

(0, 99.856)

(3, 97.232)

(5, 93.481)

(7, 96.005)

(10, 102.008)


c ≈ 100.3437

b ≈ -2.5318

a ≈ 0.2495


y ≈ 100.3437 - 2.5318 x + 0.2495 x^2



Enjoy!


Eddie



All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Sunday, September 24, 2023

HP 15C: Error Function and Lower Tail Normal Cumulative Function

 HP 15C:   Error Function and Lower Tail Normal Cumulative Function






Formulas Used 



Error Function 


erf(x) = 2 ÷ √π * ∫( e^(-t^2) dt, t = 0 to t = x)




CDF (Cumulative Normal Distribution Function)


CDF(x) = 1 ÷ √(2 * π) * ∫( e^(-t^2) dt, t = -∞ to t = x)


We can use the error function to calculate CDF(x):


0.5 - erf(x), for x < 0


0.5 + erf(x), for x ≥ 0



This program uses the integration function inside the program.  We will need 23 free memory registers to run the program, which will be terminated at the end of program execution.



Labels used:


B:  erf function (error function)

C:  lower tail normal CDF 

4:  function used for integration


(Of course, feel free to use the labels you want, just be mindful to make the appropriate adjustments)




HP 15C Program Code:  Error Function and Lower Tail Normal Cumulative Function


Steps:  35

Bytes:  40



Step # :  Key Code :  Key


001 : 42,21,12:  LBL B

002 :  0  :   0

003 : 34 : x<>y

004 : 42,20, 4 : ∫_x_y 4

005 : 43, 32 : RTN


006 : 42,21,13 : LBL C

007 : 43,30, 2 :  TEST 2 (x<0)

008 : 43, 4, 0 : SF 0

009 : 43, 16 :  ABS

010 : 0   :   0

011 : 34 :  x<>y

012 : 2   :  2

013 : 11 :  √

014 : 10 :  ÷

015 : 42,20, 4 : ∫_y_x 4

016 : 2   :  2

017 : 10 :  ÷

018 : 48 :  .

019 : 5   :  5

020 : 34 :  x<>y

021 : 43, 6, 0 :  F?0

022 : 16  : CHS

023 : 40  :  +

024 : 43, 5, 0 : CF 0

025 : 43, 32 : RTN


026 : 42,21, 4 : LBL 4

027 : 43,11 : x^2

028 : 16 :  CHS

029 : 12 :  e^x

030 :  2  :  2

031 :  20  : ×

032 :  43,26 : π

033 :  11 :  √

034 :  10 :  ÷

035 :  43,32 : RTN




Examples


x = -0.7:  erf (N/A), CDF ≈ 0.2420


x = 0.7:  erf ≈ 0.6778, CDF ≈ 0.7580


x = 1:   erf ≈ 0.8427,  CDF ≈ 0.8413


x = 2.5  erf ≈ 0.9996,  CDF ≈ 0.9938



Eddie



All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Saturday, September 23, 2023

HP 15C: Quadratic Solver Using Solve and Flags

HP 15C:  Quadratic Solver Using Solve and Flags



Introduction


This program solves the quadratic equation for real roots:


a x^2 + b x + c = 0


Instructions:


1.  Key:  a [ENTER], b [ENTER], c [ f ] [ A ]   (enter coefficients into label A)

2.  The first root is displayed. 

3.  Press [ R/S ] for the second root.


Registers Used:


R1 = a

R2 = b

R3 = c

R4 = root 1

R5 = root 2

R6 = abs((R1 + R2 + R3) ÷ 3),  I use [-R6, R6] as my initial guesses.


Flag 0: used to adjust the equation to be solved

Equation 1:  finding the first root:   (a x^2 + b x + c)  (flag 0 is set)

Equation 2:  finding the first root:   (a x^2 + b x + c) ÷ (x - root)  (flag 0 is clear, this is a depressed equation)


When a function with multiple roots is present, we can take out roots already found by dividing f(x) by (x - root).  




HP 15C Program:  Quadratic Solver

Steps:  47, Bytes:  57


Code:


Step : Key Code : Key


001 : 42,21,11 : LBL A

002 : 44, 3 : STO 3

003 : 44, 6 : STO 6

004 : 33 : R↓

005 : 44, 2 : STO 2

006 : 44,40, 6 : STO+ 6

007 : 33 : R↓

008 : 44, 1 : STO 1

009 : 44,40, 6 : STO+ 6

010 : 45, 6 : RCL 6

011 : 3 : 3

012 : 10 : ÷

013 : 43,16 : ABS

014 : 44, 6 : STO 6

015 : 43, 4, 0 : SF 0

016 : 42,21, 0 : LBL 0

017 : 45, 6 : RCL 6

018 : 16 : CHS

019 : 45, 6 : RCL 6

020 : 42,10, 1 : SOLVE 1

021 : 43, 6, 0 : F? 0

022 : 22, 2 : GTO 2

023 : 22, 3 : GTO 3

024 : 42,21, 2 : LBL 2

025 : 44, 4 : STO 4

026 : 43, 5, 0 : CF 0

027 : 22, 0 : GTO 0

028 : 42,21, 1 : LBL 1

029 : 36 : ENTER

030 : 36 : ENTER

031 : 45,20, 1 : RCL× 1

032 : 45,40, 2 : RCL+ 2

033 : 20 : ×

034 : 45,40, 3 : RCL+ 3

035 : 43, 6, 0 : F? 0

036 : 43, 32 : RTN

037 : 34 : x<>y

038 : 45, 4 : RCL 4

039 : 30 : -

040 : 10 : ÷

041 : 43,32 : RTN

042 : 42,21, 3 : LBL 3

043 : 44, 5 : STO 5

044 : 45, 4 : RCL 4

045 : 31 : R/S

046 : 34 : x<>y

047 : 43, 32 : RTN



Examples


x^2 - 5 x - 24 = 0

a = 1

b = -5

c = -24


1 [ ENTER ] 5 [CHS] [ ENTER ] 24 [ CHS ] [ f ] [ A ]

Roots:  8 [ R/S ] -3   (x = 8, x= -3)


2 x^2 + 8 x - 2 = 0

a = 2

b = 8

c = -2


2 [ ENTER ] 8 [ ENTER ] 2 [ CHS ] [ f ] [ A ]

Roots:  0.2361 [ R/S ] -4.2361  (x ≈ 0.2361, x ≈ -4.2361)


Eddie


All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Sunday, September 17, 2023

Filling the Memory of a Casio fx-4000P

Filling the Memory of a Casio fx-4000P





How many programs does it take to fill the 550 step memory?   Here are six programs that pretty much does the job.   I purposely aimed for descriptive prompts and messages. 


Spaces are added for readability.  


Here's are the six programs:


Prg 1:  Approximating the cumulative distribution function of the Normal Curve - to 3 decimal places


Mode +:  COMP,  Number of Steps: 82


"Z≥0" : ?→Z : Fix 3 : 1 - ((1+.196854 Z +.115194 Z² + .000344 Z^3 + .019527 Z^4)^ -4) ÷ 2 : Rnd : Norm : "AREA=" ◢ Ans → A


Source:  Abramowitz and Stegun, Handbook of Mathematical Functions. 1972.


Examples:


Z = 1.6

Results:  AREA = 0.945


Z = 1

Results:  AREA = 0.841


For best results, enter a positive Z.  



Prg 2:  Binomial Distribution PDF with Mean and Variance 


Mode +: COMP, Number of Steps: 86


"P(WIN)" : ?→P : "TRIALS" : ?→T : "WINS" : ?→N : "PDF=" ◢ T nCr N × P x^y N × (1-p) x^y (T-N) ◢ "MU=" ◢ T P  ◢ "VAR=" ◢ Ans (1 - P)


Note:  The combination function, nCr, is shown on the screen as a lone solid C.  I have the nCr for clarification.  


P(WIN):  probability of a successful event

TRIALS:  number of events

WINS:  number of successful events

PDF:  probability of we get the number of successful events

MU:  expected value, mean - depending on P(WIN) and TRIALS

VAR:  variance - depending on P(WIN) and TRIALS


Example:


P(WIN) = 0.7,  TRIALS = 25, WINS = 10

Results:  PDF = 1.324897424 x 10^-3, MU= 17.5, VAR= 5.25   



Prg 3:  Angles of a triangle given 3 side lengths in Degrees - Solve a SSS (side-side-side) Triangle


Mode +: COMP,  Number of Steps: 86


Deg : "A" : ?→A : "B" : ?→B : "C" : ?→C : "<A=" ◢ cos^-1 ((A² + B² - C²) ÷ (-2 B C))  → D ◢ "<B=" ◢ sin^-1(B sin D ÷ A) → E ◢ "<C=" ◢ 180-D-E→F 


Angle <A  (stored in D) is opposite of side with length A

Angle <B (stored in E) is opposite of side with length B

Angle <C (stored in F) is opposite of side with length C


Degrees mode is set in the program.  


Example:

Triangle with lengths A = 24, B = 60, C = 44

Results:  <A = 20.04997572,  <B = 58.99241697, <C = 100.9576073



Prg 4:   Free Fall with Air Resistance (from Ke!san)  

Assume coefficient is standard at k = 0.24 kg/m

(angle is not needed, hyperbolic trig does not depend on angle unit)


Site:  https://keisan.casio.com/exec/system/1231475371

(last retrieved:  February 27, 2023)


Mode +:  COMP,  Number of Steps:  88


.24 → K : 9.80665 → G : "MASS" : ?→M : "DIST" : ?→D : √(M ÷ G ÷ K) → X : "TIME=" ◢ 

X cosh^-1 (e(D K ÷ M)) → T ◢ "VEL=" ◢ X G tanh(T ÷ X) → V


SI units are assumed.


Example:

MASS = 68 kg, DIST (free fall distance) = 1874 m

Results:  TIME = 39.27745305 s, VEL (velocity at free fall) = 52.71191359 m/s



Prg 5:  Sums of 1 to n for k, k^2, k^3, and K^4


Mode +:  COMP, Number of Steps:  83


"1 TO..." : ?→N : "K =" ◢ N (N+1) ÷ 2 → S ◢ "K²=" ◢ S (2 N + 1) ÷ 3 → T ◢ 

"K◢3=" ◢ S² → U ◢ "K◢4=" ◢ T (3 N² + 3 N - 1) ÷ 5 → V 



The power character, x^y can not be used in a string or an error occurs.  The stop character, ◢, can be used.  


K:   Σ (K from K = 1 to K = N)

K²:  Σ (K^2 from K = 1 to K = N)

K◢3:  Σ (K^3 from K = 1 to K = N)

K◢4:  Σ (K^4 from K = 1 to K = N)


Example:  

N = 9

Results:

K:  45

K²:  285

K◢3:  2025

K◢4:  15333



Prg 6:  Simple Ohm's Law Wheel/Volts, Current, Resistance:  "PIE" chart


Mode +:  COMP, Number of Steps: 121


Lbl 0 : "ENT 0 TO SLV" ◢ "I" ◢ ?→I : "V" ◢ ?→ V : "R" ◢ ?→R : I=0 ⇒ Goto 1 : V=0 ⇒ Goto 2:  R=0 ⇒ Goto 3: Goto 0:  Lbl 1:  "I="  ◢ V ÷ R → I ◢ Goto 4: Lbl 2: "V=" ◢ I R → V ◢ Goto 4: Lbl 3: "R=" ◢ V ÷ I → R ◢ Lbl 4: "END"


I:  current (amps, A)

V: voltage (volts, V)

R:  resistance (ohms, Ω)


The inputs will be in this order.  Enter a zero for the variable you want to solve for.  


Examples:


Solve for I:  I = 0, V = 12, R = 3

Result:  I = 4


Solve for V:  I = 20, V = 0, R = 30

Result:  V = 600


Solve for R:  I = 17, V = 120, R = 0

Result:  R ≈ 7.05






Total Number of Programs: 6

Total Steps Used: 121 (I only have 4 left)


Eddie



All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Saturday, September 16, 2023

Integrating Absolute Value Functions

Integrating Absolute Value Functions



Calculating ∫ abs(f(x)) dx


The function abs(f(x)) can be broken into two parts, depending on the sign of f(x):


abs(f(x)) = 

{   f(x)  when f(x) > 0

{  -f(x)  when f(x) < 0



General procedure:


1.  Find the roots of f(x).  

2.  Split the integral at the roots.

3.  For areas where f(x)>0, calculate the sub-area ∫ f(x) dx.

4.  For areas where f(x)<0, calculate the sub-area ∫ -f(x) dx.

5.  Total all the sub-areas.



Let's illustrate this with examples.  Screen shots are made with the HP Prime emulator.  The absolute value function |f(x)| is graphed in blue, while f(x) is graphed in red (for illustrative purposes).  



Example 1:  ∫ |4x- 2| dx from x = 0 to x = 5





∫ |4x- 2| dx from x = 0 to x =5


The root of (4x - 2) is x = 1/2.     

When x < 1/2, (4x - 2) < 0.   

When x > 1/2, (4x - 2) > 0.


Break down the integral into:


∫ (|4x- 2| dx from x = 0 to x =5)

= ∫ ( -(4x- 2) dx from x = 0 to x =1/2) + ∫ (4x- 2 dx from x = 1/2 to x =5)

=  1/2 + 81/2

=  82/2

=  41


We can type in the entire integral into a calculator or app.  Depending on the function and the advanced engine of the calculator, the accuracy may be affected.   Calculators and apps with advanced engines include the HP Prime, Wolfram Alpha, and Desmos.   (Your mileage may vary)



Example 2:  ∫ |x^3 - 28x + 48| dx from x = 0 to x = 3





∫ |x^3 - 28x + 48| dx from x = 0 to x = 3


The roots of x^3 - 28x + 48 are at x = -6, x = 2, and x = 4.   Since the root x = 2 is the only root in the interval [0, 3], this is the one root we are working with.  


With root x = 2, 

When x < 2, x^3 - 28x + 48 > 0

When x > 2, x^3 - 28x + 48 < 0


Then:

∫ ( |x^3 - 28x + 48| dx from x = 0 to x = 3 ) 

= ∫ ( (x^3 - 28x + 48) dx from x = 0 to x = 2 ) 

+ ∫ ( -(x^3 - 28x + 48) dx from x = 2 to x = 3 ) 

= 44 + 5.75

= 49.75



Example 3:  ∫ |e^(2x) - 2| dx from x = 0 to x = 2





The root of e^(2x) - 2 is x = ln 2 ÷ 2 ≈ 0.34657


Let A = ln 2 ÷ 2, and with root x = A,

When x < A, e^(2x) - 2 < 0

When x > A, e^(2x) - 2 > 0


∫ ( |e^(2x) - 2| dx from x = 0 to x = 2 )

=  ∫ ( -(e^(2x) - 2) dx from x = 0 to x = A ) + ∫( (e^(2x) - 2) dx from x = A to x = 2)

≈ 0.193147 + 22.99222

≈ 23.18537


Where is where this method returns approximates and using different calculators and apps may not produce the same results.  


HP Prime:  23.1853693777

Desmos:  23.1853693777

Wolfram Alpha:  23.1853693776920

TI-30X Pro MathPrint:  23.18537052



Hope this technique helps, 


Eddie



All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Sunday, September 10, 2023

TI-84 Plus CE Python and TI-83 Premium CE Python Edition: Blinking Binary Strings

TI-84 Plus CE Python and TI-83 Premium CE Python Edition:  Blinking Binary Strings




Introduction


The script BINBLANK.py draws a circle's whose color changes based on a list of binary numbers:  green for 1, yellow for 0.  


The code listed is made for the TI-84 Plus CE Python and TI-83 Premium CE Python Edition.  



Program Script: BINBLANK.py



# ews 2023-06-24

# blink of lights: binary bits


from random import *

from ti_system import *

from ti_draw import *



# set subroutines

def black():

  from ti_draw import *

  set_color(0,0,0)

def white():

  from ti_draw import *

  set_color(255,255,255)

def yellow():

  from ti_draw import *

  set_color(255,255,0)

def lime():

  from ti_draw import *

  set_color(0,255,0)

def dcircle():

  from ti_draw import *

  fill_circle(150,100,50)


# using pixels, not coordinates

# initialization

blist=[]

print("Binary String Lighting")

print("Choose 1 or 2:")

print("1. generate random string")

print("2. enter a bit list")


# key press

k=0

# flag

f=0

while f==0:

  k=wait_key()

  if k==143:

    # 1 key

    n=int(input("Number of bits? "))

    blist=[randint(0,1) for i in range(n)]

    f=1

  if k==144:

    # 2 key

    blist=eval(input("List of bits: "))

    # data check

    m=len(blist)

    for i in range(m):

      if blist[i]!=0 and blist[i]!=1:

        1/0

        # force an error and terminate

    f=1


# draw the light

clear()

black()

fill_rect(0,0,319,209)

white()

dcircle()


# loop

for i in range(len(blist)):

  sleep(.1)

  if blist[i]==0:yellow()

  if blist[i]==1:lime()

  dcircle()

  black()

  sleep(.55)

  # reset

  white()

  dcircle()


# end indicator

black()

dcircle()

white()

# no line breaks in draw_text 

draw_text(50,50,"THE END: PRESS CLEAR")

show_draw()



Download the file:  https://drive.google.com/file/d/16nxEe51KEhTIX3p1E0uAzHUxCvElQTWE/view?usp=sharing


Notes:  


Three modules are used:  random, ti_system, and ti_draw.  The modules ti_system and ti_draw are Texas Instruments-specific modules.


The ti_system module allows for a get key type of function called wait_key().   wait_key() stops execution and waits for the user to press a key.   The value of the key is returned.   There are separate key combinations with [ 2nd ] or [ alpha ].  


The [ 1 ] key returns a value of 143.

The [ 2 ] key returns a value of 144. 


Here are some other key values:


[ 3 ]   145

[ 4 ]   146

[ 5 ]   147

[ 6 ]   148

[ 7 ]   149

[ 8 ]   150

[ 9 ]   151

[ 0 ]   142

[ enter ]   5

[ ← ]   2

[ → ]   1

[  ↑  ]  3

[  ↓  ]  4

[ clear ]  9



The script also checks self-entered lists to see that each entry is a 0 or 1.   If not, the program "calculates" 1 ÷ 0 to cause a program-stopping error.  


Eddie



All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Saturday, September 9, 2023

Casio fx-9750GIII: Coin Flips and Probability of Winning

Casio fx-9750GIII:  Coin Flips and Probability of Winning




Introduction


The program COINPROB answers two questions in probability:



Take a game, with the chance of winning p.    There are only wins (success) and losses (failures).  Each game is independent and has the same chance of win.


1.  What is the chance of having an amount of wins in a set amount of games?   The binomial distribution is used to answer this question.


2. How many losses must be endured before a win occurs?  For this question, we turn to the geometric distribution.


Three results are listed are:


PDF:  the probability

MEAN:  expected value

VAR:  variance


The probability of success is also listed.



Casio fx-9750GIII Program:  COINPROB


Program Code:

(most spaces are added for readability)


"EWS 2023-06-25"

.5 → P

Locate 1, 4, "P(WIN) = P"

Locate 1, 5, "P(LOSS) = 1-P"

Locate 1, 6, "GAMES = TRIALS" ◢


Lbl 0

ClrText

Menu "PROB WIN VS LOSS", "SETTINGS", S, "# WINS IN GAMES", 1,

"# LOSS BEF. WIN", 2, "EXIT", E


Lbl S

Menu "SETTINGS", "P(WIN) = 0.5", F, "SET P(WIN)", U

Lbl F

.5 → F

"P(WIN) = .5"

"P(LOSS) = .5" ◢

Goto 0


Lbl U

"P(WIN)"? → P

"P(LOSS) = "

1 - P ◢

Goto 0


Lbl 1

"# GAMES"? → N

"# WINS"? → S

N nCr S × P^S × (1 - P)^(N - S) → D

N × P → M

M × (1 - P) → V

Goto R


Lbl 2

"# LOSSES"? → F

(1 - P)^(F - 1) × P → D

P⁻¹ → M

M × (1 - P) ÷ P → V

Goto R


Lbl R

ClrText

Locate 1, 3, "PDF= "

Locate 7, 3, D

Locate 1, 4, "MEAN="

Locate 7, 4, M

Locate 1, 5, "VAR="

Locate 7, 5, V

Locate 1, 6, "P(WIN)="

Locate 9, 6, P ◢

Goto 0


Lbl E

"THANK YOU"


Note:  The bold C is the combination function (nCr).  


A typo has been corrected (see line in red).  I thank Richard Antley for pointing out my error. - 9/27/2023 


Examples



1.  A fair coin is flipped.  You win if the coin flipped is heads.   What is the probability that you win 7 out of 10 times?


P = 0.5


# GAMES?  10

# WINS?  7


Problem Type:  # WINS IN GAMES


Results:

PDF=  0.1171875

MEAN= 5

VAR= 2.5

P(WIN)= 0.5


2.  What is the probability that you flip 7 tails before flipping a head?  Assume a fair coin.


P = 0.5


Problem Type:  # LOSS BEF. WIN


# LOSSES? 7


Results:


PDF=  0.0078125

MEAN= 2

VAR= 2

P(WIN)= 0.5


3.  On any given day in Luau Town, the chance of rain is 35% each day.   On a given week, how many days are expected to be sunny?


WIN:  sunny days (because that is what we want)

P(WIN) = 1 - 35% = 65% = 0.65

We want the mean (expected value).


Problem Type:  # WINS IN GAMES


Change P(WIN) to 0.65 in settings.  


# GAMES?  7  (7 days)

# WINS? (since we are interested in the mean, we can enter any number 0 to 7)


Results:


PDF= (N/A)

MEAN= 4.55

VAR= 1.5925

P(WIN)= 0.65


A week is expected to have about 4.55 sunny days.  Expected value does not have to be an integer.


4.  Assuming the chance of rain is 35% in Luau Town, what is the chance that there are 3 rainy days before a sunny day?


WIN:  sunny days, P(WIN) = 0.65


# LOSSES?  3


Results:  


PDF= 0.079625

MEAN= 1.538461538

VAR= 0.8284023669

P(WIN)= 0.65



Eddie


All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author. 


Sunday, September 3, 2023

Python: Basic Functions with Quaternions

 Python:   Basic Functions with Quaternions   



Calculating with Quaternions



Quaternions are four-dimensional hypercomplex numbers that is written as the form:


r = a * i + b * j + c * k + d


We can also write quaternions as a four dimensional vector [a, b, c, d].


In computer graphics, quaternions are used to draw, display, move, and rotate 3D objects.  The constant term, d, is a mass-point, while the other three terms describe the location of vector point.



Some Mathematical Properties of Quaternions



The conjugate of a quaternion is similar to what a conjugate of a regular complex number would look like:


r* =  - a * i - b * j - c * k + d  


Similarly, the norm of a quaternion works in a familiar fashion:


| r | = abs(r) = √(a^2 + b^2 + c^2 + d^2)


We can multiply a quaternion by any scalar, s, as follows:


s * r = s*  a * i + s * b * j + a* r3 * k + a* r4  


To find the unit vector (unit quaternion), multiply the quaternion by the scalar 1/|r|.


Adding quaternions is communicative and associative.  However, multiplying quaternions is associative, but not communicative.  In short, when it comes to multiplying quaternions, order matters.  


For the complex coefficients:


i^2 = j^2 = k^2 = -1


i * j =  k,    j * i = -k 


j * k = i,   k * j = -i


k * i = j,  i * k = -j




List of Functions Included in QUAT.py



Quaternions are entered as vectors:


r = a * i + b * j + c * k + d -> [a, b, c, d]


qplus(r, q):   addition, r + q

qminus(r, q):  subtraction, r - q

qtimes(r,q):  multiplication,  r * q

qdivide(r,q):   division (multiply by the reciprocal):  r / q  = r * q^-1

qnorm(r):  norm or magnitude

qconj(r):  conjugate

qscale(r, a):   multiply quaternion r by a scalar a

qinv(r):  inverse or reciprocal of a quaternion,  1/r = r* / |r|

qsquare(r):  r^2


qhelp():  brings up the list of functions in QUAT plus.  


To use these functions on other Python scripts:


from QUAT import *


Important QUAT should also import the math module.  




Python:  Quaternions


This should work on all Python platforms since only the math module is used.  This was programmed on a TI-84 Plus CE Python calculator.




Python Code:  QUAT.py


# arithmetic with quaternions

# ews 2023-06-20


from math import *


# help

def qhelp():

  s="Quaternions - qhelp()"

  # add \n for new lines

  s+="\nEnter r=ai+bj+ck+d as [a,b,c,d]."

  s+="\nqplus(r,q) +  qminus(r,q) -"

  s+="\nqtimes(r,q) *  qdivide(r,q) /"

  s+="\nqnorm(r) magnitude \nqconj(r) conjugate \nqscalar(r,a) scalar mult. by a"

  s+="\nqinv(r) r**-1 qsquare(r) r**2"

  # need print to look nice

  return print(s)


# add

def qplus(r,q):

  return [r[i]+q[i] for i in range(4)]


# subtract

def qminus(r,q):

  q=[-i for i in q]

  return qplus(r,q)


# conjugate

def qconj(r):

  return [-r[i] for i in range(3)]+[r[3]]


# norm

def qnorm(r):

  s=sum([i**2 for i in r])

  return sqrt(s)


# inverse/reciprocal

def qinv(r):

  r=qconj(r)

  s=qnorm(r)

  s*=s

  return [i/s for i in r]


# scalar multiplication

def qscalar(r,a):

  return [a*i for i in r]


# multiplication

def qtimes(r,q):

  s=[0,0,0,0]

  s[0]=r[0]*q[3]+r[1]*q[2]-r[2]*q[1]+r[3]*q[0]

  s[1]=-r[0]*q[2]+r[1]*q[3]+r[2]*q[0]+r[3]*q[1]

  s[2]=r[0]*q[1]-r[1]*q[0]+r[2]*q[3]+r[3]*q[2]

  s[3]=r[3]*q[3]-r[0]*q[0]-r[1]*q[1]-r[2]*q[2]

  return s


# division

def qdivide(r,q):

  q=qinv(q)

  return qtimes(r,q)


# qsquare

def qsquare(r):

  return qtimes(r,r)


qhelp()



Download the script here:  https://drive.google.com/file/d/1iNDxw82Gm1nOVJ7jDHdnIIxA4tcaV0p3/view?usp=sharing



Source


Ron Goldman,

Understanding quaternions,

Graphical Models,

Volume 73, Issue 2,

2011,

Pages 21-49,

ISSN 1524-0703,

https://doi.org/10.1016/j.gmod.2010.10.004.

(https://www.sciencedirect.com/science/article/pii/S1524070310000172)


Science Direct requires a paid subscription or access through a university.


Alternate source:  https://www.researchgate.net/publication/220632454_Understanding_quaternions




All original content copyright, © 2011-2023.  Edward Shore.   Unauthorized use and/or unauthorized distribution for commercial purposes without express and written permission from the author is strictly prohibited.  This blog entry may be distributed for noncommercial purposes, provided that full credit is given to the author.