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.75)
  108.         alfa=0.75;
  109.      end
  110.   end
  111.   if (h2<=alfa*h1)
  112.      denoye='oui';
  113.   else
  114.      denoye='non';
  115.      x=sqrt(1-h2/h1);
  116.      beta=-2*alfa+2.6;
  117.      if (x>0.2)
  118.         KF=1-(1-x/sqrt(1-alfa))^beta;
  119.      else
  120.         KF=5*x*(1-(1-0.2/sqrt(1-alfa))^beta);
  121.      end
  122.      alfa1=1-0.14*(h2-w)/w;
  123.      if (alfa1<0.4)
  124.         alfa1=0.4;
  125.      end
  126.      if (alfa1>0.75)
  127.         alfa1=0.75;
  128.      end
  129.      if (h2<=alfa1*h1+(1-alfa1)*w)
  130.         partiel='oui';
  131.      else
  132.         partiel='non';
  133.      end
  134.   end
  135. % Seuil - Denoye ---------------------------------------------
  136.   if (surfacelibre=='oui') & (denoye=='oui')
  137.      Q=muf*L*R2G*(h1^1.5);
  138.      C=1;
  139.   end
  140. % Seuil - Noye -----------------------------------------------
  141.   if (surfacelibre=='oui') & (denoye=='non')
  142.      Q=KF*muf*L*R2G*(h1^1.5);
  143.      C=2;
  144.   end
  145. % Vanne - Denoye ---------------------------------------------
  146.   if (surfacelibre=='non') & (denoye=='oui')
  147.      Q=L*R2G*(mu*h1^1.5-mu1*(h1-w)^1.5);
  148.      C=3;
  149.   end
  150. % Vanne - Noye -----------------------------------------------
  151.   if (surfacelibre=='non') & (denoye=='non')
  152.      x1=sqrt(1-(h2-w)/(h1-w));
  153.      beta1=-2*alfa1+2.6;
  154.      if (x1>0.2)
  155.         KF1=1-(1-x1/sqrt(1-alfa1))^beta1;
  156.      else
  157.         KF1=5*x1*(1-(1-0.2/sqrt(1-alfa1))^beta1);
  158.      end
  159. %  Ennoyement partiel ---------------------------------
  160.      if (partiel=='oui')
  161.         Q=L*R2G*(KF*mu*(h1^1.5)-mu1*(h1-w)^1.5);
  162.         C=4;
  163. %  Ennoyement total -----------------------------------
  164.      else
  165.         Q=L*R2G*(KF*mu*(h1^1.5)-KF1*mu1*(h1-w)^1.5);
  166.         C=5;
  167.      end
  168.   end
  169. end
  170. % Surverse dans cas 1 et 2 ===================================
  171. if hv==0
  172.   hv=inf;
  173. end;
  174. if (T==1) | (T==2)
  175.   if (h1>w+hv)
  176.      surverse='oui';
  177.      alfa=2/3;
  178.      if (h2-w-hv<=alfa*(h1-w-hv))
  179.         denoyesurverse='oui';
  180.      else
  181.         denoyesurverse='non';
  182.      end
  183.   else
  184.      surverse='non';
  185.   end
  186.   % Surverse - Denoye -------------------------------------------
  187.   if (surverse=='oui')
  188.      if (denoyesurverse=='oui')
  189.         Q=Q+0.4*L*R2G*(h1-w-hv)^1.5;
  190.         C=C+10;
  191.      else
  192.         Q=Q+1.04*L*R2G*((h1-h2)^0.5)*(h2-w-hv);
  193.         C=C+10;
  194.      end
  195.   end
  196. end
  197. % Classique - Seuil ==========================================
  198. if (T==3)
  199.   Q=Cd*L*R2G*h1^1.5;
  200.   C=1;
  201. end
  202. % Classique - Orifice ========================================
  203. if (T==4)
  204.   if (h1>w)
  205.      Q=Cd*w*L*R2G*(h1-w/2)^0.5;
  206.      C=3;
  207.   else
  208.      Q=0;
  209.      C=0;
  210.   end
  211. end
  212. % Classique - Noyee ==========================================
  213. if (T==5)
  214.   if (h1>h2)
  215.      Q=Cd*w*L*R2G*(h1-h2)^0.5;
  216.      C=5;
  217.  else
  218.      Q=0;
  219.      C=0;
  220.   end
  221. end
  222. % Classique Cunge ===========================================
  223. if (T==6)
  224. % Conditions d'ecoulement ------------------------------------
  225.   if (h2<=2/3*h1)
  226.      denoye='oui';
  227.      if (w<=2/3*h1)
  228.         surfacelibre='non';
  229.         Q=Cd*L*R2G*w*((h1-w)^0.5);
  230.         C=3;
  231.      else
  232.         surfacelibre='oui';
  233.         Q=Cd*L*R2G/R32*h1^1.5;
  234.         C=1;
  235.      end
  236.   else
  237.      denoye='non';
  238.      if (w<=h2)
  239.         surfacelibre='non';
  240.         Q=Cd*L*R2G*w*((h1-h2)^0.5);
  241.         C=5;
  242.      else
  243.         surfacelibre='oui';
  244.         Q=Cd*L*R2G*h2*((h1-h2)^0.5);
  245.         C=2;
  246.      end
  247.   end
  248. 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