Research works undertaken in the first author’s laboratory at the University of Tokyo over the past 30 years are highlighted. Finding of the occurrence of nonlinear waves (named Free-Surface Shock Waves) in the vicinity of a ship advancing at constant speed provided the start-line for the progress of innovative technologies in the ship hull-form design. Based on these findings, a multitude of the Computational Fluid Dynamic (CFD) techniques have been developed over this period, and are highlighted in this paper. The TUMMAC code has been developed for wave problems, based on a rectangular grid system, while the WISDAM code treats both wave and viscous flow problems in the framework of a boundary-fitted grid system. These two techniques are able to cope with almost all fluid dynamical problems relating to ships, including the resistance, ship’s motion and ride-comfort issues. Consequently, the two codes have contributed significantly to the progress in the technology of ship design, and now form an integral part of the ship-designing process.

Naval architecture represents one of the most important components of transportation technology, and great progress in the subject has been made in the 20th Century. In this review paper, the progress made in the hydro-dynamical technologies in the field of naval architecture is summarized, derived mostly from the research carried out over the past 30 years in the Ship Model Basin Laboratory at the University of Tokyo.

Of principal importance is the understanding of the physical phenomena influencing the mechanisms, both of ships in particular and more general of engineering products. Without a thorough understanding of the basic physical phenomena, no useful technology can be developed within the field of design engineering. In naval hydrodynamics, there are two major physical phenomena of substantial importance: nonlinear waves (including wave breaking) and turbulent flow (including large-scale, separated flow situations). Typically, these are nonlinear and very complex phenomena. Only by understanding all the physical phenomena which occur around a ship, can the most efficient ship, with the most optimal hull configuration, be designed, leading finally to the advancement in the transportation technology.

The waves generated by ships advancing steadily in deep water had long been considered characterized by typical, linear-dispersive waves. Many wave-making theories have been proposed by mathematicians and naval architects, most of which are based on the postulate of linearity of the physical phenomena. However, resulting from a series of experiments begun in 1977 at the University of Tokyo, the occurrence of a kind of shock wave was identified in the near-field region of the ship, the characteristics of which were subsequently explained. This nonlinear wave was named the Free-Surface Shock Wave. The first part of the present article provides some details of these experimental works.

After the progress made in the elucidation of the basic physical phenomena, the techniques of Computational Fluid Dynamics (CFD) were developed to further the design of ship’s hull-forms, since theoretical fluid dynamics has definite limitations in the description of nonlinear phenomena. The governing, partial-differential equations are nonlinear, and thus need to be discretely solved using numerical techniques. In the beginning, the objective was to develop a numerical technique which can resolve the nonlinear waves generated by a ship advancing at constant speed on a straight course. Since wave resistance plays an important role in the overall resistance of the ship, quantification of the resistance due to nonlinear phenomena at sea is very important in the ultimate design of the ship, economical performance, and ultimately minimal fuel-consumption. This optimum design was achieved by numerically solving the steady-flow situation, together with the dynamics of the waves, for a ship on a steady, straight course on a calm sea. This method of approach was subsequently embodied in a code named TUMMAC. The first version of the TUMMAC code appeared in 1983: it was the first one of its kind in the world, and was distributed to the major shipbuilding companies in that year. In the TUMMAC code, a rectangular grid system was employed, and viscous effects were not taken into account (free-slip condition at wall boundaries). Subsequently, some further versions of the TUMMAC code were released, which could also deal with the wave-breaking problem.

Quite independently, CFD-type techniques were also developed in the framework of the boundary-fitted grid system approach, and employing the Reynolds-Averaged Navier-Stokes (RANS) equations of turbulence modelling; these effects were embodied in the code WISDAM. In 1987, the combined viscous-flow and wave-propagation problem of flow around a ship was solved simultaneously for the first time using the WISDAM code. Following the rapid progress in computer hardware, the technique has recently been extended to the simulation of maneuvering and motion of ships in waves. Finally, the CFD solution can also take into account the motion of the ship itself. Since arbitrary sea conditions can be realized in the numerical simulation, the motion of the ship itself can be included in the simulation, thus reproducing all motions of the ship in actual sea conditions. In the latter part of this article, this latest advancement in numerical simulation is described in detail.

It may be opportune to say that we are now working towards a virtual reality approach to ship design, based on Computational Fluid Dynamics. All attitudes, forces and moments, as well as the physical phenomena arising from the sea conditions, may be realized in terms of computer simulation derived from our application software package.

In the last two years (1975 to 1977) of a five-year period as a ship designer at the Ship Initial Design Office of IHI (Ishikawajima-Harima Heavy Industries Co., Ltd.), the first author of this article was in charge of the hull-form design of the FUTURE-32 bulk carrier (32,000 DWT); orders for 17 such vessels were received. However, the requirement of speed was very severe, and the achievement of the guaranteed trial speed written in the specifications was thought to be unrealizable. The official trial results for the first two ships vindicated the claim: they were faster than the guaranteed speed by 0.6 knots under ballast conditions. This meant that the propulsion horsepower had been over-estimated by 15%, and a 6-cylinder engine could have been installed instead of the 7-cylinder one, at considerable cost benefit.

This experience provided a very good justification for research work on non-linear ship waves. And, immediately following the delivery of first two FUTURE-32 ships, the first author of this article was invited to lead the Ship Model Basin Laboratory of the University of Tokyo as an Associate Professor, in which one of the prominent research activities was concerned with wave pattern picture analysis. Consequently, subsequent research work on non-linear ship wave started with the examination of ship wave pattern pictures. First insights were provided by an analogy to shock waves in supersonic flow.

A series of experiments was conducted. In the first of these, the non-linearity of the bow-wave phenomenon was exaggerated by use of a model ship with a round-ended bow shape, as shown in Fig. 1.^1)–3)

Fig. 1.

Free-surface shock wave (drawn with dashed lines for lucidity) around a wall-sided model with a blunt bow, at Fn = 0.12, 0.16, 0.24, 0.26 (from left to right).

It was shown that the wave pattern depended on the Froude number Fn (non-dimensional parameter of speed with reference to the ship’s length, $Fn = U/\sqrt{gL} $, where U is the ship speed, g the gravity acceleration and L the ship’s length) and the hull-form configuration, contrary to the theory of linear, dispersive waves. Other experimental results indicated that the variation of the wave-front angle of the bow wave is rather systematic as shown in Fig. 2. With an increase of apex angle (ship’s half entrance angle), the wave front angle increases, and the steepness of the wave slope obviously exceeds the limit of the linear wave system, as seen in Fig. 3.

Fig. 2.

