CAS med Sympy

Sympy er et python bibliotek som er pythons svar på Maple og arbejder lidt anderledes end almindelig python som du måske kender. Det er bedst at tænke det som et separat program løsnet fra python. Du kan derfor snildt arbejde med Sympy uden at kende ret meget til python generelt.

Numeriske beregninger og simuleringer sker med andre pakker som Numpy og Scipy. Disse to kan nemt bruges direkte SAMMEN med python. Det er nemt at integrere scipy i sympy i programmer. Scipy kan fx også integrere men ikke symbolsk.

Brug sympy til CAS! Brug andet til simulering og numeriske problemstillinger.

Sympy er genialt som CAS system både til

beregninger som hurtig lommeregner (scipy, sympy)
at løse ligninger og ligningssystemer (sympy)
hurtigt at plotte funktioner og figurer (sympy)
differentiere, integrere (scipy, sympy)
løse differentialligninger analytisk (sympy)
regne med enheder. Her i eksemplet til højre importeres de "import units as u" dvs at de kan kaldes med fx u.m, u.s, u.m/u.s osv ...

Spyder: Scipy, numpy, matplotlib

Jupyter: Sympy (scipy, numpy, matplotlib)

Installation af Jupyter har sin egen side.

Som CAS system er det klart bedst at bruge Jupyter i stedet for Spyder. Spyder virker selvfølgelig også men det er jo ikke programmer vi skal skrive, men lave hurtige beregninger og der er Jupyter Notebook bare smart og fikst.

Du åbner et nyt "notebook" ved at klippe på "New --> Python 3" i jupyter

Kort version af en god header

from sympy import * # importerer hele sympy

from sympy.physics import units as u

import scipy.constants as sc # fysiske konstanter

init_printing()

%matplotlib inline

x,y,z,t = symbols('x,y,z,t')

omega,phi,A = symbols('omega,phi,A')

k, m, n, a, b = symbols('k m n a b', integer=True)

f = Function('f') # til matematik

Et par hurtige links til uddybende sider med eksempler

Solve med sympy

Integration med sympy (og scipy)

Kraft og Arbejde sympy med fysik

Differentiation, finde maxima og rødder

Konstanter i fysik

Et simpelt eksempel med stedfunktionen:

Her er der en genial side med masser af eksempler for solve, simplify, differentiation og integration.

Du skal have Jupyter installeret før du kan bruge sympy.

Klip den grønne Header over til Jupyter for at starte med CAS og sympy!

Kopier nedenstående kode (det grønne + Stedfunktionsprogrammet) ind i et Jupyter notebook og tryk på CTRL + Enter. Læg mærke til at sympy kan plotte flere grafer "inline" men kun lave en beregning ad gangen. Du kan også se at sympy fint regner med enheder men kan IKKE plotte formler med enheder.

#Header (klip det grønne over til Jupyter)

from sympy import * # importerer hele sympy

from sympy.physics import units as u

import scipy.constants as sc # fysiske konstanter

#printe flotte formler og grafer i Jupyter

init_printing()

%matplotlib inline

#Definition af variabler og symboler

x,y,z,t = symbols('x,y,z,t')

omega,phi,A = symbols('omega,phi,A')

k, m, n, a, b = symbols('k m n a b', integer=True)

f = Function('f') # til matematik

s = Function('s') # stedfunktionen

v = Function('v') # hastighedsfunktionen

#Eksempel med stedfunktion

stal = 10*u.meter + 10*u.m/u.s*t + 0.5*sc.g*u.m/(u.s**2)*t**2

vtal = diff(stal,t)

atal = diff(vtal,t)

s = 10 + 10*t + 0.5*sc.g*t**2

v = diff(s,t)

a = diff(v,t)

plot(s,(t,0,5))

plot(v,(t,0,5))

plot(a,(t,0,5))

stal.subs(t,4*u.s)

#vtal.subs(t,4*u.s)

#atal.subs(t,4*u.s)

Her ser du sympy's beregning for stedfunktionen (stal), hastighedsfunktionen (vtal) og accelerationen (atal) for t=4s

Gennemgang af sympy's vigtigste funktioner

Se denne her supergode film så kan du det hele på 28min

