text:
\section{Conclusions}\label{sec:conc}
We have considered the problem of  quantifying rare event statistics due to internal and finite-time instabilities in general dynamical systems. Our analysis is based on the assumption that the conditional dynamics involving rare events are of low-dimensionality, although the total response of the system can still be very high dimensional.  Relying on this setup, we have formulated a conditional decomposition into low-dimensional extreme events caused by internal instabilities, and a high dimensional stochastic background. This decomposition allows for the study of the two components separately (but taking into account mutual interactions), using \emph{different} uncertainty quantification methods that (i)  take into consideration the possibly high-dimensional (broad spectrum) character of the stochastic background, and (ii) the nonlinear and unstable character of the rare events.

The adopted decomposition is in full correspondence with a partition of the phase space into a stable region, where we have no internal instabilities, and a region where non-linear instabilities  lead to extreme transitions with high probability. We quantify the statistics in the stable region using a Gaussian approximation, while the non-Gaussian distributions associated with the intermittently unstable regions of   the phase space, are inexpensively computed through order-reduction methods that take into account the strongly nonlinear character of the dynamics. The probabilistic information for the two domains is analytically synthesized through a total probability argument.

The proposed approach allows for the derivation of the statistics for any quantity of interest in a semi-analytical form, where only a few carefully selected simulations through a reduced order model are sufficient for the accurate determination of the heavy tail structure. For low-dimensional systems the developed framework allows for the derivation of fully analytical forms, while for more complex systems it provides an inexpensive computational method to determine extreme event statistics.

To demonstrate the new method we considered two systems of increasing complexity where non-trivial energy exchanges occur due to internal instabilities, leading to extreme responses. The first application was a two-degree-of-freedom nonlinear system of coupled mechanical oscillators encountered in a viariety of engineering settings. This setup leads to non-Gaussian statistics with heavy tails characterized by qualitatively different regimes. Through numerical experiments we demonstrated  that our method is able to capture   very accurately the core of the response distribution, the  exponential like heavy-tails at extreme values,  and the subsequent sub-exponential decay at very extreme values of distribution.

The second application was a prototype nonlinear envelope equation, that describes the  one-dimensional propagation of deep water waves, where   extreme waves (known as rogue or freak waves)  randomly appear due to  nonlinear focusing in a wavefield. This is an example that is particularly challenging due to its nonlinear, dispersive, and infinite dimensional character. In such case our method is very advantageous as it allows  to   separately quantify the extreme wavegroups from  the background field.      Comparisons with direct   Monte-Carlo simulation  demonstrated the effectiveness of our approach on semi-analytically and inexpensively capturing the heavy-tailed statistics   for   the distribution of   the local wave field maxima, for four different spectra of increasing  heavy-tailed statistics . We also demonstrated the value of our approach on computing the pdf of interest  for different spectrum parameters  in comparison with Monte-Carlo results where the simulation would have to be run anew.

\section*{Acknowledgments}
TPS has been supported through the Office of Naval Research grant ONR N00014-14-1-0520, the Army Research Office Young Investigator Award 66710-EG-YIP, and the DARPA grant HR0011-14-1-0060. MAM and WC have been partially supported  by the first and second grants.

\bibliographystyle{plain}

\section{Introduction}\label{sec:introduction}
Quantifying extreme or rare events is a central issue for many technological processes and natural phenomena. As extreme events, we consider  rare transient responses that  push  the system away from its statistical steady state, which often lead to catastrophic consequences. Complex systems exhibiting rare events include (i) dynamical systems found in nature, such as the occurrence of rare climate events~\cite{hense06, Frei2001, Easterling2000} and  turbulence~\cite{Majda_filter, majda1997, qi15, tong15, Majda2014a}, formation of freak water waves in the ocean \cite{Dysthe08, Olagnon2005, Kharif2003, trulsen1996, cousinsSapsis2015_JFM, cousins_sapsis}; but also (ii) dynamical systems in engineering applications involving mechanical components subjected to stochastic loads~\cite{naess_book, extreme_value09, paik03, Shadman}, ship rolling and capsizing~\cite{belenky07, Kreuzer, arnold_l}, critical events in power grids~\cite{Kishore12, Pourbeik, Susuki2012a, kundur94}, as well as chemical reactions and conformation changes in molecules~\cite{Wigner1938, Pratt1986, Bolhuis2003, E2010}.

For many systems of practical interest like those above, it has now been well established that rare transitions occur frequently enough   that they are of critical importance. These intermittent events are randomly triggered while the system evolves on the (stochastic) background attractor, and they are subsequently governed, primarily, by the spatiotemporally local and strongly nonlinear dynamics associated with finite-time instabilities. Systems with these properties  pose a significant  challenge for uncertainty quantification schemes~\cite{majda_branicki_DCDS} and in recent years a wide range of research efforts has taken place in various fields towards the quantification and short-term prediction of rare events in complex dynamical systems.

The quantification of rare events is one of the most fundamental problems in chemistry. Chemical reactions, conformation changes of molecules, and quantum tunneling are examples of rare events~\cite{Wigner1938, Pratt1986, Bolhuis2003, E2010}. These events are rare because the system has to overcome certain barriers of energetic or entropic nature in order to move from one stable state to the other. The usual setup for modeling such systems is their formulation in terms of a Langevin equation, i.e. a dynamical system with some non-quadratic potential that has multiple equilibria excited by white noise (see for example~\cite{E2010}). The goal then is to study barrier-crossing events by computing transition rates as well as shortest paths between states. For such systems the classical transition state-theory (TST)~\cite{Eyring1935,  Wigner1938} has been successful in providing the foundation for the development of computational tools that determine transition trajectories between different states. However, important limitations for transition state theory may occur when the system potential is not smooth. In this case it is essential to seek for transition tubes (i.e. ensembles of transition trajectories). The statistical framework to analyze such transition-path ensemble is known as the transition-path theory (TPT)~\cite{E.2006, Metzner2006, Metzner2009} and it has been applied successfully for applied to interesting and challenging problems in a variety of areas, for example, in chemistry, biology, and material science.

Although successful on quantifying transitions between different states, path theory can have limitations when   considering dynamical systems that exhibit rare responses due to the occurrence of intermittent instabilities (as opposed to multiple equilibria) that lead to strong energy transfers between modes, as  is the case in turbulence or nonlinear waves.  In such cases, the rare event is not the result of a transition that takes the system from one state to the other. Rather, it is the result of intermittent instabilities that can `push' the system away from its statistical steady state to a dynamical regime with a strongly transient character. The period that the system spends away from its statistical steady state attractor as well as the distance from the attractor usually takes on a continuous range of values, depending on the intensity and duration of the instability. This situation is completely different from the setup involving transitions between discrete states, for which TST and TPT theories have been developed.

Large deviations theory~\cite{varadhan84, varadhan08, dembo00, stroock84}  is a powerful method for the probabilistic quantification of extreme events in sequences of probability distributions. It has also been applied in the context of stochastic differential equations, known as Freidlin-Wentzell theory~\cite{freidlin98}, as well as for stochastic partial differential equations~\cite{chow92, sowers92, sritharan06}. In this case, the method essentially provides us with rates of convergence to probabilistic limits. For example, for a dynamical system excited by very low intensity noise, large deviations theory   gives closed form expressions bounding the probability for a big deviation of the stochastic solution from the completely unperturbed solution, i.e. a probabilistic characterization of the stochastic solution relative to the asymptotic (deterministic) limit. Despite its importance, it is not straightforward to apply this framework in order to quantify extreme events that rise out of the steady state attractor of the system due to the occurrence of intermittent instabilities, which is the problem that we are interested here.

For the probabilistic description of rare events in phenomena characterized by intermittent instabilities, the analysis is usually limited to the statistical examination of  observed statistics. For example, in ocean engineering, where it is important to analyze the probability of upcrossings and maxima for various quantities of interests (e.g. wave elevation or mechanical stresses), the standard setup involves the adoption of globally stable dynamics, for which many techniques have been developed  (see e.g.~\cite{Naess1982, Soong_Grigoriou93, naess_book}). There are numerous technical steps involved in this case that lead to elegant and useful results, but the starting point is usually the assumption of stationarity in the system response, which is not a valid hypothesis for intermittently unstable systems.

Extreme value theory~\cite{extreme_value09, Leadbetter, Thomas_extremes, Embrechts12, Galambos} is also a widely applicable method which focuses on thoroughly analyzing the extreme properties of stationary stochastic processes following various distributions. However, even in this case the analysis does not take into account any information about the unstable character of the system dynamics and is usually restricted to very specific forms of correlation functions for the response statistics~\cite{extreme_value09, Leadbetter}. To this end, it is not surprising that   for a large range of complex systems exhibiting intermittent characteristics, the Monte-Carlo method is  the only reliable computational approach to arrive at   accurate  estimates for the tails. However, for high-dimensional systems this direct approach  is usually prohibitively expensive for practical purposes.

Recently, there have been  efforts to quantify the heavy tail statistical structures of systems undergoing  transient instabilities. In \cite{mohamad2015} a probabilistic decomposition was utilized to obtain analytical approximations for the full probability density function (pdf), including the tail structure, of systems subjected to parametric (or multiplicative) noise with correlated characteristics. In~\cite{tong15} the intermittent behavior of turbulent diffusion models with a mean gradient was rigorously quantified, while in~\cite{qi15} the capacity of imperfect models to capture intermittent behavior of turbulent systems was studied.

The goal of this work is to  present a computational framework for   efficiently  quantifying  the statistical characteristics of extreme events. We focus on systems where uncertainty interacts with   system dynamics to produce~\emph{intermittent extreme events}, i.e. sporadically occurring large amplitude responses,  which give   rise to non-Gaussian statistics. Uncertainty could be due to    the initial data, parameters, or the dynamical system itself. The core idea of our probabilistic decomposition-synthesis method is the separation of   intermittent events from the background stochastic attractor,  in the spirit of the work done in~\cite{mohamad2015}, for low-dimensional systems. This decomposition  allows  us to apply~\emph{different} uncertainty quantification schemes for the two regimes (the background and the intermittent component). The background component, although potentially very high-dimensional, can be efficiently described by uncertainty quantification schemes that resolve low-order statistics. On the other hand, the intermittent component,  can be described in terms of a low-dimensional representation through a small number of  localized modes. The  probabilistic information from these two regimes  is synthesized according to a total probability decomposition argument, in order to approximate the   heavy-tailed probability distributions  for functionals of interest. Thus, the core of the approach relies on the assumption that heavy-tails are primarily due to the action of  intermittent instabilities, whereas the `main mass' of the probability density is due to the the background component.

 To illustrate the method, we apply the developed framework on two systems that exhibit intermittently extreme responses and estimate their statistical distributions. The first example is a   nonlinear system of coupled random oscillators, which serves to illustrate the various steps of the method  in a simple  prototype,  where the mechanism behind intermittent instabilities  is easy to understand  although  the statistical characteristics  of the response are highly  non-trivial. The second example, is  a nonlinear envelope equation characterizing the propagation of unidirectional water waves in  deep water. The   benefits  of the method are  highlighted in this complex example, where we are able to capture the statistics at a fraction of the cost of direct numerical simulations. In both cases, we compare our estimates for the statistical distribution with direct Monte-Carlo results, and  illustrate   the  performance of the  approach.