Free-surface shock waves around models of different half-entrance angles 15°, 10°, 5° (from left to right) for Fn = 0.277; the free-surface shock wave is traced with dashed lines.

Fig. 3.

Results of detailed wave profile measurement for a wedge model of half-entrance angle 20°, Fd = 1.1 (Froude number with reference to draft length, $Fd = U/\sqrt{gd} $, where d is the draft length).⁴⁾

This non-linear wave system appearing in the near-field of a ship advancing steadily in deep water was subsequently named the “Free-Surface Shock Wave” (FSSW). Typical characteristics are: (1) the formation of lines of discontinuity; (2) the steepness and unsteadiness of the wave front; (3) satisfaction of the shock condition; (4) the systematic change in the angle of the wave front; (5) the non-dispersive nature of the propagation; and (6) the dissipation of wave energy into momentum loss in the wake.

These characteristics were exemplified through a series of physical experiments⁴⁾ in 1979. Though linear dispersive waves were also observed in the far-field, non-linear waves are dominant in the near-field. Thus, both linear dispersive and non-linear dissipative phenomena coexist.

It was assumed that a theoretical description of this non-linear wave phenomenon would only be possible using numerical analysis. Consequently, a special computational method would need to be developed, requiring considerable development time, and access to high-speed computers. Nonetheless, though understanding of the physical phenomenon was considered to be essential in itself in the process of actual ship hull-form design. The characteristics of Free-Surface Shock Waves are governed by the Froude number based on draft, Fd. When Fd > 1.0, the FSSW becomes dominant, with the presence of a strongly interacting breaking wave. To reduce the wave front angle of the FSSW at high Froude number, the entrance (apex) angle of the ship’s bow should be made small: a decreasing angle leads to a decrease in the angle of the shock front, and consequently a reduction in the wave resistance.

This design criterion is simple and straightforward, analogous to that for a supersonic body. The idea was first applied to the hull-form of the USUKI PIONEER (26,000 DWT bulk carrier), designed in 1983; see Fig. 4.

Fig. 4.

The USUKI PIONEER, a 26,000 DWT bulk carrier; the first ship designed with a thin, long-protrudent bulb in 1983.

The bulb equipped at the bow, as seen in Fig. 4, could reduce wave resistance by 10–20% in comparison with a ship with a conventional, bulbous bow design. In the 1980’s, a bow shape of this design was communicated to the world’s shipyards. The new bulb shape is especially useful for ships operating at speeds corresponding to a Froude number (based on ship length) Fn < 0.22.

Non-linear waves similar to the FSSW generated from the bow of the ship also occur near the stern, especially for high-speed ships, as illustrated in Fig. 5, which is an aerial photograph of the SURORECHIA-MARU passenger boat (Length = 100 m, Fn = 0.30). After a series of experimental works carried out by the group, it was found that the thin, long-protrudent bulb was also effective in the reduction of the non-linear stern waves. However, due to many reasons and space restrictions, the installation of such a stern bulb is limited to be in the proximity of the waterline, as shown in Fig. 6, and the number of ships fitted with this device was limited in comparison with those fitted with the thin, long-protrudent bulb at the bow. Nonetheless, there is about 5% energy saving when this Stern-End-Bulb (SEB) is fitted, if applied to car carriers, truck ferries and container carriers for which the draft at the stern is similar.^5),6)

Fig. 5.

Aerial picture of a ship for Fn ≒ 0.3; non-linearity of waves can be observed, especially around the stern.

Fig. 6.

Photo of a Stern-End-Bulb (SEB) for a car ferry.

The research work to develop a numerical method, based on the finite-difference approach, undertaken to compute and realize the free-surface flow, including non-linear wave phenomena, around a ship on a steady course, was begun in 1979. Four years were needed for the TUMMAC-IV (Tokyo University Modified Marker And Cell method version IV) code to reach maturity.^7),8)

The algorithm is based on the MAC (Marker And Cell) method.⁹⁾ The governing equations are the Euler equations of motion and the mass-conservation equation. These are written in finite-difference form, in accordance with established CFD techniques. Viscous effects are not taken into account. A rectangular grid system is employed, and the three components of velocity, and the pressure, are defined in each cell of a fixed rectangular grid system. The free-surface is tracked by the movement of markers located on the surface. The equations, in difference form, are solved as an initial value problem in a time-marching procedure. Initially, the ship is at rest, and it is gradually accelerated by increasing the horizontal velocity of the water at the inflow boundary until a given steady-state situation is reached. Most efforts were focused on the implementation of the body-boundary condition in the framework of the rectangular grid system. Using a special technique,⁸⁾ the free-slip condition could be satisfied in the body-boundary cells for an arbitrary solid-fluid configuration.

The first version of TUMMAC-IV was completed in 1983, and computed results compared against experimental data. Computed wave patterns agreed very well with the measured values, especially in the case in which wave-breaking was rare. The comparison of wave pattern for a tanker hull in ballast condition in Fig. 7 indicates that some discrepancy remains, most likely attributable to the fact that TUMMAC-IV cannot cope with breaking waves. Nonetheless, the agreement is satisfactory, and marks a significant first step.

Fig. 7.

Wave contours around the bow of a tanker model at ballast condition: measurements (left) and computed results using the TUMMAC-IV code (right).

The TUMMAC-IV code was also applied to the diffraction wave problem: that is, the interaction of the ship’s waves with the ocean waves. The two cases of a wedge-shaped bow and a tanker bow were computed. Results revealed that the diffraction wave shows phase-dependent variation of the wave formation, as shown in Fig. 8.¹⁰⁾

Fig. 8.

Computed waves and unsteady free-surface shock wave in regular incident waves.

Because the TUMMAC-IV code allows easy grid generation, and has proven accuracy, it was distributed to the major shipbuilding companies, and subsequently employed as a useful tool in ship hull-form design. However, the wave resistance, obtained by integrating the surface pressure distribution over the hull, did not correspond fully with the experimental measurements. This was thought to be due to the viscous effect not being taken into account. But, hull-form designers soon noted that the accuracy was satisfactory in regard to the relative wave resistance for different hull designs, and was thus useful in guiding hull-form design by successive improvement. The code was subsequently extended to a version applicable to a catamaran, in which the advantage of small-wave dissipation in the far-field allowed application to the wake-wash problem, as shown in Fig. 9.

Fig. 9.

Wave contours around one hull of a catamaran, computed using TUMMAC-IV. The solid and dashed contour lines correspond to the positive and negative values, respectively, with the interval of 0.0025 L. The Froude number is Fn = 0.645.

As a consequence of this work, several ship initial-design offices of the major Japanese shipbuilding companies have employed the TUMMAC design tool for more than 20 years.

