<jupyter_start><jupyter_text>Molecule<jupyter_code>from SCFInitialGuess.utilities.dataset import Molecule

angle_HCC = 121.3/180 *np.pi
distance_CC = 1.339
distance_CH = 1.087

mol = Molecule(
    species=["C", "C", "H", "H", "H", "H"],
    positions=[
        [ 0.00000,       0.00000,        0.00000],
        [ distance_CC  ,       0.00000,        0.00000],
        [ distance_CH * np.cos(angle_HCC), distance_CH * np.sin(angle_HCC), 0.00000],
        [ distance_CH * np.cos(2*np.pi - angle_HCC), distance_CH * np.sin(2*np.pi - angle_HCC), 0.00000],
        [ distance_CC + distance_CH * np.cos(np.pi - angle_HCC), distance_CH * np.sin(np.pi - angle_HCC), 0.00000],
        [ distance_CC + distance_CH * np.cos(- np.pi + angle_HCC), distance_CH * np.sin(- np.pi + angle_HCC), 0.00000]
    ]
)
pos = np.array(mol.positions)
plt.scatter(pos[:2, 0], pos[:2, 1], label="C")
plt.scatter(pos[2:, 0], pos[2:, 1], label="H")
plt.legend()<jupyter_output><empty_output><jupyter_text>Utilities<jupyter_code>def activation(r, phi, theta, descriptor, values):
    
    radial_descriptor = descriptor.radial_descriptor
    angular_descriptor = descriptor.angular_descriptor
    
    return radial_descriptor.calculate_inverse_descriptor(
        r, 
        values[:radial_descriptor.number_of_descriptors]
    ) * angular_descriptor.calculate_inverse_descriptor(
        r,
        phi,
        theta,
        values[radial_descriptor.number_of_descriptors:]
    )<jupyter_output><empty_output><jupyter_text>Initialize Descriptors  Utilities<jupyter_code>from SCFInitialGuess.descriptors.utilities import plot_normal_model, plot_periodic_model<jupyter_output><empty_output><jupyter_text>Low Level<jupyter_code>r_model_concentrated = (
    list(np.linspace(0.2, 0.7, 11))[:-1] + \
    list(np.linspace(0.7, 2.5, 30)) + \
    list(np.linspace(2.5, 4, 11))[1:],
    [500]*10 + \
    [1000]*30 + \
    [200]*10,
)

plot_normal_model(r_model_concentrated, np.linspace(0, 5, 500))

plt.ylabel("$g(r, r_s, \eta) \quad \mathrm{ / }\quad 1 $")
plt.xlabel("$r\quad \mathrm{ / } \quad\mathrm{\AA}$")

plt.tight_layout()

#plt.savefig("/home/jo/Repos/MastersThesis/AtomicBlocks/figures/RadialModelOrigin.pdf")

phi_model_evenly = (
    np.linspace(0, 2*np.pi, 50),
    [300] * 50,
    2 * np.pi
)

plot_periodic_model(phi_model_evenly, np.linspace(0, 2 * np.pi, 500))
theta_model_evenly = (
    np.linspace(0, np.pi, 50),
    [700] * 50,
    np.pi
)

plot_periodic_model(theta_model_evenly, np.linspace(0, np.pi, 500))<jupyter_output><empty_output><jupyter_text>High Level<jupyter_code>from SCFInitialGuess.descriptors.coordinate_descriptors import \
    Gaussians, PeriodicGaussians, IndependentAngularDescriptor
from SCFInitialGuess.descriptors.cutoffs import \
    BehlerCutoff1
from SCFInitialGuess.descriptors.high_level import NonWeighted
    
    
descriptor_concentrated = NonWeighted(
    Gaussians(*r_model_concentrated),
    IndependentAngularDescriptor(
        PeriodicGaussians(*phi_model_evenly),
        PeriodicGaussians(*theta_model_evenly),
    ),
    BehlerCutoff1(5)
)
model = (
    [0.0] * 50,
    np.linspace(10, 0.1, 50)
)

descriptor_origin = NonWeighted(
    Gaussians(*model),
    IndependentAngularDescriptor(
        PeriodicGaussians(*phi_model_evenly),
        PeriodicGaussians(*theta_model_evenly),
    ),
    BehlerCutoff1(5)
)
plot_normal_model(model, np.linspace(0, 5, 200))

plt.ylabel("$g(r, r_s, \eta) \quad \mathrm{ / }\quad 1 $")
plt.xlabel("$r\quad \mathrm{ / } \quad\mathrm{\AA}$")

plt.tight_layout()

#plt.savefig("/home/jo/Repos/MastersThesis/AtomicBlocks/figures/RadialModelOrigin.pdf")

from SCFInitialGuess.descriptors.models import make_uniform

model = make_uniform(5, 5, eta_max=10, eta_min=1)

descriptor_unifrom_radial_5 = NonWeighted(
    Gaussians(*model),
    IndependentAngularDescriptor(
        PeriodicGaussians(*phi_model_evenly),
        PeriodicGaussians(*theta_model_evenly),
    ),
    BehlerCutoff1(5)
)
plot_normal_model(model, np.linspace(0, 5, 200))
from SCFInitialGuess.descriptors.models import make_uniform

