<commit_before>from scipy.interpolate import RegularGridInterpolator
from scipy import interpolate, optimize, integrate

#import multiprocessing as mp
import matplotlib.pyplot as plt
import numpy as np

import iofile
import rotations

def z_slides(x, y, z, function):
    return function([x, y, z])

def z_slides_r(y, x, z, function):
    return function([x, y, z])


def function_multi(angle, center, align, radial, ranges, function_use):
    x, y, z = ranges
    dens_rot = function_use(np.array([[[rotations.rotate_align_z([i, j, k], angle, center=center,
                                                                                   align=align,
                                                                                   radial=radial
                                                                 ) for k in z] for j in y] for i in x]))
    return {angle:dens_rot}

def function_measure(z_slide, ranges, function_use, epsabs, epsrel):
    minx, maxx, miny, maxy = ranges
    [integral, error] = integrate.dblquad(z_slides, miny, maxy, lambda x: minx, lambda x: maxx,
                                                    args=(z_slide, function_use),
                                                    epsabs=epsabs, epsrel=epsrel)

    return {z_slide: integral}


def function_measure_1d(z_coordinate, z_slide, ranges, function_use, epsabs, epsrel):
    minx, maxx, miny, maxy = ranges
    [integral, error] = integrate.quad(z_slides_r, miny, maxy,
                                                   args=(z_coordinate, z_slide, function_use),
                                                   epsabs=epsabs, epsrel=epsrel)

    return {z_slide: integral}


