From 611f46e0723aa0698816178b15867fc7fbb462cc Mon Sep 17 00:00:00 2001 From: Sebastiaan Mulders Date: Sun, 5 Nov 2017 13:49:26 +0100 Subject: [PATCH] Complete modularization of the DRC code --- DISCON/DISCON_gwin32.dll | Bin 61608 -> 61082 bytes Scripts/CompileDISCONHere.cmd | 6 + Source/Constants.f90 | 5 + Source/Controllers.f90 | 259 +++++++++++++++++ Source/DISCON.f90 | 416 ++------------------------- Source/DRC_Types.f90 | 204 ++++++------- Source/FunctionToolbox.f90 | 140 ++++++++- Source/IPC.f90 | 183 ------------ Source/Obj_win32/filters.mod | Bin 759 -> 0 bytes Source/Obj_win32/functiontoolbox.mod | Bin 747 -> 0 bytes Source/ReadParameters.f90 | 127 -------- Source/ReadSetParameters.f90 | 340 ++++++++++++++++++++++ Source/makefile | 2 +- 13 files changed, 872 insertions(+), 810 deletions(-) create mode 100644 Source/Constants.f90 create mode 100644 Source/Controllers.f90 delete mode 100644 Source/IPC.f90 delete mode 100644 Source/Obj_win32/filters.mod delete mode 100644 Source/Obj_win32/functiontoolbox.mod delete mode 100644 Source/ReadParameters.f90 create mode 100644 Source/ReadSetParameters.f90 diff --git a/DISCON/DISCON_gwin32.dll b/DISCON/DISCON_gwin32.dll index 9e37f00c08f6873aee94f77b553a7e43e11fbabc..52e897c7900aae079305b80cf2e1c572b05c0d1c 100644 GIT binary patch literal 61082 zcmeEv3w#vSz5i^oOO_;N6InEB#03`(HiS*G2>}8mY<5W?FOvX4(9LV}&?Sk<1|m`m zY!Z!el~!%J^>Ojh)@!+yTH8`f}q7YG9n1ukGxdtTNTrw6wH(OtmgkduNNOxy6)UT47q%TJK7|=%Q$|$hyof z2t^u$uun6xT=UgAdw&4^8TG=&lNpB5KO4+d< zceuGRLBs441QVwt0;0#kl7oVu2wdZ+7Rstq>uWqUxNagH+Hexi4LC!NT%md(wVfk&$0B1r&N!Sm z;S4!)1xxM>D+fF^9_c2W(Ky{4CjAm#YHddc-RGO|vju0EV?D`rakh^fNIy%AD{vh9KdR;kCQm zt#tw~LzE@+U$+Wx!HhKSm~6^2iTv&wjhDM%M*1uU97#fA0p28}caO%4UogX(#(*<` zB9w{zo(hAPnaP0r2(MbidpZnW)=UOu`n8F8TSwzr7APTKe^=o>o!CE7yJeuY=a6UC zV*)OGRp(|MS#yQ&S47<(*)!nLj~wJU-c#CHM^>NeGnxapqnHMDL`Z0ZS~DeqRe*Z; zB(B)K$@knf1X;Ee65{(;n$7(lv#EcbIobP##xtwWYVKQS*7vP6YXg0V9q;Mbp11n; zc23}O(gvVSRKPhhGSYX(yC;6dP0n4hg3oge30e=ffzuU&&@=D^$@kni2#6p&LCQb2 zrw^tKc#On;FyE}DxPY{ggAfszR0h(2gFGphypG&}PSjzb^SgoJuSP~7Gi`uN4g4IU z1_i2DX#>NpP>FlS_aY!rJgGn5Y=R^`i6@DqQ=irx*uYdYJTTzhlVteSz~6>;`&OS1 z|9ap9b@!&wV@cm%mjzo`ZbCY(?HMq1{~Y2)aNq7g2G4jRcYtd#0oZ*>W)ZIM z+e@lV6M>T$vi=F74bJ-1lR$a;Eg;mI_iY+z9eMMy`VQvC*?{DMDDX+SDEE1?+B10(eGfjfP!vny`8k?bGTk<38tK&HNj|>Qe+9!Ez+K9a?!bnDR?g8J zc&}6-J|?0fJ`vzUxQs(SXkP~P3N1h&JQzaw^HGHR=|10_1VUma66(EsvRCXD!|R1> z$^fRrXBx95uxQtM2 zdh#H$$qwR1Y8inf`NWKC-4M(D_@dNb|4MbWX612z*nsJKgy3N#4 zd5KlGnff9y;oRa`$~i@V6G-)s@=Wbrzz4U$E<58IIKtCV?xzX$!ni*PYT^zdlxmX^ z1;o8x@o(AA;R&VZD;zJO>3caB{3w?9rE^tI}^d*i_QMD)85Q&|A$uk2yQAvl!P(cYjFM zafft3s_5=u;fd~MobD%1aDbPmp@x1Wb6G2MxtyRB0-;=rK@FguMrla34~TBPh{Lnt zOySf7w3mphY+sqidv-}*r8(Pk1rceHNsb{k4%h_YQ*hUK{y?|y<93v7ACm}f+vBoU zqF|*~zXml4K#4G|wuxG;;amgvA>FWtDkr)|H0S7#`|aq9s21GuE}UEDyPaxEFpNKw zW&Z-aHoO!d)vadH!IyD4K_yttNEel0R<47U=$_xixW-pOHH#aM(L*Kb@{0t%zngPW zum7LfLS}G+qfR2(;YGkcG4W z3022zA5PpCk<*dp6C7y{1~?3?O2`Wo%f#|&=t@9Y!!x@K&%D{UbDAF@{j`tzt!BTr z->>~X!Y;+@e?`=D=%yTrX3P%FEbntPG+WlILfvVr6k+9N~jl83365 z1Bi-9TqU!s{X+xnE@f!{&;i6p%<%*MS3s{F+J%O4=g=@FySZ#^83*APf%FLU-9eC zQr+J<^!CtS;6(lLJ@gQCeD~1<=#owz-+6rJ(A;xWteBHHzJGVarv3zOi&N3MYCJN>&5Gy_AgqnUg4_(}Y}bNn+jcyCiT%2^m0Ln`)A zF=(NPKUxBkw*S<>AYbfO%}_BVSeD7WEDy4>gf2#@eysl~so8Vl-&m6S6ic@J7zx#j z0g7&3I}+uhObXTu6T4z10Dm7qjd&TRz;7s{{V?q0(4P_PT6WFKYbchUOm`2xhTr|{ z7p>FKD|AOqP@xkgRkRa5lJiT&B4#lj1D-d!Byj54d<=n35NN-Mbz~fY4l^tf z6veRUDKtn}k9lan_hbmxorf9wp%7ko-hpPI-+LCI{{C~6k$;!g|A8oKbzIN;gHb4o z5kt>s3^AStbp1yBj)~}Is%#ID>pjI(p{IeMrQ^wR?^XhNpQ2yMxMJdMlUL>OseFO`WH6(D8nar?q|1@u3Otc?KPSHMP5L9dvxw>FV$_^f|uV>ljJJV=6r^ z70?5Q36AsbNujas@RjrM_JPK)A|N5)+}j6(w|i2VIyPr9)M;)h zYlZ}~&wd12c}%OIOKb;TH4QrAT_(3{&|&Q8FhQlnrkWPhpd+b_(V1$x%9;`7dY|JO z>QYFdKFmC9f+s=>apN3F=7Hq>aF&&=CjHByy$80H0NR@!E)Z17(#N5;l3<9 zUP~mWMUr$z5)6fnfBjOHGR<&b6#Ooc7`wTYksMMi7%dwQg^)}!+(*GCk>rXb2N+2( zpf)}eLNe8G->MMB>P3=&VI;wr+jvh1NwVR-*&$Ny5lI?|L=l!7n?tDJ0au8i=R_(C zqY91OjgAm1I3Ybm)=`n_^Eo6dG`Kfj9YO^S=OWf6Kj?e7>-I3J(D>gd1gRusi-`OK zAcxe9s@YUC4Bgd1g5kIExXa>Ij!(}Q6ra}_?j4VMwf$BWpGzy8<%myD<*g>fXR^=! zCgO8dX9q$tmRFK^dLN4TY;1ScAU-{=c8X6A;&jkqGJ$^3VJ@!es^npMU1*pNL#;6z zKSS?D&_{||li}Xq>rju2PxzLa2FHQzVyIQdDvF+0lvtOk5p_4-RpX)h(r)rJK_66? zwa~}a-kOT(O+8-p1U0}I1vNP_P|MU{P3!&Cw9??kpT9|q2X6$o)H3$|1XBWBoxX!X z4GO=p20`(trvbcr^LsW!5_DchWlPaQy~>@m00hdXH=oh`Np&$m4v)bMD{Uc zbwj`L(jdujBOD#JM{EXru*>LX;^gTZN5Hp&$R;5VvFuem|Krh7Mb~vBh4(Fo~ z65)JALsAO@0Ab*?c?uB*EfsC4JPeXUz3@&pc_Fkv-MHMvs`ut_I<%7x#}WRc+P{cT zN&c0F`wG$0NM$os^rRg`%AsNJU;a8WGQ39Nze!**{5$epnDQsD4<|vc?C(SZ{T}r! z9Dab9tw*Mezl%UK$H2dr^of9Q8Gj?;ABJDU@n1rj-N=;j(+M|FlS&L-`IW0{(pD{nBLUZKLNU z^i2-a@Y~1&O#uxQ5&8L5=#r&Dq5xW-{jI@qw$^rZk6Kr|1GPu*y48djk43MFt}`Cp zX(RPth-WmU=uVL&@nN57E#i5v!z{LL*#BDGum2Hlg69o^p5=wTAXT4vs5CI|7?qSGi!GzNUmm*s*m@WO_U| z-N;Q(h>-@|Uf}gW%AX)wEBAjO9*hePcbY{{nyAN^;hU{@U&e9q2W-q%X!mFUe`KSV50TpTKLx{eF-CY`%QcVXn#Xg^ja+l-k&J6@ z61h(f48s06_dq*Q?&rFz2GyBKUqL<53;0TAdFOpUQIrRP&dcPyqzarW1 z$o8*@^*f6EE8_i*mHrh*zoW^&BGK>g_*W$P9bNtvCck5ye}&oKyPn3MD^mQu-R#2R z?`1>x6;^-m8G=B}v){t7%a%>^_czJ&e{A^-Jz`_wN__uw=JmKl-AcS2{VLL=cu*s? zewQA7*&+zbr5u#~6i3OyP@7Qdf6DKx@i;^G@dAG_Z(NT1z>`-)I-?J=o*?3_=PAeY z01w3c#DyEWq+w?P@lcg=4dVv56ne=99_XD65){xR69?l##J9 zj%M@kBAfSl8fw%P;k0=bcs@tsXf_QZn;)y$B!y$s%h<@vB}5lrC9-j<*<2Zp&3wkj zG+LV>6ZcJ`nhmYO1-jBVXMsRY$k|oY)Jd>NRgJUoi%iSKLIOF_s_*xdA{LaccvZc8M!0Mhc{U`diEQ7 z3Q(2%ODuHh!e!5)zJrGDM+|Tgx|COwX$q^)VKHmo9gyDq)Id5&u_*CXx#58pjYW{;*y*3ZX7pb zoM=XYp}UwYyJO{`ql}wT#r3b`Nl_RxK&v+6BqGJYC(N03cs{1U0sjEQjLc}fGj zAYeu{nenZl88bvP=3mQ8quS8jOQ1V2P;cU9P*aW4pa%ma^28tJ%RCMBQ+7B` zPcu$0gyEDca$2wEG%FmZYQ||G45v>qR3PiJ)SPm{ak_+Y+7pJ;3nHh}JPoyOZa7W{ zfamMo7lu=<$mwx4r+MKxZDgDdgyEDZa;j2unjel+Hsf?K45!yHG9y35sX66_60)z<%^s=Nc(pXM#zZn z5xglpPfxWPuO#|=&w_?}s%#+ho+^=ERiLLji%y?y=MZ#9EV`w^Uh159Po>8d z{sX2BUY)~zn4W4ep^ff^_e>-okpyI+FJ!60S8IVm(fK{3Dhugd=y@!BKt~qRY&BWf z4hy)24#3I6OGFC^!_a*LL3fyOJ=mLwUv8lZSNIQDQt*o2Lx!=igsUIFo~QhY^c^Ev zND{`fa5FH-!fiYajfU&Paoof>CUcGoUuB7$8q}Q1!*OyjPAOqH9mUf%>6E1AR1uET zIpF(xEnzr4E^<1Aw0{R-ghXf#Tz;PNI~)m`d##|M$Bt;NRDRzFoXT%$5|wOUFJTzE ze@M_B1*AhSJ#O&wqlM=94_I#SDjQc}%5P~n%cd~PNWxgJz4g$6bL2rvoVrC$dw3ekX=ON0 zPclxl8e=r4*&?U4YECzX<8%Y#3?8+}I*XdqP2o7v6Q-}XEexk0i=5s?>eqE+ zYiOu}4I(on!qB~955pP>|Dzso^Y@-Wj-U01w5k+K5g%R0;4~%~+J{wzk9M*97wEos z0Kfi6*8@V@-wl}?_K;Rs@f2L&v>o@t@N!Hq^HreyKAt8re}x=zKgHmv`#WedJS0E= zW99vEdchWwzYWv=q4(AFRw*Pu9q%BqKeM_QTQNV!6Ctj*VYV;uE?uK<{FqPuY0pur zHOJnGyj|h*tirR(-&3TG-78|H*FVv7149 zvpRA|F}wMkb4zD8_oD2PlgDm827|Oi!{dN}#0x=4eh6A5SMj`6z{lMz$VIY)V}2L1 zaJL&M;)_*r{pXTp{tswQI;TJ1T-lG!i~VkMecxWpZJPk8^51_Cf@}~SFM@O@ph^Dw zx8NFiEB*IBN6H4?&mRWFaNjlpW1%ogNE<+r^u^LPuav-9 z8foByc8?PPTR8*9VQkc})5ZhyV9x^}q+_t>Nq!;tAGwR&lSP8C{zW49+f)3HkokVQ z#sA1AT)?J!a7e;7^eSwqs>F+ES_NATmv!GQ9NP^#n&77dn`1lUunuNo{qL}#_{wVhVtZG< z?{=$?-mOED&z1!1V%fGau>a9+EZ*3W3}F%I#lQk@ZfH`h43~oUj#!OnHTIE63?~x9 zA--p<(z`baFY9UNYdmBW(XLE15hNHO+Gt!9XxxR3Gq&uWfjh1O-%4zG^RFoLOu%)O z_iTiJMWrVS*LVX`<=x{Z@5=8lunu#mM33BjY>UuwbLjUt?nnzi(}iurXX{b=`Oa8h z1^F`mE_fb=8B5zlvDO?XmOJR7m(k~lCpQG{ISAD_4cwazb#p*y6ba+K14N_+Y9eNi zk;M8OMk4gF_k>h3?-2)(g3NRe^+E6oMG3^ks_Q+nRM zW3v=rVk>SV@2CCmANb4KXcM@Xs)Mu!YvUWZp={qRdf#uU?)f)UKL$_{be4Tdn;Ub; z7e-$%>-fo!#y8mBs3QMEPl0Akd^d$N;h%)p{P%AqDSQPrZHC$J{|ok*4Sl*N!vD4> z!hd3Dv3|%)RoH(&`+eKrn*s(wjjuvo82S};ss6WzV9!Hj58HbZ@_vrL(dC@|Y*Y~B zKZdm|Vz`xL$iq}{L2Ebc4ccI$l81M=&5R10;83k3|47t$}E9`g8w11J6JA! zAxyBnrq;R$Y~6T^*n}po8s%D#IP|sX2k)d-4DHwa!8?z#E9OmvCcQV}LF+Nob?g>W zT5E<pPl6%R zB<2mbnJPgd&Ul}<3wXGS?dQus+2W&MK;M53^n8_SFra%$Hk1oIvkz^kcO416DW8Nq zHKq_U_nQgXA_`dqA>QBa!JWo`6J8t~A#vW*xIbm2{bkr^7Avo(LM!JTLY}l$%)mBT zu8jp&G_A+r`_cNx`IqLRjZN@Bx|tZogv1-_YYtzJxncKVWFj6|Q{!Xv$Nm#6PS7&? z9J$omkRN^x6~WN|ePZVSS89K~dx}^?$ldF^e{?ah21+UhY&QF5zL5$bb#ROVP_SpE z{P!2Y{=H;ff!~quXR+efZAUSqy1z_&7LOglPLkDV&W`#0x)&&;p9KW=>Kr>IW;gO| z2mE>vzF~VyR!;>1TY`#BG5d!dK%w&Sdjr52TAed&1%%3X_*=LLmJe?)seV<{NCdIiI34{x6ejebHx3VSep{};hmHHC-5nDM+VcW%{43nP+uMAd7N5O~+P`8YgJNGgja?cU6z}t^L!s9bv}Y6JRxf7wZg=~BNdg$%3}EK1Ywz#GLBd25^?|>l%wJK7j+=Gw zhI=W1eGgI>fbC1Z4iA{;YzL;```uQOK!eXpv$gLfzpV{rVD-MF8FC!;PXeAE`zWBA zHgFl@h4KL95XuZha;V#XlhIc~ooKK_^6u5r$PfmDs-BRMW+T-G?sypvyeWQ;+nm$s z1GL)jr54n`=}G$a-b25eNn5{{mBVII(!ZG;C{;k(pWhRO5E7%)M>WxJ%ja8{HolK1 zCw~uPvBEtZz#YW6QOC9X(-iu^4Ycm_7#Q&J2~?v9%sU zdiz1?vl|NyI|JICC$xc&^nx&L4_Yr1q-t|tgfzI{emF4g7{&BSlnuSePh5{t7g8f( zAdTO`M|!dS_QOPQ<4INmb{Y_B1K-AjY5zn*?K#v{*l&-eoq}Dtpoq`ibV=WD+3CAg)AMQHnNA&xIe4qD_(0$d7~fxfpl=XVr1$`L zD0>p2fFrG;5AAOV=qSX$eQ?CHq7QES0VrpEjOU`E?=0jYrce8nkrW>ev@^Dc`^_@BOOAL>xN(8bAXZNXtA;q_8 zWa`nk2U<5F9TL;SVtPf?#bNKB*pJ!wKrom|?7yy-txqb7$0Hxr!<|Pz!j#Ic-VgFIW#-;ov^gGPDS+X_SrQ3#va1LPVh6Q z_RIbAW5NJQ@I+zDHnsQoeiQ)Y&QlMNj{N&uX7AZZ54LA};+g*5JxOE&;bWwpOp)1- zI7dsQM~nor(X5By@qx2A9%JjdymNCPBmy*MPxPRU-9RFqT|5_~eef*y?lCe3D{g|E z9mX^Gw8;nbyYs7vo=@+1JMFDt|I@F1gfgI71U^BwTQlo;=V=8sgd`(M*~2rdNu(&0 zydO56`oX5Of#6nkRwPV6crghI2G(GxIzKp zdl(*lzl~xa-*Q`_?|G~94-yv*M}}5dNcLDFxPc%_e14$vB=urpzZA!tXx@IHZ?`8G z%9}(seUZ2^_t~_8w~_srs8k=*PC~Ki0}miGZ2-El(%()~qNhHgf)f2t@|X_jsZG$s z06ct&E%>}GDZLJS`X@@JAM#w1?p9Fj^mAJ|A z)nB0RTg0N`X-s`f*dO2;*ryJH%M&D|@8GGoSUd&eP2*YEeksDPbG+m-mlRV$C)Tno zL`>&Nk@KLA{0mH%xa%of=ct#5+QoPY)jJQ_ocEAF^!Vj+&OP2M$ zzO!iba&-lboZ8_X9L|`bPw)+yNN*i4c$8c`YmDi+w^R- zVFPFH?$`}q5YXW5Ba#XIw`+aR3Iee1reDLaD)hcP>B{?gyrKI+;CMeb8M;ri3%#NH z7r201v@oXNZdUqan!w_NN#pE0Qu`Nc`zLsJ>m$EQ+gi1Y&s_^4dJQ{xe17g_FbLh5l8(yY&9@9S3yg_VfBJ4M4 zd4<4J2wC>jj%sN5Rd(O!-Gux)y&swOH$v5X#*oS3f!hiBI5gt%IPJqoeUKYz9k_uI zeg_09ZCG!F^?|=)(7@`%$ka_@|8@^cuP>WRAL8;AnPb~0ZWyL(?|5N4r7UC3scAa_ zs9KI@7o#q%f8#I0yXakfva_{Cl($|4s^#qr-lk z2e#s|3P60RXRW?-+>ow|K;mZtDHKKSo|zbh8UueJ5fC;gOy}j+Gna@7X)SRM&qDvvAouc={hm(5>tm$SL0p%u`2S2#FhaP#!EEKa*_pRZ)N%=cev`NnZ_r#*VF?W_;fX+gFw^89}!5|KGhmk~?v?ViqF6jYU(NIxKg|S?LO8S-uE)hxR z)r0mrj;SYtwwyLA(oV&(Y@sg&e4V|7L!Y~med2b}w{56wI6d;9FKw|quRx>y}w^0mgz@$UniC}`q3^w7qgw7Tn2E$fBC7LThDm9up5AAqp@YGP#?CxxL z;nQ9Crc`q~njiP|oi4Nkl!4ku-k(|z-+EjfCiJ&8^=MtttxRnLd(bRtxuvCbO^d0m z!_`^eitoCyrh_Um{aTDY{Md96#r1v>0qPd1n0e0`zK zHrtfEuw=_O~9GXGN)tsCSq~A|E)M<`Y4>={wH9g z%H|vC7o?@~{^#-k#sfxCtY&SwpjoHaToJ3)-h~Fo1UO9~pTJ{772)K?>AexQI8#J< zoIY=N9JNOT&&8R9#|v8dLO6Fui8JDNWprFTXkzpBf_5-YyD(-g zgWtq{i;!yz$+g9R73uXb&Mn|Q9jO2|^^i*N=WuQZ{Av;I6ybC`+Kyzv$PU7Pn-HL5 zPsDbXS%|lcalGQwRv$yug(#7Ff#Iy z$mlpDJQ@#=#=@ieyqDs%_L#``HDTa+q;$A(Qn|eVST%Jn3&}U5 zsgaH>zI^-j^w%$E55{;A&(yI#BPxX(C+`=FnMvLyW_lU6vIF&mc zqkV4%jC`LB+{MIGjMJn@0lr|7x|n2}aOUHz_Kl3NxA7ugV#K>5nmva2G~smLJ2LVJ zQqC))G#cUF1|}we?!tK*`O`vV86yHB4;$mL6*y<_Q}9?wJi5SREIX(yNPhLEk&$BT za0y%Xwz%BL+atobcTtQj474yNGAcR@+zvgcUfhE-_lc1ayjEiJda3NakyPT}L0IjE z`}R^rw<3Ex&L@HU5!PapZlU$u5w|`vvVX#;N?sUK3Dh$<-M<+b`Ge>$F}`fUc(%or zMJ^5(wgqF@G`0_OAd7r>^!bsIFO>4+`eTb*uZegij8_)Ll!107&WB$Z8DWF3aP<#m za8nrSg)xynjk-<>q|-K>i91F{mf?lv80A_WsSQ`Ic`-y!KFr>UHkRsWc%7;vr-r8= zE8bGTihMPI)kV(*eT8<>ac5+&)~k!WPs5tVLDo3hW5(hWFY?)kb!2*6(=iq0?dR3Q#cFRICF3svBZpCk&!$(_UohLVjW`a+hRP5JKi>RfusqvRGf5-r3c9&J&J%^ zjZ`V8C@!fm&Qunc1nq8SzqXJ{NB$vu+d$h%wAYJ1>4-}zfH`qV_BazNJ!|3MAF_wq zrUBqTMEH^q%YcL4Ra+Rt1vYX~B6M-g4eXw^2m3%~!U|Em!DT+HfzR?{ni;~4xL!ak zpzT=plC5;E1CK?-L*ZB6w~?M(kxwyA`KtC{VO(`pT%|3px;acOTO3m#SBAU@2TQ^= zdhi9s(vvu?SZ4bQ#rDn0_R1D#i8wi`wOVwZwmA_j$gX^>F8l&1;wd-c-rv>K#b}Hs0#cyGS zt+d41h$s7ubT;YFA^$zP{gdl|uLtT+$o#irbpZX{+RY)b(|}b0b|3^Myo<2|U_&9W zm4IyqEZlJt@N+ng@8D1Q;!MZsz*&XUgYzDo5953q=L5 zMjMRta-3$IsW`K7(vgp}TD@2@$w2ZmsWhy9zLU=62sgXsln$;?>kQQ0gT2n8|1=k0)mF^BOURvYsbP2!HiZfY|L z;TCNJY{s81);6qxn{f4nc9n?rj25%ib#1o^!Y}Bm!v)z|-Rg$6_U0B39;T*_`EnFp z&ls;$eS)(Qvv8jFniePK7BMmJ?C^kJBjT+zLe6XLAS=l4g7EzixH!`n%>2PvdCnSW zD+sSd$~-$f&N_E%2W%9CH%v0nM0s&Y)afIGnmByiTd3m)b$R=x}%_QhLEc zcEy}0U*ZxM4hJ4{gFa!bJc;-OkiNmW*j`>@FA_UlkQzbzP0m#vtLxf5tjiXXK}%~0psnP$02Q{=GJ^csikEj1z`jm_n83!{y%(%|F&zh3CEc3q1$1@LQzMm;% zS+Z7SJ(Trw)>l~xGxKLIn|bHVBQxKhd1dy}?B?tZvmThWboM>7|8@4$vwuE&*X)^d zKANMQJ8^ExT=(2(=ElrRoR>O3bAHZz+x(UDPtJdT{-^W5oG+L#UWe-@S&}UomU)&t zEIpQcEcaVpvi!*suzYCw!ZJQBAuS_qQQGxsE7R)Iev7S&3m9EQ}l#z_rtg|!LWZaSQWX5kZc4ZvT_)Er&nYA+?oVjo2 z7c*0{>$88DeK7mW>?5el2_Atd+B#nDySQwAo#=Uz|N|PVt=BoZsjCaPF#kh4UYp{}-aq<#Vw! zEVo&9TP{jlmiAcMskF@W@1*ZfpOA55#!ukcIo2NQ>(i#^J?Zb&3ka(FXru-_vd*D^RJ$tKY!_b=lomd-#Wj0 z{=M^`nExh~R5hOWU|SR-WOCZIX%<*i1$%Bu>q_&bJ(~6;V&!n!yJ@G>&ZR}9U!I

