charged-shells/analysis/path_plot.py

import numpy as np
from charged_shells.rotational_path import PairRotationalPath, PathEnergyPlot
from charged_shells import expansion, patch_size
from charged_shells.parameters import ModelParams

Array = np.ndarray


zero_to_pi_half = np.linspace(0, np.pi/2, 100, endpoint=True)
pi_half_to_pi = np.linspace(np.pi/2, np.pi, 100, endpoint=True)

QuadPath = PairRotationalPath()
QuadPath.set_default_x_axis(zero_to_pi_half)
QuadPath.add_euler(beta1=np.pi/2, beta2=zero_to_pi_half)
QuadPath.add_euler(beta1=zero_to_pi_half[::-1], beta2=zero_to_pi_half[::-1])
QuadPath.add_euler(beta1=zero_to_pi_half)
QuadPath.add_euler(beta1=zero_to_pi_half[::-1], beta2=zero_to_pi_half)
QuadPath.add_euler(beta1=np.pi/2, beta2=zero_to_pi_half, alpha2=np.pi/2)
QuadPath.add_euler(beta1=np.pi/2, beta2=np.pi/2, alpha1=zero_to_pi_half[::-1])

DipolePath = PairRotationalPath()
DipolePath.set_default_x_axis(zero_to_pi_half)
DipolePath.add_euler(beta2=pi_half_to_pi[::-1])
DipolePath.add_euler(beta2=zero_to_pi_half[::-1])
DipolePath.add_euler(beta2=zero_to_pi_half, beta1=zero_to_pi_half)
DipolePath.add_euler(beta2=np.pi/2, beta1=np.pi/2, alpha2=zero_to_pi_half)
DipolePath.add_euler(beta2=np.pi/2, alpha2=np.pi/2, beta1=pi_half_to_pi)
DipolePath.add_euler(beta2=np.pi/2, beta1=pi_half_to_pi[::-1], alpha1=np.pi)
DipolePath.add_euler(beta2=zero_to_pi_half[::-1], beta1=pi_half_to_pi, alpha1=np.pi)
DipolePath.add_euler(beta2=zero_to_pi_half, beta1=pi_half_to_pi[::-1], alpha1=np.pi)
DipolePath.add_euler(beta2=pi_half_to_pi, beta1=zero_to_pi_half[::-1], alpha1=np.pi)

DipolePath2 = PairRotationalPath()
DipolePath2.set_default_x_axis(zero_to_pi_half)
DipolePath2.add_euler(beta2=pi_half_to_pi[::-1])
DipolePath2.add_euler(beta2=zero_to_pi_half[::-1])
DipolePath2.add_euler(beta2=zero_to_pi_half, beta1=zero_to_pi_half)
DipolePath2.add_euler(beta2=np.pi/2, beta1=np.pi/2, alpha2=zero_to_pi_half)
DipolePath2.add_euler(beta2=np.pi/2, alpha2=np.pi/2, beta1=pi_half_to_pi)
DipolePath2.add_euler(beta2=zero_to_pi_half[::-1], beta1=pi_half_to_pi[::-1])
DipolePath2.add_euler(beta2=zero_to_pi_half[::-1], beta1=np.pi)
DipolePath2.add_euler(beta2=zero_to_pi_half, beta1=pi_half_to_pi[::-1], alpha1=np.pi)
DipolePath2.add_euler(beta2=pi_half_to_pi, beta1=zero_to_pi_half[::-1], alpha1=np.pi)


def kappaR_dependence(kappaR: Array, abar: float, sigma_tilde=0.001, save=False):
    params = ModelParams(R=150, kappaR=kappaR)

    ex1 = expansion.MappedExpansionQuad(abar, params.kappaR, 0.001, l_max=30)
    ex2 = ex1.clone()

    # charge_profile = ex1.charge_value(theta=np.linspace(0, np.pi, 100), phi=0)
    # plt.plot(charge_profile.T)
    # plt.show()

    # expansion.plot_theta_profile(ex1, num=1000, theta_end=np.pi, phi=0)
    #
    # ps = patch_size.potential_patch_size(ex1, params, theta1=np.pi / 2,
    #                                      # match_expansion_axis_to_params=1
    #                                      )
    # print(ps)

    path_plot = PathEnergyPlot(ex1, ex2, DipolePath2, dist=2., params=params, match_expansion_axis_to_params=None)

    path_plot.plot(  # labels=[rf'$\theta$={180 * patch / np.pi:.1f}' for patch in np.array([0.5])],
        # norm_euler_angles={'beta2': np.pi},
        # save_as=Path("/home/andraz/ChargedShells/Figures/dipole_path_kappaR3.png")
    )

    # path_plot.plot(labels=[rf'$\kappa R$={kR}' for kR in kappaR],
    #                # norm_euler_angles={'beta2': np.pi},
    #                # save_as=Path("/home/andraz/ChargedShells/Figures/dipole_path_lambda5_norm2.png")
    #                )

    # path_plot.plot_sections(save_as=Path('/home/andraz/ChargedShells/Figures/dipole_path2.png'))



if __name__ == '__main__':

    # kappaR = np.array([1, 3, 10])
    params = ModelParams(R=150, kappaR=3)
    # ex1 = expansion.MappedExpansionQuad(0.4, params.kappaR, 0.001, l_max=20)
    # ex1 = expansion.Expansion24(sigma2=0.001, sigma4=0)
    # ex1 = expansion.Expansion(l_array=np.array([1]), coefs=expansion.rot_sym_expansion(np.array([1]), np.array([0.001])))
    # ex1 = expansion.GaussianCharges(omega_k=np.array([[0, 0], [np.pi, 0]]), lambda_k=np.array([5, 10]), sigma1=0.001, l_max=10)
    ex1 = expansion.GaussianCharges(omega_k=np.array([0, 0]), lambda_k=np.array([1, 5, 10]), sigma1=0.001, l_max=20)
    # ex1 = expansion.SphericalCap(np.array([[0, 0], [np.pi, 0]]), np.array([0.5]), 0.001, 20)
    ex2 = ex1.clone()

    # charge_profile = ex1.charge_value(theta=np.linspace(0, np.pi, 100), phi=0)
    # plt.plot(charge_profile.T)
    # plt.show()

    # expansion.plot_theta_profile(ex1, num=1000, theta_end=np.pi, phi=0)
    #
    # ps = patch_size.potential_patch_size(ex1, params, theta1=np.pi / 2,
    #                                      # match_expansion_axis_to_params=1
    #                                      )
    # print(ps)

    path_plot = PathEnergyPlot(ex1, ex2, DipolePath2, dist=2., params=params, match_expansion_axis_to_params=None)

    path_plot.plot(# labels=[rf'$\theta$={180 * patch / np.pi:.1f}' for patch in np.array([0.5])],
                   # norm_euler_angles={'beta2': np.pi},
                   # save_as=Path("/home/andraz/ChargedShells/Figures/dipole_path_kappaR3.png")
                   )

    # path_plot.plot(labels=[rf'$\kappa R$={kR}' for kR in kappaR],
    #                # norm_euler_angles={'beta2': np.pi},
    #                # save_as=Path("/home/andraz/ChargedShells/Figures/dipole_path_lambda5_norm2.png")
    #                )

    # path_plot.plot_sections(save_as=Path('/home/andraz/ChargedShells/Figures/dipole_path2.png'))