Wave equation

2D wave equation with Dirichlet boundary conditions

Consider the time dependent solution to 2D wave equation on annulus shaped domain

<math> \begin{align*}

   	\frac{ \partial^2 u}{\partial t^2} &= c^2 \nabla^2 u  &&\text{in } \Omega, \\
   	  u &= 0           &&\text{on } \partial \Omega_O,\\
         u &= f(t)           &&\text{on } \partial \Omega_I,

\end{align*} </math>

where $\partial\Omega_I$ denotes the inner and $\partial\Omega_O$ the outer boundary of the domain. Through the boundary condition on the inner boundary the source is introduced to the problem as a function of time

<math> \begin{align*} f(t)= u_o \sin \omega_o t. \end{align*} </math>

First the domain is constructed by subtracting a smaller circle domain from a larger one. Boundaries of the domain are populated in the same step.

// // identifier to be added to nodes on the inner boundary
int CENTRE = -10;

 // build circle domain
BallShape<Vec2d> domain({0, 0}, outer_radius);
auto discretization = domain.discretizeBoundaryWithStep(dx);

// build source domain
BallShape<Vec2d> empty({0, 0}, inner_radius);
auto discretization_empty = empty.discretizeBoundaryWithStep(dx, CENTRE);  
    
// substract the source domain
discretization -= discretization_empty;

Next the domain is populated with nodes in acordance with the desired density function.

// Lambda function for setting the density of nodes
auto fill_density = [=](const Vec2d& p) {
    double r = p.norm();
    double default_value = dx;
    double dens = default_value;
     double r1 = 15*inner_radius;
     double r2 = 0.8*outer_radius;
     if (r < r1) dens = linear(inner_radius, 0.8*default_value, r1, default_value, r );
     if (r > r2) dens = linear(r2, default_value, outer_radius, 0.8* default_value, r);
     return dens;
};

GeneralFill<Vec2d> fill;
fill.seed(fill_seed);
discretization.fill(fill, fill_density);

 // find support
FindClosest find_support(n);
discretization.findSupport(find_support);