// Série 07,  Charge et décharge de condensateur

funcprot(0); // pour ne pas indiquer une redéfinition des fonctions

function ex1()
//============
// Définition de l'équation différentielle :  UC' = - UC / RC
// Simulation d'un circuit R C, c.f. ex. 1 de la série 7
// Décharge d'un condensateur
  function y_prime=f(t, y)
  //======================
  // t représente le temps
  // y  représente la tension aux bornes du condensateur    U  en volts
  // y_prime  représente la dérivée de y
  C = 1.0e-4;  // valeur de la capacité en Farads
  R = 15e3;   // valeur de la résistance, en Ohms 
  y_prime = -y/(R*C);  // UC' = - UC / RC
  endfunction

// données :
t0    =  0;   // temps initial
tFin  =  5;   // temps final
y0    = 12;  // UC(0) = 12 [V]
N  = 500;      // nombre de points d'évaluation
at = linspace(t0, tFin, N+1);  // vecteur LIGNE

// résolution numérique en appelant la fonction  ode  de SciLab
ay = ode( y0,..     // condition initiale
          t0,..     // x_initial,
          at,..     // série des t(j) pour lesquels  y(j)  sera calculé 
          1e-4,..   // tolérance d'erreur relative
          1e-5,..   // tolérance d'erreur absolue
          f..       // fonction définissant l'équation différentielle
          );
// L'algorithme essaye d'avoir une erreur inférieure au
// maximum entre l'erreur absolue et l'erreur relative.

// graphique de la solution numérique
win0 = scf(1);
clf();
plot(at, ay(1,:), 'b');  // 'b' = blue, pour les tensions aux bornes du condensateur
gcf().figure_name="Évolution de la tension en fonction du temps. fig. n°1" // Titre de la fenêtre
title("Évolution de la tension en fonction du temps", 'fontsize', 4);
xlabel("Temps [s]", 'fontsize', 3, "color", [0, 0, 0]);
ylabel("Tension [V]", 'fontsize', 3, "color", [0, 0, 1]);
endfunction
//ex1();
// On remarque que la charge du condensateur vaut U0 / e = 4,4145 [V] après
// un temps égale à R * C = 1,5 secondes. 
// C'est bien ce qui est prédit par la théorie.

function ex2()
//============
// Définition de l'équation différentielle :  UC' = (U0 - UC) / RC   UC(t) = U0 * (1 - exp(-t/RC))
// Simulation d'un circuit R C, c.f. ex. 2 de la série 7
// Charge d'un condensateur
  function y_prime=f(t, y)
  //======================
  // t représente le temps
  // y  représente la tension aux bornes du condensateur    U  en volts
  // y_prime  représente la dérivée de y
  C = 1.0e-4;  // valeur de la capacité en Farads
  R = 15e3;   // valeur de la résistance, en Ohms 
  U0 = 12.0;  // volts
  y_prime = (U0 - y) / (R*C);  // UC' = (U0 - UC) / RC   UC(t) = U0 * (1 - exp(-t/RC))
  endfunction

// données :
t0    =  0;   // temps initial
tFin  =  5;   // temps final
y0    =  0;  // UC(0) = 0 [V]
N  = 500;      // nombre de points d'évaluation
at = linspace(t0, tFin, N+1);  // vecteur LIGNE

// résolution numérique en appelant la fonction  ode  de SciLab
ay = ode( y0,..     // condition initiale
          t0,..     // x_initial,
          at,..     // série des t(j) pour lesquels  y(j)  sera calculé 
          1e-4,..   // tolérance d'erreur relative
          1e-5,..   // tolérance d'erreur absolue
          f..       // fonction définissant l'équation différentielle
          );
// L'algorithme essaye d'avoir une erreur inférieure au
// maximum entre l'erreur absolue et l'erreur relative.

// graphique de la solution numérique
win0 = scf(2);
clf();
plot(at, ay(1,:), 'b');  // 'b' = blue, pour les tensions aux bornes du condensateur
gcf().figure_name="Évolution de la tension en fonction du temps. fig. n°2" // Titre de la fenêtre
title("Évolution de la tension en fonction du temps", 'fontsize', 4);
xlabel("Temps [s]", 'fontsize', 3, "color", [0, 0, 0]);
ylabel("Tension [V]", 'fontsize', 3, "color", [0, 0, 1]);
endfunction
//ex2();
// On remarque que la charge du condensateur vaut U0 * (1 - 1/e) = 7,585 [V] après
// un temps égale à R * C = 1,5 secondes. 
// C'est bien ce qui est prédit par la théorie.