The paper is structured as follows. In~\cref{sec:problem_setup} we describe the problem setup. In~\cref{sec:decomposition} we detail the various steps for  the proposed decomposition-synthesis method. \Cref{sec:rare_events} we    analyze  and detail the  various order reduction schemes for the statistical quantification for  rare events and the background attractor.  In~\cref{sec:osc_eg} we illustrate the method to the first example of a coupled, nonlinear system of random oscillators. In~\cref{sec:mnls}, we demonstrate the method on the second example to nonlinear  waves in  deep water.

\section{Problem setup}\label{sec:problem_setup}

Let $(\Omega,\mathcal{B},\prob)$ be a probability space, where $\Omega$ is the sample space with $\omega \in\Omega$ denoting an elementary event of the sample space, $\mathcal{B}$ the associated $\sigma$-algebra of the sample space, and $\prob$ a probability measure. (In the following, $\prob_X$ will denote the probability measure associated with a variable $X$ and $\pdf_X$ the associated  probability density function, pdf, if appropriate).

Let the \emph{dynamical system} of interest be governed by the following stochastic Partial Differential Equation (SPDE):
\begin{equation}
    \frac{\partial u(x,t)}{\partial t} = \mathcal N [ u(x,t); \omega], \quad x \in D,\;\; t \in [0,T], \; \; \omega \in \Omega,
    \label{eq:general_eq}
\end{equation}
where $\mathcal N$ is a general (nonlinear) differential operator with appropriate boundary conditions. We assume the initial state at $t=t_0$  is random and described by $u(x,t_0) = u_0(x;\omega)$. In what follows we  utilize a spatial inner product, denoted for two arbitrary functions $u(x)$ and $v(x)$ as $u \cdot v$.
Extreme events are meant in terms of a norm $\lVert \,\cdot\, \rVert$ (e.g. the spatial supreme norm). In addition, the ensemble average of a random quantity $f(\omega)$ is   denoted by  $\bar f $.

Here we are interested in determining the statistical distribution for a \emph{quantity of interest} given by a functional of the solution $u(x,t)$ or as a solution of another dynamical system subjected to $u(x,t)$:
\begin{equation}\label{eq:quantity_interest}
    q=q[u(x,t)], \quad \text{or}  \quad   \frac{dq}{dt} = \mathcal M[q, u(x,t)].
\end{equation}Examples of such quantities could be properties of inertial tracers in turbulent flows, stresses for ocean structures subjected to water waves, or strains of mechanical components subjected to parametric and/or additive excitations. The computational cost associated with the estimation of the heavy-tailed statistics for such quantities is vast given the fact that the dynamical system~\eqref{eq:general_eq}  is characterized by the occurrence of rare events, which define the heavy-tail properties for $q$. Therefore, application of direct methods, such as Monte-Carlo simulations are prohibitively expensive for these problems.

\section{The probabilistic decomposition-synthesis method}\label{sec:decomposition}

We describe a systematic method to quantify the non-Gaussian response statistics (due to rare  events caused by intermittent instabilities) in a computationally efficient manner. More specifically, in this section we  first give an overview of the steps involved for the decomposition-synthesis method; more details are provided in the following section.

\subsection*{Decomposition and assumptions}

