#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
PRESSION ATMOSPHERIQUE DANS LE MODELE ISA
Created on Sun Apr  4 11:45:19 2021

@author: Emilie Frémont
"""

## Importation des bibliothèques utiles
## ------------------------------------

import  as np                  # pour la manipulation des tableaux
import  as plt     # pour les représentations graphiques
from  import odeint  # pour la résolution des équations différentielles


## Définition des constantes du problème
## -------------------------------------

g =                             # accélération de la pesanteur (en m/s^2)
Mair =                        # masse molaire de l'air (en kg/mol)
R =                          # constante du gaz parfait (en J/K/mol)
Tsol =                          # température de l'atmosphère au niveau du sol (en K)
Psol =                       # pression de l'atmosphère au niveau du sol (en Pa)


## Définition du gradient thermique vertical selon le modèle ISA
## -------------------------------------------------------------

def kISA(z):
    """ z est l'altitude en mètres. La fonction renvoie la valeur du gradient thermique
    vertical à l'altitude z (en K/m). """
    if 0 <= z < 11e3: return ...
    elif z ....
    ....
    else: return None
    
    
## -----------------------------------------------------------------------------------------------
## A - Résolution numérique
## -----------------------------------------------------------------------------------------------

## Définition du système différentiel à résoudre
    
def systDiff(TP,z):
    """ 
    TP désigne le vecteur inconnu de dimension 2 (TP[0] : température ; TP[1] : pression) ;
    z désigne l'altitude. 
    La fonction renvoie respectivement la dérivée de la température et la dérivée de
    la pression à l'altitude z.
    """
    # Lois prévues par le modèle théorique
    dT =                # dérivée verticale de la température
    dP =      # dérivée verticale de la pression
    
    return ...

## Définition des conditions aux limites

CAL = [....]

## Définition de l'ensemble des valeurs de z pour lesquelles on cherche la solution
## numérique approchée du système différentiel précédent

z = np.linspace(0, 85e3, 10000)     # on choisit 10000 points régulièrement espacés entre 0
                                    # et 85 km d'altitude

## La résolution en elle-même
                                    
TP = odeint(systDiff, CAL, z)
T = TP[:,0]                         # extraction des valeurs de la température
P = TP[:,1]                         # extraction des valeurs de la pression

## Représentation graphique des résultats

def PisoT(z, T0 = 288):             # pour comparaison
    """ Calcule la valeur de P à l'altitude z selon le modèle isotherme. 
    Par défaut, la température est fixée à 288 K. """
    return ....

plt.figure(figsize=(13,6))
plt.subplot(1,2,1)                  # graphique de gauche
plt.title("Evolution de la température avec l'altitude")
plt.plot(z*1e-3,T, 'b-')
plt.ylim(180,300), plt.ylabel(r"$T$ (en K)")
plt.xlim(0,85), plt.xlabel(r"$z$ (en km)")
plt.grid()
plt.subplot(1,2,2)                  # graphique de droite
plt.title("Evolution de la pression avec l'altitude")
plt.plot(z*1e-3,P*1e-5, 'b-', label = "Modèle ISA")
plt.plot(z*1e-3,PisoT(z)*1e-5,'m:', label = r"Modèle isotherme ($T=T_{\rm sol}$)")
#plt.plot(z*1e-3,P*1e-3, 'b-')
plt.ylim(0,1), plt.ylabel(r"$P$ (en bar)")
plt.xlim(0,85), plt.xlabel(r"$z$ (en km)")
plt.legend(loc=0)
plt.grid(which = 'both')
plt.show()

plt.figure()                        # comparaison des modèles ISA et isotherme
plt.plot(z*1e-3, abs(P-PisoT(z))/P, 'y-')
plt.xlim(0,20), plt.xlabel(r"$z$ (en km)")
plt.ylim(0,1), plt.ylabel(r"écart relatif : $\dfrac{\vert P-P_{\rm iso}\vert}{P}$")
plt.grid()
plt.show()


## -----------------------------------------------------------------------------------------------
## B - Application au dimensionnement d'un ballon-sonde
## -----------------------------------------------------------------------------------------------
'''
## Constantes physiques supplémentaires

m =                              # masse du système (hors hélium embarqué, en kg)
MHe =                         # masse molaire de l'hélium (en kg/mol)

## Calcul du volume d'enveloppe permettant d'atteindre l'altitude maximale z

V = 

## Représentation graphique

plt.figure()
plt.plot(z*1e-3, V, 'm-')
plt.xlim(0,50), plt.xlabel(r"$z_{\max}$ (en km)")
plt.ylim(0,2000), plt.ylabel(r"$V$ (en m$^3$)")
plt.grid()
plt.show()

## Pour atteindre l'altitude z_max = 36 km...
## On recherche par dichotomie l'indice de position de l'élément du tableau de valeurs de z
## dont la valeur est la plus proche de z_max, puis on détermine le volume V correspondant

z_max =                         # altitude maximale visée (en m)
  

    
V36 =                           # volume du ballon recherché
d36 =           # diamètre correspondant
print("Volume d'enveloppe nécessaire : {:.1f} m^3".format(V36))
print("Diamètre correspondant : {:.2f} m".format(d36))


'''