Discrete Displacement Hydraulic Power Take-Off System for the Wavestar Wave Energy Converter

Abstract

The Wavestar Wave Energy Converter (WEC) is a multiple absorber concept, consisting of 20 hemisphere shaped floats attached to a single platform. The heart of the Wavestar WEC is the Power Take-Off (PTO) system, converting the wave induced motion of the floats into a steady power output to the grid. In the present work, a PTO based on a novel discrete displacement fluid power technology is explored for the Wavestar WEC. Absorption of power from the floats is performed by hydraulic cylinders, supplying power to a common fixed pressure system with accumulators for energy smoothing. The stored pressure energy is converted into electricity at a steady pace by hydraulic motors and generators. The storage, thereby, decouples the complicated process of wave power absorption from power generation. The core for enabling this PTO technology is implementing a near loss-free force control of the energy absorbing cylinders. This is achieved by using special multi-chambered cylinders, where the different chambers may be connected to the available system pressures using fast on/off valves. Resultantly, a Discrete Displacement Cylinder (DDC) is created, allowing near loss free discrete force control. This paper presents a complete PTO system for a 20 float Wavestar based on the DDC. The WEC and PTO is rigorously modeled from incident waves to the electric output to the grid. The resulting model of +600 states is simulated in different irregular seas, showing that power conversion efficiencies above 70% from input power to electrical power is achievable for all relevant sea conditions.

1. Introduction

A large group of Wave Energy Converter (WEC) concepts are based on extracting energy from ocean waves using the principle of placing buoyant bodies in the sea. As waves pass the bodies, these are forced to oscillate. This is illustrated in Figure 1, showing different version of such WECs; see [1,2] for a survey. The body movements are converted into electricity by a system called the Power Take-Off (PTO); see Figure 2a or Figure 3a.

The development of an efficient and reliable PTO system is a main challenge for WECs. The problem faced is that the movement or oscillations of bodies caused by sea waves are very slow, bidirectional and irregular. An example is given in Figure 2a, showing the resulting velocities of a 5 m diameter point absorber from the Wavestar concept. The torque, ${\tau }_{\mathrm{PTO}}$, is the load torque applied by the PTO in order to extract energy from the absorber’s pivoting motion. The graphs in Figure 2b shows the slow oscillation (≈5$\mathrm{s}$ period), resulting in a requirement of a 1 MNm bi-directional load torque for extracting an average of 27 kW. In comparison, a wind turbine operating at approximately 13 RPM would produce 1.4 MW when loading the main rotor shaft with 1 MNm. Hence, the torque density requirement for wave power is immense. Furthermore, the incoming wave power is characterized by having more than a factor of 10 between mean and peak power, as seen in Figure 2c. The peak contributes heavily to the overall production [3] and may not simply be discarded. Figure 2c also shows that the energy is grouped. Resultantly, energy storage is required to store the energy peaks and use them for maintaining production between wave groups.

Figure 1. Different embodiments of Wave Energy Converters (WECs), capturing wave energy using buoyant oscillating bodies.

Figure 1. Different embodiments of Wave Energy Converters (WECs), capturing wave energy using buoyant oscillating bodies.

Figure 2. (a) a Wavestar absorber; (b) velocities and PTO load in typical production waves; and (c), instantaneous power during two minutes production for a single absorber.

Figure 2. (a) a Wavestar absorber; (b) velocities and PTO load in typical production waves; and (c), instantaneous power during two minutes production for a single absorber.

Figure 3. (a) the definition of Power Take-Off (PTO) system; (b–d), different fluid power based transmission in PTO systems.

Figure 3. (a) the definition of Power Take-Off (PTO) system; (b–d), different fluid power based transmission in PTO systems.

To incorporate this, power absorption requires a PTO having a complex transmission with a very high gearing ratio in order to make ocean waves power a generator. Moreover, to extract adequate amounts of energy from waves, the PTO load force, $F_{\mathrm{PTO}}$, or torque, ${\tau }_{\mathrm{PTO}}$, as in Figure 2b, applied by the PTO should be controlled as a function of incident wave and body movement. Finally, to fully maximize energy extraction, the PTO should be able to operate in four quadrant mode. Four quadrant mode means that all four combination of force and velocity direction should be provided, as the PTO sometimes is required to aid the float motion (reactive control strategies). An overview of required PTO characteristics are well described in [3].

The controllability and four quadrant behavior is required to compensate for the inherit off-resonance behavior of point absorbers. Frequency-wise, point absorbers are characterized by being narrow-banded with an under-damped resonance. To correct this deficiency, the force or torque, ${\tau }_{\mathrm{PTO}}$, applied by the PTO to the absorber is controlled as a feedback of the absorber’s motion, allowing the PTO to adjust the absorber’s resonance frequency to match the wave frequency [4,5]. Optimal energy transfer from wave to an absorber is obtained when the incident wave frequency matches the resonance frequency of the absorber.

Regarding implementation of PTO systems for wave energy, several investigations have been performed on using linear generators driven directly by the absorbers’ movement. However, due to the slow linear velocities (peak linear velocity about 2 m/s), conventional permanent magnet linear generators would become very large. Typical achievable air-gap sheer-stress levels between stator and translator is about 20–25 kN/m² [6]. Resultantly, a linear generator for a Wavestar float of 5 m diameter, requiring a load force of 400 kN, would be equivalent to a generator with 16–20 m² of active air-gap surface. In [7], the weight of active magnetic material alone (copper, iron laminations, magnets and back-iron) is estimated to be 1500 kg/m². Neglecting the required support structure, the material requirements still is 24–30 tons, rendering the solution infeasible. Resultantly, effort is put into using a transmission combined with more conventional generators. Mechanical transmissions have been explored. However, these would be too massive, as a gearing ratio, which is a factor of 10 higher than a wind turbine transmission, is required, along with handling a bidirectional input.

Fluid power is a suitable technology for implementing the required transmission, as it is capable of producing high controllable forces at low velocities and easily “rectify” the bidirectional movement with compact actuators (cylinders). Unfortunately, fluid power systems are often characterized by poor efficiencies when operating at part load, which is crucial with the high ratio between peak and mean power in wave energy.

A conventional hydraulic transmission for wave energy is seen in Figure 3b, where a cylinder operates as a pump, producing a bi-directional flow, which drives a hydraulic motor. The motor adapts to the flow and rectifies the flow into a unidirectional turning of the generator. However, an optimization of such a PTO from wave-to-grid is performed in [8] for the Wavestar converter, showing an overall PTO power conversion efficiency below 65% at the optimum point, quickly dropping to 45% in smaller waves. The system in Figure 3b also has the shortcoming of not allowing energy storage/smoothing.

PTO systems with the cylinder operating as a passive pump against a steady pressure, as in Figure 3c, have been used in [9,10], where energy smoothing may be performed using hydraulic accumulators. This allows operating the hydraulic motor and generator at a fairly constant load, yielding a PTO efficiency of up to 80%. However, the cylinder is limited to providing a Coulomb-like force load, reducing the amount of absorbed energy [11].

If the active valves are used instead of the passive check valves in Figure 3c, the absorber motion may be controlled though latching to improve energy capture [12]. The latching control prolongs the natural period of the absorber motion non-linearly by locking the absorber’s movement in parts of an oscillation cycle. This may be implemented by closing the valves, blocking the cylinder and absorber motion. This may yield the same or higher energy extraction as a load force controlled as a linear damper. In [13], it is suggested that the energy capture of Figure 3c may be increased by modifying with an extra accumulator with a controllable on/off valve. By, in turn, storing and releasing energy from the accumulator, the motion amplitude of the absorber may be improved.

In Figure 3c, a PTO is shown based on a hydraulic transformer [14], capable of controlling the force of the cylinder, while having a fixed system pressure with energy smoothing accumulators. However, the part load efficiency of the hydraulic transformer is poor, as it is basically two variable displacement pump/motors back-to-back.

It would be desirable to combine the positive features of the above system characteristics, i.e., having a PTO with a common fixed pressure system with accumulators for energy smoothing, combined with an efficient force control of the cylinders. This has been investigated for the Wavestar system in [15,16], where the force control of a hydraulic cylinder is based on connecting the chambers to different fixed system pressures using an arrangement of on/off vales, as shown in Figure 4a.

In Figure 4a, a cylinder with four chambers is attached to a manifold with eight on/off valves. By either having system pressure (red) or low pressure (blue) in the chambers, 16 different pressure and chamber combinations are achievable for this particular embodiment, yielding 16 different available PTO forces, as shown in Figure 4b. This corresponds to having 16 different gears available. Thus, during a wave, the PTO force is varied discretely, as shown in Figure 4c, by shifting “gear” every half second, thereby approximating a continuous force reference. The resulting PTO force, $F_{\mathrm{PTO}}$, is generated as the sum of forces produced by the different cylinder chambers:

$$\begin{array}{c}\hfill F_{\mathrm{PTO}}=-p_{\mathrm{A}1}{A}_1+p_{\mathrm{A}2}{A}_2-p_{\mathrm{A}3}{A}_3+p_{\mathrm{A}4}{A}_4\end{array}$$

(1)

As the “gear” shifting is actually a discrete variation of the cylinder displacements using the valve system, the cylinder with manifold is referred to as a Discrete Displacement Cylinder (DDC).

Figure 4. (a) the discrete PTO system based on discrete displacement control of a multi-chambered cylinder (Discrete Displacement Cylinder (DDC)-system); (b) the forces, $F_{\mathrm{PTO}}$, the PTO may produce; and (c) how $F_{\mathrm{PTO}}$ is discretely varied using the DDC-system during a wave.

Figure 4. (a) the discrete PTO system based on discrete displacement control of a multi-chambered cylinder (Discrete Displacement Cylinder (DDC)-system); (b) the forces, $F_{\mathrm{PTO}}$, the PTO may produce; and (c) how $F_{\mathrm{PTO}}$ is discretely varied using the DDC-system during a wave.

In [15], it was shown that despite assuming infinite fast and large switching valves, a certain amount of energy will always be lost when shifting force due to the compressibility of the fluid. To assess whether a PTO system may be feasible with this unavoidable compression loss, an estimate of the PTO efficiency was calculated for the Wavestar WEC in [15]. This was performed for different wave conditions, showing that an efficiency above 90% was reachable for the DDC. In [16], it was shown that if the opening and closing time of the valves was kept less than 15 ms, the minimum loss could almost be assumed, thereby making the DDC feasible for a PTO solution.

Discrete control of a hydraulic cylinder has also been investigated for mobile hydraulics in [17], showing promising energy saving compared to, e.g., conventional load sensing systems. A PTO system with a similar approach, utilizing two asymmetric cylinders has previously been discussed in [18]. The efficiency of controlling the force of the cylinders by pressure shifting was found to be between 88% and 94%, excluding the friction of the cylinder. Currently, [19] has shown a hydraulic transmission with accumulators for power smoothing operating at 70% efficiency. A similar system as in [18] is suggested in [20]. The system of [20] is tested in a scaled version in [21], capable of applying a torque of 16 kNm. The PTO is tested in a test-rig, simulating regular waves. The overall efficiency is estimated to be from 69% to 80%. Estimates were given, as the test-rig was not fully operational.

Figure 5. (a) full-scale Wavestar C5 prototype with two floats [22]; (b) Wavestar C-concept with 20 floats total; (c) Wavestar SC-concept integration with wind turbine.

Figure 5. (a) full-scale Wavestar C5 prototype with two floats [22]; (b) Wavestar C-concept with 20 floats total; (c) Wavestar SC-concept integration with wind turbine.

The PTO system presented in this paper is based on the results in [15,16] and is a complete PTO for the 20 float Wavestar C-concept shown in Figure 5b with 5 m diameter floats. The C-concept consists of semi-submerged hemisphere-shaped floats mounted on separate arms. Each arm is mounted and rotationally supported by a common platform, which is resting on the sea bed through a number of piles. One advantage of a multiple point absorber is the increased power smoothing achieved when the waves passes through the systems. The PTO components are mainly enclosed within the main structure; only the cylinder is mounted externally. The PTO is based on a hydraulic transmission, using a fluid that is biodegradable in the marine environment. For all components inside the structure, an extra level of protection is provided, such that in the rare case of leakage, the fluid is contained inside the WEC. For the cylinder outside the machine, a sump is mounted below the cylinder, connected to automated suction lines from inside the WEC. In case of registering oil in the sump, the cylinder and float is lifted into storm protection and taken out of service until maintenance.

For storm protection and extreme seas, the Wavestar concept incorporates a jacking-system, allowing the floats and platform to be lifted out of the water. Resultantly, the PTO system and structure only have to be designed for “production” waves (which, for the C5, is 3 m significant wave height). Wavestar has in 2009 installed a two-float prototype of the 600 kW C-concept at Hanstholm at the West coast of Denmark, which has been in operation ever since [22]. The C-concept is a part of the strategy for integration of wave power into wind farms, as shown in Figure 5c, with the SC-concept (star concept). This reduces both establishment, installation and maintenance cost, while having an energy park with a higher power density and a more stable energy output.

A key property of the presented PTO concept is being scalable to future multi-megawatt Wavestar systems in Figure 5c.

2. Methods

The evaluation and presentation of the PTO is performed by first describing the layout and main components of the PTO. This description also presents the main features of PTO and explains the overall operation.

After the PTO description, the modeling of the system is attended. The modeling covers wave models, hydrodynamic-model, hydraulic system and electric system. It has not been possible to use standard component simulation software, as the system is very multi-disciplinary. Additionally, complete insight of the components models is required to avoid simulating unnecessary dynamics. Otherwise, poorly conditioned stiff systems may be obtained with slow execution times. Resultantly, modeling of all systems to the required detail level has been performed and implemented in Simulink^®. Using this approach, the model still ends up having more than 600 states, but is able to simulate about a factor of 4 slower than real time on a reasonable work station. The reasonable execution time is important, as the model presented in the paper is used for optimizing the design.

After the modeling section, the basic control of the system is presented. This consists of two parts: The Wave Power Extraction Algorithm (WPEA), which covers how a float should be loaded to extract wave power, maximizing the PTO energy production. The presented algorithm is based on reactive control methods. The second part of the control section sketches how the PTO is controlled to implement the chosen WPEA algorithm, while ensuring high efficiency and a stable power output to the grid.

Finally, simulation results for the 20 floats system in operation are given for different sea states, evaluating the performance of the PTO system. The simulations presented are of relative short duration (approximately 100) wave periods. However, based on a number of generated short waves, a selection of short wave series have been tested to give approximately the same average power extraction and PTO performance results compared to long simulation of approximately 1000 wave periods. Resultantly, these short wave series are applicable for comparing and validating design iterations.

For complete statistical background of performance evaluation, final simulations with more subsystems and more complex wave models are required. This is, however, not a part of the scope of this paper.

