---
author:
- Richard Kowar
- Otmar Scherzer
title: Photoacoustic Imaging Taking into Account Attenuation
---

Introduction
============

*Photoacoustic Imaging* is one of the recent hybrid imaging techniques, which attempts to visualize the distribution of the *electromagnetic absorption coefficient* inside a biological object.  In photoacoustic experiments, the medium is exposed to a short pulse of a relatively low frequency electromagnetic (EM) wave.  The exposed medium absorbs a fraction of the EM energy, heats up, and reacts with thermoelastic expansion.  This induces acoustic waves, which can be recorded outside the object and used to determine the electromagnetic absorption coefficient.  The combination of EM and ultrasound waves (which explains the usage of the term *hybrid*) allows one to combine high contrast in the EM absorption coefficient with high resolution of ultrasound.  The method has demonstrated great potential for biomedical applications, including functional brain imaging of animals [@WanPanKuXieSto03], soft-tissue characterization, and early stage cancer diagnostics [@KruKisMilRey00], as well as imaging of vasculature [@ZhaLauBea07].  For a general survey on biomedical applications see [@XuWan06].  In comparison with the X-Ray CT, photoacoustics is non-ionizing.  Its further advantage is that soft biological tissues display high contrasts in their ability to absorb frequency electromagnetic waves.  For instance, for radiation in the near infrared domain, as produced by a Nd:YAG laser, the absorption coefficient in human soft tissues varies in the range of $0.1/\textrm{cm}$–$0.5/\textrm{cm}$ [@ChePraWel90].  The contrast is also known to be high between healthy and cancerous cells, which makes photoacoustics a promising early cancer detection technique.  Another application arises in biology: Multispectral optoacoustic tomography technique is capable of high-resolution visualization of fluorescent proteins deep within highly light-scattering living organisms [@RazDisVinMaPerKoeNtz09].  In contrast, the current fluorescence microscopy techniques are limited to the depth of several hundred micrometers, due to intense light scattering. 

Different terms are often used to indicate different excitation sources: *Optoacoustics* refers to illumination in the visible light spectrum, *Photoacoustics* is associated with excitations in the visible and infrared range, and *Thermoacoustics* corresponds to excitations in the microwave or radio-frequency range.  In fact, the carrier frequency of the illuminating pulse is varying, which is usually not taken into account in mathematical modeling.  Since the corresponding mathematical models are equivalent, in the mathematics literature, the terms opto-, photo-, and thermoacoustics are used interchangeably.  In this article, we are addressing only the *photoacoustic tomographic technique* PAT (which is mathematically equivalent to the thermoacoustic tomography TAT). 

Various kinds of photoacoustic imaging techniques have been implemented.  One should distinguish between photoacoustic *microscopy* (PAM) and *tomography* (PAT).  In microscopy, the object is scanned pixel by pixel (or voxel by voxel).  The measured pressure data provides an image of the electromagnetic absorption coefficient [@ZhaMasStoWan06].  Tomography, on the other hand, measures pressure waves with detectors surrounding completely or partially the object.  Then the internal distribution of the absorption coefficients is reconstructed using mathematical inversion techniques (see the sections below). 

The common underlying mathematical equation of PAT is the *wave equation* for the pressure $$\label{eq:ex:wave3d}
\boxed{
\frac{1}{c_0^2} \frac{\partial^2 p}{\partial t^2}({{\vec{x}}},t ) - \nabla^2 p({{\vec{x}}},t) = \frac{d j}{d  t}(t)
      \left(  \frac{\mu_{\rm abs}({{\vec{x}}}) \beta ({{\vec{x}}}) J({{{\vec{x}}}})}{c_p({{\vec{x}}})} \right)\,,\;{{\vec{x}}}\in {{\mathbb{R}}}^3,\,t > 0 \;.}$$ Here $c_p$ denotes the specific heat capacity, $J$ is the spatial intensity distribution, $\mu_{\rm abs}$ denotes the absorption coefficient, $\beta$ denotes the thermal expansion coefficient and $c_0$ denotes the speed of sound, which is commonly assumed to be constant.  The assumption that there is no acoustic pressure before the object is illuminated at time $t  = 0$ is expressed by $$\label{eq:ex:ini3d}
\boxed{
    p({{\vec{x}}},t) = 0 \,, \qquad  {{\vec{x}}}\in {{\mathbb{R}}}^3, t < 0\;.}$$ In PAT, $j(t)$ approximates a pulse, and can be considered as a $\delta$-impulse $\delta(t)$. Introducing the shorthand notations $$\label{eq:u(x)}
\begin{aligned}
\rho({{\vec{x}}}) := \frac{\mu_{\rm abs}({{\vec{x}}}) \beta ({{\vec{x}}}) J({{\vec{x}}})}{c_p({{\vec{x}}})}\,,
\end{aligned}$$ one reduces (\[eq:ex:wave3d\]) and (\[eq:ex:ini3d\]) to $$\label{eq:ex:ivp3d}
\boxed{
\frac{1}{c_0^2} \frac{\partial^2 p}{\partial t^2} ({{\vec{x}}},t) - \nabla^2 p({{\vec{x}}},t)= 0\,,\quad  {{\vec{x}}}\in {{\mathbb{R}}}^3,  t  > 0\,,}$$ with initial values $$\label{eq:init_values}
\boxed{
   p({{\vec{x}}},0) = \rho({{\vec{x}}})\,, \quad \frac{\partial p}{\partial t}({{\vec{x}}},0) = 0 \qquad  {{\vec{x}}}\in {{\mathbb{R}}}^3\;.}$$ The quantity $\rho$ in (\[eq:ex:wave3d\]) and (\[eq:u(x)\]) is a combination of several physical parameters.  All along this paper $\rho$ should not be confused with the source term $$\label{eq:source}
\boxed{f({{\vec{x}}},t) = \frac{d j}{d  t}(t) \rho({{\vec{x}}})\,,\quad  {{\vec{x}}}\in {{\mathbb{R}}}^3,  t  > 0  \;.}$$ In PAT, some data about the pressure $p({{\vec{x}}},t)$ are measured and the main task is to reconstruct the initial pressure $\rho$ from these data.  While the excitation principle is always as described above and thus (\[eq:ex:ivp3d\]) holds, the specific type of data measured depends on the type of transducers used, and thus influences the mathematical model. 

Nowadays there is a trend to incorporate more and more modeling into photoacoustic.  In particular, taking into account locally varying *wave speed* and *attenuation*.  Even more there is a novel trend to *qualitative photoacoustics*, which is concerned with estimating physical parameters from the imaging parameter of standard photoacoustics.  In this paper we focus on attenuation correction, where we survey some recent progress.  Inversion with varying wave speed has been considered for instance in [@AgrKuc07; @HriKucNgu08], and is not further discussed here. 

The outline of this paper is as follows: First, we review existing attenuation models and discuss their causality properties, which we believe to be essential for algorithms for inversion with attenuated data.  Then, we survey causality properties of common attenuation models.  We also derive integro-differential equations which the attenuated waves are satisfying.  In addition we discuss the ill–conditionness of the inverse problem for calculating the unattenuated wave from the attenuated one. 

Attenuation {#sec:attenuation}
===========

The difficult issue of effects of and corrections for the attenuation of acoustic waves in PAT has been studied[@RivZhaAna06; @BurGruHalNusPal07; @PatGre06; @KowSchBon10], although no complete conclusion on the feasibility of these models has been reached. 

Mathematical models for describing attenuation are formulated in the frequency domain, taking into account that attenuation disperses high frequency components more rapidly over traveled distance.  Let ${\mathcal{G}}({{\vec{x}}},t)$ denote the attenuated wave which originates from an impulse ($\delta_{{{\vec{x}}},t}$-distribution) at ${{\vec{x}}}=0$ at time $t=0$. In mathematical terms ${\mathcal{G}}$ is the Green-function of attenuated wave equation.  Moreover, we denote by $$\boxed{
\label{eq:g0} {\mathcal{G}}_0({{\vec{x}}},t) =
\frac{\delta\left(t-\frac{{\left|{{\vec{x}}}\right|}}{c_0}\right)}{4\,\pi\,{\left|{{\vec{x}}}\right|}}}\index{Green function, standard wave equation}$$ the Green function of the unattenuated wave equation; that is, it is the solution of , with constant sound speed $c(x) \equiv c_0$ and initial conditions $${\mathcal{G}}_0({{\vec{x}}},0)=0 \quad\mbox{ and}\quad
   \frac{\partial {\mathcal{G}}_0}{\partial t}({{\vec{x}}},0)=\delta_{{{\vec{x}}},t}\,.$$ Common mathematical formulations of *attenuation* assume that $$\label{eq:GG0}
\boxed{
{\mathcal{F} \left\{ {\mathcal{G}}\right\}}({{\vec{x}}},\omega) =  {\exp \left(  -\beta^*
({\left|{{\vec{x}}}\right|},\omega) \right)}\, {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}}({{\vec{x}}},\omega)\,, \quad {{\vec{x}}}\in {{\mathbb{R}}}^3,\, \omega \in {{\mathbb{R}}}\;.}$$ Here ${\mathcal{F} \left\{ \cdot \right\}}$ denotes the Fourier transform with respect to time $t$ (cf.  Appendix \[sec:app\]).  Applying the inverse Fourier transform ${\mathcal{F}^{-1} \left\{ \cdot \right\}}$ to [(\[eq:GG0\])]{} gives $$\label{eq:FGG0}
\boxed{{\mathcal{G}}({{\vec{x}}},t) =   K({{\vec{x}}},t) *_t {\mathcal{G}}_0({{\vec{x}}},t)}
\qquad \mbox{($*_t$ time convolution)}$$ where $$\label{eq:kernel}
K({{\vec{x}}},t) := \frac{1}{\sqrt{2\,\pi}} {\mathcal{F}^{-1} \left\{ {\exp \left( -\beta^* ({\left|{{\vec{x}}}\right|},\cdot) \right)} \right\}}(t)\;.$$ From (\[eq:FGG0\]) and (\[eq:g0\]) it follows that $$\begin{aligned}
{\mathcal{G}}({{\vec{x}}},t) &= K({{\vec{x}}},t) *_t {\mathcal{G}}_0({{\vec{x}}},t)\\
&= \int_{{\mathbb{R}}}K({{\vec{x}}},t-\tau) \frac{\delta (\tau - \frac{{\left|{{\vec{x}}}\right|}}{c_0})}{4 \pi {\left|{{\vec{x}}}\right|}}\,d\tau \\
&= \frac{K\left({{\vec{x}}},t-\frac{{\left|{{\vec{x}}}\right|}}{c_0}\right)}{4\pi {\left|{{\vec{x}}}\right|}}\;. 
\end{aligned}$$ Consequently, $$\label{eq:Kandgreen}
\boxed{
{\mathcal{G}}({{\vec{x}}},t+{\left|{{\vec{x}}}\right|}/c_0) = K({{\vec{x}}},t)/(4\pi {\left|{{\vec{x}}}\right|})\;.}$$ Moreover, we emphasize that the Fourier transform of a real and even (real and odd) function is real and even (imaginary and odd).  Since ${\mathcal{G}}$ and ${\mathcal{G}}_0$ are real valued, $K$ must be real valued and consequently the real part $\Re(\beta^*)$ of $\beta^*$ has to be even with respect to the frequency $\omega$ and $\Im(\beta^*)$ has to be odd with respect to $\omega$. Attenuation is caused if $\Re(\beta^*)$ is positive and since then $\beta^*$ has a nonzero imaginary part due to the Kramers-Kronig relation, attenuation causes dispersion.  In the literature the following product ansatz is commonly used $$\label{eq:attenuation_law}
\boxed{
\beta^* ({\left|{{\vec{x}}}\right|},\omega) = \alpha^*(\omega)\,{\left|{{\vec{x}}}\right|} \qquad \omega\in{{\mathbb{R}}},\,{{\vec{x}}}\in{{\mathbb{R}}}^3\;.}$$ In the sequel we concentrate on these models and use the following terminology:

We call $\beta^*$ of standard form if (\[eq:attenuation\_law\]) holds.  Then the function $$\label{eq:attenuation_coefficient}
\alpha^*: {{\mathbb{R}}}\to {{\mathbb{C}}}$$ is called *standard attenuation coefficient* and $\alpha = \Re(\alpha^*)$ is called the *attenuation law*.  We also call $\beta^*$ the *attenuation coefficient*. 

From the relation (\[eq:attenuation\_law\]), it follows that $\Re (\alpha^*)$ is even, $\Im (\alpha^*)$ is odd, and $\Re (\alpha^*)>0$ (the last inequality guarantees attenuation). 

In the following we summarize common attenuation coefficients and laws: In what follows $\alpha_0$ denotes a positive parameter and $$\label{ta0} {\tilde{\alpha}}_0 = \frac{\alpha_0}{\cos
\left(\frac{\pi}{2}\gamma\right)}           \qquad\quad (0 < \gamma \not\in {{\mathbb{N}}})\,$$ is a possibly non-positive coefficient. 

-   [**Frequency Power Laws:**]{}

    -   Let $0 < \gamma \not\in {{\mathbb{N}}}$. The frequency power law *attenuation coefficient* is defined by $$\label{eq:powlaw1}
        \alpha_{pl}^*(\omega) = {\tilde{\alpha}}_0 (-{{\mathrm{i}}\,}\omega)^\gamma = {\tilde{\alpha}}_0 {\left|\omega\right|}^\gamma \left(
        \cos \left(\frac{\pi}{2} \gamma \right) - {{\mathrm{i}}\,}{\mbox{sgn}}(\omega) \sin \left(\frac{\pi}{2} \gamma \right)\right)$$ for $\omega \in {{\mathbb{R}}}$. Therefore, the *attenuation law* is given by $$\label{eq:powlaw1b}
        \boxed{
        \alpha_{pl}(\omega) = \alpha_0\,{\left|\omega\right|}^\gamma\;.}$$ These models have been considered for instance in [@Sza94; @Sza95; @WatHugBraMil00; @WaMoMi05]. 

    -   Let $\gamma=1$ and $\omega_0\neq0$, the attenuation coefficient is defined by $$\label{alphagamma1}
        \begin{aligned}
          \alpha_{pl}^*(\omega)
             := \alpha_0\,{\left|\omega\right|}
                + {{\mathrm{i}}\,}\frac{2}{\pi}\,\alpha_0\,\omega\,\log\left|\frac{\omega}{\omega_0}\right|
        \qquad \omega\in{{\mathbb{R}}}\,. 
        \end{aligned}$$ The attenuation law is $$\label{alphagamma1a}
        \boxed{
        \alpha_{pl}(\omega):=\alpha_0\,{\left|\omega\right|}\;.}$$ This model has been considered in [@Sza95; @WaMoMi05]. 

-   [[**Szabo:**]{}]{} Let $0 < \gamma\not\in{{\mathbb{N}}}$. The attenuation coefficient [^1] of Szabo’s law is defined by $$\label{alpha*szabo}
    \begin{aligned}
      \alpha_{sz}^*(\omega) &= \frac{1}{c_0}\,\sqrt{ (-{{\mathrm{i}}\,}\omega)^2
            + 2 {\tilde{\alpha}}_0 c_0 (-{{\mathrm{i}}\,}\omega)^{\gamma+1}  }
             + {{\mathrm{i}}\,}\frac{\omega}{c_0}\;. 
    \end{aligned}$$ We denote Szabo’s attenuation law by $$\boxed{
    \alpha_{sz}(\omega) := \Re (\alpha_{sz}^*(\omega))\;.}$$ For small frequencies $\alpha_{sz}(\omega)$ behaves like $\alpha_0\,|\omega|^\gamma$. This model has been considered in [@Sza94; @Sza95] where, in addition, also a model for $\gamma \in {{\mathbb{N}}}$ has been introduced. 

-   [[**Thermo-Viscous Attenuation Law:**]{} (see e.g. [@KinFreCopSan00; @Sza94]):]{} Here, for $\tau_0>0$, the attenuation coefficient is defined by $$\label{alpha*th}
    \begin{aligned}
      \alpha_{tv}^*(\omega)
         =  \frac{-{{\mathrm{i}}\,}\omega}{c_0\,\sqrt{1-{{\mathrm{i}}\,}\tau_0\,\omega}}
             +  \frac{{{\mathrm{i}}\,}\omega}{c_0}\;
    \end{aligned}$$ with attenuation law $$\label{alpha*threal}
    \begin{aligned}
    \boxed{
      \alpha_{tv}(\omega)
         =  \frac{\tau_0\, \omega^2}{ \sqrt{2}\,c_0\,\sqrt{ (1+\sqrt{1+(\tau_0\,\omega)^2})\,(1+(\tau_0\,\omega)^2)} } \;.} 
    \end{aligned}$$ For small frequencies $\alpha_{tv}(\omega)$ behaves like $\frac{\tau_0\,\omega^2}{2\,c_0}$. That is the thermo-viscous law approximates a power attenuation law with exponent $2$.

-   [[**Nachman, Smith and Waag [@NacSmiWaa90]:**]{}]{} Consider a homogeneous and isotropic fluid with density $\rho_0$ in which $N$ relaxation processes take place.  Then the attenuation coefficient of the model in [@NacSmiWaa90] reads as follows: $$\label{alpha*Nachman+}
    \begin{aligned}
      \alpha_{nsw}^*(\omega)
         =  \frac{-{{\mathrm{i}}\,}\omega}{c_0}\,\left[
                 \frac{c_0}{\tilde c_0}\, \sqrt{ \frac{1}{N}\, \sum_{m=1}^N
              \frac{1-{{\mathrm{i}}\,}\,\tilde\tau_m\,\omega}{1-{{\mathrm{i}}\,}\tau_m\,\omega} }
                -1 \right]\;. 
    \end{aligned}$$ All parameters appearing in (\[alpha\*Nachman+\]) are positive and real.  $\kappa_m$ and $\tau_m$ denote the compression modulus and the relaxation time of the $m-$th relaxation process, respectively, and $$\label{deftildec0tau0}
      \tilde c_0 := \frac{c_0}{\sqrt{1+\sum_{m=1}^N c_0^2\,\rho_0\,\kappa_m}}
    \quad\mbox{and}\quad
      \tilde \tau_m := \tau_m\,(1-N\,\tilde c_0^2\,\rho_0\,\kappa_m)\,.$$ The last two definitions imply that $$\label{proptildec0tau0}
      \frac{\tilde c_0^2}{c_0^2} = \frac{1}{N}\,\sum_{m=1}^N \frac{\tilde \tau_m}{\tau_m} \;.$$ We denote the according attenuation law by [^2] $$\boxed{
    \alpha_{nsw}(\omega) := \Re (\alpha_{nsw}^*(\omega))\;.}$$

-   [[**Greenleaf and Patch [@PatGre06]**]{}]{} consider for $\gamma \in {\{ 1,2 \}}$ the attenuation coefficient $$\alpha_{gp}^*(\omega) = \alpha_0\,{\left|\omega\right|}^\gamma\,,$$ which, since it is real, equals the attenuation law $$\label{alpha*Patch}
    \boxed{
          \alpha_{gp}(\omega) = \Re(\alpha_{gp}^*(\omega))\;.}$$

-   [**Chen and Holm [@CheHolm04]:**]{} This model describes the attenuation as a function of the absolute value of the vector-valued wave number ${\mathbf{k}}\in{{\mathbb{R}}}^3$ (instead of the frequency $\omega\in{{\mathbb{R}}}$).  Let ${\mathcal{F}}_{3D}$ denote the $3D-$Fourier transform $$\begin{aligned}
       {\mathcal{F}}_{3D} {\{  f({\mathbf{k}})  \}}({{\vec{x}}})
            = \frac{1}{\sqrt{(2\,\pi)^3}}\,\int_{{{\mathbb{R}}}^3} {\exp \left( {{\mathrm{i}}\,}{{\vec{x}}}\cdot{\mathbf{k}}\right)}\,  f({\mathbf{k}}) \,{\mbox{d}}{\mathbf{k}}\,,
    \end{aligned}$$ then the Green function of the attenuated equation is defined by $$\label{GreenChenHolm}
    \begin{aligned}
        {\mathcal{G}}({{\vec{x}}},t) = \frac{H(t)\,c_0^2}{(2\,\pi)^{3/2}}\,{\mathcal{F}}_{3D}\left\{ {\exp \left( A(\cdot)\,t \right)}\, \frac{\sin(B(\cdot)\,t)}{B(\cdot)} \right\}({{\vec{x}}})
    \end{aligned}$$ where, for given $\alpha_1>0$, $$\label{AkBk}
        A({\mathbf{k}}) := -\alpha_1\,c_0\,{\left|{\mathbf{k}}\right|}^\gamma\,,
    \qquad
        B({\mathbf{k}}) := c_0\,\sqrt{{\left|{\mathbf{k}}\right|}^2-\alpha_1^2\,{\left|{\mathbf{k}}\right|}^{2\,\gamma}}\,.$$