model = make_uniform(50, 5, eta_max=120, eta_min=40)

descriptor_unifrom_radial_50 = NonWeighted(
    Gaussians(*model),
    IndependentAngularDescriptor(
        PeriodicGaussians(*phi_model_evenly),
        PeriodicGaussians(*theta_model_evenly),
    ),
    BehlerCutoff1(5)
)
plot_normal_model(model, np.linspace(0, 5, 200))
from SCFInitialGuess.descriptors.models import make_uniform

model = make_uniform(250, 5, eta_max=500, eta_min=100)

descriptor_unifrom_radial_250 = NonWeighted(
    Gaussians(*model),
    IndependentAngularDescriptor(
        PeriodicGaussians(*phi_model_evenly),
        PeriodicGaussians(*theta_model_evenly),
    ),
    BehlerCutoff1(5)
)
plot_normal_model(model, np.linspace(0, 5, 200))<jupyter_output><empty_output><jupyter_text>Radial<jupyter_code>def plot_radial_descriptor(r, descriptor, label):
    values = descriptor.calculate_atom_descriptor(
        0,
        mol,
        descriptor.number_of_descriptors
    )
    
    n_radial = descriptor.radial_descriptor.number_of_descriptors
    radial_values = values[:n_radial]
    
    inverse_values = descriptor.radial_descriptor.calculate_inverse_descriptor(r, radial_values)
    
    plt.plot(
        r, 
        inverse_values / np.max(inverse_values),
        label=label
    )
    

r = np.linspace(0, 5, 500)

plot_radial_descriptor(r, descriptor_origin, "Origin-50")
plot_radial_descriptor(r, descriptor_concentrated, "Man-50")
plot_radial_descriptor(r, descriptor_unifrom_radial_5, "Unif-5")
plot_radial_descriptor(r, descriptor_unifrom_radial_50, "Unif-50")
plot_radial_descriptor(r, descriptor_unifrom_radial_250, "Unif-250")


plt.plot([distance_CC]*2, [0, 1.1], "--", label="CC")
plt.plot([distance_CH]*2, [0, 1.1], "--", label="CH")
plt.plot([
    np.sqrt(distance_CH**2 + distance_CC**2 - 2*distance_CH*distance_CC*np.cos(angle_HCC))
]*2, [0, 1.1], "--", label="C1H3")


plt.xlabel("r / $\AA$")
plt.ylabel("Radial Activation / 1")

plt.legend()
plt.tight_layout()

plt.savefig("/home/jo/Repos/MastersThesis/GenericDescriptors/figures/RadialActivation.pdf")<jupyter_output><empty_output><jupyter_text>Angular<jupyter_code>from SCFInitialGuess.descriptors.models import make_periodic_uniform
from SCFInitialGuess.descriptors.coordinate_descriptors import \
    Gaussians, PeriodicGaussians, IndependentAngularDescriptor
from SCFInitialGuess.descriptors.high_level import NonWeighted
from SCFInitialGuess.descriptors.cutoffs import BehlerCutoff1

model = r_model_concentrated

descriptor_pg_10 = NonWeighted(
    Gaussians(*model),
    IndependentAngularDescriptor(
        PeriodicGaussians(*make_periodic_uniform(10, 2*np.pi)),
        PeriodicGaussians(*make_periodic_uniform(50, 2*np.pi)),
    ),
    BehlerCutoff1(5)
)

descriptor_pg_50 = NonWeighted(
    Gaussians(*model),
    IndependentAngularDescriptor(
        PeriodicGaussians(*make_periodic_uniform(50, 2*np.pi)),
        PeriodicGaussians(*make_periodic_uniform(50, 2*np.pi)),
    ),
    BehlerCutoff1(5)
)
from SCFInitialGuess.descriptors.coordinate_descriptors import SPHAngularDescriptor
from SCFInitialGuess.descriptors.models import make_uniform

model = make_uniform(25, 5, eta_max=50, eta_min=10)

descriptor_sph_2 = NonWeighted(
    Gaussians(*model),
    SPHAngularDescriptor(2),
    BehlerCutoff1(5)
)

descriptor_sph_5 = NonWeighted(
    Gaussians(*model),
    SPHAngularDescriptor(5),
    BehlerCutoff1(5)
)

descriptor_sph_15 = NonWeighted(
    Gaussians(*model),
    SPHAngularDescriptor(15),
    BehlerCutoff1(5)
)
def plot_azimuthal_descriptor_pg(phi, descriptor, label):
    values = descriptor.calculate_atom_descriptor(
        0,
        mol,
        descriptor.number_of_descriptors
    )

    n_radial = descriptor.radial_descriptor.number_of_descriptors
    n_azimuthal = descriptor.angular_descriptor.azimuthal_descriptor.number_of_descriptors

    angular_values = values[n_radial:(n_radial+n_azimuthal)]

    inverse_values = descriptor.angular_descriptor.azimuthal_descriptor.calculate_inverse_descriptor(
        phi,
        #r=None,
        #phi=phi,
        #theta=theta,
        y=angular_values
    )
    
    plt.plot(
        phi / (np.pi ),
        inverse_values/np.max(inverse_values),
        label=label
    )