For the reader not interested in the detailed modeling of the PTO, the paper is constructed such that after the PTO layout section, the reader may move directly to the control section and results.

3. PTO Layout

The PTO concept is illustrated in Figure 6 for a 20 float Wavestar with 5 m diameter floats. As seen in Figure 6a, each float arm is equipped with a Discrete Displacement Cylinder (DDC)-system, consisting of a multi-chambered cylinder with integrated shifting manifold, controlling the resulting force, $F_{\mathrm{PTO}}$. These DDCs provide independent force control of the floats. All DDCs supply oil flow into a common pressure line system, consisting of three pressure lines:

• Low pressure line, $p_{\mathrm{L}}$: 10 bar–30 bar;
• Mid-(intermediate) pressure line, $p_{\mathrm{M}}$: 70 bar–170 bar;
• High pressure line, $p_{\mathrm{H}}$: 150 bar–320 bar.

These pressure lines run through the entire main tube, thereby collecting all the extracted power into common hydraulic energy storages. Hydraulic motors and generators convert the stored pressure energy to electricity. In this way, the complicated wave power absorption process and the electricity generation becomes decoupled.

3.1. Hydraulic Motors and Generators

Four main pairs of hydraulic motors and generators are placed on the high pressure line for power generation. These generators are strategically distributed across the Wavestar WEC (at float row 2, 5 and 8) to minimize pipe losses. Two parallel generators are placed at float row 5 for extra generator capacity in the center of the system during high energetic sea states. Each motor and generator pair consists of a 160 kW asynchronous generator driven by a 250 cm³/rev hydraulic motor, operating the generator at a torque of 1193 Nm at a pressure of 300 bar. The hydraulic motors are fixed displacement bent-axis motors, whose efficiency is about 95% in these operating conditions.

At each generator point, a storage battery of gas-loaded hydraulic accumulators is placed. These accumulators perform short-term energy storage for smoothing the extracted wave power. The size of the storage may be increased as desired to meet the required smoothing level.

As the pressure in the lines is only slowly varying, due to the accumulators, the torque produced by the motors is near constant, thereby driving the generators at a stable load. To adapt to the overall variation in power input from the waves during a sea state, the speed of the generators is controlled using a converter. Each set of generator has its own converter, allowing speed control of the generators, while feeding 690 V at 50 Hz to the grid. During less energetic seas, different combinations of generators and converters may be closed down, such that only the required number of generators are operating.

3.1.1. Pressure Line System and Accumulators

Power generation is only performed from the high pressure line, also having the largest storage capacity. The hydraulic motors driving the generators extract power from the high pressure line and delivers flow back to the low pressure line in a closed-circuit arrangement. The low pressure line is pressurized to avoid cavitation in the cylinders during pressure shifting in the chambers.

Figure 6. (a) an overview of the PTO concept for the Wavestar WEC in Figure 5b; (b) the DDC-system, consisting of a multi-chambered cylinder with integrated shifting manifold for discrete throttle-less force control; (c) a sketch of the shifting manifold; and (d) the illustration of the 27 available forces.

Figure 6. (a) an overview of the PTO concept for the Wavestar WEC in Figure 5b; (b) the DDC-system, consisting of a multi-chambered cylinder with integrated shifting manifold for discrete throttle-less force control; (c) a sketch of the shifting manifold; and (d) the illustration of the 27 available forces.

To filter the oil in the closed-circuit system and replenish flow leakage in the motors, a charging pump is installed, supplying flow into the low pressure line.

The mid-pressure line in the system is an intermediate pressure, which is used for increasing the efficiency of the DDCs. The compression loss associated with a pressure shift in a cylinder chamber is proportional to the square of the pressure change. Resultantly, having the option of shifting to an intermediate pressure reduces this compression loss. The extra pressure line also increases the force resolution of the DDC-system from $2^3=16$ to $3^3=27$ values.

For energy storage and stabilizing the pressure lines, the batteries of gas-loaded hydraulic piston accumulators are placed on each line. When used for energy storage, the round-trip efficiency of accumulators are ≥94% [23]. To increase utilization of the accumulators’ storage capacity, the pressure lines are allowed to vary with a minimum factor of two in pressure level. At these conditions, a 50 L piston accumulator is able to store about 23 L oil. Thus, a battery of 16 may store 375 L of high pressure oil, which may keep a 250 cm³/rev motor running at 1500 RPM for one minute. The low pressure accumulator batteries each consists of 10 × 50 L accumulators, the mid-pressure batteries of 4 × 50 L accumulators and the high pressure of 16 × 50 L.

3.1.2. DDC

The core technology for enabling the PTO concept with common pressure lines is the DDC, implementing the throttle-less force control of the cylinder. The manifold of the DDC in Figure 6a is shown in Figure 6b. The cylinder has three active chambers, where each chamber may be connected to either low, mid- or high pressure. As the piston areas of the chambers varies in size, $3^3=27$ different resulting forces, $F_{\mathrm{PTO}}$, may be generated, as illustrated in Figure 6d.

The stroke of the cylinder is 3 m. To perform the shifting of chamber pressures, a shifting manifold, as in Figure 6c, with nine fast on/off valves is used (less than 15 ms opening and closing time). The pressure drop across the valves is low, i.e., less than 3bar, at flows corresponding to a cylinder velocity of 0.5 m/s.

The manifold is directly mounted on the cylinder to reduce volume (compression losses) and to reduce transmission phenomena in the pipe connections. The DDC-system shifts force approximately every 400 ms. To supply/consume the short flow pulses required for compression or decompression of the volumes, small accumulators are mounted directly on the DDC-manifold.

An important property of the DDC is to be able to process extreme power levels with low losses, as the cylinder chambers are directly connected to the pressure lines when in a given force configuration. Only the compression losses are constant, which, thus, become more significant at lower power levels. However, efficiency is maintained by reducing the overall pressure levels, as the high DDC forces are not required in these conditions.

Another important feature of the PTO concept with the DDC is to always be able to reduce the load force independent of the float velocity and position. In this way, the WEC may always reduce the power absorption of the floats if full load capacity is reached, similar to wind turbine pitching out of the wind.

Flow to the intermediate pressure line is only supplied and extracted by the cylinders to reduce cost. Thus, the control of the floating intermediate pressure is performed collectively by all DDCs, making sure that the net-flow into the mid-pressure system is zero over time. This means that the DDC-modules actually operate as a sort of loss-less pressure amplifier. For example, the DDC first lets the float motion supply flow into the mid-pressure line. An instant later, the DDC shifts configuration, such that the mid-pressure and float motion are now combined to produce a high pressure flow. An example of such a cycle is illustrated in Figure 7a: The cylinder is supplying fluid into the mid-pressure at time instance (1); At instance (2) the DDC increases the load force, $F_{\mathrm{PTO}}$, using a configuration where the second chamber is now pumping fluid into the high pressure line, while the first chamber now consumes flow from the mid-pressure. In instance (3), the force is further increased, such that the cylinder is pumping flow into high pressure without aid from the mid-pressure.

The different configurations of pressure and chambers and the resulting cylinder forces are seen in Figure 7b.

Figure 7. (a) an example of how pressure and chamber connection may change during a wave; (b) the available cylinder forces. The numbers 1, 2 and 3 relate to which pressure is in the different cylinder areas in (a) (1 = $p_{\mathrm{L}}$, 2 = $p_{\mathrm{M}}$, 3 = $p_{\mathrm{H}}$).

Figure 7. (a) an example of how pressure and chamber connection may change during a wave; (b) the available cylinder forces. The numbers 1, 2 and 3 relate to which pressure is in the different cylinder areas in (a) (1 = $p_{\mathrm{L}}$, 2 = $p_{\mathrm{M}}$, 3 = $p_{\mathrm{H}}$).

Considering energy extraction, the DCC approximates a continuous reference with discrete values. This may affect the amount of extracted energy. However, comparing [15] (discrete control) and [8] (continuous control), the power extraction appears unaffected. The reason for the low impact of the discrete implementation is that the 35 tons of float and arm serves as an effective low pass filter, having a break frequency around 0.28 Hz. Thus, from a load control point of view, implementation of extra low pass filters, as performed in, e.g., digital power electronics, is not required.

This concludes the overall PTO description, and the PTO modeling section follows, which also further specifies and describes the PTO. To skip the modeling section and move directly to control of the PTO, see Section 7.

4. Modeling

An overview of the modeled sub-systems is seen in Figure 8. These are treated in the shown order. The derived equations are implemented in Matlab/^®Simulink^® and solved with a Runge-Kutta solver, running at a fixed step time of 0.5 ms.

Figure 8. Modeled sub-systems of the Wavestar WEC and PTO system.

Figure 8. Modeled sub-systems of the Wavestar WEC and PTO system.

5. Wave and Float Model

An irregular wave is often characterized by two quantities, the significant wave height, $H_{\mathrm{m},0}$, and the peak wave period, $T_{\mathrm{p}}$. The significant wave height, $H_{\mathrm{m},0}$, is the average wave height of the one-third highest waves, and $T_{\mathrm{p}}$ is the wave period, where most energy is concentrated. These two parameters along with an underlying assumption of the shape of the power density spectrum, $S\left(f\right)$, defines a sea state. In this study, the Pierson-Moskowitz (PM) spectrum is utilized [24], describing a fully developed sea state:

$$\begin{array}{c}\hfill S_{\mathrm{PM}}\left(f\right)=\frac{A}{f^5}{\mathrm{e}}^{\frac{-B}{f^4}}B=\frac{5}{4}\frac{1}{T_{\mathrm{p}}^4}A=\frac{1}{4}B{H}_{\mathrm{m},0}^2\end{array}$$

(2)

If the PTO may control the phase of the absorber, the PM spectrum may also be the conservative choice in regards to power extraction. The PTO often tunes the absorber resonance to match the peak period; however, the power in the PM spectrum is spread across a wide band of frequencies. In comparison, an often used spectrum, the JONSWAP (Joint North Sea Wave Project) spectrum, has the energy concentrated in a more narrow frequency band (peak enhancement factor larger than one). The PM spectrum in Equation (2) is shown in Figure 9a. The figure also shows the three sea states, which are going to be used in the simulation.

Figure 9. Wave spectra for three sea states and an example of a realization of a wave.

Figure 9. Wave spectra for three sea states and an example of a realization of a wave.

To perform time simulation, a time realization of a wave is required, complying with the spectrum. Often, the random phase method is used, where the spectrum is converted into a finite number of regular wave components, which are added with a random phase for each component. This simple method has been experienced to be adequate for estimating average power absorption [25]. However, the method does not reproduce, e.g., wave grouping correctly, according to [26]. Studies indicate that using the method, the mean length of wave groups is too short [26]. Thus, the random phase generated waves may not accurately test the PTO’s capability of, e.g., power smoothing and absorption, as this is highly dependent on the degree of wave grouping. Instead, [26] suggests that by filtering Gaussian white noise, better representation of ocean waves may be achieved, and arbitrarily long series may be generated.

The white noise method is to design a filter according to the PM-spectrum, whereby the filter shapes the power density of the white noise signal to comply with the PM-spectrum. The filers are implemented based on [11,27].

If the input process to a filter, $x\left(t\right)$, is white noise, which is defined by having a flat spectrum, $S_{\mathrm{xx}}\left(f\right)={\sigma }_{\mathrm{xx}}^2$, and the filter, $H\left(f\right)$, is designed according to the PM-spectrum in Figure 9a, $H\left(f\right)=\sqrt{S_{\mathrm{PM}}\left(f\right)}$, then:

$$\begin{array}{c}\hfill S_{\mathrm{yy}}\left(f\right)={\sigma }_{\mathrm{xx}}^2·{|\sqrt{S_{\mathrm{PM}}\left(f\right)}|}^2=S_{\mathrm{PM}}\left(f\right)\end{array}$$

(3)

where a white noise with variance, ${\sigma }_{\mathrm{xx}}^2=1$, has been used as input. Thus, by filtering white noise, $x\left(t\right)$, with the filter, $H\left(f\right)=\sqrt{S_{\mathrm{PM}}\left(f\right)}$, yields a time series, $y\left(y\right)$, having a spectrum, $S_{\mathrm{yy}}\left(f\right)$, agreeing with the PM spectrum.

The generated wave with this method is shown in Figure 9b. A spectrum analysis of the generated wave is seen in Figure 9c, showing fluctuation around the target Pierson-Moskowitz-spectrum. This is desired, as spectrum measurements of real sea waves also have a non-smooth spectrum.

5.1. Absorber Model

The Wavestar absorbers are modeled using linear wave theory, as described in, e.g., [24], yielding an adequate description in the conditions in which a WEC produces energy. The absorber is illustrated in Figure 10a, representing a single degree of freedom system expressed by its angular motion, ${\theta }_{\mathrm{arm}}$. The angular position, ${\theta }_{\mathrm{arm}}$, is defined to be zero when the float is horizontal, i.e., position at calm water. A derivation of the model for single absorber may be found in [28]. The model here is briefly restated for a single absorber with a diameter of 5 m.

Figure 10. (a) definition of variables of the Wavestar absorber; (b) some dimensions for the case study; and (c) the moment arm of the cylinder.

Figure 10. (a) definition of variables of the Wavestar absorber; (b) some dimensions for the case study; and (c) the moment arm of the cylinder.

5.1.1. Single Absorber Model

In linear wave theory, the equation of motion of the absorber is obtained by superimposing the different effects of the wave-absorber interaction:

$$\begin{array}{cc}\hfill J_{\mathrm{mech}}{\ddot{\theta }}_{\mathrm{arm}}& ={\tau }_{\mathrm{Arch}}-{\tau }_{\mathrm{G}}+{\tau }_{\mathrm{rad}}-{\tau }_{\mathrm{PTO}}+{\tau }_{\mathrm{ext}}\hfill \end{array}$$

(4)

where the torque, ${\tau }_{\mathrm{Arch}}$, is due to Archimedes force; the torque, ${\tau }_{\mathrm{G}}$, is caused by gravity; the torque, ${\tau }_{\mathrm{rad}}$, is the torque experienced on the absorber from radiating waves; ${\tau }_{\mathrm{PTO}}$ is the torque applied by the PTO system; and ${\tau }_{\mathrm{ext}}$ is the wave excitation torque. The moment of inertia term, $J_{\mathrm{mech}}$, is the total (mechanical) moment inertia in the ${\theta }_{\mathrm{arm}}$ rotation.

By identifying the terms in Equation Equation (4), the following well-known expression is obtained:

$$\begin{array}{cc}\hfill {\ddot{\theta }}_{\mathrm{arm}}& =\frac{1}{J_{\mathrm{mech}}+J_{\mathrm{add},\infty }}\left(-k_{\mathrm{res}}{\theta }_{\mathrm{arm}}-k_{\mathrm{r}}*{\dot{\theta }}_{\mathrm{arm}}-{\tau }_{\mathrm{PTO}}+{\tau }_{\mathrm{ext}}\right)\hfill \end{array}$$