-   [**In [@KowSchBon10]**]{} we proposed $$\label{eq:powlaw2}
    \begin{aligned}
       \alpha_{ksb}^*(\omega)
           = \frac{\alpha_0\,(-{{\mathrm{i}}\,}\omega)}{c_0\,\sqrt{1+(-{{\mathrm{i}}\,}\tau_0\,\omega)^{\gamma-1}}}
    \qquad\qquad (\gamma\in (1,2],\,\tau_0>0)\,,
    \end{aligned}$$\[eq:powlaw2a\] where the square root is again the complex root with positive real part. 

    Let $\gamma\in (1,2]$. Then, for small frequencies we have $$\boxed{\alpha_{ksb}(\omega) \approx \frac{\alpha_0\,\sin(\frac{\pi}{2}(\gamma-1))}{2\,c_0\,\tau_0} \,|\tau_0\,\omega|^\gamma >0\,.}$$ Thus our model behaves like a power law for small frequencies. 

Distinctive features of unattenuated wave propagation (,i.e. the solution of the standard wave equation) are *causality* and *finite wave front velocity*.  It is reasonable to assume that the attenuated wave satisfies the same distinctive properties as well.  In the following we analyze causality properties of the standard attenuation models. 

Causality
=========

In the following we present some abstract definitions and basic notations.  In the remainder ${{\vec{x}}}$ will always denote a vector in three dimensional space.  When we speak about functions, we always mean generalized functions, such as for instance distributions or tempered distributions - we recall the definitions of (tempered) distribution in the course of the paper. 

\[def:defAc\] A function $f:=f({{\vec{x}}},t)$ defined on the Euclidean space over time (i.e. in ${{\mathbb{R}}}^4$) is said to be *causal* if it satisfies $f({{\vec{x}}},t) = 0$ for $t < 0$.

Let ${\mathcal{A}}: D \to D$ be a linear operator, where $\emptyset \neq D$ is an appropriate set of functions from ${{\mathbb{R}}}^4$ to ${{\mathbb{R}}}$. In this paper we always assume that ${\mathcal{A}}$ satisfies the following properties:

-   ${\mathcal{A}}$ is *shift invariant* in space and time.  That is, for every function $f$ and every shift $L:=L({{\vec{x}}},t):=({{\vec{x}}}-{{\vec{x}}}_0,t-t_0)$, with ${{\vec{x}}}_0 \in {{\mathbb{R}}}^3$ and $t_0 \in {{\mathbb{R}}}$, it holds that $${\mathcal{A}}(f \circ L)  = ({\mathcal{A}}f) \circ L\;.$$

-   ${\mathcal{A}}$ is *rotation invariant* in space.  That is, for every function $f$ and every rotation matrix $R$, it holds that $${\mathcal{A}}(R f)  = R({\mathcal{A}}f)\;.$$

-   ${\mathcal{A}}$ is *causal*.  That is, it maps causal functions to causal functions.  From (\[eq:a\_conv\]) it follows that ${\mathcal{A}}$ is causal, if and only if the associated Green function is causal. 

\[def:defA\] The *Green function* of ${\mathcal{A}}$ is defined by $${\mathcal{G}}:={\mathcal{G}}({{\vec{x}}},t)={\mathcal{A}}\delta_{{{\vec{x}}},t} ({{\vec{x}}},t)\;.$$

The operator ${\mathcal{A}}$ is uniquely determined by ${\mathcal{G}}$ and vice versa.  This follows from the fact that $$\label{eq:a_conv}
\begin{aligned}
{\mathcal{A}}f ({{\vec{x}}}_0,t_0) &= {\mathcal{A}}\left( \int_{{\mathbb{R}}}\int_{{{\mathbb{R}}}^3} f({{\vec{x}}}_0-{{\vec{x}}},t-t_0)
\delta_{{{\vec{x}}},t} ({{\vec{x}}},t)\,d {{\vec{x}}}dt \right) \\
&= \int_{{\mathbb{R}}}\int_{{{\mathbb{R}}}^3} f({{\vec{x}}}_0-{{\vec{x}}},t-t_0) {\mathcal{G}}({{\vec{x}}},t)\,d {{\vec{x}}}dt\,. 
\end{aligned}$$

Moreover, we use the following terminology and abbreviations:

-   From the rotation invariance of ${\mathcal{A}}$ it follows that $$\label{eq:T}
          \hat{T}({{\vec{x}}}):=\sup {\{ t : {\mathcal{G}}({{\vec{x}}},\tau) = 0 \text{ for all } \tau \leq t \}}\,,$$ is rotationally symmetric, which allows us to use the shorthand notation $$\label{defT}
          T(r) = \hat{T}({{\vec{x}}}) \text{ where } r = {\left|{{\vec{x}}}\right|}\;.$$ With this notation (\[eq:T\]) can be equivalently expressed as $$\label{eq:defT(r)0a}
          {\mathcal{G}}\left({{\vec{x}}}, t + T({\left|{{\vec{x}}}\right|}) \right) = 0 \text{ for every } t<0\;.$$ In physical terms $T({\left|{{\vec{x}}}\right|})$ denotes the *travel time of a wave front* originating at position $\mathbf{0}$ at time $t=0$ and traveling to ${{\vec{x}}}$.

-   Because ${\mathcal{G}}$ is rotationally symmetric we can write $${\mathcal{G}}({{\vec{x}}},T({\left|{{\vec{x}}}\right|})) = \hat{{\mathcal{G}}}(r,T(r)) \text{ with } r = {\left|{{\vec{x}}}\right|}\;.$$ Taking the inverse function of $T$, which we denote by $S=S(t)$, we then find $${\mathcal{G}}({{\vec{x}}},T({\left|{{\vec{x}}}\right|})) = \hat{{\mathcal{G}}}(S(t),t)\,,$$

-   The *wave front* is the set $${\cal W}:={\{ ({{\vec{x}}},T({\left|{{\vec{x}}}\right|})): {{\vec{x}}}\in {{\mathbb{R}}}^3 \}}\;$$