class Calculation:

    def __init__(self,
                 electronic_density,
                 ranges,
                 order=1,
                 center=(0, 0, 0),
                 align=(0, 0, 1),
                 radial=False,
                 r_grid=(20, 20, 20),
                 ):

        self._electronic_density = None
        self._order = order
        self._center = center
        self._align = align
        self._ranges_cart = ranges
        self._radial = radial

        self._density = None
        self._total_overlap = None
        self._measure = None

        self._grid = r_grid

        self._fn_electronic_density = None
        self._fn_overlap = None
        self._fn_density2 = None
        self._fn_density = None


        print('Input cube bounds (Cartesian)')
        print(' x: {0} {1}'.format(ranges[0][0], ranges[0][-1]))
        print(' y: {0} {1}'.format(ranges[1][0], ranges[1][-1]))
        print(' z: {0} {1}'.format(ranges[2][0], ranges[2][-1]))


        self._fn_electronic_density = RegularGridInterpolator(self._ranges_cart, electronic_density, bounds_error=False, fill_value=0)

        corner = [[ranges[0][0], ranges[1][0], ranges[2][0]],
                [ranges[0][0], ranges[1][0], ranges[2][-1]],
                [ranges[0][0], ranges[1][-1], ranges[2][0]],
                [ranges[0][-1], ranges[1][0], ranges[2][0]],
                [ranges[0][-1], ranges[1][-1], ranges[2][-1]],
                [ranges[0][-1], ranges[1][-1], ranges[2][0]],
                [ranges[0][-1], ranges[1][0], ranges[2][-1]],
                [ranges[0][0], ranges[1][-1], ranges[2][-1]]
                ]

        vector = self._align

        distances = [np.linalg.norm(np.cross(vector, center-np.array(point)))/np.linalg.norm(vector) for point in corner]
        positions = [-np.dot(center-np.array(point), vector)/np.linalg.norm(vector) for point in corner]
        max_rad = np.max(distances)

        print ('Oriented cube bounds')
        if radial:
            print('Using cylindrical')

    #        r1 = np.max([ranges[0][-1]-center[0],center[0]-ranges[0][1]])
    #        r2 = np.max([ranges[1][-1]-center[1],center[1]-ranges[1][1]])
    #        r3 = np.max([ranges[2][-1]-center[2],center[2]-ranges[2][1]])

        #    max_rad = np.max([r1, r2, r3])
        #    long_c = np.linspace(np.linalg.norm(center)-max_rad,np.linalg.norm(center)+max_rad, self._radial_grid[2])

            long_c = np.linspace(np.min(positions), np.max(positions), self._grid[2])
            angle = np.linspace(0, 2 * np.pi, self._grid[1])
            r = np.linspace(0, float(max_rad), self._grid[0])

            self._ranges = [long_c, angle, r]

            print(' z:     {0} {1} ({2})'.format(self._ranges[0][0], self._ranges[0][-1], self._grid[0]))
            print(' angle: {0} {1} ({2})'.format(self._ranges[1][0], self._ranges[1][-1], self._grid[1]))
            print(' r:     {0} {1} ({2})'.format(self._ranges[2][0], self._ranges[2][-1], self._grid[2]))

        else:
            print ('Using Cartesian')
            z = np.linspace(np.min(positions), np.max(positions), self._grid[2])
            y = np.linspace(center[1] - max_rad, center[1] + max_rad, self._grid[1])
            x = np.linspace(center[0] - max_rad, center[0] + max_rad, self._grid[0])

            self._ranges = [x, y, z]

            print(' x:     {0} {1} ({2})'.format(self._ranges[0][0], self._ranges[0][-1], self._grid[0]))
            print(' y:     {0} {1} ({2})'.format(self._ranges[1][0], self._ranges[1][-1], self._grid[1]))
            print(' z:     {0} {1} ({2})'.format(self._ranges[2][0], self._ranges[2][-1], self._grid[2]))

  #          self._ranges = self._ranges_cart


    def plot_full(self, step, rotation=0):
        x, y, z = self._ranges_cart
        Y, X = np.meshgrid(y, x)
        density_slide = self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, 0],
                                                                                       rotation,
                                                                                       center=self._center,
                                                                                       align=self._align
                                                                                       ) for j in y] for i in x]))

        for val in np.arange(self._ranges[0][0]+step, self._ranges[0][-1], step):
            density_slide += self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, val],
                                                                                            rotation,
                                                                                            center=self._center,
                                                                                            align=self._align
                                                                                            ) for j in y] for i in x]))
        density_slide *=step
    #    levels = np.linspace(0, 0.5, 40)
        plt.xlabel("Bohr")
        plt.ylabel("Bohr")
        plt.contourf(X, Y, density_slide)
        plt.colorbar().set_label("Electronic density (e-/Bohr^2)")
        plt.show()


    def plot_full_rad(self, step, rotation=0):
        x, y, z = self._ranges
        Z, Y = np.meshgrid(z, y)
        density_slide = self.fn_electronic_density(np.array([[rotations.rotate_align_z([0, j, k],
                                                                                       rotation,
                                                                                       center=self._center,
                                                                                       align=self._align,
                                                                                       radial=self._radial
                                                                                       ) for k in z] for j in y]))

        for val in np.arange(self._ranges[0][0]+step, self._ranges[0][-1], step):
            density_slide += self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, val],
                                                                                            rotation,
                                                                                            center=self._center,
                                                                                            align=self._align,
                                                                                            radial=self._radial
                                                                                            ) for j in y] for i in x]))
        density_slide *=step

        plt.contourf(Y, Z, density_slide)
        if self._radial:
            plt.xlabel("Angle (Radian)")
            plt.ylabel("Radius (Bohr)")
        else:
            plt.xlabel("Bohr")
            plt.ylabel("Bohr")
        plt.colorbar().set_label("Electronic density (e-/Bohr^3)")
        plt.show()



    def plot_slide(self, val, rotation=0):
        x, y, z = self._ranges_cart

        Y, X = np.meshgrid(y, x)
        density_slide = self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, val],
                                                                                       rotation,
                                                                                       center=self._center,
                                                                                       align=self._align
                                                                                       ) for j in y] for i in x]))
    #    levels = np.linspace(0, 0.5, 40)
        plt.xlabel("Bohr")
        plt.ylabel("Bohr")
        plt.contourf(X, Y, density_slide)
        plt.colorbar().set_label("Electronic density (e-/Bohr^3)")
        plt.show()
<commit_msg>Fix coordinates range in slide plot<commit_after>from scipy.interpolate import RegularGridInterpolator
from scipy import interpolate, optimize, integrate

