error while solving the coupled ode
1 Ansicht (letzte 30 Tage)
Ältere Kommentare anzeigen
Hello,
I have been trying to solve the following ode
but a "singular jacobian error" is recurring and is very sensitive to parameters. I think it's because of the guess function issue but I'm not very sure of it. I tried many guess functions but didn't seem to work. If there is any way out or suggestion to overshoot this error, please do help.
Thanks!
4 Kommentare
Akzeptierte Antwort
Torsten
am 21 Feb. 2023
% Defining parameters
delta = 0.02; % Lower integral bound
R = 5; % Upper integral bound
theta = 0; % ArcTan(q/g)
maxPoints = 1e6; % Maximum numer of grid point used by bvpc4
initialPoints = 100; % Number of initial grid points used by bvpc4
tol = 1e-3; % Maximum allowed relative error.
L = 10;
N = 1;
n = 1;
m = 0;
g = 5;
lambda = 1;
% Boundary conditions
y0 = [0, -1, N*pi, 0];
% Initial conditions
A = @(r) (1-tanh(((L*r)/R)-(L/3)))/2;
dA = cosh(theta)*(coth(delta)-delta*csch(delta).^2);
F = @(r) (1+tanh(((L*r)/R)-(L/3)))/2;
dF = (1-delta*coth(delta))*csch(delta);
solinit = bvpinit(linspace(delta, R, initialPoints), [1 1 1 1]);
% Solves system using bvpc4
options = bvpset('RelTol', tol, 'NMax', maxPoints); % This function sets the allowed
%relative error and maximum number of grid points.
sol = bvp4c(@(r, y) heatGauge(r, y, lambda, g, m, n), @(ya, yb) bcheatGauge(ya, yb, y0),...
solinit, options);
r = linspace(delta, R, 1e4);
y = deval(sol, r);
plot(r,y(1,:),r, y(2,:))
grid on
function dy = heatGauge(r, y, lambda, g, m, n)
a = y(1);
f = y(2);
adot = y(3);
fdot = y(4);
atilde = n*(a+1.0)/r;
dy(1) = adot;
dy(2) = fdot;
dy(3) = a/r + g^2*(1+a)*(1+lambda^2*fdot^2)*sin(f)^2;
dy(4) = (-fdot/r*((2*n*adot-atilde)*lambda^2*atilde*sin(f)^2+lambda^2*r*fdot*atilde^2*sin(f)*cos(f)+1)...
+atilde^2*sin(f)*cos(f)+m^2*sin(f))/(1+lambda^2*atilde^2*sin(f)^2);
end
function res = bcheatGauge(ya, yb, y0)
res = [ya(1) - y0(1);yb(1) - y0(2);ya(2) - y0(3);yb(2) - y0(4)];
end
12 Kommentare
Torsten
am 24 Feb. 2023
These equations look a lot better than the ones you posted first.
Did you try to solve them ?
Weitere Antworten (0)
Siehe auch
Kategorien
Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!