Sympy er et CAS system med python backend. Det kan regne med symboler, generere flotte formler og nemme hurtige grafer. Den er ikke lavet til programmering men til at løse ligninger analytisk. Derfor giver det mest mening at bruge Jupyter Notebook til Sympy. Det hele virker også i Spyder, men formlerne bliver ikke flotte.

Links til tutorials og youtube

Hurtig tutorial med eksempler til CAS med sympy
Supergod film på 28 min
En anden film om sympy som guider dig igennem funktionerne
Sympy's egen dokumntation er virkelig god
Cheat Sheet

Se et par eksempler så du får et indtryk hvad sympy kan gøre:

Sympy er en bibliotek som kan løse ligninger, integraler og differentialligninger analystisk.

Er du interesseret i numeriske metoder som euler eller runge kutta så skal du bruge python på en anden måde.

Sympy virker mere som wordmat eller maple. Fx kan man lave hurtige grafer:

Sympy kan lave "flotte" formler, og grafer. Den har en latex backend som den bruger til at generere formler i.

Se her

Man skal "deklarere" sine variabler som i linjerne

from sympy import *

%matplotlib inline

init_printing()

x,y,z,t = symbols('x,y,z,t')

f = (x**2 + 1)

g = (-x**2 +15)

p1 = plot(f,g)

Eller:

omega1, phi, A = symbols('omega phi A')

Grundlæggende ting i sympy

Import af sympy

Man kan loade sympy på to meningsflydte måder

from sympy import * (arbejde KUN med sympy)
import sympy as sp (blande Sympy med andre pakker, som numpy, math, matplotlib, scipy etc)

Begge måder giver mening men man skal beslutte sig. Bruger du 2eren og importerer sympy som sp skal alle sympy funktioner kaldes med sp.funktion. Det er træls da man skal skrive sp. hele tiden, men siden der er mange funktioner i sympy som har samme navn som i fx "numpy" eller pakken "math" giver det forvirring hvis du blander pakker.

I nedenstående eksempler arbejder jeg udelukkende med sympy og derfor har jeg brugt

from sympy import *

Dvs at jeg er fri for at skrive sp. foran alle sympy kommandoer. Se forklaringen til højre for at forstå forskellen.

Spot the difference:

import sympy as sp

from sympy import *

Problemer ? --> Restart kernel

Mange problemer i sympy kan løses ved at slette alle variabler og definere dem forfra. Det gør man ved at trykke på Kernel --> Restart

Forskellen til "normal" pythonprogrammering

I python definerer man variabler som lister, arrays, booleans og strings. Bum.

Sympy omdefinerer ALLE variabler til sine egne typer. Derfor kan man IKKE bruge python funktioner som loops, print, plot osv. direkte på sympy-objekter.

Sympy arbejder med funktioner. Funktionerne er defineret i sympy's header som objekter og skal håndteres som sympy kræver det.

En god header til sympy klip den ud og ind øverst i din Jupyter notebook!

Denne header definerer x,y,z,t,omega,phi og A som symboler k, m, n defineres som heltal samt f og g er defineret som funktioner. Se også den "grønne header" øverst i filen til stedfunktioner.

from sympy import * # importerer hele sympy

init_printing() # printer flotte formler

%matplotlib inline

import scipy.constants as sc # fysiske konstanter

x,y,z,t = symbols('x,y,z,t')

omega,phi,A = symbols('omega,phi,A')

k, m, n = symbols('k m n', integer=True)

f = Function('f')

g = Function('g')

f = 4*x**2-225

print(solve(f,x))

plot(f)

Deklaration af variabler

I pyhon skal man normalt ikke deklarere variabler undtagen lister og arrays.

Sympy skal vide om et bogstav står for en funktion eller for en variabel. Det skal derfor defineres i starten af ens Jupyter notebook.

Definition af variabler:

x,y,t,omega1,phi,A = symbols('x y t omega phi A')

Du skal gøre det i den rigtige rækkefølge :-)

En god "header" for sin fil vil se ud som følgende

from sympy import *

import scipy.constants as sc

x,y,z,t = symbols('x,y,z,t')

omega,phi,A = symbols('omega,phi,A')

k, m, n = symbols('k m n', integer=True)

f = Function('f')

g = Function('g')

x y z og t er symbol variabler som sympy regner med.

k m og n er integer tal (heltal)

Definiton af en funktion