One of the most non-linear, fluid-dynamical problems is wave breaking. Since CFD technology is based either on the Euler (non-viscous) or Navier-Stokes (viscous) equations, it can, in principle, deal with any problem of high non-linearity, although the resolution of the micro-mechanical and very fast phenomena will necessitate access to top-end computing power.

The modelling of wave-breaking in 2D was accomplished with a relatively small amount of research effort. The deforming free-surface configuration is represented by a succession of piece-wise linear segments. Then, not only can the overturning motion of the wave be simulated, but also its impingement on the free-surface below, as illustrated in Fig. 10.¹¹⁾

Fig. 10.

Simulation of 2D wave-breaking ahead of a steadily advancing rectangular body.

One application of the 2D wave-breaking problem is that of a circular cylinder placed horizontally in the vicinity of a free-surface. The complex flow field of the vortices generated, strongly interacting with the free-surface wave, could be realistically simulated, as shown in Fig. 11.¹²⁾

Fig. 11.

Contours of vorticity around a circular cylinder, set horizontally, and advancing steadily, near a free-surface, as computed by TUMMAC-V.

A quite different technique was necessary to model 3D wave-breaking motions, because the water surface cannot then be treated simply by a succession of linear segments. Therefore, based on the idea of the VOF (volume-of-fluid) method, a new technique was developed to impose the free-surface boundary conditions on a 3D interface of complex topology. We called our approach the Density-Function Method. It was revealed later that a similar method, known as the Level-Set method,¹³⁾ was being developed almost simultaneously in U.S.A.¹³⁾

The first application of the TUMMAC code, incorporating the Density-Function Method applied to the free-surface, was for the flow around a vertically placed rectangular cylinder.¹⁴⁾ When the Froude number is high, the fluid flow exhibits completely non-linear behavior, with wave breaking, vortex shedding, sprays and air-entrainment; all these phenomena can be seen in Fig. 12. Since the air flow above the water surface is simultaneously computed, some features of the associated spray and air-entrainment phenomena can be observed in Fig. 13, which shows the situation in side view. When the relevant physical parameters, such as density, Reynolds number and so on, are varied, a variety of two-phase flow phenomena can appear simultaneously: e.g. oil flow on a water layer, water flow over a liquefied sand bed (scouring problem), and so forth.

Fig. 12.

Contours of wave height around a vertical, rectangular cylinder (0.1 m long, positioned in a uniform stream 0.9 s after the start of acceleration.

Fig. 13.

Side view corresponding to Fig. 12, in a longitudinal plane at one side of the pillar.

The Density Function Technique was later employed in another code, based on the finite-volume method, written to treat viscous and wave-flow problems in the framework of a boundary-fitted curvilinear coordinate system. A typical simulation result is shown in Fig. 14, which shows the bow wave generated from a semi-planing boat.¹⁵⁾ The 3D breaking waves are well reproduced, and illustrate the difference in wave formation as a function of the Froude number. Some features of the FSSW around a wedge model are also simulated with the same free-surface condition, as described in Refs. 16, 17, and the discontinuous and dissipative features are displayed in Fig. 15.

Fig. 14.

Overturning wave features simulated with the density-function method: Fn = 1.0 (above) and Fn = 0.6 (below).

Fig. 15.

Velocity vectors on the free-surface; discontinuity of flow is present on the free-surface.

The approach of imposing the body-boundary condition in a rectangular grid system, which was developed in the context of the TUMMAC code, had considerable merits with respect to the ease of the grid generation process compared to the boundary-fitted grid system. Grid generation for boundary-fitted grid system inevitably demands more resources in terms of time and manpower, while the process can be semi-automated if carried out in the framework of a rectangular grid approach. In our Laboratory, the development of the TUMMAC code was considered complete, and efforts were then focused on the development of the WISDAM code, described in the next section. However, the case for employing a rectangular grid system appears not yet to be settled. Other research groups are currently developing similar techniques, with non-slip boundary conditions, in the framework of the rectangular grid system.^18),19) These approaches are widely used within the CFD research community under the name of the Immersed Boundary Method^20),21) or the Cartesian Grid Method. ¹⁹⁾ However, a critique of these parallel approaches lies outside of the scope of the present paper.

Since the TUMMAC code used free-slip boundary conditions on body surfaces within a rectangular grid system, it was unable to cope with the viscous flow problems. To overcome this restriction, the development of a new CFD code based on the boundary-fitted grid system was initiated in the 1980’s. The new code was named WISDAM (Wave vIScous flow Difference Accurate Method). In the same way as for the TUMMAC code, a time-marching simulation method was employed instead of the steady-state method in order to simulate transient or unsteady phenomena. The finite volume method was employed, and newly adopted to satisfy the conservation properties accurately. After the development of the first version, which always requires considerable manpower commitment, the second version, WISDAM-II, could finally cope with both free-surface and viscous-flow problems simultaneously.²²⁾ Since the issue of WISDAM-II, a number of updates have been developed to simulate non-linear free-surface phenomena together with turbulent flow around a ship. The first objective of the application of WISDAM was to calculate the drag reduction of ships on a steady, straight course in calm sea, and this was almost achieved by the fifth version, WISDAM-V. In all versions of the WISDAM code, the viscous effect is considered by employing turbulence models via solution of the Reynolds-Averaged Navier-Stokes (RANS) equations.

During the 1990’s, the research work of the group was focused on the ship motion, which requires new models to be implemented to take account of the conditions of moving coordinates on moving boundaries. A simulation method for ship maneuvering was developed in the framework of the WISDAM-V code using a moving-grid technique, and including models for the propulsion by the propeller, and the rudder forces.²³⁾ The trajectory of the maneuvering motion of an oil tanker computed using this version agreed well with data from a full-scale (real ship) test. Thus, the applicability of the WISDAM simulation approach for ship maneuvering was demonstrated.

As well as leading the research work of the group for the development of a variety of simulation techniques for merchant ships, the first author was also assigned as technical director of the syndicate of the Nippon Challenge America’s Cup for the challenge to the International America's Cup Class (IACC) sailing boat race in the year 2000. In order to win such races, the team must design boats which are the fastest in the world of that class. To this purpose, we decided to develop a new CFD technology for the simulation of sailing boat performance, CFD representing the most competitive advantage of the team over other, more conventional approaches. The acceptable error in the CFD simulations for the IACC sailing boats was less than 1% in terms of the drag in the heeled condition. So, the WISDAM-VII code was developed for the realization of the ship’s motion with six degrees of freedom.^24)–26) More than 200 boats were designed and simulated using the WISDAM-VII code, and 35 model ships were built and tested in the towing tanks at the Akishima Laboratory of Mitsui Zosen Inc. and at the University of Tokyo. The computed drag agreed satisfactorily with the experimental data, especially in terms of ranking of the drag. After the America’s Cup race in 2000, the WISDAM-VIII code was developed to predict the performance of a semi-planing boat with a transom stern in unsteady motion.¹⁵⁾

