Is it possible to use non-constant Neumann boundary conditions with the parabolic pde solver?
4 次查看(过去 30 天)
显示 更早的评论
I am looking to solve the 2D heat equation T = T(r, theta, t) on a circle. I have a non-constant Neumann BC on the outside of the circle (r = a), where I have a heat flux as a function of theta (-k*dT/dx = q(theta) at r = a).
I know that I can decompose theta in to atan(y/x) -- in general though, I am unclear on how to use either the pdebound or assemb functions to prescribe these boundary conditions.
Thank you for your help!
0 个评论
采纳的回答
Bill Greene
2012-8-31
Hi,
You are definitely on the right track in assuming that creating a pdebound function is a good way to define this BC. I created a simple example below that assumes you have an inward heat flux of 5*sin(theta), defines a pdebound function for this BC (I called it boundaryConditions), and then uses that to solve the heat equation on a circle. Hopefully this example will get you past some of the sticky issues.
Regards,
Bill
function transientHeatCircle( )
radius = 2;
gdm=[1 0 0 radius]';
g = decsg(gdm, 'C1', ('C1')');
[p,e,t]=initmesh(g);
c = 1; d = 2; a = 0; f = 0;
b = @boundaryConditions;
u0=0; tlist=0:.001:1;
u=parabolic(u0, tlist, b,p,e,t,c,a,f,d);
figure; pdeplot(p, e, t, 'xydata', u(:,end), 'mesh', 'on'); axis equal;
title 'Final Temperature Distribution'
n=4; %grid point at theta=pi/2
figure; plot(tlist, u(n,:)); grid on;
title(sprintf('Temperature at (%3.1f,%3.1f) as a Function of Time', ...
p(1,n), p(2,n)));
end
function [ q, g, h, r ] = boundaryConditions( p, e, u, time )
N = 1;
ne = size(e,2);
q = zeros(N^2, ne);
% calculate coordinates of edge mid-points
xy1 = p(:,e(1,:));
xy2 = p(:,e(2,:));
xyMidEdge = (xy1+xy2)/2;
g = 5*sin(atan2(xyMidEdge(2,:),xyMidEdge(1,:)));
h = zeros(N^2, 2*ne);
r = zeros(N, 2*ne);
end
2 个评论
Bill Greene
2012-9-2
My code handled the simple case where you have a single BC expression that applies to all edges on the boundary. A more general version of a boundary function might look something like this.
function [ q, g, h, r ] = boundFunc( p, e, u, time )
N = 1; % number of PDEs in the system
ne = size(e,2);
q = zeros(N^2, ne);
g = zeros(N, ne);
h = zeros(N^2, 2*ne);
r = zeros(N, 2*ne);
for i=1:ne
ei = e(5,i);
if(ei == 1)
% apply appropriate BCs on edge 1
else if(ei == 3)
% apply appropriate BCs on edge 3
else if(ei == ...)
% apply appropriate BCs on this edge
end
end
end
更多回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Eigenvalue Problems 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!