\documentclass[11pt,twoside]{article}\makeatletter

\IfFileExists{xcolor.sty}%
  {\RequirePackage{xcolor}}%
  {\RequirePackage{color}}
\usepackage{colortbl}
\usepackage{wrapfig}
\usepackage{ifxetex}
\ifxetex
  \usepackage{fontspec}
  \usepackage{xunicode}
  \catcode`⃥=\active \def⃥{\textbackslash}
  \catcode`❴=\active \def❴{\{}
  \catcode`❵=\active \def❵{\}}
  \def\textJapanese{\fontspec{Noto Sans CJK JP}}
  \def\textChinese{\fontspec{Noto Sans CJK SC}}
  \def\textKorean{\fontspec{Noto Sans CJK KR}}
  \setmonofont{DejaVu Sans Mono}
  
\else
  \IfFileExists{utf8x.def}%
   {\usepackage[utf8x]{inputenc}
      \PrerenderUnicode{–}
    }%
   {\usepackage[utf8]{inputenc}}
  \usepackage[english]{babel}
  \usepackage[T1]{fontenc}
  \usepackage{float}
  \usepackage[]{ucs}
  \uc@dclc{8421}{default}{\textbackslash }
  \uc@dclc{10100}{default}{\{}
  \uc@dclc{10101}{default}{\}}
  \uc@dclc{8491}{default}{\AA{}}
  \uc@dclc{8239}{default}{\,}
  \uc@dclc{20154}{default}{ }
  \uc@dclc{10148}{default}{>}
  \def\textschwa{\rotatebox{-90}{e}}
  \def\textJapanese{}
  \def\textChinese{}
  \IfFileExists{tipa.sty}{\usepackage{tipa}}{}
\fi
\def\exampleFont{\ttfamily\small}
\DeclareTextSymbol{\textpi}{OML}{25}
\usepackage{relsize}
\RequirePackage{array}
\def\@testpach{\@chclass
 \ifnum \@lastchclass=6 \@ne \@chnum \@ne \else
  \ifnum \@lastchclass=7 5 \else
   \ifnum \@lastchclass=8 \tw@ \else
    \ifnum \@lastchclass=9 \thr@@
   \else \z@
   \ifnum \@lastchclass = 10 \else
   \edef\@nextchar{\expandafter\string\@nextchar}%
   \@chnum
   \if \@nextchar c\z@ \else
    \if \@nextchar l\@ne \else
     \if \@nextchar r\tw@ \else
   \z@ \@chclass
   \if\@nextchar |\@ne \else
    \if \@nextchar !6 \else
     \if \@nextchar @7 \else
      \if \@nextchar (8 \else
       \if \@nextchar )9 \else
  10
  \@chnum
  \if \@nextchar m\thr@@\else
   \if \@nextchar p4 \else
    \if \@nextchar b5 \else
   \z@ \@chclass \z@ \@preamerr \z@ \fi \fi \fi \fi
   \fi \fi  \fi  \fi  \fi  \fi  \fi \fi \fi \fi \fi \fi}
\gdef\arraybackslash{\let\\=\@arraycr}
\def\@textsubscript#1{{\m@th\ensuremath{_{\mbox{\fontsize\sf@size\z@#1}}}}}
\def\Panel#1#2#3#4{\multicolumn{#3}{){\columncolor{#2}}#4}{#1}}
\def\abbr{}
\def\corr{}
\def\expan{}
\def\gap{}
\def\orig{}
\def\reg{}
\def\ref{}
\def\sic{}
\def\persName{}\def\name{}
\def\placeName{}
\def\orgName{}
\def\textcal#1{{\fontspec{Lucida Calligraphy}#1}}
\def\textgothic#1{{\fontspec{Lucida Blackletter}#1}}
\def\textlarge#1{{\large #1}}
\def\textoverbar#1{\ensuremath{\overline{#1}}}
\def\textquoted#1{‘#1’}
\def\textsmall#1{{\small #1}}
\def\textsubscript#1{\@textsubscript{\selectfont#1}}
\def\textxi{\ensuremath{\xi}}
\def\titlem{\itshape}
\newenvironment{biblfree}{}{\ifvmode\par\fi }
\newenvironment{bibl}{}{}
\newenvironment{byline}{\vskip6pt\itshape\fontsize{16pt}{18pt}\selectfont}{\par }
\newenvironment{citbibl}{}{\ifvmode\par\fi }
\newenvironment{docAuthor}{\ifvmode\vskip4pt\fontsize{16pt}{18pt}\selectfont\fi\itshape}{\ifvmode\par\fi }
\newenvironment{docDate}{}{\ifvmode\par\fi }
\newenvironment{docImprint}{\vskip 6pt}{\ifvmode\par\fi }
\newenvironment{docTitle}{\vskip6pt\bfseries\fontsize{22pt}{25pt}\selectfont}{\par }
\newenvironment{msHead}{\vskip 6pt}{\par}
\newenvironment{msItem}{\vskip 6pt}{\par}
\newenvironment{rubric}{}{}
\newenvironment{titlePart}{}{\par }

\newcolumntype{L}[1]{){\raggedright\arraybackslash}p{#1}}
\newcolumntype{C}[1]{){\centering\arraybackslash}p{#1}}
\newcolumntype{R}[1]{){\raggedleft\arraybackslash}p{#1}}
\newcolumntype{P}[1]{){\arraybackslash}p{#1}}
\newcolumntype{B}[1]{){\arraybackslash}b{#1}}
\newcolumntype{M}[1]{){\arraybackslash}m{#1}}
\definecolor{label}{gray}{0.75}
\def\unusedattribute#1{\sout{\textcolor{label}{#1}}}
\DeclareRobustCommand*{\xref}{\hyper@normalise\xref@}
\def\xref@#1#2{\hyper@linkurl{#2}{#1}}
\begingroup
\catcode`\_=\active
\gdef_#1{\ensuremath{\sb{\mathrm{#1}}}}
\endgroup
\mathcode`\_=\string"8000
\catcode`\_=12\relax

\usepackage[a4paper,twoside,lmargin=1in,rmargin=1in,tmargin=1in,bmargin=1in,marginparwidth=0.75in]{geometry}
\usepackage{framed}

\definecolor{shadecolor}{gray}{0.95}
\usepackage{longtable}
\usepackage[normalem]{ulem}
\usepackage{fancyvrb}
\usepackage{fancyhdr}
\usepackage{graphicx}
\usepackage{marginnote}

\renewcommand{\@cite}[1]{#1}


\renewcommand*{\marginfont}{\itshape\footnotesize}

\def\Gin@extensions{.pdf,.png,.jpg,.mps,.tif}

  \pagestyle{fancy}

\usepackage[pdftitle={A Numerical Solution for the Coexisting Field of Surface and Internal Solitary Waves},
 pdfauthor={}]{hyperref}
\hyperbaseurl{}

	 \paperwidth210mm
	 \paperheight297mm
              
\def\@pnumwidth{1.55em}
\def\@tocrmarg {2.55em}
\def\@dotsep{4.5}
\setcounter{tocdepth}{3}
\clubpenalty=8000
\emergencystretch 3em
\hbadness=4000
\hyphenpenalty=400
\pretolerance=750
\tolerance=2000
\vbadness=4000
\widowpenalty=10000

\renewcommand\section{\@startsection {section}{1}{\z@}%
     {-1.75ex \@plus -0.5ex \@minus -.2ex}%
     {0.5ex \@plus .2ex}%
     {\reset@font\Large\bfseries}}
\renewcommand\subsection{\@startsection{subsection}{2}{\z@}%
     {-1.75ex\@plus -0.5ex \@minus- .2ex}%
     {0.5ex \@plus .2ex}%
     {\reset@font\Large}}
\renewcommand\subsubsection{\@startsection{subsubsection}{3}{\z@}%
     {-1.5ex\@plus -0.35ex \@minus -.2ex}%
     {0.5ex \@plus .2ex}%
     {\reset@font\large}}
\renewcommand\paragraph{\@startsection{paragraph}{4}{\z@}%
     {-1ex \@plus-0.35ex \@minus -0.2ex}%
     {0.5ex \@plus .2ex}%
     {\reset@font\normalsize}}
\renewcommand\subparagraph{\@startsection{subparagraph}{5}{\parindent}%
     {1.5ex \@plus1ex \@minus .2ex}%
     {-1em}%
     {\reset@font\normalsize\bfseries}}


\def\l@section#1#2{\addpenalty{\@secpenalty} \addvspace{1.0em plus 1pt}
 \@tempdima 1.5em \begingroup
 \parindent \z@ \rightskip \@pnumwidth 
 \parfillskip -\@pnumwidth 
 \bfseries \leavevmode #1\hfil \hbox to\@pnumwidth{\hss #2}\par
 \endgroup}
\def\l@subsection{\@dottedtocline{2}{1.5em}{2.3em}}
\def\l@subsubsection{\@dottedtocline{3}{3.8em}{3.2em}}
\def\l@paragraph{\@dottedtocline{4}{7.0em}{4.1em}}
\def\l@subparagraph{\@dottedtocline{5}{10em}{5em}}
\@ifundefined{c@section}{\newcounter{section}}{}
\@ifundefined{c@chapter}{\newcounter{chapter}}{}
\newif\if@mainmatter 
\@mainmattertrue
\def\chaptername{Chapter}
\def\frontmatter{%
  \pagenumbering{roman}
  \def\thechapter{\@roman\c@chapter}
  \def\theHchapter{\roman{chapter}}
  \def\thesection{\@roman\c@section}
  \def\theHsection{\roman{section}}
  \def\@chapapp{}%
}
\def\mainmatter{%
  \cleardoublepage
  \def\thechapter{\@arabic\c@chapter}
  \setcounter{chapter}{0}
  \setcounter{section}{0}
  \pagenumbering{arabic}
  \setcounter{secnumdepth}{6}
  \def\@chapapp{\chaptername}%
  \def\theHchapter{\arabic{chapter}}
  \def\thesection{\@arabic\c@section}
  \def\theHsection{\arabic{section}}
}
\def\backmatter{%
  \cleardoublepage
  \setcounter{chapter}{0}
  \setcounter{section}{0}
  \setcounter{secnumdepth}{2}
  \def\@chapapp{\appendixname}%
  \def\thechapter{\@Alph\c@chapter}
  \def\theHchapter{\Alph{chapter}}
  \appendix
}
\newenvironment{bibitemlist}[1]{%
   \list{\@biblabel{\@arabic\c@enumiv}}%
       {\settowidth\labelwidth{\@biblabel{#1}}%
        \leftmargin\labelwidth
        \advance\leftmargin\labelsep
        \@openbib@code
        \usecounter{enumiv}%
        \let\p@enumiv\@empty
        \renewcommand\theenumiv{\@arabic\c@enumiv}%
	}%
  \sloppy
  \clubpenalty4000
  \@clubpenalty \clubpenalty
  \widowpenalty4000%
  \sfcode`\.\@m}%
  {\def\@noitemerr
    {\@latex@warning{Empty `bibitemlist' environment}}%
    \endlist}

\def\tableofcontents{\section*{\contentsname}\@starttoc{toc}}
\parskip0pt
\parindent1em
\def\Panel#1#2#3#4{\multicolumn{#3}{){\columncolor{#2}}#4}{#1}}
\newenvironment{reflist}{%
  \begin{raggedright}\begin{list}{}
  {%
   \setlength{\topsep}{0pt}%
   \setlength{\rightmargin}{0.25in}%
   \setlength{\itemsep}{0pt}%
   \setlength{\itemindent}{0pt}%
   \setlength{\parskip}{0pt}%
   \setlength{\parsep}{2pt}%
   \def\makelabel##1{\itshape ##1}}%
  }
  {\end{list}\end{raggedright}}
\newenvironment{sansreflist}{%
  \begin{raggedright}\begin{list}{}
  {%
   \setlength{\topsep}{0pt}%
   \setlength{\rightmargin}{0.25in}%
   \setlength{\itemindent}{0pt}%
   \setlength{\parskip}{0pt}%
   \setlength{\itemsep}{0pt}%
   \setlength{\parsep}{2pt}%
   \def\makelabel##1{\upshape ##1}}%
  }
  {\end{list}\end{raggedright}}
\newenvironment{specHead}[2]%
 {\vspace{20pt}\hrule\vspace{10pt}%
  \phantomsection\label{#1}\markright{#2}%

  \pdfbookmark[2]{#2}{#1}%
  \hspace{-0.75in}{\bfseries\fontsize{16pt}{18pt}\selectfont#2}%
  }{}
      \def\TheFullDate{2020-01-15 (revised: 15 January 2020)}
\def\TheID{\makeatother }
\def\TheDate{2020-01-15}
\title{A Numerical Solution for the Coexisting Field of Surface and Internal Solitary Waves}
\author{}\makeatletter 
\makeatletter
\newcommand*{\cleartoleftpage}{%
  \clearpage
    \if@twoside
    \ifodd\c@page
      \hbox{}\newpage
      \if@twocolumn
        \hbox{}\newpage
      \fi
    \fi
  \fi
}
\makeatother
\makeatletter
\thispagestyle{empty}
\markright{\@title}\markboth{\@title}{\@author}
\renewcommand\small{\@setfontsize\small{9pt}{11pt}\abovedisplayskip 8.5\p@ plus3\p@ minus4\p@
\belowdisplayskip \abovedisplayskip
\abovedisplayshortskip \z@ plus2\p@
\belowdisplayshortskip 4\p@ plus2\p@ minus2\p@
\def\@listi{\leftmargin\leftmargini
               \topsep 2\p@ plus1\p@ minus1\p@
               \parsep 2\p@ plus\p@ minus\p@
               \itemsep 1pt}
}
\makeatother
\fvset{frame=single,numberblanklines=false,xleftmargin=5mm,xrightmargin=5mm}
\fancyhf{} 
\setlength{\headheight}{14pt}
\fancyhead[LE]{\bfseries\leftmark} 
\fancyhead[RO]{\bfseries\rightmark} 
\fancyfoot[RO]{}
\fancyfoot[CO]{\thepage}
\fancyfoot[LO]{\TheID}
\fancyfoot[LE]{}
\fancyfoot[CE]{\thepage}
\fancyfoot[RE]{\TheID}
\hypersetup{citebordercolor=0.75 0.75 0.75,linkbordercolor=0.75 0.75 0.75,urlbordercolor=0.75 0.75 0.75,bookmarksnumbered=true}
\fancypagestyle{plain}{\fancyhead{}\renewcommand{\headrulewidth}{0pt}}

\date{}
\usepackage{authblk}

\providecommand{\keywords}[1]
{
\footnotesize
  \textbf{\textit{Index terms---}} #1
}

\usepackage{graphicx,xcolor}
\definecolor{GJBlue}{HTML}{273B81}
\definecolor{GJLightBlue}{HTML}{0A9DD9}
\definecolor{GJMediumGrey}{HTML}{6D6E70}
\definecolor{GJLightGrey}{HTML}{929497} 

\renewenvironment{abstract}{%
   \setlength{\parindent}{0pt}\raggedright
   \textcolor{GJMediumGrey}{\rule{\textwidth}{2pt}}
   \vskip16pt
   \textcolor{GJBlue}{\large\bfseries\abstractname\space}
}{%   
   \vskip8pt
   \textcolor{GJMediumGrey}{\rule{\textwidth}{2pt}}
   \vskip16pt
}

\usepackage[absolute,overlay]{textpos}

\makeatother 
      \usepackage{lineno}
      \linenumbers
      
\begin{document}

             \author[1]{Taro  Kakinuma}

             \author[2]{Kei  Yamashita}

             \affil[1]{  Kagoshima University}

\renewcommand\Authands{ and }

\date{\small \em Received: 6 December 2019 Accepted: 5 January 2020 Published: 15 January 2020}

\maketitle


\begin{abstract}
        


The numerical solutions for the coexisting fields of surface and internal solitary waves have been obtained, where the set of nonlinear equations based on the variational principle for steady waves are solved using the Newton- Raphson method. The relative phase velocity of surface-mode solitary waves is smaller in the coexisting fields of surface and internal solitary waves than in the cases without the coexistence of internal waves. The relative phase velocity of internal-mode solitary waves is also smaller in the coexisting fields of surface and internal solitary waves than in the cases without surface waves. The interfacial position of an internalmode internal solitary wave in a coexisting field of surface and internal waves can exceed the critical level determined in the corresponding case without a surface wave. The wave height ratio between internal-mode surface and internal solitary waves is smaller than the corresponding linear shallow water wave solution, and the difference increases, as the relative wave height of internal-mode internal solitary waves is increased.

\end{abstract}


\keywords{solitary wave, internal wave, free surface, nonlinear wave equation, numerical solution}

\begin{textblock*}{18cm}(1cm,1cm) % {block width} (coords) 
\textcolor{GJBlue}{\LARGE Global Journals \LaTeX\ JournalKaleidoscope\texttrademark}
\end{textblock*}

\begin{textblock*}{18cm}(1.4cm,1.5cm) % {block width} (coords) 
\textcolor{GJBlue}{\footnotesize \\ Artificial Intelligence formulated this projection for compatibility purposes from the original article published at Global Journals. However, this technology is currently in beta. \emph{Therefore, kindly ignore odd layouts, missed formulae, text, tables, or figures.}}
\end{textblock*}


\let\tabcellsep& 	 	 		 
\section[{Introduction}]{Introduction}\par
urface and internal waves coexist in the ocean with stratification development. The behaviors of waves in such coexisting fields of surface and internal waves show more complicated characteristics than those which exist individually. For instance, the traveling time for a distant tsunami is delayed due to the influence of density stratification in the ocean, according to the theoretical analyses for linear waves 1), 2) . Fructus and Grue 3) used a pressure field for two-layer fluids sandwiched by two fixed horizontal plates, to obtain the surface waves caused by large-amplitude internal waves. A coexisting field of surface and internal waves can be established even in nearshore zones, where surface long waves have great influence on sediment motion and coastal structures as an external force, and conversely, internal waves may greatly affect the coastal environment through water salinity and temperature. Surface and internal waves, however, have often been studied individually: especially, the nonlinear characteristics of surface and internal waves have been investigated independently by e.g. Longuet-Higgins and Fenton \hyperref[b3]{4)} and Choi and Camassa \hyperref[b4]{5)} . In also the research by Fructus and Grue \hyperref[b2]{3)} mentioned above, the interaction between surface and internal waves has not been considered.\par
In the present study, solitary wave solutions for coexistence fields of surface and internal waves have been numerically calculated using the set of nonlinear wave equations based on the variational principle \hyperref[b5]{6)} for two-layer fluids with a free water surface, to examine the characteristics of surface and internal solitary waves, where the phases of both the steady surface and internal solitary waves are assumed to be the same, with a surface mode or an internal mode. 
\section[{II.}]{II.} 
\section[{Fundamental Equations}]{Fundamental Equations}\par
The motion in two-layer inviscid and incompressible fluids is assumed to be irrotational. The upper and lower layers are called the first and second layers, respectively, and the fluids in each layer do not mix even in motion. The velocity potential ? i in the i-th layer (i = 1 or 2) is expanded into the power series of derivation process of nonlinear surface wave (1)   where Ni is the number of terms and fi,?i is the weightings of the power series.\par
By applying the variational principle, the nonlinear surface/internal wave equations \hyperref[b5]{6)} are obtained as follows:\par
The upper and 1 st layer S ( ) ( )? ? ? ? = ? = 1 0 , , , , i i i i N i i z t f t z ? ? ? ? x x ( ) ? ? ( ) 0 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 , 1 1 1 1 1 1 1 , 1 1 1 1 1 = ? ? + ? ? ? ? + + + ? ? ? ? ? ? + ? + + + + + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? f f t t , 0 2 1 2 1 1 1 1 1 1 1 1 1 1 1 , 1 , 1 2 1 1 , 1 , 1 , 1 = + + ? ? + ? ? ? + + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? g f f f f t f ,\textbf{(2)}\par
(3)© 2020 Global Journals 1 Year 2020 
\section[{Global}]{Global}\par
Journal of Researches in Engineering ( ) Volume Xx X Issue III V ersion I E Abstract-The numerical solutions for the coexisting fields of surface and internal solitary waves have been obtained, where the set of nonlinear equations based on the variational principle for steady waves are solved using the Newton-Raphson method. The relative phase velocity of surface-mode solitary waves is smaller in the coexisting fields of surface and internal solitary waves than in the cases without the coexistence of internal waves. The relative phase velocity of internal-mode solitary waves is also smaller in the coexisting fields of surface and internal solitary waves than in the cases without surface waves. The interfacial position of an internalmode internal solitary wave in a coexisting field of surface and internal waves can exceed the critical level determined in the corresponding case without a surface wave. The wave height ratio between internal-mode surface and internal solitary waves is smaller than the corresponding linear shallow water wave solution, and the difference increases, as the relative wave height of internal-mode internal solitary waves is increased. vertical position z, in the manner similar to that for the equations \hyperref[b6]{7)} , as , Japan. e-mail: kyamashita@irides.tohoku.ac.jp (4)   The lower and 2 nd layer (5) In this study, we focus on solitary waves, such that the number of terms for the expanded velocity potential expressed by Eq. (  {\ref 1}) is three for both upper and lower layers, i.e., N 1 = N 2 = N = 3, based on the accuracy verification \hyperref[b7]{8)} for the surface and internal solitary waves obtained using the fundamental equations.(6) 0 2 1 2 1 1 , 1 , 1 2 1 1 , 1 , 1 , 1 1 1 1 1 1 1 1 1 1 1 = + + + ? ? + ? ? ? + + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? p g f f f f t f , ( ) ? ? ( ) 0 1 1 1 2 2 2 2 2 2 2 2 2 2 2 , 2 1 1 2 2 2 2 , 2 1 1 2 2 = ? ? + ? ? ? ? + + + ? ? ? + ? + + + + + ? ? ? ? ? ? = ? + + + + ? ? + ? ? ? + + gh p g f f f f t f ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? 
\section[{a) Determinant in the Newton-Raphson method}]{a) Determinant in the Newton-Raphson method}\par
For the propagation of nonlinear surface/internal waves, the fundamental differential equations, i.e., Eqs.\par
(2), (3), (5), and (7), are transformed to finite difference equations, which are solved using an implicit scheme \hyperref[b8]{9)} . In the present study, numerical solutions for surface/internal solitary waves are obtained using the method introduced by Yamashita and Kakinuma \hyperref[b7]{8)} , where the Newton-Raphson method is applied to solve the fundamental equations for steady waves in a coexisting fields of surface and internal waves. We substitute the advection equation ?F/?t = ?C ?F/?x into the time derivative terms of Eqs. (2), (3), (5), and (7), and then solve the resulting nonlinear wave equations for steady waves traveling in the direction of the x-axis, where C is the phase velocity of the waves, and the physical quantity F is the water surface displacement ?, the interface displacement ?, and the weightings of the expanded velocity potential f i,? . In this method, an arbitrary phase velocity C is given, and these unknown physical quantities for a steady wave with phase velocity C are evaluated using the Newton-Raphson method. Note that in the resulting equations for steady waves, the physical quantities F are functions of only x, for the time derivative terms are eliminated.\par
For the discretization in the Newton-Raphson method, the second-order central finite difference is used for spatial differentiation. The computational domain is the region of 1 ? m ? M, where m is grid point number. The grid points of m = 0 and m = M + 1 are virtual grid points for the central finite difference at the lateral boundaries.\par
( )0 2 1 2 1 1 2 1 2 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 , 1 , 1 2 1 1 , 1 , 1 , 1 2 1 1 2 1 , 2 , 2 2 2 2 , 2 , 2 , 2 = ? ? ? + ? ? ? ? ? ? + ? ? ? + ? ? ? ? ? ? ? ? ? + + ? ? + ? ? ? + + ? + + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? f f f f t f h g f f f f t f (7)\par
.\par
The method to solve the determinant J? = D, which represents the simultaneous difference equations obtained by the discretization above, is the Gaussian elimination method, with partial pivots of high computational stability, where J = J(m) (m = (1, 2, ??? , M)) is the\par
Jacobian matrix, and ? = ?(m) is a column vector composed of the difference ?F between the numerical solution F at the k th and that at the (k + 1) th iterative calculations for convergence.\par
where ?, ?, b, p, h 1 , and ? i are the water surface displacement, interface displacement, seabed position, pressure at the interface, the upper-layer thickness in still water, and fluid density of the i-th layer, respectively. The fluid density ? i is constant in each layer. The It should be noted that the sum rule of product is used for the subscripts ? i , ? i , and ? i : for example, ? 1 in the first term on the left-hand side of Eq. (  {\ref 3}) is the power of ?.\par
From Eqs. (  {\ref 4}) and (  {\ref 6}), p is eliminated to obtain the following equation:\par
The number of elements of the Jacobian matrix J is \{(2 + 2 N) M\} \hyperref[b1]{2} . For example, if the number of grids in the computational domain is 2,500, the total number of elements is about 400 million, such that it is not efficient to store the Jacobian matrix J in one array, from the viewpoint of memory capacity. Therefore, considering that the Jacobian matrix J is a band matrix, we secure only both the elements required to the pivot operation and those of the Jacobian matrix J corresponding to ?(m) = (?f 1,? , ?? , ?? , ?f 2,? ) m for one computational grid point, such that the Jacobian matrix Consequently, the number of elements has been reduced to around 640,000, and the calculation efficiency could be improved significantly.J is composed of (2 + 2 N) × 4 (2 + 2 N) × M. horizontal partial differential operator ? is (?/?x, ?/?y), 
\section[{b) Initial values in the Newton-Raphson method}]{b) Initial values in the Newton-Raphson method}\par
The initial values in the Newton-Raphson method are the surface and interface profiles, as well as the velocity potential, obtained through the KdV theory for small amplitude solitary waves. In two-layer fluids, there are two types of solitary waves with different restoring forces: solitary waves with a surface-wave mode due to gravity, as sketched in Fig. \hyperref[fig_1]{1}, and solitary waves with an internal-wave mode owing to the effective gravity between the two layers, as illustrated in Fig. \hyperref[fig_2]{2}. For the former, the initial values in the Newton-Raphson method are the KdV solutions for a one-layer fluid, and for the latter, those are the KdV solutions for two-layer fluids, the upper surface of which contacts with a fixed horizontal plate.   
\section[{c) Lateral boundary conditions for approximating solitary waves in the finite domain}]{c) Lateral boundary conditions for approximating solitary waves in the finite domain}\par
Solitary waves have the property that the horizontal gradient dF/dx of the physical quantity asymptotically approaches zero at a distance in the horizontal direction. In the numerical calculation, however, the target domain is a finite region, such that the property should be described using boundary conditions. First, as a boundary condition of the calculation using the central finite difference, we assumed dF/dx = 0 for the physical quantities F at the virtual grid points, i.e., m = 0 and m = M + 1, and then the calculation diverged immediately. Second, although we extrapolated the physical quantities F at the virtual grid points m = 0 and M + 1 using the first-or secondorder approximation, the calculation also diverged. displacement ? near the boundaries oscillates without asymptotically approaching zero toward the boundary, which means that in order to obtain stable solutions, it is necessary to suppress such oscillation and express that dF/dx approaches zero toward the boundaries. Finally, we adopted F 0 = F 1 and F M+1 = F M , which means that the gradient of physical quantities in the virtual regions adjacent to the boundaries is assumed to be zero, although it does not mean dF/dx = 0 at the boundaries. For example, dF/dx at the boundary m = 1 is expressed as (F 2 ? F 1 )/ 2?x, which has the same sign as dF/dx at the position m ? 1.5, and the absolute value is 1/2 of dF/dx at the position m ? 1.5, such that the oscillation due to sign reversal around the boundaries is suppressed, and the property of solitary waves, where dF/dx approaches zero toward the boundary, is approximately expressed.\par
The illustration in Fig. \hyperref[fig_1]{1} is our schematic for surface and internal solitary waves with a surface-wave mode, where the still water depth h = h 1 + h 2 is uniform, and the thickness of the upper layer h 1 is 0.2h in still water. By applying the final method described above, we obtain numerical solutions for surface-mode solitary waves, where the phases of both surface-mode surface and internal solitary waves are assumed to be the same as shown in Fig. \hyperref[fig_1]{1}.\par
The density ratio of the lower and upper layers, ? 2 /? 1 , is 1.02, which is close to the density ratio of seawater and freshwater. The total length of the calculation domain, L, is 100.0h, and the grid width in the x direction, ?x, is 0.05h. Shown in Fig. \hyperref[fig_5]{3} are the numerical results for the water surface profiles of the surface-mode surface solitary waves, where the horizontal and vertical axes indicate horizontal distance from the position of the wave-profile peak and the ratio of surface displacement from the still water level to still water depth h. The ratio of the wave height of the surface solitary waves to still water depth, a s /h, is 0.1, 0.3, and 0.5. Comparing the water surface profiles of the surface solitary waves for the one-layer fluid indicated by the black solid lines and those for the two-layer fluids drawn with the red broken lines, a significant difference is not observed between the two, although the latter is slightly sharpened. 
\section[{IV.}]{IV.} 
\section[{Surface-Mode Solitary Waves}]{Surface-Mode Solitary Waves}\par
Figure  {\ref 4} shows the relationship between the relative representative wavelength of surface-mode surface solitary waves, ? s /h, and the ratio of wave height to still water depth, a s /h, where the red solid line shows the numerical solution for the two-layer fluids, and the black solid and broken lines show the numerical solution and the KdV solution for the one-layer fluid, respectively.\par
The representative wavelength ? s of surface solitary waves is defined by (8)   Fig.  {\ref 4}: Relationship between the relative representative wavelength of surface-mode surface solitary waves, ? s /h, and the ratio of wave height to still water depth, a s /h, where ? s is defined by Eq. (  {\ref 8}); h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.02. Fig.  {\ref 4} indicates that the relative representative wavelength decreases, as the ratio of wave height to still water depth, a s /h, is increased. Although the representative wavelength for the two-layer fluids is slightly shorter than that for the one-layer fluid, there is almost no difference between the two. The representative wavelength from the KdV theory for the one-layer fluid is shorter than those through the numerical calculation for the one-layer fluid and the two-layer fluids, for the wavelength by the KdV theory decreases as the wave height is increased, satisfying the assumption that O(a s /h) = O((h/? s ) 2 ). Conversely, in the derivation process of the set of fundamental equations \hyperref[b5]{6)} , no assumptions are made regarding both the ratio of wave height to water depth and the ratio of water depth to wavelength, when the number of the expansion terms for velocity potential, N, is infinity.   
\section[{Global Journal of Researches in Engineering ( ) Volume Xx X Issue III V ersion I}]{Global Journal of Researches in Engineering ( ) Volume Xx X Issue III V ersion I}\par
Shown in Fig.  {\ref 6} is the ratio of the wave height a i of surface-mode internal solitary waves to the wave height a s of surface-mode surface solitary waves, for the two-layer fluids. Although the numerical solution of a i /a s is close to 0.8, which is the value through the linear theoretical solution for small-amplitude surface solitary waves, the difference between the value of a i /a s through the numerical calculation and that from the linear theory increases, a s /h is increased. Fig.  {\ref 6}: Relationship between the wave height ratio a i /a s and the ratio of wave height to still water depth, a s /h, where a i and a s are the wave height of surface-mode internal and surface solitary waves, respectively; h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.02.\par
In the following cases, the density ratio of lower and upper layers, ? 2 /? 1 , is 1.20. In the numerical calculation, the total length of the calculation domain, L, is 50.0h, and the grid width in the x direction, ?x, is 0.02h. Figure \hyperref[fig_7]{7} shows the numerical results for the water surface profiles of surface-mode surface solitary waves, where the ratio of the wave height of surface solitary waves to still water depth, a s /h, is 0.5. The distance between the front and back surfaces of the wave profile at each height of the surface-mode surface solitary wave is shorter in the coexisting field of surface and internal waves than in the case without internal waves. is smaller than the KdV solution, and the numerical solution of C/C s,0 is smaller for the two-layer fluids than for the one-layer fluid. The difference ?(C/C s,0 ) between the numerical solution for the one-layer fluid and that for the two-layer fluids is 2.0 × 10 ?3 , 1.8 × 10 ?3 , 1.4 × 10 ?3 , and 6.0 × 10 ?5 , when a s /h = 0.1, 0.3, 0.5, and 0.6, respectively, where ?(C/C s,0 ) decreases as a s /h is increased.\par
Figure \hyperref[fig_6]{5} shows the relationship between the relative phase velocity C/C s,0 and the ratio of wave height to still water depth, a s /h, for the surface-mode surface solitary waves, where C s,0 = ?ð??"ð??"? is the phase velocity of linear shallow water waves for a one-layer fluid. In Fig. \hyperref[fig_6]{5}, the red solid line shows the numerical solution for the two-layer fluids, and the black solid and broken lines indicate the numerical solution and the KdV solution, respectively, for the one-layer fluid. The relative phase velocity C/C s,0 through the numerical calculation  Figure  {\ref 8} indicates the relationship between the relative representative wavelength of surface-mode surface solitary waves, ? s /h, and the ratio of their wave height to water depth, a s /h, where the thick and thin lines show the numerical solutions for the two-layer fluids and for the one-layer fluid, respectively. The representative wavelength ? s is defined by Eq. (  {\ref 8}). The relative representative wavelength decreases, as the ratio of wave height to still water depth, a s /h, is increased, as in the case shown in Fig.  {\ref 4}.\par
Fig.  {\ref 8}: Relationship between the relative representative wavelength of surface-mode surface solitary waves, ? s /h, and the ratio of their wave height to water depth, a s /h, where ? s is defined by Eq. (  {\ref 8}); h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.20.\par
Figure  {\ref 9} shows the relationship between the relative phase velocity C/C s,0 and the ratio of wave height to water depth, a s /h, for surface-mode surface solitary waves, where C s,0 = ?ð??"ð??"? is the phase velocity of linear shallow water waves for a one-layer fluid. The numerical solution for relative phase velocity C/C s,0 is smaller for the two-layer fluids than for the one-layer fluid, which is the same as in the case shown in Fig. \hyperref[fig_6]{5}. 
\section[{Global}]{Global}\par
Journal of Researches in Engineering ( ) Volume Xx X Issue III V ersion I E Fig.  {\ref 9}: Relationship between the relative phase velocity C/C s,0 and the ratio of wave height to water depth, a s /h, for surface-mode surface solitary waves, where C s,0 = ?ð??"ð??"? is the phase velocity of linear shallow water waves for a onelayer fluid; h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.20. Shown in Fig. \hyperref[fig_1]{10} is the ratio of the wave height a i of surface-mode internal solitary waves to the wave height a s of surface-mode surface solitary waves, for the two-layer fluids. The surface-mode wave height ratio a i /a s decreases, as the relative wave height of surfacemode surface solitary waves, a s /h, is increased, as in the case shown in Fig.  {\ref 6}. Fig. \hyperref[fig_1]{10}: Relationship between the wave height ratio a i /a s and the ratio of wave height to still water depth, a s /h, where a i and a s are the wave height of surface-mode internal and surface solitary waves, respectively; h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.20. 
\section[{V.}]{V.} 
\section[{Internal-Mode Solitary Waves}]{Internal-Mode Solitary Waves}\par
Illustrated in Fig. \hyperref[fig_2]{2} are internal-mode surface and internal solitary waves, where the still water depth h is uniform, and the thickness of the upper layer h 1 is 0.2h in still water. By applying the same method, the numerical solutions for internal-mode solitary waves are obtained, where the phases of both internal-mode surface and internal solitary waves are assumed to be the same as shown in Fig. \hyperref[fig_2]{2}. The total length of the calculation domain, L, is 25.0h, and the grid width in the x direction, ?x, is 0.005h. First, the density ratio of the lower and upper layers, ? 2 /? 1 , is 1.02.\par
The numerical solutions for the interface profiles of internal-mode internal solitary waves are shown in Fig. \hyperref[fig_10]{11}. The red lines indicate the interface profiles for the coexisting field of both surface and internal solitary waves, where the ratio of wave height to upper-layer thickness in still water, a i /h 1 , is 0.15, 0.5, and 1.0, as well as 1.493, which is the maximum value obtained by numerical calculation. On the other hand, the black line shows the numerical solution for the interface profile of the internal solitary wave with the obtained maximum wave height, where the upper surface is in contact with a fixed horizontal plate. In the absence of a free water surface, the downward convex interface of stable internal waves cannot appear below the height of (z c + h 1 )/h 1 = ?1.488, which is called the critical level \hyperref[b9]{10)} . Figure \hyperref[fig_10]{11}, however, indicates that (? min + h 1 )/h 1 = ?1.493, such that the interfacial minimum position ? min can exceed the critical level, when the free water surface coexists.  Fig. \hyperref[fig_9]{12}: Relationship between the relative representative wavelength ? i /h 1 and the ratio of wave height to upper layer thickness in still water, a i /h 1 , for internal-mode internal solitary waves, where the representative wavelength ? i is defined by Eq. (  {\ref 9}); h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.02. the representative wavelength of internal-mode solitary waves in the coexistence field of surface and internal waves is shorter than that for the case without the coexistence of surface waves. These numerical solutions are larger than the corresponding KdV solution, which is similar to surface-mode surface solitary waves shown in Fig.  {\ref 4}. 
\section[{Global}]{Global}\par
Figure \hyperref[fig_5]{13} shows the relative phase velocity C i /C i,0 of internal-mode solitary waves, where is the phase velocity of linear internal shallow water waves without the coexistence of surface waves. As indicated in Fig. \hyperref[fig_5]{13}, the relative phase velocity C i /C i,0 decreases in the coexistence field of surface and internal waves than in the case without the coexistence of surface waves, where the difference between the two decreases as a i /h 1 is increased, as for the case of surface-mode solitary waves shown in Fig. \hyperref[fig_6]{5}. . As shown in Fig. \hyperref[fig_9]{12}, the numerical solution for Fig. \hyperref[fig_5]{13}: Relationship between the relative phase velocity C i /C i,0 and the ratio of wave height to upper-layer thickness in still water, a i /h 1 , for internal-mode internal solitary waves, where C i,0 is the phase velocity of linear internal shallow water waves without the coexistence of surface waves; h 2 /h 1 = 4.0 and ? 2 /? 1 = 1.02. internal-mode surface solitary waves to that of internalmode internal solitary waves, a s /a i . The wave height ratio a s /a i decreases, as the ratio a i /h 1 is increased. Conversely, the wave height ratio a s /a i from the linear shallow water wave theory for the coexisting field of surface and internal waves does not depend on the ratio a i /h 1 , for a s /a i = (1? ? 1 /? 2 ) h 2 /h = 0.016. Second, we compare the numerical solutions for two cases, where the density ratio of the lower and upper layers, ? 2 /? 1 , is 1.02 and 1.20. Figure \hyperref[fig_6]{15} shows the relative representative wavelength ? i /h 1 for internalmode internal solitary waves, where ? i is defined by Eq. (9). As shown in Fig. \hyperref[fig_6]{15}, although the representative wavelength ? i of internal-mode internal solitary waves in the coexisting field of surface and internal waves is larger in the case where ? 2 /? 1 = 1.02 than in the case where ? 2 /? 1 = 1.20, when a i /h 1 is relatively small, the opposite is true, when a i /h 1 is relatively large. Shown in Fig. \hyperref[fig_11]{14} is the ratio of wave height of   \hyperref[fig_6]{15}: Relationship between the relative representative wavelength ? i /h 1 and the ratio of wave height to upper layer thickness in still water, a i /h 1 , for internal-mode internal solitary waves, where the representative wavelength ? i is defined by Eq. (  {\ref 9}), and h 2 /h 1 = 4.0.( ) ( ) 2 1 1 2 2 1 1 2 0 , i / h h h gh C ? ? ? ? + ? = i 2 2 1 i a dx h L L ? ? + = ? ?\par
Figure \hyperref[fig_1]{16} shows the relative phase velocity C i /C i,0 of internal-mode solitary waves, where is the phase velocity of linear internal shallow water waves without the coexistence of surface waves. The relative phase velocity C i /C i,0 is larger when ? 2 /? 1 = 1.02 than when ? 2 /? 1 = 1.20.\par
Fig. \hyperref[fig_1]{16}: Relationship between the relative phase velocity C i /C i,0 and the ratio of wave height to upper-layer thickness in still water, a i /h 1 , for internal-mode internal solitary waves, where C i,0 is the phase velocity of linear internal shallow water waves without the coexistence of surface waves, and h 2 /h 1 = 4.0.\par
Shown in Fig. \hyperref[fig_7]{17} are the ratios of wave height, a s /a i , where a s and a i are the wave height of internalmode surface and internal solitary waves, respectively. The wave height ratio a s /a i is larger when ? 2 /? 1 = 1.20 than when ? 2 /? 1 = 1.02. The numerical solutions for wave height ratio a s /a i decrease, as the relative wave height a i /h 1 is increased, although that through the linear shallow water wave theory for the coexisting field of surface and internal waves does not depend on the ( ) ( )2 1 1 2 2 1 1 2 0 , i / h h h gh C ? ? ? ? + ? =\par
relative wave height a i /h 1 , for a s /a i = [(? 2 /? 1 ) ? 1] / [(? 2 /? 1 )/(h 2 /h 1 ) + 1], such that a s /a i ? 0.154 when ? 2 /? 1 = 1.20, and a s /a i ? 0.016 when ? 2 /? 1 = 1.02. 
\section[{Global}]{Global}\par
Journal of Researches in Engineering ( ) Volume Xx X Issue III V ersion I E Fig. \hyperref[fig_7]{17}: Relationship between wave height ratio a s /a i and relative wave height a i /h 1 , where a s and a i are the wave height of internal-mode surface and internal solitary waves, respectively, and h 2 /h 1 = 4.0. 
\section[{VI.}]{VI.} 
\section[{Conclusions}]{Conclusions}\par
The numerical solutions for the solitary waves in the coexisting fields of surface and internal waves were obtained for the two-layer fluids with a free water surface, where the phases of both the steady surface and internal solitary waves were assumed to be the same, with a surface mode or an internal mode. The set of nonlinear equations based on the variational principle for steady waves were solved using the Newton-Raphson method.\par
The relative phase velocity of surface-mode solitary waves was smaller in the coexisting fields of surface and internal waves than in the cases without the coexistence of internal waves. The difference in the relative phase velocity between the two decreased, as the relative wave height of surface-mode surface solitary waves was increased.\par
The relative phase velocity of internal-mode solitary waves was also smaller in the coexisting fields of surface and internal waves than in the cases without the coexistence of surface waves. The difference in the relative phase velocity between the two decreased, as the relative wave height of internal-mode internal solitary waves was increased.\par
The interfacial position of the internal-mode internal solitary waves in the coexisting fields of surface and internal waves exceeded the critical level determined in the cases without the coexistence of surface waves.\par
The wave height ratio between internal-mode surface and internal solitary waves was smaller than the corresponding linear shallow water wave solution, and the difference increased, as the relative wave height of internal-mode internal solitary waves was increased.\begin{figure}[htbp]
\noindent\textbf{}\includegraphics[]{image-2.png}
\caption{\label{fig_0}}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{1}\includegraphics[]{image-3.png}
\caption{\label{fig_1}Fig. 1 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{2}\includegraphics[]{image-4.png}
\caption{\label{fig_2}Fig. 2 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{}\includegraphics[]{image-5.png}
\caption{\label{fig_3}}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{3}\includegraphics[]{image-6.png}
\caption{\label{fig_5}Fig. 3 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{5}\includegraphics[]{image-7.png}
\caption{\label{fig_6}Fig. 5 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{7}\includegraphics[]{image-8.png}
\caption{\label{fig_7}Fig. 7 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{}\includegraphics[]{image-9.png}
\caption{\label{fig_8}}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{12}\includegraphics[]{image-10.png}
\caption{\label{fig_9}Figure 12}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{11}\includegraphics[]{image-11.png}
\caption{\label{fig_10}Fig. 11 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{14}\includegraphics[]{image-12.png}
\caption{\label{fig_11}Fig. 14 :}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{}\includegraphics[]{image-13.png}
\caption{\label{fig_12}Fig.}\end{figure}
 \begin{figure}[htbp]
\noindent\textbf{} \par 
\begin{longtable}{P{0.03207547169811321\textwidth}P{0.8179245283018868\textwidth}}
III.\tabcellsep Calculation Method for Steady\\
\tabcellsep Wave Solutions in a Coexisting Field\\
\tabcellsep of Surface and Internal\\
\tabcellsep Solitary Waves\end{longtable} \par
  {\small\itshape [Note: GlobalJournal of Researches in Engineering ( ) Volume Xx X Issue III V ersion I]} 
\caption{\label{tab_0}}\end{figure}
 			\footnote{© 2020 Global Journals} 		 		\backmatter  			  				\begin{bibitemlist}{1}
\bibitem[Kakinuma (ed.) ()]{b5}\label{b5} 	 		\textit{A set of fully nonlinear equations for surface and internal gravity waves},  		 			T Kakinuma 		.  		Coastal Engineering V (ed.)  		2001. WIT Press. p. .  	 
\bibitem[Tsai et al. ()]{b0}\label{b0} 	 		‘Estimating the effect of earth elasticity and variable water density on tsunami speeds’.  		 			V C Tsai 		,  		 			J.-P Ampuero 		,  		 			H Kanamori 		,  		 			D J Stevenson 		.  	 	 		\textit{Geophys. Res. Lett}  		2013. 40 p. .  	 
\bibitem[Choi and Camassa ()]{b4}\label{b4} 	 		‘Fully nonlinear internal waves in a two-fluid system’.  		 			W Choi 		,  		 			R Camassa 		.  	 	 		\textit{J. Fluid Mech}  		1999. 396 p. .  	 
\bibitem[Fructus and Grue ()]{b2}\label{b2} 	 		‘Fully nonlinear solitary waves in a layered stratified fluid’.  		 			D Fructus 		,  		 			J Grue 		.  	 	 		\textit{J. Fluid Mech}  		2004. 505 p. .  	 
\bibitem[Nakayama and Kakinuma ()]{b8}\label{b8} 	 		‘Internal waves in a two-layer system using fully nonlinear internal-wave equations’.  		 			K Nakayama 		,  		 			T Kakinuma 		.  	 	 		\textit{Int. J. Numer. Meth. Fluids}  		2010. 62 p. .  	 
\bibitem[Funakoshi and Oikawa ()]{b9}\label{b9} 	 		‘Long internal waves of large amplitude in a two-layer fluid’.  		 			M Funakoshi 		,  		 			M Oikawa 		.  	 	 		\textit{J. Phys. Soc. Jpn}  		1986. 55 p. .  	 
\bibitem[Longuet-Higgins and Fenton ()]{b3}\label{b3} 	 		‘On the mass, momentum, energy and circulation of a solitary wave’.  		 			M S Longuet-Higgins 		,  		 			J D Fenton 		.  	 	 		\textit{Proc. R. Soc. Lond}  		1974. II  (1623)  p. .  	 
\bibitem[Yamashita and Kakinuma ()]{b7}\label{b7} 	 		‘Properties of surface and internal solitary waves’.  		 			K Yamashita 		,  		 			T Kakinuma 		.  	 	 		\textit{ASCE, waves. 45, 15 pages},  				 			Coastal Eng 		 (ed.)  		2015.  	 
\bibitem[Isobe ()]{b6}\label{b6} 	 		‘Time-dependent mild-slope equations for random waves’.  		 			M Isobe 		.  	 	 		\textit{Coastal Engineering},  				1994. 1995. ASCE. p. .  	 
\bibitem[Watada ()]{b1}\label{b1} 	 		‘Tsunami speed variations in densitystratified compressible global oceans’.  		 			S Watada 		.  	 	 		\textit{Geophys. Res. Lett}  		2013. 40 p. .  	 
\end{bibitemlist}
 			 		 	 
\end{document}
