Révision 272 pobysoPythonSage/src/pobyso.py
| pobyso.py (revision 272) | ||
|---|---|---|
| 521 | 521 |
precSo = pobyso_constant_from_int(precSa) |
| 522 | 522 |
pobyso_set_prec_so_so(precSo) |
| 523 | 523 |
sollya_lib_clear_obj(precSo) |
| 524 |
sollyaPrecChanged = True |
|
| 524 |
sollyaPrecChanged = True |
|
| 525 |
## Free variable stuff. |
|
| 526 |
freeVariableNameChanged = False |
|
| 527 |
polyFreeVariableNameSa = \ |
|
| 528 |
str(polySa.variables()[0]) |
|
| 529 |
currentFreeVariableNameSa = pobyso_get_free_variable_name_so_sa() |
|
| 530 |
if polyFreeVariableNameSa != currentFreeVariableNameSa: |
|
| 531 |
#print "Free variable names do not match.", polyFreeVariableNameSa |
|
| 532 |
sollya_lib_name_free_variable(polyFreeVariableNameSa) |
|
| 533 |
freeVariableNameChanged = True |
|
| 525 | 534 |
## Get exponents and coefficients. |
| 526 | 535 |
exponentsSa = polySa.exponents() |
| 527 | 536 |
coefficientsSa = polySa.coefficients() |
| ... | ... | |
| 546 | 555 |
polySo = sollya_lib_build_function_add(polySo, polyTermSo) |
| 547 | 556 |
if sollyaPrecChanged: |
| 548 | 557 |
pobyso_set_prec_so_so(initialSollyaPrecSo) |
| 549 |
sollya_lib_clear_obj(initialSollyaPrecSo) |
|
| 558 |
sollya_lib_clear_obj(initialSollyaPrecSo) |
|
| 559 |
## Do not set back the free variable name in Sollya to its initial value: |
|
| 560 |
# it will change it back, in the Sollya polynomial, to what it was in the |
|
| 561 |
# first place. |
|
| 550 | 562 |
return polySo |
| 551 | 563 |
# End pobyso_float_poly_sa_so |
| 552 | 564 |
|
| ... | ... | |
| 1569 | 1581 |
def pobyso_rat_poly_sa_so(polySa, precSa = None): |
| 1570 | 1582 |
""" |
| 1571 | 1583 |
Create a Sollya polynomial from a Sage rational polynomial. |
| 1584 |
We first convert the rational polynomial into a floating-point |
|
| 1585 |
polynomial. |
|
| 1572 | 1586 |
""" |
| 1573 | 1587 |
## TODO: filter arguments. |
| 1574 | 1588 |
## Precision. If no precision is given, use the current precision |
| ... | ... | |
| 1581 | 1595 |
P_RRR = RRR[polySa.variables()[0]] |
| 1582 | 1596 |
polyFloatSa = P_RRR(polySa) |
| 1583 | 1597 |
## Make sure no precision is provided: pobyso_float_poly_sa_so will |
| 1584 |
# recover it all by itself and not make an extra conversion.
|
|
| 1598 |
# recover it all by itself and will not make any extra conversion.
|
|
| 1585 | 1599 |
return pobyso_float_poly_sa_so(polyFloatSa) |
| 1586 | 1600 |
|
| 1587 | 1601 |
# End pobyso_rat_poly_sa_so |
| ... | ... | |
| 1807 | 1821 |
currentApproxErrorSo, |
| 1808 | 1822 |
approxAccurSo, |
| 1809 | 1823 |
debug=False): |
| 1824 |
""" |
|
| 1825 |
From an input approximation polynomial, compute an output one with |
|
| 1826 |
smaller coefficients and yet yields a sufficient approximation accuracy. |
|
| 1827 |
""" |
|
| 1810 | 1828 |
if debug: |
| 1811 | 1829 |
print "Input arguments:" |
| 1812 | 1830 |
print "Polynomial: ", ; pobyso_autoprint(polySo) |
| ... | ... | |
| 1822 | 1840 |
# no possible gain. |
| 1823 | 1841 |
if currentApproxErrorSa >= approxAccurSa / 2: |
| 1824 | 1842 |
#### Do not return the initial argument but copies: caller may free |
| 1825 |
# as inutile after call. |
|
| 1843 |
# the former as inutile after call.
|
|
| 1826 | 1844 |
return (sollya_lib_copy_obj(polySo), |
| 1827 | 1845 |
sollya_lib_copy_obj(currentApproxErrorSo)) |
| 1828 | 1846 |
degreeSa = pobyso_polynomial_degree_so_sa(polySo) |
| ... | ... | |
| 1833 | 1851 |
print "intervalSa :", intervalSa.str(style='brackets') |
| 1834 | 1852 |
print "currentApproxErrorSa :", currentApproxErrorSa |
| 1835 | 1853 |
print "approxAccurSa :", approxAccurSa |
| 1836 |
### Start with a 0 value expression. |
|
| 1837 | 1854 |
radiusSa = intervalSa.absolute_diameter() / 2 |
| 1838 | 1855 |
if debug: |
| 1839 | 1856 |
print "log2(radius):", RR(radiusSa).log2() |
| 1840 | 1857 |
iterIndex = 0 |
| 1858 |
## Build the "shaved" polynomial. |
|
| 1841 | 1859 |
while True: |
| 1860 |
### Start with a 0 value expression. |
|
| 1842 | 1861 |
resPolySo = pobyso_constant_0_sa_so() |
| 1843 | 1862 |
roundedPolyApproxAccurSa = approxAccurSa / 2 |
| 1844 | 1863 |
currentRadiusPowerSa = 1 |
| ... | ... | |
| 1851 | 1870 |
else: |
| 1852 | 1871 |
roundingPowerSa = \ |
| 1853 | 1872 |
floor(((currentRadiusPowerSa/roundedPolyApproxAccurSa)*(degreeSa+1)).log2()) |
| 1854 |
## Under extreme conditions the above formulas can evaluate under 2, which is the
|
|
| 1855 |
# minimal precision of an MPFR number. |
|
| 1873 |
## Under extreme conditions the above formulas can evaluate under 2, |
|
| 1874 |
# which is the minimal precision of an MPFR number.
|
|
| 1856 | 1875 |
if roundingPowerSa < 2: |
| 1857 | 1876 |
roundingPowerSa = 2 |
| 1858 | 1877 |
if debug: |
| ... | ... | |
| 1901 | 1920 |
# free the former as inutile after call. |
| 1902 | 1921 |
return (sollya_lib_copy_obj(polySo), |
| 1903 | 1922 |
sollya_lib_copy_obj(currentApproxErrorSo)) |
| 1904 |
else: # Round 0, got round 1
|
|
| 1923 |
else: # Round 0 (agressive rounding), got round 1 (proved rounding)
|
|
| 1905 | 1924 |
sollya_lib_clear_obj(resPolySo) |
| 1906 | 1925 |
sollya_lib_clear_obj(infNormSo) |
| 1907 | 1926 |
iterIndex += 1 |
| ... | ... | |
| 2018 | 2037 |
""" |
| 2019 | 2038 |
Compute the Taylor expansion without the variable change |
| 2020 | 2039 |
x -> x-intervalCenter. |
| 2040 |
If errorTypeSo is None, absolute is used. |
|
| 2041 |
If sollyaPrecSo is None, Sollya internal precision is not changed. |
|
| 2021 | 2042 |
""" |
| 2022 | 2043 |
# Change internal Sollya precision, if needed. |
| 2023 | 2044 |
(initialSollyaPrecSo, initialSollyaPrecSa) = pobyso_get_prec_so_so_sa() |
| ... | ... | |
| 2037 | 2058 |
taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, |
| 2038 | 2059 |
intervalCenterSo, |
| 2039 | 2060 |
rangeSo, errorTypeSo, None) |
| 2040 |
# taylorFormListSaSo is a Python list of Sollya objects references that
|
|
| 2041 |
# are copies of the elements of taylorFormSo. |
|
| 2061 |
# Object taylorFormListSaSo is a Python list of Sollya objects references
|
|
| 2062 |
# that are copies of the elements of taylorFormSo.
|
|
| 2042 | 2063 |
# pobyso_get_list_elements_so_so clears taylorFormSo. |
| 2043 | 2064 |
(taylorFormListSaSo, numElementsSa, isEndEllipticSa) = \ |
| 2044 | 2065 |
pobyso_get_list_elements_so_so(taylorFormSo) |
| 2066 |
## Copy needed here since polySo will be returned and taylorFormListSaSo |
|
| 2067 |
# will be cleared. |
|
| 2045 | 2068 |
polySo = sollya_lib_copy_obj(taylorFormListSaSo[0]) |
| 2046 | 2069 |
#print "Num elements:", numElementsSa |
| 2047 | 2070 |
sollya_lib_clear_obj(taylorFormSo) |
| 2048 |
#polySo = taylorFormListSaSo[0] |
|
| 2049 |
#errorRangeSo = sollya_lib_copy_obj(taylorFormListSaSo[2]) |
|
| 2050 |
errorRangeSo = taylorFormListSaSo[2] |
|
| 2051 | 2071 |
# No copy_obj needed here: a new objects are created. |
| 2052 |
maxErrorSo = sollya_lib_sup(errorRangeSo) |
|
| 2053 |
minErrorSo = sollya_lib_inf(errorRangeSo) |
|
| 2072 |
maxErrorSo = sollya_lib_sup(taylorFormListSaSo[2]) |
|
| 2073 |
minErrorSo = sollya_lib_inf(taylorFormListSaSo[2]) |
|
| 2074 |
# List taylorFormListSaSo is not needed anymore. |
|
| 2075 |
pobyso_clear_taylorform_sa_so(taylorFormListSaSo) |
|
| 2054 | 2076 |
absMaxErrorSo = sollya_lib_abs(maxErrorSo) |
| 2055 | 2077 |
absMinErrorSo = sollya_lib_abs(minErrorSo) |
| 2056 | 2078 |
sollya_lib_clear_obj(maxErrorSo) |
| 2057 | 2079 |
sollya_lib_clear_obj(minErrorSo) |
| 2058 | 2080 |
absMaxErrorSa = pobyso_get_constant_as_rn_so_sa(absMaxErrorSo) |
| 2059 | 2081 |
absMinErrorSa = pobyso_get_constant_as_rn_so_sa(absMinErrorSo) |
| 2060 |
# If changed, reset the Sollya working precision. |
|
| 2082 |
# |
|
| 2083 |
if errorTypeIsNone: |
|
| 2084 |
sollya_lib_clear_obj(errorTypeSo) |
|
| 2085 |
## If changed, reset the Sollya working precision. |
|
| 2061 | 2086 |
if sollyaPrecChanged: |
| 2062 | 2087 |
sollya_lib_set_prec(initialSollyaPrecSo) |
| 2063 | 2088 |
sollya_lib_clear_obj(initialSollyaPrecSo) |
| 2064 |
if errorTypeIsNone: |
|
| 2065 |
sollya_lib_clear_obj(errorTypeSo) |
|
| 2066 |
pobyso_clear_taylorform_sa_so(taylorFormListSaSo) |
|
| 2089 |
## According to what error is the largest, return the errors. |
|
| 2067 | 2090 |
if absMaxErrorSa > absMinErrorSa: |
| 2068 | 2091 |
sollya_lib_clear_obj(absMinErrorSo) |
| 2069 | 2092 |
return (polySo, intervalCenterSo, absMaxErrorSo) |
| ... | ... | |
| 2099 | 2122 |
taylorFormSo = sollya_lib_taylorform(functionSo, degreeSo, \ |
| 2100 | 2123 |
intervalCenterSo, \ |
| 2101 | 2124 |
rangeSo, errorTypeSo, None) |
| 2102 |
# taylorFormListSaSo is a Python list of Sollya objects references that
|
|
| 2103 |
# are copies of the elements of taylorFormSo. |
|
| 2125 |
# Object taylorFormListSaSo is a Python list of Sollya objects references
|
|
| 2126 |
# that are copies of the elements of taylorFormSo.
|
|
| 2104 | 2127 |
# pobyso_get_list_elements_so_so clears taylorFormSo. |
| 2105 |
(taylorFormListSo, numElements, isEndElliptic) = \ |
|
| 2128 |
(taylorFormListSaSo, numElements, isEndElliptic) = \
|
|
| 2106 | 2129 |
pobyso_get_list_elements_so_so(taylorFormSo) |
| 2107 |
polySo = taylorFormListSo[0] |
|
| 2108 |
errorRangeSo = taylorFormListSo[2] |
|
| 2109 |
maxErrorSo = sollya_lib_sup(errorRangeSo) |
|
| 2110 |
minErrorSo = sollya_lib_inf(errorRangeSo) |
|
| 2130 |
sollya_lib_clear_obj(taylorFormSo) |
|
| 2131 |
polySo = taylorFormListSaSo[0] |
|
| 2132 |
## Maximum error computation with taylorFormListSaSo[2], a range |
|
| 2133 |
# holding the actual error. Bounds can be negative. |
|
| 2134 |
maxErrorSo = sollya_lib_sup(taylorFormListSaSo[2]) |
|
| 2135 |
minErrorSo = sollya_lib_inf(taylorFormListSaSo[2]) |
|
| 2111 | 2136 |
absMaxErrorSo = sollya_lib_abs(maxErrorSo) |
| 2112 | 2137 |
absMinErrorSo = sollya_lib_abs(minErrorSo) |
| 2113 | 2138 |
sollya_lib_clear_obj(maxErrorSo) |
| ... | ... | |
| 2118 | 2143 |
sollya_lib_build_function_free_variable(),\ |
| 2119 | 2144 |
sollya_lib_copy_obj(intervalCenterSo)) |
| 2120 | 2145 |
polyVarChangedSo = sollya_lib_evaluate(polySo, changeVarExpSo) |
| 2121 |
sollya_lib_clear_obj(polySo) |
|
| 2146 |
# List taylorFormListSaSo is not needed anymore. |
|
| 2147 |
pobyso_clear_taylorform_sa_so(taylorFormListSaSo) |
|
| 2122 | 2148 |
sollya_lib_clear_obj(changeVarExpSo) |
| 2123 | 2149 |
# If changed, reset the Sollya working precision. |
| 2124 | 2150 |
if sollyaPrecChanged: |
| ... | ... | |
| 2126 | 2152 |
sollya_lib_clear_obj(initialSollyaPrecSo) |
| 2127 | 2153 |
if errorTypeIsNone: |
| 2128 | 2154 |
sollya_lib_clear_obj(errorTypeSo) |
| 2129 |
sollya_lib_clear_obj(taylorFormSo) |
|
| 2130 | 2155 |
# Do not clear maxErrorSo. |
| 2131 | 2156 |
if absMaxErrorSa > absMinErrorSa: |
| 2132 | 2157 |
sollya_lib_clear_obj(absMinErrorSo) |
Formats disponibles : Unified diff