#import multiprocessing as mp
import matplotlib.pyplot as plt
import numpy as np

import iofile
import rotations

def z_slides(x, y, z, function):
    return function([x, y, z])

def z_slides_r(y, x, z, function):
    return function([x, y, z])


def function_multi(angle, center, align, radial, ranges, function_use):
    x, y, z = ranges
    dens_rot = function_use(np.array([[[rotations.rotate_align_z([i, j, k], angle, center=center,
                                                                                   align=align,
                                                                                   radial=radial
                                                                 ) for k in z] for j in y] for i in x]))
    return {angle:dens_rot}

def function_measure(z_slide, ranges, function_use, epsabs, epsrel):
    minx, maxx, miny, maxy = ranges
    [integral, error] = integrate.dblquad(z_slides, miny, maxy, lambda x: minx, lambda x: maxx,
                                                    args=(z_slide, function_use),
                                                    epsabs=epsabs, epsrel=epsrel)

    return {z_slide: integral}


def function_measure_1d(z_coordinate, z_slide, ranges, function_use, epsabs, epsrel):
    minx, maxx, miny, maxy = ranges
    [integral, error] = integrate.quad(z_slides_r, miny, maxy,
                                                   args=(z_coordinate, z_slide, function_use),
                                                   epsabs=epsabs, epsrel=epsrel)

    return {z_slide: integral}


