SUBROUTINE CAMIL(U,Y,PARA)
PARAMETER(NBCM=$NBCM,NBUM=$NBUM,NSING=$NTO,NBIEF=$IB)
PARAMETER(NPOIN=$NPOIN,NSECT=$NSC)
CHARACTER CFIC1*3
DIMENSION Y(*),U(*),PARA(*)
C-----index of the current controller (INURE),
C-----current number of calls of the controller n (INUAP(n))
COMMON/P1/IFIC1,INURE,INUAP(NBCM),IFIC2
COMMON/P2/INU,IFIC3(2),FIC1(3)
COMMON/REGU/IFIC4,CFIC1(NBUM),IFIC5(NBUM,2),NSUG(NBUM)
1 ,IFIC6(NBUM),NKSU(NBUM),FIC2(NBUM,2)
COMMON/REG2/NUD(NBCM),IFIC7(NBCM),BFIC1(NBCM)
COMMON/NGEO/IFIC8(NBIEF,2),IP1(NSECT),IFIC9(NSECT,2)
COMMON/RGEO/IFIC10,FIC3(NSECT),YGEO(NPOIN),FIC4(NPOIN,3)
1 ,FIC5(NSECT,3)
COMMON/RSIN/FIC6(NSING),ZS0(NSING),FIC7(NSING,3)
C-----Equations de Goussard
C-----Affectations
C-----PARA(1) Modele (D parmis LISTD) de vanne AMIL
D=PARA(1)
Dec=PARA(2)
R=PARA(3)
HA=PARA(4)
C-----On met un filtre d'ordre 1 : S/E = (1-a)Z-1 / 1-aZ-1 (FILTRE = a)
FILTRE=PARA(5)
C-----ZD : hauteur du radier % cote de fond
ZD=ZS0(NKSU(NUD(INURE)))-YGEO(IP1(NSUG(NUD(INURE))))
C-----Calcul pour cette methode
alfmax=asin(0.45*D/R)-asin(0.2205*D/R)
alfb=HA/R
alfb=min(alfb,1.)
alfb=max(alfb,-1.)
alfb=acos(alfb)
alfc=alfb+alfmax
UMAX=HA-R*cos(alfc)
C-----Y(1) est en tirant d'eau (mode TYS obligatoire sur le fichier .REG)
Z=Y(1)-HA-ZD
C-----Initialisation de Zold au premier appel
IF(INUAP(INURE).EQ.1) PARA(10)=Z
Z=FILTRE*PARA(10)+(1.0-FILTRE)*Z
C-----On stocke Zold filtré dans PARA(10) pour le filtre
PARA(10)=Z
C-----Regular setting
if (Z.LE.-Dec) then
alfa=0
c=0
elseif (z.GE.0.) then
alfa=alfmax
c=1
else
alfa=alfmax-asin(-Z/Dec*sin(alfmax))
c=(1 - (0.45*D/R)**2)**0.5 * sin(alfa)
1 / (0.2295*D/R) + 1.9608 * (1 - cos(alfa))
endif
C-----Calcul de la commande
U(1)=c*UMAX
C-----Si plusieurs U (SIMO) on duplique la meme ouverture
C-----(interessant dans le cas de plusieurs vannes identiques en parallele
C en mode V par exemple)
IF(INU.GT.1) THEN
DO 10 I=2,INU
U(I)=U(1)
10 CONTINUE
ENDIF
C-----Calcul du debit
C En réglage classique Zmax=0, en réglage alternatif Zmax=Dec
C Kj=1 + 2.2*zmax/D - 1.6*(zmax/D)**2
C Kq=1 + 2.82*zmax/D - 34.8*(zmax/D)**2 + 164*(zmax/D)**3
C x=Y(1)-z2
C if (x.LT.0) then
C Q=0.
C elseif (x.LE.(0.45*(1 - 0.7211*c**0.3)*Kj*D)) then
C Q=0.7096*c*Kq*D**2.103*x**0.397
C else
C Q=0.3113*c*Kq*Kj**0.397*D**2.5
C endif
END