ASKSAGE: Sage Q&A Forum - RSS feedhttps://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Mon, 02 Jul 2018 10:11:25 +0200Solving a power series equationhttps://ask.sagemath.org/question/42778/solving-a-power-series-equation/I'd like to solve the equation `f(t) == 1` where `f` is a power series:
var('n,t')
a1 = float(4 + sqrt(8))
b1 = float(3 + sqrt(8))
c1 = b1
an = 4*n*(n-1)+a1*n
an1 = 4*n*(n+1) + a1*(n+1)
anr = b1 + 2*n*(4*n+a1)
bn = 1+2*n*(n-1)+n*(b1-1)
bn1 = 1+2*n*(n+1)+(n+1)*(b1-1)
bnr = c1 + 2*n*(2*n+b1-1)
def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum(an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t ,n,1,1000).n()
I've tried simply plugging in values for t and I estimate that if f(t) = 1 then t is approximately 1.51. I'd like to possibly do something like
solve(f(t) == 1,t)
but this appears to lock sage up and I have to interrupt the process. The true f function is actually an infinite series, but I am truncating it to the first 1000 terms. I know there are methods of solving series equations such as the method of regula falsi. Does Sage have anything like this? Thanks!Fri, 29 Jun 2018 19:52:53 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/Comment by slelievre for <p>I'd like to solve the equation <code>f(t) == 1</code> where <code>f</code> is a power series:</p>
<pre><code>var('n,t')
a1 = float(4 + sqrt(8))
b1 = float(3 + sqrt(8))
c1 = b1
an = 4*n*(n-1)+a1*n
an1 = 4*n*(n+1) + a1*(n+1)
anr = b1 + 2*n*(4*n+a1)
bn = 1+2*n*(n-1)+n*(b1-1)
bn1 = 1+2*n*(n+1)+(n+1)*(b1-1)
bnr = c1 + 2*n*(2*n+b1-1)
def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum(an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t ,n,1,1000).n()
</code></pre>
<p>I've tried simply plugging in values for t and I estimate that if f(t) = 1 then t is approximately 1.51. I'd like to possibly do something like</p>
<pre><code>solve(f(t) == 1,t)
</code></pre>
<p>but this appears to lock sage up and I have to interrupt the process. The true f function is actually an infinite series, but I am truncating it to the first 1000 terms. I know there are methods of solving series equations such as the method of regula falsi. Does Sage have anything like this? Thanks!</p>
https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42807#post-id-42807If an answer solves your question, please accept it by clicking the
"accept" button (the one with a check mark, at the top left of the answer).
This will mark the question as solved in the Ask Sage's list of questions,
as well as in lists of questions related to a particular query or keyword.Mon, 02 Jul 2018 10:11:25 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42807#post-id-42807Answer by rburing for <p>I'd like to solve the equation <code>f(t) == 1</code> where <code>f</code> is a power series:</p>
<pre><code>var('n,t')
a1 = float(4 + sqrt(8))
b1 = float(3 + sqrt(8))
c1 = b1
an = 4*n*(n-1)+a1*n
an1 = 4*n*(n+1) + a1*(n+1)
anr = b1 + 2*n*(4*n+a1)
bn = 1+2*n*(n-1)+n*(b1-1)
bn1 = 1+2*n*(n+1)+(n+1)*(b1-1)
bnr = c1 + 2*n*(2*n+b1-1)
def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum(an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t ,n,1,1000).n()
</code></pre>
<p>I've tried simply plugging in values for t and I estimate that if f(t) = 1 then t is approximately 1.51. I'd like to possibly do something like</p>
<pre><code>solve(f(t) == 1,t)
</code></pre>
<p>but this appears to lock sage up and I have to interrupt the process. The true f function is actually an infinite series, but I am truncating it to the first 1000 terms. I know there are methods of solving series equations such as the method of regula falsi. Does Sage have anything like this? Thanks!</p>
https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?answer=42781#post-id-42781What your code for `f(t)` does is first take the sum symbolically up to 1000 terms, and then evaluate numerically.
What you want to do instead is calculate the terms numerically and sum them.
For example, like this:
def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum((an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t).subs(n=k).n() for k in xrange(1,1000))
The method `solve` is used to find exact solutions. For approximate solutions you can use e.g. `find_root`:
find_root(lambda t: f(t) - 1, 1.48, 1.52)
It yields `1.5111899890344738`, and `f(1.5111899890344738)` is approximately `1.00000000000000`.Fri, 29 Jun 2018 20:29:51 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?answer=42781#post-id-42781Comment by Daniel L for <p>What your code for <code>f(t)</code> does is first take the sum symbolically up to 1000 terms, and then evaluate numerically.</p>
<p>What you want to do instead is calculate the terms numerically and sum them.</p>
<p>For example, like this:</p>
<pre><code>def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum((an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t).subs(n=k).n() for k in xrange(1,1000))
</code></pre>
<p>The method <code>solve</code> is used to find exact solutions. For approximate solutions you can use e.g. <code>find_root</code>:</p>
<pre><code>find_root(lambda t: f(t) - 1, 1.48, 1.52)
</code></pre>
<p>It yields <code>1.5111899890344738</code>, and <code>f(1.5111899890344738)</code> is approximately <code>1.00000000000000</code>.</p>
https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42794#post-id-42794When I run this code I get the following traceback error: ValueError: rtol too small (4.5e-16 < 8.88178e-16)
I'm not sure what this means or why it's happening.Sat, 30 Jun 2018 22:58:02 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42794#post-id-42794Comment by rburing for <p>What your code for <code>f(t)</code> does is first take the sum symbolically up to 1000 terms, and then evaluate numerically.</p>
<p>What you want to do instead is calculate the terms numerically and sum them.</p>
<p>For example, like this:</p>
<pre><code>def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum((an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t).subs(n=k).n() for k in xrange(1,1000))
</code></pre>
<p>The method <code>solve</code> is used to find exact solutions. For approximate solutions you can use e.g. <code>find_root</code>:</p>
<pre><code>find_root(lambda t: f(t) - 1, 1.48, 1.52)
</code></pre>
<p>It yields <code>1.5111899890344738</code>, and <code>f(1.5111899890344738)</code> is approximately <code>1.00000000000000</code>.</p>
https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42795#post-id-42795That seems to be [a bug](https://ask.sagemath.org/question/37540/valueerror-rtol-too-small/) in some version(s) of Sage. You can work around it by specifying the `rtol` keyword argument to `find_root` explicitly, or using the `partial` method in the thread I just linked, or (I would hope) by upgrading to the latest version of SageMath.Sat, 30 Jun 2018 23:25:18 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42795#post-id-42795Comment by Daniel L for <p>What your code for <code>f(t)</code> does is first take the sum symbolically up to 1000 terms, and then evaluate numerically.</p>
<p>What you want to do instead is calculate the terms numerically and sum them.</p>
<p>For example, like this:</p>
<pre><code>def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum((an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t).subs(n=k).n() for k in xrange(1,1000))
</code></pre>
<p>The method <code>solve</code> is used to find exact solutions. For approximate solutions you can use e.g. <code>find_root</code>:</p>
<pre><code>find_root(lambda t: f(t) - 1, 1.48, 1.52)
</code></pre>
<p>It yields <code>1.5111899890344738</code>, and <code>f(1.5111899890344738)</code> is approximately <code>1.00000000000000</code>.</p>
https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42796#post-id-42796The workaround on the above thread works well. Thanks for the help!Sun, 01 Jul 2018 00:13:48 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42796#post-id-42796Comment by rburing for <p>What your code for <code>f(t)</code> does is first take the sum symbolically up to 1000 terms, and then evaluate numerically.</p>
<p>What you want to do instead is calculate the terms numerically and sum them.</p>
<p>For example, like this:</p>
<pre><code>def f(t):
return (1+b1)^-t + 3/2*a1^-t + 1/2*(2*b1)^-t + b1^-t + 2*sum((an^-t +an1^-t + anr^-t + 2*bn^-t +2*bn1^-t + bnr^-t).subs(n=k).n() for k in xrange(1,1000))
</code></pre>
<p>The method <code>solve</code> is used to find exact solutions. For approximate solutions you can use e.g. <code>find_root</code>:</p>
<pre><code>find_root(lambda t: f(t) - 1, 1.48, 1.52)
</code></pre>
<p>It yields <code>1.5111899890344738</code>, and <code>f(1.5111899890344738)</code> is approximately <code>1.00000000000000</code>.</p>
https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42801#post-id-42801You're welcome. You can mark this answer as correct.Sun, 01 Jul 2018 12:46:26 +0200https://ask.sagemath.org/question/42778/solving-a-power-series-equation/?comment=42801#post-id-42801