Direkte definition

f_1 = (1/x) + (x * sin(x) - 1)/(x**2 - 1)

Ligningsdefinition

expr = Eq(y,(1/x) + (x * sin(x) - 1)/(x**2 - 1))

f_1 og expr er "det samme" de representerer funktionen. en variabel af typen "object" som sympy definerer. Disse representationer kan give lidt problemer. Se længere nede afsnit "Hvad er Eq" til at forstå hvorfor dte er svært.

f_1 er dog forskellig fra expr idet det "kun" viser højre siden på funktionen. Den kan differentieres fx direkte ved at skrive diff(f_1,x)

Se næste afsnit og se hvordan du bruger det.

Funktioner

De fleste funktioner giver sig selv.

Trigonometrisk: sin cos tan

Inverse arcus: asin acos atan acot

Hyberboliske: sinh cosh tanh coth

Invers hyberbolisk asinh acosh atanh acoth

Exponential funktionen: exp(x)

Kvadratrod: sqrt(x)

Log(a) til basis b: log(a, b)

Ln(a): log(a)
Numerisk værdi: abs(a)

Differential og Integralregning

Solve

HINT: her er det smart at bruge Eq definitionen på funktionen da det ofte er en hel ligning
y =(x-2*y)+2 man vil løse og ikke bare en funktion som (x-2*y)+2

For at komme lidt mere i dybden med selve solve funktionen klik her (intern side).

For at løse en ligning kan man bruge følgende

simplify((1/x) + (x * sin(x) - 1)/(x**2 - 1))

Eller definere funktionen først og så vise den

f_1 = Eq(y,(1/x) + (x * sin(x) - 1)/(x**2 - 1))

f_1

For at erstatte x med et tal fx 3 bruger man subs og hvis man vil have talværdien og ikke den symbolske løsning bruger man evalf(4) som angiver resultatet med 4 decimaler

f_1 = Eq(y,(1/x) + (x * sin(x) - 1)/(x**2 - 1))

simplify(f_1)

f_1.subs(x,3).evalf(8)

Forvirring! Hvorfor bruger man både f_1 og y?

Definerer man en funktion med sympy så sætter den automatisk alle led på den ene side og =0 på den anden. f_1 skal ses nu som en representation af funktionen nedenunder. med 0 og -y og hele pivtøjet.

y er en symbolsk variabel og en del af formlen

f_1 eller expr er python variabel af typen objekt og indeholder selve formlen

Meget kort vejledning

Header:

from sympy import *

x,y=symbols('x y')

definér funktionen f_1 = Eq(y,(x-2*y)+2)
solve(f_1,x) --> x = [3𝑦−2]
Beregne y for x=3 : f_1.subs(x,3) --> [𝑦 = 5 - 2𝑦]
Solve solve(f_1.subs(x,3).evalf(),y) giver (y=) 1.66666666667
solve løser ligningen mht til y fordi der står ,y til sidst. subs erstatter x med et 3tal. evalf() laver resultatet fra brøk om til decimaltal

Definition:

Og løsning med den numeriske løsning til y(3).

Skriver du funktionen direkte ser den sådan ud:

Cheet Sheet til at solve funktioner

Forenkle en funktion:
simplify((1/x) + (x * sin(x) - 1)/(x**2 - 1))
Solve ligning f_1 mht x
solve(f_1,x)
Finde rødder i et polynomium
solve(x**2+2*x-4,x)
Finde numeriske rødder med evalf(4) for 4 decimaler
solve(f_1.evalf(),x)eller
solve(N(f_1),x)
Beregne en funktion med konstanter defineret.
Husk at scipy pakken har en masse fysiske konstanter .
a=5; b=2;c=14
solve(N(a*x**2+b*x+c),x)

Integration og differentiation

HINT: lad være med at bruge Eq definitionen på funktioner, skriv dem direkte!

Her er der en side om integration med scipy og sympy med mange eksempler.

Hurtig differentiation

from sympy import *

x = symbols('x')

diff(exp(-x**2),x)

Hurtig integration

from sympy import *

x = symbols('x')

integrate(exp(-x**2),(x,-oo,oo))

integration og differentiation er straight forward.

Man skrive funktionen direkte i diff

diff((x**2) + (2*x**3 -4*x),x)

Eller definere en funktion

