Below are the first 10 and last 10 pages of uncorrected machine-read text (when available) of this chapter, followed by the top 30 algorithmically extracted key phrases from the chapter as a whole.

Intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text on the opening pages of each chapter.
Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.

Do not use for reproduction, copying, pasting, or reading; exclusively for search engines.

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
SESSION 1
WAVY/FREE-SURFACE FLOW: PANEL METHODS 1

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
This page in the original is blank.

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fully Nonlinear Water Wave Computations Using the Desingularized Method
R.F.Beck, Y.Cao, and T.-H.Lee
(University of Michigan, USA)
ABSTRACT
The use of Euler-Lagrange time stepping methods to solve numerically fully nonlinear marine hydrodynamic problems is discussed. The mixed boundary value problem that arises at each time step is solved using a desingularized approach. In this approach the singularities generating the flow field are outside the fluid domain. This allows the singularity distribution to be replaced by simple isolated singularities with a resultant reduction in computational complexity and computer time.
Various examples of the use of the method are presented including two-dimensional water wave problems, the added mass and damping due to sinusoidal oscillations of two-dimensional and axial symmetric bodies, and the wave pattern and wave resistance of a Wigley hull moving at constant forward speed. The results show excellent agreement with other published computations and good agreement with experiments.
INTRODUCTION
In recent years, there has been significant progress in solving marine hydrodynamic problems with fully nonlinear free surface boundary conditions. While some of the work has incorporated fluid viscosity and surface tension, the majority has assumed that the fluid is incompressible and inviscid and that the flow is irrotational. These assumptions lead to a boundary value problem for which the Laplace equation is valid everywhere in the fluid domain. This does not imply that viscous effects are unimportant and in certain situations they may dominate. However, the inviscid problem is an order of magnitude easier to solve and thus has been the basis of much of the research. Throughout this paper it is assumed that the fluid is incompressible and inviscid and that the flow is irrotational.
Fully nonlinear computations can be performed by a variety of methods. For steady forward speed, an iterative procedure can be used to satisfy the nonlinear free surface boundary conditions. In the work of Jensen, et al. (1989), and more recently Kim and Lucas (1990, 1992) and Raven (1992), a series of linearized boundary value problems based on the solution to the previous iteration and satisfied on the deformed free surface of that iteration are solved. The iteration is continued until convergence to the complete, nonlinear solution is obtained. Bai et al. (1992) also use an iterative procedure but employ a localized finite-element method to solve the boundary value problem at each iteration. For steady problems the iteration procedure may converge to the fully nonlinear solution faster than a time stepping method.
Longuet-Higgins and Cokelet (1976) first introduced the mixed Euler-Lagrange method for solving two-dimensional fully nonlinear water wave problems. This time stepping procedure requires two major tasks at each time step: the linear field equation is solved in an Eulerian frame; then the fully nonlinear boundary conditions are used to track individual Lagrangian points on the free surface to update their position and potential values. Variations of this method have been applied to a wide variety of two- and three-dimensional problems. Most of the variations involve how the linear field equation is solved at each time step. Among the researchers who applied the method to two-dimensional problems are Faltinsen (1977), Vinje and Brevig (1981), Baker et al. (1982) and more recently Grosenbaugh and Yeung (1988), Cointe et al. (1990) and Saubestre (1990). In three dimensions, the computations become much more difficult because of the large number of unknowns that are required. Results have been obtained by a number of researchers, including Lin et al. (1984), Dommermuth and Yue (1987), Kang and Gong (1990), Zhou and Gu (1990), Cao (1991), Cao et al.

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
(1990, 1991, 1992), and Lee (1992).
The successful implementation of an Euler-Lagrange method to solve fully nonlinear water wave problems requires a fast and accurate method to solve the mixed boundary value problem that results at each time step and a stable time stepping scheme. In our research, the mixed boundary value problem is solved using a desingularized boundary integral method. Similar to conventional boundary integral methods, it reformulates the boundary value problem into a boundary integral equation. The difference is that the desingularized method separates the integration and control surfaces, resulting in nonsingular integrals. The solution is constructed by integrating a distribution of fundamental singularities over a surface (the integration surface) outside the fluid domain. The integral equation for the unknown distribution is obtained by satisfying the boundary conditions on the surface (the control surface) surrounding the fluid domain.
A variety of problems have been solved by the desingularized approach. Webster (1975) was probably the first to apply the technique to panel methods. He used triangular patches of linearly distributed sources “submerged” within the body surface to study the steady flow past an arbitrary three-dimensional body. Schultz and Hong (1989) showed the effectiveness and accuracy of the desingularized method for two-dimensional potential flow problems. Cao et al. (1991) gave convergence rates and error limits for simple three-dimensional flows including a source-sink pair traveling below a free surface.
The desingularized method has been successfully applied to several free surface problems. Cao et al. (1993a) have investigated the formation of solitons propagating ahead of a disturbance moving in shallow water near the critical Froude number. The wave resistance, lift force and pitch moment acting on a submerged spheroid traveling beneath the free surface were computed by Bertram et al. (1991) and compared to linear theory results. Using an iterative technique, the wave resistance was computed by Jensen et al. (1986, 1989) and Raven (1992). The desingularized method has been applied to ship motion problems by Bertram (1990) and Lee (1992).
There are several computational advantages to the desingularized method. Because of the desingularization, the kernels of the integral equation are no longer singular and special care is not required to evaluate integrals over the panels. Simple numerical quadratures can be used to greatly reduce the computational effort, particularly by avoiding transcendental functions. In fact, for the source distribution method the distributed sources may be replaced by simple isolated sources. Higher order singularities such as dipoles can easily be incorporated. Isolated Rankine sources also allow the direct computation of the induced velocities on the free surface without further numerical integration or differentiation. The resulting code does not require any special logic and is thus easily vectorized. We have not yet installed the code on a parallel machine but the algorithm is straight forward and should not cause any difficulties. At present, the method is O(N2), but by using multipole expansions it could be reduced to an O(N lnN) method.
In the next section, the desingularized approach using isolated sources to solve fully nonlinear marine hydrodynamic problems will be discussed. The solution is developed in the time domain starting from rest. Following the theoretical development, some of the details of the numerical methods will be presented. Finally, numerical results will be presented and compared to experimental results where available.
FULLY NONLINEAR PROBLEM FORMULATION
An ideal, incompressible fluid is assumed and surface tension is neglected. The problem is started from rest so that the flow remains irrotational. This implies the existence of a velocity potential such that the fluid velocity is given by its gradient and the governing equation in the fluid domain is the Laplace equation.
A coordinate system Oxyz translating in the negative x direction relative to a space fixed frame is used. The time dependent velocity of translation is given by Uo(t). The Oxyz axis system is chosen such that the z=0 plane corresponds to the calm water level and z is positive upwards. The total velocity potential of the flow can then be expressed as
1
where is the perturbation potential. Both Φ and satisfy the Laplace equation
(2)
Boundary conditions must be applied on all surfaces surrounding the fluid domain: the free surface (SF), the body surface (SH), the bottom (SB) and the surrounding surface at infinity (S∞). A kinematic body boundary condition is applied on the instantaneous position of the body wetted surface:
(3)

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
where n=(n1,n2,n3) is the unit normal vector to the surface (out of the fluid domain) and VH is the velocity of the body relative to the moving coordinate system. The subscripts 1,2,3 refer to the x,y, and z axis directions respectively. There is also a kinematic condition applied on the bottom:
(4)
where VB is the velocity of the bottom relative to the Oxyz system. For an infinitely deep ocean equation (4) reduces to
(5)
Finite depth will increase the computational time because of the additional unknowns necessary to meet the bottom boundary condition but there is no increase in computational difficulty. In fact, the flatness of the bottom is immaterial. The only overhead relative to a flat bottom is computing the necessary geometrical parameters of a nonflat bottom. This contrasts with the typical Green function approach where a finite depth Green function is significantly harder to compute than an infinite depth Green function and a nonflat bottom can not in general be accommodated.
On the instantaneous free surface both the kinematic and dynamic conditions must be satisfied. The kinematic condition is
(6)
where z=η (x,y,t) is the free surface elevation. The dynamic condition requires that the pressure everywhere on the free surface equals the ambient pressure, Pa. The ambient pressure is assumed known and may be a function of space and time. Normally it would be set equal to zero. Using Bernoulli 's equation the dynamic condition becomes:
(7)
where ρ is the fluid density and g the gravitational acceleration.
Because we are solving an initial value problem with no incident waves, the fluid disturbance must vanish at infinity:
(8)
In addition, the initial values of the potential and free surface elevation must be specified such that
(9)
In the Euler-Lagrange method a time stepping procedure is used in which a boundary value problem is solved at each time step. At each step, the value of the potential is given on the free surface (a Dirichlet boundary condition) and the value of the normal derivative of the potential (a Neumann boundary condition) is known on the body surface and bottom surface. The potential and its normal derivative are updated at the end of each time step. The free surface potential and elevation are determined by integrating with respect to time (or time marching) the free surface boundary conditions. In our calculations, a Runge-Kutta-Fehlberg technique is used to do the time stepping. The body and bottom boundary conditions are prescribed for the forced motion problem. In the case of a freely floating body, the equations of motion must be integrated with respect to time in a manner similar to the free surface conditions.
On the free surface, the kinematic condition is used to time step the free surface elevation and the dynamic condition is used to update the potential. Many different approaches are possible to time march the free surface boundary conditions. The most common is a material node approach in which the nodes or collocation points follow the individual fluid particles. Another technique is to prescribe the horizontal movement of the node but allow the node to follow the vertical displacement of the free surface. The prescribed movement may be zero such that the node locations remain fixed in the x-y plane. Depending on the problem, one of the techniques may be easier to apply than the others.
It is convenient to rewrite the free surface boundary conditions, equations (6) and (7), in terms of the time derivative of a point moving with a prescribed velocity v relative to the Oxyz coordinate system. By adding to both sides of (6) and to both sides of (7) and after some algebraic manipulation, the kinematic and dynamic conditions can be put in the form
(10)
and
(11)

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
where
(12)
is the time derivative following the moving node. This is similar to the usual material derivative of fluid mechanics except the velocity is the given v rather than the fluid velocity.
If v is set equal to the node follows a prescribed path with velocity (U(t), V(t)) in the x-y plane and moves vertically with the free surface. Setting results in the x-y locations of the nodes remaining fixed in the Oxyz coordinate system and equations (10) and (11) reduce to:
(13)
and
(14)
In the material node approach, the velocity is set equal to the fluid velocity such that , resulting in
(15)
and
(16)
where XF=(xF(t),yF(t),zF(t)) is the position vector of a fluid particle on the free surface and
(17)
is the usual material derivative.
The form of the free surface boundary conditions given by the above equations allows the value of the elevation and potential to be stepped forward in time. The left hand sides of equations (10)–(16) are the derivatives with respect to time of the potential and wave elevation moving with the node. The quantities on the right hand side are all known at each time step; the spatial gradient of the potential can be determined analytically after solving the boundary value problem and the wave elevation is known. The difficulty is the gradient of the free surface elevation ( in equations (10) and (13)) that must be evaluated numerically. This leads to increased computer time and numerical inaccuracies. However, this term is only needed in the prescribed horizontal node movement approach. In the material node approach no extra derivatives need to be evaluated and this probably explains why this is the approach most often used. With material nodes one must always be concerned that the nodes do not penetrate the body surface between time steps since they are unconstrained. In zero forward speed problems, material nodes or fixed nodes seem to be the most appropriate. In problems with forward speed, the material node approach has difficulties near the body because nodes tend to pile up near the bow and stern stagnation regions. The prescribed horizontal node movement approach does not have this difficulty since the node movement is constrained. An appropriate choice of v is one which parallels the body waterline and is close to . In this case, the contribution of the term to the right hand side of (10) will be small and numerical inaccuracies will be minimized. Consequently, fast, simple numerical derivatives can be used to evaluate the term.
At each time step a mixed boundary value problem must be solved; the potential is given on the free surface and the normal derivative of the potential is known on the body surface and the bottom. In terms of a desingularized source distribution, the potential at any point in the fluid domain is given by
(18)
where Ω is the integration surface outside the fluid domain.
Applying the boundary conditions, the integral equations that must be solved to determine the unknown source strengths σ(xs) are
(19)

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
and
(20)
where xs=a point on the integration surface, Ω
xc=a point on the real boundary
= the given potential value at xc
Гd=surface on which is given
= the given normal velocity at xc
Гn= surface on which is given
The hydrodynamic forces acting on the body due to a prescribed motion are found by integrating the pressure over the instantaneous wetted surface. The generalized force acting on the body in the jth direction is thus given by:
(21)
where nj is the generalized unit normal into the hull defined as
(n1, n2, n3) = n
(n4, n5, n6) = r×n
(22)
n=unit normal to body surface (out of fluid)
=body axis system
and j=1,2,3 corresponds to the directions of the axis respectively.
The pressure in the moving coordinate system is given by Bernoulli 's equation:
(23)
where is the time derivation of the potential following a moving node on the body and v is the velocity of the node relative to the Oxyz system.
NUMERICAL TECHNIQUES
In the usual manner, the integrals may be discretized to form a system of linear equations to be solved at each time step. In the desingularized method, the source distribution is outside the fluid domain so that the source points never coincide with
Figure 1: Schematic of source and node locations
the collocation or node points and the integrals are nonsingular. In addition, because of the desingularization we can use simple isolated sources, rather than a distribution and obtain the equivalent accuracy. This greatly reduces the complexity of the form of the influence coefficients that make up the elements of the kernel matrix. As shown in figure 1, nodes (the collocation points) are distributed on the free surface and body surface. The isolated sources are distributed a small distance above each of the nodes. The nondimensional desingularized distance is given by
(24)
where ℓd and ν are constants to be chosen by the user and Dm is a measure of the local mesh size (typically the square root of the panel area in three-dimensional problems). Cao et al. (1991) found values of ℓd=1.0 and ν=.5 to be about optimum. The accuracy and convergence of the solutions are relatively insensitive to the choices of ℓd and ν.
At each time step the mixed boundary value problem that must be solved is defined by equations (19) and (20). Using an array of isolated simple sources above the free surface and inside the body surface to construct the velocity potential, equation (18) becomes
(25)
where NF is the number of the isolated sources above the free surface,
is the location of the jth source and
is the strength of the jth source at

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Similarly,
NB, and are the number, location and strength of the sources inside the body surface.
The integral equations (19) and (20) are satisfied at the nodes on the free surface and body surface such that
(26)
with N=NF+NB being the total number of unknowns. Equations (26) can also be written in matrix form as
A · Σ=B (27)
or
(28)
where
i=1…NF, j=1…NF
i=1…NF, j=(NF+1)…N
i=1(NF+1)…N, j=1…NF
i=(NF+1)…N,
j=(NF+1)…N
i=1…NF
i=(NF+1)…N
i=1…NF
i=1(NF+1)…N
where and are the node points on the free surface and body surface respectively.
Equations (27) are a system of (NF+NB) equations for the unknown source strength Σ. Once Σ is known, the velocity potential and fluid velocity can be analytically calculated anywhere in the fluid domain.
Equations (27) can be solved either by a direct or iterative solver depending on the size of the matrix. For most of the two-dimensional problems considered in this paper, we use a LU decomposition solver. For the three-dimensional problems, an iterative solver GMRES is used. For problems where different grid spacings are used on different parts of the domain boundary, solving equation (27) as a whole system may not be efficient. A domain-decomposition technique is often better. In this technique Σ(1) and Σ(2) are determined separately through an iterative procedure. Rewrite equation (28) as
A(11) · Σ(1)=B(1)−A(12) · Σ(2) (29)
A(22) · Σ(2)=B(2)−A(21)· Σ(1) (30)
Σ(1) is determined by solving (29) with Σ(2) known (or guessed), then Σ(1) is substituted into (30) to determine Σ(2). The new Σ(2) is then used in (29) to update Σ(1). The procedure repeats until both equations are satisfied to the given error tolerance.
After Σ is determined, the fluid velocity on the free surface can be computed so that the right hand sides of the free surface conditions (10) and (11) can be found. The free surface elevation and the potential are then updated by a time stepping procedure. In the free surface updating, if material nodes are used, the spatial derivatives of the free surface elevation, are not required. The use of generalized nodes moving with prescribed velocity v requires the spatial derivatives of η and central-differencing or cubic splines are used. The time-stepping of the free surface is done by using a fourth-order Runge-Kutta-Fehlberg scheme.

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
For forced body motion problems, the force acting on the body does not enter the time-stepping procedure. Therefore, the calculation of the forces acting on the body can be post-processed. Bernoulli 's equation (23) is used to compute the pressure on the body and the pressure is then integrated over the body surface to determine the forces. The pressure calculation requires the time-derivative of the potential on the body surface which is determined using central-differencing.
NUMERICAL RESULTS
In this section we shall present numerical results for a variety of cases including surface piercing bodies and a two-dimensional wave tank. The wave tank is useful to illustrate several properties of the calculations. Two-dimensional results are presented for a rectangular box in heave. Three-dimensional results are shown for axially symmetric cylinders and a Wigley hull.
Numerical Wave Tank
The two-dimensional wave tank consists of side walls, a wavemaker at one end and an absorbing beach at the other and a bottom at depth h. We have used several types of wavemakers including a paddle, a plunging wedge, and a variable pressure distribution on the free surface. The most effective beaches (cf. Cao, et al. 1993b) seem to be an energy absorbing layer on the free surface proportional to either or .
Figure 2 shows the waves generated by a pneumatic sinusoidal wavemaker started from rest. The tank has walls at each end and is .625 wavelengths deep. Two lines are shown: the dash line was computed using the traditional material nodes and the solid line is for fixed horizontal nodes. As can be seen the agreement between the two methods of calculation is excellent.
Figure 3 presents the comparison between two and three dimensional calculations. In this case, the wave tank is 12 units long, 4 units deep and the wavelength is 7.5 units. For the three-dimensional calculations the tank was 6 units wide. For the two-dimensional computations the fundamental singularity is ln(r) while in the three-dimensional calculations it is 1/r. As can be seen, the agreement is again very good; it is not perfect because of the panelization of the walls in the three-dimensional calculations. As more panels are used on the walls, the agreement continues to improve. For figure 3, only 8 panels were used at each section to divide the wall from the free surface to the bottom and the tank was less than one wavelength wide. The required number of panels needed to represent the walls is important because for a given total number of panels, the more panels on the walls the less will be available to represent the free surface and the body.
The conservation of energy in the wave tank with a pneumatic wavemaker and a beach proportional to is shown in figure 4. The solid line represents the energy build up with respect to time as the waves are made. The pneumatic wavemaker is forced to oscillate through one cycle with a period of . Thus, we see that the energy in the tank, ET , builds up until and then remains constant until the waves start to be absorbed by the beach. EA represents the energy absorbed by the beach and it starts to build up at approximately . The sum of ET+EA is the energy in the tank plus that absorbed by the beach. As expected the total energy remains constant. At the end of the run some energy remains in the tank because of wave reflection from the beach. The energy in the tank will decrease further if the simulation is allowed to run longer, but it would take a long time for the reflected waves to travel down to the wall at the wavemaker end and return back to the beach. The ratio of the energy remaining in the tank to the total energy is a measure of the reflection from the beach. This particular beach has not been optimized and it appears to absorb about 90% of the energy. One of the problems of fully nonlinear calculations is to find a means of measuring the effectiveness of a beach. The usual reflection coefficient is only valid for linear waves. We have tried several different measures, but none have been very satisfactory.
Dommermuth et al. (1988) showed excellent agreement between their fully nonlinear mixed Euler-Lagrange computations and experimental results. The experiments and numerical calculations were for a piston-type wavemaker programmed to produce wave breaking at mid-tank. To allow comparisons, the Fourier amplitudes and phases of the time history of the wavemaker's velocity are given in their paper. Figure 5 compares the computed wave elevation at five different locations with those presented in Dommermuth et al. (1988, figure 2). To within the accuracy of the graphics, the agreement is excellent. For our computations, the wave tank was 20 units long and 1 unit deep. A absorbing beach is placed at the far end. 501 nodes were placed on the free surface with 26 nodes on the piston wavemaker and 26 nodes on the far wall. The bottom was accounted for by using image sources for each desingularized source. Material nodes were used and the time step size was . A fourth-order Runga-Kutta-Fehlberg time stepping scheme is used. The running time on a CRAY-YMP was 1.6 hours.

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Figure 2: Comparison of fixed horizontal nodes versus material nodes.
Figure 3: Comparison of two-dimensional versus three-dimensional computations
Figure 4: Energy in a wave tank as a function of time
Figure 5: Free surface elevation at various distances from a piston-type wavemaker as a function of time.
Present calculations ——;
Dommermuth et al. (1988)
experiments---;
calculations •

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Heaving two-dimensional rectangle
Two-dimensional computations have been performed on a rectangle in heave and sway; details of the calculations may be found in Lee (1992). The convergence of the numerical results and agreement with experiments are similar for both heave and sway. Because of space limitations only a sampling of the heaving results will be presented here.
Three node spacings on the body have been tried. The first is simple equal spacing between the nodes. The second is a cosine spacing in which the nodes are more closely packed around the sharp corner and the free surface. The third is a semi-cosine spacing that retains the cosine spacing on the bottom but uses only half a cosine on the sides. This keeps the closely packed spacing around the corner but increases the panel size near the free surface. The semi-cosine spacing is useful because for stability reasons the panels on the free surface near the body need to be approximately the same size as the body panels near the waterline. Thus, larger free surface panels can be used for a given number of body panels. Consistent with the findings of Cao et al. (1992) for a Karman-Trefftz airfoil, no node is placed at the corner in any of the spacing schemes in order to avoid a spike in the tangential velocity. Calculations have shown little variation in the numerical results with node spacing schemes.
Figure 6 shows the typical wave profiles generated by a heaving rectangle in infinite water depth. The beam of the rectangle, B, is two and the draft, H, is one. The nondimensional heaving frequency is and the motion amplitude is 10% of the draft (a/H=.1). Semi-cosine spacing of is used on the body and approximately 30 material nodes per wavelength were used on the free surface. Rather than a beach in the far field, Lee (1992) used an algebraically increasing panel size starting at four wavelengths from the body. The increasing panels are carried out for 80 more wavelengths with the last panel being many wavelengths long. This far field panel arrangement introduces significant numerical damping and greatly delays the significantly attenuated reflected waves from entering the near field of interest. As can be seen in the figure, there is essentially no reflected waves at the edge of the inner region for the time interval shown.
Figure 7 shows the convergence of the numerical calculations to the experimental results of Vughts (1968). Vughts conducted experiments with an oscillating model and compared the results to linear theory calculations using both conformal mapping and Lewis forms. The fully nonlinear calculations converge to a value between Vugts' experimental and numerical results for a nondimensional frequency of π/3.
The comparisons between Vughts' results for both the added mass and damping coefficients over the entire frequency range are shown in figures 8 and 9. There figures are for a nondimensional amplitude of .15, the largest amplitude tested. The results for smaller amplitudes are similar. Analogous to the experiments, our results show very little sensitivity to motion amplitude. The variation of the added mass and damping coefficients with amplitude fall approximately within the circles used to plot the fully nonlinear results in figures 8 and 9. As can be seen from the figures, the fully nonlinear computations agree well with the experiments and with linear theory.
Heaving axially-symmetric bodies
The flow due to a heaving three-dimensional cylinder has also been studied, the details may be found in Lee (1992). Axially-symmetric problems ca n be reduced to a one-dimensional integral equation by employing ring sources such as used by Dommermuth and Yue (1987). Because the desingularized technique is being developed for general body shapes, ring sources were not used in our calculations. Instead, the calculations were performed using multiple planes of symmetry. The first calculations were done using two planes of symmetry so that the unknown sources were distributed over one quarter of the domain. More computational efficiency is gained by increasing the number of planes of symmetry. Figure 10 shows the waves produced at different periods for a circular cylinder (R/H=1.0) in heave. Two sets of results are shown corresponding to using two planes of symmetry and 60 planes of symmetry. The actual panelization of the free surface is very different for the two cases since in the two-plane case the number of nodes is increased as the distance from the origin increases and in the other it is constant. In both cases the number of nodes on the cylinder and at the waterline on the free surface are the same. The results are so close that only one line is apparent in the figure for each period shown.
A more accurate measure of the difference between the predicted wave elevation is given in table I. Table I presents the RMS difference between the wave elevations produced by the two calculations for the different periods. The RMS difference is a measure of the difference in predicted wave profiles at a total of 120 fixed wave probes. The probes were uniformly distributed over the computational domain of interest. Cubic spline interpolation is used to determine the wave elevation at each probe from the values of the wave elevation at the computational

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.5 Wave making resistance coefficient Cw of point doublet
Fig.6(a) Wave profiles due to circular cylinder
Fig.6(b) Wave profiles due to circular cylinder
Fig.7 Wave making resistance coefficient Cw of circular cylinder
Fig.8 Pressure distribution on circular cylinder surface
Fig.9 Comparison of 2-D panel method and SQCM

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.10 Wave profiles due to hydrofoil (NACA0005, h/c=0.956, α=5º)
Fig.11 Lift and drag with analytical solution
Fig.12 Pressure distribution on hydrofoil surface

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.13(a) Comparison of source distribution between HSA and LSA
Fig.13(b) Comparison of source distribution between HSA and LSA
Fig.14(a) Comparison of wave profiles between HSA and LSA
Fig.14(b) Comparison of wave profiles between HSA and LSA
Fig.15(a) Change of interference factor versus stagger
Fig.15(b) Change of interference factor versus stagger

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.16(a) Wave profiles versus stagger
Fig.16(b) Wave profiles versus stagger

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.17(a) Change of interference factor versus Fn
Fig.17(b) Change of interference factor versus Fn
Fig.18 Wave profiles versus Fn

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.19 Pressure coef. on 3-D mono-foil surface in unbounded flow
Fig.20 Pressure coef. on 3-D tandem-foil surfaces in unbounded flow
Fig.21 interference factor of lift and drag against stagger
Fig.22(a) Comparison of wave profiles and patterns between LSA and HSA (Coarse mesh)

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.22(b) Comparison of wave profiles and patterns between LSA and HSA (Coarse mesh)
Fig.23 Comparison of wave profiles and patterns between LSA and HSA (Fine mesh)
Fig.24 Lift coef. and drag coef. of mono-hydrofoil

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.25(a) Comparison of wave profiles and patterns between LSA and HSA (tandem foil, s/c=6.0)
Fig.25(b) Comparison of wave profiles and patterns between LSA and HSA (tandem foil, s/c=6.0)
Fig.25(c) Comparison of wave profiles and patterns between LSA and HSA (tandem foil, s/c=6.0)

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
Fig.26(a) Comparison of CL (fore foil)
Fig.26(b) Comparison of CL (aft foil)
Fig.27(a) Comparison of CD (fore foil)
Fig.27(b) Comparison of CD (aft foil)
Fig.28(a) Comparison of CL (tandem foil)
Fig.28(b) Comparison of CD (tandem foil)

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
DISCUSSION
by Dr. Spyros Kinnas, MIT
The authors should be congratulated for their extensive study of the Rankine Source Method. I have two questions concerning the so-called SQCM:
What is the theoretical basis of combining a QCM (an inherently linearized formulation) to a source method with the sources distributed on the exact hydrofoil surface? Is this a consistent formulation?
Can their method predict the correct effect of thickness on the hydrofoil loading? (See for example JSR paper on “A General Theory for the Coupling Between Thickness and Loading,” S.Kinnas, 1992). In other words, is their method predicting the correct value of the derivative of the lift L with respect to the τmax (dL/dτmax) while the camber and angle of attack remains the same?
Author's Reply
We can select an arbitrary form and location of distributing vortices around the foil. Then we use those of the QCM formally and satisfy the Kutta's condition. Though there is no assumption of a thin wing, we think this is consistent.
In Figure A1 we show the SQCM results of the lift variation due to thickness comparing with our panel method and the QCM for the 2-D NACA symmetric wing sections. We can see the thickness increases lift as well as the panel method and the discusser's method.
DISCUSSION
by Kazu-hiro Mori, Hiroshima University
Can we not expect to have the same results even if the expressions of the basic solution are not the same? This can be supported by the results of Cw shown in Fig. 7. The differences in wave elevation shown in Fig. 4 are suspected to have appeared from numerical reasons; the elevation is obtained by the finite difference of the velocity potential with respect to x.
Author's Reply
Thank you for your discussion.
We agree with the discusser in thinking that both results should agree with each other on an ideal numerical calculation. However, we have tried to use finer mesh on still water surface using 2-D circular cylinder problem and have also checked them using double precision (REAL$*8) and four times precision (REAL$*16) options, but it is confirmed that the depression becomes smaller but never agree with the HSA solution. We may say that both results never agree at high speed.
In order to understand this phenomenon we show calculated results in Figure A2. The condition is corresponding to Figure 4(b) of the text, that is, one of the 2-D circular cylinder problem with a point doublet. We inversely obtained the source distribution on the still water surface from the wave profile which is obtained by the HSA. This source distribution differs only near the steep peaks from results of the LSA. This difference causes the peculiar depression above the doublet.
We also show other two results by the ordinary LSA; one is the usual case in which calculated region on the still water surface(−15.0<x <15.0) is divided into about 150 pieces and the other is the 25 times finer mesh case. Even if we use very fine mesh, we can not follow the steep peaks obtained above.
Then, we think the HSA is very effective in order to avoid this numerical difficulty in the high speed range.
DISCUSSION
by Dr. H.Raven, MARIN
Your paper suggests that you have tested two different linearizatons: LSA, with the basic flow being a double body flow, and HSA, with the basic flow being an ‘inverse image flow.' You then conclude that HSA is better because it more closely agrees with analytic solutions. But these analytic solutions are linearized with respect to uniform flow, and will never tell you which linearization is better. To get the best

