I need to solve the standard acoustic wave equation in 2D: u_tt = c(x,y)^2 (u_xx + u_yy - f) with absorbing boundary conditions. Have written code in c for that, which uses the PML approach for the absorbing boundaries. Now I need to embed my solver in matlab, to obtain a matlab function which recieves grid data, c^2, f, initial conditions, and returns the solution u. For this I embarked on writing a mex interface - and this turns out to be a heavier task then anticipated. 

So my question is: is it feasible (performance wise) to write a solver for a 2D wave equation directly in matlab code? How is its performance expected to compare with c/fortran code for the same task? Does any one have good/bad experience in using matlab for similar PDEs (without resorting to c/fortran) ?

Thanks!

I should add - 

My c code uses a finite difference scheme with 4th order differences for approximating spatial derivatives. The grids I need to work with range from 200 to 1000 points in each direction, i.e. a total of between 200^2=4x10^4 to 10^6 grid points... 




"Ilan" wrote in message <ijtrk4$3ii$1@fred.mathworks.com>...
> I need to solve the standard acoustic wave equation in 2D: u_tt = c(x,y)^2 (u_xx + u_yy - f) with absorbing boundary conditions. Have written code in c for that, which uses the PML approach for the absorbing boundaries. Now I need to embed my solver in matlab, to obtain a matlab function which recieves grid data, c^2, f, initial conditions, and returns the solution u. For this I embarked on writing a mex interface - and this turns out to be a heavier task then anticipated. 
> 
> So my question is: is it feasible (performance wise) to write a solver for a 2D wave equation directly in matlab code? How is its performance expected to compare with c/fortran code for the same task? Does any one have good/bad experience in using matlab for similar PDEs (without resorting to c/fortran) ?
> 
> Thanks!

On 2/21/2011 6:08 AM, Ilan wrote:
> I need to solve the standard acoustic wave equation in 2D: u_tt = c(x,y)^2 (u_xx + u_yy - f) with
>  absorbing boundary conditions. Have written code in c for that, which uses the PML approach for the
>  absorbing boundaries. Now I need to embed my solver in matlab, to obtain a matlab function which
>  recieves grid data, c^2, f, initial conditions, and returns the solution u. For this I embarked
>  on writing a mex interface - and this turns out to be a heavier task then anticipated.
>

If you already wrote the solver in C, why do you want to run it inside Matlab?


> So my question is: is it feasible (performance wise) to write a solver for a 2D wave
>  equation directly in matlab code? How is its performance expected to compare with c/fortran
>  code for the same task? Does any one have good/bad experience in using matlab for similar PDEs
>  (without resorting to c/fortran) ?
>
> Thanks!

What is the scale of your problem? how large the grid? is this a school project?
real project? performance is all relative.

For school projects, I never have a problem with performance, using Matlab 2010a,
new PC, but then again, the pde solvers I write all run on small size grid (512 by 512 is
largest I tried). 2-3 minutes to run. iterative solver. much faster using sparse
direct solver, depends on the problem.

--Nasser

Nasser - thanks for your quick reply.

> If you already wrote the solver in C, why do you want to run it inside Matlab?
> 

I am researching inverse problems, this solver is just part of a larger process of finding c^2 if we have partial knowledge on the solution.

> 
> What is the scale of your problem? how large the grid? is this a school project?
> real project? performance is all relative.
> 

As I wrote in the ps above scale is between 200 by 200, to 1000 by 1000.

> For school projects, I never have a problem with performance, using Matlab 2010a,
> new PC, but then again, the pde solvers I write all run on small size grid (512 by 512 is
> largest I tried). 2-3 minutes to run. iterative solver. much faster using sparse
> direct solver, depends on the problem.
> 

So your running times with the sparse direct solver are around 10-20 sec? Thats slightly faster than what my c code would give for a 512 by 512 problem (computer with i5 processor, but the absorbing boundaries more or less double the amount of work)

-Ilan
> --Nasser

On 2/21/2011 7:55 AM, Ilan wrote:

>
>> For school projects, I never have a problem with performance, using Matlab 2010a,
>> new PC, but then again, the pde solvers I write all run on small size grid (512 by 512 is
>> largest I tried). 2-3 minutes to run. iterative solver. much faster using sparse
>> direct solver, depends on the problem.
>>
>
> So your running times with the sparse direct solver are around 10-20 sec?
>  Thats slightly faster than what my c code would give for a 512 by 512 problem
>  (computer with i5 processor, but the absorbing boundaries more or less double the amount of work)
>

here is some numbers for one HW, this is for solving poisson pde on 2D, Dirichlet B.C.
using Matlab 2010a, just to give you an idea. on i7 intel CPU, 4 cores, 8 GB RAM.

(these are cputime() )

N=128, sparse A\b time=0.25 seconds, iterative solver SOR=0.6 seconds
N=256,sparse A\b time=1.5 seconds, iterative solver SOR=5 seconds
N=512,  sparse A\b time=5.5 seconds, iterative solver SOR=50 seconds
N=1024, sparse A\b time=27 seconds, iterative solver SOR=530 seconds

--Nasser

> here is some numbers for one HW, this is for solving poisson pde on 2D, Dirichlet B.C.
> using Matlab 2010a, just to give you an idea. on i7 intel CPU, 4 cores, 8 GB RAM.
> 
> (these are cputime() )
> 
> N=128, sparse A\b time=0.25 seconds, iterative solver SOR=0.6 seconds
> N=256,sparse A\b time=1.5 seconds, iterative solver SOR=5 seconds
> N=512,  sparse A\b time=5.5 seconds, iterative solver SOR=50 seconds
> N=1024, sparse A\b time=27 seconds, iterative solver SOR=530 seconds
> 
> --Nasser

Thanks for the data. 

The equation I reffered to is a bit different, it is u_tt = c^2( u_xx - f), the one you give data for is u_xx = f  (Poisson's eq). So its not directly comparable. But certainly the answer can be given only by such timing data. I will hopefuly try to make such tables myself and post them here (but this is going to take at least a few weeks). 

Thanks again!

-Ilan