Design Optimization of Hydropower Generators
Design Optimization of Hydropower Generators
Generators
This Master’s thesis is written as the final assignment of the five year long Master of Science program
at the department of Electrical Power Engineering at the Norwegian University of Science and
Technology. The report is centered around optimization of the synchronous generator where a
genetic algorithm has been implemented in MATLAB to minimize the cost of the machine. The
final result is intended to be used for education and research purposes.
I want to thank my supervisor, Arne Nysveen, who has been very helpful with clearing up aspects
about the assignment and, in particular, the synchronous generator. I also want to thank PhD
candidate Erlend Engevik who has made himself available throughout my assignment, whether
concering MATLAB syntax or questions about the genetic algorithm. Furthermore, the work done
on this assignment would not have been possible without the design calculation program made by
former students Aleksander Lundseng and Ivar Vikan.
I want to thank my fellow Master’s students at room E320; Anders Espeseth, Andreas Bock and
Petter Christiansen. They have made each day of writing a bit more interesting and fun. Lastly, I
want to thank my good friends Kailin Jones and Duncan Scott for assisting me with English spelling
and grammar.
Date
Signature
i
ii
Abstract
For education and research purposes, a design calculation program for hydropower generators has
been made in MATLAB. The user can change design parameters to optimize the generator’s cost
and efficiency. However, manual optimization is a slow and tedious process and the user might
make mistakes which will lead to incorrect results. By implementing an optimization in a computer
program, this process will be made faster and more reliable. In this assignment, a genetic algorithm
is implemented in MATLAB to optimize the synchronous generator with respect to the total cost of
the machine, which includes cost of materials and machine-losses. The purpose of the optimization
was to find the lowest machine-cost by varying chosen geometrical parameters. The variables were
selected during conversations with supervisor Arne Nysveen and Ph.D. candidate Erlend Engevik.
MATLAB has embedded functions that can perform genetic algorithms when provided with objec-
tive function, number of variables and constraints. The design calculation program mentioned has
been used actively and made the foundation for the objective function to the implemented optimiza-
tion. It has been adjusted to meet the demands of objective functions for optimization problems
in MATLAB. Furthermore, limits to the machine-design have been implemented as nonlinear con-
straints. These are gathered in a script under the function name nonlcons. For every set of variable
values, the algorithm call this script to check if the constraints are met.
To test the implemented optimization, the price of load-dependent losses was increased to determine
whether this made the algorithm reduce these losses. Furthermore, the power factor was increased
in the ambition this would make the machine cheaper. Lastly, the synchronous reactance was
increased. For high values of the synchronous reactance the air gap drops, which in turn led the
total magnetization to reduce, and consequently the total cost of the machine decreased also. Ten
optimizations were performed for each case in anticipation of variance in results from the genetic
algorithm. The results were then gathered in Excel in order to comparatively examine how the
algorithm had performed.
Results show that the optimization performs as expected. The load-dependent losses dropped by
8.1% when the cost of these increased from 15000 to 20000 kroner per kW. The cost of constant
losses were at 30000 kroner per kW for both cases. They did not change considerably. When the
power factor was increased from 0.8 to 1.0, the total cost of the machine dropped by 4.4%. This
is because of a reduced magnetization of the machine that reduces the required field current and
consequently the copper losses. For the case of a varying synchronous reactance, the magnetization
of the machine reduced by 9.3% as the reactance moved from 1.2pu to 4.9pu. Consequently the
total machine-cost decreased by 7.6%.
Test results from one optimization to the next also showed large variation between ideal solutions.
It is believed that a larger population could bring down this variation. However, the genetic algo-
rithm includes randomly generated mutation which for complex problems necessarily lead to varying
results. Furthermore, the ratio between cost of losses and materials, set by the user affects machine-
design to a great extent. Increasing the cost of losses while keeping the material costs constant causes
the use of copper in the machine to increase significantly. The constraints added ensure that the
algorithm designs reliable generators that fulfill demands of the Norwegian grid. In consideration of
this impact on machine-design, it is important that the user can adjust these constraints according
to the problem at hand.
iii
iv
Sammendrag
v
vi
NORGES TEKNISK-NATURVITENSKAPELIGE UNIVERSITET
NTNU
MASTEROPPGAVE
Fag : ELKRAFTTEKNIKK
Oppgavens tekst:
For education and research purposes, we are developing a design calculation program for hydropower
machines. So far, it uses manual optimization such that it is the user that finds the optimum design by
changing design parameters. This is time-consuming and may lead to wrong results for inexperienced
uses. The program calculates the machine rating and parameters using analytical formulas.
The purpose with this project is to establish a calculation tool that implements optimization algorithms on
electrical machines. The objective will be to establish a MATLAB program that can perform a Genetic
Algorithm on an optimization problem involving a synchronous machine. With such a tool, one can
change a parameter or constraint value from one optimization to the next to get an understanding of how
the parameter affects the optimal design.
Further details of the work to be discussed with the supervisors during the project period.
vii
Oppgaven gitt : 15. januar 2016
Oppgaven revidert: : 19. mai 2016
Besvarelsen leveres innen : 10. juni 2016
Besvarelsen levert :
Utført ved (institusjon, bedrift) : Inst. for elkraftteknikk/NTNU
Kandidatens veileder : Erlend Engevik, NTNU
Faglærer : Professor Arne Nysveen
Arne Nysveen
faglærer
viii
Table of contents
Preface i
Abstract iii
Sammendrag v
Figures xi
Tables xii
1 Introduction 1
2 Theory 2
2.1 Synchronous generators . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
2.1.1 Construction and design calculations . . . . . . . . . . . . . . . . . . . . . . . 2
2.2 Optimization problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2.2.1 Bounds and constraints . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
2.2.2 Solving an optimization problem . . . . . . . . . . . . . . . . . . . . . . . . . 14
2.3 Genetic algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.3.1 Initializing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.3.2 Crossover and mutation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.3.3 Final generations and termination . . . . . . . . . . . . . . . . . . . . . . . . 21
3 Earlier work 23
4 Method 24
4.1 Optimization methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
4.2 GenProg . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
4.3 Modifying GenProg . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
4.4 Objective function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
4.5 Selecting variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
4.6 Bounds, constraints and penalty . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
4.7 Performing the algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
4.8 Testing of the implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
5 Results 37
5.1 Changing the price of losses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
5.2 Changing the power factor . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
5.3 Changing the constraint for synchronous reactance . . . . . . . . . . . . . . . . . . . 40
6 Discussion 42
7 Conclusion 48
7.1 Further work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
Referances 50
Appendix B, MATLAB-scripts 52
Appendix C, Optimization results 54
x
Figures
1 Rotor and stator in synchronous generators [1] . . . . . . . . . . . . . . . . . . . . . 2
2 Lamination in stator [2] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
3 Roebel bar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
4 Damper windings [2] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
5 Rotor parameters [2] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
6 Peaks, MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
7 Feasible region for an optimization problem [3] . . . . . . . . . . . . . . . . . . . . . 14
8 Evolved antenna ST5-3-10 [4] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
9 Initial population . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
10 Fully performed Genetic Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
11 Types of children for the next generation [4] . . . . . . . . . . . . . . . . . . . . . . . 19
12 Evolution in the genetic algorithm [5] . . . . . . . . . . . . . . . . . . . . . . . . . . 19
13 Genetic algorithm when crossover fraction is 1 . . . . . . . . . . . . . . . . . . . . . 20
14 Genetic algorithm when crossover fraction is 0 . . . . . . . . . . . . . . . . . . . . . 21
15 Required input values, GenProg . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
16 Genkalk collapsed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
17 nonlcons collapsed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
18 Optimization application, MATLAB . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
19 Spacing between adjacent poles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
20 Code section that implements a penalty for high cooling air flow speed . . . . . . . . 34
21 Examples on how to implement option settings . . . . . . . . . . . . . . . . . . . . . 35
22 Flowchart of the implemented genetic algorithm . . . . . . . . . . . . . . . . . . . . . 37
23 Optimal set of variable values from ten optimization for the normal case scenario . . 43
24 load-dependent and constant losses in when changing the price of losses . . . . . . . 44
25 Total cost when changing the power factor . . . . . . . . . . . . . . . . . . . . . . . . 45
26 Cost of losses when changing the power factor . . . . . . . . . . . . . . . . . . . . . . 46
27 Total cost when the synchronous reactance increases . . . . . . . . . . . . . . . . . . 47
28 Total magnetization when the synchronous reactance increases . . . . . . . . . . . . 47
29 cost . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
30 nonlcons collapsed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
31 Changing cost of losses, losses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
32 Changing cost of losses, geometrical parameters . . . . . . . . . . . . . . . . . . . . . 55
33 Changing power factor, costs . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56
34 Changing, geometrical parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
35 Changing constraint for synchronous reactance, costs . . . . . . . . . . . . . . . . . . 58
36 Changing constraint for synchronous reactance, geometrical parameters . . . . . . . 59
xi
Tables
1 Losses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
2 Geometrical parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3 Cost of materials and losses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
4 Geometrical parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
5 Cost of materials and losses . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
6 Geometrical parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
7 Magnetization and current . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41
8 Generator specifications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
xii
1 Introduction
The first electromagnetic generator was made by Michael Faraday who developed the operating
principle in the years of 1831-1832. The generator, called the Faraday Disk, produced a small DC
voltage and was very inefficient. The first generator to be used in industry, the Dynamo, was made
in 1844. It used electromagnetic induction to produce direct current with commutators. Later, the
AC alternator, first two phase and then three phase, made its entry into the market. Today, the syn-
chronous generator is by far the most common machine used for electricity production in the world.
Though the main principles remained unaltered, the machine has evolved since its inception and
small improvements increasing stability and efficiency continue to be made. For teaching purposes,
a MATLAB program that performs the design calculations on synchronous generators from a set of
input values has been made at NTNU. The program was made by former students during work on
their Master’s theses. To get an even better understanding of the complexity of the machine, it is
in NTNU’s interest to research if the generator design can be optimized with the use of computers.
In this context, the following project was suggested: ”The purpose of this project is to establish a
calculation tool that implements optimization algorithms on electrical machines. The objective is
to establish a MATLAB program that can perform a genetic algorithm on an optimization problem
involving a synchronous machine. With such a tool, one can change a parameter or constraint value
from one optimization to the next to get an understanding of how the parameter affects the optimal
design.”
Thus, the main objective of the assignment is to develop a MATLAB script that implements a
genetic algorithm to optimize the design of the synchronous generator. The design script made by
previous students is provided as a modelling tool, however the focus of this project will be directed
more towards the optimization of the generator, not the design calculation itself. The process of
designing a synchronous generator is very complex. Around 100 parameters have to be decided
on, and the design involves both mechanical, electromagnetic and thermal calculations. A detailed
understanding is outside the scope of this report.
This assignment will first go through necessary theory on the geometry and design calculations
of the synchronous generator, before previous research made by other students is presented. The
methodology section will then explain the work done to make the script, and how it will be used to
inform the design of the synchronous generator. Results from the various tests are then presented,
before being analysed in the discussion section. Lastly, the conclusion will reiterate the general
findings of the study.
1
2 Theory
Synchronous generators are preferred in large power plants and produce the majority of the electricity
in the world. While the induction generator can only consume reactive power, the synchronous
generator has the capability to both produce and consume reactive power. In fact, an induction
generator has to be connected to an external source of reactive power to maintain its magnetic field
in the stator. As the grid is inductive by nature the ability of a synchronous generator to produce
reactive power is a great advantage. However, the induction generators simple construction makes it
a good option for windmills and other supplementary power sources connected to an existing power
system. As for large hydropower plants, the synchronous generator can therefore be seen as a more
advantageous choice.
A synchronous generator consists of a stator and a rotor. The rotor is locked to the synchronous
speed of the generator. The generator has a separate field circuit with a DC field current. This
current produces a magnetic field that rotates at the same speed as the rotor. The magnetic field
induces an AC voltage in the stator winding, in turn producing a lagging current flow when connected
to a load. In some machines the rotor magnetic field is produced from permanent magnets, however
that will not be discussed further in this report. Synchronous generators have either salient or non-
salient poles. A non-salient pole means that the pole is embedded within the surface of the rotor,
while a salient pole protrudes radially from the shaft of the rotor. Figure 1 shows the stator and
rotor of synchronous generators with salient and non-salient poles. For machines with more than 4
poles, salient poles are the most common. Only these will be discussed further in this report.
When designing a synchronous generator, some dimension criteria has to be determined. The cus-
tomer will decide on the number of generators, their power rating Sn , and power factor pf , also
2
referred to as cos(φ). The frequency f is given from the grid, which is 50Hz in Europe. The customer
will also decide on an appropriate synchronous speed ns measured in rounds per minute (rpm). The
number of poles Np can be calculated from Equation 1:
120f
Np = (1)
ns
The output voltage is usually not set by the customer. In Europe most generators are connected to
a transformer, so it makes little difference for the customer what the output voltage is. Having the
output voltage as a free parameter is also very convenient for the generator constructor.
The cost of losses in the machine are of great concern to the customer and therefore important when
deciding on the material to be used in the stator core. The customer usually specifies the price of
load-dependent and constant losses. It is in their interest to buy a machine with minimal losses
while still meeting the dimension criteria. The cost of constant losses are commonly higher than
load-dependent losses since the machine will not run at full load at all times.
The synchronous reactance Xd is of significance for the static stability of the generator and an upper
limit is often specified by the customer. The transient reactance Xd0 and sub transient reactance Xd00
are less commonly set by the customer. In addition to demands from the customer, the generator
constructor has to meet the formal requirements for rotating electrical machines, as described in
IEC 60034-1 [6]. If demands from the customer are less strict than the standard, the standard
requirements should be followed rather than demands from the customer. Listed are some of the
types of requirements in the IEC 60034-1:
• Thermal tests are to be performed on the generator. This is to ensure that the winding will
not overheat, which can lead to a breakdown of the insulation. Limits of temperature increase
during testing are given in the standard.
• The harmonics in the generator should be minimized. Limits to the Total Harmonic Distortion
(THD) are given in the standard.
• The generator should be able to run normally for thirty seconds with a current 1.5 times
the nominal current. It should also withstand a short circuit current of twenty-one-times the
nominal current [7].
In addition, the system operator often specifies requirements to the moment of inertia to ensure
stability in the grid. This will be discussed further in section 4.5.
When the requirements are set one can start calculating the geometrical parameters of the generator.
The inner diameter, Di , is measured from the shaft center to the stator, including the air gap between
rotor and stator. It can be set by the customer or calculated from empirical formulas. The gross
iron length, lb , is also calculated from empirical formulas if it is not given. From Di and lb , the
utilization factor C can be calculated using:
Sn [kV A]
C= (2)
Di2 ·lb ·ns
where Sn is apparant power in kVA. If the gross iron length is not yet known, an empirical formula
for the utilization factor freqently used is:
3
Next, the number of stator slot, Qs , is decided. The stator is often made out of sections of laminated
steel put together to form a circle. The circles are placed after each other to the desired length, with
air gaps in between for cooling. The gross iron length lb includes air gaps between the steel plates.
A typical section is shown in figure 2.
Qs has to be dividable by the number of phases and by 2 to obtain balance in the machine. Often,
the number of slots per pole is calculated using:
Qs F ·N
= (4)
Np F ·U
where F is the largest common multiplier. N is the number of slots in which a winding repeats itself
and the number of slots per pole pair. U is a factor depending on F and N. The number of parallel
circuits pnr in the machine can be chosen as a factor of F. If F is 8, then pnr can be 1, 2, 4 or 8.
Further, the slot pitch τs can be calculated as:
π·Di
τs = (5)
Qs
The slot pitch is the arc length from one slot to the next and can be divided in slot width bs and
tooth width bt . Usually, the ratio bs /bt is within 0.55-0.75. The pole pitch τp is also calculated as:
π·Di
τp = (6)
Np
An estimated armature loading As can be calculated:
C · 60 · km · 10
As = fj
[At/cm] (7)
π2 · kf · fw · 1+bs /bt · Bt
where
• km is the ratio between maximal and average flux density, often set as 1.515
• kf is constant equal to 1.11
• fw is the winding factor, in average around 0.925
• fj is the ratio between between net and gross iron length, which is usually within 0.80-0.85
• Bt is the maximum tooth flux density when the generator is running at no load
From the armature loading, the nominal current can be calculated:
As · τs · pnr
In = (8)
2
4
The nominal voltage is then calculated as:
Sn
Un = √ (9)
3 · In
W
= (0.135 − 0.003 · U )·(1 + 0.002(∆T − 60))
cm2
(10)
where U is measured in kV and ∆T is allowed
temperature rise in Celsius. The insulation
thicknesses are decided from the nominal volt-
age. The cobber width in each slot can now be
calculated from knowing both slot width and in-
sulation thickness. The actual armature loading
can be calculated as:
W 2.1·10−6 ·Ic2
= (11) (b) Intertwining of a Roebel bar [8]
cm2 Acu,s (bs + hs − hk )
where Figure 3: Roebel bar
• Ic is the current in each bar, Ic = In /pnr,
• Acu,s is the cobber area in each bar
5
• hk is the yoke height
Iteration on the slot height is normally performed to find an appropriate value that satisfies the
armature loading limit. The current density should also be calculated to check if this is too high.
Normally this will be between 2 [A/mm2 ] and 5 [A/mm2 ].
The next objective is to find the outer diameter. First, the distribution factor Kd is calculated:
sin m·β
2
Kd = (12)
m· sin β2
where
• m is the number of slots per pole per phase
• β is the angular displacement between slots, β = 180°·Np /Qs
Next, the pitch factor Kp should be chosen to minimize the 5th and 7th harmonics. This is done
most efficiently when the coil span is around 0.80-0.85 times the pole pitch. This ratio is called
relative pole pitch and will be denoted as y:
π
Kp = cos(y· ) (13)
2
In this report, skewing of the coils will not be considered, and so the skewing factor can be omitted
when calculating the total winding factor fw :
fw = Kd ·Kp (14)
where Ns is the number of turns in series per phase and Ef is the induced phase voltage when the
generator is running at no load. If not yet decided, the net iron length ln can now be calculated
from first finding the air gap length δ between rotor and stator. An empirical formula often used is:
As
δ = γ·τp · (16)
Bδ
where γ is an empirical factor and Bδ is the maximal flux density in the air gap at no load, from
now assumed to be 0.95T. The net iron length is found from:
φm ·km
ln = (17)
0.95· τbst ·Bt ·τp
where the number 0.95 is an iron fill factor. Next, the yoke height hk can be calculated as:
φm
hk = (18)
2·0.95·ln ·Bd
where Bd is the yoke flux. Lastly, the outer diameter can be calculated:
Dy = Di + 2·(hs + hk ) (19)
6
2.1.1.3 Armature reaction and reactances
When current is flowing in the stator windings, the stator magnetic field will produce armature
reaction voltages in the stator. The armature reaction is modelled by the armature reaction reac-
tance Xad . For salient pole machines the armature reaction acts differently for the direct-axis and
quadrature-axis. This means that we get different armature reactances for the two axis. It is desir-
able to adjust the air gap in the machine so that the direct axis synchronous reactance is around
1.1pu. First, the general armature reaction reactance is found. To calculate this, the maximum
ampere turns because of armature reactance at rated load is found:
2.7 · fw · In · Ns
Fa = (20)
Np
Xd = Xad + Xl (25)
When the stator calculation are done, one can start calculating the parameters in the rotor. First,
the damper windings are designed, which are important for the stability of the generator. Figure 4
shows how the damper windings are placed in the rotor. The damper windings are shorted in both
ends so that it forms a cage. When the rotor is pulled into synchronism, voltage will be induced in
the windings, leading to a current flow that sets up a magnetic field. This field connects to the stator
field producing rotor torque. When the rotor reaches synchronous speed, no voltage will be induced
since the rotor has the same speed as the magnetic field from the stator windings. No current will
then flow in the damper windings. However, if the load is changing quickly, the rotor speed might
drop or increase momentarily. Again current in the damper windings will set up a magnetic field
7
to pull the rotor back into synchronism. To reduce harmonics, the damper winding span τr should
not be the same length as the pole pitch. Normal values for τr are 0.82 - 0.88 or 1.15 times the pole
pitch.
Further, rotor parameters as the pole shoe height and width are decided. The pole shoe height is
often set as 5cm, while the width is chosen to fit with the pole pitch.
The magnetic calculations are then done. First, the magnetic voltage drops are calculated. The
magnetic voltage drop over the air gap between stator and rotor is found from:
2Bδ
Um,δ = ·δme (26)
πµ0
where µ0 is permeability in free space. The magnetic voltage drops in slots and yoke are found by
integrating over the magnetic H-field. After finding the voltage drops, one can calculate the leakage
inductances in the machine. The total leakage inductance Lσ omitting the skew leakage inductance
is:
Lσ = Lδ + Lu + Ld + Lw (27)
where
• Lδ is the air gap leakage inductance
• Lu is the leakage inductance in the slot
• Ld is the stator tooth leakage inductance
• Lw is the leakage inductance in the end windings
The total magnetization of the machine is then calculated as:
where Um,tot is the total magnetic voltage drop and 1.1 is included to account for saturation in the
machine. After the magnetic calculations are done, the field winding can be designed. An allowed
field loading can be calculated from an empirical formula:
W πDi n
= (0.22 + 0.0055 )[1 + 0.016(∆T − 60)] (29)
cm2 60
where ∆T is the allowed temperature rise. Some parameters now have to be set:
• end plate height hkr
• insulation thickness between conductors in the field winding bif and between field windings
and the pole bi
8
• temporary field winding width bcuf
• temporary end plate thickness Le
An empirical ratio between the average and maximum winding length is given as:
Lf md
= 0.935 (30)
Lf mx
From this the winding height can be calculated:
v
u 2.1·10−6 · Lf md
u
Lf mx
hf [cm] = · ΘmN (31)
t
W
cm 2 · 0.85·b cuf
To find the number of conductors in the field winding, the formula for field voltage Vf can be
rearranged:
Vf ·(Af nf )
nf = (34)
2.1·10−6 Lf md ·Np ·ΘmN
Further, the conductor height is found:
hf − bif (nf − 1)
hcuf = (35)
nf
9
A new field winding width can now be calculated:
Af n f
bcuf = ·1.03 (36)
nf ∗ hcuf
where the factor 1.03 is added to compensate for rounded corners on the conductors. Af can now
be found:
π
Af = hcuf [bcuf − (1 − )hcuf ] (37)
4
The actual field loading is then calculated from:
W −6 Θ2mN Lf md
= 2.1·10 · (38)
cm2 (Af nf )hf Lf mx
which should be less than the allowed field loading to prevent the field winding from overheating.
Lastly the current density in the field winding is found:
ΘmN
Sf = (39)
Af n f
2.1.1.5 Losses
The machine’s geometrical parameters are now set, and the losses can now be estimated. Copper
losses are related to heat produced by the current in stator and rotor windings. In the stator
windings it should also be accounted for skin effects since AC-current is flowing here. When Rac is
the AC-resistance in the stator, the stator copper losses are found from:
The field winding only carries DC-current so no skin effect needs to be taken into account:
Next, the iron losses are calculated. For the stator yoke, the iron losses are found from:
2
Bys
PF e,y = kF e,y ·P10 · ·mF e (43)
1T
where
• kF e,y is a factor implemented to include the effect of saturation in the yoke
• P10 is specific loss for steel in [W/kg]
• Bys is yoke flux density
10
• mF e is the stator mass excluding teeth and copper
Similarly, iron losses in the stator teeth can be calculated:
2
Bts
PF e,t = kF e,t ·P10 · ·mds (44)
1T
where
• kF e,t is a factor implemented to include the effect of saturation in the teeth
• Bts is the average tooth flux density
• mds is the mass of the stator teeth
Bts can be found by calculating the average of the maximum and minimum tooth flux density.
The mechanical losses from moving parts in the machine are now found. These losses are mainly
from friction, fans and bearings. The fan and bearing losses are difficult to calculate and should be
measured for every machine. The losses can be estimated from an empirical formula:
p
Pf r ≈ kv ·D3 ·n2s · lb ·10−5 [kW ] (45)
If the generator has an external DC-generator to supply the field circuit, the magnetizing losses can
be set as 15% of the rotor copper losses. For the case with static magnetization 7% of the rotor
copper losses is sufficient.
The total losses can now be added together:
11
2.2 Optimization problems
The objective of an optimization problem is to find the best possible solution from all feasible
solutions. Since the aim often is to decrease the cost of production, most cases of optimization
problems are minimization problems, although maximization and minimization are two versions of
the same. Generally, we would write optimization problems in this form:
minimize f (x)
x
where f(x) is the function we want to minimize and x is the variable or a vector of several variables.
If a solver that finds the minimum of a problem were to be used on a maximization problem, we can
go about easily this by writing the problem as:
minimize − f (x)
x
12
2.2.1 Bounds and constraints
Often one are only interested in finding the best solution within a certain area of the solution space.
In some cases one often has an idea of where the global optima is located before the algorithm
is performed. We can then create bounds in which the algorithm will search inside of. This will
decrease the search area and consequently reduce the time needed to find the optimal solution. If
one were to find the global minima for the peaks-function in figure 6, it would be wise to create
bounds to reduce the search area. A lower bound at -3 and an upper bound at 3 for both the x-
and y-variable would work well in this case. However, it is important not to create too narrow of
bounds that the solution we are looking for is left outside the search area for the algorithm. The
lower and upper bound will later be denoted by lb and ub. A common way to specify bounds is to
treat lb and ub as vectors of length n, where n is the number of variables. One can then write our
bounds as [lb;ub].
In other cases, some of the variables may have other requirements that must be taken into account.
For example, if one of the variables in figure 6 was distance or time, then only positive values would
make sense in a physical world. To meet this condition we create a lower bound at 0 for that variable.
In some cases this will leave out the optimal solution for the unbound problem, but as we are only
interested in physically possible solutions, it will give us the best feasible solution. Solutions outside
bounds or constraints are called non-feasible solutions.
Some variables may only take integer values. For example, if a firm produces tables and chairs and
wants to find the optimal amount of units produced to maximize revenue, only integer values would
be functional. The genetic algorithm in MATLAB works well with integer value problems although
it has difficulties with simultaneously dealing with integer and equality constraints.
In more complex problems we have other constraints that must be taken into account. These can
be linear and nonlinear equalities and inequalities. The linear constraints are written as:
A · x ≤ b,
Aeq · x = beq
where A is a matrix m · n with m number of linear inequalities and n variables, b is a vector of
length m, Aeq is a matrix p · n with p number of linear equalities and beq is a vector of length p.
Often, constraints have to be implemented as less than or equal constraints. We can implement a
greater than or equal constraint by multiplying with -1 on each side of the inequality. Meaning, if
a constraint Ai · x ≥ bi is to be implemented we can rewrite the constraint as −Ai · x ≤ −bi . The
nonlinear constraints are written as follows:
ci (x) ≤ 0, i = 1, . . . , q
ceqi (x) = 0, i = q + 1, . . . , qt
where c(x) is the nonlinear inequalities, ceq(x) is the nonlinear equalities, q is the number of non-
linear inequalities and qt is the total number of nonlinear constraints.
13
We call the set of all possible solutions that satisfy all constraints the feasible region or a solution
space to a problem. Figure 7 shows the feasible region of a problem with one variable x and three
linear constraints. If we were to maximize y, the best possible solution in the feasible region would
be within the circle illustrated.
The first step to solve an optimization problem is to define the problem with its objective function,
bounds and constraints. Secondly, an optimization method should be selected based on the problem
type. Sophisticated methods might work better or even be essential on complicated problems, but are
likely to be slow on less complicated problems. For each optimization problem the best technique
is the one that finds the global optima in the shortest amount of time. In its simplest form an
optimization problem could look like this:
which would not even require an optimization method to be performed. The optimal solution is
clearly x = 1 which gives an optima of f (1) = 0. For a slightly more complex problem like:
one can find the derivative of the objective function and its zeros. Still there is no need for sophis-
ticated optimization methods. As problems get more complex with several variables and linear and
nonlinear constraints, the use of optimizing methods becomes more and more essential. Listed below
are some techniques commonly used
• Simplex algorithm
• Combinatorial optimization
14
• Newton’s method
• Local search
• Pattern search
• Evolutionary algorithms
• Simulated annealing
• Particle swarm optimization
Some of these would take too long to do by hand and has to be implemented by a computer to
reduce the computational time. Simulated annealing and evolutionary algorithms both perform
a high number of function evaluations that make them very difficult to use without a computer.
However, it makes them thorough and they are commonly used on problems with non-differentiable
objective functions.
Lastly, the optimization technique should be performed on the problem. For a great deal of these
techniques, they will not perform properly unless they are adjusted correctly to the problem. For
the genetic algorithm, this is explained in detail in section 2.3.2. Sometimes it may be necessary to
run the optimization several times as it might not give the same result every time. This is often the
case for techniques that use random or semi-random numbers in its algorithm.
15
2.3 Genetic algorithm
A genetic algorithm is one out of many techniques used to solve optimization problems. Its name
comes from Charles Darwin’s theory of natural selection, as it involves a crossover, inheritance,
selection, and mutation. Genetic algorithms are a part of Evolutionary algorithms that must also
follow the principles of: reproduction, natural selection and diversity. From a population of several
individuals or function evaluations, a genetic algorithm will pick the most fit individuals to move
onto and reproduce the next generation. In a minimization problem, the most fit individuals would
have a gene or variable composition that leads to the lowest function evaluation. Meaning, the
most fit individual would have the lowest fitness value. Because of the mutation and high number of
function evaluations in each generation, the algorithm works well on complex optimization problems.
As in every optimizing method, there is no guarantee that a genetic algorithm will find the global
optima. It is wise to run the algorithm several times as it will not give the same answer ev-
ery time. This comes from the randomly performed mutations. The high number of function
evaluations also makes the algorithm slow on simpler problems, where less sophisticated meth-
ods are better off. However, the evolutionary process sometimes leads to unimaginable solutions.
An example of an unlikely solution found
using genetic algorithms is an antenna cre-
ated by NASA for the ST5 mission in
2006, illustrated in figure 8. The mis-
sion was created to identify, develop, build,
and test innovative technologies and concepts
for use in future missions. NASA stated
the ST5-3-10 antenna was 100% compliant
with the mission antenna performance require-
ments, and was confirmed by testing the pro-
totype antenna in an anechoic test cham-
ber at NASA Goddard Space Flight Center
[9].
2.3.1 Initializing
The algorithm begins by randomly selecting a population within the bounds of the problem. It is
possible to set an initial population range in MATLAB. This will cause the initial population to
take values within a specific range. If we know before the algorithm is performed where the best
solution is likely to be, we can possibly decrease the run time by setting an initial population range
16
around this point. Figure 9 shows the peaks-function in a 2D-plot and the first generation in the
genetic algorithm. The color shading implies different function values, blue and yellow meaning low
and high values as in figure 6. The red dots are individuals that makes up the population in this
first generation. In figure 9a, the initial population bound is [-1 -3; 1 -1]. In this case, the initial
population is selected randomly within these bounds, while in case (b), the population is randomly
selected from the full solution space. The initial population bounds are selected to include the global
minima. Clearly, case (a) has more individuals closer to the global optima than case (b).
By defining suitable bounds for the initial population we tell the algorithm to evaluate points close
to where we think the optimal solution is. This will make the algorithm concentrate more on this
area than the rest of the solution space and possibly reduce the number of generations and function
evaluations performed. Caution must be used when setting bounds to ensure that the optimal
solution is included. Otherwise, the algorithm may terminate in a different place than the optimal
solution.
17
Figure 10 shows the optimization problem after
the algorithm has terminated in three different
cases. Case (a) begins with an initial popula-
tion located around the global optimum as in
figure 9a. The global optima is found and the
algorithm terminates. It is worth noting that
in case (a),areas outside the initial population
bounds are evaluated as well. The initial pop-
ulation bounds do not apply after initializing,
and so the algorithm starts looking for solu-
tions within lb and ub. However, in this case
these individuals are not likely to move on to
the next generation as there are more fit indi-
viduals closer to the global optima. (a) With initial population bounds
Case (b) has no initial population bounds. Al-
though the algorithm finds the global optima,
points from all over the solution space are evalu-
ated. Many of these points are far away from the
global optima, and the algorithm spends some
time finding the relevant area. The run times for
case (a) and (b) were approximately the same,
reason being the simplicity of this problem. De-
spite the lack of initial population bounds, the
algorithm very quickly finds where to look for
the global optima. For more complex problems,
significant time can be saved by defining initial
population bounds.
In case (c), the initial population bounds are [-3
(b) No initial population bounds
1; -1 3]. The algorithm quickly moves out of its
start point and towards a local optima. How-
ever, it never makes it to the global optima, and
instead it get stuck and terminates at the local
optima. As the algorithm finds this area of low
fitness values, it starts evaluating many points
close by, all of which have lower fitness values
than the surrounding areas. A lucky mutation
could help the algorithm find the global optima,
but as the mutation is randomly generated that
will not happen every time.
In conclusion, applying well defined initial pop-
ulation bounds can possibly decrease the run
time of an optimizing problem. It also increases
the algorithm’s chances of finding the global op- (c) Poorly placed initial population bounds
tima. It is not given that the algorithm will find
the global optima every time. Case (a) has a Figure 10: Fully performed Genetic Algorithm
greater chance of success than case (b). Lastly,
poorly placed bounds increases the chances of
the algorithm getting stuck in local optimas.
18
2.3.2 Crossover and mutation
The three different types of children are all important to bring the algorithm forward towards
the global optimum. Elite children make sure we bring the best known genotypes to the next
generation. This secures that the the algorithm never moves backwards to a less optimal solution.
Crossovers are important to investigate if a recombination of genes from parents with low fitness
value can create even lower fitness values. Mutations adds to the diversity of the population. They
can potentially create superior children, but their greatest feature is to reduce the chance of the
algorithm terminating at a local optima.
19
There are several ways to adjust the performance of the genetic algorithm. By increasing the
crossover fraction, less mutations are performed, and so we get a less diverse population. If the
crossover fraction is set to 1 there will be no mutation whatsoever, and the algorithm is thereby
unable to create new genes. Only genes from the initial population can be recombined. With initial
population bounds, a crossover fraction of 1 would disable the population from ever leaving those
bounds. On the other hand, a crossover fraction of 0 means there are only mutation children besides
the elite children. In this case, the algorithm will struggle to carry good genes forward from one
generation to the next, causing difficulty with finding a minima. The next figures are reproduced
from Mathworks’ homepage and are good illustrations of what happens when the crossover fraction
is set either too low or too high. The objective function f (x) in question is defined as:
In figure 13 the crossover fraction is set to 1. The topmost plot shows fitness values of the population
as they change from one generation to the next. Blue represents the mean fitness value of the
population while black is the best fitness value. The bottom plot shows the average distance between
individuals in each generation and so it expresses the diversity of the population. The algorithm is
also run with a crossover fraction of 0.8, although the plots for this case are not given here. With
a crossover fraction of 0.8, the best fitness when the algorithm terminated was 0.08. However, in
figure 13 the best fitness value is 0.66. The algorithm is not able to generate new genes as there are
no mutations, and the diversity of the population is quickly dropping. Already from generation 19
all the individuals in the population have the same genes. This makes the algorithm unable to find
the global minima.
20
Figure 14: Genetic algorithm when crossover fraction is 0
[10]
The case showed in Figure 14 on the other hand has a crossover fraction of 0, and so all the individuals
except the elite children are mutation children. The initial population bounds in this problem is set
to -1 and 1 for all the variables, although the problem in general is unbound. This means that after
the first generation the variables can take any value. This may lead to high fitness values, while all
the fitness values in the first generation have to be 10 or less. This explains why the mean fitness
seems to be around 5 in the first generation before it increases drastically. The initial population
bounds are chosen to make sure there is a low fitness value in the first generation to illustrate the
following: while the best individual moves on to the next generations, there are no improvements
after the first generation. The mutations only happens to the other individuals, which improve
slowly as can be seen in the top plot. Since we have no crossovers, the improved genes are never
recombined with each other to generate better individuals. Other adjustments are also implemented
in MATLAB. Among them are shrink, which controls how much the average amount of mutations
decreases from one generation to the next. An effect of this is seen in the bottom plot, where the
average distance between individuals is decreasing. By reducing the amount of mutations performed,
the diversity of the next generations decreases.
21
Several termination criteria can be set to tell the algorithm when to stop. In MATLAB, the listed
stopping criteria are implemented. All of these have default values, though they can be altered to
fit the problem at hand.
• Generations - maximum number of generations to be performed
• Time limit - maximum run time
• Fitness limit - a satisfactory fitness value
• Stall generations - maximum number of generations performed with a relative change in the
fitness value function less than Function tolerance
• Stall time limit - maximum run time without any improvement of the best fitness value
• Stall test - which can be either average or geometric weighed. This generates the stall
condition for Stall generations.
• Function tolerance - satisfactory average change in fitness values
The algorithm terminates when one of these conditions is met. However, in some cases this happens
while the best fitness value is still improving, and so the parameters should be adjusted to prevent
early terminations.
Whenever a function evaluation is performed the algorithm checks if that value is feasible. However,
MATLAB allows small violations of the constraints. The Nonlinear constrain tolerance, which
is also adjustable, sets an upper bound to violations of the nonlinear constraints. Violations of the
linear constraints have to be smaller than the square root of this tolerance for the solution to be
feasible.
22
3 Earlier work
Former students at NTNU completed research on the construction of synchronous generators and
their optimization. The work that has been done by these students provides knowledge beneficial to
solve the tasks given in this report. Two Master’s theses in particular have been used throughout
the work of this report. Their most relevant results are presented here.
In 2010, Aleksander Lundseng and Ivar Vikan completed a Master’s thesis [7] on generator design
for hydropower station upgrade. The underlying motivation for the assignment has its reason in
that Europe is facing a period of upgrading and modernization of hydropower generators. This is
explained by aging generators or changing needs. The objective of the thesis was to implement the
construction of synchronous generators in a data program. The program was to be tested against
actual machines or other calculation programs to verify its accuracy. A second part of the thesis
was to analyse alternative solutions for generators during modernization of a hydropower plant.
As a preparatory study, Lundseng and Vikan analysed how the synchronous generator is constructed
and made a collection of all the calculations that are done to finalize the machine. Their work on
the Master’s thesis culminated in a MATLAB script called GenProg that does the full calculation
for the synchronous generator. The inputs to the program are parameters the customer usually
provides the manufacturer. Amongst these are power rating and power factor, speed of rotation,
number of poles, and constraints on temperature rise and synchronous reactance. The script creates
a machine with these start parameters and sends outputs to an Excel spread-sheet. Parameters in
the output file include but are not limited to:
• voltage and current ratings
• resistances and reactances
• thermal characteristics and cooling
• magnetic parameters
• geometrical parameters
• losses
Further study concluded that GenProg gave accurate results. This was after comparing the outputs
with real life values.
The other Master’s thesis that has been an important source of information for the work in this
report is Erlend Engevik’s research on optimal design of tidal power generators from spring 2014 [11].
In this report, genetic algorithms (GA) and particle swarm optimization (PSO) were used to reduce
the cost of permanent magnet synchronous generators. While in general this study dealt with tidal
power, the understanding of how these techniques perform on generator optimization is relevant
both for tidal power and traditional hydropower generators. Engevik’s work involved comparing
GA and PSO, and hybrid systems in which GA and PSO techniques were combined separately with
gradient-based algorithms. The algorithm that performed best and gave the lowest cost was hybrid
GA. Going from a pure GA to a hybrid GA led to a cost reduction of 31.2%. Engevik concluded
that stochastic algorithms like GA and PSO in general are good at finding the optimal region of a
design space, but that a hybrid version is superior at finding the actual optima. Furthermore, the
genetic algorithm performed better than the particle swarm optimization.
23
4 Method
The aim of the work done with this thesis is to implement an optimization technique on a synchronous
generator to reduce its cost. Choosing an appropriate technique and a program to perform it were
essential to succeed in getting reliable results. Next, it was important to set up the problem so
that it reflects the real life problem in a good way. This includes choosing the objective function,
optimization variables and constraints. Lastly, to reduce the run time of the optimization, the
implementation had to be done as streamlined as possible. The second part of the thesis is to test
the optimization by changing appropriate machine parameters.
The first decision made when solving the problem in this assignment was that MATLAB was to
be used to perform the optimization. From former students Aleksander Lundseng and Ivar Vikan,
it already existed a MATLAB script that performed the calculations on a synchronous generator
when dimensioning machine-parameters were decided. This made it very convenient to continue
with MATLAB and use the script that was made. A second reason for using MATLAB was personal
experience with the program. Selecting an unfamiliar program would demand time on learning
how to use it, which would give less time to solve the assignment. Thirdly, MATLAB provides a
well-working optimization toolbox called Global Optimization Toolbox that has been used by former
students at NTNU with good results.
As described in the theory section about optimization techniques, there exist a great number of
possible optimization techniques. Many of these are implemented in the optimization toolbox in
MATLAB. MATLAB only needs information about the problem to be solved and what type of
technique to be used. Erlend Engevik compared genetic algorithms with particle swarm optimization
in his Master’s thesis [11] and found that genetic algorithms worked better when the optimal design
of tidal power generators were to be found. It is through discussions with Engevik that the genetic
algorithm was chosen as optimization technique. No research was intended to be done on comparing
different optimization techniques. This concludes the process in deciding implementation program
and optimization technique and was to a great extent determined before the assignment was handed
out.
4.2 GenProg
GenProg had in Aleksander Lundseng’s and Ivar Vikan’s research given good results when comparing
the output values with real life machines. It was decided that GenProg was to be used in the work
with this assignment because it did the design calculations thoroughly. Considerable amounts of
time and effort were spent on becoming familiar with the script. This was important to be able
to use GenProg for something it was not originally meant to do. The script is trying to design a
functional machine based on iterations and empirical formulas from a set of inputs that are extracted
from an Excel spread-sheet. Figure 15 shows a screen shot of how the required values are listed in
Excel. These are the values that GenProg needs to run correctly. If one of these is left out the script
might crash during the calculations.
24
Figure 15: Required input values, GenProg
The actual values given in the figure are strictly illustrative and can be altered to fit the machine that
the user wants to examine. Other than the required values, the Excel spread-sheet include Optional
values, Slot dimensions and Pole dimensions. These are all parameters that the user can set if he
wants, and so they are all practically optional parameters. If any of these are specified, GenProg
will treat that parameter as a constant. If not specified, the script will calculate an appropriate
value for the parameter based on iteration and empiricism.
GenProg first extracts all the specified values from Excel. Then, in the following order, the script
calculates
• inner diameter and iron length
• number of slots
• nominal voltage, armature loading and nominal current
• stator parameters
• rotor parameters
• leakage inductances
25
• losses
• necessary cooling
• moment of inertia and total weight
Lastly, GenProg writes the calculated data to an Excel file called Output and terminates.
It quickly became clear when reading through GenProg that it needed to be altered to be useful in an
optimization problem. MATLAB spends considerably amounts of time reading and writing to Excel,
which caused the run time to be around ten seconds. Keeping in mind that solving an optimization
problem involves calling on the objective function a high number of times, it was important to reduce
the run time of GenProg. The first change was therefore to comment out the writing to Excel. This
does not change the script in any way except stopping it from creating an output file. Also, the
reading from Excel was commented out. Though this creates an obvious problem. The script now
has no information about the machine in question and will crash during calculations. To solve this,
instead of reading from Excel the values were simply declared at the beginning of the script.
Secondly, two-way communication between the user and the program had to be removed. For
example, GenProg calculates possible Qs -values write these numbers to an Excel file. It then waits
for the user to type the desired number of slots. This would not work for an optimization technique
because it would be too time consuming. As discussed before, the run time has to be low because of
the high number of objective function evaluations. This was first omitted by giving the number of
slots a value at the beginning of the script. Later, as explained in section 4.5 this was again changed.
When starting to modify GenProg, the focus was on making the script as efficient as possible. After
the changes mentioned above were made, the script still gave the exact same calculated data as
before, but the run time had now reduced to well under a second. To keep GenProg unchanged, the
code was copied and pasted into a new script called Genkalk as it was clear that more changes had
to be done.
When setting up an optimization problem, it is common to start with the objective function. In
this case, the objective function to be minimized is the cost of the generator. Further, the cost can
be divided into two parts, the cost of materials and the cost of losses when the machine is running.
The material costs depend on the price of the materials [kr/kg] and the amount of materials [kg]
needed. For the synchronous generator in question, the rotor is made from steel, the stator is
made from laminated steel and the windings are made out of copper. The prices of these materials
varies with manufacturer, diameter and amount of details. In Engevik’s Master’s thesis the cost
of laminations was e 4/kg based on price of M235-50A for approximated 30 pieces of 2m diameter
winch laminations. The price of copper wire was e 11/kg based on quotation from Dahrentråd on a
5mm*2mm wire with MV insulation and the cost of steel was e 6/kg. These prices were adapted for
the work with this assignment. To get the cost of the machine in Norwegian currency, the exchange
rate between Euro and Norwegian krone was noted March 7th. 2016 at 9.3583. Both the price of
materials and the exchange rate are easily changed in the MATLAB script. Genkalk calculates the
weight of all the major components in the machine, which are the stator core, stator winding, rotor
core and rotor winding.
The second part are the prices of losses. These will vary widely with machine size and kWh-price.
However, the machine only runs at full load parts of the time, and so the the constant losses should
26
be more expencive than the load-dependent losses. An estimated price of constant losses is 30000
kroner per kW. These are the iron losses, no load losses in the rotor, the friction and windage losses
and the magnetizing losses in the magnetizing machine. The load-dependent losses are the copper
losses in rotor and stator and the additional losses. These losses can be estimated to be around
20000 kroner per kW. The following objective function was implemented in MATLAB:
Ctot = Cloaddep · Ploaddep + Cconst · Pconst + Csteel · Grotor + Clam · Gstator + Ccopper · Gwire
where
• Cloaddep is the price of load-dependent losses in kroner per kW
• Ploaddep is the load-dependent losses in kW
• Cconst is the price of constant dependent losses in kroner per kW
• Pconst is the constant losses in kW
• Csteel is the price steel in kroner per kg
• Grotor is the rotor weight in kg
• Clam is the price laminations in kroner per kg
• Gstator is the stator weight in kg
• Ccopper is the price of copper wire in kroner per kg
• Gwire is the combined wire weight in kg
To optimize the synchronous generator, some parameters have to vary while others are kept constant.
It is of interest to analyse how the cost of the machine vary with some parameters while the defining
machine parameters stay the same. The machine should for example not change its power rating
when a change is done. Other parameters that are set to be constant are the power factor, the
pole number and the synchronous speed. These are values that the customer would decide and
they to great extent define how the machine performs. We want to find the cheapest machine that
still provides what the customer expects. During conversations with Arne Nysveen and Erlend
Engevik some parameters were discussed as potential variables. These were the number of slots,
inner diameter and length, slot dimensions and pole core dimensions. The diameter and length of
the machine are both important parameters as they determine the size of the machine. With a
constant power rating it was desirable to examine how the algorithm dimensioned the machine so
that the cost was minimized. It was also interesting to look at how the slots adjusted to get the
cheapest generator. For a high number of slots, the slot width becomes narrower. This again leads
to a higher current density in the slot and possibly overheating. As an option the slot height can be
increased, which would decrease the current density again. However, this increases the number of
strands stacked on top of each other that in turn increases the AC-resistance in the stator windings.
Number of slots, slot width and slot height are closely linked together and was therefore chosen as
variables. Lastly, the pole core height and pole core width were added to the list of variables. The
pole shoe height is often set to 5cm as a well working value based on empiricism. The pole core
dimensions were therefore selected as variables instead.
Genkalk was now transformed from a script to a function by declaring the function Genkalk.
27
Figure 16: Genkalk collapsed
Here, most of the script is collapsed to illustrate the function declaration. Genkalk takes in an input
x that is a vector containing all the variables. For every function call to Genkalk with a set of
variables the function returns the machine-cost. The genetic algorithm picks a set of variables and
evaluates that set by how it performs when Genkalk sends its output back to the algorithm. It is
therefore important that when a set of variables is sent to Genkalk, the script gives one and only
one output back to the algorithm. Furthermore, Genkalk should give the same output every time it
is called on with the same input. However, the same output may be produced from more than one
set of inputs.
The variables were now examined further. First, the number of poles Qs can only take integer values.
Furthermore, it has to be divisible by 3 for the phase balance to be achieved and by 2 for symmetry.
Instead of keeping the number of slots as a variable, Qs /6 was selected as the variable and called
Qss . Qss can also only take integer values. However, it can take all positive integer values, while Qs
only could take every 6th integer number. Integer variables are easily implemented for the genetic
algorithm in MATLAB. The second variable, the inner diameter Di , is continous and was therefore
kept as it was. This was also the case for the slot width and pole core width. The iron length on the
other hand can not take all possible values. The stator is made from laminated plates with cooling
ducts in between, both with a specific width. This implies that the gross iron length has to satisfy
Equation 49:
lb = (nv + 1)bcs + nv · bv (49)
where nv is the number of cooling ducts, bsc is the length of the laminated plate and bv is the
cooling duct length. nv was selected as variable instead of the gross iron length. As for Qs , nv can
only take integer values. Similarly, the slot height and pole core height both need to provide space
for respectively the stator and rotor windings. The stator winding consists of a number of strands
stacked on top of each other as shown in figure 3a in the theory section. The slot height hs needs
to be large enough to fit the Roebel bars:
ndl was chosen as variable instead of the slot height and can only take integer values. In the same
way, the rotor winding is made up of a number of conductors as can be seen in figure 5. The
conductor height hcuf was adopted from the work of Aleksander Lundseng and Ivar Vikan that had
given a value of 1.249mm. The pole core height can be calculated as:
28
Then when the pole core height was set, a new while-loop calculated the necessary copper area and
the width of the field winding. This structure with while-loops makes the script slower and should
be avoided if possible when optimizing. To eliminate the while-loop, the field winding width was
introduced as a variable. Lastly, during testing it seemed as the synchronous reactance was limited
by the script in some way. A high Xd makes the machine less stable, but is cheaper to design.
However, it did not take as high values as expected. The synchronous reactance is closely linked
to the air gap and increases when the air gap length decreases. In Genprog the minimum air gap
length δ0 is either extracted from the excel spread-sheet as an input or calculated from empirical
formulas. During the transformation to Genkalk, the empirical formula were chosen to calculate δ0 .
This empirical approximation limited the air gap length and made it impossible for the synchronous
reactance to vary freely. As it was an option to extract δ0 from the excel spread-sheet directly, δ0
should also be able to be selected as a variable and not be calculated from these formulas. Later
testing verified that the synchronous reactance increased and varied freely when the air gap was a
variable.
In total nine variables were selected with Qss , nv , ndl and nf being integer variables:
x1 Qss
x2 Di
x3 nv
x4 bu
x= x5 = ndl
x6 bpk
x7 nf
x8 bcuf
x9 δ0
Genkalk was now tested against GenProg and gave as expected similar results. From an output
Excel file from GenProg, the eight variables were copied and sent to Genkalk as input. Genkalk
gave almost the exact same machine-design as GenProg, which was what Genkalk was intended to
do. One while-loop was left as it was in Genkalk. The iteration over necessary cooling air flow was
performed at the end of the script and did not affect the geometry of the machine in any way. For
each iteration the air flow was increased until the temperature rise for all the different nodes in the
machine were within the maximum temperature bounds. The calculations were performed to be sure
that the cooling air flow prevented the machine from overheating. As there are no costs related to
temperature rise in the machines in this assignment, it is satisfactory to know that the machine will
not overheat. Although there are costs related to losses when current is flowing through a resistance
and generates heat, these are already accounted for.
The dimensioning parameters of the generator were now set. Although the intention is to optimize
the generator, some parameters need to be set to initialize the optimization. Some of the most
important parameters are listed below. A full list can be found in the appendix.
• Sn = 50MVA
• cos φ = 0.9
• Np = 16
• ns = 375rpm
which are typical values for a medium sized hydropower generator.
29
A lower and upper bound were now set for all the variables. To be sure that the optimal solution
were within these bounds, they were set fairly wide. In GenProg the lower and upper limit of number
of slots is the circumference divided by 0.12 and 0.04 respectively. An inner diameter between 3 and
5 meters is expected for a hydropower synchronous generator at 50MVA. For a generator with an
inner diameter of 3m the bounds become 79 and 236 slots, while with a 5m inner diameter they are
131 and 393 slots. The bounds for Qs were set to 78 and 390 that led to a lower bound at 13 and an
upper bound at 65 for Qss . For the inner diameter a lower bound of 2m was decided to be sure not
to leave the optimal solution outside the bounds. The peripheral speed on the rotor increases with
an increasing diameter, and the construction becomes difficult if the speed is too high. A common
limit on the peripheral speed is 150m/s for salient pole synchronous generator. They are designed
to withstand stress at runaway speed, which differs with different turbines. For a Pelton turbine
the runaway speed is around 1.85pu while it can be from 2.0pu to 2.2pu for a Francis turbine [12],
though these numbers differ from different sources and machine sizes. A runaway speed of 1.8pu
was chosen. From the runaway speed and maximal peripheral speed, the maximum inner diameter
can be calculated:
150 · 60
Di,max = (52)
π · nr
With nr equal to 1.8pu = 675pu the maximum inner diameter becomes 4.244m. This was chosen as
the upper limit for the inner diameter. For the rest of the variables, bounds were set from trial and
error.
The next step was to include constraints associated with maximum current density, synchronous
reactance and flux density. It was clear that at least some of these constraints were nonlinear. For
the genetic algorithm in MATLAB to deal with nonlinear constraints, it has to be provided with a
function that tells the algorithm if a set of variables obeys the nonlinear constraints or not. The
function and script nonlcons was created to keep track of the nonlinear constraints. Figure 17 shows
the function declaration of nonlcons. Most of the script is collapsed.
As can be seen, the function takes in the variable x and returns an array with two components, c
and ceq . When the genetic algorithm runs, it will perform a high number of function evaluations.
For every set of variables, the algorithm also checks if it obeys the nonlinear constraints. So the
input x is the set of variables discussed in the previous section. The output c is a vector where
each component represents a nonlinear constraint. The constraints have to be based on less than
statements. For example, if the following constraints were to be implemented:
x1 x2 < 50
2x21 − 3x2 x3 > 15
x1 x2 − 50 < 0
− (2x21 − 3x2 x3 − 15) < 0
30
and represented in the script as:
c(1) = x1 x2 − 50
c(2) = −(2x21 − 3x2 x3 − 15)
If all the components in c are negative for a set of variables, the constraints are satisfied. The vector
ceq contains the nonlinear equalities. However, the genetic algorithm does not handle both integer
value variables and nonlinear equalities. Therefore, it has to be set as an empty vector in nonlcons.
The nonlinear constraints would be very complicated and chaotic if they were to be expressed only
by the inputs x1 -x9 . There are a high number of calculations done going from the inputs to values as
reactances and flux density. Furthermore, the script nonlcons does not have access to the parameters
in Genkalk. Consequently the code in Genkalk was copied to nonlcons. After deciding what constants
were needed, the excess code was removed to improve efficiency. More efficient ways of accessing
these parameters were not examined.
The constraints could now be implemented. First, the synchronous reactance Xd in a synchronous
generator should be limited to maintain static stability for the generator. A common upper limit for
the synchronous reactance for hydropower generators is 1.2pu. Additionally, the transient reactance
should not be above 0.4pu. The sub transient reactance is important for limiting the short circuit
currents and to prevent failures from spreading to other parts of the power system. A lower limit of
0.15pu was implemented for the sub transient reactance.
The current density needed to be limited to prevent overheating in both stator and rotor. GenProg
itself makes sure that sufficient cooling is present for each machine. However, run time can be
reduced by ruling out some of the machines with high current density before GenProg is called.
Only an upper limit was set for the current density in rotor and stator. Common values are from 2
to 5A/mm2 . An upper limit of 5A/mm2 was implemented for the stator and rotor current density
Ss and Sf .
To avoid too much magnetic saturation in the machine that causes high losses, constraints were
implemented in nonlcons. High flux densities also cause noise that should be limited. Norm values
give guidelines for the flux density in different parts of the machine. The magnetic flux density is
higher in the stator teeth than in the stator yoke because the area that the magnetic flux spread
over is smaller in the teeth than in the yoke. Similarly, the flux density is higher in the pole core
than in the rotor yoke. The constraints set were therefore lower for the rotor and stator yoke than
in the pole core and stator teeth. For the bottom parts of the stator teeth the maximum allowable
flux density set was 1.7T. The flux density is largest closer to the air gap at the bottom parts of
the stator teeth because of leakage flux. In the same way the flux density is higher at the bottom
parts of the pole core. Here, a constraint of 1.6T was set. For the rotor yoke the flux density is
assumed constant. A limit of 1.3T was set for the rotor yoke. The stator yoke flux density is in
Genkalk calculated from Bymx, which is the maximum yoke flux density. This parameter is set as
1.3T, and consequently the stator yoke flux density will always be 1.3T. Hence no constraint for the
stator yoke flux density is needed.
For large generators connected to the grid, the system operator often specifies restrictions to its
moment of inertia. This is to maintain enough swing mass in the power system, which ensure
system stability. Statnett, the Norwegian system operator says in FIKS 2012 [13] that the ratio
between the time constants for the moment of inertia TM and the water way Tw should be more
than 4. Tw is usually set to 1s. TM is calculated as two times the inertia constant H:
2
1 Jωm
H= (53)
2 Sn [kV A]
where J is the moment of inertia and ωm is the angular velocity in mechanical radians per second.
31
Although this is not an absolute requirement but more of a guideline, a constraint was implemented
to make sure that TM /Tw > 4.
The optimization was continuously tested in the optimization application in MATLAB. The screen
shot illustrated in figure 18 shows the application. The fitness function is provided along with the
number of variables, variable bounds, nonlinear constraint function and what variables are integer
variables. If not specified, the algorithm options are set to default values. During testing, some ma-
chine constructions that the algorithm suggested were not physically possible or the algorithm would
crash during the calculations. Occasionally, machine-designs resulted in an imaginary machine-cost,
which the genetic algorithm could not deal with, reason being that MATLAB tried calculating
the logarithm of negative numbers. It proved that the genetic algorithm would sometimes suggest
variable sets that made some geometrical parameters negative. Often, the slot width bu would be
larger than the slot pitch τs , and so the slot tooth width bd became negative. To prevent this a
nonlinear constraint saying that the slot width had to be less than 0.8 times the coil pitch was im-
plemented. Furthermore, a constraint had to be implemented to make sure that adjacent poles with
field windings were not overlapping. The parameter polklaring was calculated as the space between
the bottommost part of adjacent field windings, and a constraint saying that this parameter had to
be greater than 0.002m was implemented in nonlcons.
32
Figure 19: Spacing between adjacent poles
taumn is the clearance between adjacent poles and should also be greater than zero. It is very
unlikely that taumn is negative if polklaring is greater than zero. However, a constraint saying that
taumn have to be greater than 0.002m was implemented to be sure this will not happen. In total
12 nonlinear constraints were implemented. They are listed below:
• Xd < 1.2pu
• Xdt < 0.4pu
• Xdtt > 0.15pu
• Ss < 5A/mm2
• Sf < 5A/mm2
• Bdmax < 1.7T
• Bdrmx < 1.6T
• Byr < 1.3T
• TM
Tw >4
• bu < 0.8 · τs
• polklaring > 0.002m
• taumn > 0.002m
They all have to satisfied for the algorithm to certify a set of variables as feasible. However, small
violations of the constraints are by default allowed. The Constraint tolerance is an upper limit of
how big these violations can be. The default value of the constraint tolerance is 0.001.
Lastly, the cooling air flow had to be limited. During testing, some machines ended up with a very
high cooling need. A common limit for the cooling air flow speed vim is around 16m/s. Setting
this as a constraint was avoided because of the iterative process of calculating the cooling air flow.
The iteration would for some variable sets get stuck because the cooling need was simply too high.
Another solutions was to implement a penalty for high air flow speed values. A common way to limit
the solution space for the genetic algorithm is to penalize the fitness function whenever something
undesirable occurs. In this case, it was desirable to increase the cost of the machine when the cooling
33
air flow speed became more than 16m/s. The code section in figure 20 shows how the penalty was
implemented. It was placed in the while-loop that calculates the cooling need.
Figure 20: Code section that implements a penalty for high cooling air flow speed
If vim becomes more than 16m/s, a penalty value is calculated. The iteration will continue until
either the cooling need is meet or vim becomes more than 30m/s. The binary variable aqth1 is
then set to 1 instead of 0 that tells MATLAB to exit the while-loop. Lastly, the penalty value is
added to the objective function. The penalty value is calculated to greatly penalize the objective
function as soon as vim becomes more than 16m/s. A logarithmic function is used to adjust down
the penalization of high cooling air flow speeds. This makes the plotting of function evaluations
easier.
Furthermore, it was not enough to makes sure that Genkalk would not crash from imaginary num-
bers. nonlcons contains a lot of the same code as Genkalk and would crash similarly by the reasoning
made earlier. An if -statement was implemented saying that only when certain conditions were met,
the calculations of the constraint parameters would be done. If not, the calculations would not be
performed and the nonlinear constraints were set as not satisfied.
A fourth script cost was made to contain the information about the optimization problem and
to start the genetic algorithm. This script creates a struct called problem with all the information
MATLAB needs to perform the optimization. A struct or a structure array is a data type that groups
data were each field can contain any type of data. Besides setting the fitness function, number of
variables, bounds and constraints, the optimization options are also set here. The options are set with
the optimoptions-function in MATLAB and are collected in a field called options in the structure
problem. options is itself a structure with different fields for all the option settings. The following
code example shows how the options are specified in cost. The first statement sets the crossover
fraction to 0.6. To not cancel out option settings that may already be set, option is added to the
statement before setting the crossover fraction. The second statement tells MATLAB to plot the
best fitness values for each generation. It also makes MATLAB display the iteration. This include
but is not limited to the generation number, function call count, the best and mean fitness value
and the number of stall generations.
34
Figure 21: Examples on how to implement option settings
The options are structured in groups. Some of these are listed below, although the full list can be
found on-line at Mathworks’ homepage [14].
• plot options
• population options
• fitness scaling options
• reproduction options
• mutation options
• crossover options
• stopping criteria options
Initial population bounds were specified in the option settings to reduce the run time of the opti-
mization. The bounds were set after running the optimization a considerable number of times to be
sure that the optimal solution were within these bounds. cost is fully given in appendix B.
Besides implementing the optimization, it was desirable to test and analyse how changes in the
constant parameters affected the optimal solution. This is a good way of examining the correlation
between a parameter and the machine-design. First, the relationship between cost of losses and the
optimal solution were to be analysed. Normally the load-dependent losses are less costly than than
the constant losses. To figure out what happens to the generator design when these two are equally
expencive, a test was designed. The optimization were to be run ten times for three different cases.
In case 1, the price of constant losses was 30000 kroner per kW, while the price load-dependent
losses was 15000 kroner per kW. For case 2, the constant losses cost 30000 kroner per kW, while
the load-dependent losses cost 20000 kroner per kW. In case 3 both prices were going to be 25000
kroner per kW. The losses in these three cases were then to be compared. The test should reveal if
a change in the price changes the composition of these losses.
The second test to be performed was to change the power factor cos φ. A typical value for the power
factor is 0.8. In the test the optimization would be run ten times with the power factor at 0.8, 0.9
and finally 1.0. The requirement that synchronous generators have to be able to generate reactive
power makes the generator more expencive. Changing the power factor from 0.8 to 1.0 should result
in the optimization algorithm finding cheaper machines, Furthermore, it is desirable to examine if
the geometrical parameters changes. Lastly, the optimization were to be performed without the
constraint for moment of inertia to see if the algorithm suggested even cheaper designs. The power
factor would then be 1.0.
Thirdly, the synchronous reactance Xd is limited in nonlcons to be under 1.2pu. This is to maintain
the machine stability. Often, this requirement is not set for smaller machines. They can then be
made cheaper and smaller because the air gap can be smaller, however it also makes the machine
more unstable. It is desirable to test if the machine-design changes when this requirement changes.
Ten optimizations each were to be performed with Xd < 1.2, Xd < 2.0pu and with Xd unlimited.
35
The transient and sub transient reactance were not to be constrained in the cases of Xd < 2.0pu
and when Xd was unlimited.
For all the tests, a normal case scenario is included. The normal case is when all constraints are
included, the power factor is 0.9 and the cost of load-dependent and constant losses is respectively
20000 kroner per kW and 30000 kroner per kW.
36
5 Results
The main part of this project was to develop a MATLAB script that implements a genetic algorithm
to optimize the design of the synchronous generator. In total three scripts were made. How these
scripts were implemented was explained thoroughly in the methodology section. A closer look at how
they relate to each other is gone through here. Figure 22 shows a flowchart of how the optimization
works.
cost
Optimization
problem is defined
Initial population
MATLAB de-
New population nonlcons
c>0 fines penalty
c<0
Mutation Genkalk
Ranking of
Crossover
fitness values
Selection Stopping
Stop
of parents no criteria met? yes
37
Cost contains the information about the optimization problem and performs the genetic algorithm
on the problem. After the problem has been defined in cost, the script calls on the embedded
MATLAB-function ga. It randomly selects the initial population from within the initial population
bounds. Every individual is then checked up against the nonlinear constraints in nonlcons. If
any of components in the output c is larger than the constraint tolerance, MATLAB penalize that
individual. The individuals that satisfy the constraints are evaluated in Genkalk. All the individuals
are then ranked after their fitness value. If one of the stopping criteria explained in section 2.3.3 are
met, the algorithm terminates. Otherwise, parents are selected from the current generation. A new
population is then created through crossovers and mutations. The algorithm will continue until one
of the stopping criteria is met or the user terminates it manually.
Table 1 and 2 show the results from the first test where the prices of losses were varied. The given
values are average values from ten runs with standard deviation in parentheses.
Table 1: Losses
The losses are given in table 1. For the first case, Ploaddep is 8.7% higher than for case 2, and 13%
higher than for case 3. The standard deviation is 7.2kW for the first case, 8.3kW for case 2 and
6.5kW for the last case. This is respectively 2.1%, 2.6% and 2.1% of the average values. Similarly,
the constant losses in case 1 are 0.15% lower than in case 2 and 1.5% lower than in case 3. Case 1
has the highest total losses at 683.91kW, being 4.2% higher than in case 2 and 5.5% than in case 3.
The standard deviations are respectively 2.1%, 2.4% and 2.3% for case 1, 2 and 3. The full results
are given in appendix C.
Cloaddep Cconst Qs Di [m] lb [m] bu [mm] hs [mm] bpk [mm] hpk [mm]
15000 30000 178.2 3.70 1.26 22.4 108 363 185
(14.21) (0.08) (0.09) (2.23) (4.12) (11.1) (11.96)
20000 30000 183.6 3.61 1.42 22.0 117 345 191
(12.06) (0.12) (0.20) (0.84) (5.72) (29.0) (7.28)
25000 25000 181.8 3.67 1.36 22.4 130 355 193
(16.11) (0.12) (0.17) (1.45) (4.91) (20.4) (10.71)
Table 2 shows the inner diameter, gross iron length, slot width and height and pole core width
and height for the three cases. These values represent most of the geometrical parameters that
were selected as variables. As in the previous table, these values are also average values from ten
optimizations with standard deviation in parentheses.
38
5.2 Changing the power factor
For the second test, the power factor cos φ was varied. It was set to 0.8 in case 1, 0.9 in case 2 and
1.0 in case 3. The constraint for the moment of inertia was removed for case 4 (*), where the power
factor was 1.0.
Table 3 shows the cost of losses, laminations, steel and wire in addition to the total cost of the
machine for the four cases. In table 4, some of the major geometrical parameters are shown. In both
tables, the values given are the average values from ten optimizations with the standard deviation in
parentheses. The cost of losses is 6.8% lower for case 3 than case 1 and 5.1% lower than for case 2.
Similarly, the cost of laminations and steel are lowest for the case with cos φ equal to 1.0. However,
case 3 has the highest wire cost, being 27.2% higher than in case 1 and 13.2% higher than in case 2.
The total costs are again lowest in case 3, and at 20907 kNOK (thousand Norwegian kroner), it is
4.4% lower than for case 1 and 3.5% lower than for case 2. For the geometrical parameters the slot
height is 8.7% higher for case 3 than for case 1 while the standard deviations are 3.9% and 6.3%
respectively for case 1 and 3. The results are similar for the pole core height where in case 3 it is
229mm on average, 28% higher than i case 1. The pole core width on the other hand is largest for
case 1, 9.7% higher than in case 3. Here, the standard deviations are 16.9mm and 13.8mm for case
1 and 3, or 4.7% and 4.2% of the average values. Case 4 has the lowest total cost at 19960 kNOK
or 4.5% lower than in case 3. In fact, the loss cost and all the material costs are lower in case 4
than 3. Of the geometrical parameters, the diameter has dropped from 3.61m to 3.32m, a decrease
of 8.0%. The slot height has increased to 131.3mm, whereas both the pole core height and width
has decreased.
39
5.3 Changing the constraint for synchronous reactance
In the third test, the constraint for the synchronous reactance was changed. In case 1, the upper limit
is at its normal value at 1.2pu. For case 2, Xd is limited to 2.0pu while it is unconstrained in case
3. For case 2 and 3 the transient and sub transient reactance are not constrained. Again, the tables
show the average values of ten optimizations with the standard deviation in parentheses. In addition
to the costs, Xd and the nominal current In are included in table 5. The synchronous reactance is
for case 1 and 2 pushed towards the limit as can be seen in the table. For the unconstrained case,
Xd ends up at 4.9pu with a standard deviation of 0.2pu. This case also has the largest nominal
current at 5509A, 22.8% more than in case 1 and 8.1% more than in case 2. The costs are lowest
for case 3, except the wire cost, which for case 3 is 1785 kNOK, 9.8% higher than for case 1. The
total cost of the machine is 21673 kNOK, 20565 kNOK and 20021 kNOK for case 1, 2 and 3. Here
the standard deviations are 0.98%, 0.53% and 0.30% of the average values.
Table 6 again shows the geometrical parameters representing the selected variables. The slot width
for case 3 is 15% higher than for case 1 with the standard deviation being 0.8mm and 1.54mm or
3.6% and 6.1% of the average values for case 1 and 3 respectively. The pole core height is 14% higher
in case 3 than in case 1 and 5.3% higher than for case 2. Here the standard deviations are 7.3mm,
12.8mm and 7.3mm for case 1, 2 and 3. In percentage of the average values this equal 3.8%, 6.2%
and 3.3%.
In table 7, the magnetic voltage drop over the air gap Um,δ along with the total voltage drop Um,tot
and total magnetization of the machine ΘmN are given. Also field and stator current with their
current densities are shown. The magnetic voltage drop over the air gap reduces with an increased
synchronous reactance. In the unconstrained case Um,δ is only 27.4% of the equivalent value in case
1, while the total magnetization drops with 9.3% from case 1 to case 3. Both the current density in
the rotor and stator circuit reduces when Xd increases. However, the nominal current is increasing
at the same time.
40
Table 7: Magnetization and current
41
6 Discussion
The main purpose of this assignment was to implement an optimization technique for the syn-
chronous generator using geometrical parameters as variables. Further, it was desirable to test the
optimization by varying some major machine parameters or constraints. This was important to
examine the implemented optimization and to research if it gave reliable results. The discussion
part of this assignment is therefore going to answer both how well the results can be trusted as well
as what they mean.
There are several reasons why MATLAB was selected as optimization program. First, the script
GenProg was written in MATLAB-code and was essential for the work in this assignment. Secondly,
knowledge about MATLAB made it a convenient program to use instead of getting to know a new
program. Lastly, the work of Erlend Engevik on optimization of tidal power generators was done in
MATLAB and had shown good results. MathWorks provides sufficient information on-line on how
optimization problems can be implemented. This has been important during the work on this as-
signment. Information is easily accessible and well described. In addition to explaining syntax and
options for the implementations, more thorough information on how the optimization algorithms
works is given. This has made MATLAB a suitable choice for this assignment. Furthermore, the
optimization application makes the implementation as convenient as it can get, only demanding a
fitness function and the number of variables. As for any program the debugging can be tedious
and challenging. Although MATLAB uses embedded functions for optimization problems, knowl-
edge about what these functions do is not crucial when debugging. On the other hand, important
requirements to the objective function are that:
• all inputs give real outputs.
• the output is a number and not a vector, matrix etc..
Furthermore, when a nonlinear function is provided:
• the output vectors c and ceq need to be defined.
• c should not change size for different inputs.
• ceq have to be empty when integer value variables are present.
It is essential that Genkalk work as intended and do the machine calculation correctly. During
the Master’s theses of Aleksander Lundseng and Ivar Vikan, GenProg was tested against real life
machines and gave satisfactory results. Still, errors have been found and corrected as a consequence
of the work on this assignment. Errors can still be present in the script and this should be kept in
mind when concluding from the results. A second source of misleading results is the implementation
of variables in Genkalk. Some of the variables selected are possible to calculate from empirical
formulas. This means that there are two possible ways of selecting some of the parameters for the
same machine. However, this is also the case in GenProg where the input Excel file contains both
required and optional values. If not specified, the optional parameters are calculated from empirical
formulas. Results have shown no sign that the implementation of variables has caused invalid results.
The genetic algorithm was chosen based mostly on the work of Erlend Engevik that had shown
a superior performance to the particle swarm optimization. Both techniques are suitable for the
problem at hand, reason being the complexity of the problem with nonlinear constraints and integer
value variables. The simulations showed that the solutions to the optimization varied from each
simulation. For the normal case scenario, the best set of variable values from ten optimizations are
shown, and are plotted as a percentage of the average value in figure 23. It is clear that the best
solution for each optimization varies significantly. Variable 3 is the number of cooling ducts and
42
varies most in this plot with the highest value being 32.9% larger than the average value and 66.7%
larger than the lowest value.
Figure 23: Optimal set of variable values from ten optimization for the normal case scenario
The varying results can be explained by the population being too low, which does not allow the
algorithm to search the whole solution space. Increasing the population size will presumably give
less diversified results, but will significantly increase the run time of the optimization. Also, to
further narrow down the result spread when dealing with genetic algorithms, a hybrid solver is often
implemented. However, in this case hybrid solvers can not be used because the variable set includes
integral value variables. Although the algorithm gives varying results, the aim of the optimization
is not to find the same result every time, but to find low fitness values. By running the optimization
several times, the chances of finding the optimal solution increases. The genetic algorithm is in fact
often used to find several different solutions and then try them out in practice to see which one
actually performs the best. This was the case for NASA when using the genetic algorithm to design
antennas for the ST5-mission. The antenna in figure 8 was just one out of many designs that the
algorithm suggested. Therefore, for complex problems, varying results are not indicative of unreliable
results. However, the high standard deviation indicate that only performing ten optimizations for
each case may not be enough.
In the first test the cost of load-dependent and constant losses were varied. The test results showed
a significant decrease in the load-dependent losses when the price of the these was increased, which
was expected. Figure 24 shows load-dependent and constant losses for the three cases, case 1 having
the lowest price for load-dependent losses at 15000 kroner per kW. The values plotted are the average
values of the ten optimizations performed for each case. The standard deviation is included as error
bars. The decrease in load-dependent losses is clear, while the increase in constant losses moving
from case 1 to 3 is minimal. in fact, this increase can not be ruled out as a random variation due to
the relatively high standard deviation.
The significant decrease in load-dependent losses while the constant losses stay fairly constant can
be explained by their ability to change. The load-dependent losses consists of copper losses in stator,
full load copper losses in rotor and the additional losses. The copper losses are able to change largely
by adjusting the nominal current and the copper area in the machine. Figure 31 and 32 in appendix
C shows the full results from test 1. The DC copper losses in the stator winding decreases drastically
43
from case 1 to case 3. It is calculated as three times the DC copper losses in each phase, which in
turn is calculated as the DC resistance in each phase, RDC , times the nominal current squared. It
is clear that this component of the load-dependent losses can vary largely when the nominal current
changes.
Figure 24: load-dependent and constant losses in when changing the price of losses
In addition, the slot height increases from case 1 to case 3 as can be seen in table 2, which in turn
increases the copper area in the stator windings and by this decreases the DC-resistance. There
is also a smaller decrease in the full load copper losses in the rotor windings. This coincides well
with the increase of the pole core height, making more space for the field winding which reduces
the current density and thus the copper losses. However, the standard deviation for both the pole
core height and the full load rotor copper losses is high and indicate that more tests should be
performed to conclude that this trend is in fact representative for what happens when the price of
losses is changed. The constant losses on the other hand consists of magnetizing losses, iron losses,
no load copper losses in rotor and friction and windage losses, which do not have the same ability
to change. The magnetizing losses are small, around 5kw, and close to constant for all tests done.
The iron losses are proportional to the weight of lamination in the stator. It also varies with the
magnetic flux density in the stator, but that will stay close to constant. The iron losses increases as
the load-dependent losses become more expensive, which implies that the amount of lamination in
the machine is increasing. On the other hand, the no load copper losses in the rotor is decreasing at
the same time. These losses are proportional to the field current and follows the same trend as the
full load rotor copper losses as can be seen in figure 31. The decrease of the no load copper losses to
some degree compensate for the increase of iron losses. Thus there are two reasons why the constant
losses do not change as much as the load-dependent losses. First, the constant losses do not have
the same ability to change as the load-dependent losses. Secondly, the no load copper losses in the
rotor decreases at the same time as the iron losses increases. The friction and windage losses stay
fairly constant for all cases.
Furthermore, the total losses are larger for case 1 than for case 2 and 3. As the price of losses
changes, how the algorithm priorities to minimize the losses will also change. The change from case
1 to case 2 makes this very clear. The constant losses are equally expencive for both cases while the
44
load-dependent losses are 15000 kroner per kW in case 1 and 20000 kroner per kW in case 2. Since
the cost of losses increases from case 1 to case 2, the algorithm will try to decrease the total losses
at expense of the total cost of materials. This also explains some of the reason why the iron losses
increases from case 1 to case 2. As the loss costs increases, the algorithm tries to minimize the losses
by increasing the copper area in the machine and reducing current density. Similarly, the weight of
the stator core is allowed to increase because the cost of materials now are relatively lower compared
to the cost of losses. There is also an increase of the iron length from case 1 to case 2, which increase
the machine weight including the weight of lamination. The standard deviation for the iron length
is too large to say if this a qualitative increase. Nevertheless, as the amount of lamination increases,
the iron losses increases too.
The second test was intended to reveal the relationship between power factor and machine-design.
From the results it was clear that cheaper machines could be designed when lowering power factor
requirements, as can be seen in figure 25.
Again, the standard deviation is given by the error bars in the figure. The decrease in total cost is
related to the decrease of the cost of losses shown in figure 26. As a consequence of a higher power
factor, the total magnetization in the machine can be reduced that in turn demands a lower field
current. This leads to a lower current density and consequently reduced rotor copper losses.
From table 3, it should be noted that the cost of copper wire is rising with an increased power
factor while the cost of lamination and steel decrease. The material costs are proportional to the
weight, and so an increased cost indicate an equally increased mass. As the total magnetization in
the machine decreases, the flux density in rotor and stator will decrease. Thus, the cross section
area of rotor and stator can be reduced to lower the cost of materials. This will increase the flux
densities again, although they are limited by the nonlinear constraints. Since the machine weight is
reduced, the moment of inertia drops as well which in some cases causes the constraint TM /Tw > 4
not to be met. The algorithm therefore compensate by adding more copper to the machine. This is
an unfortunate solution as copper is expensive. For a real life problem, a better solution would be
to make the rotor ring larger.
45
Figure 26: Cost of losses when changing the power factor
A first step to improve the optimization would be to implement two prices for steel in the rotor.
One price for steel used in the poles, and a lower price for steel in the rotor ring. In case 4, this
constraint was removed. As a consequence, the algorithm reduces the diameter to design a more
compact machine as can be seen in table 4. As in case 3, the use of lamination and rotor steel
are lower than in case 1 and 2. However, copper is no longer used to compensate for the reduced
moment of inertia. Instead the machine becomes smaller and cheaper than in the other cases. The
slot height increases to keep the stator copper losses low. With a smaller diameter there is less
space for slots in the stator, so the algorithm compensate by adding copper to each slot, and at the
same time lowering the number of slots. Similarly the pole core height and width decreases with the
reduced diameter.
Lastly, the third test was constructed to examine how the requirement of the synchronous reactance
defines the machine. Testing showed that there is a significant relationship between the synchronous
reactance and the total cost of the machine as can be seen in figure 27. The synchronous reactance
is given in per unit value. Of the component of the total cost, the only one to increase when Xd
increased was the copper cost. The synchronous reactance is closely linked to the air gap between
rotor and stator. When Xd decreases, this is because the air gap becomes smaller. A smaller air
gap reduces the leakage flux and the magnetic voltage drop over the air gap, and consequently the
total magnetization of the machine drops. This is illustrated in figure 28. Similarly to test 2, the
reduced magnetization led to a lower flux density that again made the algorithm reduce the use of
lamination and steel. The earlier prediction that a smaller air gap lead to a smaller machine does
not seem to be the case here. Table 6 states that the inner diameter stays fairly constant when Xd
is varied. The iron length on the other hand reduces from 1.42m to 1.24m going from Xd equal to
1.2pu to Xd equal to 2.0pu, though the length does not have the same impact on the moment of
inertia as the diameter. Since the constraint related to the moment of inertia is active, the copper
area in rotor and stator increases and the diameter is not allowed to shrink as it was predicted to.
46
Figure 27: Total cost when the synchronous reactance increases
47
7 Conclusion
Performing genetic algorithms in MATLAB is well established because of embedded functions and the
optimization application which simplifies their implementation. For the purpose of this assignment,
the MATLAB-script GenProg was important in creating the objective function to be minimized.
The script was rewritten to meet MATLAB’s demands for objective functions, and to minimize
the optimization run time. To include requirements of the machine-design the script nonlcons
contains nonlinear constraints that corresponds to the requirements. Lastly, a final script cost was
designed to contain optimization information and options, and to begin the algorithm. Default
values were used for parameters as crossover fraction and population size. Results show that the
optimal solution varies for each simulation. It is believed that increasing the population size will
significantly reduce the result diversification. However, this will also lead to increased run time. It
is not possible to include a hybrid solver because MATLAB cannot simultaneously deal with integer
value variables and hybrid solvers. Although the results are diverse this is not necessarily indicative
of a malfunctioning optimization, rather more likely indicates that the problem is very complex.
The testing of the optimization gave results as expected. Increasing the cost of load-dependent
losses led to their decrease, however the material cost became greater. By changing the power factor
from 0.8 to 1.0 the total cost of the machine went down 4.4%, which is a result of lower required
magnetization. This was also the case when increasing the synchronous reactance. The results
show that the algorithm has the ability to adjust the machine-design and find other solutions when
parameters are altered by the user. However, more than ten optimizations should be performed for
each case in order to provide results more credibility.
The genetic algorithm is a fitting optimization technique for synchronous generators due to proven
performance on problems with high complexity and integer value variables. It is important, however,
to adjust optimization options to fit the problem at hand. Doing so reduces the run time and increases
the chances of the algorithm finding the best possible solution. In particular, the population size
must be large enough to enable the algorithm to search the whole solution space. Limiting this
size will ensure the run time is not increased too much. The user of the implemented optimization
should also be aware of some parameters and constraints that affect the machine-design greatly.
The ratio between price of losses and price of materials impacts the optimal solution to a large
degree. Increasing the price of losses leads to a significant increase of copper in the machine. The
constraints implemented limit machine-design, and can be omitted if necessary. They are chosen
to ensure reliable machines that fulfill requirements of generators in the Norwegian grid. However,
the requirement related to the moment of inertia is a guideline suggested by the Norwegian system
operator Statnett SF. This constraint to a great extent decides the machine-design, and particularly
affects the inner diameter.
To further test the validity of the optimization, suggested designs should be compared with real
life generators. This would help verify that the algorithm is able to create functional machines.
Comparisons may also reveal if there are considerably flaws to the design. Additional constraints
may then have to be implemented. Furthermore, there may still be errors in the machine calculations
performed in Genkalk. One way of testing this could be to send input values to Genkalk that represent
a real life machine. The full machine-design can be viewed in the output Excel spread-sheet, which
should then match those of the real life machine in question.
The prices of materials and losses are crucial for the machine-design. One step toward improving
optimization would be to implement several prices on materials used in different parts of the machine.
48
In particular should the price of steel used in poles and rotor ring be different. Lastly, optimization
options should be adjusted as this has not been in examined thoroughly in this assignment.
49
References
50
Appendix A, Generator specifications
51
Appendix B, MATLAB-scripts
52
Figure 30: nonlcons collapsed
53
Appendix C, Optimization results
54
Figure 32: Changing cost of losses, geometrical parameters
55
56
Figure 33: Changing power factor, costs
57
Figure 34: Changing, geometrical parameters
Figure 35: Changing constraint for synchronous reactance, costs
58
59
Figure 36: Changing constraint for synchronous reactance, geometrical parameters