1

I'm using the bisection method from the scipy.optimize package within a for loop. The idea is to get a value of "sig" with the bisection method for each element (value) in the "eps_komp" vector. I've coded this much:

import numpy as np
import scipy.optimize as optimize

K=300
n = 0.43
E = 210000
Rm = 700
sig_a = []
RO_K = 300
RO_n = 0.43

eps_komp =     [0.00012893048999999997,
0.018839115269999998,
0.01230539995,
0.022996934109999999,
-0.0037319012899999999,
0.023293921169999999,
0.0036927752099999997,
0.020621037629999998,
0.0063656587500000002,
0.020324050569999998,
-0.0025439530500000001,
0.018542128209999998,
0.01230539995,
0.019730076449999998,
0.0045837363899999999,
0.015275270549999997,
-0.0040288883499999999,
0.021215011749999999,
-0.0031379271699999997,
0.023590908229999999]

def eps_f(i):
    return eps_komp[i]

for j in range(len(eps_komp)):
    eps_komp_j = eps_f(j)
    if j <= len(eps_komp):
        def func(sig):
            return eps_komp_j - sig/E - (sig/RO_K)**(1/RO_n) 
        sig_a.append(optimize.bisect(func, 0, Rm))
    else:
        break

print(sig_a)

Now if I change the the value of "j" in eps_f(j) to 0:

eps_komp_j = eps_f(0)

it works, and so it does for all other values that I insert by hand, but if I keep it like it is in the for loop, the "j" value doesnt change automatically and I get an error:

f(a) and f(b) must have different signs

Has anyone a clue what is the problem and how could this be solved?

Regards,

L

P.S. I did post another topic on this problem yesterday, but I wasnt very specific with the problem and got negative feedback. However, I do need to solve this today, so I was forced to post it again, however I did manage to get a bit further with the code then I did in the earlier post, so it isn't a repost...

mcluka
  • 285
  • 1
  • 7
  • 17
  • 2
    When you get an error, please paste the *full* traceback, not just the final message. – Wayne Werner Apr 05 '16 at 15:13
  • I can't believe that `j` isn't being updated by the for-loop correctly. What happens if you insert `print(j, eps_komp_j)` before `sig_a.append(optimize.bisect...`? – Alex Hall Apr 05 '16 at 15:26
  • Funny, cause when I run this code I get `NameError: name 'Emod' is not defined` so you might want to look at how to post a [Minimal, Complete, Verifiable Example](http://stackoverflow.com/help/mcve) – Wayne Werner Apr 05 '16 at 15:32
  • oops, my bad, Emod=E, this code is just a part of a bigger code so I only posted some parts. E is defined withint this code, you can just change the Emod to E within the function. I updated the code now. – mcluka Apr 05 '16 at 15:37

1 Answers1

1

If you read the docs you'll find that:

Basic bisection routine to find a zero of the function f between the arguments a and b. f(a) and f(b) cannot have the same signs. Slow but sure.

In your code:

    def func(sig):
        return eps_komp_j - sig/Emod - (sig/RO_K)**(1/RO_n) 
    sig_a.append(optimize.bisect(func, 0, Rm))

You're passing it func(0) and func(700).

By replacing the optimize.bisect line with print(func(0), func(700)) I get the following output:

0.00012893048999999997 -7.177181168628421
0.018839115269999998 -7.158470983848421
0.01230539995 -7.165004699168421
0.02299693411 -7.15431316500842
-0.00373190129 -7.1810420004084206
0.02329392117 -7.154016177948421
0.0036927752099999997 -7.173617323908421
0.02062103763 -7.156689061488421
0.00636565875 -7.17094444036842
0.02032405057 -7.156986048548421
-0.00254395305 -7.17985405216842
0.018542128209999998 -7.15876797090842
0.01230539995 -7.165004699168421
0.019730076449999998 -7.157580022668421
0.00458373639 -7.172726362728421
0.015275270549999997 -7.162034828568421
-0.00402888835 -7.181338987468421
0.02121501175 -7.156095087368421
-0.0031379271699999997 -7.1804480262884205
0.02359090823 -7.153719190888421

Note the multiple pairs that have the same signs. optimize.bisect can't handle those. I don't know what you're trying to accomplish, but this is the wrong approach.

Wayne Werner
  • 49,299
  • 29
  • 200
  • 290
  • if you run the code that I just updated and assign the "j" an actual value here: eps_komp_j = eps_f(j) it works, so this works right for sure, the only thing that isnt working properly is the for loop... – mcluka Apr 05 '16 at 15:43
  • What do you mean "it works"? I have `eps_komp_j = eps_f(j)` (which, at least at the moment does nothing useful) the code that produced these outputs. `optimize.bisect` will never work in this particular configuration. – Wayne Werner Apr 05 '16 at 15:57
  • well if you insert `eps_komp_j = eps_f(3)` for example, the vector will be filled with the right value, btu only for the 4th (j=3) element of vector eps_komp... – mcluka Apr 05 '16 at 15:59
  • That's because it produces `0.02299693411 -7.15431316500842` for the values that are passed to `bisect`. Did you read the documentation, or my answer at all? If there's a point that I made you do not understand I'd be happy to clarify, but assuming that people know exactly what your thinking is a dangerous assumption to make. – Wayne Werner Apr 05 '16 at 16:01
  • ah, I'm sorry, I see the point that you are trying to make, there are cases where both are negative so I go pass 0... Well I don't know what else I could do... – mcluka Apr 05 '16 at 16:23
  • What is it that you're *trying* to do? `bisect` is for finding the zero of a function. Is that what you're trying to do? If so, you can't find the zero of a function that never crosses zero. And if you can't prove that your function is crossing zero, how can scipy guarantee that it won't try to find a zero that doesn't exist? – Wayne Werner Apr 05 '16 at 16:34
  • Ok so I managed to solve the problem, I just made all eps_komp values positive. If so, it will always pass the zero and I will get the value. This is also physical correct, because I'm making a stress-strain diagram based on the Ramberg-Osgood equation and if strain is negative (compression) I now used its absolute values... I do need to confess, I'm new to python so ofcourse I have some basic problems, but as far as mathematically, physically or numerically goes, I know very well what I'm doing and what I would like to achieve. Thanks for your help Wayne! – mcluka Apr 05 '16 at 16:46