function ex3()
//============
// Définition de l'équation différentielle :  
// Simulation d'un circuit R C, c.f. ex. 2 de la série 7
// Charge d'un condensateur
      
  function y_prime=f(t,y)
  //=====================
  // t = le temps
  // y(1) = tension aux bornes du condensateur au temps  t.
  // Lors de la décharge : y' = -y / (R*C)
  // Lors de la charge :   y' = (U0 - y) / (R*C)
  // On charge si y < Umin
  // On décharge si y > Umax
  global fCharge;  // pour dire d'utiliser la variable fCharge déclarée à l'ext. de la fonction.
  
  C = 1.0e-4;  // valeur de la capacité en Farads
  R = 15e3;   // valeur de la résistance, en Ohms 
  U0 = 12.0;  // volts
  Umin =  1;  // volts
  Umax = 11;  // volts
  if (y < Umin) fCharge = %t; end;
  if (y > Umax) fCharge = %f; end;
  if (fCharge) y_prime = (U0 - y) / (R*C); // cas de charge
  else         y_prime = ( 0 - y) / (R*C); // cas de décharge
  end;  
  endfunction

// données :
fCharge = %t;
t0    =  0;   // temps initial
tFin  = 20;   // temps final
y0    =  0;  // UC(0) = 0 [V]
N  = 500;      // nombre de points d'évaluation
at = linspace(t0, tFin, N+1);  // vecteur LIGNE

// résolution numérique en appelant la fonction  ode  de SciLab
ay = ode( y0,..     // condition initiale
          t0,..     // x_initial,
          at,..     // série des t(j) pour lesquels  y(j)  sera calculé 
          1e-4,..   // tolérance d'erreur relative
          1e-5,..   // tolérance d'erreur absolue
          f..       // fonction définissant l'équation différentielle
          );
// L'algorithme essaye d'avoir une erreur inférieure au
// maximum entre l'erreur absolue et l'erreur relative.

// graphique de la solution numérique
win0 = scf(3);
clf();
plot(at, ay(1,:), 'b');  // 'b' = blue, pour les tensions aux bornes du condensateur
gcf().figure_name="Évolution de la tension en fonction du temps. fig. n°3" // Titre de la fenêtre
title("Évolution de la tension en fonction du temps", 'fontsize', 4);
xlabel("Temps [s]", 'fontsize', 3, "color", [0, 0, 0]);
ylabel("Tension [V]", 'fontsize', 3, "color", [0, 0, 1]);
endfunction
//ex3();

function ex4()
//============
// Définition de l'équation différentielle :  
// Simulation d'un circuit R L C, 
  function y_prime=f(x, y)
  //======================
  // x représente le temps
  // y(1)  représente la tension aux bornes du condensateur    UC  en volts
  // y(2)  représente l'intensité du courant traversant la self  IL  en ampères
  // y_prime(1)  représente la dérivée de y(1)
  // y_prime(2)  représente la dérivée de y(2)
  // y et y_prime sont des vecteurs colonne
  C = 1.0e-6;  // valeur de la capacité en Farads
  L = 4;  // valeur de la self (inductance) en Henri
  R = 5e4;   // valeur de la résistance, en Ohms 
  UC = y(1);
  IL = y(2)
  y_prime = [ ( (1/C) * (IL - UC/R) ); ..  // UC' = (1/C) * (IL - UC/R)
               -(1/L) * UC];    // IL' = -(1/L) * UC
  // On peut aussi écrire
  //y_prime(1) = (1/C) * (y(2) - y(1)/R);  // UC' = (1/C) * (IL - UC/R)
  //y_prime(2) = -(1/L) * y(1);    // IL' = -(1/L) * UC
  endfunction

// données :
t0    =  0;    // temps initial
tFin  = 0.4;   // temps final
UC0   = 12;    // UC(0) = 12 [V]
IL0   =  0;    // IL(0) = 0 [A]
y0    = [UC0; IL0]; // vecteur COLONNE
N  = 500;      // nombre de points d'évaluation
at = linspace(t0, tFin, N+1);  // vecteur LIGNE

// résolution numérique en appelant la fonction  ode  de SciLab
ay = ode( y0,..     // condition initiale
          t0,..     // x_initial,
          at,..     // série des t(j) pour lesquels  y(j)  sera calculé 
          1e-4,..   // tolérance d'erreur relative
          1e-5,..   // tolérance d'erreur absolue
          f..       // fonction définissant l'équation différentielle
          );
// L'algorithme essaye d'avoir une erreur inférieure au
// maximum entre l'erreur absolue et l'erreur relative.