(5)

where $k_{\mathrm{res}}$ is the spring coefficient representing the hydro-static restoring torque; ${\tau }_{\mathrm{res}}$, which is a linearization of the Archimedes force; ${\tau }_{\mathrm{Arch}}$, and gravity, ${\tau }_{\mathrm{g}}$. In Equation (4), the radiation term, ${\tau }_{\mathrm{rad}}$, may be expressed as:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{rad}}=-J_{\mathrm{add},\infty }{\ddot{\theta }}_{\mathrm{arm}}-{\int }_{-\infty }^{\infty }{k}_{\mathrm{r}}(t-\tau ){\dot{\theta }}_{\mathrm{arm}}\left(\tau \right)\mathrm{d}\tau \end{array}$$

(6)

where the function, $k_{\mathrm{r}}\left(t\right)$, is the radiation-force impulse-response function and $J_{\mathrm{add},\infty }$ is the added inertia at infinite high frequency.

The wave excitation torque, ${\tau }_{\mathrm{ext}}$, is the torque an incoming wave applies to a float held fixed, including the contribution from the float diffracting the wave. A filter for calculating the wave torque based on the wave, ${\eta }_{\mathrm{w}}$, may be obtained given as an impulse response function, $h_{\mathrm{ext}}\left(t\right)$. The exciting wave torque may then be found by the convolution:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{ext}}={\int }_0^{\infty }{h}_{\mathrm{ext}}(t-\tau ){\eta }_{\mathrm{w}}\left(\tau \right)\mathrm{d}\tau \end{array}$$

(7)

The impulse response, $h_{\mathrm{ext}}$, is non-zero for $t<0$, rendering the filter non-causal. This mean that the current excitation force depends on the future incident waves. This is partly due to the fact that the waves hit the float before reaching the reference point in the center of the float. Furthermore, the defined wave height is not the direct cause to the wave excitation torque, but is a quantity defined model-wise. Thus, both are caused by some unknown process, and therefore, their relation is not forced to be causal.

To solve this in the simulation, the wave excitation torque is pre-computed before simulation by discrete convolution:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{ext}}(n·T_{\mathrm{s}})=T_{\mathrm{s}}\sum _{k=-n_{-}}^{k=n_{+}}{\eta }_{\mathrm{w}}\left(k{T}_{\mathrm{s}}\right)h_{\mathrm{ext}}\left((n-k)T_{\mathrm{s}}\right)\end{array}$$

(8)

where the impulse response is approximated with a finite impulse response, requiring $n_{-}$ samples of future wave knowledge. As the excitation torque is assumed independent of the float motion, but only dependent on the incident wave height, the pre-computation of ${\tau }_{\mathrm{ext}}$ is allowed. As the entire wave signal, ${\eta }_{\mathrm{w}}$, is available beforehand, the non-causality is not a problem.

5.1.2. Multi-Absorber Model

In an array of closely spaced absorbers, such as the Wavestar, the absorbers interact. The absorbers interact in two ways:

• By diffracting the incident wave;
• By radiating waves, interacting with neighboring absorbers.

The diffraction part may be obtained by placing 20 floats in the given array and, then, identifying the force filter, $h_{\mathrm{ext}}$, for each float. The resulting diffraction effect experienced by the individual floats will be included in their force filter for calculating the exciting wave torque for each float:

$$\begin{array}{c}\hfill {\mathit{\tau }}_{\mathrm{ext}}={\int }_{-\infty }^{\infty }{\mathit{h}}_{\mathrm{ext}}(t-\tau ){\eta }_{\mathrm{w}}\left(\tau \right)\mathbf{\text{d}}\tau ,{\mathit{\tau }}_{\mathrm{ext}}={[{\tau }_{\mathrm{ext},1}\cdots {\tau }_{\mathrm{ext},20}]}^T,{\mathit{h}}_{\mathrm{ext}}\left(t\right)={[h_{\mathrm{ext},1}\left(t\right)\cdots h_{\mathrm{ext},20}\left(t\right)]}^T\end{array}$$

(9)

Considering the radiation effect, the equation of motion for each float, before given as Equation (5), now becomes coupled:

$$\begin{array}{cc}\hfill {\ddot{\mathit{\theta }}}_{\mathrm{arm}}& =\frac{1}{J_{\mathrm{mech}}+J_{\mathrm{add},\infty }}\left(-{\mathit{K}}_{\mathrm{res}}{\mathit{\theta }}_{\mathrm{arm}}-{\int }_0^{\infty }{\mathit{K}}_{\mathbf{\text{r}}}(t-\tau ){\dot{\mathit{\theta }}}_{\mathrm{arm}}\left(\tau \right)\mathbf{\text{d}}\tau -{\mathit{\tau }}_{\mathrm{PTO}}+{\mathit{\tau }}_{\mathrm{ext}}\right)\hfill \end{array}$$

(10)

as the matrix, ${\mathit{K}}_{\mathbf{\text{r}}}\left(t\right)$, containing the radiation impulse responses is non-diagonal:

$$\begin{array}{c}\hfill {\mathit{\omega }}_{\mathrm{arm}}=\left[\begin{array}{c}{\omega }_{\mathrm{arm},1}\\ {\omega }_{\mathrm{arm},2}\\ ⋮\\ {\omega }_{\mathrm{arm},20}\end{array}\right],{\mathit{K}}_{\mathrm{res}}=\left[\begin{array}{cccc}{k}_{\mathrm{res}}& & & 0\\ & k_{\mathrm{res}}& & \\ & & \ddots & \\ 0& & & k_{\mathrm{res}}\end{array}\right],{\mathit{K}}_{\mathrm{r}}\left(t\right)=\left[\begin{array}{cccc}{k}_{\mathrm{r}1,1}\left(t\right)& k_{\mathrm{r}1,2}\left(t\right)& \cdots & k_{\mathrm{r}1,20}\left(t\right)\\ k_{\mathrm{r}2,1}\left(t\right)& k_{\mathrm{r}2,2}\left(t\right)& \cdots & k_{\mathrm{r}2,20}\left(t\right)\\ ⋮& & \ddots & ⋮\\ k_{\mathrm{r}20,1}\left(t\right)& k_{\mathrm{r}20,2}\left(t\right)& \cdots & k_{\mathrm{r}20,20}\left(t\right)\end{array}\right]\end{array}$$

(11)

e.g., the impulse response, $k_{\mathrm{r}2,20}\left(t\right)$, describes the cross-coupling from the velocity of float number 20 to force on float number 2. Note that $k_{\mathrm{r}\mathrm{x},\mathrm{y}}\left(t\right)==k_{\mathrm{r}\mathrm{y},\mathrm{x}}\left(t\right)$, rendering ${\mathit{K}}_{\mathbf{\text{r}}}\left(t\right)$ symmetric.

With the focus on validating the PTO, the radiation cross terms are neglected in this paper, i.e., $k_{\mathrm{r}\mathrm{x},\mathrm{y}}=0$ for ($x\ne y$). Inclusion of the cross terms mostly affects how the absorbers should be controlled to further increase power absorption. Thus, the influence of the float dynamics on the PTO and the cylinder control may be modeled with acceptable accuracy by using a single absorber hydrodynamic model. However, it is still required to have the effect of the floats operating out of phase due to their distributed locations. This is obtained by applying individual force filters, as mentioned previously.

5.1.3. Model Parameters

To obtain the hydro-dynamic model parameters, WAMIT has been applied, which is a program for computing wave loads and motions of structures in waves.

The force excitation filters are shown in Figure 11 for a wave direction of ${\theta }_{\mathrm{w}}=0^{\circ }$ and ${\theta }_{\mathrm{w}}=-{45}^{\circ }$. The wave “measurement” point is shown in Figure 11c, along with a definition of float numbering and wave angle, ${\theta }_{\mathrm{w}}$. The filters show that for ${\theta }_{\mathrm{w}}=0^{\circ }$, the wave reaches the floats farthest back (numbers 10 and 20), approximately 7 s later after exciting float number 1. From the plot, it is also seen that there is some amount of shadow effect, as the magnitude of the wave excitation filter reduces for each float. For the ${\theta }_{\mathrm{w}}=0^{\circ }$, the input for the float is pair-wise equal, whereas this is not true for the wave angle, ${\theta }_{\mathrm{w}}=-{45}^{\circ }$.

Note that the excitation torque is independent of the float movement, i.e., the excitation torque may be calculated “off-line” by the convolution in Equation (7).

From WAMIT, all the impulse-response functions, $k_{\mathrm{r}\mathrm{x},\mathrm{y}}\left(t\right)$, describing the radiation-force are also obtained, including the cross terms ($x\ne y$). However, as theses are neglected, only the single absorber radiation impulse response is used; see Figure 12a. To avoid performing the convolution term, ${\tau }_{\mathrm{rad},\mathrm{kr}}=k_{\mathrm{r}}\left(t\right)\times {\omega }_{\mathrm{arm}}\left(t\right)$, the term is approximated as a system of Ordinary Differential Equations (ODE). This is performed using Prony’s method [28]. The used fifth order approximation is shown in Figure 12a.

Figure 11. The wave excitation force filter impulse responses, $h_{\mathrm{ext}}$, for different floats. (a) is for wave direction, ${\eta }_{\mathrm{w}}=0^{\circ }$; (b) is for wave direction, ${\eta }_{\mathrm{w}}=-{45}^{\circ }$; and (c), the wave “measurement” point, float numbering and wave angle, ${\theta }_{\mathrm{w}}$, are defined.

Figure 11. The wave excitation force filter impulse responses, $h_{\mathrm{ext}}$, for different floats. (a) is for wave direction, ${\eta }_{\mathrm{w}}=0^{\circ }$; (b) is for wave direction, ${\eta }_{\mathrm{w}}=-{45}^{\circ }$; and (c), the wave “measurement” point, float numbering and wave angle, ${\theta }_{\mathrm{w}}$, are defined.

Figure 12. (a) the impulse response function, $k_{\mathrm{r}}\left(t\right)$, being approximated; and (b) Bode diagram of Equation (12), ${\theta }_{\mathrm{arm}}\left(s\right)/{\tau }_{\mathrm{ext}}\left(s\right)$.

Figure 12. (a) the impulse response function, $k_{\mathrm{r}}\left(t\right)$, being approximated; and (b) Bode diagram of Equation (12), ${\theta }_{\mathrm{arm}}\left(s\right)/{\tau }_{\mathrm{ext}}\left(s\right)$.

The absorber dynamics may now be written as the following transfer function:

$$\begin{array}{c}\hfill \frac{{\theta }_{\mathrm{arm}}\left(s\right)}{{\tau }_{\mathrm{ext}}\left(s\right)-{\tau }_{\mathrm{PTO}}\left(s\right)}=\frac{1}{(J_{\mathrm{mech}}+J_{\mathrm{add},\infty })s^2+K_{\mathrm{r}}\left(s\right)·s+k_{\mathrm{res}}}\end{array}$$

(12)

A Bode-diagram of Equation (12) is given in Figure 12b with the parameters in Table 1. The Bode diagram shows that the system has a resonance peak at $f_{\mathrm{r}}=0.285\mathrm{Hz}$, corresponding to a period of 3.51 s. Hence, the Wavestar point absorber will have optimal absorption at a wave period of 3.51 s. At Wavestar C5 production sites, this period corresponds to the sea states with the shortest peak-period, in which a substantial part of the yearly wave energy is concentrated.

Table 1. Parameter values for the Wavestar C5 absorber.

Parameter	Symbol	Value	Unit
Inertia of arm and float (with ballast water)	$J_{\mathrm{mech}}$	2.45 $\times {10}^6$	[kgm2]
Hydrostatic restoring torque coefficient	$k_{\mathrm{res}}$	14.0 $\times {10}^6$	[Nm/rad]
Added-inertia $J_{\mathrm{add}}\left(\omega \right)$ for $\omega \to \infty$	$J_{\mathrm{add},\infty }$	1.32 $\times {10}^6$	[kgm2]
Transfer-function coefficients for $K_{\mathrm{r}}\left(s\right)$:
$\left(b_0,b_1,\cdots ,b_5\right)=\left(0.0001,0.0144,0.624,8.16,13.1,1.44\right)\times {10}^6$
$\left(a_0,a_1,\cdots ,a_5\right)=\left(0.0010,0.0906,1.67,6.31,13.3,9.18\right)$

Table 1. Parameter values for the Wavestar C5 absorber.

Parameter	Symbol	Value	Unit
Inertia of arm and float (with ballast water)	$J_{\mathrm{mech}}$	2.45 $\times {10}^6$	[kgm2]
Hydrostatic restoring torque coefficient	$k_{\mathrm{res}}$	14.0 $\times {10}^6$	[Nm/rad]
Added-inertia $J_{\mathrm{add}}\left(\omega \right)$ for $\omega \to \infty$	$J_{\mathrm{add},\infty }$	1.32 $\times {10}^6$	[kgm2]
Transfer-function coefficients for $K_{\mathrm{r}}\left(s\right)$:
$\left(b_0,b_1,\cdots ,b_5\right)=\left(0.0001,0.0144,0.624,8.16,13.1,1.44\right)\times {10}^6$
$\left(a_0,a_1,\cdots ,a_5\right)=\left(0.0010,0.0906,1.67,6.31,13.3,9.18\right)$

6. Absorption System

The PTO concept is modeled as displayed in Figure 8. First, the DDC-module is developed, which connects to the pressure line system. Secondly, the pressure line system is described, which consists of transmission pipe lines and accumulator batteries.

6.1. Discrete Displacement Cylinder—DDC

The DDC-module is shown in Figure 13a, consisting of a cylinder, a manifold with valves, accumulators and pipe-lines.

Figure 13. (a) diagram of DDC-system, consisting of cylinder, manifold with nine on/off valves, accumulators and pipe lines. Fitting are indicated, where ξ denotes the fitting resistance coefficient; (b) illustration of hose model; and (c) illustration for combined valve and hose model.

Figure 13. (a) diagram of DDC-system, consisting of cylinder, manifold with nine on/off valves, accumulators and pipe lines. Fitting are indicated, where ξ denotes the fitting resistance coefficient; (b) illustration of hose model; and (c) illustration for combined valve and hose model.

6.1.1. Cylinder

The pressure dynamics of the cylinder chambers are described by the flow continuity equation:

$$\begin{array}{cc}\hfill {\dot{p}}_{\mathrm{A}1}& =\frac{{\beta }_{\mathrm{eff}}\left(p_{\mathrm{A}1}\right)}{A_1(x_{\mathrm{c},\max }-x_{\mathrm{c}})+V_{0,\mathrm{A}1}}\left(Q_{\mathrm{A}1}+v_{\mathrm{c}}{A}_1\right)\hfill \end{array}$$