We assume that all rare event states due to internal instabilities, defined by the condition $\lVert u\rVert>\zeta,$ with $\zeta$ being the rare event threshold,   `live' in a low dimensional subspace $V_s$. We then decompose the response of the system as:
\begin{equation}\label{eq:decomp}
    u(x,t)=u_b(x,t)+u_r(x,t), \quad \text{with }\;   u_r=\Pi_{V_s}[u],\quad \text{if }\;\lVert u\rVert>\zeta,  \quad\text{and }\;  u_b=u-u_r,
\end{equation}
where $\Pi_{V_s}$ denotes the linear projection to the subspace $V_s$.
A similar decomposition onto a subspace of interest and a background component, but without taking into account the conditioning on rare events, has been utilized successfully for the uncertainty quantification and filtering of turbulent systems~\cite{sapsis_majda_mqgdo, Majda2014}. This conditional decomposition  will allow for the study of the two  components separately (but taking into account mutual interactions), using
\emph{different} uncertainty quantification methods that (i)  take into consideration
the possibly high-dimensional (broad spectrum) character of the stochastic background,
and (ii) the nonlinear and unstable character of   rare events.

\textit{In this work we are interested in rare events that occur  due to transient instabilities}. To this end, we denote as $R_e$ the set of all background states that trigger an instability that   lead to a rare event. In~\cref{fig:phaseDecomp} we demonstrate the adopted decomposition in a three-dimensional system where the rare event subspace is defined by the linear span of $u_3$.
\begin{figure}[tb]
    \centering
    \includegraphics{phasespace}
    \caption{The conditional decomposition~\eqref{eq:decomp} partitions the system response only when a rare event occurs due to an instability. This happens when the state of the system enters the instability region $R_e$. In this example the subspace  associated with rare events due to instabilities is $V_s=\spn\{u_3\}$.}
    \label{fig:phaseDecomp}
\end{figure}

The application of this decomposition onto a stochastic background and rare events relies on the  following assumptions:
\begin{enumerate}
\item The existence of intermittent events have negligible effects on the statistical characteristics of the stochastic attractor and can be ignored when analyzing  the  background state $u_b$.
\item Rare events are statistically independent from each other.
\item Rare events are characterized by low-dimensional dynamics.
\end{enumerate}

The first assumption  allows for the application of closure models or representation
methods that can deal with the high dimensional character of the stochastic background attractor. It expresses the property
that the rare events, although of large magnitude, are localized in time and space and can  induce only negligible modifications to the statistics of the background state.
The second assumption   follows from the rare character of   extreme events.
As for the third assumption,  related to the  low-dimensionality of the dynamics
of rare events,  follows naturally from their spatially or temporally localized
character. We emphasize  that all these assumptions do not imply any restrictions on the dimensionality of the stochastic background state.

\subsection*{Analysis of the various regimes}
The analysis of the two regimes will consist of the following steps:
    \begin{enumerate}
        \item \textbf{Order-reduction} in the subspace $V_s$ in order to model the rare event dynamics, expressed through $u_r.$  Then using the approximation   $u(x,t)\simeq u_r(x,t)$ we will compute the conditional pdf $\pdf(q\mid \lVert u \rVert>\zeta,u_{b}\in R_{e})$, under the condition that  an extreme event occurs due to an internal instability in $R_e$.
        \item \textbf{Quantification of the instability region} $R_e$ using the reduced-order model, by analyzing the conditions that lead to a rare event.
        \item \textbf{Description of the  background dynamics}, expressed through the statistics of $u_b$, which is not influenced by any internal instabilities in $R_e$.  Thus, when the response is dominated only by the background dynamics, we   have  $u(x,t)= u_b(x,t)$ and the pdf for the quantity of interest is given by  $\pdf(q\mid u=u_{b}).$
        \item \textbf{Probability for rare events due to internal instabilities} $\prob(\lVert u \rVert>\zeta,u_{b}\in R_{e})$, which quantifies the total time/space that the response spends in the rare event regime due to the occurrence of instabilities.
\end{enumerate}

\subsection*{Probabilistic Synthesis}

The next step of our technique is to probabilistically synthesize the information obtained from the previous analysis.  Using a total probability argument, in the spirit of~\cite{mohamad2015}, we  obtain the statistics for the quantity of interest  $q$   by
\begin{equation}\label{eq:probLaw}
\pdf(q) =  \pdf(q \mid {\left\Vert u \right\Vert>\zeta,u_{b}\in R_{e}} ) \prob_{r} +  \pdf(q \mid u=u_{b})\left( 1-\prob_{r}\right),
\end{equation}
where $\prob_{r}= \prob(\lVert u \rVert>\zeta, u_{b}\in R_{e}) $ is the probability of a rare event due to an instability.  The first term expresses the contribution of   rare events due to internal instabilities and is the heavy-tailed part of the distribution for $q$. The second term expresses the contribution of the background state and is the main    main probability mass in the pdf for $q$.

This total probability decomposition separates the full response  into the
conditionally extreme response and the conditionally background response,
weighted by their appropriate probabilities. The decomposition~\eqref{eq:probLaw}
separates  statistical quantities according to the total probability law through conditioning on  \emph{dynamical regimes}. In this manner,   our approach connects the statistical quantities that we are interested in with important dynamical regimes that determine the dominant statistical features (e.g. a Gaussian core due to the background state and exponential like heavy-tails due to  intermittent bursts). An outline of all the steps involved is presented in~\cref{fig:methodFlowChart}.
\begin{figure}[tb]
    \centering
    \includegraphics{chart}
    \caption{Outline of the steps  of the decomposition-synthesis method.}
    \label{fig:methodFlowChart}
\end{figure}

For low-dimensional systems, this decomposition
can   lead to   analytical results for the probability distribution of
the response, which would otherwise be impossible to obtain from the dynamical
system and properties of the noise. For more complex systems, the primary
benefits are computational. We can resolve tail statistics with far fewer simulations
than direct Monte-Carlo simulation would require, since \textit{the decomposition allows for the evaluation of the tails by targeted simulations of rare events}, as opposed to the
random sampling of Monte-Carlo simulations.

We emphasize that the expensive computation in the proposed algorithm is the conditional statistics for rare events $\pdf(q \mid {\lVert u \rVert>\zeta,u_{b}\in R_{e}} ).$  Once this has been obtained, it is simple to compute  rare event statistics for different configurations of the background state, since the other (non-rare) quantities are easy to compute multiple times. This contrasts with the Monte-Carlo method, where
sampling must be performed anew when noise and/or system parameters are changed.
Furthermore, with the proposed method the statistics of arbitrary functionals
of the response can be easily obtained with minimal additional computational expense, whereas Monte-Carlo simulations may be prohibitively expensive in this case.

\section{Detailed description of the analysis steps }\label{sec:rare_events}

We now provide the details for the analysis of each component starting from the reduced-order description of rare events and continuing with the statistical quantification of the background dynamics and the probability of rare events.

\subsection{Order reduction of the rare event dynamics}

The first  step for the effective description of the dynamics of rare events is the selection of a reduced-order basis or a reduced-order subspace denoted as $V_s$. These  events are spatially localized structures,   to this end we  employ  localized basis around a neighborhood of an arbitrary point $x_c$ at time $t_c$:
\begin{equation}
\hat {  v }_{ n } ( x-{ x }_{ c },t-t_c ) ,\quad n=1,\ldots,s,
\end{equation}
where $s$ is the dimension of the subspace. Note that the subspace  explicitly depends on the spatiotemporal location $(x_c, t_c)$ of the rare event, which is arbitrary. There are numerous methods with variable complexity that   can be utilized to choose or compute the subspace of rare events. The simplest choice in this case is a steady localized basis, such as Gabor basis or wavelets~\cite{cousins_sapsis}. However, if there is important spatial translation of the rare event during its lifetime then it may be beneficial to utilize adaptive methods for the evolution of the basis elements    $\hat{v}_{ n }$, such as the dynamically orthogonal field equations~\cite{SapsisLermusiaux09, sapsis11a}.

Apart from the selection of an appropriate set of basis functions, there is a variety of options for performing the order-reduction of the dynamics. Here we   discuss two different approaches for obtaining the local dynamics, which we   later illustrate through two specific problems.

\subsubsection{Projection of the dynamics around a zero background state}\label{proj1bk}
The simplest strategy to study the reduced-order dynamics of rare events is to perform a Galerkin projection ignoring the background state $u_b(x,t)$ for the evolution of the rare events inside the reduced-order subspace $V_s$. \textit{Such an assumption is valid if the background state plays a role only for the initial triggering of a nonlinear  instability}, while its small magnitude (compared with the intense local  responses associated with the   nonlinear dynamics in $V_{s}$) has very small influence on the reduced-order dynamics of the rare event.

Therefore, with the proposed order-reduction strategy we study rare events assuming that their evolution is isolated from the stochastic background. The coupling with the background state is introduced only through the background conditions at $t=t_c$, which is the initial condition for the reduced-order dynamics within $V_s$.
Based on this setup, we express the intermittent event in terms of the localized basis that span  the subspace $V_s$:
\begin{equation}\label{eq:projec2}
u_{ r }(x,t)=\sum _{ n=1 }^{ s } { a }_{ n } ( t )  \hat{v}_n(x -  x_c, t - t_c).
\end{equation}
Projecting the local (in the sense of $x_c$) dynamics of the full nonlinear equation on the subspace $V_s$, we obtain the following local dynamical system:
\begin{equation}\label{ic_proj1}
        \begin{aligned}
             \dot a_i &=\mathcal N \biggl[\,  \sum _{n=1}^{s}{a_n \hat v_n } ;\omega \biggr] \cdot {\hat v_i},\quad &i=1,\ldots,s, \\
            a_{i0} &=   (u_b\cdot\hat{v}_i )|_{ t=t_c },   &i=1,\ldots,s,
        \end{aligned}
\end{equation}where we again emphasize that the reduced-order system contains information about the background state only through the initial condition at $t=t_{c}$.
Through this reduced-order description we  can obtain the conditions (in terms of $a_{i0}$ and $\omega$) that define $R_{e}$, the domain of attraction to rare events.

\subsubsection{Projection of the dynamics around the stochastic background attractor}\label{red2}
A more accurate approach for quantifying localized instabilities is the formulation of reduced-order models around the background attractor, instead of a zero background state. This can be critical if the background state   not only triggers intermittent instabilities but also determines their evolution forward in time.

Here we  utilize the assumption that the presence of the localized instability $u_{r}$ does not have   important influence on the evolution of the background state $u_{b}$. We then express the full solution as in decomposition (\ref{eq:decomp}):
\begin{equation}
u(x,t)=\Pi _{V_s^{\bot}} [u(x,t) ]+\sum _{ n=1 }^{ s }{  a
_n } ( t )  \hat v_n (  x-x_c ,t-t_c  ),\\

\end{equation}
where the projection operator  onto the orthogonal complement is defined as:
\begin{equation}
\Pi_{V_s^\bot}[ u ]\equiv u-\sum_{ n=1 }^{ s }{ (u \cdot\ \hat{v} _{ n })  \hat{ v }_{ n }}.
\end{equation}
Using this setup we project the original equation on the basis elements $\hat{ v } _{ i }$ to obtain:
\begin{equation}
{ \dot  a   }_{ i } =\mathcal N\biggl[\Pi _{V_s^{\bot}} [u]+ \sum_{n=1}^s{  a_n  \hat  v_n } ;\omega \biggr] \cdot \hat
v_i - \Pi_{V_s^{\bot}}\biggl[\frac{\partial u}{\partial t}\biggr] \cdot\hat
{ v } _{ i }, \quad
i=1,\ldots,s,
\end{equation}
where the last term on the right hand side vanishes identically. Therefore, the localized dynamical system takes the form:
\begin{equation}
\label{ROE2}
\begin{aligned}
  \dot{a}_i &= \mathcal N\biggl[\Pi _{V^{\bot}_{s}} [u] +
\sum_{n=1 }^{ s }{ {a_n} \hat v_n }; \omega \biggr] \cdot \hat {v_i },\quad  &i=1,\ldots,s,\\
a_{i0}&=   (u_b \cdot \hat v_i  )|_{t=t_c},    &i=1,\ldots,s.
\end{aligned}
\end{equation}
Clearly, if we set the background state to zero, the   dynamical system above reduces to the formulation derived previously. The description of the background state is obtained by projecting the full equation in $V_s^\perp$ and taking into account that the evolution of $u_r$ has negligible effect on the projected dynamics of the background component. However, the infinite dimensional character of the dynamics in the orthogonal complement, $V_s^\perp$, makes it impractical to utilize the full equations for $\Pi _{V^{\bot}_{s}} [u]$. To this end, the reduced-order dynamical system (\ref{ROE2})  should be seen
as a starting point, where appropriate finite-dimensional truncations of
$\Pi _{V^{\bot}_{s}} [u]$ can be utilized that   capture the essential
effects of the stochastic background on  the evolution of the  localized
rare
events. Such an approach will be demonstrated in the problem involving water waves.

\subsection{Quantification of the instability region}\label{sec:domain_attraction}

To simplify the analysis we
  consider the case where the reduced-order dynamics have no important dependence
on $\omega$. We define $R_e$ as the set of background states for which the
reduced-order dynamics in $V_s$ posses at least one finite-time Lyapunov exponent that is positive.
More specifically, let the flow map $\phi^{t}_{t_{c}}:V_s\rightarrow V_s$
that maps any point $\mathbf a_0=\Pi _{V_{s}} [u]$ to its position at time time $t$ under
the effect of the dynamical flow (\ref{ic_proj1}):
\begin{equation}
\phi^{t}_{t_c}:V_s\rightarrow V_s,\quad \mathbf a_0 \mapsto \mathbf a(t,t_c,\mathbf
a_0;u_{b}).
\end{equation}
Note that this dynamical flow may also depend on the background state  $u_b \in V_s^\perp$. For each initial condition $\mathbf a_0$ we define the maximum finite-time
Lyapounov exponent over some finite-interval~$\tau$:
\begin{equation}
\lambda_{\tau}=\max _{i=1,\ldots,s}\log l_{i}\Bigl( \bigl[ \nabla_{\mathbf{a}_0}\phi^{t_c+\tau}_{t_{c}}
\bigr]^{*}\nabla_{\mathbf a_0}\phi^{t_c+\tau}_{t_{c}}\Bigr),
\end{equation}
where $  l_{i}$ denotes the eigenvalue of the Cauchy-Green tensor, which
is by definition symmetric and positive-definite. The domain of attraction
to rare events $R_e$  is then   defined as the set of background states
for which we have expansion in the reduced-order subspace:
\begin{equation}
R_{e}= \bigl\{ \mathbf a_0 \in V_{s},u_b \in V_s^\perp    \bigm\vert  \lambda_{\tau} ( \mathbf{a}_0
;u_{b}) >0 \bigr\}.
\end{equation}

\subsection{Description of the  background dynamics}
Here we discuss methods for the representation of the statistics of the background stochastic attractor. This refers to the part of the decomposition (\ref{eq:decomp})  $u_b$  for which instabilities have no role. There are numerous ways to approach the problem and here we review some of the techniques that can be used.

\subsubsection{Gaussian closure }
Based on our setup,  the stochastic background $u_b$ does not contain  intermittent events due to instabilities; as a consequence, it is reasonable to assume that its statistics can be approximated by a Gaussian distribution. This assumption can be the starting point for the application of closure schemes.

For systems  that are characterized by a stable mean state $\overbar  u $, the finite variance of the steady state attractor is caused by the   external stochastic excitation (see~\cite{sapsis_majda_qgdo} for more details). In this case partial linearization of the dynamics or a Gaussian closure of the infinite system of moment equations (see e.g.~\cite{epstein69})    is  an effective option to capture the conditional statistics of the system in the state $R_e^{c}$, where we have no rare events occurring. For a such situation the governing equations can be linearized to give
\begin{equation}
\frac{\partial u_b}{\partial t} = \mathcal N [ \overbar u;\omega]+ \mathcal L_{\overbar{u}}  u_b, \quad
 \omega \in \Omega,
\end{equation}
where the first term on the right-hand side contains deterministic and stochastic noise terms that are independent of the state $u_b$, while the second term denotes the linearization of the system,
\begin{equation}
 \mathcal L_ {\overbar u} = \frac{\delta \mathcal N}{\delta u}\bigg|_{\overbar u}
\end{equation}
around the mean $\overbar u,$  where all the eigenvalues have a negative real part (stable mean state). For such system it is straightforward to formulate the  second order moment equations and obtain an expression for the Gaussian measure that characterizes the statistics of the background attractor.

Such an approach will not be effective if the system under consideration has persistent instabilities that lead to nonlinear energy transfers between modes, i.e. the mean state is not stable. For such cases other methods may be used to obtain  representations of the background attractor statistics, such mean stochastic models~\cite{sapsis_majda_tur} or quasi-linear Gaussian closures~\cite{sapsis_majda_mqg,sapsis_majda_tur}.

\subsubsection{Gaussian representations based on conditional Monte-Carlo simulation}\label{sec:gausconditional}

A more direct approach involves the numerical simulation of the system and the conditional sampling of the  second-order statistics of the background state using the representation (\ref{eq:decomp}).
It is also often the case that the spectrum of the full random field has been measured or estimated. Such an approach is typical, for example, in water waves or other geophysical systems.

\subsubsection{Analytical approach using the Fokker-Planck-Kolmogorov equation}\label{sec:decomp_analyticalfpk}

For systems having background states that can be modeled by  dynamical systems of a  special form, analytical  descriptions of their stationary probability measure may be available. In particular, for systems excited by white noise one can formulate the corresponding Fokker-Planck-Kolmogorov (FPK) equation defining the evolution of the state pdf~\cite{Sobczyk91,Soong_Grigoriou93}. For the special case of vibrational systems possessing a Hamiltonian structure perturbed by  additive and/or parametric white noise under linear damping, for special conditions, the stationary measure has an explicit form   in terms of the Hamiltonian of the system~\cite{soize88,soize94}. Furthermore, in~\cite{caughey82,zhu90,wang00} an analytical approach is utilized for determining the stationary pdf of more generic vibrational systems, where the steady FPK equation is solved by splitting to simpler partial differential equations.

To demonstrate the analytical approach, we consider the case of a background state described by a collection of decoupled vibrational modes, each governed by an equation of the following form:
\begin{equation}
\ddot u + \delta \dot u + \nabla V(u) =\sigma\dot W(t;\omega).
\end{equation}
Then, the pdf in the statistical steady state  will be given in terms of the Hamiltonian \(\mathbb{H}(u,\dot u)\) of the system as~\cite{Sobczyk91}:
\begin{equation}
\pdf_{u\dot u}(u,\dot u)=C\exp\left(-\frac{2\delta}{\sigma^2}\mathbb{H}(u,\dot u)\right),
\end{equation}
 where $\mathbb{H}(u,\dot u)=\frac{1}{2}\dot u^2 + V(u)$, and $C$ is a normalization
constant. We   utilize this approach in the first example involving a nonlinear vibrational system.

\subsection{Probability for rare events $\prob(\lVert u \rVert>\zeta,u_{b}\in R_e)$}

Having determined the conditional statistics of the background dynamics, the final step is to determine the probability of   rare transitions, i.e. the likelihood of  $\lVert u \rVert>\zeta$ when an instability occurs, i.e. $u_{b}\in
R_e$.  This will be defined for an arbitrary point in space, $x_0$, through the integral:
\begin{equation}
\prob(\lVert u \rVert>\zeta,u_{b}\in
R_e)=\frac{1}{T}\int\limits_{t \in T} \mathbbm{1}(\lVert u(x_{0},t) \rVert>\zeta,u_{b}\in R_e)\, dt,
\end{equation} where  $ \mathbbm{1}$ is the indicator function. This integral measures the duration of rare events compared with the overall time interval.

Note that the above probability is not directly equal to the probability of the background state crossing into  the instability region $R_e$, which is a condition just for the occurrence of a rare event and does not contain information regarding the duration of the rare event, which can last even if the system background has moved outside the instability region.

The probability that we are interested in will be found using information obtained from the reduced-order model
developed previously. In particular, using the reduced-order model we will obtain the temporal extent $T_{e}(u)$  for each rare event that corresponds to any unstable point $u \in R_e$.  Note that for states $u$ not associated with instabilities $T_{e}(u)=0$.

Then the
probability of rare transitions due to instabilities will be approximated by
\begin{equation}\label{eq:rare_trans}
\mathbb{P}(\left\Vert u \right\Vert>\zeta,u_{b}\in
R_e)=\frac{1}{T}\int\limits_{u } \! T_{e}\left(u\right)\pdf_{u_b}(u)\,du ,
\end{equation}
where $\pdf_{u_b}=\pdf(u \ \mid u=u_{b})$ is the conditional pdf of the background state, which contains  no information on rare events.

\section{Unidirectional nonlinear  water waves}\label{sec:mnls}

The second application   involves the problem of local extremes for nonlinear dispersive water waves. In particular, our goal  is to quantify the non-Gaussian, heavy-tailed distribution  of the local wave field maxima  in a nonlinear envelope equation characterizing the propagation of unidirectional water waves. This example has many ocean engineering applications, since quantifying  extreme water waves   is critical   for ocean structures and naval operations  due their catastrophic consequences. Indeed, extreme waves, termed  freak or rogue waves, have caused considerable damage to ships, oil rigs, and human life ~\cite{Dysthe08, Olagnon2005, Kharif2003, trulsen1996}. A wave is termed a freak if its crest-to-trough height exceeds twice the significant wave height $H_s$, with $H_s$ being equal to \emph{four times} the standard deviation of the surface elevation (\cref{fig:exampleFreakWave}). Therefore, such  waves are extreme responses that `live' in the tails of the distribution of the wave elevation.
\begin{figure}[b]
        \centering
        \includegraphics{wavefield.pdf}
        \caption{Example large amplitude freak wave formed due to nonlinear interactions.}
        \label{fig:exampleFreakWave}