def plot_azimuthal_descriptor_sph(phi, descriptor, label):
    values = descriptor.calculate_atom_descriptor(
        0,
        mol,
        descriptor.number_of_descriptors
    )

    n_radial = descriptor.radial_descriptor.number_of_descriptors

    angular_values = values[n_radial:]

    inverse_values = descriptor.angular_descriptor.calculate_inverse_descriptor(
        r=None,
        phi=phi,
        theta=theta,
        y=angular_values
    )
    
    plt.plot(
        phi / (np.pi ),
        inverse_values/np.max(inverse_values),
        label=label
    )

theta = np.pi/ 2

phi = np.linspace(-np.pi*0.5, 2.5 * np.pi, 500)

plot_azimuthal_descriptor_pg(phi, descriptor_pg_10, "PG-10")
plot_azimuthal_descriptor_pg(phi, descriptor_pg_50, "PG-50")

plot_azimuthal_descriptor_sph(phi, descriptor_sph_5, "SPH-5")
plot_azimuthal_descriptor_sph(phi, descriptor_sph_15, "SPH-15")

plt.plot([121.4/180]*2, [0, 1.1], "--", label="HCC")
plt.plot([(121.4 + 117.2)/180]*2, [0, 1.1], "--", label="HCC + HCH")

plt.xlabel(r"$\varphi$ / $\pi$")
plt.ylabel("Azimuthal Activation / 1")

#plt.legend(bbox_to_anchor=(0.3, 1),
#          bbox_transform=plt.gcf().transFigure)
plt.legend()
plt.tight_layout()

plt.savefig("/home/jcartus/Repos/MastersThesis/GenericDescriptors/figures/AzimuthalActivation.pdf")<jupyter_output><empty_output><jupyter_text>Polar<jupyter_code>n_polar = descriptor.angular_descriptor.polar_descriptor.number_of_descriptors
theta = np.linspace(0, np.pi, 500)


plt.plot(
    theta / (np.pi),
    descriptor_concentrated.angular_descriptor.polar_descriptor.calculate_inverse_descriptor(
        theta,
        values[(n_radial + n_azimuthal):]
    )
)<jupyter_output><empty_output><jupyter_text>Surface<jupyter_code>descriptor = descriptor_concentrated
values = descriptor.calculate_atom_descriptor(
    0,
    mol,
    descriptor.number_of_descriptors
)
R, Phi = np.meshgrid(r, np.linspace(- 1.5 * np.pi, 1.5 * np.pi, 500))

theta = np.pi / 2

z = activation(
    R.flatten(),
    Phi.flatten(),
    theta,
    descriptor_concentrated,
    values
)


plt.contour(
    Phi / np.pi,
    R,
    z.reshape(*R.shape),
    locator=mpl.ticker.LogLocator()
)
plt.colorbar()
plt.ylabel("R / 1")
plt.xlabel("$\phi$ / 1/$\phi$")
R, Phi = np.meshgrid(
    np.linspace(0, 3, 500), 
    np.linspace(0, 2 * np.pi, 500)
)

theta = np.pi / 2



z = activation(
    R.flatten(),
    Phi.flatten(),
    theta,
    descriptor,
    values
)
th = 1e-10
z[z<=th] = th
plt.rc("figure", figsize=(6,4))
fig, ax = plt.subplots()

plt.contourf(
    (R * np.cos(Phi)),
    (R * np.sin(Phi)),
    z.reshape(*R.shape),
    locator=mpl.ticker.LogLocator()
)



plt.xlabel("x / $\AA$")
plt.ylabel("y / $\AA$")

plt.colorbar(label="Total Activation / 1")

plt.scatter(pos[2:, 0], pos[2:, 1], marker="x", color="r", s=60, label="H")
plt.scatter(pos[:2, 0], pos[:2, 1], marker="x", color="g", label="C")
#plt.scatter(pos[:2, 0], pos[:2, 1], marker="x", color="orange", label="C")


plt.legend()

plt.gca().set_aspect('equal', adjustable='box')
plt.draw()

#ax = plt.gca()
#ax.set_ylabel=("y / $\AA$")
#ax.set_xlabel=("x / $\AA$")

plt.tight_layout()


plt.savefig("/home/jo/Repos/MastersThesis/AtomicBlocks/figures/XYPlaneActivation.pdf")<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1377396131494492, 'fraction_numerical': 0.03469511666232978, 'mean_word_length': 4.226654578422484, 'pattern_counts': {'":': 0, '<': 38, '<?xml version=': 0, '>': 37, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '259731', 'license': 'MIT', 'n_tokens_mistral': 4670, 'n_tokens_neox': 4340, 'n_words': 603, 'path': 'thesis/notebooks/Descriptors/Descriptors2.ipynb', 'repo_name': 'jcartus/SCFInitialGuess', 'size': '11529'}