#BEGIN OUTFILE1 # before write_ruby top matter ALWAYS = 1 INFO = 2 DEBUGL = 3 DEBUGMASSIVE = 4 MAX_UNCHANGED = 10 ATS_MAX_TERMS = 40 # before generate init omniout const $glob_iolevel = INFO # before write_ats library and user def block #BEGIN ATS LIBRARY BLOCK # Begin Function number 2 def omniout_str(iolevel,str) if ($glob_iolevel >= iolevel) then # if number 1 puts str end# end if 1 end # End Function number 2 # Begin Function number 3 def omniout_str_noeol(iolevel,str) if ($glob_iolevel >= iolevel) then # if number 1 printf("%s", str) end# end if 1 end # End Function number 3 # Begin Function number 4 def omniout_labstr(iolevel,label,str) if ($glob_iolevel >= iolevel) then # if number 1 puts label + str end# end if 1 end # End Function number 4 # Begin Function number 5 def omniout_float(iolevel,prelabel,prelen,value,vallen,postlabel) if ($glob_iolevel >= iolevel) then # if number 1 puts prelabel.ljust(30) + value.to_s + postlabel end# end if 1 end # End Function number 5 # Begin Function number 6 def omniout_int(iolevel,prelabel,prelen,value,vallen,postlabel) if ($glob_iolevel >= iolevel) then # if number 1 puts prelabel.ljust(32) + value.to_s + postlabel end# end if 1 end # End Function number 6 # Begin Function number 7 def omniout_float_arr(iolevel,prelabel,elemnt,prelen,value,vallen,postlabel) if ($glob_iolevel >= iolevel) then # if number 1 puts (prelabel.ljust(12) + " " + elemnt.to_s + " " + value.to_s + " " + postlabel) end# end if 1 end # End Function number 7 # Begin Function number 8 def dump_series(iolevel,dump_label,series_name,arr_series,numb) if ($glob_iolevel >= iolevel) then # if number 1 i = 1 while (i <= numb) do puts dump_label + series_name + "[" + i.to_s + "]" + arr_series[i].to_s i += 1 end end end # End Function number 8 # Begin Function number 9 def dump_series_2(iolevel,dump_label,series_name2,arr_series2,numb,subnum,arr_x) if ($glob_iolevel >= iolevel) then # if number 2 sub = 1; while (sub <= subnum) do i = 1; while (i <= numb) do puts dump_label + series_name2 + "[" + sub.to_s + "," + i.to_s + "]" + arr_series2[sub,i].to_s end# end do number 0 sub += 1; end# end do number -1 end# end if 2 end # End Function number 9 # Begin Function number 10 def cs_info(iolevel,str) if ($glob_iolevel >= iolevel) then # if number 2 puts "cs_info " + str + " $glob_correct_start_flag = " , $glob_correct_start_flag.to_s + "$glob_h := " + $glob_h + "$glob_reached_optimal_h := " + $glob_reached_optimal_h end# end if 2 end # End Function number 10 # Begin Function number 11 def logitem_time(fd,secs_in) fd.printf("") if (secs_in >= 0) then # if number 2 years_int = int_trunc(secs_in / $glob_sec_in_year) sec_temp = (secs_in % $glob_sec_in_year) days_int = int_trunc(sec_temp / $glob_sec_in_day) sec_temp = (sec_temp % $glob_sec_in_day) hours_int = int_trunc(sec_temp / $glob_sec_in_hour) sec_temp = (sec_temp % $glob_sec_in_hour) minutes_int = int_trunc(sec_temp / $glob_sec_in_minute) sec_int = (sec_temp % $glob_sec_in_minute) if (years_int > 0) then # if number 3 fd.printf(years_int.to_s + " Years " + days_int.to_s + " Days " + hours_int.to_s + " Hours " + minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds") elsif (days_int > 0) then # if number 4 fd.printf(days_int.to_s + " Days " + hours_int.to_s + " Hours " + minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds") elsif (hours_int > 0) then # if number 5 fd.printf(hours_int.to_s + " Hours " + minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds") elsif (minutes_int > 0) then # if number 6 fd.printf(minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds") else fd.printf(sec_int.to_s + " Seconds") end# end if 6 else fd.printf(" 0.0 Seconds") end# end if 5 fd.printf("") end # End Function number 11 # Begin Function number 12 def omniout_timestr(secs_in) if (secs_in >= 0) then # if number 5 years_int = int_trunc(secs_in / $glob_sec_in_year) sec_temp = (secs_in % $glob_sec_in_year) days_int = int_trunc(sec_temp / $glob_sec_in_day) sec_temp = (sec_temp % $glob_sec_in_day) hours_int = int_trunc(sec_temp / $glob_sec_in_hour) sec_temp = (sec_temp % $glob_sec_in_hour) minutes_int = int_trunc(sec_temp / $glob_sec_in_minute) sec_int = (sec_temp % $glob_sec_in_minute) if (years_int > 0) then # if number 6 puts years_int.to_s + " Years " + days_int.to_s + " Days " + hours_int.to_s + " Hours " + minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds" elsif (days_int > 0) then # if number 7 puts days_int.to_s + " Days " + hours_int.to_s + " Hours " + minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds" elsif (hours_int > 0) then # if number 8 puts hours_int.to_s + " Hours " + minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds" elsif (minutes_int > 0) then # if number 9 puts minutes_int.to_s + " Minutes " + sec_int.to_s + " Seconds" else puts sec_int.to_s + " Seconds" end# end if 9 else puts " 0.0 Seconds" end# end if 8 end # End Function number 12 # Begin Function number 13 def zero_ats_ar(arr_a) iii = 1 while (iii <= ATS_MAX_TERMS) do # do number -1 arr_a [iii] = $glob__0 iii = iii + 1 end# end do number -1 end # End Function number 13 # Begin Function number 14 def ats(mmm_ats,arr_a,arr_b,jjj_ats) ret_ats = $glob__0 if (jjj_ats <= mmm_ats) then # if number 8 ma_ats = mmm_ats + 1 iii_ats = jjj_ats while (iii_ats <= mmm_ats) do # do number -1 lll_ats = ma_ats - iii_ats if ((lll_ats <= ATS_MAX_TERMS and (iii_ats <= ATS_MAX_TERMS) )) then # if number 9 ret_ats = ret_ats + c(arr_a[iii_ats])*c(arr_b[lll_ats]) end# end if 9 iii_ats = iii_ats + 1 end# end do number -1 end# end if 8 return ( ret_ats) end # End Function number 14 # Begin Function number 15 def att(mmm_att,arr_aa,arr_bb,jjj_att) ret_att = $glob__0 if (jjj_att < mmm_att) then # if number 8 ma_att = mmm_att + 2 iii_att = jjj_att while ((iii_att < mmm_att) and (iii_att <= ATS_MAX_TERMS) ) do # do number -1 lll_att = ma_att - iii_att al_att = (lll_att - 1) if ((lll_att <= ATS_MAX_TERMS and (iii_att <= ATS_MAX_TERMS) )) then # if number 9 ret_att = ret_att + c(arr_aa[iii_att])*c(arr_bb[lll_att])* c(al_att) end# end if 9 iii_att = iii_att + 1 end# end do number -1 ret_att = ret_att / c(mmm_att) end# end if 8 return ( ret_att) end # End Function number 15 # Begin Function number 16 def logditto(file) file.printf("") file.printf("ditto") file.printf("") end # End Function number 16 # Begin Function number 17 def logitem_integer(file,n) file.printf("") file.printf("%d",n) file.printf("") end # End Function number 17 # Begin Function number 18 def logitem_str(file,str) file.printf("") file.printf(str.to_s) file.printf("") end # End Function number 18 # Begin Function number 19 def logitem_good_digits(file,rel_error) file.printf("") file.printf("%d",$glob_min_good_digits) file.printf("") end # End Function number 19 # Begin Function number 20 def log_revs(file,revs) file.printf(revs.to_s) end # End Function number 20 # Begin Function number 21 def logitem_float(file,x) file.printf("") file.printf(x.to_s) file.printf("") end # End Function number 21 # Begin Function number 22 def logitem_h_reason(file) file.printf("") if ($glob_h_reason == 1) then # if number 8 file.printf("Max H") elsif ($glob_h_reason == 2) then # if number 9 file.printf("Display Interval") elsif ($glob_h_reason == 3) then # if number 10 file.printf("Optimal") elsif ($glob_h_reason == 4) then # if number 11 file.printf("Pole Accuracy") elsif ($glob_h_reason == 5) then # if number 12 file.printf("Min H (Pole)") elsif ($glob_h_reason == 6) then # if number 13 file.printf("Pole") elsif ($glob_h_reason == 7) then # if number 14 file.printf("Opt Iter") else file.printf("Impossible") end# end if 14 file.printf("") end # End Function number 22 # Begin Function number 23 def logstart(file) file.printf("") end # End Function number 23 # Begin Function number 24 def logend(file) file.printf("") end # End Function number 24 # Begin Function number 25 def chk_data() errflag = false if ($glob_max_iter < 2) then # if number 14 omniout_str(ALWAYS,"Illegal max_iter") errflag = true end# end if 14 if (errflag) then # if number 14 end# end if 14 end # End Function number 25 # Begin Function number 26 def comp_expect_sec(t_end2,t_start2,t2,clock_sec2) ms2 = c(clock_sec2) sub1 = c(t_end2-t_start2) sub2 = c(t2-t_start2) if (sub1 == $glob__0) then # if number 14 sec_left = $glob__0 else if (sub2 > $glob__0) then # if number 15 rrr = (sub1/sub2) sec_left = rrr * c(ms2) - c(ms2) else sec_left = $glob__0 end# end if 15 end# end if 14 return ( sec_left) end # End Function number 26 # Begin Function number 27 def comp_percent(t_end2,t_start2, t2) sub1 = (t_end2-t_start2) sub2 = (t2-t_start2) if (sub2 > $glob_small_float) then # if number 14 rrr = ($glob__100*sub2)/sub1 else rrr = 0.0 end# end if 14 return ( rrr) end # End Function number 27 # Begin Function number 28 def comp_rad_from_ratio(term1,term2,last_no) #TOP TWO TERM RADIUS ANALYSIS if (float_abs(term2) > $glob__0) then # if number 14 ret = float_abs(term1 * $glob_h / term2) else ret = $glob_larger_float end# end if 14 return ( ret) #BOTTOM TWO TERM RADIUS ANALYSIS end # End Function number 28 # Begin Function number 29 def comp_ord_from_ratio(term1,term2,last_no) #TOP TWO TERM ORDER ANALYSIS if (float_abs(term2) > $glob__0) then # if number 14 ret = $glob__1 + float_abs(term2) * c(last_no) * ln(float_abs(term1 * $glob_h / term2))/ln(c(last_no)) else ret = $glob_larger_float end# end if 14 return ( ret) #BOTTOM TWO TERM ORDER ANALYSIS end # End Function number 29 # Begin Function number 30 def c(in_val) #To Force Conversion when needed ret = in_val return ( ret) #End Conversion end # End Function number 30 # Begin Function number 31 def comp_rad_from_three_terms(term1,term2,term3,last_no) #TOP THREE TERM RADIUS ANALYSIS temp = float_abs(term2*term2*c(last_no)+$glob__m2*term2*term2-term1*term3*c(last_no)+term1*term3) if (float_abs(temp) > $glob__0) then # if number 14 ret = float_abs((term2*$glob_h*term1)/(temp)) else ret = $glob_larger_float end# end if 14 return ( ret) #BOTTOM THREE TERM RADIUS ANALYSIS end # End Function number 31 # Begin Function number 32 def comp_ord_from_three_terms(term1,term2,term3,last_no) #TOP THREE TERM ORDER ANALYSIS ret = float_abs(($glob__4*term1*term3*c(last_no)-$glob__3*term1*term3-$glob__4*term2*term2*c(last_no)+$glob__4*term2*term2+term2*term2*c(last_no*last_no)-term1*term3*c(last_no*last_no))/(term2*term2*c(last_no)-$glob__2*term2*term2-term1*term3*c(last_no)+term1*term3)) return ( ret) #TOP THREE TERM ORDER ANALYSIS end # End Function number 32 # Begin Function number 33 def comp_rad_from_six_terms(term1,term2,term3,term4,term5,term6,last_no) #TOP SIX TERM RADIUS ANALYSIS if ((term5 != $glob__0) and (term4 != $glob__0) and (term3 != $glob__0) and (term2 != $glob__0) and (term1 != $glob__0)) then # if number 14 rm0 = term6/term5 rm1 = term5/term4 rm2 = term4/term3 rm3 = term3/term2 rm4 = term2/term1 nr1 = c(last_no-1)*rm0 - $glob__2*c(last_no-2)*rm1 + c(last_no-3)*rm2 nr2 = c(last_no-2)*rm1 - $glob__2*c(last_no-3)*rm2 + c(last_no-4)*rm3 dr1 = $glob__m1/rm1 + $glob__2/rm2 - $glob__1/rm3 dr2 = $glob__m1/rm2 + $glob__2/rm3 - $glob__1/rm4 ds1 = $glob__3/rm1 - $glob__8/rm2 + $glob__5/rm3 ds2 = $glob__3/rm2 - $glob__8/rm3 + $glob__5/rm4 if ((float_abs(nr1 * dr2 - nr2 * dr1) == $glob__0) or (float_abs(dr1) == $glob__0)) then # if number 15 rad_c = $glob_larger_float ord_no = $glob_larger_float else if (float_abs(nr1*dr2 - nr2 * dr1) > $glob__0) then # if number 16 rcs = ((ds1*dr2 - ds2*dr1 +dr1*dr2)/(nr1*dr2 - nr2 * dr1)) #(Manuels) rcs = (ds1*dr2 - ds2*dr1)/(nr1*dr2 - nr2 * dr1) ord_no = (rcs*nr1 - ds1)/($glob__2*dr1) -c(last_no)/$glob__2 if (float_abs(rcs) != $glob__0) then # if number 17 if (rcs > $glob__0) then # if number 18 rad_c = sqrt(rcs) * float_abs($glob_h) else rad_c = $glob_larger_float ord_no = $glob_larger_float end# end if 18 else rad_c = $glob_larger_float ord_no = $glob_larger_float end# end if 17 else rad_c = $glob_larger_float ord_no = $glob_larger_float end# end if 16 end# end if 15 else rad_c = $glob_larger_float ord_no = $glob_larger_float end# end if 14 $glob_six_term_ord_save = ord_no return ( rad_c) #BOTTOM SIX TERM RADIUS ANALYSIS end # End Function number 33 # Begin Function number 34 def comp_ord_from_six_terms(term1,term2,term3,term4,term5,term6,last_no) #TOP SIX TERM ORDER ANALYSIS #TOP SAVED FROM SIX TERM RADIUS ANALYSIS return ( $glob_six_term_ord_save) #BOTTOM SIX TERM ORDER ANALYSIS end # End Function number 34 # Begin Function number 35 def factorial_2(nnn) if (nnn > 0) then # if number 14 return ( nnn * factorial_2(nnn - 1)) else return ( 1.0) end# end if 14 end # End Function number 35 # Begin Function number 36 def factorial_1(nnn) if (nnn <= ATS_MAX_TERMS) then # if number 14 if ($array_fact_1[nnn] == 0) then # if number 15 ret = factorial_2(nnn) $array_fact_1[nnn] = ret else ret = $array_fact_1[nnn] end# end if 15 else ret = factorial_2(nnn) end# end if 14 return ( ret) end # End Function number 36 # Begin Function number 37 def factorial_3(mmm,nnn) if ((nnn <= ATS_MAX_TERMS) and (mmm <= ATS_MAX_TERMS)) then # if number 14 if ($array_fact_2[mmm][nnn] == 0) then # if number 15 ret = factorial_1(mmm)/factorial_1(nnn) $array_fact_2[mmm][nnn] = ret else ret = $array_fact_2[mmm][nnn] end# end if 15 else ret = factorial_2(mmm)/factorial_2(nnn) end# end if 14 return ( ret) end # End Function number 37 # Begin Function number 38 def convfloat(mmm) return(mmm.to_f) end # End Function number 38 # Begin Function number 39 def elapsed_time_seconds() t = Time.now return(t.to_i) end # End Function number 39 # Begin Function number 40 def expt(x,y) if ((x <= 0.0) and (y < 0.0)) then # if number 14 puts "expt error x = " + x.to_s + "y = " + y.to_s end# end if 14 return(x**y) end # End Function number 40 # Begin Function number 41 def Si(x) return(0.0) end # End Function number 41 # Begin Function number 42 def Ci(x) return(0.0) end # End Function number 42 # Begin Function number 43 def float_abs(x) return(x.abs) end # End Function number 43 # Begin Function number 44 def int_abs(x) return(x.abs) end # End Function number 44 # Begin Function number 45 def exp(x) return(Math.exp(x)) end # End Function number 45 # Begin Function number 46 def int_trunc(x) return(x.to_i) end # End Function number 46 # Begin Function number 47 def floor(x) return(x.floor) end # End Function number 47 # Begin Function number 48 def sin(x) return(Math.sin(x)) end # End Function number 48 # Begin Function number 49 def neg(x) return(-x) end # End Function number 49 # Begin Function number 50 def cos(x) return(Math.cos(x)) end # End Function number 50 # Begin Function number 51 def tan(x) return(Math.tan(x)) end # End Function number 51 # Begin Function number 52 def arccos(x) return(Math.acos(x)) end # End Function number 52 # Begin Function number 53 def arccosh(x) return(Math.acosh(x)) end # End Function number 53 # Begin Function number 54 def arcsin(x) return(Math.asin(x)) end # End Function number 54 # Begin Function number 55 def arcsinh(x) return(Math.asinh(x)) end # End Function number 55 # Begin Function number 56 def arctan(x) return(Math.atan(x)) end # End Function number 56 # Begin Function number 57 def arctanh(x) return(Math.atanh(x)) end # End Function number 57 # Begin Function number 58 def cosh(x) return(Math.cosh(x)) end # End Function number 58 # Begin Function number 59 def erf(x) return(Math.erf(x)) end # End Function number 59 # Begin Function number 60 def log(x) return(Math.log(x)) end # End Function number 60 # Begin Function number 61 def ln(x) return(Math.log(x)) end # End Function number 61 # Begin Function number 62 def log10(x) return(Math.log10(x)) end # End Function number 62 # Begin Function number 63 def sinh(x) return(Math.sinh(x)) end # End Function number 63 # Begin Function number 64 def tanh(x) return(Math.tanh(x)) end # End Function number 64 # Begin Function number 65 def sqrt(x) return(Math.sqrt(x)) end # End Function number 65 # Begin Function number 66 def array2d(op3,op4) i = 0 x1 = Array.new(op3.to_i + 1) while i <= op3.to_i + 1 do # do number -1 x1[i] = Array.new(op4.to_i + 1) i += 1 end# end do number -1 return x1 end # End Function number 66 # Begin Function number 67 def estimated_needed_step_error(x_start,x_end,estimated_h,estimated_answer) omniout_float(ALWAYS,"$glob_desired_digits_correct",32,$glob_desired_digits_correct,32,"") desired_abs_gbl_error = expt($glob__10,c( -$glob_desired_digits_correct)) * c(float_abs(c(estimated_answer))) omniout_float(ALWAYS,"estimated_h",32,estimated_h,32,"") omniout_float(ALWAYS,"estimated_answer",32,estimated_answer,32,"") omniout_float(ALWAYS,"desired_abs_gbl_error",32,desired_abs_gbl_error,32,"") range = (x_end - x_start) omniout_float(ALWAYS,"range",32,range,32,"") estimated_steps = range / estimated_h omniout_float(ALWAYS,"estimated_steps",32,estimated_steps,32,"") step_error = (c(float_abs(desired_abs_gbl_error) /sqrt(c( estimated_steps))/c(ATS_MAX_TERMS))) omniout_float(ALWAYS,"step_error",32,step_error,32,"") return ( (step_error)) end # End Function number 67 #END ATS LIBRARY BLOCK #BEGIN USER FUNCTION BLOCK #BEGIN BLOCK 3 #BEGIN USER DEF BLOCK def exact_soln_y (x) x = c(x) return(c(5.0) * ln(c(0.1) * c(x) + c(0.2)) * ( c(0.1) * c(x) + c(0.2)) - c(0.5) * c(x) - c(1.0)) end #END USER DEF BLOCK #END BLOCK 3 #END USER FUNCTION BLOCK # before write_aux functions # Begin Function number 2 def display_poles() if (($glob_type_given_pole == 1) or ($glob_type_given_pole == 2)) then # if number 1 rad_given = sqrt(($array_x[1] - $array_given_rad_poles[1][1]) * ($array_x[1] - $array_given_rad_poles[1][1]) + $array_given_rad_poles[1][2] * $array_given_rad_poles[1][2]) omniout_float(ALWAYS,"Radius of convergence (given) for eq 1 ",4,rad_given,4," ") omniout_float(ALWAYS,"Order of pole (given) ",4,$array_given_ord_poles[1][1],4," ") if (rad_given < $glob_least_given_sing) then # if number 2 $glob_least_given_sing = rad_given end# end if 2 elsif ($glob_type_given_pole == 3) then # if number 2 omniout_str(ALWAYS,"NO POLE (given) for Equation 1") elsif ($glob_type_given_pole == 5) then # if number 3 omniout_str(ALWAYS,"SOME POLE (given) for Equation 1") else omniout_str(ALWAYS,"NO INFO (given) for Equation 1") end# end if 3 if ($array_rad_test_poles[1][1] < $glob_large_float) then # if number 3 omniout_float(ALWAYS,"Radius of convergence (ratio test) for eq 1 ",4,$array_rad_test_poles[1][1],4," ") if ($array_rad_test_poles[1][1]< $glob_least_ratio_sing) then # if number 4 $glob_least_ratio_sing = $array_rad_test_poles[1][1] end# end if 4 omniout_float(ALWAYS,"Order of pole (ratio test) ",4, $array_ord_test_poles[1][1],4," ") else omniout_str(ALWAYS,"NO POLE (ratio test) for Equation 1") end# end if 3 if (($array_rad_test_poles[1][2] > $glob__small) and ($array_rad_test_poles[1][2] < $glob_large_float)) then # if number 3 omniout_float(ALWAYS,"Radius of convergence (three term test) for eq 1 ",4,$array_rad_test_poles[1][2],4," ") if ($array_rad_test_poles[1][2]< $glob_least_3_sing) then # if number 4 $glob_least_3_sing = $array_rad_test_poles[1][2] end# end if 4 omniout_float(ALWAYS,"Order of pole (three term test) ",4, $array_ord_test_poles[1][2],4," ") else omniout_str(ALWAYS,"NO REAL POLE (three term test) for Equation 1") end# end if 3 if (($array_rad_test_poles[1][3] > $glob__small) and ($array_rad_test_poles[1][3] < $glob_large_float)) then # if number 3 omniout_float(ALWAYS,"Radius of convergence (six term test) for eq 1 ",4,$array_rad_test_poles[1][3],4," ") if ($array_rad_test_poles[1][3]< $glob_least_6_sing) then # if number 4 $glob_least_6_sing = $array_rad_test_poles[1][3] end# end if 4 omniout_float(ALWAYS,"Order of pole (six term test) ",4, $array_ord_test_poles[1][3],4," ") else omniout_str(ALWAYS,"NO COMPLEX POLE (six term test) for Equation 1") end# end if 3 end # End Function number 2 # Begin Function number 3 def my_check_sign( x0 ,xf) if (xf > x0) then # if number 3 ret = $glob__1 else ret = $glob__m1 end# end if 3 return ( ret) end # End Function number 3 # Begin Function number 4 def est_size_answer() min_size = $glob_estimated_size_answer if (float_abs($array_y[1]) < min_size) then # if number 3 min_size = float_abs($array_y[1]) omniout_float(ALWAYS,"min_size",32,min_size,32,"") end# end if 3 if (min_size < $glob__1) then # if number 3 min_size = $glob__1 omniout_float(ALWAYS,"min_size",32,min_size,32,"") end# end if 3 return ( min_size) end # End Function number 4 # Begin Function number 5 def test_suggested_h() max_estimated_step_error = $glob__small no_terms = ATS_MAX_TERMS hn_div_ho = $glob__0_5 hn_div_ho_2 = $glob__0_25 hn_div_ho_3 = $glob__0_125 omniout_float(ALWAYS,"hn_div_ho",32,hn_div_ho,32,"") omniout_float(ALWAYS,"hn_div_ho_2",32,hn_div_ho_2,32,"") omniout_float(ALWAYS,"hn_div_ho_3",32,hn_div_ho_3,32,"") est_tmp = float_abs($array_y[no_terms-3] + $array_y[no_terms - 2] * hn_div_ho + $array_y[no_terms - 1] * hn_div_ho_2 + $array_y[no_terms] * hn_div_ho_3) if (est_tmp >= max_estimated_step_error) then # if number 3 max_estimated_step_error = est_tmp end# end if 3 omniout_float(ALWAYS,"max_estimated_step_error",32,max_estimated_step_error,32,"") return ( max_estimated_step_error) end # End Function number 5 # Begin Function number 6 def track_estimated_error() no_terms = ATS_MAX_TERMS hn_div_ho = $glob__0_5 hn_div_ho_2 = $glob__0_25 hn_div_ho_3 = $glob__0_125 est_tmp = c(float_abs($array_y[no_terms-3])) + c(float_abs($array_y[no_terms - 2])) * c(hn_div_ho) + c(float_abs($array_y[no_terms - 1])) * c(hn_div_ho_2) + c(float_abs($array_y[no_terms])) * c(hn_div_ho_3) if ($glob_prec * c(float_abs($array_y[1])) > c(est_tmp)) then # if number 3 est_tmp = c($glob_prec) * c(float_abs($array_y[1])) end# end if 3 if (c(est_tmp) >= c($array_max_est_error[1])) then # if number 3 $array_max_est_error[1] = c(est_tmp) end# end if 3 end # End Function number 6 # Begin Function number 7 def reached_interval() if (($glob_check_sign * $array_x[1]) >= ($glob_check_sign * $glob_next_display - $glob_h/$glob__10)) then # if number 3 ret = true else ret = false end# end if 3 return(ret) end # End Function number 7 # Begin Function number 8 def display_alot(iter) #TOP DISPLAY ALOT if (reached_interval()) then # if number 3 if (iter >= 0) then # if number 4 ind_var = $array_x[1] omniout_float(ALWAYS,"x[1] ",33,ind_var,20," ") closed_form_val_y = exact_soln_y(ind_var) omniout_float(ALWAYS,"y[1] (closed_form) ",33,closed_form_val_y,20," ") term_no = 1 numeric_val = $array_y[term_no] abserr = float_abs(numeric_val - closed_form_val_y) omniout_float(ALWAYS,"y[1] (numeric) ",33,numeric_val,20," ") if (c(float_abs(closed_form_val_y)) > c($glob_prec)) then # if number 5 relerr = abserr*$glob__100/float_abs(closed_form_val_y) if (c(relerr) > c($glob_prec)) then # if number 6 $glob_good_digits = -int_trunc(log10(c(relerr))) + 3 else $glob_good_digits = 16 end# end if 6 else relerr = $glob__m1 $glob_good_digits = -16 end# end if 5 if ($glob_good_digits < $glob_min_good_digits) then # if number 5 $glob_min_good_digits = $glob_good_digits end# end if 5 if ($glob_apfp_est_good_digits < $glob_min_apfp_est_good_digits) then # if number 5 $glob_min_apfp_est_good_digits = $glob_apfp_est_good_digits end# end if 5 if (c(float_abs(numeric_val)) > c($glob_prec)) then # if number 5 est_rel_err = c($array_max_est_error[1])*$glob__100 * c(sqrt($glob_iter))*c(24)*c(ATS_MAX_TERMS)/c(float_abs(numeric_val)) if (est_rel_err > $glob_prec) then # if number 6 $glob_est_digits = -int_trunc(log10(est_rel_err)) + 3 else $glob_est_digits = 16 end# end if 6 else relerr = $glob__m1 $glob_est_digits = -16 end# end if 5 $array_est_digits[1] = $glob_est_digits if ($glob_iter == 1) then # if number 5 $array_1st_rel_error[1] = relerr else $array_last_rel_error[1] = relerr end# end if 5 $array_est_rel_error[1] = est_rel_err omniout_float(ALWAYS,"absolute error ",4,abserr,20," ") omniout_float(ALWAYS,"relative error ",4,relerr,20,"%") omniout_int(INFO,"Desired digits ",32,$glob_desired_digits_correct,4," ") omniout_int(INFO,"Estimated correct digits ",32,$glob_est_digits,4," ") omniout_int(INFO,"Correct digits ",32,$glob_good_digits,4," ") omniout_float(ALWAYS,"h ",4,$glob_h,20," ") end# end if 4 #BOTTOM DISPLAY ALOT end# end if 3 end # End Function number 8 # Begin Function number 9 def prog_report(x_start,x_end) #TOP PROGRESS REPORT clock_sec1 = elapsed_time_seconds() total_clock_sec = (clock_sec1) - ($glob_orig_start_sec) $glob_clock_sec = (clock_sec1) - ($glob_clock_start_sec) left_sec = ($glob_max_sec) + ($glob_orig_start_sec) - (clock_sec1) expect_sec = comp_expect_sec((x_end),(x_start),($array_x[1]) + ($glob_h) ,( clock_sec1) - ($glob_orig_start_sec)) opt_clock_sec = ( clock_sec1) - ($glob_optimal_clock_start_sec) $glob_optimal_expect_sec = comp_expect_sec((x_end),(x_start),($array_x[1]) +( $glob_h) ,( opt_clock_sec)) $glob_total_exp_sec = $glob_optimal_expect_sec + c(total_clock_sec) percent_done = comp_percent((x_end),(x_start),($array_x[1]) + ($glob_h)) $glob_percent_done = percent_done omniout_str_noeol(INFO,"Total Elapsed Time ") omniout_timestr((total_clock_sec)) omniout_str_noeol(INFO,"Elapsed Time(since restart) ") omniout_timestr(($glob_clock_sec)) if (c(percent_done) < $glob__100) then # if number 3 omniout_str_noeol(INFO,"Expected Time Remaining ") omniout_timestr((expect_sec)) omniout_str_noeol(INFO,"Optimized Time Remaining ") omniout_timestr(($glob_optimal_expect_sec)) omniout_str_noeol(INFO,"Expected Total Time ") omniout_timestr(($glob_total_exp_sec)) end# end if 3 omniout_str_noeol(INFO,"Time to Timeout ") omniout_timestr((left_sec)) omniout_float(INFO, "Percent Done ",33,percent_done,4,"%") #BOTTOM PROGRESS REPORT end # End Function number 9 # Begin Function number 10 def check_for_pole() #TOP CHECK FOR POLE tmp_rad = $glob_larger_float prev_tmp_rad = $glob_larger_float tmp_ratio = $glob_larger_float rad_c = $glob_larger_float $array_rad_test_poles[1][1] = $glob_larger_float $array_ord_test_poles[1][1] = $glob_larger_float found_sing = 1 last_no = ATS_MAX_TERMS - 1 - 10 cnt = 0 while (last_no < ATS_MAX_TERMS-3 and found_sing == 1) do # do number 1 tmp_rad = comp_rad_from_ratio($array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) if (float_abs(prev_tmp_rad) > $glob__0) then # if number 3 tmp_ratio = tmp_rad / prev_tmp_rad else tmp_ratio = $glob_large_float end# end if 3 if ((cnt > 0 ) and (tmp_ratio < $glob_upper_ratio_limit) and (tmp_ratio > $glob_lower_ratio_limit)) then # if number 3 rad_c = tmp_rad elsif (cnt == 0) then # if number 4 rad_c = tmp_rad elsif (cnt > 0) then # if number 5 found_sing = 0 end# end if 5 prev_tmp_rad = tmp_rad cnt = cnt + 1 last_no = last_no + 1 end# end do number 1 if (found_sing == 1) then # if number 5 if (rad_c < $array_rad_test_poles[1][1]) then # if number 6 $array_rad_test_poles[1][1] = rad_c last_no = last_no - 1 tmp_ord = comp_ord_from_ratio($array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) $array_rad_test_poles[1][1] = rad_c $array_ord_test_poles[1][1] = tmp_ord end# end if 6 end# end if 5 #BOTTOM general radius test1 tmp_rad = $glob_larger_float prev_tmp_rad = $glob_larger_float tmp_ratio = $glob_larger_float rad_c = $glob_larger_float $array_rad_test_poles[1][2] = $glob_larger_float $array_ord_test_poles[1][2] = $glob_larger_float found_sing = 1 last_no = ATS_MAX_TERMS - 1 - 10 cnt = 0 while (last_no < ATS_MAX_TERMS-4 and found_sing == 1) do # do number 1 tmp_rad = comp_rad_from_three_terms($array_y_higher[1][last_no-2],$array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) if (float_abs(prev_tmp_rad) > $glob__0) then # if number 5 tmp_ratio = tmp_rad / prev_tmp_rad else tmp_ratio = $glob_large_float end# end if 5 if ((cnt > 0 ) and (tmp_ratio < $glob_upper_ratio_limit) and (tmp_ratio > $glob_lower_ratio_limit)) then # if number 5 rad_c = tmp_rad elsif (cnt == 0) then # if number 6 rad_c = tmp_rad elsif (cnt > 0) then # if number 7 found_sing = 0 end# end if 7 prev_tmp_rad = tmp_rad cnt = cnt + 1 last_no = last_no + 1 end# end do number 1 if (found_sing == 1) then # if number 7 if (rad_c < $array_rad_test_poles[1][2]) then # if number 8 $array_rad_test_poles[1][2] = rad_c last_no = last_no - 1 tmp_ord = comp_ord_from_three_terms($array_y_higher[1][last_no-2],$array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) $array_rad_test_poles[1][2] = rad_c if (rad_c < $glob_min_pole_est) then # if number 9 $glob_min_pole_est = rad_c end# end if 9 $array_ord_test_poles[1][2] = tmp_ord end# end if 8 end# end if 7 #BOTTOM general radius test1 tmp_rad = $glob_larger_float prev_tmp_rad = $glob_larger_float tmp_ratio = $glob_larger_float rad_c = $glob_larger_float $array_rad_test_poles[1][3] = $glob_larger_float $array_ord_test_poles[1][3] = $glob_larger_float found_sing = 1 last_no = ATS_MAX_TERMS - 1 - 10 cnt = 0 while (last_no < ATS_MAX_TERMS-7 and found_sing == 1) do # do number 1 tmp_rad = comp_rad_from_six_terms($array_y_higher[1][last_no-5],$array_y_higher[1][last_no-4],$array_y_higher[1][last_no-3],$array_y_higher[1][last_no-2],$array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) if (float_abs(prev_tmp_rad) > $glob__0) then # if number 7 tmp_ratio = tmp_rad / prev_tmp_rad else tmp_ratio = $glob_large_float end# end if 7 if ((cnt > 0 ) and (tmp_ratio < $glob_upper_ratio_limit) and (tmp_ratio > $glob_lower_ratio_limit)) then # if number 7 rad_c = tmp_rad elsif (cnt == 0) then # if number 8 rad_c = tmp_rad elsif (cnt > 0) then # if number 9 found_sing = 0 end# end if 9 prev_tmp_rad = tmp_rad cnt = cnt + 1 last_no = last_no + 1 end# end do number 1 if (found_sing == 1) then # if number 9 if (rad_c < $array_rad_test_poles[1][3]) then # if number 10 $array_rad_test_poles[1][3] = rad_c last_no = last_no - 1 tmp_ord = comp_ord_from_six_terms($array_y_higher[1][last_no-5],$array_y_higher[1][last_no-4],$array_y_higher[1][last_no-3],$array_y_higher[1][last_no-2],$array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) $array_rad_test_poles[1][3] = rad_c if (rad_c < $glob_min_pole_est) then # if number 11 $glob_min_pole_est = rad_c end# end if 11 $array_ord_test_poles[1][3] = tmp_ord end# end if 10 end# end if 9 #BOTTOM general radius test1 #START ADJUST ALL SERIES if (float_abs($glob_min_pole_est) * $glob_ratio_of_radius < float_abs($glob_h)) then # if number 9 h_new = $glob_check_sign * $glob_min_pole_est * $glob_ratio_of_radius omniout_str(ALWAYS,"SETTING H FOR POLE") $glob_h_reason = 6 if ($glob_check_sign * $glob_min_h > $glob_check_sign * h_new) then # if number 10 omniout_str(ALWAYS,"SETTING H FOR MIN H") h_new = $glob_min_h $glob_h_reason = 5 end# end if 10 term = 1 ratio = 1.0 while (term <= ATS_MAX_TERMS) do # do number 1 $array_y[term] = $array_y[term]* ratio $array_y_higher[1][term] = $array_y_higher[1][term]* ratio $array_x[term] = $array_x[term]* ratio ratio = ratio * h_new / float_abs($glob_h) term = term + 1 end# end do number 1 $glob_h = h_new end# end if 9 #BOTTOM ADJUST ALL SERIES if (reached_interval()) then # if number 9 display_poles() end# end if 9 end # End Function number 10 # Begin Function number 11 def atomall() #TOP ATOMALL #END OUTFILE1 #BEGIN OUTFILE2 #END OUTFILE2 #BEGIN ATOMHDR1 #emit pre mult CONST - LINEAR $eq_no = 1 i = 1 $array_tmp1[1] = $array_const_0D1[1] * $array_x[1] #emit pre add LINEAR - CONST $eq_no = 1 i = 1 $array_tmp2[1] = $array_tmp1[1] + $array_const_0D2[1] #emit pre sqrt 1 $eq_no = 1 $array_tmp3[1] = sqrt($array_tmp2[1]) #emit pre ln 1 FULL $eq_no = 1 $array_tmp4[1] = ln($array_tmp3[1]) #emit pre add CONST FULL $eq_no = 1 i = 1 $array_tmp5[1] = $array_const_0D0[1] + $array_tmp4[1] #emit pre assign xxx $eq_no = 1 i = 1 $min_hdrs = 5 if ( not $array_y_set_initial[1][2]) then # if number 1 if (1 <= ATS_MAX_TERMS) then # if number 2 temporary = c($array_tmp5[1]) * (expt(($glob_h) , c(1))) * c(factorial_3(0,1)) if (2 <= ATS_MAX_TERMS) then # if number 3 $array_y[2] = temporary $array_y_higher[1][2] = temporary end# end if 3 temporary = c(temporary) / c($glob_h) * c(1) $array_y_higher[2][1] = c(temporary) end# end if 2 end# end if 1 kkk = 2 #END ATOMHDR1 #BEGIN ATOMHDR2 #emit pre mult CONST - LINEAR $eq_no = 1 i = 2 $array_tmp1[2] = $array_const_0D1[1] * $array_x[2] #emit pre add LINEAR - CONST $eq_no = 1 i = 2 $array_tmp2[2] = $array_tmp1[2] #emit pre sqrt 2 $eq_no = 1 $array_tmp3[2] = $array_tmp2[2] / $array_tmp3[1]/$glob__2 #emit pre ln 2 FULL $eq_no = 1 $array_tmp4[2] = $array_tmp3[2] / $array_tmp3[1] #emit pre add CONST FULL $eq_no = 1 i = 2 $array_tmp5[2] = $array_tmp4[2] #emit pre assign xxx $eq_no = 1 i = 2 $min_hdrs = 5 if ( not $array_y_set_initial[1][3]) then # if number 1 if (2 <= ATS_MAX_TERMS) then # if number 2 temporary = c($array_tmp5[2]) * (expt(($glob_h) , c(1))) * c(factorial_3(1,2)) if (3 <= ATS_MAX_TERMS) then # if number 3 $array_y[3] = temporary $array_y_higher[1][3] = temporary end# end if 3 temporary = c(temporary) / c($glob_h) * c(2) $array_y_higher[2][2] = c(temporary) end# end if 2 end# end if 1 kkk = 3 #END ATOMHDR2 #BEGIN ATOMHDR3 #emit pre sqrt ID_LINEAR iii = 3 $eq_no = 1 $array_tmp3[3] = 0 $array_tmp3[3] = neg(ats(3,$array_tmp3,$array_tmp3,2)) / $array_tmp3[1] /$glob__2 #emit pre ln ID_FULL iii = 3 $eq_no = 1 #emit pre ln 3 $eq_no = 1 $array_tmp4[3] = ( $array_tmp3[3] - att(2,$array_tmp3,$array_tmp4,2) ) / $array_tmp3[1] #emit pre add CONST FULL $eq_no = 1 i = 3 $array_tmp5[3] = $array_tmp4[3] #emit pre assign xxx $eq_no = 1 i = 3 $min_hdrs = 5 if ( not $array_y_set_initial[1][4]) then # if number 1 if (3 <= ATS_MAX_TERMS) then # if number 2 temporary = c($array_tmp5[3]) * (expt(($glob_h) , c(1))) * c(factorial_3(2,3)) if (4 <= ATS_MAX_TERMS) then # if number 3 $array_y[4] = temporary $array_y_higher[1][4] = temporary end# end if 3 temporary = c(temporary) / c($glob_h) * c(3) $array_y_higher[2][3] = c(temporary) end# end if 2 end# end if 1 kkk = 4 #END ATOMHDR3 #BEGIN ATOMHDR4 #emit pre sqrt ID_LINEAR iii = 4 $eq_no = 1 $array_tmp3[4] = 0 $array_tmp3[4] = neg(ats(4,$array_tmp3,$array_tmp3,2)) / $array_tmp3[1] /$glob__2 #emit pre ln ID_FULL iii = 4 $eq_no = 1 #emit pre ln 4 $eq_no = 1 $array_tmp4[4] = ( $array_tmp3[4] - att(3,$array_tmp3,$array_tmp4,2) ) / $array_tmp3[1] #emit pre add CONST FULL $eq_no = 1 i = 4 $array_tmp5[4] = $array_tmp4[4] #emit pre assign xxx $eq_no = 1 i = 4 $min_hdrs = 5 if ( not $array_y_set_initial[1][5]) then # if number 1 if (4 <= ATS_MAX_TERMS) then # if number 2 temporary = c($array_tmp5[4]) * (expt(($glob_h) , c(1))) * c(factorial_3(3,4)) if (5 <= ATS_MAX_TERMS) then # if number 3 $array_y[5] = temporary $array_y_higher[1][5] = temporary end# end if 3 temporary = c(temporary) / c($glob_h) * c(4) $array_y_higher[2][4] = c(temporary) end# end if 2 end# end if 1 kkk = 5 #END ATOMHDR4 #BEGIN ATOMHDR5 #emit pre sqrt ID_LINEAR iii = 5 $eq_no = 1 $array_tmp3[5] = 0 $array_tmp3[5] = neg(ats(5,$array_tmp3,$array_tmp3,2)) / $array_tmp3[1] /$glob__2 #emit pre ln ID_FULL iii = 5 $eq_no = 1 #emit pre ln 5 $eq_no = 1 $array_tmp4[5] = ( $array_tmp3[5] - att(4,$array_tmp3,$array_tmp4,2) ) / $array_tmp3[1] #emit pre add CONST FULL $eq_no = 1 i = 5 $array_tmp5[5] = $array_tmp4[5] #emit pre assign xxx $eq_no = 1 i = 5 $min_hdrs = 5 if ( not $array_y_set_initial[1][6]) then # if number 1 if (5 <= ATS_MAX_TERMS) then # if number 2 temporary = c($array_tmp5[5]) * (expt(($glob_h) , c(1))) * c(factorial_3(4,5)) if (6 <= ATS_MAX_TERMS) then # if number 3 $array_y[6] = temporary $array_y_higher[1][6] = temporary end# end if 3 temporary = c(temporary) / c($glob_h) * c(5) $array_y_higher[2][5] = c(temporary) end# end if 2 end# end if 1 kkk = 6 #END ATOMHDR5 #BEGIN OUTFILE3 #Top Atomall While Loop-- outfile3 while (kkk <= ATS_MAX_TERMS) do # do number 1 #END OUTFILE3 #BEGIN OUTFILE4 #emit sqrt LINEAR $eq_no = 1 $array_tmp3[kkk] = 0 $array_tmp3[kkk] = neg(ats(kkk,$array_tmp3,$array_tmp3,2)) /$array_tmp3[1] / $glob__2 #emit ln FULL $eq_no = 1 $array_tmp4[kkk] = ($array_tmp3[kkk] - att(kkk-1,$array_tmp3,$array_tmp4,2))/$array_tmp3[1] #emit NOT FULL - FULL add $eq_no = 1 $array_tmp5[kkk] = $array_tmp4[kkk] #emit assign $eq_no = 1 order_d = 1 if (kkk + order_d <= ATS_MAX_TERMS) then # if number 1 if ( not $array_y_set_initial[1][kkk + order_d]) then # if number 2 temporary = c($array_tmp5[kkk]) * expt(($glob_h) , c(order_d)) * c(factorial_3((kkk - 1),(kkk + order_d - 1))) $array_y[kkk + order_d] = c(temporary) $array_y_higher[1][kkk + order_d] = c(temporary) term = kkk + order_d - 1 adj2 = kkk + order_d - 1 adj3 = 2 while ((term >= 1) and (term <= ATS_MAX_TERMS) and (adj3 < order_d + 1)) do # do number 1 if (adj3 <= order_d + 1) then # if number 3 if (adj2 > 0) then # if number 4 temporary = c(temporary) / c($glob_h) * c(adj2) else temporary = c(temporary) end# end if 4 $array_y_higher[adj3][term] = c(temporary) end# end if 3 term = term - 1 adj2 = adj2 - 1 adj3 = adj3 + 1 end# end do number 1 end# end if 2 end# end if 1 kkk = kkk + 1 end# end do number 1 #BOTTOM ATOMALL #END OUTFILE4 #BEGIN OUTFILE5 #BOTTOM ATOMALL ??? end # End Function number 12 #END OUTFILE5 # before write_ruby_main_top # Begin Function number 12 def main() #BEGIN OUTFILEMAIN # before first input block #BEGIN FIRST INPUT BLOCK #BEGIN BLOCK 1 #BEGIN FIRST INPUT BLOCK # Digits:=32; ELIMINATED in preodein.rb max_terms=40 #END BLOCK 1 #END FIRST INPUT BLOCK #START OF INITS AFTER INPUT BLOCK $glob_html_log = true #END OF INITS AFTER INPUT BLOCK # before generate arrays $array_y_init= Array.new(ATS_MAX_TERMS + 1) $array_norms= Array.new(ATS_MAX_TERMS + 1) $array_fact_1= Array.new(ATS_MAX_TERMS + 1) $array_1st_rel_error= Array.new(2 + 1) $array_last_rel_error= Array.new(2 + 1) $array_est_rel_error= Array.new(2 + 1) $array_max_est_error= Array.new(2 + 1) $array_type_pole= Array.new(2 + 1) $array_type_real_pole= Array.new(2 + 1) $array_type_complex_pole= Array.new(2 + 1) $array_est_digits= Array.new(2 + 1) $array_y= Array.new(ATS_MAX_TERMS + 1) $array_x= Array.new(ATS_MAX_TERMS + 1) $array_tmp0= Array.new(ATS_MAX_TERMS + 1) $array_tmp1= Array.new(ATS_MAX_TERMS + 1) $array_tmp2= Array.new(ATS_MAX_TERMS + 1) $array_tmp3= Array.new(ATS_MAX_TERMS + 1) $array_tmp4= Array.new(ATS_MAX_TERMS + 1) $array_tmp5= Array.new(ATS_MAX_TERMS + 1) $array_m1= Array.new(ATS_MAX_TERMS + 1) $array_y_higher = array2d(2 + 1 ,ATS_MAX_TERMS + 1) $array_y_higher_work = array2d(2 + 1 ,ATS_MAX_TERMS + 1) $array_y_higher_work2 = array2d(2 + 1 ,ATS_MAX_TERMS + 1) $array_y_set_initial = array2d(2 + 1 ,ATS_MAX_TERMS + 1) $array_given_rad_poles = array2d(2 + 1 ,3 + 1) $array_given_ord_poles = array2d(2 + 1 ,3 + 1) $array_rad_test_poles = array2d(2 + 1 ,4 + 1) $array_ord_test_poles = array2d(2 + 1 ,4 + 1) $array_fact_2 = array2d(ATS_MAX_TERMS + 1 ,ATS_MAX_TERMS + 1) # before generate constants # before generate globals definition #Top Generate Globals Definition $glob__small=c(0.1e-50) $glob_small_float=c(0.1e-50) $glob_smallish_float=c(0.1e-60) $glob_large_float=c(1.0e100) $glob_larger_float=c(1.1e100) $glob__m2=c(-2) $glob__m1=c(-1) $glob__0=c(0) $glob__1=c(1) $glob__2=c(2) $glob__3=c(3) $glob__4=c(4) $glob__5=c(5) $glob__8=c(8) $glob__10=c(10) $glob__100=c(100) $glob__pi=c(0.0) $glob__0_5=c(0.5) $glob__0_8=c(0.8) $glob__m0_8=c(-0.8) $glob__0_25=c(0.25) $glob__0_125=c(0.125) $glob_prec=c(1.0e-16) $glob_est_digits=1 $glob_check_sign=c(1.0) $glob_desired_digits_correct=c(8.0) $glob_max_estimated_step_error=c(0.0) $glob_ratio_of_radius=c(0.1) $glob_percent_done=c(0.0) $glob_subiter_method=3 $glob_total_exp_sec=c(0.1) $glob_optimal_expect_sec=c(0.1) $glob_estimated_size_answer=c(100.0) $glob_html_log=true $glob_min_good_digits=99999 $glob_good_digits=0 $glob_min_apfp_est_good_digits=99999 $glob_apfp_est_good_digits=0 $glob_max_opt_iter=10 $glob_dump=false $glob_djd_debug=true $glob_display_flag=true $glob_djd_debug2=true $glob_h_reason=0 $glob_sec_in_minute=60 $glob_min_in_hour=60 $glob_hours_in_day=24 $glob_days_in_year=365 $glob_sec_in_hour=3600 $glob_sec_in_day=86400 $glob_sec_in_year=31536000 $glob_almost_1=c(0.9990) $glob_clock_sec=c(0.0) $glob_clock_start_sec=c(0.0) $glob_not_yet_finished=true $glob_initial_pass=true $glob_not_yet_start_msg=true $glob_reached_optimal_h=false $glob_optimal_done=false $glob_disp_incr=c(0.1) $glob_h=c(0.1) $glob_diff_rc_fm=c(0.1) $glob_diff_rc_fmm1=c(0.1) $glob_diff_rc_fmm2=c(0.1) $glob_diff_ord_fm=c(0.1) $glob_diff_ord_fmm1=c(0.1) $glob_diff_ord_fmm2=c(0.1) $glob_six_term_ord_save=c(0.1) $glob_guess_error_rc=c(0.1) $glob_guess_error_ord=c(0.1) $glob_max_h=c(0.1) $glob_min_h=c(0.000001) $glob_type_given_pole=0 $glob_least_given_sing=c(9.9e200) $glob_least_ratio_sing=c(9.9e200) $glob_least_3_sing=c(9.9e100) $glob_least_6_sing=c(9.9e100) $glob_last_good_h=c(0.1) $glob_optimize=false $glob_look_poles=false $glob_display_interval=c(0.1) $glob_next_display=c(0.0) $glob_dump_closed_form=false $glob_abserr=c(0.1e-10) $glob_relerr=c(0.1e-10) $glob_min_pole_est=c(0.1e+10) $glob_max_hours=c(0.0) $glob_max_iter=1000 $glob_max_rel_trunc_err=c(0.1e-10) $glob_max_trunc_err=c(0.1e-10) $glob_no_eqs=0 $glob_optimal_clock_start_sec=c(0.0) $glob_optimal_start=c(0.0) $glob_upper_ratio_limit=c(1.0001) $glob_lower_ratio_limit=c(0.9999) $glob_unchanged_h_cnt=0 $glob_warned=false $glob_warned2=false $glob_max_sec=c(10000.0) $glob_orig_start_sec=c(0.0) $glob_start=0 $glob_iter=0 $glob_normmax=c(0.0) $glob_max_minutes=c(0.0) #Bottom Generate Globals Deninition # before generate const definition $array_const_1= Array.new(ATS_MAX_TERMS + 1) $array_const_0D0= Array.new(ATS_MAX_TERMS + 1) $array_const_0D1= Array.new(ATS_MAX_TERMS + 1) $array_const_0D2= Array.new(ATS_MAX_TERMS + 1) # before arrays initialized term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_y_init[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_norms[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_fact_1[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_1st_rel_error[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_last_rel_error[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_est_rel_error[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_max_est_error[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_type_pole[term] = 0 term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_type_real_pole[term] = 0 term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_type_complex_pole[term] = 0 term = term + 1 end# end do number 1 term = 1 while (term <= 2) do # do number 1 $array_est_digits[term] = 0 term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_y[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_x[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_tmp0[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_tmp1[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_tmp2[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_tmp3[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_tmp4[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_tmp5[term] = c(0.0) term = term + 1 end# end do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 1 $array_m1[term] = c(0.0) term = term + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 2 $array_y_higher[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 2 $array_y_higher_work[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 2 $array_y_higher_work2[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 2 $array_y_set_initial[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= 3) do # do number 2 $array_given_rad_poles[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= 3) do # do number 2 $array_given_ord_poles[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= 4) do # do number 2 $array_rad_test_poles[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=2) do # do number 1 term = 1 while (term <= 4) do # do number 2 $array_ord_test_poles[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 ord = 1 while (ord <=ATS_MAX_TERMS) do # do number 1 term = 1 while (term <= ATS_MAX_TERMS) do # do number 2 $array_fact_2[ord][term] = c(0.0) term = term + 1 end# end do number 2 ord = ord + 1 end# end do number 1 # before symbols initialized #BEGIN SYMBOLS INITIALIZATED zero_ats_ar($array_y) zero_ats_ar($array_x) zero_ats_ar($array_tmp0) zero_ats_ar($array_tmp1) zero_ats_ar($array_tmp2) zero_ats_ar($array_tmp3) zero_ats_ar($array_tmp4) zero_ats_ar($array_tmp5) zero_ats_ar($array_m1) zero_ats_ar($array_const_1) $array_const_1[1] = c(1) zero_ats_ar($array_const_0D0) $array_const_0D0[1] = c(0.0) zero_ats_ar($array_const_0D1) $array_const_0D1[1] = c(0.1) zero_ats_ar($array_const_0D2) $array_const_0D2[1] = c(0.2) zero_ats_ar($array_m1) $array_m1[1] = $glob__m1 #END SYMBOLS INITIALIZATED # before generate factorials init #Initing Factorial Tables iiif = 0 while (iiif <= ATS_MAX_TERMS) do # do number 1 jjjf = 0 while (jjjf <= ATS_MAX_TERMS) do # do number 2 $array_fact_1[iiif] = 0 $array_fact_2[iiif][jjjf] = 0 jjjf = jjjf + 1 end# end do number 2 iiif = iiif + 1 end# end do number 1 #Done Initing Factorial Table # before generate constants assign $glob_iolevel = 5 $glob_yes_pole = 4 $glob_no_pole = 3 $glob_not_given = 0 $glob_no_sing_tests = 4 $glob_ratio_test = 1 $glob_three_term_test = 2 $glob_six_term_test = 3 $glob_log_10 = log(c(10.0)) # before generate globals assign $glob__small = c(0.1e-50) $glob_small_float = c(0.1e-50) $glob_smallish_float = c(0.1e-60) $glob_large_float = c(1.0e100) $glob_larger_float = c(1.1e100) $glob__m2 = c(-2) $glob__m1 = c(-1) $glob__0 = c(0) $glob__1 = c(1) $glob__2 = c(2) $glob__3 = c(3) $glob__4 = c(4) $glob__5 = c(5) $glob__8 = c(8) $glob__10 = c(10) $glob__100 = c(100) $glob__pi = c(0.0) $glob__0_5 = c(0.5) $glob__0_8 = c(0.8) $glob__m0_8 = c(-0.8) $glob__0_25 = c(0.25) $glob__0_125 = c(0.125) $glob_prec = c(1.0e-16) $glob_est_digits = 1 $glob_check_sign = c(1.0) $glob_desired_digits_correct = c(8.0) $glob_max_estimated_step_error = c(0.0) $glob_ratio_of_radius = c(0.1) $glob_percent_done = c(0.0) $glob_subiter_method = 3 $glob_total_exp_sec = c(0.1) $glob_optimal_expect_sec = c(0.1) $glob_estimated_size_answer = c(100.0) $glob_html_log = true $glob_min_good_digits = 99999 $glob_good_digits = 0 $glob_min_apfp_est_good_digits = 99999 $glob_apfp_est_good_digits = 0 $glob_max_opt_iter = 10 $glob_dump = false $glob_djd_debug = true $glob_display_flag = true $glob_djd_debug2 = true $glob_h_reason = 0 $glob_sec_in_minute = 60 $glob_min_in_hour = 60 $glob_hours_in_day = 24 $glob_days_in_year = 365 $glob_sec_in_hour = 3600 $glob_sec_in_day = 86400 $glob_sec_in_year = 31536000 $glob_almost_1 = c(0.9990) $glob_clock_sec = c(0.0) $glob_clock_start_sec = c(0.0) $glob_not_yet_finished = true $glob_initial_pass = true $glob_not_yet_start_msg = true $glob_reached_optimal_h = false $glob_optimal_done = false $glob_disp_incr = c(0.1) $glob_h = c(0.1) $glob_diff_rc_fm = c(0.1) $glob_diff_rc_fmm1 = c(0.1) $glob_diff_rc_fmm2 = c(0.1) $glob_diff_ord_fm = c(0.1) $glob_diff_ord_fmm1 = c(0.1) $glob_diff_ord_fmm2 = c(0.1) $glob_six_term_ord_save = c(0.1) $glob_guess_error_rc = c(0.1) $glob_guess_error_ord = c(0.1) $glob_max_h = c(0.1) $glob_min_h = c(0.000001) $glob_type_given_pole = 0 $glob_least_given_sing = c(9.9e200) $glob_least_ratio_sing = c(9.9e200) $glob_least_3_sing = c(9.9e100) $glob_least_6_sing = c(9.9e100) $glob_last_good_h = c(0.1) $glob_optimize = false $glob_look_poles = false $glob_display_interval = c(0.1) $glob_next_display = c(0.0) $glob_dump_closed_form = false $glob_abserr = c(0.1e-10) $glob_relerr = c(0.1e-10) $glob_min_pole_est = c(0.1e+10) $glob_max_hours = c(0.0) $glob_max_iter = 1000 $glob_max_rel_trunc_err = c(0.1e-10) $glob_max_trunc_err = c(0.1e-10) $glob_no_eqs = 0 $glob_optimal_clock_start_sec = c(0.0) $glob_optimal_start = c(0.0) $glob_upper_ratio_limit = c(1.0001) $glob_lower_ratio_limit = c(0.9999) $glob_unchanged_h_cnt = 0 $glob_warned = false $glob_warned2 = false $glob_max_sec = c(10000.0) $glob_orig_start_sec = c(0.0) $glob_start = 0 $glob_iter = 0 $glob_normmax = c(0.0) $glob_max_minutes = c(0.0) # before generate set diff initial $array_y_set_initial[1][1] = true $array_y_set_initial[1][2] = false $array_y_set_initial[1][3] = false $array_y_set_initial[1][4] = false $array_y_set_initial[1][5] = false $array_y_set_initial[1][6] = false $array_y_set_initial[1][7] = false $array_y_set_initial[1][8] = false $array_y_set_initial[1][9] = false $array_y_set_initial[1][10] = false $array_y_set_initial[1][11] = false $array_y_set_initial[1][12] = false $array_y_set_initial[1][13] = false $array_y_set_initial[1][14] = false $array_y_set_initial[1][15] = false $array_y_set_initial[1][16] = false $array_y_set_initial[1][17] = false $array_y_set_initial[1][18] = false $array_y_set_initial[1][19] = false $array_y_set_initial[1][20] = false $array_y_set_initial[1][21] = false $array_y_set_initial[1][22] = false $array_y_set_initial[1][23] = false $array_y_set_initial[1][24] = false $array_y_set_initial[1][25] = false $array_y_set_initial[1][26] = false $array_y_set_initial[1][27] = false $array_y_set_initial[1][28] = false $array_y_set_initial[1][29] = false $array_y_set_initial[1][30] = false $array_y_set_initial[1][31] = false $array_y_set_initial[1][32] = false $array_y_set_initial[1][33] = false $array_y_set_initial[1][34] = false $array_y_set_initial[1][35] = false $array_y_set_initial[1][36] = false $array_y_set_initial[1][37] = false $array_y_set_initial[1][38] = false $array_y_set_initial[1][39] = false $array_y_set_initial[1][40] = false # set default block #Write Set Defaults $glob_orig_start_sec = elapsed_time_seconds() $glob_display_flag = true $glob_no_eqs = 1 $glob_iter = -1 opt_iter = -1 $glob_max_iter = 50000 $glob_max_hours = (0.0) $glob_max_minutes = (15.0) omniout_str(ALWAYS,"##############ECHO OF PROBLEM#################") omniout_str(ALWAYS,"##############temp/ln_sqrtpostode.ode#################") omniout_str(ALWAYS,"diff ( y , x , 1 ) = ln ( sqrt ( 0.1 * x + 0.2 ) ) ; ") omniout_str(ALWAYS,"!") omniout_str(ALWAYS,"#BEGIN FIRST INPUT BLOCK") omniout_str(ALWAYS,"# Digits:=32; ELIMINATED in preodein.rb") omniout_str(ALWAYS,"max_terms=40 ") omniout_str(ALWAYS,"!") omniout_str(ALWAYS,"#END FIRST INPUT BLOCK") omniout_str(ALWAYS,"#BEGIN SECOND INPUT BLOCK") omniout_str(ALWAYS,"x_start=c(10.0) ") omniout_str(ALWAYS,"x_end=c(11.0) ") omniout_str(ALWAYS,"$array_y_init[0 + 1] = exact_soln_y(x_start)") omniout_str(ALWAYS,"$glob_look_poles=true ") omniout_str(ALWAYS,"$glob_min_h=0.001 ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"$glob_type_given_pole=1 ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"$array_given_rad_poles[1][1]=c(-2.0) ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"$array_given_rad_poles[1][2]=c(0.0) ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"$array_given_ord_poles[1][1]=c(0.5) ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"$array_given_ord_poles[1][2]=c(0.0) ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"#END SECOND INPUT BLOCK") omniout_str(ALWAYS,"#BEGIN OVERRIDE BLOCK") omniout_str(ALWAYS,"$glob_desired_digits_correct=8 ") omniout_str(ALWAYS,"$glob_max_minutes=(3.0) ") omniout_str(ALWAYS,"$glob_subiter_method=3 ") omniout_str(ALWAYS,"$glob_max_iter=100000 ") omniout_str(ALWAYS,"$glob_upper_ratio_limit=c(1.000001) ") omniout_str(ALWAYS,"$glob_lower_ratio_limit=c(0.999999) ") omniout_str(ALWAYS,"$glob_look_poles=false ") omniout_str(ALWAYS,"$glob_h=c(0.001) ") omniout_str(ALWAYS,"$glob_display_interval=c(0.01) ") omniout_str(ALWAYS,"#END OVERRIDE BLOCK") omniout_str(ALWAYS,"!") omniout_str(ALWAYS,"#BEGIN USER DEF BLOCK") omniout_str(ALWAYS,"def exact_soln_y (x)") omniout_str(ALWAYS,"x = c(x)") omniout_str(ALWAYS,"") omniout_str(ALWAYS,"return(c(5.0) * ln(c(0.1) * c(x) + c(0.2)) * ( c(0.1) * c(x) + c(0.2)) - c(0.5) * c(x) - c(1.0)) ") omniout_str(ALWAYS,"end") omniout_str(ALWAYS,"#END USER DEF BLOCK") omniout_str(ALWAYS,"#######END OF ECHO OF PROBLEM#################") $glob_unchanged_h_cnt = 0 $glob_warned = false $glob_warned2 = false $glob_small_float = $glob__0 $glob_smallish_float = $glob__0 $glob_large_float = c(1.0e100) $glob_larger_float = c( 1.1e100) $glob_almost_1 = c( 0.99) # before second block #TOP SECOND INPUT BLOCK #BEGIN SECOND INPUT BLOCK #BEGIN BLOCK 2 #END FIRST INPUT BLOCK #BEGIN SECOND INPUT BLOCK x_start=c(10.0) x_end=c(11.0) $array_y_init[0 + 1] = exact_soln_y(x_start) $glob_look_poles=true $glob_min_h=0.001 # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb $glob_type_given_pole=1 # ELIMINATED in preodein.rb $array_given_rad_poles[1][1]=c(-2.0) # ELIMINATED in preodein.rb $array_given_rad_poles[1][2]=c(0.0) # ELIMINATED in preodein.rb $array_given_ord_poles[1][1]=c(0.5) # ELIMINATED in preodein.rb $array_given_ord_poles[1][2]=c(0.0) # ELIMINATED in preodein.rb #END SECOND INPUT BLOCK #BEGIN OVERRIDE BLOCK $glob_desired_digits_correct=8 $glob_max_minutes=(3.0) $glob_subiter_method=3 $glob_max_iter=100000 $glob_upper_ratio_limit=c(1.000001) $glob_lower_ratio_limit=c(0.999999) $glob_look_poles=false $glob_h=c(0.001) $glob_display_interval=c(0.01) #END OVERRIDE BLOCK #END BLOCK 2 #END SECOND INPUT BLOCK #BEGIN INITS AFTER SECOND INPUT BLOCK $glob_last_good_h = $glob_h $glob_max_sec = (60.0) * ($glob_max_minutes) + (3600.0) * ($glob_max_hours) # after second input block $glob_check_sign = c(my_check_sign(x_start,x_end)) $glob__pi = arccos($glob__m1) $glob_prec = 1.0e-16 if ($glob_optimize) then # if number 9 #BEGIN OPTIMIZE CODE omniout_str(ALWAYS,"START of Optimize") #Start Series -- INITIALIZE FOR OPTIMIZE found_h = false $glob_min_pole_est = $glob_larger_float last_min_pole_est = $glob_larger_float $glob_least_given_sing = $glob_larger_float $glob_least_ratio_sing = $glob_larger_float $glob_least_3_sing = $glob_larger_float $glob_least_6_sing = $glob_larger_float $glob_min_h = float_abs($glob_min_h) * $glob_check_sign $glob_max_h = float_abs($glob_max_h) * $glob_check_sign $glob_h = float_abs($glob_min_h) * $glob_check_sign $glob_display_interval = c((float_abs(c($glob_display_interval))) * ($glob_check_sign)) display_max = c(x_end) - c(x_start)/$glob__10 if (($glob_display_interval) > (display_max)) then # if number 10 $glob_display_interval = c(display_max) end# end if 10 chk_data() min_value = $glob_larger_float est_answer = est_size_answer() opt_iter = 1 est_needed_step_err = estimated_needed_step_error(x_start,x_end,$glob_h,est_answer) omniout_float(ALWAYS,"est_needed_step_err",32,est_needed_step_err,16,"") estimated_step_error = $glob_small_float while ((opt_iter <= 100) and ( not found_h)) do # do number 1 omniout_int(ALWAYS,"opt_iter",32,opt_iter,4,"") $array_x[1] = c(x_start) $array_x[2] = c($glob_h) $glob_next_display = c(x_start) order_diff = 1 #Start Series $array_y term_no = 1 while (term_no <= order_diff) do # do number 2 $array_y[term_no] = $array_y_init[term_no] * expt($glob_h , c(term_no - 1)) / c(factorial_1(term_no - 1)) term_no = term_no + 1 end# end do number 2 rows = order_diff r_order = 1 while (r_order <= rows) do # do number 2 term_no = 1 while (term_no <= (rows - r_order + 1)) do # do number 3 it = term_no + r_order - 1 if (term_no < ATS_MAX_TERMS) then # if number 10 $array_y_higher[r_order][term_no] = $array_y_init[it]* expt($glob_h , c(term_no - 1)) / (c(factorial_1(term_no - 1))) end# end if 10 term_no = term_no + 1 end# end do number 3 r_order = r_order + 1 end# end do number 2 atomall() if ($glob_check_sign * $glob_min_h >= $glob_check_sign * $glob_h) then # if number 10 omniout_str(ALWAYS,"SETTING H FOR MIN H") $glob_h = $glob_check_sign * float_abs($glob_min_h) $glob_h_reason = 1 found_h = true end# end if 10 if ($glob_check_sign * $glob_display_interval <= $glob_check_sign * $glob_h) then # if number 10 omniout_str(ALWAYS,"SETTING H FOR DISPLAY INTERVAL") $glob_h_reason = 2 $glob_h = $glob_display_interval found_h = true end# end if 10 if ($glob_look_poles) then # if number 10 check_for_pole() end# end if 10 if ( not found_h) then # if number 10 est_answer = est_size_answer() est_needed_step_err = estimated_needed_step_error(x_start,x_end,$glob_h,est_answer) omniout_float(ALWAYS,"est_needed_step_err",32,est_needed_step_err,16,"") estimated_step_error = test_suggested_h() omniout_float(ALWAYS,"estimated_step_error",32,estimated_step_error,32,"") if (estimated_step_error < est_needed_step_err) then # if number 11 omniout_str(ALWAYS,"Double H and LOOP") $glob_h = $glob_h*$glob__2 else omniout_str(ALWAYS,"Found H for OPTIMAL") found_h = true $glob_h_reason = 3 $glob_h = $glob_h/$glob__2 end# end if 11 end# end if 10 opt_iter = opt_iter + 1 end# end do number 1 if (( not found_h) and (opt_iter == 1)) then # if number 10 omniout_str(ALWAYS,"Beginning $glob_h too large.") found_h = false end# end if 10 if ($glob_check_sign * $glob_max_h <= $glob_check_sign * $glob_h) then # if number 10 omniout_str(ALWAYS,"SETTING H FOR MAX H") $glob_h = $glob_check_sign * float_abs($glob_max_h) $glob_h_reason = 1 found_h = true end# end if 10 else found_h = true $glob_h = $glob_h * $glob_check_sign end# end if 9 #END OPTIMIZE CODE if ($glob_html_log) then # if number 9 html_log_file = File.new("entry.html","w") end# end if 9 #BEGIN SOLUTION CODE if (found_h) then # if number 9 omniout_str(ALWAYS,"START of Soultion") #Start Series -- INITIALIZE FOR SOLUTION $array_x[1] = c(x_start) $array_x[2] = c($glob_h) $glob_next_display = c(x_start) $glob_min_pole_est = $glob_larger_float $glob_least_given_sing = $glob_larger_float $glob_least_ratio_sing = $glob_larger_float $glob_least_3_sing = $glob_larger_float $glob_least_6_sing = $glob_larger_float order_diff = 1 #Start Series $array_y term_no = 1 while (term_no <= order_diff) do # do number 1 $array_y[term_no] = $array_y_init[term_no] * expt($glob_h , c(term_no - 1)) / c(factorial_1(term_no - 1)) term_no = term_no + 1 end# end do number 1 rows = order_diff r_order = 1 while (r_order <= rows) do # do number 1 term_no = 1 while (term_no <= (rows - r_order + 1)) do # do number 2 it = term_no + r_order - 1 if (term_no < ATS_MAX_TERMS) then # if number 10 $array_y_higher[r_order][term_no] = $array_y_init[it]* expt($glob_h , c(term_no - 1)) / (c(factorial_1(term_no - 1))) end# end if 10 term_no = term_no + 1 end# end do number 2 r_order = r_order + 1 end# end do number 1 current_iter = 1 $glob_clock_start_sec = elapsed_time_seconds() $glob_clock_sec = elapsed_time_seconds() $glob_iter = 0 omniout_str(DEBUGL," ") $glob_reached_optimal_h = true $glob_optimal_clock_start_sec = elapsed_time_seconds() while (($glob_iter < $glob_max_iter) and ($glob_check_sign * $array_x[1] < $glob_check_sign * x_end ) and ((($glob_clock_sec) - ($glob_orig_start_sec)) < ($glob_max_sec))) do # do number 1 #left paren 0001C if (reached_interval()) then # if number 10 omniout_str(INFO," ") omniout_str(INFO,"TOP MAIN SOLVE Loop") end# end if 10 $glob_iter = $glob_iter + 1 $glob_clock_sec = elapsed_time_seconds() track_estimated_error() atomall() track_estimated_error() display_alot(current_iter) if ($glob_look_poles) then # if number 10 check_for_pole() end# end if 10 if (reached_interval()) then # if number 10 $glob_next_display = $glob_next_display + $glob_display_interval end# end if 10 $array_x[1] = $array_x[1] + $glob_h $array_x[2] = $glob_h #Jump Series $array_y order_diff = 2 #START PART 1 SUM AND ADJUST #START SUM AND ADJUST EQ =1 #sum_and_adjust $array_y #BEFORE ADJUST SUBSERIES EQ =1 ord = 2 calc_term = 1 #adjust_subseries$array_y iii = ATS_MAX_TERMS while (iii >= calc_term) do # do number 2 $array_y_higher_work[2][iii] = $array_y_higher[2][iii] / expt($glob_h , c(calc_term - 1)) / c(factorial_3(iii - calc_term , iii - 1)) iii = iii - 1 end# end do number 2 #AFTER ADJUST SUBSERIES EQ =1 #BEFORE SUM SUBSERIES EQ =1 temp_sum = $glob__0 ord = 2 calc_term = 1 #sum_subseries$array_y iii = ATS_MAX_TERMS while (iii >= calc_term) do # do number 2 temp_sum = temp_sum + $array_y_higher_work[ord][iii] iii = iii - 1 end# end do number 2 $array_y_higher_work2[ord][calc_term] = temp_sum * expt($glob_h , c(calc_term - 1)) / c(factorial_1(calc_term - 1)) #AFTER SUM SUBSERIES EQ =1 #BEFORE ADJUST SUBSERIES EQ =1 ord = 1 calc_term = 2 #adjust_subseries$array_y iii = ATS_MAX_TERMS while (iii >= calc_term) do # do number 2 $array_y_higher_work[1][iii] = $array_y_higher[1][iii] / expt($glob_h , c(calc_term - 1)) / c(factorial_3(iii - calc_term , iii - 1)) iii = iii - 1 end# end do number 2 #AFTER ADJUST SUBSERIES EQ =1 #BEFORE SUM SUBSERIES EQ =1 temp_sum = $glob__0 ord = 1 calc_term = 2 #sum_subseries$array_y iii = ATS_MAX_TERMS while (iii >= calc_term) do # do number 2 temp_sum = temp_sum + $array_y_higher_work[ord][iii] iii = iii - 1 end# end do number 2 $array_y_higher_work2[ord][calc_term] = temp_sum * expt($glob_h , c(calc_term - 1)) / c(factorial_1(calc_term - 1)) #AFTER SUM SUBSERIES EQ =1 #BEFORE ADJUST SUBSERIES EQ =1 ord = 1 calc_term = 1 #adjust_subseries$array_y iii = ATS_MAX_TERMS while (iii >= calc_term) do # do number 2 $array_y_higher_work[1][iii] = $array_y_higher[1][iii] / expt($glob_h , c(calc_term - 1)) / c(factorial_3(iii - calc_term , iii - 1)) iii = iii - 1 end# end do number 2 #AFTER ADJUST SUBSERIES EQ =1 #BEFORE SUM SUBSERIES EQ =1 temp_sum = $glob__0 ord = 1 calc_term = 1 #sum_subseries$array_y iii = ATS_MAX_TERMS while (iii >= calc_term) do # do number 2 temp_sum = temp_sum + $array_y_higher_work[ord][iii] iii = iii - 1 end# end do number 2 $array_y_higher_work2[ord][calc_term] = temp_sum * expt($glob_h , c(calc_term - 1)) / c(factorial_1(calc_term - 1)) #AFTER SUM SUBSERIES EQ =1 #END SUM AND ADJUST EQ =1 #END PART 1 #START PART 2 MOVE TERMS to REGULAR Array term_no = ATS_MAX_TERMS while (term_no >= 1) do # do number 2 $array_y[term_no] = $array_y_higher_work2[1][term_no] ord = 1 while (ord <= order_diff) do # do number 3 $array_y_higher[ord][term_no] = $array_y_higher_work2[ord][term_no] ord = ord + 1 end# end do number 3 term_no = term_no - 1 end# end do number 2 #END PART 2 HEVE MOVED TERMS to REGULAR Array end# end do number 1#right paren 0001C omniout_str(ALWAYS,"Finished!") if ($glob_iter >= $glob_max_iter) then # if number 10 omniout_str(ALWAYS,"Maximum Iterations Reached before Solution Completed!") end# end if 10 if (elapsed_time_seconds() - ($glob_orig_start_sec) >= ($glob_max_sec )) then # if number 10 omniout_str(ALWAYS,"Maximum Time Reached before Solution Completed!") end# end if 10 $glob_clock_sec = elapsed_time_seconds() omniout_str(INFO,"diff ( y , x , 1 ) = ln ( sqrt ( 0.1 * x + 0.2 ) ) ; ") omniout_int(INFO,"Iterations ",32,$glob_iter,4," ") prog_report(x_start,x_end) if ($glob_html_log) then # if number 10 logstart(html_log_file) logitem_str(html_log_file,"2015-04-08T21:32:58-05:00") logitem_str(html_log_file,"Ruby") logitem_str(html_log_file,"ln_sqrt") logitem_str(html_log_file,"diff ( y , x , 1 ) = ln ( sqrt ( 0.1 * x + 0.2 ) ) ; ") logitem_float(html_log_file,x_start) logitem_float(html_log_file,x_end) logitem_float(html_log_file,$array_x[1]) logitem_float(html_log_file,$glob_h) logitem_h_reason(html_log_file) logitem_str(html_log_file,"16") logitem_float(html_log_file,$glob_desired_digits_correct) if ($array_est_digits[1] != -16) then # if number 11 logitem_integer(html_log_file,$array_est_digits[1]) else logitem_str(html_log_file,"Unknown") end# end if 11 if ($glob_min_good_digits != -16) then # if number 11 logitem_integer(html_log_file,$glob_min_good_digits) else logitem_str(html_log_file,"Unknown") end# end if 11 if ($glob_good_digits != -16) then # if number 11 logitem_integer(html_log_file,$glob_good_digits) else logitem_str(html_log_file,"Unknown") end# end if 11 logitem_str(html_log_file,"NA") logitem_str(html_log_file,"NA") logitem_integer(html_log_file,ATS_MAX_TERMS) if ($glob_type_given_pole == 0) then # if number 11 logitem_str(html_log_file,"Not Given") logitem_str(html_log_file,"NA") elsif ($glob_type_given_pole == 4) then # if number 12 logitem_str(html_log_file,"No Solution") logitem_str(html_log_file,"NA") elsif ($glob_type_given_pole == 5) then # if number 13 logitem_str(html_log_file,"Some Pole") logitem_str(html_log_file,"????") elsif ($glob_type_given_pole == 3) then # if number 14 logitem_str(html_log_file,"No Pole") logitem_str(html_log_file,"NA") elsif ($glob_type_given_pole == 1) then # if number 15 logitem_str(html_log_file,"Real Sing") logitem_float(html_log_file,$glob_least_given_sing) elsif ($glob_type_given_pole == 2) then # if number 16 logitem_str(html_log_file,"Complex Sing") logitem_float(html_log_file,$glob_least_given_sing) end# end if 16 if ($glob_least_ratio_sing < $glob_large_float) then # if number 16 logitem_float(html_log_file,$glob_least_ratio_sing) else logitem_str(html_log_file,"NONE") end# end if 16 if ($glob_least_3_sing < $glob_large_float) then # if number 16 logitem_float(html_log_file,$glob_least_3_sing) else logitem_str(html_log_file,"NONE") end# end if 16 if ($glob_least_6_sing < $glob_large_float) then # if number 16 logitem_float(html_log_file,$glob_least_6_sing) else logitem_str(html_log_file,"NONE") end# end if 16 logitem_integer(html_log_file,$glob_iter) logitem_time(html_log_file,($glob_clock_sec)) if (c($glob_percent_done) < $glob__100) then # if number 16 logitem_time(html_log_file,($glob_total_exp_sec)) 0 else logitem_str(html_log_file,"Done") 0 end# end if 16 log_revs(html_log_file," 305.fix.digits.div.zero ") logitem_str(html_log_file,"ln_sqrt diffeq.rb") logitem_str(html_log_file,"ln_sqrt Ruby results") logitem_str(html_log_file,"Poor Accuracy") logend(html_log_file) end# end if 15 if ($glob_html_log) then # if number 15 html_log_file.close; end# end if 15 end# end if 14 #END OUTFILEMAIN end # End Function number 12 main()