# Path: pyqha/constants.py
# RY_KBAR          = 10.0 * AU_GPA / 2.0
# 
# K_BOLTZMANN_RY   = K_BOLTZMANN_SI / RYDBERG_SI
# 
# KB1=1.0/K_BOLTZMANN_RY/RY_TO_CMM1 # inverse Boltzmann constant in cm^{-1}/K
# 
# Path: pyqha/fitfreqgrun.py
# def freqmingrun(afreq, min0, nq, modes, ibrav, typefreq):
#     """
#     This function calculates the frequencies and the Gruneisen parameters
#     from the fitted polynomials coeffients (one
#     for each q point and mode) at the minimun point *min0* given in input. 
#     *afreq* is a :math:`nq*modes` numpy matrix containing the fitted polynomial coefficients.
#     It can be obtained from :py:func:`fitfreqxx`.
# 
#     It returns a :math:`nq*modes` matrix, each element [i,j] being the fitted frequency 
#     In addition, it returns a :math:`nq*modes*6` with the Gruneisein parameters.
#     Each element [i,j,k] is the the Gruneisein parameter at *nq=i*, *mode=j* and direction
#     *k* (for example, in hex systems *k=0* is *a* direction, *k=2* is *c* direction, others are zero)
# 
#     Note that the Gruneisein parameters are not multiplied for the lattice parameters. 
#     """
#     #nq = 10 # for testing quickly
#     f = np.zeros((nq,modes))
#     grun = np.zeros((6,nq,modes))   # initialize the Gruneisein parameters matrix
#     for i in range(0,nq):
#         for j in range(0,modes):
#             if typefreq=="quadratic":
#                 f[i,j] = fquadratic(min0,afreq[i][j],ibrav=4)
#                 if f[i,j]<1E-5:
#                     grun[:,i,j] = 0.0
#                 else:
#                     grun[:,i,j] = fquadratic_der(min0,afreq[i][j],ibrav)/f[i,j]*min0
#             elif typefreq=="quartic":
#                 f[i,j] = fquartic(min0,afreq[i][j],ibrav=4) 
#                 if f[i,j]<1E-5:
#                     grun[:,i,j] = 0.0
#                 else:
#                     grun[:,i,j] = fquartic_der(min0,afreq[i][j],ibrav)/f[i,j]*min0
#     
#     return f, grun 
# 
# Path: pyqha/properties_anis.py
# def compute_volume(celldms,ibrav=4):
#     """
#     Compute the volume given the *celldms*. Only for ibrav=4 for now, else 
#     returns 0.
#     """
#     if ibrav==4:
#         return 0.866025404*celldms[0]*celldms[0]*celldms[2]
# 
#     return 0
# 
# Path: pyqha/fitC.py
# def fS(aS,mintemp,typeCx):
#     """
#     An auxiliary function returning the elastic compliances 6x6 tensor at the
#     set of lattice parameters given in input as *mintemp*. These should be the
#     lattice parameters at a given temperature obtained from the free energy
#     minimization, so that S(T) can be obtained.
#     Before calling this function, the polynomial coefficients resulting from 
#     fitting the elastic compliances over a grid of lattice parameters, i.e. over
#     different geometries, must be obtained and passed as input in *aS*. 
#     *typeCx* defines what kind of polynomial to use for fitting ("quadratic" or
#     "quartic")
#     """
#     S = np.zeros((6,6))
#     for i in range(0,6):
#        for j in range(0,6):
#            if typeCx=="quadratic":
#                S[i,j] = fquadratic(mintemp,aS[i,j],ibrav=4)
#            elif typeCx=="quartic":
#                S[i,j] = fquartic(mintemp,aS[i,j],ibrav=4)  
#     return S


