HP-35s SOLVE with integration

[Rick314]

I have a program on my HP-35s to solve the problem in the original post. Thank you to those that replied, especially Thomas Klemm -- I recognized your HP-15C program as being compatible with the JRPN 15C simulator, ran it, and enjoyed commenting it (below) to understand it. Notably:

• Solving f(x)*constant = 0 gives the same result as solving f(x) = 0.
• Register arithmetic.
• X intercept of a line through (x1, y1) and (x2, y2).
• Regula falsi, that I have used but did not know its name or history.
  

Code:

# Thomas Klemm, https://www.hpmuseum.org/forum/thread-21638.html
# Registers:
# 0: n     Bessel order
# 1: x
# 2: a     = x1
# 3: b     = x2
# 4: f(a)  = y1 = Jn(x1)*pi
# 5: f(b)  = y2 = Jn(x2)*pi
# x0 = (b*f(a)-a*f(b))/(f(a)-f(b)) = (x2*y1 - x1*y2)/(y1 - y2)
#    = where a line through (x1, y1), (x2, y2) crosses the x axis.
# A null in Jn(x) is also a null in Jn(x)*pi, what is found below.

   001 { 42 21 11 } f LBL A   # main
   002 {    44  3 } STO 3     # b
   003 {       33 } Rv
   004 {    44  2 } STO 2     # a
   005 {    32 12 } GSB B
   006 {    44  4 } STO 4     # f(a)
   007 {    45  3 } RCL 3
   008 {    32 12 } GSB B
   009 {    44  5 } STO 5     # f(b)

   010 { 42 21  0 } f LBL 0   # find x0
   011 {    45  4 } RCL 4     # f(a)
   012 { 45 20  3 } RCL * 3   # b*f(a)
   013 {    45  5 } RCL 5     # f(b)
   014 { 45 20  2 } RCL * 2   # a*f(b)
   015 {       30 } -         # b*f(a)-a*f(b)
   016 {    45  4 } RCL 4     # f(a)
   017 { 45 30  5 } RCL - 5   # f(a)-f(b)
   018 {       10 } /         # (b*f(a)-a*f(b))/(f(a)-f(b)) = x0
   019 {    32 12 } GSB B     # find Jn(x0)*pi
   020 {    43 34 } g RND     # round to display resolution
   021 {    43 20 } g x==0
   022 {    22  2 } GTO 2     # if f(x0) == 0, DONE
   023 {    43 36 } g LSTx    # full-resolution f(x0), nearing 0
   024 { 45 20  4 } RCL * 4   # f(a)*f(x0), x0 going to a or b?
   025 { 43 30  1 } g TEST 1  # x>0?
   026 {    22  1 } GTO 1     # if x0->a, branch ahead,
   027 {    43 36 } g LSTx    # else x0->b...
   028 {    44  5 } STO 5     # f(x0) -> f(b)
   029 {    45  1 } RCL 1     # x0
   030 {    44  3 } STO 3     # x0 -> b
   031 {    22  0 } GTO 0     # branch up for another iteration

   032 { 42 21  1 } f LBL 1
   033 {    43 36 } g LSTx    # f(x0)
   034 {    44  4 } STO 4     # f(x0) -> f(a)
   035 {    45  1 } RCL 1     # x0
   036 {    44  2 } STO 2     # x0 -> a
   037 {    22  0 } GTO 0     # branch up for another iteration

   038 { 42 21  2 } f LBL 2   # DONE
   039 {    45  1 } RCL 1     # x0 final value
   040 {    43 32 } g RTN     # DONE

   041 { 42 21 12 } f LBL B   # integrate f(t) 0 to pi dt
   042 {    44  1 } STO 1     # x
   043 {        0 } 0
   044 {    43 26 } g PI
   045 { 42 20 13 } f integrate C
   046 {    43 32 } g RTN

   047 { 42 21 13 } f LBL C   # f(t), Ry=t, Rx=t
   048 { 45 20  0 } RCL * 0   # n*t
   049 {       34 } X<=>Y     # t
   050 {       23 } SIN       # sin(t)
   051 { 45 20  1 } RCL * 1   # x*sin(t)
   052 {       30 } -         # n*t - sin(t)
   053 {       24 } COS       # cos(n*t - sin(t))
   054 {    43 32 } g RTN

# End.