Mesh FunctionsĀ¶
MeshFunction
arithmeticsĀ¶
Given several MeshFunction
s (e.g. exact solution and discrete solution, i.e. a GlobalFunction
and a FEFunction
) we want to be able to compute the sum, difference, product, ... of these MeshFunction
s as a new MeshFunction
. This requires to implement the corresponding arithmetics for MeshFunction
s.
- Overload arithmetic operations like
__add__
,__sub__
,__mul__
,__div__
,__invert__
,__neg__
,__pos__
forMeshFunction
s (including tests, demo notebook) - Use them to compute the difference between exact and discrete solution
- Provide a gradient operator for
MeshFunction
s (non-implemented for generalMeshFunction
s but implemented forFEFunction
s)
InĀ [1]:
import import_hack
from methodsnm.fes import *
from methodsnm.forms import *
from methodsnm.formint import *
from methodsnm.mesh_1d import *
from methodsnm.mesh_2d import *
from methodsnm.visualize import *
import math
from math import pi, cos, sin
from scipy.sparse.linalg import spsolve
InĀ [2]:
mesh = Mesh1D((0,1),10)
mesh2d = StructuredRectangleMesh(5, 5)
Rechnen mit Mesh FunctionsĀ¶
InĀ [3]:
h=GlobalFunction(lambda x:x,mesh)
f2d=GlobalFunction(lambda x:2*x[0],mesh2d)
c2d=ConstantFunction(3,mesh2d)
c=ConstantFunction(3,mesh)
c2=ConstantFunction(2,mesh)
g=c+h*h-c/c2
f=(g)/g
g2d=f2d-c2d
Die Operatoren, die die Funktionen miteinander verbinden kƶnnen mit term
angezeigt werden.
InĀ [4]:
DrawFunction1D(g)
print(g.term())
DrawFunction2D(g2d)
print(g2d.term())
DrawFunction1D(f)
print(f.term())
Um die Funktion zu bekommen die sich hinter einem Bustaben verbirgt kann get_fct
benutzt werden
InĀ [5]:
DrawFunction1D(f.get_fct("f0"))
DrawFunction1D(f.get_fct("f1"))
des weiteren kƶnnen - und ~ (1/f) benutzt werden
InĀ [6]:
DrawFunction1D(-g)
print((-g).term())
DrawFunction2D(~g2d)
print(((~g2d)).term())
Mit gradient
kann der Gradient einer Funktion berechnet werden
InĀ [7]:
print(f.gradient().evaluate(np.array([0.5]),mesh.trafo(0)))
print(f2d.gradient().evaluate(np.array([0.5,0.5]),mesh2d.trafo(0)))
Unterschied zwichen $u_h$ und $u_{\text{exact}}$Ā¶
InĀ [8]:
def uhue1d(n):
mesh = Mesh1D((0,1),n)
fes = P1_Segments_Space(mesh)
ue=GlobalFunction(lambda x: math.cos(x*math.pi), mesh = mesh)
lf = LinearForm(fes)
uef = GlobalFunction(lambda x:math.cos(x*math.pi)+math.cos(x*math.pi)*math.pi**2,mesh = mesh)
lf += SourceIntegral(uef)
lf.assemble()
blf = BilinearForm(fes)
c = GlobalFunction(lambda x:1, mesh = mesh)
blf += LaplaceIntegral(c)
blf += MassIntegral(c)
blf.assemble()
uh = FEFunction(fes)
uh.vector = spsolve(blf.matrix, lf.vector)
return (uh,ue)
def uhue2d(n):
mesh = StructuredRectangleMesh(n, n)
fes = P1_Triangle_Space(mesh)
ue = GlobalFunction(lambda x: cos(pi*x[0]), mesh = mesh)
blf = BilinearForm(fes)
c = GlobalFunction(lambda x: 1, mesh = mesh)
blf += LaplaceIntegral(c)
blf += MassIntegral(c)
blf.assemble()
lf = LinearForm(fes)
f = GlobalFunction(lambda x: (1 + pi**2) * cos(pi*x[0]), mesh = mesh)
lf += SourceIntegral(f)
lf.assemble()
uh = FEFunction(fes)
uh.vector = spsolve(blf.matrix, lf.vector)
return uh,ue,mesh
InĀ [9]:
uh,ue = uhue1d(30)
dif=uh-ue
DrawFunction1D(ue)
DrawFunction1D(uh)
DrawFunction1D(dif)
InĀ [10]:
uh,ue,_ = uhue2d(3)
dif=ue-uh
DrawFunction2D(ue)
DrawFunction2D(uh)
DrawFunction2D(dif)
DrawFunction2D(ue,contour= True)
DrawFunction2D(uh,contour= True)
DrawFunction2D(dif,contour= True)