#! /usr/bin/env python # def r8_ci ( x ): #*****************************************************************************80 # ## R8_CI evaluates the cosine integral Ci of an R8 argument. # # Discussion: # # The cosine integral is defined by # # CI(X) = - integral ( X <= T < Infinity ) ( cos ( T ) ) / T dT # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # Original FORTRAN77 version by Wayne Fullerton. # Python version by John Burkardt. # # Reference: # # Wayne Fullerton, # Portable Special Function Routines, # in Portability of Numerical Software, # edited by Wayne Cowell, # Lecture Notes in Computer Science, Volume 57, # Springer 1977, # ISBN: 978-3-540-08446-4, # LC: QA297.W65. # # Parameters: # # Input, real X, the argument. # # Output, real VALUE, the cosine integral Ci evaluated at X. # import numpy as np from r8_csevl import r8_csevl from r8_inits import r8_inits from machine import r8_mach from sys import exit cics = np.array ( [ \ -0.34004281856055363156281076633129873, \ -1.03302166401177456807159271040163751, \ 0.19388222659917082876715874606081709, \ -0.01918260436019865893946346270175301, \ 0.00110789252584784967184098099266118, \ -0.00004157234558247208803840231814601, \ 0.00000109278524300228715295578966285, \ -0.00000002123285954183465219601280329, \ 0.00000000031733482164348544865129873, \ -0.00000000000376141547987683699381798, \ 0.00000000000003622653488483964336956, \ -0.00000000000000028911528493651852433, \ 0.00000000000000000194327860676494420, \ -0.00000000000000000001115183182650184, \ 0.00000000000000000000005527858887706, \ -0.00000000000000000000000023907013943, \ 0.00000000000000000000000000091001612, \ -0.00000000000000000000000000000307233, \ 0.00000000000000000000000000000000926 ] ) nci = r8_inits ( cics, 19, 0.1 * r8_mach ( 3 ) ) xsml = np.sqrt ( r8_mach ( 3 ) ) if ( x <= 0.0 ): print ( '' ) print ( 'R8_CI - Fatal error!' ) print ( ' X <= 0.0.' ) exit ( 'R8_CI - Fatal error!' ) elif ( x <= xsml ): y = - 1.0 value = np.log ( x ) - 0.5 + r8_csevl ( y, cics, nci ) elif ( x <= 4.0 ): y = ( x * x - 8.0 ) * 0.125 value = np.log ( x ) - 0.5 + r8_csevl ( y, cics, nci ) else: f, g = r8_sifg ( x ) sinx = np.sin ( x ) value = f * sinx - g * np.cos ( x ) return value def r8_ci_test ( ): #*****************************************************************************80 # ## R8_CI_TEST tests R8_CI. # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # John Burkardt # import platform from ci_values import ci_values print ( '' ) print ( 'R8_CI_TEST:' ) print ( ' Python version: %s' % ( platform.python_version ( ) ) ) print ( ' R8_CI evaluates the cosine integral.' ) print ( '' ) print ( ' X CI(X) R8_CI(X) Diff' ) print ( '' ) n_data = 0 while ( True ): n_data, x, fx1 = ci_values ( n_data ) if ( n_data == 0 ): break fx2 = r8_ci ( x ) print ( ' %14.4f %14.6g %14.6g %14.6g' % ( x, fx1, fx2, abs ( fx1 - fx2 ) ) ) # # Terminate. # print ( '' ) print ( 'R8_CI_TEST:' ) print ( ' Normal end of execution.' ) return def r8_cin ( x ): #*****************************************************************************80 # ## R8_CIN evaluates the alternate cosine integral Cin of an R8 argument. # # Discussion: # # CIN(X) = gamma + log(X) # + integral ( 0 <= T <= X ) ( cos ( T ) - 1 ) / T dT # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # Original FORTRAN77 version by Wayne Fullerton. # Python version by John Burkardt. # # Reference: # # Wayne Fullerton, # Portable Special Function Routines, # in Portability of Numerical Software, # edited by Wayne Cowell, # Lecture Notes in Computer Science, Volume 57, # Springer 1977, # ISBN: 978-3-540-08446-4, # LC: QA297.W65. # # Parameters: # # Input, real X, the argument. # # Output, real VALUE, the cosine integral Cin evaluated at X. # import numpy as np from r8_csevl import r8_csevl from r8_inits import r8_inits from machine import r8_mach eul = 0.57721566490153286060651209008240 cincs = np.array ( [ \ 0.37074501750909688741654801228564992, \ -0.05893574896364446831956864397363697, \ 0.00538189642113569124048745326203340, \ -0.00029860052841962135319594906563410, \ 0.00001095572575321620077031054467306, \ -0.00000028405454877346630491727187731, \ 0.00000000546973994875384912457861806, \ -0.00000000008124187461318157083277452, \ 0.00000000000095868593117706609013181, \ -0.00000000000000920266004392351031377, \ 0.00000000000000007325887999017895024, \ -0.00000000000000000049143726675842909, \ 0.00000000000000000000281577746753902, \ -0.00000000000000000000001393986788501, \ 0.00000000000000000000000006022485646, \ -0.00000000000000000000000000022904717, \ 0.00000000000000000000000000000077273, \ -0.00000000000000000000000000000000233 ] ) ncin = r8_inits ( cincs, 18, 0.1 * r8_mach ( 3 ) ) xmin = np.sqrt ( r8_mach ( 1 ) ) absx = abs ( x ) if ( absx <= xmin ): value = 0.0 elif ( absx <= 4.0 ): value = r8_csevl ( ( x * x - 8.0 ) * 0.125, cincs, ncin ) * x * x else: f, g = r8_sifg ( absx ) sinx = np.sin ( absx ) value = - f * sinx + g * np.cos ( absx ) + np.log ( absx ) + eul return value def r8_cin_test ( ): #*****************************************************************************80 # ## R8_CIN_TEST tests R8_CIN. # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # John Burkardt # import platform from cin_values import cin_values print ( '' ) print ( 'R8_CIN_TEST:' ) print ( ' Python version: %s' % ( platform.python_version ( ) ) ) print ( ' R8_CIN evaluates the alternate hyperbolic cosine integral.' ) print ( '' ) print ( ' X CIN(X) R8_CIN(X) Diff' ) print ( '' ) n_data = 0 while ( True ): n_data, x, fx1 = cin_values ( n_data ) if ( n_data == 0 ): break fx2 = r8_cin ( x ) print ( ' %14.4f %14.6g %14.6g %14.6g' % ( x, fx1, fx2, abs ( fx1 - fx2 ) ) ) # # Terminate. # print ( '' ) print ( 'R8_CIN_TEST:' ) print ( ' Normal end of execution.' ) return def r8_si ( x ): #*****************************************************************************80 # ## R8_SI evaluates the sine integral Si of an R8 argument. # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # Original FORTRAN77 version by Wayne Fullerton. # Python version by John Burkardt. # # Reference: # # Wayne Fullerton, # Portable Special Function Routines, # in Portability of Numerical Software, # edited by Wayne Cowell, # Lecture Notes in Computer Science, Volume 57, # Springer 1977, # ISBN: 978-3-540-08446-4, # LC: QA297.W65. # # Parameters: # # Input, real X, the argument. # # Output, real VALUE, the sine integral Si evaluated at X. # import numpy as np from r8_csevl import r8_csevl from r8_inits import r8_inits from machine import r8_mach sics = np.array ( [ \ -0.1315646598184841928904275173000457, \ -0.2776578526973601892048287660157299, \ 0.0354414054866659179749135464710086, \ -0.0025631631447933977658752788361530, \ 0.0001162365390497009281264921482985, \ -0.0000035904327241606042670004347148, \ 0.0000000802342123705710162308652976, \ -0.0000000013562997692540250649931846, \ 0.0000000000179440721599736775567759, \ -0.0000000000001908387343087145490737, \ 0.0000000000000016669989586824330853, \ -0.0000000000000000121730988368503042, \ 0.0000000000000000000754181866993865, \ -0.0000000000000000000004014178842446, \ 0.0000000000000000000000018553690716, \ -0.0000000000000000000000000075166966, \ 0.0000000000000000000000000000269113, \ -0.0000000000000000000000000000000858 ] ) nsi = r8_inits ( sics, 18, 0.1 * r8_mach ( 3 ) ) xsml = np.sqrt ( r8_mach ( 3 ) ) absx = abs ( x ) if ( absx < xsml ): value = x elif ( absx <= 4.0 ): value = x * ( 0.75 + r8_csevl ( ( x * x - 8.0 ) * 0.125, sics, nsi ) ) else: f, g = r8_sifg ( absx ) cosx = np.cos ( absx ) value = 0.5 * np.pi - f * cosx - g * np.sin ( x ) if ( x < 0.0 ): value = - value return value def r8_si_test ( ): #*****************************************************************************80 # ## R8_SI_TEST tests R8_SI. # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # John Burkardt # import platform from si_values import si_values print ( '' ) print ( 'R8_SI_TEST:' ) print ( ' Python version: %s' % ( platform.python_version ( ) ) ) print ( ' R8_SI evaluates the sine integral.' ) print ( '' ) print ( ' X SI(X) R8_SI(X) Diff' ) print ( '' ) n_data = 0 while ( True ): n_data, x, fx1 = si_values ( n_data ) if ( n_data == 0 ): break fx2 = r8_si ( x ) print ( ' %14.4g %14.6g %14.6g %14.6g' % ( x, fx1, fx2, abs ( fx1 - fx2 ) ) ) # # Terminate. # print ( '' ) print ( 'R8_SI_TEST:' ) print ( ' Normal end of execution.' ) return def r8_sifg ( x ): #*****************************************************************************80 # ## R8_SIFG is a utility routine. # # Licensing: # # This code is distributed under the GNU LGPL license. # # Modified: # # 24 April 2016 # # Author: # # Original FORTRAN77 version by Wayne Fullerton. # Python version by John Burkardt. # # Reference: # # Wayne Fullerton, # Portable Special Function Routines, # in Portability of Numerical Software, # edited by Wayne Cowell, # Lecture Notes in Computer Science, Volume 57, # Springer 1977, # ISBN: 978-3-540-08446-4, # LC: QA297.W65. # # Parameters: # # Input, real X, the argument. # # Output, real F, G. # import numpy as np from r8_csevl import r8_csevl from r8_inits import r8_inits from machine import r8_mach from sys import exit f1cs = np.array ( [ \ -0.1191081969051363610348201965828918, \ -0.0247823144996236247590074150823133, \ 0.0011910281453357821268120363054457, \ -0.0000927027714388561748308600360706, \ 0.0000093373141568270996868204582766, \ -0.0000011058287820557143938979426306, \ 0.0000001464772071460162169336550799, \ -0.0000000210694496287689532601227548, \ 0.0000000032293492366848236382857374, \ -0.0000000005206529617529375828014986, \ 0.0000000000874878884570278750268316, \ -0.0000000000152176187056123668294574, \ 0.0000000000027257192405419573900583, \ -0.0000000000005007053075968556290255, \ 0.0000000000000940240902726068511779, \ -0.0000000000000180014444791803678336, \ 0.0000000000000035062621432741785826, \ -0.0000000000000006935282926769149709, \ 0.0000000000000001390925136454216568, \ -0.0000000000000000282486885074170585, \ 0.0000000000000000058031305693579081, \ -0.0000000000000000012046901573375820, \ 0.0000000000000000002525052443655940, \ -0.0000000000000000000533980268805594, \ 0.0000000000000000000113855786274122, \ -0.0000000000000000000024462861505259, \ 0.0000000000000000000005293659320439, \ -0.0000000000000000000001153184940277, \ 0.0000000000000000000000252786568318, \ -0.0000000000000000000000055738645378, \ 0.0000000000000000000000012358245621, \ -0.0000000000000000000000002754350842, \ 0.0000000000000000000000000616906808, \ -0.0000000000000000000000000138817443, \ 0.0000000000000000000000000031375329, \ -0.0000000000000000000000000007121249, \ 0.0000000000000000000000000001622778, \ -0.0000000000000000000000000000371206, \ 0.0000000000000000000000000000085221, \ -0.0000000000000000000000000000019633, \ 0.0000000000000000000000000000004538, \ -0.0000000000000000000000000000001052, \ 0.0000000000000000000000000000000245 ] ) f2cs = np.array ( [ \ -0.03484092538970132330836049733745577, \ -0.01668422056779596873246786312278676, \ 0.00067529012412377385045207859239727, \ -0.00005350666225447013628785577557429, \ 0.00000626934217790075267050759431626, \ -0.00000095266388019916680677790414293, \ 0.00000017456292242509880425504427666, \ -0.00000003687954030653093307097646628, \ 0.00000000872026777051395264075816938, \ -0.00000000226019703919738748530423167, \ 0.00000000063246249765250612520444877, \ -0.00000000018889118884717869240911480, \ 0.00000000005967746729997813372620472, \ -0.00000000001980443117372239011196007, \ 0.00000000000686413954772103383713264, \ -0.00000000000247310193070199106074890, \ 0.00000000000092263594549941404196042, \ -0.00000000000035523634999261784497297, \ 0.00000000000014076049625351591461820, \ -0.00000000000005726228499747652794311, \ 0.00000000000002386537545413171810106, \ -0.00000000000001017141890764597142232, \ 0.00000000000000442594531078364424968, \ -0.00000000000000196344933049189761979, \ 0.00000000000000088688748314810461024, \ -0.00000000000000040743345027311546948, \ 0.00000000000000019016837215675339859, \ -0.00000000000000009009707297478042442, \ 0.00000000000000004329211274095668667, \ -0.00000000000000002108144465322479526, \ 0.00000000000000001039637907026452274, \ -0.00000000000000000518891007948931936, \ 0.00000000000000000261955324869899371, \ -0.00000000000000000133690399951301570, \ 0.00000000000000000068941057702931664, \ -0.00000000000000000035905362610437250, \ 0.00000000000000000018878077255791706, \ -0.00000000000000000010016125265594380, \ 0.00000000000000000005360725691578228, \ -0.00000000000000000002893198974944827, \ 0.00000000000000000001574065100202625, \ -0.00000000000000000000863027106431206, \ 0.00000000000000000000476715602862288, \ -0.00000000000000000000265222739998504, \ 0.00000000000000000000148582865063866, \ -0.00000000000000000000083797235923135, \ 0.00000000000000000000047565916422711, \ -0.00000000000000000000027169073353112, \ 0.00000000000000000000015612738881686, \ -0.00000000000000000000009024555078347, \ 0.00000000000000000000005246097049119, \ -0.00000000000000000000003066450818697, \ 0.00000000000000000000001801996250957, \ -0.00000000000000000000001064443050752, \ 0.00000000000000000000000631942158881, \ -0.00000000000000000000000377013812246, \ 0.00000000000000000000000225997542918, \ -0.00000000000000000000000136100844814, \ 0.00000000000000000000000082333232003, \ -0.00000000000000000000000050025986091, \ 0.00000000000000000000000030526245684, \ -0.00000000000000000000000018705164021, \ 0.00000000000000000000000011508404393, \ -0.00000000000000000000000007108714611, \ 0.00000000000000000000000004408065533, \ -0.00000000000000000000000002743760867, \ 0.00000000000000000000000001714144851, \ -0.00000000000000000000000001074768860, \ 0.00000000000000000000000000676259777, \ -0.00000000000000000000000000426981348, \ 0.00000000000000000000000000270500637, \ -0.00000000000000000000000000171933331, \ 0.00000000000000000000000000109636138, \ -0.00000000000000000000000000070132573, \ 0.00000000000000000000000000045001784, \ -0.00000000000000000000000000028963835, \ 0.00000000000000000000000000018697009, \ -0.00000000000000000000000000012104646, \ 0.00000000000000000000000000007859065, \ -0.00000000000000000000000000005116867, \ 0.00000000000000000000000000003340627, \ -0.00000000000000000000000000002186851, \ 0.00000000000000000000000000001435340, \ -0.00000000000000000000000000000944523, \ 0.00000000000000000000000000000623117, \ -0.00000000000000000000000000000412101, \ 0.00000000000000000000000000000273208, \ -0.00000000000000000000000000000181558, \ 0.00000000000000000000000000000120934, \ -0.00000000000000000000000000000080737, \ 0.00000000000000000000000000000054022, \ -0.00000000000000000000000000000036227, \ 0.00000000000000000000000000000024348, \ -0.00000000000000000000000000000016401, \ 0.00000000000000000000000000000011074, \ -0.00000000000000000000000000000007497, \ 0.00000000000000000000000000000005091, \ -0.00000000000000000000000000000003470, \ 0.00000000000000000000000000000002377 ] ) g1cs = np.array ( [ \ -0.3040578798253495954499726682091083, \ -0.0566890984597120587731339156118269, \ 0.0039046158173275643919984071554082, \ -0.0003746075959202260618619339867489, \ 0.0000435431556559843679552220840065, \ -0.0000057417294453025046561970723475, \ 0.0000008282552104502629741937616492, \ -0.0000001278245892594642727883913223, \ 0.0000000207978352948687884439257529, \ -0.0000000035313205921990798042032682, \ 0.0000000006210824236308951068631449, \ -0.0000000001125215474446292649336987, \ 0.0000000000209088917684421605267019, \ -0.0000000000039715831737681727689158, \ 0.0000000000007690431314272089939005, \ -0.0000000000001514696742731613519826, \ 0.0000000000000302892146552359684119, \ -0.0000000000000061399703834708825400, \ 0.0000000000000012600605829510933553, \ -0.0000000000000002615029250939483683, \ 0.0000000000000000548278844891796821, \ -0.0000000000000000116038182129526571, \ 0.0000000000000000024771654107129795, \ -0.0000000000000000005330672753223389, \ 0.0000000000000000001155666075598465, \ -0.0000000000000000000252280547744957, \ 0.0000000000000000000055429038550786, \ -0.0000000000000000000012252208421297, \ 0.0000000000000000000002723664318684, \ -0.0000000000000000000000608707831422, \ 0.0000000000000000000000136724874476, \ -0.0000000000000000000000030856626806, \ 0.0000000000000000000000006995212319, \ -0.0000000000000000000000001592587569, \ 0.0000000000000000000000000364051056, \ -0.0000000000000000000000000083539465, \ 0.0000000000000000000000000019240303, \ -0.0000000000000000000000000004446816, \ 0.0000000000000000000000000001031182, \ -0.0000000000000000000000000000239887, \ 0.0000000000000000000000000000055976, \ -0.0000000000000000000000000000013100, \ 0.0000000000000000000000000000003074, \ -0.0000000000000000000000000000000723 ] ) g2cs = np.array ( [ \ -0.1211802894731646263541834046858267, \ -0.0316761386394950286701407923505610, \ 0.0013383199778862680163819429492182, \ -0.0000895511011392252425531905069518, \ 0.0000079155562961718213115249467924, \ -0.0000008438793322241520181418982080, \ 0.0000001029980425677530146647227274, \ -0.0000000139295750605183835795834444, \ 0.0000000020422703959875980400677594, \ -0.0000000003196534694206427035434752, \ 0.0000000000528147832657267698615312, \ -0.0000000000091339554672671033735289, \ 0.0000000000016426251238967760444819, \ -0.0000000000003055897039322660002410, \ 0.0000000000000585655825785779717892, \ -0.0000000000000115229197730940120563, \ 0.0000000000000023209469119988537310, \ -0.0000000000000004774355834177535025, \ 0.0000000000000001000996765800180573, \ -0.0000000000000000213533778082256704, \ 0.0000000000000000046277190777367671, \ -0.0000000000000000010175807410227657, \ 0.0000000000000000002267657399884672, \ -0.0000000000000000000511630776076426, \ 0.0000000000000000000116767014913108, \ -0.0000000000000000000026935427672470, \ 0.0000000000000000000006275665841146, \ -0.0000000000000000000001475880557531, \ 0.0000000000000000000000350145314739, \ -0.0000000000000000000000083757732152, \ 0.0000000000000000000000020191815152, \ -0.0000000000000000000000004903567705, \ 0.0000000000000000000000001199123348, \ -0.0000000000000000000000000295170610, \ 0.0000000000000000000000000073113112, \ -0.0000000000000000000000000018217843, \ 0.0000000000000000000000000004565148, \ -0.0000000000000000000000000001150151, \ 0.0000000000000000000000000000291267, \ -0.0000000000000000000000000000074125, \ 0.0000000000000000000000000000018953, \ -0.0000000000000000000000000000004868, \ 0.0000000000000000000000000000001256, \ -0.0000000000000000000000000000000325 ] ) g3cs = np.array ( [ \ -0.0280574367809472928402815264335299, \ -0.0137271597162236975409100508089556, \ 0.0002894032638760296027448941273751, \ -0.0000114129239391197145908743622517, \ 0.0000006813965590726242997720207302, \ -0.0000000547952289604652363669058052, \ 0.0000000055207429918212529109406521, \ -0.0000000006641464199322920022491428, \ 0.0000000000922373663487041108564960, \ -0.0000000000144299088886682862611718, \ 0.0000000000024963904892030710248705, \ -0.0000000000004708240675875244722971, \ 0.0000000000000957217659216759988140, \ -0.0000000000000207889966095809030537, \ 0.0000000000000047875099970877431627, \ -0.0000000000000011619070583377173759, \ 0.0000000000000002956508969267836974, \ -0.0000000000000000785294988256492025, \ 0.0000000000000000216922264368256612, \ -0.0000000000000000062113515831676342, \ 0.0000000000000000018384568838450977, \ -0.0000000000000000005610887482137276, \ 0.0000000000000000001761862805280062, \ -0.0000000000000000000568111050541451, \ 0.0000000000000000000187786279582313, \ -0.0000000000000000000063531694151124, \ 0.0000000000000000000021968802368238, \ -0.0000000000000000000007754666550395, \ 0.0000000000000000000002791018356581, \ -0.0000000000000000000001023178525247, \ 0.0000000000000000000000381693403919, \ -0.0000000000000000000000144767895606, \ 0.0000000000000000000000055779512634, \ -0.0000000000000000000000021817239071, \ 0.0000000000000000000000008656646309, \ -0.0000000000000000000000003482157895, \ 0.0000000000000000000000001419188130, \ -0.0000000000000000000000000585714314, \ 0.0000000000000000000000000244660482, \ -0.0000000000000000000000000103387099, \ 0.0000000000000000000000000044177299, \ -0.0000000000000000000000000019080079, \ 0.0000000000000000000000000008326038, \ -0.0000000000000000000000000003669553, \ 0.0000000000000000000000000001632875, \ -0.0000000000000000000000000000733357, \ 0.0000000000000000000000000000332327, \ -0.0000000000000000000000000000151906, \ 0.0000000000000000000000000000070020, \ -0.0000000000000000000000000000032539, \ 0.0000000000000000000000000000015240, \ -0.0000000000000000000000000000007193, \ 0.0000000000000000000000000000003420, \ -0.0000000000000000000000000000001638, \ 0.0000000000000000000000000000000790, \ -0.0000000000000000000000000000000383 ] ) tol = 0.1 * r8_mach ( 3 ) nf1 = r8_inits ( f1cs, 43, tol ) nf2 = r8_inits ( f2cs, 99, tol ) ng1 = r8_inits ( g1cs, 44, tol ) ng2 = r8_inits ( g2cs, 44, tol ) ng3 = r8_inits ( g3cs, 56, tol ) xbig = np.sqrt ( 1.0 / r8_mach ( 3 ) ) xmaxf = np.exp ( min ( - np.log ( r8_mach ( 1 ) ), \ np.log ( r8_mach ( 2 ) ) ) - 0.01 ) xmaxg = 1.0 / np.sqrt ( r8_mach ( 1 ) ) xbnd = np.sqrt ( 50.0 ) xbndg = np.sqrt ( 200.0 ) if ( x < 4.0 ): print ( '' ) print ( 'R8_SIFG - Fatal error!' ) print ( ' Approximation invalid for X < 4.' ) exit ( 'R8_SIFG - Fatal error!' ) elif ( x <= xbnd ): f = ( 1.0 + r8_csevl ( ( 1.0 / x / x - 0.04125 ) / 0.02125, f1cs, nf1 ) ) / x g = ( 1.0 + r8_csevl ( ( 1.0 / x / x - 0.04125 ) / 0.02125, g1cs, ng1 ) ) / x / x elif ( x <= xbig ): f = ( 1.0 + r8_csevl ( 100. / x / x - 1.0, f2cs, nf2 ) ) / x if ( x <= xbndg ): g = ( 1.0 + r8_csevl ( ( 10000.0 / x / x - 125.0 ) / 75.0, g2cs, ng2 ) ) / x / x else: g = ( 1.0 + r8_csevl ( 400.0 / x / x - 1.0, g3cs, ng3 ) ) / x / x else: if ( x < xmaxf ): f = 1.0 / x else: f = 0.0 if ( x < xmaxg ): g = 1.0 / x / x else: g = 0.0 return f, g if ( __name__ == '__main__' ): from timestamp import timestamp timestamp ( ) r8_ci_test ( ) r8_cin_test ( ) r8_si_test ( ) timestamp ( )