Using a given array as an initial guess for bvp4c using bvpinit

7 visualizaciones (últimos 30 días)
Marcus Rosales
Marcus Rosales el 14 de Mayo de 2024
Comentada: Torsten el 15 de Mayo de 2024
I tried using the solution, G, to the "free" version of a given ODE (given in code with a=0) to initialize the bvp4c solver. Given a mesh xint and a function guess(x) (made by interp1, G and xint), I used the following:
solinit = bvpinit(xint, @guess);
The details are in the attached code (forgot to include: xl=-50, xr=2). The free case works with @guess replaced with [0 0].
I keep getting the error: "Index exceeds matrix dimensions" when the solver attempts the ode solution. I have tried a lot to get this to work, but things keep failing. Any clues on how I am probbaly screwing up "solinit"?
function [Sxint, x] =bvp4c_mathworks_half_VAV_2(xl,xr, N, theta, r0, Guess)
global tt
global r
global G
global xx
xint = linspace(xl,xr,N); %BC at [-inf, inf] replaced with [xl, xr]
xx = xint;
tt = theta; %parameters: winding number and vortex radius, say pi/6, 16
r = r0;
G = zeros(1,length(Guess)); %What I want to initialize with
solinit = bvpinit(xint, @guess);
sol = bvp4c(@myode,@mybc,solinit);
Sxint = deval(sol,xint);
end
function dy = myode(t,y)
size(y)
size(t)
a = 0*10^(-5); %a reduced coupling to pdw times uniform gap here
dy(1,1) = y(2);
dy(2,1) = 1/2*(1-exp(2*t)*cos(2*y(1)))*sin(2*y(1))+a*a*exp(2*t)*cos(f(t))*cos(2*y(1));
end
function phase_diff = f(t)
global tt;
global r;
dd = 80;
phase_diff = pi*tanh(r*exp(t)/dd)-tt;
end
function res = mybc(ya,yb)
res = [ya(1)
yb(1)-pi/4];
end
function [y1, y2] = guess(x)
global G;
global xx;
ee = 10^(-10);
y1 = interp1(xx,G,x).';
y2 = (interp1(xx,G,x+ee).'-y1)/ee;
%y2 = gradient(G).';
end
  3 comentarios
Marcus Rosales
Marcus Rosales el 15 de Mayo de 2024
I originally used the gradient here. Will this also give me some issues? I figured the end point might run into some issues since there is nothing to interpolate with (should be zero though).
Torsten
Torsten el 15 de Mayo de 2024
If you can supply G, don't you have access to dG/dx ?
But the bad approximation for dG/dx is not the reason for the error message. As said: If xx and G have the same length, your code should work.

Iniciar sesión para comentar.

Respuestas (0)

Categorías

Más información sobre Numerical Integration and Differential Equations en Help Center y File Exchange.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by