C(MFd46a4c0kgR`R{x9&r6R<%ne+`*K<$rnIKx=C14F12P zy>13p2=j#&Gx+1e3_P(mG&jy@&dAQ1QFj|n1*TRAGn%l3ZpJMxtbdt7RI6x73Zt#< zY^iEmrK7pA1&ekF7!eWG z(S&IOlM<@eVb-J3#C{lBtd?EZw{~LLo)DoKM~}Zk}TZa?z(=J53Omup}ybUxYo$`?7#)^0$;?FIBMmt+Ou>47g`1~vr zg79azUj%Pg)!d92iZ~|&u!0kSF9{IiY*^J=k5aK~NLYk6iWWzCZfkQ1(XF%`vtEeI zVd~=!|Aol8(jBQ7IZwKSr9xyL=^Ep#uer_L+}Pw1v@03c7$+?RTT3%E+8Z_WW}F#S z5s3nO3pQv#<5|l}L3@+N25^kb;@g~s%;65kAx0wj4kL(as_AH2UE_9o1l@IHSX7gz zrq)SNK{w-9U^!0v0)8)g06(wbjH4d_-#}vE$gW!jT^4QR9H;#c{9Z(+eSr&x%JH$V zPUPxrYel)k@w(ss2sk{$b*l{CaR}fDjXoF=Eb8&6n`nNiwar=MZWOenss73Uy0e6q zR&mm;6+LtU>NL#@3Yvs=SDUk~-PO?C$q>DZkx zE}^r9E$OXyQmzn{u$nb8{7#>MeG(p~+jtZ_^wLM_uz&y;xPKqXIIRQa7x^VFV_;gF z3;Q*gNgcwt1sXyW!K-EXoeZuGnvM&s{cXHWh>hC47VLEt{!x2AhBE{m`6cQla*0;g z(uz6T#ug!JFO`K}h@1&*nGaAq1Z&j5cH+^=mkBy4eub!A#6ZyTo56iVCg@2PU9zBy z6e5Z4Vl=aLSZ>(bj;YTMPlu4Gi8}Bjy!Fyy9ikBRDp$(wsu7}IqiZcHeRG4yiDk$_ z)awM%ixBifsFhI%2^Onc)#>(NvkK;ZquzKP?BY2D%MXR9-xG>iSKqw48C5UpFy$xf zR<*9~fE(TZxL|xb%-frfDZI35d0^CtkgBMIf@)U6)GEXZk=M$k( zXNdSYojc80S5w#I5=;9df^OHj(-qK<3HpMLDvt{b;GMOYX1_&<`h>s(2n{Ntnif}U zX9ufbgo7PIQGccUJvzK_XaNxgg%I^Og6-2au$VYYHwSbLpckS(qnm@e1{9_c^*P-f z23jrnd_gxux`w(Y+?=DEqq>HUTQTP!^(Eb%m>|yNM+ur?*n5)hX+B@P4@~Gl6HwC* z=dZ#fInzUYbjNTQCE!D1(a_ zLcEB8Ei>Ve3=GMr5h7nj<<^t8sCKa$uN5Mn#x(;{fH7c%5c#46q(~Jzq)7UM-{Z93 z0h5^KuIwz)x(ON$78sAyQUk%^)WV3FMh$@{GISaTB5LcK(YEkrKI#_D4S0Y_r~OzcIOya5lXCzvZBID z;#46-$~h9PqbP476H`%{S5zcuWsKs|rB#A1*s zTXC5)uO#2OpggYxNEzwqS%ON3!^pA<`txH{l#-Ga!uU~Wc|}V?LLT6q%^tedIDJw56 zDqXn5?p#tq*>IU%h?BGJWreD2gPcwNDzfK+CSC@XEvYDQ7BgQ92{M@Y7oo_JD_kst z^Gla3kiZjUu!LP&CKzR?gbe6J8CsT?@2n^&Ew6MsmXz2k3rkCcOO(0g_R1yYB^AP@ zGT2eNq};h^2@1QYu(+^NNR+csT$Ri6%Asqy-C;+hkY^UyOYDw98)90xOhzrTm5E_= zxeP?%<-?Qtg^uOI6*8EdRaj!nx7&my8Hy6NFSM6C^BthDvjV6Pu9Wjj%PO6D3zrC! zWKiDH(!zZ5yVFrx?##>2FBh(ob1RBUmpKbj?)WLrNT-qzOk-qc3zOv>&KBs-(lW$9 z@wTyem?Gnd#ZpmO4xgYDu9mYXf3WLNLqtth_O;ct$9OZUDj@OohRpr?ccJ($^xqQ z!hAU^uV`7`@(QN|off*dJRl&bgT?U6Vi_T@`g_U-}A!o8Cqd3jE zv@p+!h+HV=l$2J&yJdNW<$35RgaSE78UfH4gZJfyA(^6%*R$3@SftJ=DCA>`#d4Nt z5V=n%lEFMwi%P`~vO*}9^C&p0Yy}G!2qj7&EV2ltGK}ot5l1bTP$uW_5EG*j-2-J` zKRVlyR|H8IO0kA{L2;o_E~AKJ6~t$LA?pB9G}Mw+$SB+oEHZ>j87vj$^`bMD$XR*$ zi#W0c57q07V}CL2pOU>`ho z(cw{|gN-BkFq7sNM-4MS*Ti92JY$_#j*%F3M#}p#a=bXE6j2x`djX~vXw1fdRa8^> zz>(cApbrCx<3t8vI_tR3lOR$YWAkA#Q%u%b8$n9Nbb)#BaSV@kLezu^GWil(iQ{o% z4<;7fxLi)`hIBbO`l<$3OC1o?B5ntVaoTeHPQtYEdr93o+ z1PnvS!hA|15^R)2CM=*tn@~|msV*U(rQ;INy#P>w@v5z;bio2UqJabS3D@(B@!(wI z6mOyv3Nb3L;yE!1MfQbxXgU~&*o4r|4dX>N$Sk8wj1F;Afha}eK=+0O2YM56d>)@b zeV(JR9Dzfb3EY$_);}`V7c&fYIU&Km6!M9J-5V1YaGn(y-C?{d=UkFdTr4VjsiaGg zEHOl(An-C#t03TV5kMnqChZECUnz@73ds)=U71jvSFza1)MP}HIDlPW#dKtsrUVQk ziLKLyV5&e5FTA2?uYs(p!Q54yo{%Ub%Q`Bno@aHe zZnrwsgf5oL#P}sKF+${i;zmIVcB?vQwsp1B)rHwATBBB5i>Ry*ia||m)D!VUc?OBi zo3cs0sK4oN#ot3qKto}}cmjh$78->%1*{P7ED)wLQlS-qh!6lC0V@la6eD_Nq@ujV zcuXj%L@did5fK1{whC;;yq^p#DMc9OA>d_@4a1#uJdnwd5_>^j5eC8ILtqtn!dt$e zyfA;E9gkfyZb4qfGCKxwF(I%8r6pJ(RwTrRz$yxhFs5O;UlamfTvCMkA+xJkyqs)w z<}IztOiM>qk+W!^RfNYT8M4q$4SzxiWN{Xrlw<_-O6Y@x2_c}x*?2q(!6?pmdt(_&$hDZbvP-zK~Uc@txGFM_zUX zQ1tx8i2ChFgjS1>Mqou~ggz4bw7FefCJCUGUKYg%vm{7j5Tbai6t8vSMG)9w3@AjH zBArfYLf9$jlcC5?oXnsRA@b+Aj&WkHzO9bFO@WD|T4DSt4WR&B=U&SpA2LW(13nqV zUPeVmup< zvWeHwodI|5fAB9eVH)j!;Ls!p@#7LS7i%K*SLmneCt&got)m84F$9T6SVowIn3#FE zMR;8!6i&c|`~;xK#LSMF6QjYc?jj*lPd9}VvayY3oB{ZR7IjGwj{ZfXDV#vq<0nHz zv^AJD|FWFRP$VKiS2%$KVkbk?46wUsqUPGmIBA@Yi8ovW4)I#1L*WETi-csbOaj5h zOgo(fpCD=n&ZM0&H9mY-yake7R-DyeU&6 zG+b#JijGtA*TAh_vqD>xqCAx~k#)1cK1vsc&vg zZ6ZJ65^oOqb)$=3*z#+7d&;lqjVPKs+~CD4y-=FY-dokQcsc;XD@$5Ofh)A%bTu6; zC8^6Kd3?>^!|i4TK~plH_VL-zg?LuoA#w&`u4t-d2!mEf48`Q68-~vpbR=yCW+R7E+&$+y(dAKu}RRG@f7O8^ui zHT1diOGF4yQfu+7C>|7{ks5R$qmERza)j@${AvUtu0cyeRUb=(Jhk|aU`uM9v$ghS zA|qZQ4PXzJI;Ds%8))lwn#dYnqiKA1#A;r4?(s|c5y1wc@pR_VhA}> z>+v-o@t(=jprVdCB-?{JYr=OPzXshtoWvngQ;%QbC?BcqvM{nTQd1{qhFtB$^=7;! z<2a%FO*vBP&2?-?i|ilC^VF4c&)eu-I%|{Y`;2r*-WXv1A?XaP_>sIp7I1GUve}K~ z4HNEoB&cx+LRM@IT5+`sNz$3L<2vJ;Q_e~(oZ+8NBZB`IIv>0MofijnCM{Ex#SF*+ zoZQcd!!Uj(sxNSEK@syhNRX?6xe2b9A-RB*K_3}X56A#cx$e;4xG{lRhV&DTg6gM) zb5SUt*9iwFD`gTc{E2Xo%8*llSnzigW$CfV+!BPSOD8*)w>`vEfSf=BEJNs9>t*QP<+3;e zNG_&&WE}cyw_7mil_AmOs*6J*lL4s)9~mbbko9O}WX@%P44{-{oJK$>eq_kE2?EYC z75#PH^%#=NINJa@ffxERWFH{8aJCHj6Chj0g+e|BqzqQboY`j{Kq6C3MbT9Qk3*z5 zR?-E?+N_dQq>^=TG-qr4Xh_yM6xP@IaEnS- zflAha(VVTZqaj&`FF@8)Dp~m|S^GwFw(3VivbJ4-tZgb;xhh$EMsv1mM?jLM1C%C2Oln7OfpuSeiu* zFNciX9tKph%qm$=sbmFIvI6LWzfRT}l`NA=){`n(Csnd&EdF(}&Z=b5UtUxACv%HR z)(MrYtrs9mC>q^=C8}gStdez9C5!f_e4Sl-l`Nx5)@GHgA(gDI^T-0SuoHbwJdj z;HTy?B)UcgpkAsi51R;_>1588G=g+{Rl9LZfmewT6Qa^2_ z+R7`6IBxB-u+#Ti{Zn%Nw5HYj^AEj@IHkEMl z@))17a0qV#yb)g#>Yp?~HYp(GfLKZ;9LZzbm=8(%rEu)k#4LD79y6LGkKrQB(MZy- zAIP3c~d?l+s{b0Yp#f; zitr&O?nWK7m!<(Q?w|1e4t9OD5Kcbx;TSmXc3S(=&~_disVq1S=|i*H)auEH^epT_ zSHpEyx8a1n2e!9Wl2^j(6O_ep!rLxB`DS**bIMkNW;qp*%=IdqR29xT6^>bj)1|_> zT!rIN;l!(O+Eh4t6^>hl^A#$P(yk^I&L=9IdKJ#;(KuNrM?*4Ks;FoVS}DCsg>zsu zj&(Ci`)jB&%T#<`7|kcMNQLvX3a3DY^RNmhUxjmz3MW^E<5A({sBo%OIN2&3hYH84 z!bw-*SX4MuR5&Rr9HR;+S%q_MG$hkJ8j^K#G{jo3D$@`uv1g>OhVz~=uLsh_VsFh%3il6_=UXR`-1$EZs z;7%z}^{zI|=#-SxqC#-HNM~w6m3#>2V&K#V?IJI!adraQcAnlrY$4N=fv2p47cl-w zC;OCC7x0v1J*=X70E#R4%oI?r3W(J-8YlC_XvoZ!K`Mw4t^$KSA`Su_$9zD}NPP-} zlmpT&`iJA(3J9M&MZd#vHUUBh-Eqjj0%8%PjzfM1h+8ZRtg_&r8x^)IWc>~}QX9po z-UQ@Hg@$K|vjP%_hYzXs<5W`t;WJ&5hV%|?lY;6-Kz;Z$xwq}Go^HUM&39QSg_ zgMdh_7OW6Iu{xx%Vmokh6}GLMjb`K$u)(_S`G*uvd8WQM5^6f`o96PMIn8kNQLkaOXp2Ml%@O`Abhq^@>l{5 zW)&qkjUZyJsUPH!Ws%Jd?%8Cw>ImS7qnQv1WLSOka>j8LUTAKnOv zfIpxv%W@M$;jxVbQP}=HK%~B#bN&e+TZ5K5#m^pWQD`^-96rx4*?t0$lR+xdKOHsjj2aMS0g_bD{Ilpg7D$dPazAWsFUC@&un zsXyiAUI$2(0%tWK>lBb~K*-55=N|x)tHAjsAYBT`Za|VjD97qCK={09&@Mm*6y7l4 z!F?a7s72-)UJZz{ea;5NtdQjZM4B;Uy$n82!RsK{h}UO#B}Ogbvw1Zve7Sq2b>E$yU@s1M28u16|(jL(iDU^`KQl?g9_W<1kNT!d7TwG3xY{# z)?gWu02ZT;i!g;?0!*&&mjiMjD2veZ0oknJTn&g8WO9Ge1xUQYhYteMrr`4!Abc(k zF{Yt+?`nOX2M$4H{~Q8@4jJ-iK%P@@)?h7?G}Fc7Ck_xgw1(rPXlBt+%DTS zACNsLK8h-Cdlev5E9Ek+1B8M>j?4}~jw)nr24qM9`86P$6n53bWPal_v)#bR6bb^F$KEP1BCr|G&-^s(@mtH#W1z_tT%*^CvcLp$i853Em1&B&OQ zp@yduJQKRMsL->sRbYnAtqCutkK&4pUSSRH*FvY z^NUcMsK+BsOWRbmsz~y@a$6=y`m&p?A+;u>| zv0`wTh8eDpE0R=OyrJrr8r{%|EJ$WiD%>*I$_@%hUbk|f>k;4-jP7e z{c{{qiycU@@F&S3M!6gSKv|uE{9osK^@!^YPr!VvZV4Q;wTJpMH3HQ z0RORTu*14l!WKy6Peen{-I*%HQ5$k@H-G~vOSi2qJY)Cisp3|B@e54p?iY6yh#t}G z6;#KH>aCn%sVK}bmiwr{C6N#E}vRrK0*hYaf*nGZk4+}9x}E~?2Kue>=&CJj-|X!2A|w)ag2#T){HR1i_;;!DdnSdl^TUi*@s;hOhm7k z9QrAV;n4?VXMyv46J z+os*3BNXN|SK3ym=W2&Wp#=uxSRh6&$s=elUS6Jc*RW$$a@0@#Qy93g`e%`0={d?iw zbI$Ml&Tl)vv)pseos?WZAcPBopu{~gA_!X%(v{1e=YDh{JNokHqlL|*U!Sm5QT+OZ zs>O{R+P3!AhW46fZEa0UORGb>$gXYgY|%EhX!FY|wau+{_LQ+>$LK`TQG@tms-!mDgtGv*-<2yY9878Fh*P4F(SjQUFh>xER8)}AhLk!AJWm1Y+!cSr&Vl~T zR}y6I9B_!~TcFeRIds~-6}lwnw+cr_uUXf-LZ|Lspi_GGmY~eZsGj$FcXy8Cd<=t- zCc!%emwQ9{4tN@COnvYshm)HC=f@ps%f1c)F!O2Tg#=!GhSjyNJe zl&@2gUx3@lez5S|RDq(y;3wsh)sgEd5M>zb{IO^F`;igwGz@aCo=>ZxpFr)(FgV-_ zk+^00Uj_t%C-mj(wBV#B_Qa87>Q(AI+flA(P;=j)b60}qxxuepJNuWM3;)N!2kOoN z|Es$9BW6J|aMxg~XA>cvRDy2zvtTcRdv|&|c*dj3K`zBOVD~2IM7X+lH;Fb`1m47u z)sG5oFxFR&SjtncM?s}-&%j{o$R9VcbypPukHZ!g-32+f(Cnm8-lsj=2 zi&)%=Q&^ftA;(Jq9*C($BFe$_++8*o zvGSfTcnB0Pow!hR4g$}g&MuM83X#qxIUQ0BMPdTvSrbHOk4Q%+(y5cvp&>#>CqIbJ zK9SDJSzNCKIh~6`(K!Qrf4%mLbnX-B9ONNTy!N>y6rIgMbPkGiibXnWCgxzvyU!_&LNS`wwc^MXL#tZS7InShk)mApQ9ojhe&6$oDPkYGI{O` zqH|oNGexA+B&TzEC_2SKbUY%Ri&A?!I`5!_61_#I^thllR~ML z46%VkgIGYd5lyyxJ*6?HN`@Si5UNX?JUWXzhq6 zX*U`96j~9=>Mr%%+Wl@Uw6lN^X_pQ0dEKP}Om(-rm}y4{n(o^PNZJ)6Ii$+PlWR8% zF|J*{HrRoqjo;%usZigXnaNbjFKx3gmRMLeZHNL}#@~ z=fxaupYM3+uUB>`I){PhZ=W?HolcR?7C9Y2eZzK%`SG(z#zw2mPCzPG=CE z^&*|IBApUB9gG5UI>|wF9uw)j0xu!^sN{4oR?6vo1blzJHi~q9B+_}2hyM1#%tKD+ z!5})DL^`P=olZF&Ow8nTs)Oh}BhvW`S`q1$ET@ABu$<1=AUc~xI$K0KA0g~ZWl-O= z?~aIbjy`ioVf5JwLNxkZ0$e`&e4~PB7=3Prfie0Jf~Na&;6mQs)O|<}sRr@nqt711 z6eO`~mp8TzvEb2%*56?Dz?Z6YNDUxHHWqq&Y2r!N`RjCUrHD}4eFBlrF}B_YH8$)A zA<}LhaJhC@;$Ns|{;(a`QZbfihIJ5l^bN>yb}kJ3^Uy@yc_m zcRUPE83mchDf_(I{T^V_?rZe<(2fu^-J1wV+G&yWuhjV=q#jbm0gG#&fEd?4-W!WS zjI{Uc!gK;uLn%mXJu7&==RPl9ytD*#LkPYehZoQ5PdcI#CDJMO=6mTZ2t{X#7tc#4 ziRkP~=hjww^SyL#3`OTS$M^T|DNHYs&Ngqpm(ES0=se)X^U9+qI%biM1EE{BX(jUs z|JtbsHP=`gF95}hDRU;ariC>;P;jegjC8ALgmdq8?Q!kJN*;p`0OHzn1ZfLPub{N5 zi|`V20rA9x=UuyKb?%0nc8crQTdAgp=bVyOeR*2cAooA-JcI089Ju;lS`DutV@;FhcMo$|ATJ#1T9Nd<5O#i@?n>pXcS?0ZPDf_dp+<*sXd73c2)% zN4W(&)o@ymd%J79`!$yuX{GCP29I-Z?-P*ss?(!zMfC`e4YnTn*!}8}4{?9)eifjF zj{r36-s?_ax$Zayp@sMUZC6oaajQ0iS2f#~i2`YDUlq1y+M;^C>Quc%zTP)0?wVOK zeer!WqMX|mJ%iLmLV`}VnZgjldiDaM)~r2=WS>&!7Tkk}kG_BOq2Lay@!X9 z(SGA$7yUZ{gu^Z}iCZ`^iXpq4yP~i&bK}m#ZelWpF?k-HxzL}Gl<`*YpH^Jz4ks2T zm{6L&7WB&~d3=EOk3=rp4|_=9+a!U#_K}Oe4WhgElJH|-Mt4Mj^lnJ+g7hB94TTQ7 zMag^j^f`1&zCgG~^Yre?r;t{=m8*VG_{#!?Z}JMi`AV{8?;cKcyQcdS+E~eC>h~q+ zP+8BD?N4Y9yGX=uGhWqnzX=p~I5Ty-=kUAky?xklbAGJob!nd7p?Ug`z1t_d0a8!? ztk0}-EBoBapTX^>INh&_a=MNk`S!>O_ueC~9@)NAw4G`fEW&C@7i&O(wv0KR_Q{kEFPa{^#CJN$j`9K^T2svCh&LgAKoEq|B^p zyH8Qo#=@Tx`e(rH(D?{2p$5?hZDYaU^DNN*0YMe?sX}*Gg;o)#Lm91GwFhnl*X@mU znE=qa2jLarh$XSuUe_Q?O>ymY?Sp?r9@*!94ds=t9q1^xyLPyDx^^8Ic72NQ71vJJ zXRbHgN}UUi>wevJk^zc!>LahX;bw>Sx?ewX483cePI+X>u`5d&CZHv88faOFG9xU2u1Zb^JAcEGXx6O0*zP!b zKByTFo9ewbox*LpkJ;26P8N;p>H0mh;L$BC$lc6>^}j_h=-$q3IS4`t7@{#x#vUmd zlNk;Bqs*i*lPmHCeF~da(MaOjLz5(!Z?F3eD9UK<11CoNHAINl`%}m)j5%t7?Cwh5 zhsrwgnp?G>+WvOeG1o^hqWZ`#klcG@=b?QVl8ztQeq_5V>mXGt_VADF-B~}-2YE;r zB`2q3^xfQ0jJcAtMh}46-8-BNL;3HwSHj)mV5kIUs4+;zlc5|9eF`Y3?k{JKcfar~ zndbh5Pjk~CD@ydQK}D~E86)P%-`wFhK=keI9dMdK*PH0(UO#di&kv6L4ISPw8b&#Z zI7X1F-Jgw;XDT_jy(fNweAFimt2&4{#aow2ye{`k?ITWVAA?Vkm_6D}Ea-ZI1xGg` zpms4x-pyM_yss(~yzN5EVyp&W?|ob&*GQ0LrQhMP5aXGf7bKBo}WPTu>651_)N4u&t;M{GasU{!f?~Y^b7H zaZZdXdhhTxRH=MCJJ0&eeH`TC?qB~a1aNh#Pp7XYU6D7y6mU7X-VR>R3m(j^Q9!r3c}Df zNswHG7FPMix=ggqv0o6UD&1mDcP-GIcSK>l9*r0VMM=H!OuZ_d#g(sH01?bk)m50l z)yVW*R61ay<_Q zmE0nh?8~$8!aakW{T+!YJj;1n;l82RaT(&p&eLJ;8_FFQx^JjDF%Cnh#_JD!Jaz7B z!xWR+JI;@|O4TS+PP7nK5pAfHubYRM*N7ODRLmP-D0d>ojsU_;#|$|963y0Ve$dyg zTMZ}bznx~(F)ma9PC-TCK$6lqsP0AY)?XNdK;56OTilC9i!pe~5X-=hxhNDT7K-78 zl)#SfFQjRDRDWR{E9)Ep3o6}&(kP-m1U?88s3^!ZAs&B}c(2gKAs%y_=SO)E_Wg_w zto5xb#jM+^_x2RST6oFn-QL@ikEl42?%mznQ$;z0n(ZFt z_M=KqGxiiG?}8w@-l5)}a^&^8IP%W^ypi7R;$*qss_E^|B|^Mdny*;TEPx>fDf0At zh$d>WC*D^OHPo{JIZ#&^>RITG)p=u!#aPdKL%<=dCK00P*-ERuUkwRV8r#AM?)v@# zJT)C#7*v!hi2>iGcQVb&5ePspJY*f57KX8Sb(~JrtUWbR60UCx zsps@04vosM(l!qjDmyzI+Mz-l{p^dYOvqa;s>s%KZs{@mO^?L>c}OMD#rI?tfQhrd|7JuU4YgpK$nCFPz>V z1I4)e#RQU!tvK!F=ofj>xpC;7XE`(FeN29CO6N9Kke0C|mF|A)5~9V{tvGAV-48L) zI>O@Qmt54>Pb~mUdiL@<(DWo+fK|_}K8k$i0^mh>k~jS!fwOmp>?( zVbP%_FV;55gM}WMb`xuCp;=Q062BK9*FDQk9XkpX~`u1Xyi`F>3pp5|a zWzk8uo)OL3-(#7%kM$x$Jp-tklXs)>_1uTp3HT`dS=rY^T@N|(Dm0Sbo!h^?{Oh+T z?CkGQ_I%ZQs#CQUjhS`Ub8-HocZd_iprG!5a5V~|F~ogMMA1#tw4Qo>oLL3sC{GU& z#0t%en~1}fjfCFw44&UIDK;a9t_Xb-5fUKjP|q5~ssZy%pt6}d&!<;28o|<^WaYL2 zJ8U&+G>8}oNVMS6ezabkfF$~(pAZ66R-*{y^85&J3h(D(Fpt0dJv`=%JPJh~4~RT2 z7I_#DI$5ovl_rMe)E|Ii=;}&W`@&s-BVF)4m$%y`O5%^YUa6W6|E-V(nmc zhw9}6G7ahlyDnr)s*_$!-93K+NvfAEtX^I`-+Gw?GM-*4Yu5E|nMb@U!5fre4x@NEm&T*vVKA-Fy|g~@(n3R!)ik5^Ho5|;>1aZrn#Or+s^Om! zyoW?qAm_xT;{3jkf%?l8&ePSs_fxhbk%}`%h7qagEEN5Ud8oqiCIK1%2%+;AUO6uW zQ|^&W4DQ{5I$tBo$y^q4>O?tDT*hd7!ja#XN@IP`SDMv2*baSaUwqGdn$;JvgxZ~p zw_cs+9<00J1v5?7uRE@Ap1wkJ_ZajBsCRWAn|_cl~CPm?2BqzeCaw3hX8{%pjuMsL~h{mG79)aga zSao)fU%+#;Rls(u%t*OZv}@!}dxYH z(LnF*x?@E1M6u!jD0S-)SxNYilnA(4FVF13DU#2OX9!*NrhqwI&QRAO79*e8ir5Cy zhlcWpQ9r0+sy060Yn=3{i(5szTkgcwpbvIxs!?LA1{7>h0qpto&+wJRWkDEPKa4jS z#)=>eiyy{Q45KRuqudXpk70Ck4DCX(ZbiS2Cn+Hpt5H6|Kq^L*I*M{#lSFS;i*hKL z#Q3-XB;ASZsoKa8r!dGGJUyMmew9wjYr0kAX?WO&T)gA{3LUmH%Y$NJn(lQ-&@{M- z%8zejR2+D2n1=_9y!J&{9Lo*;{q#og(ngT$Wdr_Jv{j+^!AwA3ra5-+=lm_Bcei_{ z%eh~o$!87j%Q38wO+<@Dk>fdsXz*bj z0NkG5-&iN$c^cbGZrhoREwdnzre`V5(3y0s7vlOV+ygF14KD>D_rPN;zK-JSHnR8w zEWU~2>z=_=(|xIi76fQ3rhZ_6uz~I66(epGI8+_*kF_aAU)kSQifeoX;#cNZ>ZW(-m z9bD+?LJZu5-kvtZ`q^`+ryFsqjiH_uV*XH%12MYbT}sb(5GQTkA}0JD8({u+(~Dv5 zS;`;UK63fjf0{7Vzf=)4sa>f>_D=8U7}vwQfEXsp3jPh4G@X48q?|9(8~+HtZvsXi zTL~a8(?O=MzZhA)Lm8j;?zyv?dVF-JN}9ut_vZF>#v2Y`p+#KZ(>$lVXtnY}ET3k) z+W8r&5k~T$$Tye(^xMB+TZumW-cu(Y_ZyX9(FCclVx!qVKHhzB{0GK6iO+v?yju=F z?(uw}tAj(jRs6WaYS^7t{AnZX@VgIboMG7IKJ0d@)= z&o`(wt5n$M_Y&qL-K?{9p1xeOwgYdt9LSeLo?l{zYkMvq8a1cVmWMSFM^Q^9)l2(g`tT9;jX&4|(tt5+M0cx#GtVAYO1;lbxE!9G-W-LTOEG<4uhgkGwI>c(R zH$dcgN8=LvZmYNewMUXfT!8uwEjx&YCbyTC0#@11F4AS1wOXG@Y)$M1dkQ(hWy}dG z-~fFHTTAj38uG6q~|1h`(djP zYfPF|Ps{0%F#D*mS3J@%6J!Eco>s1smvx?oUw#4Icth^nLe0K^p3eD+Zb&H4vsDb` zD(&@^Z9}>0I&#fi?wYyWHFLRZ=5p7}<*u2_T~q8^Sc}8p@hi+g+^RM|7g#waP#kIz zjAsdWvgVM>o5MfM)4+JLzR|a=rdhe|WQZ!hQgg2>GPp|2t9TWk!t5b%)otZHp@>Y- zti2x{g@5PGQJ-59>88^;%Qw0rlr&dr)_PdYQk{yOX)9|-S9(QiZKk%Ot|>0ZT5=R?=7qd{7sfR} zG_)fJL>nqha6q)7Lai2}p=+fD8+ET=tGbdT8^kQct#bKGwwhj60~b=ArzdOHJ{u&N z(=%BnR4IC!!l>HCWVIM>wK#6Icy6@>ZZ$2pnvPp7$!oRA{#Mg5*$(;3HbO58!Lm)% ztTp>&WBv5Y+-eht6s#hx%C^Og{=FjW4r7?PtZ~d-=oQIa)&y>@+qBGF=&2|Ai}d3Q z$y~?%MY=PTNRu>c2QfVj96Ba>rCDxw(9nS@Y(^i=>|~8$cCyAXJ6YqIovaDWPF9ks zb+=U~4j+<4&cR7khi5QJajQK3LVbsKLYbYCG;6Q*2_^1aBzfzQ300^@71p8-QH2<+ zW7x1pqcn!K!gw~U6(+D@tx(H`wL+bDSo8C(B$BI+$#v4d16dYIE_nT0XvqG?dS3E+ z#e7=<22^Ex)nWu{d&go_BlmjmsK%?2d%cx=y;b!3>$unZjd9HDUyUP+`7ejB{L~Y| zV(|JoKFLIj!Rzx{959Enyrz8;n&NyLCIT4J?Ee(6Phqn8ugkAo97;A!ao$n~cYhCH zDV$9IDcf}LU&&vYCMQCEw2Mjv3`1K`c32H;THm49KUXTvi-%uloxj>$-sWceZEh94 zl8KYrkFoZMH}C!}v~1vBw~c0P4!+$g&h^>bJUFUJHR7Cv(wOAZ3#sSW5tK^Yi0Ii$ zkuJfdo}^OaH=z^+r6--D+Ux!Z6Is_+yTaVZ@I5jXh+Jcdl}p98{`dyqF}I5OrMQl? zax5meu2kj_$FP1#v=V7gariVQhSGeW54C7kwB1RsjGFG*A7NsC3VXuryd|4~Q$9vT zo*z&7SyFyr$Ay)J?mr-(2J}_)Y3-5C)%)8H@sS;4`tUluZg}_W<=bK1dx`!_;`>8T zgD3jd|04QHnpKfBposKc68)s%3$gs`qc4~56QX|)wz-*p^86Pl{|!39jrW%j@~x)) zI4S=@Hmt!u?@9Klru+riw`b+qR3^Z_6Dj|3y!f1j|6lO}IQ+%g_1{MM7Aar0o!C9) zL*!qB6YG6)*{-Lu2fvne3uM=1DtC~|u|)!|dH+0Z|E))d@Vghb*_aTabI0sU&3e$> zDAt2w5=JiuxpAxq&5dV0XfDP?ych({zOKT2`s{Pc8^GVMy7@#*yE(AM!q81LGK|2i!416E*>RxjrGN){@0Ye35rdl3BtU*Sxo z6}sVv*s7&(n)U~K9lh>cF*K^Yk;YHuVhoQM8ycGM;#C+`Q|qA7w8qhjaRaYAG);=v zw%c#TG_BL#;i&H|)$T5wD9(eZpQBm0bMOl4%o8v{9ex-sSsWpKBR6{!EmStPl$C%N z>Ycw@(K}TPMr?Yug1^zojTcdjQ3);DId~J-+VUJlLxsAQ$~GBCHgZ3K>C~{BUiw*= zAE1ms#~EOf72}FZ<;?a6NS8BTT-*E3KP;KJ|T4uUM=ZM%1$NLHhfd! z&-!JDGzrGE9I8U&@(F!X1EvX`oOF@MA926dm#qF!0bPWODVmSU2-@$J6|k) zGqU+w!tm2F8Z{zPfe&dSb2;kGq}SRc=O;-+!hBn6GriU})D)vfEibI;z-uk5{9Ham z#cM6iP;>bVHJ8s&bNLK4m(Ngx`qv@6%?F-!SQ!^+Gv28Nv94@t1?#1(i~UDg7lZZi z11rjUtT$KpNg4V{uT- z{vBnHM|%h?eOk=EkFxt@*;;Y#*-6>TFM>NvOANo&CW*LVA$|;uozT7%nEcLN38Vs>UsdU- zPgeF}(4>_)p%-)?5JGoh?K1n?E#`M)9 z+Z!0bydEEx?R4alup*h>a3plSc}mYQWN#3K>SfYND0ynneaJKnLN-?WTj`DB#1~Xi z;(ibDV>0yNlM)E70UmzI&-;6XlF%#Y-D13?`T@s9ePJ0pJI9s6H}M}}%q6l-IDy66tF9NGT+GPF;44RqA2GD# z1kW*gJx9d(O|IkD&To{C4NMnKfZCDXpHMg_BuOBh;~-UFl})gUXAoobhSSXd9g_cV zow5d|R1uVnm`q9M>4-ayvJVTi&eM@t-b6-UD&^oPQqR4VO4L*8Jx$MDlp1SB>IDmd z_G-Fc0_X+8?#olNowYnrcJ7S2>stb9oV`?JT;J`={@)4$u8G`bV2^)zb2k4XNg`Tfs{V13Vl+eY2> zD*Yz$z%^eMA!H`3p4oYL^0>;xs+aBH>=4o1*!BzHUm_@z3j_QHr|>_n`X> z!YMdv{4pKmzMEX+(6C#Ti<-ucFo4+QT(0gM1A;*|g%3Gd-UHiVa<{lq(WW$wgjbfPGoD_^=52hjmwd&x_mJhAP>D z6d;J^p4(BCnx0NPU6V+l;hPa<^}%fV(SOit)+C`o5CTUxng?}@_0=VYE5vY_7&^qT zO$?jFaIqNHiQxhk?Y732`c^>e?4jU{z|F~qG+y@F zCAJE&bY4?UM~A&b5Gw37b=VdHDbValsn5_0d9BzjYHez=xBF!@wzM|YVXLnWhAC=s z*c-4K)Vc&4wN0%{wM!Scp2*j1JWssT z^Te~9A6~loJn_=c56^O*@>zaBzQ@{rlmj#7VuE-l4n@eo&HP)Id(^FqrC;EzZz@dcXqhcn7QqZ zaN`#3qT95UDdj2J5;z)K(-iHsH7&_=8_|@S;lC}~xe(3XUI$F+8pDOGZovts_Krr! zZB%Pj_S(g0ZmkWsP0?0?VW}N#+nZ`y>N=!y2*=E6L;-DUXS;T0qhnU*qU5%k20OJf zCVJzdPHgZZ7SlC=3GIZ99 z4tk|fehF<6geKxGokDnYiiZ|f&n60RD(Hk6!wBX?nRLu6bS}w(vZA75*#Kp^09CZl zfaq-6Ks$#VbT;R#`E(G9PC?NTsI$?at$`oXmvCsa!Do{uJ8dvf-SG0rgS>Ri5H%FCG@+MK289WQr>D%1YVzn^L!^p`xI? z*qT?6WStC~q1hE<*ld~u*{T|w?b^|(!@Sn!Dy??3Hv=|<_dwAgAiOdF36<^#m~Sib zPqw$VI!emQ%99M}9!E3(FJBr!ai+bcs{K|bQ&NY$cCvPKQhQAu>r7c}Y6pRm#RrH; z*`t$k2y1lGTppA*BP3OgNmJ>XoMfG1NlltK+Pb8@a;~*JY2uV3rBm)nQz4*gLf7r)ugjheC6l)CRIBQ{8e)OY71W zZCi)Ev#u5Lvgr{`R`7?K9IdsjOY29so1xn{fGQ3$pS#j6M2juDTZ!-KabT8kR*|io;f2v%;o%RB8w3071pe-5<*Z1iC!q)P!a6hK zwddy}BlU!Us|h!Cw7PX8BMWK3nMpc=_T5wutCZI0IH!uL?hj+4>ZJfct!YU?Y} zN^4};*3jTbLx4-DZb*jpxZ{2~GBSn+J6r^#dWrduk&#;nfQw*za5n+QdKe~eT z(g8*`To2sYL=*NQe@h4#6l1s*XA^HN?tI*P9vB&6yNM!ReAs`9Xa$i(XEE-B4~~pH zh>+6?D;op1Z37W4K&x?AJ%kMdKVETR9+5^Kqp=A%3m*=kF_UO41C6ukVMZC^zX#<~ z=>-(>1-6MUMfm4XDWE~*N1lxy`3=$A`Ut+#!8l2xPAYqy*iUvu zC(MF6(Fq07+B~0qNFVaOdBFbA5}~sruV-l?2}Snc-iq?r zDPuOPfz7OujSS&N#Agr*NP9MYNmsh#_Kl2OOEd!PDhPD|k(~L+KPBdi{Yjp;KbaZ5 za9(s(Ui8Ao5Z%P=$hzoqaWs zkl45%{W3!7PC|MLZawaWY48Vxc(S)hTLV_IOh&fo%ORBa@9Te?1?rCaB;JJGUEMK3 zSmy^j3D`Wq_W8kt{}O~XfVuo&3jkXOSg7kb;Ae2heSkmZi`$615cfRX4%};SKY;rQ z+}m)|wGZLDxR2vLiTexO-{V$)h`t^7CAf9CQ*c{w)0K~Kp%`C=@IKs|aPPtG!7Utz z47ihU=i^?0dl~Nga6g0lAnr4`6QFZ0?ghB%`rj`9n*x?yvG4-zI#om}z$uV+9K9KCReIi48-#q_wWI$vz9q#ZC05 znXCjWZ*Q%&^2G;1I1*OW0&X=;SStu-kP=>OuUR4k3x@FG)|xtzFR0$HtZM22;W;%; zop#|xrMUV5J`QLg2ro+zwjP4{9BoozFQ&!53>z+q3Y)#Xy#+hx8x*$2)ieT8W*U+P2#S;dd14 zu!FZ!wWPkSy|KlCN%*9*ZjNC{JIi)jzhG;?%8jGFrp1P}bZl7JIvk+a0DmhB^X0X6 zkQQWjLHL;;TwERUX8!zadA1rzi)m1}k7kF%R@>Ct0UZV5tw4Bdo4rLi5C~&Ca&OCF zwmMjvXt)CNMEl{ONC>!8=y*ULDQ;`Bx7c8UW*FBeD;cPx9$S$$?d>(F0mWJKY+k;q zP&l=Vu>)C4o0?v{F`@8KrL^%#dPSZiU!oEn3I`TzLiyOUsf%-oCUi7ef@uF z0fn%=H39d*mj_#wIBfp^r4jyrRXPT9U4vQx!(Z8(Ur+-U7+fZq^-QYAl zY1mr8c~F4KVNQPb~C@0rx*E6nBQTg>aszcU{-e``)k%}Z@d?MvO7 z`gW?2c4^wIwA!?Rv?tSCX`iQQ)6>##PVYCWZPGOJv}x1wr`?SpCI(=V7FKYil#Dbv%Y&zQb+`cJ06JpJA25t$P+%QA1y{CVbYGvCkrCiAMS z8Ck7a_hkJ(>p<2>)+O09vlnIGoBc%gpRzy8zC34IPD9ShoTqbM&G|Y<&~{fd zcfI~@{UiD}^&jgmG9()o8tyU-7=C70Z+OD6*)VAM(hzCPHrk9WMu+hZ<7(qij1L*# zhRrTDrI~6?cbJ|uZ8v>sy3jnye2qEPoMSFBSDWuKuQ$JK{=yud8kssd)s&i*Iy3c_ zRA=h%QeR8`Hg#NDUfK<5eQ6J+{b$LFVntFQ=~_vKb8K6^xf%i zq`#B?LHft(pQeAGKFTuIa*1W4WwIsBl5HulR9LRJ)L2?9ODsRK+yi_5!t!g&Z!B9Z zt1`A_yqWQC#<7ggG8(2$oNmr6$Xt;50BkcoYgX1x+4k(F?2ha;*>O3EIalT+=alBm z%c;$2$@xvrmYi)lJ9ECwp#^C@?4r@f>nG^1(P!wb`da;R{ayMs`k(6mq(7oRsh?m- zGGrSH40jlM3~LNOHS9DD8!j*=8g<4LW2&*#IL}yXyxI7;@i)eU#uLWRjM1j?rfgHb zsl{}^X~dKWU+OeJZO%_EPMwo_Q)*pmSL&8jPwJ`EGpXvdn6z=}x2Jcd-wo~GOaCbS z1oZzp{kwFPCDNj`JZafudBfte9I+g?d~7*oIc@pga$$xsqafp!j5{)Zobg!3Ug-Z# zhA?f)w47-b@Shda?w$6?v@O&2OnYzIC({(u5or$O^?acWzNlXWUkJ9 zA@i+FPi9V5Th?!~c4U2%6`gI(ekA+3?BVQ@?5doHa~_A+ya4;Wn)6Q1;hYb1hI9T3 zAEH0HD&z`68|)RXAB~ncUVkNQX3|gB=jmtbtMxbP>-5d~E`5)_57zq+{iFJ)^jq{V z>R-{np?^p3(f?KdrTzkg1=gHps4^@v+-z88xD(#@3&U2!ZnWO_4SzQL)$o-;VH_>C z;cJX0<1}M|vBX$yTwrW4wi=fj?=;?R>^I(Te9-u)@k!(FjL#bfjk}Hej0cQ|jmM3r zjNcif(6-}EiP+>WGTBT`CWmR2X|3r$P0yOPnf93eV*183%6x&@U@kOYZ|*n$%KW1F z9rItzs?@~P)YOvH+SH|J>o29Yr(KwSY5EoE*QA@$r>7UBm!wyv-qQ(C- z{lWBK!F#u$&A*cFN*=iHS)XNfWnYr>9ND{^Ee8li@WdsC z7~>}6hsIBh<4tDMG}CPKT1!oLn(j8;YxHe(&^gnGa{alzBKaENenmZk8?Uj;sgaKd++qc{l5WtdAl0H(6oX zW3n&Ez9f5M_SM}lCEvaijqK`Xp7`)Ap|&i-xozU*V@o#y7$<#gxVlk;fK=^SC9 zAUp)<=<^c(aW;V+--zt>-3xXw^(=tlqaF#7ms4Z93)7(O&yXfzmSLi-w1 zgQ?AQk7=Fh*XWnFnqD^@F&#I3Vu~|gX`W)vG|xc&-fI4-`GVAmsaL0dpBjsP+mUv6 z+9PQP(%wz`AnoI{&(gk03rio9enI*r=@ZkhPB*4cOP`T`ZF&W4wFv!sd;0S9mFbVA z3vGfB1sNBb@aW?4fege?@wig2Mo*$(`wA3~(#IhE@9Ups0a`B)W5kr~skaZdUhPMm z*1^Ms$MNL(ep&yH{vWWwqtHRDe~dmDvr+%o<4#GL%Kvw?*G|O=aeg>wDqk*`iWRH+ z#)helCQI7X+S_O=Af-~6x){gRrru)5dC#d-s+o>4v9gOg8=LBURFYYa$(UT*-rCWT z+`^7za&(hXWFg>)iq=}G#Z1b?d!9gzqhJ{6REXJ~-HV~>-Vb!xi%N`6ZO3C7N zt(`c|E`%vY(ejrNo{c&Nc`u%l7DhYk1+C)15oD=sZH{(9aggFj-~hd#c!$JLs%(pD zuwbb8BW0;ow)&=;20`&B%8Nn@2hjwDi&8Nt23AB=;U+*FigwoH$e5sbm(uYno1?j{ zz6DQbYpz))gstUh9k+GZ;6Lp`*sBCk z6XtCwL)&8oc415_eWp++gl96jkivfDwzP_i?yZ=G#G*~pzJQ>JZMV1C+S={) zjmvBX8y&0^!u(-o4h#2(SvX7?+tdn1Vks3%38P{=TiC}gbvDWsB4U@YPKKw{vG`8L z!DJhaiice4aFrmOMg-}I6Hu7H~Vh0frR6H@Xhl&Ym;zd#Bbm2lcmAe4lY%PwcwzgyIqr=gG7wCw6 zzrt?13TIb^h}XGLP4*fg;th%`(dZlN9X5QPBSgGO5VZ)wID}Riv7ca3s^-on2R>H9 z-cQ6^&x2eHhv3Mn5b-vl=u~x$OB&JiA`Vi1lB&6NNe9gEM^ZJDQp=HIy;`oS3C)YX zed>TMTvR%r5bQSCv9VQ%a1&@2E3>ExHHoj5BHkldG4eO= z41j)4&}}MeJa!zIwk^Uo+bu%G7X%)JYfu%{wAfoaJ6HoF9DE@b@g?Q&Qenff1%=2d zgov*RwntUZeBv}E_NnSoUWoV`CHAZ8QJF%-HuUxGfTi1xbg?IB2wqc!5zqgu-gblvh%2vzF%DW>i>9fn+jbWbq0e z3M0)dsLzd2Sw=!u3ZsM2ti{!Stb{SX?1}h(Jw}oBUiY<2hK06o*{w9`M?r(S-BABgG$JN zj`u;!t@*afS!ESfw!-StysDzIQsKhD+=_y#>Wb1z;UXWnu&lbmc5O8(ySS*Ns7i?U zWudyN=2|NtYehj}0X&5)Gozrips*+pJ}q4ALoLoL7v1O*9}tz74@>4370wqf^?}J) zMWuQ91$jb(4~iNtm|0L^%P&NU0#*T)LZUCfth~x*omnkh<^x&hlojQZ-EDc^^Ql0$w)OJXwQ;u5hnO@I9Z_E%F5yYL_3f9!$cpBSS^)R6|f0v z;R;_Cnw7Pgm~%WRStRSG&EKAEG^I;&8a6bQ>JuF4m#41~?760Qn_ z%qn8$DlCE!Nr73Fep!i1 zliaZaLZic_TCAkL1x64|=Rhf}C@Wz|=8#A*?3{v%N;Jz#A=QV$d-VW!Ls!O;Q3+{2 zO!5?75jKB}^pF@{0v2C3YZB#^v#JYc2pK-GNK@)uZ3XiR@~T<)G0lg>TZMVXeDtE! z9~8}lHK+UM&oic`ovoTP&w&_N(^&zS2BXDWv)R5JvHzP}R5gp*y^!O}vKG&^&abo; zV$i}6X9WVBTFCVwRF|TwFQ#T$Tu>@leHoRriVCX?!VDke`T}aS^Qtgv+c>Y%d?C-5 zk2YOW278E&kUAS7-Vx0x&v_~sm6UgA}vgUuuP zG?UgBgQl6^D5CMR?u@jx0y8lhi~`f;WO#8-DWWitWd@cOXwJrf^Qfipi6cwTppR6= zc_ITanN?io%fM2cWAkY-6U-;GGK_?ZTn6Ia*D<`>2@&JMNaYLZB%Z^D-yFc9JGa^MA?$Wg7^SShb3A_g1E&V@iyq?Cdb;XH{5o80X&HA-Ar{jR(lDQb zu-H5b!eeJppp307qEHo^&%#l$7+wIV#C$ccxNOFZ0(b)lsAI3=k1##Ox09U0=?Hu_YyzsM1 z1|yn^FN-Ex#dvMFiP?ZttAFM+K&%7 zBOjs~@JbW}GU9b82xP)*P!PzB*Bw7l6{B@YFvy4(9Us>9G~qMgEyo9geEn#Dvjoxu2D5Y{3tp&?)C&q=^hk zO`Yi}c(2s@jv^k3K9_a5*sDOxn>3mbu$Lq52eHNNi0bQ#SQ8lCNxgdb#g_s_%{NBm zyCY#rC4Ni^C&I(j;gF}T$-Yc00hH2aQH(cBg2a0v0lP{uN*gwUKn`<2Axs-?vq=lW zHbI>PLAK*&3gJ&rKZ|&z4Qusnwe*`!SV&qVj6R_t6o6}+mUGCb3=&b_*;30kqsCl4 z0*p~g6*d4bQvVh4Sj1>w#zyZucqK;vZ4<^!9R)Ce`TjsX0pqDaWb|lN_))GPXLQCW^Y!` z#hLg6SCJw>RWyzRq9%aVRFE4xUUAjMylAwFvDaJ)3NcD1L(w=%ir6HOOhADPm~<)$ zK2DSllu5ccza*3$??*;O;$3(V$dE)AjTh+ekj9&ZnZ!dChCbl`-i5C%u>W6#r$X^< z^jJI<*lVIs02Sxn0rxusAPC+w08~~`h0!U3@JAGii4*6~);DC$pRb>qJ#k`~0^}=! z5+;-ogJ6dK2MDAi0=}aM=s}x0v3GNBv5*^tNaqfIxDb&1Pc@vofFN`@>Ka>97L%P2 z#m*s*H`r;zmd9!Ll*edC6x|&X*m$K4rDV2mRny|=01R7}bdCZs^xqUq1_{bQQaWaU z4OKyCSk#nEXA>w7Nhz%6(3QeYj_|LxzaN2%E6|fr)5p;u$0GcWWJ^k|t##4ORE$b9 zI*bNByJHtx6ty5+x8tru)CcJWWW6La!cf=401mM85Tv`L7cWIJUxsp+G5B36b@=r< zG0k`>&`^UeNq4WziqM0{?@?|KZlVybsKb*e`mU6AA2ZT2Tv6-G^owmrybil$9LN7V zc8-*Cb1nO#&1WCU@-&ok%j5S7ZB;mg!hfKX%1PE3WcDHMoSOervIc3ut>K@`;9%rygH?bm1T7!r7C=r(bXc)IK+4g8xSa{~G$4bxNeI#(-|-@) zMd*XP2FODEabsVplYs2^LL|@RGLvU+gk%~b=Dno)(gMOLE2eOEUJTTQ*hhJ_2!$iI z!Dm2(*H`Gayd880&!fK9K@(nIE|ehK`7WuxsN>}IRoK{s6HmT2$#J!iMUOvoeXjbN zPD%sES38w}EJR)SAh!TAh?hzqWF;V4)SwUY3qZOq_pgiR34#{`AIPv zEf^(qFDm+>m$Kp{>LBov`i(sx0#cwub2s@#et7&m{ICslxF6C8MEdwH$q#z~;(mzV z%Va+s>BF%#w$<8@)LtQ441OKzx)e+vga2@g5k@lE&z!S|G^(sPA z520uP#GDWeNm~HdUD&uup#w-6kc*V$g>1Dh&%Xobu zBzUdHRFjMt{~;FT-m6))rUfQ;86 z887<&;)iTFOU5fs#*5xu0yH(#8%zKsjppKK;{{}4Eqte(@zRETH<#ul9Eav4n4*Wc zehE4Y&q0TZJO!Kq2?xK!gl)9O7U8HDBJtvMh}Q!sg{iaOC7=I*Xz+YPoW1XoW>z6q zgoxxI>5f-%-&^h5fE@KgBzXh|A$dMTDYB;TlHOSNczMa+0#foSs&w=6>F6_CH46_R zkBDQX(q!SY;}lzo`&n@uO&z-HufLvVd^m{>I&hMpD}`{ZRz~M+I8@{)*#{G8pMAau zWPbqUa?DxUu9I*ii=FaX>>|uBuL$5TS@@w1cZMd-R|g-H|#~VbP zu?&YiTYfc)IMK)oSKNvx82eq4?s`dg#npbfJSw$!?@DPP35F!qyHYE@)rPb0S7-$x=P>9xSCmz3iR|CxGGeNip(`kIHP zAG*G%mUw+#8qz{!A1!dWec}<5qfwnwP4YUWPebL`c7P6f-x;xd6m-$IBK3fA`-C2M z@KsX_;dC<_o&~3=oz8aDKXVQosV+DU)#n`;h-6$W}a%t@%0Fv4z!%30htdQa8WH`%YIG4z995S338BUuFM=ir? zlHq)h1{A2*Vj0dCGMqXY&dFe$wBx~$)CICq@5*rI$#C`s@&jd-l}e>vGyo@U zTQHr}Vj0d8GMrg5oCjn$`7)d}GMro)jzfl%DZ`m3!?DP43S~HE8IDngqnF`Kl;KQ~ z;l#;sl4Lk%f+4B8U`X2WV2HU+R;PZY#F6n8l;R_!-@d24r4GL@ zaXk%)cWh#iZGa5S<-J769XfG^YBj^5UV;@1c>FU*2U-{P%B;vlDv+UUFEYgJl!SZB zs~5>kJ!l{u!g(J!k9qYXE6H)T0~(K^itKle9@3bE{DAhAdLHv1Bk31dY8^xffXLsH z%m-zq?h7uJE(E`SnYF<bjV^1l)v}gN?^c}GVc^h3Ys(z+ zB|$Jf@Ij&pUC{bLk^mv!_u=F-9RAkE>CBNq7Rw;Z0MUBij)_YjAbf>NE6Vc_Abf3y zAY6;(Y$izLIIrI#x;VfZJE&9ufT|SF-0wT3{E)Tu2tq z5fqu{ANTG{Q$cK$ zseni$EVq0qAoBunY5-Xg0BHq;jO?Sl5|G>gob`Y#3xNCvkR%lJ`RX1(_)4jy`!PU- z08PIEWDkJUBXbF(ql5eBD*(|2@G<}*t)a0|20zHa+aSn@bfyDGnhA5t)qrRNY`6#z zX-)z<;!oIn0=)Nj;P}R6qTdS$T|WD417uHtgdYI1EP&THfE*0KnTTbKt^kOh7BXngK>0ix zYY58AD3W{v%66pt7ofze;9F20Ms0ouGGQWxUjJ}Ibaqocg70gGk zQr&k!S;MDtH}8YORa;tl1`4OBjF6v#@+unr8&H0aC@(;Wp|JI{bhFHo)EU&{Vy2qEweYBKE}@< z*ei~U9#BLVJ4o}YV!Dscc&K}{>oaTc3Gc@H6g(@CrKXVZp~&*8&09A&H)Kc*31G2X z2zZv95K`}AnO`EesksH-lqXQU@5?XF`FWo`m?K0SlVU!f_UBi^bU3mdfGpoqa12)0 zYIi#KHlc0s)RDVQtfbV~Vtv^3d@%`!8cmjY!CkH0>27Sb3?M|7GyK1Kr>Q^daVa!> zgDyFr}r}rQaoe1G1h~ zns0Ac$kJN1sUJ9>0zD+J4uj#lUGPiA0>kbG7*>bxc43u-%+_7lFxGp#Xi}PU51kV( zq*inQ#_+U}52fE)ua|zSIslt*(*E(BQu@t%sr1dGwmTbls?14+DDceQ70MC93=|HN zVky_XOA;{rwZZ?mNMeI+FuFFq!2s0*`4P*00W{L*TM$U)0`Y0bJg`xy%W6_853xX1IMvg#~AK^Ma(43+5ig2 zTSzwMWOO~Egb$A1f4py4dBX4`W2XH8C1D;w%}7@x1x4kuzK!fp#@<$26L~f_P}Cir zC$j}OIIt8X{0N$5nxS>_@Gxtxu2A{Qlb)mQ&Vsp(?d6wCbdi(L9E!V z7FA35Vd56Kcpoy_IggH?o=Np9HuHA{k%&uqkc$^dX|wsrW(i}B&L(05DY|F;0~Dpm zmZ(KY$CFy4274k>=ZkEHs&UI4yP`TJ(r9?kRH|=;GO3{;T1cN@o{hjk1@3}6m2^9d z8n#j}DU%wqVnSW936SyrX|Vy5A7PhzyHC(N20v?Z1q@VeQA|TLk+wK5>9OM8gxO!>uo1i(XouT z_BiV;@av`X0@cgpW_WoQHAaDT7JI9`o3wL>Sp!s>gJj#m4qT-m^O_#QYLj-gH^w4Y z57JpClwBaTobdiZ8pix!&eh-$2nDct=(E-CexbMk#d&bxX`N^qFkxbTx0+#;$*hMu zo-AK17RF}JL!(wvoMG1^^)SGl6Xy~`r-3iAK5Uk*k62?Tt~o%_5vOyrQ*beJzjOh^ ztB8RVvn*IBQ8X1Tqa7sZBEj5^9_Gxezq%*t2m>#t?%c!d3e4TJ(n-8*USn<`Y#NW? z+ePvQ(Z5{!)vz8#QI?pkb3Ge;Qn$s~Jiyd*9RqkkOS0F5a;nT?r%;x4P=%$7*!vZO z@;T~w{H_$NPl)WNAQj0!aK>^As4g&Rluwol{i?JHEe314S*j13B|c&xbF`MwKiuGb z5KX`9$^N=BjCuvVlL&~PbLKkm$#Q1PdGoz(1ake=KK6b#iJ9RBsoqO`MHk}C0b{~_ zRy0RHzjGSDV?4SrNA?W2i_>wAq{tQfg}rXI9%ZjK1H(5eP#p&HNDOeRG6^}JQzaxe zPZur$dB%phWZ%ut=Wg3wdAQNpmEY_jQ?2nKth%L!*f5B*N-od8C=;&-n?w<~JP#7y z;|Jr@yf{SxAnrOWmD9J?131@0(ElTJWA_@!|SuESR#c7?u;j9h_5U5Jdz#{ zrCZ}^Z*_zy4U?X2Dw|I8*5J+y*!ap>@q53a=o+taS`5Ap$sGHulU!~JaP8ikvXZts q8tQ+^4Tl^93LpIwb9ujuz9y_{)H3BXr?J@@9})nul +rmdir /s /q Obj_win32 + cd %DISCONSourceDir%\Source mingw32-make.exe +del /f /s /q Obj_win32 1>nul +rmdir /s /q Obj_win32 + pause diff --git a/Source/Constants.f90 b/Source/Constants.f90 new file mode 100644 index 00000000..b1cb78cf --- /dev/null +++ b/Source/Constants.f90 @@ -0,0 +1,5 @@ +MODULE Constants + REAL(4), PARAMETER :: RPS2RPM = 9.5492966 ! Factor to convert radians per second to revolutions per minute. + REAL(4), PARAMETER :: R2D = 57.295780 ! Factor to convert radians to degrees. + REAL(4), PARAMETER :: PI = 3.14159265359 ! Mathematical constant pi +END MODULE Constants \ No newline at end of file diff --git a/Source/Controllers.f90 b/Source/Controllers.f90 new file mode 100644 index 00000000..2f61e493 --- /dev/null +++ b/Source/Controllers.f90 @@ -0,0 +1,259 @@ +MODULE Controllers + + USE, INTRINSIC :: ISO_C_Binding + USE FunctionToolbox + USE Filters + + IMPLICIT NONE + +CONTAINS + SUBROUTINE PitchControl(avrSWAP, CntrPar, LocalVar, objInst) + + USE DRC_Types, ONLY : ControlParameters, LocalVariables, ObjectInstances + + ! Local Variables: + REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + INTEGER(4) :: K ! Loops through blades. + + TYPE(ControlParameters), INTENT(INOUT) :: CntrPar + TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + TYPE(ObjectInstances), INTENT(INOUT) :: objInst + + !.............................................................................................................................. + ! Pitch control + !.............................................................................................................................. + ! Set the pitch override to yes + avrSWAP(55) = 0.0 ! Pitch override: 0=yes + + IF (CntrPar%VS_ControlMode == 0 .AND. LocalVar%GenTrq >= CntrPar%PC_RtTq99) THEN + LocalVar%PC_MaxPitVar = CntrPar%PC_MaxPit + ELSEIF (CntrPar%VS_ControlMode == 1 .AND. LocalVar%GenTrqAr >= CntrPar%VS_GenTrqArSatMax*0.99) THEN + LocalVar%PC_MaxPitVar = CntrPar%PC_MaxPit + ELSE + LocalVar%PC_MaxPitVar = CntrPar%PC_SetPnt + END IF + + ! Compute the gain scheduling correction factor based on the previously + ! commanded pitch angle for blade 1: + LocalVar%PC_KP = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_kp, LocalVar%PC_PitComT) + LocalVar%PC_KI = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_ki, LocalVar%PC_PitComT) + LocalVar%PC_KD = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_kd, LocalVar%PC_PitComT) + LocalVar%PC_TF = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_tf, LocalVar%PC_PitComT) + + ! Compute the current speed error and its integral w.r.t. time; saturate the + ! integral term using the pitch angle limits: + LocalVar%PC_SpdErr = CntrPar%PC_RefSpd - LocalVar%GenSpeedF ! Speed error + LocalVar%PC_PwrErr = CntrPar%VS_RtPwr - LocalVar%VS_GenPwr ! Power error + LocalVar%Y_MErr = LocalVar%Y_M + CntrPar%Y_MErrSet ! Yaw-alignment error + + ! Compute the pitch commands associated with the proportional and integral + ! gains: + LocalVar%PC_PitComT = PIController(LocalVar%PC_SpdErr, LocalVar%PC_KP, LocalVar%PC_KI, CntrPar%PC_SetPnt, LocalVar%PC_MaxPitVar, LocalVar%DT, CntrPar%PC_SetPnt, .FALSE., 2) ! + DFController(LocalVar%PC_SpdErr, LocalVar%PC_KD, LocalVar%PC_TF, LocalVar%DT, 1) + IF (CntrPar%VS_ControlMode == 1) THEN + LocalVar%PC_PitComT = LocalVar%PC_PitComT + PIController(LocalVar%PC_PwrErr, -4.0E-09, -4.0E-09, CntrPar%PC_SetPnt, LocalVar%PC_MaxPitVar, LocalVar%DT, CntrPar%PC_SetPnt, .FALSE., 5) + END IF + + ! Individual pitch control + IF ((CntrPar%IPC_ControlMode == 1) .OR. (CntrPar%Y_ControlMode == 2)) THEN + CALL IPC(CntrPar, LocalVar, objInst) + ELSE + LocalVar%IPC_PitComF = 0.0 ! THIS IS AN ARRAY!! + END IF + + ! Combine and saturate all pitch commands: + DO K = 1,LocalVar%NumBl ! Loop through all blades, add IPC contribution and limit pitch rate + LocalVar%PC_PitComT_IPC(K) = LocalVar%PC_PitComT + LocalVar%IPC_PitComF(K) ! Add the individual pitch command + LocalVar%PC_PitComT_IPC(K) = saturate(LocalVar%PC_PitComT_IPC(K), CntrPar%PC_MinPit, CntrPar%PC_MaxPit) ! Saturate the overall command using the pitch angle limits + + ! PitCom(K) = ratelimit(LocalVar%PC_PitComT_IPC(K), LocalVar%BlPitch(K), PC_MinRat, PC_MaxRat, LocalVar%DT) ! Saturate the overall command of blade K using the pitch rate limit + LocalVar%PitCom(K) = saturate(LocalVar%PC_PitComT_IPC(K), CntrPar%PC_MinPit, CntrPar%PC_MaxPit) ! Saturate the overall command using the pitch angle limits + LocalVar%PitCom(K) = LPFilter(LocalVar%PitCom(K), LocalVar%DT, CntrPar%CornerFreq, LocalVar%iStatus, .FALSE., objInst%instLPF) + END DO + + ! Command the pitch demanded from the last + ! call to the controller (See Appendix A of Bladed User's Guide): + avrSWAP(42) = LocalVar%PitCom(1) ! Use the command angles of all blades if using individual pitch + avrSWAP(43) = LocalVar%PitCom(2) ! " + avrSWAP(44) = LocalVar%PitCom(3) ! " + avrSWAP(45) = LocalVar%PitCom(1) ! Use the command angle of blade 1 if using collective pitch + END SUBROUTINE PitchControl + + SUBROUTINE VariableSpeedControl(avrSWAP, CntrPar, LocalVar, objInst) + + USE DRC_Types, ONLY : ControlParameters, LocalVariables, ObjectInstances + + REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + + TYPE(ControlParameters), INTENT(INOUT) :: CntrPar + TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + TYPE(ObjectInstances), INTENT(INOUT) :: objInst + + !.............................................................................................................................. + ! VARIABLE-SPEED TORQUE CONTROL: + !.............................................................................................................................. + avrSWAP(35) = 1.0 ! Generator contactor status: 1=main (high speed) variable-speed generator + avrSWAP(56) = 0.0 ! Torque override: 0=yes + + ! Filter the HSS (generator) speed measurement: + ! Apply Low-Pass Filter + LocalVar%GenSpeedF = SecLPFilter(LocalVar%GenSpeed, LocalVar%DT, CntrPar%CornerFreq, 0.7, LocalVar%iStatus, .FALSE., objInst%instSecLPF) ! This is the first instance of a second order LPFilter + + ! Compute the generator torque, which depends on which region we are in: + + LocalVar%VS_SpdErrAr = CntrPar%VS_RtSpd - LocalVar%GenSpeedF ! Current speed error - Above-rated PI-control + LocalVar%VS_SpdErrBr = CntrPar%VS_MinOM - LocalVar%GenSpeedF ! Current speed error - Below-rated PI-control + IF (LocalVar%PC_PitComT >= CntrPar%VS_Rgn3MP) THEN ! We are in region 3 + LocalVar%GenTrqAr = PIController(LocalVar%VS_SpdErrAr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_Rgn2MaxTq, CntrPar%VS_GenTrqArSatMax, LocalVar%DT, CntrPar%VS_GenTrqArSatMax, .TRUE., 1) + LocalVar%GenTrqBr = PIController(LocalVar%VS_SpdErrBr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_MinTq, CntrPar%VS_Rgn2MinTq, LocalVar%DT, CntrPar%VS_Rgn2MinTq, .TRUE., 4) + IF (CntrPar%VS_ControlMode == 1) THEN ! Constant power tracking + LocalVar%GenTrq = CntrPar%VS_RtPwr/LocalVar%GenSpeedF + ELSE ! Constant torque tracking + LocalVar%GenTrq = CntrPar%VS_RtTq + END IF + ELSE + LocalVar%GenTrqAr = PIController(LocalVar%VS_SpdErrAr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_Rgn2MaxTq, CntrPar%VS_GenTrqArSatMax, LocalVar%DT, CntrPar%VS_Rgn2MaxTq, .FALSE., 1) + LocalVar%GenTrqBr = PIController(LocalVar%VS_SpdErrBr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_MinTq, CntrPar%VS_Rgn2MinTq, LocalVar%DT, CntrPar%VS_Rgn2MinTq, .FALSE., 4) + IF (LocalVar%GenTrqAr >= CntrPar%VS_Rgn2MaxTq*1.01) THEN + LocalVar%GenTrq = LocalVar%GenTrqAr + CONTINUE + ELSEIF (LocalVar%GenTrqBr <= CntrPar%VS_Rgn2MinTq*0.99) THEN ! We are in region 1 1/2 + LocalVar%GenTrq = LocalVar%GenTrqBr + CONTINUE + ELSEIF (LocalVar%GenSpeedF < CntrPar%VS_MaxOM) THEN ! We are in region 2 - optimal torque is proportional to the square of the generator speed + LocalVar%GenTrq = CntrPar%VS_Rgn2K*LocalVar%GenSpeedF*LocalVar%GenSpeedF + ELSE ! We are in region 2 1/2 - simple induction generator transition region + LocalVar%GenTrq = CntrPar%VS_Rgn2MaxTq + END IF + END IF + + ! Saturate the commanded torque using the maximum torque limit: + LocalVar%GenTrq = MIN(LocalVar%GenTrq, CntrPar%VS_MaxTq) ! Saturate the command using the maximum torque limit + + ! Saturate the commanded torque using the torque rate limit: + IF (LocalVar%iStatus == 0) LocalVar%VS_LastGenTrq = LocalVar%GenTrq ! Initialize the value of LocalVar%VS_LastGenTrq on the first pass only + LocalVar%GenTrq = ratelimit(LocalVar%GenTrq, LocalVar%VS_LastGenTrq, -CntrPar%VS_MaxRat, CntrPar%VS_MaxRat, LocalVar%DT) ! Saturate the command using the torque rate limit + + ! Reset the value of LocalVar%VS_LastGenTrq to the current values: + LocalVar%VS_LastGenTrq = LocalVar%GenTrq + + ! Set the generator contactor status, avrSWAP(35), to main (high speed) + ! variable-speed generator, the torque override to yes, and command the + ! generator torque (See Appendix A of Bladed User's Guide): + avrSWAP(47) = LocalVar%VS_LastGenTrq ! Demanded generator torque + END SUBROUTINE VariableSpeedControl + + SUBROUTINE YawRateControl(avrSWAP, CntrPar, LocalVar, objInst) + + USE DRC_Types, ONLY : ControlParameters, LocalVariables, ObjectInstances + + REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + + TYPE(ControlParameters), INTENT(INOUT) :: CntrPar + TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + TYPE(ObjectInstances), INTENT(INOUT) :: objInst + + !.............................................................................................................................. + ! Yaw control + !.............................................................................................................................. + + IF (CntrPar%Y_ControlMode == 1) THEN + avrSWAP(29) = 0 ! Yaw control parameter: 0 = yaw rate control + IF (LocalVar%Time >= LocalVar%Y_YawEndT) THEN ! Check if the turbine is currently yawing + avrSWAP(48) = 0.0 ! Set yaw rate to zero + + LocalVar%Y_ErrLPFFast = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPFast, LocalVar%iStatus, .FALSE., objInst%instLPF) ! Fast low pass filtered yaw error with a frequency of 1 + LocalVar%Y_ErrLPFSlow = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPSlow, LocalVar%iStatus, .FALSE., objInst%instLPF) ! Slow low pass filtered yaw error with a frequency of 1/60 + + LocalVar%Y_AccErr = LocalVar%Y_AccErr + LocalVar%DT*SIGN(LocalVar%Y_ErrLPFFast**2, LocalVar%Y_ErrLPFFast) ! Integral of the fast low pass filtered yaw error + + IF (ABS(LocalVar%Y_AccErr) >= CntrPar%Y_ErrThresh) THEN ! Check if accumulated error surpasses the threshold + LocalVar%Y_YawEndT = ABS(LocalVar%Y_ErrLPFSlow/CntrPar%Y_Rate) + LocalVar%Time ! Yaw to compensate for the slow low pass filtered error + END IF + ELSE + avrSWAP(48) = SIGN(CntrPar%Y_Rate, LocalVar%Y_MErr) ! Set yaw rate to predefined yaw rate, the sign of the error is copied to the rate + LocalVar%Y_ErrLPFFast = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPFast, LocalVar%iStatus, .TRUE., objInst%instLPF) ! Fast low pass filtered yaw error with a frequency of 1 + LocalVar%Y_ErrLPFSlow = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPSlow, LocalVar%iStatus, .TRUE., objInst%instLPF) ! Slow low pass filtered yaw error with a frequency of 1/60 + LocalVar%Y_AccErr = 0.0 ! " + END IF + END IF + END SUBROUTINE YawRateControl + + SUBROUTINE IPC(CntrPar, LocalVar, objInst) + !------------------------------------------------------------------------------------------------------------------------------- + ! Individual pitch control subroutine + ! + ! Variable declaration and initialization + !------------------------------------------------------------------------------------------------------------------------------ + USE DRC_Types, ONLY : ControlParameters, LocalVariables, ObjectInstances + + ! Local variables + REAL(4) :: PitComIPC(3) + INTEGER(4) :: K ! Integer used to loop through turbine blades + REAL(4) :: axisTilt, axisYaw, axisYawF ! Direct axis and quadrature axis outputted by Coleman transform + REAL(4), SAVE :: IntAxisTilt, IntAxisYaw ! Integral of the direct axis and quadrature axis + REAL(4) :: IntAxisYawIPC ! IPC contribution with yaw-by-IPC component + REAL(4) :: Y_MErrF, Y_MErrF_IPC ! Unfiltered and filtered yaw alignment error [rad] + REAL(4) :: PitComIPC_woYaw(3) + + TYPE(ControlParameters), INTENT(INOUT) :: CntrPar + TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + TYPE(ObjectInstances), INTENT(INOUT) :: objInst + + !------------------------------------------------------------------------------------------------------------------------------ + ! Body + !------------------------------------------------------------------------------------------------------------------------------ + ! Calculates the commanded pitch angles. + ! NOTE: if it is required for this subroutine to be used multiple times (for 1p and 2p IPC for example), the saved variables + ! IntAxisTilt and IntAxisYaw need to be modified so that they support multiple instances (see LPFilter in the Filters module). + !------------------------------------------------------------------------------------------------------------------------------ + ! Filter rootMOOPs with notch filter + + !DO K = 1,LocalVar%NumBl + ! Instances 1-3 of the Notch Filter are reserved for this routine. + ! rootMOOPF(K) = LocalVar%rootMOOP(K) ! Notch filter currently not in use + !END DO + + ! Initialization + ! Set integrals to be 0 in the first time step + + IF(LocalVar%iStatus==0) THEN + IntAxisTilt = 0.0 + IntAxisYaw = 0.0 + END IF + + ! Pass rootMOOPs through the Coleman transform to get the direct and quadrature axis + CALL ColemanTransform(LocalVar%rootMOOP, LocalVar%Azimuth, axisTilt, axisYaw) + + ! High-pass filter the MBC yaw component and filter yaw alignment error, and compute the yaw-by-IPC contribution + IF (CntrPar%Y_ControlMode == 2) THEN + axisYawF = HPFilter(axisYaw, LocalVar%DT, CntrPar%IPC_omegaHP, LocalVar%iStatus, .FALSE., objInst%instHPF) + Y_MErrF = SecLPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%IPC_omegaLP, CntrPar%IPC_zetaLP, LocalVar%iStatus, .FALSE., objInst%instSecLPF) + Y_MErrF_IPC = PIController(Y_MErrF, CntrPar%Y_IPC_KP(1), CntrPar%Y_IPC_KI(1), -100.0, 100.0, LocalVar%DT, 0.0, .FALSE., 3) + ELSE + axisYawF = axisYaw + Y_MErrF = 0.0 + END IF + + ! Integrate the signal and multiply with the IPC gain + IF (CntrPar%IPC_ControlMode == 1) THEN + IntAxisTilt = IntAxisTilt + LocalVar%DT * CntrPar%IPC_KI * axisTilt + IntAxisYaw = IntAxisYaw + LocalVar%DT * CntrPar%IPC_KI * axisYawF + ELSE + IntAxisTilt = 0.0 + IntAxisYaw = 0.0 + END IF + + ! Add the yaw-by-IPC contribution + IntAxisYawIPC = IntAxisYaw + Y_MErrF_IPC + + ! Pass direct and quadrature axis through the inverse Coleman transform to get the commanded pitch angles + CALL ColemanTransformInverse(IntAxisTilt, IntAxisYawIPC, LocalVar%Azimuth, CntrPar%IPC_phi, PitComIPC) + + ! Filter PitComIPC with second order low pass filter + DO K = 1,LocalVar%NumBl + ! Instances 1-3 of the Second order Low-Pass Filter are reserved for this routine. + ! LocalVar%IPC_PitComF(K) = SecLPFilter(PitComIPC(K), LocalVar%DT, CntrPar%IPC_omegaLP, CntrPar%IPC_zetaLP, LocalVar%iStatus, K) + LocalVar%IPC_PitComF(K) = PitComIPC(K) + END DO + END SUBROUTINE IPC +END MODULE Controllers \ No newline at end of file diff --git a/Source/DISCON.f90 b/Source/DISCON.f90 index 8ba5a47b..d23dbe9d 100644 --- a/Source/DISCON.f90 +++ b/Source/DISCON.f90 @@ -1,24 +1,14 @@ !======================================================================= -! SUBROUTINE DISCON (avrSWAP, from_SC, to_SC, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME='DISCON') -SUBROUTINE DISCON (avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME='DISCON') +! SUBROUTINE DISCON(avrSWAP, from_SC, to_SC, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME='DISCON') +SUBROUTINE DISCON(avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAME='DISCON') +! DO NOT REMOVE or MODIFY LINES starting with "!DEC$" or "!GCC$" +! !DEC$ specifies attributes for IVF and !GCC$ specifies attributes for gfortran !DEC$ ATTRIBUTES DLLEXPORT :: DISCON - ! 06/09/2017 - - ! This Bladed-style DLL controller is used to implement a variable-speed - ! generator-torque controller, PI collective blade pitch controller, individual pitch - ! controller and yaw controller for the NREL Offshore 5MW baseline wind turbine. - ! This routine was extended by S.P. Mulders, J. Hoorneman and J. Govers of TU Delft. - ! The routine is based on the routine as written by J. Jonkman of NREL/NWTC. - - ! DO NOT REMOVE or MODIFY LINES starting with "!DEC$" or "!GCC$" - ! !DEC$ specifies attributes for IVF and !GCC$ specifies attributes for gfortran - USE, INTRINSIC :: ISO_C_Binding -USE :: ReadParameters -USE :: FunctionToolbox -USE :: Filters -USE DRC_Types, ONLY : ObjectInstances +USE :: DRC_Types +USE :: ReadSetParameters +USE :: Controllers IMPLICIT NONE #ifndef IMPLICIT_DLLEXPORT @@ -29,7 +19,7 @@ SUBROUTINE DISCON (avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAM ! Variable declaration and initialization !------------------------------------------------------------------------------------------------------------------------------ - ! Passed Variables: +! Passed Variables: !REAL(C_FLOAT), INTENT(IN) :: from_SC(*) ! DATA from the supercontroller !REAL(C_FLOAT), INTENT(INOUT) :: to_SC(*) ! DATA to the supercontroller @@ -38,392 +28,26 @@ SUBROUTINE DISCON (avrSWAP, aviFAIL, accINFILE, avcOUTNAME, avcMSG) BIND (C, NAM CHARACTER(KIND=C_CHAR), INTENT(IN) :: accINFILE(NINT(avrSWAP(50))) ! The name of the parameter input file CHARACTER(KIND=C_CHAR), INTENT(IN) :: avcOUTNAME(NINT(avrSWAP(51))) ! OUTNAME (Simulation RootName) CHARACTER(KIND=C_CHAR), INTENT(INOUT) :: avcMSG(NINT(avrSWAP(49))) ! MESSAGE (Message from DLL to simulation code [ErrMsg]) The message which will be displayed by the calling program if aviFAIL <> 0. +CHARACTER(SIZE(avcOUTNAME)-1) :: RootName ! a Fortran version of the input C string (not considered an array here) [subtract 1 for the C null-character] +CHARACTER(SIZE(avcMSG)-1) :: ErrMsg ! a Fortran version of the C string argument (not considered an array here) [subtract 1 for the C null-character] - ! Types - -TYPE(ObjectInstances) :: objInst - - ! Local Variables: - -REAL(4), SAVE :: LastGenTrq ! Commanded electrical generator torque the last time the controller was called, [Nm]. -REAL(4), SAVE :: PitComT ! Total command pitch based on the sum of the proportional and integral terms, [rad]. -REAL(4), SAVE :: PitComT_IPC(3) ! Total command pitch based on the sum of the proportional and integral terms, including IPC term [rad]. -REAL(4), SAVE :: Y_AccErr ! Accumulated yaw error [rad]. -REAL(4), SAVE :: Y_YawEndT ! Yaw end time, [s]. Indicates the time up until which yaw is active with a fixed rate -REAL(4), SAVE :: testValue ! TestValue - -INTEGER(4) :: I ! Generic index. -INTEGER(4) :: instLPF ! Instance counter for all first-order low-pass filters -INTEGER(4) :: K ! Loops through blades. -INTEGER(4), PARAMETER :: UnDb = 85 ! I/O unit for the debugging information -INTEGER(4), PARAMETER :: UnDb2 = 86 ! I/O unit for the debugging information - -LOGICAL(1), PARAMETER :: DbgOut = .FALSE. ! Flag to indicate whether to output debugging information - -CHARACTER(1), PARAMETER :: Tab = CHAR(9) ! The tab character. -CHARACTER(25), PARAMETER :: FmtDat = "(F8.3,99('"//Tab//"',ES10.3E2,:)) " ! The format of the debugging data - -CHARACTER(SIZE(avcOUTNAME)-1) :: RootName ! a Fortran version of the input C string (not considered an array here) [subtract 1 for the C null-character] -CHARACTER(SIZE(avcMSG)-1) :: ErrMsg ! a Fortran version of the C string argument (not considered an array here) [subtract 1 for the C null-character] - - - - ! Read avrSWAP array into derived types/variables -CALL ReadAvrSWAP(avrSWAP) - - ! Initialize aviFAIL to 0: -aviFAIL = 0 - - ! Initialize all filter instance counters at 1 -objInst%instLPF = 1 -objInst%instSecLPF = 1 -objInst%instHPF = 1 -objInst%instNotchSlopes = 1 -objInst%instNotch = 1 - - ! Read any External Controller Parameters specified in the User Interface - ! and initialize variables: -IF (LocalVar%iStatus == 0) THEN ! .TRUE. if we're on the first call to the DLL - - ! Inform users that we are using this user-defined routine: - aviFAIL = 1 - ErrMsg = ' '//NEW_LINE('A')// & - 'Running the Delft Research Controller (DRC) '//NEW_LINE('A')// & - 'A wind turbine controller for use in the scientific field '//NEW_LINE('A')// & - 'Written by S.P. Mulders, Jan-Willem van Wingerden '//NEW_LINE('A')// & - 'Delft University of Technology, The Netherlands '//NEW_LINE('A')// & - 'Visit our GitHub-page to contribute to this project: '//NEW_LINE('A')// & - 'https://github.com/TUDelft-DataDrivenControl ' - - CALL ReadControlParameterFileSub() - - ! Initialize testValue (debugging variable) - ! testValue = 0.4 - - ! Initialize the SAVEd variables: - ! NOTE: LastGenTrq, though SAVEd, is initialized in the torque controller - ! below for simplicity, not here. - LocalVar%PitCom = LocalVar%BlPitch ! This will ensure that the variable speed controller picks the correct control region and the pitch controller picks the correct gain on the first call - Y_AccErr = 0.0 ! This will ensure that the accumulated yaw error starts at zero - Y_YawEndT = -1.0 ! This will ensure that the initial yaw end time is lower than the actual time to prevent initial yawing - - !.............................................................................................................................. - ! Check validity of input parameters: - !.............................................................................................................................. - - IF (CntrPar%CornerFreq <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'CornerFreq must be greater than zero.' - ENDIF - - IF (LocalVar%DT <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'DT must be greater than zero.' - ENDIF - - IF (CntrPar%VS_CtInSp < 0.0) THEN - aviFAIL = -1 - ErrMsg = 'VS_CtInSp must not be negative.' - ENDIF - - IF (CntrPar%VS_MinOM <= CntrPar%VS_CtInSp) THEN - aviFAIL = -1 - ErrMsg = 'VS_MinOM must be greater than VS_CtInSp.' - ENDIF - - IF (CntrPar%VS_MaxRat <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'VS_MaxRat must be greater than zero.' - ENDIF - - IF (CntrPar%VS_RtTq < 0.0) THEN - aviFAIL = -1 - ErrMsg = 'VS_RtTw must not be negative.' - ENDIF - - IF (CntrPar%VS_Rgn2K < 0.0) THEN - aviFAIL = -1 - ErrMsg = 'VS_Rgn2K must not be negative.' - ENDIF - - IF (CntrPar%VS_MaxTq < CntrPar%VS_RtTq) THEN - aviFAIL = -1 - ErrMsg = 'VS_RtTq must not be greater than VS_MaxTq.' - ENDIF - - IF (CntrPar%VS_KP(1) > 0.0) THEN - aviFAIL = -1 - ErrMsg = 'VS_KP must be greater than zero.' - ENDIF - - IF (CntrPar%VS_KI(1) > 0.0) THEN - aviFAIL = -1 - ErrMsg = 'VS_KI must be greater than zero.' - ENDIF - - IF (CntrPar%PC_RefSpd <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'PC_RefSpd must be greater than zero.' - ENDIF - - IF (CntrPar%PC_MaxRat <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'PC_MaxRat must be greater than zero.' - ENDIF - - IF (CntrPar%PC_MinPit >= CntrPar%PC_MaxPit) THEN - aviFAIL = -1 - ErrMsg = 'PC_MinPit must be less than PC_MaxPit.' - ENDIF - - IF (CntrPar%IPC_KI <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'IPC_KI must be greater than zero.' - ENDIF - - IF (CntrPar%IPC_omegaLP <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'IPC_omegaLP must be greater than zero.' - ENDIF - - IF (CntrPar%IPC_omegaNotch <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'IPC_omegaNotch must be greater than zero.' - ENDIF - - IF (CntrPar%IPC_phi <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'IPC_phi must be greater than zero.' - ENDIF - - IF (CntrPar%IPC_zetaLP <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'IPC_zetaLP must be greater than zero.' - ENDIF - - IF (CntrPar%IPC_zetaNotch <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'IPC_zetaNotch must be greater than zero.' - ENDIF - - IF (CntrPar%Y_ErrThresh <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'Y_ErrThresh must be greater than zero.' - ENDIF - - IF (CntrPar%Y_Rate <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'CntrPar%Y_Rate must be greater than zero.' - ENDIF - - IF (CntrPar%Y_omegaLPFast <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'Y_omegaLPFast must be greater than zero.' - ENDIF - - IF (CntrPar%Y_omegaLPSlow <= 0.0) THEN - aviFAIL = -1 - ErrMsg = 'Y_omegaLPSlow must be greater than zero.' - ENDIF - - !.............................................................................................................................. - ! Initializing debug file - !.............................................................................................................................. - - ! If we're debugging, open the debug file and write the header: - IF (CntrPar%LoggingLevel > 0) THEN - OPEN (UnDb, FILE=TRIM(RootName)//'.dbg', STATUS='REPLACE') - WRITE (UnDb,'(A)') ' LocalVar%Time ' //Tab//'PitComT ' //Tab//'LocalVar%PC_SpdErr ' //Tab//'LocalVar%PC_KP ' //Tab//'LocalVar%PC_KI ' //Tab//'LocalVar%Y_M ' //Tab//'LocalVar%rootMOOP(1) '//Tab//'VS_RtPwr '//Tab//'LocalVar%GenTrq' - WRITE (UnDb,'(A)') ' (sec) ' //Tab//'(rad) ' //Tab//'(rad/s) '//Tab//'(-) ' //Tab//'(-) ' //Tab//'(rad) ' //Tab//'(?) ' //Tab//'(W) '//Tab//'(Nm) ' - END IF - - IF (CntrPar%LoggingLevel > 1) THEN - OPEN(UnDb2, FILE=TRIM(RootName)//'.dbg2', STATUS='REPLACE') - WRITE(UnDb2,'(/////)') - WRITE(UnDb2,'(A,85("'//Tab//'AvrSWAP(",I2,")"))') 'LocalVar%Time ', (i,i=1,85) - WRITE(UnDb2,'(A,85("'//Tab//'(-)"))') '(s)' - END IF - -ENDIF +TYPE(ControlParameters), SAVE :: CntrPar +TYPE(LocalVariables), SAVE :: LocalVar +TYPE(ObjectInstances), SAVE :: objInst !------------------------------------------------------------------------------------------------------------------------------ ! Main control calculations !------------------------------------------------------------------------------------------------------------------------------ - +! Read avrSWAP array into derived types/variables +CALL ReadAvrSWAP(avrSWAP, LocalVar) +CALL SetParameters(avrSWAP, aviFAIL, ErrMsg, SIZE(avcMSG), CntrPar, LocalVar, objInst) IF ((LocalVar%iStatus >= 0) .AND. (aviFAIL >= 0)) THEN ! Only compute control calculations if no error has occurred and we are not on the last time step - ! Abort if the user has not requested a pitch angle actuator (See Appendix A - ! of Bladed User's Guide): - IF (NINT(avrSWAP(10)) /= 0) THEN ! .TRUE. if a pitch angle actuator hasn't been requested - aviFAIL = -1 - ErrMsg = 'Pitch angle actuator not requested.' - ENDIF - - ! Set unused outputs to zero (See Appendix A of Bladed User's Guide): - avrSWAP(36) = 0.0 ! Shaft brake status: 0=off - avrSWAP(41) = 0.0 ! Demanded yaw actuator torque - avrSWAP(46) = 0.0 ! Demanded pitch rate (Collective pitch) - avrSWAP(65) = 0.0 ! Number of variables returned for logging - avrSWAP(72) = 0.0 ! Generator start-up resistance - avrSWAP(79) = 0.0 ! Request for loads: 0=none - avrSWAP(80) = 0.0 ! Variable slip current status - avrSWAP(81) = 0.0 ! Variable slip current demand - - ! Filter the HSS (generator) speed measurement: - ! Apply Low-Pass Filter - LocalVar%GenSpeedF = SecLPFilter(LocalVar%GenSpeed, LocalVar%DT, CntrPar%CornerFreq, 0.7, LocalVar%iStatus, .FALSE., objInst%instSecLPF) ! This is the first instance of a second order LPFilter - - ! Calculate yaw-alignment error - LocalVar%Y_MErr = LocalVar%Y_M + CntrPar%Y_MErrSet - - !.............................................................................................................................. - ! VARIABLE-SPEED TORQUE CONTROL: - !.............................................................................................................................. - - ! Compute the generator torque, which depends on which region we are in: - - LocalVar%VS_SpdErrAr = CntrPar%VS_RtSpd - LocalVar%GenSpeedF ! Current speed error - Above-rated PI-control - LocalVar%VS_SpdErrBr = CntrPar%VS_MinOM - LocalVar%GenSpeedF ! Current speed error - Below-rated PI-control - IF (PitComT >= CntrPar%VS_Rgn3MP) THEN ! We are in region 3 - LocalVar%GenTrqAr = PIController(LocalVar%VS_SpdErrAr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_Rgn2MaxTq, CntrPar%VS_GenTrqArSatMax, LocalVar%DT, CntrPar%VS_GenTrqArSatMax, .TRUE., 1) - LocalVar%GenTrqBr = PIController(LocalVar%VS_SpdErrBr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_MinTq, CntrPar%VS_Rgn2MinTq, LocalVar%DT, CntrPar%VS_Rgn2MinTq, .TRUE., 4) - IF (CntrPar%VS_ControlMode == 1) THEN ! Constant power tracking - LocalVar%GenTrq = CntrPar%VS_RtPwr/LocalVar%GenSpeedF - ELSE ! Constant torque tracking - LocalVar%GenTrq = CntrPar%VS_RtTq - END IF - ELSE - LocalVar%GenTrqAr = PIController(LocalVar%VS_SpdErrAr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_Rgn2MaxTq, CntrPar%VS_GenTrqArSatMax, LocalVar%DT, CntrPar%VS_Rgn2MaxTq, .FALSE., 1) - LocalVar%GenTrqBr = PIController(LocalVar%VS_SpdErrBr, CntrPar%VS_KP(1), CntrPar%VS_KI(1), CntrPar%VS_MinTq, CntrPar%VS_Rgn2MinTq, LocalVar%DT, CntrPar%VS_Rgn2MinTq, .FALSE., 4) - IF (LocalVar%GenTrqAr >= CntrPar%VS_Rgn2MaxTq*1.01) THEN - LocalVar%GenTrq = LocalVar%GenTrqAr - CONTINUE - ELSEIF (LocalVar%GenTrqBr <= CntrPar%VS_Rgn2MinTq*0.99) THEN ! We are in region 1 1/2 - LocalVar%GenTrq = LocalVar%GenTrqBr - CONTINUE - ELSEIF (LocalVar%GenSpeedF < CntrPar%VS_MaxOM) THEN ! We are in region 2 - optimal torque is proportional to the square of the generator speed - LocalVar%GenTrq = CntrPar%VS_Rgn2K*LocalVar%GenSpeedF*LocalVar%GenSpeedF - ELSE ! We are in region 2 1/2 - simple induction generator transition region - LocalVar%GenTrq = CntrPar%VS_Rgn2MaxTq - END IF - END IF - - ! Saturate the commanded torque using the maximum torque limit: - - LocalVar%GenTrq = MIN(LocalVar%GenTrq, CntrPar%VS_MaxTq) ! Saturate the command using the maximum torque limit - - ! Saturate the commanded torque using the torque rate limit: - IF (LocalVar%iStatus == 0) LastGenTrq = LocalVar%GenTrq ! Initialize the value of LastGenTrq on the first pass only - LocalVar%GenTrq = ratelimit(LocalVar%GenTrq, LastGenTrq, -CntrPar%VS_MaxRat, CntrPar%VS_MaxRat, LocalVar%DT) ! Saturate the command using the torque rate limit - - ! Reset the value of LastGenTrq to the current values: - LastGenTrq = LocalVar%GenTrq - - ! Set the generator contactor status, avrSWAP(35), to main (high speed) - ! variable-speed generator, the torque override to yes, and command the - ! generator torque (See Appendix A of Bladed User's Guide): - avrSWAP(35) = 1.0 ! Generator contactor status: 1=main (high speed) variable-speed generator - avrSWAP(56) = 0.0 ! Torque override: 0=yes - avrSWAP(47) = LastGenTrq ! Demanded generator torque - - !.............................................................................................................................. - ! Pitch control - !.............................................................................................................................. - - IF (CntrPar%VS_ControlMode == 0 .AND. LocalVar%GenTrq >= CntrPar%PC_RtTq99) THEN - LocalVar%PC_MaxPitVar = CntrPar%PC_MaxPit - ELSEIF (CntrPar%VS_ControlMode == 1 .AND. LocalVar%GenTrqAr >= CntrPar%VS_GenTrqArSatMax*0.99) THEN - LocalVar%PC_MaxPitVar = CntrPar%PC_MaxPit - ELSE - LocalVar%PC_MaxPitVar = CntrPar%PC_SetPnt - END IF - - ! Compute the gain scheduling correction factor based on the previously - ! commanded pitch angle for blade 1: - LocalVar%PC_KP = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_kp, PitComT) - LocalVar%PC_KI = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_ki, PitComT) - LocalVar%PC_KD = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_kd, PitComT) - LocalVar%PC_TF = interp1d(CntrPar%PC_GS_angles, CntrPar%PC_GS_tf, PitComT) - - ! Compute the current speed error and its integral w.r.t. time; saturate the - ! integral term using the pitch angle limits: - LocalVar%PC_SpdErr = CntrPar%PC_RefSpd - LocalVar%GenSpeedF ! Speed error - LocalVar%PC_PwrErr = CntrPar%VS_RtPwr - LocalVar%VS_GenPwr ! Power error - - ! Compute the pitch commands associated with the proportional and integral - ! gains: - ! PitComT = NotchFilter(PC_SpdErr, LocalVar%DT, 1.59, 0.01, 0.2, LocalVar%iStatus, 1) - PitComT = PIController(LocalVar%PC_SpdErr, LocalVar%PC_KP, LocalVar%PC_KI, CntrPar%PC_SetPnt, LocalVar%PC_MaxPitVar, LocalVar%DT, CntrPar%PC_SetPnt, .FALSE., 2) ! + DFController(LocalVar%PC_SpdErr, LocalVar%PC_KD, LocalVar%PC_TF, LocalVar%DT, 1) - IF (CntrPar%VS_ControlMode == 1) THEN - PitComT = PitComT + PIController(LocalVar%PC_PwrErr, -4.0E-09, -4.0E-09, CntrPar%PC_SetPnt, LocalVar%PC_MaxPitVar, LocalVar%DT, CntrPar%PC_SetPnt, .FALSE., 5) - END IF - - ! Individual pitch control - IF ((CntrPar%IPC_ControlMode == 1) .OR. (CntrPar%Y_ControlMode == 2)) THEN - CALL IPC(LocalVar, CntrPar, objInst) - ELSE - LocalVar%IPC_PitComF = 0.0 ! THIS IS AN ARRAY!! - END IF - - ! Combine and saturate all pitch commands: - DO K = 1,LocalVar%NumBl ! Loop through all blades, add IPC contribution and limit pitch rate - PitComT_IPC(K) = PitComT + LocalVar%IPC_PitComF(K) ! Add the individual pitch command - PitComT_IPC(K) = saturate(PitComT_IPC(K), CntrPar%PC_MinPit, CntrPar%PC_MaxPit) ! Saturate the overall command using the pitch angle limits - - ! PitCom(K) = ratelimit(PitComT_IPC(K), LocalVar%BlPitch(K), PC_MinRat, PC_MaxRat, LocalVar%DT) ! Saturate the overall command of blade K using the pitch rate limit - LocalVar%PitCom(K) = saturate(PitComT_IPC(K), CntrPar%PC_MinPit, CntrPar%PC_MaxPit) ! Saturate the overall command using the pitch angle limits - LocalVar%PitCom(K) = LPFilter(LocalVar%PitCom(K), LocalVar%DT, CntrPar%CornerFreq, LocalVar%iStatus, .FALSE., objInst%instLPF) - END DO - - ! Set the pitch override to yes and command the pitch demanded from the last - ! call to the controller (See Appendix A of Bladed User's Guide): - - avrSWAP(55) = 0.0 ! Pitch override: 0=yes - - avrSWAP(42) = LocalVar%PitCom(1) ! Use the command angles of all blades if using individual pitch - avrSWAP(43) = LocalVar%PitCom(2) ! " - avrSWAP(44) = LocalVar%PitCom(3) ! " - - avrSWAP(45) = LocalVar%PitCom(1) ! Use the command angle of blade 1 if using collective pitch - - !.............................................................................................................................. - ! Yaw control - !.............................................................................................................................. - - IF (CntrPar%Y_ControlMode == 1) THEN - avrSWAP(29) = 0 ! Yaw control parameter: 0 = yaw rate control - IF (LocalVar%Time >= Y_YawEndT) THEN ! Check if the turbine is currently yawing - avrSWAP(48) = 0.0 ! Set yaw rate to zero - - LocalVar%Y_ErrLPFFast = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPFast, LocalVar%iStatus, .FALSE., objInst%instLPF) ! Fast low pass filtered yaw error with a frequency of 1 - LocalVar%Y_ErrLPFSlow = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPSlow, LocalVar%iStatus, .FALSE., objInst%instLPF) ! Slow low pass filtered yaw error with a frequency of 1/60 - - Y_AccErr = Y_AccErr + LocalVar%DT*SIGN(LocalVar%Y_ErrLPFFast**2, LocalVar%Y_ErrLPFFast) ! Integral of the fast low pass filtered yaw error - - IF (ABS(Y_AccErr) >= CntrPar%Y_ErrThresh) THEN ! Check if accumulated error surpasses the threshold - Y_YawEndT = ABS(LocalVar%Y_ErrLPFSlow/CntrPar%Y_Rate) + LocalVar%Time ! Yaw to compensate for the slow low pass filtered error - END IF - ELSE - avrSWAP(48) = SIGN(CntrPar%Y_Rate, LocalVar%Y_MErr) ! Set yaw rate to predefined yaw rate, the sign of the error is copied to the rate - LocalVar%Y_ErrLPFFast = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPFast, LocalVar%iStatus, .TRUE., objInst%instLPF) ! Fast low pass filtered yaw error with a frequency of 1 - LocalVar%Y_ErrLPFSlow = LPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%Y_omegaLPSlow, LocalVar%iStatus, .TRUE., objInst%instLPF) ! Slow low pass filtered yaw error with a frequency of 1/60 - Y_AccErr = 0.0 ! " - END IF - END IF - !.............................................................................................................................. - ! Output debugging information if requested: - IF (CntrPar%LoggingLevel > 0) THEN - WRITE (UnDb,FmtDat) LocalVar%Time, PitComT, LocalVar%PC_SpdErr, LocalVar%PC_KP, LocalVar%PC_KI, LocalVar%Y_MErr, LocalVar%rootMOOP(1), CntrPar%VS_RtPwr, LocalVar%GenTrq - END IF - IF (CntrPar%LoggingLevel > 1) THEN - WRITE (UnDb2,FmtDat) LocalVar%Time, avrSWAP(1:85) - END IF - - -ENDIF + CALL VariableSpeedControl(avrSWAP, CntrPar, LocalVar, objInst) + CALL PitchControl(avrSWAP, CntrPar, LocalVar, objInst) +! CALL YawRateControl(avrSWAP, CntrPar, LocalVar, objInst) +END IF avcMSG = TRANSFER(TRIM(ErrMsg)//C_NULL_CHAR, avcMSG, SIZE(avcMSG)) - RETURN - END SUBROUTINE DISCON diff --git a/Source/DRC_Types.f90 b/Source/DRC_Types.f90 index 83ff3770..9e04e5d5 100644 --- a/Source/DRC_Types.f90 +++ b/Source/DRC_Types.f90 @@ -2,109 +2,109 @@ MODULE DRC_Types IMPLICIT NONE -REAL(4), PARAMETER :: RPS2RPM = 9.5492966 ! Factor to convert radians per second to revolutions per minute. -REAL(4), PARAMETER :: R2D = 57.295780 ! Factor to convert radians to degrees. - - TYPE, PUBLIC :: ControlParameters - REAL(4) :: CornerFreq ! Corner frequency (-3dB point) in the first-order low-pass filter, [rad/s] - INTEGER(4) :: LoggingLevel ! 0 = write no debug files, 1 = write standard output .dbg-file, 2 = write standard output .dbg-file and complete avrSWAP-array .dbg2-file - REAL(4) :: IPC_KI ! Integral gain for the individual pitch controller, [-]. 8E-10 - INTEGER(4) :: IPC_ControlMode ! Turn Individual Pitch Control (IPC) for fatigue load reductions (pitch contribution) on = 1/off = 0 - REAL(4) :: IPC_omegaHP ! High-pass filter cut-in frequency used to separate yaw-by-IPC contribution from blade load reduction contribution, [rad/s]. - REAL(4) :: IPC_omegaLP ! Low-pass filter corner frequency for the individual pitch controller, [rad/s]. - REAL(4) :: IPC_omegaNotch ! Notch filter corner frequency for the individual pitch controller, [rad/s]. - REAL(4) :: IPC_phi ! Phase offset added to the azimuth angle for the individual pitch controller, [rad]. - REAL(4) :: IPC_zetaHP ! High-pass filter damping value, [-]. - REAL(4) :: IPC_zetaLP ! Low-pass filter damping factor for the individual pitch controller, [-]. - REAL(4) :: IPC_zetaNotch ! Notch filter damping factor for the individual pitch controller, [-]. - INTEGER(4) :: PC_GS_n ! Amount of gain-scheduling table entries - REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_angles ! Gain-schedule table: pitch angles - REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_kp ! Gain-schedule table: pitch controller kp gains - REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_ki ! Gain-schedule table: pitch controller ki gains - REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_kd ! Gain-schedule table: pitch controller kd gains - REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_tf ! Gain-schedule table: pitch controller tf gains (derivative filter) - REAL(4) :: PC_MaxPit ! Maximum physical pitch limit, [rad]. - REAL(4) :: PC_MinPit ! Minimum physical pitch limit, [rad]. - REAL(4) :: PC_MaxRat ! Maximum pitch rate (in absolute value) in pitch controller, [rad/s]. - REAL(4) :: PC_MinRat ! Minimum pitch rate (in absolute value) in pitch controller, [rad/s]. - REAL(4) :: PC_RefSpd ! Desired (reference) HSS speed for pitch controller, [rad/s]. - REAL(4) :: PC_SetPnt ! Record 5: Below-rated pitch angle set-point (deg) [used only with Bladed Interface] - REAL(4) :: PC_Switch ! Angle above lowest minimum pitch angle for switch [rad] - INTEGER(4) :: VS_ControlMode ! Generator torque control mode in above rated conditions, 0 = constant torque / 1 = constant power - REAL(4) :: VS_CtInSp ! Transitional generator speed (HSS side) between regions 1 and 1 1/2, [rad/s]. - REAL(4) :: VS_GenTrqArSatMax ! Above rated generator torque PI control saturation, [Nm] -- 212900 - REAL(4) :: VS_MaxOM ! Optimal mode maximum speed, [rad/s]. - REAL(4) :: VS_MaxRat ! Maximum torque rate (in absolute value) in torque controller, [Nm/s]. - REAL(4) :: VS_MaxTq ! Maximum generator torque in Region 3 (HSS side), [Nm]. -- chosen to be 10% above VS_RtTq - REAL(4) :: VS_MinTq ! Minimum generator (HSS side), [Nm]. - REAL(4) :: VS_MinOM ! Optimal mode minimum speed, [rad/s] - REAL(4) :: VS_Rgn2K ! Generator torque constant in Region 2 (HSS side), N-m/(rad/s)^2 - REAL(4) :: VS_RtPwr ! Wind turbine rated power [W] - REAL(4) :: VS_RtTq ! Rated torque, [Nm]. - REAL(4) :: VS_RtSpd ! Rated generator speed [rad/s] - INTEGER(4) :: VS_n ! Number of controller gains - REAL(4), DIMENSION(:), ALLOCATABLE :: VS_KP ! Proportional gain for generator PI torque controller, used in the transitional 2.5 region - REAL(4), DIMENSION(:), ALLOCATABLE :: VS_KI ! Integral gain for generator PI torque controller, used in the transitional 2.5 region - INTEGER(4) :: Y_ControlMode ! Yaw control mode: (0 = no yaw control, 1 = yaw rate control, 2 = yaw-by-IPC) - REAL(4) :: Y_ErrThresh ! Error threshold [rad]. Turbine begins to yaw when it passes this. (104.71975512) -- 1.745329252 - INTEGER(4) :: Y_IPC_n ! Number of controller gains (yaw-by-IPC) - REAL(4), DIMENSION(:), ALLOCATABLE :: Y_IPC_KP ! Yaw-by-IPC proportional controller gain Kp - REAL(4), DIMENSION(:), ALLOCATABLE :: Y_IPC_KI ! Yaw-by-IPC integral controller gain Ki - REAL(4) :: Y_MErrSet ! Yaw alignment error, setpoint [rad] - REAL(4) :: Y_omegaLPFast ! Corner frequency fast low pass filter, 1.0 [Hz] - REAL(4) :: Y_omegaLPSlow ! Corner frequency slow low pass filter, 1/60 [Hz] - REAL(4) :: Y_Rate ! Yaw rate [rad/s] - - REAL(4) :: PC_RtTq99 ! 99% of the rated torque value, using for switching between pitch and torque control, [Nm]. - REAL(4) :: VS_Rgn2MaxTq ! Maximum torque at the end of the below-rated region 2, [Nm] - REAL(4) :: VS_Rgn2MinTq ! Minimum torque at the beginning of the below-rated region 2, [Nm] - REAL(4) :: VS_Rgn3MP ! Minimum pitch angle at which the torque is computed as if we are in region 3 regardless of the generator speed, [rad]. - END TYPE ControlParameters - - TYPE, PUBLIC :: LocalVariables - ! From avrSWAP - INTEGER(4) :: iStatus - REAL(4) :: Time - REAL(4) :: DT - REAL(4) :: VS_GenPwr - REAL(4) :: GenSpeed - REAL(4) :: Y_M - REAL(4) :: HorWindV - REAL(4) :: rootMOOP(3) - REAL(4) :: BlPitch(3) - REAL(4) :: Azimuth - INTEGER(4) :: NumBl - - ! Internal controller variables - REAL(4) :: GenSpeedF ! Filtered HSS (generator) speed [rad/s]. - REAL(4) :: GenTrq ! Electrical generator torque, [Nm]. - REAL(4) :: GenTrqAr ! Electrical generator torque, for above-rated PI-control [Nm]. - REAL(4) :: GenTrqBr ! Electrical generator torque, for below-rated PI-control [Nm]. - REAL(4) :: IPC_PitComF(3) ! Commanded pitch of each blade as calculated by the individual pitch controller, F stands for low pass filtered, [rad]. - REAL(4) :: PC_KP ! Proportional gain for pitch controller at rated pitch (zero), [s]. - REAL(4) :: PC_KI ! Integral gain for pitch controller at rated pitch (zero), [-]. - REAL(4) :: PC_KD ! Differential gain for pitch controller at rated pitch (zero), [-]. - REAL(4) :: PC_TF ! First-order filter parameter for derivative action - REAL(4) :: PC_MaxPitVar ! Maximum pitch setting in pitch controller (variable) [rad]. - REAL(4) :: PC_PwrErr ! Power error with respect to rated power [W] - REAL(4) :: PC_SpdErr ! Current speed error (pitch control) [rad/s]. - REAL(4) :: PitCom(3) ! Commanded pitch of each blade the last time the controller was called, [rad]. - REAL(4) :: VS_SpdErrAr ! Current speed error (generator torque control) [rad/s]. - REAL(4) :: VS_SpdErrBr ! Current speed error (generator torque control) [rad/s]. - REAL(4) :: Y_ErrLPFFast ! Filtered yaw error by fast low pass filter [rad]. - REAL(4) :: Y_ErrLPFSlow ! Filtered yaw error by slow low pass filter [rad]. - REAL(4) :: Y_MErr ! Measured yaw error, measured + setpoint [rad] - END TYPE LocalVariables +TYPE, PUBLIC :: ControlParameters + REAL(4) :: CornerFreq ! Corner frequency (-3dB point) in the first-order low-pass filter, [rad/s] + INTEGER(4) :: LoggingLevel ! 0 = write no debug files, 1 = write standard output .dbg-file, 2 = write standard output .dbg-file and complete avrSWAP-array .dbg2-file + REAL(4) :: IPC_KI ! Integral gain for the individual pitch controller, [-]. 8E-10 + INTEGER(4) :: IPC_ControlMode ! Turn Individual Pitch Control (IPC) for fatigue load reductions (pitch contribution) on = 1/off = 0 + REAL(4) :: IPC_omegaHP ! High-pass filter cut-in frequency used to separate yaw-by-IPC contribution from blade load reduction contribution, [rad/s]. + REAL(4) :: IPC_omegaLP ! Low-pass filter corner frequency for the individual pitch controller, [rad/s]. + REAL(4) :: IPC_omegaNotch ! Notch filter corner frequency for the individual pitch controller, [rad/s]. + REAL(4) :: IPC_phi ! Phase offset added to the azimuth angle for the individual pitch controller, [rad]. + REAL(4) :: IPC_zetaHP ! High-pass filter damping value, [-]. + REAL(4) :: IPC_zetaLP ! Low-pass filter damping factor for the individual pitch controller, [-]. + REAL(4) :: IPC_zetaNotch ! Notch filter damping factor for the individual pitch controller, [-]. + INTEGER(4) :: PC_GS_n ! Amount of gain-scheduling table entries + REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_angles ! Gain-schedule table: pitch angles + REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_kp ! Gain-schedule table: pitch controller kp gains + REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_ki ! Gain-schedule table: pitch controller ki gains + REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_kd ! Gain-schedule table: pitch controller kd gains + REAL(4), DIMENSION(:), ALLOCATABLE :: PC_GS_tf ! Gain-schedule table: pitch controller tf gains (derivative filter) + REAL(4) :: PC_MaxPit ! Maximum physical pitch limit, [rad]. + REAL(4) :: PC_MinPit ! Minimum physical pitch limit, [rad]. + REAL(4) :: PC_MaxRat ! Maximum pitch rate (in absolute value) in pitch controller, [rad/s]. + REAL(4) :: PC_MinRat ! Minimum pitch rate (in absolute value) in pitch controller, [rad/s]. + REAL(4) :: PC_RefSpd ! Desired (reference) HSS speed for pitch controller, [rad/s]. + REAL(4) :: PC_SetPnt ! Record 5: Below-rated pitch angle set-point (deg) [used only with Bladed Interface] + REAL(4) :: PC_Switch ! Angle above lowest minimum pitch angle for switch [rad] + INTEGER(4) :: VS_ControlMode ! Generator torque control mode in above rated conditions, 0 = constant torque / 1 = constant power + REAL(4) :: VS_CtInSp ! Transitional generator speed (HSS side) between regions 1 and 1 1/2, [rad/s]. + REAL(4) :: VS_GenTrqArSatMax ! Above rated generator torque PI control saturation, [Nm] -- 212900 + REAL(4) :: VS_MaxOM ! Optimal mode maximum speed, [rad/s]. + REAL(4) :: VS_MaxRat ! Maximum torque rate (in absolute value) in torque controller, [Nm/s]. + REAL(4) :: VS_MaxTq ! Maximum generator torque in Region 3 (HSS side), [Nm]. -- chosen to be 10% above VS_RtTq + REAL(4) :: VS_MinTq ! Minimum generator (HSS side), [Nm]. + REAL(4) :: VS_MinOM ! Optimal mode minimum speed, [rad/s] + REAL(4) :: VS_Rgn2K ! Generator torque constant in Region 2 (HSS side), N-m/(rad/s)^2 + REAL(4) :: VS_RtPwr ! Wind turbine rated power [W] + REAL(4) :: VS_RtTq ! Rated torque, [Nm]. + REAL(4) :: VS_RtSpd ! Rated generator speed [rad/s] + INTEGER(4) :: VS_n ! Number of controller gains + REAL(4), DIMENSION(:), ALLOCATABLE :: VS_KP ! Proportional gain for generator PI torque controller, used in the transitional 2.5 region + REAL(4), DIMENSION(:), ALLOCATABLE :: VS_KI ! Integral gain for generator PI torque controller, used in the transitional 2.5 region + INTEGER(4) :: Y_ControlMode ! Yaw control mode: (0 = no yaw control, 1 = yaw rate control, 2 = yaw-by-IPC) + REAL(4) :: Y_ErrThresh ! Error threshold [rad]. Turbine begins to yaw when it passes this. (104.71975512) -- 1.745329252 + INTEGER(4) :: Y_IPC_n ! Number of controller gains (yaw-by-IPC) + REAL(4), DIMENSION(:), ALLOCATABLE :: Y_IPC_KP ! Yaw-by-IPC proportional controller gain Kp + REAL(4), DIMENSION(:), ALLOCATABLE :: Y_IPC_KI ! Yaw-by-IPC integral controller gain Ki + REAL(4) :: Y_MErrSet ! Yaw alignment error, setpoint [rad] + REAL(4) :: Y_omegaLPFast ! Corner frequency fast low pass filter, 1.0 [Hz] + REAL(4) :: Y_omegaLPSlow ! Corner frequency slow low pass filter, 1/60 [Hz] + REAL(4) :: Y_Rate ! Yaw rate [rad/s] - TYPE, PUBLIC :: ObjectInstances - INTEGER(4) :: instLPF - INTEGER(4) :: instSecLPF - INTEGER(4) :: instHPF - INTEGER(4) :: instNotchSlopes - INTEGER(4) :: instNotch - END TYPE ObjectInstances + REAL(4) :: PC_RtTq99 ! 99% of the rated torque value, using for switching between pitch and torque control, [Nm]. + REAL(4) :: VS_Rgn2MaxTq ! Maximum torque at the end of the below-rated region 2, [Nm] + REAL(4) :: VS_Rgn2MinTq ! Minimum torque at the beginning of the below-rated region 2, [Nm] + REAL(4) :: VS_Rgn3MP ! Minimum pitch angle at which the torque is computed as if we are in region 3 regardless of the generator speed, [rad]. +END TYPE ControlParameters + +TYPE, PUBLIC :: LocalVariables + ! From avrSWAP + INTEGER(4) :: iStatus + REAL(4) :: Time + REAL(4) :: DT + REAL(4) :: VS_GenPwr + REAL(4) :: GenSpeed + REAL(4) :: Y_M + REAL(4) :: HorWindV + REAL(4) :: rootMOOP(3) + REAL(4) :: BlPitch(3) + REAL(4) :: Azimuth + INTEGER(4) :: NumBl - !TYPE, PUBLIC :: PersistentVariables - !END TYPE PersistentVariables + ! Internal controller variables + REAL(4) :: GenSpeedF ! Filtered HSS (generator) speed [rad/s]. + REAL(4) :: GenTrq ! Electrical generator torque, [Nm]. + REAL(4) :: GenTrqAr ! Electrical generator torque, for above-rated PI-control [Nm]. + REAL(4) :: GenTrqBr ! Electrical generator torque, for below-rated PI-control [Nm]. + INTEGER(4) :: GlobalState ! Current global state to determine the behavior of the different controllers [-]. + REAL(4) :: IPC_PitComF(3) ! Commanded pitch of each blade as calculated by the individual pitch controller, F stands for low pass filtered, [rad]. + REAL(4) :: PC_KP ! Proportional gain for pitch controller at rated pitch (zero), [s]. + REAL(4) :: PC_KI ! Integral gain for pitch controller at rated pitch (zero), [-]. + REAL(4) :: PC_KD ! Differential gain for pitch controller at rated pitch (zero), [-]. + REAL(4) :: PC_TF ! First-order filter parameter for derivative action + REAL(4) :: PC_MaxPitVar ! Maximum pitch setting in pitch controller (variable) [rad]. + REAL(4) :: PC_PitComT ! Total command pitch based on the sum of the proportional and integral terms, [rad]. + REAL(4) :: PC_PitComT_IPC(3) ! Total command pitch based on the sum of the proportional and integral terms, including IPC term [rad]. + REAL(4) :: PC_PwrErr ! Power error with respect to rated power [W] + REAL(4) :: PC_SpdErr ! Current speed error (pitch control) [rad/s]. + REAL(4) :: PitCom(3) ! Commanded pitch of each blade the last time the controller was called, [rad]. + INTEGER(4) :: TestType ! Test variable, no use + REAL(4) :: VS_LastGenTrq ! Commanded electrical generator torque the last time the controller was called, [Nm]. + REAL(4) :: VS_SpdErrAr ! Current speed error (generator torque control) [rad/s]. + REAL(4) :: VS_SpdErrBr ! Current speed error (generator torque control) [rad/s]. + REAL(4) :: Y_AccErr ! Accumulated yaw error [rad]. + REAL(4) :: Y_ErrLPFFast ! Filtered yaw error by fast low pass filter [rad]. + REAL(4) :: Y_ErrLPFSlow ! Filtered yaw error by slow low pass filter [rad]. + REAL(4) :: Y_MErr ! Measured yaw error, measured + setpoint [rad]. + REAL(4) :: Y_YawEndT ! Yaw end time, [s]. Indicates the time up until which yaw is active with a fixed rate +END TYPE LocalVariables +TYPE, PUBLIC :: ObjectInstances + INTEGER(4) :: instLPF + INTEGER(4) :: instSecLPF + INTEGER(4) :: instHPF + INTEGER(4) :: instNotchSlopes + INTEGER(4) :: instNotch +END TYPE ObjectInstances END MODULE DRC_Types \ No newline at end of file diff --git a/Source/FunctionToolbox.f90 b/Source/FunctionToolbox.f90 index 0d845624..f2d75770 100644 --- a/Source/FunctionToolbox.f90 +++ b/Source/FunctionToolbox.f90 @@ -1,6 +1,7 @@ ! This module contains basic functions MODULE FunctionToolbox +USE Constants IMPLICIT NONE CONTAINS @@ -180,4 +181,141 @@ END FUNCTION DFController ! !END FUNCTION PRBSgen !------------------------------------------------------------------------------------------------------------------------------- -END MODULE FunctionToolbox + ! Stata machine, determines the state of the wind turbine to determine the corresponding control actions + ! States: + ! - 1, idling, wind ad rotor speed too low for start-up: set pitch to vane position and torque to minimum + ! - 2, start-up mode, set pitch demand to start-up pitch angle for maximum aerodynamic torque and torque demand to minimum + ! - 3, start-up2normal + ! - 4, Region 1.5 operation, torque control to keep the rotor at cut-in speed towards the Cp-max operational curve + ! - 5, Region 2, operation, maximum rotor power efficiency (Cp-max) tracking, keep TSR constant at a fixed fine-pitch angle + ! - 6, Region 2.5, transition between below and above-rated operating conditions (near-rated region) using PI torque control + ! - 7, Region 3, above-rated operation using pitch control + !REAL FUNCTION StateMachine(LocalVar, CntrPar) + ! + ! USE DRC_Types, ONLY : LocalVariables, ControlParameters + ! + ! IMPLICIT NONE + ! + ! ! Inputs + ! TYPE(ControlParameters), INTENT(IN) :: CntrPar + ! + ! ! Inputs/outputs + ! TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + ! + ! ! Local + ! + ! + !END FUNCTION StateMachine + !------------------------------------------------------------------------------------------------------------------------------- + SUBROUTINE Debug(LocalVar, CntrPar, avrSWAP, RootName, size_avcOUTNAME) + USE, INTRINSIC :: ISO_C_Binding + USE DRC_Types, ONLY : LocalVariables, ControlParameters + + IMPLICIT NONE + + TYPE(ControlParameters), INTENT(IN) :: CntrPar + TYPE(LocalVariables), INTENT(IN) :: LocalVar + + INTEGER(4), INTENT(IN) :: size_avcOUTNAME + INTEGER(4) :: I ! Generic index. + CHARACTER(1), PARAMETER :: Tab = CHAR(9) ! The tab character. + CHARACTER(25), PARAMETER :: FmtDat = "(F8.3,99('"//Tab//"',ES10.3E2,:)) " ! The format of the debugging data + INTEGER(4), PARAMETER :: UnDb = 85 ! I/O unit for the debugging information + INTEGER(4), PARAMETER :: UnDb2 = 86 ! I/O unit for the debugging information, avrSWAP + REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + CHARACTER(size_avcOUTNAME-1), INTENT(IN) :: RootName ! a Fortran version of the input C string (not considered an array here) [subtract 1 for the C null-character] + + !.............................................................................................................................. + ! Initializing debug file + !.............................................................................................................................. + IF (LocalVar%iStatus == 0) THEN ! .TRUE. if we're on the first call to the DLL + ! If we're debugging, open the debug file and write the header: + IF (CntrPar%LoggingLevel > 0) THEN + OPEN (UnDb, FILE=TRIM(RootName)//'.dbg', STATUS='REPLACE') + WRITE (UnDb,'(A)') ' LocalVar%Time ' //Tab//'LocalVar%PC_PitComT ' //Tab//'LocalVar%PC_SpdErr ' //Tab//'LocalVar%PC_KP ' //Tab//'LocalVar%PC_KI ' //Tab//'LocalVar%Y_M ' //Tab//'LocalVar%rootMOOP(1) '//Tab//'VS_RtPwr '//Tab//'LocalVar%GenTrq' + WRITE (UnDb,'(A)') ' (sec) ' //Tab//'(rad) ' //Tab//'(rad/s) '//Tab//'(-) ' //Tab//'(-) ' //Tab//'(rad) ' //Tab//'(?) ' //Tab//'(W) '//Tab//'(Nm) ' + END IF + + IF (CntrPar%LoggingLevel > 1) THEN + OPEN(UnDb2, FILE=TRIM(RootName)//'.dbg2', STATUS='REPLACE') + WRITE(UnDb2,'(/////)') + WRITE(UnDb2,'(A,85("'//Tab//'AvrSWAP(",I2,")"))') 'LocalVar%Time ', (i,i=1,85) + WRITE(UnDb2,'(A,85("'//Tab//'(-)"))') '(s)' + END IF + ELSE + !.............................................................................................................................. + ! Output debugging information if requested: + IF (CntrPar%LoggingLevel > 0) THEN + WRITE (UnDb,FmtDat) LocalVar%Time, LocalVar%PC_PitComT, LocalVar%PC_SpdErr, LocalVar%PC_KP, LocalVar%PC_KI, LocalVar%Y_MErr, LocalVar%rootMOOP(1), CntrPar%VS_RtPwr, LocalVar%GenTrq + END IF + + IF (CntrPar%LoggingLevel > 1) THEN + WRITE (UnDb2,FmtDat) LocalVar%Time, avrSWAP(1:85) + END IF + END IF + + IF (MODULO(LocalVar%Time, 10.0) == 0.0) THEN + !LocalVar%TestType = LocalVar%TestType + 10 + !PRINT *, LocalVar%TestType + END IF + END SUBROUTINE Debug + !------------------------------------------------------------------------------------------------------------------------------- + !The Coleman or d-q axis transformation transforms the root out of plane bending moments of each turbine blade + !to a direct axis and a quadrature axis + SUBROUTINE ColemanTransform(rootMOOP, aziAngle, axisTilt, axisYaw) + !............................................................................................................................... + + IMPLICIT NONE + + ! Inputs + + REAL(4), INTENT(IN) :: rootMOOP(3) ! Root out of plane bending moments of each blade + REAL(4), INTENT(IN) :: aziAngle ! Rotor azimuth angle + + ! Outputs + + REAL(4), INTENT(OUT) :: axisTilt, axisYaw ! Direct axis and quadrature axis outputted by this transform + + ! Local + + REAL(4), PARAMETER :: phi2 = 2.0/3.0*PI ! Phase difference from first to second blade + REAL(4), PARAMETER :: phi3 = 4.0/3.0*PI ! Phase difference from first to third blade + + ! Body + + axisTilt = 2.0/3.0 * (cos(aziAngle)*rootMOOP(1) + cos(aziAngle+phi2)*rootMOOP(2) + cos(aziAngle+phi3)*rootMOOP(3)) + axisYaw = 2.0/3.0 * (sin(aziAngle)*rootMOOP(1) + sin(aziAngle+phi2)*rootMOOP(2) + sin(aziAngle+phi3)*rootMOOP(3)) + + END SUBROUTINE ColemanTransform + !------------------------------------------------------------------------------------------------------------------------------- + !The inverse Coleman or d-q axis transformation transforms the direct axis and quadrature axis + !back to root out of plane bending moments of each turbine blade + SUBROUTINE ColemanTransformInverse(axisTilt, axisYaw, aziAngle, phi, PitComIPC) + !............................................................................................................................... + + IMPLICIT NONE + + ! Inputs + + REAL(4), INTENT(IN) :: axisTilt, axisYaw ! Direct axis and quadrature axis + REAL(4), INTENT(IN) :: aziAngle ! Rotor azimuth angle + REAL(4), INTENT(IN) :: phi ! Phase shift added to the azimuth angle + + ! Outputs + + REAL(4), INTENT(OUT) :: PitComIPC (3) ! Root out of plane bending moments of each blade + + ! Local + + REAL(4), PARAMETER :: phi2 = 2.0/3.0*PI ! Phase difference from first to second blade + REAL(4), PARAMETER :: phi3 = 4.0/3.0*PI ! Phase difference from first to third blade + + ! Body + + PitComIPC(1) = cos(aziAngle+phi)*axisTilt + sin(aziAngle+phi)*axisYaw + PitComIPC(2) = cos(aziAngle+phi+phi2)*axisTilt + sin(aziAngle+phi+phi2)*axisYaw + PitComIPC(3) = cos(aziAngle+phi+phi3)*axisTilt + sin(aziAngle+phi+phi3)*axisYaw + + END SUBROUTINE ColemanTransformInverse + !------------------------------------------------------------------------------------------------------------------------------- +END MODULE FunctionToolbox \ No newline at end of file diff --git a/Source/IPC.f90 b/Source/IPC.f90 deleted file mode 100644 index b0cb5aa0..00000000 --- a/Source/IPC.f90 +++ /dev/null @@ -1,183 +0,0 @@ -!------------------------------------------------------------------------------------------------------------------------------- -! Individual pitch control subroutine -!SUBROUTINE IPC(rootMOOP, aziAngle, phi, Y_MErr, DT, KInter, Y_IPC_KP, Y_IPC_KI, omegaHP, omegaLP, omegaNotch, zetaHP, zetaLP, zetaNotch, iStatus, IPC_ControlMode, Y_ControlMode, NumBl, PitComIPCF, objInst) -SUBROUTINE IPC(LocalVar, CntrPar, objInst) -!............................................................................................................................... - - USE :: FunctionToolbox - USE :: Filters - USE DRC_Types, ONLY : LocalVariables, ControlParameters, ObjectInstances - - IMPLICIT NONE - - !------------------------------------------------------------------------------------------------------------------------------ - ! Variable declaration and initialization - !------------------------------------------------------------------------------------------------------------------------------ - - ! Inputs - - !REAL(4), INTENT(IN) :: LocalVar%Azimuth ! Rotor azimuth angle - !REAL(4), INTENT(IN) :: LocalVar%DT ! Time step - !REAL(4), INTENT(IN) :: CntrPar%IPC_KI ! Gain for the integrator - !REAL(4), INTENT(IN) :: CntrPar%IPC_omegaHP ! High-pass filter cut-in frequency - !REAL(4), INTENT(IN) :: CntrPar%IPC_omegaLP ! Low-pass filter cut-off frequency - !REAL(4), INTENT(IN) :: CntrPar%IPC_omegaNotch ! Notch filter frequency - !REAL(4), INTENT(IN) :: CntrPar%IPC_phi ! Phase offset added to the azimuth angle - !REAL(4), INTENT(IN) :: LocalVariables%rootMOOP(3) ! Root out of plane bending moments of each blade - !REAL(4), INTENT(IN) :: LocalVariables%Y_MErr ! Yaw alignment error, measured [rad] - !REAL(4), INTENT(IN) :: CntrPar%IPC_zetaHP ! High-pass filter damping value - !REAL(4), INTENT(IN) :: CntrPar%IPC_zetaLP ! Low-pass filter damping value - !REAL(4), INTENT(IN) :: CntrPar%IPC_zetaNotch ! Notch filter damping value - !INTEGER(4), INTENT(IN) :: LocalVariables%iStatus ! A status flag set by the simulation as follows: 0 if this is the first call, 1 for all subsequent time steps, -1 if this is the final call at the end of the simulation. - !INTEGER(4), INTENT(IN) :: LocalVariables%NumBl ! Number of turbine blades - !INTEGER(4), INTENT(IN) :: CntrPar%Y_ControlMode ! Yaw control mode - !REAL(4), INTENT(IN) :: CntrPar%Y_IPC_KP ! Yaw-by-IPC proportional controller gain Kp - !REAL(4), INTENT(IN) :: CntrPar%Y_IPC_KI ! Yaw-by-IPC integral controller gain Ki - !INTEGER(4), INTENT(IN) :: CntrPar%IPC_ControlMode ! Turn Individual Pitch Control (IPC) for fatigue load reductions (pitch contribution) on = 1/off = 0 - - ! Outputs - - !REAL(4), INTENT(OUT) :: LocalVar%IPC_PitComF(3) ! Filtered pitch angle of each rotor blade - - ! Inputs/outputs - - TYPE(LocalVariables), INTENT(INOUT) :: LocalVar - TYPE(ControlParameters), INTENT(IN) :: CntrPar - TYPE(ObjectInstances), INTENT(INOUT) :: objInst - - ! Local variables - - REAL(4), PARAMETER :: PI = 3.14159265359 ! Mathematical constant pi - REAL(4) :: PitComIPC(3) - INTEGER(4) :: K ! Integer used to loop through turbine blades - REAL(4) :: axisTilt, axisYaw, axisYawF ! Direct axis and quadrature axis outputted by Coleman transform - REAL(4), SAVE :: IntAxisTilt, IntAxisYaw ! Integral of the direct axis and quadrature axis - REAL(4) :: IntAxisYawIPC ! IPC contribution with yaw-by-IPC component - REAL(4) :: Y_MErrF, Y_MErrF_IPC ! Unfiltered and filtered yaw alignment error [rad] - REAL(4) :: PitComIPC_woYaw(3) - - !------------------------------------------------------------------------------------------------------------------------------ - ! Body - !------------------------------------------------------------------------------------------------------------------------------ - ! Calculates the commanded pitch angles. - ! NOTE: if it is required for this subroutine to be used multiple times (for 1p and 2p IPC for example), the saved variables - ! IntAxisTilt and IntAxisYaw need to be modified so that they support multiple instances (see LPFilter in the Filters module). - !------------------------------------------------------------------------------------------------------------------------------ - ! Filter rootMOOPs with notch filter - - !DO K = 1,LocalVar%NumBl - ! Instances 1-3 of the Notch Filter are reserved for this routine. - ! rootMOOPF(K) = LocalVar%rootMOOP(K) ! Notch filter currently not in use - !END DO - - ! Initialization - ! Set integrals to be 0 in the first time step - - IF(LocalVar%iStatus==0) THEN - IntAxisTilt = 0.0 - IntAxisYaw = 0.0 - END IF - - ! Body - ! Pass rootMOOPs through the Coleman transform to get the direct and quadrature axis - - CALL ColemanTransform(LocalVar%rootMOOP, LocalVar%Azimuth, axisTilt, axisYaw) - - ! High-pass filter the MBC yaw component and filter yaw alignment error, and compute the yaw-by-IPC contribution - - IF (CntrPar%Y_ControlMode == 2) THEN - axisYawF = HPFilter(axisYaw, LocalVar%DT, CntrPar%IPC_omegaHP, LocalVar%iStatus, .FALSE., objInst%instHPF) - Y_MErrF = SecLPFilter(LocalVar%Y_MErr, LocalVar%DT, CntrPar%IPC_omegaLP, CntrPar%IPC_zetaLP, LocalVar%iStatus, .FALSE., objInst%instSecLPF) - Y_MErrF_IPC = PIController(Y_MErrF, CntrPar%Y_IPC_KP(1), CntrPar%Y_IPC_KI(1), -100.0, 100.0, LocalVar%DT, 0.0, .FALSE., 3) - ELSE - axisYawF = axisYaw - Y_MErrF = 0.0 - END IF - - ! Integrate the signal and multiply with the IPC gain - IF (CntrPar%IPC_ControlMode == 1) THEN - IntAxisTilt = IntAxisTilt + LocalVar%DT * CntrPar%IPC_KI * axisTilt - IntAxisYaw = IntAxisYaw + LocalVar%DT * CntrPar%IPC_KI * axisYawF - ELSE - IntAxisTilt = 0.0 - IntAxisYaw = 0.0 - END IF - - ! Add the yaw-by-IPC contribution - - IntAxisYawIPC = IntAxisYaw + Y_MErrF_IPC - - ! Pass direct and quadrature axis through the inverse Coleman transform to get the commanded pitch angles - - CALL ColemanTransformInverse(IntAxisTilt, IntAxisYawIPC, LocalVar%Azimuth, CntrPar%IPC_phi, PitComIPC) - - ! Filter PitComIPC with second order low pass filter - - DO K = 1,LocalVar%NumBl - ! Instances 1-3 of the Second order Low-Pass Filter are reserved for this routine. - ! LocalVar%IPC_PitComF(K) = SecLPFilter(PitComIPC(K), LocalVar%DT, CntrPar%IPC_omegaLP, CntrPar%IPC_zetaLP, LocalVar%iStatus, K) - LocalVar%IPC_PitComF(K) = PitComIPC(K) - END DO - -CONTAINS - - !------------------------------------------------------------------------------------------------------------------------------- - !The Coleman or d-q axis transformation transforms the root out of plane bending moments of each turbine blade - !to a direct axis and a quadrature axis - SUBROUTINE ColemanTransform(rootMOOP, aziAngle, axisTilt, axisYaw) - !............................................................................................................................... - - IMPLICIT NONE - - ! Inputs - - REAL(4), INTENT(IN) :: rootMOOP(3) ! Root out of plane bending moments of each blade - REAL(4), INTENT(IN) :: aziAngle ! Rotor azimuth angle - - ! Outputs - - REAL(4), INTENT(OUT) :: axisTilt, axisYaw ! Direct axis and quadrature axis outputted by this transform - - ! Local - - REAL(4), PARAMETER :: phi2 = 2.0/3.0*PI ! Phase difference from first to second blade - REAL(4), PARAMETER :: phi3 = 4.0/3.0*PI ! Phase difference from first to third blade - - ! Body - - axisTilt = 2.0/3.0 * (cos(aziAngle)*rootMOOP(1) + cos(aziAngle+phi2)*rootMOOP(2) + cos(aziAngle+phi3)*rootMOOP(3)) - axisYaw = 2.0/3.0 * (sin(aziAngle)*rootMOOP(1) + sin(aziAngle+phi2)*rootMOOP(2) + sin(aziAngle+phi3)*rootMOOP(3)) - - END SUBROUTINE ColemanTransform - !------------------------------------------------------------------------------------------------------------------------------- - !The inverse Coleman or d-q axis transformation transforms the direct axis and quadrature axis - !back to root out of plane bending moments of each turbine blade - SUBROUTINE ColemanTransformInverse(axisTilt, axisYaw, aziAngle, phi, PitComIPC) - !............................................................................................................................... - - IMPLICIT NONE - - ! Inputs - - REAL(4), INTENT(IN) :: axisTilt, axisYaw ! Direct axis and quadrature axis - REAL(4), INTENT(IN) :: aziAngle ! Rotor azimuth angle - REAL(4), INTENT(IN) :: phi ! Phase shift added to the azimuth angle - - ! Outputs - - REAL(4), INTENT(OUT) :: PitComIPC (3) ! Root out of plane bending moments of each blade - - ! Local - - REAL(4), PARAMETER :: phi2 = 2.0/3.0*PI ! Phase difference from first to second blade - REAL(4), PARAMETER :: phi3 = 4.0/3.0*PI ! Phase difference from first to third blade - - ! Body - - PitComIPC(1) = cos(aziAngle+phi)*axisTilt + sin(aziAngle+phi)*axisYaw - PitComIPC(2) = cos(aziAngle+phi+phi2)*axisTilt + sin(aziAngle+phi+phi2)*axisYaw - PitComIPC(3) = cos(aziAngle+phi+phi3)*axisTilt + sin(aziAngle+phi+phi3)*axisYaw - - END SUBROUTINE ColemanTransformInverse - !------------------------------------------------------------------------------------------------------------------------------- -END SUBROUTINE IPC diff --git a/Source/Obj_win32/filters.mod b/Source/Obj_win32/filters.mod deleted file mode 100644 index 3b84261c572bf48a26ec1502afbdabe845b1bd20..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 759 zcmVOw=bO)stDu8D#L z?C&2yj7CqvY`vUg?SJ)yAO5F>znc52 z*WtG?m|h3u-#bj0m7;Q@d|by8r5smaZ|~-qa2PM=ilGIjkekoe;SEnkS>ep))1}IS zhqXUlEdsxS3g?q0#6O~PP~a3NMiy*1aN&Vg6>y@UQ3~Vu5p0wdaFl}zPz%%sbwJ(n z!b&I=(Ss2u`Szt3B`GM;m6DZLlua5u@IY(^;R~4_Sq?@V( zW1FOLl4&zru|q0Leh$K(A^ws6zvq*`0hB0p3^9TIZl%!pv!o8H-^CS@0>GIi3nGKgC!=%5lD) zRAplnNVrUcWW}ccZ1KPv4Et8w@*aosCyf5R@i2Nczvmb8|M~pVGA3~u2PS}3Gi-Lq zUjwoVZ=3+EE-5wE0M-Q70#1Q5;2gLBUIkvu!D}iwuLDWYL{Yd1)0z*`_F z2nK@Z*;`jhg#fbl(~YDpDMe1PwLfl?S9DPWQ3uh;?nFX`w+bXkk~qo8@W4s^o8el< zXY`e|x6O-Qg3fp_c*sl7`IV!JhZJE<#rgC~)mOamf-2@wIxAhZfps@6U9~BT)%rtp zWpegRcQ0LxUs1<0M7C2u^$UH*-U*;}tbu8dWWMiF-bylA0gcgUJTN7YA{17Gj-0tT0&rODZw=d?EAIwYnFCjC+`Q zSiUYkNy#EyKKaomxL95L>l{5+XQ=LPZ>7%F8e8I;V|fVa4Rq1YNw5yC4mq<`-|1Ox zMZFlm22eYEyj-P!!Se-oVz>aab@+4XM+F4#(#py=^!Ot1*SFFaxz~5w%(BIEtKUfc z!0GSSz1VRgOZe3$J*aTF!7eNA`e}79QN3|m-OE##FE5T>U0~OEroE!jQd@K?8{os1 dSEI7GRJ%zkTR5d@!{$y@v}YbRtuE6M001Y5aohj^ diff --git a/Source/ReadParameters.f90 b/Source/ReadParameters.f90 deleted file mode 100644 index abb579d5..00000000 --- a/Source/ReadParameters.f90 +++ /dev/null @@ -1,127 +0,0 @@ -MODULE ReadParameters - - USE DRC_Types - IMPLICIT NONE - - TYPE(ControlParameters), SAVE :: CntrPar - TYPE(LocalVariables) :: LocalVar - -CONTAINS - !.............................................................................................................................. - ! Read all constant control parameters from ControllerParameters.in parameter file - !.............................................................................................................................. - SUBROUTINE ReadControlParameterFileSub() - - INTEGER(4), PARAMETER :: UnControllerParameters = 89 - OPEN(unit=UnControllerParameters, file='ControllerParameters.in', status='old', action='read') - - !------------------- GENERAL CONSTANTS ------------------- - READ(UnControllerParameters, *) CntrPar%CornerFreq - READ(UnControllerParameters, *) CntrPar%LoggingLevel - - !------------------- IPC CONSTANTS ----------------------- - READ(UnControllerParameters, *) CntrPar%IPC_KI - READ(UnControllerParameters, *) CntrPar%IPC_ControlMode - READ(UnControllerParameters, *) CntrPar%IPC_omegaHP - READ(UnControllerParameters, *) CntrPar%IPC_omegaLP - READ(UnControllerParameters, *) CntrPar%IPC_omegaNotch - READ(UnControllerParameters, *) CntrPar%IPC_phi - READ(UnControllerParameters, *) CntrPar%IPC_zetaHP - READ(UnControllerParameters, *) CntrPar%IPC_zetaLP - READ(UnControllerParameters, *) CntrPar%IPC_zetaNotch - - !------------------- PITCH CONSTANTS ----------------------- - READ(UnControllerParameters, *) CntrPar%PC_GS_n - - ALLOCATE(CntrPar%PC_GS_angles(CntrPar%PC_GS_n)) - READ(UnControllerParameters,*) CntrPar%PC_GS_angles - - ALLOCATE(CntrPar%PC_GS_kp(CntrPar%PC_GS_n)) - READ(UnControllerParameters,*) CntrPar%PC_GS_kp - - ALLOCATE(CntrPar%PC_GS_ki(CntrPar%PC_GS_n)) - READ(UnControllerParameters,*) CntrPar%PC_GS_ki - - ALLOCATE(CntrPar%PC_GS_kd(CntrPar%PC_GS_n)) - READ(UnControllerParameters,*) CntrPar%PC_GS_kd - - ALLOCATE(CntrPar%PC_GS_tf(CntrPar%PC_GS_n)) - READ(UnControllerParameters,*) CntrPar%PC_GS_tf - - READ(UnControllerParameters, *) CntrPar%PC_MaxPit - READ(UnControllerParameters, *) CntrPar%PC_MinPit - READ(UnControllerParameters, *) CntrPar%PC_MaxRat - READ(UnControllerParameters, *) CntrPar%PC_MinRat - READ(UnControllerParameters, *) CntrPar%PC_RefSpd - READ(UnControllerParameters, *) CntrPar%PC_SetPnt - READ(UnControllerParameters, *) CntrPar%PC_Switch - - !------------------- TORQUE CONSTANTS ----------------------- - READ(UnControllerParameters, *) CntrPar%VS_ControlMode - READ(UnControllerParameters, *) CntrPar%VS_CtInSp - READ(UnControllerParameters, *) CntrPar%VS_GenTrqArSatMax - READ(UnControllerParameters, *) CntrPar%VS_MaxOM - READ(UnControllerParameters, *) CntrPar%VS_MaxRat - READ(UnControllerParameters, *) CntrPar%VS_MaxTq - READ(UnControllerParameters, *) CntrPar%VS_MinTq - READ(UnControllerParameters, *) CntrPar%VS_MinOM - READ(UnControllerParameters, *) CntrPar%VS_Rgn2K - READ(UnControllerParameters, *) CntrPar%VS_RtPwr - READ(UnControllerParameters, *) CntrPar%VS_RtTq - READ(UnControllerParameters, *) CntrPar%VS_RtSpd - READ(UnControllerParameters, *) CntrPar%VS_n - - ALLOCATE(CntrPar%VS_KP(CntrPar%VS_n)) - READ(UnControllerParameters,*) CntrPar%VS_KP - - ALLOCATE(CntrPar%VS_KI(CntrPar%VS_n)) - READ(UnControllerParameters,*) CntrPar%VS_KI - - !------------------- YAW CONSTANTS ----------------------- - READ(UnControllerParameters, *) CntrPar%Y_ControlMode - READ(UnControllerParameters, *) CntrPar%Y_ErrThresh - READ(UnControllerParameters, *) CntrPar%Y_IPC_n - - ALLOCATE(CntrPar%Y_IPC_KP(CntrPar%Y_IPC_n)) - READ(UnControllerParameters,*) CntrPar%Y_IPC_KP - - ALLOCATE(CntrPar%Y_IPC_KI(CntrPar%Y_IPC_n)) - READ(UnControllerParameters,*) CntrPar%Y_IPC_KI - - READ(UnControllerParameters, *) CntrPar%Y_MErrSet - READ(UnControllerParameters, *) CntrPar%Y_omegaLPFast - READ(UnControllerParameters, *) CntrPar%Y_omegaLPSlow - READ(UnControllerParameters, *) CntrPar%Y_Rate - - !------------------- CALCULATED CONSTANTS ----------------------- - CntrPar%PC_RtTq99 = CntrPar%VS_RtTq*0.99 - CntrPar%VS_Rgn2MinTq = CntrPar%VS_Rgn2K*CntrPar%VS_MinOM**2 - CntrPar%VS_Rgn2MaxTq = CntrPar%VS_Rgn2K*CntrPar%VS_MaxOM**2 - CntrPar%VS_Rgn3MP = CntrPar%PC_SetPnt + CntrPar%PC_Switch - - CLOSE(UnControllerParameters) - END SUBROUTINE ReadControlParameterFileSub - - SUBROUTINE ReadAvrSWAP(avrSWAP) - - USE, INTRINSIC :: ISO_C_Binding - REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. - - ! Load variables from calling program (See Appendix A of Bladed User's Guide): - LocalVar%iStatus = NINT(avrSWAP(1)) - LocalVar%Time = avrSWAP(2) - LocalVar%DT = avrSWAP(3) - LocalVar%BlPitch(1) = avrSWAP(4) - LocalVar%VS_GenPwr = avrSWAP(15) - LocalVar%GenSpeed = avrSWAP(20) - LocalVar%Y_M = avrSWAP(24) - LocalVar%HorWindV = avrSWAP(27) - LocalVar%rootMOOP(1) = avrSWAP(30) - LocalVar%rootMOOP(2) = avrSWAP(31) - LocalVar%rootMOOP(3) = avrSWAP(32) - LocalVar%BlPitch(2) = avrSWAP(33) - LocalVar%BlPitch(3) = avrSWAP(34) - LocalVar%Azimuth = avrSWAP(60) - LocalVar%NumBl = NINT(avrSWAP(61)) - END SUBROUTINE ReadAvrSWAP -END MODULE ReadParameters \ No newline at end of file diff --git a/Source/ReadSetParameters.f90 b/Source/ReadSetParameters.f90 new file mode 100644 index 00000000..fed826af --- /dev/null +++ b/Source/ReadSetParameters.f90 @@ -0,0 +1,340 @@ +MODULE ReadSetParameters + + USE, INTRINSIC :: ISO_C_Binding + IMPLICIT NONE + +CONTAINS + !.............................................................................................................................. + ! Read all constant control parameters from ControllerParameters.in parameter file + !.............................................................................................................................. + SUBROUTINE ReadControlParameterFileSub(CntrPar) + USE DRC_Types, ONLY : ControlParameters + + INTEGER(4), PARAMETER :: UnControllerParameters = 89 + TYPE(ControlParameters), INTENT(INOUT) :: CntrPar + + OPEN(unit=UnControllerParameters, file='ControllerParameters.in', status='old', action='read') + + !------------------- GENERAL CONSTANTS ------------------- + READ(UnControllerParameters, *) CntrPar%CornerFreq + READ(UnControllerParameters, *) CntrPar%LoggingLevel + + !------------------- IPC CONSTANTS ----------------------- + READ(UnControllerParameters, *) CntrPar%IPC_KI + READ(UnControllerParameters, *) CntrPar%IPC_ControlMode + READ(UnControllerParameters, *) CntrPar%IPC_omegaHP + READ(UnControllerParameters, *) CntrPar%IPC_omegaLP + READ(UnControllerParameters, *) CntrPar%IPC_omegaNotch + READ(UnControllerParameters, *) CntrPar%IPC_phi + READ(UnControllerParameters, *) CntrPar%IPC_zetaHP + READ(UnControllerParameters, *) CntrPar%IPC_zetaLP + READ(UnControllerParameters, *) CntrPar%IPC_zetaNotch + + !------------------- PITCH CONSTANTS ----------------------- + READ(UnControllerParameters, *) CntrPar%PC_GS_n + + ALLOCATE(CntrPar%PC_GS_angles(CntrPar%PC_GS_n)) + READ(UnControllerParameters,*) CntrPar%PC_GS_angles + + ALLOCATE(CntrPar%PC_GS_kp(CntrPar%PC_GS_n)) + READ(UnControllerParameters,*) CntrPar%PC_GS_kp + + ALLOCATE(CntrPar%PC_GS_ki(CntrPar%PC_GS_n)) + READ(UnControllerParameters,*) CntrPar%PC_GS_ki + + ALLOCATE(CntrPar%PC_GS_kd(CntrPar%PC_GS_n)) + READ(UnControllerParameters,*) CntrPar%PC_GS_kd + + ALLOCATE(CntrPar%PC_GS_tf(CntrPar%PC_GS_n)) + READ(UnControllerParameters,*) CntrPar%PC_GS_tf + + READ(UnControllerParameters, *) CntrPar%PC_MaxPit + READ(UnControllerParameters, *) CntrPar%PC_MinPit + READ(UnControllerParameters, *) CntrPar%PC_MaxRat + READ(UnControllerParameters, *) CntrPar%PC_MinRat + READ(UnControllerParameters, *) CntrPar%PC_RefSpd + READ(UnControllerParameters, *) CntrPar%PC_SetPnt + READ(UnControllerParameters, *) CntrPar%PC_Switch + + !------------------- TORQUE CONSTANTS ----------------------- + READ(UnControllerParameters, *) CntrPar%VS_ControlMode + READ(UnControllerParameters, *) CntrPar%VS_CtInSp + READ(UnControllerParameters, *) CntrPar%VS_GenTrqArSatMax + READ(UnControllerParameters, *) CntrPar%VS_MaxOM + READ(UnControllerParameters, *) CntrPar%VS_MaxRat + READ(UnControllerParameters, *) CntrPar%VS_MaxTq + READ(UnControllerParameters, *) CntrPar%VS_MinTq + READ(UnControllerParameters, *) CntrPar%VS_MinOM + READ(UnControllerParameters, *) CntrPar%VS_Rgn2K + READ(UnControllerParameters, *) CntrPar%VS_RtPwr + READ(UnControllerParameters, *) CntrPar%VS_RtTq + READ(UnControllerParameters, *) CntrPar%VS_RtSpd + READ(UnControllerParameters, *) CntrPar%VS_n + + ALLOCATE(CntrPar%VS_KP(CntrPar%VS_n)) + READ(UnControllerParameters,*) CntrPar%VS_KP + + ALLOCATE(CntrPar%VS_KI(CntrPar%VS_n)) + READ(UnControllerParameters,*) CntrPar%VS_KI + + !------------------- YAW CONSTANTS ----------------------- + READ(UnControllerParameters, *) CntrPar%Y_ControlMode + READ(UnControllerParameters, *) CntrPar%Y_ErrThresh + READ(UnControllerParameters, *) CntrPar%Y_IPC_n + + ALLOCATE(CntrPar%Y_IPC_KP(CntrPar%Y_IPC_n)) + READ(UnControllerParameters,*) CntrPar%Y_IPC_KP + + ALLOCATE(CntrPar%Y_IPC_KI(CntrPar%Y_IPC_n)) + READ(UnControllerParameters,*) CntrPar%Y_IPC_KI + + READ(UnControllerParameters, *) CntrPar%Y_MErrSet + READ(UnControllerParameters, *) CntrPar%Y_omegaLPFast + READ(UnControllerParameters, *) CntrPar%Y_omegaLPSlow + READ(UnControllerParameters, *) CntrPar%Y_Rate + + !------------------- CALCULATED CONSTANTS ----------------------- + CntrPar%PC_RtTq99 = CntrPar%VS_RtTq*0.99 + CntrPar%VS_Rgn2MinTq = CntrPar%VS_Rgn2K*CntrPar%VS_MinOM**2 + CntrPar%VS_Rgn2MaxTq = CntrPar%VS_Rgn2K*CntrPar%VS_MaxOM**2 + CntrPar%VS_Rgn3MP = CntrPar%PC_SetPnt + CntrPar%PC_Switch + + CLOSE(UnControllerParameters) + END SUBROUTINE ReadControlParameterFileSub + + SUBROUTINE ReadAvrSWAP(avrSWAP, LocalVar) + USE DRC_Types, ONLY : LocalVariables + + REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + + ! Load variables from calling program (See Appendix A of Bladed User's Guide): + LocalVar%iStatus = NINT(avrSWAP(1)) + LocalVar%Time = avrSWAP(2) + LocalVar%DT = avrSWAP(3) + LocalVar%BlPitch(1) = avrSWAP(4) + LocalVar%VS_GenPwr = avrSWAP(15) + LocalVar%GenSpeed = avrSWAP(20) + LocalVar%Y_M = avrSWAP(24) + LocalVar%HorWindV = avrSWAP(27) + LocalVar%rootMOOP(1) = avrSWAP(30) + LocalVar%rootMOOP(2) = avrSWAP(31) + LocalVar%rootMOOP(3) = avrSWAP(32) + LocalVar%BlPitch(2) = avrSWAP(33) + LocalVar%BlPitch(3) = avrSWAP(34) + LocalVar%Azimuth = avrSWAP(60) + LocalVar%NumBl = NINT(avrSWAP(61)) + END SUBROUTINE ReadAvrSWAP + + SUBROUTINE Assert(LocalVar, CntrPar, avrSWAP, aviFAIL, ErrMsg, size_avcMSG) + + USE, INTRINSIC :: ISO_C_Binding + USE DRC_Types, ONLY : LocalVariables, ControlParameters + + IMPLICIT NONE + + ! Inputs + TYPE(ControlParameters), INTENT(IN) :: CntrPar + TYPE(LocalVariables), INTENT(IN) :: LocalVar + INTEGER(4), INTENT(IN) :: size_avcMSG + REAL(C_FLOAT), INTENT(IN) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + + ! Outputs + INTEGER(C_INT), INTENT(OUT) :: aviFAIL ! A flag used to indicate the success of this DLL call set as follows: 0 if the DLL call was successful, >0 if the DLL call was successful but cMessage should be issued as a warning messsage, <0 if the DLL call was unsuccessful or for any other reason the simulation is to be stopped at this point with cMessage as the error message. + CHARACTER(size_avcMSG-1), INTENT(OUT) :: ErrMsg ! a Fortran version of the C string argument (not considered an array here) [subtract 1 for the C null-character] + + ! Local + + !.............................................................................................................................. + ! Check validity of input parameters: + !.............................................................................................................................. + + IF (CntrPar%CornerFreq <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'CornerFreq must be greater than zero.' + ENDIF + + IF (LocalVar%DT <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'DT must be greater than zero.' + ENDIF + + IF (CntrPar%VS_CtInSp < 0.0) THEN + aviFAIL = -1 + ErrMsg = 'VS_CtInSp must not be negative.' + ENDIF + + IF (CntrPar%VS_MinOM <= CntrPar%VS_CtInSp) THEN + aviFAIL = -1 + ErrMsg = 'VS_MinOM must be greater than VS_CtInSp.' + ENDIF + + IF (CntrPar%VS_MaxRat <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'VS_MaxRat must be greater than zero.' + ENDIF + + IF (CntrPar%VS_RtTq < 0.0) THEN + aviFAIL = -1 + ErrMsg = 'VS_RtTw must not be negative.' + ENDIF + + IF (CntrPar%VS_Rgn2K < 0.0) THEN + aviFAIL = -1 + ErrMsg = 'VS_Rgn2K must not be negative.' + ENDIF + + IF (CntrPar%VS_MaxTq < CntrPar%VS_RtTq) THEN + aviFAIL = -1 + ErrMsg = 'VS_RtTq must not be greater than VS_MaxTq.' + ENDIF + + IF (CntrPar%VS_KP(1) > 0.0) THEN + aviFAIL = -1 + ErrMsg = 'VS_KP must be greater than zero.' + ENDIF + + IF (CntrPar%VS_KI(1) > 0.0) THEN + aviFAIL = -1 + ErrMsg = 'VS_KI must be greater than zero.' + ENDIF + + IF (CntrPar%PC_RefSpd <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'PC_RefSpd must be greater than zero.' + ENDIF + + IF (CntrPar%PC_MaxRat <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'PC_MaxRat must be greater than zero.' + ENDIF + + IF (CntrPar%PC_MinPit >= CntrPar%PC_MaxPit) THEN + aviFAIL = -1 + ErrMsg = 'PC_MinPit must be less than PC_MaxPit.' + ENDIF + + IF (CntrPar%IPC_KI <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'IPC_KI must be greater than zero.' + ENDIF + + IF (CntrPar%IPC_omegaLP <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'IPC_omegaLP must be greater than zero.' + ENDIF + + IF (CntrPar%IPC_omegaNotch <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'IPC_omegaNotch must be greater than zero.' + ENDIF + + IF (CntrPar%IPC_phi <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'IPC_phi must be greater than zero.' + ENDIF + + IF (CntrPar%IPC_zetaLP <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'IPC_zetaLP must be greater than zero.' + ENDIF + + IF (CntrPar%IPC_zetaNotch <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'IPC_zetaNotch must be greater than zero.' + ENDIF + + IF (CntrPar%Y_ErrThresh <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'Y_ErrThresh must be greater than zero.' + ENDIF + + IF (CntrPar%Y_Rate <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'CntrPar%Y_Rate must be greater than zero.' + ENDIF + + IF (CntrPar%Y_omegaLPFast <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'Y_omegaLPFast must be greater than zero.' + ENDIF + + IF (CntrPar%Y_omegaLPSlow <= 0.0) THEN + aviFAIL = -1 + ErrMsg = 'Y_omegaLPSlow must be greater than zero.' + ENDIF + + ! Abort if the user has not requested a pitch angle actuator (See Appendix A + ! of Bladed User's Guide): + IF (NINT(avrSWAP(10)) /= 0) THEN ! .TRUE. if a pitch angle actuator hasn't been requested + aviFAIL = -1 + ErrMsg = 'Pitch angle actuator not requested.' + ENDIF + + END SUBROUTINE Assert + + SUBROUTINE SetParameters(avrSWAP, aviFAIL, ErrMsg, size_avcMSG, CntrPar, LocalVar, objInst) + + USE DRC_Types, ONLY : ControlParameters, LocalVariables, ObjectInstances + + INTEGER(4), INTENT(IN) :: size_avcMSG + TYPE(ControlParameters), INTENT(INOUT) :: CntrPar + TYPE(LocalVariables), INTENT(INOUT) :: LocalVar + TYPE(ObjectInstances), INTENT(INOUT) :: objInst + + REAL(C_FLOAT), INTENT(INOUT) :: avrSWAP(*) ! The swap array, used to pass data to, and receive data from, the DLL controller. + INTEGER(C_INT), INTENT(OUT) :: aviFAIL ! A flag used to indicate the success of this DLL call set as follows: 0 if the DLL call was successful, >0 if the DLL call was successful but cMessage should be issued as a warning messsage, <0 if the DLL call was unsuccessful or for any other reason the simulation is to be stopped at this point with cMessage as the error message. + CHARACTER(size_avcMSG-1), INTENT(OUT) :: ErrMsg ! a Fortran version of the C string argument (not considered an array here) [subtract 1 for the C null-character] + + ! Set aviFAIL to 0 in each iteration: + aviFAIL = 0 + + ! Initialize all filter instance counters at 1 + objInst%instLPF = 1 + objInst%instSecLPF = 1 + objInst%instHPF = 1 + objInst%instNotchSlopes = 1 + objInst%instNotch = 1 + + ! Set unused outputs to zero (See Appendix A of Bladed User's Guide): + avrSWAP(36) = 0.0 ! Shaft brake status: 0=off + avrSWAP(41) = 0.0 ! Demanded yaw actuator torque + avrSWAP(46) = 0.0 ! Demanded pitch rate (Collective pitch) + avrSWAP(65) = 0.0 ! Number of variables returned for logging + avrSWAP(72) = 0.0 ! Generator start-up resistance + avrSWAP(79) = 0.0 ! Request for loads: 0=none + avrSWAP(80) = 0.0 ! Variable slip current status + avrSWAP(81) = 0.0 ! Variable slip current demand + + ! Read any External Controller Parameters specified in the User Interface + ! and initialize variables: + IF (LocalVar%iStatus == 0) THEN ! .TRUE. if we're on the first call to the DLL + + ! Inform users that we are using this user-defined routine: + aviFAIL = 1 + ErrMsg = ' '//NEW_LINE('A')// & + 'Running the Delft Research Controller (DRC) '//NEW_LINE('A')// & + 'A wind turbine controller for use in the scientific field '//NEW_LINE('A')// & + 'Written by S.P. Mulders, Jan-Willem van Wingerden '//NEW_LINE('A')// & + 'Delft University of Technology, The Netherlands '//NEW_LINE('A')// & + 'Visit our GitHub-page to contribute to this project: '//NEW_LINE('A')// & + 'https://github.com/TUDelft-DataDrivenControl ' + + CALL ReadControlParameterFileSub(CntrPar) + + ! Initialize testValue (debugging variable) + LocalVar%TestType = 0 + + ! Initialize the SAVEd variables: + ! NOTE: LocalVar%VS_LastGenTrq, though SAVEd, is initialized in the torque controller + ! below for simplicity, not here. + LocalVar%PitCom = LocalVar%BlPitch ! This will ensure that the variable speed controller picks the correct control region and the pitch controller picks the correct gain on the first call + LocalVar%Y_AccErr = 0.0 ! This will ensure that the accumulated yaw error starts at zero + LocalVar%Y_YawEndT = -1.0 ! This will ensure that the initial yaw end time is lower than the actual time to prevent initial yawing + + ! Check validity of input parameters: + CALL Assert(LocalVar, CntrPar, avrSWAP, aviFAIL, ErrMsg, size_avcMSG) + + ENDIF + END SUBROUTINE SetParameters +END MODULE ReadSetParameters \ No newline at end of file diff --git a/Source/makefile b/Source/makefile index 05f008b8..aa3b2c85 100644 --- a/Source/makefile +++ b/Source/makefile @@ -19,7 +19,7 @@ BITS = 32 DLL_DIR = . -SOURCE_FILE = DRC_Types.f90 ReadParameters.f90 FunctionToolbox.f90 Filters.f90 IPC.f90 DISCON.f90 +SOURCE_FILE = Constants.f90 DRC_Types.f90 FunctionToolbox.f90 Filters.f90 ReadSetParameters.f90 Controllers.f90 DISCON.f90 # Name of compiler to use and flags to use.