SPH kernels#

This example shows the all the SPH kernels

 8 import matplotlib.pyplot as plt
 9 import numpy as np
10
11 import shamrock

utilities & check integral == 1

17 def compute_integ_3d(q, W):
18     if hasattr(np, "trapezoid"):
19         integrate_func = getattr(np, "trapezoid")
20     else:
21         integrate_func = getattr(np, "trapz")
22     return integrate_func(4 * np.pi * q**2 * W, q)
23
24
25 def plot_test_sph_kernel(q, f, df, W, dW, title, ax):
26     integral_result = compute_integ_3d(q, W)
27     assert np.abs(integral_result - 1) < 1e-6, (
28         "3D integration of 4\pi q^2 W(q) is not 1, kernel: " + title
29     )
30
31     ax[0].plot(q, W, label=f"$W_{{{title}}}(q)$")
32     ax[1].plot(q, dW, label=f"$dW_{{{title}}}(q)$")
33
34
35 q = np.linspace(0, 4, 1000)
/__w/Shamrock/Shamrock/doc/sphinx/examples/sph/run_show_all_sph_kernels.py:28: SyntaxWarning: invalid escape sequence '\p'
  "3D integration of 4\pi q^2 W(q) is not 1, kernel: " + title

Cubic splines

40 fig, axs = plt.subplots(1, 2, figsize=(10, 5))
41
42 f_M4 = [shamrock.math.sphkernel.M4_f(x) for x in q]
43 w3d_M4 = [shamrock.math.sphkernel.M4_W3d(x, 1) for x in q]
44 df_M4 = [shamrock.math.sphkernel.M4_df(x) for x in q]
45 dW3d_M4 = [shamrock.math.sphkernel.M4_dW3d(x, 1) for x in q]
46 plot_test_sph_kernel(q, f_M4, df_M4, w3d_M4, dW3d_M4, "M4", axs)
47
48 f_M6 = [shamrock.math.sphkernel.M6_f(x) for x in q]
49 w3d_M6 = [shamrock.math.sphkernel.M6_W3d(x, 1) for x in q]
50 df_M6 = [shamrock.math.sphkernel.M6_df(x) for x in q]
51 dW3d_M6 = [shamrock.math.sphkernel.M6_dW3d(x, 1) for x in q]
52 plot_test_sph_kernel(q, f_M6, df_M6, w3d_M6, dW3d_M6, "M6", axs)
53
54 f_M8 = [shamrock.math.sphkernel.M8_f(x) for x in q]
55 w3d_M8 = [shamrock.math.sphkernel.M8_W3d(x, 1) for x in q]
56 df_M8 = [shamrock.math.sphkernel.M8_df(x) for x in q]
57 dW3d_M8 = [shamrock.math.sphkernel.M8_dW3d(x, 1) for x in q]
58 plot_test_sph_kernel(q, f_M8, df_M8, w3d_M8, dW3d_M8, "M8", axs)
59
60 axs[0].legend()
61 axs[1].legend()
62 axs[0].set_xlabel(r"$q$")
63 axs[1].set_xlabel(r"$q$")
64 plt.show()
run show all sph kernels

Wendland kernels

70 fig, axs = plt.subplots(1, 2, figsize=(10, 5))
71
72 f_C2 = [shamrock.math.sphkernel.C2_f(x) for x in q]
73 w3d_C2 = [shamrock.math.sphkernel.C2_W3d(x, 1) for x in q]
74 df_C2 = [shamrock.math.sphkernel.C2_df(x) for x in q]
75 dW3d_C2 = [shamrock.math.sphkernel.C2_dW3d(x, 1) for x in q]
76 plot_test_sph_kernel(q, f_C2, df_C2, w3d_C2, dW3d_C2, "C2", axs)
77
78 f_C4 = [shamrock.math.sphkernel.C4_f(x) for x in q]
79 w3d_C4 = [shamrock.math.sphkernel.C4_W3d(x, 1) for x in q]
80 df_C4 = [shamrock.math.sphkernel.C4_df(x) for x in q]
81 dW3d_C4 = [shamrock.math.sphkernel.C4_dW3d(x, 1) for x in q]
82 plot_test_sph_kernel(q, f_C4, df_C4, w3d_C4, dW3d_C4, "C4", axs)
83
84 f_C6 = [shamrock.math.sphkernel.C6_f(x) for x in q]
85 w3d_C6 = [shamrock.math.sphkernel.C6_W3d(x, 1) for x in q]
86 df_C6 = [shamrock.math.sphkernel.C6_df(x) for x in q]
87 dW3d_C6 = [shamrock.math.sphkernel.C6_dW3d(x, 1) for x in q]
88 plot_test_sph_kernel(q, f_C6, df_C6, w3d_C6, dW3d_C6, "C6", axs)
89
90 axs[0].legend()
91 axs[1].legend()
92 axs[0].set_xlabel(r"$q$")
93 axs[1].set_xlabel(r"$q$")
94 plt.show()
run show all sph kernels

Double hump kernels