\end{figure}

We consider waves on the surface of a fluid of infinite depth and work with approximate equations that govern the dynamics of the wave envelope. In particular, the evolution of a unidirectional, narrow-banded wave field  is described well by the modified Nonlinear Schrodinger (MNLS) equation of Dysthe~\cite{dysthe1979},  a high order approximation of the fully nonlinear model:
\begin{equation}        \label{eq:MNLS}
  \frac{\partial u}{\partial t} +\frac{1}{2}\frac{\partial u}{\partial x} + \frac{i}{8} \frac{\partial^2 u}{\partial x^2} - \frac{1}{16} \frac{\partial^3 u}{\partial x^3} + \frac{i}{2} |u|^2 u + \frac{3}{2} |u|^2 \frac{\partial u}{\partial x} + \frac{1}{4} u^2 \frac{\partial u^*}{\partial x} + i u \frac{\partial \phi}{\partial x} \bigg|_{z=0} = 0,
\end{equation}
where $x$ is space, $t$ is time, and $u$ is the envelope of the modulated carrier wave. The velocity potential $\phi$ at the surface may be expressed explicitly in terms of $u$, giving $\partial \phi/\partial x |_{z=0} = - \mathcal F^{-1}\bigl[ \abs{k} \mathcal F\bigl[\abs{u}^2\bigr] \bigr]/2$, where $\mathcal F$ is the Fourier transform. \Cref{eq:MNLS} has been nondimensionalized with $x = k_0 \tilde x$, $t = \omega_0 \tilde t$, $u = k_0 \tilde u$, where $\tilde x$, $\tilde t$, $\tilde u$  are physical space, time, and envelope variables, with $k_0$ the dominant spatial frequency of the surface and $\omega_0  = \sqrt{ g k_0 }$.  To leading order,  the   nondimensionalized surface elevation around the undistributed level is given by  $\eta(x,t) = \real[ u(x,t) e^{i\left(x-  t\right)} ]$. The local maxima of the surface elevation is described by the modulus of the envelope $|u|$, which is the quantity of interest in this example.

We consider a wave field with with initial characteristics given by
\begin{equation}
        u(x,0) = \sum_{k = -N/2+1}^{N/2} \sqrt{2 \Delta_k F(k\Delta_k)} e^{i(\omega_k x + \xi_k)}, \quad  F(k) = \frac{\epsilon^2}{\sigma \sqrt{2\pi}} e^{\frac{-k^2}{2\sigma^2}}
        \label{eq:gaussSpec}
\end{equation}
where $\xi_k$ are independent, uniformly distributed random phases between $0$ and $2\pi$.

If the Benjamin-Feir Index (BFI), the ratio of steepness to bandwidth $\epsilon / \sigma$ is large enough, we have important probability for the occurrence of rare events due to nonlinear focusing. Such nonlinear focusing is triggered by the energy localization over a specific region, which is the result of
phase differences between the various harmonics~\cite{cousinsSapsis2015_JFM} in the stochastic background $u_b$. These relative phases continuously change primarily due to the effect of linear dispersion and as the BFI increases there is a higher probability for them to result in important energy localization for $u_b$ and a subsequent nonlinear focusing event $u_r$.

It has been shown~\cite{cousinsSapsis2015_PRE} that for unidirectional waves the occurrence of nonlinear focusing of an arbitrary  wavegroup (formed due to linear dispersion) in $u_b$ is controlled by the wavegroup amplitude $A$ and its lengthscale $L$. This fact leads us to the adoption of the wavegroup characteristics  $A$ and $L$ as a way to parameterize  the instability region, i.e. the region where nonlinear focusing occurs.

\subsection{Rare events subspace and  instability region $R_e$}\label{sec:mnls_rareventspace}

 We  describe the  dynamics of focusing wavegroups through a reduced
order subspace, $V_s$, obtained using a proper orthogonal decomposition~\cite{Holmes_et_al96} of
the focusing wavegroups under MNLS dynamics. The proper orthogonal decomposition is appropriate  in this case  as  the order reduction scheme for focusing wavegroups,  since we do not have important spatial translations of the focusing waves and this allows us to capture the dynamics with just a few modes.

To capture the variations in the dynamics for different wavegroup lengthscales and amplitudes, we choose $n = 8$ sets of simulations of the MNLS for various $(A,L)$  that undergo nonlinear focusing. For each simulation we take snapshots $u_{s}(x,t_k; A_i,L_i)$ in time, where for each simulation we ensure that the snapshots are capturing the dominant focusing action. Stacking the snapshots for all the   simulations   gives the matrix $X$
\begin{equation}
        X = \bigl[ u_s(x,t_1; A_1,L_1) \;\; u_s(x,t_2;A_2,L_2)  \;\;\ldots  \;\; u_s(x,t_m; A_n,L_n)\bigr],
\end{equation}
where $m$ is total number of snapshots. We use the method of snapshots~\cite{Sirovich87} to determine the orthonormal POD modes, by solving the $m\times m$ eigenvalue problem $X^T X U = U \Lambda$, for the modes $V = X U \Lambda^{-1/2}$.

By this  procedure we obtain the local basis  $\hat v_{n}(x,t)$ and represent  $u_r$ as
\begin{equation}\label{eq:mnls_pod}
u_{ r }(x,t)=\sum _{n=1}^{s} {a}_{n} (t)  \hat{v}_n(x -  x_c, t - t_c),
\end{equation}
which is conditional upon a background state at $(x_c,t_c)$. We used the projection of the full MNLS around a zero background state as described in~\cref{proj1bk}. We found that a projection upon $s = 14$ modes approximated the dynamics of $V_s$ well across a range of  $L$ and  $A$ (fig. \ref{fig:comp}).
 \begin{figure}[tb]
        \centering
        \includegraphics[width=0.9\textwidth]{comparison.png}
        \caption{Simulation of an extreme wave group for  $A = 0.11,$ $L = 11$, comparing the exact MNLS (left) and the reduced order model (right) with $14$ modes.}
\label{fig:comp}
\end{figure}

With this reduced order description of the dynamics for  the intermittent component $u_r$, we  investigate the evolution of various wavegroups   $(A,L)$ in order to quantify the conditions for the occurrence of rare transitions,   which will determine  $R_e$. More specifically, we consider initial wavegroups of the form
\begin{equation*}
u_{r}(x,t_c)=\Pi _{V_s} [A\sech( x/L ) ].
\end{equation*}
We   compute the value of the first spatiotemporal local maximum  of $u_r(x,t)$ for a range of $(A,L)$ to investigate group dynamics under different lengthscales  $L$ and amplitudes $A$. In the left pane of~\cref{fig:focusing} we display the value of this map  divided by the initial amplitude $A$. This wave group amplification factor is  $1$ for defocusing groups and greater than $1$ for groups that undergo nonlinear focusing. Importantly,  by this  map we can partition the space (right pane of~\cref{fig:focusing}) and subscribe the region where   $u_\text{max}(A,L)/A > 1$ as the critical conditions that cause  rare transitions, since for these states  $\lambda_T > 0$   and $u_r(x,t)$  (conditional on the background state $(A,L)$) is a focusing wave packet.  This gives the set
\begin{equation}
R_e = \biggl\{ (A,L)   \biggm\vert \frac{1}{A} u_\text{max}(A,L) >1 \biggr\},
\end{equation}
which   parameterizes the instability region in terms of just  two  parameters~$(A,L)$.
\begin{figure}[htb]
\centering
\includegraphics[width=0.8\textwidth]{ALMapPODpartition.png}
\caption{\textit{Left: }The wave group amplification factor $\dfrac{1}{A}
u_\text{max}(A,L)$ computed using  the reduced-order dynamical system for
the rare events  $u_r(x,t)$. \textit{Right:} Partition of $(A,L)$ plane into
focusing and defocusing regions according to the reduced order model.}
\label{fig:focusing}
\end{figure}

\subsection{Conditional statistics of rare events}
\label{sec:extremeResponse}