f_2=(x**2) + (2*x**3 -4*x)

f_3=diff(f_2)

Integration er det samme

integrate(f_3,x)

Se screenshots fra Jupyter her til højre.

Et bestemt integral løses ved

integrate(f_3,(x,a,b))

Husk de ekstra parenteser

Her en funktion som afhænger af flere variabler og som bliver differentieret 6 gange.

Du kan taste funktionen direkte i parentesen

Eller definere en funktion fx f_2 og differentiere den

Integralet af e^(-x) fra 0 til uendelig skrives på denne måde:

Nice:

Plotte funktioner i sympy (flere funktioner klik her)

HINT: Det er Mathplotlib som laver "flotte" grafer. Sympy er kun smart hvis du accepterer dens layout, farver osv .. hvis du vil lave flotte figurer hvor du styrer text, grid, farver osv så skal du ud af sympy og tilbage til almindelig python.

Husk at starte dit sheet med

from sympy import *

%matplotlib inline

init_printing()

x,y,z,t = symbols('x,y,z,t')

init_printing() funktionen connecter til printeren i dit system.

Du kan plotte to eller flere funktioner adskildt af komma og i (funktion,(x fra -6 til 6)),(anden funktion,(x+7...))

Lige som med differentiationen kan du plotte direkte eller definere en funktion her er f_3 defineret som:

plot((x**2, (x, -6, 6)), (x+7, (x, -5, 5)),title="Min funktion",xlabel="x-akse",ylabel="y-akse",axis_center=(-6,0),line_color='blue')

Hurtig plot

from sympy import *

%matplotlib inline

var('x')

plot(sin(x),(x,0, 2*pi))

Hurtig 3d - plot

from sympy import *

%matplotlib inline

var('x y')

plotting.plot3d((exp(-(x**2+y**2))),(x,-3,3),(y,-3,3))

Løsning af differentialligninger metode 1 og metode 2

Grundlæggende ting

Start med at definere din header:

from sympy import *

init_printing()

x = symbols('x')

f = Function('f')(x)

f # vil printe funktionen f(x)

f.diff(x) #denne linje tager f'(x)

f.diff(x,x) #vil være den ener f''(x)

f.diff(x,x,x) #vil være den ener f'''(x)

Metode 1 radioaktiv henfald

Vi vil løse denne differentialligning:

I sympy vil det se ud som følgende:

Her er koden:

from sympy import *

t,k = symbols('t,k')

N = Function('N')(t)

diffeq = Eq(N.diff(t) + k * N,0)

display(diffeq)

dsolve(diffeq, N)

Metode 1 definition via Eq statement.

diffeq = Eq(f.diff(x,x)-5*f,0)

Ovenstående definerer differentialligningen. Læg mærke til at ,0 til sidst betyder at ligningen sættes lige med 0

Løsningen er nu straight forward ved at bruge funktionen dsolve(diffeq,f)

Metode 2 uden Eq

Her definerer man bare ligningen direkte hvor man sætter ligningen til 0 på den ene side. Som man kan se er bette metoder lige gode!

from sympy import *

t,k = symbols('t,k')

N = Function('N')(t)

diffeq = N.diff(t) + k * N

display(diffeq)

dsolve(diffeq, N)

Harmonisk oscillator med metode 2

Her bruges Eq definitionen IKKE!

Differentialligningen til den harmoniske oscillator:

Er meget tæt relateret til opgaven ovenover, med det lille forskel at der er et + mellem de to led hvis de står på samme side. Det gør at løsningen bliver kompleks.

Se først hvordan man erstatter x med y og t med omega.

Løsningen indeholder to e-funktionen her -iωt i eksponenten. Det skal omskrives til noget der ligner det her:

Koden er her

from sympy import *

t, omega = symbols('t omega')

y = Function('y')(t)

z = Function('z')(t)

diffeq = y.diff(t,t)+omega**2*y

display(diffeq)

dsolve(diffeq,y)

Eulers formler kan bruges til at forstå sammenhængen mellem den imaginære løsning og cosinus funktionen.

Løsningen fra differentialligningen snildt kan omkrives til en cosinus.

Funfact:

Python er god til at arbejde med komplekse tal. Man kan plotte løsningen til venstre med følgende lille program som IKKE bruger sympy:

import numpy as np

