You are now following this question
- You will see updates in your followed content feed.
- You may receive emails, depending on your communication preferences.
Open boundary in wave equation with “applyBoundaryCondition”
8 views (last 30 days)
Show older comments
Yaser Dehghan
on 14 Aug 2022
I am implementing a simulation of the wave equation using "solvepde” function in MATLAB. the wave is travelling to right from edge E1 and propagate off forever the edge E3 without being disturbed. how can I specify open boundary (full absorbing/nonreflecting boundary) condition with “applyBoundaryCondition” command?
I am aware there are many academic papers discussing nonreflecting/absorbing boundary conditions but most seem to focus on analytic solutions. On open boundary we have:
I cannot figure out how to implement nonreflecting boundaries numerically in my simulation. This is the code of incident wave on edge E1:
fun1=@(location,state) 2*sin(2*pi/T*state.time)
applyBoundaryCondition(model,'dirichlet','Edge',1,'u',fun1)
Does anyone know how can I adjust the “applyBoundaryCondition” to have open boundaries on edge E3?
19 Comments
Torsten
on 14 Aug 2022
As far as I know, you only have the choice between Dirichlet, Neumann and Mixed for setting up boundary conditions. None of these reflects an "Open Boundary".
Yaser
on 14 Aug 2022
as showen as in figure on open boundary(edge E3) we have du/dx=du/dt (if c=1). how does Dirichlet, Neumann conditin reflect this condition in boundary?
Torsten
on 14 Aug 2022
As I said, MATLAB does not offer an appropriate boundary condition type to set du/dt = c*du/dx at a boundary.
Use a tool especially designed to cope with hyperbolic PDEs, e.g. CLAWPACK available under
Yaser
on 15 Aug 2022
there are sum commans like diff and gradient which in my opinion can reflect the du/dt, but i get some errors as follow:
fun=@(location,state) diff(state.u,state.time)
applyBoundaryCondition(model,'neumann','Edge',3,'g',fun,'q',0);
error: Difference order N must be a positive integer scalar.
Or
fun=@(location,state) gradient(state.u)/gradient(state.time)
applyBoundaryCondition(model,'neumann','Edge',3,'g',fun,'q',0);
error: Solution failed to reach the requested end time.
Torsten
on 15 Aug 2022
If it says in the documentation that g can depend on x,t and u, then you will have to accept this.
Yaser
on 16 Aug 2022
how this error can be correct?
Warning: Failure at t=0.000000e+00. Unable to meet integration tolerances without reducing
the step size below the smallest value allowed (7.905050e-323) at time t.
> In ode15s (line 661)
In pde.EquationModel/solveTimeDependent (line 101)
In pde.PDEModel/solvepde (line 57)
In for_dr (line 44)
Error using pde.EquationModel/solveTimeDependent
Solution failed to reach the requested end time.
Error in pde.PDEModel/solvepde (line 57)
[u,dudt] = self.solveTimeDependent(coefstruct, u0, ut0, tlist, ...
Error in for_dr (line 44)
result = solvepde(model,tlist);
Torsten
on 16 Aug 2022
If the code works if you change the boundary condition at the exit to a supported type, you know what the reason for the time integration failure is ...
Yaser
on 19 Aug 2022
Edited: Yaser
on 19 Aug 2022
clc
clear all
close all
m=1
m = 1
d=0
d = 0
c=1
c = 1
a=0
a = 0
f=0
f = 0
% ...........
Domain=[3 4 0 0 10 10 0 12 12 0]';
ns=[];
sf=[];
geom = decsg(Domain)%decsg(Domain,sf,ns)
geom = 7×4
2 2 2 2
0 0 10 10
0 10 10 0
0 12 12 0
12 12 0 0
0 0 0 0
1 1 1 1
model = createpde; % Set model geometry
geometryFromEdges(model,geom)
ans =
AnalyticGeometry with properties:
NumCells: 0
NumFaces: 1
NumEdges: 4
NumVertices: 4
Vertices: [4×2 double]
generateMesh(model,'Hmax',.5,'Hmin',.2);hold on
pdemesh(model)
pdegplot(model,'EdgeLabels','on'),
specifyCoefficients(model,'m',m,'d',d,'c',c,'a',a,'f',f);
T=2
T = 2
fun1=@(location,state) 1*sin( +2*pi/T*state.time)
fun1 = function_handle with value:
@(location,state)1*sin(+2*pi/T*state.time)
bc=applyBoundaryCondition(model,'dirichlet','Edge',1,'u',fun1,'Vectorized','off'); %
% % % ==============================================================================================
fun4=@(location,state) -gradient(state.u)./gradient(location.x)
fun4 = function_handle with value:
@(location,state)-gradient(state.u)./gradient(location.x)
applyBoundaryCondition(model,'neumann','Edge',3,'g',fun4,'q',0);
xlabel x
ylabel y
% initial conditions:
u0=0.1
u0 = 0.1000
ut0=0
ut0 = 0
setInitialConditions(model,u0,ut0);
% initial:
tlist=linspace(0,5);
model.SolverOptions.ReportStatistics ='on';
result = solvepde(model,tlist);
Warning: Failure at t=0.000000e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (7.905050e-323) at time t.
0 successful steps
615 failed attempts
2462 function evaluations
1 partial derivatives
616 LU decompositions
2461 solutions of linear systems
Error using pde.EquationModel/solveTimeDependent
Solution failed to reach the requested end time.
Solution failed to reach the requested end time.
Error in pde.PDEModel/solvepde (line 57)
[u,dudt] = self.solveTimeDependent(coefstruct, u0, ut0, tlist, ...
u = result.NodalSolution;
Yaser
on 21 Aug 2022
Edited: Yaser
on 21 Aug 2022
according to matlab help:
Note that applyBoundaryCondition uses the default Neumann boundary condition with g = 0 and q = 0 for equations for which you do not explicitly specify a boundary condition.
So no boundary condition at edge2 and edge4 wont have any problem.
as mentioned earlier we want to apply open/non reflective boundary condition at edge 4.
Torsten
on 21 Aug 2022
Edited: Torsten
on 21 Aug 2022
as mentioned earlier we want to apply open/non reflective boundary condition at edge 4.
Edge 4 ? I thought it was edge 3.
And as I mentionned earlier: This is not possible with the PDE toolbox.
Or why do you think the time integrator cannot solve your problem ? It's the boundary condition at edge 3 that causes the problems.
Yaser
on 21 Aug 2022
Sory i had mistake. edge 3 is true. we want to apply open/non reflective boundary condition at edge 3.
Further, we dont use Pde toolbox. we are using solvepde that can solve time dependent problem which PDE toolbox cannot. in matlab help we have:
Torsten
on 21 Aug 2022
Edited: Torsten
on 21 Aug 2022
"solvepde" is part of the PDE toolbox.
Look at the left-hand side of the page
to see this.
I do not doubt that the PDE toolbox can solve the wave equation, but not with a non-reflecting boundary condition.
Try a different boundary condition at edge 3 - if "solvepde" works with the changed setting, you know at least the reason why the integrator fails.
Answers (0)
See Also
Categories
Find more on Eigenvalue Problems in Help Center and File Exchange
Tags
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!An Error Occurred
Unable to complete the action because of changes made to the page. Reload the page to see its updated state.
Select a Web Site
Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .
You can also select a web site from the following list
How to Get Best Site Performance
Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.
Americas
- América Latina (Español)
- Canada (English)
- United States (English)
Europe
- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)
- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)
Asia Pacific
- Australia (English)
- India (English)
- New Zealand (English)
- 中国
- 日本Japanese (日本語)
- 한국Korean (한국어)