Software Package on Integrated Nonlinear Dynamic Modeling and
Contents
1. 150 200 250 0 50 100 150 200 250 300 Time t min Figure 2 6 Electric input power and regenerative power 32 The temperature profiles for various nodes taken from different areas of the motor are shown below The highest temperatures pertain to the nodes closest to the copper windings The ambient temperature was held at 22 C This flight profile showed only a 0 4 C temperature increase over 5 minute period The results are saved in temperature xlsx Temperature Distribution 22 4 22 35 4 223 E pe e il TT B m A E O 2225 gt LA ue eg 3 222 y LL m y pu p 5 E e O M Pus a gt di c Jf Ly A a i 5 2215 J gt p Pa a p j y eU gt d hi gt Pu a 22 1 j P d Po p pa e e lupe AT di di gt E ii p E uu e e 22 05 dd g y coy 7 m 7 n Nd Y wat E uu gi eee ee adi AD a S 22 Time t min Figure 2 7 Temperature profile of the motor 33 Part Ill Source Code foc m Nonlinear Dynamic Modeling and Field Oriented Control of Permanent Magnet PM Motor o9 o9 o Written in SI or MKS Unit System Ao oe Authors David Woodburn Dr Lei Zhou Dr Thomas X Wu o oe o o clear all close all clc motorParameter o Load EMA mission profile Mission xlsread missionProfile xlsx Load i2. plot tN 60 tauLN ylabel Load Torque it tau_L rm N m xlabel Time it t xrm min title Load Torque figure Name Magnetic Torque plot ENR 60 tauMN g 5 xlsSwrite tauM xlsx tNR tauMN 3 ylabel Magnetic Torque it tau_M rm N m xlabel Time Yitit ikes min title Magnetic Torque Current Tlyurc Heme current plot tNR 60 idN tNR 60 iqgN xlswrite idlqg xlsx ENR 1idN igNl ylabel Currents it i_d rm A and it i_q rm A xlabeli Time XititlXrm min title dq Currents legend XiEb L dj Xacrix gi Voltage figure Name Voltage plot tNR 60 udN tNR 60 uqN xlswrite uduq xlsx tNR udN uqN v Xlawrite udcdug xlsx tHR UdH UgdH for high resolution very slow ylabel Voltages it u_d rm V and it u_q rm V xlabel Time Xitit irm min title dq Voltages legend iibim HI Mich apk Copper Loss figure Name Copper Loss plot ENR 60 PcuN 1000 xlswrite Ponu xlsk tNR PCUN ylabel Power Loss it P_ cu rm kW xlabel Time Xit it iiem min title Power Loss in Windings 43 Oo Electric Power rigure Name Eleoctrio POwer xlsSwrite Pin xlsX tNB PainN l j PinNN ENN qs nO 1 hold on for n l length tNR 1 PinNN PinNN PinN n tNN tNN tNR n if PinN n PinN n 1 lt O PinNN PinNN 0 tadd t
3. Tth 13 Rth 14 Ith 8 TEh l13 7yREB I0 7CEh 13 dtTbhs 1 40 id idP ig GP dir de d dP die dig de digP dt Ld LdP Lo LGF theta me theta meP omega me omega meP alpha me alpha meP CHR ntR t udH ntR ud uqH ntR ug theta_meR ntR theta_me omega meR ntR omega me tauLR ntR tauLP nck rs Xs ER ER end end hi res simulation Store time ENR nt t Store motor stroke theta meN nt theta me 5 rad omega meN nt omega me rad s Oo 5 Store currents idN nt id A iqN nt iq A Oo Store inductances LdN nt Ld H LqN nt Lg H Oo Store voltages udN nt ud V ugN nt uq 5 V Oo Store resistance RsN nt Rs ohm Oo oLOre IOLOF Lorque tauMN nt tauMP N m Oo Store powers W 4 PinN nt 3 2 iq uq id ud 3 Input power PcuN nt 3 2 iq 2 id 2 Rs 3 Copper loss for whole motor PfricN nt tauFricP omega me 2 p Oo SLOre LOmperaLures LCAN ne e TEO Oo end End time stepping Unpad profile arrays EN CNA NEL theta meSN theta meSN 2 Nt 1 tauLN tauLN 2 Nt 1 tR ntR end theta meR ntR end omega meR ntR end tauLR ntR end disp Dynamic Modeling and Field Oriented Control of PM Motor Calculate stroke speed and acceleration xN t
4. th2Lh 2 Midpoint between theta me limits ER R Th2R EhO i lerodt tho th3R CERAR En0 Llercseat tho Chill Ch thzLh i ilrSat tho Chol tho th2L 1 roat FP THU Oo Set threshold to show limits in plots limThreshold 0 8 35 o NERO round NtR 1 05 ER zeros 1 NtRO uaR zeros 1 NtRO ubR zeros 1 NtRO0 tauLR zeros 1 NtRO theta meR zeros 1 NtR0 omega meR zeros 1 NtR0 o Initialize hi res arrays o Add 5 buffer Hi res time record s Hi res phase A voltage record V Hi res phase B voltage record V Hi res load torque record N m Hi res theta_me record rad Hi res omega_me record rad s oO o9 oO Ae ae o Initialize low res arrays tNR zeros 1 Nt theta meN zeros 1 Nt o Actual theta me rad theta meN 1 theta meSN 1 omega meN zeros 1 Nt idN zeros 1 Nt iqN zeros 1 Nt phiN zeros 1 Nt o Actual omega me rad s Direct current A Quadrature current A Torque angle rad o oO oP LAN zeros 1 Nt LdN 1 Ld0 5 Direct inductance H LqN zeros 1 Nt LqN 1 Lq0 5 Quadrature inductance H udN zeros 1 Nt Direct voltage V ugN zeros 1 Nt Quadrature voltage V RsN zeros 1 Nt RsN 1 Rs Resistance ohm tauMN zeros 1 Nt Machine torque N m PinN zeros 1 Nt 5 Power in W PcuN zeros 1 Nt Copper loss W Pfri
5. 2 929048419 2400420585304 248439321887 2905950657 2 771095598 2199172065535 2 1106524174 2 090325 02901l 2 659 1174 ZO OYU Z2blZISUISS Z 0921424T1721 74 2572607 460 26005902515 2 9360498402 2 519900421 0 001 gt H GSN 4612083399397 4459269 655 44 276907416 3919619295 349109909594 2209915783 2 1403690271 2445203741 XL 999 al 199252454 1 5023051455 140045900442 1 9608691069 1 4707298535 1386369063 1 31L3170018 1249202936 1 1924251298 1 141458078 1 0958965961 1 0546L090985 1 01726986601 0 99509395 0 951930145 0 923114329 0 896516968 0 971994934 0 584926359 025009849653 0808445183 0789983231 0001 H Define thermal parameters Rtn 9000153995 0 007697 Oso 9 001996 16 07 042085717 2 155097 10 977 2 159914 a A 5 0337 37 40 2267 2615 gt Thermal Resistance o C W Crh L39e09r UF 92525 Lovells OF 596 055 l4 9099 T 095 6 725 OI0DITOsX OF 5 15 DL 52 56 3 Thermal Capacitance J o C Nth length Rth Troom 22 Room Temperature o C Tth ones length Rth 1 Troom IS 0 stator loss IW 0 winding loss IBB 0 rear bearing loss IBF front bearing loss IM 0 magnet loss IR7a 0 IR14 0 NL length iL iLMin min iL iLMax max iL o Nex Get current min and mex o Total Coupling Ratio Define gear train coupling ratio 49 8728 0 0254 rad m 45 errDrive 0 75
6. LqsNi nbq wLa LasNiniaqrL gt TH oo The above is equivalent to using interpl but much faster LdP interpl iL LdSN abs iqP LAF anterpl 1L L9SN abs 192 oo oo Get next magnetic and air torques p 2 is used to get cummulative torque not to convert to electrical units tauMP 3 p 4 igP lambdaPM idP LdP LqP N m tauFricP c omega meP 2 p N m Oo Get next angular acceleration alpha meP p 2 I tauMP tauLP tauFricP 4 DIrad s 2 o Comes from tauM tauL I alpha m alpha me alpha m p 2 Get next angular speed and angle omega meP omega me alpha dt alpha me beta dt alpha meP dtE rad s theta meP theta me alpha dt omega me beta dt omega meP dtE rad theta meP is updated here because some parameters might depend on the value of theta me oo Get next currents Current should be updated last since it is used for the convergence test didP dt 1 LdP Rs idP omega meP LgP igP ud A s digP dt 1 LqP Rs iqP omega meP LdP idP uq omega meP lambdaPmM A78 ide 1d alpha dt did dt beta dt oidP dt dtE A iqPOld iqP igP ig aloha de gt di dL beta di digi dt cil Bi A Check for electrical convergence if abs iqP iqPOld iLim lt epsilon break end end end nConverge 39 Oo Add convergence iterations NConverge NConverge nGConverge Update heat gene
7. angle multiply p 2 rad Given when stroke is O Table 2 4 Simulation Limits These parameters are given motor torque limit voltage limit m s 17 Excel file missionProfile xlsx provides data for mission profile The mission profile data is assumed to have been provided already and is not subject for discussion in this manual The input parameters in mission profile is summarized in Table 5 Table 2 5 Mission Profile Parameters may come from aerodynamic measurement or simulation Mission Profile Explanations Units Parameters xSN array for stroke profile FLN array for load force profile Due to the fact that there can be various mission profiles and motor designs it 1s important to ensure that your file names in the FOC folder match the file names as represented in the code above There are also input parameters used to set up simulation These parameters are summarized in Table 2 6 Table 2 6 Simulation Parameters These parameters are given spe samples per cycle of highest frequency tauTh thermal time constant s e __ proportional parameter for PI contol d __ dt forward time stepping weight for implicit algorithm in electrical loop MEN epsilon error for electrical convergence loop Ncvg maximum number of convergence iterations If the experience of the user 1s adequate in the design of the control scheme there are few parameters that can be changed to decrease error in the
8. dt R Rig Wo U R RQU R Roly t 20 12 0 21 8 R R dt U 3 U Oz U dU C IT t 27 A E te yg chi Q2 Buc UM curis qe 26 qu pg 23 Ra KR R R dt U Uy U U dU e e 24 E xz cw WU Q4 U U U U dU A e TO E ut 25 R R dt U U U U U U dU 6 7 ae 6 8 a 9 6 ae C 6 0 26 Roa Ro R dt U U U U dU 7 12 6 7 Co 7 Iy 27 R Ro dt U U U U dU 12 l4 14 4 tC 4 Io 28 Ri Ra dt Because the time constant about 10 ms for the thermal component is much larger than for the electrical component in our case the thermal component was calculated less often thereby increasing the speed of the simulation V References 1 Shi K L Chan T F Wong Y K Ho S L Modeling and Simulation of the Three phase Induction Motor Using Simulink Int J Elect Enging Educ Vol 36 pp 163 172 Manchester U P Great Britain 1999 2 Lipo Thomas A Consoli Alfio Modeling and Simulation of Induction Motors with Saturable Leakage Reactances IEEE Transactions on Industry Applications V ol IA 20 No 1 pp 180 189 1984 13 10 11 12 13 14 15 16 Soe Nyein N Yee Thet T H Aung S S Dynamic Modeling and Simulation of Three phase Small Power Induction Motor World Academy of Science Engineering and Technology Vol 42 No 79 pp 421 424 2008 Otto Jens Dynamic Simulation of Electro
9. for input motor parameters and mission profiles motorParameters m contains motor parameters specific to a particular motor design as well as calculated parameters If you are to change the name of this folder or make a new motor design the new name of the folder must be input into the main program The code is displayed in Part IV but these parameters show the motor parameters for this particular example These values can be changed to accommodate other motor EMA designs Tables 2 1 2 4 summarize the input parameters in motorParameters m Table 2 1 Motor Electrical Parameters Motor Explanations Source Electrical Parameters Rs room phase resistance at room temperature O FEM or Measurement o p numberofpoles Given inertia moment of rotor friction coefficient Measurement corresponding to current array i Measurement corresponding to current array 1 Measurement array of quadrature current 16 Table 2 2 Motor Thermal Parameters Motor Thermal Explanations Units Source Parameters winding loss Simulation Electrical Simulation Electrical Simulation Mechanical Simulation Mechanical Simulation Electrical Simulation E FEM B front bearing loss ILE oo RR IR7a conduction heat loss to gear box axle IR14 conduction heat loss to gear box case IC IS IBF IM Tj lt Table 2 3 Drive Train Parameters Parameters effDrive drivetrainefficiency Given theta melntial initial rotor
10. o Drive train efficiency theta meInitial pi 6 Rotor initial mechanical angle p 2 when RATE 505 001580297 2 D51 1 927 7 2 puu 1 5 51 0 927 5 Oo Define simulation limits xMinLim 0 xMaxLim 4 05 0 0254 tauMLim 1 9256 PLim 688 uLim 165 iLim 19 2 didtLim 5 10 45 vLim 0 086 omega meLim vLim Ncr p 2 o o Ad o o9 o o o stroke min m Stroke max m Motor torque limit N m Rated output power W Maximum instantaneous voltage allowed per Maximum instantaneous current allowed per Maximum current change rate allowed Maximum speed allowed m s stroke 0 phase V phase A A s Maximum mechanical angular speed p 2 rad s 46
11. the green line The dotted green lines represent current limits of the motor The current results are saved Currents d A and ly A I 1 dq Currents Sell O EAT N aga sa XS yoy W any a Mtns bal pv y J AA ul L ow BA E WIP te e PS Ji EDU Y np fe l la ni Y ys AN f AN y MW PY Y i Ne 1 5 2 2 5 3 3 5 4 Time t min Figure 2 3 Direct and quadrature currents 28 4 5 4 The direct and quadrature voltages are shown here in Figure 2 4 and are similar to the currents The direct voltage should remain around 0 and the quadrature voltage should contribute to the requirements of the selected flight profile dq Voltages 80 u d 60 u q 40 gt gt 20 A O c 49 LI h gt 0 ha i het m 1 peers ri em i aa Mae w MT b PUR e Aer fh AGO ceat NORINCO PHA RS O S 20 E av O gt 40 60 E 80 0 0 5 1 1 5 2 2 5 3 3 5 4 4 5 9 Time t min Figure 2 4 a Direct and quadrature voltages The ug voltage does not vary too much and generally stays around zero as seen in the zoomed section of the graph The voltage results are saved in Excel file uduq xlsx 29 Volt
12. AFRL RZ WP TR 2011 2005 SOFTWARE PACKAGE ON INTEGRATED NONLINEAR DYNAMIC MODELING AND FIELD ORIENTED CONTROL FOC OF PERMANENT MAGNET PM MOTOR FOR HIGH PERFORMANCE ELECTROMECHANICAL ACTUATORS EMAs Quinn Leland Mechanical Energy Conversion Branch Power Division Thomas X Wu Louis Chow David Woodburn Lei Zhou Jared Bindl Yang Hu and Wendell Brokaw University of Central Florida JANUARY 2011 Interim Report STINFO COPY AIR FORCE RESEARCH LABORATORY PROPULSION DIRECTORATE WRIGHT PATTERSON AIR FORCE BASE OH 45433 7251 AIR FORCE MATERIEL COMMAND UNITED STATES AIR FORCE NOTICE AND SIGNATURE PAGE Using Government drawings specifications or other data included in this document for any purpose other than Government procurement does not in any way obligate the U S Government The fact that the Government formulated or supplied the drawings specifications or other data does not license the holder or any other person or corporation or convey any rights or permission to manufacture use or sell any patented invention that may relate to them This report was cleared for public release by the USAF 88 Air Base Wing 88 ABW Public Affairs Office PAO and is available to the general public including foreign nationals Copies may be obtained from the Defense Technical Information Center DTIC http www dtic mil AFRL RZ WP TR 2011 2005 HAS BEEN REVIEWED AND IS APPROVED FOR PUBLICATION IN ACCORDANCE WITH THE ASSIG
13. Base OH 45433 7251 AGENCY REPORT NUMBER S Air Force Materiel Command AFRL RZ WP TR 2011 2005 United States Air Force 12 DISTRIBUTION AVAILABILITY STATEMENT Approved for public release distribution unlimited 13 SUPPLEMENTARY NOTES PAO Case Number 88ABW 2010 6498 Clearance Date 13 December 2010 Report contains color 14 ABSTRACT The development of all electric aircraft is a high priority in the avionics community Current aircraft use a combination of hydraulic pneumatic and electric systems for flight control However the expectation for future airplanes is a single electric system using electromechanical actuators EMAs Such a system would reduce the cost to build operate and maintain aircraft It would also make aircraft lighter more reliable safer and more easily reconfigurable reducing the turnaround for new technology One of the greatest hurdles to replacing all hydraulic actuators with EMAs is heat generation a consequence of the absence of cooling hydraulic fluid Accurately quantifying the heat generated is complicated by the highly transient and localized nature of the power demands of an EMA s motor an especially significant issue in aircraft Thus accurate modeling must be dynamic 15 SUBJECT TERMS integrated nonlinear dynamic field oriented control permanent magnet EMA electromechanical actuator lumped element model inductance flux linkage all electric aircraft coenergy eddy current three phase c
14. NED DISTRIBUTION STATEMENT signature signature QUINN LELAND JACK VONDRELL Engineer Chief Mechanical Energy Conversion Branch Mechanical Energy Conversion Branch This report is published in the interest of scientific and technical information exchange and its publication does not constitute the Government s approval or disapproval of its ideas or findings Disseminated copies will show signature stamped or typed above the signature blocks ridet ria aate tata tt ra t ttt tall raa ara aa aa ora aM aM aaa MA al aga aaa REPORT DOCUMENTATION PAGE The public reporting burden for this collection of information is estimated to average 1 hour per response inciuding the time for reviewing instructions searching existing data sources gathering and maintaining the data needed and completing and reviewing the collection of information Send comments regarding this burden estimate or any other aspect of this collection of information including suggestions for reducing this burden to Department of Defense Washington Headquarters Services Directorate for Information Operations and Reports 0704 0188 1215 Jefferson Davis Highway Suite 1204 Arlington VA 22202 4302 Respondents should be aware that notwithstanding any other provision of law no person shail be subject to any penalty for failing to comply with a collection of information if it does not display a currently valid OMB controi number PLEASE DO NOT RETURN YOU
15. NR n ENR n 1 tNR n abs PinN n abs PinN n abs P1nN n 1 tNN tNN tadd if PinN n gt O plot tNN 60 PinNN 5 else plot tNN 60 PinNN end PinNN 0 tNN tadd E end end PinNN PinNN PinN end tNN tNN tNR end if PinN end gt O plot tNN 60 PinNN b else plore ENN 60 Pi0NN e end ylabel Power it P rm W xlabel Time it t rm min title Electrical Power hold BILL o 5 Thermal figure Name Thermal DIOE ENR OO0 TEnhN LIT 2 3 45 0 7 9 9 10 12 15 44 292 35 xlswrite temperature xlsx tNR TthN ylabel Temperature it T rm C xlabel Time it t rm min title Temperature Distribution Loc 44 Part IV Motor Parameter Code motorParameter m o Rs room 0 775 p 10 L 115 2 10 1 6 7 c 0 00001 lambdaPM gs o Ad o9 o ce o Build direct inductance Ld and Define primary motor parameters This example is for a Danaher s motor Phase resistance at room temperature ohm Number of poles Inertia moment kg m 2 rad Friction Coefficient kg m 2 s rad PM flux amplitude Wh quadrature inductance Lq Load current and inductance arrays iL 0 10 300 LASN 4 469299574 4 4556 155906 Av 054450601 4 2014 1 6515 4a 0952212439 3895196891 3741345461 299904698 2424624592095 3 94295182235 25290596076 3 15910712 3 098657234242 34093991034 22 970607945 7
16. R FORM TO THE ABOVE ADDRESS 1 REPORT DATE DD MM YY 2 REPORT TYPE 3 DATES COVERED From To January 2011 Interim 22 January 2008 04 June 2010 4 TITLE AND SUBTITLE 5a CONTRACT NUMBER SOFTWARE PACKAGE ON INTEGRATED NONLINEAR DYNAMIC FA8650 09 2 2940 MODELING AND FIELD ORIENTED CONTROL FOC OF PERMANENT MAGNET PM MOTOR FOR HIGH PERFORMANCE 5c PROGRAM ELEMENT NUMBER ELECTROMECHANICAL ACTUATORS EMAs Quinn Leland Power Division Mechanical Energy Conversion Branch 3145 5e TASK NUMBER 5f WORK UNIT NUMBER AFRL RZPG Thomas X Wu Louis Chow David Woodburn Lei Zhou Jared Bindl 8 PERFORMING ORGANIZATION REPORT NUMBER Yang Hu and Wendell Brokaw University of Central Florida 7 PERFORMING ORGANIZATION NAME S AND ADDRESS ES Power Division Mechanical Energy Conversion Branch AFRL RZPG 4000 Central Florida Blvd Air Force Research Laboratory Propulsion Directorate Orlando FL 32816 Wright Patterson Air Force Base OH 45433 72581 Air Force Materiel Command University of Dayton United States Air Force 444 E 2 St Dayton OH 45402 9 SPONSORING MONITORING AGENCY NAME S AND ADDRESS ES University of Central Florida SPONSORING MONITORING AGENCY ACRONYM S AFRL RZPG Air Force Research Laboratory Propulsion Directorate SPONSORING MONITORING Wright Patterson Air Force
17. ability for the motor to follow the flight profile data Users can change the proportional constant or k value and the integral constant or k value of the control 18 22 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 231 238 239 240 33 34 35 36 37 38 39 40 41 42 43 44 45 dn 38 39 Set PI coefficients 40 En 0 1 41 A s rad 4 Get igs tauM iq 3 p 4 lambdaPM id Ld Lq N m A igS iq 2 p tauM ig If kKp efmega me dtE alpha me AP met t me dtE A Set ids id5 Of 3 id should almost always be zero A Rate limit ig5 reduces transients digdtS5 iqS iq dtE Predicted current rate A s diqdtS sign diqdt5 min didtLim abs diqdt5 Cap diqdt A s ig5 ig digdtS dtE Rebuilt igS A Saturation limit abs of iqs5 igs sign 1aS abs 1q5 lt 13 abs 15 abs iq3 i1 abs 15 iLim 1 rSat 2 4 iLim rSat abs 1q5 gt 3 13 iLim 3 A omega meP omega me Up reQnega_ me Two parameters alpha dt and beta dt located in the main code are weighted averages used in calculating next step values Generally the beta value is set to 5 but can be altered at the user s discretion for algorithm convergence Another parameter that can be changed to help improve design is spc or samples per frequency as shown below D Set up simu
18. ages u d V and Ug V dq Voltages 60 00 d M j IM N Ay P 0 9 0 95 1 1 05 1 1 1 15 Time t min Figure 2 4 b Zoomed in Section of Direct and Quadrature Voltages 30 5 The copper winding inside the motor contain the highest source of heat The transient heat analysis is represented here with peak values near 12 Watts which is much higher than the steady state power loss The results are saved in Excel file Pcu xlsx Power Loss in Windings 0 016 I 0 014 J 0 012 J kW 0 01 J cu 0 008 a 0 006 E Power Loss P 0 004 0 002 l AM MW nota Ps s y nny ang J yk x uw 0 pl L LU HMM 0 0 5 1 1 5 2 2 5 3 3 5 4 4 5 5 Time t min Figure 2 5 Power loss in the motor windings 3 6 Of great interest to power systems engineers is the power in positive values and out negative values of the motor The power out is often called regenerative power Here we show the power in and out of the motor in Figure 2 6 The results are saved in Excel file Pin xlsx Electrical Power 250 200 150 50 le ee e r e re Ms Sc ooo o Roo a RC Power P W O 100
19. ar inductance L Much of the motor modeling research done in the past 1 and 8 and even recent work 9 use constant parameter values in simulations although some did use nonlinear inductances 2 and 7 and 10 Since the BH curve is nonlinear and H is a function of current inductance is also a nonlinear function of current For our model we used FEM to obtain L and L as functions of i At this point we are neglecting changes in i since our control algorithm already seeks to maintain a zero value for i The inductance curves we do have were calculated in ANSYS using the magnetization curve for our particular soft magnetic material and the geometry of our motor The analysis was done for zero direct current and for a wide range of quadrature currents Figure 1 5 Inductance as a Function of Current Inductance Ls La mH 200 200 100 0 100 200 300 Current i A Figure 1 5 Direct and quadrature inductance as a function of current Time Stepping The size of variable arrays affects the speed of a simulation Using prerecorded input profiles high resolution is necessary to accurately capture transient moments A uniform high resolution would result in a very large profile To keep our array sizes small but maintain high resolution where needed we used a line simplification algorithm to reduce the number of points where little activity would occur This compression technique allowed us to reduce
20. ation of hydraulic pneumatic and electric systems for flight control However the expectation for future airplanes is a single electric system using electromechanical actuators EMAs Such a system would reduce the cost to build operate and maintain aircraft It would also make aircraft lighter more reliable safer and more easily reconfigurable reducing the turnaround for new technology One of the greatest hurdles to replacing all hydraulic actuators with EMAs is heat generation a consequence of the absence of cooling hydraulic fluid Accurately quantifying the heat generated is complicated by the highly transient and localized nature of the power demands of an EMA s motor an especially significant issue in aircraft Thus accurate modeling must be dynamic The parameters that characterize the motor are nonlinear functions In steady state analysis these nonlinearities can simply be averaged out and often the parameters are considered to be constant 1 But during transient analysis these nonlinearities can become critical 2 Consequently time averaged linear models 3 are inadequate for thermal analysis and management designs So proper treatment of the problem requires dynamic nonlinear analysis The motor parameters can be determined by careful experimentation 2 3 but this is a time intensive process and can be expensive especially when considering multiple potential designs for a motor Jens Otto suggested using a reduced
21. ay foractual stroke m dN array for direct voltage V dN Oo array for quadrature voltage V o idN if array fordirectcurrent A O e O PouN array for instantaneous copper loss W_ PinN _arrayforinstantaneousinputpower__ W TthN 2D array for temperature at different nodes 22 VI Expected Results The simulation has fully run when the following statistics are displayed in the command window to include Mean winding power loss max motor torque max actuator force max velocity max acceleration min electrical time step max electrical time step and the total elapsed time for the program to run Dynamic Modeling and Field Oriented Contro Mean winding power loss 1 8149 W Max motor torque 2 1564 Nm Max actuator force 2873 9049 N Max velocity 0 092486 m s Max acceleration 2 685 m s Min electrical time step 0 0005 s Max electrical time step 0 00398539 3 Elapsed time is 300 858010 seconds Also figures are plotted are a result of the flight mission profile used specific to this manual and will vary based on other motor parameters and flight mission data This profile was a five minute section of a full flight profile Acceptable values will be specific each set of parameters or profile data however some generalizations can be made about parameters such as power loss 23 The stroke profile shows that movement of the shaft in the EMA The green line r
22. cN zeros 1 Nt v Friction loss W TthN zeros Nth Nt TthN 1 Tth Temperature C o od o9 o Al o9 oO o9 o9 o o o Ao o A low resolution low res profile is analyzed one segment step at a time breaking each segment into many sub segements that have a high resolution Interpolation is used to accomplish this and the number and size of the time steps are based on an average hi res time step Size based on the maximum electrical speed of the motor During playback there is only one low res segment Nt 2 and all the profile points are hi resolution At that point the hi res time step size and the number of time steps are based on the pre recorded hi res time array Any technique used for interpolating omega on the low res scale will be yield different results than on the hi res scale just because of the scale even if the same technique is used However using a perfect N point derivative PNPD should yield fairly accurate results without averaging For each low res segment simulate 36 for nt Z NEC Oo Calculate number and size of steps for this segment dt tN nt 1 tN nt 5 Update low res segment time step size Calculate variable hi res time step based on dthetaE dt dtheta me theta meSN nt 1 theta meSN nt dtheta meAbs abs dtheta me omegaSAbs spc dtheta meAbs dt spc samples per cycle omegaSAbs omegaSAbs omegaSAbs omegaSMin omegaSMin omegaSAbs om
23. ctive features of this model is that it is able to incorporate various motor designs The software manual assumes that the mission profile is provided to the user For users familiar with control coding will be able to manipulate some constants in the control helping to tune for their specific motor designs c User Privileges The main advantages to the user is to put in their own motor parameters for alternate designs and test whether the design is capable of providing proper control and also find out the thermal performance The motor parameter code located in Appendix A contains a variety of values that must be determined by the user if an alternate design is used Values such as phase resistance number of magnetic poles and the moment of inertia must be determined Motor limitations are also set in this file to ensure that these limits are not exceeded during operation The inductance values as well as the thermal resistance and capacitance values can be obtained from FEM models If a new motor model is to be used that maintains the non linear nature of the simulation the new inductance values thermal resistance and capacitance should be recalculated using FEM Although determining these values is outside the scope of 15 this manual the user is able to manipulate any of the data inside the motor simulation to suit the needs of their design II Input Parameters There are the two file names motorParameter and missionProfile xlsx
24. dtBar mean diff tN Note that dt might not be constant tN tN 1 dtBar tN tN end dtBar tN end 2 dtBar theta meSN theta meSN 1 theta meSN theta meSN end theta meSN end tauLN tauLN 1 tauLN tauLN end tauLN end gt Initialize states Rs Rs room Initial phase winding resistance ohm id 0 Phase d current A ig U Phase q current A did dt 0 Rate of change of id A s didq dt 0 o Rate of change of iq A s o Ld at zero current LG Lt Zero current LdO LdSN 1 LqU LgqgsN I o Ld Ld0 Initial Ld Lq Lq0 Initial Lq ud 0 Phase d voltage V ug 0 Phase q voltage V theta me theta meSN 1 Start up Mechanical angle p 2 rad omega me 0 Mechanical angular frequency p 2 rad s alpha me 0 Mechanical angular acceleration p 2 rad s 2 Os 0 Heat source flux from windings O zeros Nth 1 Heat flux NConverge 0 Counter for number of convergence iterations Oo Saturation curve points roat 0 10 Saturation ratio o il iLim 1 rSat Beginning of current saturation A 13 iLim 1 rSat Ending of current saturation A ul uLim 1 rSat Beginning of voltage saturation V o s lim 1 reset Ending of voltage saturation V th2L xMinLim Ncr p 2 theta meInitial theta me max limit th2R xMaxLim Ncr p 2 theta melnitial theta me min limit tho thZR
25. e 6 shaft body AAA A Cl z 7 shaft front end R17 8 shaft back end cover plate R4a 9 magnet R14a 10 air gap ISO 11 connection gear t CD toss 12 front bearing Coppe R1 PS O R4 oy 13 back bearing IBF t loo ora NNV 3 TVW e EE VW 14 front cover Ram air Bearing C1 Stator i 15 moving mechanism loss Shaft EN C3 Tc Surface R7 Frontend A 16 wing surface R6a IM T 2 IM mo 17 gear box case R6 I agnet loss RO T R IG t IF t D Gear Sed box loss Wing ose R11 R15 Surface R16 2s 9 10 D pin body C6 Magnet C Air gap WO 2 E e A R8 1 T E NA di e AN Connection Jy e 16 1 a 8 Air flow outside ring Windage R7a gear R6b Tai C15 Shaft Back end C11 Motor assembly Figure 1 7 Lumped node model of EMA 2 Dynamic Thermal Equations After all the thermal resistor and capacitor values of the lumped node network model of the motor in Figure 5 are known the same values of the thermal resistors and capacitors can be used to simulate the temperature response of the motor parts with any combination of heat losses and boundary conditions For this standard resistor capacitor electrical network a set of first order ordinary differential equations can be written as U si dU C ZA 18 R dt c t 18 U U U U U U dU AAA APP SEIT 19 R Rio R
26. e run button at the top of the editor window E Editor C Users WUTX Desktop Working foc m File Edt Tet Go Cell Tools Debug Desktop Window Help MT ELLET BA ho tia x st ot Nonlinear Dynamic Modeling and Field Oriented Control eye Be E BS of Permanent Magnet PM Motor Written in 5I or MES Unit System NN HE Li Ri P wo dP wh dn dP 2 The following menu will pop up to add the foc m file to the MATLAB path MATLAB Editor File C Users WUTX Desktop Working foc m is not found in the current directory or on the MATLAB path To run this file you can either change the MATLAB current directory or add its directory to the MATLAB path AddtoPath Cancel Help e Click Change Directory Different versions of MATLAB will show different popup boxes when the directory is not specified Simply ensure that all the files needed to run the program are located in the same folder before changing the path 3 The simulation will begin running The runtime will vary based on the system 21 V There will be a statistical report and a total of 9 graphs displayed as well once the simulation is finished Output Parameters The Output parameters are summarized in Tables 2 7 and 2 8 for reporting statistics and plot figures Table 2 7 Output Parameters to Report Statistics Table 2 8 Output Parameters in Figures NR array for outputtime sequence s N arr
27. eat Flow 1 f m Electrical Conductvity W C m Thermal Conductivity R p Resistance R pC W Thermal Resistance F Capacitance C E Thermal Capacitance 9 AX R conductive 2 r 16 1 R convective h As AV C p V C 17 Equations 16 and 17 show how to calculate of the R and C values for simple one dimensional geometry Accurate estimates for R and C values are not possible for complicated geometries such as an electric motor The traditional lumped node thermal network is based on a forward direction modeling process In this process with all the material geometry and construction information provided the thermal resistance RO and capacitance CO are calculated based on empirical and theoretical relations For example when the winding wire type winding structure potting material slot width and teeth thickness are given the thermal resistance from winding to stator and from winding to rotor can be calculated Because this type of approach requires large number of empirical relations to achieve high accuracy a reverse modeling process is introduced in this paper A detailed 3 D solid model of a target motor is constructed and imported to an FEA software like ANSYS 16 to perform steady state and transient simulation With these results we can estimate and select the values for the thermal resistances and capacitances The advantage of this method is that we can evaluate the fidelity of the lumped node model w
28. egaSMin omegaSAbs omegaSAbs omegaSAbs lt omegaSMax omegaSMax omegaSAbs gt omegaSMax dtE 2 pi omegaSAbs NtE ceil dt dtE 1 Update number of sub hi res time points dtE dE NEE I S Update hi res segments time step size Count the total number of simulation points NtETotal NtETotal NtE 1 Reset micro loop s energy PcuBar 0 This is used to find the average power in a micro loop which is then used for the thermal component o o Simulate this low res segment in hi res Inr nth Z22NECE oo Q D Ct j i ae O G Ct 109 cr O 3 O ct O K 109 ll 3 G i C Ct Eis O ES o Variables ending in S represent desired values Variables ending in P represent next step values Variables without these are this step s values de oe o The process of getting the inputs to the motor simulation transitions from this time step to the next The inputs belong to the resolved convergent state of the next time step Therefore udP uqP tauLP thetaEP and omegaEP all belong to the next time step The motor simulation accepts certain of those inputs as given values and then resolves the remaining states of the motor such that those inputs could be true Ad o9 oe de o9 o9 oe o Determine index in Nt scale from index in NtE scale wt ntE NtE Progress in this nt segment 0 to 1 5 Use weighting to get t thetaS and omegaS t tN nt 1 wt tN n
29. epresents the actual movement of the actuator based on the mission profile The blue line 1s the desired movement of the EMA It is hard to see the blue because the green is covering indicating good control The desired stroke is in excel file missionProfile xlsx The actual stroke is saved in strokeActual xlsx Mechanical Stroke Following 2 6 Xdesired 2 4 4 actual 22 L El Stroke x in ly h FM Pen u D mes Ro V Ma M es PN 0 8 l O 0 5 1 1 5 2 2 5 3 3 5 4 4 5 5 Time t min Figure 2 1 a EMA stroke To further see the accuracy of the stroke following plot a zoomed section of Figure 1 is shown below 24 Mechanical Stroke Following 2 4 m desired Xactual 2 2 Stroke x in 0 8 0 35 0 4 0 45 0 5 0 55 0 6 0 65 Time f min Figure 2 1 b Zoomed section of stroke 2 In the following Figure 2 2 a and b the load force from mission profile and load torque converted by F Load T No TT dee where Ner is coupling ratio and Marive 1S drive train efficiency In Figure 2 2 c the magnetic torque is shown Comparing Figure 2 2 b and c we found they are close This is because the moment of inertia of the motor we used is small The magnetic t
30. heta meN theta meInitial 2 p Ncr 5 m vN diff xN diff tNR 2 m s NT length tNR tNRmid tNR 1 NT 1 tNR 2 NT 2 amp aN diff vN diff tNRmid m 8 2 o Calculate and Report statistics PcuMean mean PcuN mean winding power loss W tauMMax max abs tauMN 5 N m tauLMax max abs tauLN 5 N m FMax tauLMax Ncr effDrive 5 N vMax max abs vN m s aMax max abs aN m8 2 dtEMin min diff tHR s dtEMax max diff tHR s Spaces disp spaces Mean winding power loss num2str PcuMean W disp spaces Max motor torque num2str tauMMax Nm disp spaces Max actuator force num2str FMax N disp spaces Max velocity num2str vMax m s disp spaces Max acceleration num2str aMax m s 2 disp spaces Min electrical time step num2str dtEMin s disp spaces Max electrical time step num2str dtEMax s o Stroke 42 figure Name Stroke pLOL EN 60 Xx9N 0 0254 ENR 00 xN 0 0254 3 xlswrite strokeActunal xlsx tNR xN ylabel Stroke it x rm in xlabel Time XAititliem min title Mechanical Stroke Following legend Xitix desired xXitix actual Torque figure Name Load Force plot tN 60 FLN ylabel Load Force it F_ load rm N xlabel Time it t rm min title Load Force figure Name Load Torque
31. ith a real motor add or reduce nodes to increase the model accuracy and efficiency and estimate the maximum error between the node temperature and maximum temperature in real motor This procedure can also be extended to model the gear box motor driver and drive train and even include the aircraft wing surfaces and frames The thermal network can also be incorporated into a multi physics model to simulate the electrical thermal and mechanical performance of the whole EMA and its supporting structure Once the proper thermal resistances and capacitances are selected this simulation engine can be used with various time dependent boundary conditions during the whole mission duration including the air temperature aircraft speed altitude and sunlight 10 The computational requirement of such a simulation is negligible compared to the FEA simulation Due to the symmetry the geometry of the model shown in Figure 1 1 is divided into a quarter section and imported into ANSYS for the FEA simulation The nodes are numbered in Figure 1 6 These numbers are also the same as those in the lumped node network shown in Figures 1 7 Air gap Node 10 Magnet Node 9 Rotor core Node 6 Case Node 4 Case surface Node 5 Stator iron Node 3 Copper winding Node 1 Figure 1 6 Thermal node locations on motor 11 1 copper Detailed Lumped node model Of EMA transient 3 stator 4 case and back cover 5 case surface
32. lation parameters H ipc 10 t Samples per cycle of highest frequency aulh 0 0 t Thermal time constant s t Set PI coefficients kp 0 1 ki 0 0001 A s rad Set differential time stepping weights 19 HI Software Startup Once all the proper m files are verified click on the Start Menu and open the MATLAB software through the program listing ME e File Bi Editor C Users WUTX Desktop Working foc m _ Once MATLAB is open go to File and click on Open Go to through the above mentioned directory path to find the FOC folder Click on the FOC folder and then click on the file name FOC m The following code should be displayed in the MATLAB editor screen Tet Go Cell Tool Debug Desktop Window Help OGE Saee Ma P E BBM BB Sick Bese 10 ui x off of 0 i cO 1 Ch CF to Ri FP Kk Pp FP FP PR PR FP F HK Pp i ro 5 Oe C PRO p oc 2 Nonlinear Dynamic Modelin i of Permanent Magnet EM 3 t Written in SI or MES Unit t 2 Authors 3 David Woodburn 3 Dr Lei Zhou t Dr Thomas X Wu clear all close all clc motorParameter 20 q and Field Oriented Control Motor System IV Software Execution and Procedure 1 In order to run the simulation you can either click in the editor window and press the F5 shortcut key or click th
33. mechanical Systems using ANSYS and CASPOC 2002 International ANSYS Conference ANSYS 2002 Roshen Waseem Iron Loss Model for PM Synchronous Motors in Transportation 2005 IEEE Conference on Vehicle Power and Propulsion pp 4 ISBN 0 7803 9280 9 2005 Dolinar D Weerdt R De Freeman E M Calculation of Two axis Induction Motor Model Parameters Using Finite Elements IEEE Transactions on Energy Conversion Vol 12 No 2 pp 133 142 1997 Topcu E E Kamis Z Yuksel I Simplified numerical solution of electromechanical systems by look up tables Mechatronics Vol 18 No 10 pp 559 565 Elsevier 2008 Mohamed M A Nagrial M H Modelling and Simulation of Vector controlled Reluctance Motors Drive System International Conference on Simulation No 457 pp 380 384 ISBN 0 85296 709 8 1998 Sun Fengchun Li Jian Sun Liging Zhai Li Cguo Fen Modeling and Simulation of Vector Control AC Motor Used by Electric Vehicle Journal of Asian Electric Vehicles Vol 3 No 1 pp 669 672 Asian Electric Vehicle Society 2005 Demerdash N A Gillott D H A New Approach for Determination of Eddy Current and Flux Penetration in Nonlinear Ferromagnetic Materials IEEE Trans MAG 10 pp 682 685 1974 Merzouki R Cadiou J C Estimation of backlash phenomenon in the electromechanical actuator Control Engineering Practice Vol 13 No 8 pp 973 983 Elsevier 2004 Mellor P H Robert
34. nto mission profile EN Missrion 1 time sequence XSN Mission 2 stroke sequence FLN Mission 23 load force sequence RON XSN 0 0254 If stroke is in inch needs to convert to meter FLN FLN 4 4482 If load force is in lbf needs to covert to Newton o o oe o9 oe oP oe theta meSN xSN Ncr p 2 theta meInitial tauLN FLN Ncr effDrive 5 Load torque N m spc 10 Samples per cycle of highest frequency tauTh 0 01 Thermal time constant s Set PI coefficients kp 0 1 ki 0 00015 L s rad l o Oo Set differential time stepping weights alpha dt 0 5 beta dt 1 alpha dt Set motor simulation convergence parameters epsilon 0 0001 5 Convergence limit o Ncvg 50 Maximum number of convergence iterations T tN end tN 1 T could have been redefined 34 Nt length tN Initial length of profile arrays nfMax 1 Order of highest frequency harmonic of interest NtR round spc nfMax omega meLim 2 pi T High res time steps o dtEO T NtR 1 Mean hi res time step omegaSMin 2 pi 0 002 Min sampling frequency rad sample s omegaSMax 2 pi dtE0 Max sampling frequency rad sample s Mek 2 Hi res index for recording NtETotal 1 Counter for the total number of simulation points trh 0 Last time for thermal calculation s Pad arrays for interpolation referencing Lengths are Nt 3
35. onvergence direct quadrature reference frame 16 SECURITY CLASSIFICATION OF 17 LIMITATION 18 NUMBER a REPORT b ABSTRACT c THIS PAGE OF ABSTRACT OF PAGES Unclassified Unclassified Unclassified SAR 54 19a NAME OF RESPONSIBLE PERSON Monitor Quinn Leland 19b TELEPHONE NUMBER nclude Area Code N A Standard Form 298 Rev 8 98 Prescribed by ANSI Std 239 18 Table of Contents Part Ls Technical MODO eius A A A MM Ye e 1 PM e e AR E l IL Simulation Flow Chart 3 MAN OO a 4 l Diam ca f at OS ias A 4 Za Tne SUS OI PA e I OEA 6 3 Pica Onenied Control TOC sr ro tao 7 NAAA q e 9 l hueso Me E E 9 2 Diname Thermal EOUOUOLDS aaa 12 i SU e RET 13 Ign Hs Users Manda cee ceive cite 15 MB o USO MNT e o o 15 a Purpose Ol UC SOR ANG sii RN UE UE 15 b POEL Wy ate RII TH aklama 15 C IU o PAVI EET 15 W ANT Pa AC YY EE 16 M ORES SUD los 20 AE a io dud Procedure MEAN Bee ya O ya e ais 20 Vo Output Parame LEE SE ga re ai lr n a PPEEEPEPU lr g l RU yu lr g l caga 2 RDA m 22 iii Part TIL source Code LOC iet oeste caia ParttlIV MotorParameterCode motorparameters m IV Part Technical Manual I Introduction The development of all electric aircraft is a high priority in the avionics community Current aircraft use a combin
36. order coenergy model 4 which yielded results comparable to FEM results Likewise Roshen considered more involved empirical formulas which accounted for excess eddy current losses 5 Through FEM modeling the nonlinear parameters such as self and mutual inductances of the three phase windings can be determined 6 7 Though dynamic FEM can be very accurate 4 such detailed simulation is very slow Fortunately FEM is only necessary to quantify the parameters it 1s not necessary for simulating the motor over lengthy mission profiles Instead a nonlinear lumped element model NL LEM can use the parameters from the FEM model and then dynamically simulate the motor just as accurately as the FEM model but at a much lower computational cost To our knowledge we are the first to address the nonlinear dynamic modeling of a permanent magnet motor and describe both the control and thermal performance of the motor in following highly transient mission profiles Figure 1 1 is the 3 D model of a typical Permanent Magnet Synchronous Machine PMSM servo motor This design features a 12 slot stator and a 10 pole rotor Figure 1 2 In the following we will discuss the procedures for electrical and thermal modeling Figure 1 1 3 D model of the PMSM motor design Case Copper winding Iron Figure 1 2 Cross section of motor and rotor II Simulation Flow Chart Before describing
37. orque and load torque are related by do dt TM TTL The magnetic torque is saved in Excel file tauM xlsx 25 Load Torque t N m 3000 Load Force 2000 1000 F 1000 2000 3000 O o URP A nut 0 5 2 2 Time f min 5 3 3 5 Figure 2 2 a Load Force Load Torque 4 5 1 5 A m A oves LT A 0 5 espe ho V Y TAN NM T n JY FAL rod eme A n y y M M W i J A de NI M NI WY i ATUN hn f Nihat 1h N JU pe Y AA QAM v sh ae eN Y Nun 2 2 5 Time t min 26 3 3 5 Figure 2 2 b Load Torque Aue gt 2 9 1 5 gt o o o Magnetic Torque TIN m O zu Oi 2 Magnetic Torque 0 5 1 1 5 2 2 5 3 3 5 4 Time ft min Figure 2 2 c Magnetic Force in Excel file idiq xlsx 21 4 5 The dq currents provided to the motor are represented in Figure 2 3 Notice that the current ig nearly constant around 0 and doesn t exceed the dotted blue line limitations This is a proper reading as direct current does not contribute to the rotation of the rotor Only i contributes to the torque as seen by the variation in
38. que is the mechanical angular acceleration multiplied by the number of pole pairs p 2 I is the rotor s moment of inertia and c 1s the rotor s coefficient of friction from windage and bearings where u is input voltage 7 is direct current i is quadrature current R is phase resistance L is direct inductance L 1s quadrature inductance p 1s the number of poles is the mechanical frequency multiplied by the number of pole pairs p 2 A py 18 the flux linkage from the permanent magnet 7 is the motor torque generated by the magnetic fields 7 is the load torque is the mechanical angular acceleration multiplied by the number of pole pairs p 2 I is the rotor s moment of inertia and c 1s the rotor s coefficient of friction from windage and bearings The motor losses are related to the motor parameters such as R c L and L Since losses have an effect on motor behavior they should be modeled through the motor parameters not calculated 1n post processing as 1s often done 5 The copper loss can be calculated as E i rd i R 5 or equivalently af 5 ti kR 6 for balanced 3 phases In either case note that the resistance 1s the parameter directly associated with the power loss in the windings Likewise windage and bearing losses are directly associated with the coefficient of friction c All iron losses hysteresis classical eddy and excess eddy losses can be associated with the nonline
39. rated Pou 3 2 iq 2 id 2 Rs Next step source heat W PcuBar PcuBar Pcu NtE t Tih gt Eauih Get next step winding resistance Rs RS room dtIh t tIh tTh t Calculate temperatures Oo IC Pcu 4 3 winding loss only quarter portion is used Oo Calculate next step temperatures TthP Tth l4 IR1A4 Tth 4 Tth 14 BRth 3 2 Tth l4 4BEH l3 Cth b4 dtlrhs TthP 1 Tth 1 IC Tth 3 Tth 1 Rth 1 Cth 1 dtTh TEHP 2 Tth 2 TthP 3 Tth 3 IS Tth 4 Tth eee Ith i0 Tth 3 Bth 12 du UM ie Rtn 1 7 Cen A dera Tehp 4 Tth d TERh 4 TEh 1 qo 3 lone anaes bo Rth 14 Pee TERh 4 JRER Z TEh 5 TEh A4 REh A CEh A dtTh TENPTI5 mq HA REh 4 TEh 2 REh A4 Reh 5 TthP 6 Tth 6 Tth 9 Tth 6 Rth 6 i Tth 7 Tth 6 Rth 7 Tth 8 Tth 6 Rth 8 Cth 6 dtTh TChP 7 Tth 7 IR7a Tth 6 Tth 7 Rth a Lo Tth 12 Tth 7 REh 9 Cth 7 dtTh TthP 8 Tth 8 Tth 6 Tth 8 REh 8 Tth 13 Tth 8 Rth 10 Cth 6 dtTh TthP 9 Tth 9 IM Tth 6 Tth 9 Rth 6 Tth 10 Tth 9 7Rth 11 Ceh 9 deTh TthP 10 Tth 3 Rth 11 Tth 9 Rth 12 Rth 11 Rth 12 IW Rth 11 Rth 12 TthP 12 Tth 12 IBF Tth 14 Tth 12 Rth 13 Tth 7 Tth 12 Rth 9 Cth 12 dtTh TthP 13 Tth 13 IBB Tth 4
40. s D Turner D R Lumped parameter thermal model for electrical machines of TEFC design IEEE Proc B Vol 138 No5 Sept 1991 DiGerlando A Vistoilo I Thermal network of induction motors for steady state and transient operations analysis ICEM 1994 Paris Motor CAD v3 1 7 Motor Design Ltd www motor design com Y K Chin D A Staton Transient thermal analysis using both lumped circuit approach and finite element method of a permanent magnet traction motor IEEE AFRICON 2004 ANSYS V12 0 ANSYS Inc www ansys com 14 I Part Il User s Manual Introduction Purpose of the Software This software provides integrated non linear dynamic modeling including both electrical model and thermal model together with a field oriented control scheme Different motor configurations can be entered into the software The simulation provides the user with all the necessary motor information such as mission profile following motor torque motor current phase voltage input power power losses in the windings and a thermal profile This comprehensive and detailed look at motor control together with heat generation and transfer will provide solid foundation for users to design highly efficient EMAs Software Capability Package of a nonlinear dynamic modeling of a permanent magnet motor that describes both the control and thermal performance of the motor in following highly transient mission profiles One of the most attra
41. sonable simplification since for a round rotor 4 is independent of i should normally be zero and the i factor is divided out u pen ti L L 15 q For the non round rotor case our direct interpolation method could be used to include the dependence on theta Even after the desired control values were calculated the reality of limits current and voltage needed to be included Instead of directly capping those values we realized that the system converged more quickly and ran more smoothly if the functions used were all smooth 11 For our case we used a segmented curve of first and second order polynomials IV Thermal Model 1 Model Setup A lumped node thermal network is to represent the temperature of every solid part of the EMA The thermal resistances and capacitances between the nodes can be treated as electrical resistances and capacitances Hence the temperature of every node can be solved as the voltage in this equivalent network Table 1 1 This approach is well developed in motor design industry 12 13 Some commercial motor design software package has already included the thermal network simulation i e Motor CAD 14 Some studies have been made with FEA analysis and experimental testing have shown that such an approach is valid 15 Table 1 1 The analogy of equivalent thermal circuit Electrical Circuit Thermal Circuit V V Voltage PE Temperature I Al Current Flow Q i H
42. t 1 wt s tauLP tauLN nt l wt tauLN nt 1 wt 5 N m theta mes theta meSN nt l wt theta meSN ntt1 wt 5 rad den EN erimez EN besle nelik s omega me3 theta_meSN nt nt 2 theta_meSN nt 1 nt 1 dt3 gt rad s omega meS wt lt 0 5 omega me3 1 0 5 wt 37 omega me3 2 0 5 wt Lus wt 0 5 omega me3 2 wt omega me3 3 wt 0 5 gt Drag Get theta and omega errors eTheta me theta meS theta me Mechanical rad eOmega me omega meS omega me Mechanical rad s Get igs tauM iq 3 p 4 lambdaPM id Ld Lq N m A iqS iq 2 p tauM ig I kp eOmega me dtE alpha me C kp eOmega me ki eTheta me dtE 5 A Set ids idS 0 id should almost always be zero A Oo Rate limit iqS reduces transients diqdtS iqS iq dtE Predicted current rate A s digdis sSxon digdt5 muintldidtLuim abs diqdgts 3 5 C 1qS iq digdtS dtH Rebuilt igs A ap diqdt A s 5 saturation limit bs or 195 19S srgm iigs i abs igs r3 abs ni1qo5 abs uqo o ril i abs 1qs 1Lim l YSat 2 M 1Lim rsat abs 195 gt 13 iLim A omega meP omega me kp eOmega me Oo Calculate inductances based on desired currents mLq NL 1 abs iqS iLMin iLMax iLMin 1 Unbounded real wLq mod mLq 1 Real number 0 1 nLq mLq wLq Unbo
43. the details of the simulation we give the overall flow chart in Figure 1 3 as follows Input Motor Parameters Input EMA Mission Profile Set up Simulation Parameters Low Resolution Loop Start Calculate High Resolution Time Step High Resolution Loop Start Controller Ug Motor Electrical Simulation Motor Thermal Simulation High Resolution Loop End Low Resolution Loop End Post Processing Plot figures etc Figure 1 3 Flow chart for simulation For the motor electrical simulation block the flow chart is given in Figure 1 4 Motor Electrical Loop Start Calculate 1 1 6 0 Yes Motor Electrical Loop End Figure 1 4 Flow chart for motor electrical simulation III Electrical Model 1 Dynamical Equations The motor dynamics are modeled by four primary dynamical equations in direct quadrature reference frame dg0 di u Ri tL Ope L 1 1 me qq di u R i tL uro lala 0 Ap 2 fy i F i i ti L L 3 2 qom a co 4 p where u 1s direct input voltage u 1s quadrature input voltage i 1s direct current i 18 quadrature current R is phase resistance L is direct inductance L 1s quadrature inductance p is the number of poles is the mechanical frequency multiplied by the number of pole pairs p 2 Apy is the flux linkage from the permanent magnet 7 is the motor torque generated by the magnetic fields 7 is the load tor
44. the profiles to less than 4 of their original size 6 That said the electrical time constant is often much smaller than the time step size of the compressed profiles Since the frequency of the electrical component of the simulation 1s directly related to the first derivative of the input angular position profile for the motor we can dynamically scale the electrical time step non uniformly Therefore we had two loops one macro loop that stepped through the elements of the profile and a nested micro loop with very small but variable time steps for the electrical equations The result of this apparent complication was a many fold speed improvement This is a tremendous advantage of not relying on MATLAB built in ODE solvers which do not have the advantage of knowing how big of a time step to take but must rather continually search for it Field Oriented Control FOC A thorough analysis of motor heat generation should include how the motor is driven Our control equations are derived from the primary simulation equations 1 through 4 Starting with equation 4 and noting that in 3 i can be factored out we get e a dem 7 L p which can be rearranged to be T e c0 gt 8 At this point we reinterpret amp which is do dt tobe e dt where e is the error in with a tuning coefficient k in front This gives me 2 T Ik tco E p dt i nm 9 Tu ig
45. unded integer Las L wbq LasN ng wLa LasN nbartlls H Las l1 wLq LqgSN nLq wLqg LgSN nLq 1 5 H The above is equivalent to using interpl but much faster Las interpl iL LbdSN abs ig5 LqS interpl iL LgSsN abs 19S oP o9 o Calculate new voltages uds Rs idsS Ld5 ridS id dtE omega meP Lg5 igs5 2 V ugs Rs i0S LgS igS iq dtE omega meP LdS i1dS Oo omega meP lambdaPM 5 V Saturation limit abs of voltages ud sSxgn udS abs udS ua abs udS abs udS gt ul lt abs udS uLim 1l rSat 2 4 uLim rSat abs uds gt u3 uLim V ug sion 1as abs uqgs lt us abs 1qs abs ugqS gt ul a abs uqS uLim 1 rSat 2 A4 uLim rSat abs ugS gt u3 uLim Vv Seek convergence for one hi res time step 38 Assume next values from present values Simulation should not assume any knowledge of S values idf id did dt dtEb A Lor ig s dig di SL A alpha meP alpha me Angular acceleration rad s 2 omega meP omega me alpha me dtE Oo Run dynamical equations until convergence for nConverge 1 Ncvg Oo Get next inductances by look up table mLq NL 1 abs iqP iLMin iLMax iLMin 1 Unbounded real wLq mod mLq 1 Real number 0 1 nig mhLq wha Unbounded integer LAP l wLdq LbdSNinLbdq wLa LasN aLa L LH LOP l whq
46. where i 7 is the desired current for the next time step and Qu is the desired next step angular speed equal to w k e Using 4 to substitute for 7 in 9 we arrive at 10 Tu i 10 2 P ME REO P At Ty e s On an lie 11 gp At 7 This will track the desired motor velocity well but there will be some rotor angle drift over time To counter this drift error we added a theta error term ke dt Qs Cg i 1 K zn f a 12 2 du ue e Jen 7 dt n The k and k coefficients are the only PI values that need to be tuned in the control algorithm For maximum efficiency it is generally desired that the direct current be zero since any flux generated in the direct axis would not contribute to torque However direct current does not necessarily directly translate to direct flux As the following equation shows direct current can contribute to useful torque if the direct and quadrature inductances are not equal as given in 3 Most of the time however the desired direct current set point on should be zero With the desired currents determined the desired voltages can be found M u R i L Ong Ly 13 T u R i L Pu L i FO Apy gt 14 where L and L are based on the desired current values It should be noted that while many of the variables are updated based upon the desired currents 7 was not updated for the control equations This is a rea
Download Pdf Manuals
Related Search