To determine   $\pdf_{\abs{u}}(q \mid  \lVert u \rVert > \zeta,u_b \in R_e) $ we express it in terms of the wavegroup parameters, i.e. in the form    $\pdf_{\abs{u}}(q \mid\lVert u \rVert > \zeta , \, A, L\in R_e)$, where $A$ and $L$ denote the parameters  of  an arbitrary wavegroup formed in the stochastic background $u_b$ and lead to a rare event defined by $\zeta = \sigma_0$, the standard deviation of the initial spectrum. In the rare event regime we approximate $u\simeq u_r$. The estimation of this quantity involves simulations of $u_r$ through the reduced-order model. We perform the simulation until the wavegroup  relaxes back to the background state, necessarily a short time simulation. We also sample spatial points for which the response has important magnitude, in practice this means that we consider points $x \in [-2L,2L]$. This is done for all unstable wavegroup parameters $(A,L)$ for which the probability of occurrence $\pdf_{AL}(A,L)$ within the background dynamics is finite. A visual demonstration of this procedure is displayed in~\cref{fig:extGroupVis}.
\begin{figure}[tb]
        \centering
        \includegraphics[width=0.8\textwidth]{patches_example.pdf}
        \caption{Procedure for computing the  conditional statistics for rare events.
 First we  determine the probability  of     wavegroups, $\pdf_{AL}(A,L)$, in $R_e$ for a given
spectrum. Then we  compute the pdf for the quantity of interest for each focusing group, using the reduced order model.}
        \label{fig:extGroupVis}
\end{figure}

We explicitly compute $\pdf_{AL}(A,L)$ by generating many random fields according
to~\eqref{eq:gaussSpec}. This random field represents typical realizations of the background state $u_b$, where the dominant effect is linear dispersion, which continuously mixes the phases between harmonics. We note that we ignore  the evolution of the spectrum due to weak nonlinearities (typically the spectrum tends to broaden~\cite{dysthe2003, xiao13})  and we use the initial spectrum as the  spectrum of the background $u_b$. This  simplification, which could be avoided by sampling the spectrum of the steady state using a relatively short simulation, does not cause any serious discrepancies to our final results. However, for more complex cases, such as two-dimensional nonlinear waves, where the spectrum is continuously evolving, the   transient character of the background statistics may have to be taken into consideration. Our framework can support such situation without any modifications.

For each random field realization, we apply a group detection procedure (described in detail in~\cite{cousinsSapsis2015_JFM}), which returns a set of the groups in the field along with the amplitude and length scale of each group. This is a fast computation since it only requires generating and analyzing random realizations out of a given spectrum. After generating many random fields and computing the groups, we have a set of samples of $(A,L)$; subsequently we estimate the joint density $\pdf_{AL}(A,L)$, scaling the amplitudes in the resulting pdf by their corresponding lengthscale.

Example group densities for two different spectra are displayed in~\cref{fig:groupDensity}.  In each figure we overlay the $R_e$ boundary.   Notice that in the case displayed on the right panel, the spectral width is larger, meaning that groups are narrower.  The effect of this change is to \emph{reduce} the number of focusing groups, and thus reduce the number of extreme waves.  This is consistent with the Benjamin-Feir Index ($2\sqrt{2} \epsilon/\sigma$), used in the nonlinear water-waves community to evaluate the  likelihood for  rare events,  which is reduced by increasing $\sigma$.
\begin{figure}[tb]
        \centering
        \includegraphics[width=\textwidth]{probAL.png}
        \caption{Group density $\pdf_{AL}(A,L)$ for Gaussian spectra with
$\epsilon=0.05$, $\sigma = 0.1$ (left) and $\epsilon=0.05$, $\sigma=0.2$
(right). The red curve is the $R_e$ boundary.}
        \label{fig:groupDensity}
\end{figure}

\subsection{Conditional statistics of background dynamics}\label{sec:mnls_regularcomponenet}

The next ingredient is to  estimate   the  conditional statistics for the background dynamics $\pdf_{\abs{u}}(q \mid u=u_{b})$. We will have only background components if and only if the occurring wavegroups belong in $R_e^c$. As a result, we have that  $\pdf_{\abs{u}}(q \mid u=u_{b})=\pdf_{\abs{u}}(q \mid A, L\in R_e^c).$ Moreover, in this regime we have predominantly linear dynamics and therefore the statistics of the background state can be approximated by a Gaussian distribution. More accurate representations that take into account weak nonlinearities~\cite{Tayfun1980} may also  be utilized to improve the accuracy of the main core of the distribution.

For linear waves $\eta(t)$  with a  Gaussian distribution their envelope is Rayleigh distributed. Therefore, the conditional distribution for the quantity of interest in the background state is given by
\begin{equation}\label{eq:stableregimedistribution}
 \pdf_{\abs{u}}(q\mid R_e^{\,c}) = \frac{q}{\sigma_s^2} e^{-q^2/ 2 \sigma_s^2},
\end{equation}
where $\sigma_s^2$ is the variance   of   \emph{non-focusing} wave groups, computed by taking  the variance of wave  groups in $R_e^{c}$. This is estimated   along side the computation described in~\cref{sec:mnls_likelihood} for determining $\pdf_{AL}(A,L)$, by   using the scale selection algorithm  to   identify the  wave  groups in $R_e^{c}$
 and then by taking the  variance of this set of wave packets.

\subsection{Probability for rare events}\label{sec:mnls_likelihood}

The final component  is the computation of the total probability to have a rare event in an arbitrary spatial location. We assume that the system is ergodic, and we use~\eqref{eq:rare_trans} in  the following spatial form
\begin{equation}
\mathbb{P}(\left\Vert u \right\Vert>\zeta,u_{b}\in
R_e)=\frac{1}{X_{D}}\int\limits_{u } \! X_{e}\left(u\right)\pdf_{u_b}(u)\,du
,
\end{equation}
where $X_{e}\left(u\right)$ is the spatial extend of each wavegroup associated with nonlinear focusing, and $X_D$ is the spatial size of the domain. Taking into account that the pdf $\pdf_{AL}$ already has been weighted with respect to the spatial extend of the wavegroups the last formula takes the form
\begin{align}
        \mathbb{P}_{r}=\mathbb{P}(\left\Vert u \right\Vert>\zeta,u_{b}\in
R_e)   = \iint\limits_{R_e}\pdf_{AL}(A, L\mid R_e) \, dA dL.
        \label{eq:decompExtremeGroup}
\end{align}
This probability, of course, depends on the particular choice of the spectral properties of the field. Increasing  the spectral width, for example, will shift the distribution of $\rho_{AL}$ to the left, reducing in this way the total likelihood for the occurrence of a rare event (\cref{fig:groupDensity}).

\subsection{Synthesis and comparison to Monte-Carlo simulations}
\label{sec:envelopeSynthesis}

We have now determined all the components required to compute the heavy tailed    distribution for $\abs{u}$ according to  the decomposition-synthesis procedure. We have the following final result, by the total probability  probability argument:
\begin{equation}\label{eq:mnls_decompfinalresult}
        \pdf_{|u|}(q) = \pdf_{|u|}(q \mid R_e^{c}) ( 1-\mathbb{P}_r ) + \prob_{r}\iint\limits_{R_e}
\pdf_{|u|}(q \mid A,L) \pdf_{AL}(A,L)\, dAdL.
\end{equation}

