Mesh FunctionsĀ¶

MeshFunction arithmeticsĀ¶

Given several MeshFunctions (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 MeshFunctions as a new MeshFunction. This requires to implement the corresponding arithmetics for MeshFunctions.

  • Overload arithmetic operations like __add__, __sub__, __mul__, __div__, __invert__, __neg__, __pos__ for MeshFunctions (including tests, demo notebook)
  • Use them to compute the difference between exact and discrete solution
  • Provide a gradient operator for MeshFunctions (non-implemented for general MeshFunctions but implemented for FEFunctions)
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
source module for methodsNM imported.
InĀ [2]:
mesh = Mesh1D((0,1),10)
mesh2d = StructuredRectangleMesh(5, 5)

Rechnen mit Mesh FunctionsĀ¶

Implementation in meshfct.py

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())
No description has been provided for this image
(f0+(f1*f1))-(f0/f2)
No description has been provided for this image
f0-f1
No description has been provided for this image
((f0+(f1*f1))-(f0/f2))/((f0+(f1*f1))-(f0/f2))

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"))
No description has been provided for this image
No description has been provided for this image

des weiteren kƶnnen - und ~ (1/f) benutzt werden

InĀ [6]:
DrawFunction1D(-g)
print((-g).term())
DrawFunction2D(~g2d)
print(((~g2d)).term())
No description has been provided for this image
(-((f0+(f1*f1))-(f0/f2))
No description has been provided for this image
1/((f0-f1)

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)))
[0.]
[2.00000001 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)
Warning: Using NumPy Gauss rules!
Warning: Using numerical differentiation for deriv evaluation in <class 'methodsnm.fe_1d.P1_Segment_FE'> object.
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
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)
Warning: Using numerical differentiation for deriv evaluation in <class 'methodsnm.fe_2d.P1_Triangle_FE'> object.