Ask Your Question
1

How to solve numericaly with arbitrary precision

asked 2020-03-09 10:19:57 +0100

JonasA gravatar image

updated 2020-03-09 17:19:42 +0100

Hello,

I want to solve a numerical equation for which I can only access by a lambda function, the M matrix is to big for the determinant to be computed on the symbolic ring.

f = lambda om: RRR(M.subs({omega:om}).change_ring(RRR).det())

But , find-root solve in the built in float type of Python which is lacking the precision I need. Is there a way in sage to solve numericaly with arbitrary precision?

Thank you

Regards

(edited with example more in line with what I am trying to achieve)

sage: x=var('x')
sage: M= Matrix(SR,[[cos(x),cosh(x)],[sin(x),sinh(x)]])
sage: RRR = RealField(200)
sage: f = lambda om: M.subs({x:om}).change_ring(RRR).det()
sage: find_root(f,1,2)
edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
4

answered 2020-03-09 14:47:59 +0100

vdelecroix gravatar image

updated 2020-03-11 10:21:05 +0100

This depends very much on the type of equation you are considering. For polynomial equations

sage: x = polygen(ZZ)
sage: (x**5 - x - 1).roots(RealField(128))
[1.1673039782614186842560458998548421807]

If you want a more specific answer, ask a more precise question.

EDITED:

Sadly, Sage has no rountines for more general arbitrary general root finding. Here is a straightforward implementation of the dichotomy

def dichotomy(f, left, right, n):
    sl = f(left).sign()
    sr = f(right).sign()
    if sl == sr:
        raise ValueError
    for _ in range(n):
        middle = (left + right) / 2
        sm = f(middle).sign()
        if sm == 0:
            return sm
        elif sm == sl:
            left = middle
        else:
            right = middle
    return (left,right)

which you can use as

sage: f = lambda x: x.sin()
sage: R = RealField(256)
sage: dichotomy(f, R(3), R(3.5), 10)
(3.141113281250000000000000000000000000000000000000000000000000000000000000000,
 3.141601562500000000000000000000000000000000000000000000000000000000000000000)
sage: dichotomy(f, R(3), R(3.5), 256)
(3.141592653589793238462643383279502884197169399375105820974944592307816406286,
 3.141592653589793238462643383279502884197169399375105820974944592307816406286)

EDITED BIS: For a non-trivial example involving functions with cos/sin/cosh/sinh

sage: f(x) = -2*cos(x)**2 + 5 * sin(x) + cos(3*x) * cosh(x) - sinh(2*x - 1)
sage: g = fast_callable(f, vars=[x], domain=R)
sage: dichotomy(g, R(-0.3), R(0), 100)
(-0.09451007784336178406072003155719521150577190161782642363759963932690055798957,
 -0.09451007784336178406072003155695855323420559807620290506801480345799859605904)

To find the appropriate initial conditions for f you need to either carry on some analysis on paper or look at its graph via

sage: plot(f)
edit flag offensive delete link more

Comments

Thank you for answering. I have a matrix of hyperbolic and trigonometric function of a variable, I want to know for which values of the variable the matrix is singular. It is analog to the physics problem of knowing the natural frequencies of an oscillator. Although my matrix is not very dense, it is dense enough so that the determinant can't be computed in the symbolic ring. I have to use a lambda function to substitute de variable by its value, change the matrix to the real ring then compute the determinant (the computation of the determinant is too slow if I stay in the symbolic ring). Which is the lambda function I write on my first message with RRR a custom Real Field with 200 bits precision.

A code example will follow in the next message

JonasA gravatar imageJonasA ( 2020-03-09 17:11:11 +0100 )edit

So

sage: x=var('x')
sage: M= Matrix(SR,[[cos(x),cosh(x)],[sin(x),sinh(x)]])
sage: RRR = RealField(200)

Then

 sage: f = lambda om: RRR(M.subs({x:om}).change_ring(RRR).det())

gives f a very long runtime, find_root don't converge on a reasonable time and

sage: f = lambda om: M.subs({x:om}).change_ring(RRR).det()
sage: find_root(f,1,2)

while f has a much faster runtime, throw an exception

JonasA gravatar imageJonasA ( 2020-03-09 17:13:21 +0100 )edit

Dichotomy solve the problem although my function is too much unsmooth to be used so I will stop trying. My function was taking those value: (f(37) was returning 0, but when I used a better resolution it was more 1e67)

sage: f(36.99)
-1.0936253623915059621862511135588106826765847154466062182129e99
sage: f(36.999)
1.0936253623915059621862511135588106826765847154466062182129e99
sage: f(36.99999)
-5.4681268119575298109312555677940534133829235772330310910644e98
sage: f(36.9999999)
5.4681268119575298109312555677940534133829235772330310910644e98
sage: f(36.999999999)
-1.6404380435872589432793766703382160240148770731699093273193e99
JonasA gravatar imageJonasA ( 2020-03-12 14:07:01 +0100 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

2 followers

Stats

Asked: 2020-03-09 10:19:57 +0100

Seen: 825 times

Last updated: Mar 11 '20