In~\cref{fig:exDecompEnv} we demonstrate the decomposition-synthesis method on four different spectra alongside the `true' density  from Monte-Carlo simulations. To compute the Monte-Carlo statistics we use $10^{4}$ realizations of the MNLS equation on a spatial domain $256 \pi$ and a total sampled duration of $1000$ time units, discarding the first $500$ time units for cases with  $\sigma=0.10$ case and $500$ time units for the cases with $\sigma=0.20$ case. Overall, we see that the  decomposition-synthesis  approach gives good approximations to the true density.   For the most heavy-tailed cases, we observe that large waves occur many orders of magnitude more frequently than predicted by linear dynamics (Rayleigh distribution).  Our  method compares favorably with the Monte-Carlo results   at a fraction of the computational cost.

We note that the computational savings of  the decomposition-synthesis  method is due to several important advantages over  the direct brute force Monte-Carlo method. First and very importantly, we do not have to wait for extreme waves to emerge from the background wave field; moreover, we do not have to wait until a sufficient number of extreme waves occur in order to obtain reliable tail statistics. This is because we simulate  extreme wave groups  directly by carefully choosing the initial/background conditions that trigger them. The fact that we use  a low dimensional  reduced order model to evolve  extreme wave groups, makes applying our estimation procedure very inexpensive  for computing  the  extreme components.  And finally, the distribution of the background dynamics was obtained analytically.

We also emphasize that  once we performed the decomposition-synthesis procedure, computing the distribution for different background spectra is  extremely  cheap, since the only quantity that needs to be recomputed is the distribution of various groups $\pdf_{AL}(A,L)$ and $\sigma_s^2$, which is easy and fast to determine since it does not require the simulation of the original SPDE or the solution of the reduced-order model.
\begin{figure}[htb]
        \centering
                \includegraphics{envDecomp_POD_groupV_nModes14_updatepatchalg_Fine_eps0pt025_sig0pt200.pdf}
                \includegraphics{envDecomp_POD_groupV_nModes14_updatepatchalg_Fine_eps0pt025_sig0pt100.pdf}
                \includegraphics{envDecomp_POD_groupV_nModes14_updatepatchalg_Fine_eps0pt050_sig0pt200.pdf}
                \includegraphics{envDecomp_POD_groupV_nModes14_updatepatchalg_Fine_eps0pt050_sig0pt100.pdf}
        \caption{Probability distribution function approximation by the decomposition-synthesis method for the MNLS equation~\eqref{eq:MNLS}.  In each case, the blue curve denotes Monte Carlo simulation of the MNLS and the red line represents the approximation by the decomposition-synthesis procedure.  The vertical dashed line denotes $4\epsilon$ (i.e. 4 standard deviations) and $x$-axis plotted to $8\epsilon$.  We also display the Rayleigh distribution   of $|u|$  (green dashed line). Ordered in increasing BFI regimes: top left $\text{BFI} = 0.35 $, top right $\text{BFI} = 0.71$, bottom left $\text{BFI} = 0.71 $, bottom right $\text{BFI} = 1.41 $.}
        \label{fig:exDecompEnv}
\end{figure}

\section{Nonlinear  system of coupled  oscillators}\label{sec:osc_eg}
We first illustrate the decomposition-synthesis method on   a  two-degree-of-freedom system composed of a linear oscillator  coupled quadratically to a nonlinear oscillator with cubic stiffness,   both under the action of white-noise excitation. The system is given by
\begin{equation}\label{eq:osc_system}
        \begin{aligned}
        &\dd{x}{t} +  c_x \d{x}{t} + \bigl(k_x + a_x y(t) \bigr)x  + s x^3  = \sigma_{f_{x}} \dot W_x , \\
        &\dd{y}{t} + c_y \d{y}{t} + k_y y  = \sigma_{f_{y}} \dot W_y ,
        \end{aligned}
\end{equation}
where $x,y$ are the two state variables, $c_x,c_y > 0$ are the damping parameters, $k_x,k_y>0$ are   the linear stiffness parameters, and $s\neq0$ is the nonlinearity parameter. The constants $\sigma_{f_x},\sigma_{f_y} > 0$ are the noise intensities, and $W_x,W_y$ are two independent scalar Wiener processes. With $a_x\neq 0$, the system (\ref{eq:osc_system}) possesses regimes with  intermittent instabilities  that lead to fat-tailed equilibrium pdfs in the variable $x$, whereas the variable $y$ converges to a  Gaussian   pdf (see~\cref{fig:oscil_real}). This system is introduced such that  $y$  can trigger an intermittent response in $x$  through a large deviation from its mean value.

System~\eqref{eq:osc_system}      is prototypical of the action of intermittency in      more complex systems   where similar interactions are at play between system modes, e.g. internal instabilities   associated with transfers of energy between modes in turbulence, buckling of beams and plates under stochastic excitation, ship rolling under parametric stochastic resonance, just to mention a few. Although  system~\eqref{eq:osc_system} is   low-dimensional, quantifying the stationary pdf for the variable $x$ is  challenging. For one, the  finite time-correlated parametric excitation term  due to  $y$   on the variable $x$ (which is the mechanism behind  transient instabilities), makes application of the Fokker-Planck equation, unpractical and computationally prohibitive (the resulting equations become high-dimensional). Moreover,  there is a  nonlinear restoring force ($s\neq0$) that has to be taken into consideration; this nonlinear  term   has a significant impact on the pdf of $x$.

Besides its prototypical character, part of the motivation behind the current example is to apply our proposed method to a system where the mechanism behind rare events is transparent. As a result, we do  not need to apply the  order-reduction schemes   described in~\cref{sec:rare_events} for the rare event component (this will be done in the next application).  This will allow for the demonstration of the   decomposition-synthesis procedure   in a   manner that can be adapted to different problems that share similar characteristics.

\begin{figure}[htb]
        \centering
        \includegraphics{realization}
        \caption{Realization of the intermittent variable $x(t)$ (top, left) with large amplitude   bursts and a fat-tailed equilibrium pdf (top, right), alongside the exciting signal  $y(t)$ (bottom, left) with a Gaussian response pdf (bottom, right).}
        \label{fig:oscil_real}
\end{figure}

\subsection{Decomposition and instability region $R_e$}
In this case the decomposition step is trivial since all the rare events involve the variable $x$, thus $V_s=\spn\{ (1,0)^{T} \}$. In particular, intermittency in the variable $x$ is a consequence of $\kappa(t) \coloneqq k_x + a_x y(t)$ switching signs from positive to negative values, during which $x(t)$ transitions from its regular  response to a domain where the likelihood of an instability is high or guaranteed.  This switching in $\kappa(t)$ is the triggering mechanism behind localized instabilities in the variable $x(t)$. Therefore, we define the instability region\begin{equation}
R_{e}= \{ y  \mid  \kappa =k_x + a_x y<0 \}.
\end{equation}

\subsection{Conditional statistics of rare events}\label{subsec:osc_domattraction}

The intensity of the rare events depends both of the average magnitude of $\kappa$ when this becomes negative, but also on the duration of the   downcrossing. Therefore, we will choose the rare event description presented in section \ref{red2}, where the background state is taken into account for the evolution of the rare events. We choose to parameterize the instability region $R_e$ with the   average amplitude $\alpha$ and average duration $\xi$ of each downcrossing event of $\kappa$, which expresses the background state information.

We compute a few realization of the background state \begin{equation}
        \dd{y}{t} + c_y \d{y}{t} + k_y y  = \sigma_{f_y} \dot W_y,
\end{equation}
and for each realization  we identify all the regions where $\kappa(t)  = k_x + a_x y(t) < 0$. For each of  these events, that is, a zero downcrossing followed by a subsequent zero upcrossing, we take the duration of time that $\kappa<0$  as  $\xi$ and assign a characteristic amplitude $\alpha$ by taking the  mean of $\kappa$ over the duration that $\kappa<0$ . Performing this procedure for all realizations gives us a set of samples of $(\alpha,\xi)$. Since we are interested in the probability at a given temporal location  we have an event $\kappa < 0$ with duration $\xi$ and amplitude $\alpha$, we scale the amplitudes in the resulting histogram  by their corresponding durations $\xi$ in order to correctly weight the samples. This gives $\pdf_{\alpha\xi}(\alpha,\xi \mid R_e)$: the probability of finding a rare event characterized by an excitation duration $\xi$ and amplitude $\alpha$. We display this pdf in~\cref{fig:oscil_PAL}.
\begin{figure}[tb]
        \centering
        \includegraphics{PAL.png}
        \caption{Joint density $\pdf_{\alpha \xi}(\alpha,\xi)$ (left), alongside the marginalized densities for   amplitude $\alpha$ (top right) and length scale $\xi$ (bottom right) for $c_y=0.5,$ $k_y=1$, $\sigma_{f_y}=2.5$.}
        \label{fig:oscil_PAL}
\end{figure}

Next, we proceed with the computation of the conditional statistics for rare events, $\pdf_x(x \mid \alpha,\xi)$.
 This will be based on a few simulation of the equation determining $x$ taking into account the background state (expressed through $(\alpha,\xi)$). More specifically, for a given $\alpha_i,\xi_i$ we evolve the rare event dynamics
 according to
\begin{equation}\label{eq:osc_extremerun}
        \dd{x}{t} +  c_x \d{x}{t} + \tilde \kappa(t;\alpha_i,\xi_i) x  + s x^3
= \sigma_{f_x} \dot W_x,
\end{equation}
where $\tilde\kappa(t;\alpha_i,\xi_i)$ is a  local representation of the function  $\kappa(t) = k_x + a_x y(t)$  in  the critical regime with
the same mean amplitude $\alpha_i$ (over $0\leq t\leq \xi_i$)  and duration
$\xi_i$. For simplicity, we choose the following function to approximate $\kappa$:
\begin{equation}
        \tilde\kappa(t,\alpha_i,\xi_i) =
        \begin{dcases}
                \frac{\pi\alpha_i}{2} \sin( \pi t /\xi_i), \quad &\text{for }
0 \le t \le \xi_i,\\
                k_x,\quad &\text{otherwise}.
        \end{dcases}
\end{equation}
Initial conditions   for $x$ are chosen from the background state,  which we will describe later. In this regime,   we are only interested
in extreme responses for which $\left\vert x(t) \right\vert>\zeta=\sigma_{b}$,
where $\sigma_{b}^2=\overbar{({x}|_{\kappa>0})^2}$ is the variance of the background stochastic
state. To perform the computation we simulate an instability
according to (\ref{eq:osc_extremerun})   for
a length of time that is long enough to ensure that   $x(t)$ returns back to the
 background regime, where $\left\vert x(t) \right\vert<\zeta$,  and  sample only the points before the response relaxes
back, i.e. the points for which $\left\vert x(t) \right\vert>\zeta$.
Following
this procedure we obtain the quantity $\pdf_x(x \mid \alpha,\xi)$.

\subsection{Conditional statistics of background dynamics}
\label{sec:osc_backgroundattractor}
Next we determine the statistical characteristics of the background attractor
$\pdf_x(x \mid R_e^{c})$. In this regime, by definition,  we have no rare
events, and we  choose to derive the statistical distribution    analytically.
This quantity could also be obtained more directly  through  the conditional
Monte-Carlo approach described in~\cref{sec:gausconditional}  by constraining
the sampling to the regime $\kappa=k_x + a_x y(t) > 0$.

Since  $\kappa > 0$ in  $R_e^{c}$, we choose to approximate $\kappa(t)$  by its conditionally average value  $\overbar{\kappa}|_{\kappa>0}$,
\begin{equation}\label{eq:osc_x_aved}
         \dd{x}{t} +  c_x \d{x}{t} + \overbar{\kappa}|_{\kappa>0} x  + s
x^3  = \sigma_{f_{x}} \dot W_x.
\end{equation}
The steady state pdf can be found
by   solution of the associated  FPK equation:
\begin{equation}
        \pdf_{x\dot x}(q,\dot q \mid R_e^{c}) = C \exp\biggl( - \frac{2c_x}{\sigma_{f_x}^2}
\mathbb H(q,\dot q) \biggr), \quad\text{where } \;\;  \mathbb H(q,\dot q) = \frac{1}{2}\dot
q^2 + \frac{1}{2} \overbar{\kappa}|_{\kappa>0} q^2 + \frac{1}{4}s q^4,
\end{equation}
and $C$ is a normalization constant. Marginalizing out $\dot x$ in the equation
above, gives the  conditional pdf for $x(t)$ in $R_e^{c}$:
\begin{equation}\label{eq:osc_stablepdf}
        \pdf_x(q \mid R_e^{c}) = C \exp\biggl[ - \frac{2c_x}{\sigma_{f_x}^2}
\biggl(\frac{1}{2} \overbar{\kappa}|_{\kappa>0} q^2 + \frac{1}{4}s q^4\biggr)
\biggr],
\end{equation}
where $C$ is the normalization constant.

To determine $\overbar{\kappa}|_{\kappa>0}$ we utilize the steady state pdf
of $y(t)$  which is    Gaussian distributed:
\begin{equation}\label{eq:yresponsepdf}
        \pdf_y(q) =  \sqrt{\frac{c_y k_y }{\pi \sigma_{f_y}^2}} \exp\biggl( -
\frac{c_y k_y}{\sigma_{f_y}^2} q^2\biggr).
\end{equation}
with variance $\sigma_y^2 = \sigma_{f_y}^2/(2 c_y  k_y)$. Therefore, we obtain
\begin{equation}
        \overbar{\kappa}|_{\kappa>0} = k_x +   a_x\sigma \frac{\phi\bigl(-\frac{k_x}{
a_x\sigma_y} \bigr)}{1-\Phi\bigl(-\frac{k_x}{ a_x\sigma_y}\bigr)},
\end{equation}
where $\phi(\,\cdot\,)$ is the normal probability distribution function and
$\Phi(\,\cdot\,)$ the normal cumulative distribution function.

\subsection{Probability for rare events}\label{sec:osc_liklihood}

We now proceed to compute the probability of rare events due to instabilities in $R_e$. This can be found by analyzing the duration of the rare transitions using~\eqref{eq:osc_extremerun}. Denoting by  $T_e(\alpha_i,\xi_i)$ the duration of the rare event corresponding to an instability in $R_e$ of average magnitude $\alpha_i$ and duration $\xi_i$, according to~\eqref{eq:rare_trans},  we have that
\begin{equation}
\label{prob_dur_osc}
\mathbb{P}(x>\zeta,y\in R_e)=\frac{1}{T} \int T_{e}\left(\alpha, \xi \right) \pdf_{\alpha \xi}(\alpha, \xi \mid \kappa < 0) \mathbb{P} (\kappa<0)\, d\alpha d\xi,
\end{equation}
where
\begin{equation}
        \prob(\kappa<0) = \prob(y < -k_x/a_x) = \Phi \biggl( -\frac{k_x}{ a_x\sigma}\biggr).
\end{equation}
The probability~\eqref{prob_dur_osc} can  be approximated using the analytical argument presented~\cite{mohamad2015}, where the typical temporal duration of
the growth and relaxation phase of a rare event are examined. More specifically,   consider a single  representative extreme response with an average growth $\overbar{\Lambda}_+$ and decay rate $\overbar{\Lambda}_-$. During the growth phase we have that
\begin{equation}\label{eq:osc_typicalgrowth}
        u_p = u_0 e^{\overbar{\Lambda}_+ T_{\kappa<0}},
\end{equation}
where $T_{\kappa<0}$ is the temporal duration of the growth event and $u_p$ is the peak value of the response and $u_0$ an arbitrary initial condition. Similarly, over the decay phase we have
\begin{equation}\label{eq:osc_typicaldecay}
        u_0 = u_p e^{-\overbar{\Lambda} - T_\text{decay}},
\end{equation}
and by connecting these~\cref{eq:osc_typicalgrowth,eq:osc_typicaldecay},
\begin{equation}
        \frac{T_\text{decay}}{T_{\kappa<0}} = \frac{\overbar{\Lambda}_+}{\overbar{\Lambda}_-}.
\end{equation}
The  average duration of an extreme transition  is   then given by  $T_{e} = (1 + \overbar{\Lambda}_+/{\overbar{\Lambda}_-}) T_{\kappa<0}$. Therefore, by dividing over  the total time  duration $T$ we have
\begin{equation}
        \mathbb{P}(x>\zeta,y\in
R_e)\simeq\  \biggl(1 + \frac{\overbar{\Lambda}_+}{\overbar{\Lambda}_-}\biggr) \prob(\kappa<0).
\end{equation}To leading order (neglecting the nonlinear restoring term) $\Lambda_+  = \sqrt{- \kappa|_{\kappa<0}},$ since dissipation and additive forcing have a small role during the growth phase of an extreme response. Therefore, the average growth exponent is   $\overbar{\Lambda}_+ = \expec\bigl(\sqrt{-\kappa|_{\kappa<0}}\bigr),$ which is straightforward  to compute using $\pdf_{\alpha\xi}$. Similarly, we find that during the decay phase is $\Lambda_- = c/2$, and thus  $\overbar{\Lambda}_- = c/2$.

\subsection{Synthesis and comparison to Monte-Carlo simulations}

We have now determined  all the components required to construct the distribution of the response according to the decomposition-synthesis method. We synthesize  the results of the previous sections by the total probability law:
\begin{equation}\label{eq:osc_decompfinalresult}
\pdf_x(q) = \pdf_x(q \mid R_e^{c})\left( 1-\mathbb{P}_r \right) + \prob_{r}\iint \pdf_x(q \mid \alpha,\xi) \pdf_{\alpha \xi}(\alpha,\xi)\, d\alpha d\xi,
\end{equation}
where $\mathbb{P}_r=\mathbb{P}(x>\zeta,y\in R_e)$.

In~\cref{fig:osc_pdfresults} we compare the decomposition-synthesis results for the equilibrium pdf of $x(t)$ alongside the `true' density from Monte Carlo simulations,  for increasing values of the cubic stiffness parameter for a hardening spring $s>0$. Monte-Carlo results are computed using $5000$ realizations of (\ref{eq:osc_system}), integrated using the Euler–Maruyama method with time step $\Delta t = 0.002$ to $T = 500$ time units, discarding the first $t = 60$ data  to ensure a statistical steady state. Overall, we have very good quantitative agreement for both the tails and the core of the distribution  between the decomposition-synthesis approach  and the   Monte Carlo results.  For larger values of the nonlinearity parameter $s$, tail values are suppressed due to larger restoring forces. We see that this effect is also accurately captured for increasing $s$ values in the decomposition-synthesis approach. Here we emphasize the non-uniform decay character of
the tails which can still be captured very accurately. In particular we have favorable comparison for all three   qualitatively different  regimes: the core of the   distribution, the  exponential like heavy-tails at extreme values,  and the subsequent sub-exponential decay at very extreme values of distribution (where nonlinearity is very important).
\begin{figure}[bt]
        \centering
        \includegraphics{{{response-pdf_xParams-c=1.00_m=4.00_s=0.40_sigx=0.75_yParams-b=0.50_k=1.00_sigy=2.50}}}
        \includegraphics{{{response-pdf_xParams-c=1.00_m=4.00_s=0.60_sigx=0.75_yParams-b=0.50_k=1.00_sigy=2.50}}}
        \includegraphics{{{response-pdf_xParams-c=1.00_m=4.00_s=0.80_sigx=0.75_yParams-b=0.50_k=1.00_sigy=2.50}}}
        \includegraphics{{{response-pdf_xParams-c=1.00_m=4.00_s=1.00_sigx=0.75_yParams-b=0.50_k=1.00_sigy=2.50}}}
        \caption{Probability distribution function approximation with the decomposition-synthesis method for the variable $x(t)$ for the nonlinear system of coupled oscillators~\eqref{eq:osc_system}. The blue line denotes the Monte Carlo simulation and the red line denotes
the approximation  by the decomposition-synthesis procedure.  The vertical dashed lines denotes $4$ standard deviations of $x(t)$. Parameters  are:  $c_x =1,$ $c_y=0.5,$ $k_x = 4,$ $k_y=1$, $\sigma_x = 0.75,$ $\sigma_y= 2.5$.}
        \label{fig:osc_pdfresults}
