Equation and devices calibration (gates and weirs)

Appendix
- Available translations on this page : français -

We provide below the scripts of 2 Matlab .m functions (Qouvrage.m and calage_f.m) and a script file .m (calage_r.m) usefull for calibrating a discharge coefficient from mesurement datas using the equations programmed in SIC.

Devices modelisation (gates and weirs) (Script to copy in a file Qouvrage.m)

  1. function [Q,C,w] = Qouvrage(T,h1,h2,w,hv,L,Cd)
  2. % T=1 : Seuil Orifice avec surverse éventuelle (SIC)
  3. %   2 : Seuil Vanne avec surverse éventuelle (SIC)
  4. %   3 : Classique Seuil
  5. %   4 : Classique Orifice
  6. %   5 : Classique Noyee
  7. %   6 : Cunge avec differentes conditions d'ecoulement
  8. % h1  : cote amont (par rapport au radier)
  9. % h2  : cote aval (par rapport au radier)
  10. % w   : ouverture (par rapport au radier)
  11. % hv  : hauteur de la vanne (pour le calcul de la surverse éventuelle)
  12. % L   : largeur
  13. % Cd  : coefficient de debit
  14. % Q   : debit
  15. % C   : code couleur correspondant aux conditions d'ecoulement
  16. %       C=0; debit nul
  17. %       C=1; surface libre denoye
  18. %       C=2; surface libre noye
  19. %       C=3; charge denoye
  20. %       C=4; charge noye partiel
  21. %       C=5; charge noye total
  22. %       C=C+10; si surverse
  23. % Initialisation ---------------------------------------------
  24. R2G = sqrt(2*9.81);
  25. R32 = 3*sqrt(3)/2;
  26. % Tests ------------------------------------------------------
  27. Q=0;
  28. C=0;
  29. if (w==0 & hv==0)
  30.   Q=0;
  31.   return
  32. end
  33. if (h2>h1)
  34.   Q=0;
  35.   return
  36. end
  37. if (h1==0)
  38.   Q=0;
  39.   return
  40. end
  41. % Seuil - Orifice ===========================================
  42. if (T==1 & w>0)
  43. % Conditions d'ecoulement ------------------------------------
  44.   if (h1<=w)
  45.      surfacelibre='oui';
  46.   else
  47.      surfacelibre='non';
  48.   end
  49.   alfa=2/3;
  50.   if (h2<=alfa*h1)
  51.      denoye='oui';
  52.   else
  53.      denoye='non';
  54.      if (h2<=2/3*h1+w/3)
  55.         partiel='oui';
  56.      else
  57.         partiel='non';
  58.      end
  59.   end
  60. % Seuil - Denoye ---------------------------------------------
  61.   if (surfacelibre=='oui') & (denoye=='oui')
  62.      Q=Cd*L*R2G*(h1^1.5);
  63.      C=1;
  64.   end
  65. % Seuil - Noye -----------------------------------------------
  66.   if (surfacelibre=='oui') & (denoye=='non')
  67.      Cs=R32*Cd;
  68.      Q=Cs*L*R2G*((h1-h2)^0.5)*h2;
  69.     C=2;
  70.   end
  71. % Orifice - Denoye -------------------------------------------
  72.   if (surfacelibre=='non') & (denoye=='oui')
  73.      Q=Cd*L*R2G*(h1^1.5-(h1-w)^1.5);
  74.      C=3;
  75.   end
  76. % Orifice - Noye ---------------------------------------------
  77.   if (surfacelibre=='non') & (denoye=='non')
  78. %  Ennoyement partiel ---------------------------------
  79.      if (partiel=='oui')
  80.         Q=Cd*L*R2G*(R32*h2*(h1-h2)^0.5-(h1-w)^1.5);
  81.         C=4;
  82. %  Ennoyement total -----------------------------------
  83.      else
  84.         Cs=R32*Cd;
  85.         Q=Cs*L*R2G*((h1-h2)^0.5)*w;
  86.        C=5;
  87.      end
  88.   end
  89. end
  90. % Seuil - Vanne =============================================
  91. if (T==2 & w>0)
  92. % Calcul de parametres ---------------------------------------
  93.   mu0=2/3*Cd;
  94. % Conditions d'ecoulement ------------------------------------
  95.   if (h1<=w)
  96.      surfacelibre='oui';
  97.      muf=mu0-0.08;
  98.      alfa=0.75;
  99.   else
  100.      surfacelibre='non';
  101.      mu =mu0-0.08/(h1/w);
  102.      mu1=mu0-0.08/(h1/w-1);
  103.      alfa=1-0.14*h2/w;
  104.      if (alfa<0.4)
  105.         alfa=0.4;
  106.      end
  107.      if (alfa>0    alfa=0.75;
  108.      end
  109.   end
  110.   if (h2<=alfa*h1)
  111.      denoye='oui';
  112.   else
  113.      denoye='non';
  114.      x=sqrt(1-h2/h1);
  115.      beta=-2*alfa+2.6;
  116.      if (x>0.2)
  117.         KF=1-(1-x/sqrt(1-alfa))^beta;
  118.      else
  119.         KF=5*x*(1-(1-0.2/sqrt(1-alfa))^beta);
  120.      end
  121.      alfa1=1-0.14*(h2-w)/w;
  122.      if (alfa1<0.4)
  123.         alfa1=0.4;
  124.      end
  125.      if (alfa1>0.75)
  126.         alfa1=0.75;
  127.      end
  128.      if (h2<=alfa1*h1+(1-alfa1)*w)
  129.         partiel='oui';
  130.      else
  131.         partiel='non';
  132.      end
  133.   end
  134. % Seuil - Denoye ---------------------------------------------
  135.   if (surfacelibre=='oui') & (denoye=='oui')
  136.      Q=muf*L*R2G*(h1^1.5);
  137.      C=1;
  138.   end
  139. % Seuil - Noye -----------------------------------------------
  140.   if (surfacelibre=='oui') & (denoye=='non')
  141.      Q=KF*muf*L*R2G*(h1^1.5);
  142.      C=2;
  143.   end
  144. % Vanne - Denoye ---------------------------------------------
  145.   if (surfacelibre=='non') & (denoye=='oui')
  146.      Q=L*R2G*(mu*h1^1.5-mu1*(h1-w)^1.5);
  147.      C=3;
  148.   end
  149. % Vanne - Noye -----------------------------------------------
  150.   if (surfacelibre=='non') & (denoye=='non')
  151.      x1=sqrt(1-(h2-w)/(h1-w));
  152.      beta1=-2*alfa1+2.6;
  153.      if (x1>0.2)
  154.         KF1=1-(1-x1/sqrt(1-alfa1))^beta1;
  155.      else
  156.         KF1=5*x1*(1-(1-0.2/sqrt(1-alfa1))^beta1);
  157.      end
  158. %  Ennoyement partiel ---------------------------------
  159.      if (partiel=='oui')
  160.         Q=L*R2G*(KF*mu*(h1^1.5)-mu1*(h1-w)^1.5);
  161.         C=4;
  162. %  Ennoyement total -----------------------------------
  163.      else
  164.         Q=L*R2G*(KF*mu*(h1^1.5)-KF1*mu1*(h1-w)^1.5);
  165.         C=5;
  166.      end
  167.   end
  168. end
  169. % Surverse dans cas 1 et 2 ===================================
  170. if hv==0
  171.   hv=inf;
  172. end;
  173. if (T==1) | (T==2)
  174.   if (h1>w+hv)
  175.      surverse='oui';
  176.      alfa=2/3;
  177.      if (h2-w-hv<=alfa*(h1-w-hv))
  178.         denoyesurverse='oui';
  179.      else
  180.         denoyesurverse='non';
  181.      end
  182.   else
  183.      surverse='non';
  184.   end
  185.   % Surverse - Denoye -------------------------------------------
  186.   if (surverse=='oui')
  187.      if (denoyesurverse=='oui')
  188.         Q=Q+0.4*L*R2G*(h1-w-hv)^1.5;
  189.         C=C+10;
  190.      else
  191.         Q=Q+1.04*L*R2G*((h1-h2)^0.5)*(h2-w-hv);
  192.         C=C+10;
  193.      end
  194.   end
  195. end
  196. % Classique - Seuil ==========================================
  197. if (T==3)
  198.   Q=Cd*L*R2G*h1^1.5;
  199.   C=1;
  200. end
  201. % Classique - Orifice ========================================
  202. if (T==4)
  203.   if (h1>w)
  204.      Q=Cd*w*L*R2G*(h1-w/2)^0.5;
  205.      C=3;
  206.   else
  207.      Q=0;
  208.      C=0;
  209.   end
  210. end
  211. % Classique - Noyee ==========================================
  212. if (T==5)
  213.   if (h1>h2)
  214.      Q=Cd*w*L*R2G*(h1-h2)^0.5;
  215.      C=5;
  216.  else
  217.      Q=0;
  218.      C=0;
  219.   end
  220. end
  221. % Classique Cunge ===========================================
  222. if (T==6)
  223. % Conditions d'ecoulement ------------------------------------
  224.   if (h2<=2/3*h1)
  225.      denoye='oui';
  226.      if (w<=2/3*h1)
  227.         surfacelibre='non';
  228.         Q=Cd*L*R2G*w*((h1-w)^0.5);
  229.         C=3;
  230.      else
  231.         surfacelibre='oui';
  232.         Q=Cd*L*R2G/R32*h1^1.5;
  233.         C=1;
  234.      end
  235.   else
  236.      denoye='non';
  237.      if (w<=h2)
  238.         surfacelibre='non';
  239.         Q=Cd*L*R2G*w*((h1-h2)^0.5);
  240.         C=5;
  241.      else
  242.         surfacelibre='oui';
  243.         Q=Cd*L*R2G*h2*((h1-h2)^0.5);
  244.         C=2;
  245.      end
  246.   end
  247. end

Download

Intermediate function used for calibration (script to copy in a f.m file)

  1. function Q = f(Cd)
  2. % fonction utilisee pour le calage dans calage_r.m
  3. % T=1 : Seuil Orifice (SIC)
  4. %   2 : Seuil Vanne (SIC)
  5. %   3 : Classique Seuil
  6. %   4 : Classique Orifice
  7. %   5 : Classique Noyee
  8. % h1  : cote amont (par rapport au radier)
  9. % h2  : cote aval (par rapport au radier)
  10. % w   : ouverture
  11. % hv  : hauteur de la vanne (pour le calcul de la surverse éventuelle)
  12. % L   : largeur
  13. % Cd  : coefficient de debit
  14. % Q   : debit calcule par Qouvrage
  15. % Qm  : debit mesure pour le calage
  16. % C   : code de couleur correspondant aux conditions d'ecoulement
  17. % Initialisation
  18. T = 2;
  19. hr = 41.95
  20. h1 = 43.28 - hr;
  21. h2 = 42.77 - hr;
  22. hv = 0;
  23. w = 0.66;
  24. L = 1.39;
  25. Qm = 1.909;
  26. % Calcul
  27. [Q,C,w] = Qouvrage(T,h1,h2,w,hv,L,Cd);
  28. Q = Q - Qm;

Download

Calibration of discharge coefficient (script to copy in a file calage_r.m)

  1. % Fichier .M utilisé pour caler Cd
  2. % fzero est une fonction MatLab qui trouve le zéro
  3. % d'une fonction proche d'une valeur donnée
  4. Cd = fzero('calage_f',0.5)

Download

Once these codes have been copied in the .m file as above, you have to run the .m script (calage_r.m) for calculating the discharge coefficient Cd from measures h1, h2, Qm and the description of the device (h1, hv, w, L) using the equations of the evices programmed in SIC (several choices are proposed).

Search

dans le siteSpip

Navigation