Add existing code

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