Nonlinear Bessel-ish Function Need Help Solving
5 次查看(过去 30 天)
显示 更早的评论
Hi there, I'm looking to solve a nonlinear differential equation that is very similar to a bessel function, but with an added T^4 kick.
where g is a constant, and which simplifies using product rule to
I tried solving it using bvp4c with the boundary conditions before I knew it was similar to a bessel function. (I used this approach for a very similar problem with no issues)
where T sub L is just a positive number specified in the code below. However, bvp4c told me there was an issue with a singular jacobian matrix. Code listed below. Is there something I'm doing wrong here, or is this because bessel functions cause problems for this sort of analysis?
I'm wondering if there's a better way to solve this differential equation using modified bessel functions, because I noticed my equations were similar to those in this youtube tutorial for fin convection heat transfer, and this paper, although I cannot find either of these seperate equations sourced anywhere else, so if there's names for those variations of bessel equations, I would appreciate if you could let me know how to search them up. Regardless, would it be possible to use a shooting algorithm to guess and check values for T, since my rightmost term is T^3 away from being a normal bessel problem? How would I best implement this using Matlab's bessel functions?
Thanks
% Pass variables into functions. Need global variables, although I know
% this is bad practice, I do not currently know how to work around this
global Tb;
global Tr;
global g1;
% Constants
k = 1300; % [W/m/K]
Sigma_Boltzmann = 5.670374419*10^(-8); % [W/(m^2*K^4)] Stefan Boltzmann constant
a = 0.0255/1000;
emissivity = 1;
b = 0.0001; %[m]
Tb = 300; %[K]
% Parameters
L = 0.4; %[m] length of radiator
Tr = 157.74801710; %[K] estimated temp at x0=0, want to iterate this until you get dT/dx0 = 0 at x0=0
N = 40; % Number of points in mesh
g1 = Sigma_Boltzmann*emissivity*L/a/k; % NOTE NO 2 IN FRONT BECAUSE ASSUME RADIATES ONE DIRECTION ONLY
% Boundary Value Problem
x0mesh = linspace(0,L,N);
solinit = bvpinit(x0mesh, @guess);
sol = bvp4c(@bvpfcn, @bcfcn, solinit);
x0vec = sol.x;
Tvec = sol.y(1,:);
Tpvec = sol.yp(1,:);
Tppvec = sol.yp(2,:);
format long
dTdx0atL = Tpvec(end)
qdotstart = k*a*b*Tpvec(1)
figure(1),clf
plot(x0vec, Tvec, '-o')
xlim([0,L])
ylim([0 Tb])
xlabel('x0 distance (m)'), ylabel('T temperature of radiator (K)')
title('Temperature vs Distance Plot')
set(1,'Position',[200 200 800 500])
% FUNCTIONS
function dTdx0 = bvpfcn(x0,T) % equation to solve
global g1;
dTdx0 = zeros(2,1);
dTdx0 = [T(2);+g1*T(1).^4/x0 - T(2)/x0];
end
%--------------------------------
function res = bcfcn(TaR,TaL) % boundary conditions at x0=0 and x0=L respectively
global Tb;
global Tr;
res = [TaR(1)-Tr; TaL(1)-Tb];
end
%--------------------------------
function g = guess(x0) % initial guess for y and y'
g = [-250*x0+300;5*x0];
end
%--------------------------------
1 个评论
David Goodmanson
2020-8-25
Hi Kyle,
I'm afraid the nonlinearity is going to do in any kind of useful relation to bessel functions. Case in point, the function (9*g*x)^(-1/3) is a solution, although not one that satisfies the boundary conditions.
回答(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!