Freiflächen in 3d plotten mit malplotlib

[stefanO]

Hallo,

ich möchte 3-dimensionale Formen wie Kegel, Zylinder usw. mit Matplotlib darstellen. Im konkreten Bsp. geht es um die MohrCoulomb yield function in 3D. Diese ist eine Pyramide mit einer unregelmäßigen hexagonalen Grundfläche - siehe Bild

Diese Yield Surface ist durch 6 Fläche begrenzt und Ich habe nun versucht die formeln in python zu definieren und dann mit matplotlib darzustellen - leider kommt hier absolut nicht das heraus was ich haben möchte.....

Ich schaffe es auch nicht den Kegel so darzustellen dass er die entlang der hydrostatischen Achse im dreidimensionalen Raum liegt.

leider finde ich online nur sehr wenig zu diesem spezifischen Thema, deshalb hoffe ich das mir hier jemand weiterhelfen kann.
Vielen Dank schon mal im Voraus!

Mohr Coulomb Yield Surface:





der code den ich für Mohr Coulomb geschrieben habe lautet wie folgt:

from mpl_toolkits.mplot3d import axes3d
from matplotlib.collections import PolyCollection
from mpl_toolkits.mplot3d.art3d import Poly3DCollection as p3dc
import matplotlib.pyplot as plt
import numpy as np

phi = 30/180*np.pi
c = 0

x = np.linspace(0,1000,5)
y = z = x*((1-np.sin(phi))/(1+np.sin(phi))) - (2*c*np.cos(phi))/(1+np.sin(phi))


#print(x,y,z)

def f1(x,z):
return (x-z)-(x+z)*np.sin(phi)-2*c*np.cos(phi)

X1, Z1 = np.meshgrid(x, z)
Y1 = f1(X1, Z1)

def f2(x,y):
return (x-y)-(x+y)*np.sin(phi)-2*c*np.cos(phi)

X2, Y2 = np.meshgrid(x, y)
Z2 = f2(X2,Y2)

def f3(x,y):
return (y-x)-(y+x)*np.sin(phi)-2*c*np.cos(phi)

X3, Y3 = np.meshgrid(x, y)
Z3 = f3(X3,Y3)

def f4(y,z):
return (y-z)-(y+z)*np.sin(phi)-2*c*np.cos(phi)

Y4, Z4 = np.meshgrid(y,z)
X4 = f4(Y4,Z4)

def f5(y,z):
return (z-y)-(z+y)*np.sin(phi)-2*c*np.cos(phi)

Y5, Z5 = np.meshgrid(y,z)
X5 = f5(Y5, Z5)

def f6(x,z):
return (z-x)-(z+x)*np.sin(phi)-2*c*np.cos(phi)

X6, Z6 = np.meshgrid(x,z)
Y6 = f5(X6, Z6)


fig = plt.figure()
ax = plt.axes(projection='3d')
#ax.contour3D(X, Y, Z, 30, cmap='binary')
ax.plot_surface(X1,Y1,Z1)
ax.plot_wireframe(X2,Y2,Z2)
ax.plot_wireframe(X3,Y3,Z3)
ax.plot_wireframe(X4,Y4,Z4)
ax.plot_wireframe(X5,Y5,Z5)
ax.plot_wireframe(X6,Y6,Z6)
ax.set_xlabel('Sigma 1')
ax.set_ylabel('Sigma 3')
ax.set_zlabel('Sigma 3');

plt.show()

[Sirius3]

Deine ganzen f-Funktionen sind identisch, das kann man also zusammenfassen:

Code: Alles auswählen

def f1(x,z):
    return (x-z)-(x+z)*np.sin(phi)-2*c*np.cos(phi)

X1, Z1 = np.meshgrid(x, z)
Y1 = f1(X1, Z1)
X2, Y2 = np.meshgrid(x, y)
Z2 = f1(X2, Y2)
X3, Y3 = X2, Y2
Z3 = f1(Y3, X3)
Y4, Z4 = np.meshgrid(y, z)
X4 = f1(Y4, Z4)
Y5, Z5 = Y4, Z4
X5 = f1(Z5, Y5)
X6, Z6 = X1, Z1
Y6 = f1(Z6, X6)

Ich habe jetzt nicht alles nachgerechnet, aber Deine Formel scheint nicht zu stimmen. Du hast ja nur 6 Ebenengleichungen, statt Dreiecke zu berechnen.

[stefanO]

Hallo,

danke für deine Antwort. Die Formeln habe ich aus einem Paper zu dem Thema. Die mathematische Herleitung habe ich noch nicht überprüft.

Aber um meine Frage vielleicht etwas allgemeiner zu formulieren:

Es gibt noch weitere Bruchmodelle welche die Form eines Zylinders oder Kegel haben welcher um eine beliebige Achse im Raum gedreht ist. Bis jetzt habe ich es geschafft den Zylinder bzw. Kegel entlang einer Hauptachse zu modellieren, bei der Drehung stehe ich allerdings an.

Aussehen sollte es am Ende so:



bzw. so:




Das ist mein Code für den Zylinder:

import math as math
import matplotlib.pylab as plt
import numpy as np

from mpl_toolkits.mplot3d import Axes3D

def f(x, y):
return np.sqrt(x ** 2 + y ** 2)

x = np.linspace(-10, 10, 100)
y = np.linspace(-10, 10, 100)

X, Y = np.meshgrid(x, y)
Z = f(X, Y)

fig = plt.figure()
ax = plt.axes(projection='3d')
#ax.contour3D(X, Y, Z, 30, cmap='binary')
ax.plot_surface(X,Y,Z)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z');

plt.show()

[einfachTobi]

`math` wird importiert, aber nicht verwendet. `as` verwendet man zum Umbenennen des Imports - das machst du hier aber gar nicht. Die Punkte kannst du einfach mit einer entsprechenden Rotationsmatrix drehen. Drei Drehungen um jeweils eine Achse durch Multiplikation der Punkte-Vektoren mit sowas sollte funktionieren:

Code: Alles auswählen

alpha = 30
rotation_x = np.array([[1, 0, 0], [0, np.cos(alpha), -np.sin(alpha)], [0, np.sin(alpha), np.cos(alpha)]])
rotation_y = np.array([[np.cos(alpha), 0, np.sin(alpha)], [0, 1, 0], [-np.sin(alpha), 0, np.cos(alpha]])

Für Z analog.