05.01.2020

<< | Home | >>

[ Lukas ]

Die Delaunay - Tetraedisierung wurde nun implementiert und funktioniert mehr oder weniger zuverlässig. Das Programm wählt nach dem Kugelkriterium die richtigen Tetraeder aus und fügt diese dem Viewport hinzu.

Kugelkriterium Es wird über eine zufällige Punktmenge iteriert und jede mögliche Kombination für vier verschiedene Kerne überprüft. Vier Punkte bestimmen (mathematisch bewiesen!) eine Kugel eindeutig, somit testen wir für alle Kombinationen (äquivalent zum zweidimensionalen Fall), ob sich weitere Kerne innerhalb der aufgespannten Kugel befinden. Ist dies der Fall, handelt es sich bei der zu überprüfenden Kombination NICHT um ein Tetraeder der Delaunay - Tetraedisierung (Begriff von uns eingeführt) - andernfalls handelt es sich um eines der gesuchten Volumina.

Außerdem wurde der lästige Fehler „LinalgError Singular Matrix“ durch Ausschluss der Determinante = 0 ausgemerzt.

import bpy
from random import random, randint
from math import *
import numpy as np
 
 
 
#   KLASSEN
 
class Sphere(object):
    pos = np.array([0,0,0])
    radius = 0
    epsilon = 0.001
 
    def __init__(self, pos, radius):
        self.pos = pos
        self.radius = radius
 
    def punktInKugel(self, p):
        if (abstandPunkte(self.pos, p) < self.radius - self.epsilon):
            print(abstandPunkte(self.pos, p), self.radius)
            return True
        else:
            return False
 
class Tetraeder(object):
    verts = []
    faces = []
    sp = np.array([0,0,0])
    name = "tetraeder"
    umkugel = Sphere(np.array([0,0,0]), 0)
 
    def __init__(self, verts, name):
        self.name = name
        self.verts = verts
        self.faces = [(0,1,2),(1,2,3),(2,3,0),(3,0,1)]
        if (self.schwerpunkt()[0] != None):
            self.sp = self.schwerpunkt()
            self.umkugel = Sphere(self.sp,abstandPunkte(self.sp,self.verts[0]))
        else:
            self.sp = np.array([None,None,None])
 
    def schwerpunkt(self):
        v1 = np.subtract(self.verts[1],self.verts[0])
        v2 = np.subtract(self.verts[2],self.verts[0])
        v3 = np.subtract(self.verts[3],self.verts[0])
 
        sp = tetraSchwerpunkt(self.verts[0],v1,v2,v3)
        return sp
 
    def selfRender(self, ker):
        frei = True
        for m in range(0, len(ker)):
            if (self.umkugel.punktInKugel(ker[m]) == True):
                frei = False
                print("shit", ker[m])
        if (frei == True):
            generateObj("Tetraeder", t.verts, t.faces, (0,0,0))
            #zeichnePunkt(t.sp, "Schwerpunkt")
 
class Gerade(object):
    ov = (0,0,0)
    rv = (0,0,0)
 
    def __init__(self,ov,rv):
        self.ov = ov
        self.rv = rv
 
    def punktAn(self,k):
        punkt = np.add(self.ov,(self.rv[0]*k,self.rv[1]*k,self.rv[2]*k))
        return punkt
 
 
 
#   FUNKTIONEN
 
def schnittTest(g1,g2):
    a = [[g2.rv[0],-g1.rv[0]],[g2.rv[1],-g1.rv[1]]]
    b = [[g1.ov[0]-g2.ov[0]],[g1.ov[1]-g2.ov[1]]]
 
    print("Det :",np.linalg.det(a))
    if (np.linalg.det(a) != 0):
        return np.linalg.solve(a,b)[0][0]
    else: return None
 
def schnittPunkt(g1,g2):
    z = schnittTest(g1,g2)
    if (z != None):
        return g2.punktAn(z)
    else: return np.array([None,None,None])
 