import numpy as np
import math
from .constants import RY_KBAR, K_BOLTZMANN_RY, KB1
from .fitfreqgrun import freqmingrun
from .properties_anis import compute_volume
from .fitC import fS
from multiprocessing import Pool
    from grunc import c_qvc  # This is the same routine c_qv implemented in C to speed it up
    |                                       | the lattice parameters at the minimun free energy at each temperature|  
    |                                       | . *S* is the elastic compliances tensor at 0 K calculated at *min0*. |
    |                                       | *weights* are the weights of the phonon frequencies.  *afreq*,       |
    |                                       | *fittypefreq* are the coefficients of the fitted polynomials for the |
    |                                       | phonon frequencies and the type of polynomial. They are used here to |
    |                                       | compute the average frequencies and Gruneisen parameters at each     |
    |                                       | temperature.                                                         |
    +---------------------------------------+----------------------------------------------------------------------+
    |*weights,afreq,fittypefreq,aS,fittypeS*| Use temperarature-dependent (0 K) volume, elastic compliances,       |
    |                                       | average frequencies and Gruneisen parameters. The volume is obtained |
    |                                       | from *minT*, a numpy array with the lattice parameters at the minimun|   
    |                                       | free energy at each temperature. *weights* are the weights of the    |  
    |                                       | phonon frequencies. *afreq*,                                         |
    |                                       | *fittypefreq* are the coefficients of the fitted polynomials for the |
    |                                       | phonon frequencies and the type of polynomial. They are used here to |
    |                                       | the average frequencies and Gruneisen parameters at each temperature.|
    |                                       | *aS,fittypeS* are the coefficients of the fitted polynomials for the |
    |                                       | elastic compliances and the polynomial type. They are used here to   |
    |                                       | compute the elastic compliances at each temperature as in the quasi- |
    |                                       | static approximation.                                                |
    +---------------------------------------+----------------------------------------------------------------------+

    By default, all the input parameters in the above table are ==None and if nothing
    is given in input the function will return None without computing anything.
    
    *nproc* is the number of processes to be run in parallel.
    
    """
    
    if (min0!=None):        # if the lattice parameters at the minimun 0 K energy are given

# Path: pyqha/constants.py
# RY_KBAR          = 10.0 * AU_GPA / 2.0
# 
# K_BOLTZMANN_RY   = K_BOLTZMANN_SI / RYDBERG_SI
# 
# KB1=1.0/K_BOLTZMANN_RY/RY_TO_CMM1 # inverse Boltzmann constant in cm^{-1}/K
# 
# Path: pyqha/fitfreqgrun.py
# def freqmingrun(afreq, min0, nq, modes, ibrav, typefreq):
#     """
#     This function calculates the frequencies and the Gruneisen parameters
#     from the fitted polynomials coeffients (one
#     for each q point and mode) at the minimun point *min0* given in input. 
#     *afreq* is a :math:`nq*modes` numpy matrix containing the fitted polynomial coefficients.
#     It can be obtained from :py:func:`fitfreqxx`.
# 
#     It returns a :math:`nq*modes` matrix, each element [i,j] being the fitted frequency 
#     In addition, it returns a :math:`nq*modes*6` with the Gruneisein parameters.
#     Each element [i,j,k] is the the Gruneisein parameter at *nq=i*, *mode=j* and direction
#     *k* (for example, in hex systems *k=0* is *a* direction, *k=2* is *c* direction, others are zero)
# 
#     Note that the Gruneisein parameters are not multiplied for the lattice parameters. 
#     """
#     #nq = 10 # for testing quickly
#     f = np.zeros((nq,modes))
#     grun = np.zeros((6,nq,modes))   # initialize the Gruneisein parameters matrix
#     for i in range(0,nq):
#         for j in range(0,modes):
#             if typefreq=="quadratic":
#                 f[i,j] = fquadratic(min0,afreq[i][j],ibrav=4)
#                 if f[i,j]<1E-5:
#                     grun[:,i,j] = 0.0
#                 else:
#                     grun[:,i,j] = fquadratic_der(min0,afreq[i][j],ibrav)/f[i,j]*min0
#             elif typefreq=="quartic":
#                 f[i,j] = fquartic(min0,afreq[i][j],ibrav=4) 
#                 if f[i,j]<1E-5:
#                     grun[:,i,j] = 0.0
#                 else:
#                     grun[:,i,j] = fquartic_der(min0,afreq[i][j],ibrav)/f[i,j]*min0
#     
#     return f, grun 
# 
# Path: pyqha/properties_anis.py
# def compute_volume(celldms,ibrav=4):
#     """
#     Compute the volume given the *celldms*. Only for ibrav=4 for now, else 
#     returns 0.
#     """
#     if ibrav==4:
#         return 0.866025404*celldms[0]*celldms[0]*celldms[2]
# 
#     return 0
# 
# Path: pyqha/fitC.py
# def fS(aS,mintemp,typeCx):
#     """
#     An auxiliary function returning the elastic compliances 6x6 tensor at the
#     set of lattice parameters given in input as *mintemp*. These should be the
#     lattice parameters at a given temperature obtained from the free energy
#     minimization, so that S(T) can be obtained.
#     Before calling this function, the polynomial coefficients resulting from 
#     fitting the elastic compliances over a grid of lattice parameters, i.e. over
#     different geometries, must be obtained and passed as input in *aS*. 
#     *typeCx* defines what kind of polynomial to use for fitting ("quadratic" or
#     "quartic")
#     """
#     S = np.zeros((6,6))
#     for i in range(0,6):
#        for j in range(0,6):
#            if typeCx=="quadratic":
#                S[i,j] = fquadratic(mintemp,aS[i,j],ibrav=4)
#            elif typeCx=="quartic":
#                S[i,j] = fquartic(mintemp,aS[i,j],ibrav=4)  
#     return S