OCR for page 1

Proceedings of the Sixth International Conference on Numerical Ship Hydrodynamics
agreement with these analytic solutions you should impose a Kelvin condition.
On the other hand, (6) is not a useful double-body linearization but is equivalent to a Kelvin condition, could you clarify this?
Author's Reply
Thank you for your discussion.
The Eq. (6) is exactly the Kelvin condition which is derived substituting to the linearized free surface condition , but is not double-body linearization, i.e.,
in case of normal image.
However, the difference between the LSA and the HSA is not one of linearization but one of the equation expressed in the text. That is, in the case of the LSA, Eq. (6) becomes from and it becomes also from in the case of the HSA.
DISCUSSION
by Dr. Henry T.Wang, Naval Research Laboratory
The authors are to be commended for bringing out the High Speed Approximation (HSA) as an alternative to the Low Speed Approximation (LSA), which is usually used in ship hydrodynamics. The authors' figures indicate that the wave profile and force coefficient results for LSA and HSA are similar at low Fn and diverge at high Fn. Have the authors considered other criteria, which may indicate a priori the approximate point of divergence of the two theories? For example, could one such criterion be related to the vertical velocity at the free surface (which is identically zero for LSA) reaching a certain fraction of the hydrofoil velocity? Also, can the authors comment on the conditions in which HSA, rather than the presently used LSA, would be more appropriate as the basic flow for surface ship problems?
Author's Reply
Thank you for your discussion.
It seems to be very useful to use such an index. In fact we have tried to search for a criteria, but we have found nothing indicating those matters generally. As the LSA is effective at a somewhat high speed and the HSA is also useful in a low speed problem to a degree, it is relative to use the LSA or the HSA. The point is that the speed range should be respected in which both results agree, and when they begin to diverge numerically, we should choose the HSA or the LSA.
We wish to consider the surface ship problem in the near future.