\end{figure}

We point out that computing the response pdf via the decomposition-synthesis for different background parameters $c_x,s,\sigma_x$ parameters is extremely cheap since   $\pdf_{\alpha \xi}(\alpha, \xi),$ which involves the rare event dynamics, remains fixed.  Also, if we are interested in the response pdf for different $k_x,a_x$,  this also has minimal computational cost because  we can store rare event realizations of $y$ and use them to determine $\pdf_{\alpha \xi}(\alpha, \xi)$, as required. Moreover, the computation to determine $\pdf_{\alpha \xi}(\alpha, \xi)$  is easy and fast, since we do not need a large number of realizations  to give good results when used in (\ref{eq:osc_decompfinalresult}) and  only requires the simulation of the background variable $y(t)$.

redpajama_set_name: RedPajamaArXiv


text:
Drinks Trade
Negroni with love from Melbourne
A new player in the premium premixed cocktail market, Vero Cocktails, with their flagship Vero Negroni are introducing the tradition of the Italian Aperitivo into Australian homes.
"An aperitivo is best enjoyed shared, so rather than the usual single serve pre-made cocktails we decided to create a 700ml bottle that serves nine glasses – perfect for dinner parties and gatherings," James says.
Now, the Aussie-made, Italian-inspired Vero Negroni is being enjoyed across the country, available in over 50 stores across Victoria and NSW and online at Dan Murphy's and the Vero Cocktails website.
Vero Cocktails' Negroni with love from Melbourne
Sign up to Vero
December 2021 Competition
Under the Liquor Control Reform Act 1998 it is an offence: to supply alcohol to a person under the age of 18 years (penalty exceeds $19,000). For a person under the age of 18 years to purchase or receive liquor (penalty exceeds $800).
Liquor license number 33775866
@vero_cocktailsInstagram
© 2022, Vero Cocktails

redpajama_set_name: RedPajamaCommonCrawl


The HTC Vive. Image source: HTC.
The concept of virtual reality headsets is hardly a new one; they've popped up in various forms over the years, but have always ultimately failed to gain traction. Over the last several years, the concept of VR has, once again, become quite popular. Whether virtual reality sticks this time or ultimately winds up being another fad is anyone's guess.
However, for those interested in companies (and their stocks) that could benefit if this round of virtual reality hits it big, here are three that spring to mind.
Right now there are two major VR platforms that seem to have people's attention: Rift from Facebook-owned Oculus and the Vive from struggling smartphone maker HTC .
If Oculus were a stand-alone publicly traded company, it might be an interesting way to play a potential VR boom. However, it's owned by Facebook, already a large, successful company trading at a lofty valuation. Even if the Oculus Rift is successful (by some reasonable definition of success), it's not clear that this will mean much to Facebook's top and bottom lines.
HTC, on the other hand, seems to be a much more interesting way to play a potential boom in VR headset sales. In addition to what seem to be quite a few voices in the tech press claiming that the HTC Vive is a superior platform to the Oculus Rift (though there are some that say the opposite), HTC's business could be positively, and meaningfully, impacted if VR takes off since its market capitalization is roughly equivalent to $2.53 billion(remember that Facebook bought Oculus for $2 billion).
That said, it's important to keep in mind that betting on the proverbial little guy is inherently riskier; Oculus, backed by Facebook, is likely to be able to outspend HTC on hardware development, software ecosystem development, and marketing.
It's almost a common-sense notion that doing virtual reality properly is going to require immense amounts of graphics processing power. Graphics specialist NVIDIA is the largest vendor of high-performance graphics processors for personal computers.
Since headsets like the HTC Vive and the Oculus Rift require very high performance graphics cards in order to play virtual reality games at acceptable performance levels, it would seem that this plays right into the graphics specialist's strengths.
The company seems to be taking virtual reality quite seriously, too. Its recently released Pascal graphics architecture implements specialized circuitry that helps to accelerate virtual reality workloads. Additionally, NVIDIA announced VR Works, which the company describes as a "comprehensive suite of APIs, libraries, and engines that enable application and headset developers to create amazing virtual reality experiences."
Should virtual reality ultimately take off, it could lead to both unit growth in the company's core graphics processor business as well as average selling price growth (as customers buy more powerful graphics processors in order to take advantage of virtual reality applications).
Looking to get in on the action in virtual reality, console maker Sony is planning to launch its own VR headset, known as the PlayStation VR, later this year. Sony has the advantage of being able to market its new headset to a fairly large PlayStation 4 installed base. It also controls both the hardware and software ecosystem around PlayStation, which could allow the company, and developers for the PlayStation 4 platform, to deliver relatively good, polished experiences.
If PlayStation VR is well received, Sony could benefit in three ways. First, the appeal of PlayStation VR could draw gamers to its platform over competing game consoles, driving console unit sales.
Next, Sony could benefit from game developers potentially releasing "special" virtual reality-capable titles (either remixes of pre-existing titles or totally new titles), driving PlayStation game sales, ultimately resulting in increased royalty revenue to Sony from incremental game sales.
And, finally, assuming that Sony is selling its PlayStation VR headsets for a positive gross profit, it stands benefit directly from sales of the PlayStation VR headsets.
The article 3 Stocks That Could Benefit From a Virtual Reality Boom originally appeared on Fool.com.
Ashraf Eassa has no position in any stocks mentioned. The Motley Fool owns shares of and recommends Facebook and NVIDIA. Try any of our Foolish newsletter services free for 30 days. We Fools may not all hold the same opinions, but we all believe that considering a diverse range of insights makes us better investors. The Motley Fool has a disclosure policy.


