diff options
| author | Nao Pross <np@0hm.ch> | 2022-04-17 23:55:32 +0200 |
|---|---|---|
| committer | Nao Pross <np@0hm.ch> | 2022-04-17 23:55:32 +0200 |
| commit | 85d45cb2a89143fee1f4ac0e188c574b963645c1 (patch) | |
| tree | 4d0857b49e23e664e05a4c2395fa4afc96099cc7 | |
| download | FourierOnS2-85d45cb2a89143fee1f4ac0e188c574b963645c1.tar.gz FourierOnS2-85d45cb2a89143fee1f4ac0e188c574b963645c1.zip | |
Add existing code
Diffstat (limited to '')
| -rw-r--r-- | manim/harmonics.py | 67 | ||||
| -rw-r--r-- | sph_harm.py | 71 | ||||
| -rw-r--r-- | sph_harm_sandbox.py | 32 |
3 files changed, 170 insertions, 0 deletions
diff --git a/manim/harmonics.py b/manim/harmonics.py new file mode 100644 index 0000000..a908e25 --- /dev/null +++ b/manim/harmonics.py @@ -0,0 +1,67 @@ +from manim import * +from manim.mobject.opengl.opengl_surface import OpenGLSurface +from manim.mobject.opengl.opengl_three_dimensions import OpenGLSurfaceMesh + + +import numpy as np +from scipy.special import sph_harm + +# configure style +config.background_color = '#202020' +config.tex_template.add_to_preamble( + r"\usepackage[p,osf]{scholax}" + r"\usepackage{amsmath}" + r"\usepackage[scaled=1.075,ncf,vvarbb]{newtxmath}" +) + +class FourierOnS1(ThreeDScene): + @staticmethod + def harmonic_cos(m, n, theta, phi): + # get only the cos() part + r = np.real(sph_harm(m, n, theta, phi)) * 5 # Y_n^m + + # coordinates transformation + x = r * np.cos(theta) * np.sin(phi) + y = r * np.sin(phi) * np.sin(theta) + z = r * np.cos(phi) + + return [x, y, z] + + @staticmethod + def harmonic_surf(m, n): + surf = Surface(lambda u, v: FourierOnS1.harmonic_cos(m, n, u, v), + u_range=[0, 2 * PI], v_range=[-PI, PI], resolution=[50, 50], + fill_opacity=.8, checkerboard_colors=[BLUE_D, BLUE_E]) + + return surf + + def construct(self): + axes = ThreeDAxes(x_range=[-1, 1], y_range=[-1, 1]) + self.add(axes) + + self.set_camera_orientation(theta=30 * DEGREES, phi=75 * DEGREES) + self.begin_ambient_camera_rotation(rate=.1) + + last_surf = None + last_label = None + for n in range(0, 5): + for m in range(0, n): + surf = FourierOnS1.harmonic_surf(m, n) + label = MathTex(r"Y_{%d}^{%d}(\theta, \phi)" % (n, m)) + + if last_surf is not None: + self.play(ReplacementTransform(last_surf, surf)) + + if last_label is not None: + self.remove(last_label) + + self.play(Create(surf)) + self.add_fixed_in_frame_mobjects(label) + label.to_corner(UL) + + last_surf = surf + last_label = label + + self.wait() + + self.wait(5) diff --git a/sph_harm.py b/sph_harm.py new file mode 100644 index 0000000..1d6227b --- /dev/null +++ b/sph_harm.py @@ -0,0 +1,71 @@ +# -*- coding: utf-8 -*-
+"""
+date of creation: Fri Mar 25 16:20:46 2022
+-------------------------------------------------
+@author : Manuel Cattaneo
+@contact: cattaneo.manuel@hotmail.com
+-------------------------------------------------
+
+Description:
+ Program used to play around with the spherical harmonics.
+
+Modules:
+"""
+import numpy as np
+# plot stuff
+import matplotlib.pyplot as plt
+import mpl_toolkits.mplot3d.axes3d as axes3d # to plot 3d
+import matplotlib.gridspec as gs # for grid layout
+# spherical harmonics
+from scipy.special import sph_harm
+"""
+------------------------------------------------------------------------------
+"""
+
+if __name__ == '__main__':
+ # grid spec
+ row = 3
+ col = 2*row-1
+ grid_spec = gs.GridSpec(row,col)
+ fig = plt.figure()
+
+ # creating spherical domain
+ N = 50
+ theta, phi = np.meshgrid(np.linspace(0, 2*np.pi, N), # theta in [0, 2pi]
+ np.linspace(0, np.pi, N)) # phi in [0, pi]
+
+ import matplotlib
+ # m < n
+ for y_grid in range(0,row):
+ for i,x_grid in enumerate(range(col//2-y_grid, col-col//2+y_grid)):
+ # real(...) to obtain only the cos(.) part. This is done in order
+ # to have the same functions as the ones at p.518 chapter 12.2.
+ m = y_grid
+ n = m+i+1
+ r = np.real( sph_harm(m, n, theta, phi) ) # Y_n^m
+ # coordinates transformation
+ x = r*np.cos(theta)*np.sin(phi)
+ y = r*np.sin(phi) *np.sin(theta)
+ z = r *np.cos(phi)
+
+ # adding subplot
+ ax = fig.add_subplot(grid_spec[y_grid,x_grid], projection='3d')
+ ax.set_axis_off()
+ ax.set_title(rf'$Y_{n}^{m}(\theta, \phi)$')
+
+ # colors
+ color_dim = r
+ color_min, color_max = color_dim.min(), color_dim.max()
+ norm = matplotlib.colors.Normalize(color_min,color_max)
+ m = plt.cm.ScalarMappable(norm, cmap='jet')
+ m.set_array([])
+ fc = m.to_rgba(color_dim)
+
+ # plot
+ plot = ax.plot_surface(x, y, z, facecolors = fc,
+ rstride=1, cstride=1,
+ # cmap=plt.get_cmap('jet'),
+ linewidth=0, antialiased=False, alpha=0.5)
+
+
+ plt.show()
\ No newline at end of file diff --git a/sph_harm_sandbox.py b/sph_harm_sandbox.py new file mode 100644 index 0000000..51828bc --- /dev/null +++ b/sph_harm_sandbox.py @@ -0,0 +1,32 @@ +# -*- coding: utf-8 -*-
+"""
+date of creation: Fri Mar 10 03:20:46 2022
+-------------------------------------------------
+@author : Manuel Cattaneo
+@contact: cattaneo.manuel@hotmail.com
+-------------------------------------------------
+
+Description:
+ Programma per capire come cazzo funziona sta SHExpandDHC(.) di merda.
+
+Modules:
+"""
+from scipy.special import sph_harm
+from pyshtools.expand import SHExpandDHC, SHExpandDH
+
+import numpy as np
+"""
+------------------------------------------------------------------------------
+"""
+
+_theta = np.linspace(0, np.pi , 1001)[1:] # theta in ]0, pi]
+_phi = np.linspace(0, 2*np.pi, 1001)[:-1] # phi in [0, 2pi[
+
+theta, phi = np.meshgrid(_theta, _phi)
+
+# n >! m
+m = 2
+n = 4
+y = sph_harm(m, n, theta=phi, phi=theta) # angles inverted due to convention
+
+cilm = np.abs(SHExpandDHC(y))
\ No newline at end of file |