Over the last decade, new versions of the WISDAM code have been developed to compute a variety of motion problems, both for conventional merchant ships and high-speed passenger boats, with different hull configurations. For the fulfillment of the free-surface condition in the body-fixed grid system, the Density Function Method is introduced to cope with the large amplitude motion of ships. The non-linear free-surface phenomena are computed with the density function, and a computational grid fixed to the body (hull) allows the ship’s motion to be simulated without deforming the grid. The method was embodied in an improved version of the WISDAM-V code, which we re-named WISDAM-V motion, and successfully applied to a variety of ship motion simulations, e.g. in Refs. 27–29.

Although the WISDAM-V motion code can simulate ship motions in a variety of regular wave conditions, there are some difficulties in computing ship motions of large-amplitude, and in evaluating the added resistance (resistance increase due to the interaction with the incoming ocean waves). In addition, the WISDAM-V motion code was unable to simulate the motions of high-speed ships with practical hull forms, because it could not deal with complicated hull configurations, such as multi-hull forms. Thus, two new versions of the code, named WISDAM-X and WISDAM-XI, were developed to remove these limitations.

The WISDAM-X code was developed for the rigorous prediction of ship motion, and the estimation of the added resistance for conventional ships subjected to arbitrary wave conditions.^30),31) To model the interaction of a ship with incident waves, and to compute the resultant ship’s motions, an overlapping grid system is employed. The overall computational domain is divided into two solution domains. In each, the computational grid is constructed independently. The incident waves are generated by specifying the fluid velocity components and wave height explicitly at the inflow boundary of the (outer) solution domain. The motion of the ship is simultaneously solved by coupling the equations of the ship’s motion with the flow computation, thus providing the hydro-dynamical forces and moments. In the following Sections 7–8, some application cases of the WISDAM-X code are presented.

The WISDAM-XI method was specially developed for unsteady ship motion applied to high-speed, multi-hull vessels.³²⁾ For the treatment of the complicated configuration of a multi-hull vessel, a multi-block grid system is employed. Details are given in Section 9.

The WISDAM-X code is rather versatile, and can be applied to a variety of problems. One application is the prediction of the magnitude of the hull surface pressure caused by the interaction of the ship’s motion with the prevailing wave motion. Sometimes ships may receive damage as a consequence of this interaction, especially at the bow, a phenomenon called “slamming” in rough seas. Here, examination is made of the degree of accuracy of the computed pressures on the hull surface evaluated using WISDAM-X, and comparing with the experimental results. An oil tanker model SR221C is chosen for the simulation, since it has a typical hull form of today’s large oil tankers. The calculations were performed using an overlapping grid system consisting of inner and outer grids, as shown in Fig. 16. Simulations were made under head-wave conditions when the wave-length-to-ship-length ratio λ/L = 0.4–2.0, where λ is the wave-length and L is the ship length. Comparisons were made with the measured pressures at the locations on the hull surface shown in Fig. 17.

Fig. 16.

Close-up view of computational grid for the SR221C model.

Fig. 17.

Locations of pressure measurements on the SR221C model.

Time histories of the hull surface pressure contours are shown in Fig. 18 for the case of λ/L = 1.0. The non-dimensional incident wave amplitude ζ_A/L = 0.01. In this figure, the distributions of non-dimensional pressure ϕ ≡ P/(ρU²) + ρgz/(ρU²) on the hull surface (where P is the pressure and z is the vertical coordinate with origin z = 0 denoting the calm water level) are shown at intervals of 1/4 of the encounter wave period T_e together with the time histories of pressure at the measurement location P1. The pressure measured at the locations P1–9 is non-dimensionalized by ρgζ_A. From these figures, it can be clearly seen that large pressures are generated on the bow flare part of the model ship, supposed to be due to the occurrence of weak bow flare slamming phenomena.

Fig. 18.

Time histories of computed bow pressures and pressure contours on a SR221C tanker, λ/L = 1.0 and ζ_A/L = 0.01, where Te denotes encounter wave period. The black circles in the right figures indicate the instances of the left figures.

Comparisons were also made of the time histories of the pressures on the hull surface; results are shown in Fig. 19 for the cases λ/L = 0.6 and 1.2. As can be seen in these figures, computed results (shown by bold lines) are in good agreement with the experimental data (shown by dotted lines). The simulation results reproduce quite well the shape of the pressure variation with time: from a triangular shape to a more rounded shape with increase of the wave-length of the incident waves.

Fig. 19.

Comparisons of time variation of surface pressures on a SR221C tanker in regular head waves for λ/L = 0.6 (left) and 1.2 (right).

Comparisons of the amplitudes of surface pressure with the measured data are shown in Fig. 20 as a function of the non-dimensional incident wave-length λ/L. The agreement of the calculated results with the experimental data is quite good for all the points of pressure measurement, except for the case λ/L = 0.4.

Fig. 20.

Comparison of bow pressure amplitudes on a SR221C tanker at Fn = 0.150 in regular head waves of ζ_A/L = 0.01.

It is to be noted that the computed pressure agrees well with the experimental data at the fore end of the hull (points P1–P3), and that the magnitude of the surface pressure increases significantly in the region near the bow, contributing significantly to the added resistance and the local wave load magnitude on the hull. From these findings, we assert that the present CFD method is useful for practical applications. By using the pressure prediction method based on the present CFD simulations, high-performance ships of smaller resistance, and subsequently low fuel-consumption, together with sufficient structural strength against wave-induced forces, can be designed.

Another application of the WISDAM-X code is to compute the detailed structure and characteristics of the diffraction flow-field generated around a ship’s hull in the presence of ocean waves coming from different directions. Ships are meant to be operated in the presence of ocean waves, so the interaction phenomena of the hull with the waves, causing additional forces and moments, is of paramount interest with respect to the performance, safety and economic operation of the ship.

CFD simulations were conducted for the oil-tanker model SR221C in regular waves, over a range of wave directions from 180° (head wave) to 90° (beam wave) at Fn = 0.150. The length and amplitude of the incident waves were kept constant for all the cases at λ/L = 0.5 and ζ_A/L = 0.01. In these simulations, the ship is subject to heave, pitch, roll and surge motions.

