質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう!

ただいまの
回答率

91.02%

  • Python

    5540questions

    Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

  • NumPy

    320questions

    NumPyはPythonのプログラミング言語の科学的と数学的なコンピューティングに関する拡張モジュールです。

定義した関数がうまく使えない。quadpack.error: Supplied function does not return a valid float.

受付中

回答 1

投稿 編集

  • 評価
  • クリップ 0
  • VIEW 219

youseiranbu7

score 0

前提・実現したいこと

積分を含むある関数を定義して、使おうとしているのですが、具体的な値を代入するとうまくいきません。

最終的にやりたいことは、以下の関数m_model_00(x)のグラフをxy平面上に書くことです。

発生している問題・エラーメッセージ

quadpack.error: Supplied function does not return a valid float.

該当のソースコード

python

import math
import numpy as np
import scipy.optimize
from scipy import integrate
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
import scipy

def dL(z,omega_m,omega_l,H0,c):
def integrand(x):
return 1. / np.sqrt((1.+x)**2.*(1.+x*omega_m)-x*(2.+x)*omega_l)

if omega_m + omega_l ==1:
i=integrate.quad(integrand,0.,z)
return c*(1.+z)/H0*i[0]
elif omega_m + omega_l > 1:
i=integrate.quad(integrand,0.,z)
return  c*(1.+z)/(H0*np.sqrt(abs(1.-omega_m-omega_l)))*np.sin(np.sqrt(abs(1.-omega_m-omega_l))*i[0])
else:
i=integrate.quad(integrand,0.,z)
return  c*(1.+z)/(H0*np.sqrt(abs(1.-omega_m-omega_l)))*np.sinh(np.sqrt(abs(1.-omega_m-omega_l))*i[0])

dL = numpy.vectorize(dL, excluded=("omega_m","omega_l","H0","c"))

def m_model(z,omega_m,omega_l,H0,MB,):
return MB + 5.0*np.log(dL(z,omega_m,omega_l,H0,1.))+25.0

def m_model_00(x):
return m_model(x,0.,0.,0.7,-1.8,1.)

以下、ターミナルで実行した結果です。

execfile("problem_set_04_1b.py")
m_model_00(0.4)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "problem_set_04_1b.py", line 44, in m_model_00
return m_model(x,0.,0.,0.7,-1.8)
File "problem_set_04_1b.py", line 31, in m_model
return MB + 5.0*np.log(dL(z,omega_m,omega_l,H0))+25.0
File "problem_set_04_1b.py", line 26, in dL
return  c*(1.+z)/(H0*np.sqrt(abs(1.-omega_m-omega_l)))*np.sinh(np.sqrt(abs(1.-omega_m-omega_l))*integrate.quad(integrand,0.,z))
File "/System/Library/Frameworks/Python.framework/Versions/2.7/Extras/lib/python/scipy/integrate/quadpack.py", line 252, in quad
retval = _quad(func,a,b,args,full_output,epsabs,epsrel,limit,points)
File "/System/Library/Frameworks/Python.framework/Versions/2.7/Extras/lib/python/scipy/integrate/quadpack.py", line 317, in _quad
return _quadpack._qagse(func,a,b,args,full_output,epsabs,epsrel,limit)
quadpack.error: Supplied function does not return a valid float.

補足情報(言語/FW/ツール等のバージョンなど)