100 fig, axs = plt.subplots(1, 2, figsize=(10, 5))
101
102 f_M4DH = [shamrock.math.sphkernel.M4DH_f(x) for x in q]
103 w3d_M4DH = [shamrock.math.sphkernel.M4DH_W3d(x, 1) for x in q]
104 df_M4DH = [shamrock.math.sphkernel.M4DH_df(x) for x in q]
105 dW3d_M4DH = [shamrock.math.sphkernel.M4DH_dW3d(x, 1) for x in q]
106 plot_test_sph_kernel(q, f_M4DH, df_M4DH, w3d_M4DH, dW3d_M4DH, "M4DH", axs)
107
108 f_M4DH3 = [shamrock.math.sphkernel.M4DH3_f(x) for x in q]
109 w3d_M4DH3 = [shamrock.math.sphkernel.M4DH3_W3d(x, 1) for x in q]
110 df_M4DH3 = [shamrock.math.sphkernel.M4DH3_df(x) for x in q]
111 dW3d_M4DH3 = [shamrock.math.sphkernel.M4DH3_dW3d(x, 1) for x in q]
112 plot_test_sph_kernel(q, f_M4DH3, df_M4DH3, w3d_M4DH3, dW3d_M4DH3, "M4DH3", axs)
113
114 f_M4DH5 = [shamrock.math.sphkernel.M4DH5_f(x) for x in q]
115 w3d_M4DH5 = [shamrock.math.sphkernel.M4DH5_W3d(x, 1) for x in q]
116 df_M4DH5 = [shamrock.math.sphkernel.M4DH5_df(x) for x in q]
117 dW3d_M4DH5 = [shamrock.math.sphkernel.M4DH5_dW3d(x, 1) for x in q]
118 plot_test_sph_kernel(q, f_M4DH5, df_M4DH5, w3d_M4DH5, dW3d_M4DH5, "M4DH5", axs)
119
120 f_M4DH7 = [shamrock.math.sphkernel.M4DH7_f(x) for x in q]
121 w3d_M4DH7 = [shamrock.math.sphkernel.M4DH7_W3d(x, 1) for x in q]
122 df_M4DH7 = [shamrock.math.sphkernel.M4DH7_df(x) for x in q]
123 dW3d_M4DH7 = [shamrock.math.sphkernel.M4DH7_dW3d(x, 1) for x in q]
124 plot_test_sph_kernel(q, f_M4DH7, df_M4DH7, w3d_M4DH7, dW3d_M4DH7, "M4DH7", axs)
125
126 axs[0].legend()
127 axs[1].legend()
128 axs[0].set_xlabel(r"$q$")
129 axs[1].set_xlabel(r"$q$")
130 plt.show()
run show all sph kernels

shifted cubic splines kernels

136 fig, axs = plt.subplots(1, 2, figsize=(10, 5))
137
138 f_M4Shift2 = [shamrock.math.sphkernel.M4Shift2_f(x) for x in q]
139 w3d_M4Shift2 = [shamrock.math.sphkernel.M4Shift2_W3d(x, 1) for x in q]
140 df_M4Shift2 = [shamrock.math.sphkernel.M4Shift2_df(x) for x in q]
141 dW3d_M4Shift2 = [shamrock.math.sphkernel.M4Shift2_dW3d(x, 1) for x in q]
142 plot_test_sph_kernel(q, f_M4Shift2, df_M4Shift2, w3d_M4Shift2, dW3d_M4Shift2, "M4Shift2", axs)
143
144 f_M4Shift4 = [shamrock.math.sphkernel.M4Shift4_f(x) for x in q]
145 w3d_M4Shift4 = [shamrock.math.sphkernel.M4Shift4_W3d(x, 1) for x in q]
146 df_M4Shift4 = [shamrock.math.sphkernel.M4Shift4_df(x) for x in q]
147 dW3d_M4Shift4 = [shamrock.math.sphkernel.M4Shift4_dW3d(x, 1) for x in q]
148 plot_test_sph_kernel(q, f_M4Shift4, df_M4Shift4, w3d_M4Shift4, dW3d_M4Shift4, "M4Shift4", axs)
149
150 f_M4Shift8 = [shamrock.math.sphkernel.M4Shift8_f(x) for x in q]
151 w3d_M4Shift8 = [shamrock.math.sphkernel.M4Shift8_W3d(x, 1) for x in q]
152 df_M4Shift8 = [shamrock.math.sphkernel.M4Shift8_df(x) for x in q]
153 dW3d_M4Shift8 = [shamrock.math.sphkernel.M4Shift8_dW3d(x, 1) for x in q]
154 plot_test_sph_kernel(q, f_M4Shift8, df_M4Shift8, w3d_M4Shift8, dW3d_M4Shift8, "M4Shift8", axs)
155
156 f_M4Shift16 = [shamrock.math.sphkernel.M4Shift16_f(x) for x in q]
157 w3d_M4Shift16 = [shamrock.math.sphkernel.M4Shift16_W3d(x, 1) for x in q]
158 df_M4Shift16 = [shamrock.math.sphkernel.M4Shift16_df(x) for x in q]
159 dW3d_M4Shift16 = [shamrock.math.sphkernel.M4Shift16_dW3d(x, 1) for x in q]
160 plot_test_sph_kernel(q, f_M4Shift16, df_M4Shift16, w3d_M4Shift16, dW3d_M4Shift16, "M4Shift16", axs)
161
162 axs[0].legend()
163 axs[1].legend()
164 axs[0].set_xlabel(r"$q$")
165 axs[1].set_xlabel(r"$q$")
166 plt.show()
run show all sph kernels

Total running time of the script: (0 minutes 3.808 seconds)

Estimated memory usage: 159 MB

Gallery generated by Sphinx-Gallery