John Alan Wreford, 65 of Berrien Springs went to be with his Lord and Savior on Monday, July 30, 2018 at his home in the care of his family.
John was born Saturday, February 14, 1953 in Gatooma, Rhodesia, now known as Kadoma, Zimbabwe, a British settlement in South Africa, the son of the late Frank Alan and Mary (Porter) Wreford. On August 5, 1989 he married Tammy Sue Herr. John was a member of Berrien Center Bible Church and a media tech at Life Action Ministries.
Survivors include his wife of 28 years, Tammy; a son: Mark Alan Wreford; a brother: Peter Wreford; sisters: Margaret (Nigel) Tayor, Ann (Jeff) Kuper, Veronica (Jonathan) Johnson-Butcher and Patricia ((Mick) Stanton.
John was also preceded in death by his first wife, Helen (Nicoll) Wreford in 1982 and a sister, Elizabeth VanBuren.
Arrangements have been entrusted to Pike Funeral and Cremation Services, 9191 Red Arrow Hwy, Bridgman.


The year was 1966 and I was a newly assigned machine gunner in an infantry unit that was protecting a high mountain pass that was between Da Nang and Phu Bai, Vietnam. The pass is called Hai Van Pass, which means Place of the Clouds. It is located 30 km north of Da Nang. The Marine Corps name of the hill that my unit was on was Hill 724. It was a dangerous and treacherous area that was of high strategic value. There was only one narrow highway that wound up steep cliffs (Highway 1) that reached the peak at Hai Van. At the base of the hill was a small village called Namo. My story starts in that village.
On a patrol through the village, my squad noticed a small wooden cage outside of a hooch. We were new to Vietnam and unaware of the customs of the people. We were always alert for booby traps and we moved closer to the cage. Inside the cage was a small, brown puppy that brightened up when he saw us approach his cage. He began wagging his tail and wanting to be held. It was love at first sight for all of us. Don't tell anybody, but most Marines are real softies when it comes to animals. The owner of the hooch came out and I decided to offer to buy his puppy. I just didn't like his living conditions. The old villager started his trying to up the price and after much debate, he decided on a price. He then pulled out a sharp knife and proceeded to open the cage.
All of us were shocked that he pulled a knife out and we all raised our weapons to him. He looked very frightened and assured us that he meant us no harm. He explained in gestures that he was just going to prepare the dog for us. We were confused and then he made the motion of drawing his knife across his neck to tell us that he was going to kill the dog and dress it out for us to Chop-Chop. I can only tell you that the anger level went up 100 notches at that suggestion. We then realized that the dog was being raised to eat by the villager. I unbuttoned my flak jacket and reached in and took the dog from his hands. I threw the money on the ground and placed the pup next to my heart so I could carry him up the long torturous climb back to our base camp.
The first event that happened was all of us knew that it was going to be a challenge to get our superiors to allow us to keep him. It was a miserable walk back to our hill. Our fears were baseless because we had a wonderful Captain, named Capt. Silva, and he allowed us to keep the puppy as a mascot. We loved him before, but we really loved the Captain after that. Next we had the heated debate on what the dogs name should be. It was not an easy process since marines are very bull-headed and strong willed. After much serious and highly intelligent discussion we arrived at the unique name of Brown Dog. You have to know Marines to appreciate this.
Brown Dog was the darling of the hill and only one Marine hated him and that was OK, because we all hated him too. We all decided that if Brown Dog didn't like him that there had to be a damn good reason. Later he proved to be a coward in combat and was removed from our hill before he had an Accident. Brown Dog was very happy that he left. I lost a stripe because I caught this guy kicking at Brown Dog and I explained to him not to do that again. I explained a little too harsh and the 1st Sgt. explained to me the error of my ways.
Brown Dog had a ritual of his life on the hill. At night we were in bunkers staring out into a dark, fog filled jungle. We were issued Seismic Listening Devices which consisted of probes planted in the ground in front of our bunkers. We had a small console inside that had earphones. We could hear footsteps approaching or animals moving. We got pretty good with the device. Brown Dog would make the perimeter of our positions and visit every bunker to check on His marines. He was always a welcome visitor and he spent all of his non patrolling time in his daddys bunker, MINE! He also went on patrols with us and had an intense hatred for the Vietnamese. He would growl and really act up when he would see or smell one.
On April 1, 1967, we were dug in and the fog was pea soup thick that night. I was in the machine gun bunker and we were really spooked. About 3 am Brown dog shot up and went on full alert. I rang the field phone and informed the Command bunker that Brown Dog had alerted. Our Lt., ( Naval Academy IDIOT) advised me to trust my Seismic Device and not a damn dog!
My bunker was the forward bunker and the most vulnerable. I looked at my machine gun crew and whispered that I was going to disobey the Lt.s stupid order. They all nodded as I prepared the pop-up flare to shoot into the sky. I popped the flare and Lo and Behold we had Beaucoup Gooks in the wire! All hell broke loose and it was a very violent battle that night. It was up close and personal fighting and many people on both sides were killed and wounded. Brown Dog was hit by shrapnel but continued to fight the enemy. I saw him attack the leg of a NVA before he was zapped.
The attack failed and for what seemed like an eternity, we waited for the sun to come up. There were dead people in the wire, burning, and moaning out in pain in the darkness. I held Brown Dog in my arms and awaited the medivac helicopters that were coming to help our wounded. I really thought Brown Dog was going to die in my arms.
When the choppers came I handed my baby to the door gunner and asked him to get him medical help and told him that we would all be dead if it wasn't for Brown Dog. Marine to Marine, he gave me his word and I watched the dust-off. My Capt. pulled me aside and told me that he was glad that I never followed orders. He rubbed my head and told me that he was going to call to headquarters and tell them the story and he ASSURED me that Brown Dog was going to get the best of care.
Brown Dog had lost a lot of blood so I really didn't have much hope. The next day we received word that Brown Dog had been taken to Army Vets and they had saved his life. He had over 100 stitches, needed blood and antibiotics but he was going to live. He was the hero to all of us. He received a canine award for heroism and we promoted him to Cpl. I also got my stripe back much to the chagrin of the Lt.
Cpl. Brown Dog returned to a heros welcome as we had a full formation to welcome him home. I left the hill June 9th, 1967 and he was in very good hands with the Marines on the hill. I tried to take him home but that was impossible. I really hurt having to leave him on the hill, but I knew the Marines would take excellent care of him.
I have always thought it was so ironic that once he was going to be eaten by the Vietnamese and in turn he caused the death of so damn many of them!
All of my friends still have his picture and we all know that we would not be here today if it were not for a little, mixed breed dog named Brown Dog. He will be in our hearts until we die and a part of our souls forever.
When we assemble for our reunions, we always toast Brown Dog.
Rest in peace little Warrior and wait for Your Marines to join you. We will always be Semper Fidelis to you and your memory. A Salute and three cheers for the finest Marine on our hill!
The above story came through the internet to my Email today, 4-25-06. It is a tale that at least half of the warriors alive today can relate to and I am sure could tell us a similar experience with a stray dog who became Mascot of a tent, foxhole, squad of troops, platoon, company, or battalion. During my 24 years it happened more than once.
Tramp was a three legged brown dog and the Mascot of VMA(AW)-533 in Vietnam. The troops smuggled Tramp into MCAS Iwakuni, Japan when the squadron rotated out of Vietnam around 1970. I arrived as Commanding Officer of the squadron in 1971 and Tramp treated me like any other FNG boot that came aboard. I had to show myself worthy of being a Night Hawk in 533. Tramp was married to a four legged dog named Hubert who produced pups for Marine Unit Mascots. Good work by Tramp on one hind leg!
John McNamee provided these pictures of Tramp and Hubert taken in 1972 after 533 began flying combat missions from the Rose Garden in Nam Phong,Thailand. Tramp accepted me just as he did every FNG. In fact, when I separated from the squadron in June 72, Coy Cline and his administrators drafted a Mast citation which Tramp signed. I have kept that citation all these many years because it probably is the only such document in existence.
pleasant journey I enjoy to that lofty aerie where Tramp awaits us all.


Same day, same lens, messing with flash.
I may unconsciously be compensating for the down time the blog will suffer while I'm in Cuba next week... But what the heck. It's two for one today — although this one's a bit late.
Same production Fujinon 55-200mm zoom from earlier this morning, this time outdoors with a bare SB-900 speedlight on a light stand. I'm using a set of Cactus V5 triggers with the X-Pro1, something I purchased recently. They're fussy to get going and the on/off/mode switch feels like it was designed for lilliputians... But they do fire every time. I'm doing almost a day for night here, with the help of trees and an overcast sky, shooting 1/125 sec at f/8, ISO 200. No flash on the last two images.
Vancouver fashion photographer (and fellow canadian X Photographer) Kale JF commented on this morning's post via Twitter, saying how much he loved this lens in spite of his preference for primes. I agree. I had sworn off zooms and had made an oath with myself, god and the tooth fairy that should I ever need one it would never, ever be a variable aperture lens. Ever.
I'm such a liar. More to come.
Morning, sitting on the couch, playing with a zoom.


Show your gratitude with personal and professional style!
Did a friend do you a favor you could only hope to repay? Did you just get out of a meeting with an important connection at work? It's time to thank the people who gave you their time by taking the time to show them your appreciation, in the form of a thank-you note. Whether you're following up after an interview or letting a family member know you appreciate them, here are six steps to writing a thank-you note that will get your gratitude across.
Here is where you get personal, and where your thank-you note will truly shine. In order to make it clear that you really do appreciate whatever it is this person did for you, and that you aren't just going through the motions of thanking them, tell them how what they did made you feel. Do you feel happy? Appreciated? Supported? Encouraged? Loved? People like to feel that they've had a positive impact on you, so specificity is once again key here.
This is a must-do for letters in general, and especially with thank-you notes. This person did something for you, and you want good things to happen to them, too! Phrases like "I wish you the best of luck with…" or "I hope you are having a wonderful day" are good, and it's always nice to tell someone that you hope your paths will cross again in the future.
I recommend typing up your thank-you note first, so that spellcheck can help you out and you can more easily make draft edits. But after you have the thank-you note worded the way you want it, it's time to hand-write it! As communication becomes increasingly digital, a hand-written note stands out as a sign that you put thought, time, and effort into expressing your gratitude. Go to a store and buy a nice card and envelope. Use a pen to write the note, and take your time while writing to avoid mistakes!
If you're going to see the person again in the next couple days, feel free to hand the card to them in person. Otherwise, look up their address and send it in the mail. If you put the thank-you card in the mail the day after the person did the thing you're thanking them for, it will arrive just in time to give them a pleasant reminder!
Congratulations, you just wrote a terrific thank-you note!
Need an example of steps one through four in action? Here is a thank-you note for a hypothetical situation in which my mother took me to a restaurant and we talked about my future.
I wanted to thank you for treating me to lunch at One World the other day. It was great to get to spend some time with you; it's been too long, and your advice and support about my future as always left me feeling loved and supported! I hope you are having a wonderful day, and I look forward to seeing you again soon.
Did I respond too quickly? Ugh, auto-correct! Why is he taking so long to respond?
1. Over analyzing punctuation and phrasing.
I hear Ye Old Cafe has an awesome lunch menu!
2. Predicting a rejection and assuming the worst.
Great, he hates me! He thinks I'm a total weirdo and is probably mocking my very existence right now.
4. Those three little dots of dread.
Wow, your response time is impeccable... NOT! Just say what you need to say!
5. Assuming the worst when someone doesn't respond.
She has probably been attacked by zombies...and I was too slow to save her. Oh god! What if she's still being attacked? What do I do?
6. Feeling like a bother when you text first.
Hey! Oh dang, I'm probably annoying her...I take it back!
7. Trying to decipher the exact meaning of excess letters.
"Funnyyy!" OK, three y's, that means he thinks I'm actually funny? No, he's definitely mocking me.
8. Deciding on a context appropriate emoji.
9. Immediately regretting a text and wishing there was a way to undo it.
"LOL, you're sooooooo funny :)" OH GOD NO, that sounded way too eager! ABORT MISSION!
10. Wondering what you did wrong when someone is online but ignores your text.
Great, that status was probably about me...she could at least say it to my face!
11. The fear of misinterpreting a text.
He didn't use a smiley face...that means he's mad at me! Or is he just busy? Or maybe he just didn't see it...should I send it again?