Three time-history maps of computed wave-height (ζ) contours are shown in Fig. 21. Note that ζ is the wave height, and ζ_A is the amplitude of incident waves. These are shown at intervals of 1/4 of the encounter wave period T_e. The time-related variation of the formation of the diffraction wave is clearly discernible over one cycle of the wave encounter. From the figures, it is seen that diffraction of the incident waves become significant as the wave direction changes from head to beam. The interaction process can be clearly seen in the case of the oblique head wave condition at χ = 120° (where χ is the angle of the incident wave to the ship’s advancing direction). On the weather-side of the hull, wave reflection on the hull intensifies with decrease of angle of the incident wave, and very steep waves, with maximum heights more than two times that of the incident waves, are generated around the bow, and propagate towards the stern quartering direction.

Fig. 21.

Time histories of computed wave-height (ζ) contour maps around a SR221C tanker for Fn = 0.150 in waves of λ/L = 0.5, χ = 120°, 150°, 180°, and at ζ_A/L = 0.01. The interval between the contours is 0.001 L; contours of positive value are drawn in solid lines and those of negative value are drawn in dotted lines.

Close-up views of the evolution of the wave-height contours near the bow of the ship are shown in Fig. 22, for the same instants as those given in Fig. 21. It may be noted that the elevated wave formation and the depressed wave formation in the vicinity of the bow show very similar characteristics. For the case of an oblique head incident wave at χ = 120°, the contour maps at 1/16 T_e and 5/16 T_e resemble closely those at 9/16 T_e and 13/16 T_e respectively, when the sign of the wave elevation becomes reversed. This feature is also noted for the case of an oblique head wave at χ = 150°, as well as for the head wave case, χ = 180°, although the similarity then becomes somewhat more obscure. This implies that the location of the maximum slope, whether it be positive or negative, does not vary greatly. This feature of wave formation is very similar to those observed for unsteady FSSWs for a series of wedge models advancing in regular head waves,¹⁰⁾ and exemplifies the occurrence of unsteady FSSWs in oblique wave conditions.

Fig. 22.

Close-up view of time histories of computed wave-height (ζ) contour maps around a SR221C tanker for Fn = 0.150 in waves of λ/L = 0.5, χ = 120°, 150°, 180°, and at ζ_A/L = 0.01. The interval between the contours is 0.001 L; contours of positive value are drawn in solid lines and those of negative value are drawn in dotted lines.

Three sets of time-histories of the hull-surface pressure contours on the weather-side of the hull are shown in Fig. 23 for the cases χ = 120°, 150° and 180°. The instantaneous pressure distributions are given at intervals of 1/4 of the encounter wave period T_e, and at the same instants as those in Fig. 21. The time variations of the hull surface pressures are very noticeable over one cycle of the wave encounter. From the figures, it can be clearly seen that quite high values of pressure are generated on the bow flare part of the model ship under oblique wave incidence conditions. The generation of such large pressures is closely related to the deformation of the incident waves, as illustrated in the case of χ = 120° at 9/16 T_e.

Fig. 23.

Time histories of computed pressure distributions on the weather-side surface of a SR221C tanker hull for Fn = 0.150 in waves of λ/L = 0.5, χ = 120°, 150°, 180°, and at ζ_A/L = 0.01. The interval between the contours is 0.02; pressure contours of positive value are drawn in solid lines and those of negative value are drawn in dotted lines.

Since diffraction flow phenomena under oblique wave conditions cause an increase in added resistance, and is generally greater than those under other wave conditions, accurate evaluation of the phenomena by means of numerical simulation methods, as performed using the present version of the WISDAM-X code is invaluable in the design of superior hull forms, with low resistance, and consequently low fuel-consumption capability, under a variety of ocean wave conditions.

The motivation for the development of the WISDAM-XI code is based on the needs of the ship industry, especially the demand for a fast sea-transportation system, a goal that has increased substantially in importance over the past two decades. About 50% of currently operating fast ferries are catamarans, and they are satisfactorily making services. However, some shortcomings have been recognized in catamarans, such as the high acceleration force due to the large transverse restoring force with excessive transverse metacenter height (TKM). With this background, the focus of the group’s activities shifted to the trimaran concept, of which the TKM is between those of the monohull and the catamaran. In addition, several experimental studies and numerical resistance tests have revealed the superior property of the low-wave making resistance of the trimaran.

Accurate estimation of the wave-induced motion and loads on a trimaran vessel is essential for the design of fast ships of this type, combined with satisfactory ride comfort and structural strength. A full-scale measurement test was conducted for an actual ship: RV Triton,^33),34) constructed by Austal Ships in 2005.^35),36) The trimaran is 126 m long, and is capable of a speed of up to 40 knots with a deadweight of 500 tons. Several tank tests, and numerical simulations, were carried out to assess its seakeeping performance.^37)–40) However, these studies could not address the issue of large amplitude ship motion, including the slamming phenomenon, in which the bottom of the ship’s bow appears above water level in rough weather conditions. In order to prevent accidents, and damage to the vessel, it was necessary to predict the wave-induced motion of the ship, and the wave-impact loads on the hull, in a more quantitative way. Furthermore, it is essential that the ride-comfort property is evaluated quantitatively for such high-speed ships, and with a variety of hull configurations.

Therefore, further development of the WISDAM-XI code was undertaken to produce a numerical tool for predicting resistance and seakeeping performance, including the ride-comfort factor, for multi-hull vessels, by way of a CFD approach. The numerical method employed was based on the WISDAM-V code,²⁹⁾ extended to the overlapping grid system previously incorporated in the WISDAM-X version of the code,^31),41) as described in the previous section. For the simulation of multi-hull vessels, WISDAM-V motion was modified to incorporate a multi-block grid system approach, in order to deal with a multi-hull vessel configuration. The new code was named WISDAM-XI.

The governing equations are the three-dimensional, time-dependent, incompressible Reynolds-Averaged Navier-Stokes (RANS) equations, incorporating the Baldwin-Lomax turbulence model.⁴²⁾ The density-function method²⁹⁾ was used for modelling the free-surface deformation. The finite-volume method was used for the space discretization, and a boundary-fitted curvilinear coordinate system employed, with a staggered variable arrangement.⁹⁾ For the treatment of the ship motion, a body-fixed coordinate system was employed. The grid system fixed to the ship moved in accordance with the ship’s motion, as shown schematically in Fig. 24. The trajectory and attitude of the ship were defined in terms of a space-fixed coordinate system. The incident waves are assumed to be sinusoidal, in infinitely deep water, and generated by the velocities and wave heights given at the inflow boundary.

Fig. 24.

Moving, body-fixed coordinate system.

