Difference between revisions of "Wave equation"
From Medusa: Coordinate Free Mehless Method implementation
(Created page with "TODO: Jure MB") |
|||
| Line 1: | Line 1: | ||
| − | + | ==2D wave equation with dirtchlet 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. | ||
| + | |||
| + | <syntaxhighlight lang="c++" line> | ||
| + | // // 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; | ||
| + | </syntaxhighlight> | ||
| + | |||
| + | |||
| + | Next the domain is populated with nodes in acordance with the desired density function. | ||
| + | |||
| + | |||
| + | <syntaxhighlight lang="c++" line> | ||
| + | // 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); | ||
| + | </syntaxhighlight> | ||
Revision as of 19:10, 7 May 2019
2D wave equation with dirtchlet boundary conditions
Consider the time dependent solution to 2D wave equation on annulus shaped domain
\( \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*} \)
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
\( \begin{align*} f(t)= u_o \sin \omega_o t. \end{align*} \)
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.
1 // // identifier to be added to nodes on the inner boundary
2 int CENTRE = -10;
3
4 // build circle domain
5 BallShape<Vec2d> domain({0, 0}, outer_radius);
6 auto discretization = domain.discretizeBoundaryWithStep(dx);
7
8 // build source domain
9 BallShape<Vec2d> empty({0, 0}, inner_radius);
10 auto discretization_empty = empty.discretizeBoundaryWithStep(dx, CENTRE);
11
12 // substract the source domain
13 discretization -= discretization_empty;
Next the domain is populated with nodes in acordance with the desired density function.
1 // Lambda function for setting the density of nodes
2 auto fill_density = [=](const Vec2d& p) {
3 double r = p.norm();
4 double default_value = dx;
5 double dens = default_value;
6 double r1 = 15*inner_radius;
7 double r2 = 0.8*outer_radius;
8 if (r < r1) dens = linear(inner_radius, 0.8*default_value, r1, default_value, r );
9 if (r > r2) dens = linear(r2, default_value, outer_radius, 0.8* default_value, r);
10 return dens;
11 };
12
13 GeneralFill<Vec2d> fill;
14 fill.seed(fill_seed);
15 discretization.fill(fill, fill_density);
16
17 // find support
18 FindClosest find_support(n);
19 discretization.findSupport(find_support);
