icon 2018/12/03 [Sympy]3次関数と解グラフ化(Matplotlib編)
高次方程式の一般式でなく、3次関数(ax^3 + bx^2 + cx + d)でグラフ化してみる。
*定数は任意に入力
a=0.1, b=0.1, c=-20, d=-30
x = [-20, 20]

sympy のインストール


pip3 install sympy

3次方程式と解


def function(a, b, c, d):
    x = Symbol('x')
    return a*x**3 + b*x**2 + c*x + d
 
def answer(a, b, c, d):
    return solve(function(a, b, c, d))
Symbol() で方程式の変数を作成し、solve で方程式の解を計算する

実数解に虚数部が含まれる


解が実数解でも虚数部(誤差)がでるため、matplotlib の座標として扱えない。
[-10.4646576744031 + 0.e-20*I, -0.0999101526839223 + 0.e-19*I, 9.56456782708699 + 0.e-20*I]

as_coeff_Add() で配列から実数部分だけ抜き出す(解の誤差を無視)
*as_coeff_Add() は、SymPy Modules Reference参照
def convertFloat(value):
    value = (S(value).as_coeff_Add())
    if (type(value[0]) == Float):
        return value[0]

Y座標の作成


式を取得して、subs() で xに値を代入する
def plots(a, b, c, d, min_x, max_x, step):
    x = np.arange(min_x, max_x, step)
    y = [function(a, b, c, d).subs(Symbol('x'), value) for value in x]
    return x, y

サンプル


- 定数a, b, c, d を標準入力
- グラフ用のX軸 max, min を標準入力
import matplotlib.pyplot as plt
import numpy as np
from sympy import *
import sys
 
def convertFloat(value):
    value = (S(value).as_coeff_Add())
    if (type(value[0]) == Float):
        return value[0]
        
def function(a, b, c, d):
    x = Symbol('x')
    return a*x**3 + b*x**2 + c*x + d
 
def answer(a, b, c, d):
    return solve(function(a, b, c, d))

def plots(a, b, c, d, min_x, max_x, step):
    x = np.arange(min_x, max_x, step)
    y = [function(a, b, c, d).subs(Symbol('x'), value) for value in x]
    return x, y
 
step = 0.1
a, b, c, d = map(float, input('Please input a, b, c, d.\n').split())
min_x, max_x = map(float, input('Please input min x, max x.\n').split())
 
if (min_x >= max_x): sys.exit('invalid min x max x!')
 
#fig, ax = plt.subplots()
for value in answer(a, b, c, d):
    value = convertFloat(value)
    #ax.annotate(value, (value, 0))
    plt.scatter([value], [0], label = 'x = %s' % value)

x, y = plots(a, b, c, d, min_x, max_x, step)
plt.plot(x, y)
plt.grid(color='gray')
plt.title("%sx^3 + %sx^2 + %sx + %s" % (a, b, c, d))
plt.legend()
plt.show()
a = 0 にすれば2次関数の解にも使える

このサイトについて

HTML5 & CSS3化しつつあるので、現在IEには対応してません。
できれば、Google Chromeやら Apple SafariのWebKit系をお勧めします。

DBからプログラムまで一応全て自作なので、バグってたらすいません。
実験でFlash版(Flex版)を先に作りましたが、ちょっと停止してます。

プロフィール

新宿近辺でSE & プログラマーしてます。
Webアプリの開発・設計とか、最近はiPhoneとか奮闘してます。
デザインはさっぱりです。

音楽は、昔からCubase打ち込み人間で、そっちの方が経歴は長いですが、最近はやる暇がないです。。。

今は、Gon's Privates ってバンドのキーボードやってます。
単発的に、なんちゃってジャズ系のライブもやってます。

名古屋生まれなのでドラゴンズ好きです。

Info && SNS

Gmail

 yohei.yoshikawa@gmail.com

Twitter

 http://twitter.com/yoo_yoo_yoo

あんまつぶやきませんが、一応技術系メインで使ってます。情報交換はこちらへ

FaceBook

 http://www.facebook.com/#!/profile.php?id=1439130626

海外の知り合いがいないので閑散としてます。

mixi

 http://mixi.jp/show_profile.pl?id=230072

音楽仲間とかはこっちメインでやってます。興味があればこちらへ