より詳細な情報

  • 気になる質問をクリップする

    クリップした質問は、後からいつでもマイページで確認できます。

    またクリップした質問に回答があった際、通知やメールを受け取ることができます。

    クリップを取り消します

  • 良い質問の評価を上げる

    以下のような質問は評価を上げましょう

    • 質問内容が明確
    • 自分も答えを知りたい
    • 質問者以外のユーザにも役立つ

    評価が高い質問は、TOPページの「注目」タブのフィードに表示されやすくなります。

    質問の評価を上げたことを取り消します

  • 評価を下げられる数の上限に達しました

    評価を下げることができません

    • 1日5回まで評価を下げられます
    • 1日に1ユーザに対して2回まで評価を下げられます

    質問の評価を下げる

    teratailでは下記のような質問を「具体的に困っていることがない質問」、「サイトポリシーに違反する質問」と定義し、推奨していません。

    • プログラミングに関係のない質問
    • やってほしいことだけを記載した丸投げの質問
    • 問題・課題が含まれていない質問
    • 意図的に内容が抹消された質問
    • 広告と受け取られるような投稿

    評価が下がると、TOPページの「アクティブ」「注目」タブのフィードに表示されにくくなります。

    質問の評価を下げたことを取り消します

    この機能は開放されていません

    評価を下げる条件を満たしてません

    評価を下げる理由を選択してください

    詳細な説明はこちら

    上記に当てはまらず、質問内容が明確になっていない質問には「情報の追加・修正依頼」機能からコメントをしてください。

    質問の評価を下げる機能の利用条件

    この機能を利用するためには、以下の事項を行う必要があります。

質問への追記・修正、ベストアンサー選択の依頼

  • LouiS0616

    2017/11/14 02:27

    インデントが崩れて読めません。https://teratail.com/help#about-markdown こちらに従って、プレビューを確認しながらコードブロックを設定してください。

    キャンセル

回答 1

0

エラーの
File "problem_set_04_1b.py", line 26, in dL
return  c*(1.+z)/(H0*np.sqrt(abs(1.-omega_m-omega_l)))*np.sinh(np.sqrt(abs(1.-omega_m-omega_l))*integrate.quad(integrand,0.,z))
の部分に相当するものが上に貼り付けてあるコードに含まれてないのは気のせいでしょうか?


追記:
不思議な質問ですね。
貼り付けてあるコードは普通に動きます。(一部numpy->npのミスがありますが)
正しく修正した後に過去のエラーを勘違いしていませんか?
i=integrate.quad(integrand,0.,z)
.....*i[0]
とご自身で訂正されています?

import math
import numpy as np
import scipy.optimize
from scipy import integrate
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt
import scipy

def dL(z,omega_m,omega_l,H0,c):
    def integrand(x):
        return 1. / np.sqrt((1.+x)**2.*(1.+x*omega_m)-x*(2.+x)*omega_l)

    if omega_m + omega_l ==1:
        i=integrate.quad(integrand,0.,z)
        return c*(1.+z)/H0*i[0]
    elif omega_m + omega_l > 1:
        i=integrate.quad(integrand,0.,z)
        return  c*(1.+z)/(H0*np.sqrt(abs(1.-omega_m-omega_l)))*np.sin(np.sqrt(abs(1.-omega_m-omega_l))*i[0])
    else:
        i=integrate.quad(integrand,0.,z)
        return  c*(1.+z)/(H0*np.sqrt(abs(1.-omega_m-omega_l)))*np.sinh(np.sqrt(abs(1.-omega_m-omega_l))*i[0])

dL = np.vectorize(dL, excluded=("omega_m","omega_l","H0","c"))

def m_model(z,omega_m,omega_l,H0,MB,sth):
    return MB + 5.0*np.log(dL(z,omega_m,omega_l,H0,sth))+25.0

def m_model_00(x):
    return m_model(x,0.,0.,0.7,-1.8,1.)

print(m_model_00(0.4))

投稿

編集

  • 回答の評価を上げる

    以下のような回答は評価を上げましょう

    • 正しい回答
    • わかりやすい回答
    • ためになる回答

    評価が高い回答ほどページの上位に表示されます。

  • 回答の評価を下げる

    下記のような回答は推奨されていません。

    • 間違っている回答
    • 質問の回答になっていない投稿
    • スパムや攻撃的な表現を用いた投稿

    評価を下げる際はその理由を明確に伝え、適切な回答に修正してもらいましょう。

15分調べてもわからないことは、teratailで質問しよう!

  • ただいまの回答率 91.02%
  • 質問をまとめることで、思考を整理して素早く解決
  • テンプレート機能で、簡単に質問をまとめられる

関連した質問

同じタグがついた質問を見る

  • Python

    5540questions

    Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

  • NumPy

    320questions

    NumPyはPythonのプログラミング言語の科学的と数学的なコンピューティングに関する拡張モジュールです。