The validation of the WISDAM-XI code was first carried out in terms of the resistance problem. Three types of trimaran, with different side hull positions with respect to the main hull, were used for the validation tests. These revealed that the computed values of drag, trim, sinkage and the transverse wave profile agreed well with the measurements taken at the towing tank facility at the University of Tokyo.⁴³⁾

The second validation case involved the roll motion in beam waves under stationary conditions. Three types of ship — a monohull, a catamaran and a trimaran — were employed, as illustrated in Fig. 25. The principal particulars of the three designs (abbreviated as Mono, Cat and Tri, respectively) are shown in Table 1, where Lpp is the ship’s length, B the breadth, d the draft, and GM the metacentric height. The displacements of these ships (in tons) are also given in the table. The design speed of Mono is 30 knots, while that of Cat and Tri is 35 knots for practical designs, each with respect to the nominal engine power.

Fig. 25.

Hull shapes used for seakeeping performance simulation: Mono (left), Cat (middle) and Tri (right).

In this second validation case, the ship models are subject to four modes of motion – sway, heave, roll and pitch – but restricted in other modes of motion. The amplitude of the incident wave is 1.0 × 10⁻² Lpp for all the cases. The wave-length is set longer than half the ship’s length. In order to simulate wave propagation with a sufficient degree of accuracy, forty grids are allotted for one-wave length, and six grids are allotted in the vertical direction for the wave height. The total number of grid points is 1.4 × 10⁶ for Mono, 1.6 × 10⁶ for Cat and 2.9 × 10⁶ for Tri. The minimum grid spacing on the hull is set at 1.0 × 10⁻³ Lpp. Figure 26 (left) shows the whole computational domain for Tri using the multi-block grid generation technique; and the grids on the hull are shown in Fig. 26 (right).

Fig. 26.

Multi-block grid system for the trimaran model (left), and grid distribution on the hull and the waterline (right).

The computed roll amplitude is compared with measurements in Fig. 27 as the function of the non-dimensionalized wave length λ/Lpp. The roll amplitude ϕ_roll is non-dimensionalized by the wave number k and the wave amplitude ζ_A. It is noted that computed results agree well with the measured values, especially for the cases of Mono and Tri, while the roll amplitude is slightly overestimated.

Fig. 27.

Comparison of roll motion amplitude between computation and measurement.

Simulations were carried out for the Mono, Cat and Tri designs for the purpose of elucidating the seakeeping performance of multi-hull vessels. Ships were subject to sway, heave, roll and pitch motions, but were restricted in other modes of motion. The forward speed of the vessel was set at 30 knots for Mono, and 35 knots for Cat and Tri, which correspond to Froude numbers based on ship length (Fn) of 0.40, 0.44 and 0.40, respectively. The condition of the incident wave is listed in Table 2. The regular incident wave was introduced from five directions: head wave (χ = 180°), oblique head wave (χ = 150°), beam wave (χ = 90°), oblique follow wave (χ = 30°), and follow wave (χ = 0°).

Time-sequential images and time-histories of the ship’s motion for the Mono design are shown in Fig. 28. The non-dimensionalized time tU/Lpp is used for the horizontal axis, where t is time and U is the forward speed (m/s). In the case of an incident wave at 150°, a relatively small roll amplitude of 0.25° is recorded, and the observed mean heel angle is about −0.75°. The results for the Cat and Tri designs are shown in Figs. 29 and 30, respectively. Ship motion of large amplitude, and the noticeable deformation of the free surface, are to be noted in these results. Computed non-linear free surface features around the bulbous bow are shown in Fig. 28 (top), and it is to be noted that the bottom of the bow appears out of the water in Fig. 30 (top). The time history of the roll motion for the Cat design is shown in Fig. 29 (bottom), where the features of the (non-simple) harmonic motion may be observed, mostly due to the coupling effects of a number of motions.

Fig. 28.

Time-sequential drawing of the ship motion at every 1/4 encounter period (top), time history of incident wave elevation at the center of gravity, heave motion (middle), roll and pitch angle (bottom) for Mono design for waves of χ = 150° and λ/Lpp = 0.98.

Fig. 29.

Time-sequential drawing of the ship motion at every 1/4 encounter period (top), time history of incident wave elevation at the center of gravity, heave motion (middle), roll and pitch angle (bottom) for Cat design for waves of χ = 150° and λ/Lpp = 0.88.

Fig. 30.

Time-sequential drawing of the ship motion at every 1/4 encounter period (top), time history of incident wave elevation at the center of gravity, heave motion (middle), roll and pitch angle (bottom) for Tri design for waves of χ = 150° and λ/Lpp = 0.72.

In Figs. 31, 32 and 33, comparisons are made of the heave, roll, pitch amplitudes and the vertical acceleration for the three types of ship design. Only selected results for the cases of oblique head wave (χ = 150°) and head wave (χ = 180°) are reported here, while all the results are reported elsewhere.⁴³⁾ From the viewpoint of a comparison of ride-comfort performance, results are shown in dimensional units. Under head or oblique head wave conditions, the amplitudes of the heave and pitch motions of the Cat design, and the vertical acceleration, are significantly larger than those for the Mono and Tri designs, especially when the wave length is long (150 m). The vertical acceleration of the Cat design exceeds 0.2 g, while that of the Tri design is less than 0.1 g. Under follow (χ = 0°) or oblique follow (χ = 30°) wave conditions, there is little difference observed for the three vessel types. From these results, it may be safe to assert that the trimaran has some advantages over other two hull shapes with respect to the (heave, pitch and roll) motions and acceleration characteristics.

Fig. 31.

Amplitude of heave motion for an incident wave of ζ_A = 1.5 m, for χ = 150° (top) and χ = 180° (bottom).

Fig. 32.

Amplitude of pitch motion for an incident wave of ζ_A = 1.5 m, for χ = 150° (top) and χ = 180° (bottom).

Fig. 33.

Amplitude of vertical acceleration at the center of gravity for an incident wave of ζ_A = 1.5 m, for χ = 150° (top) and χ = 180° (bottom).

Validation of the WISDAM-XI code for the ship’s motion for the trimaran under incident waves is also reported in Ref. 32. When the wave conditions of the total ocean area of the ship’s service are taken into account, it would appear that a multi-hull vessel with higher ride-comfort and lower fuel-consumption is preferable over other design concepts, and can be designed with the help of the described CFD simulation method.

