WHOI 3D parabolic equation code for Matlab.
WHOI rev 11. CARPE3D.  (Cartesian 3D PE)

This is an underwater acoustic propagation simulation code. Cartesian coordinates, 
FFT method used ('Tappert style, 1974,1977'). 

When publishing or presenting results made with
this, please cite this document:

Duda, T. F., Initial results from a Cartesian three-dimensional 
parabolic equation acoustical propagation code, WHOI Tech. Rept.,
WHOI-2006-041, 20 pages, 2006.  http://hdl.handle.net/1912/1428/

Tim Duda, WHOI,  14 July 2023

-----------------------------------------------
Rev 11 notes: This updates initially released 2009 code by removing a few unused 
variables and  uninformative comments. A few new comments are inserted. 

This code uses a fluid seabed approximation, as do most underwater acoustic PE codes.
It also imposes uniform density everywhere, even the seabed. Therefore the procedure 
of smoothing the medium discontinuities, to prevent aliasing by Fourier transforms 
to the wavenumber domain at each step, is not needed. This means that the sound field 
pressure is solved for, not reduced pressure, which factors in the density.

The script New_runtestcase.m is set up to simulate 200 Hz in a wedge environment. Here,
the parameters are set up, a function (CARPE3D) is called that does a few more things, 
and finally the function that is the actual wide-angle WAPE propagator is executed by 
CARPE3D.

The PE steps in the x direction to compute the field psi(x,y,z), and the wedge gets
deeper in the +y direction. The z unit vector points upward in some output plots, to
make a right-handed coordinate system. Unfortiately this code was built from older 
PE codes that have z increasing downwards, and some inputs are specified that way, 
things like seabed depth, sound source depth, depths for input sound speeds.

The main propagator is in the function propWAPErev11.m. The wide-angle propagator is 
used. This is a wide-angle (Thomson and Chapman, 1983; Feit and Fleck, 1978). Some of
the grid and code set up coding lines are from the narrow-angle Cartesian 3D PE 
fortran codes (standard PE) used by Martin and Flatte', 1988. 

The water sound speed 3D file format in the example setup is flexible. If range-dependent
bathymetry and/or sound-speed is desired, these fields are updated by executing an 'envfile'
at a range-setup interval that is given by a parameter. The name of the file to execute is an 
input parameter. In this case, the 'envfile' needs data from another file and that is another 
parameter. In the wedge example, the sound speed in water (cw) is in slices in the transverse 
direction-depth (y,z) slices, for a number of stepping positions (x) starting at the source 
plane, x=0. The format has cw in 121x81 vertical slice. The range-dependent example has small 
sinusoidal oscillations in x of isoveolcity depths, so that isovelocity surfaces of constant 
cw are like a corrugated metal roof with the crests aligned in the y direction

dat=load('layered302rd.mat')
dat = 
  struct with fields:

      cw: [301×121×81 double]
       x: [1×301 double]
       y: [121×81 double]
      y1: [1×121 double]
    zout: [301×121×81 double]

This version has a flag to run with range-independent cw. This runs a little faster because the
cw is only interpolated onto the compute (y,z) grid once.


References:

Feit, M. D., and Fleck, J. A. , Jr. (1978). “Light propagation in graded-index 
   fibers,” Appl. Opt. 17, 3990–3998.

Martin, J. M. , and Flatté, S. M. (1988). “Intensity images and statistics 
  from numerical simulation of wave propagation in 3-D random media,” 
  Appl. Opt. 27, 2111–2126.

Tappert, F. D. (1974). “Parabolic equation method in underwater acoustics,”
   J. Acoust. Soc. Am. 55, S34.

Tappert, F. D. (1977). “The parabolic equation method,” in Wave Propagation
   and Underwater Acoustics, edited by J. B. Keller and J. Papadakis, 
   Lecture Notes in Physics Vol. 70 (Springer, New York), pp. 224–286.

Thomson, D. J. , and Chapman, N. R. (1983). “A wide-angle split-step 
   algorithm for the parabolic equation,” J. Acoust. Soc. Am. 74, 1848–1854.




---Prior version notes follow ----------------------------------------------------------------

WHOI 3D parabolic equation code for Matlab 7. Release version 1,
WHOI rev 10.

Acoustic propagation simulation code.

Tim Duda  23 June 2009
-----------------------------------------------

rev 10 notes: Start by looking at or running runtestcase.m . This simple
example reads in data to generate a y-z grid of refractive index in
a domain of fixed water depth, replicating the domain each and
every x step (The solution is marched in the x direction). The seafloor
sound speed and attenuation have simple structure.

When publishing or presenting results made with
this, please cite this document:

Duda, T. F., Initial results from a Cartesian three-dimensional 
parabolic equation acoustical propagation code, WHOI Tech. Rept.,
WHOI-2006-041, 20 pages, 2006.  http://hdl.handle.net/1912/1428/

This code runs the Thomson-Chapman two-step split-step algorithm
with a standard starter. Density is uniform everywhere, typical
and virtually required for split-step codes.

The grid is generated in units scaled with the wavelength, i.e. the
code will run for any frequency. 

The major hurdle when running this code is providing 2-d grids
of refractive index every step without flooding memory or excessively 
accessing the disk drive and slowing down the operation. You will see
that the mechanism for doing this is intertwined with the propagation 
loop (as it must be), but is handled by repeatedly calling an 
environment file. Optimizing the operation of the enviromnent file is key.

The example uses 1024 x 2048 fft's and will run fine on a windows PC with
2GB RAM.