def dreieckSchwerpunkt(pos,a,b):
    print("DreieckSchwerpunkt")
 
    g1 = Gerade(pos,a)
    g2 = Gerade(pos,b)
 
    c = kreuzProdukt(a,b)
 
    d = kreuzProdukt(a,c)
    e = kreuzProdukt(b,c)
 
    g3 = Gerade(g1.punktAn(.5), np.divide(d, 40))
    g4 = Gerade(g2.punktAn(.5), np.divide(e, 40))
 
    schnitt = schnittPunkt(g3,g4)
    #eichnePunkt(schnitt, "Dreieckschwerpunkt")
    return schnitt
 
def tetraSchwerpunkt(pos, v1, v2, v3):
    print("tetraSchwerpunkt")
 
    sp1 = dreieckSchwerpunkt(pos, v1, v2)
    sp2 = dreieckSchwerpunkt(pos, v1, v3)
 
    print(pos, v1, v2, v3, sp1, sp2)
 
    if (sp1[0] != None and sp2[0] != None):
        g1 = Gerade(sp1, kreuzProdukt(v1, v2))
        g2 = Gerade(sp2, kreuzProdukt(v1, v3))
 
        schnitt = schnittPunkt(g1, g2)
        print("Tetraschwerpunkt :",schnitt)
        return schnitt
    else: return np.array([None,None,None])
 
def abstandPunkte(p1, p2):
    a = sqrt(((p1[0]-p2[0])**2) + ((p1[1]-p2[1])**2) + ((p1[2]-p2[2])**2))
    return a
 
def zeichnePunkt(pos, name = "Punkt"):
    verts = [(0.1,0.1,0.1),(-0.1,0.1,0.1),(-0.1,-0.1,0.1),(0.1,-0.1,0.1)
            ,(0.1,0.1,-0.1),(-0.1,0.1,-0.1),(-0.1,-0.1,-0.1),(0.1,-0.1,-0.1)]
    faces = [(0,1,2,3),(4,5,6,7),(0,1,4,5),(1,2,4,5),(2,3,6,7),(3,0,7,4)]
    generateObj(name,verts,faces,pos)
 
def kreuzProdukt(a,b):
    c = (a[1]*b[2]-a[2]*b[1] , a[2]*b[0]-a[0]*b[2] , a[0]*b[1]-a[1]*b[0])
    return c
 
def genKerne(anzahl, x_max, y_max, z_max):
    ker = []
 
    for n in range(0, anzahl):
        x = randint(-x_max, x_max)*((random()/(1+random())))
        y = randint(-y_max, y_max)*((random()/(1+random())))
        z = randint(-z_max, z_max)*((random()/(1+random())))
        if (len(ker) != 0):
            neu = True
            for m in range(0, len(ker)):
                if (x == ker[m][0] and y == ker[m][1] and z == ker[m][2]):
                    neu = False
            if (neu == True):
                ker.append([x,y,z])
                zeichnePunkt(ker[n])
            else:
                n = n-1
        else:
            ker.append([x,y,z])
            zeichnePunkt(ker[n])
    return ker
 
def generateObj(name, verts, faces, loc):
    mymesh = bpy.data.meshes.new(name)
    myobject = bpy.data.objects.new(name,mymesh)
 
    myobject.location = loc #(0,0,0) #bpy.context.scene.cursor.location
    bpy.context.collection.objects.link(myobject)
 
    mymesh.from_pydata(verts,[],faces)
    mymesh.update(calc_edges=True)
 
 
 
 
 
#   CODE
 
x_max = 10
y_max = 10
z_max = 10
ker = genKerne(20, x_max, y_max, z_max)
 
for i in range(0, len(ker)):
    for j in range(i+1, len(ker)):
        for k in range(j+1, len(ker)):
            for l in range(k+1, len(ker)):
                print(str(i),str(j),str(k),str(l))
 
                t = Tetraeder([ker[i], ker[j], ker[k], ker[l]],"Tetraeder")
                t.selfRender(ker)