(13)

$$\begin{array}{cc}\hfill {\dot{p}}_{\mathrm{A}2}& =\frac{{\beta }_{\mathrm{eff}}\left(p_{\mathrm{A}2}\right)}{A_2{x}_{\mathrm{c}}+V_{0,\mathrm{A}2}}\left(Q_{\mathrm{A}2}-v_{\mathrm{c}}{A}_2\right)\hfill \end{array}$$

(14)

$$\begin{array}{cc}\hfill {\dot{p}}_{\mathrm{A}3}& =\frac{{\beta }_{\mathrm{eff}}\left(p_{\mathrm{A}3}\right)}{A_3(x_{\mathrm{c},\max }-x_{\mathrm{c}})+V_{0,\mathrm{A}3}}\left(Q_{\mathrm{A}3}-v_{\mathrm{c}}{A}_3\right)\hfill \end{array}$$

(15)

where $V_{0,1}$ and $V_{0,2}$ are the volumes of the connecting hoses and ${\beta }_{\mathrm{eff}}$ is the effective bulk modulus of the fluid as a function of pressure. The used value of the effective bulk modulus is shown in Figure 14b. The cylinder force is calculated as:

$$\begin{array}{c}\hfill F_{\mathrm{c}}=-p_{\mathrm{A}1}A1+p_{\mathrm{A}2}A2-p_{\mathrm{A}3}A3-F_{\mathrm{fric},\mathrm{c}};F_{\mathrm{fric},\mathrm{c}}=\left\{\begin{array}{ccc}& \tanh \left(a{v}_{\mathrm{c}}\right)\left|F_{\mathrm{c}}\right|(1-{\eta }_{\mathrm{c}})\hfill & \hfill ;F_{\mathrm{c}}{v}_{\mathrm{c}}>0\\ & \tanh \left(a{v}_{\mathrm{c}}\right)\left|F_{\mathrm{c}}\right|\left(\frac{1}{{\eta }_{\mathrm{c}}}-1\right)\hfill & \hfill ;F_{\mathrm{c}}{v}_{\mathrm{c}}\le 0\end{array}\right.\end{array}$$

(16)

where the cylinder friction, $F_{\mathrm{fric},\mathrm{c}}$, is defined, such that the cylinder has a constant efficiency of ${\eta }_{\mathrm{c}}$. The function, $\tanh \left(a{v}_{\mathrm{c}}\right)$, is used instead of $\mathrm{sign}\left(v_{\mathrm{c}}\right)$ to avoid discontinuity, where a adjusts the steepness around zero velocity.

Figure 14. (a) modeling a part of the pressure line system; and (b) the used bulk modulus model.

Figure 14. (a) modeling a part of the pressure line system; and (b) the used bulk modulus model.

The moment arm, $d_{\mathrm{a}}$, at which the cylinder, $F_{\mathrm{c}}$, applies ${\tau }_{\mathrm{PTO}}$ is found from Figure 10c to be:

$$\begin{array}{cc}\hfill d_{\mathrm{a}}& =\frac{a_2{b}_2\sin ({\theta }_{\mathrm{a}}-{\alpha }_0)}{(x_{\mathrm{c}}+c_2)}\hfill \end{array}$$

(17)

where:

$$\begin{array}{c}\hfill x_{\mathrm{c}}=-c_2+\sqrt{-2a_2{b}_2\cos \left({\theta }_{\mathrm{arm}}-{\alpha }_0\right)+(a_2^2+b_2^2)}\end{array}$$

(18)

The function is shown in Figure 10c.

6.1.2. Transmission Line Model

Considering the pipe connection between manifold and cylinder, the fast on/off valves that are utilized may produce steep pressure transients, exciting the transmission line dynamics. In other terms, the time required to accelerate the mass of the fluid within the hoses may not be neglected.

The applied pipe-line or transmission line model is to view the fluid within the line as a single mass column, which is accelerated by the pressure difference at the line inlet and outlet; see Figure 13b. The dynamics of the line flow velocity, $v_{\mathrm{l}}$, may then be described using the momentum equation:

$$\begin{array}{c}\hfill {\dot{v}}_{\mathrm{l}}=\frac{(p_{\mathrm{l}2}-p_{\mathrm{l}1})A_{\mathrm{l}}-p_{\mathrm{f}}{A}_{\mathrm{l}}}{l_{\mathrm{l}}\rho A_{\mathrm{l}}}\Rightarrow {\dot{Q}}_{\mathrm{l}}=\frac{(p_{\mathrm{l}2}-p_{\mathrm{l}1})A_{\mathrm{l}}-p_{\mathrm{f}}\left(Q_{\mathrm{l}}\right)A_{\mathrm{l}}}{l_{\mathrm{l}}\rho }\end{array}$$

(19)

where $p_{\mathrm{f}}\left(Q_{\mathrm{l}}\right)$ is the pressure drop due to the total line resistance (fittings, orifices, hose). Pressure drop in a straight pipe/hose, $p_{\lambda }$, may be described as:

$$\begin{array}{c}\hfill p_{\lambda }=\underset{\mathrm{turbulent}}{\underbrace{\frac{0.3164}{{\mathit{Re}}^{0.25}}\frac{l_{\mathrm{l}}}{d_{\mathrm{l}}}\rho \frac{1}{2}{\left(\frac{Q_{\mathrm{l}}}{\frac{1}{4}{d}_{\mathrm{l}}^2\pi }\right)}^2}}\underset{\approx 0\mathrm{for}\mathit{Re}<2200}{\underbrace{\left(\frac{1}{2}+\frac{1}{2}\tanh \left(\frac{2300-\mathit{Re}}{100}\right)\right)}}+\underset{\mathrm{laminar}}{\underbrace{\frac{128\nu \rho l_{\mathrm{l}}{Q}_{\mathrm{l}}}{\pi d_{\mathrm{l}}^4}}}\underset{\approx 0\mathrm{for}\mathit{Re}>2400}{\underbrace{\left(\frac{1}{2}+\frac{1}{2}\tanh \left(\frac{-2300+\mathit{Re}}{100}\right)\right)}}\end{array}$$

(20)

where $l_{\mathrm{l}}$ is pipe length, $d_{\mathrm{l}}$ is the lines inner diameter and ν is the kinematic viscosity of the fluid. The decision of flow regime is formed using Reynold number: Re:

$$\begin{array}{c}\hfill \mathit{Re}=\frac{v_{\mathrm{l}}{d}_{\mathrm{l}}}{\nu }\end{array}$$

(21)

The use of the hyperbolic-tangent expression in Equation (20) is for creating a continues transition between linear and turbulent pressure drop for increased robustness of the simulation. $\mathit{Re}=2300$ is used as the transition number. The pressure drops, $p_{\xi },$ of fitting are described as:

$$\begin{array}{c}\hfill p_{\xi }=\xi \frac{\rho }{2}{Q}_{\mathrm{l}}^2{\left(\frac{1}{\frac{1}{4}{d}_{\mathrm{l}}^2\pi }\right)}^2\end{array}$$

(22)

where ξ if a friction coefficient for a given fitting type. Thus, for a line with n fittings and k line pieces, the pressure drop, $p_{\mathrm{f}}\left(Q_{\mathrm{l}}\right)$, due to the total line resistance, may be written as:

$$\begin{array}{c}\hfill p_{\mathrm{f}}\left(Q_{\mathrm{l}}\right)=p_{\lambda ,1}\left(Q_{\mathrm{l}}\right)+\cdots +p_{\lambda ,k}\left(Q_{\mathrm{l}}\right)+p_{\xi ,1}\left(Q_{\mathrm{l}}\right)+\cdots +p_{\xi ,n}\left(Q_{\mathrm{l}}\right)\end{array}$$

(23)

6.1.3. Manifold with Valves and Cylinder Pipe Connection

The valve flows, $Q_{\mathrm{v}\mathrm{xx}}$, are described by the orifice equation, here illustrated for valve $vH3$ (see Figure 13c):

$$\begin{array}{cc}\hfill Q_{\mathrm{vH}3}& =\mathrm{sign}(p_{\mathrm{MH}}-p_{\text{v,A3}})C_{\mathrm{d}}{A}_{\mathrm{v}}\left({\alpha }_{\mathrm{v}}\right)\sqrt{\frac{2}{\rho }|p_{\mathrm{MH}}-p_{\text{v,A3}}|}\hfill \end{array}$$

(24)

where $C_{\mathrm{d}}=0.65$ is the discharge coefficient; ρ is the fluid density and $A_{\mathrm{v}}\left({\alpha }_{\mathrm{v}}\right)$ is the valve opening area. The opening and closing of the valve is described as a ramp function, which is also illustrated in Figure 15a:

$$\begin{array}{cc}\hfill A\left({\alpha }_{\mathrm{v}}\right)& =\alpha A_{\text{o}}\hfill \end{array}$$

(25)

$$\begin{array}{cc}\hfill {\dot{\alpha }}_{\mathrm{v}}\left(t\right)& =\left\{\begin{array}{cc}\hfill \frac{1}{t_{\mathrm{v}}}& ;u_{\mathrm{v}}<1\wedge u_{\text{v3H}}=1\hfill \\ \hfill -\frac{1}{t_{\mathrm{v}}}& ;u_{\mathrm{v}}>0\wedge u_{\text{v3H}}=0\hfill \\ \hfill 0& ;\text{else}\hfill \end{array}\right.\hfill \end{array}$$

(26)

where $A_{\text{o}}$ is the valve opening area when fully opened. The areas, $A_{\text{o}}$, are chosen, such that the pressure drop across a fully open valve is 5 bar at the steady state flow occurring when the cylinder displaces with 0.5 m/s. Thus, the actual value of $C_{\mathrm{d}}$ is not important in this study, as it is compensated for in the chosen $A_{\text{o}}$.

Figure 15. (a) the valve model; and (b) the accumulator model.

Figure 15. (a) the valve model; and (b) the accumulator model.

To avoid simulating the dynamics of the pressure downstream of the valve, $p_{\text{v,A3}}$, which is a small volume with fast dynamics, an algebraic relation between $p_{\text{v,A3}}$ and $p_{\mathrm{MH}}$ is assumed, which is given by Equation (24). This yields the following combined hose and valve model to describe the flow, ${\dot{Q}}_{\mathrm{vH}3}$:

$$\begin{array}{cc}\hfill \mathbf{\text{if}}\left|\alpha \right|>0.05& \end{array}$$

$$\begin{array}{cc}\hfill {\dot{Q}}_{\mathrm{vH}3}& =\frac{(p_{\text{v,A3}}-p_{\mathrm{A}3})A_{\text{l,A3}}-p_{\mathrm{f}}\left(Q_{\mathrm{vH}3}\right)A_{\text{l,A3}}}{l_{\mathrm{A}3}\rho }\hfill \end{array}$$

(27)

$$\begin{array}{cc}\hfill p_{\text{v,A3}}& =p_{\mathrm{MH}}-Q_{\mathrm{vH}3}^2{\left(\frac{1}{C_{\mathrm{d}}{A}_{\mathrm{v}}\left({\alpha }_{\mathrm{v}}\right)}\right)}^2\frac{\rho }{2}\hfill \end{array}$$

(28)

$$\begin{array}{cc}\hfill p_{\mathrm{f}}\left(Q_{\mathrm{A}3}\right)& =p_{\lambda ,A3}\left(Q_{\mathrm{A}3}\right)+p_{\xi \mathrm{M}}\left(Q_{\mathrm{A}3}\right)+p_{\xi \mathrm{c}}\left(Q_{\mathrm{vH}3}\right)\hfill \end{array}$$

(29)

$$\begin{array}{cc}\hfill \text{else}& \end{array}$$

$$\begin{array}{cc}\hfill Q_{\mathrm{vH}3}& =0\hfill \end{array}$$

(30)

The above “if”-condition is to avoid dividing by zero when the opening area goes to zero, i.e., the pressure resistance goes to infinity. The cylinder flows are calculated as $Q_{\mathrm{A}1}=Q_{\text{vH1}}+Q_{\text{vM1}}+Q_{\text{vL1}}$. The fitting resistance, ${\xi }_{\mathrm{M}}$, shown in Figure 13a, is added to take into account internal pressure drops in the manifold. The fittings, ${\xi }_{\mathrm{c}}$, are the cylinder inlet loss, and ${\xi }_{\mathrm{A}1}$ is an extra loss due to an extra bend in this particular pipe connection.

6.1.4. Manifold Accumulators

The operating region of charging and discharging of the gas accumulators will be a combination of adiabatic processes for the fast shifting transients and between adiabatic and isothermal processes for the slower energy smoothing cycles. Resultantly, modeling of the gas temperature and its energy exchange with the surroundings is required to get an adequate estimate on the accumulator losses.

The energy balance of the gas in the accumulator may be written as [29]:

$$\begin{array}{c}\hfill m_{\mathrm{g}}{c}_{\mathrm{v}}\dot{T}=h{A}_{\mathrm{w}}\left(T_{\mathrm{w}}-T\right)-p_{\mathrm{a}}{\dot{v}}_{\mathrm{g}}{m}_{\mathrm{g}}\end{array}$$

(31)

and the standard gas law:

$$\begin{array}{c}\hfill p_{\mathrm{a}}{v}_{\mathrm{g}}=RT\end{array}$$

(32)

where:

h :	Heat exchange coefficient between gas and environment	$A_{\mathrm{w}}$ :	Wall area
$m_{\mathrm{g}}$ :	Gas mass	$T_{\mathrm{w}}$ :	Wall temperature
$v_{\mathrm{g}}$ :	Gas specific volume ($V_{\mathrm{g}}/m_{\mathrm{g}}$)	$V_{\mathrm{g}}$ :	Gas volume
T :	Gas temperature	$p_{\mathrm{a}}$ :	Gas pressure
$c_{\mathrm{v}}$ :	Gas specific heat at constant volume	R :	Ideal gas constant

To simplify the expression, the thermal time constant, ${\tau }_{\mathrm{a}}$, is defined [29]:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{a}}=\frac{m_{\mathrm{g}}{c}_{\mathrm{v}}}{h{A}_{\mathrm{w}}}\end{array}$$

(33)

yielding the equations:

$$\begin{array}{cc}\hfill \dot{T}& =\frac{1}{{\tau }_{\mathrm{a}}}\left(T_{\mathrm{w}}-T\right)-\frac{RT}{c_{\mathrm{v}}{v}_{\mathrm{g}}}{\dot{v}}_{\mathrm{g}}\hfill \end{array}$$

(34)

$$\begin{array}{cc}\hfill {\dot{p}}_{\mathrm{a}}& =\frac{R}{v_{\mathrm{g}}}\dot{T}-\frac{RT}{v_{\mathrm{g}}^2}{\dot{v}}_{\mathrm{g}}\hfill \end{array}$$

(35)

To improve simulation speed, the external volume seen by the accumulator is combined into the dynamics of the accumulator; see Figure 15b. The volume and accumulator dynamics are combined by letting $Q_{\text{V}}=Q_{\mathrm{a}}$ and by adding the external volume, $V_{\mathrm{ext}}$, to the accumulator oil volume. This yields the following flow continuity equation:

$$\begin{array}{cc}\hfill {\dot{p}}_{\mathrm{a}}& =\frac{{\beta }_{\mathrm{eff}}}{V_{\text{a0}}-V_{\mathrm{g}}+V_{\mathrm{ext}}}(Q_{\mathrm{a}}+{\dot{V}}_{\mathrm{g}})\hfill \end{array}$$

(36)

where $V_{\text{a0}}$ is the size of the accumulator; ${\beta }_{\mathrm{eff}}$ is the bulk modulus of fluid in $V_{\mathrm{ext}}$ and $V_{\mathrm{g}}=V_{\text{a0}}$ when $p_{\mathrm{a}}$ is below the pre-charge pressure, $p_{\text{a0}}$.

Rewriting Equation (36), Equation (34) and Equation (35) yields:

$$\begin{array}{cc}\hfill {\dot{V}}_{\mathrm{g}}& =-Q_{\mathrm{a}}+{\dot{p}}_{\mathrm{a}}\frac{V_{\text{a0}}-V_{\mathrm{g}}+V_{\mathrm{ext}}}{{\beta }_{\mathrm{eff}}}\hfill \end{array}$$

(37)

$$\begin{array}{cc}\hfill {\dot{V}}_{\mathrm{g}}& =-{\dot{p}}_{\mathrm{a}}\frac{V_{\mathrm{g}}}{p_{\mathrm{a}}}+\frac{V_{\mathrm{g}}}{T}\dot{T}\hfill \end{array}$$

(38)

$$\begin{array}{cc}\hfill \dot{T}& =\frac{1}{{\tau }_{\mathrm{a}}}\left(T_{\mathrm{w}}-T\right)-\frac{RT}{c_{\mathrm{v}}{V}_{\mathrm{g}}}{\dot{V}}_{\mathrm{g}}\hfill \end{array}$$

(39)

$$\begin{array}{c}\hfill \Downarrow \end{array}$$

$$\begin{array}{cc}\hfill {\dot{p}}_{\mathrm{a}}& =\frac{Q_{\mathrm{a}}+\frac{1}{1+\frac{R}{c_{\mathrm{v}}}}\frac{V_{\mathrm{g}}}{T}\frac{1}{{\tau }_{\mathrm{a}}}\left(T_{\mathrm{w}}-T\right)}{\frac{V_{\text{a0}}-V_{\mathrm{g}}+V_{\mathrm{ext}}}{{\beta }_{\mathrm{eff}}}+\frac{1}{1+\frac{R}{c_{\mathrm{v}}}}\frac{V_{\mathrm{g}}}{p_{\mathrm{a}}}}\hfill \end{array}$$

(40)

Thus, solved in the order, Equation (40), Equation (37) and Equation (39), these equations yield the state equations of the accumulator, including thermal losses. The initial conditions of the accumulator state are related through the standard gas law:

$$\begin{array}{cc}\hfill V_{\mathrm{g}}& =p_{\text{a0}}{V}_{\text{a0}}\frac{1}{T_0}\frac{T}{p_{\mathrm{a}}}\hfill \end{array}$$

(41)

where $p_{\text{a0}}$ is the pre-charge pressure at temperature, $T_0$. The pressure, $p_{\text{V}}$, in the volume seen by the components connected to the accumulator may then be found by the expression for the inlet orifice with opening, $A_{\mathrm{a}}$:

$$\begin{array}{c}\hfill p_{\text{V}}=p_{\mathrm{a}}+Q_{\mathrm{a}}^2{\left(\frac{1}{C_{\mathrm{d}}{A}_{\mathrm{a}}}\right)}^2\frac{\rho }{2}\end{array}$$

(42)

The inlet flow to the accumulators is given as, $Q_{\text{aL}}=Q_{\text{vL1}}+Q_{\text{vML2}}+Q_{\text{vL3}}-Q_{\text{Lh}}$. The external volume, $V_{\mathrm{ext}}$, is set to half the volume of the connecting transmission lines. Depending on the lines, the volume, $V_{\mathrm{ext}}$, is between 3–6 L.

6.1.5. Hoses from Manifold to Pressure Lines

The hoses are modeled as in Section 6.1.2. Here shown for $Q_{\text{ML}}$:

$$\begin{array}{c}\hfill {\dot{Q}}_{\text{ML}}=\frac{(p_{\text{vL}}-p_{\text{ML}})A_{\text{ML}}-p_{\mathrm{f}}\left(Q_{\text{ML}}\right)A_{\text{ML}}}{l_{\text{ML}}\rho }\end{array}$$

(43)

The model parameters for the DDC-system are summarized in Table 2. The hose connection, $l_{\text{MM}}$, length is 3 m for float 1–10 and 5 m for float 11–20. This is performed to spread out transmission line resonance frequencies, yielding a more steady operating pressure line system.

Table 2. Values used for the float module, where “Par.” and “Val.” abbreviate parameter and value, respectively. * hose connection length is 3 m for float 1–10 and 5 m for float 11–20.

Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.
$A_1$ [cm2]:	111	${\eta }_{\mathrm{c}}$ [-]:	0.97	ρ [kg/m3]:	900	$d_{\mathrm{A}1}$ [in]:	1.5	$l_{\mathrm{A}1}$ [m]:	2.0
$A_2$ [cm2]:	196	$x_{\mathrm{c},\max }$ [m]:	3.0	ν [m2s]:	26 $\times {10}^{-6}$	$d_{\mathrm{A}2}$ [in]:	1.5	$l_{\mathrm{A}2}$ [m]:	0.40
$A_3$ [cm2]:	72	${\xi }_{\mathrm{c}}$ [-]:	1.0	${\xi }_{\mathrm{A}1}$ [-]:	1.2	$d_{\mathrm{A}3}$ [in]:	1.5	$l_{\mathrm{A}3}$ [m]:	0.40
$A_{\mathrm{a}}$ [cm2]:	7.9	$A_{\text{o1}}$ [cm2]:	2.8	$A_{\text{o2}}$ [cm2]:	7.7	$A_{\text{o3}}$ [cm2]:	2.8	$t_{\mathrm{v}}$ [ms]:	12
$p_{\text{0ML}}$ [bar]:	10	$p_{\text{0MM}}$ [bar]:	70	$p_{\text{0MH}}$ [bar]:	145			$V_{\text{AM}}$ [L]:	3.6
$T_{\text{W}}$ [° C]:	50	$T_{\text{0}}$ [° C]:	50	R [J/kg/K]:	276	$c_{\mathrm{v}}$ [J/kg/K]:	760	${\tau }_{\text{aH}}$ [s]	15
$d_{\text{ML}}$ [in]:	2.0	$l_{\text{ML}}$ [m]:	3, 5 *	$d_{\text{MM}}$ [in]:	2.0	$l_{\text{MM}}$ [m]:	3, 5 *	${\tau }_{\text{aM}}$ [s]	8
$d_{\mathrm{MH}}$ [in]:	2.0	$l_{\mathrm{MH}}$ [m]:	3, 5 *	${\xi }_{\text{Min}}$ [-]:	0.6	${\xi }_{\mathrm{M}}$ [-]:	1.3	${\tau }_{\text{aL}}$ [s]	4.5

Table 2. Values used for the float module, where “Par.” and “Val.” abbreviate parameter and value, respectively. * hose connection length is 3 m for float 1–10 and 5 m for float 11–20.

Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.
$A_1$ [cm2]:	111	${\eta }_{\mathrm{c}}$ [-]:	0.97	ρ [kg/m3]:	900	$d_{\mathrm{A}1}$ [in]:	1.5	$l_{\mathrm{A}1}$ [m]:	2.0
$A_2$ [cm2]:	196	$x_{\mathrm{c},\max }$ [m]:	3.0	ν [m2s]:	26 $\times {10}^{-6}$	$d_{\mathrm{A}2}$ [in]:	1.5	$l_{\mathrm{A}2}$ [m]:	0.40
$A_3$ [cm2]:	72	${\xi }_{\mathrm{c}}$ [-]:	1.0	${\xi }_{\mathrm{A}1}$ [-]:	1.2	$d_{\mathrm{A}3}$ [in]:	1.5	$l_{\mathrm{A}3}$ [m]:	0.40
$A_{\mathrm{a}}$ [cm2]:	7.9	$A_{\text{o1}}$ [cm2]:	2.8	$A_{\text{o2}}$ [cm2]:	7.7	$A_{\text{o3}}$ [cm2]:	2.8	$t_{\mathrm{v}}$ [ms]:	12
$p_{\text{0ML}}$ [bar]:	10	$p_{\text{0MM}}$ [bar]:	70	$p_{\text{0MH}}$ [bar]:	145			$V_{\text{AM}}$ [L]:	3.6
$T_{\text{W}}$ [° C]:	50	$T_{\text{0}}$ [° C]:	50	R [J/kg/K]:	276	$c_{\mathrm{v}}$ [J/kg/K]:	760	${\tau }_{\text{aH}}$ [s]	15
$d_{\text{ML}}$ [in]:	2.0	$l_{\text{ML}}$ [m]:	3, 5 *	$d_{\text{MM}}$ [in]:	2.0	$l_{\text{MM}}$ [m]:	3, 5 *	${\tau }_{\text{aM}}$ [s]	8
$d_{\mathrm{MH}}$ [in]:	2.0	$l_{\mathrm{MH}}$ [m]:	3, 5 *	${\xi }_{\text{Min}}$ [-]:	0.6	${\xi }_{\mathrm{M}}$ [-]:	1.3	${\tau }_{\text{aL}}$ [s]	4.5

6.2. Pressure Line System

The pressure line system is shown in Figure 16. Accumulator batteries are placed in three locations in the system. The size of accumulator batteries is seen in Figure 14a. The three high pressure accumulator batteries each consist of 16 × 50 L piston accumulators for energy storage. Each intermediate pressure accumulator battery consist of 4 × 50 L, and the low pressure batteries consist of 10 × 50 L.

Figure 16. Pressure line system layout.

Figure 16. Pressure line system layout.

To show how the pressure line system is modeled, a small section is modeled, see Figure 14a, where an example is given for a high pressure line. The same line model as in Section 6.1.2. is used:

$$\begin{array}{c}\hfill {\dot{Q}}_{\text{12H}}=\frac{(p_{\text{H1}}-p_{\text{H2}})A_{\mathrm{H}}-p_{\mathrm{f}}\left(Q_{\text{12H}}\right)A_{\mathrm{H}}}{l_{\mathrm{H}}\rho }\end{array}$$

(44)

To describe the fusion of the flows from float 1 and float 11 at the T-fitting, a pressure at volume, $V_{\text{T}}=\frac{1}{2}{l}_{\mathrm{H}}{A}_{\text{H1}}+l_{\mathrm{MH}}{A}_{\mathrm{MH}}$, is modeled as:

$$\begin{array}{c}\hfill {\dot{p}}_{\text{H1}}=\frac{{\beta }_{\mathrm{eff}}\left(p_{\text{H1}}\right)}{V_{\text{T}}}\left(Q_{\text{MH1}}+Q_{\text{MH11}}-Q_{\text{12H}}\right)\end{array}$$

(45)

As the accumulators are placed close together, they may model-wise be combined into one equivalent accumulator, here shown for the high pressure battery:

$$\begin{array}{cc}\hfill {\dot{p}}_{\text{aH1}}& =\frac{Q_{\text{aH1}}+\frac{1}{1+\frac{R}{c_{\mathrm{v}}}}\frac{V_{\mathrm{g}}}{T}\frac{1}{\tau }\left(T_{\mathrm{w}}-T\right)}{\frac{16·V_{\text{a0H}}-V_{\mathrm{g}}+V_{\mathrm{ext}}}{{\beta }_{\mathrm{eff}}}+\frac{1}{1+\frac{R}{c_{\mathrm{v}}}}\frac{V_{\mathrm{g}}}{p_{\mathrm{a}}}}\hfill \end{array}$$

(46)

$$\begin{array}{cc}\hfill {\dot{V}}_{\mathrm{g}}& =-Q_{\text{aH1}}+{\dot{p}}_{\text{aH1}}\frac{V_{\text{a0}}-V_{\mathrm{g}}+V_{\mathrm{ext}}}{{\beta }_{\mathrm{eff}}}\hfill \end{array}$$

(47)

$$\begin{array}{cc}\hfill \dot{T}& =\frac{1}{\tau }\left(T_{\mathrm{w}}-T\right)-\frac{RT}{c_{\mathrm{v}}{V}_{\mathrm{g}}}{\dot{V}}_{\mathrm{g}}\hfill \end{array}$$

(48)

$$\begin{array}{cc}\hfill V_{\mathrm{g}}& =p_{\text{aH0}}·16V_{\text{a0H}}\frac{1}{T_0}\frac{T}{p_{\text{aH1}}}\hfill \end{array}$$

(49)

$$\begin{array}{cc}\hfill p_{\text{H2}}& =p_{\text{aH1}}+Q_{\text{aH1}}^2{\left(\frac{1}{C_{\mathrm{d}}·16A_{\text{A}}}\right)}^2\frac{\rho }{2}\hfill \end{array}$$

(50)

$$\begin{array}{cc}\hfill Q_{\text{aH1}}& =Q_{\text{12H}}-Q_{\text{M1B}}-Q_{\text{23H}}+Q_{\text{MH2}}+Q_{\text{MH12}}-Q_{\text{F}}\hfill \end{array}$$

(51)

The model parameters for the pressure line system are summarized in Table 3.

Table 3. Parameter values used for the pressure line system, where “Par.” and “Val.” abbreviate parameter and value, respectively.

Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.	Val.	Par.
$A_{\text{A}}$ [cm2]:	11.4	$d_{\text{P}}$ [in]:	2.0	$l_{\text{P}}$ [m]:	6	${\xi }_{\text{P1}}$ [-]:	1.14	${\xi }_{\text{P2}}$ [-]:	1.14
$p_{\text{0L}}$ [bar]:	10	$p_{\text{0M}}$ [bar]:	70	$p_{\mathrm{MH}}$ [bar]:	145	$V_{\text{0H}}$ [L]:	50	$V_{\text{0M}}$ [L]:	50
${\tau }_{\text{AL}}$ [s]:	23	${\tau }_{\text{AM}}$ [s]:	34	${\tau }_{\text{AH}}$ [s]:	50	$V_{\text{0L}}$ [L]:	50

Table 3. Parameter values used for the pressure line system, where “Par.” and “Val.” abbreviate parameter and value, respectively.

Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.	Val.	Par.
$A_{\text{A}}$ [cm2]:	11.4	$d_{\text{P}}$ [in]:	2.0	$l_{\text{P}}$ [m]:	6	${\xi }_{\text{P1}}$ [-]:	1.14	${\xi }_{\text{P2}}$ [-]:	1.14
$p_{\text{0L}}$ [bar]:	10	$p_{\text{0M}}$ [bar]:	70	$p_{\mathrm{MH}}$ [bar]:	145	$V_{\text{0H}}$ [L]:	50	$V_{\text{0M}}$ [L]:	50
${\tau }_{\text{AL}}$ [s]:	23	${\tau }_{\text{AM}}$ [s]:	34	${\tau }_{\text{AH}}$ [s]:	50	$V_{\text{0L}}$ [L]:	50

6.3. Generation System

In this section, the four identical modules of a motor, generator, inverter and grid inverter are described and modeled.

6.3.1. Hydraulic Motor

The hydraulic motors in the system are all 250 cc fixed displacement closed-circuit bent-axis motors. The model in this simulation is based on the efficiency plots for a Sauer-Danfoss 250 cc Series 51-1 bent-axis motor at full displacement. The efficiency plot is seen in Figure 17a. Using the Schloesser formula [30] for friction and flow loss, the following expression gives motor torque and flow inlet:

$$\begin{array}{cc}\hfill Q_{\mathrm{M}}({\omega }_{\mathrm{M}},\Delta p)& =D_{\omega }{\omega }_{\mathrm{M}}-\Delta p{C}_{Q1}\hfill \end{array}$$

(52)

$$\begin{array}{cc}\hfill {\tau }_{\mathrm{M}}({\omega }_{\mathrm{M}},\Delta p)& =D_{\omega }\Delta p-\left(C_{\tau 1}+C_{\tau 2}\Delta p+C_{\tau 3}{\omega }_{\mathrm{M}}+C_{\tau 4}{\omega }_{\mathrm{M}}^2\right)\hfill \end{array}$$

(53)

The efficiency of the fitted model is seen in Figure 17b. Bent-axis motors at frame size 250 cc are allowed to operate at up to 2200 RPM continuously.

Figure 17. (a) efficiency plots (from datasheet) of a Sauer-Danfoss 250 cc Series 51-1 bent-axis motor at full displacement; (b) the efficiency of the modelled motor; and (c) the efficiency of a generator when inverter control is configured for optimal efficiency.

Figure 17. (a) efficiency plots (from datasheet) of a Sauer-Danfoss 250 cc Series 51-1 bent-axis motor at full displacement; (b) the efficiency of the modelled motor; and (c) the efficiency of a generator when inverter control is configured for optimal efficiency.

To replenish the fluid leaked by the motor and to give the necessary flushing and filtering, a charge pump is installed in the system. The pump is set to be $20\%$ of the total motor displacement; thus, $D_{\text{charge}}=200\text{cc}$. The speed of the pump is controlled by a motor with an inverter. The power consumed for flushing is calculated as:

$$\begin{array}{c}\hfill P_{\text{flush}}=D_{\omega \text{charge}}{p}_{\mathrm{L}}{\omega }_{\text{charge}}\frac{1}{0.95·0.95·0.97}\end{array}$$

(54)

where it is assumed that the flushing pump has an efficiency of 95%, the electric motor an efficiency of 95% and the converter an efficiency of 97%.

6.3.2. Generator and Inverter

The generator setup consists of an asynchronous generator and a converter for a grid connection and to enable variable speed control. The input to the generator is the hydraulic motor torque, and the output of the power system is the angular velocity of the generator, ${\omega }_{\text{Gn}}$.

The transient behavior of the electrical part of the generator is neglected. This is assumed acceptable, as the inverter-generator closed loop dynamics are fast compared to the slow load variation. An equivalent circuit for a phase of a three-phase Δ-connected induction motor is seen in Figure 18b, where γ denotes the motor slip. The slip is defined as:

$$\begin{array}{c}\hfill \gamma =1-\frac{n_{\text{pp}}{\omega }_{\text{Gn}}}{{\omega }_{\text{V}}}\end{array}$$

(55)

where $n_{\text{pp}}$ is the number of pole pairs and ${\omega }_{\text{V}}$ is the frequency of the supply voltage. The resistor, $R_2\frac{1-\gamma }{\gamma }$, represents the mechanical power input to the generator, i.e., ${\tau }_{\text{Gn}}{\omega }_{\text{Gn}}=I_2{R}_2\frac{1-\gamma }{\gamma }$; for reference, see e.g., [31].

As the current, $I_2\gg I_{\mathrm{M}}$, the generator torque is given as:

$$\begin{array}{cc}\hfill {\tau }_{\text{Gn}}& =\frac{3n_{\text{pp}}{R}_2}{\gamma {\omega }_{\text{V}}}\frac{V_{\text{RMS}}^2}{{\left(R_1+R_2+\frac{1-\gamma }{\gamma }{R}_2\right)}^2+{\left({\omega }_{\text{V}}(L_1+L_2)\right)}^2}\hfill \end{array}$$

(56)

where $V_{\text{RMS}}$ is the RMS (Root-Mean-Square)-value of the line-to-line voltage.

The steady-state phase current, $I_{\text{P}}$, of the generator is given as:

$$\begin{array}{c}\hfill I_{\text{P}}=\frac{V_{\text{P}}}{|H_{\text{Gn}}\left(j{\omega }_{\text{V}}\right)|},H_{\text{Gn}}\left(j{\omega }_V\right)=\frac{V_{\mathrm{p}}}{I_{\text{P}}}={\left.\frac{Z_2\left(s\right)Z_{\mathrm{M}}\left(s\right)}{Z_2\left(s\right)+Z_{\mathrm{M}}\left(s\right)}+Z_1\left(s\right)\right|}_{s=j{\omega }_{\text{V}}}\end{array}$$

(57)

where $Z_1=R_1+L_{1}s$; $Z_2=\frac{R_2}{\gamma }+L_{2}s$; and $Z_{\mathrm{M}}=\frac{R_{\text{Fe}}s{L}_{\mathrm{M}}}{R_{\text{Fe}}+s{L}_{\mathrm{M}}}$.

The electrical output power of the generator is three-times the power per phase:

$$\begin{array}{c}\hfill P_{\text{Gn,out}}=3V_{\text{P}}{I}_{\text{P}}\cos (\angle H_{\text{Gn}}\left(j{\omega }_{\text{V}}\right))=\sqrt{3}{V}_{\mathrm{L}}{I}_{\mathrm{L}}\cos (\angle H_{\text{Gn}}\left(j{\omega }_{\text{V}}\right))\end{array}$$

(58)

As the chosen hydraulic motors have their best efficiency in the range of 500–2000 RPM, a four-poled generator is used. The chosen generator in the simulation is a asynchronous high-efficiency four-pole 160 kW generator, where the parameters of the equivalent circuit have been obtained. The efficiency result is seen in Figure 18c, where the efficiency is plotted as a function of load at 50 Hz and 400 V.

Figure 18. (a) Delta-connected induction motor; (b) per phase equivalent circuit; and (c) efficiency as a function of load relative to rated torque at 50 Hz and 400 V.

Figure 18. (a) Delta-connected induction motor; (b) per phase equivalent circuit; and (c) efficiency as a function of load relative to rated torque at 50 Hz and 400 V.

The angular velocity of the generator is given by:

$$\begin{array}{c}\hfill {\dot{\omega }}_{\text{Gn}}=\frac{1}{J_{\text{Gn}}+J_{\mathrm{M}}}({\tau }_{\mathrm{M}}-{\tau }_{\text{Gn}})\end{array}$$

(59)

where $J_{\text{Gn}}$ and $J_{\mathrm{M}}$ are the moment of inertia of the generator and the hydraulic motor, respectively.

The torque of the generator is controlled by an inverter. To control the torque of the generator, the expressions for finding the appropriate voltage and voltage frequency for the generator is implemented in the inverter; see Figure 19. The control of the inverter is set to optimize the generator efficiency, i.e., the voltage is reduced and the slip increased at lower speeds to reduce losses. Using the control, the generator efficiency seen in Figure 17c is obtained, widening the high efficiency region of the generator.

The inverter loses are modeled as an efficiency map ${\eta }_{\text{inv}}({\omega }_{\text{V}},P_{\text{inv}})$. Based on the data-sheet for a 160 kW Danfoss VLT-Drive, a look-up table of the inverter efficiency in Figure 19b has been constructed. Finally, the inverter to grid is assumed to have a constant efficiency of ${\eta }_{\text{inv,grid}}=97\%$. The parameters are summarized in Table 4.

Figure 19. (a) torque control of the generator and inverter model; and (b) the efficiency map used for the inverter.

Figure 19. (a) torque control of the generator and inverter model; and (b) the efficiency map used for the inverter.

Table 4. Parameter values used for the hydraulic motor, generator and inverter system, where “Par.” and “Val.” abbreviate parameter and value, respectively.

Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.
$C_{\tau 1}$ [Nm] :	11.22	$C_{\tau 2}$ [Nm/Pa]:	0.17$·{10}^{-6}$	$C_{\tau 3}$ [Nms]:	0.0085	$C_{\tau 4}$ [Nms2]:	0.68$·{10}^{-3}$
$C_{Q1}$ [m3/s/Pa]:	5.4$·{10}^{-12}$	$D_{\omega }$ [m3/rad]:	39.79$·{10}^{-6}$	$J_{\mathrm{M}}$ [kgm2] :	0.048

Table 4. Parameter values used for the hydraulic motor, generator and inverter system, where “Par.” and “Val.” abbreviate parameter and value, respectively.

Par.	Val.	Par.	Val.	Par.	Val.	Par.	Val.
$C_{\tau 1}$ [Nm] :	11.22	$C_{\tau 2}$ [Nm/Pa]:	0.17$·{10}^{-6}$	$C_{\tau 3}$ [Nms]:	0.0085	$C_{\tau 4}$ [Nms2]:	0.68$·{10}^{-3}$
$C_{Q1}$ [m3/s/Pa]:	5.4$·{10}^{-12}$	$D_{\omega }$ [m3/rad]:	39.79$·{10}^{-6}$	$J_{\mathrm{M}}$ [kgm2] :	0.048

7. PTO Control

The PTO control may be divided into two parts, as shown in Figure 20. The first part is the float control, Figure 20a, which controls a DDC-system, compromising between maintaining efficiency and proper power absorption. A Wave Power Extraction Algorithm (WPEA) generates the force reference for optimizing power extraction, while taking into account the PTO efficiency.

The second part is the system control, Figure 20b, which handles the overall control of pressure lines and power generation, i.e., the numbers of active generators and their speed references. The system control may also manipulate the float control to put more or less flow into the mid-pressure line (control parameter, φ) or to reduce the power absorption if the pressure lines are saturated.

7.1. Wave Power Extraction Algorithm

The power extracted by the absorbers from the wave, $P_{\mathrm{ext}}$, is the product of the applied PTO torque, ${\tau }_{\mathrm{PTO}}$, and arm velocity, ${\dot{\theta }}_{\mathrm{arm}}$ (or -$F_{\mathrm{c}}·v_{\mathrm{c}}$). To maximize the extracted power, the float velocity should be in phase with the exiting wave torque, ${\tau }_{\mathrm{ext}}$, i.e., the natural frequency of the float and arm should match the incoming wave. However, frequency-wise point absorbers are narrow-banded with an under-damped resonance frequency, as seen in Figure 12b. Accordingly, point absorbers are prone to operate off-resonance, as the wave period varies wave to wave. On a larger time scale, the average wave period also varies from sea state to sea state, as seen for the wave spectra in Figure 9.

Figure 20. (a) float control (20 parallel system), handling the manifold control to track the generated force reference to extract wave power; and (b) system control for controlling generation and pressure lines.

Figure 20. (a) float control (20 parallel system), handling the manifold control to track the generated force reference to extract wave power; and (b) system control for controlling generation and pressure lines.

To improve the behavior of point absorbers, it is well established that control of the PTO load force may be used, increasing the amount of energy extracted from waves [4,5]. Basically, the torque, ${\tau }_{\mathrm{PTO}}$, applied to the absorber is controlled as a feedback of the absorber’s motion as:

$$\begin{array}{c}\hfill {\tau }_{\mathrm{PTO}}=b_{\mathrm{PTO}}{\dot{\theta }}_{\mathrm{arm}}+k_{\mathrm{PTO}}{\theta }_{\mathrm{arm}}\end{array}$$

(60)

allowing the PTO to adjust the absorber’s resonance frequency to match the wave frequency. The term, $k_{\mathrm{PTO}}{\theta }_{\mathrm{arm}}$, is the reactive term and imposes the PTO to implement a virtual spring element, changing the external experienced dynamics of the absorber.

Optimally, the feedback ensures that the velocity of the absorber is always in phase with the wave excitation force and the PTO damping matching the absorbers’ hydrodynamic damping [4]. However, this requires future knowledge of ${\eta }_{\mathrm{w}}$ for some time, $t_{\text{future}}$, rendering the control non-causal. The required wave prediction horizon to remove the non-causal part is approximately the same as the settling time of the system [4].

A causal and more robust control described for Wavestar in [28] is to continuously gather statistically information of the current sea state over a window of 100 or more waves. The absorbers narrow frequency response is then widened and adjusted to a fixed resonance frequency, yielding the best average abortion in the current sea state. This may be implemented by having fixed control coefficients in the feedback law Equation (60). This is a suboptimal approach, but has been applied with good performance results, as reported in [32].

Optimal adjustment of the absorbers frequency, as the above, requires the PTO to transfer energy to the absorber (assisting its movement) in parts of an oscillation cycle. Thus, the above WPEAs are termed reactive strategies and require a PTO capable of four-quadrant behavior to supply the reactive power. If reactive control is not offered by the PTO, the reactive terms of the control are removed and damping control remains, i.e., for linear damping control only the coefficient, $b_{\mathrm{PTO}}$, is non-zero in Equation (60). Linear damping is only optimal when the wave frequency and the absorber’s resonance frequency match by other means.

For the Wavestar concept, it has been verified from prototype measurement [22] that increased production is obtained through reactive control; see Figure 21a. Each point is an average of extracted power over a 10 min measurement, where both linear damping and reactive control have been tested. The measured average power results are given as a function of significant wave height. Consequently, some of the scatter is due to the different wave periods. The applied control is as in Equation (60), where the coefficients are fixed for a given sea state. The system continuously gather statistically information of the current sea state over a window of 100 or more waves and, then, choose the coefficients yielding the best average abortion in the current sea state.

Figure 21. (a) average extracted power measurements for a single absorber of the prototype in Figure 5a [22]. Each point is a 10 min measurement; (b) comparison of power extraction of linear damping and reaction as a function of efficiency in a sea state with $H_{\mathrm{m},0}=1.75\text{m}$ and $T_{\mathrm{p}}=4.5\mathrm{s}$ [11]; and (c), yearly production of a single absorber at the prototype site [28] as a function of PTO torque limitation and efficiency [11].

Figure 21. (a) average extracted power measurements for a single absorber of the prototype in Figure 5a [22]. Each point is a 10 min measurement; (b) comparison of power extraction of linear damping and reaction as a function of efficiency in a sea state with $H_{\mathrm{m},0}=1.75\text{m}$ and $T_{\mathrm{p}}=4.5\mathrm{s}$ [11]; and (c), yearly production of a single absorber at the prototype site [28] as a function of PTO torque limitation and efficiency [11].

When optimizing reactive control schemes, it is important to take into account PTO power conversion efficiency, as power is lost in each conversion. Hence, having reactive power oscillating between float and PTO consumes energy. The control in Equation (60) should always optimize the power output to the grid and not the amount of extracted power.

In [8,28], it is shown how the PTO force constraints and conversion efficiency may be included into the control design. The results are a map of control parameters, $b_{\mathrm{PTO}}$ and $k_{\mathrm{PTO}}$, as a function of sea state and PTO efficiency. For detailed information on the methods, see [11].

Figure 22. (a) model for reactive control optimization [8]; and (b) control parameters for two sea states as a function of PTO power conversion efficiency, ${\eta }_{\mathrm{PTO}}$.

Figure 22. (a) model for reactive control optimization [8]; and (b) control parameters for two sea states as a function of PTO power conversion efficiency, ${\eta }_{\mathrm{PTO}}$.

The optimal parameters in [8,28] are found by numerically optimizing the output of the simulation model shown in Figure 22. The model includes the absorber dynamics, and the input to the model is an irregular sea wave, ${\eta }_{\mathrm{w}}$. The PTO is characterized by a force limitation and an efficiency. The PTO power output may then be formulated as a function of the extracted power, $P_{\mathrm{ext}}\left(t\right)$:

$$\begin{array}{c}\hfill P_{\text{out}}\left(t\right)=\left\{\begin{array}{cc}{P}_{\mathrm{ext}}\left(t\right){\eta }_{\mathrm{PTO}}\hfill & \hfill ;P_{\mathrm{ext}}\left(t\right)>0\\ P_{\mathrm{ext}}\left(t\right)\frac{1}{{\eta }_{\mathrm{PTO}}}\hfill & \hfill ;P_{\mathrm{ext}}\left(t\right)\le 0\end{array}\right.\end{array}$$

(61)

The optimal control parameters for Equation (60) is then found as the arguments maximizing the energy output of the PTO:

$$\begin{array}{cc}\hfill \underset{k_{\mathrm{PTO}},b_{\mathrm{PTO}}}{\text{arg}\text{max}}& {\int }_0^{t_{\mathrm{f}}}{P}_{\text{out}}\left(t\right)\mathrm{d}t\hfill \end{array}$$

(62)

The optimal values of $b_{\mathrm{PTO}}$ and $k_{\mathrm{PTO}}$ for two sea states are shown in Figure 22b. These results are from [11]. Here, it is seen as expected, that when the PTO conversion efficiency, ${\eta }_{\mathrm{PTO}}$, is reduced, the reactive part ($k_{\mathrm{PTO}}$) reduces and the resistive part increases. To further justify the reactive control, a comparison of reactive control and linear damping is shown in Figure 21b [11]. It shows that despite having a power conversion efficiency of, e.g., 70%, 50%, more power is produced to the grid compared to linear damping.

The torque requirement of the PTO is chosen to be $1\text{MNm}$, which is a compromise made according to Figure 21c, as the PTO efficiency is expected to be about 70%. In Figure 21c, the yearly production for a single absorber has been computed in [11] by calculating the average power output for all sea states and, then, assuming a yearly wave distribution.

To summarize, given the sea state and PTO information, $(H_{\mathrm{m},0},T_{\mathrm{p}},{\eta }_{\mathrm{PTO}})$, the force reference to the cylinder is calculated as:

$$\begin{array}{c}\hfill F_{\text{ref}}=\gamma \left(b_{\mathrm{PTO}}(H_{\mathrm{m},0},T_{\mathrm{p}},{\eta }_{\mathrm{PTO}}){\dot{\theta }}_{\mathrm{arm}}+k_{\mathrm{PTO}}(H_{\mathrm{m},0},T_{\mathrm{p}},{\eta }_{\mathrm{PTO}}){\theta }_{\mathrm{arm}}\right)\frac{1}{d_{\mathrm{a}}\left({\theta }_{\mathrm{arm}}\right)}\end{array}$$

(63)

where $d_{\mathrm{a}}$ is the moment arm and γ is a coefficient, which the system control may use to reduce the power absorption if the pressure lines are saturated. To reduce power absorption, a method could also have been to reduce the efficiency input to the WPEA algorithm, forcing the WPEA to move towards linear damping, thereby reducing power absorption and system stress.

The actual PTO efficiency is not independent of the WPEA defined in Equation (63). Adjusting the reactive power may either improve or reduce component performances, depending on whether the increased power level brings the components closer to their rated/optimum power level. To this account, an initial efficiency estimate is given to find the WPEA, and through simulation of the PTO, a new efficiency is obtained, which is used to update the WPEA. A number of iterations have been performed to find the proper settings. This is also discussed in [8] for the Wavestar converter.

The force control is going to be implemented discretely, i.e., the continuous reference generated by Equation (63) is subjected to a quantification. To investigate the influence, comparison of a continuous control and a discrete control has been performed.

The continuous control is simulated as in Figure 22a, with the optimum parameters, a PTO efficiency of 70% and a maximum torque of ±1 MNm. The simulation of the discrete version is implemented likewise with the same WPEA parameters, but the PTO torque is quantified before being applied to the float. The available torque values besides zero are set to 9 evenly spaced values between $\pm 1\text{MNm}$. Moreover, a limit is set on the shift rate, limiting the force to shift at maximum one time per 300 ms. Thus, after a shift has occurred, the force value is locked for 300 ms, whereafter the quantification process may choose the force closest to the continuous reference. This avoids the risk of jittering and PWM (Pulse Width Modulated)-like behavior and keeps the number of shifts at a reasonable level. The results are seen in Figure 23, showing the two power matrices of average produced power. As seen, the difference is relatively low despite the relatively rough quantization. Thus, the force approximation of the DCC having 27 force values will be adequate.

Figure 23. Power matrices of average produced power in [kW]. The matrices are shown for a continuous torque control and a discrete control with $\pm 4$ force values.

Figure 23. Power matrices of average produced power in [kW]. The matrices are shown for a continuous torque control and a discrete control with $\pm 4$ force values.

The reason for the low impact of the discrete variation is that the 35 tons of the float and arm system serve as an effective low pass filter, having a break frequency around 0.28 Hz, c.f. the earlier calculated natural period of 3.5 s. Thus, from a load control point of view, implementation of extra dampers/low pass filters is not required.

Another aspect is the compliance of the mechanical structure when performing the rapid load force changes; however, these issues are not treated in this paper. One aspect to notice, though, is that the force steps are relatively small, e.g., a single step from zero to maximum is never performed.

7.2. DDC Control—Force Shifting Algorithm (FSA)

The purpose of the Force Shifting Algorithm (FSA), controlling the DDC, is to choose the appropriate force level to approximate the reference, $F_{\text{ref}}$, generated by the WPEA. However, always shifting to the force closest to the continuous reference, $F_{\text{ref}}$, may not energy-wise be optimal, as some force-shifts are more energy-expensive than others, due to the compression loss, $E_{\beta }$. The minimum compression loss of a volume when shifting $\Delta p$, independent of the pipe line dynamics [33], is given as [16]:

$$\begin{array}{c}\hfill E_{\beta }=\frac{1}{2}\Delta p^2\frac{V}{\beta }\end{array}$$

(64)

This loss is illustrated in Figure 24a for a system consisting of four system pressures and a symmetric cylinder. The figure shows the compression loss in [kJ] from shifting between forces. The forces are sorted according to size. As illustrated, if the current force is force number 9 and an increase in force is desired, it is cheaper to shift to force number 11 instead of 10. Likewise, from number 11, it is cheaper to shift to number 14 and skip 12 and 13. Thus, to avoid doing very expensive force shifts, a more suitable strategy is to calculate the energy expense of possible force-shifts and make a compromise between tracking and energy-cost. This is the basis of the FSAstrategy developed in [15] and which is also used in this study. Note that the loss matrix changes with pressure and cylinder position.

Figure 24. (a) illustration of the shifting loses due to compression; (b) illustration of the force shifting algorithm; and (c) the opening and closing procedure of the valves during shifting.

Figure 24. (a) illustration of the shifting loses due to compression; (b) illustration of the force shifting algorithm; and (c) the opening and closing procedure of the valves during shifting.

The shifting loss from shifting from force x to y is found summing the losses Equation (64) for pressure changes in the individual volumes:

$$\begin{array}{c}\hfill E_{\text{shift}}(x,y)=\frac{1}{2}{(p_{\text{A1,old}}-p_{\text{A1,new}})}^2\frac{V_{\mathrm{A}1}\left(x_{\mathrm{c}}\right)}{\beta }+\frac{1}{2}{(p_{\text{A2,old}}-p_{\text{A2,new}})}^2\frac{V_{\mathrm{A}2}\left(x_{\mathrm{c}}\right)}{\beta }+\frac{1}{2}{(p_{\text{A3,old}}-p_{\text{A3,new}})}^2\frac{V_{\mathrm{A}3}\left(x_{\mathrm{c}}\right)}{\beta }\end{array}$$

(65)

Note that the pressures, $p_{\text{Ax,old}}=p_{\text{Ax,new}}$, may be equal, as the pressure is not necessarily shifted in all chambers.

To allow compromising between force tracking and performing less expensive force shifts, a maximum allowed tracking error, $F_{\text{b}}$, is defined, meaning that $F_{\mathrm{c}}$ must stay within a band of $\pm F_{\text{b}}$ about $F_{\text{ref}}$. However, within the band, the FSA may choose the force steps with the lowest shift cost. A fixed time limit, $T_{\text{min}}$, on how frequently shifting is allowed is also added to reduce tracking cost.

The FSA algorithm is illustrated in Figure 24b. The values, $k_{-}$ and $k_{+}$, are the numbers of the forces, which are the cheapest to shift to within the band, $\pm F_{\text{b}}$. As seen, these are continuously updated. When the force reference comes closer to one of the forces, $k_{-}$ or $k_{+}$, a shift is performed, and a new lock down period of $T_{\text{min}}$ is initiated.

Based on the optimization procedure in [15], it has been found that $F_{\text{b}}=150\text{kN}$ and $T_{\text{min}}=350\text{ms}$ for optimizing the power output of the Wavestar C5.

Thus, when a force shift is initiated, the control sends out a matrix, ${\mathit{u}}_{\mathrm{c}}$, of control values for the nine on/off valves, corresponding to the desired pressure configuration. Regarding valve timing, it was found in [16] that a small amount of overlap between opening and closing of the valves for a single volume was desirable; thus, a 3 ms overlap is used. As the valves have 12 ms opening and closing time, the signal to the opening valve is delayed 9 ms. This is illustrated in Figure 24c.

To help in controlling the pressure in the mid-pressure line, the FSA continuously identifies the force combinations that would currently supply or consume flow from the mid-pressure line. For these configurations, an “artificial” energy loss, φ, may be added in Equation (65) by the system control to either penalize supplying or consuming flow from the mid-pressure line.

7.3. System Control

The purpose of the system control is to:

• Avoid the high pressure accumulator storage from depletion or saturation;
• Keep the mid-pressure line floating between high and low pressure;
• Ensure as steady a power production as possible, while satisfying 1 and 2;
• Choose the proper number of generators for a given sea state;
• Reduce power absorption when full load capacity is reached.

To maintain a stable power production, the generators are set to initially produce the expected average power, $p_{\text{avg,expected}}$, in the current sea state. An initial guess is given based on the current sea state when starting production, whereafter a moving average is used based on the absorbed power over a window of 5 min. The number of active generators, $k_{\text{GN}}$, is then chosen, such that generation capacity is roughly $p_{\text{avg,expected}}$ plus 30%. All active generators are operated at the same speed ${\omega }_{\text{GN}}$.

To avoid the high pressure accumulator storage from depleting or saturating, the power generation is increased or decreased based on the pressure in the accumulators. This is performed through the coefficient, ψ. The speed reference is then given as ${\omega }_{\text{GN}}$:

$$\begin{array}{c}\hfill {\omega }_{\text{GN,ref}}=\frac{p_{\text{avg,expected}}\frac{1}{{\eta }_{\text{to-Gn}}}\psi }{p_{\mathrm{H}}{k}_{\text{GN}}{D}_{\mathrm{M}}}\end{array}$$

(66)

where $D_{\mathrm{M}}$ is the total active motor displacement in [m³/rad]. As $p_{\text{avg,expected}}$ is the absorber power, the efficiency from cylinder to power out of the generator, ${\eta }_{\text{to-Gn}}$, is required to calculate the average generator power. The maximum allowed speed is set according to 180 kW per generator, which is 15% overload. The lowest speed is set to 400 RPM.

The high pressures are set to be between 150 bar and 300 bar; thus, the map in Figure 25 for ψ is used. Finally, if pressure is still reaching 300 bar, the value, γ, is manipulated, such that the float reduces power absorption; see Figure 25. This is similar to turbines pitching out of the wind when rated production is reached.

Figure 25. Maps for the different system control parameters, ψ, φ and γ.

Figure 25. Maps for the different system control parameters, ψ, φ and γ.

The mid-pressure is to be floating between high and low pressure for optimizing the DDC efficiency. Hence, the “reference”, $p_{\text{M,ref}}$, for the mid-pressure line is:

$$\begin{array}{c}\hfill p_{\text{M,ref}}=\frac{p_{\mathrm{H}}+p_{\mathrm{L}}}{2}\end{array}$$

(67)

Based on this, the penalty value, φ, is set according to Figure 25.

8. Results and Discussion

The modeled PTO system is evaluated for the three sea states, 1, 2 and 3, which were defined in Figure 9. These correspond to small, medium and high production cases, respectively. To evaluate the PTO, a test wave has been created using the described white noise method. The first 100 periods of the wave have been checked to have a representative energy content complying to the PM-spectrum. Additionally, the test wave contains both high intensity peaks and calm periods between wave groups. To compare the simulations, this wave is reused for all three sea states by scaling it in time and height correspondingly. The results are shown for an incoming wave angle of ${\theta }_{\mathrm{w}}=0^{\circ }$. The PTO performance is nearly independent of wave direction, as the accumulators handle the wave-to-wave power fluctuation.

First, the PTO is overall evaluated for the three sea states. Afterwards, the cylinder force tracking, and control performance is inspected.

8.1. Overall System Performance

In Figure 26, approximately 6 min have been simulated of the complete Wavestar C5 for sea state 1. The power absorption is the instantaneous mechanical power input to the 20 cylinders. For sea state 1 with an average power production of about 120 kW, two generators are active, namely, generators 1 and 3, such that pipe losses are reduced.

As seen in Figure 26, the PTO system is capable of smoothing out the varying power absorption using the accumulators, whereby the power output to the grid is fairly stable. Note that at around times, $t=50\mathrm{s}$ and $t=150\mathrm{s}$, long periods of almost no absorbed power are experienced, but the system is able to maintain output. The results show that all available “peak power” of the waves is extracted, as the high pressure, and, thereby, the accumulators’ storage, never saturates [$p_{\mathrm{H}}\left(t\right)<$320 bar]. The control system stabilizes the mid-pressure between high and low pressure as desired. Note that at small waves, the high pressure level is kept a bit lower to reduce compression losses, as the highest cylinder forces are not required during smaller waves.

Figure 26. (a) PTO simulation results for Wavestar C5 for sea state 1, $H_{\mathrm{m},0}$ = 1.00 m, $T_{\mathrm{p}}$ = 4.5 s; and (b) PTO simulation results for Wavestar C5 for sea state 2, $H_{\mathrm{m},0}$ = 1.75 m, $T_{\mathrm{p}}$ = 5.5 s.

Figure 26. (a) PTO simulation results for Wavestar C5 for sea state 1, $H_{\mathrm{m},0}$ = 1.00 m, $T_{\mathrm{p}}$ = 4.5 s; and (b) PTO simulation results for Wavestar C5 for sea state 2, $H_{\mathrm{m},0}$ = 1.75 m, $T_{\mathrm{p}}$ = 5.5 s.

The generators are operated at the same speed, which is seen in the bottom of Figure 26. The speed is slowly varying between 850 to 1,200 RPM; likewise, the generator input torque is also stable around 800 Nm. Thus, the operating conditions for hydraulic motors, generators and inverters are very favorable. Hence, the PTO strategy of using the hydraulic accumulators and DDCs effectively separates the generation from the complicated and highly fluctuating wave power absorption process. Looking at the power absorption of a single float, the factor between peak and average power is 12, whereas the fluctuation of output power is within $\pm 50\%$ of mean output.

In the simulation for sea state 2 in Figure 26, all four generators are active. Compared to the PTO behavior in sea sate 1, the storage capacity is now almost fully used, with the high pressure varying between 160 bar and 300 bar. This means that the full storage capacity of the accumulator battery is being utilized. The system completely absorbs the 10 s long 2 MW peak, while maintaining a steady production around 400 kW.

Moving to sea state 3 (Figure 27), the system is operating at full load most of the time. As seen from the pressure, the storage saturates; however, the system is kept below the maximum of 320 bar by reducing the load force on the absorber. Thus, no extra energy has to be dissipated internally in the system using, e.g., relief valves, although these are installed. Looking at the power absorption of a single float, the mean power is 45 kW and the peak power is 600 kW, yielding a factor of 12 between mean and peak power.

Figure 27. PTO simulation results for Wavestar C5 for sea state 3, $H_{\mathrm{m},0}$ = 2.50 m, $T_{\mathrm{p}}$ = 6.5 s.

Figure 27. PTO simulation results for Wavestar C5 for sea state 3, $H_{\mathrm{m},0}$ = 2.50 m, $T_{\mathrm{p}}$ = 6.5 s.

To evaluate how efficiently the PTO converts the absorbed power, the power input and output of the different sub systems have been integrated during the simulation to identify the sub-system losses.

$$\begin{array}{cc}\hfill P_{\text{in,avg}}=\frac{1}{t_{\text{end}}}{\int }_0^{t_{\text{end}}}{P}_{\text{in}}\left(t\right)\mathrm{d}t,& P_{\text{out,avg}}=\frac{1}{t_{\text{end}}}{\int }_0^{t_{\text{end}}}{P}_{\text{out}}\left(t\right)\mathrm{d}t,\eta =\frac{P_{\text{out,avg}}}{P_{\text{in,avg}}}\hfill \end{array}$$

(68)

The system is divided as shown in Figure 28a, where the flushing and charging power consumption is also taken into account. The results are summarized in Table 5 for the different sea state, showing the average power inputs and outputs.

Figure 28. (a) the definition of average power input and output of the different subsystems; and (b) how the difference in stored energy is calculated.

Figure 28. (a) the definition of average power input and output of the different subsystems; and (b) how the difference in stored energy is calculated.

Table 5. Efficiencies of subsystems in Figure 28. SS1, SS2 and SS3 are abbreviations of the three tested sea states (1, 2, 3).

Table 5. Efficiencies of subsystems in Figure 28. SS1, SS2 and SS3 are abbreviations of the three tested sea states (1, 2, 3).

To calculate efficiency accurately, it is required to take into account that the stored pressure energy in accumulators and oil volumes is not the same at the initial state and at simulation end. Resultantly, a post-simulation is performed on the pressure line system and manifold using the end states. Here, how much energy is required to fill or drain the pressure lines to get to the initial state is simulated, which is illustrated in Figure 28b. The difference in energy is then mapped into an equivalent average power input, $P_{\Delta \text{init}}$, to the pressure line system.

In Table 5, the total PTO efficiency, ${\eta }_{\text{tot}}$, is defined as the ratio between total absorbed energy and the energy to grid, taking into account the difference in the stored energy:

$$\begin{array}{c}\hfill {\eta }_{\text{tot}}=\frac{E_{\text{grid}}}{E_{\text{in}}+E_{\Delta \text{init}}\frac{1}{{\eta }_{\text{DDC}}}}=\frac{P_{\text{grid}}}{P_{\text{in}}+P_{\Delta \text{init}}\frac{1}{{\eta }_{\text{DDC}}}}\end{array}$$

(69)

In the above expression, the difference in stored energy, $P_{\Delta \text{init}}$, is scaled with $\frac{1}{{\eta }_{\text{DDC}}}$, such that it seems that it has been absorbed from the waves.

From Table 5, it is seen that the total PTO efficiency, ${\eta }_{\text{tot}}$, is 70%, 73%, 71% for the three sea states, respectively. Thus, despite the varying power levels and a factor of 10 between the mean and peak power, the PTO system is able to maintain efficiency.

The motors and generator all operate around 94% and 95%, as expected in the operating conditions. Most important in the results is that the DDC-system is able to maintain an overall efficiency of ${\eta }_{\text{DDC}}=90\%$, which is the core technology for enabling the PTO concept. The actual power conversion efficiency is actually higher, as reactive control is being performed, making the average processed power larger than the average input power. This may be viewed from the single float power absorption, which is periodically negative. The average negative power, $P_{\text{ext,−}}$, and average positive, $P_{\text{ext,+}}$, power have been calculated for float 1 in the three sea states:

• Sea state 1: $P_{\text{ext,−}}=0.57\text{kW}$ and $P_{\text{ext,+}}=7.5\text{kW}$: Ratio of reactive and real power: $c=\frac{P_{\text{ext,+}}}{P_{\text{ext,−}}}=0.075$;
• Sea state 2: $P_{\text{ext,−}}=2.5\text{kW}$ and $P_{\text{ext,+}}=19.5\text{kW}$: Ratio of reactive and real power: $c=\frac{P_{\text{ext,+}}}{P_{\text{ext,−}}}=0.088$;
• Sea state 3: $P_{\text{ext,−}}=3.1\text{kW}$ and $P_{\text{ext,+}}=49.0\text{kW}$: Ratio of reactive and real power: $c=\frac{P_{\text{ext,+}}}{P_{\text{ext,−}}}=0.063$.

Assuming equal power conversion efficiency in both direction, the relation of the actual power conversion efficiency of the DDC-system, ${\eta }_{\text{act}}$, and the ‘measured’ average efficiency, ${\eta }_{\text{DDC}}$, is given as:

$$\begin{array}{c}\hfill {\eta }_{\text{DDC}}=\frac{P_{\text{ext,+}}{\eta }_{\text{act}}-\frac{1}{{\eta }_{\text{act}}}{P}_{\text{ext,-}}}{P_{\text{ext,+}}-P_{\text{ext,−}}}\Rightarrow {\eta }_{\text{act}}=\frac{{\eta }_{\text{DDC}}(1-c)+\sqrt{{\eta }_{\text{DDC}}^2{(1-c)}^2+4c}}{2}\end{array}$$

(70)

The actual power conversion efficiency for the manifold and cylinder is, thus, $90\%$, $92\%$ and $91\%$ for the sea states, 1, 2 and 3, respectively.

8.2. Control System Performance

In Figure 29, the operation of the cylinder and manifold of float number 1 is shown. The control tracks the force reference with the discrete steps adequately and is excellent at not performing extra unnecessary shifts, when the variation in force reference is small, for example, at time 160 s and at time 166 s.

Looking at the pressures within the three chambers, the shifting algorithm optimizes efficiency by mostly shifting one pressure at a time and using the mid-pressure as an intermediate level. The chamber, $A2$, is the largest chamber and, thereby, the most energy expensive to shift. This is seen in the control output, as $A2$ is the chamber that experiences the lowest number of pressure shifts, whereas the smallest chamber, $A3$, experiences the most pressure shifts.

In the zoom in on the pressure, the shifting transients in pressure may be seen. The fast oscillatory behavior in the transients are due to the fact that the flow may oscillate between the cylinder chamber and manifold accumulator, due to the pipe line between the valve and cylinder. These are most pronounced for the $A1$ chamber, as this chamber has the longest pipe connection, thereby increasing the effect of pipe “inductance”. The slower oscillation in pressures $p_{\mathrm{H}}$ and $p_{\mathrm{M}}$ are likewise caused by flow oscillating between the small manifold accumulators and the large storage accumulator batteries.

Figure 29. Force tracking of cylinder.

Figure 29. Force tracking of cylinder.

In Figure 30, the system control performance is shown. The mid-pressure reference, $p_{\text{M,ref}}$, is set between low and high pressure, and the mid-pressure is effectively kept floating around this reference by manipulating the force shifting algorithms. This was implemented by giving a penalty for shifting to force configuration, which will take the mid-pressure further away from the reference. The penalizing “energy” input, φ, to the force shifting algorithms is seen in the figure.

Figure 30. Controller performance and gas temperature of accumulators.

Figure 30. Controller performance and gas temperature of accumulators.

In Figure 30, the temperature variation in one of the high pressure storage accumulators is seen. The process is not isothermal, as the temperature is highly varying, and the process is not purely adiabatic either, as the temperature and pressure do not always return to the same initial conditions. Thus, the inclusion of heat transfer between accumulator and environment in the model is justified.

Finally, Figure 30 shows how the system control uses the γ input to downscale the aggressiveness of the float absorption, such that high pressure stays below 320 bar.

9. Conclusions

The paper has demonstrated a PTO-system for the Wavestar WEC capable of solving the following challenges:

• Handling peak power input, which is a factor of 12 higher than mean power, while maintaining component efficiency—the DDC maintains above 90% in these conditions, and the remaining components are above 94% in efficiency, each;
• Full controllability of load force on absorbers and four quadrant mode—the DDC offers 27 force steps and four quadrant mode;
• Incorporating a short-term storage for supplying reactive power—reactive power is only processed by DDCs and accumulators;
• Incorporating an efficient energy storage for power smoothing—generators are operated independent of the wave absorption; energy storage for operating generators at 1500 RPM for one minute;
• Maintaining PTO efficiency in small waves when operating at 15% of full load capacity—total efficiency maintained above 70% in all sea states;
• Being able to reduce power absorption when full load capacity is reached—the DDCs reduce absorption when the WEC reaches full load;
• Being scalable to future multi-MW systems.

Based on rigorously modeling the system from wave-to-wire, the overall efficiency of the PTO was found to beyond 70% in all sea conditions, while providing excellent power absorption.

The utilized power absorption algorithms (WPEA) were based on reactive control methods, tuning the resonance of the system. The algorithms were tuned to take into account the efficiency of the PTO to maximize electrical power production.

The WPEA determined load force reference was applied using the specially developed DDC-modules mounted on each float, consisting of discrete controlled multi-chambered cylinders. The DDC modules directly convert and store the absorbed wave power as high pressure energy in the accumulators with a loss less than 10%.

The implemented PTO system had an extra intermediate pressure line for improved efficiency of the DCC. It was shown that by overall system control, the net flow to the intermediate line could be kept at zero; hence, extra pumps and motors for supporting this line are not required. The simulations verified that when reaching full load capacity, the system would reduce the power absorption using the DDCs, such that full load is sustained, but no extra energy has to be dissipated internally in the system.

Regarding scalability, the DDC-modules are fully scalable to be increased to larger systems. Currently, commercial valves with the required transient opening and closing times (≤15 ms) are available with only 2% loss at $1.8\text{MW}$ and with a peak power level of more than 4 MW. This is more than required for, e.g., a 6 MW Wavestar system with 20 floats. Furthermore, valves may be easily used in parallel, thereby also increasing redundancy. Regarding hydraulic motors, commercial 1000 cc high speed motors are available, producing 750 kW at 300 bar and 1500 RPM. Within the wind turbine industry, hydraulic transmission are also being investigated, leading to development of fast multi-MW motors with high efficiency; see, e.g., Digital Displacement^® technology. These variable displacement also shows high part load efficiency, leading to the option of operating at a fixed generator speed of the 1500 RPM, thereby potentially removing the power converters.

The size of the storage may be increased as desired. The storage size is a cost optimization problem between power smoothness and accumulator cost. Increasing the storage does not reduce the overall efficiency, as the round-trip efficiency of accumulators is around 97%. Instead, increasing storage may increase efficiency and durability, as it narrows/stabilizes the operating region of the remaining PTO components, thereby making them operate near their optimum point at a constant load.

Looking at potential improvements, the control of the DDC-modules presented in this paper is a simplified version, not taking into account, e.g., cylinder velocity and simultaneous shifting of multiple-chambers. The system control is also a very simplified version, which does not use the energy storage at full potential to stabilize and increase energy production. Additionally, further optimization on the pressure line network may be performed, reducing pipe losses and improving transient behavior. Improved hydraulic motor-efficiencies are also obtainable with commercial available components (96% efficiency) and coming digital displacement motors. Using the rigid-body assumption, it was through simulation verified that the discrete force control yielded approximately the same energy extraction performance as a continuous control. However, regarding compliance of the mechanical system, the presented work did not cover this issue. For example, when using the DDC, there is a risk of exciting the structural eigen-frequencies. Thus, when applying the described DDC, a detailed dynamical mapping of the mechanical system should be performed as a part of the DDC design.

With the potential improvements, the PTO concept has been assessed to be able to reach about 80% efficiency from mechanical input to electrical output.