-   The *wave front speed* $V$ is the variation of the location of the wave front as a function of time.  That is, $$\label{eq:v}
          V(t)=\frac{{\mbox{d}}S}{{\mbox{d}}t}(t)= \left.  \frac{1}{T'(r)}\right|_{r=S(t)}\,.$$ Here $T'$ denotes the derivative with respect to the radial component $r$.

-   We say that ${\mathcal{A}}$ has a *finite speed of propagation* if there exists a constant $\hat{c}_0$ such that $$\label{eq:defT(r)}
          0 < \left(T'(r)\right)^{-1} \leq  \hat{c}_0 < \infty\;.$$ In this case it follows from (\[eq:v\]) that the wave front velocity satisfies $$\label{eq:v2}
          V(t) \leq \hat{c}_0 < \infty\;.$$

-   We call an operator ${\mathcal{A}}$ strongly causal, if it is causal and satisfies the finite propagation speed property. 

The following lemma addresses the case of attenuation coefficients of standard form and gives examples of strongly causal operators ${\mathcal{A}}$.

\[lemm:cconst\] Let $\beta^*({\left|{{\vec{x}}}\right|},\omega) = \alpha^*(\omega) {\left|{{\vec{x}}}\right|}$ be of the standard form (\[eq:attenuation\_law\]) and ${\mathcal{A}}$ (\[eq:a\_conv\]) be the operator defined by the Green function ${\mathcal{G}}$, which is defined in (\[eq:FGG0\]).  Then ${\mathcal{A}}$ is strongly causal if and only if for every ${{\vec{x}}}\in {{\mathbb{R}}}^3$ the function $$t \to \frac{1}{\sqrt{2\pi}} {\mathcal{F}^{-1} \left\{ {\exp \left( - \alpha^*(\omega){\left|{{\vec{x}}}\right|} \right)} \right\}}\,,$$ defined in (\[eq:kernel\]), is causal. 

We assume that ${\mathcal{A}}$ is strongly causal.  It follows from [@KowSchBon10 Theorem 3.1] that there exists a constant $c$, which is smaller than or equal to the wave speed $c_0$ from (\[eq:ex:ivp3d\]), which satisfies $T({\left|{{\vec{x}}}\right|}) = \frac{{\left|{{\vec{x}}}\right|}}{c}$ for all ${{\vec{x}}}\in {{\mathbb{R}}}^3$. Using the definitions of the travel time $T({\left|x\right|})$ and (\[eq:kernel\]), it follows from (\[eq:Kandgreen\]) that $t\to K({{\vec{x}}},t)$ is causal. 

Now, for every ${{\vec{x}}}\in {{\mathbb{R}}}^3$ let $K$ be causal.  Then from (\[eq:Kandgreen\]) it follows that $t \to {\mathcal{G}}\left({{\vec{x}}},t+{\left|{{\vec{x}}}\right|}/c_0)\right)$ is causal.  Since $T({\left|{{\vec{x}}}\right|})$ denotes the largest positive time period for which $t \to {\mathcal{G}}\left({{\vec{x}}},t+T({\left|{{\vec{x}}}\right|})\right)$ is causal, we have for $r>0$: $$\label{eq:lem1Tc0}
     0 <  \frac{r}{c_0} \leq   T(r) = \int_0^{r} \frac{1}{V(s)}\;.$$ Here $V$ is parameterized with respect to the distance $s$ at time $t$ of the wave front from its origin.  As shown in the proof of [@KowSchBon10 Theorem 3.1], the fact that $\beta^*$ is of standard form together with (\[eq:lem1Tc0\]) implies that there exist a constant $c$ such that $T(r) = r/c$ for all $r >0$. But then from  (\[eq:lem1Tc0\]) it follows $0<r/c_0\leq r/c<\infty$ and consequently $0 < c \leq c_0 < \infty$.

Finally we explain the above notation for the standard wave equation:

In the case of the standard wave equation the wave front is the support of the Green function ${\mathcal{G}}_0$, the wave front velocity is $c_0$, and $T({\left|{{\vec{x}}}\right|})=\frac{{\left|{{\vec{x}}}\right|}}{c_0}$ denotes the travel time of the wave front. 

Strong Causality of Attenuation Laws {#sec:n}
====================================

In this section we analyze causality properties of attenuation laws.  We split the section into two parts, where the first concerns numerical studies to determine the kernel function $K$, defined in [(\[eq:kernel\])]{}, and the second part contains analytical investigations. 

In Figures \[fig:powlaw\], \[fig:szabo\] and \[fig:thviscous\] we represent the attenuation kernels according to power, Szabo’s, and the thermo-viscous law. 

![ Simulation of $K({{\vec{x}}},t)$ for the frequency power law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. In the first example $\gamma < 1$ and thus the function is causal.  For all other cases it is non causal.[]{data-label="fig:powlaw"}](PowLaw01.eps "fig:"){height="4.0cm"} ![ Simulation of $K({{\vec{x}}},t)$ for the frequency power law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. In the first example $\gamma < 1$ and thus the function is causal.  For all other cases it is non causal.[]{data-label="fig:powlaw"}](PowLaw02.eps "fig:"){height="4.0cm"}\
![ Simulation of $K({{\vec{x}}},t)$ for the frequency power law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. In the first example $\gamma < 1$ and thus the function is causal.  For all other cases it is non causal.[]{data-label="fig:powlaw"}](PowLaw03.eps "fig:"){height="4.0cm"} ![ Simulation of $K({{\vec{x}}},t)$ for the frequency power law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. In the first example $\gamma < 1$ and thus the function is causal.  For all other cases it is non causal.[]{data-label="fig:powlaw"}](PowLaw04.eps "fig:"){height="4.0cm"}

![ Simulation of $K({{\vec{x}}},t)$ for Szabo’s frequency law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. []{data-label="fig:szabo"}](Szabo01.eps "fig:"){height="4.0cm"} ![ Simulation of $K({{\vec{x}}},t)$ for Szabo’s frequency law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. []{data-label="fig:szabo"}](Szabo02.eps "fig:"){height="4.0cm"}\
![ Simulation of $K({{\vec{x}}},t)$ for Szabo’s frequency law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. []{data-label="fig:szabo"}](Szabo03.eps "fig:"){height="4.0cm"} ![ Simulation of $K({{\vec{x}}},t)$ for Szabo’s frequency law with $(\gamma,\alpha_0)\in
{\{ (0.5,0.1581),\,(1.5,0.0316),\,(2.7,0.0071),\,(3.3,0.0027) \}}$, $c_0=1$ and ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. []{data-label="fig:szabo"}](Szabo04.eps "fig:"){height="4.0cm"}

![[* Left:*]{} $K({{\vec{x}}},t)$ defined by the complex thermo-viscous attenuation law with $\tau_0=10^{-5}$, $c_0=1$ and fixed ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. [*Right:*]{} The proposed law (\[eq:powlaw2\]) for $\gamma=2$ with $\alpha_1=1$, $\tau_0=10^{-5}$, $c_0=1$ and fixed ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$ is causal. []{ data-label="fig:thviscous"}](ThViscous01.eps "fig:"){height="4.0cm"} ![[* Left:*]{} $K({{\vec{x}}},t)$ defined by the complex thermo-viscous attenuation law with $\tau_0=10^{-5}$, $c_0=1$ and fixed ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$. [*Right:*]{} The proposed law (\[eq:powlaw2\]) for $\gamma=2$ with $\alpha_1=1$, $\tau_0=10^{-5}$, $c_0=1$ and fixed ${\left|{{\vec{x}}}\right|}=\frac{1}{4}$ is causal. []{ data-label="fig:thviscous"}](CThViscous01.eps "fig:"){height="4.0cm"}

The figures already indicate that power laws with index greater than $1$ violate causality.  In the following we support these computational studies by analytical considerations.  Thereby we make use of distribution theory, which we recall first.  Generally speaking *Distributions* are generalized functions:

\[def:distribution\] We use the abbreviations:

-   ${\mathcal{D}}:=C_0^\infty({{\mathbb{R}}},{{\mathbb{C}}})$ is the space of infinitely often differentiable functions from ${{\mathbb{R}}}$ to ${{\mathbb{C}}}$ which have compact support. 

-   ${\mathcal{S}}$ is the space of infinitely often differentiable functions from ${{\mathbb{R}}}$ to ${{\mathbb{C}}}$ which are *rapidly decreasing*.  A function $f: {{\mathbb{R}}}\to {{\mathbb{C}}}$ is rapidly decreasing if for all $i,j \in {{\mathbb{N}}}_0$ $${\left|x\right|}^i {\left|f^{(j)}(x)\right|} \to 0 \text{ for } {\left|x\right|} \to \infty\;.$$

-   ${\mathcal{S}}$ is a locally convex space (see [@Yos95] for a definition) with the topology induced by the family of semi-norms $$p_{P,j}(f) = \sup_{x \in {{\mathbb{R}}}} {\left|P(x) f^{(j)}(x)\right|} ,$$ where $P$ is a polynomial and $j \in {{\mathbb{N}}}_0$. The topology on a locally convex set is defined as follows: $U \subseteq {\mathcal{S}}$ is open, if for every $f \in U$ there exists ${\varepsilon}> 0$ and a finite non-empty set $J'$ of polynomials and a finite set of indices $K'$ such that $$\bigcap_{P \in J', k \in K'} {\{ g \in {\mathcal{S}}: p_{P,k}(g-f) < {\varepsilon}\}} \subseteq U\;.$$

-   The space of *tempered distributions*, ${\mathcal{S}}'$, is the space of linear continuous functionals and ${\mathcal{S}}$.

-   A functional $L: {\mathcal{S}}\to {{\mathbb{C}}}$ is continuous if there exists a constant $C > 0$ and a seminorm $p_{P,j}$ such that $${\left|L u\right|} \leq C p_{P,j}(u)\,, \text{ for every } u \in {\mathcal{S}}$$ (see [@Yos95 Sect.  I.6, Thm.  1]

In the following we give some examples of tempered distributions and review some of their properties.  The examples are taken from [@Yos95 Sec6.2, Ex.  3] and [@DauLio02_2 Remark 6]. 

[**(Examples of Tempered Distributions)**]{} \[ex:distr\_ex\]

-   Let $1 \leq p \leq \infty$ and $f \in L^p({{\mathbb{R}}},{{\mathbb{C}}})$, then the linear operator $T\phi = \int_{{{\mathbb{R}}}} f(x)\phi(x)\,dx$ is a tempered distribution.  In the following we identify $f$ and $T$, and this clarifies the terminology $f \in {\mathcal{S}}'$ later on. 

-   ${\mathcal{S}}\subseteq {\mathcal{S}}'$ - thereby already the above relation between functions and tempered distributions is used. 

-   Distributions with compact support are tempered distributions.  For instance the $\delta$-Distribution is a tempered distribution. 

-   Polynomials are tempered distributions. 

-   The functions $f$ of $L_{loc}^1({{\mathbb{R}}})$ which are uniformly bounded by a polynomial for ${\left|x\right|}$ sufficiently large, are tempered distributions. [ ^3]

\[le:distr\_ex\]

-   The pointwise limit $f:{{\mathbb{R}}}\to {{\mathbb{C}}}$ of a sequence of functions ${\{ f_n:{{\mathbb{R}}}\to {{\mathbb{C}}}\}} \subseteq {\mathcal{S}}'$, is again a tempered distribution. 

-   Let $f \in {\mathcal{S}}'$, then ${\mathcal{F} \left\{ f \right\}} \in {\mathcal{S}}'$ and ${\mathcal{F}^{-1} \left\{ f \right\}} \in {\mathcal{S}}'$.

In the following we review Theorem 4 on p294 ff from [@DauLio02_5] which characterized when a generalized function $f\in{\mathcal{S}}'({{\mathbb{R}}})$ is causal, that is, when ${\mbox{supp}}(f) \subseteq [0,\infty)$. Below we use the following notation $${{\mathbb{C}}}_{\varepsilon}:={\{ z\in{{\mathbb{C}}}: \Im(z)\geq {\varepsilon}\}}\;.$$

(Theorem 4 on p294 ff in [@DauLio02_5]) \[th:lion\] Let $f \in {\mathcal{S}}'({{\mathbb{R}}})$. Then $f$ is causal [^4] if and only if

1.   \[it1\_Lions\] There exists a function $F:{{\mathbb{C}}}_0:={\{ \xi + {{\mathrm{i}}\,}\eta : \eta \geq 0 \}}\to{{\mathbb{C}}}$, which is holomorphic in the interior $\mathring{{{\mathbb{C}}}}_0:={\{ \xi + {{\mathrm{i}}\,}\eta : \eta > 0 \}}$. [^5]

2.   \[it2\_Lions\] For all fixed $\eta>0$ and $\xi\in{{\mathbb{R}}}$, $F(\xi + {{\mathrm{i}}\,}\eta)$ is a tempered distribution with respect to the variable $\xi$ and for $\eta\to 0$ $F(\xi + {{\mathrm{i}}\,}\eta)$ is convergent (with respect to the weak topology on ${\mathcal{S}}'$) to ${\mathcal{F} \left\{ f \right\}}:{{\mathbb{R}}}\to{{\mathbb{C}}}$. [^6]

3.   \[it3\_Lions\] For every ${\varepsilon}> 0$, there exists a polynomial $P$ such that $${\left|F(z)\right|} \leq P({\left|z\right|})  {\qquad \text{for} \qquad} z \in {{\mathbb{C}}}_{\varepsilon}\,.$$

\[re:holomophic\] The definition of the Fourier transform in this chapter has a different sign as in [@DauLio02_5] and consequently also ${{\mathbb{C}}}_0$ denotes the upper half plane and not the lower half plane as in [@DauLio02_5]. 

For analyzing attenuation laws, we use the following corollary, which is derived from Theorem \[th:lion\]. 

\[co:lion\] Let $\alpha^*: {{\mathbb{R}}}\to {{\mathbb{C}}}$ be continuous and let there exist a holomorphic extension to ${{\mathbb{C}}}_0$, which for the sake of simplicity of notation is again denoted be $\alpha^*$. In addition, let $\alpha^* (\xi+{{\mathrm{i}}\,}\eta) \to \alpha^* (\xi)$ for $\eta \to 0$ pointwise.  We denote by $\alpha: {{\mathbb{C}}}_0 \to {{\mathbb{C}}}$ the real part of $\alpha^*$. [^7] Moreover, we assume that there exists a constant $C$ such that $$\label{estimate:rapidly}
\alpha(\omega) \geq C \text{ for all }\omega \in {{\mathbb{R}}}\;.$$

1.   If in addition $$\label{estimate:polynomial}
    \alpha(z) \geq C \text{ for all } z \in {{\mathbb{C}}}_0\;.$$ Then, for every ${{\vec{x}}}\in {{\mathbb{R}}}^3$, the function $$\label{31}
    t \to K({{\vec{x}}},t) := \frac{1}{\sqrt{2\pi}}{\mathcal{F}^{-1} \left\{ {\exp \left( -\alpha^*(\cdot){\left|{{\vec{x}}}\right|} \right)} \right\}}(t)$$ is causal. 

2.   On the other hand, if there exists $C_1>0$, $\mu>0$ and $C_2 \in {{\mathbb{R}}}$ and a sequence ${\{ z_n \}}$ in $\mathring{{{\mathbb{C}}}}_0$ such that $$\label{estimate:polynomial_negativ}
    \alpha(z_n) \leq - C_1 {\left|z_n\right|}^\mu - C_2\,,$$ then $K$ violates causality. 

Let ${{\vec{x}}}\in {{\mathbb{R}}}^3$ fixed.  We apply Theorem \[th:lion\] to $f(\cdot)=K({{\vec{x}}},\cdot)$. Therefore, we have $$\sqrt{2\pi} {\mathcal{F} \left\{ f \right\}}(\omega) = {\exp \left(  -\alpha^*(\omega) {\left|{{\vec{x}}}\right|} \right)}\;.$$ Under the assumption (\[estimate:rapidly\]), taking into account that $\alpha_{pl}^*$ is continuous, ${\mathcal{F} \left\{ f \right\}}$ is in $L_{loc}^1({{\mathbb{R}}})$ and bounded by a constant polynomial, thus in ${\mathcal{S}}'$ (cf.  Example \[le:distr\_ex\]), and consequently, according to Lemma \[le:distr\_ex\], $f \in {\mathcal{S}}'$. Therefore, the general assumption of Theorem \[th:lion\] is satisfied. 

The function $z \in {{\mathbb{C}}}_0 \mapsto F(z):={\exp \left(  -\alpha^*(z) {\left|{{\vec{x}}}\right|} \right)}$ is an extension of $2\pi {\mathcal{F} \left\{ f \right\}}(\omega)$, which is holomorphic in $\mathring{{{\mathbb{C}}}}_0$. Thus Item \[it1\_Lions\] of Theorem \[th:lion\] holds. 

-   For proving the first assertion, it follows from (\[estimate:polynomial\]) that for $z = \xi + {{\mathrm{i}}\,}\eta$ with $\eta \geq 0$, $$\label{eq:new}
    {\left|F(z)\right|} \leq {\exp \left( -C{\left|{{\vec{x}}}\right|} \right)}\,,$$ which, in particular, shows that for all $\eta \geq 0$, the functions $\xi \to F(\xi+{{\mathrm{i}}\,}\eta)$ is a tempered distribution (cf.  Lemma \[le:distr\_ex\]).  Hence Item \[it1\_Lions\] of Theorem \[th:lion\] holds. 

    Moreover, since by assumption $\alpha^*(\xi+{{\mathrm{i}}\,}\eta) \to \alpha^*(\xi)$ for $\eta \to 0$ pointwise, $F(\xi+{{\mathrm{i}}\,}\eta)$ converges to $F(\xi)=2\pi {\mathcal{F} \left\{ f \right\}}(\xi)$ pointwise.  Because the limit is a tempered distribution and the convergence is with respect to the weak topology ${\mathcal{S}}'$ (cf.  Lemma \[le:distr\_ex\]).  Hence Item \[it2\_Lions\] of Theorem \[th:lion\] holds. 

    Moreover, from (\[eq:new\]) it follows that ${\left|F(z)\right|}$ is bounded by a constant polynomial.  Hence Item (\[it3\_Lions\]) of Theorem \[th:lion\] holds and therefore Theorem \[th:lion\] guarantees that $t \mapsto K({{\vec{x}}},t)$ is causal. 

-   For the second case, Item \[it3\_Lions\] of Theorem \[th:lion\] is violated.  Consequently, $K({{\vec{x}}},\cdot)$ is not causal. 

Power Laws {#power-laws .unnumbered}
----------

\[coro:powlaw1\] Let $0 < \gamma \in {{\mathbb{R}}}$, be not an odd number, and $\omega \in {{\mathbb{R}}}\mapsto \alpha_{pl}^*(\omega)={\tilde{\alpha}}_0 (-{{\mathrm{i}}\,}\omega)^\gamma$ be the power law attenuation coefficient from (\[eq:powlaw1\]), with ${\tilde{\alpha}}_0=\alpha_0/\cos \left(\frac{\pi}{2}\gamma\right)$ as in (\[ta0\]).  Then, the function $K$, defined in (\[31\]), is causal if and only if $\gamma\in (0,1)$.

Let ${{\vec{x}}}\in {{\mathbb{R}}}^3$ be fixed.  The function $z \in {{\mathbb{C}}}\mapsto \alpha_{pl}^*(z)={\tilde{\alpha}}_0 (-{{\mathrm{i}}\,}z)^\gamma$ is the holomorphic extension of $\omega \in {{\mathbb{R}}}\mapsto \alpha_{pl}^*(\omega)$. We prove or disprove causality by using Corollary \[co:lion\]. 

For $z={\left|z\right|} {\exp \left( {{\mathrm{i}}\,}\phi \right)}$ it follows from (\[eq:powlaw1\]) that $$\Re((-{{\mathrm{i}}\,}z)^\gamma)
        = \Re({\left|z\right|}^\gamma\,{\exp \left( {{\mathrm{i}}\,}\gamma\,(\phi-\pi/2) \right)})
        = {\left|z\right|}^\gamma\,\cos(\gamma(\phi-\pi/2))\,.$$ This implies that $$\alpha_{pl}(z) =  {\tilde{\alpha}}_0 \Re((-{{\mathrm{i}}\,}z)^\gamma)  =  {\tilde{\alpha}}_0 \cos(\gamma(\phi-\pi/2)){\left|z\right|}^\gamma\;.$$ In particular, if $z = \omega \in {{\mathbb{R}}}$, then $\phi$ is either $0$ or $\pi$. Taking into account the definition of ${\tilde{\alpha}}_0$ and that the $\cos$-function is symmetric around the origin, it follows that $$\label{relalphapowlaw1a}
\alpha_{pl}(\omega) = \alpha_0 {\left|w\right|}^\gamma \geq 0\;.$$ Thus (\[estimate:rapidly\]) holds. 

-   Let $\gamma \in (0,1)$: Every $z={\left|z\right|}\,{\exp \left( {{\mathrm{i}}\,}\phi \right)}\in {{\mathbb{C}}}_0$ satisfies $\phi\in [0,\pi]$. Consequently $\gamma(\phi-\pi/2) \in [-\pi/2,\pi/2]$ and thus $\cos(\gamma(\phi-\pi/2))$ is uniformly non-negative.  Even more for $\gamma \in (0,1)$ the coefficient ${\tilde{\alpha}}_0$, defined in (\[ta0\]), is positive.  In summary, we have that there exists a constant $C_1 \geq 0$ such that $$\label{relalphapowlaw1b}
    \alpha_{pl}(z) \geq C_1 {\left|z\right|}^\gamma \geq 0 {\qquad \text{ for } \qquad} z\in {{\mathbb{C}}}_0\;.$$ Thus (\[estimate:polynomial\]) holds and application of Corollary \[co:lion\] shows that $K$ is causal. 

-   Let $\gamma\in (1,3)\cup (5,7)\cup\cdots$. Then ${\tilde{\alpha}}_0<0$. The sequence $$\label{eq:34a}
    {\{ z_n := n {\exp \left( {{\mathrm{i}}\,}\pi/2 \right)} = {{\mathrm{i}}\,}n \}}_{n \in {{\mathbb{N}}}}$$ consists of elements of $\mathring{{{\mathbb{C}}}}_0$ and satisfies assumption (\[estimate:polynomial\_negativ\]), that is, $$\label{eq:azn1}
    \alpha_{pl}(z_n) = \underbrace{{\tilde{\alpha}}_0}_{<0} {\left|z_n\right|}^\gamma\;.$$ Application of Corollary \[co:lion\] shows that $K$ is not causal. 

-   Let $\gamma\in (3,5)\cup (7,9)\cup\cdots$. We fix some $0<\delta<\pi/2$, and define $$\phi := \left(1+\frac{1}{\gamma} \right)\frac{\pi}{2}  + \frac{\delta}{\gamma}\;.$$ The sequence $$\label{eq:34b}
    {\{ z_n:=n\,{\exp \left( {{\mathrm{i}}\,}\phi \right)} \}}$$ consists of elements of $\mathring{{{\mathbb{C}}}}_0$. Under the above assumptions, it follows that ${\tilde{\alpha}}_0 > 0$ and therefore $$\label{eq:azn2}
    \alpha_{pl} (z_n) = \underbrace{{\tilde{\alpha}}_0}_{>0} \underbrace{\cos(\pi/2+\delta)}_{<0} {\left|z_n\right|}^\gamma\;.$$ Thus form Corollary \[co:lion\] the assertion follows. 

In the following we analyze the following family of variants of power laws: $$\label{szabo_model}
\alpha_{pl+}^*(\omega) = {\tilde{\alpha}}_0 (-{{\mathrm{i}}\,}\omega)^\gamma  + \alpha_1(-{{\mathrm{i}}\,}\omega)\,,$$ which have been considered in [@Sza95; @WatHugBraMil00]. 

\[le:pow\_law\_extended\] Let $0 < \gamma \notin {{\mathbb{N}}}$ and $\alpha_{pl+}^*$ as defined in (\[szabo\_model\]).  Moreover, let $K$ be as in (\[31\]).  Then, if $\gamma >1$, $K$ is not causal.  For $\gamma \in (0,1)$ $K$ is causal if and only if $\alpha_1\in [0,\infty)$.

The holomorphic extension of $\omega \in {{\mathbb{R}}}\to \alpha_{pl+}^*(\omega)$ is the function $$\alpha_{pl+}^*(z) = {\tilde{\alpha}}_0 (-{{\mathrm{i}}\,}z)^\gamma  + \alpha_1(-{{\mathrm{i}}\,}z)$$ and consequently $$\alpha_{pl+}(z)
= {\tilde{\alpha}}_0 {\left|z\right|}^\gamma \left(1+\frac{\alpha_1}{{\tilde{\alpha}}_0} {\left|z\right|}^{1-\gamma} \right)
\cos \left(\gamma \left(\phi-\frac{\pi}{2}\right)\right)\;.$$

-   For $\gamma > 1$ we have that $1+\frac{\alpha_1}{{\tilde{\alpha}}_0} {\left|z\right|}^{1-\gamma} \to 1$ for ${\left|z\right|} \to \infty$. Let ${\{ z_n \}}$ as defined in (\[eq:34a\]) or (\[eq:34b\]).  Then, since for both sequences ${\left|z_n\right|} \to \infty$, it follows from (\[eq:azn1\]), (\[eq:azn2\]) that the according sequences ${\{ z_n \}}$ satisfy (\[estimate:polynomial\_negativ\]) for $n$ sufficiently large, respectively.  Thus $K$ is not causal. 

-   For $\gamma \in (0,1)$ and $\alpha_1 \geq 0$ the assertion follows already from the fact that $\alpha_{pl+}(\omega) \geq \alpha_{pl}(\omega)$ and that the later already satisfies (\[estimate:polynomial\]).  Thus $K$ is causal. 

-   Let $\gamma \in (0,1)$ and $\alpha_1 < 0$. Then for some $0 < \delta < -\alpha_1$ fixed, we can find a constant $C_2$ such that for all $z \in {{\mathbb{C}}}$ $${\tilde{\alpha}}_0 {\left|z\right|}^\gamma + \alpha_1 {\left|z\right|} \leq \underbrace{(\alpha_1+\delta)}_{< 0} {\left|z\right|} - C_2\;.$$ Consequently, for ${\{ z_n={{\mathrm{i}}\,}n \}}$, we have $$\alpha_{pl+}(z_n) \leq (\alpha_1+\delta) {\left|n\right|} - C_2\;.$$ which shows (\[estimate:polynomial\_negativ\]).  Thus $K$ is not causal. 

Powerlaw with $\gamma = 1$ {#powerlaw-with-gamma-1 .unnumbered}
--------------------------

\[thcaus02\] Let $\alpha_{pl}^*$ be as in defined in (\[alphagamma1\]).  Then the function $K$, defined in (\[31\]) is not causal. 

First, we prove that $$\begin{aligned}
z \in {{\mathbb{C}}}_0 \mapsto \hat{\alpha}_{pl}^*(z) := \alpha_0 z
        + {{\mathrm{i}}\,}\frac{2 \alpha_0}{\pi} z \log \left(\frac{z}{\omega_0}\right)
\end{aligned}$$ is the holomorphic extension of $\omega \to \alpha_{pl}^*(\omega)$. This assertion follows from the facts $$\begin{aligned}
\hat{\alpha}_{pl}^*(\omega)& =\alpha_{pl}^*(\omega)  \text{ for } & \omega>0\,,\\
\lim_{\eta \to 0+} {{\mathrm{i}}\,}\frac{2}{\pi} \log \left(\frac{\omega + {{\mathrm{i}}\,}\eta}{\omega_0}\right)
      & = {{\mathrm{i}}\,}\frac{2}{\pi}\,\log\left|\frac{\omega}{\omega_0}\right|
        - 2 \text{ for } & \omega<0\;. 
\end{aligned}$$ Since $$\alpha_{pl}(\omega) = \Re( \alpha_{pl}^*(\omega)) = \alpha_0\,{\left|\omega\right|} \geq 0\,,$$ Corollary \[co:lion\] is applicable.  For the elements of the sequence ${\{ z_n:={{\mathrm{i}}\,}n \}}_{n \in {{\mathbb{N}}}}$ in ${{\mathbb{C}}}_0$ $$\hat{\alpha}_{pl}^*(z_n)
       = - \frac{2 \alpha_0}{\pi} n \log \left(\frac{n}{\omega_0}\right)$$ is real and therefore equals $\alpha_{pl}(z_n)$ and thus (\[estimate:polynomial\_negativ\]) holds.  Thus Corollary \[co:lion\] gives the assertion. 

Szabo’s Model: {#szabos-model .unnumbered}
--------------

\[theo:szabo\] For $\alpha_0>0$ let $\alpha_{sz}^*$ be the coefficient of Szabo’s model (\[alpha\*szabo\]).  Then, for $\gamma\in (0,1)$, the function $K$ (\[eq:kernel\]) is a causal function and for $\gamma>1$ with $\gamma\not\in{{\mathbb{N}}}$, $K$ violates causality. 

Without loss of generality we assume that $c_0=1$. The holomorphic extension of $\alpha_{sz}^*:{{\mathbb{R}}}\to {{\mathbb{C}}}$ from (\[alpha\*szabo\]) is $$z\in {{\mathbb{C}}}_0 \to \alpha_{sz}^*(z) = (-{{\mathrm{i}}\,}z)
  \left[ \sqrt{1 + 2 {\tilde{\alpha}}_0 (-{{\mathrm{i}}\,}z)^{\gamma-1}} - 1\right]\,.$$

First, we make some general manipulations which can be used in several ways: Let $z = \xi+{{\mathrm{i}}\,}\eta \in {{\mathbb{C}}}_0$. We use the polar representation $$z = {\left|z\right|} {\exp \left( {{\mathrm{i}}\,}\phi \right)}\,, \quad \phi := \phi(z) \in \left[0,\pi\right]\;.$$ Then $$\alpha_{sz}^*(z) = {\left|z\right|} {\exp \left( {{\mathrm{i}}\,}(\phi-\pi/2) \right)} \Psi(z)\,,$$ where $$\label{eq:delta}
\begin{aligned}
\Psi(z) &:= \sqrt{\hat{\Psi}(z)} - 1\,,\\
\hat{\Psi}(z) &:= 1 + 2 {\tilde{\alpha}}_0 {\left|z\right|}^{\gamma-1} {\exp \left( {{\mathrm{i}}\,}\delta \right)}\,,\\
\delta &:= \delta(z):=(\phi-\pi/2)(\gamma-1)\;. 
\end{aligned}$$ With this notation we have $$\label{eq:hatpsi}
\begin{aligned}
\Re (\hat{\Psi}(z)) &= 1 + 2 {\tilde{\alpha}}_0 \cos(\delta) {\left|z\right|}^{\gamma-1}\,,\\
\Im (\hat{\Psi}(z)) &= 2 {\tilde{\alpha}}_0 \sin(\delta) {\left|z\right|}^{\gamma-1}\,,\\
{\left|\hat{\Psi}(z)\right|} &= {\left|z\right|}^{\gamma-1} \sqrt{(1 + 2 {\tilde{\alpha}}_0 \cos(\delta))^2+4 {\tilde{\alpha}}_0^2 \sin(\delta)^2}\;. 
\end{aligned}$$ Representing $\hat{\Psi}$ in polar coordinates, $$\hat{\Psi}(z) = {\left|\hat{\Psi}(z)\right|} {\exp \left( {{\mathrm{i}}\,}\theta(z) \right)}\,,$$ we get $$\label{theta}
\begin{aligned}
     \sqrt{\hat{\Psi}(z)} &= \sqrt{{\left|\hat{\Psi}(z)\right|}} {\exp \left( {{\mathrm{i}}\,}\theta(z)/2 \right)}
\qquad\quad \text{ with }\\
      \theta (z) &= \arctan (\Im (\hat{\Psi}(z))/ \Re (\hat{\Psi}(z))))  \qquad \theta\in (-\pi,\pi) \;. 
\end{aligned}$$ Note, that $\sqrt{\hat{\Psi}(z)}$ is the complex root with non-negative real part, which meets the general assumption of the paper.  Moreover, we have $$\label{szz}
\alpha_{sz}(z) = \Re (\alpha_{sz}^*(z))
          = \eta \Re (\Psi(z)) + \xi \Im (\Psi (z))\;.$$

First, we prove that $\Re (\Psi(z)) \geq 0$: We use the elementary inequality $$\cos (\theta(z)) \leq \cos^2 (\theta(z)/2)\,,$$ and $\cos (\theta(z)/2)\geq 0$ which imply that $$\label{eq:real}
\begin{aligned}
    \Re \left(\sqrt{\hat{\Psi}(z)}\right)
      &= \sqrt{{\left|\hat{\Psi}(z)\right|}} \cos (\theta(z)/2) \\
      &=  \sqrt{{\left|\hat{\Psi}(z)\right|}\, \cos^2 (\theta(z)/2) }\\
      &\geq  \sqrt{\Re (\hat{\Psi}(z))}\;. 
\end{aligned}$$

-   Now, let $\gamma\in (0,1)$. Since $\Re (\hat{\Psi}(z)) \geq 1$ for $\gamma\in (0,1)$, it follows that for all $z \in {{\mathbb{C}}}_0$ $$\eta \Re (\Psi(z))
               = \eta \Re \left(\sqrt{\hat{\Psi}(z)}\right) - \eta
              \geq  \eta \sqrt{\Re (\hat{\Psi}(z))} - \eta
    \geq 0\,.$$ Thus $\eta \Re (\Psi(z)) \geq 0$.

    Now we show that $z \to \xi \Im (\Psi (z))$ is uniformly bounded from below by $0$ in ${{\mathbb{C}}}_0$. Thus according to (\[szz\]) $\alpha_{sz}$ is uniformly bounded from below, and thus from Corollary \[co:lion\], it follows that $t \to K({{\vec{x}}},t)$ is causal. 

    Using the definition of $\theta$, (\[theta\]), and the facts that $\delta \in [0,(1-\gamma)\pi/2]$ for $\phi \in [0,\pi/2]$ and $\delta \in [(\gamma-1)\pi/2,0)$ for $\phi \in (\pi/2,\pi]$ it follows from the monotonicity of $\tan$ on $(-\pi,\pi)$ that $$\theta(z) = \arctan\left(\frac{2 {\tilde{\alpha}}_0 \sin(\delta) {\left|z\right|}^{\gamma-1}}{1 + 2 {\tilde{\alpha}}_0 \cos(\delta) {\left|z\right|}^{\gamma-1}}\right) \in
    \left\{ \begin{array}{rcl}
    ~[0,\delta] & \text{ for all } & \phi \in [0,\pi/2]\,,\\
    ~[\delta,0) & \text{ for all } & \phi \in (\pi/2,\pi]
    \end{array}
    \right.$$ Now, noting that ${\mbox{sgn}}(\xi) = {\mbox{sgn}}\left(\sin (\theta(z)/2)\right)$ it follows that $$\label{eq:imaginary}
    \xi \Im (\Psi(z)) = \xi \sqrt{{\left|\hat{\Psi}(z)\right|}} \sin (\theta(z)/2) \geq 0\,.$$ Thus the assertion follows from Corollary \[co:lion\]. 

-   Assume that $\gamma  > 1$. Let $z=\xi+{{\mathrm{i}}\,}\eta$ with $\eta = 0$. Since the square root in [(\[eq:delta\])]{} is such that $\Re(\Psi(z))>0$ for $z\in{{\mathbb{C}}}_0$, property (\[propsqrt\]) in the Appendix implies $\xi \Im (\Psi(z))>0$ for $z=\xi$ and hence $\alpha_{sz}(z=\xi) \geq 0$. Thus (\[estimate:rapidly\]) holds and we can apply Corollary \[co:lion\]. 

    -   Let $\gamma\in (1,3)\cup (5,7)\cup\cdots$, which implies that $\cos\left(\gamma\,\pi/2\right)<0$, and consequently, ${\tilde{\alpha}}_0 < 0$. For sufficiently large $n$ the elements of the sequence ${\{ z_n:={{\mathrm{i}}\,}n \}}$ satisfy $$\alpha_{sz}(z_n) = n\,\Re\left(\sqrt{1 - 2{\left|{\tilde{\alpha}}_0\right|}\, n^{\gamma-1}}-1\right) \leq -n\,,$$ which shows that (\[estimate:polynomial\_negativ\]) holds with $\mu=1$, $C_1=1/2$ and $C_2=0$, and hence the assertion follows from Corollary \[co:lion\]. 

    -   Let $\gamma\in (3,5)\cup (7,9)\cup\cdots$, which implies that $\cos\left(\gamma\,\pi/2\right)>0$, and consequently, ${\tilde{\alpha}}_0>0$. Now, let $z_n:=n {\exp \left( {{\mathrm{i}}\,}\phi \right)}$ with $$\phi:=  \frac{\pi}{\gamma-1} + \frac{\pi}{2}\,.$$ Since $\gamma>3$, we have $\phi\in (\pi/2,\pi)$, and therefor $\Re(z_n)=n\,\cos(\phi)<0$ and $\Im(z_n)=n\,\sin(\phi)>0$. Moreover, $$\hat \Psi(z_n) = 1-2\,{\left|{\tilde{\alpha}}_0\right|}\,n^{\gamma-1}\,$$ and thus for sufficiently large $n$ we have $\Re(\Psi(z_n))=-1$ and\
        $\Im\left(\sqrt{\hat\Psi(z_n)}\right)\geq C_1\,n^{(\gamma-1)/2}$ for some constant $C_1>0$. Hence it follows that $$\begin{aligned}
        \alpha_{sz}(z_n)
                 &=  -\Im(z_n) + \Re(z_n)\, \Im(\sqrt{\hat\Psi(z_n)})\\
                 &= - n{\left|\sin(\phi)\right|} - n{\left|\cos(\phi)\right|}\, C_1\,n^{(\gamma-1)/2}
        \end{aligned}$$ and therefore (\[estimate:polynomial\_negativ\]) holds.  Thus from Corollary \[co:lion\] the assertion follows. 

Thermo-Viscous Attenuation Law {#thermo-viscous-attenuation-law .unnumbered}
------------------------------

Let $c_0,\tau_0>0$ and let $\alpha_{tv}^*$ as defined in (\[alpha\*th\]).  Then the kernel function $K$ violates causality. 

Since the function $z\in {{\mathbb{C}}}_0 \to 1-{{\mathrm{i}}\,}\tau_0\,z$ does not vanish, the function $$z \in {{\mathbb{C}}}_0 \to \alpha_{tv}^*(z)
     =  \frac{-{{\mathrm{i}}\,}z}{c_0\,\sqrt{1-{{\mathrm{i}}\,}\tau_0 z}} +  \frac{{{\mathrm{i}}\,}z}{c_0}$$ is the holomorphic extension of $\omega \in {{\mathbb{R}}}\to \alpha_{tv}^*(\omega)$. That (\[estimate:rapidly\]) holds follows from the identity (\[alpha\*threal\]).  For the sequence ${\{ z_n  \}}_{n\in{{\mathbb{N}}}} := {\{ {{\mathrm{i}}\,}n \}}_{n\in{{\mathbb{N}}}}$ we get for $n$ sufficiently large $$\alpha_{tv}^*(z_n) = \frac{n}{c_0}\,\left[ \frac{1}{\sqrt{1+\tau_0\,n}} -1 \right] \leq - \frac{1}{2\,c_0} n$$ Thus (\[estimate:polynomial\_negativ\]) holds. 

The second part of Corollary \[co:lion\] implies that $K$ is not causal. 

Model of Nachman, Smith and Waag {#model-of-nachman-smith-and-waag .unnumbered}
--------------------------------

Let $\alpha_{nsw}^*$ as in (\[alpha\*Nachman+\]).  If $$\label{asstau}
  \tilde \tau_m <  \tau_m  \qquad \mbox{ for all} \qquad m\in\{1,\,\ldots\,,N\}\,,$$ then the kernel function $K$ is causal. 

Since for all $z \in {{\mathbb{C}}}_0$ and all $1 \leq m \leq N$, $1-{{\mathrm{i}}\,}\tau_m\,z$ does not vanish, $$\begin{aligned}
z \in {{\mathbb{C}}}_0 \to  \alpha_{nsw}^*(z)
     =  \frac{-{{\mathrm{i}}\,}z}{c_0}\,\left[
             \frac{c_0}{\tilde c_0}\, \sqrt{ \frac{1}{N}\, \sum_{m=1}^N
          \frac{1-{{\mathrm{i}}\,}\,\tilde\tau_m\,z}{1-{{\mathrm{i}}\,}\tau_m\,z} }
            -1 \right]
\end{aligned}$$ is the holomorphic extension of $\omega \in {{\mathbb{R}}}\to \alpha_{nsw}^*(\omega)$.

We use a similar notation as in Proposition \[theo:szabo\].  $$z = \xi + {{\mathrm{i}}\,}\eta = {\left|z\right|} {\exp \left( {{\mathrm{i}}\,}\phi \right)} \in {{\mathbb{C}}}_0\,,$$ with some $\phi \in [0,\pi]$. $$\alpha_{nsw}^*(z) = \frac{{\left|z\right|}}{c_0} {\exp \left( {{\mathrm{i}}\,}(\phi-\pi/2) \right)} \,(\Psi(z)-1)\,,$$ where $$\Psi(z) = \sqrt{\sum_{m=1}^N \hat{\Psi}_m(z)}\, \quad\mbox{ with }\quad
  \hat{\Psi}_m(z) :=  \frac{1}{N}\, \frac{c_0^2}{\tilde c_0^2}\,
          \frac{1-{{\mathrm{i}}\,}\,\tilde\tau_m\,z}{1-{{\mathrm{i}}\,}\tau_m\,z}  \,,$$

In the following we show that for all $z\in {{\mathbb{C}}}_0$ $$\label{z2A1z1A2}
c_0 \alpha_{nsw}(z) = \eta (\Re (\Psi)(z)-1) + \xi \Im (\Psi)(z) >0\,,$$ which means that (\[estimate:polynomial\]) holds.  Then, according to Corollary \[co:lion\] the function $t \to K({{\vec{x}}},t)$ is causal. 

As in the proof of Proposition \[theo:szabo\] we prove $\eta (\Re (\Psi)(z)-1)>0$ and $\xi \Im (\Psi)(z) >0$.

-   Taking into account (\[proptildec0tau0\]) we define $$\label{s}
       s := \frac{1}{N}\, \frac{c_0^2}{\tilde c_0^2}
                 = \left(\sum_{m=1}^N \frac{\tilde\tau_m}{\tau_m}\right)^{-1}\,.$$ Using this notation, we get $$\begin{aligned}
      \hat{\Psi}_m(z) &=
              s\,\frac{(1+\tilde\tau_m\,\tau_m\,{\left|z\right|}^2) +{{\mathrm{i}}\,}\,(\tau_m\,\bar{z} -\tilde\tau_m\,z)}
               {|1-{{\mathrm{i}}\,}\,\tau_m\,z|^2}  \\
           &=   s\,\frac{\tilde\tau_m}{\tau_m}\,
          \frac{(\tau_m/\tilde \tau_m + \tau_m^2\,{\left|z\right|}^2 + (\tau_m^2/\tilde \tau_m+\tau_m)\,\eta)
                       +{{\mathrm{i}}\,}\,(\tau_m^2/\tilde \tau_m -\tau_m)\,\xi}
               {1+\tau_m^2\,{\left|z\right|}^2 + 2\,\tau_m\,\eta} \,,
    \end{aligned}$$ Because $\tau_m/\tilde \tau_m>1$, by assumption (\[asstau\]), it follows that for all $z \in {{\mathbb{C}}}_0$ $$\Re(\hat{\Psi}_m(z)) > s\,\frac{\tilde\tau_m}{\tau_m}\;.$$ Consequently, by using the definition of $s$, (\[s\]), it follows that $$\Re\left(\sum_{m=1}^N \hat{\Psi}_m(z)\right) > s\, \sum_{m=1}^N\frac{\tilde\tau_m}{\tau_m} = 1.$$ Now, using (\[eq:real\]) it follows $$\Re(\Psi(z))
          =    \Re\left(\sqrt{\sum_{m=1}^N \Psi_m(z)}\right)
         \geq  \sqrt{\Re\left(\sum_{m=1}^N \Psi_m(z)\right)} > 1$$ and consequently $\eta (\Re (\Psi(z))-1)>0$.

-   We have $$\begin{aligned}
      \Im(\hat{\Psi}_m(z)) =  s\,\frac{\tilde\tau_m}{\tau_m}\,
          \frac{(\tau_m^2/\tilde \tau_m -\tau_m)\,\xi}
               {1+\tau_m^2\,{\left|z\right|}^2 + 2\,\tau_m\,\eta}
    \end{aligned}$$ together with the assumption (\[asstau\]), which state that $\tau_m/\tilde \tau_m>1$, it follows that ${\mbox{sgn}}(\Im(\hat{\Psi}_m(z)))={\mbox{sgn}}(\xi)$. According to our assumption, we take that complex root, such that the real part of the argument is non-negative which together with property (\[propsqrt\]) in the Appendix implies $${\mbox{sgn}}\left( \Im\left(\sqrt{\sum_{m=1}^N \hat{\Psi}_m(z)} \right) \right)
         = {\mbox{sgn}}\left( \Im\left( \sum_{m=1}^N \hat{\Psi}_m(z)  \right) \right) \,.$$ Therefore, $$\begin{aligned}
    {\mbox{sgn}}(\Im(\Psi(z))) &= {\mbox{sgn}}\left(\Im\left(\sqrt{\sum_{m=1}^N \hat{\Psi}_m(z)}\}\right)\right)\\
                   &= {\mbox{sgn}}\left(\Im\left(\sum_{m=1}^N \hat{\Psi}_m(z)\right)\right)
                   = {\mbox{sgn}}(\xi)\;. 
    \end{aligned}$$ This shows the assertion. 

Our Model {#our-model .unnumbered}
---------

\[thcaus01\] For $\alpha_0,\,\tau_0>0$ and $\gamma\in (1,2]$ let $\alpha_{ksb}^*$ be defined as in (\[eq:powlaw2\]).  Then $K$, as defined in (\[31\]), is causal. 

The function $\hat{z}\to 1+(-{{\mathrm{i}}\,}\tau_0 \hat{z})$ does not vanish in ${{\mathbb{C}}}_0$. Thus the holomorphic extension of $\omega \to \alpha_{ksb}^*(\omega)$ is given by $$\hat{z} \in {{\mathbb{C}}}_0 \to \alpha_{ksb}^*(\hat{z})
          = \frac{\alpha_0(-{{\mathrm{i}}\,}\hat{z})}{c_0\sqrt{1+(-{{\mathrm{i}}\,}\tau_0 \hat{z})^{\gamma-1}}}\;.$$ In the following let $\hat{z} \in {{\mathbb{C}}}_0$. For proving (\[estimate:polynomial\]) we make a variable transformation $$\alpha_{ksb}^*(\hat{z})
    =  \frac{-{{\mathrm{i}}\,}\tau_0 \hat{z}}{\sqrt{1+(-{{\mathrm{i}}\,}\tau_0 \hat{z})^{\gamma-1}}}
    =  \frac{\alpha_0}{\tau_0 c_0} \frac{-{{\mathrm{i}}\,}z}{\sqrt{1+(-{{\mathrm{i}}\,}z^{\gamma-1})}}\,,$$ and define $$\Psi(z) = \frac{1}{\sqrt{\hat{\Psi}(z)}} \quad\text{ and }\quad
       \hat{\Psi}(z)= 1+(-{{\mathrm{i}}\,}\,z)^{\gamma-1}\;.$$ Then, with this notation, in order to prove causality of $K$, it suffices to prove that for all $z \in {{\mathbb{C}}}_0$ $$\label{positive}
\frac{\tau_0 c_0}{\alpha_0}\alpha_{ksb}(\hat{z}) = \eta \Re (\Psi(z)) + \xi \Im (\Psi(z)) \geq 0\;.$$ As in the proof of Proposition \[theo:szabo\] we show that both terms $\eta \Re (\Psi(z))$ and $\xi \Im (\Psi(z))$ are non-negative, and then from Corollary \[co:lion\] the assertion follows. 

In order to prove (\[positive\]) we note that the function $\hat{\Psi}$ here is the same as in (\[eq:delta\]) in the proof of Proposition \[theo:szabo\] when ${\tilde{\alpha}}_0$ is set to $1/2$. Thus we can already rely on the series of manipulations for $\hat{\Psi}$ developed in the proof of Proposition \[theo:szabo\]. 

-   Since $\eta \geq 0$ it suffices to show that $\Re (\Psi(z)) \geq 0$. We note that for a complex number $a+{{\mathrm{i}}\,}b$ $$\Re \left(\frac{1}{a+{{\mathrm{i}}\,}b}\right)
               = \Re \left(\frac{a- {{\mathrm{i}}\,}b}{a^2 + b^2}\right)
               = \frac{1}{a^2 + b^2} \Re(a+{{\mathrm{i}}\,}b)\;.$$ Taking into account the definition of $\Psi$ it therefore suffices to show that $\Re \left(\sqrt{\hat{\Psi}(z)}\right) \geq 0$ in ${{\mathbb{C}}}_0$. Since $\Re \left(\hat{\Psi}(z)\right) \geq 0$ in ${{\mathbb{C}}}_0$ for $\gamma\in (1,2]$, it follows that $\Re \left(\sqrt{\hat{\Psi}(z)}\right) \geq 0$ in ${{\mathbb{C}}}_0$.

-   Now, using that $$\Im \left(\frac{1}{a+{{\mathrm{i}}\,}b}\right) = \Im \left(\frac{a- {{\mathrm{i}}\,}b}{a^2 + b^2}\right) = -\frac{1}{a^2 + b^2} \Im(a+{{\mathrm{i}}\,}b)\,,$$ it suffices to show that $-\xi \Im (\sqrt{\hat{\Psi}(z)}) \geq 0$ for proving that $\xi \Im (\Psi(z)) \geq 0$. The proof is along the lines as the analogous part in Proposition \[theo:szabo\] by taking into account that here $\gamma \in (1,2)$ (in Proposition \[theo:szabo\] $\gamma \in (0,1)$).  In this case we have now that sign of $\delta$ is exactly opposite as in the proof of Proposition \[theo:szabo\], which in turn gives that $\Im (\Psi(z))$ has the opposite sign as well, and consequently $-\xi \Im (\sqrt{\hat{\Psi}(z)}) \geq 0$. Thus the assertion follows from Corollary \[co:lion\]. 

In experiments it has been discovered that several biological tissues satisfy a frequency power law (\[eq:powlaw1b\]) with exponent $\gamma \in (1,2)$ (cf.  [@Web00; @BurRoiBauPal010]).  However, as it has been shown in Theorem \[coro:powlaw1\], such models are not causal.  Our proposed model approximates the frequency power law for small frequencies, which is actually the range where it has been experimentally validated.  So, our proposed model, is valid in the actual range of experimentally measured data and extrapolates the measured data in a causal way.  Figure \[fig:comp\] shows a comparison of $\alpha_{pl}$ and $\alpha_{ksb}$ in an experimental frequency range. 

![ For $\gamma=1.5$: Comparison $\omega \to \alpha_{ksb}(\omega)$ (as defined in (\[eq:powlaw2\])) where $\alpha_0 := 2c_0\tau_0/{\left|\cos (\frac{\pi}{2} \gamma)\right|}$ (dashed line) and the power law $\alpha_{pl}(\omega) = {\left|\tau_0\,\omega\right|}^\gamma$ (as defined in (\[eq:powlaw1b\])).  For liquids: $\tau_0=10^{-6}\,MHz$ (left picture) and for gases: $\tau_0=10^{-4}\,MHz$ (right picture) (cf.  [@KinFreCopSan00]).  Experiments for determining the power law coefficient are performed in the range $0-60\,MHz$ (cf.  e.g. [@Sza95]), which is the basis for the range of the represented data. []{ data-label="fig:comp"}](compL.eps "fig:"){width="45.00000%"} ![ For $\gamma=1.5$: Comparison $\omega \to \alpha_{ksb}(\omega)$ (as defined in (\[eq:powlaw2\])) where $\alpha_0 := 2c_0\tau_0/{\left|\cos (\frac{\pi}{2} \gamma)\right|}$ (dashed line) and the power law $\alpha_{pl}(\omega) = {\left|\tau_0\,\omega\right|}^\gamma$ (as defined in (\[eq:powlaw1b\])).  For liquids: $\tau_0=10^{-6}\,MHz$ (left picture) and for gases: $\tau_0=10^{-4}\,MHz$ (right picture) (cf.  [@KinFreCopSan00]).  Experiments for determining the power law coefficient are performed in the range $0-60\,MHz$ (cf.  e.g. [@Sza95]), which is the basis for the range of the represented data. []{ data-label="fig:comp"}](compR.eps "fig:"){width="45.00000%"}

Model of Greenleaf and Patch {#model-of-greenleaf-and-patch .unnumbered}
----------------------------

\[theo:gp\] For $\alpha_0>0$ let $\alpha_{gp}^*$ be defined as in (\[alpha\*Patch\]) with the specified values $\gamma\in{\{ 1,2 \}}$. Then $K$, as defined in (\[31\]), is not causal. 

For the two specified models we have $\alpha_{gp1}(\omega)= a_0 {\left|\omega\right|} >0$ and $\alpha_{gp2}(\omega) = a_0 \omega^2$. The respective holomorphic extensions are given by (cf.  Proof of Theorem \[thcaus02\] and Theorem \[coro:powlaw1\]) $$\begin{aligned}
z \in {{\mathbb{C}}}_0 \mapsto \hat{\alpha}_{gp1}^*(z) := a_0 z
        + {{\mathrm{i}}\,}\frac{2 a_0}{\pi} z \log \left(\frac{z}{\omega_0}\right) \qquad (\omega_0\neq0, \mbox{fixed})
\end{aligned}$$ and $$\begin{aligned}
 z \in {{\mathbb{C}}}_0 \mapsto \hat{\alpha}_{gp2}^*(z) := a_0\,(-{{\mathrm{i}}\,}\,z)^2  \,. 
\end{aligned}$$ The assertion for $\gamma=1$ follows as in the proof of Theorem \[thcaus02\] and the second assertion follows from Theorem \[coro:powlaw1\] for $\gamma=2$.

Model of Chen and Holm {#model-of-chen-and-holm .unnumbered}
----------------------

\[th:chenholm\] Let $0<\alpha_1<1$, $\gamma\in (0,2)$ and ${\mathcal{G}}$ as in (\[GreenChenHolm\]).  Then there does not exist a constant $c>0$ such that $$\label{CausRef}
  {\mbox{supp}}( {\mathcal{G}}\left(\cdot,t \right) ) \subseteq B_{c\,t}({{\mathbf{0}}})
              \qquad \mbox{ for } \qquad t>0\,,$$ i.e. for each $c>0$ the function $t\mapsto {\mathcal{G}}({{\vec{x}}},t+{\left|{{\vec{x}}}\right|}/c)$ is not causal. 

Let $t>0$ be fixed.  Assume that ${{\vec{x}}}\mapsto {\mathcal{G}}({{\vec{x}}},t)$ has support in $ B_{c\,t}({{\mathbf{0}}})$ for some $c>0$. Then according to the Paley-Wiener-Schwartz Theorem (Cf.  [@GasWit99; @Hoe03]) the map ${\mathbf{k}}\mapsto {\mathcal{F}^{-1} \left\{ {\mathcal{G}}\right\}}({\mathbf{k}},t)$ is infinitely differentiable.  We show that this is not possible.  According to (\[GreenChenHolm\]) and (\[AkBk\]), we have $$\begin{aligned}
  {\mathcal{F}}_{3D}^{-1}\left\{ {\mathcal{G}}\right\} ({\mathbf{k}},t)
    = \frac{H(t)\,c_0^2}{(2\,\pi)^{3/2}}\,{\exp \left( A({\mathbf{k}})\,t \right)}\, \frac{\sin(B({\mathbf{k}})\,t)}{B({\mathbf{k}})}
\end{aligned}$$ with $$A({\mathbf{k}}) := -\alpha_1\,c_0\,{\left|{\mathbf{k}}\right|}^\gamma\,,
\qquad
    B({\mathbf{k}}) := c_0\,\sqrt{{\left|{\mathbf{k}}\right|}^2-\alpha_1^2\,{\left|{\mathbf{k}}\right|}^{2\,\gamma}}\,.$$ Since $\gamma \in (0,2)$, the function ${\mathbf{k}}\mapsto {\exp \left( A({\mathbf{k}})\,t \right)}$ is not infinitely often differentiable at ${\mathbf{k}}={{\mathbf{0}}}$ and since the holomorphic function $\frac{\sin(B({\mathbf{k}})\,t)}{B({\mathbf{k}})}$ does not vanish at ${\mathbf{k}}={{\mathbf{0}}}$, it follows that ${\mathbf{k}}\mapsto {\mathcal{F}}_{3D}^{-1}\{{\mathcal{G}}\}({\mathbf{k}},t)$ is not infinitely often differentiable at ${\mathbf{k}}={{\mathbf{0}}}$. Consequently, ${{\vec{x}}}\mapsto {\mathcal{G}}({{\vec{x}}},t)$ cannot have compact support, which concludes the proof. 

Integro-Differential Equations Describing Attenuation {#sec:integro}
=====================================================

In the following we derive the integro-differential equations for the attenuated pressure ${p_{\text{att}}}$ for various attenuation laws.  Thereby, we first derive equations which the according attenuated Green functions ${\mathcal{G}}$ (cf. ( \[eq:FGG0\])) are satisfying, and then, by convolution, we derive the equations for ${p_{\text{att}}}$. The integro-differential equations are general in the sense, that they apply to arbitrary source terms $f$, and in particular to the source term $f$ (\[eq:source\]) of the forward problem of photoacoustic imaging with attenuated waves. 

For this purpose, we rewrite $\nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}\right\}}$ by using its definition (\[eq:FGG0\]), i.e. ${\mathcal{G}}= K *_t {\mathcal{G}}_0 $, and the product differentiation rule, which gives $$\label{eq:laplG}
\begin{aligned}
&\frac{1}{\sqrt{2\,\pi}}\,\nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}\right\}}\\
= & \nabla^2 {\mathcal{F} \left\{ K \right\}} \cdot {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}} +
2 \nabla {\mathcal{F} \left\{ K \right\}} \cdot \nabla {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}}
+ {\mathcal{F} \left\{ K \right\}} \cdot \nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}}. 
\end{aligned}$$ To evaluate the expression on the right hand side, we calculate $\nabla {\mathcal{F} \left\{ K \right\}}$ and $\nabla^2 {\mathcal{F} \left\{ K \right\}}$. From (\[eq:kernel\]), it follows that $$\label{eq:defg2b}
\begin{aligned}
  \nabla {\mathcal{F} \left\{ K \right\}} = -\beta^*{}' \cdot {\mathcal{F} \left\{ K \right\}} \cdot {\mbox{sgn}}\,,
\end{aligned}$$ where $\beta^*{}'$ denotes the derivative of $\beta^*(r,\omega)$ (cf. ( \[eq:kernel\])) with respect to $r$. This together with the formula (\[eq:der\_sgn\]) in the Appendix implies that $$\label{eq:defg2c}
\begin{aligned}
& \nabla^2 {\mathcal{F} \left\{ K \right\}}\\
         =&- \nabla \cdot \left( \beta^*{}' \cdot {\mathcal{F} \left\{ K \right\}} \cdot {\mbox{sgn}}\right)\\
         =& - (\nabla \cdot {\mbox{sgn}}) \cdot \beta^*{}' \cdot {\mathcal{F} \left\{ K \right\}}
            - ({\mbox{sgn}}\cdot \nabla \beta^*{}') \cdot {\mathcal{F} \left\{ K \right\}}
            - ({\mbox{sgn}}\cdot \nabla {\mathcal{F} \left\{ K \right\}}) \cdot \beta^*{}' \\
         =& \left[ -\frac{2}{{\left|{{\vec{x}}}\right|}} \cdot\beta^*{}'
                   - \beta^*{}''
                   +  \left(\beta^*{}'\right)^2\right] \cdot {\mathcal{F} \left\{ K \right\}}. 
\end{aligned}$$ Inserting [(\[eq:defg2b\])]{} and [(\[eq:defg2c\])]{} into [(\[eq:laplG\])]{} and using again the identity ${\mathcal{G}}= K *_t {\mathcal{G}}_0$, shows that $$\label{eq:laplG2}
 \begin{aligned}
  & \frac{1}{\sqrt{2\,\pi}}\,\nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}\right\}} \\
  &\qquad = \frac{1}{\sqrt{2\,\pi}}\,\left[ -\frac{2}{{\left|{{\vec{x}}}\right|}} \cdot \beta^*{}' - \beta^*{}''
      +  \left(\beta^*{}'\right)^2\right]\cdot{\mathcal{F} \left\{ {\mathcal{G}}\right\}}\\
   &\qquad\quad -2 \beta^*{}' \cdot {\mathcal{F} \left\{ K \right\}} \cdot ({\mbox{sgn}}\cdot \nabla {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}})
      + {\mathcal{F} \left\{ K \right\}} \cdot \nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}}. 
 \end{aligned}$$ From this identity, together with the two following properties of ${\mathcal{G}}_0$, $$\label{eq:g1a}
\nabla{\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}}
= \left[\frac{{{\mathrm{i}}\,}\omega}{c_0} -\frac{1}{{\left|{{\vec{x}}}\right|}}\right] \cdot {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}} \cdot {\mbox{sgn}}\,,$$ and $$\label{eq:g1b}
  \nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}} + \frac{\omega^2}{c_0^2}{\mathcal{F} \left\{ {\mathcal{G}}_0 \right\}}
      = -\frac{1}{\sqrt{2\pi}}\delta_{{\vec{x}}}\,,$$ it follows that $$\label{eq:laplG3}
 \begin{aligned}
   \nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}\right\}}
      =& \left[ -\frac{2}{{\left|{{\vec{x}}}\right|}} \cdot \beta^*{}'
                   - \beta^*{}''
                   +  \left(\beta^*{}'\right)^2\right] \cdot {\mathcal{F} \left\{ {\mathcal{G}}\right\}}\\
       & - 2 \left[\frac{{{\mathrm{i}}\,}\omega}{c_0} -\frac{1}{{\left|{{\vec{x}}}\right|}}\right] \cdot
                   \beta^*{}' \cdot {\mathcal{F} \left\{ {\mathcal{G}}\right\}}\\
       &-\frac{\omega^2}{c_0^2} \cdot {\mathcal{F} \left\{ {\mathcal{G}}\right\}}\\
       &-{\mathcal{F} \left\{ K \right\}} \cdot \delta_{{\vec{x}}}\;. 
 \end{aligned}$$ Inserting the identity ${\mathcal{F} \left\{ K \right\}}({{\vec{x}}},\omega) \cdot \delta_{{\vec{x}}}= {\mathcal{F} \left\{ K \right\}}({\bf 0},\omega) \cdot \delta_{{\vec{x}}}$ in [(\[eq:laplG3\])]{} gives the *Helmholtz equation* $$\label{eq:helmholtz}
\begin{aligned}
~& \nabla^2 {\mathcal{F} \left\{ {\mathcal{G}}\right\}}
       -\left[\beta^*{}' + \frac{(-{{\mathrm{i}}\,}\, \omega)}{c_0}\right]^2 \cdot{\mathcal{F} \left\{ {\mathcal{G}}\right\}}\\
 =& -\beta^*{}'' \cdot {\mathcal{F} \left\{ {\mathcal{G}}\right\}}
           - {\mathcal{F} \left\{ K \right\}}({\bf 0},\cdot)\cdot \delta_{{\vec{x}}}\\
 =& -\beta^*{}'' \cdot {\mathcal{F} \left\{ {\mathcal{G}}\right\}}
           -\frac{1}{\sqrt{2\pi}} {\exp \left( -\beta^*({\bf 0},\omega) \right)} \cdot \delta_{{\vec{x}}}\,. 
\end{aligned}$$ To reformulate [(\[eq:helmholtz\])]{} in space–time coordinates, we introduce two convolution operators: $$\label{eq:defD*D*'}
   D_*f := K_* *_t f {\qquad \text{and} \qquad}
   D_*'f :=K_*' *_t f ,$$ where the kernels $K_*$ and $K_*'$ are given by $$\label{eq:defKK'}
   K_*:=K_*({{\vec{x}}},t) := K_*({\left|{{\vec{x}}}\right|},t) {\quad \text{and} \quad}
   K_*(r,t)  := \frac{1}{\sqrt{2\pi}} {\mathcal{F}^{-1} \left\{ \beta^*{}' \right\}}(r,t)$$ and $$\label{eq:defKK'2}
  K_*':=K_*'({{\vec{x}}},t) := K_*'({\left|{{\vec{x}}}\right|},t)
{\quad \text{and} \quad}
    K_*'(r,t)  = \frac{1}{\sqrt{2\pi}} {\mathcal{F}^{-1} \left\{ \beta^*{}'' \right\}}(r,t).$$ Using these operators and applying the inverse Fourier transform to [(\[eq:helmholtz\])]{} gives $$\label{eq:waveeq+}
  \nabla^2 {\mathcal{G}}-\left[D_*
         + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2 {\mathcal{G}}=  -D_*' {\mathcal{G}}- K({\bf 0},\cdot) \delta_{{\vec{x}}}\;.$$ In the case that $\beta^*({\left|{{\vec{x}}}\right|},\omega) = \alpha^*(\omega){\left|{{\vec{x}}}\right|}$ is of standard form (\[eq:attenuation\_coefficient\]), it follows that $$\label{eq:defKK'+}
\begin{aligned}
K_*(t) = \frac{1}{\sqrt{2\pi}} {\mathcal{F}^{-1} \left\{ \alpha^* \right\}}(t) \qquad\text{ and } \qquad
K_*'\equiv 0\;. 
\end{aligned}$$

For a general source term $f$, we denote the attenuated wave by ${p_{\text{att}}}$. That is $${p_{\text{att}}}:= {p_{\text{att}}}({{\vec{x}}},t) = {\mathcal{G}}*_{{{\vec{x}}},t} f =: {\mathcal{A}}f\,,$$ where ${\mathcal{A}}$ is the convolution operator according to the Green function ${\mathcal{G}}$. This then shows that ${p_{\text{att}}}$ satisfies the integro-differential equation $$\label{eq:waveeq+2}
\boxed{
\nabla^2 {p_{\text{att}}}-\frac{1}{c_0^2} \frac{\partial^2 {p_{\text{att}}}}{\partial t^2} =  - {\mathcal{A}}_s f\;,}$$ where ${\mathcal{A}}_s$ denotes the space–time convolution operator with kernel $$\label{eq:defKs}
\begin{aligned}
   K_s:=K_s({{\vec{x}}},t) :=  -({\mathcal{B}}{\mathcal{G}})({{\vec{x}}},t)
           +(D_*' {\mathcal{G}})({{\vec{x}}},t)
          + K({\bf 0},t) \cdot \delta_{{\vec{x}}}({{\vec{x}}})
\end{aligned}$$ and $$\label{eq:B}
    {\mathcal{B}}:= D_*^2 + \frac{2}{c_0} D_*\frac{\partial}{\partial t}\;.$$ Equation [(\[eq:waveeq+2\])]{} is called *pressure wave equation with attenuation coefficient* $\beta^*$. We emphasize that $\beta^*$ determines the operators $D_*$ and $D'_*$ which in turn determine the operator ${\mathcal{A}}_s$, which in turn determines ${p_{\text{att}}}$ - this reveals the dependence of ${p_{\text{att}}}$ from $\beta_*$.

\[rema:alpha\*2\] Let $\beta^*(r,\omega)=\alpha^*(\omega)r$ be the standard attenuation model (cf. ( \[eq:attenuation\_law\])).  Assuming that the associated kernel $K$ (cf. ( \[31\])) is causal, it follows that $${\left|\nabla K\right|} = \frac{1}{\sqrt{2\pi}} {\left|{\mathcal{F}^{-1} \left\{ \alpha^* \cdot {\exp \left( -\alpha^* {\left|{{\vec{x}}}\right|} \right)} \right\}}\right|}\;.$$ Using some sequence ${\{ {{\vec{x}}}_n \}}$ satisfying ${{\vec{x}}}_n \neq {\bf 0}$ and ${{\vec{x}}}_n \to {\bf 0}$ shows that $$\lim_{n \to \infty} {\left|\nabla K\right|}({{\vec{x}}}_n,t) = \frac{1}{\sqrt{2\pi}} {\left|{\mathcal{F}^{-1} \left\{ \alpha^* \right\}}(t)\right|} \underbrace{=}_{(\ref{eq:defKK'+})} {\left|K_*(t)\right|}\;.$$ Due to the causality of $K$ the left hand side is zero for $t<0$, and thus $K_*$ is also causal. 

Because the convolution of causal distributions is well-defined, the operator $D_*$ is well-defined on all causal distributions.  Moreover, since $K_*'=0$, it follows that $D_*'\equiv 0$. Using that $K_*$ depends only on $t$ it follows that $$(D_* {\mathcal{G}}) *_{{{\vec{x}}},t} f =  [K_* *_t {\mathcal{G}}] *_{{{\vec{x}}},t} f
       =  K_* *_t [{\mathcal{G}}*_{{{\vec{x}}},t} f] = D_*({\mathcal{G}}*_{{{\vec{x}}},t} f).$$ Convolving each term in [(\[eq:waveeq+\])]{} with a function $f$, using the previous identity and that $D_*' \equiv 0$, it follows that $$\label{eq:waveeq+3}
  \nabla^2 {p_{\text{att}}}-\left[D_* + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2  {p_{\text{att}}}= - f \;$$ where $$\label{eq:D*sp}
   D_* \cdot = \frac{1}{\sqrt{2\,\pi}}\,{\mathcal{F} \left\{ \alpha^* \right\}}(t) *_t \cdot \,.$$

In the following we derive the common forms of the wave equation models corresponding to the various attenuation models listed in Section \[sec:attenuation\]. 

Power Laws {#power-laws-1 .unnumbered}
----------

-   Let $0 < \gamma \not\in {{\mathbb{N}}}$ and $0 < \alpha_0$. We note that the *Riemann-Liouville fractional derivative* with respect to time, denote by $D_t^\gamma$ (see [@KilSriTru06; @Pod99]), is defined in the Fourier domain by $$\label{eq:defDtga}
            {\mathcal{F} \left\{ D_t^{\gamma}f \right\}} = (-{{\mathrm{i}}\,}\omega)^\gamma {\mathcal{F} \left\{ f \right\}}\,,$$ and satisfies $$\label{eq:propDtgamma}
        D_t^{2\gamma}f = D_t^\gamma D_t^{\gamma}f {\qquad \text{and} \qquad}
        \frac{\partial}{\partial t} D_t^{\gamma}f = D_t^\gamma \frac{\partial}{\partial t} f = D_t^{\gamma+1}f.$$ From this together with (\[eq:powlaw1\]) and (\[eq:D\*sp\]), we infer $$\label{deftildealpha0}
            D_* = \tilde\alpha_0\, D_t^\gamma  \qquad \mbox{ with }\qquad
       \tilde \alpha_0 := \frac{\alpha_0}{\cos(\pi\,\gamma/2)}$$ and thus wave equation (\[eq:waveeq+3\]) reads as follows $$\label{eq:standwaveeq}
    \boxed{
    \begin{aligned}
      \nabla^2 {p_{\text{att}}}&-\left[\tilde\alpha_0 \,D_t^\gamma
             + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2  {p_{\text{att}}}= -f  \;. 
    \end{aligned}}$$

-   Let $\gamma=1$, then for the frequency power law (\[alphagamma1\]) it follows from the Fourier transform table I in [@Lig64] $$\begin{aligned}
             D_* &= -\frac{4\,\alpha_0}{2\,\pi} \left[\,
                        \frac{H(t)}{t^2} - (\log|\omega_0|)\,\delta'_t
                          \right] *_t \\
                 &= -\frac{4\,\alpha_0}{2\,\pi}\,  \frac{H(t)}{t^2} *_t \;\;
                    + \frac{4\,\alpha_0}{\sqrt{2\,\pi}}\,(\log|\omega_0|)\,\frac{\partial }{\partial t}\,. 
    \end{aligned}$$

Szabo’s Attenuation Law: {#szabos-attenuation-law .unnumbered}
------------------------

Let $0 < \alpha_0$ and $0 < \gamma\not\in{{\mathbb{N}}}$. From (\[alpha\*szabo\]) and (\[eq:defDtga\]), we get $$\left[D_* + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2
            = \frac{1}{c_0^2}\frac{\partial^2 }{\partial t^2}
       + \frac{2\,\tilde\alpha_0}{c_0}\,\frac{\partial }{\partial t}\,D_t^\gamma$$ and thus wave equation (\[eq:waveeq+3\]) reads as follows $$\label{szaboseq}
\begin{aligned}
\boxed{
    \nabla^2 {p_{\text{att}}}-\frac{1}{c_0^2}\frac{\partial^2 {p_{\text{att}}}}{\partial t^2}
       - \frac{2\,\tilde\alpha_0}{c_0}\,\frac{\partial }{\partial t}\,D_t^\gamma {p_{\text{att}}}= -f({{\vec{x}}},t)\,.} 
\end{aligned}$$

Thermo-Viscous Attenuation Law: {#thermo-viscous-attenuation-law-1 .unnumbered}
-------------------------------

From (\[alpha\*th\]) we get $$\left(\mbox{Id} + \tau_0\,\frac{\partial}{\partial t}\right)\,
          \left[D_* + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2
            = \frac{1}{c_0^2} \frac{\partial^2}{\partial t^2}$$ and thus (\[eq:waveeq+3\]) becomes $$\label{thviscwaveeq}
\begin{aligned}
\boxed{
    \left(\mbox{Id} + \tau_0\,\frac{\partial}{\partial t}\right)\,\nabla^2 {p_{\text{att}}}- \frac{1}{c_0^2}\frac{\partial^2 {p_{\text{att}}}}{\partial t^2}
            = -\left(\mbox{Id} + \tau_0\,\frac{\partial}{\partial t}\right)\,f\,.} 
\end{aligned}$$ This equation is called the *thermo-viscous wave equation*. 

Nachman, Smith and Waag [@NacSmiWaa90]: {#nachman-smith-and-waag .unnumbered}
---------------------------------------

We carry out the details only for one relaxation process. 

-   $N=1$: From (\[alpha\*Nachman+\]) we get $$\left(\alpha^*(\omega) + \frac{(-{{\mathrm{i}}\,}\,\omega)}{c_0} \right)^2
              = \frac{(-{{\mathrm{i}}\,}\,\omega)^2}{\tilde c_0^2}\,\frac{1-{{\mathrm{i}}\,}\,\tilde \tau_1\,\omega}{1-{{\mathrm{i}}\,}\,\tau_1\,\omega}$$ which implies $$\left(\mbox{Id} + \tau_1\,\frac{\partial}{\partial t}\right)\,
              \left[D_* + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2
                = \left(\mbox{Id} + \tilde\tau_1\,\frac{\partial}{\partial t}\right)\,
                            \frac{1}{\tilde c_0^2} \frac{\partial^2}{\partial t^2}\,.$$ Thus (\[eq:waveeq+3\]) reads as follows $$\label{Nachmanwaveeq}
    \boxed{
    \begin{aligned}
        \left(\mbox{Id} + \tau_1\,\frac{\partial}{\partial t}\right)\,\nabla^2 {p_{\text{att}}}&
                - \frac{1}{\tilde c_0^2}\left(\mbox{Id} + \tilde\tau_1\,\frac{\partial}{\partial t}\right)\,
                          \frac{\partial^2 {p_{\text{att}}}}{\partial t^2}\\
                &= -\left(\mbox{Id} + \tau_1\,\frac{\partial}{\partial t}\right)\,f\,. 
    \end{aligned}}$$ If the term with $\tilde\tau_1=0$ is dropped and $\tilde c_0$ is replaced by $c_0$, then we obtain the thermo-viscous wave equation (\[thviscwaveeq\]). 

-   $N>1$: For the general case we refer to equation (26) in [@NacSmiWaa90]. 

Greenleaf and Patch [@PatGre06]: {#greenleaf-and-patch .unnumbered}
--------------------------------

-   For $\gamma=2$ the attenuation coefficient equals to $$\alpha^*(\omega) = \alpha_0\,\omega^2 = \tilde\alpha_0\,(-{{\mathrm{i}}\,}\omega)^2\,,$$ where $\tilde\alpha_0$ is defined as in (\[deftildealpha0\]) and thus $$D_*  = -\alpha_0\,\frac{\partial^2 }{\partial t^2} = -\alpha_0\,D_t^ 2\,,$$ which gives wave equation (\[eq:standwaveeq\]) with $\gamma=2$.

-   For $\gamma=1$ we have $$\alpha^*(\omega) = \alpha_0\,|\omega| = \alpha_0\,(-{{\mathrm{i}}\,}\omega)\,{{\mathrm{i}}\,}{\mbox{sgn}}(\omega)$$ and thus $$D_* = \alpha_0\,\mathbf{D}^{-1}
            = -\alpha_0\,\frac{\partial}{\partial t} {\mathcal{H}}\,,$$ where $\mathbf{D}^{-1}$ and ${\mathcal{H}}$ denote the *Riesz fractional differentiation operator* and the Hilbert transform (cf.  Appendix), respectively.  Therefore the wave equation reads as follows $$\boxed{
    \begin{aligned}
      \nabla^2 {p_{\text{att}}}&-\left[\alpha_0\,\mathbf{D}^{-1}
             + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2  {p_{\text{att}}}= -f  \;. 
    \end{aligned}}$$

Chen and Holm [@CheHolm04]: {#chen-and-holm .unnumbered}
---------------------------

Let $\gamma \in (0,2)$. The Green function defined by (\[GreenChenHolm\]) satisfies the Helmholtz equation $$\label{Holmeq02}
\begin{aligned}
~ &\frac{\partial^2 {\mathcal{F} \left\{ {\mathcal{G}}\right\}}}{\partial t^2}({\mathbf{k}},t)
   + 2\,\alpha_1\,c_0\,{\left|{\mathbf{k}}\right|}^\gamma \,\frac{\partial {\mathcal{F} \left\{ {\mathcal{G}}\right\}}}{\partial t}({\mathbf{k}},t)
   +  c_0^2\,{\left|{\mathbf{k}}\right|}^2 \,{\mathcal{F} \left\{ {\mathcal{G}}\right\}}({\mathbf{k}},t)\\
= & \frac{c_0^2}{(2\,\pi)^{3/2}} \,\delta(t)
\end{aligned}$$ for $t\in{{\mathbb{R}}}$ and ${\mathbf{k}}\in{{\mathbb{R}}}^3$. Since the fractional Laplacian for a rotational symmetric function $f$ and $\gamma\in (0,2)$ is defined by (cf.  Definition (2.10.1) in ([@KilSriTru06]) $$\begin{aligned}
   & \left(-\nabla^2\right)^{\gamma/2} f({{\vec{x}}})\\
:= &\frac{1}{\sqrt{(2\,\pi)^3}}\,\int_{{{\mathbb{R}}}^3} {\exp \left( {{\vec{x}}}\cdot{\mathbf{k}}\right)}\,
                          \left[ {\left|{\mathbf{k}}\right|}^\gamma\,
     \frac{1}{\sqrt{(2\,\pi)^3}}\,\int_{{{\mathbb{R}}}^3} {\exp \left( -{{\vec{x}}}\cdot{\mathbf{k}}\right)}\, f({{\vec{x}}})\,{\mbox{d}}{{\vec{x}}}\right]
             \,{\mbox{d}}{\mathbf{k}}\\
= & {\mathcal{F}}_{3D}\left\{ {\left|{\mathbf{k}}\right|}^\gamma\,{\mathcal{F}}_{3D}^{-1}\{ f\}({\mathbf{k}}) \right\}({{\vec{x}}})\,,
\end{aligned}$$ we obtain the following wave equation for ${p_{\text{att}}}:={\mathcal{G}}\ast_{{{\vec{x}}},t}f$ $$\label{Holmeq01}
\boxed{
\begin{aligned}
    \nabla^2 {p_{\text{att}}}-\frac{1}{c_0^2}\frac{\partial^2 {p_{\text{att}}}}{\partial t^2}
       - \frac{2\,\alpha_1}{c_0}\,\frac{\partial }{\partial t}\,\left(-\nabla^2\right)^{\gamma/2} {p_{\text{att}}}= -f({{\vec{x}}},t)\;. 
\end{aligned}}$$ We note that Chen and Holm used instead of $2\,\alpha_1/c_0$ the term $2\,\alpha_1/c_0^{1-\gamma}$.

Our Model [@KowSchBon10]: {#our-model-1 .unnumbered}
-------------------------

From (\[eq:powlaw2\]) we get $$\left(\mbox{Id} + \tau_0^{\gamma-1}\,D_t^{\gamma-1}\right)\,
          \left[D_* + \frac{1}{c_0} \frac{\partial}{\partial t}\right]^2
            = \frac{1}{c_0^2} \frac{\partial^2}{\partial t^2} \,
            \left(\alpha_0\,\mbox{Id} + L^{1/2}\right)^2$$ where the time convolution operator $L^{1/2}$ is the convolution operator with kernel $$l(t):=L^{1/2} (\delta_t)
   = \frac{1}{\sqrt{2\,\pi}}\, {\mathcal{F}}\left\{ \sqrt{1+(-{{\mathrm{i}}\,}\tau_0\,\omega)^{\gamma-1}} \right\}\;.$$ Consequently, $L:=(L^{1/2})^2=\mbox{Id} + \tau_0\,D_t^{\gamma-1}$ and (\[eq:waveeq+3\]) can be rewritten as follows $$\label{causthviscwaveeq}
\boxed{
\begin{aligned}
    \left(\mbox{Id} + \tau_0^{\gamma-1}\,D_t^{\gamma-1}\right)\,\nabla^2 {p_{\text{att}}}&
            - \frac{1}{c_0^2}\,\left(\alpha_0\,\mbox{Id} + L^{1/2}\right)^2
                       \,\frac{\partial^2 {p_{\text{att}}}}{\partial t^2}\\
            &= -\left(\mbox{Id} + \tau_0^{\gamma-1}\,D_t^{\gamma-1}\right)\,f\,. 
\end{aligned}}$$

Pressure Relation
=================

In this section we derive the relation between ${p_{\text{att}}}$ and $p_0$ when the source term is of the form (\[eq:source\]).  This chapter is a special instance of Section \[sec:integro\].  However, utilizing the special structure of the source term different formulas can be derived. 

Attenuation is defined as a multiplicative law (in the frequency domain) relating the amplitudes of an attenuated and an unattenuated wave initialized by a delta impulse.  Here we are concerned in deriving the convolution relation between the solution $p_0$ of [(\[eq:ex:wave3d\])]{} (or equivalently of  [(\[eq:ex:ivp3d\])]{} and [(\[eq:init\_values\])]{}) and the attenuated wave function ${p_{\text{att}}}$, which, according to (\[eq:FGG0\]) and (\[eq:a\_conv\]), is given by $$\label{eq:defpatt}
 {p_{\text{att}}}= {\mathcal{G}}*_{{{\vec{x}}},t} f  = (K *_t {\mathcal{G}}_0) *_{{{\vec{x}}},t} f\,,$$ with $f$ from (\[eq:source\]).  Using (\[eq:g0\]) and the rotational symmetry of $K$, it follows that $$\begin{aligned}
  &  \qquad  \qquad \qquad \qquad  \qquad \qquad{p_{\text{att}}}({{\vec{x}}},t) \\
       &= \int_{{\mathbb{R}}}\int_{{{\mathbb{R}}}^3} (K \ast_{t}{\mathcal{G}}_0)({{\vec{x}}}-{{\vec{x}}}',t-t'')
                    \rho({{\vec{x}}}')\,d{{\vec{x}}}' \frac{\partial \delta_t}{\partial t}(t'')\,dt''\\
       &= \int_{{\mathbb{R}}}\int_{{{\mathbb{R}}}^3} \int_{{\mathbb{R}}}K({{\vec{x}}}-{{\vec{x}}}',t-t'-t''){\mathcal{G}}_0({{\vec{x}}}-{{\vec{x}}}',t')\,dt'
                    \rho({{\vec{x}}}')\,d{{\vec{x}}}' \frac{\partial \delta_t}{\partial t}(t'')\,dt''\\
       &= \int_{{\mathbb{R}}}\int_{{{\mathbb{R}}}^3} \int_{{\mathbb{R}}}\frac{\partial }{\partial t}K({{\vec{x}}}-{{\vec{x}}}',t-t'-t''){\mathcal{G}}_0({{\vec{x}}}-{{\vec{x}}}',t')
                    \rho({{\vec{x}}}') \delta_t (t'')\,dt'\,d{{\vec{x}}}' \,dt''\\
    &= \int_{{{\mathbb{R}}}^3} \int_{{\mathbb{R}}}\frac{\partial }{\partial t} K({\left|{{\vec{x}}}-{{\vec{x}}}'\right|},t-t')\,
    \frac{\delta_t(t'-{\left|{{\vec{x}}}-{{\vec{x}}}'\right|}/c_0)}{4\pi {\left|{{\vec{x}}}-{{\vec{x}}}'\right|}} \,\rho({{\vec{x}}}')\,{\mbox{d}}t'\,{\mbox{d}}{{\vec{x}}}'\\
    &= \int_{{{\mathbb{R}}}^3}\frac{\partial }{\partial t} K({\left|{{\vec{x}}}-{{\vec{x}}}'\right|},t-{\left|{{\vec{x}}}-{{\vec{x}}}'\right|}/c_0)\,\frac{\rho({{\vec{x}}}')}{4\pi {\left|{{\vec{x}}}-{{\vec{x}}}'\right|}} \,{\mbox{d}}{{\vec{x}}}'\,.\\
\end{aligned}$$

Using the representation ${{\vec{x}}}'-{{\vec{x}}}= r' \vec{s}$ with $r' \geq 0$ and $\vec{s} \in S^2$, it follows that $$\label{eq:pattrep01}
{p_{\text{att}}}({{\vec{x}}},t) = \frac{1}{4 \pi}\, \int_0^\infty \frac{\partial }{\partial t} K(r',t-r'/c_0) \,r'\,
                                \int_{S^2}  \rho({{\vec{x}}}+r'\vec{s})\,d\vec{s}\,dr'\;.$$ Moreover, $$\begin{aligned}
p_0({{\vec{x}}},t)
&= \frac{\partial }{\partial t} \int_{{{\mathbb{R}}}^3} \frac{\delta_t(t-{\left|{{\vec{x}}}-{{\vec{x}}}'\right|}/c_0)}{4\pi {\left|{{\vec{x}}}-{{\vec{x}}}'\right|}} \rho({{\vec{x}}}')\,{\mbox{d}}{{\vec{x}}}'\\
&= \frac{\partial }{\partial t} \int_0^\infty r'^2 \int_{S^2}  \frac{\delta_t(t-r'/c_0)}{4\pi r'} \rho({{\vec{x}}}+ r' \vec{s})\,d\vec{s}\,dr'\\
&= \frac{\partial }{\partial t} \int_0^\infty r' \frac{\delta_t(t-r'/c_0)}{4\pi} \int_{S^2} \rho({{\vec{x}}}+ r' \vec{s})\,d\vec{s}\,dr'\\
&= \frac{\partial }{\partial t} \int_0^\infty \frac{c_0^2 r''}{4\pi} \delta_t(t-r'') \int_{S^2} \rho({{\vec{x}}}+ c_0 r'' \vec{s})\,d\vec{s}\,d r''\\
&= \frac{\partial }{\partial t} \left(\frac{c_0^2 t}{4\pi} \int_{S^2} \rho({{\vec{x}}}+ (c_0t) \vec{s})\,d\vec{s} \right)\;. 
\end{aligned}$$ This gives $$\label{eq:summ}
\begin{aligned}
  r'\, \int_{S^2}  \rho({{\vec{x}}}+r'\vec{s})\,d \vec{s}
  =  \frac{4\,\pi}{c_0} \int_0^{r'/c_0} p_0({{\vec{x}}},t')\,{\mbox{d}}t'\;. 
\end{aligned}$$ Now, denoting $$F(t,r') := \int_0^{r'}  \frac{\partial}{\partial t} K(r'',t-r''/c_0) \,dr''\,,\qquad
G(r') := \int_0^{r'/c_0} p_0({{\vec{x}}},r'')\,dr''\,$$ and $$F(t,\infty) = \lim_{r' \to \infty} F(t,r')$$ it follows from (\[eq:pattrep01\]) and (\[eq:summ\]) and the fact that $p_0({{\vec{x}}},0) = 0$ that $$\label{eq:pattrep02a}
\begin{aligned}
{p_{\text{att}}}({{\vec{x}}},t) =& \frac{1}{c_0} \int_0^\infty F'(t,r') G(r')\,dr'\\
=& - \frac{1}{c_0} \int_0^\infty F(t,r') \underbrace{G'(r')}_{=p_0({{\vec{x}}},r'/c_0)/c_0}\,dr' + \frac{1}{c_0} \left.  F(t,r') G(r') \right|_{r'=0}^\infty\\
=& \frac{1}{c_0^2} \left( F(t,\infty) \int_0^\infty p_0({{\vec{x}}},r'/c_0)\,dr'- \int_0^\infty F(t,r') p_0({{\vec{x}}},r'/c_0)\,dr'\right)\\
=& \frac{1}{c_0} \left( F(t,\infty) \int_0^\infty p_0({{\vec{x}}},t')\,dt' - \int_0^\infty  F(t,c_0 t') p_0({{\vec{x}}},t')\,d t'\right)\\
=:& \int_0^\infty {\mathcal{M}}(t,t') p_0({{\vec{x}}},t')\,dt'\,,
\end{aligned}$$ where $$\label{eq:mos}
\begin{aligned}
  {\mathcal{M}}(t,t') :&= \frac{1}{c_0} \,(F(t,\infty)-F(t,c_0 t')) \,. 
\end{aligned}$$ In the following we derive an equivalent representation of ${\mathcal{M}}$ in terms of the attenuation coefficient,under the assumption that the attenuation coefficient $\alpha^*$ is such that $K$ is causal.  From [(\[eq:kernel\])]{}, [(\[eq:attenuation\_law\])]{} and Item \[item:Fdelta\] in the Appendix, it follows that $$\begin{aligned}
  K(r',t-r'/c_0)
    = \frac{1}{\sqrt{2\,\pi}}\,  {\mathcal{F}^{-1} \left\{  {\exp \left( -\alpha^*(\omega)\,r'+{{\mathrm{i}}\,}\,\frac{\omega}{c_0}\,r' \right)}   \right\}}(t) \,. 
\end{aligned}$$ which implies $$\label{eq:pattrep02b}
\begin{aligned}
~ & F(t,c_0\,t')\\
    =& \int_0^{c_0\,t'} \frac{\partial }{\partial t} K(r',t-r'/c_0)\,{\mbox{d}}r' \\
    =& \frac{1}{\sqrt{2\,\pi}}\,  {\mathcal{F}^{-1} \left\{ \frac{-{{\mathrm{i}}\,}\,\omega\,\left[ {\exp \left( -\alpha^*(\omega)\,c_0\,t'+{{\mathrm{i}}\,}\,\omega\,t' \right)} -1 \right]}
                                                 {-\alpha^*(\omega)+{{\mathrm{i}}\,}\,\omega/c_0}\,  \right\}}(t)  \\
    =& \frac{1}{\sqrt{2\,\pi}}\, {\mathcal{F}^{-1} \left\{ \frac{-{{\mathrm{i}}\,}\,\omega\,}{-\alpha^*(\omega)+{{\mathrm{i}}\,}\,\omega/c_0}\,  \right\}}(t)  \ast_t
          \left[ K(c_0\,t',t-t')  - \delta_t(t)\right]\,. 
\end{aligned}$$ Since $K$ is causal, it satisfies $K(c_0\,t',t-t')   = 0$ for $t < t'$, and therefore $$F(t,\infty) = \lim_{t'\to \infty} F(t,c_0\,t')
          =  \frac{1}{\sqrt{2\,\pi}}\, {\mathcal{F}^{-1} \left\{ \frac{-{{\mathrm{i}}\,}\,\omega\,}{-\alpha^*(\omega)+{{\mathrm{i}}\,}\,\omega/c_0}\,  \right\}}(t)\,.$$ Hence [(\[eq:pattrep02b\])]{} can be written as follows: $$\begin{aligned}
  F(t,c_0\,t')
    &= F(t,\infty)  \ast_t
          K(c_0\,t',t-t')  - F(t,\infty)\,. 
\end{aligned}$$ and therefore [(\[eq:mos\])]{} simplifies to $$\label{eq:defM}
\begin{aligned}
   {\mathcal{M}}(t,t')
      &=   - \frac{1}{c_0}\,F(t,\infty)  \ast_t K(c_0\,t',t-t') \\
      & =\frac{1}{\sqrt{2\,\pi}}\,  {\mathcal{F}^{-1} \left\{ \frac{{{\mathrm{i}}\,}\,\omega \,{\exp \left( -\alpha^*(\omega)\,c_0\,t'+{{\mathrm{i}}\,}\,\omega\,t' \right)}}
                                                  {-\alpha^*(\omega)\,c_0+{{\mathrm{i}}\,}\,\omega}\,  \right\}}(t) \,. 
\end{aligned}$$ Note that ${\mathcal{M}}(t,0)= - \frac{1}{c_0}\,F(t,\infty)$. The following lemma shows that if $K$ is causal $$\label{eq:propM}
         {\mathcal{M}}(t,t') = 0 \qquad \mbox{ if } \qquad 0< t < t'\,$$ and therefore the upper limit of integration in the last term [(\[eq:pattrep02a\])]{} can be replaced by $t$. This means that the set of attenuated pressure values $${\{ {p_{\text{att}}}({{\vec{x}}},s) : 0\leq s\leq t \}}$$ depend only on the unattenuated pressure values $${\{ p_0({{\vec{x}}},s) : 0\leq s\leq t \}}\;.$$

\[lemm02\] Let $K$ from [(\[eq:kernel\])]{} be causal with $\beta^*(r,\omega):=\alpha^*(\omega)\,r$. Moreover, assume that $\alpha^*(\omega) \neq {{\mathrm{i}}\,}\,\omega/c_0$ for every $\omega\in{{\mathbb{R}}}$, and let ${\mathcal{M}}$ be as defined in [(\[eq:defM\])]{}.  Then,

-   the function $t \to {\mathcal{M}}(t,0)$ is causal. 

-   For every $t'>0$ $$\label{propM1}
        {\mathcal{M}}(t,t')=0 \text{ for all } t < |t'|\,.$$

<!-- -- >

-   In order to prove causality of $t \to {\mathcal{M}}(t,0)$ we verify the three assumptions of Theorem \[th:lion\] for the tempered distribution $${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,0) = \frac{1}{\sqrt{2\,\pi}}\,\frac{\omega }{k(\omega)\,c_0}
     \qquad\mbox{with}\qquad  k(\omega):= {{\mathrm{i}}\,}\,\alpha^*(\omega) + \frac{\omega}{c_0}\,.$$ Since $K$ is causal, as has been shown in Remark \[rema:alpha\*2\], also the function $$t\mapsto K_*(t)=\frac{1}{\sqrt{2\pi}} {\mathcal{F}^{-1} \left\{ \alpha^* \right\}}(t)$$ is causal.  Now, using Theorem \[th:lion\], it follows that

    1.   $\alpha_*$ is holomorphic in $\mathring{{{\mathbb{C}}}}_0$,

    2.   $\alpha^*(\xi+{{\mathrm{i}}\,}\,\eta)\to \alpha^*(\xi)$ for $\eta\to 0$ in ${\mathcal{S}}'$ and

    3.   for each $ \epsilon>0$ there exists a polynomial $P$ such that ${\left|\alpha^*(z)\right|}\leq P({\left|z\right|})$ for $z\in{{\mathbb{C}}}_\epsilon$.

    Since $\alpha_*(\omega)\neq {{\mathrm{i}}\,}\,\omega/c_0$ for all $\omega\in {{\mathbb{R}}}$ together with in [@BelWoh66 Theorem 2.7] it follows that $z \to \alpha_*(z)$ is unique holomorphic extension to ${{\mathbb{C}}}_0$ and therefore $z \to \alpha_*(z)$ cannot be identical to $z \to {{\mathrm{i}}\,}\,z/c_0$, the holomorphic extension of ${{\mathrm{i}}\,}\,\omega/c_0$. $\alpha_*(z)\neq{{\mathrm{i}}\,}\,z/c_0$ for $z\in {{\mathbb{C}}}_0$ implies that $k$ has no zeros and hence $z/k(z)$ is holomorphic on $\mathring{{{\mathbb{C}}}}_0$. This shows that Item \[it1\_Lions\] in Theorem \[th:lion\] is satisfied for $t \to {\mathcal{M}}(t,0)$.

    Since $1/(k(\xi+{{\mathrm{i}}\,}\,\eta)\,k(\xi))$ is bounded and $k(\xi+{{\mathrm{i}}\,}\,\eta)\to k(\xi)$ for $\eta\to 0$ in $S'$, it follows that for all $\psi\in {\mathcal{S}}$ $$\lim_{\eta\to 0} \int_{{\mathbb{R}}}\left[ \frac{k(\xi)-k(\xi+{{\mathrm{i}}\,}\,\eta)}{k(\xi+{{\mathrm{i}}\,}\,\eta)\,k(\xi)}\right]\,\psi(\xi)\,{\mbox{d}}\xi \to 0$$ i.e. $1/k(\xi+{{\mathrm{i}}\,}\,\eta)\to 1/k(\xi)$ for $\eta\to 0$ in $S'$. Hence ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\xi+{{\mathrm{i}}\,}\,\eta,0)\to {\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\xi,0)$ for $\eta\to 0$ in ${\mathcal{S}}'$, which shows Item \[it2\_Lions\] in Theorem \[th:lion\]. 

    Since $k(z)$ does not vanish on ${{\mathbb{C}}}_0$ and ${\left|k(z)\right|}$ is bounded by a polynomial in ${\left|z\right|}$ for $z\in\mathring{{{\mathbb{C}}}}_0$, it follows that ${\left|z/k(z)\right|}$ is bounded by a polynomial in ${\left|z\right|}$. Hence Item \[it3\_Lions\] in Theorem \[th:lion\] is satisfied and consequently $t \to {\mathcal{M}}(t,0)$ is causal. 

-   Property (\[propM1\]) is satisfied if $$\label{MMM}
    \begin{aligned}
        {\mathcal{M}}(t + {\left|t'\right|},t')=0 \qquad \mbox{ for } \quad t < 0\,. 
    \end{aligned}$$ From [(\[eq:defM\])]{} and $$K(c_0\,t',t)={\mathcal{F}}^{-1}\left\{ \frac{e^{-\alpha^*(\omega)\,c_0\,t'}}{\sqrt{2\,\pi}}\,  \right\}(t)\,,$$ it follows that $${\mathcal{M}}(t + {\left|t'\right|},t')
              = {\mathcal{M}}(t,0) \ast_t K(c_0\,t',t)\,.$$ Since $t\mapsto{\mathcal{M}}(t,0)$ and $t\mapsto K(c_0\,t',t)$ are causal, their convolution is also causal (cf.  Item \[item:conv\] in the Appendix).  This proves property (\[MMM\]) and concludes the proof. 

Assume that the attenuation coefficient is given by $$\alpha^*(\omega) = {{\mathrm{i}}\,}\,\omega/c_0  \qquad \mbox{ for }\qquad \omega\in{{\mathbb{R}}}\,  .$$ Then $$K(r,t) = \frac{1}{\sqrt{2\,\pi}}\,{\mathcal{F}^{-1} \left\{ {\exp \left( -{{\mathrm{i}}\,}\,\omega\,r/c_0 \right)} \right\}}(t)  = \delta(t+r/c_0)$$ which implies together with [(\[eq:FGG0\])]{} and [(\[eq:g0\])]{} that $${\mathcal{G}}(r,t) = \frac{\delta(t)}{4\,\pi\,r} \,.$$ But this function does not correspond to the intuition of an attenuated wave, which is manifested by the convolution equation (\[eq:FGG0\]), which should give a smooth decay of frequency components over travel distance.  With this Green function ${\mathcal{G}}$ the input impulse collapses immediately and consequently, in this case, the assumption $\alpha^*(\omega)\neq{{\mathrm{i}}\,}\,\omega/c_0$ in Theorem \[lemm02\] reflects physical reality. 

Solution of the Integral Equation
=================================

The inverse problem of photoacoustics with attenuated waves reduces to solving the integral equation (\[eq:pattrep02a\]) for $p_0$, and to the standard photoacoustical inverse problem, which consists in calculating the initial pressure $\rho$ in the wave equation (\[eq:ex:ivp3d\]) from measurements of $p_0({{\vec{x}}},t)$ over time on a manifold surrounding the object of interest.  The standard photoacoustical imaging problem is not discussed here further, but we focus on the the integral equation (\[eq:pattrep02a\]). 

In the following we investigate the ill–conditionness of the integral equation (\[eq:pattrep02a\]), where the kernel ${\mathcal{M}}$ is given from the attenuation law [(\[eq:powlaw2\])]{} with $\gamma \in (1,2]$. In this case the model is causal and the parameter range $\gamma \in (1,2]$ is relevant for biological imaging. 

In order to estimate the ill–conditionness of the integral equation (\[eq:pattrep02a\]) it is rewritten in Fourier domain: $$\label{eq:Mfrequency}
    {\mathcal{F} \left\{ {p_{\text{att}}}\right\}}({{\vec{x}}}_0,\omega)
         =  \int_{{\mathbb{R}}}{\mathcal{F} \left\{ {\mathcal{M}}\right\}} (\omega,t')\,p_0({{\vec{x}}}_0,t')\,{\mbox{d}}t'\;.$$ After discretization the ill-conditionness of this equation is reflected by the decay rate of the singular values of the matrix ${\mathcal{F} \left\{ {\mathcal{M}}\right\}} (\omega,t')$ at certain discrete frequencies and time instances. 

We consider simple test examples of attenuation coefficients $\rho$ (as in (\[eq:u(x)\])), which are characteristic functions of balls with center at the origin and radiii $R$. For these examples we investigate the dependence of the ill–conditionedness of (\[eq:Mfrequency\]) on the radius $R$ and the location ${{\vec{x}}}_0$. For applications in photoacoustic imaging ${{\vec{x}}}_0$ would be the location of a detector outside of the object of interest, to be imaged.  Then, by solving the integral equation (\[eq:Mfrequency\]) $p_0$ can be calculated, and in turn, the absorption energy $\rho$ can be reconstructed with standard backprojection formulas.  Since $${\mbox{supp}}(p_0({{\vec{x}}}_0,\cdot))
          = \left[({\left|{{\vec{x}}}_0\right|}-R)/c_0,({\left|{{\vec{x}}}_0\right|}+R)/c_0\right]\,,$$ the integral equation [(\[eq:Mfrequency\])]{} can be rewritten as $$\label{eq:Mfrequency*}
    {\mathcal{F} \left\{ {p_{\text{att}}}\right\}}({{\vec{x}}}_0,\omega)
         =  \int_{(R_0-R)/c_0}^\infty  {\mathcal{F} \left\{ {\mathcal{M}}\right\}} (\omega,t')\,p_0({{\vec{x}}}_0,t')\,{\mbox{d}}t'\,,$$ where $R_0 = {\left|{{\vec{x}}}\right|}_0$.

In the following we analyze the integral equation (\[eq:Mfrequency\*\]) in terms of the two parameters $R$ and $R_0={\left|{{\vec{x}}}_0\right|}$. This gives a clue on the effect of attenuation in terms of the size of the object and the distance of the location ${{\vec{x}}}_0$ to the simple object.  In order to show the effect of attenuation on the single frequencies, we make a singular value decomposition of the kernel of the integral equation (\[eq:Mfrequency\*\]). 

\[exam:distance2\] For small frequencies the attenuation law of castor oil, which behaves very similar to biological soft tissue, is approximately a power law with exponent $\gamma=1.66$ and $\hat{\alpha}_0 \approx 4\cdot 10^{-2}\,\frac{1}{cm\,(MHz)^\gamma}$, i.e. $$\alpha_{pl} (\omega) \approx  4\cdot 10^{-2}\cdot \omega^{1.66} \cdot cm^{-1} \qquad \mbox{($\omega$ in $MHz$)}.$$ The sound speed of castor oil is $1490\cdot\frac{m}{s}$ at $25$ degree Celsius.  In units of $cm$ and $MHz$ we have $$c_0\approx 0.15\cdot cm\cdot MHz\,.$$ Since  [(\[eq:powlaw2\])]{} approximates the power law (cf.  Figure \[fig:comp\]) it follows that $$\hat{\alpha}_0\,{\left|\omega\right|}^\gamma \approx
    \frac{\alpha_0\,\sin(\frac{\pi}{2}\,(\gamma-1))}{2\,c_0\,\tau_0}\,{\left|\tau_0\,\omega\right|}^\gamma$$ and consequently the coefficients of  [(\[eq:powlaw2\])]{} satisfy $$\alpha_0 \approx \frac{2\,c_0\,\hat{\alpha}_0}
                        {\tau_0^{(\gamma-1)}\,\sin(\frac{\pi}{2}\,(\gamma-1))}
                 \approx 6\,.$$ We note that the relaxation time is $\tau\approx 10^{-4}\,\frac{1}{MHz}$ for liquids (cf.  [@KinFreCopSan00]). 

For the calculation of the singular value decomposition of the discretized kernel of the integral equation (\[eq:Mfrequency\*\]) we used a frequency range\
$\omega \in [-80,80] MHz$ and step size $\Delta \omega=\frac{2\,\pi}{N-1}\, MHz$ with $N=2^9$. The time interval has been set to $[0,\frac{2\,\pi}{\Delta \omega}] MHz^{-1}$ and a step size $\Delta t=\frac{2\,\pi}{80}\,MHz^{-1}$ was used. 

The upper left picture in Fig.  \[fig:distance2\] visualizes the discretized kernel of the integral equation (\[eq:Mfrequency\*\]) for $R_0=R$, i.e., when ${{\vec{x}}}_0$ is directly on the surface of the object of interest.  The upper right picture shows the singular values of the discretized kernel in a logarithmic scale.  Two properties of the singular values become apparent:

1.   For large indices the decay rate is exponential, which can be seen from the linear decay in the logarithmic scale. 

2.   Secondly, there is a range of indices, where the singular values do not decay that rapidly.  As a consequence, for solving the integral equation this means that the Fourier coefficients of $p_0$ according to the first block of singular values can be determined in a stable manner. 

For increasing distance $L={\left|{{\vec{x}}}_0\right|}-R$ of ${{\vec{x}}}_0$ to the object the singular values of the discretization of the integral equation (\[eq:Mfrequency\*\]) show a drastically more exponentially decay rate for increasing $L$ (see bottom right picture Fig.  \[fig:distance2\]).  This means that if the object is further away from ${{\vec{x}}}_0$ attenuation is more drastically, and solution of the integral equation is more unstable.  We analyze the dependence of the number of largest singular values from $L$. For this purpose we denote by $n_{cut}$ the index of the singular value that is about $0.1\%$ of the maximal singular value.  For the numerical solution of (\[eq:Mfrequency\*\]) it means that if we make a truncated singular valued decomposition with only $n_{cut}$ singular values, the error amplification can be bounded by a factor $1000$. The dependence of $n_{cut}$ on $L$ is shown in the lower left picture of Fig.  \[fig:distance2\].  The picture reveals that for increasing distance (from about $2cm$) only about four Fourier modes of $p_0$ are significant when a maximal error amplification of a factor $1000$ is required.  This reveals that in general the solution of the integral equation (\[eq:Mfrequency\*\]) is significantly ill–posed and worse if the data recording is far away from the object. 

![ Case: $\gamma=1.66$ (castor oil).  The upper left and right pictures visualize the kernel ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ and its singular values for $L:=R_0-R=0cm$. The lower right and left pictures visualize ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ for the detector distances $L=0\cdot cm$, $L=2\cdot cm$ and $L=4\cdot cm$ and the respective indices $n_{cut}$ for which the singular values are about $0.1$ per cent of the maximal singular value.[]{data-label="fig:distance2"}](ExSingVal_2_00.eps){width="\textwidth"}

\[exam:distance1\] An analogous numerical example as in Example \[exam:distance2\] for the case $\gamma=1.1$ is presented in Fig.  \[fig:distance1\].  From the lower left picture of Fig.  \[fig:distance1\], we see that if the distance is about $2\cdot cm$ from the boundary of the object, then $17$ singular values are available for the numerical estimation. 

![ Case: $\gamma=1.1$. The upper left and right pictures visualize the kernel ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ and its singular values.  The lower right and left pictures visualize ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ for the detector distances $L=0\cdot cm$, $L=2\cdot cm$ and $L=4\cdot cm$ and the respective indices $n_{cut}$ for which the singular values are about $0.1$ per cent of the maximal singular value.[]{data-label="fig:distance1"}](ExSingVal_1_00.eps){width="\textwidth"}

\[exam:distance3\] An analogous numerical example as in Example \[exam:distance2\] for the case $\gamma=2$ is presented in Fig.  \[fig:distance3\].  From the lower left picture of Fig.  \[fig:distance3\], we see that if the distance is about $2\cdot cm$ from the boundary of the object, then only $4$ singular values are available for the numerical estimation. 

![ Case: $\gamma=1.1$. The upper left and right pictures visualize the kernel ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ and its singular values.  The lower right and left pictures visualize ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ for the detector distances $L=0\cdot cm$, $L=2\cdot cm$ and $L=4\cdot cm$ and the respective indices $n_{cut}$ for which the singular values are about $0.1$ per cent of the maximal singular value.[]{data-label="fig:distance3"}](ExSingVal_3_00.eps){width="\textwidth"}

\[exam:distance5\] An analogous numerical example as in Example \[exam:distance2\] for the frequency power law $$\alpha_*(\omega) = \alpha_0^{pl}\cdot (-{{\mathrm{i}}\,}\,\omega)^{0.66}
  \qquad\quad (\mbox{$\alpha_0^{pl}$ as in Example~\ref{exam:distance2}})$$ is presented in Fig.  \[fig:distance5\].  From the lower left picture of this figure, we see that if the distance is about $2\cdot cm$ from the boundary of the object, then $77$ singular values are available for the numerical estimation.  If the distance is about $4cm$, then $46$ singular values are available. 

![ Case: $\gamma=0.66$. The upper left and right pictures visualize the kernel ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ and its singular values.  The lower right and left pictures visualize ${\mathcal{F} \left\{ {\mathcal{M}}\right\}}(\omega,t')$ for the detector distances $L=0\cdot cm$, $L=2\cdot cm$ and $L=4\cdot cm$ and the respective indices $n_{cut}$ for which the singular values are about $0.1$ per cent of the maximal singular value.[]{data-label="fig:distance5"}](ExSingVal_5_00.eps){width="\textwidth"}

Comparing all numerical examples shows that the larger $\gamma$ (stronger attenuation), the more rapidly decrease the singular values. 

Appendix: Nomenclature and Elementary Facts {#sec:app}
===========================================

[Sets:]{} $B_R$ denotes the open ball with center at $\mathbf{0}$ and radius $R$. $S^n \subseteq {{\mathbb{R}}}^n$ denotes the $n$-dimensional unit sphere. 

[Real and Complex Numbers:]{} ${{\mathbb{C}}}$ denotes the space of complex numbers, ${{\mathbb{R}}}$ the space of reals.  For a complex number $c=a+{{\mathrm{i}}\,}b$ $a=\Re{(c)}$, $b=\Im{(c)}$ denote the real and imaginary parts, respectively. 

For a complex number $c$ we denote by ${\left|c\right|}$ the absolute value and by $\phi \in (-\pi,\pi]$ the argument.  That is $$c = {\left|c\right|} {\exp \left( i \phi \right)}\;.$$ As a consequence, when $w=r\,{\exp \left( {{\mathrm{i}}\,}\phi \right)}$ then $$\label{defpower}
        w^\gamma = {\exp \left( \gamma\,\left(\text{log}(r)+{{\mathrm{i}}\,}\phi\right) \right)}\;.$$ Consequently $w^\gamma$ has absolute value $r^\gamma$ and the argument is $\gamma \phi$ modulo $2\pi$. In this paper all power functions are defined on ${{\mathbb{C}}}\backslash {{\mathbb{R}}}_-$. We note that $$\label{propsqrt}
         w\in{{\mathbb{C}}}\backslash {{\mathbb{R}}}_-   \quad\mbox{and}\quad \Re(\sqrt{w})>0  \qquad \Rightarrow \qquad
         \Im(w)\,\Im(\sqrt{w})\geq 0  \,.$$

[Differential Operators:]{} $\nabla$ denotes the gradient.  $\nabla \cdot$ denotes divergence, and $\nabla^2$ denotes the Laplacian. 

[Product:]{} When we write $\cdot$ between two functions, then it means a pointwise product, it can be a scaler product or if the functions are vector valued an inner product.  The product between a function and a number is not explicitly stated. 

[Composition:]{} The composition of operators ${\mathcal{A}}$ and ${\mathcal{B}}$ is written as ${\mathcal{A}}{\mathcal{B}}$.

[Special functions:]{}

-   The *signum* function is defined by $${\mbox{sgn}}:= {\mbox{sgn}}({{\vec{x}}}) := \frac{{{\vec{x}}}}{{\left|{{\vec{x}}}\right|}}\;.$$ In ${{\mathbb{R}}}^3$ it satisfies $$\label{eq:der_sgn} \nabla \cdot {\mbox{sgn}}= \frac{2}{{\left|{{\vec{x}}}\right|}}\;.$$

-   The *Heaviside* function $${H}:= {H}(t) := \left\{ \begin{array}{rcl}
                          0 & \text{ for } & t < 0\\
                          1 & \text{ for } & t > 0\\
                         \end{array} \right.$$ satisfies $${H}:= \frac{1}{2} (1+{\mbox{sgn}})\;.$$

-   The $\delta$-distribution is the derivative of the Heaviside function at $0$ and is denoted by $\delta_t := \delta_t(t)$. In our terminology $\delta_t$ denotes a *one*-dimensional distribution.  Sometimes, if the context is clear, we will omit the subscript at the $\delta$-distributions. 

-   The three dimensional $\delta$-distribution $\delta_{{\vec{x}}}$ is the tensor product of the three one-dimensional distributions $\delta_{x_i}$, $i=1,2,3$. Moreover, $$\label{eq:hatdelta} \delta_{{{\vec{x}}},t} := \delta_{{{\vec{x}}},t}({{\vec{x}}},t) = \delta_{{\vec{x}}}\cdot \delta_t,$$ is a four dimensional distribution in space and time.  If we do not add a subscript $\delta$ denotes a one-dimensional $\delta$-distribution. 

-   $\chi_\Omega$ denotes the characteristic set of $\Omega$, i.e., it attains the value $1$ in $\Omega$ and is zero else. 

[Properties related to functions:]{} ${\mbox{supp}}(g)$ denote the *support* of the function $g$, that is the closure of the set of points, where $g$ does not vanish. 

[Derivative with respect to radial components:]{} We use the notation $$r:=r({{\vec{x}}}) = {\left|{{\vec{x}}}\right|},$$ and denote the derivative of a function $f$, which is only dependent on the radial component ${\left|{{\vec{x}}}\right|}$, with respect to $r$ (i.e., with respect to ${\left|{{\vec{x}}}\right|}$) by $\cdot'$.

Let $\beta = \beta(r)$, then it is also identified with the function $\beta = \beta ({\left|{{\vec{x}}}\right|})$ and therefore $$\nabla \beta = \frac{{{\vec{x}}}}{{\left|{{\vec{x}}}\right|}} \beta'\;.$$

[Convolutions:]{} Three different types of convolutions are considered: $*_t$ and $*_\omega$ denote *convolutions* with respect to time and frequency, respectively.  Let $f$, $\hat{f}$, $g$ and $\hat{g}$ be functions defined on the real line with complex values.  Then $$\begin{aligned}
 &f *_t g := \int_{{\mathbb{R}}}f(t-t')g(t') d t' , \quad\quad
 &\hat{f} *_\omega \hat{g}
        := \int_{{\mathbb{R}}}\hat{f}(\omega-\omega')\hat{g}(\omega') d \omega'. 
\end{aligned}$$ $*_{{{\vec{x}}},t}$ denotes space–time convolution and is defined as follows: Let $f,g$ be functions defined on the Euclidean space ${{\mathbb{R}}}^3$ with complex values, then $$f *_{{{\vec{x}}},t} g := \int_{{{\mathbb{R}}}^3} \int_{{\mathbb{R}}}f({{\vec{x}}}- {{\vec{x}}}',t-t')g({{\vec{x}}}',t') d{{\vec{x}}}' d t'\;.$$

[Fourier transform:]{} For more background we refer to [@Lig64; @Tit48; @Pap62; @Yos95; @Hoe03].  All along this paper ${\mathcal{F} \left\{ \cdot \right\}}$ denotes the *Fourier transformation* with respect to $t$, and the *inverse Fourier transform* ${\mathcal{F}^{-1} \left\{ \cdot \right\}}$ is with respect to $\omega$. In this paper we use the following definitions of the transforms: $$\begin{aligned}
{\mathcal{F} \left\{ f \right\}}(\omega) &= \frac{1}{\sqrt{2\pi}} \int_{{\mathbb{R}}}{\exp \left( {{\mathrm{i}}\,}\omega t \right)} f(t) d t\,,\\
{\mathcal{F}^{-1} \left\{ \hat{f} \right\}}(t) &= \frac{1}{\sqrt{2\pi}} \int_{{\mathbb{R}}}{\exp \left( -{{\mathrm{i}}\,}\omega t \right)} \hat{f}(\omega) d\omega\;. 
\end{aligned}$$ The Fourier transform and its inverse have the following properties:

1.   \[item:derF\] $${\mathcal{F} \left\{ \frac{\partial }{\partial t}f \right\}} (\omega) =  (-{{\mathrm{i}}\,}\omega){\mathcal{F} \left\{ f \right\}}(\omega)\;.$$

2.   $$\begin{aligned}
          {\mathcal{F} \left\{ f \cdot g \right\}} & =\frac{1}{\sqrt{2\pi}} {\mathcal{F} \left\{ f \right\}} *_\omega {\mathcal{F} \left\{ g \right\}} \text{ and }\\
          {\mathcal{F} \left\{ f \right\}} \cdot {\mathcal{F} \left\{ g \right\}} &= \frac{1}{\sqrt{2\pi}} {\mathcal{F} \left\{ f *_t g \right\}},\\
          {\mathcal{F}^{-1} \left\{ \hat{f} \cdot \hat{g} \right\}}
                 &= \frac{1}{\sqrt{2\pi}}{\mathcal{F}^{-1} \left\{ \hat{f} \right\}} *_t {\mathcal{F}^{-1} \left\{ \hat{g} \right\}} \text{ and }\\
          {\mathcal{F}^{-1} \left\{ \hat{f} \right\}} \cdot {\mathcal{F}^{-1} \left\{ \hat{g} \right\}}
                 &= \frac{1}{\sqrt{2\pi}} {\mathcal{F}^{-1} \left\{ \hat{f} *_\omega \hat{g} \right\}}\;. 
          \end{aligned}$$

3.   For $a \in {{\mathbb{R}}}$ $${\mathcal{F} \left\{ f(t-a) \right\}}(\omega) = {\exp \left( i a \omega \right)} \cdot{\mathcal{F} \left\{ f(t) \right\}}(\omega)$$

4.   \[item:Fdelta\] The $\delta$-distribution at $a \in {{\mathbb{R}}}$ satisfies $$\delta_t(t-a)=\frac{1}{\sqrt{2\pi}}{\mathcal{F}^{-1} \left\{ \exp({{\mathrm{i}}\,}a \omega) \right\}}(t)\;.$$

5.   \[item:even\] Let $f$ be real and even, odd respectively, then ${\mathcal{F} \left\{ f \right\}}$ is real and even, imaginary and odd, respectively. 

6.   \[item:temp\] The Fourier transformation of a tempered distribution is a tempered distribution. 

7.   \[item:conv\] Let $\tau_1,\,\tau_2\in {{\mathbb{R}}}$. If $f_1$ and $f_2$ are two distributions with support in $[\tau_1,\infty)$ and $[\tau_2,\infty)$, respectively, then $f_1*f_2$ is well-defined and (cf.   [@Hoe03]) $$\label{propIM}
                \mbox{supp} (f_1 * f_2)
               \subseteq  \mbox{supp} (f_1) + \mbox{supp} (f_2) \subseteq [\tau_1+\tau_2,\infty)\,.$$

[The Hilbert transform for $L^2-$functions is defined by]{} $${\mathcal{H} \left\{ f \right\}} (t) = \frac{1}{\pi} \Xint-_{{\mathbb{R}}}\frac{f(s)}{t-s}ds\;,$$ where $\Xint-_{{\mathbb{R}}}f(s) ds$ denotes the Cauchy principal value of $\int_{{\mathbb{R}}}f(s) ds$.

A more general definition of the Hilbert transform can be found in [@BelWoh66].  The Hilbert transform satisfies

-   ${\mathcal{H} \left\{ {\mathcal{F} \left\{ f \right\}} \right\}}(\omega) = -i{\mathcal{F} \left\{ {\mbox{sgn}}f \right\}}(\omega)$,

-   ${\mathcal{H} \left\{ {\mathcal{H} \left\{ f \right\}} \right\}} = - f$.

From the first of these properties the Kramers-Kronig relation can be formally derived as follows.  Since $f(t)$ is a causal function if and only if $ f = {H}\cdot f$ and ${H}=(1+{\mbox{sgn}})/2$, it follows that ${\mathcal{F} \left\{ f \right\}} =  [{\mathcal{F} \left\{ f \right\}} +{{\mathrm{i}}\,}{\mathcal{H} \left\{ {\mathcal{F} \left\{ f \right\}} \right\}}]/2 $, which is equivalent to ${\mathcal{F} \left\{ f \right\}} = i{\mathcal{H} \left\{ {\mathcal{F} \left\{ f \right\}} \right\}}$, i.e. $$\Re({\mathcal{F} \left\{ f \right\}}) =  -\Im({\mathcal{H} \left\{ {\mathcal{F} \left\{ f \right\}} \right\}})  {\qquad \text{and} \qquad}
      \Im({\mathcal{F} \left\{ f \right\}}) =  \Re({\mathcal{H} \left\{ {\mathcal{F} \left\{ f \right\}} \right\}}) .$$

[The inverse Laplace transform]{} of $f$ is defined by $${\mathcal{L}^{-1} \left\{ f \right\}}(t) = \left\{ \begin{array}{ccl}
0 &\text{ for }& t < 0,\\
\frac{1}{2 \pi i} \int_{\gamma - i\infty}^{\gamma + {{\mathrm{i}}\,}\infty} {\exp \left( st \right)} f(s)ds, &\text{ for }& t > 0,
\end{array}\right.$$ where $\gamma$ is appropriately chosen. 

The inverse Laplace transform satisfies (see e.g. [@Heu91]) $$\label{eq:pr1} {\mathcal{L}^{-1} \left\{ h(s-a) \right\}}(t) = {\exp \left( a t \right)} {\mathcal{L}^{-1} \left\{ h(s) \right\}}(t) \text{ for all } a, t \in {{\mathbb{R}}}$$ and $$\label{eq:pr2}
{\mathcal{L}^{-1} \left\{ s^{-r} \right\}}(t) = \frac{H(t) t^{r-1}}{\Gamma(r)}\qquad (r>0)\;.$$

Acknowledgement {#acknowledgement .unnumbered}
===============

This work has been supported by the Austrian Science Fund (FWF) within the national research network Photoacoustic Imaging in Biology and Medicine, project S10505-N20. 

[\[KLR73\]]{}

Agranovsky, M. and Kuchment, P.: , vol.23, No.  5, 2089–2102, 2007. 

Beltrami, E. J. and Wohlers, M. R.: .  Academic Press, New York and London, 1966. 

Burgholzer, P. and Bauer-Marschallinger, J. and Grün, H. and Haltmeier, M. and Paltauf, G.: Temporal back-projection algorithms for photoacoustic tomography with integrating line detectors. [* Inverse Probl.*]{},  23(6):65-80, 2007. 

Burgholzer, P. and Grün, H. and Haltmeier, M. and Nuster, R. and Paltauf, G.: Compensation of acoustic attenuation for high-resolution photoacoustic imaging with line detectors.  In A.A. Oraevsky and L.V. Wang, editors, [*Photons Plus Ultrasound: Imaging and Sensing 2007: The Eighth Conference on Biomedical Thermoacoustics, Optoacoustics, and Acousto-optics*]{}, volume 6437 of [ *Proceedings of SPIE*]{}, page 643724.  SPIE, 2007. 

Burgholzer, P., Roitner, H., Bauer-Marschallinger, J., Paltauf, G.: Image Reconstruction in Photoacoustic Tomography Using Integrating Detectors Accounting for Frequency-Dependent Attenuation.  Volume 7564, page 75640O. Proc.  SPIE, 2010. 

Chen, W and Holm, S.: Fractional Laplacian time-space models for linear and nonlinear lossy media exhibiting arbitrary frequency power-law dependency.  115 (4), April 2004. 

Cheong, W. F. and Prahl, S. A. and Welch, A. J.: A review of the optical properties of biological tissues. ,  26(12):2166–2185, 1990. 

Dautray, R. and Lions, J.-L.: Mathematical Analysis and Numerical Methods for Science and Technology.  Volume 1. [* Springer-Verlag*]{}, New York, 2000. 

Dautray, R. and Lions, J.-L.: Mathematical Analysis and Numerical Methods for Science and Technology.  Volume 2. [* Springer-Verlag*]{}, New York, 2000. 

Dautray, R. and Lions, J.-L.: Mathematical Analysis and Numerical Methods for Science and Technology.  Volume 5. [* Springer-Verlag*]{}, New York, 2000. 

Finch, D. and Patch, S. and Rakesh: Determining a function from its mean values over a family of spheres. [* Siam J. Math.  Anal.*]{}  Vol.  35, No.  5, pp.  1213-1240. 

Gasquet, C and Witomski, P: Fourier Analysis and Applications. [* Springer-Verlag*]{}, New York, 1999. 

Gusev, V. E. and Karabutov, A. A.: Laser Optoacoustics. [* American Institute of Physics*]{}, New York, 1993. 

Haltmeier, M. and Scherzer, O. and Burgholzer, P. and Paltauf, P.: Thermoacoustic imaging with large planar receivers. [* Inverse Probl.*]{},  20(5):1663-1673, 2004. 

Hanyga, A. and Seredynska, M.: Power-law attenuation in acoustic and isotropic anelastic media. [* Geophys.  J. Int*]{}, 155:830-838, 2003. 

H. Heuser. .  Teubner, Stuttgart, second edition, 1991. 

Hörmander, L.: The Analysis of Linear Partial Differential Operators I. [*Springer Verlag*]{}, New York, 2nd edition, 2003. 

Hristova, Y. and Kuchment, P. and Nguyen, L.: Reconstruction and time reversal in thermoacoustic tomography in acoustically homogeneous and inhomogeneous media. [* Inverse Problems*]{}, 24(5):055006 (25pp), 2008. 

John, F.: Partial Differential Equations. [* Springer Verlag*]{}, New York, 1982. 

Kilbas, A.A. and Srivastava, H.M. and Trujillo, J.J.: , volume 204 of [*North-Holland Mathematics Studies*]{}.  Elsevier Science B.V., Amsterdam, 2006. 

Kinsler, L. E., Frey, A. R., Coppens, A. B., Sanders, J. V.: Fundamentals of Acoustics. ,  New York, 2000. 

Kowar, R .  and Scherzer, O. and Bonnefond, X.: to apper in: 2010, DOI: 10.1002/mma.1344

Kruger, R. A. and Kiser, W. L. and Miller, K. D. and Reynolds, H. E.: Thermoacoustic [CT]{}: imaging principles.  Proc.  SPIE 2000, vol.  3916, 150–159. 

Ku, G. and Wang, X. and Stoica, G. and Wang, L. V.: Multiple-bandwidth photoacoustic tomography. [* Phys.  Med.  Biol.*]{},  49:1329–1338, 2004. 

Kunyansky, L. A.: Explicit inversion formulae for the spherical mean Radon transform. [* Inverse Probl.*]{},  23,373-383, 2007. 

Kuchment, P., Kunyansky, L. A.: Mathematics of thermoacoustic and photoacoustic tomography. [* European J. Appl.  Math.*]{},  19:191–224, 2008. 

Landau, L. D. and E.M. Lifschitz, E. M.: Lehrbuch der theoretischen Physik, Band VII: Elastizitätstheorie. [* Akademie Verlag*]{}, Berlin, 1991. 

Lighthill, M.J.: .  Cambridge University Press, London, 1964. 

Nachman, A. I. and Smith, J. F., III and Waag, R. C.: An equation for acoustic propagation in inhomogeneous media with relaxation losses. [* J. Acoust.  Soc.  Am.*]{}  88 (3), Sept. 1990. 

Oraevsky, A. and Wang, L.V., editors: , volume 6437 of [*Proceedings of SPIE*]{}, 2007. 

Patch, S. K. and Scherzer, O.: Special section on photo- and thermoacoustic imaging. [* Inverse Probl.*]{},  23:S1–S122, 2007. 

Patch, S. K. and Greenleaf, A.: Equations governing waves with attenuation according to power law. [* Technical report*]{}, Department of Physics, University of Wisconsin-Milwaukee, 2006. 

Papoulis, A.: .  McGraw-Hill, New York, 1962. 

Podlubny, I.: , volume 198 of [ *Mathematics in Science and Engineering*]{}.  Academic Press Inc., San Diego, CA, 1999. 

Razansky, D. and Distel, M. and Vinegoni, C. and Ma, R. and Perrimon, N. and Köster, R. W. and Ntziachristos, V.: Multispectral opto-acoustic tomography of deep-seated fluorescent proteins in vivo, [*Nature Photonics*]{} 3: 412-417, 2009. 

La Rivi[é]{}re, P. J. and Zhang, J. and Anastasio, M. A.: Image reconstruction in optoacoustic tomography for dispersive acoustic media. [* Opt.  Letters*]{}, 31(6):781–783, 2006. 

Scherzer, O. and Grossauer, H. and Lenzen, F. and Grasmair, M. and Haltmeier, M.: Variational Methods in Inmaging. ,  New York, 2009. 

Sushilov, N. V. and Cobbold, R. S. C.: Frequency-domain wave equation and its time-domain solution in attenuating media. [* Journal of the Acoustical Society of America*]{}, 115:1431–1436, 2005. 

Szabo, T. L.: Time domain wave equations for lossy media obeying a frequency power law. [* J. Acoust.  Soc.  Amer.*]{},  96:491–500, 1994. 

Szabo, T. L.: Causal theories and data for acoustic attenuation obeying a frequency power law. [* J. Acoust.  Soc.  Amer.*]{},  97:14–24, 1995. 

Tam, A. C.: Applications of photoacoustic sensing techniques. [* Rev. Modern Phys.*]{},  58(2):381–431, 1986. 

Titchmarch, E. C.: .  Clarendon Press, Oxford, 1948. 

Waters, K. R. and Hughes, M. S. and Brandenburger, G. H. and Miller, J. G.: On a time-domain representation of the [K]{}ramers-[K]{}r[ö]{}nig dispersion relation. ,  108(5):2114–2119, 2000. 

Waters, K.R. and Mobely, J. and Miller, J. G.: Causality-Imposed ([K]{}ramers-[K]{}r[ö]{}nig) Relationships Between Attenuation and Dispersion. ,  vol.  52, no.  5, May 2005. 

Webb, S., editor: , Bristol, Philadelphia, 2000.  reprint of the 1988 edition. 

Wang, L. V.: Prospects of photoacoustic tomography. [* Med.  Phys.*]{},  35(12):5758–5767, 2008. 

Wang, X. D. and Pang, Y. J. and Ku, G. and Xie, X. Y. and Stoica, G. and Wang, L. V.: 2003, vol.  21, No.  7, 803–806

Xu, Y. and Feng, D. and Wang, L. V.: Exact Frequency-Domain Reconstruction for Thermoacoustic Tomography - I: Planar Geometry. [* IEEE Trans.  Med.  Imag.*]{},  Vol.  21, N0.  7, July 2002. 

Xu, Y. and Xu, M. and Wang, L. V.: Exact Frequency-Domain Reconstruction for Thermoacoustic Tomography - II: Cylindrical Geometry. [* IEEE Trans.  Med.  Imag.*]{},  Vol.  21, N0.  7, July 2002. 

Xu, M. and Xu, Y. and Wang, L. V.: Time-Domain Reconstruction Algorithms and Numerical Simulation for Thermoacoustic Tomography in Various Geometries. [* IEEE Trans.  Biomed.  Eng.*]{},  Vol.  50, N0.  9, Sept. 2003. 

Xu, Y. and Wang, L. V. and Ambartsoumian, G. and Kuchment, P.: Reconstructions in limited-view thermoacoustic tomography. [* Med.  Phys.*]{},  31 (4), April 2004. 

Xu, M. and Wang, L. V.: Universal back-projection algorithm for photoacoustic computed tomography. [* Phys.  Rev. E 71*]{}, 2005.  Article ID 016706. 

Xu, M. and Wang, L. V.: Photoacoustic imaging in biomedicine. [* Rev. Sci.  Instruments*]{}, 77(4):1–22, 2006.  Article ID 041101. 

Yosida, K.: Functional analysis.  Springer-Verlag, Berlin, Heidelberg, New York 1995, 5th edition

Zhang, E.Z. and Laufer, J. and Beard, P.: Three-dimensional photoacoustic imaging of vascular anatomy in small animals using an optical detection system.  In [*[@OraWan07]*]{}, 2007. 

Zhang, H. and Maslov, K. and Stoika, G. and Wang.  V.L.: Functional photoacoustic microscopy for high-resolution and noninvasive in vivo imaging. ,  24:848 – 851, 2006. 

[^1]: In this paper the root of a complex number is always the one with non-negative real part. 

[^2]: In [@NacSmiWaa90] they use the notion $c_\infty$ for $c_0$ and $c$ for $\tilde c_0$.

[^3]: A function $f$ is an element of $L_{loc}^1({{\mathbb{R}}})$ if it is in $L^1$ on every compact set. 

[^4]: In Theorem 4 on p294 ff in [@DauLio02_5] the assumption that $f$ is strongly causal is expressed by $f \in \mathring{{\mathcal{D}}}_+$, which is the set of distributions with support in $[0,+\infty)$.

[^5]: A function $F:{{\mathbb{C}}}_0 \to{{\mathbb{C}}}$ is holomorphic in $\mathring{{{\mathbb{C}}}}_0$ if it is complex differentiable in $\mathring{{{\mathbb{C}}}}_0$. Sometimes the functions are also refered to as analytic or regular functions or conformal maps. 

[^6]: A function $F:{{\mathbb{C}}}_0\to{{\mathbb{C}}}$ which satisfies Items \[it1\_Lions\],\[it2\_Lions\] of Theorem \[th:lion\] is called *holomorphic extension* of ${\mathcal{F} \left\{ f \right\}}$.

[^7]: $\alpha: {{\mathbb{C}}}_0 \to {{\mathbb{C}}}$ extends the function $\omega \in {{\mathbb{R}}}\to \alpha(\omega)$ but is not an holomorphic extension.