class Calculation:

    def __init__(self,
                 electronic_density,
                 ranges,
                 order=1,
                 center=(0, 0, 0),
                 align=(0, 0, 1),
                 radial=False,
                 r_grid=(20, 20, 20),
                 ):

        self._electronic_density = None
        self._order = order
        self._center = center
        self._align = align
        self._ranges_cart = ranges
        self._radial = radial

        self._density = None
        self._total_overlap = None
        self._measure = None

        self._grid = r_grid

        self._fn_electronic_density = None
        self._fn_overlap = None
        self._fn_density2 = None
        self._fn_density = None


        print('Input cube bounds (Cartesian)')
        print(' x: {0} {1}'.format(ranges[0][0], ranges[0][-1]))
        print(' y: {0} {1}'.format(ranges[1][0], ranges[1][-1]))
        print(' z: {0} {1}'.format(ranges[2][0], ranges[2][-1]))


        self._fn_electronic_density = RegularGridInterpolator(self._ranges_cart, electronic_density, bounds_error=False, fill_value=0)

        corner = [[ranges[0][0], ranges[1][0], ranges[2][0]],
                [ranges[0][0], ranges[1][0], ranges[2][-1]],
                [ranges[0][0], ranges[1][-1], ranges[2][0]],
                [ranges[0][-1], ranges[1][0], ranges[2][0]],
                [ranges[0][-1], ranges[1][-1], ranges[2][-1]],
                [ranges[0][-1], ranges[1][-1], ranges[2][0]],
                [ranges[0][-1], ranges[1][0], ranges[2][-1]],
                [ranges[0][0], ranges[1][-1], ranges[2][-1]]
                ]

        vector = self._align

        distances = [np.linalg.norm(np.cross(vector, center-np.array(point)))/np.linalg.norm(vector) for point in corner]
        positions = [-np.dot(center-np.array(point), vector)/np.linalg.norm(vector) for point in corner]
        max_rad = np.max(distances)

        print ('Oriented cube bounds')
        if radial:
            print('Using cylindrical')

    #        r1 = np.max([ranges[0][-1]-center[0],center[0]-ranges[0][1]])
    #        r2 = np.max([ranges[1][-1]-center[1],center[1]-ranges[1][1]])
    #        r3 = np.max([ranges[2][-1]-center[2],center[2]-ranges[2][1]])

        #    max_rad = np.max([r1, r2, r3])
        #    long_c = np.linspace(np.linalg.norm(center)-max_rad,np.linalg.norm(center)+max_rad, self._radial_grid[2])

            long_c = np.linspace(np.min(positions), np.max(positions), self._grid[2])
            angle = np.linspace(0, 2 * np.pi, self._grid[1])
            r = np.linspace(0, float(max_rad), self._grid[0])

            self._ranges = [long_c, angle, r]

            print(' z:     {0} {1} ({2})'.format(self._ranges[0][0], self._ranges[0][-1], self._grid[0]))
            print(' angle: {0} {1} ({2})'.format(self._ranges[1][0], self._ranges[1][-1], self._grid[1]))
            print(' r:     {0} {1} ({2})'.format(self._ranges[2][0], self._ranges[2][-1], self._grid[2]))

        else:
            print ('Using Cartesian')
            z = np.linspace(np.min(positions), np.max(positions), self._grid[2])
            y = np.linspace(center[1] - max_rad, center[1] + max_rad, self._grid[1])
            x = np.linspace(center[0] - max_rad, center[0] + max_rad, self._grid[0])

            self._ranges = [x, y, z]

            print(' x:     {0} {1} ({2})'.format(self._ranges[0][0], self._ranges[0][-1], self._grid[0]))
            print(' y:     {0} {1} ({2})'.format(self._ranges[1][0], self._ranges[1][-1], self._grid[1]))
            print(' z:     {0} {1} ({2})'.format(self._ranges[2][0], self._ranges[2][-1], self._grid[2]))

  #          self._ranges = self._ranges_cart


    def plot_full(self, step, rotation=0):
        x, y, z = self._ranges
        Y, X = np.meshgrid(y, x)
        density_slide = self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, 0],
                                                                                       rotation,
                                                                                       center=self._center,
                                                                                       align=self._align
                                                                                       ) for j in y] for i in x]))

        for val in np.arange(self._ranges[0][0]+step, self._ranges[0][-1], step):
            density_slide += self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, val],
                                                                                            rotation,
                                                                                            center=self._center,
                                                                                            align=self._align
                                                                                            ) for j in y] for i in x]))
        density_slide *=step
    #    levels = np.linspace(0, 0.5, 40)
        plt.xlabel("Bohr")
        plt.ylabel("Bohr")
        plt.contourf(X, Y, density_slide)
        plt.colorbar().set_label("Electronic density (e-/Bohr^2)")
        plt.show()


    def plot_full_rad(self, step, rotation=0):
        x, y, z = self._ranges
        Z, Y = np.meshgrid(z, y)
        density_slide = self.fn_electronic_density(np.array([[rotations.rotate_align_z([0, j, k],
                                                                                       rotation,
                                                                                       center=self._center,
                                                                                       align=self._align,
                                                                                       radial=self._radial
                                                                                       ) for k in z] for j in y]))

        for val in np.arange(self._ranges[0][0]+step, self._ranges[0][-1], step):
            density_slide += self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, val],
                                                                                            rotation,
                                                                                            center=self._center,
                                                                                            align=self._align,
                                                                                            radial=self._radial
                                                                                            ) for j in y] for i in x]))
        density_slide *=step

        plt.contourf(Y, Z, density_slide)
        if self._radial:
            plt.xlabel("Angle (Radian)")
            plt.ylabel("Radius (Bohr)")
        else:
            plt.xlabel("Bohr")
            plt.ylabel("Bohr")
        plt.colorbar().set_label("Electronic density (e-/Bohr^3)")
        plt.show()



    def plot_slide(self, val, rotation=0):
        x, y, z = self._ranges

        Y, X = np.meshgrid(y, x)
        density_slide = self.fn_electronic_density(np.array([[rotations.rotate_align_z([i, j, val],
                                                                                       rotation,
                                                                                       center=self._center,
                                                                                       align=self._align
                                                                                       ) for j in y] for i in x]))
    #    levels = np.linspace(0, 0.5, 40)
        plt.xlabel("Bohr")
        plt.ylabel("Bohr")
        plt.contourf(X, Y, density_slide)
        plt.colorbar().set_label("Electronic density (e-/Bohr^3)")
        plt.show()

{'fraction_non_alphanumeric': 0.12713097713097712, 'fraction_numerical': 0.020166320166320167, 'mean_word_length': 1.2019913023575188, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 210, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2844825', 'n_tokens_mistral': 5352, 'n_tokens_neox': 5183, 'n_words': 1045}