iccsa-20-wind

main.tex (27566B)

---

 \documentclass[runningheads]{llncs}
 
 \usepackage{amsmath}
 \usepackage{booktabs}
 \usepackage{graphicx}
 \usepackage{url}
 \usepackage{textcomp}
 
 \newcommand{\VectorR}[1]{\left[\begin{array}{r}#1\end{array}\right]}
 \newcommand{\VectorL}[1]{\left[\begin{array}{l}#1\end{array}\right]}
 \newcommand{\Length}[1]{\big|#1\big|}
 
 \begin{document}
 
 \title{Virtual Testbed: Simulation of air flow around ship hull and its effect
     on ship motions\thanks{Supported by Saint Petersburg State University (grants
     no.~51129371 and~51129725) and Council for grants of the President of the Russian
     Federation (grant no.~\mbox{MK-383.2020.9}).}}
 \author{%
     Anton Gavrikov\orcidID{0000-0003-2128-8368} \and\\
     Alexander Degtyarev\orcidID{0000-0003-0967-2949} \and\\
     Denis Egorov\orcidID{0000-0001-5478-2798}\and\\
     Ivan Gankevich\textsuperscript{*}\orcidID{0000-0001-7067-6928} \and\\
     Artemii Grigorev\orcidID{0000-0003-0501-7656} \and\\
     Vasily Khramushin\orcidID{0000-0002-3357-169X} \and\\
     Ivan Petriakov\orcidID{0000-0001-5835-9313}
 }
 
 \titlerunning{Simulation of air flow around ship hull and its effect on ship motions}
 \authorrunning{A.\,Gavrikov et al.}
 
 \institute{Saint Petersburg State University\\
     7-9 Universitetskaya Emb., St Petersburg 199034, Russia\\
 \email{st047437@student.spbu.ru},\\
 \email{a.degtyarev@spbu.ru},\\
 \email{st047824@student.spbu.ru},\\
 \email{i.gankevich@spbu.ru},\\
 \email{st016177@student.spbu.ru},\\
 \email{v.khramushin@spbu.ru},\\
 \email{st049350@student.spbu.ru}\\
 \url{https://spbu.ru/}}
 
 \maketitle
 
 \begin{abstract}
 
 Strong wind causes heavy load on the ship in a seaway bending and pushing it in
 the direction of the wind. In this paper we investigate how wind can be
 simulated in the framework of Virtual testbed~--- a near real-time ship motion
 simulator. We propose simple model that describes air flow around ship hull
 with constant initial speed and direction which is based on the law of
 reflection. On the boundary the model reduces to the known model for potential
 flow around a cylinder, and near the boundary they are not equivalent, but
 close enough to visualise the effect of the hull on the flow. Then we apply
 this model to simulate air flow around real-world ship hull and conclude that
 for any real-world situation ship roll angle and ship speed caused by the wind
 is small to not cause capsizing, but large enough to be considered in onboard
 intelligent systems that determine real roll, pitch and yaw angles during ship
 operation and similar applications.
 
 \keywords{%
     wind field
 \and law of reflection
 \and flow around cylinder
 \and uniform translational motion
 \and OpenMP
 \and OpenCL
 \and GPGPU.
 }
 \end{abstract}
 
 \section{Introduction}
 
 Ship motion simulation studies focus on interaction between the ship and ocean
 waves~--- a physical phenomena that gives the largest contribution to
 oscillatory motion~--- however, intelligent onboard systems require taking
 other forces into account. One of the basic functionality of such a system is
 determination of initial static ship stability parameters (roll angle, pitch
 angle and draught) from the recordings of various ship motion parameters, such
 as instantaneous roll, pitch and yaw angles, and their first and second
 instantaneous derivatives (e.g.~angular velocity and angular acceleration).
 During ship operation these initial static ship stability parameters deviate
 from the original values as a result of moving cargo between compartments,
 damaging the hull, compartment flooding etc. These effects are especially
 severe for fishing and military vessels, but can occur with any vessel
 operating in extreme conditions.
 
 Intelligent onboard system needs large amount of synchronous recordings of ship
 motions parameters to operate, mainly angular displacement, velocity and
 acceleration and draught, but these parameters depend on the shape of the ship
 hull and obtaining them in model tests is complicated, let alone field tests.
 Field tests are too expensive to perform and do not allow to simulate
 particular phenomena such as compartment flooding.  Model tests are too
 time-consuming for such a task and there is no reliable way to obtain all the
 derivatives for a particular parameter: sensors measure one particular
 derivative and all other derivatives have to be calculated by numerical
 differentiation or integration, and integration has low accuracy for time
 series of measurements~\cite{kok2017sensors}. The simplest way to obtain those
 parameters is to simulate ship motion on the computer and save all the
 parameters in the file for future analysis.
 
 Arguably, the largest contribution to ship motion besides ocean waves is given
 by wind forces: air has lesser density than water, but air motion acts on the
 area of ship hull which may be greater than underwater area due to ship
 superstructure. Steady wind may produce non-nought roll angle, and thus have to
 be taken into account when determining initial static ship stability
 parameters.  In this paper we investigate how wind velocity field can be
 simulated on the boundary and near the boundary of the ship hull.  We derive a
 simple mathematical model for uniform translational motion of the air on the
 above-water boundary of the ship hull.  Then we generalise this model to
 calculate wind velocity near the boundary still taking into account the shape
 of the above-water part of the ship hull.  Finally, we measure the effect of
 wind velocity on the ship roll angle and carry out computational performance
 analysis of our programme.
 
 \section{Related work}
 
 Studies on the effect of the wind on ship motions mostly focus on capsizing
 probability~\cite{bulian2004roll,paroka2006wind}, whereas our work focuses on
 the direct effect of wind on ship motions and on how to incorporate ship roll
 angle change due to wind in onboard intelligent systems. As a result, we do not
 use probabilistic methods, but we use direct simulation of air flow.
 
 Similar simulations can be performed in a wind tunnel~\cite{andersen2013wind}
 but in the case of onboard intelligent systems we need to gather a lot of
 statistical data to be able to tune the system for each ship hull shape.
 Performing large number of simulations in a wind tunnel is time-consuming and a
 computer programme is the most efficient option for this task.
 
 \section{Methods}
 
 \subsection{Analytic representation of wind velocity field}
 
 Air motion without turbulence can be decomposed into two components:
 translational motion~--- air particles travel in the same direction with
 constant velocity, and circular motion~--- air particles travel on a circle.
 Translational motion describe sea breeze, that occurs on the shore on the
 sunrise and after the sunset. Circular motion describe storms such as typhoons
 and hurricanes. Translational motion is a particular case of circular motion
 when the radius of the circle is infinite. Due to the fact that the scale of
 circular motion is much larger than the size of a typical ship hull we consider
 only translational motion in this paper.
 
 Since there is no rotational component, air flow is described by equations for
 irrotational inviscid incompressible fluid. In this context fluid velocity
 \(\vec\upsilon\) is determined as a vector gradient \(\vec\nabla\) of scalar
 velocity potential \(\phi\) and continuity equation and equation of motion are
 written as
 \begin{equation}
     \label{eq-governing}
     \begin{aligned}
         & \Delta\phi = 0; \qquad \vec\upsilon=\vec\nabla\phi; \\
         & \rho\frac{\partial\phi}{\partial{}t} +
           \frac{1}{2}\rho\Length{\vec\nabla\phi}^2 +
           p + \rho g z = p_0.
     \end{aligned}
 \end{equation}
 Here \(p_0\) is atmospheric pressure, \(g\) is gravitational acceleration,
 \(\rho\) is air density, \(p\) is pressure. We seek solutions to this system
 of equations for velocity potential \(\phi\). Continuity equation restricts
 the type of the function that can be used as the solution, and equation of
 motion gives the pressure for a particular velocity potential value.
 
 Ship hull boundary is defined by a parametric surface \(\vec{S}\) and surface
 normal \(\vec{n}\):
 \begin{equation*}
 \vec{S}=\vec{S}\left(a,b,t\right)
 \qquad
 a,b\in{}A=[0,1];
 \qquad
 \vec{n}=\frac{\partial\vec S}{\partial a} \times \frac{\partial\vec S}{\partial b}
 \end{equation*}
 The simplest parametric surface is infinite plane with constant normal.
 The computer model of a real ship hull is composed of many triangular panels
 with different areas and different orientations that approximate continuous
 surface. On the boundary the projection of wind velocity on the surface normal
 is nought:
 \begin{equation}
     \label{eq-boundary}
     \vec\nabla\phi\cdot\vec{n} = 0;
     \qquad
     \vec{r} = \vec{S}.
 \end{equation}
 
 The solutions to the governing system of equations differ in how boundary is
 incorporated into them: in our model the boundary is taken into account by adding
 velocity of a reflected air particle in the solution. Velocity
 \(\vec\upsilon_r\) of the particle that is reflected from the surface with
 surface normal \(\vec{n}\) is given by the law of reflection
 (fig.~\ref{fig-law-of-reflection}):
 \begin{equation}
     \label{eq-reflected}
     \vec\upsilon_r = \vec\upsilon - 2\left(\vec\upsilon\cdot\vec{n}\right)\vec{n}.
 \end{equation}
 When we add velocity of incident and reflected air particles we get a vector
 that is parallel to the boundary. As we move away from the boundary its impact
 on the velocity decays quadratically with the distance. The reason for using
 the law of reflection to describe
 air flow on the boundary is that the corresponding solution reduces to the known
 analytic solution for the potential flow around a cylinder
 (see~sec.~\ref{sec-cylinder}). Quadratic decay term is borrowed from this known
 solution.
 
 \begin{figure}
     \centering
     \includegraphics{build/inkscape/law-of-reflection.eps}
     \caption{The law of reflection diagram for incident and reflected air particle
     velocities \(\vec\upsilon\) and \(\vec\upsilon_r\) and surface normal
     \(\vec{n}\).\label{fig-law-of-reflection}}
 \end{figure}
 
 In the following subsections we describe the solution that we obtained for the
 velocity field \emph{on} the boundary and \emph{near} the boundary.
 
 \subsection{Uniform translational motion on the static body surface}
 
 On the surface we neglect the impact of neighbouring panels on the velocity
 field on the ground that the real ship hull surface is smooth,
 i.e.~neighbouring panels have approximately the same normals. This assumption
 does not hold for aft and bow of some ships, and, as a result, velocity field
 near these features has stream lines with sharp edges. We consider this effect
 negligible for the determination of roll angle caused by the wind, since the
 area of panels that distort wind field is small compared to the area of all
 other panels. 
 
 We seek solutions to the governing system of equations~\eqref{eq-governing} with
 boundary condition~\eqref{eq-boundary} of the form
 \begin{equation*}
 \phi = \vec\upsilon\cdot\vec{r}
 + C \left(\vec\upsilon_r\cdot\vec{r}\right);
 \qquad
 \vec{r}=\left(x,y,z\right),
 \end{equation*}
 Here \(\vec{r}\) is spatial coordinate, \(C\) is the coefficient, and
 \(\vec\upsilon_r\) is velocity of reflected air particle defined
 in~\eqref{eq-reflected}. This solution is independent for each panel.
 Plugging the solution into boundary condition~\eqref{eq-boundary} gives
 \begin{equation*}
 \left(\vec\upsilon + C\vec\upsilon_r\right)\cdot\vec{n} = 0,
 \end{equation*}
 hence
 \begin{equation*}
 C = -\frac{ \vec\upsilon\cdot\vec{n} }{ \vec\upsilon_r\cdot\vec{n} } = 1
 \end{equation*}
 and velocity is written simply as
 \begin{equation}
     \label{eq-solution-on-the-boundary}
     \vec\nabla\phi = \vec\upsilon + \vec\upsilon_r.
 \end{equation}
 
 This solution satisfies continuity equation. It gives velocity only at the
 centre of each ship hull panel, but this is sufficient to calculate pressure
 and force moments acting on the ship hull.
 
 \subsection{Uniform translational motion near the static body surface}
 
 Near the surface there are no neighbouring panels, the impact of which we can
 neglect, instead we add reflected particle velocities for all the panels and
 decay the velocity quadratically with the distance to the panel. Here we can
 neglect panels surface normals of which has large angles with the wind
 direction for efficiency, but they do not blow up the solution.
 
 We seek solutions of the form
 \begin{equation*}
 \phi = \vec\upsilon\cdot\vec{r}
 + \iint\limits_{a,b\,\in{}A}
 C \frac{\vec\upsilon_r\cdot\vec{r}}{1+\Length{\vec{r}-\vec{S}}^2}
 da\,db,
 \end{equation*}
 where \(\Length{\cdot}\) denotes vector length.  Plugging the solution into boundary
 condition and assuming that neighbouring panels do not affect each other (this allows
 removing the integral) gives the same coefficient \(C=1\), but velocity vector is written
 differently as
 \begin{equation}
 \label{eq-solution-near-the-boundary}
 \vec\nabla\phi =
 \vec\upsilon +
 \iint\limits_{a,b\,\in{}A}
 \left(
     \frac{1}{s} \vec\upsilon_r
     - \frac{2}{s^2} \left(\vec\upsilon_r\cdot\vec{r}\right) \left(\vec{r}-\vec{S}\right)
 \right)
 da\,db;
 \qquad
 s = 1+\Length{\vec{r}-\vec{S}}^2.
 \end{equation}
 Besides the term for reflected air particle velocity that decays quadratically with
 the distance to the panel, there is a term that decays quaternary with the distance and
 that can be neglected.
 
 This solution reduces to the solution on the boundary when \(\vec{r}=\vec{S}\)
 and takes into account impact of each panel on the velocity direction which decays
 quadratically with the distance to the panel.
 
 \section{Results}
 
 \subsection{Verification of the solution on the example of potential flow around a cylinder}
 \label{sec-cylinder}
 
 Potential flow around a cylinder in two dimensions is described by the
 following well-known formula:
 \begin{equation*}
     \phi\left(r,\theta\right) = U r \left( 1 + \frac{R^2}{r^2} \right) \cos\theta.
 \end{equation*}
 Here \(r\) and \(\theta\) are polar coordinates, \(R\) is cylinder radius and \(U\)
 is \(x\) component of the velocity. Cylinder is placed at the origin. To prove that our
 solution on the boundary~\eqref{eq-solution-on-the-boundary} reduces to this solution
 we write it in Cartesian form using polar coordinate identities
 \begin{equation*}
     r = \sqrt{x^2+y^2}; \qquad \theta = \arccos{\frac{x}{\sqrt{x^2+y^2}}}.
 \end{equation*}
 Then in Cartesian coordinates the solution is written as
 \begin{equation*}
     \phi(x,y) = U x \left(1 + \frac{R^2}{x^2+y^2}\right)
 \end{equation*}
 and the velocity is written as
 \begin{equation}
     \label{eq-solution-cylinder}
     \vec\nabla\phi = \VectorL{U \left(R^2 (y^2-x^2) + (x^2+y^2)^2 \right) \\
     -2 R^2 U x y} / (x^2+y^2)^2.
 \end{equation}
 On the boundary \(x^2+y^2=R^2\) and the velocity is written as
 \begin{equation*}
     \vec\nabla\phi = \frac{2 U}{R^2}\VectorR{y^2 \\ -x y\phantom{^2}}.
 \end{equation*}
 Now if we write surface normal as \(\vec{n}=(x/R,y/R)\) and let
 \(\vec\upsilon=(U,0)\), our solution \eqref{eq-solution-on-the-boundary}
 quite surprisingly reduces to the same expression.
 
 To reduce solution near the boundary~\eqref{eq-solution-near-the-boundary} to
 the solution for potential flow around a cylinder, we let
 \(s=\Length{\vec{r}}^2/\Length{\vec{S}}^2\) (here \(\vec{r}\) is the radius
 vector in Cartesian coordinates). Then the solution is written as
 \begin{equation*}
     \vec\nabla\phi = \vec\upsilon + \frac{1}{s} \vec\upsilon_r
 \end{equation*}
 and reduces to general form of the solution for potential flow around a
 cylinder given in~\eqref{eq-solution-cylinder}.
 
 We compared both solutions (see fig.~\ref{fig-verification}) for a cylinder
 with \(R=27.4429\), and discovered that maximum distance between velocity
 fields produced from both of them is 11\% of the maximum velocity. Our formula
 shows slightly smaller decay near the boundary, but this problem cannot be
 solved by introducing coefficients. Perhaps, comparing our solution to the
 solution produced by CFD methods may shed light on which one is closer to the
 reality. On the ship hull boundary the solutions are equivalent.
 
 \begin{figure}
     \centering
     \includegraphics{build/gnuplot/verification.eps}
     \caption{Wind velocity fields for air flow around a cylinder: left~--- known
     solution~\eqref{eq-solution-cylinder}, right~--- our
     solution~\eqref{eq-solution-near-the-boundary}.\label{fig-verification}}
 \end{figure}
 
 \subsection{Ship roll angle and velocity}
 
 Wind causes non-nought force moment and force acting on a ship that bend and
 move the ship in the direction of the wind. This effect is stronger for the
 ships with large hull areas exposed to the wind, like fully-loaded
 containerships, but other ships are also affected.
 
 We measured how wind speed affects Aurora's transversal velocity and roll
 angle.  For that purpose we made wind blow directly in the starboard of the
 ship and varied wind speed. We stopped the experiment after 60 seconds and
 measured maximum roll angle and maximum transversal velocity.  We have found
 that in order to produce 1\textdegree{} static roll angle we need wind speed of
 \(\approx{}35\) m/s (fig.~\ref{fig-roll-velocity}), and wind with that speed
 makes the ship move in transversal direction with the speed of \(\approx{}0.2\)
 m/s. We expect these numbers to be smaller for smaller ships.
 
 We have found that the law of reflection~\eqref{eq-reflected} in its original
 form does not allow to calculate the effect of wind on the symmetric ship hull:
 it happens because the pressure on the leeward and windward side of the ship is
 the same~--- which is not the case for real-world phenomena where the pressure
 is different due to turbulence. To overcome this problem we introduce a
 coefficient \(\alpha\) that controls reflection ratio:
 \begin{equation*}
     \vec\upsilon_r = \vec\upsilon - 2 \alpha \left(\vec\upsilon\cdot\vec{n}\right)\vec{n}.
 \end{equation*}
 When \(\alpha=1\) this formula equals~\eqref{eq-reflected}, when \(\alpha=0\) there
 is no reflection and the wind velocity does not change its direction near the
 ship hull.  In our tests we used \(\alpha=0.5\). A better solutions would be to
 incorporate turbulence in the model which is one of the directions of future research.
 
 \begin{figure}
     \centering
     \includegraphics{build/gnuplot/velocity.eps}
     \caption{Dependence of Aurora's roll angle on wind speed (left) and
     dependence of Aurora's transversal velocity on wind speed
     (right).\label{fig-roll-velocity}}
 \end{figure}
 
 \subsection{Computational performance analysis}
 
 We implemented solutions~\eqref{eq-solution-on-the-boundary}
 and~\eqref{eq-solution-near-the-boundary} in Virtual testbed wind solver.
 Virtual testbed is a programme for workstations that simulates ship motions in
 extreme conditions and physical phenomena that causes them: ocean waves, wind,
 compartment flooding etc.
 
 We performed benchmarks for three ships: Diogen, Aurora and MICW.  Diogen is a
 small-size fishing vessel, Aurora is mid-size cruiser and MICW is a large-size
 ship with small moment of inertia for the current waterline
 (fig.~\ref{fig-ships}). Ships in our database do not have superstructures, they
 have only hulls and compartments. The main difference between them that affects
 benchmarks is the number of panels into which the hull is decomposed. These numbers
 and sizes are shown in tab.~\ref{tab-ships}.
 
 \begin{table}
 	\centering
 	\caption{Parameters of ship hulls that were used in the
 	benchmarks.\label{tab-ships}}
 	\begin{tabular}{lrrr}
 		\toprule
 					  & Diogen & Aurora & MICW  \\
 		\midrule               
 		Length, m     &   60   & 126.5  &   260 \\
 		Beam, m       &   15   &  16.8  &    32 \\
 		Depth, m      &   15   &  14.5  &    31 \\
 		No. of panels & 4346   &  6335  &  9252 \\
 		\bottomrule
 	\end{tabular}
 \end{table}
 
 Benchmarks were performed using three workstations: DarkwingDuck, GPUlab,
 Capybara. DarkwingDuck is a laptop, GPUlab is a desktop workstation, and
 Capybara is a desktop with professional graphical accelerator server-grade
 processor (tab.~\ref{tab-config}). All of the workstations are equipped with
 graphical accelerators that allow to greatly increase their performance.
 
 Wind solver was written for both OpenMP and OpenCL to make use of graphical
 accelerator available on most modern workstations. The solvers use single
 precision floating point numbers. Benchmark results are presented in
 tab.~\ref{tab-benchmark}.
 
 \begin{table}
 	\centering
 	\caption{Hardware configurations for benchmarks. For all benchmarks we
 	used GCC version 9.1.0 compiler and optimisation flags \texttt{-O3
 	-march=native}.\label{tab-config}}
 	\begin{tabular}{lllrr}
 		\toprule
 		     &     &     & \multicolumn{2}{c}{GPU GFLOPS} \\
 		\cmidrule(r){4-5}
 		Node           & CPU              & GPU            & Single & Double \\
 		\midrule
         DarkwingDuck   & Intel i7-3630QM  & NVIDIA GT740M  &  622   &     \\
 		GPUlab         & AMD FX-8370      & NVIDIA GTX1060 & 4375   & 137 \\
 		Capybara       & Intel E5-2630 v4 & NVIDIA P5000   & 8873   & 277 \\
 		\bottomrule
 	\end{tabular}
 \end{table}
 
 \begin{table}
     \centering
     \caption{Performance benchmarks results. Numbers represent average time in milliseconds
     that is needed to compute wind field on the ship hull and near the ship
     hull.\label{tab-benchmark}}
     \begin{tabular}{l@{\hskip 3mm}r@{\hskip 2mm}r@{\hskip 3mm}r@{\hskip 2mm}r@{\hskip 3mm}r@{\hskip 2mm}r}
 		\toprule
 				 & \multicolumn{2}{l}{Diogen}
 				 & \multicolumn{2}{l}{Aurora}
 				 & \multicolumn{2}{l}{MICW} \\
 		\cmidrule(r){2-3}
 		\cmidrule(r){4-5}
 		\cmidrule(r){6-7}
 		Node     & MP & CL 
                  & MP & CL 
                  & MP & CL  \\
 		\midrule
 		DarkwingDuck  
                  &      114  &    14.00
                  &      164  &    16
                  &      314  &    29
                  \\
 		GPUlab  
 		         &       62  &     1.35            
                  &       90  &     2
 		         &      175  &     3
                  \\
 		Capybara &       33  &     0.87            
 		         &       48  &     1
 		         &       99  &     6
                  \\
 		\bottomrule
 	\end{tabular}
 \end{table}
 
 \begin{figure}
 	\centering
 	\includegraphics{build/ships/diogen.eps}\hfill
 	\includegraphics{build/ships/aurora.eps}\hfill
 	\includegraphics{build/ships/micw.eps}
 	\caption{Diogen, Aurora and MICW three-dimensional ship hull
 	models.\label{fig-ships}}
 \end{figure}
 
 \section{Discussion}
 
 Solution on the boundary \eqref{eq-solution-on-the-boundary} provides simple
 explanation of areas with the highest and lowest pressure for potential flow
 around a cylinder.  At left-most and right-most points on the cylinder boundary
 (fig.~\ref{fig-verification}) velocity is nought because incident and reflected
 particle velocities have opposite directions and cancel each other out. At
 top-most and bottom-most points incident and reflected particle velocities have
 the same direction and total velocity is two times larger than the velocity of
 the flow.
 
 In order to be compatible with the surface of any object, solution near the
 boundary \eqref{eq-solution-near-the-boundary} uses different term \(s\) than
 the solution for potential flow around a cylinder which makes reflected
 velocity term reach maximum value for the point on the boundary, and for the
 point near the boundary the solution includes reflected velocity vectors for
 each panel. 
 
 Wind speed of 35 m/s that we obtained in ship roll angle experiments matches
 hurricane with 14 m high waves on the Beaufort scale which would affect the
 ship more severely than the wind. The main reason for such a large value is
 that ships in our database do not have superstructures and hence the surface
 area affected by the wind is much smaller than in reality. Nevertheless, in
 intelligent onboard systems even small variations in static roll at less severe
 weather conditions have to be considered for the correct operation of the
 system.
 
 The introduction of the coefficient \(\alpha\) that controls reflection ratio
 is the simplest way of taking turbulence into account. \(\alpha<1\) increases
 the wind speed on the leeward side of the ship as a result of wind ``going
 around'' the ship hull. More sophisticated turbulence model would give more
 accurate results.
 
 Performance benchmarks showed that performance of both OpenMP and OpenCL
 solvers increases from the least powerful (DarkwingDuck) workstation to the
 most powerful one (Capybara). In addition to this, performance of OpenCL is
 always better than of OpenMP by a factor of 16--58. Our solver uses explicit
 analytic formula to compute wind field and to compute wind field at each point
 iterates over all panels of the ship. All panels are stored and accessed
 sequentially and all points of wind field are stored sequentially and accessed
 in parallel, which makes the solver easy to implement in OpenCL and allows to
 achieve high performance on the graphical accelerator. Finally, as expected
 performance increases and the ratio between OpenCL and OpenMP performance
 decreases from the large-size to small-size ship.  The only exception from the
 above-mentioned observations is the performance of OpenCL on Capybara for MICW
 hull. This behaviour requires further investigation. 
 
 
 \section{Conclusion}
 
 This paper proposes a new simple mathematical model for wind field around the
 ship hull.  On the ship hull boundary this model is equivalent to the known
 formula for potential flow around a cylinder. Near the boundary this model is
 close to this formula, but has slightly smaller decay. In both cases the model
 satisfies boundary conditions and continuity equation (conservation of mass),
 which makes it suitable for physical simulations. The main advantage of the
 model is its simplicity, the use of Cartesian coordinates and its applicability
 to bodies of any form, not just cylinders.
 
 We applied this model to simulate ship motions under the effect of wind with
 constant speed and direction (and in the absence of all other effects except
 buoyancy force), and discovered that to get static roll angle of 1\textdegree{}
 we need wind speed of a hurricane (12 on the Beaufort scale). Also, simulation
 of ship motions due to wind is not possible without taking into account
 turbulence, but in the absence of turbulence model we used the simple
 coefficient that controls reflection ratio to adapt the model for this kind of
 simulation.
 
 From the computational standpoint the proposed model shows high performance on
 modern processors as well as graphical accelerators due to linear memory access
 pattern and absence of synchronisation and data transfer between parallel
 processes.  Using any up-to-date workstation is enough to perform real-world
 simulations.
 
 Future work is to include stratification and circular motion in the model.
 Stratification~--- an increase of wind speed with height~--- is known phenomena
 in atmosphere which affects wind field around tall ships, and thus may improve
 accuracy of our solver.  The motivation behind including circular motion is to
 better understand air motion around object and the fact that linear motion is a
 special case of it when the circle radius is infinite. Another possible
 direction of future work is to use more advanced turbulence model.
 
 \subsubsection*{Acknowledgements.}
 Research work is supported by Saint Petersburg State University (grants
 no.~51129371 and~51129725) and Council for grants of the President of the Russian Federation
 (grant no.~MK-383.2020.9).
 
 \bibliographystyle{splncs04}
 \bibliography{references}
 
 \end{document}