Thorough and rigorous understanding of the physical phenomena involved constitutes a sound basis for all advanced technologies. Before the existence of the non-linear waves named Free-Surface Shock Waves around ships was confirmed in 1979, ship hull-form design had been accomplished mostly in terms of physical experiments, which required considerable time, money and manpower input. Although linear theories for ship waves were abundant at that time, they were almost useless for hull-form design with the major shipbuilding companies around the world, because the waves had been typically modelled as linear dispersive ones, and the importance of their non-linear characteristics had not been fully recognized.

Once ship waves were understood to be a combination of linear and non-linear waves, the only available analysis method was detailed numerical computation, which is capable of treating non-linear waves and turbulent flow conditions. With this recognition, the development of the Computational Fluid Dynamic (CFD) techniques accelerated at the University of Tokyo, aided by the rapid progress in computer hardware. First, the new technology was applied to ships pursuing a steady, straight course in calm seas. But, after the first success the CFD technology, the method was extended to almost all fluid-dynamical problems, including simulation of the ship motion in the presence of incident waves. Since then, the technology has been communicated to the major shipbuilding companies over a period of more than 20 years, since 1983, and has contributed to the progress in the design for ships of high performance and low fuel-consumption.

The principal activities of the Ship Model Basin Laboratory at the University of Tokyo have been highlighted in this article. No attempt has been made at a complete review of the new technology of CFD in ship design around the world. This task is left to a future article in this series.

Hideaki Miyata was born in Matsuyama, Japan in 1948. He graduated from the Department of Naval Architecture, Faculty of Engineering of the University of Tokyo in 1970 and continued his study at the graduate school as a Masters student, and graduated in 1972. From 1972 to 1977, he worked as a ship designer at Ishikawajima-Harima Heavy Industry (IHI) Co. Ltd. He was invited to the Department of Naval Architecture, the University of Tokyo as an Assistant in 1977, and since then progressed as a Lecturer (1978), Associate Professor (1979–1994), and Professor (1994–2012). He retired from the university in 2012, and started work for Social System Design Co. Ltd. as its President. His major research field is naval hydrodynamics and the development of ship-design technology. He has also worked on the development of new ship systems, including the sailing boat for the America’s Cup Yacht Race, when he was Technical Director and Chief Designer from 1993 to 2000. He was assigned to be a Fellow of the World Innovation Foundation in 2002. He was awarded the Imperial Prize and Japan Academy Prize in 2011.

Hideo Orihara was born in Tokyo, Japan in 1973. He graduated from the Department of Naval Architecture and Ocean Engineering, Faculty of Engineering of Osaka University in 1995. He continued his studies at the graduate school of the University of Tokyo as a Masters student, and graduated from the Department of Naval Architecture and Ocean Engineering in 1997. After graduation, he worked for Hitachi Zosen Co. as a research engineer from 1997 to 2002. Then, he worked in Universal Shipbuilding Co. as a research engineer from 2002 to 2012. During this period, he studied as a part-time student on the doctor course in the same department as the master course at the University of Tokyo, and received a Ph. D. degree in 2005. Since 2013, he has been working for Japan Marine United Co. as a research engineer. His major research field is naval hydrodynamics and the development of environmentally friendly ship technology.

Yohei Sato was born in Yokohama, Japan in 1972. He graduated from the Department of Biomaterial Sciences, Faculty of Agriculture of the University of Tokyo in 1995. After, he worked for Tokio Marine Fire Insurance Co., Ltd., and studied as a Masters student at the Department of Naval Architecture and Ocean Engineering, Faculty of Engineering of the University of Tokyo, graduating in 1999. After graduation, he worked for Asahi Glass Co., Ltd. as a research engineer for three years. From 2002 to 2009, he worked for the Center for CFD Research at the National Maritime Research Institute. During this period, he studied as a part-time student on a doctoral course at the University of Tokyo, and received his Ph. D. degree in 2006. Since 2009, he has been working in the Thermal Hydraulics Laboratory of the Paul Scherrer Institute (Switzerland) as a scientist. His major research field is Computational Fluid Dynamics, especially for free-surface flows and phase-change phenomena.

Finally, the first author would like to express sincere thanks to the large number of students who worked diligently on the development of the many versions of the TUMMAC and WISDAM CFD codes produced at the University of Tokyo over the last 30 years.