import matplotlib.pyplot as plt

T=5

omega = 2*np.pi/T

c1 = 0.5

c2 = 0.5

im = 1j #man definerer komplexe tal med j!

x = np.linspace(0,10,1000)

y = c1*np.exp(-im*omega*x)+c2*np.exp(im*omega*x)

plt.plot(x,y)

Som giver løsningen:

Python skriver også følgende fordi den imaginære del af funktionen of course ikke kan plottes:

ComplexWarning: Casting complex values to real discards the imaginary part

Løsning af 2 ligninger med to ubekendte

Her skal objekterne defineres som Eq. Er løsningen entydig formaterer sympy den som typen "dictionary" ellers som liste

from sympy import *

x,y = symbols('x,y')

exp1=Eq(2*x+y-5,0)

exp2=Eq(x-y-1,0)

loesning=solve((exp1,exp2),(x,y))

print('Løsningen for x = {0} og y = {1}'.format(loesning[x],loesning[y]))

Nedenunder er en løsning for et mere kompliceret system med 4 mulige løsninger for x og y

exp1= Eq(2*x**2 + 4*y**2 - 16*x*y-4,0)

exp2= Eq(4*x + x*y,0)

Variabel type Dictionary:

mindict={'første':42,'anden':43,'tredie':'Knud'}

print(mindict['anden'])

print(mindict['tredie']

Knud

Løsning af tre ligninger med tre ubekendte

lig1= Eq(venstreside , højreside)

lig1 er IKKE en funktion i sigselv, men et navn for funktionen.

Løsning med sympy:

lig1 = Eq(x+y+z,0)

lig2 = Eq(2*x-y-z,10)

lig3 = Eq(y+2*z,5)

solve((lig1,lig2,lig3),(x,y,z))

Hvad er "Eq"

Dette afsnit er lidt nørdet og kun nødvendigt hvis man vil forstå Eq vs direkte definition af funktionen

Prøv at kikke på koden til højre.

diffeq1 = Eq(y.diff(t,t)-omega**2*y,y)

diffeq2 = Derivative(y, (t,2))-omega**2*y

Begge funktioner definerer den samme differential ligning.

diffeq1 er enVariable defineret af sympy. Det er en "funktion" og grundlæggende af typen "double" dvs at det er liste med to elementer.

diffeq1.args[0] og diffeq1.args[1]

elementerne kan kun tilgås med funktionen .args.

Elementerne i funktionerne diffeq1 og diffeq2 ser sådan ud

diffeq1[0] = -omega**2*y(t) + Derivative(y(t), (t, 2))

diffeq1[1] = y(t)

diffeq2[0] = -omega**2*y(t)

diffeq2[1] = Derivative(y(t), (t, 2))

Så man kan ikke bare arbejde med "diffeq1 og 2" som almindelige python variabler. De er mere af alt funktioner som kun sympy kan læse. Så prøv ikke at plotte dem med matplotlib eller lingende .. det failer.

du kan dog tilgå de to elementer og arbejde med dem fx sådan her:

print('diffeq1[0] = ',diffeq1.args[0])

print('diffeq1[1] = ',diffeq1.args[1])

som giver:

diffeq1[0] = -omega**2*y(t) + Derivative(y(t), (t, 2))

diffeq1[1] = y(t)

Funfact :

diffeq1 = Eq(Derivative(y, (t,2))-omega**2*y,y)

er det samme som

diffeq1 = Eq(y.diff(t,t)-omega**2*y,y)

så

diffeq2 = Derivative(y, (t,2))-omega**2*y)

diffeq3 = y.diff(t,t)-omega**2*y

er også det samme. Vi bruger diff, men som du ser spytter programmet "Derivative" ud.

from sympy import *

t, omega = symbols('t omega')

y = Function('y')(t)

diffeq1 = Eq(y.diff(t,t)-omega**2*y,y)

diffeq2 = y.diff(t,t)-omega**2*y

display("diffeq1",diffeq1)

display("diffeq2",diffeq2)

print('diffeq1[0] = ',diffeq1.args[0])

print('diffeq1[1] = ',diffeq1.args[1])

print('diffeq2[0] = ',diffeq2.args[0])

print('diffeq2[1] = ',diffeq2.args[1])

Page updated

Google Sites

Report abuse