This example is an application from wave scattering theory. You send a plane wave into a wave guide that has a sharp bend the size of which is comparable to the wave-length. Due to the reflections on the walls of the wave-guide, the plane wave is scattered. Say you are interested to see how the shape of the bending affects the solution at the exit of the wave-guide. You compare 3 different shapes by computing the energy of the scattered waves at the exit.
You will see how working with the Graphical User Interface (GUI), which is specialized for PDE solving, is combined with the command-line mode, to make full use of MATLAB's general programming environment.
A wave guide
Note: the equation has complex solutions due to the boundary conditions.
Three geometries that bend a planar wave
The first two domains are polygons, and the third is defined by the set formula P1+(C1-SQ1).
Set the corresponding boundary conditions and export the information from the GUI to your MATLAB workspace. Call the geometry and boundary matrices g1, b1, g2, b2, g3, b3 denoting the first, second and third domain.
The PDE Toolbox offers you the possibility to choose the environment that suits you best at different stages of your work. In this case you may want to split your work as follows:
You want to have 8 - 10 grid-points per wavelength If you put the triangle size of the initial mesh to 0.2 (a wavelength), then 3 uniform refinements produce the desired density of mesh-points.
[p,e,t] = initmesh(g3,'Hmax',0.2); for j = 1:3 [p,e,t] = refinemesh(g3,p,e,t); endThe p, e, t matrices (points, edges, triangles) describe the triangular mesh. They are created by the function
u = assempde(b3,p,e,t,1,(2*pi*i*5)^2,0); pdeplot(pde,'xydata',u)The pdeplot command produces the following plot:
Scattering in a bended wave-guide.
The vector u is a discrete approximation of the solution to
-div( grad(u) ) - 100*pi^2*u=0.
The entries of u are the values of the solution corresponding to the mesh-points given in p. Find the nodes at the exit, y == 0.8, and then the corresponding values in the solution:
index = find( p(2,:)==0.8 ); x = p(1,index); v = u(index);
x contains the grid-points at the exit, while v contains the complex values of the solution. Because the original mesh is unstructured, you must sort p (and v )
[x,index] = sort(x); v = v(index); plot(x,abs(v));
Absolute values of the solutions at the exit.
At this stage, x and v are ordinary vectors containing grid-points and grid-values. The whole array of MATLAB facilities are available to process them. Compute the energy of the scattered waves at the exit and compare them with the energy of the ingoing wave (the ingoing energy is 0.4). The comparison is done by the following code, using the trapezoidal rule:
N = length(x); w = abs(v).^2; E = sum(0.5*(w(2:N)+w(1:N-1)).*diff(x))/0.4;The L-shape bending recovers only 11.02% of the energy, compared to 80.42% and 89.53% in the plane mirror and curved mirror cases.