• 1)   Inui,  T.,  Kajitani,  H. and  Miyata,  H. (1979) Experimental investigations on the wave making in the near-field of ships. J. Kansai Soc. Naval Arch. Japan 173, 95–107.
• 2)   Miyata,  H.,  Inui,  T. and  Kajitani,  H. (1980) Free surface shock waves around ships and their effects on ship resistance. J. Soc. Naval Arch. Japan 147, 1–9.
• 3)   Kawamura,  N.,  Kajitani,  H.,  Miyata,  H. and  Tsuchiya,  Y. (1980) Experimental investigation on the resistance component due to free surface shock waves on series ship models. J. Kansai Soc. Naval Arch. Japan 179, 45–55.
• 4)  Miyata, H. and Inui, T. (1984) Nonlinear ship waves. In Advances in Applied Mechanics, (eds. Hutchinson, W.J. and Wu, Y.T.). Academic Press, New York, pp. 215–288.
• 5)   Miyata,  H.,  Tsuchiya,  Y. and  Inui,  T. (1980) Resistance reduction by stern-end-bulb: First report. J. Soc. Naval Arch. Japan 148, 10–16.
• 6)   Miyata,  H.,  Tsuchiya,  Y. and  Inui,  T. (1981) Resistance reduction by stern-end-bulb: Second report. J. Soc. Naval Arch. Japan 149, 1–10.
• 7)   Miyata,  H. and  Nishimura,  S. (1985) Finite-difference simulation of nonlinear ship waves. J. Fluid Mech. 157, 327–357.
• 8)   Miyata,  H.,  Nishimura,  S. and  Masuko,  A. (1985) Finite difference simulation of nonlinear waves generated by ships of arbitrary three-dimensional configuration. J. Comput. Phys. 60, 391–436.
• 9)   Harlow,  F.H. and  Welch,  J.E. (1965) Numerical calculation of time-dependent viscous incompressible flow of fluid with free surface. Phys. Fluids 8, 2182–2189.
• 10)   Miyata,  H.,  Kanai,  M.,  Yoshiyasu,  N. and  Furuno,  Y. (1990) Diffraction waves about an advancing wedge model in deep water. J. Ship Research 34, 105–122.
• 11)   Miyata,  H. (1986) Finite-difference simulation of breaking waves. J. Comput. Phys. 65, 179–214.
• 12)   Akimoto,  H.,  Sugihara,  M. and  Miyata,  H. (1991) Vortex motions and forces about a horizontal cylinder advancing beneath the waves. J. Soc. Naval Arch. Japan 170, 253–263.
• 13)   Sussman,  M.,  Smereka,  P. and  Osher,  S. (1994) A level set approach for computing solutions to incompressible two-phase flow. J. Comput. Phys. 114, 146–159.
• 14)   Miyata,  H.,  Katsumata,  M.,  Lee,  Y.G. and  Kajitani,  H. (1988) A finite-difference simulation method for strongly interacting two-layer flow. J. Soc. Naval Arch. Japan 163, 1–16.
• 15)   Orihara,  H. and  Miyata,  H. (2000) Numerical simulation method for flows about a semi-planing boat with a transom stern. J. Ship Research 44, 170–185.
• 16)   Kanai,  A. and  Miyata,  H. (1995) Elucidation of the structure of free surface shock waves about a wedge model by finite-difference method. J. Soc. Naval Arch. Japan 177, 147–159.
• 17)   Kanai,  A. and  Miyata,  H. (1996) Numerical analysis of structure of free-surface shock wave about a wedge model. J. Ship Research 40, 278–287.
• 18)   Udaykumar,  H.S.,  Shyy,  W. and  Rao,  M.M. (1996) Elafint: A mixed Eulerian–Lagrangian method for fluid flows with complex and moving boundaries. Int. J. Numer. Methods Fluids 22, 691–712.
• 19)   Ye,  T.,  Mittal,  R.,  Udaykumar,  H.S. and  Shyy,  W. (1999) An accurate Cartesian grid method for viscous incompressible flows with complex immersed boundaries. J. Comput. Phys. 156, 209–240.
• 20)   Yang,  J. and  Stern,  F. (2009) Sharp interface immersed-boundary/level-set method for wave–body interactions. J. Comput. Phys. 228, 6590–6616.
• 21)   Yang,  J. and  Stern,  F. (2012) A simple and efficient direct forcing immersed boundary framework for fluid–structure interactions. J. Comput. Phys. 231, 5029–5061.
• 22)   Miyata,  H.,  Sato,  T. and  Baba,  N. (1987) Difference solution of a viscous flow with free-surface wave about an advancing ship. J. Comput. Phys. 72, 393–421.
• 23)  Omori, T. (1994) Ship maneuvering simulations by a computational fluid dynamical method. Ph.D. thesis, University of Tokyo (in Japanese).
• 24)   Akimoto,  H. (1997) Finite-volume simulation of the flow around a sailing boat with unsteady motion. J. Soc. Naval Arch. Japan 181, 35–44.
• 25)   Akimoto,  H. and  Miyata,  H. (2002) Finite-volume simulation method to predict the performance of a sailing boat. J. Mar. Sci. Technol. 7, 31–42.
• 26)  Miyata, H. (2002) Towards virtual reality by computational physics. In 24th Symp. Naval Hydro., Fukuoka, Japan.
• 27)  Sato, T., Izumi, K. and Miyata, H. (1998) Numerical Simulation of Maneuvering Motion. In 22nd Symp. Naval Hydro., Washington D.C., U.S.A.
• 28)   Takada,  N. and  Miyata,  H. (1998) CFD Simulation of 3-dimentional motion of a vehicle with movable wings: Application to the keel of a racing yacht. J. Soc. Naval Arch. Japan 184, 37–45.
• 29)   Sato,  Y.,  Miyata,  H. and  Sato,  T. (1999) CFD simulation of 3-dimensional motion of a ship in waves: application to an advancing ship in regular heading waves. J. Mar. Sci. Technol. 4, 108–116.
• 30)  Orihara, H. and Miyata, H. (2003) A numerical method for arbitrary ship motions in arbitrary wave conditions using overlapping grid system. In 8th Int. Conf. Numer. Ship Hydro., Busan, Korea.
• 31)   Orihara,  H. and  Miyata,  H. (2003) Evaluation of added resistance in regular incident waves by computational fluid dynamics motion simulation using an overlapping grid system. J. Mar. Sci. Technol. 8, 47–60.
• 32)  Sato, Y., Uzawa, K. and Miyata, H. (2007) Validation of motion prediction method for trimaran vessels. In 9th Int. Conf. Numer. Ship Hydro., Ann Arbor, U.S.A.
• 33)  Baldwin, C., Kiddy, J., Salter, T., Chen, P. and Niemczuk, J. (2002) Fiber optic structural health monitoring system: Rough sea trials of the RV Triton. In Ocean’s 2002 Conference and Exhibition, Mississippi, U.S.A, pp. 1807–1814.
• 34)  Renilson, M.R., Scrace, R.J., Johnson, M. and Richardson, C. (2004) Trials to measure the hydrodynamic performance of RV Triton. In Int. Conf. Design & Operation of Trimaran Ships, London, UK, pp. 5–18.
• 35)   Armstrong,  N.A. (2004) Coming soon to a port near you: The 126m trimaran. Naval Architect. 66–78.
• 36)   Anderson,  G. (2005) Benchijigua Express: World’s first ro-pax trimaran in service. Naval Architect. 42–48.
• 37)   Fang,  M.-C. and  Chen,  T.-Y. (2008) A parametric study of wave loads on trimaran ships traveling in waves. Ocean Eng. 35, 749–762.
• 38)   Maki,  K.J.,  Doctors,  L.J.,  Scher,  R.M.,  Wilson,  W.M.,  Rhee,  S.H.,  Troesch,  A.W. and  Beck,  R.F. (2009) Conceptual design and hydrodynamic analysis of a high-speed sealift adjustable-length trimaran. Trans. - Soc. Nav. Archit. Mar. Eng. 116, 16–39.
• 39)   Jia,  J.B.,  Zong,  Z. and  Shi,  H.Q. (2009) Model experiments of a trimaran with transom stern. Int. Shipbuild. Prog. 56, 119–133.
• 40)   Wu,  C.S.,  Zhou,  D.C.,  Gao,  L. and  Miao,  Q.M. (2011) CFD computation of ship motions and added resistance for a high speed trimaran in regular head waves. Int. J. Naval Archi. Ocean Eng. 3, 105–110.
• 41)   Orihara,  H. (2011) Comparison of CFD simulations with experimental data for a tanker model advancing in waves. Int. J. Naval Archi. Ocean Eng. 3, 1–8.
• 42)  Baldwin, B. and Lomax, H. (1978) Thin-layer approximation and algebraic model for separated turbulent flows. In 16th Aerospace Sciences Meeting, American Institute of Aeronautics and Astronautics, AIAA Paper 78–257.
• 43)  Sato, Y., Orihara, H. and Miyata, H. (2006) Practical application of two CFD codes for ship motions in arbitrary waves. In 26th Symp. Naval Hydro., Rome, Italy.