# Path: pyqha/alphagruneisenp.py
import numpy as np
import math
from .constants import RY_KBAR, K_BOLTZMANN_RY, KB1
from .fitfreqgrun import freqmingrun
from .properties_anis import compute_volume
from .fitC import fS
from multiprocessing import Pool
    from grunc import c_qvc  # This is the same routine c_qv implemented in C to speed it up

    |                                       | the lattice parameters at the minimun free energy at each temperature|  
    |                                       | . *S* is the elastic compliances tensor at 0 K calculated at *min0*. |
    |                                       | *weights* are the weights of the phonon frequencies.  *afreq*,       |
    |                                       | *fittypefreq* are the coefficients of the fitted polynomials for the |
    |                                       | phonon frequencies and the type of polynomial. They are used here to |
    |                                       | compute the average frequencies and Gruneisen parameters at each     |
    |                                       | temperature.                                                         |
    +---------------------------------------+----------------------------------------------------------------------+
    |*weights,afreq,fittypefreq,aS,fittypeS*| Use temperarature-dependent (0 K) volume, elastic compliances,       |
    |                                       | average frequencies and Gruneisen parameters. The volume is obtained |
    |                                       | from *minT*, a numpy array with the lattice parameters at the minimun|   
    |                                       | free energy at each temperature. *weights* are the weights of the    |  
    |                                       | phonon frequencies. *afreq*,                                         |
    |                                       | *fittypefreq* are the coefficients of the fitted polynomials for the |
    |                                       | phonon frequencies and the type of polynomial. They are used here to |
    |                                       | the average frequencies and Gruneisen parameters at each temperature.|
    |                                       | *aS,fittypeS* are the coefficients of the fitted polynomials for the |
    |                                       | elastic compliances and the polynomial type. They are used here to   |
    |                                       | compute the elastic compliances at each temperature as in the quasi- |
    |                                       | static approximation.                                                |
    +---------------------------------------+----------------------------------------------------------------------+

    By default, all the input parameters in the above table are ==None and if nothing
    is given in input the function will return None without computing anything.
    
    *nproc* is the number of processes to be run in parallel.
    
    """
    
    if (min0!=None):        # if the lattice parameters at the minimun 0 K energy are given

        V=compute_volume(min0,ibrav)  # compute eq. volume at 0 K