// graphique de la solution numérique
win0 = scf(4);
clf();
plot(at, ay(1,:), 'b');  // 'b' = blue, pour les tensions aux bornes du condensateur
plot(at, 1000*ay(2,:), 'r');  // 'r' = red , pour les intensités à travers la self
gcf().figure_name="Évolution de la tension et du courant en fonction du temps. fig. n°4" // Titre de la fenêtre
title("Tension [V] en bleue et Courant [mA] en rouge", 'fontsize', 4);
xlabel("Temps [s]", 'fontsize', 3, "color", [0, 0, 0]);
ylabel("Tension [V] et Courant [mA]", 'fontsize', 3, "color", [0, 0, 1]);
endfunction
//ex4();

function ex5(decal)
//=================
// Simulation d'un oscillateur  LC, entretenu
// decal est le décalage de la fréquence d'entretien, relativement
// a la fréquence de résonnance.
// decal = 1  => fréquence de résonnance.
// decal = 1 + 0.1 => fréquence 10% plus haut que le fréquence de résonnance.
// decal = 1 - 0.15 => fréquence 15% plus bas que la fréquence de résonnance.

  // Définition de l'équation différentielle :  
  // Simulation d'un circuit R L C.
  function y_prime=f(t, y)
  //======================
  // t représente le temps
  // y(1)  représente la tension aux bornes du condensateur    UC  en volts
  // y(2)  représente l'intensité du courant traversant la self  IL  en ampères
  // y_prime(1)  représente la dérivée de y(1)
  // y_prime(2)  représente la dérivée de y(2)
  // y et y_prime sont des vecteurs colonne
  C = 1.0e-6;  // valeur de la capacité en Farads
  L = 4;  // valeur de la self (inductance) en Henri
  R = 5e4;   // valeur de la résistance, en ohms 

  w0 = 1 / sqrt(L*C);  // pulsation de résonnance  [1/s]
  if (t == 0) then 
    disp("frequ_réson.=" + string(w0/2/%pi)); 
    disp("w0 =" + string(w0)); 
    end
  global w;
  w = decal*w0;  // décalage par rapport à la pulsation de résonnance.
  w = 500;
  // fréquence de résonnance = 1 / (2 * %pi * sqrt(L*C))  [Hz]
  // Période correspondante = 2 * %pi * sqrt(L*C)  [s]
//  U0 = 0.02*cos(w*t);  // générateur de tension
  U0 = 0.02*sign(sin(w*t));  // générateur de tension  
  UC = y(1);
  IL = y(2)
  y_prime = [ ( (1/C) * (IL - UC/R) ); // UC' = (1/C) * (IL - UC/R)
                (1/L) * (U0 - UC)];    // IL' = (1/L) * (U0 - UC)
  // On peut aussi écrire
  //y_prime(1) = (1/C) * (y(2) - y(1)/R); // dérivée de la tension Uc
  //y_prime(2) = (1/L) * (U0 - y(1));     // dérivée du courant traversant la self  IL.
  endfunction

global w

// données :
t0    =  0;   // temps initial
tFin  = 0.4;   // temps final
UC0   =  0;    // UC(0) = 0 [V]
IL0   =  0;    // IL(0) = 0 [A]
y0    = [UC0; IL0]; // vecteur COLONNE
N  = 5000;      // nombre de points d'évaluation
at = linspace(t0, tFin, N+1);  // vecteur LIGNE

// résolution numérique en appelant la fonction  ode  de SciLab
ay = ode( y0,..     // condition initiale
          t0,..     // x_initial,
          at,..     // série des t(j) pour lesquels  y(j)  sera calculé 
          1e-4,..   // tolérance d'erreur relative
          1e-5,..   // tolérance d'erreur absolue
          f..       // fonction définissant l'équation différentielle
          );
// L'algorithme essaye d'avoir une erreur inférieure au
// maximum entre l'erreur absolue et l'erreur relative.

// graphique de la solution numérique
win0 = scf(5);
clf();
plot(at, ay(1,:), 'b');  // 'b' = blue, pour les tensions aux bornes du condensateur
plot(at, 100*ay(2,:), 'r');  // 'r' = red , pour les intensités à travers la self
plot(at, 0.02*sign(sin(w*at)), 'g');  
gcf().figure_name="Évolution de la tension et du courant en fonction du temps. fig. n°5" // Titre de la fenêtre
title("Tension [V] en bleue et Courant [mA] en rouge", 'fontsize', 4);
xlabel("Temps [s]", 'fontsize', 3, "color", [0, 0, 0]);
ylabel("Tension [V] et Courant [mA]", 'fontsize', 3, "color", [0, 0, 1]);
endfunction
//ex5(0.99);
