#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 1 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("") secs_in = secs_in.to_i 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) secs_in = secs_in.to_i 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 2 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 2 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 3 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 3 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 4 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("") if (rel_error != -1.0) then # if number 8 if (rel_error > $glob_prec) then # if number 9 good_digits = 3-int_trunc(log10(rel_error)) file.printf("%d",good_digits) else good_digits = 16 file.printf("%d",good_digits) end# end if 9 else file.printf("Unknown") end# end if 8 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 ( not term2.maybe_zero) 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 ( not term2.maybe_zero) 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 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 ( not temp.maybe_zero) 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 30 # Begin Function number 31 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 31 # Begin Function number 32 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 ( not (nr1*dr2 - nr2 * dr1).maybe_zero) 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 32 # Begin Function number 33 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 33 # Begin Function number 34 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 34 # Begin Function number 35 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 35 # Begin Function number 36 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 36 # Begin Function number 37 def convfloat(mmm) return(mmm.to_f) end # End Function number 37 # Begin Function number 38 def elapsed_time_seconds() t = Time.now return(t.to_i) end # End Function number 38 # Begin Function number 39 def expt(x,y) if ((c(x) <= $glob__0) and (c(y) < $glob__0)) then # if number 14 puts "expt error x = " + x.to_s + "y = " + y.to_s end# end if 14 return(c(x)**c(y)) end # End Function number 39 # Begin Function number 40 def Si(x) return(0.0) end # End Function number 40 # Begin Function number 41 def Ci(x) return(0.0) end # End Function number 41 # Begin Function number 42 def float_abs(x) return(c(x).abs) end # End Function number 42 # Begin Function number 43 def int_abs(x) return(c(x).abs) end # End Function number 43 # Begin Function number 44 def exp(x) return(c(x).exp) end # End Function number 44 # Begin Function number 45 def int_trunc(x) return(c(x).to_i) end # End Function number 45 # Begin Function number 46 def floor(x) return(c(x).floor) end # End Function number 46 # Begin Function number 47 def sin(x) return(c(x).sin) end # End Function number 47 # Begin Function number 48 def cos(x) return(c(x).cos) end # End Function number 48 # Begin Function number 49 def tan(x) return(c(x).tan) end # End Function number 49 # Begin Function number 50 def arccos(x) return(c(x).acos) end # End Function number 50 # Begin Function number 51 def arccosh(x) return(c(x).acosh) end # End Function number 51 # Begin Function number 52 def arcsin(x) return(c(x).asin) end # End Function number 52 # Begin Function number 53 def arcsinh(x) return(c(x).asinh) end # End Function number 53 # Begin Function number 54 def arctan(x) return(c(x).atan) end # End Function number 54 # Begin Function number 55 def arctanh(x) return(c(x).atanh) end # End Function number 55 # Begin Function number 56 def cosh(x) return(c(x).cosh) end # End Function number 56 # Begin Function number 57 def erf(x) return(c(x).erf) end # End Function number 57 # Begin Function number 58 def log(x) return(c(x).log) end # End Function number 58 # Begin Function number 59 def ln(x) return(c(x).log) end # End Function number 59 # Begin Function number 60 def log10(x) return(c(x).log10) end # End Function number 60 # Begin Function number 61 def sinh(x) return(c(x).sinh) end # End Function number 61 # Begin Function number 62 def tanh(x) return(c(x).tanh) end # End Function number 62 # Begin Function number 63 def sqrt(x) return(c(x).sqrt) end # End Function number 63 # Begin Function number 64 def array2d(op3,op4) i = 0 x1 = Array.new(op3.to_i + 1) while i <= op3.to_i + 1 do # do number 4 x1[i] = Array.new(op4.to_i + 1) i += 1 end# end do number 5 return x1 end # End Function number 64 # Begin Function number 65 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 65 #END ATS LIBRARY BLOCK #BEGIN USER FUNCTION BLOCK #BEGIN BLOCK 3 #BEGIN USER DEF BLOCK def exact_soln_y (x) x = c(x) return($glob__m1 * cos(x) * c(2.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 = DESIRED_DIGITS end# end if 6 else relerr = $glob__m1 $glob_good_digits = -16 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(19)*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 = 32 end# end if 6 else relerr = $glob__m1 $glob_est_digits = -16 end# end if 5 $array_est_digits[1] = $glob_est_digits $array_good_digits[1] = $glob_good_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 2 tmp_rad = comp_rad_from_ratio($array_y_higher[1][last_no-1],$array_y_higher[1][last_no],last_no) if not prev_tmp_rad.maybe_zero 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 3 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 3 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 not prev_tmp_rad.maybe_zero 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 4 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 4 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 not prev_tmp_rad.maybe_zero 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 5 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 = $glob__1 while (term <= ATS_MAX_TERMS) do # do number 5 $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 6 $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 sin 1 $eq_no = 1 $array_tmp1[1] = sin($array_x[1]) $array_tmp1_g[1] = cos($array_x[1]) #emit pre mult CONST FULL $eq_no = 1 i = 1 $array_tmp2[1] = $array_const_2D0[1] * $array_tmp1[1] #emit pre add CONST FULL $eq_no = 1 i = 1 $array_tmp3[1] = $array_const_0D0[1] + $array_tmp2[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_tmp3[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(11) $array_y_higher[2][1] = c(temporary) end# end if 2 end# end if 1 kkk = 2 #END ATOMHDR1 #BEGIN ATOMHDR2 #emit pre sin ID_LINEAR iii = 2 $eq_no = 1 $array_tmp1[2] = $array_tmp1_g[1] * $array_x[2] / c(1) $array_tmp1_g[2] = $glob__m1*$array_tmp1[1] * $array_x[2] / c(1) #emit pre mult CONST FULL $eq_no = 1 i = 2 $array_tmp2[2] = $array_const_2D0[1] * $array_tmp1[2] #emit pre add CONST FULL $eq_no = 1 i = 2 $array_tmp3[2] = $array_tmp2[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_tmp3[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(21) $array_y_higher[2][2] = c(temporary) end# end if 2 end# end if 1 kkk = 3 #END ATOMHDR2 #BEGIN ATOMHDR3 #emit pre sin ID_LINEAR iii = 3 $eq_no = 1 $array_tmp1[3] = $array_tmp1_g[2] * $array_x[2] / c(2) $array_tmp1_g[3] = $glob__m1*$array_tmp1[2] * $array_x[2] / c(2) #emit pre mult CONST FULL $eq_no = 1 i = 3 $array_tmp2[3] = $array_const_2D0[1] * $array_tmp1[3] #emit pre add CONST FULL $eq_no = 1 i = 3 $array_tmp3[3] = $array_tmp2[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_tmp3[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(31) $array_y_higher[2][3] = c(temporary) end# end if 2 end# end if 1 kkk = 4 #END ATOMHDR3 #BEGIN ATOMHDR4 #emit pre sin ID_LINEAR iii = 4 $eq_no = 1 $array_tmp1[4] = $array_tmp1_g[3] * $array_x[2] / c(3) $array_tmp1_g[4] = $glob__m1*$array_tmp1[3] * $array_x[2] / c(3) #emit pre mult CONST FULL $eq_no = 1 i = 4 $array_tmp2[4] = $array_const_2D0[1] * $array_tmp1[4] #emit pre add CONST FULL $eq_no = 1 i = 4 $array_tmp3[4] = $array_tmp2[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_tmp3[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(41) $array_y_higher[2][4] = c(temporary) end# end if 2 end# end if 1 kkk = 5 #END ATOMHDR4 #BEGIN ATOMHDR5 #emit pre sin ID_LINEAR iii = 5 $eq_no = 1 $array_tmp1[5] = $array_tmp1_g[4] * $array_x[2] / c(4) $array_tmp1_g[5] = $glob__m1*$array_tmp1[4] * $array_x[2] / c(4) #emit pre mult CONST FULL $eq_no = 1 i = 5 $array_tmp2[5] = $array_const_2D0[1] * $array_tmp1[5] #emit pre add CONST FULL $eq_no = 1 i = 5 $array_tmp3[5] = $array_tmp2[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_tmp3[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(51) $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 62 #END OUTFILE3 #BEGIN OUTFILE4 #emit sin LINEAR $eq_no = 1 $array_tmp1[kkk] = $array_tmp1_g[kkk - 1] * $array_x[2] / c(kkk - 1) $array_tmp1_g[kkk] = $glob__m1*$array_tmp1[kkk - 1] * $array_x[2] / c(kkk - 1) #emit mult CONST FULL $eq_no = 1 i = 1 $array_tmp2[kkk] = $array_const_2D0[1] * $array_tmp1[kkk] #emit NOT FULL - FULL add $eq_no = 1 $array_tmp3[kkk] = $array_tmp2[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_tmp3[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 2 end# end if 2 end# end if 1 kkk = kkk + 1 end# end do number 63 #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 require_relative "../rapfp/ApConfig" require_relative "../rapfp/Ap" rconst = $ApConst.new rconst.one.setrconst(rconst) rconst.init2 rconst.one.setrconst(rconst) $RConst = rconst cconst = ApcConst.new(rconst) cconst.one.setcconst(cconst) cconst.one.setrconst(rconst) $CConst = cconst # 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_good_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_g= 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_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=ap_in("0.1e-20+/-1.0e-100") $glob_small_float=ap_in("0.1e-20+/-1.0e-100") $glob_smallish_float=ap_in("0.1e-24+/-1.0e-100") $glob_large_float=ap_in("1.0e100+/-1.0e90") $glob_larger_float=ap_in("1.1e100+/-1.1e90") $glob__m2=c0(-2) $glob__m1=c0(-1) $glob__0=c0(0) $glob__1=c0(1) $glob__2=c0(2) $glob__3=c0(3) $glob__4=c0(4) $glob__5=c0(5) $glob__8=c0(8) $glob__10=c0(10) $glob__100=c0(100) $glob__0_5=c(0.5) $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_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_2D0= Array.new(ATS_MAX_TERMS + 1) # before arrays initialized term = 1 while (term <= ATS_MAX_TERMS) do # do number 6 $array_y_init[term] = c(0.0) term = term + 1 end# end do number 7 term = 1 while (term <= ATS_MAX_TERMS) do # do number 7 $array_norms[term] = c(0.0) term = term + 1 end# end do number 8 term = 1 while (term <= ATS_MAX_TERMS) do # do number 8 $array_fact_1[term] = c(0.0) term = term + 1 end# end do number 9 term = 1 while (term <= 2) do # do number 9 $array_1st_rel_error[term] = c(0.0) term = term + 1 end# end do number 10 term = 1 while (term <= 2) do # do number 10 $array_last_rel_error[term] = c(0.0) term = term + 1 end# end do number 11 term = 1 while (term <= 2) do # do number 11 $array_est_rel_error[term] = c(0.0) term = term + 1 end# end do number 12 term = 1 while (term <= 2) do # do number 12 $array_max_est_error[term] = c(0.0) term = term + 1 end# end do number 13 term = 1 while (term <= 2) do # do number 13 $array_type_pole[term] = 0 term = term + 1 end# end do number 14 term = 1 while (term <= 2) do # do number 14 $array_type_real_pole[term] = 0 term = term + 1 end# end do number 15 term = 1 while (term <= 2) do # do number 15 $array_type_complex_pole[term] = 0 term = term + 1 end# end do number 16 term = 1 while (term <= 2) do # do number 16 $array_est_digits[term] = 0 term = term + 1 end# end do number 17 term = 1 while (term <= 2) do # do number 17 $array_good_digits[term] = 0 term = term + 1 end# end do number 18 term = 1 while (term <= ATS_MAX_TERMS) do # do number 18 $array_y[term] = c(0.0) term = term + 1 end# end do number 19 term = 1 while (term <= ATS_MAX_TERMS) do # do number 19 $array_x[term] = c(0.0) term = term + 1 end# end do number 20 term = 1 while (term <= ATS_MAX_TERMS) do # do number 20 $array_tmp0[term] = c(0.0) term = term + 1 end# end do number 21 term = 1 while (term <= ATS_MAX_TERMS) do # do number 21 $array_tmp1_g[term] = c(0.0) term = term + 1 end# end do number 22 term = 1 while (term <= ATS_MAX_TERMS) do # do number 22 $array_tmp1[term] = c(0.0) term = term + 1 end# end do number 23 term = 1 while (term <= ATS_MAX_TERMS) do # do number 23 $array_tmp2[term] = c(0.0) term = term + 1 end# end do number 24 term = 1 while (term <= ATS_MAX_TERMS) do # do number 24 $array_tmp3[term] = c(0.0) term = term + 1 end# end do number 25 term = 1 while (term <= ATS_MAX_TERMS) do # do number 25 $array_m1[term] = c(0.0) term = term + 1 end# end do number 26 ord = 1 while (ord <=2) do # do number 26 term = 1 while (term <= ATS_MAX_TERMS) do # do number 27 $array_y_higher[ord][term] = c(0.0) term = term + 1 end# end do number 28 ord = ord + 1 end# end do number 28 ord = 1 while (ord <=2) do # do number 28 term = 1 while (term <= ATS_MAX_TERMS) do # do number 29 $array_y_higher_work[ord][term] = c(0.0) term = term + 1 end# end do number 30 ord = ord + 1 end# end do number 30 ord = 1 while (ord <=2) do # do number 30 term = 1 while (term <= ATS_MAX_TERMS) do # do number 31 $array_y_higher_work2[ord][term] = c(0.0) term = term + 1 end# end do number 32 ord = ord + 1 end# end do number 32 ord = 1 while (ord <=2) do # do number 32 term = 1 while (term <= ATS_MAX_TERMS) do # do number 33 $array_y_set_initial[ord][term] = c(0.0) term = term + 1 end# end do number 34 ord = ord + 1 end# end do number 34 ord = 1 while (ord <=2) do # do number 34 term = 1 while (term <= 3) do # do number 35 $array_given_rad_poles[ord][term] = c(0.0) term = term + 1 end# end do number 36 ord = ord + 1 end# end do number 36 ord = 1 while (ord <=2) do # do number 36 term = 1 while (term <= 3) do # do number 37 $array_given_ord_poles[ord][term] = c(0.0) term = term + 1 end# end do number 38 ord = ord + 1 end# end do number 38 ord = 1 while (ord <=2) do # do number 38 term = 1 while (term <= 4) do # do number 39 $array_rad_test_poles[ord][term] = c(0.0) term = term + 1 end# end do number 40 ord = ord + 1 end# end do number 40 ord = 1 while (ord <=2) do # do number 40 term = 1 while (term <= 4) do # do number 41 $array_ord_test_poles[ord][term] = c(0.0) term = term + 1 end# end do number 42 ord = ord + 1 end# end do number 42 ord = 1 while (ord <=ATS_MAX_TERMS) do # do number 42 term = 1 while (term <= ATS_MAX_TERMS) do # do number 43 $array_fact_2[ord][term] = c(0.0) term = term + 1 end# end do number 44 ord = ord + 1 end# end do number 44 # 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_g) zero_ats_ar($array_tmp1) zero_ats_ar($array_tmp2) zero_ats_ar($array_tmp3) 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_2D0) $array_const_2D0[1] = c(2.0) 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 44 jjjf = 0 while (jjjf <= ATS_MAX_TERMS) do # do number 45 $array_fact_1[iiif] = 0 $array_fact_2[iiif][jjjf] = 0 jjjf = jjjf + 1 end# end do number 46 iiif = iiif + 1 end# end do number 46 #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 = ap_in("0.1e-20+/-1.0e-100") $glob_small_float = ap_in("0.1e-20+/-1.0e-100") $glob_smallish_float = ap_in("0.1e-24+/-1.0e-100") $glob_large_float = ap_in("1.0e100+/-1.0e90") $glob_larger_float = ap_in("1.1e100+/-1.1e90") $glob__m2 = c0(-2) $glob__m1 = c0(-1) $glob__0 = c0(0) $glob__1 = c0(1) $glob__2 = c0(2) $glob__3 = c0(3) $glob__4 = c0(4) $glob__5 = c0(5) $glob__8 = c0(8) $glob__10 = c0(10) $glob__100 = c0(100) $glob__0_5 = c(0.5) $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_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/mult_c_sinpostode.ode#################") omniout_str(ALWAYS,"diff ( y , x , 1 ) = 2.0 * sin ( x ) ; ") 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(0.1) ") omniout_str(ALWAYS,"x_end=c(5.0) ") omniout_str(ALWAYS,"$array_y_init[0 + 1] = exact_soln_y(x_start)") omniout_str(ALWAYS,"$glob_look_poles=true ") 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=3 ") 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=(20.0) ") omniout_str(ALWAYS,"$glob_subiter_method=3 ") omniout_str(ALWAYS,"$glob_max_iter=1000 ") omniout_str(ALWAYS,"$glob_upper_ratio_limit=c(1.11) ") omniout_str(ALWAYS,"$glob_lower_ratio_limit=c(0.99) ") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") omniout_str(ALWAYS,"# ELIMINATED in preodein.rb") 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,"return($glob__m1 * cos(x) * c(2.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(0.1) x_end=c(5.0) $array_y_init[0 + 1] = exact_soln_y(x_start) $glob_look_poles=true # 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=3 # ELIMINATED in preodein.rb #END SECOND INPUT BLOCK #BEGIN OVERRIDE BLOCK $glob_desired_digits_correct=8 $glob_max_minutes=(20.0) $glob_subiter_method=3 $glob_max_iter=1000 $glob_upper_ratio_limit=c(1.11) $glob_lower_ratio_limit=c(0.99) # ELIMINATED in preodein.rb # ELIMINATED in preodein.rb #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)) 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 $glob_prec = expt($glob__10,c(-NUM_DIGITS)) 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 46 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 47 $array_y[term_no] = c($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 48 rows = order_diff r_order = 1 while (r_order <= rows) do # do number 48 term_no = 1 while (term_no <= (rows - r_order + 1)) do # do number 49 it = term_no + r_order - 1 if (term_no < ATS_MAX_TERMS) then # if number 10 $array_y_higher[r_order][term_no] = c($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 50 r_order = r_order + 1 end# end do number 50 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 check_for_pole() 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 50 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 end# end if 9 found_h = true $glob_h = ap_in("0.05+/-0.1e-100") #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 50 $array_y[term_no] = c($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 51 rows = order_diff r_order = 1 while (r_order <= rows) do # do number 51 term_no = 1 while (term_no <= (rows - r_order + 1)) do # do number 52 it = term_no + r_order - 1 if (term_no < ATS_MAX_TERMS) then # if number 10 $array_y_higher[r_order][term_no] = c($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 53 r_order = r_order + 1 end# end do number 53 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 53 #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) check_for_pole() 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 54 $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 55 #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 55 temp_sum = temp_sum + $array_y_higher_work[ord][iii] iii = iii - 1 end# end do number 56 $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 56 $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 57 #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 57 temp_sum = temp_sum + $array_y_higher_work[ord][iii] iii = iii - 1 end# end do number 58 $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 58 $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 59 #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 59 temp_sum = temp_sum + $array_y_higher_work[ord][iii] iii = iii - 1 end# end do number 60 $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 60 $array_y[term_no] = $array_y_higher_work2[1][term_no] ord = 1 while (ord <= order_diff) do # do number 61 $array_y_higher[ord][term_no] = $array_y_higher_work2[ord][term_no] ord = ord + 1 end# end do number 62 term_no = term_no - 1 end# end do number 62 #END PART 2 HEVE MOVED TERMS to REGULAR Array end# end do number 62#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 ) = 2.0 * sin ( x ) ; ") 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-03-25T04:04:01-05:00") logitem_str(html_log_file,"Ruby Apfp") logitem_str(html_log_file,"mult_c_sin") logitem_str(html_log_file,"diff ( y , x , 1 ) = 2.0 * sin ( x ) ; ") 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,32) 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 ($array_good_digits[1] != -16) then # if number 11 logitem_integer(html_log_file,$array_good_digits[1]) else logitem_str(html_log_file,"Unknown") end# end if 11 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," 277 ") logitem_str(html_log_file,"mult_c_sin diffeq.rapfp.rb") logitem_str(html_log_file,"mult_c_sin Ruby Apfp results") logitem_str(html_log_file,"OK") 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()