forked from RoberAgro/nurbspy
-
Notifications
You must be signed in to change notification settings - Fork 0
/
demo_nurbs_surface_revolution.py
118 lines (85 loc) · 4.85 KB
/
demo_nurbs_surface_revolution.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
""" Example showing the creation of an extrusion NURBS surface """
# -------------------------------------------------------------------------------------------------------------------- #
# Importing packages
# -------------------------------------------------------------------------------------------------------------------- #
import numpy as np
import nurbspy as nrb
import matplotlib.pyplot as plt
# -------------------------------------------------------------------------------------------------------------------- #
# Revolution surface example
# -------------------------------------------------------------------------------------------------------------------- #
# Define the array of control points
P = np.zeros((3,5))
P[:, 0] = [0.20, 0.00, 0.00]
P[:, 1] = [0.50, 0.00, 0.25]
P[:, 2] = [0.55, 0.00, 0.50]
P[:, 3] = [0.45, 0.00, 0.75]
P[:, 4] = [0.30, 0.00, 1.00]
# Define the array of control point weights
W = np.asarray([1, 2, 3, 2 ,1])
# Create the generatrix NURBS curve
nurbsGeneratrix = nrb.NurbsCurve(control_points=P, weights=W)
# Set the a point to define the axis of revolution
axis_point = np.asarray([0.0, 0.0, 0.0])
# Set a direction to define the axis of revolution (needs not be unitary)
axis_direction = np.asarray([0.2, -0.2, 1.0])
# Set the revolution angle
theta_start, theta_end = 0.00, 2*np.pi
# Create and plot the NURBS surface
nurbsRevolutionSurface = nrb.NurbsSurfaceRevolution(nurbsGeneratrix, axis_point, axis_direction, theta_start, theta_end).NurbsSurface
nurbsRevolutionSurface.plot(surface=True, control_points=True, boundary=True, normals=False)
# -------------------------------------------------------------------------------------------------------------------- #
# Revolution surface example: A sphere
# -------------------------------------------------------------------------------------------------------------------- #
# Create the generatrix NURBS curve
O = np.asarray([0.00, 0.00, 0.00]) # Circle center
X = np.asarray([1.00, 0.00, 0.00]) # Abscissa direction (negative to see normals pointing outwards)
Y = np.asarray([0.00, 0.00, 1.00]) # Ordinate direction
R = 0.5
theta_start, theta_end = np.pi/2, 3/2*np.pi
nurbsGeneratrix = nrb.CircularArc(O, X, Y, R, theta_start, theta_end).NurbsCurve
# Set the a point to define the axis of revolution
axis_point = np.asarray([0.0, 0.0, 0.0])
# Set a direction to define the axis of revolution (needs not be unitary)
axis_direction = np.asarray([0.0, 0.0, 1.0])
# Set the revolution angle
theta_start, theta_end = np.pi/2, 2*np.pi
# Create and plot the NURBS surface
nurbsRevolutionSurface = nrb.NurbsSurfaceRevolution(nurbsGeneratrix, axis_point, axis_direction, theta_start, theta_end).NurbsSurface
fig, ax = nurbsRevolutionSurface.plot(surface=True, control_points=False, boundary=False, normals=False)
# Plot isoparametric curves
nurbsRevolutionSurface.plot_isocurve_u(fig, ax, np.linspace(0, 1, 10), linewidth=0.75)
nurbsRevolutionSurface.plot_isocurve_v(fig, ax, np.linspace(0, 1, 10), linewidth=0.75)
# Validate the computation of the mean and gaussian curvatures
Nu, Nv, h = 50, 50, 1e-6
u = np.linspace(0+h, 1-h, Nu) # Small offset to avoid the poles
v = np.linspace(0+h, 1-h, Nv) # Small offset to avoid the poles
[uu, vv] = np.meshgrid(u, v, indexing='ij')
u, v = uu.flatten(), vv.flatten()
K_mean, K_gauss = nurbsRevolutionSurface.get_curvature(u, v)
error_mean = np.sum(1/(Nu*Nv) * (K_mean - 1/R)**2)**(1/2)
error_gauss = np.sum(1/(Nu*Nv) * (K_gauss - 1/R**2)**2)**(1/2)
print('The two-norm error of the mean curvature computation is : ', error_mean)
print('The two-norm error of the gaussian curvature computation is : ', error_gauss)
# -------------------------------------------------------------------------------------------------------------------- #
# Revolution surface example: A torus
# -------------------------------------------------------------------------------------------------------------------- #
# Create the generatrix NURBS curve
O = np.asarray([1.00, 0.00, 0.00]) # Circle center
X = np.asarray([1.00, 0.00, 0.00]) # Abscissa direction
Y = np.asarray([0.00, 0.00, 1.00]) # Ordinate direction
R = 0.5
theta_start, theta_end = 0, 2*np.pi
nurbsGeneratrix = nrb.CircularArc(O, X, Y, R, theta_start, theta_end).NurbsCurve
# Set the a point to define the axis of revolution
axis_point = np.asarray([0.0, 0.0, 0.0])
# Set a direction to define the axis of revolution (needs not be unitary)
axis_direction = np.asarray([0.0, 0.0, 1.0])
# Set the revolution angle
theta_start, theta_end = 0.00, 2*np.pi
# Create and plot the NURBS surface
nurbsRevolutionSurface = nrb.NurbsSurfaceRevolution(nurbsGeneratrix, axis_point, axis_direction, theta_start, theta_end).NurbsSurface
nurbsRevolutionSurface.plot(surface=True, surface_color='gaussian_curvature', colorbar=True,
control_points=False, boundary=True, normals=False)
# Show the figures
plt.show()